package edu.stanford.rsl.conrad.geometry.splines;
import java.util.ArrayList;
import java.util.Iterator;
import edu.stanford.rsl.conrad.geometry.AbstractShape;
import edu.stanford.rsl.conrad.geometry.shapes.simple.PointND;
import edu.stanford.rsl.conrad.utils.Configuration;
import edu.stanford.rsl.conrad.utils.RegKeys;
public class MotionDefectTimeVariantSurfaceBSpline extends
TimeVariantSurfaceBSpline {
/**
*
*/
private static final long serialVersionUID = -1439361309222709769L;
private PointND min, max;
private double samplingStep = 0.1;
private double phaseShift = 0.0;
private ArrayList<PointND> defectPoints = new ArrayList<PointND>();
private ArrayList<PointND> interpolationPoints = new ArrayList<PointND>();
public MotionDefectTimeVariantSurfaceBSpline(
ArrayList<SurfaceBSpline> splines, PointND min, PointND max) {
super(splines);
this.min = new PointND(0,0,0);
this.max = new PointND(0,0,0);
// Make sure that min is always smaller than max.
for (int i =0 ; i < 3; i++){
this.min.set(i, Math.min(min.get(i), max.get(i)));
this.max.set(i, Math.max(min.get(i), max.get(i)));
}
SurfaceBSpline first = (SurfaceBSpline)timeVariantShapes.get(0);
for (double v = 0; v < first.vSplines.length; v++){
ArrayList<PointND> controlPointsU = first.vSplines[(int)v].getControlPoints();
for (double u = 0; u < controlPointsU.size(); u++){
boolean inside = true;
for (int i =0 ; i < 3; i++){
double coord = controlPointsU.get((int)u).get(i);
if (!(coord >min.get(i) && coord < max.get(i))){
inside = false;
break;
}
}
if (inside) {
defectPoints.add(new PointND(
v / first.vSplines.length,
u / controlPointsU.size()));
interpolationPoints.add(new PointND(
v / first.vSplines.length,
u / controlPointsU.size(),
1));
} else {
interpolationPoints.add(new PointND(
v / first.vSplines.length,
u / controlPointsU.size(),
0));
}
}
}
double sampleBasis = 0.1;
String key = Configuration.getGlobalConfiguration().getRegistry().get(RegKeys.XCAT_HEART_MOTION_DEFECT_FLEXIBILITY);
if (key!= null) {
sampleBasis = Double.parseDouble(key);
}
setSamplingStep(sampleBasis);
key = Configuration.getGlobalConfiguration().getRegistry().get(RegKeys.XCAT_HEART_MOTION_DEFECT_PHASE_SHIFT);
if (key!= null) {
phaseShift = Double.parseDouble(key);
}
}
public MotionDefectTimeVariantSurfaceBSpline(MotionDefectTimeVariantSurfaceBSpline mds) {
super(mds);
min = (mds.min != null) ? mds.min.clone() : null;
max = (mds.max != null) ? mds.max.clone() : null;
samplingStep = mds.samplingStep;
phaseShift = mds.phaseShift;
if (mds.defectPoints != null){
Iterator<PointND> it = mds.defectPoints.iterator();
defectPoints = new ArrayList<PointND>();
while (it.hasNext()) {
PointND p = it.next();
defectPoints.add((p!=null) ? p.clone() : null);
}
}
else{
defectPoints = null;
}
if (mds.interpolationPoints != null){
Iterator<PointND> it = mds.interpolationPoints.iterator();
interpolationPoints = new ArrayList<PointND>();
while (it.hasNext()) {
PointND p = it.next();
interpolationPoints.add((p!=null) ? p.clone() : null);
}
}
else{
interpolationPoints = null;
}
}
private boolean insideDefect(PointND pointUV){
for (int i=0; i < defectPoints.size(); i++){
double distance = pointUV.euclideanDistance(defectPoints.get(i));
if (distance < samplingStep) return true;
}
return false;
}
private double interpolateWeight(PointND pointUV){
double weightH = 0;
double interpolationWeight = 0;
for (int i= 0; i < interpolationPoints.size(); i++){
PointND localPoint = interpolationPoints.get(i);
double distance = Math.sqrt(
Math.pow(localPoint.get(0) - pointUV.get(0), 2) +
Math.pow(localPoint.get(1) - pointUV.get(1), 2));
double localWeight = Math.exp(-0.5 * Math.pow(distance/samplingStep,2));
interpolationWeight += interpolationPoints.get(i).get(2) * localWeight;
weightH += localWeight;
}
return interpolationWeight/weightH;
}
public PointND evaluate(double u, double v, double t){
double [] p = new double [timeVariantShapes.get(0).getControlPoints().get(0).getDimension()];
double internal = ((t * timeSpline.getKnots().length) -3.0);
int numPts = timeSpline.getControlPoints().size();
double step = internal- Math.floor(internal);
PointND result = null;
if (insideDefect(new PointND(u,v))){
PointND pointNotMoving = ((SurfaceBSpline)timeVariantShapes.get(0)).evaluate(u, v);
if (phaseShift != 0){
t-=phaseShift;
t = t%1.0;
internal = ((t * timeSpline.getKnots().length) -3.0);
double [] weights = UniformCubicBSpline.getWeights(step);
double [] p2 = new double [timeVariantShapes.get(0).getControlPoints().get(0).getDimension()];
for (int i=0;i<4;i++){
double [] loc = (internal+i < 0)?
((SurfaceBSpline)timeVariantShapes.get(0)).evaluate(u, v).getCoordinates()
: (internal+i>=numPts)?
((SurfaceBSpline)timeVariantShapes.get(numPts-1)).evaluate(u, v).getCoordinates()
: ((SurfaceBSpline)timeVariantShapes.get((int) (internal+i))).evaluate(u, v).getCoordinates();
for (int j = 0; j < loc.length; j++){
p2[j] += (loc[j] * weights[i]);
}
}
pointNotMoving = new PointND(p2);
}
// Unbewegte Punkt p_u
double weight = interpolateWeight(new PointND(u,v));
double [] weights = UniformCubicBSpline.getWeights(step);
for (int i=0;i<4;i++){
double [] loc = (internal+i < 0)?
((SurfaceBSpline)timeVariantShapes.get(0)).evaluate(u, v).getCoordinates()
: (internal+i>=numPts)?
((SurfaceBSpline)timeVariantShapes.get(numPts-1)).evaluate(u, v).getCoordinates()
: ((SurfaceBSpline)timeVariantShapes.get((int) (internal+i))).evaluate(u, v).getCoordinates();
for (int j = 0; j < loc.length; j++){
p[j] += (loc[j] * weights[i]);
}
}
PointND pointMoving = new PointND(p);
pointMoving.getAbstractVector().multiplyBy(1-weight);
pointNotMoving .getAbstractVector().multiplyBy(weight);
pointNotMoving.getAbstractVector().add(pointMoving.getAbstractVector());
result = pointNotMoving;
} else {
// Bewegte Punkt p_b
double [] weights = UniformCubicBSpline.getWeights(step);
for (int i=0;i<4;i++){
double [] loc = (internal+i < 0)?
((SurfaceBSpline)timeVariantShapes.get(0)).evaluate(u, v).getCoordinates()
: (internal+i>=numPts)?
((SurfaceBSpline)timeVariantShapes.get(numPts-1)).evaluate(u, v).getCoordinates()
: ((SurfaceBSpline)timeVariantShapes.get((int) (internal+i))).evaluate(u, v).getCoordinates();
for (int j = 0; j < loc.length; j++){
p[j] += (loc[j] * weights[i]);
}
}
PointND pointMoving = new PointND(p);
result = pointMoving;
}
return result;
}
public void setSamplingStep(double samplingStep) {
this.samplingStep = samplingStep;
}
public double getSamplingStep() {
return samplingStep;
}
@Override
public AbstractShape clone() {
return new MotionDefectTimeVariantSurfaceBSpline(this);
}
}
/*
* Copyright (C) 2010-2014 Andreas Maier
* CONRAD is developed as an Open Source project under the GNU General Public License (GPL).
*/