Mercurial > dive4elements > river

comparison flys-artifacts/src/main/java/de/intevation/flys/artifacts/geom/Polygon2D.java @ 3812:f788d2d901d6

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
merged flys-artifacts/pre2.6-2011-12-05
author Thomas Arendsen Hein <thomas@intevation.de>
date Fri, 28 Sep 2012 12:14:53 +0200
parents 6f83d9d434f2
children 5642a83420f2
comparison
equal deleted inserted replaced
3808:5fab0fe3c445 3812:f788d2d901d6
1 package de.intevation.flys.geom;
2
3 import java.io.Serializable;
4
5 import java.awt.Shape;
6
7 import java.awt.geom.Path2D;
8 import java.awt.geom.Point2D;
9
10 import java.util.ArrayList;
11 import java.util.List;
12 import java.util.Arrays;
13 import java.util.Comparator;
14 import java.util.Collections;
15
16 import de.intevation.flys.artifacts.math.Linear;
17
18 import static de.intevation.flys.geom.VectorUtils.X;
19 import static de.intevation.flys.geom.VectorUtils.Y;
20 import static de.intevation.flys.geom.VectorUtils.EPSILON;
21
22 public class Polygon2D
23 implements Serializable
24 {
25 public static final Comparator<Point2D> POINT_X_CMP =
26 new Comparator<Point2D>() {
27 public int compare(Point2D a, Point2D b) {
28 double d = X(a) - X(b);
29 if (d < 0d) return -1;
30 return d > 0d ? + 1 : 0;
31 }
32 };
33
34 public static final Comparator<Point2D []> FIRST_POINT_X =
35 new Comparator<Point2D []>() {
36 public int compare(Point2D [] a, Point2D [] b) {
37 if (a.length == 0) return -1; // should not happen.
38 if (b.length == 0) return +1; // should not happen.
39 double d = X(a[0]) - Y(b[0]);
40 if (d < 0d) return -1;
41 return d > 0d ? + 1 : 0;
42 }
43 };
44
45 protected List<Point2D> points;
46
47 public Polygon2D() {
48 points = new ArrayList<Point2D>();
49 }
50
51 public Polygon2D(List<Point2D> points) {
52 this.points = points;
53 }
54
55 public void add(double x, double y) {
56 points.add(new Point2D.Double(x, y));
57 }
58
59 protected static boolean addCheck(Point2D p, List<Point2D> points) {
60 switch (points.size()) {
61 case 0:
62 points.add(p);
63 return true;
64 case 1:
65 if (VectorUtils.epsilonEquals(points.get(0), p)) {
66 return false;
67 }
68 points.add(p);
69 return true;
70 default:
71 int L = points.size()-1;
72 Point2D last = points.get(L);
73 if (VectorUtils.epsilonEquals(last, p)) {
74 return false;
75 }
76 Point2D before = points.get(L-1);
77 if (VectorUtils.collinear(before, last, p)) {
78 points.set(L, p);
79 }
80 else {
81 points.add(p);
82 }
83 return true;
84 }
85 }
86
87 public boolean addCheck(Point2D p) {
88 return addCheck(p, points);
89 }
90
91 public void addReversed(List<Point2D> other) {
92 for (int i = other.size()-1; i >= 0; --i) {
93 addCheck(other.get(i));
94 }
95 }
96
97 public double area() {
98 double area = 0d;
99
100 for (int i = 0, N = points.size(); i < N; ++i) {
101 int j = (i + 1) % N;
102 Point2D pi = points.get(i);
103 Point2D pj = points.get(j);
104 area += X(pi)*Y(pj);
105 area -= X(pj)*Y(pi);
106 }
107
108 return 0.5d*area;
109 }
110
111 public Shape toShape() {
112 Path2D.Double path = new Path2D.Double();
113
114 int N = points.size();
115
116 if (N > 0) {
117 Point2D p = points.get(0);
118 path.moveTo(X(p), Y(p));
119 for (int i = 1; i < N; ++i) {
120 p = points.get(i);
121 path.lineTo(X(p), Y(p));
122 }
123 path.closePath();
124 }
125
126 return path;
127 }
128
129 protected static List<Point2D []> splitByNaNs(
130 double [] xs,
131 double [] ys
132 ) {
133 List<Point2D []> parts = new ArrayList<Point2D []>();
134
135 List<Point2D> tmp = new ArrayList<Point2D>(xs.length);
136
137 for (int i = 0; i < xs.length; ++i) {
138 double x = xs[i];
139 double y = ys[i];
140
141 if (Double.isNaN(x) || Double.isNaN(y)) {
142 if (!tmp.isEmpty()) {
143 Point2D [] part = new Point2D[tmp.size()];
144 parts.add(tmp.toArray(part));
145 tmp.clear();
146 }
147 }
148 else {
149 tmp.add(new Point2D.Double(x, y));
150 }
151 }
152
153 if (!tmp.isEmpty()) {
154 Point2D [] part = new Point2D[tmp.size()];
155 parts.add(tmp.toArray(part));
156 }
157
158 return parts;
159 }
160
161 protected static boolean removeNoneIntersecting(
162 List<Point2D []> As,
163 List<Point2D []> Bs
164 ) {
165 int B = Bs.size()-1;
166 OUTER: for (int i = 0; i < As.size();) {
167 Point2D [] a = As.get(i);
168 int lo = 0, hi = B;
169 while (lo <= hi) {
170 int mid = (lo + hi) >> 1;
171 Point2D [] b = Bs.get(mid);
172 if (X(a[0]) > X(b[b.length-1])) lo = mid+1;
173 else if (X(a[a.length-1]) < X(b[0])) hi = mid-1;
174 else {
175 // found: keep
176 ++i;
177 continue OUTER;
178 }
179 }
180 // not found: remove
181 As.remove(i);
182 }
183
184 return As.isEmpty();
185 }
186
187 protected static void buildPolygons(
188 Point2D [] as,
189 Point2D [] bs,
190 Point2D [][] rest,
191 List<Polygon2D> positives,
192 List<Polygon2D> negatives
193 ) {
194 List<Point2D> apoints = new ArrayList<Point2D>();
195 List<Point2D> bpoints = new ArrayList<Point2D>();
196
197 double ax = X(as[0]);
198 double bx = X(bs[0]);
199
200 int ai = 1;
201 int bi = 1;
202
203 boolean intersect = false;
204
205 if (ax == bx) {
206 apoints.add(as[0]);
207 bpoints.add(bs[0]);
208 }
209 else if (ax > bx) {
210 apoints.add(as[0]);
211 double bx1;
212 while ((bx1 = X(bs[bi])) < ax) ++bi;
213 if (bx1 == ax) {
214 bpoints.add(bs[bi]);
215 }
216 else { // need to calculate start b point.
217 intersect = true;
218 double by1 = Linear.linear(
219 ax,
220 X(bs[bi-1]), bx1,
221 Y(bs[bi-1]), Y(bs[bi]));
222
223 bpoints.add(new Point2D.Double(ax, by1));
224 }
225 }
226 else { // bx > ax: Symmetric case
227 bpoints.add(bs[0]);
228 double ax1;
229 while ((ax1 = X(as[ai])) < bx) ++ai;
230 if (ax1 == bx) {
231 apoints.add(as[ai]);
232 }
233 else { // need to calculate start b point.
234 intersect = true;
235 double ay1 = Linear.linear(
236 bx,
237 X(as[ai-1]), ax1,
238 Y(as[ai-1]), Y(as[ai]));
239
240 apoints.add(new Point2D.Double(bx, ay1));
241 }
242 }
243
244 // now we have a point in each list, decide if neg/pos.
245 boolean neg = Y(bpoints.get(0)) > Y(apoints.get(0));
246
247 // Continue with inner points
248
249 Line line = new Line();
250
251 while (ai < as.length && bi < bs.length) {
252 double xan = X(as[ai]);
253 double xbn = X(bs[bi]);
254 if (xan == xbn) {
255 double yan = Y(as[ai]);
256 double ybn = Y(bs[ai]);
257
258 if (neg) {
259 if (yan > ybn) { // intersection
260 Point2D ip = VectorUtils.intersection(
261 apoints.get(apoints.size()-1), as[ai],
262 bpoints.get(bpoints.size()-1), bs[bi]);
263 addCheck(ip, apoints);
264 addCheck(ip, bpoints);
265 Polygon2D p = new Polygon2D(
266 new ArrayList<Point2D>(apoints));
267 p.addReversed(bpoints);
268 negatives.add(p);
269 apoints.clear();
270 bpoints.clear();
271 apoints.add(ip);
272 bpoints.add(ip);
273 neg = !neg;
274 }
275 else { // no intersection
276 addCheck(as[ai], apoints);
277 addCheck(bs[bi], bpoints);
278 }
279 }
280 else { // not neg
281 if (ybn > yan) { // intersection
282 Point2D ip = VectorUtils.intersection(
283 apoints.get(apoints.size()-1), as[ai],
284 bpoints.get(bpoints.size()-1), bs[bi]);
285 addCheck(ip, apoints);
286 addCheck(ip, bpoints);
287 Polygon2D p = new Polygon2D(
288 new ArrayList<Point2D>(apoints));
289 p.addReversed(bpoints);
290 positives.add(p);
291 apoints.clear();
292 bpoints.clear();
293 apoints.add(ip);
294 bpoints.add(ip);
295 neg = !neg;
296 }
297 else { // no intersection
298 addCheck(as[ai], apoints);
299 addCheck(bs[bi], bpoints);
300 }
301 }
302 ++ai;
303 ++bi;
304 }
305 else if (xan > xbn) {
306 line.set(apoints.get(apoints.size()-1), as[ai]);
307 double dir = neg ? -1d: 1d; // XXX: correct sign?
308 while (bi < bs.length
309 && X(bs[bi]) < xan
310 && line.eval(bs[bi])*dir > EPSILON)
311 ++bi;
312 if (bi == bs.length) {
313 // b run out of points
314 // calculate Y(last_a, as[ai]) for X(bs[bi-1])
315 double ay1 = Linear.linear(
316 X(bs[bi-1]),
317 X(apoints.get(apoints.size()-1)), X(as[ai]),
318 Y(apoints.get(apoints.size()-1)), Y(as[ai]));
319 addCheck(new Point2D.Double(X(bs[bi-1]), ay1), apoints);
320 addCheck(bs[bi-1], bpoints);
321 Polygon2D p = new Polygon2D(
322 new ArrayList<Point2D>(apoints));
323 p.addReversed(bpoints);
324 apoints.clear();
325 bpoints.clear();
326 (neg ? negatives : positives).add(p);
327 break;
328 }
329 else {
330 // TODO: intersect line and/or X(bs[bi]) >= xan?
331 }
332 }
333 else { // xbn > xan
334 line.set(bpoints.get(bpoints.size()-1), bs[bi]);
335 // TODO: continue symmetric
336 }
337 }
338
339 // TODO: Continue with closing segment
340 }
341
342 public static final class Line {
343
344 private double a;
345 private double b;
346 private double c;
347
348 public Line() {
349 }
350
351 public Line(Point2D p1, Point2D p2) {
352 set(p1, p2);
353 }
354
355 public void set(Point2D p1, Point2D p2) {
356 Point2D p3 =
357 VectorUtils.normalize(
358 VectorUtils.sub(p1, p2));
359
360 Point2D n = VectorUtils.ortho(p3);
361
362 a = X(n);
363 b = Y(n);
364
365 // a*x + b*y + c = 0
366 // c = -a*x -b*y
367
368 c = -a*X(p1) - b*Y(p1);
369 }
370
371 public double eval(Point2D p) {
372 return a*X(p) + b*Y(p) + c;
373 }
374 }
375
376 public static void createPolygons(
377 double [] xAs, double [] yAs,
378 double [] xBs, double [] yBs,
379 List<Polygon2D> positives,
380 List<Polygon2D> negatives
381 ) {
382 if (xAs.length == 0 || xBs.length == 0) {
383 return;
384 }
385
386 List<Point2D []> splAs = splitByNaNs(xAs, yAs);
387 List<Point2D []> splBs = splitByNaNs(xBs, yBs);
388
389 // They feeded us with NaNs only.
390 if (splAs.isEmpty() || splBs.isEmpty()) {
391 return;
392 }
393
394 // Sort each part by x to ensure that the first
395 // is the smallest.
396 for (Point2D [] splA: splAs) {
397 Arrays.sort(splA, POINT_X_CMP);
398 }
399
400 for (Point2D [] splB: splBs) {
401 Arrays.sort(splB, POINT_X_CMP);
402 }
403
404 // Now sort all parts by there first elements.
405 // Should be good enough to find overlapping regions.
406 Collections.sort(splAs, FIRST_POINT_X);
407 Collections.sort(splBs, FIRST_POINT_X);
408
409 // Check if the two series intersect at all.
410 // If no then there will be no area between them.
411
412 Point2D [] p1 = splAs.get(0);
413 Point2D [] p2 = splBs.get(splBs.size()-1);
414
415 // Sort out the ranges that are not intersecting
416 // the ranges in the other series.
417 // We are going to merge them anyway
418 // so this is not strictly required.
419 // Keep it to recude cases.
420 if (removeNoneIntersecting(splAs, splBs)
421 || removeNoneIntersecting(splBs, splAs)
422 ) {
423 // They do not intersect at all.
424 return;
425 }
426
427 // TODO: Intersect/split the two series parts.
428 }
429 }
430 // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :
http://dive4elements.wald.intevation.org