Mercurial > dive4elements > river
view flys-artifacts/src/main/java/de/intevation/flys/artifacts/geom/Polygon2D.java @ 2230:59af81364eb1
Improved the 'historical discharge' calculation: implemented findValueForW().
flys-artifacts/trunk@3872 c6561f87-3c4e-4783-a992-168aeb5c3f6f
author | Ingo Weinzierl <ingo.weinzierl@intevation.de> |
---|---|
date | Wed, 01 Feb 2012 15:41:11 +0000 |
parents | 6f83d9d434f2 |
children | 5642a83420f2 |
line wrap: on
line source
package de.intevation.flys.geom; import java.io.Serializable; import java.awt.Shape; import java.awt.geom.Path2D; import java.awt.geom.Point2D; import java.util.ArrayList; import java.util.List; import java.util.Arrays; import java.util.Comparator; import java.util.Collections; import de.intevation.flys.artifacts.math.Linear; import static de.intevation.flys.geom.VectorUtils.X; import static de.intevation.flys.geom.VectorUtils.Y; import static de.intevation.flys.geom.VectorUtils.EPSILON; public class Polygon2D implements Serializable { public static final Comparator<Point2D> POINT_X_CMP = new Comparator<Point2D>() { public int compare(Point2D a, Point2D b) { double d = X(a) - X(b); if (d < 0d) return -1; return d > 0d ? + 1 : 0; } }; public static final Comparator<Point2D []> FIRST_POINT_X = new Comparator<Point2D []>() { public int compare(Point2D [] a, Point2D [] b) { if (a.length == 0) return -1; // should not happen. if (b.length == 0) return +1; // should not happen. double d = X(a[0]) - Y(b[0]); if (d < 0d) return -1; return d > 0d ? + 1 : 0; } }; protected List<Point2D> points; public Polygon2D() { points = new ArrayList<Point2D>(); } public Polygon2D(List<Point2D> points) { this.points = points; } public void add(double x, double y) { points.add(new Point2D.Double(x, y)); } protected static boolean addCheck(Point2D p, List<Point2D> points) { switch (points.size()) { case 0: points.add(p); return true; case 1: if (VectorUtils.epsilonEquals(points.get(0), p)) { return false; } points.add(p); return true; default: int L = points.size()-1; Point2D last = points.get(L); if (VectorUtils.epsilonEquals(last, p)) { return false; } Point2D before = points.get(L-1); if (VectorUtils.collinear(before, last, p)) { points.set(L, p); } else { points.add(p); } return true; } } public boolean addCheck(Point2D p) { return addCheck(p, points); } public void addReversed(List<Point2D> other) { for (int i = other.size()-1; i >= 0; --i) { addCheck(other.get(i)); } } public double area() { double area = 0d; for (int i = 0, N = points.size(); i < N; ++i) { int j = (i + 1) % N; Point2D pi = points.get(i); Point2D pj = points.get(j); area += X(pi)*Y(pj); area -= X(pj)*Y(pi); } return 0.5d*area; } public Shape toShape() { Path2D.Double path = new Path2D.Double(); int N = points.size(); if (N > 0) { Point2D p = points.get(0); path.moveTo(X(p), Y(p)); for (int i = 1; i < N; ++i) { p = points.get(i); path.lineTo(X(p), Y(p)); } path.closePath(); } return path; } protected static List<Point2D []> splitByNaNs( double [] xs, double [] ys ) { List<Point2D []> parts = new ArrayList<Point2D []>(); List<Point2D> tmp = new ArrayList<Point2D>(xs.length); for (int i = 0; i < xs.length; ++i) { double x = xs[i]; double y = ys[i]; if (Double.isNaN(x) || Double.isNaN(y)) { if (!tmp.isEmpty()) { Point2D [] part = new Point2D[tmp.size()]; parts.add(tmp.toArray(part)); tmp.clear(); } } else { tmp.add(new Point2D.Double(x, y)); } } if (!tmp.isEmpty()) { Point2D [] part = new Point2D[tmp.size()]; parts.add(tmp.toArray(part)); } return parts; } protected static boolean removeNoneIntersecting( List<Point2D []> As, List<Point2D []> Bs ) { int B = Bs.size()-1; OUTER: for (int i = 0; i < As.size();) { Point2D [] a = As.get(i); int lo = 0, hi = B; while (lo <= hi) { int mid = (lo + hi) >> 1; Point2D [] b = Bs.get(mid); if (X(a[0]) > X(b[b.length-1])) lo = mid+1; else if (X(a[a.length-1]) < X(b[0])) hi = mid-1; else { // found: keep ++i; continue OUTER; } } // not found: remove As.remove(i); } return As.isEmpty(); } protected static void buildPolygons( Point2D [] as, Point2D [] bs, Point2D [][] rest, List<Polygon2D> positives, List<Polygon2D> negatives ) { List<Point2D> apoints = new ArrayList<Point2D>(); List<Point2D> bpoints = new ArrayList<Point2D>(); double ax = X(as[0]); double bx = X(bs[0]); int ai = 1; int bi = 1; boolean intersect = false; if (ax == bx) { apoints.add(as[0]); bpoints.add(bs[0]); } else if (ax > bx) { apoints.add(as[0]); double bx1; while ((bx1 = X(bs[bi])) < ax) ++bi; if (bx1 == ax) { bpoints.add(bs[bi]); } else { // need to calculate start b point. intersect = true; double by1 = Linear.linear( ax, X(bs[bi-1]), bx1, Y(bs[bi-1]), Y(bs[bi])); bpoints.add(new Point2D.Double(ax, by1)); } } else { // bx > ax: Symmetric case bpoints.add(bs[0]); double ax1; while ((ax1 = X(as[ai])) < bx) ++ai; if (ax1 == bx) { apoints.add(as[ai]); } else { // need to calculate start b point. intersect = true; double ay1 = Linear.linear( bx, X(as[ai-1]), ax1, Y(as[ai-1]), Y(as[ai])); apoints.add(new Point2D.Double(bx, ay1)); } } // now we have a point in each list, decide if neg/pos. boolean neg = Y(bpoints.get(0)) > Y(apoints.get(0)); // Continue with inner points Line line = new Line(); while (ai < as.length && bi < bs.length) { double xan = X(as[ai]); double xbn = X(bs[bi]); if (xan == xbn) { double yan = Y(as[ai]); double ybn = Y(bs[ai]); if (neg) { if (yan > ybn) { // intersection Point2D ip = VectorUtils.intersection( apoints.get(apoints.size()-1), as[ai], bpoints.get(bpoints.size()-1), bs[bi]); addCheck(ip, apoints); addCheck(ip, bpoints); Polygon2D p = new Polygon2D( new ArrayList<Point2D>(apoints)); p.addReversed(bpoints); negatives.add(p); apoints.clear(); bpoints.clear(); apoints.add(ip); bpoints.add(ip); neg = !neg; } else { // no intersection addCheck(as[ai], apoints); addCheck(bs[bi], bpoints); } } else { // not neg if (ybn > yan) { // intersection Point2D ip = VectorUtils.intersection( apoints.get(apoints.size()-1), as[ai], bpoints.get(bpoints.size()-1), bs[bi]); addCheck(ip, apoints); addCheck(ip, bpoints); Polygon2D p = new Polygon2D( new ArrayList<Point2D>(apoints)); p.addReversed(bpoints); positives.add(p); apoints.clear(); bpoints.clear(); apoints.add(ip); bpoints.add(ip); neg = !neg; } else { // no intersection addCheck(as[ai], apoints); addCheck(bs[bi], bpoints); } } ++ai; ++bi; } else if (xan > xbn) { line.set(apoints.get(apoints.size()-1), as[ai]); double dir = neg ? -1d: 1d; // XXX: correct sign? while (bi < bs.length && X(bs[bi]) < xan && line.eval(bs[bi])*dir > EPSILON) ++bi; if (bi == bs.length) { // b run out of points // calculate Y(last_a, as[ai]) for X(bs[bi-1]) double ay1 = Linear.linear( X(bs[bi-1]), X(apoints.get(apoints.size()-1)), X(as[ai]), Y(apoints.get(apoints.size()-1)), Y(as[ai])); addCheck(new Point2D.Double(X(bs[bi-1]), ay1), apoints); addCheck(bs[bi-1], bpoints); Polygon2D p = new Polygon2D( new ArrayList<Point2D>(apoints)); p.addReversed(bpoints); apoints.clear(); bpoints.clear(); (neg ? negatives : positives).add(p); break; } else { // TODO: intersect line and/or X(bs[bi]) >= xan? } } else { // xbn > xan line.set(bpoints.get(bpoints.size()-1), bs[bi]); // TODO: continue symmetric } } // TODO: Continue with closing segment } public static final class Line { private double a; private double b; private double c; public Line() { } public Line(Point2D p1, Point2D p2) { set(p1, p2); } public void set(Point2D p1, Point2D p2) { Point2D p3 = VectorUtils.normalize( VectorUtils.sub(p1, p2)); Point2D n = VectorUtils.ortho(p3); a = X(n); b = Y(n); // a*x + b*y + c = 0 // c = -a*x -b*y c = -a*X(p1) - b*Y(p1); } public double eval(Point2D p) { return a*X(p) + b*Y(p) + c; } } public static void createPolygons( double [] xAs, double [] yAs, double [] xBs, double [] yBs, List<Polygon2D> positives, List<Polygon2D> negatives ) { if (xAs.length == 0 || xBs.length == 0) { return; } List<Point2D []> splAs = splitByNaNs(xAs, yAs); List<Point2D []> splBs = splitByNaNs(xBs, yBs); // They feeded us with NaNs only. if (splAs.isEmpty() || splBs.isEmpty()) { return; } // Sort each part by x to ensure that the first // is the smallest. for (Point2D [] splA: splAs) { Arrays.sort(splA, POINT_X_CMP); } for (Point2D [] splB: splBs) { Arrays.sort(splB, POINT_X_CMP); } // Now sort all parts by there first elements. // Should be good enough to find overlapping regions. Collections.sort(splAs, FIRST_POINT_X); Collections.sort(splBs, FIRST_POINT_X); // Check if the two series intersect at all. // If no then there will be no area between them. Point2D [] p1 = splAs.get(0); Point2D [] p2 = splBs.get(splBs.size()-1); // Sort out the ranges that are not intersecting // the ranges in the other series. // We are going to merge them anyway // so this is not strictly required. // Keep it to recude cases. if (removeNoneIntersecting(splAs, splBs) || removeNoneIntersecting(splBs, splAs) ) { // They do not intersect at all. return; } // TODO: Intersect/split the two series parts. } } // vim:set ts=4 sw=4 si et sta sts=4 fenc=utf8 :