# HG changeset patch # User Sascha L. Teichmann # Date 1263465245 0 # Node ID ab29e4ff2fdae10041c28ed51e9835466c869289 # Parent a6a33ef35809cf4a85539a4fe1318a61fcff5ac6 Added area interpolation needed for "Horizontalschnitt" gnv-artifacts/trunk@540 c6561f87-3c4e-4783-a992-168aeb5c3f6f diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/ChangeLog --- a/gnv-artifacts/ChangeLog Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/ChangeLog Thu Jan 14 10:34:05 2010 +0000 @@ -1,3 +1,28 @@ +2010-01-13 Sascha L. Teichmann + + * src/test/java/de/intevation/gnv/artifacts/util/ShapeFileWriterTestCase.java, + src/main/java/de/intevation/gnv/state/profile/verticalcrosssection/VerticalCrossSectionOutputState.java: + Removed needless imports. + + * src/main/java/de/intevation/gnv/math/AreaInterpolation.java: New. Interpolates + area for a given bounding box, taking gaps and DEM into account. + Not very fast. Use bilinear interpolation to match the "Profilschnitt". + Possible TODOs: + - speed up by assuming the grid is not sparse. + - use higher interpolation methods. + + * src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContext.java, + src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContextFactory.java: + Added configuration for ground interpolation. + + * src/main/java/de/intevation/gnv/state/profile/horizontalcrosssection/HorizontalCrossSectionMeshOutputState.java: + Add helper functions to access configuration. + + * src/main/java/de/intevation/gnv/math/Interpolation2D.java: Simplified Code. + + * src/main/java/de/intevation/gnv/math/L1Comparator.java: add setReference() + method. + 2010-01-13 Ingo Weinzierl * src/main/java/de/intevation/gnv/artifacts/GNVArtifactBase.java: 'advance' diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContext.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContext.java Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContext.java Thu Jan 14 10:34:05 2010 +0000 @@ -46,7 +46,13 @@ DEFAULT_HORIZONTAL_CROSS_SECTION_PROFILE_SAMPLES = Integer.valueOf(250); public static final Integer - DEFAULT_HORIZONTAL_CROSS_SECTION_SAMPLES = Integer.valueOf(250); + DEFAULT_HORIZONTAL_CROSS_SECTION_SAMPLES = Integer.valueOf(1024); + + public static final String HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION_KEY = + "gnv.horizontal.cross.section.ground.interpolation"; + + public static final String DEFAULT_HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION = + "bilinear"; public static final String PALETTES_KEY = "gnv.color.palettes"; diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContextFactory.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContextFactory.java Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/artifacts/context/GNVArtifactContextFactory.java Thu Jan 14 10:34:05 2010 +0000 @@ -72,6 +72,9 @@ public final static String HORIZONTAL_CROSS_SECTION_PROFILE_SAMPLES = "/artifact-database/gnv/horizontal-cross-section-profile/samples/@number"; + public final static String HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION = + "/artifact-database/gnv/horizontal-cross-section/ground/@interpolation"; + public final static String HORIZONTAL_CROSS_SECTION_SAMPLES = "/artifact-database/gnv/horizontal-cross-section/samples/@number"; @@ -323,6 +326,33 @@ configureHorizontalCrossSectionSamples(config, context); configureHorizontalCrossSectionResultShapeFilePath(config, context); + configureHorizontalCrossSectionGroundInterpolation(config, context); + } + + protected void configureHorizontalCrossSectionGroundInterpolation( + Document config, + GNVArtifactContext context + ) + { + log.info( + "configuration of horizontal cross section ground interpolation"); + + String interpolation = Config.getStringXPath( + config, + HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION); + + if (interpolation == null + || (interpolation = interpolation.trim()).length() == 0) { + interpolation = GNVArtifactContext + .DEFAULT_HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION; + } + + log.info("ground interpolation: " + interpolation); + + context.put( + GNVArtifactContext + .HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION_KEY, + interpolation); } protected void configureHorizontalCrossSectionResultShapeFilePath( diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java Thu Jan 14 10:34:05 2010 +0000 @@ -0,0 +1,156 @@ +package de.intevation.gnv.math; + +import java.io.Serializable; + +import java.util.ArrayList; +import java.util.Arrays; +import java.util.List; +import java.util.HashMap; +import java.util.Collections; + +import java.awt.Dimension; + +import com.vividsolutions.jts.geom.Coordinate; +import com.vividsolutions.jts.geom.Envelope; + +import com.vividsolutions.jts.index.quadtree.Quadtree; + +import org.apache.log4j.Logger; + +/** + * @author Sascha L. Teichmann (sascha.teichmann@intevation.de) + */ +public class AreaInterpolation +implements Serializable +{ + private static Logger log = Logger.getLogger(AreaInterpolation.class); + + protected double [] raster; + + protected int width; + protected int height; + + public AreaInterpolation() { + } + + public int getWidth() { + return width; + } + + public int getHeight() { + return height; + } + + public double [] getRaster() { + return raster; + } + + public boolean interpolate( + List points, + Envelope boundingBox, + Dimension samples, + XYDepth depth + ) { + if (points == null || points.isEmpty()) { + return false; + } + + double [] buffer = Interpolation2D.calculateBuffer(points); + double dxMax = buffer[0]; + double dyMax = buffer[1]; + + Quadtree spatialIndex = new Quadtree(); + + for (int i = points.size()-1; i >= 0; --i) { + Point2d p = points.get(i); + spatialIndex.insert(p.envelope(), p); + } + + int W = samples.width; + int H = samples.height; + + double cellWidth = boundingBox.getWidth() / W; + double cellHeight = boundingBox.getHeight() / H; + + Envelope queryBuffer = new Envelope(); + Point2d [] neighbors = new Point2d[4]; + Coordinate center = new Coordinate(); + L1Comparator invL1 = new L1Comparator(center); + + double [] raster = new double[W*H]; + Arrays.fill(raster, Double.NaN); + + int row = 0; + for (int j = 0; j < H; ++j, row += W) { + double y = j*cellHeight + 0.5d*cellHeight; + double x = 0.5d*cellWidth; + for (int i = 0; i < W; ++i, x += cellWidth) { + queryBuffer.init( + x - dxMax, x + dxMax, + y - dyMax, y + dyMax); + + List potential = spatialIndex.query(queryBuffer); + + if (potential.isEmpty()) { + continue; + } + + center.x = x; center.y = y; + Collections.sort(potential, invL1); + + neighbors[0] = neighbors[1] = + neighbors[2] = neighbors[3] = null; + + int mask = 0; + // reversed order is essential here! + for (int k = potential.size()-1; k >= 0; --k) { + Point2d n = (Point2d)potential.get(k); + int code = n.x > center.x ? 1 : 0; + if (n.y > center.y) code |= 2; + neighbors[code] = n; + mask |= 1 << code; + } + + int numNeighbors = Integer.bitCount(mask); + + // only interpolate if we have all four neighbors + // we do not have any gaps and we are below surface. + if (numNeighbors == 4 + && !neighbors[0].hasIGap(neighbors[1]) + && !neighbors[1].hasJGap(neighbors[3]) + && !neighbors[3].hasIGap(neighbors[2]) + && !neighbors[2].hasJGap(neighbors[0]) + && depth.depth(center) <= 0d + ) { + double z1 = Interpolation2D.interpolate( + neighbors[0].x, neighbors[0].z, + neighbors[1].x, neighbors[1].z, + center.x); + double z2 = Interpolation2D.interpolate( + neighbors[2].x, neighbors[2].z, + neighbors[3].x, neighbors[3].z, + center.x); + double y1 = Interpolation2D.interpolate( + neighbors[0].x, neighbors[0].y, + neighbors[1].x, neighbors[1].y, + center.x); + double y2 = Interpolation2D.interpolate( + neighbors[2].x, neighbors[2].y, + neighbors[3].x, neighbors[3].y, + center.x); + raster[row+i] = Interpolation2D.interpolate( + y1, z1, + y2, z2, + center.y); + } + } + } + + this.raster = raster; + this.width = W; + this.height = H; + + return true; + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 : diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java Thu Jan 14 10:34:05 2010 +0000 @@ -135,6 +135,8 @@ int missedInterpolations = 0; int interpolations = 0; + L1Comparator invL1 = new L1Comparator(center); + for (double p = from; p <= to; p += dP) { if (!linearToMap.locate(p, center)) { continue; @@ -145,7 +147,6 @@ List potential = spatialIndex.query(queryBuffer); - L1Comparator invL1 = new L1Comparator(center); Collections.sort(potential, invL1); neighbors[0] = neighbors[1] = diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/src/main/java/de/intevation/gnv/math/L1Comparator.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/L1Comparator.java Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/L1Comparator.java Thu Jan 14 10:34:05 2010 +0000 @@ -12,10 +12,17 @@ { private Coordinate ref; + public L1Comparator() { + } + public L1Comparator(Coordinate ref) { this.ref = ref; } + public void setReference(Coordinate ref) { + this.ref = ref; + } + public int compare(Object a, Object b) { Coordinate pa = (Coordinate)a; Coordinate pb = (Coordinate)b; diff -r a6a33ef35809 -r ab29e4ff2fda 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 Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontalcrosssection/HorizontalCrossSectionMeshOutputState.java Thu Jan 14 10:34:05 2010 +0000 @@ -3,10 +3,12 @@ */ package de.intevation.gnv.state.profile.horizontalcrosssection; +import java.io.File; import java.io.IOException; import java.io.OutputStream; import java.io.OutputStreamWriter; import java.io.UnsupportedEncodingException; + import java.util.Arrays; import java.util.Collection; import java.util.Iterator; @@ -16,42 +18,53 @@ import org.apache.log4j.Logger; import org.jfree.chart.ChartTheme; + import org.w3c.dom.Node; import com.vividsolutions.jts.geom.Coordinate; import com.vividsolutions.jts.geom.Envelope; import com.vividsolutions.jts.geom.Polygon; -import com.vividsolutions.jts.io.ParseException; -import com.vividsolutions.jts.io.WKTReader; import au.com.bytecode.opencsv.CSVWriter; + import de.intevation.gnv.artifacts.cache.CacheFactory; import de.intevation.gnv.artifacts.context.GNVArtifactContext; + import de.intevation.gnv.chart.Chart; import de.intevation.gnv.chart.ChartLabels; + import de.intevation.gnv.geobackend.base.Result; + import de.intevation.gnv.geobackend.base.query.QueryExecutor; import de.intevation.gnv.geobackend.base.query.QueryExecutorFactory; import de.intevation.gnv.geobackend.base.query.exception.QueryException; + import de.intevation.gnv.state.InputData; + import de.intevation.gnv.state.exception.StateException; + import de.intevation.gnv.state.timeseries.TimeSeriesOutputState; + import de.intevation.gnv.statistics.Statistics; + import de.intevation.gnv.utils.StringUtils; import de.intevation.gnv.utils.WKTUtils; import de.intevation.artifactdatabase.Config; + import de.intevation.artifacts.CallContext; +import de.intevation.gnv.geobackend.sde.datasources.RasterObject; + /** - * @author Tim Englich - * + * @author Tim Englich (tim.englich@intevation.de) + * @author Sascha L. Teichmann (sascha.teichmann@intevation.de) */ public class HorizontalCrossSectionMeshOutputState - extends TimeSeriesOutputState { - +extends TimeSeriesOutputState +{ private static Logger log = Logger - .getLogger(HorizontalCrossSectionMeshOutputState.class); + .getLogger(HorizontalCrossSectionMeshOutputState.class); /** * The UID of this Class @@ -95,23 +108,6 @@ log.info("Chart not in cache yet."); log.warn("This sort of chart is not implemented yet."); - /* TODO Implement a special chart for this sort of charts. - chart = new HorizontalProfileChart( - chartLables, - chartTheme, - parameters, - measurements, - result, - dates, - locale - ); - chart.generateChart(); - - if (CACHE_CHART) { - log.info("Put chart into cache."); - purifyChart(chart, uuid); - } - */ InputData meshPolygon = inputData.get("mesh_polygon"); String meshPolygonWkt = null; @@ -226,11 +222,12 @@ * java.util.Collection) */ @Override - protected void createCSV(OutputStream outputStream, - Collection chartResult) - throws UnsupportedEncodingException, - IOException, - StateException { + protected void createCSV( + OutputStream outputStream, + Collection chartResult + ) throws UnsupportedEncodingException, IOException, StateException + { + /* if (chartResult != null) { try { CSVWriter writer = new CSVWriter(new OutputStreamWriter( @@ -260,6 +257,7 @@ throw new StateException( "No Data given for generating an CSV-File."); } + */ } /** @@ -282,4 +280,24 @@ : GNVArtifactContext.DEFAULT_HORIZONTAL_CROSS_SECTION_SAMPLES; } + private static File shapefileDirectory(CallContext callContext) { + GNVArtifactContext context = + (GNVArtifactContext)callContext.globalContext(); + File dir = (File)context.get( + GNVArtifactContext.HORIZONTAL_CROSS_SECTION_RESULT_SHAPEFILE_PATH_KEY); + return dir != null + ? dir + : GNVArtifactContext.DEFAULT_HORIZONTAL_CROSS_SECTION_PROFILE_SHAPEFILE_PATH; + } + + private static int getGroundInterpolation(CallContext callContext) { + GNVArtifactContext context = + (GNVArtifactContext)callContext.globalContext(); + + String interpolation = (String)context.get( + GNVArtifactContext.HORIZONTAL_CROSS_SECTION_GROUND_INTERPOLATION_KEY); + + return RasterObject.getInterpolationType(interpolation); + } } +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 : diff -r a6a33ef35809 -r ab29e4ff2fda 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 Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/verticalcrosssection/VerticalCrossSectionOutputState.java Thu Jan 14 10:34:05 2010 +0000 @@ -57,7 +57,6 @@ import de.intevation.gnv.math.IJKey; import de.intevation.gnv.math.LinearMetrics; import de.intevation.gnv.math.Interpolation3D; -import de.intevation.gnv.math.ConstantXYDepth; import de.intevation.gnv.state.InputData; diff -r a6a33ef35809 -r ab29e4ff2fda gnv-artifacts/src/test/java/de/intevation/gnv/artifacts/util/ShapeFileWriterTestCase.java --- a/gnv-artifacts/src/test/java/de/intevation/gnv/artifacts/util/ShapeFileWriterTestCase.java Wed Jan 13 23:10:56 2010 +0000 +++ b/gnv-artifacts/src/test/java/de/intevation/gnv/artifacts/util/ShapeFileWriterTestCase.java Thu Jan 14 10:34:05 2010 +0000 @@ -13,7 +13,6 @@ import org.apache.log4j.Logger; import com.vividsolutions.jts.geom.MultiLineString; -import com.vividsolutions.jts.io.ParseException; import com.vividsolutions.jts.io.WKTReader; import de.intevation.gnv.utils.Pair;