# HG changeset patch
# User Sascha L. Teichmann <sascha.teichmann@intevation.de>
# Date 1264443919 0
# Node ID 44415ae01ddb57caddf6874d25f5bc7e6a715fd9
# Parent  f598702b2a10a2c18bf10fc7269105aa3d89d5c5
Fixed issue gnv/issue159


gnv-artifacts/trunk@624 c6561f87-3c4e-4783-a992-168aeb5c3f6f

diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/ChangeLog
--- a/gnv-artifacts/ChangeLog	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/ChangeLog	Mon Jan 25 18:25:19 2010 +0000
@@ -1,3 +1,46 @@
+2010-01-25	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
+
+	Fix/Workaround for gnv/issue159
+
+	* src/main/java/de/intevation/gnv/state/profile/horizontalcrosssection/HorizontalCrossSectionMeshOutputState.java,
+	  src/main/java/de/intevation/gnv/state/profile/verticalcrosssection/VerticalCrossSectionOutputState.java,
+	  src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileMeshCrossOutputState.java:
+	  Deactivated the index buffer by default because applying it gives to
+	  less data points in some circumstances. When processing the 
+	  "horizontalen Schnittprofile", the "Profilschnitte" and the
+	  "Horizontalschnitte" the whole set of mesh cells for
+	  the given parameters are retrieved. This results in much
+	  more traffic from database (~ 10x) but seems not to be much slower.
+	  The memory consumption increases accordingly.
+
+	  To turn the index buffer on again use the boolean system properties:
+
+		  "gnv.horizontal.profile.mesh.cross.index.buffer"
+		  "gnv.vertical.cross.section.index.buffer"
+		  "gnv.horizontal.profile.mesh.cross.index.buffer"
+
+	* src/main/java/de/intevation/gnv/math/Interpolation2D.java: Because much
+	  more data is processed a more sophisticated point culling strategy is
+	  implemented now. 
+	  
+	  If the number of points are higher than a given threshold (system property 
+	  "gnv.interpolation2d.cull.point.threshold", default: 1000) a culling is 
+	  performed. Two cases:
+
+	  a) given a path P ("horizontales Schnittprofil", "Profilschnitt")
+	     Calculates a bounding box B for P and continues with (b).
+
+	  b) a bounding box B ("Horizontalschnitt", and above after (a))
+	     Calculates a bounding box C for the input points. If the
+		 area of B is greater than 80% of the area of C no culling is
+		 performed. Else B is extended by 10% in each direction.
+		 If the new area of B is smaller than 10% of the area of C
+		 B is adjusted to have at least 10% the size of C.
+
+	  All input points that are not inside B are culled. Afterwards the
+	  i/j grid cell units are created and pushed into a spatial index
+	  as described below.
+
 2010-01-25  Ingo Weinzierl <ingo.weinzierl@intevation.de>
 
 	  Issue157
diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java
--- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java	Mon Jan 25 18:25:19 2010 +0000
@@ -22,9 +22,6 @@
 {
     private static Logger log = Logger.getLogger(AreaInterpolation.class);
 
-    public static final int CULL_POINT_THRESHOLD = Integer.getInteger(
-        "gnv.areainterpolation.cull.point.threshold", 1000);
-
     protected double [] raster;
 
     protected int width;
@@ -61,17 +58,11 @@
             return false;
         }
 
