comparison gnv-artifacts/src/main/java/de/intevation/gnv/math/Interpolation2D.java @ 365:f66088a43ecc

Added horizontal crossprofile charts to chart pallet. Fixed some bugs before interpolation. gnv-artifacts/trunk@440 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author Ingo Weinzierl <ingo.weinzierl@intevation.de>
date Wed, 16 Dec 2009 19:29:05 +0000
parents aec85d00d82c
children 086e3af38b96
comparison
equal deleted inserted replaced
364:2413273f1c13 365:f66088a43ecc
1 package de.intevation.gnv.math; 1 package de.intevation.gnv.math;
2 2
3 import java.util.ArrayList;
3 import java.util.List; 4 import java.util.List;
5 import java.util.HashMap;
4 import java.util.Collections; 6 import java.util.Collections;
5 7
6 import com.vividsolutions.jts.geom.Coordinate; 8 import com.vividsolutions.jts.geom.Coordinate;
7 import com.vividsolutions.jts.geom.Envelope; 9 import com.vividsolutions.jts.geom.Envelope;
8 10
9 import com.vividsolutions.jts.index.quadtree.Quadtree; 11 import com.vividsolutions.jts.index.quadtree.Quadtree;
12
13 import org.apache.log4j.Logger;
10 14
11 /** 15 /**
12 * @author Sascha L. Teichmann 16 * @author Sascha L. Teichmann
13 */ 17 */
14 public final class Interpolation2D 18 public final class Interpolation2D
15 { 19 {
20 private static Logger log = Logger.getLogger(Interpolation2D.class);
21
16 public interface Consumer { 22 public interface Consumer {
17 void interpolated(Coordinate point); 23 void interpolated(Coordinate point);
18 } // interface Consumer 24 } // interface Consumer
19 25
20 private Interpolation2D() { 26 private Interpolation2D() {
30 Consumer consumer 36 Consumer consumer
31 ) { 37 ) {
32 int N = path.size(); 38 int N = path.size();
33 int M = points.size(); 39 int M = points.size();
34 40
41 log.debug("Size of path: " + N);
42 log.debug("Size of points: " + M);
43
35 if (M < 1 || N < 2) { // nothing to do 44 if (M < 1 || N < 2) { // nothing to do
36 return; 45 return;
37 } 46 }
38 // figure out max delta(p[i].x, p[i-1].x) 47
39 Collections.sort(points, Point2d.X_COMPARATOR); 48 HashMap<Integer, ArrayList<Point2d>> map = new HashMap<Integer, ArrayList<Point2d>>();
49
50 for (int k = M-1; k >= 0; --k) {
51 Point2d p = points.get(k);
52
53 ArrayList<Point2d> list = map.get(p.j);
54
55 if (list == null) {
56 map.put(p.j, list = new ArrayList<Point2d>());
57 }
58 list.add(p);
59 }
60
40 double dxMax = -Double.MAX_VALUE; 61 double dxMax = -Double.MAX_VALUE;
41 for (int i = 1; i < M; ++i) { 62
42 double dx = Math.abs(path.get(i).x - path.get(i-1).x); 63 for (ArrayList<Point2d> v: map.values()) {
43 if (dx > dxMax) { 64 Collections.sort(v, Point2d.X_COMPARATOR);
44 dxMax = dx; 65 for (int i = 1, L = v.size(); i < L; ++i) {
45 } 66 double dx = Math.abs(v.get(i).x - v.get(i-1).x);
46 } 67 if (dx > dxMax) {
47 68 dxMax = dx;
48 dxMax = dxMax*0.5d + 1e-5d; 69 }
49 70 }
50 // figure out max delta(p[i].y, p[i-1].y) 71 }
51 Collections.sort(path, Point2d.X_COMPARATOR); 72
73 dxMax = dxMax + 1e-5d;
74
75 map.clear();
76
77 for (int k = M-1; k >= 0; --k) {
78 Point2d p = points.get(k);
79
80 ArrayList<Point2d> list = map.get(p.i);
81
82 if (list == null) {
83 map.put(p.i, list = new ArrayList<Point2d>());
84 }
85 list.add(p);
86 }
87
52 double dyMax = -Double.MAX_VALUE; 88 double dyMax = -Double.MAX_VALUE;
53 for (int i = 1; i < M; ++i) { 89
54 double dy = Math.abs(path.get(i).y - path.get(i-1).y); 90 for (ArrayList<Point2d> v: map.values()) {
55 if (dy > dyMax) { 91 Collections.sort(v, Point2d.Y_COMPARATOR);
56 dyMax = dy; 92 for (int i = 1, L = v.size(); i < L; ++i) {
57 } 93 double dy = Math.abs(v.get(i).y - v.get(i-1).y);
58 } 94 if (dy > dyMax) {
59 95 dyMax = dy;
60 dyMax = dyMax*0.5d + 1e-5d; 96 }
97 }
98 }
99
100 dyMax = dyMax + 1e-5d;
101
102 map = null;
103
104 log.debug("buffer size: " + dxMax + " / " + dyMax);
61 105
62 // put into spatial index to speed up finding neighbors. 106 // put into spatial index to speed up finding neighbors.
63 Quadtree spatialIndex = new Quadtree(); 107 Quadtree spatialIndex = new Quadtree();
64 108
65 for (int i = 0; i < M; ++i) { 109 for (int i = 0; i < M; ++i) {
76 120
77 Envelope queryBuffer = new Envelope(); 121 Envelope queryBuffer = new Envelope();
78 122
79 Point2d [] neighbors = new Point2d[4]; 123 Point2d [] neighbors = new Point2d[4];
80 124
81 for (double p = to; p <= from; p += dP) { 125 int missedInterpolations = 0;
126 int interpolations = 0;
127
128 for (double p = from; p <= to; p += dP) {
82 if (!linearToMap.locate(p, center)) { 129 if (!linearToMap.locate(p, center)) {
83 continue; 130 continue;
84 } 131 }
85 queryBuffer.init( 132 queryBuffer.init(
86 center.x - dxMax, center.x + dxMax, 133 center.x - dxMax, center.x + dxMax,
139 center.z = interpolate( 186 center.z = interpolate(
140 y1, z1, 187 y1, z1,
141 y2, z2, 188 y2, z2,
142 center.y); 189 center.y);
143 consumer.interpolated(center); 190 consumer.interpolated(center);
144 } 191 ++interpolations;
145 } 192 }
193 else {
194 ++missedInterpolations;
195 }
196 }
197
198 log.debug("interpolations: " + interpolations + " / " + missedInterpolations);
146 } 199 }
147 200
148 public static final double interpolate( 201 public static final double interpolate(
149 double x1, double y1, 202 double x1, double y1,
150 double x2, double y2, 203 double x2, double y2,

http://dive4elements.wald.intevation.org