From 6cf0e2a3e712fa87c5c39b01d085ec8ff477089e Mon Sep 17 00:00:00 2001 From: Karl Wright Date: Wed, 6 Apr 2016 09:22:11 -0400 Subject: [PATCH] LUCENE-7173: Add complex test logic for creating nested polygons --- .../lucene/spatial3d/TestGeo3DPoint.java | 209 ++++++++++++++++++ 1 file changed, 209 insertions(+) diff --git a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/TestGeo3DPoint.java b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/TestGeo3DPoint.java index 3caf039a2d0..01f5e830645 100644 --- a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/TestGeo3DPoint.java +++ b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/TestGeo3DPoint.java @@ -46,6 +46,9 @@ import org.apache.lucene.spatial3d.geom.GeoPolygonFactory; import org.apache.lucene.spatial3d.geom.GeoShape; import org.apache.lucene.spatial3d.geom.PlanetModel; import org.apache.lucene.spatial3d.geom.XYZBounds; +import org.apache.lucene.spatial3d.geom.SidedPlane; +import org.apache.lucene.spatial3d.geom.Plane; +import org.apache.lucene.spatial3d.geom.GeoPolygon; import org.apache.lucene.index.DirectoryReader; import org.apache.lucene.index.IndexReader; import org.apache.lucene.index.IndexWriter; @@ -60,6 +63,7 @@ import org.apache.lucene.search.IndexSearcher; import org.apache.lucene.search.Query; import org.apache.lucene.search.SimpleCollector; import org.apache.lucene.store.Directory; +import org.apache.lucene.geo.Polygon; import org.apache.lucene.util.FixedBitSet; import org.apache.lucene.util.IOUtils; import org.apache.lucene.util.LuceneTestCase; @@ -821,4 +825,209 @@ public class TestGeo3DPoint extends LuceneTestCase { assertFalse(q.equals(Geo3DPoint.newShapeQuery("point", shape2))); } + + protected static double MINIMUM_EDGE_ANGLE = Math.toRadians(5.0); + protected static double MINIMUM_ARC_ANGLE = Math.toRadians(1.0); + + /** Cook up a random Polygon that makes sense, with possible nested polygon within. + * This is part of testing more complex polygons with nested holes. Picking random points + * doesn't do it because it's almost impossible to come up with nested ones of the proper + * clockwise/counterclockwise rotation that way. + */ + protected Polygon makePoly(final PlanetModel pm, final GeoPoint pole, final boolean clockwiseDesired) { + // Polygon edges will be arranged around the provided pole, and holes will each have a pole selected within the parent + // polygon. + final int pointCount = TestUtil.nextInt(random(), 3, 10); + // The point angles we pick next. The only requirement is that they are not all on one side of the pole. + // We arrange that by picking the next point within what's left of the remaining angle, but never more than 180 degrees, + // and never less than what is needed to insure that the remaining point choices are less than 180 degrees always. + // These are all picked in the context of the pole, + final double[] angles = new double[pointCount]; + final double[] arcDistance = new double[pointCount]; + double accumulatedAngle = 0.0; + for (int i = 0; i < pointCount; i++) { + final int remainingEdgeCount = pointCount - i; + final double remainingAngle = 2.0 * Math.PI - accumulatedAngle; + if (remainingEdgeCount == 1) { + angles[i] = remainingAngle; + } else { + // The maximum angle is 180 degrees, or what's left when you give a minimal amount to each edge. + double maximumAngle = remainingAngle - (remainingEdgeCount-1) * MINIMUM_EDGE_ANGLE; + if (maximumAngle > Math.PI) { + maximumAngle = Math.PI; + } + // The minimum angle is MINIMUM_EDGE_ANGLE, or enough to be sure nobody afterwards needs more than + // 180 degrees. And since we have three points to start with, we already know that. + final double minimumAngle = MINIMUM_EDGE_ANGLE; + // Pick the angle + final double angle = random().nextDouble() * (maximumAngle - minimumAngle) + minimumAngle; + angles[i] = angle; + accumulatedAngle += angle; + } + // Pick the arc distance randomly + arcDistance[i] = random().nextDouble() * (Math.PI * 0.5 - MINIMUM_ARC_ANGLE) + MINIMUM_ARC_ANGLE; + } + if (clockwiseDesired) { + // Reverse the signs + for (int i = 0; i < pointCount; i++) { + angles[i] = -angles[i]; + } + } + + // Now, use the pole's information plus angles and arcs to create GeoPoints in the right order. + final List polyPoints = convertToPoints(pm, pole, angles, arcDistance); + + // Create the geo3d polygon, so we can test out our poles. + GeoPolygon poly = GeoPolygonFactory.makeGeoPolygon(pm, polyPoints, null); + + // Next, do some holes. No more than 2 of these. The poles for holes must always be within the polygon, so we're + // going to use Geo3D to help us select those given the points we just made. + + final int holeCount = TestUtil.nextInt(random(), 0, 2); + + final Polygon[] holes = new Polygon[holeCount]; + + for (int i = 0; i < holeCount; i++) { + // Choose a pole. The poly has to be within the polygon, but it also cannot be on the polygon edge. + // We try indefinitely to find a good pole... + while (true) { + final GeoPoint poleChoice = new GeoPoint(pm, toRadians(randomLat()), toRadians(randomLon())); + if (!poly.isWithin(poleChoice)) { + continue; + } + // We have a pole within the polygon. Now try 100 times to build a polygon that does not intersect the outside ring. + // After that we give up and pick a new pole. + boolean foundOne = false; + for (int j = 0; j < 100; j++) { + final Polygon insidePoly = makePoly(pm, poleChoice, !clockwiseDesired); + // Verify that the inside polygon is OK. If not, discard and repeat. + if (!verifyPolygon(pm, insidePoly, poly)) { + continue; + } + holes[i] = insidePoly; + foundOne = true; + } + if (foundOne) { + break; + } + } + } + + // Finally, build the polygon and return it + final double[] lats = new double[polyPoints.size() + 1]; + final double[] lons = new double[polyPoints.size() + 1]; + + for (int i = 0; i < polyPoints.size(); i++) { + lats[i] = polyPoints.get(i).getLatitude() * 180.0 / Math.PI; + lons[i] = polyPoints.get(i).getLongitude() * 180.0 / Math.PI; + } + lats[polyPoints.size()] = lats[0]; + lons[polyPoints.size()] = lons[0]; + return new Polygon(lats, lons, holes); + } + + protected static List convertToPoints(final PlanetModel pm, final GeoPoint pole, final double[] angles, final double[] arcDistances) { + // To do the point rotations, we need the sine and cosine of the pole latitude and longitude. Get it here for performance. + final double sinLatitude = Math.sin(pole.getLatitude()); + final double cosLatitude = Math.cos(pole.getLatitude()); + final double sinLongitude = Math.sin(pole.getLongitude()); + final double cosLongitude = Math.cos(pole.getLongitude()); + final List rval = new ArrayList<>(); + for (int i = 0; i < angles.length; i++) { + rval.add(createPoint(pm, angles[i], arcDistances[i], sinLatitude, cosLatitude, sinLongitude, cosLongitude)); + } + return rval; + } + + protected static GeoPoint createPoint(final PlanetModel pm, + final double angle, + final double arcDistance, + final double sinLatitude, + final double cosLatitude, + final double sinLongitude, + final double cosLongitude) { + // From the angle and arc distance, convert to (x,y,z) in unit space. + // We want the perspective to be looking down the x axis. The "angle" measurement is thus in the Y-Z plane. + // The arcdistance is in X. + final double x = Math.cos(arcDistance); + final double yzScale = Math.sin(arcDistance); + final double y = Math.cos(angle) * yzScale; + final double z = Math.sin(angle) * yzScale; + // Now, rotate coordinates so that we shift everything from pole = x-axis to actual coordinates. + // This transformation should take the point (1,0,0) and transform it to the pole's actual (x,y,z) coordinates. + // Coordinate rotation formula: + // x1 = x0 cos T - y0 sin T + // y1 = x0 sin T + y0 cos T + // We're in essence undoing the following transformation (from GeoPolygonFactory): + // x1 = x0 cos az - y0 sin az + // y1 = x0 sin az + y0 cos az + // z1 = z0 + // x2 = x1 cos al - z1 sin al + // y2 = y1 + // z2 = x1 sin al + z1 cos al + // So, we reverse the order of the transformations, AND we transform backwards. + // Transforming backwards means using these identities: sin(-angle) = -sin(angle), cos(-angle) = cos(angle) + // So: + // x1 = x0 cos al + z0 sin al + // y1 = y0 + // z1 = - x0 sin al + z0 cos al + // x2 = x1 cos az + y1 sin az + // y2 = - x1 sin az + y1 cos az + // z2 = z1 + final double x1 = x * cosLatitude + z * sinLatitude; + final double y1 = y; + final double z1 = - x * sinLatitude + z * cosLatitude; + final double x2 = x1 * cosLongitude + y1 * sinLongitude; + final double y2 = - x1 * sinLongitude + y1 * cosLongitude; + final double z2 = z1; + + // Scale final (x,y,z) to land on planet surface + // Equation of ellipsoid: x^2 / a^2 + y^2 / b^2 + z^2 / c^2 - 1 = 0 + // Use a parameterization, e.g. x = t * x2, y = t * y2, z = t * z2, and find t. + // t^2 ( x2^2 / a^2 + y2^2 / b^2 + z2^2 / c^2 ) = 1 + // t = +/- sqrt( 1 / ( x2^2 / a^2 + y2^2 / b^2 + z2^2 / c^2 ) ) + // We want the + variant because we're scaling in the same direction as the original vector. + final double t = Math.sqrt( 1.0 / (x2 * x2 * pm.inverseAbSquared + y2 * y2 * pm.inverseAbSquared + z2 * z2 * pm.inverseCSquared)); + return new GeoPoint(x2 * t, y2 * t, z2 * t); + } + + protected static boolean verifyPolygon(final PlanetModel pm, final Polygon polygon, final GeoPolygon outsidePolygon) { + // Each point in the new poly should be inside the outside poly, and each edge should not intersect the outside poly edge + final double[] lats = polygon.getPolyLats(); + final double[] lons = polygon.getPolyLons(); + final List polyPoints = new ArrayList<>(lats.length-1); + for (int i = 0; i < lats.length - 1; i++) { + final GeoPoint newPoint = new GeoPoint(pm, Math.toRadians(lats[i]), Math.toRadians(lons[i])); + if (!outsidePolygon.isWithin(newPoint)) { + return false; + } + polyPoints.add(newPoint); + } + // We don't need to construct the world to find intersections -- just the bordering planes. + for (int planeIndex = 0; planeIndex < polyPoints.size(); planeIndex++) { + final GeoPoint startPoint = polyPoints.get(planeIndex); + final GeoPoint endPoint = polyPoints.get(legalIndex(planeIndex + 1, polyPoints.size())); + final GeoPoint beforeStartPoint = polyPoints.get(legalIndex(planeIndex - 1, polyPoints.size())); + final GeoPoint afterEndPoint = polyPoints.get(legalIndex(planeIndex + 2, polyPoints.size())); + final SidedPlane beforePlane = new SidedPlane(endPoint, beforeStartPoint, startPoint); + final SidedPlane afterPlane = new SidedPlane(startPoint, endPoint, afterEndPoint); + final Plane plane = new Plane(startPoint, endPoint); + + // Check for intersections!! + if (outsidePolygon.intersects(plane, null, beforePlane, afterPlane)) { + return false; + } + } + return true; + } + + protected static int legalIndex(int index, int size) { + if (index > size) { + index -= size; + } + if (index < 0) { + index += size; + } + return index; + } }