-        Envelope relevantArea = null;
-        
-        if (points.size() > CULL_POINT_THRESHOLD) {
-            relevantArea = new Envelope(boundingBox);
-            relevantArea.expandBy(
-                0.05d*boundingBox.getWidth(),
-                0.05d*boundingBox.getHeight());
-        }
-
         List<GridCell> cells = GridCell.pointsToGridCells(
-            points, relevantArea);
+            points,
+            Interpolation2D.relevantArea(
+                boundingBox,
+                points));
 
         if (cells.isEmpty()) {
             log.warn("no cells to interpolate");
diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java
--- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java	Mon Jan 25 18:25:19 2010 +0000
@@ -28,9 +28,88 @@
     private Interpolation2D() {
     }
 
-    public static Envelope pathBoundingBox(
+    public static final Envelope relevantArea(
         List<? extends Coordinate> path,
-        double extra
+        List<? extends Coordinate> points
+    ) {
+        return relevantArea(path, points, CULL_POINT_THRESHOLD);
+    }
+
+    public static final Envelope relevantArea(
+        Envelope                   pathBBox,
+        List<? extends Coordinate> points
+    ) {
+        return relevantArea(pathBBox, points, CULL_POINT_THRESHOLD);
+    }
+
+    public static final Envelope relevantArea(
+        Envelope                   pathBBox,
+        List<? extends Coordinate> points,
+        int                        threshold
+    ) {
+        return points.size() < threshold
+            ? null
+            : relevantArea(
+                pathBBox,
+                pointsBoundingBox(points));
+    }
+
+    public static final Envelope relevantArea(
+        List<? extends Coordinate> path,
+        List<? extends Coordinate> points,
+        int                        threshold
+    ) {
+        return points.size() < threshold || path.isEmpty()
+            ? null
+            : relevantArea(
+                pointsBoundingBox(path),
+                pointsBoundingBox(points));
+    }
+
+    public static final Envelope relevantArea(
+        Envelope pathBBox,
+        Envelope pointsBBox
+    ) {
+        double pathArea   = pathBBox.getWidth()*pathBBox.getHeight();
+        double pointsArea = pointsBBox.getWidth()*pointsBBox.getHeight();
+
+        if (pathArea > 0.8d*pointsArea) { return null; }
+
+        double nArea = 1.44d * pathArea;
+        if (nArea < 0.1d*pointsArea) nArea = 0.1d*pointsArea;
+        double w = pathBBox.getWidth();
+        double h = pathBBox.getHeight();
+        double [] d = solveQuadratic(1d, w+h, pathArea - nArea);
+
+        if (d == null) { return null; }
+
+        double extra = pos(d);
+
+        pathBBox.expandBy(extra);
+
+        return pathBBox;
+    }
+
+    public static final double [] solveQuadratic(
+        double a, double b, double c
+    ) {
+        double d = b*b - 4d*a*c;
+        if (d < 0d) { return null; }
+
+        d = Math.sqrt(d);
+        a = 1d/(2d*a);
+        b = -b;
+
+        return new double [] { a*(b + d), a*(b - d) };
+    }
+
+    public static final double pos(double [] x) {
+        return x[0] >= 0 ? x[0] : x[1];
+    }
+
+
+    public static Envelope pointsBoundingBox(
+        List<? extends Coordinate> path
     ) {
         int N = path.size();
         Envelope area = new Envelope(path.get(N-1));
@@ -39,10 +118,6 @@
             area.expandToInclude(path.get(i));
         }
 
-        area.expandBy(
-            extra*area.getWidth(),
-            extra*area.getHeight());
-
         return area;
     }
 
@@ -69,12 +144,8 @@
             return;
         }
 
