# HG changeset patch
# User Sascha L. Teichmann <sascha.teichmann@intevation.de>
# Date 1263053530 0
# Node ID ccd976fc0f7be0c509b611916158834d8fbce22b
# Parent  23d5cc37dd5b62acee64c873d4a362b23b46d42c
Implemented bilinear interpolation on raster tiles.


geo-backend/trunk@520 c6561f87-3c4e-4783-a992-168aeb5c3f6f

diff -r 23d5cc37dd5b -r ccd976fc0f7b geo-backend/ChangeLog
--- a/geo-backend/ChangeLog	Sat Jan 09 12:41:26 2010 +0000
+++ b/geo-backend/ChangeLog	Sat Jan 09 16:12:10 2010 +0000
@@ -1,3 +1,8 @@
+2009-01-09	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
+
+	* src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java:
+	  Implemented bilinear interpolation.
+
 2009-01-09	Sascha L. Teichmann	<sascha.teichmann@intevation.de>
 
 	* src/main/java/de/intevation/gnv/geobackend/sde/datasources/ArcSDEStatement.java:
diff -r 23d5cc37dd5b -r ccd976fc0f7b geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java
--- a/geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java	Sat Jan 09 12:41:26 2010 +0000
+++ b/geo-backend/src/main/java/de/intevation/gnv/geobackend/sde/datasources/RasterObject.java	Sat Jan 09 16:12:10 2010 +0000
@@ -89,6 +89,71 @@
         return getValue(coordinate, NEAREST_NEIGHBOR);
     }
 
+    public final double interpolateBilinear(double px, double py) {
+
+        if (px < 0.5d) { // left border
+            if (py < 0.5d) { // upper left edge
+                return rasterData[0];
+            }
+            if (py > tileHeight-0.5d) { // lower left edge
+                return rasterData[rasterData.length-tileWidth];
+            }
+            // center left
+            int    i = (int)(py -= 0.5d);
+            double t = py - i;
+            i *= tileWidth;
+            return (1d-t)*rasterData[i] + t*rasterData[i+tileWidth];
+        }
+
+        if (px > tileWidth-0.5d) { // right border
+            if (py < 0.5d) { // upper right edge
+                return rasterData[tileWidth-1];
+            }
+            if (py > tileHeight-0.5d) { // lower right edge
+                return rasterData[rasterData.length-1];
+            }
+            // center right
+            int    i = (int)(py -= 0.5d);
+            double t = py - i;
+            i = i*tileWidth + tileWidth-1;
+            return (1d-t)*rasterData[i] + t*rasterData[i+tileWidth];
+        }
+
+        if (py < 0.5d) { // top border
+            int    i = (int)(px -= 0.5d);
+            double t = px - i;
+            return (1d-t)*rasterData[i] + t*rasterData[i+1];
+        }
+
+        if (py > tileHeight-0.5d) { // bottom border
+            int    i = (int)(px -= 0.5d);
+            double t = px - i;
+            i += rasterData.length-tileWidth;
+            return (1d-t)*rasterData[i] + t*rasterData[i+1];
+        }
+
+        // center
+        int    i  = (int)(py -= 0.5d);
+        int    j  = (int)(px -= 0.5d);
+        double ti = py - i;
+        double tj = px - j;
+
+        int idx = i*tileWidth + j;
+
+        double v1j1 = rasterData[idx];
+        double v2j1 = rasterData[idx+1];
+
+        idx += tileWidth;
+
+        double v1j2 = rasterData[idx];
+        double v2j2 = rasterData[idx+1];
+
+        double v1 = (1d-tj)*v1j1 + tj*v2j1;
+        double v2 = (1d-tj)*v1j2 + tj*v2j2;
+
+        return (1d-ti)*v1 + ti*v2;
+    }
+
     public double getValue(Coordinate coordinate, int interpolationType) {
 
         double px = mx*coordinate.x + bx;
@@ -98,11 +163,9 @@
             return Double.NaN;
         }
 
-        /*
         if (interpolationType == BILINEAR) {
-            // TODO: implement me!
+            return interpolateBilinear(px, py);
         }
-        */
 
         int posX = Math.min(tileWidth-1,  (int)Math.round(px));
         int posY = Math.min(tileHeight-1, (int)Math.round(py));