Merge pull request #150 from tomrink/master

Trajectory: make sure we get streamlines with careful selection of control parameters.

core/src/visad/TrajectoryManager.java

@@ -307,8 +307,10 @@ public class TrajectoryManager {
         uInterp = new NoneInterpolator(numSpatialPts);
         vInterp = new NoneInterpolator(numSpatialPts);
         wInterp = new NoneInterpolator(numSpatialPts);
-        trajDoIntrp = false;
-        numIntrpPts = 1;
+      }
+      
+      if (!trajDoIntrp) {
+         numIntrpPts = 1;
       }
       
       values0 = null;
@@ -552,7 +554,7 @@ public class TrajectoryManager {
              intrpW = mean(intrpW, intrpW_1);
            }
 
-           for (int t=0; t<numTrajectories; t++) {
+           for (int t=0; t<numTrajectories; t++) { //Euler method if NOT HySplit
              Trajectory traj = trajectories.get(t);
              traj.currentTimeIndex = direction*i;
              traj.currentTime = direction*times[i];

core/src/visad/TrajectoryParams.java

@@ -84,6 +84,7 @@ public class TrajectoryParams {
   boolean terrainFollowEnabled = true;
   boolean trcrStreamingEnabled = false;
   boolean saveTracerLocations = false;
+  boolean trajDoIntrp = true;
   
   int trajForm = LINE;
   float cylWidth = 0.00014f;
@@ -133,6 +134,7 @@ public class TrajectoryParams {
     this.trcrStreamingEnabled = params.getTracerStreamingEnabled();
     this.saveTracerLocations = params.getSaveTracerLocations();
     this.timeStepScaleFactor = params.getTimeStepScaleFactor();
+    this.trajDoIntrp = params.getTrajDoIntrp();
   }
 
   public TrajectoryParams(double trajVisibilityTimeWindow, double trajRefreshInterval, int numIntrpPts, int startSkip, SmoothParams smoothParams) {
@@ -386,6 +388,14 @@ public class TrajectoryParams {
      this.timeStepScaleFactor = fac;
   }
   
+  public boolean getTrajDoIntrp() {
+     return this.trajDoIntrp;
+  }
+  
+  public void setTrajDoIntrp(boolean yesno) {
+     this.trajDoIntrp = yesno;
+  }
+  
   public boolean equals(Object obj) {
     if (obj == null || !(obj instanceof TrajectoryParams)) {
       return false;
@@ -446,6 +456,9 @@ public class TrajectoryParams {
       else  if (this.timeStepScaleFactor != trajParams.timeStepScaleFactor) {
          return false;
       }
+      else  if (this.trajDoIntrp != trajParams.trajDoIntrp) {
+         return false;
+      }      
     }
     return true;
   }

core/src/visad/XTrackScanLatLonSet.java

+/*
+ * To change this license header, choose License Headers in Project Properties.
+ * To change this template file, choose Tools | Templates
+ * and open the template in the editor.
+ */
+package visad;
+
+/**
+ *
+ * @author rink
+ */
+public class XTrackScanLatLonSet extends GriddedLatLonSet {
+   
+   GriddedLatLonSet[] scanSets;
+   int numOfScans;
+   int linesPerScan;
+   
+   int lastSetIdx;
+   
+   public XTrackScanLatLonSet(MathType type, float[][] samples, int lengthX, int lengthY, int linesPerScan) throws VisADException {
+      super(type, samples, lengthX, lengthY, null, null, null, false, false);
+      
+      this.linesPerScan = linesPerScan;
+      numOfScans = lengthY/linesPerScan;
+      scanSets = new GriddedLatLonSet[numOfScans];
+      int scanLen = linesPerScan*lengthX;
+      
+      float[] scanLonArray = new float[linesPerScan*lengthX];
+      float[] scanLatArray = new float[linesPerScan*lengthX];
+      
+      for (int k=0; k<numOfScans; k++) {
+         System.arraycopy(lons, k*scanLen, scanLonArray, 0, scanLen);
+         System.arraycopy(lats, k*scanLen, scanLatArray, 0, scanLen);
+         
+         scanSets[k] = new GriddedLatLonSet(RealTupleType.SpatialEarth2DTuple, new float[][] {scanLonArray, scanLatArray}, lengthX, linesPerScan);
+      }
+   }
+   
+   public float[][] valueToGrid(float[][] value) throws VisADException {
+      return valueToGrid(value, null);
+   }
+   
+   public float[][] valueToGrid(float[][] value, int[] guess) throws VisADException {
+    if (value.length < DomainDimension) {
+      throw new SetException("Gridded2DSet.valueToGrid: value dimension " +
+                             value.length + " not equal to Domain dimension " +
+                             DomainDimension);
+    }
+    if (ManifoldDimension < 2) {
+      throw new SetException("Gridded2DSet.valueToGrid: Manifold dimension " +
+                             "must be 2, not " + ManifoldDimension);
+    }
+    if (Length > 1 && (Lengths[0] < 2 || Lengths[1] < 2)) {
+      throw new SetException("Gridded2DSet.valueToGrid: requires all grid " +
+                             "dimensions to be > 1");
+    }
+    if (guess != null && guess.length != 2) {
+      throw new SetException("Gridded2DSet.valueToGrid: guess length " 
+          + guess.length + " must equal 2");
+    }
+    
+    
+    int length = Math.min(value[0].length, value[1].length);
+    float[][] grid = new float[ManifoldDimension][length];
+
+    // (gx, gy) is the current grid box guess
+    int gx = (LengthX-1)/2;
+    int gy = (TrackLen-1)/2;
+    
+    if (guess != null && guess[0] >= 0 && guess[1] >= 0) {
+      gx = guess[0];
+      gy = guess[1];
+    }
+    else if (lgxy[0] != -1 && lgxy[1] != -1) {
+       gx = lgxy[0];
+       gy = lgxy[1];
+    }
+
+
+    int setIdx = gy/linesPerScan;
+    //int gy_set = gy % linesPerScan;
+    lastSetIdx = setIdx;
+    
+    int idx = gy*LengthX+gx;
+    if (Float.isNaN(lons[idx]) || Float.isNaN(lats[idx]) || (Math.abs(lats[idx]) > 90)) {
+       throw new VisADException("initial grid box guess cannot be invalid or missing");
+    }
+    
+    float[][] lonlat = new float[2][1];
+    int[] guess_set = new int[2];
+    
+    for (int i=0; i<length; i++) {
+
+      if (Length == 1) {
+        if (Float.isNaN(value[0][i]) || Float.isNaN(value[1][i])) {
+           grid[0][i] = grid[1][i] = Float.NaN;
+        } else {
+           grid[0][i] = 0;
+           grid[1][i] = 0;
+        }
+        continue;
+      }
+
+      // test for missing
+      if ( (i != 0) && grid[0][i-1] != grid[0][i-1] ) {
+        // use last valid
+        if (lgxy[0] != -1 && lgxy[1] != -1) {
+          gx = lgxy[0];
+          gy = lgxy[1];
+        }
+      }
+      
+      // if the loop doesn't find the answer, the result should be NaN
+      grid[0][i] = grid[1][i] = Float.NaN;
+      
+      float targetLon = value[lonI][i];
+      float targetLat = value[latI][i];
+      lonlat[lonI][0] = targetLon;
+      lonlat[latI][0] = targetLat;
+      float[][] gxgy;
+      
+      while (lastSetIdx >= 0 && lastSetIdx < numOfScans) {
+         
+        GriddedLatLonSet scanSet = scanSets[lastSetIdx];
+        gxgy = scanSet.valueToGrid(lonlat, guess_set);
+      
+        if (gxgy[0][0] == Float.NaN || gxgy[1][0] == Float.NaN) {
+           int gx_set = scanSet.lgxy[0];
+           int gy_set = scanSet.lgxy[1];
+           if (gy_set == 0) {
+              lastSetIdx -= 1;
+              guess_set[0] = gx_set;
+              guess_set[1] = linesPerScan-1;
+           }
+           else if (gy_set == linesPerScan-1) {
+              lastSetIdx += 1;
+              guess_set[0] = gx_set;
+              guess_set[1] = 0;
+           }
+           else if (gx_set == 0 || gx_set == LengthX-1) {
+              break;
+           }
+        }
+        else {
+           grid[0][i] = gxgy[0][0];
+           grid[1][i] = (lastSetIdx-1)*linesPerScan + gxgy[1][0];
+           lgxy[0] = scanSet.lgxy[0];
+           lgxy[1] = scanSet.lgxy[1] + (lastSetIdx-1)*linesPerScan;
+           break;
+        }
+      }
+      
+    }
+    
+    return grid;
+  }
+   
+}

core/src/visad/java3d/ShadowFunctionOrSetTypeJ3D.java

@@ -1432,10 +1432,12 @@ System.out.println("Texture.BASE_LEVEL_LINEAR = " + Texture.BASE_LEVEL_LINEAR);
     
     // We don't really have a time interval defined beyond the next to last point
     int computeLength = dataDomainLength-1;
-    if (trajParams.getInterpolationMethod() == TrajectoryParams.InterpolationMethod.None) {
-       if (trajParams.getMethod() == TrajectoryParams.Method.RK4) { // need 3 time steps
-          computeLength = dataDomainLength-2;
-       }
+    if (trajParams.getMethod() == TrajectoryParams.Method.RK4) { // need 3 time steps
+      computeLength = dataDomainLength-2;
+    }
+    else if (trajParams.getMethod() == TrajectoryParams.Method.Euler && 
+            ((trajParams.getInterpolationMethod() == TrajectoryParams.InterpolationMethod.None) || !trajParams.getTrajDoIntrp())) {
+      computeLength = dataDomainLength;
     }
     
     for (int k=0; k<computeLength; k++) {