# HG changeset patch # User Sascha L. Teichmann # Date 1264184531 0 # Node ID d9d933e06875715f9060c1d4d351fc282de3bcd4 # Parent ca5048e4e515413c9dbbcc36421287a8ab84e86e Fixed gnv/issue153 gnv-artifacts/trunk@608 c6561f87-3c4e-4783-a992-168aeb5c3f6f diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/ChangeLog --- a/gnv-artifacts/ChangeLog Fri Jan 22 15:45:59 2010 +0000 +++ b/gnv-artifacts/ChangeLog Fri Jan 22 18:22:11 2010 +0000 @@ -1,3 +1,33 @@ +2010-01-22 Sascha L. Teichmann + + * src/main/java/de/intevation/gnv/math/GridCell.java: New. + A 4-tupel of neighbored points in the mesh. It is valid + to interpolate in this area. + + * src/main/java/de/intevation/gnv/math/AreaInterpolation.java: + The algorithm how neighbored points in the mesh are determined + has changed. Now all incoming points are tiled into GridCells. + If there are gaps in i,j the corresponding tile is omited. + These tiles are stored in an R tree. To lookup a point in + world coordinates the spatial index is queried. If no result + is found the point is ignore as a gap. If a fitting grid cell + is found the interpolation in done between the four points + of that cell is performed. Special gap checking is not needed any + longer. This fixes gnv/issue153 because there are no assumptions + about axis aligned points any more. + + * src/main/java/de/intevation/gnv/math/Interpolation2D.java: Used + euclid distance to estimate spatial buffer size. TODO: Remove + this code when adjusting the "Profilschnitte" to the same logic + as used in "Horzontalschnitte" now. + + * src/main/java/de/intevation/gnv/math/LinearFunction.java: Added + author's email. + + * src/main/java/de/intevation/gnv/chart/TimeSeriesChart.java, + src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileOutputState.java: + Cleanup imports. + 2010-01-22 Ingo Weinzierl * src/main/java/de/intevation/gnv/artifacts/GNVArtifactBase.java: Removed @@ -27,6 +57,7 @@ Added the possibility to configure the different width of the Meshes. 2010-01-22 Tim Englich + * src/test/java/de/intevation/gnv/artifacts/TestArtifactDatabase.java (serviceNamesAndDescriptions),(process): Fixed Compiler-Error after adding further Methods to the Interface ArtifactDatabase. diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/src/main/java/de/intevation/gnv/chart/TimeSeriesChart.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/chart/TimeSeriesChart.java Fri Jan 22 15:45:59 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/chart/TimeSeriesChart.java Fri Jan 22 18:22:11 2010 +0000 @@ -1,7 +1,13 @@ package de.intevation.gnv.chart; -import java.text.DateFormat; -import java.text.SimpleDateFormat; +import de.intevation.gnv.artifacts.ressource.RessourceFactory; + +import de.intevation.gnv.geobackend.base.Result; + +import de.intevation.gnv.state.describedata.KeyValueDescibeData; + +import de.intevation.gnv.timeseries.gap.TimeGap; + import java.util.Collection; import java.util.Date; import java.util.HashMap; @@ -11,25 +17,23 @@ import org.apache.log4j.Logger; +import org.jfree.chart.ChartFactory; import org.jfree.chart.ChartTheme; -import org.jfree.chart.ChartFactory; + import org.jfree.chart.axis.Axis; import org.jfree.chart.axis.DateAxis; -import org.jfree.chart.plot.XYPlot; + import org.jfree.chart.plot.PlotOrientation; +import org.jfree.chart.plot.XYPlot; + import org.jfree.data.general.Series; + +import org.jfree.data.time.Minute; import org.jfree.data.time.TimeSeries; -import org.jfree.data.time.Minute; import org.jfree.data.time.TimeSeriesCollection; -import de.intevation.gnv.artifacts.ressource.RessourceFactory; -import de.intevation.gnv.geobackend.base.Result; -import de.intevation.gnv.state.describedata.KeyValueDescibeData; -import de.intevation.gnv.timeseries.gap.TimeGap; - - /** - * @author Ingo Weinzierl + * @author Ingo Weinzierl (ingo.weinzierl@intevation.de) */ public class TimeSeriesChart extends AbstractXYLineChart diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java Fri Jan 22 15:45:59 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/AreaInterpolation.java Fri Jan 22 18:22:11 2010 +0000 @@ -3,15 +3,13 @@ import com.vividsolutions.jts.geom.Coordinate; import com.vividsolutions.jts.geom.Envelope; -import com.vividsolutions.jts.index.quadtree.Quadtree; +import com.vividsolutions.jts.index.strtree.STRtree; import java.awt.Dimension; import java.io.Serializable; -import java.util.ArrayList; import java.util.Arrays; -import java.util.Collections; import java.util.List; import org.apache.log4j.Logger; @@ -44,6 +42,9 @@ return raster; } + public static final double EPS = 1e-6d; + + public boolean interpolate( List points, Envelope boundingBox, @@ -53,9 +54,14 @@ boolean debug = log.isDebugEnabled(); if (points == null || points.isEmpty()) { - if (debug) { - log.debug("no points to interpolate"); - } + log.warn("no points to interpolate"); + return false; + } + + List cells = GridCell.pointsToGridCells(points); + + if (cells.isEmpty()) { + log.warn("no cells to interpolate"); return false; } @@ -64,47 +70,25 @@ double cellWidth = boundingBox.getWidth() / W; double cellHeight = boundingBox.getHeight() / H; - Envelope gridEnvelope = new Envelope(boundingBox); - - gridEnvelope.expandBy(2d*cellWidth, 2d*cellWidth); - - ArrayListrelevantPoints = - new ArrayList(); - for (int i = points.size()-1; i >= 0; --i) { - if (gridEnvelope.contains(points.get(i))) { - relevantPoints.add(points.get(i)); - } - } + STRtree spatialIndex = new STRtree(); - double [] buffer = Interpolation2D - .calculateBufferRemoveOutlier(relevantPoints); - double dxMax = buffer[0]; - double dyMax = buffer[1]; - - Quadtree spatialIndex = new Quadtree(); - - for (int i = relevantPoints.size()-1; i >= 0; --i) { - Point2d p = relevantPoints.get(i); - spatialIndex.insert(p.envelope(), p); + for (GridCell cell: cells) { + spatialIndex.insert(cell.getEnvelope(), cell); } if (debug) { - log.debug("points: " + relevantPoints.size()); log.debug("width: " + boundingBox.getWidth()); log.debug("height: " + boundingBox.getHeight()); log.debug("sample width: " + W); log.debug("sample height: " + H); log.debug("cell width: " + cellWidth); log.debug("cell height: " + cellHeight); - log.debug("buffer x: " + dxMax); - log.debug("buffer y: " + dyMax); } - Envelope queryBuffer = new Envelope(); - Point2d [] neighbors = new Point2d[4]; - Coordinate center = new Coordinate(); - L1Comparator invL1 = new L1Comparator(center); + Envelope queryBuffer = new Envelope(); + Coordinate center = new Coordinate(); + GridCell.CellFinder finder = new GridCell.CellFinder(); double [] raster = new double[W*H]; Arrays.fill(raster, Double.NaN); @@ -119,64 +103,38 @@ double y = j*cellHeight + 0.5d*cellHeight + minY; double x = 0.5d*cellWidth + minX; for (int i = 0; i < W; ++i, x += cellWidth) { - queryBuffer.init( - x - dxMax, x + dxMax, - y - dyMax, y + dyMax); - List potential = spatialIndex.query(queryBuffer); + queryBuffer.init(x - EPS, x + EPS, y - EPS, y + EPS); + center.x = x; center.y = y; + finder.prepare(center); + spatialIndex.query(queryBuffer, finder); - if (potential.size() < 4) { + GridCell found = finder.found; + + if (found == null) { 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); - } + double z1 = Interpolation2D.interpolate( + found.p1.x, found.p1.z, + found.p2.x, found.p2.z, + center.x); + double z2 = Interpolation2D.interpolate( + found.p3.x, found.p3.z, + found.p4.x, found.p4.z, + center.x); + double y1 = Interpolation2D.interpolate( + found.p1.x, found.p1.y, + found.p2.x, found.p2.y, + center.x); + double y2 = Interpolation2D.interpolate( + found.p3.x, found.p3.y, + found.p4.x, found.p4.y, + center.x); + raster[row+i] = Interpolation2D.interpolate( + y1, z1, + y2, z2, + center.y); } } diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/src/main/java/de/intevation/gnv/math/GridCell.java --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/GridCell.java Fri Jan 22 18:22:11 2010 +0000 @@ -0,0 +1,134 @@ +package de.intevation.gnv.math; + +import com.vividsolutions.jts.geom.Coordinate; +import com.vividsolutions.jts.geom.Envelope; +import com.vividsolutions.jts.geom.Geometry; +import com.vividsolutions.jts.geom.GeometryFactory; +import com.vividsolutions.jts.geom.LinearRing; +import com.vividsolutions.jts.geom.Point; +import com.vividsolutions.jts.geom.Polygon; + +import com.vividsolutions.jts.index.ItemVisitor; + +import java.io.Serializable; + +import java.util.ArrayList; +import java.util.HashMap; +import java.util.List; + +/** + * @author Sascha L. Teichmann (sascha.teichmann@intevation.de) + */ +public class GridCell +implements Serializable +{ + public static final class CellFinder + implements ItemVisitor + { + public GridCell found; + + protected Point point; + + public CellFinder() { + } + + public void prepare(Coordinate center) { + found = null; + point = GEOMETRY_FACTORY.createPoint(center); + } + + public void visitItem(Object item) { + if (found == null) { + GridCell cell = (GridCell)item; + if (cell.contains(point)) { + found = cell; + } + } + } + } // class CellFinder + + public Point2d p1; + public Point2d p2; + public Point2d p3; + public Point2d p4; + + protected Polygon polygon; + + public static final GeometryFactory GEOMETRY_FACTORY + = new GeometryFactory(); + + public GridCell() { + } + + public GridCell(Point2d p1, Point2d p2, Point2d p3, Point2d p4) { + this.p1 = p1; + this.p2 = p2; + this.p3 = p3; + this.p4 = p4; + createPolygon(); + } + + protected void createPolygon() { + LinearRing shell = GEOMETRY_FACTORY + .createLinearRing(new Coordinate [] { p1, p2, p3, p4, p1 }); + polygon = GEOMETRY_FACTORY.createPolygon(shell, null); + } + + public Envelope getEnvelope() { + return polygon.getEnvelopeInternal(); + } + + public boolean contains(Geometry coord) { + return polygon.contains(coord); + } + + public static List pointsToGridCells( + List points + ) { + int minI = Integer.MAX_VALUE; + int maxI = Integer.MIN_VALUE; + int minJ = Integer.MAX_VALUE; + int maxJ = Integer.MIN_VALUE; + + HashMap> rows = + new HashMap>(); + + for (Point2d p: points) { + + if (p.i < minI) minI = p.i; + if (p.i > maxI) maxI = p.i; + if (p.j < minJ) minJ = p.j; + if (p.j > maxJ) maxJ = p.j; + + HashMap row = rows.get(p.i); + + if (row == null) { + rows.put(p.i, row = new HashMap()); + } + + row.put(p.j, p); + } + + ArrayList cells = new ArrayList(points.size()); + + for (int i = minI + 1; i <= maxI; ++i) { + HashMap row1 = rows.get(i-1); + HashMap row2 = rows.get(i); + if (row1 != null && row2 != null) { + for (int j = minJ + 1; j < maxJ; ++j) { + Point2d p1 = row1.get(j-1); + Point2d p2 = row1.get(j); + Point2d p3 = row2.get(j); + Point2d p4 = row2.get(j-1); + + if (p1 != null && p2 != null && p3 != null && p4 != null) { + cells.add(new GridCell(p1, p2, p3, p4)); + } + } // for all columns in row + } + } // for all rows + + return cells; + } +} +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 : diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java Fri Jan 22 15:45:59 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java Fri Jan 22 18:22:11 2010 +0000 @@ -72,7 +72,8 @@ for (ArrayList v: jMap.values()) { Collections.sort(v, Point2d.Y_COMPARATOR); for (int i = 1, L = v.size(); i < L; ++i) { - double dy = Math.abs(v.get(i).y - v.get(i-1).y); + // double dy = Math.abs(v.get(i).y - v.get(i-1).y); + double dy = v.get(i).distance(v.get(i-1)); deltas.add(dy); mean.increment(dy); sd .increment(dy); @@ -104,7 +105,8 @@ for (ArrayList v: iMap.values()) { Collections.sort(v, Point2d.X_COMPARATOR); for (int i = 1, L = v.size(); i < L; ++i) { - double dx = Math.abs(v.get(i).x - v.get(i-1).x); + //double dx = Math.abs(v.get(i).x - v.get(i-1).x); + double dx = v.get(i).distance(v.get(i-1)); deltas.add(dx); mean.increment(dx); sd .increment(dx); diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/src/main/java/de/intevation/gnv/math/LinearFunction.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/math/LinearFunction.java Fri Jan 22 15:45:59 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/math/LinearFunction.java Fri Jan 22 18:22:11 2010 +0000 @@ -7,7 +7,7 @@ import org.apache.commons.math.analysis.UnivariateRealFunction; /** - * @author Sascha L. Teichmann + * @author Sascha L. Teichmann (sascha.teichmann@intevation.de) */ public class LinearFunction implements ParametricRealFunction @@ -54,4 +54,4 @@ return new double [] { x, 1f }; } } -// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8: +// vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 : diff -r ca5048e4e515 -r d9d933e06875 gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileOutputState.java --- a/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileOutputState.java Fri Jan 22 15:45:59 2010 +0000 +++ b/gnv-artifacts/src/main/java/de/intevation/gnv/state/profile/horizontal/HorizontalProfileOutputState.java Fri Jan 22 18:22:11 2010 +0000 @@ -3,14 +3,41 @@ */ package de.intevation.gnv.state.profile.horizontal; -import com.vividsolutions.jts.geom.Point; +import com.vividsolutions.jts.io.ParseException; import com.vividsolutions.jts.io.WKTReader; -import com.vividsolutions.jts.io.ParseException; + +import de.intevation.artifacts.CallContext; + +import de.intevation.gnv.chart.Chart; +import de.intevation.gnv.chart.ChartLabels; +import de.intevation.gnv.chart.HorizontalProfileChart; + +import de.intevation.gnv.exports.DefaultExport; +import de.intevation.gnv.exports.DefaultProfile; + +import de.intevation.gnv.exports.Export.Profile; + +import de.intevation.gnv.exports.ShapeDataCollector; + +import de.intevation.gnv.geobackend.base.Result; + +import de.intevation.gnv.state.describedata.KeyValueDescibeData; + +import de.intevation.gnv.state.exception.StateException; + +import de.intevation.gnv.state.timeseries.TimeSeriesOutputState; + +import de.intevation.gnv.statistics.HorizontalProfileStatistics; +import de.intevation.gnv.statistics.Statistics; + +import de.intevation.gnv.utils.WKTUtils; import java.io.IOException; import java.io.OutputStream; import java.io.UnsupportedEncodingException; + import java.text.SimpleDateFormat; + import java.util.Collection; import java.util.Date; import java.util.Iterator; @@ -20,26 +47,8 @@ import org.jfree.chart.ChartTheme; -import de.intevation.artifacts.CallContext; - -import de.intevation.gnv.chart.Chart; -import de.intevation.gnv.chart.ChartLabels; -import de.intevation.gnv.chart.HorizontalProfileChart; -import de.intevation.gnv.exports.DefaultExport; -import de.intevation.gnv.exports.ShapeDataCollector; -import de.intevation.gnv.exports.DefaultProfile; -import de.intevation.gnv.exports.Export.Profile; -import de.intevation.gnv.geobackend.base.Result; -import de.intevation.gnv.state.describedata.KeyValueDescibeData; -import de.intevation.gnv.state.exception.StateException; -import de.intevation.gnv.state.timeseries.TimeSeriesOutputState; -import de.intevation.gnv.statistics.HorizontalProfileStatistics; -import de.intevation.gnv.statistics.Statistics; -import de.intevation.gnv.utils.WKTUtils; - /** - * @author Tim Englich - * + * @author Tim Englich (tim.englich@intevation.de) */ public class HorizontalProfileOutputState extends TimeSeriesOutputState