From 267400485e5d683dce4989968875d26069657230 Mon Sep 17 00:00:00 2001 From: Nick Knize Date: Thu, 31 Dec 2015 16:08:25 +0000 Subject: [PATCH] LUCENE-6908: Add space segmentation for handling irregular rectangle accuracy at the poles. git-svn-id: https://svn.apache.org/repos/asf/lucene/dev/trunk@1722448 13f79535-47bb-0310-9956-ffa450edef68 --- .../lucene/util/GeoProjectionUtils.java | 16 +++++ .../apache/lucene/util/GeoRelationUtils.java | 71 ++++++++++++++----- 2 files changed, 70 insertions(+), 17 deletions(-) diff --git a/lucene/sandbox/src/java/org/apache/lucene/util/GeoProjectionUtils.java b/lucene/sandbox/src/java/org/apache/lucene/util/GeoProjectionUtils.java index 0c797fd969f..4ff80f2e123 100644 --- a/lucene/sandbox/src/java/org/apache/lucene/util/GeoProjectionUtils.java +++ b/lucene/sandbox/src/java/org/apache/lucene/util/GeoProjectionUtils.java @@ -419,4 +419,20 @@ public class GeoProjectionUtils { return pt; } + /** + * Finds the bearing (in degrees) between 2 geo points (lon, lat) using great circle arc + * @param lon1 first point longitude in degrees + * @param lat1 first point latitude in degrees + * @param lon2 second point longitude in degrees + * @param lat2 second point latitude in degrees + * @return the bearing (in degrees) between the two provided points + */ + public static double bearingGreatCircle(double lon1, double lat1, double lon2, double lat2) { + double dLon = (lon2 - lon1) * TO_RADIANS; + lat2 *= TO_RADIANS; + lat1 *= TO_RADIANS; + double y = SloppyMath.sin(dLon) * SloppyMath.cos(lat2); + double x = SloppyMath.cos(lat1) * SloppyMath.sin(lat2) - SloppyMath.sin(lat1) * SloppyMath.cos(lat2) * SloppyMath.cos(dLon); + return Math.atan2(y, x) * TO_DEGREES; + } } diff --git a/lucene/sandbox/src/java/org/apache/lucene/util/GeoRelationUtils.java b/lucene/sandbox/src/java/org/apache/lucene/util/GeoRelationUtils.java index 6c4d8ed2860..40b46b0293c 100644 --- a/lucene/sandbox/src/java/org/apache/lucene/util/GeoRelationUtils.java +++ b/lucene/sandbox/src/java/org/apache/lucene/util/GeoRelationUtils.java @@ -200,35 +200,72 @@ public class GeoRelationUtils { || SloppyMath.haversin(centerLat, centerLon, rMinY, rMaxX)*1000.0 <= radiusMeters; } + /** + * Compute whether any of the 4 corners of the rectangle (defined by min/max X/Y) are outside the circle (defined + * by centerLon, centerLat, radiusMeters) + * + * Note: exotic rectangles at the poles (e.g., those whose lon/lat distance ratios greatly deviate from 1) can not + * be determined by using distance alone. For this reason the approx flag may be set to false, in which case the + * space will be further divided to more accurately compute whether the rectangle crosses the circle + */ private static boolean rectAnyCornersOutsideCircle(final double rMinX, final double rMinY, final double rMaxX, final double rMaxY, final double centerLon, final double centerLat, final double radiusMeters, final boolean approx) { if (approx == true) { return rectAnyCornersOutsideCircleSloppy(rMinX, rMinY, rMaxX, rMaxY, centerLon, centerLat, radiusMeters); } - double w = Math.abs(rMaxX - rMinX); - if (w <= 90.0) { + // if span is less than 70 degrees we can approximate using distance alone + if (Math.abs(rMaxX - rMinX) <= 70.0) { return GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, rMinX) > radiusMeters || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, rMinX) > radiusMeters || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, rMaxX) > radiusMeters || GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, rMaxX) > radiusMeters; } - // partition - w /= 4; - final double p1 = rMinX + w; - final double p2 = p1 + w; - final double p3 = p2 + w; + return rectCrossesOblateCircle(centerLon, centerLat, radiusMeters, rMinX, rMinY, rMaxX, rMaxY); + } - return GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, rMinX) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, rMinX) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, p1) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, p1) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, p2) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, p2) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, p3) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, p3) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxY, rMaxX) > radiusMeters - || GeoDistanceUtils.haversin(centerLat, centerLon, rMinY, rMaxX) > radiusMeters; + /** + * Compute whether the rectangle (defined by min/max Lon/Lat) crosses a potentially oblate circle + * + * TODO benchmark for replacing existing rectCrossesCircle. + */ + public static boolean rectCrossesOblateCircle(double centerLon, double centerLat, double radiusMeters, double rMinLon, double rMinLat, double rMaxLon, double rMaxLat) { + double w = Math.abs(rMaxLon - rMinLon); + final int segs = (int)Math.ceil(w / 45.0); + w /= segs; + short i = 1; + double p1 = rMinLon; + double maxLon, midLon; + double[] pt = new double[2]; + + do { + maxLon = (i == segs) ? rMaxLon : p1 + w; + + final double d1, d2; + // short-circuit if we find a corner outside the circle + if ( (d1 = GeoDistanceUtils.haversin(centerLat, centerLon, rMinLat, p1)) > radiusMeters + || (d2 = GeoDistanceUtils.haversin(centerLat, centerLon, rMinLat, maxLon)) > radiusMeters + || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxLat, p1) > radiusMeters + || GeoDistanceUtils.haversin(centerLat, centerLon, rMaxLat, maxLon) > radiusMeters) { + return true; + } + + // else we treat as an oblate circle by slicing the longitude space and checking the azimuthal range + // OPTIMIZATION: this is only executed for latitude values "closeTo" the poles (e.g., 88.0 > lat < -88.0) + if ( (rMaxLat > 88.0 || rMinLat < -88.0) + && (pt = GeoProjectionUtils.pointFromLonLatBearingGreatCircle(p1, rMinLat, + GeoProjectionUtils.bearingGreatCircle(p1, rMinLat, p1, rMaxLat), radiusMeters - d1, pt))[1] < rMinLat || pt[1] < rMaxLat + || (pt = GeoProjectionUtils.pointFromLonLatBearingGreatCircle(maxLon, rMinLat, + GeoProjectionUtils.bearingGreatCircle(maxLon, rMinLat, maxLon, rMaxLat), radiusMeters - d2, pt))[1] < rMinLat || pt[1] < rMaxLat + || (pt = GeoProjectionUtils.pointFromLonLatBearingGreatCircle(maxLon, rMinLat, + GeoProjectionUtils.bearingGreatCircle(maxLon, rMinLat, (midLon = p1 + 0.5*(maxLon - p1)), rMaxLat), + radiusMeters - GeoDistanceUtils.haversin(centerLat, centerLon, rMinLat, midLon), pt))[1] < rMinLat + || pt[1] < rMaxLat == false ) { + return true; + } + p1 += w; + } while (++i <= segs); + return false; } private static boolean rectAnyCornersOutsideCircleSloppy(final double rMinX, final double rMinY, final double rMaxX, final double rMaxY,