-        Envelope relevantArea = M > CULL_POINT_THRESHOLD
-			? pathBoundingBox(path, 0.05d)
-			: null;
-        
         List<GridCell> cells = GridCell.pointsToGridCells(
-            points, relevantArea);
+            points, relevantArea(path, points));
 
         if (cells.isEmpty()) {
             log.warn("no cells found");
diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation3D.java
--- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation3D.java	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation3D.java	Mon Jan 25 18:25:19 2010 +0000
@@ -22,9 +22,6 @@
 {
     private static Logger log = Logger.getLogger(Interpolation3D.class);
 
-    public static final int CULL_POINT_THRESHOLD = Integer.getInteger(
-        "gnv.interpolation3d.cull.point.threshold", 1000);
-
     public static final int DEFAULT_WIDTH  = 1024;
     public static final int DEFAULT_HEIGHT =  768;
 
@@ -109,12 +106,8 @@
             return false;
         }
 
-        Envelope relevantArea = M > CULL_POINT_THRESHOLD
-			? Interpolation2D.pathBoundingBox(path, 0.05d)
-			: null;
-        
         List<GridCell> cells = GridCell.pointsToGridCells(
-            points, relevantArea);
+            points, Interpolation2D.relevantArea(path, points));
 
         if (cells.isEmpty()) {
             log.warn("no cells found");
diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileMeshCrossOutputState.java
--- a/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileMeshCrossOutputState.java	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileMeshCrossOutputState.java	Mon Jan 25 18:25:19 2010 +0000
@@ -1,6 +1,3 @@
-/**
- *
- */
 package de.intevation.gnv.state.profile.horizontal;
 
 import com.vividsolutions.jts.geom.Coordinate;
@@ -54,20 +51,20 @@
  * @author Ingo Weinzierl      (iweinzierl@intevation.de)
  * @author Sascha L. Teichmann (sascha.teichmann@intevation.de)
  */
-public class HorizontalProfileMeshCrossOutputState extends
-                                                  HorizontalProfileOutputState {
+public class HorizontalProfileMeshCrossOutputState
+extends      HorizontalProfileOutputState 
+{
 
-    /**
-     *
-     */
+    public static final boolean USE_INDEX_BUFFER =
+        Boolean.getBoolean("gnv.horizontal.profile.mesh.cross.index.buffer");
+
     private static final long serialVersionUID = 2205958041745637263L;
     
     /**
      * the logger, used to log exceptions and additonaly information
      */
-    private static Logger log = 
-                  Logger.getLogger(HorizontalProfileMeshCrossOutputState.class);
-    
+    private static Logger log = Logger.getLogger(
+        HorizontalProfileMeshCrossOutputState.class);
 
     private String ijkQueryID = null;
    
@@ -178,12 +175,13 @@
                 }
 
                 try {
-                    String additionWhere =
-                        WKTUtils.worldCoordinatesToIndex(
+                    String additionWhere = USE_INDEX_BUFFER
+                        ? WKTUtils.worldCoordinatesToIndex(
                             coords,
                             result,
                             meshId.getValue(),
-                            ijkQueryID);
+                            ijkQueryID)
+                        : WKTUtils.TRUE_EXPRESSION;
 
                     String[] addedFilterValues = StringUtils.append(
                         generateFilterValuesFromInputData(),
diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontalcrosssection/HorizontalCrossSectionMeshOutputState.java
--- a/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontalcrosssection/HorizontalCrossSectionMeshOutputState.java	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontalcrosssection/HorizontalCrossSectionMeshOutputState.java	Mon Jan 25 18:25:19 2010 +0000
@@ -85,6 +85,9 @@
      */
     private static final long serialVersionUID = 3233620652465061860L;
     
+    public static final boolean USE_INDEX_BUFFER =
+        Boolean.getBoolean("gnv.horizontal.cross.section.mesh.index.buffer");
+
     private String ijkQueryID;
 
     private Boolean shapeFileLock = new Boolean(true);
@@ -371,17 +374,24 @@
         try {
             Envelope env  = p.getEnvelopeInternal();
 
-            Coordinate [] coords = new Coordinate [] {
-                new Coordinate(env.getMinX(), env.getMinY()),
-                new Coordinate(env.getMinX(), env.getMaxY()),
-                new Coordinate(env.getMaxX(), env.getMaxY()),
-                new Coordinate(env.getMaxX(), env.getMinY()) };
+            String additionWhere;
 
-            String additionWhere =
-                WKTUtils.worldEnvelopeCoordinatesToIndex(
-                    coords,
-                    meshId.getValue(),
-                    ijkQueryID);
+            if (USE_INDEX_BUFFER) {
+                Coordinate [] coords = new Coordinate [] {
+                    new Coordinate(env.getMinX(), env.getMinY()),
+                    new Coordinate(env.getMinX(), env.getMaxY()),
+                    new Coordinate(env.getMaxX(), env.getMaxY()),
+                    new Coordinate(env.getMaxX(), env.getMinY()) };
+
+                additionWhere =
+                    WKTUtils.worldEnvelopeCoordinatesToIndex(
+                        coords,
+                        meshId.getValue(),
+                        ijkQueryID);
+            }
+            else {
+                additionWhere = WKTUtils.TRUE_EXPRESSION;
+            }
 
             String[] addedFilterValues = StringUtils.append(
                 generateFilterValuesFromInputData(),
diff -r f598702b2a10 -r 44415ae01ddb gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/verticalcrosssection/VerticalCrossSectionOutputState.java
--- a/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/verticalcrosssection/VerticalCrossSectionOutputState.java	Mon Jan 25 17:06:01 2010 +0000
+++ b/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/verticalcrosssection/VerticalCrossSectionOutputState.java	Mon Jan 25 18:25:19 2010 +0000
@@ -93,6 +93,9 @@
 
     public static final Integer GROUND_FILL_INDEX = Integer.valueOf(-2);
 
+    public static final boolean USE_INDEX_BUFFER =
+        Boolean.getBoolean("gnv.vertical.cross.section.index.buffer");
+
     public static final String[] ATTRIBUTE_LIST = {
         "SHAPE",
         "Z",
@@ -191,11 +194,13 @@
         }
 
         try {
-            String additionWhere = WKTUtils.worldCoordinatesToIndex(
-                coords,
-                result,
-                meshId.getValue(),
-                ijkQueryID);
+            String additionWhere = USE_INDEX_BUFFER
+                ? WKTUtils.worldCoordinatesToIndex(
+                    coords,
+                    result,
+                    meshId.getValue(),
+                    ijkQueryID)
+                : WKTUtils.TRUE_EXPRESSION;
 
             String[] addedFilterValues = StringUtils.append(
                 generateFilterValuesFromInputData(),