LUCENE-7127: remove epsilon-based testing from lucene/spatial, fix distance bugs.

This commit is contained in:
Robert Muir 2016-03-22 16:53:25 -04:00
parent c5da271b9d
commit ee1aca8643
13 changed files with 234 additions and 1094 deletions

View File

@ -29,6 +29,10 @@ Optimizations
* LUCENE-7115: Speed up FieldCache.CacheEntry toString by setting initial
StringBuilder capacity (Gregory Chanan)
Bug Fixes
* LUCENE-7127: Fix corner case bugs in GeoPointDistanceQuery. (Robert Muir)
======================= Lucene 6.0.0 =======================
System Requirements

View File

@ -16,29 +16,11 @@
*/
package org.apache.lucene.document;
import java.io.IOException;
import java.util.BitSet;
import org.apache.lucene.codecs.FilterCodec;
import org.apache.lucene.codecs.PointsFormat;
import org.apache.lucene.codecs.PointsReader;
import org.apache.lucene.codecs.PointsWriter;
import org.apache.lucene.codecs.lucene60.Lucene60PointsReader;
import org.apache.lucene.codecs.lucene60.Lucene60PointsWriter;
import org.apache.lucene.index.IndexReader;
import org.apache.lucene.index.IndexWriterConfig;
import org.apache.lucene.index.RandomIndexWriter;
import org.apache.lucene.index.SegmentReadState;
import org.apache.lucene.index.SegmentWriteState;
import org.apache.lucene.search.IndexSearcher;
import org.apache.lucene.search.ScoreDoc;
import org.apache.lucene.search.Sort;
import org.apache.lucene.search.TopDocs;
import org.apache.lucene.store.Directory;
import org.apache.lucene.util.LuceneTestCase;
import org.apache.lucene.util.SloppyMath;
import org.apache.lucene.util.TestUtil;
import org.apache.lucene.util.bkd.BKDWriter;
/** Simple tests for {@link LatLonPoint#newDistanceQuery} */
public class TestLatLonPointDistanceQuery extends LuceneTestCase {
@ -97,94 +79,4 @@ public class TestLatLonPointDistanceQuery extends LuceneTestCase {
assertTrue(expected.getMessage().contains("radiusMeters"));
assertTrue(expected.getMessage().contains("is invalid"));
}
/** Run a few iterations with just 10 docs, hopefully easy to debug */
public void testRandom() throws Exception {
for (int iters = 0; iters < 100; iters++) {
doRandomTest(10, 100);
}
}
/** Runs with thousands of docs */
@Nightly
public void testRandomHuge() throws Exception {
for (int iters = 0; iters < 10; iters++) {
doRandomTest(2000, 100);
}
}
private void doRandomTest(int numDocs, int numQueries) throws IOException {
Directory dir = newDirectory();
IndexWriterConfig iwc = newIndexWriterConfig();
int pointsInLeaf = 2 + random().nextInt(4);
iwc.setCodec(new FilterCodec("Lucene60", TestUtil.getDefaultCodec()) {
@Override
public PointsFormat pointsFormat() {
return new PointsFormat() {
@Override
public PointsWriter fieldsWriter(SegmentWriteState writeState) throws IOException {
return new Lucene60PointsWriter(writeState, pointsInLeaf, BKDWriter.DEFAULT_MAX_MB_SORT_IN_HEAP);
}
@Override
public PointsReader fieldsReader(SegmentReadState readState) throws IOException {
return new Lucene60PointsReader(readState);
}
};
}
});
RandomIndexWriter writer = new RandomIndexWriter(random(), dir, iwc);
for (int i = 0; i < numDocs; i++) {
double latRaw = -90 + 180.0 * random().nextDouble();
double lonRaw = -180 + 360.0 * random().nextDouble();
// pre-normalize up front, so we can just use quantized value for testing and do simple exact comparisons
double lat = LatLonPoint.decodeLatitude(LatLonPoint.encodeLatitude(latRaw));
double lon = LatLonPoint.decodeLongitude(LatLonPoint.encodeLongitude(lonRaw));
Document doc = new Document();
doc.add(new LatLonPoint("field", lat, lon));
doc.add(new StoredField("lat", lat));
doc.add(new StoredField("lon", lon));
writer.addDocument(doc);
}
IndexReader reader = writer.getReader();
IndexSearcher searcher = newSearcher(reader);
for (int i = 0; i < numQueries; i++) {
double lat = -90 + 180.0 * random().nextDouble();
double lon = -180 + 360.0 * random().nextDouble();
double radius = 50000000 * random().nextDouble();
BitSet expected = new BitSet();
for (int doc = 0; doc < reader.maxDoc(); doc++) {
double docLatitude = reader.document(doc).getField("lat").numericValue().doubleValue();
double docLongitude = reader.document(doc).getField("lon").numericValue().doubleValue();
double distance = SloppyMath.haversinMeters(lat, lon, docLatitude, docLongitude);
if (distance <= radius) {
expected.set(doc);
}
}
TopDocs topDocs = searcher.search(LatLonPoint.newDistanceQuery("field", lat, lon, radius), reader.maxDoc(), Sort.INDEXORDER);
BitSet actual = new BitSet();
for (ScoreDoc doc : topDocs.scoreDocs) {
actual.set(doc.doc);
}
try {
assertEquals(expected, actual);
} catch (AssertionError e) {
for (int doc = 0; doc < reader.maxDoc(); doc++) {
double docLatitude = reader.document(doc).getField("lat").numericValue().doubleValue();
double docLongitude = reader.document(doc).getField("lon").numericValue().doubleValue();
double distance = SloppyMath.haversinMeters(lat, lon, docLatitude, docLongitude);
System.out.println("" + doc + ": (" + docLatitude + "," + docLongitude + "), distance=" + distance);
}
throw e;
}
}
reader.close();
writer.close();
dir.close();
}
}

View File

@ -20,8 +20,6 @@ import org.apache.lucene.document.Document;
import org.apache.lucene.document.LatLonPoint;
import org.apache.lucene.spatial.util.BaseGeoPointTestCase;
import org.apache.lucene.spatial.util.GeoRect;
import org.apache.lucene.spatial.util.GeoRelationUtils;
import org.apache.lucene.util.SloppyMath;
public class TestLatLonPointQueries extends BaseGeoPointTestCase {
@ -50,19 +48,6 @@ public class TestLatLonPointQueries extends BaseGeoPointTestCase {
return LatLonPoint.newPolygonQuery(FIELD_NAME, lats, lons);
}
@Override
protected Boolean rectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
assert Double.isNaN(pointLat) == false;
if (rect.minLon < rect.maxLon) {
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, rect.maxLon);
} else {
// Rect crosses dateline:
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, -180.0, rect.maxLon)
|| GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, 180.0);
}
}
@Override
protected double quantizeLat(double latRaw) {
return LatLonPoint.decodeLatitude(LatLonPoint.encodeLatitude(latRaw));
@ -72,34 +57,4 @@ public class TestLatLonPointQueries extends BaseGeoPointTestCase {
protected double quantizeLon(double lonRaw) {
return LatLonPoint.decodeLongitude(LatLonPoint.encodeLongitude(lonRaw));
}
@Override
protected Boolean polyRectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
// TODO write better random polygon tests
// note: logic must be slightly different than rectContainsPoint, to satisfy
// insideness for cases exactly on boundaries.
assert Double.isNaN(pointLat) == false;
assert rect.crossesDateline() == false;
double polyLats[] = new double[] { rect.minLat, rect.maxLat, rect.maxLat, rect.minLat, rect.minLat };
double polyLons[] = new double[] { rect.minLon, rect.minLon, rect.maxLon, rect.maxLon, rect.minLon };
// TODO: separately test this method is 100% correct, here treat it like a black box (like haversin)
return GeoRelationUtils.pointInPolygon(polyLats, polyLons, pointLat, pointLon);
}
@Override
protected Boolean circleContainsPoint(double centerLat, double centerLon, double radiusMeters, double pointLat, double pointLon) {
double distanceMeters = SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon);
boolean result = distanceMeters <= radiusMeters;
//System.out.println(" shouldMatch? centerLon=" + centerLon + " centerLat=" + centerLat + " pointLon=" + pointLon + " pointLat=" + pointLat + " result=" + result + " distanceMeters=" + (distanceKM * 1000));
return result;
}
@Override
protected Boolean distanceRangeContainsPoint(double centerLat, double centerLon, double minRadiusMeters, double radiusMeters, double pointLat, double pointLon) {
final double d = SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon);
return d >= minRadiusMeters && d <= radiusMeters;
}
}

View File

@ -19,7 +19,6 @@ package org.apache.lucene.spatial.geopoint.search;
import org.apache.lucene.search.MultiTermQuery;
import org.apache.lucene.spatial.geopoint.document.GeoPointField.TermEncoding;
import org.apache.lucene.spatial.util.GeoRect;
import org.apache.lucene.spatial.util.GeoRelationUtils;
import org.apache.lucene.util.SloppyMath;
/** Package private implementation for the public facing GeoPointDistanceQuery delegate class.
@ -54,15 +53,28 @@ final class GeoPointDistanceQueryImpl extends GeoPointInBBoxQueryImpl {
@Override
protected boolean cellCrosses(final double minLat, final double maxLat, final double minLon, final double maxLon) {
return GeoRelationUtils.rectCrossesCircle(minLat, maxLat, minLon, maxLon,
distanceQuery.centerLat, centerLon, distanceQuery.radiusMeters, true);
if (maxLat < GeoPointDistanceQueryImpl.this.minLat ||
maxLon < GeoPointDistanceQueryImpl.this.minLon ||
minLat > GeoPointDistanceQueryImpl.this.maxLat ||
minLon > GeoPointDistanceQueryImpl.this.maxLon) {
return false;
} else {
return true;
}
}
@Override
protected boolean cellWithin(final double minLat, final double maxLat, final double minLon, final double maxLon) {
return GeoRelationUtils.rectWithinCircle(minLat, maxLat, minLon, maxLon,
distanceQuery.centerLat, centerLon,
distanceQuery.radiusMeters, true);
// TODO: we call cellCrosses because of how the termsEnum logic works, helps us avoid some haversin() calls here.
if (cellCrosses(minLat, maxLat, minLon, maxLon) && maxLon - centerLon < 90 && centerLon - minLon < 90 &&
SloppyMath.haversinMeters(distanceQuery.centerLat, centerLon, minLat, minLon) <= distanceQuery.radiusMeters &&
SloppyMath.haversinMeters(distanceQuery.centerLat, centerLon, minLat, maxLon) <= distanceQuery.radiusMeters &&
SloppyMath.haversinMeters(distanceQuery.centerLat, centerLon, maxLat, minLon) <= distanceQuery.radiusMeters &&
SloppyMath.haversinMeters(distanceQuery.centerLat, centerLon, maxLat, maxLon) <= distanceQuery.radiusMeters) {
// we are fully enclosed, collect everything within this subtree
return true;
}
return false;
}
@Override

View File

@ -18,16 +18,12 @@ package org.apache.lucene.spatial.util;
import org.apache.lucene.util.SloppyMath;
import static org.apache.lucene.util.SloppyMath.TO_RADIANS;
/**
* Reusable geo-spatial distance utility methods.
*
* @lucene.experimental
*/
public class GeoDistanceUtils {
/** error threshold for point-distance queries (in percent) NOTE: Guideline from USGS is 0.005 **/
public static final double DISTANCE_PCT_ERR = 0.005;
// No instance:
private GeoDistanceUtils() {
@ -58,7 +54,7 @@ public class GeoDistanceUtils {
lat = StrictMath.toRadians(lat);
// get the diameter at the latitude
final double diameter = 2 * GeoProjectionUtils.SEMIMAJOR_AXIS;
final double diameter = 2 * GeoUtils.SEMIMAJOR_AXIS;
// compute inverse haversine
double a = StrictMath.sin(distance/diameter);
@ -125,7 +121,7 @@ public class GeoDistanceUtils {
*/
public static double distanceToDegreesLat(double lat, double distance) {
// get the diameter at the latitude
final double diameter = 2 * GeoProjectionUtils.SEMIMAJOR_AXIS;
final double diameter = 2 * GeoUtils.SEMIMAJOR_AXIS;
// compute inverse haversine
double a = StrictMath.sin(distance/diameter);

View File

@ -84,12 +84,6 @@ public final class GeoEncodingUtils {
return (val / LAT_SCALE) + MIN_LAT_INCL;
}
/** Compare two position values within a {@link GeoEncodingUtils#TOLERANCE} factor */
public static double compare(final double v1, final double v2) {
final double delta = v1-v2;
return Math.abs(delta) <= TOLERANCE ? 0 : delta;
}
/** Convert a geocoded morton long into a prefix coded geo term */
public static void geoCodedToPrefixCoded(long hash, int shift, BytesRefBuilder bytes) {
geoCodedToPrefixCodedBytes(hash, shift, bytes);

View File

@ -1,158 +0,0 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.lucene.spatial.util;
import static java.lang.StrictMath.sqrt;
import static org.apache.lucene.util.SloppyMath.asin;
import static org.apache.lucene.util.SloppyMath.cos;
import static org.apache.lucene.util.SloppyMath.TO_DEGREES;
import static org.apache.lucene.util.SloppyMath.TO_RADIANS;
import static org.apache.lucene.spatial.util.GeoUtils.MAX_LAT_INCL;
import static org.apache.lucene.spatial.util.GeoUtils.MAX_LON_INCL;
import static org.apache.lucene.spatial.util.GeoUtils.MIN_LAT_INCL;
import static org.apache.lucene.spatial.util.GeoUtils.MIN_LON_INCL;
import static org.apache.lucene.spatial.util.GeoUtils.PIO2;
import static org.apache.lucene.spatial.util.GeoUtils.normalizeLat;
import static org.apache.lucene.spatial.util.GeoUtils.normalizeLon;
import static org.apache.lucene.spatial.util.GeoUtils.sloppySin;
import static org.apache.lucene.spatial.util.GeoUtils.sloppyTan;
/**
* Reusable geo-spatial projection utility methods.
*
* @lucene.experimental
*/
public class GeoProjectionUtils {
// WGS84 earth-ellipsoid parameters
/** major (a) axis in meters */
public static final double SEMIMAJOR_AXIS = 6_378_137; // [m]
/** earth flattening factor (f) */
public static final double FLATTENING = 1.0/298.257223563;
/** minor (b) axis in meters */
public static final double SEMIMINOR_AXIS = SEMIMAJOR_AXIS * (1.0 - FLATTENING); //6_356_752.31420; // [m]
/** first eccentricity (e) */
public static final double ECCENTRICITY = sqrt((2.0 - FLATTENING) * FLATTENING);
/** major axis squared (a2) */
public static final double SEMIMAJOR_AXIS2 = SEMIMAJOR_AXIS * SEMIMAJOR_AXIS;
/** minor axis squared (b2) */
public static final double SEMIMINOR_AXIS2 = SEMIMINOR_AXIS * SEMIMINOR_AXIS;
private static final double E2 = (SEMIMAJOR_AXIS2 - SEMIMINOR_AXIS2)/(SEMIMAJOR_AXIS2);
private static final double EP2 = (SEMIMAJOR_AXIS2 - SEMIMINOR_AXIS2)/(SEMIMINOR_AXIS2);
/** min longitude value in radians */
public static final double MIN_LON_RADIANS = TO_RADIANS * MIN_LON_INCL;
/** min latitude value in radians */
public static final double MIN_LAT_RADIANS = TO_RADIANS * MIN_LAT_INCL;
/** max longitude value in radians */
public static final double MAX_LON_RADIANS = TO_RADIANS * MAX_LON_INCL;
/** max latitude value in radians */
public static final double MAX_LAT_RADIANS = TO_RADIANS * MAX_LAT_INCL;
// No instance:
private GeoProjectionUtils() {
}
/**
* Converts from geodesic lat lon alt to geocentric earth-centered earth-fixed
* @param lat geodesic latitude
* @param lon geodesic longitude
* @param alt geodesic altitude
* @param ecf reusable earth-centered earth-fixed result
* @return either a new ecef array or the reusable ecf parameter
*/
public static final double[] llaToECF(double lat, double lon, double alt, double[] ecf) {
lon = TO_RADIANS * lon;
lat = TO_RADIANS * lat;
final double sl = sloppySin(lat);
final double s2 = sl*sl;
final double cl = cos(lat);
if (ecf == null) {
ecf = new double[3];
}
if (lat < -PIO2 && lat > -1.001D * PIO2) {
lat = -PIO2;
} else if (lat > PIO2 && lat < 1.001D * PIO2) {
lat = PIO2;
}
assert (lat >= -PIO2) || (lat <= PIO2);
if (lon > StrictMath.PI) {
lon -= (2*StrictMath.PI);
}
final double rn = SEMIMAJOR_AXIS / StrictMath.sqrt(1.0D - E2 * s2);
ecf[0] = (rn+alt) * cl * cos(lon);
ecf[1] = (rn+alt) * cl * sloppySin(lon);
ecf[2] = ((rn*(1.0-E2))+alt)*sl;
return ecf;
}
/**
* Finds a point along a bearing from a given lat,lon geolocation using great circle arc
*
* @param lat origin latitude in degrees
* @param lon origin longitude in degrees
* @param bearing azimuthal bearing in degrees
* @param dist distance in meters
* @param pt resulting point
* @return the point along a bearing at a given distance in meters
*/
public static final double[] pointFromLonLatBearingGreatCircle(double lat, double lon, double bearing, double dist, double[] pt) {
if (pt == null) {
pt = new double[2];
}
lat *= TO_RADIANS;
lon *= TO_RADIANS;
bearing *= TO_RADIANS;
final double cLat = cos(lat);
final double sLat = sloppySin(lat);
final double sinDoR = sloppySin(dist / GeoProjectionUtils.SEMIMAJOR_AXIS);
final double cosDoR = cos(dist / GeoProjectionUtils.SEMIMAJOR_AXIS);
pt[1] = asin(sLat*cosDoR + cLat * sinDoR * cos(bearing));
pt[0] = TO_DEGREES * (lon + Math.atan2(sloppySin(bearing) * sinDoR * cLat, cosDoR - sLat * sloppySin(pt[1])));
pt[1] *= TO_DEGREES;
return pt;
}
/**
* Finds the bearing (in degrees) between 2 geo points (lat, lon) using great circle arc
* @param lat1 first point latitude in degrees
* @param lon1 first point longitude in degrees
* @param lat2 second point latitude in degrees
* @param lon2 second point longitude in degrees
* @return the bearing (in degrees) between the two provided points
*/
public static double bearingGreatCircle(double lat1, double lon1, double lat2, double lon2) {
double dLon = (lon2 - lon1) * TO_RADIANS;
lat2 *= TO_RADIANS;
lat1 *= TO_RADIANS;
double y = sloppySin(dLon) * cos(lat2);
double x = cos(lat1) * sloppySin(lat2) - sloppySin(lat1) * cos(lat2) * cos(dLon);
return Math.atan2(y, x) * TO_DEGREES;
}
}

View File

@ -16,8 +16,6 @@
*/
package org.apache.lucene.spatial.util;
import org.apache.lucene.util.SloppyMath;
/**
* Reusable geo-relation utility methods
*/
@ -277,252 +275,4 @@ public class GeoRelationUtils {
return !(rectCrossesPolyApprox(rMinLat, rMaxLat, rMinLon, rMaxLon, shapeLats, shapeLons, sMinLat, sMaxLat, sMinLon, sMaxLon)
|| !pointInPolygon(shapeLats, shapeLons, rMinLat, rMinLon));
}
/////////////////////////
// Circle relations
/////////////////////////
private static boolean rectAnyCornersInCircle(final double rMinLat, final double rMaxLat, final double rMinLon,
final double rMaxLon, final double centerLat, final double centerLon,
final double radiusMeters, final boolean approx) {
if (approx == true) {
return rectAnyCornersInCircleSloppy(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters);
}
double w = Math.abs(rMaxLon - rMinLon);
if (w <= 90.0) {
return SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMinLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMinLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMaxLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMaxLon) <= radiusMeters;
}
// partition
w /= 4;
final double p1 = rMinLon + w;
final double p2 = p1 + w;
final double p3 = p2 + w;
return SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMinLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMinLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, p1) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, p1) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, p2) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, p2) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, p3) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, p3) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMaxLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMaxLon) <= radiusMeters;
}
private static boolean rectAnyCornersInCircleSloppy(final double rMinLat, final double rMaxLat, final double rMinLon, final double rMaxLon,
final double centerLat, final double centerLon, final double radiusMeters) {
return SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMinLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMinLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMaxLon) <= radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMaxLon) <= 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 lat/lon 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 rMinLat, final double rMaxLat, final double rMinLon,
final double rMaxLon, final double centerLat, final double centerLon,
final double radiusMeters, final boolean approx) {
if (approx == true) {
return rectAnyCornersOutsideCircleSloppy(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters);
}
// if span is less than 70 degrees we can approximate using distance alone
if (Math.abs(rMaxLon - rMinLon) <= 70.0) {
return SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMinLon) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMinLon) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMaxLon) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMaxLon) > radiusMeters;
}
return rectCrossesOblateCircle(centerLat, centerLon,
radiusMeters,
rMinLat, rMaxLat,
rMinLon, rMaxLon);
}
/**
* Compute whether the rectangle (defined by min/max Lat/Lon) crosses a potentially oblate circle
*
* TODO benchmark for replacing existing rectCrossesCircle.
*/
private static boolean rectCrossesOblateCircle(double centerLat, double centerLon,
double radiusMeters,
double rMinLat, double rMaxLat,
double rMinLon, double rMaxLon) {
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 = SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, p1)) > radiusMeters
|| (d2 = SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, maxLon)) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, p1) > radiusMeters
|| SloppyMath.haversinMeters(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(rMinLat, p1,
GeoProjectionUtils.bearingGreatCircle(rMinLat, p1, rMaxLat, p1), radiusMeters - d1, pt))[1] < rMinLat || pt[1] < rMaxLat
|| (pt = GeoProjectionUtils.pointFromLonLatBearingGreatCircle(rMinLat, maxLon,
GeoProjectionUtils.bearingGreatCircle(rMinLat, maxLon, rMaxLat, maxLon), radiusMeters - d2, pt))[1] < rMinLat || pt[1] < rMaxLat
|| (pt = GeoProjectionUtils.pointFromLonLatBearingGreatCircle(rMinLat, maxLon,
GeoProjectionUtils.bearingGreatCircle(rMinLat, maxLon, rMaxLat, (midLon = p1 + 0.5*(maxLon - p1))),
radiusMeters - SloppyMath.haversinMeters(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 rMinLat, final double rMaxLat, final double rMinLon, final double rMaxLon,
final double centerLat, final double centerLon, final double radiusMeters) {
return SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMinLon) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMinLon) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMaxLat, rMaxLon) > radiusMeters
|| SloppyMath.haversinMeters(centerLat, centerLon, rMinLat, rMaxLon) > radiusMeters;
}
/**
* Computes whether a rectangle is within a circle. Note: approx == true will be faster but less precise and may
* fail on large rectangles
*/
public static boolean rectWithinCircle(final double rMinLat, final double rMaxLat, final double rMinLon, final double rMaxLon,
final double centerLat, final double centerLon, final double radiusMeters,
final boolean approx) {
return rectAnyCornersOutsideCircle(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters, approx) == false;
}
/**
* Computes whether a rectangle crosses a circle. Note: approx == true will be faster but less precise and may
* fail on large rectangles
*
* <p>NOTE: this is basic method that does not handle dateline or pole crossing. Unwrapping must be done before
* calling this method.
*/
public static boolean rectCrossesCircle(final double rMinLat, final double rMaxLat, final double rMinLon, final double rMaxLon,
final double centerLat, final double centerLon, final double radiusMeters,
final boolean approx) {
if (approx == true) {
if (rectAnyCornersInCircle(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters, approx)) {
return true;
}
} else {
if (rectAnyCornersInCircle(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters, approx) &&
rectAnyCornersOutsideCircle(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters, approx)) {
return true;
}
}
if (isClosestPointOnRectWithinRange(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, radiusMeters, approx)) {
return true;
}
return false;
}
private static boolean isClosestPointOnRectWithinRange(final double rMinLat, final double rMaxLat,
final double rMinLon, final double rMaxLon,
final double centerLat, final double centerLon,
final double radiusMeters,
final boolean approx) {
double[] closestPt = {0, 0};
GeoDistanceUtils.closestPointOnBBox(rMinLat, rMaxLat, rMinLon, rMaxLon, centerLat, centerLon, closestPt);
boolean haverShortCut = SloppyMath.haversinMeters(centerLat, centerLon, closestPt[0], closestPt[1]) <= radiusMeters;
if (approx == true || haverShortCut == true) {
return haverShortCut;
}
double lon1 = rMinLon;
double lon2 = rMaxLon;
double lat1 = rMinLat;
double lat2 = rMaxLat;
if (closestPt[1] == rMinLon || closestPt[1] == rMaxLon) {
lon1 = closestPt[1];
lon2 = lon1;
} else if (closestPt[0] == rMinLat || closestPt[0] == rMaxLat) {
lat1 = closestPt[0];
lat2 = lat1;
}
return lineCrossesSphere(lat1, lon1, 0,
lat2, lon2, 0,
centerLat, centerLon, 0,
radiusMeters);
}
/**
* Computes whether or a 3dimensional line segment intersects or crosses a sphere
*
* @param lon1 longitudinal location of the line segment start point (in degrees)
* @param lat1 latitudinal location of the line segment start point (in degrees)
* @param alt1 altitude of the line segment start point (in degrees)
* @param lon2 longitudinal location of the line segment end point (in degrees)
* @param lat2 latitudinal location of the line segment end point (in degrees)
* @param alt2 altitude of the line segment end point (in degrees)
* @param centerLon longitudinal location of center search point (in degrees)
* @param centerLat latitudinal location of center search point (in degrees)
* @param centerAlt altitude of the center point (in meters)
* @param radiusMeters search sphere radius (in meters)
* @return whether the provided line segment is a secant of the
*/
private static boolean lineCrossesSphere(double lat1, double lon1, double alt1,
double lat2, double lon2, double alt2,
double centerLat, double centerLon, double centerAlt,
double radiusMeters) {
// convert to cartesian 3d (in meters)
double[] ecf1 = GeoProjectionUtils.llaToECF(lat1, lon1, alt1, null);
double[] ecf2 = GeoProjectionUtils.llaToECF(lat2, lon2, alt2, null);
double[] cntr = GeoProjectionUtils.llaToECF(centerLat, centerLon, centerAlt, null);
// convert radius from arc radius to cartesian radius
double[] oneEighty = GeoProjectionUtils.pointFromLonLatBearingGreatCircle(centerLat, centerLon, 180.0d, radiusMeters, new double[3]);
GeoProjectionUtils.llaToECF(oneEighty[1], oneEighty[0], 0, oneEighty);
radiusMeters = GeoDistanceUtils.linearDistance(oneEighty, cntr);// Math.sqrt(oneEighty[0]*cntr[0] + oneEighty[1]*cntr[1] + oneEighty[2]*cntr[2]);
final double dX = ecf2[0] - ecf1[0];
final double dY = ecf2[1] - ecf1[1];
final double dZ = ecf2[2] - ecf1[2];
final double fX = ecf1[0] - cntr[0];
final double fY = ecf1[1] - cntr[1];
final double fZ = ecf1[2] - cntr[2];
final double a = dX*dX + dY*dY + dZ*dZ;
final double b = 2 * (fX*dX + fY*dY + fZ*dZ);
final double c = (fX*fX + fY*fY + fZ*fZ) - (radiusMeters*radiusMeters);
double discrim = (b*b)-(4*a*c);
if (discrim < 0) {
return false;
}
discrim = StrictMath.sqrt(discrim);
final double a2 = 2*a;
final double t1 = (-b - discrim)/a2;
final double t2 = (-b + discrim)/a2;
if ( (t1 < 0 || t1 > 1) ) {
return !(t2 < 0 || t2 > 1);
}
return true;
}
}

View File

@ -16,10 +16,6 @@
*/
package org.apache.lucene.spatial.util;
import java.util.ArrayList;
import org.apache.lucene.util.SloppyMath;
import static java.lang.Math.max;
import static java.lang.Math.min;
import static java.lang.Math.PI;
@ -30,11 +26,6 @@ import static org.apache.lucene.util.SloppyMath.cos;
import static org.apache.lucene.util.SloppyMath.TO_DEGREES;
import static org.apache.lucene.util.SloppyMath.TO_RADIANS;
import static org.apache.lucene.spatial.util.GeoEncodingUtils.TOLERANCE;
import static org.apache.lucene.spatial.util.GeoProjectionUtils.MAX_LAT_RADIANS;
import static org.apache.lucene.spatial.util.GeoProjectionUtils.MAX_LON_RADIANS;
import static org.apache.lucene.spatial.util.GeoProjectionUtils.MIN_LAT_RADIANS;
import static org.apache.lucene.spatial.util.GeoProjectionUtils.MIN_LON_RADIANS;
import static org.apache.lucene.spatial.util.GeoProjectionUtils.SEMIMAJOR_AXIS;
/**
* Basic reusable geo-spatial utility methods
@ -53,6 +44,19 @@ public final class GeoUtils {
/** Maximum latitude value. */
public static final double MAX_LAT_INCL = 90.0D;
/** min longitude value in radians */
public static final double MIN_LON_RADIANS = TO_RADIANS * MIN_LON_INCL;
/** min latitude value in radians */
public static final double MIN_LAT_RADIANS = TO_RADIANS * MIN_LAT_INCL;
/** max longitude value in radians */
public static final double MAX_LON_RADIANS = TO_RADIANS * MAX_LON_INCL;
/** max latitude value in radians */
public static final double MAX_LAT_RADIANS = TO_RADIANS * MAX_LAT_INCL;
// WGS84 earth-ellipsoid parameters
/** major (a) axis in meters */
public static final double SEMIMAJOR_AXIS = 6_378_137; // [m]
// No instance:
private GeoUtils() {
@ -153,7 +157,7 @@ public final class GeoUtils {
// some sloppyish stuff, do we really need this to be done in a sloppy way?
// unless it is performance sensitive, we should try to remove.
static final double PIO2 = Math.PI / 2D;
private static final double PIO2 = Math.PI / 2D;
/**
* Returns the trigonometric sine of an angle converted as a cos operation.
@ -169,26 +173,8 @@ public final class GeoUtils {
* @see Math#sin(double)
*/
// TODO: deprecate/remove this? at least its no longer public.
static double sloppySin(double a) {
private static double sloppySin(double a) {
return cos(a - PIO2);
}
/**
* Returns the trigonometric tangent of an angle converted in terms of a sin/cos operation.
* <p>
* Note that this is probably not quite right (untested)
* <p>
* Special cases:
* <ul>
* <li>If the argument is {@code NaN} or an infinity, then the result is {@code NaN}.
* </ul>
* @param a an angle, in radians.
* @return the tangent of the argument.
* @see Math#sin(double) aand Math#cos(double)
*/
// TODO: deprecate/remove this? at least its no longer public.
static double sloppyTan(double a) {
return sloppySin(a) / cos(a);
}
}

View File

@ -22,23 +22,30 @@ import org.apache.lucene.spatial.util.GeoEncodingUtils;
import org.apache.lucene.spatial.geopoint.document.GeoPointField;
import org.apache.lucene.spatial.geopoint.document.GeoPointField.TermEncoding;
import org.apache.lucene.spatial.util.BaseGeoPointTestCase;
import org.apache.lucene.spatial.util.GeoDistanceUtils;
import org.apache.lucene.spatial.util.GeoRect;
import org.apache.lucene.spatial.util.GeoRelationUtils;
import org.apache.lucene.util.SloppyMath;
import static org.apache.lucene.spatial.util.GeoDistanceUtils.DISTANCE_PCT_ERR;
/**
* random testing for GeoPoint query logic
*
* @lucene.experimental
*/
public class TestGeoPointQuery extends BaseGeoPointTestCase {
@Override
protected boolean forceSmall() {
return false;
protected double maxRadius(double latitude, double longitude) {
// TODO: clean this up
return GeoDistanceUtils.maxRadialDistanceMeters(latitude, longitude);
}
@Override
protected double quantizeLat(double lat) {
return GeoEncodingUtils.mortonUnhashLat(GeoEncodingUtils.mortonHash(lat, 0));
}
@Override
protected double quantizeLon(double lon) {
return GeoEncodingUtils.mortonUnhashLon(GeoEncodingUtils.mortonHash(0, lon));
}
@Override
@ -58,7 +65,9 @@ public class TestGeoPointQuery extends BaseGeoPointTestCase {
@Override
protected Query newDistanceRangeQuery(String field, double centerLat, double centerLon, double minRadiusMeters, double radiusMeters) {
return new GeoPointDistanceRangeQuery(field, TermEncoding.PREFIX, centerLat, centerLon, minRadiusMeters, radiusMeters);
// LUCENE-7126: currently not valid for multi-valued documents, because it rewrites to a BooleanQuery!
// return new GeoPointDistanceRangeQuery(field, TermEncoding.PREFIX, centerLat, centerLon, minRadiusMeters, radiusMeters);
return null;
}
@Override
@ -66,63 +75,4 @@ public class TestGeoPointQuery extends BaseGeoPointTestCase {
return new GeoPointInPolygonQuery(field, TermEncoding.PREFIX, lats, lons);
}
@Override
protected Boolean rectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
if (GeoEncodingUtils.compare(pointLon, rect.minLon) == 0.0 ||
GeoEncodingUtils.compare(pointLon, rect.maxLon) == 0.0 ||
GeoEncodingUtils.compare(pointLat, rect.minLat) == 0.0 ||
GeoEncodingUtils.compare(pointLat, rect.maxLat) == 0.0) {
// Point is very close to rect boundary
return null;
}
if (rect.minLon < rect.maxLon) {
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, rect.maxLon);
} else {
// Rect crosses dateline:
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, -180.0, rect.maxLon)
|| GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, 180.0);
}
}
@Override
protected Boolean polyRectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
return rectContainsPoint(rect, pointLat, pointLon);
}
@Override
protected Boolean circleContainsPoint(double centerLat, double centerLon, double radiusMeters, double pointLat, double pointLon) {
if (radiusQueryCanBeWrong(centerLat, centerLon, pointLon, pointLat, radiusMeters)) {
return null;
} else {
return SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon) <= radiusMeters;
}
}
@Override
protected Boolean distanceRangeContainsPoint(double centerLat, double centerLon, double minRadiusMeters, double radiusMeters, double pointLat, double pointLon) {
if (radiusQueryCanBeWrong(centerLat, centerLon, pointLon, pointLat, minRadiusMeters)
|| radiusQueryCanBeWrong(centerLat, centerLon, pointLon, pointLat, radiusMeters)) {
return null;
} else {
final double d = SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon);
return d >= minRadiusMeters && d <= radiusMeters;
}
}
private static boolean radiusQueryCanBeWrong(double centerLat, double centerLon, double ptLon, double ptLat,
final double radius) {
final long hashedCntr = GeoEncodingUtils.mortonHash(centerLat, centerLon);
centerLon = GeoEncodingUtils.mortonUnhashLon(hashedCntr);
centerLat = GeoEncodingUtils.mortonUnhashLat(hashedCntr);
final long hashedPt = GeoEncodingUtils.mortonHash(ptLat, ptLon);
ptLon = GeoEncodingUtils.mortonUnhashLon(hashedPt);
ptLat = GeoEncodingUtils.mortonUnhashLat(hashedPt);
double ptDistance = SloppyMath.haversinMeters(centerLat, centerLon, ptLat, ptLon);
double delta = StrictMath.abs(ptDistance - radius);
// if its within the distance error then it can be wrong
return delta < (ptDistance*DISTANCE_PCT_ERR);
}
}

View File

@ -22,11 +22,8 @@ import org.apache.lucene.spatial.util.GeoEncodingUtils;
import org.apache.lucene.spatial.geopoint.document.GeoPointField;
import org.apache.lucene.spatial.geopoint.document.GeoPointField.TermEncoding;
import org.apache.lucene.spatial.util.BaseGeoPointTestCase;
import org.apache.lucene.spatial.util.GeoDistanceUtils;
import org.apache.lucene.spatial.util.GeoRect;
import org.apache.lucene.spatial.util.GeoRelationUtils;
import org.apache.lucene.util.SloppyMath;
import static org.apache.lucene.spatial.util.GeoDistanceUtils.DISTANCE_PCT_ERR;
/**
* random testing for GeoPoint query logic (with deprecated numeric encoding)
@ -36,8 +33,19 @@ import static org.apache.lucene.spatial.util.GeoDistanceUtils.DISTANCE_PCT_ERR;
public class TestLegacyGeoPointQuery extends BaseGeoPointTestCase {
@Override
protected boolean forceSmall() {
return false;
protected double maxRadius(double latitude, double longitude) {
// TODO: clean this up
return GeoDistanceUtils.maxRadialDistanceMeters(latitude, longitude);
}
@Override
protected double quantizeLat(double lat) {
return GeoEncodingUtils.mortonUnhashLat(GeoEncodingUtils.mortonHash(lat, 0));
}
@Override
protected double quantizeLon(double lon) {
return GeoEncodingUtils.mortonUnhashLon(GeoEncodingUtils.mortonHash(0, lon));
}
@Override
@ -57,7 +65,9 @@ public class TestLegacyGeoPointQuery extends BaseGeoPointTestCase {
@Override
protected Query newDistanceRangeQuery(String field, double centerLat, double centerLon, double minRadiusMeters, double radiusMeters) {
return new GeoPointDistanceRangeQuery(field, TermEncoding.NUMERIC, centerLat, centerLon, minRadiusMeters, radiusMeters);
// LUCENE-7126: currently not valid for multi-valued documents, because it rewrites to a BooleanQuery!
// return new GeoPointDistanceRangeQuery(field, TermEncoding.NUMERIC, centerLat, centerLon, minRadiusMeters, radiusMeters);
return null;
}
@Override
@ -65,63 +75,14 @@ public class TestLegacyGeoPointQuery extends BaseGeoPointTestCase {
return new GeoPointInPolygonQuery(field, TermEncoding.NUMERIC, lats, lons);
}
// legacy encoding is too slow somehow for this random test, its not up to the task.
@Override
protected Boolean rectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
if (GeoEncodingUtils.compare(pointLon, rect.minLon) == 0.0 ||
GeoEncodingUtils.compare(pointLon, rect.maxLon) == 0.0 ||
GeoEncodingUtils.compare(pointLat, rect.minLat) == 0.0 ||
GeoEncodingUtils.compare(pointLat, rect.maxLat) == 0.0) {
// Point is very close to rect boundary
return null;
}
if (rect.minLon < rect.maxLon) {
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, rect.maxLon);
} else {
// Rect crosses dateline:
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, -180.0, rect.maxLon)
|| GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, 180.0);
}
public void testRandomDistance() throws Exception {
assumeTrue("legacy encoding is too slow/hangs on this test", false);
}
@Override
protected Boolean polyRectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
return rectContainsPoint(rect, pointLat, pointLon);
}
@Override
protected Boolean circleContainsPoint(double centerLat, double centerLon, double radiusMeters, double pointLat, double pointLon) {
if (radiusQueryCanBeWrong(centerLat, centerLon, pointLon, pointLat, radiusMeters)) {
return null;
} else {
return SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon) <= radiusMeters;
}
}
@Override
protected Boolean distanceRangeContainsPoint(double centerLat, double centerLon, double minRadiusMeters, double radiusMeters, double pointLat, double pointLon) {
if (radiusQueryCanBeWrong(centerLat, centerLon, pointLon, pointLat, minRadiusMeters)
|| radiusQueryCanBeWrong(centerLat, centerLon, pointLon, pointLat, radiusMeters)) {
return null;
} else {
final double d = SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon);
return d >= minRadiusMeters && d <= radiusMeters;
}
}
private static boolean radiusQueryCanBeWrong(double centerLat, double centerLon, double ptLon, double ptLat,
final double radius) {
final long hashedCntr = GeoEncodingUtils.mortonHash(centerLat, centerLon);
centerLon = GeoEncodingUtils.mortonUnhashLon(hashedCntr);
centerLat = GeoEncodingUtils.mortonUnhashLat(hashedCntr);
final long hashedPt = GeoEncodingUtils.mortonHash(ptLat, ptLon);
ptLon = GeoEncodingUtils.mortonUnhashLon(hashedPt);
ptLat = GeoEncodingUtils.mortonUnhashLat(hashedPt);
double ptDistance = SloppyMath.haversinMeters(centerLat, centerLon, ptLat, ptLon);
double delta = StrictMath.abs(ptDistance - radius);
// if its within the distance error then it can be wrong
return delta < (ptDistance*DISTANCE_PCT_ERR);
public void testRandomDistanceHuge() throws Exception {
assumeTrue("legacy encoding is too slow/hangs on this test", false);
}
}

View File

@ -21,6 +21,7 @@ import java.text.DecimalFormat;
import java.text.DecimalFormatSymbols;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.BitSet;
import java.util.HashSet;
import java.util.List;
import java.util.Locale;
@ -28,9 +29,16 @@ import java.util.Set;
import java.util.concurrent.CountDownLatch;
import java.util.concurrent.atomic.AtomicBoolean;
import org.apache.lucene.codecs.FilterCodec;
import org.apache.lucene.codecs.PointsFormat;
import org.apache.lucene.codecs.PointsReader;
import org.apache.lucene.codecs.PointsWriter;
import org.apache.lucene.codecs.lucene60.Lucene60PointsReader;
import org.apache.lucene.codecs.lucene60.Lucene60PointsWriter;
import org.apache.lucene.document.Document;
import org.apache.lucene.document.Field;
import org.apache.lucene.document.NumericDocValuesField;
import org.apache.lucene.document.StoredField;
import org.apache.lucene.index.DirectoryReader;
import org.apache.lucene.index.IndexReader;
import org.apache.lucene.index.IndexWriter;
@ -39,16 +47,22 @@ import org.apache.lucene.index.LeafReaderContext;
import org.apache.lucene.index.MultiDocValues;
import org.apache.lucene.index.NumericDocValues;
import org.apache.lucene.index.RandomIndexWriter;
import org.apache.lucene.index.SegmentReadState;
import org.apache.lucene.index.SegmentWriteState;
import org.apache.lucene.index.Term;
import org.apache.lucene.search.IndexSearcher;
import org.apache.lucene.search.Query;
import org.apache.lucene.search.ScoreDoc;
import org.apache.lucene.search.SimpleCollector;
import org.apache.lucene.search.Sort;
import org.apache.lucene.search.TopDocs;
import org.apache.lucene.store.Directory;
import org.apache.lucene.util.FixedBitSet;
import org.apache.lucene.util.IOUtils;
import org.apache.lucene.util.LuceneTestCase;
import org.apache.lucene.util.SloppyMath;
import org.apache.lucene.util.TestUtil;
import org.apache.lucene.util.bkd.BKDWriter;
import org.junit.BeforeClass;
// TODO: cutover TestGeoUtils too?
@ -72,17 +86,8 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
originLat = GeoUtils.normalizeLat(GeoUtils.MIN_LAT_INCL + latRange + (GeoUtils.MAX_LAT_INCL - GeoUtils.MIN_LAT_INCL - 2 * latRange) * random().nextDouble());
}
/** Return true when testing on a non-small region may be too slow (GeoPoint*Query) */
protected boolean forceSmall() {
return false;
}
// A particularly tricky adversary for BKD tree:
public void testSamePointManyTimes() throws Exception {
// For GeoPointQuery, only run this test nightly:
assumeTrue("GeoPoint*Query is too slow otherwise", TEST_NIGHTLY || forceSmall() == false);
int numPoints = atLeast(1000);
boolean small = random().nextBoolean();
@ -100,12 +105,8 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
}
public void testAllLatEqual() throws Exception {
// For GeoPointQuery, only run this test nightly:
assumeTrue("GeoPoint*Query is too slow otherwise", TEST_NIGHTLY || forceSmall() == false);
int numPoints = atLeast(10000);
boolean small = forceSmall() || random().nextBoolean();
boolean small = random().nextBoolean();
double lat = randomLat(small);
double[] lats = new double[numPoints];
double[] lons = new double[numPoints];
@ -151,12 +152,8 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
}
public void testAllLonEqual() throws Exception {
// For GeoPointQuery, only run this test nightly:
assumeTrue("GeoPoint*Query is too slow otherwise", TEST_NIGHTLY || forceSmall() == false);
int numPoints = atLeast(10000);
boolean small = forceSmall() || random().nextBoolean();
boolean small = random().nextBoolean();
double theLon = randomLon(small);
double[] lats = new double[numPoints];
double[] lons = new double[numPoints];
@ -204,10 +201,6 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
}
public void testMultiValued() throws Exception {
// For GeoPointQuery, only run this test nightly:
assumeTrue("GeoPoint*Query is too slow otherwise", TEST_NIGHTLY || forceSmall() == false);
int numPoints = atLeast(10000);
// Every doc has 2 points:
double[] lats = new double[2*numPoints];
@ -289,19 +282,10 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
double latDoc2 = lats[2*docID+1];
double lonDoc2 = lons[2*docID+1];
Boolean result1 = rectContainsPoint(rect, latDoc1, lonDoc1);
if (result1 == null) {
// borderline case: cannot test
continue;
}
boolean result1 = rectContainsPoint(rect, latDoc1, lonDoc1);
boolean result2 = rectContainsPoint(rect, latDoc2, lonDoc2);
Boolean result2 = rectContainsPoint(rect, latDoc2, lonDoc2);
if (result2 == null) {
// borderline case: cannot test
continue;
}
boolean expected = result1 == Boolean.TRUE || result2 == Boolean.TRUE;
boolean expected = result1 || result2;
if (hits.get(docID) != expected) {
String id = s.doc(docID).get("id");
@ -447,6 +431,10 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
protected double quantizeLon(double lon) {
return lon;
}
protected double maxRadius(double latitude, double longitude) {
return 50000000D; // bigger than earth, shouldnt matter
}
protected GeoRect randomRect(boolean small, boolean canCrossDateLine) {
double lat0 = randomLat(small);
@ -482,16 +470,44 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
protected abstract Query newPolygonQuery(String field, double[] lats, double[] lons);
/** Returns null if it's borderline case */
protected abstract Boolean rectContainsPoint(GeoRect rect, double pointLat, double pointLon);
static final boolean rectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
assert Double.isNaN(pointLat) == false;
/** Returns null if it's borderline case */
protected abstract Boolean polyRectContainsPoint(GeoRect rect, double pointLat, double pointLon);
if (rect.minLon < rect.maxLon) {
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, rect.maxLon);
} else {
// Rect crosses dateline:
return GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, -180.0, rect.maxLon)
|| GeoRelationUtils.pointInRectPrecise(pointLat, pointLon, rect.minLat, rect.maxLat, rect.minLon, 180.0);
}
}
static final boolean polyRectContainsPoint(GeoRect rect, double pointLat, double pointLon) {
// TODO write better random polygon tests
// note: logic must be slightly different than rectContainsPoint, to satisfy
// insideness for cases exactly on boundaries.
assert Double.isNaN(pointLat) == false;
assert rect.crossesDateline() == false;
double polyLats[] = new double[] { rect.minLat, rect.maxLat, rect.maxLat, rect.minLat, rect.minLat };
double polyLons[] = new double[] { rect.minLon, rect.minLon, rect.maxLon, rect.maxLon, rect.minLon };
/** Returns null if it's borderline case */
protected abstract Boolean circleContainsPoint(double centerLat, double centerLon, double radiusMeters, double pointLat, double pointLon);
// TODO: separately test this method is 100% correct, here treat it like a black box (like haversin)
return GeoRelationUtils.pointInPolygon(polyLats, polyLons, pointLat, pointLon);
}
protected abstract Boolean distanceRangeContainsPoint(double centerLat, double centerLon, double minRadiusMeters, double radiusMeters, double pointLat, double pointLon);
static final boolean circleContainsPoint(double centerLat, double centerLon, double radiusMeters, double pointLat, double pointLon) {
double distanceMeters = SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon);
boolean result = distanceMeters <= radiusMeters;
//System.out.println(" shouldMatch? centerLon=" + centerLon + " centerLat=" + centerLat + " pointLon=" + pointLon + " pointLat=" + pointLat + " result=" + result + " distanceMeters=" + (distanceKM * 1000));
return result;
}
static final boolean distanceRangeContainsPoint(double centerLat, double centerLon, double minRadiusMeters, double radiusMeters, double pointLat, double pointLon) {
final double d = SloppyMath.haversinMeters(centerLat, centerLon, pointLat, pointLon);
return d >= minRadiusMeters && d <= radiusMeters;
}
private static abstract class VerifyHits {
@ -525,7 +541,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
for(int docID=0;docID<maxDoc;docID++) {
int id = (int) docIDToID.get(docID);
Boolean expected;
boolean expected;
if (deleted.contains(id)) {
expected = false;
} else if (Double.isNaN(lats[id])) {
@ -534,8 +550,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
expected = shouldMatch(lats[id], lons[id]);
}
// null means it's a borderline case which is allowed to be wrong:
if (expected != null && hits.get(docID) != expected) {
if (hits.get(docID) != expected) {
// Print only one failed hit; add a true || in here to see all failures:
if (failFast == false || failed.getAndSet(true) == false) {
@ -568,7 +583,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
/** Return true if we definitely should match, false if we definitely
* should not match, and null if it's a borderline case which might
* go either way. */
protected abstract Boolean shouldMatch(double lat, double lon);
protected abstract boolean shouldMatch(double lat, double lon);
protected abstract void describe(int docID, double lat, double lon);
}
@ -666,7 +681,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
verifyHits = new VerifyHits() {
@Override
protected Boolean shouldMatch(double pointLat, double pointLon) {
protected boolean shouldMatch(double pointLat, double pointLon) {
return rectContainsPoint(rect, pointLat, pointLon);
}
@Override
@ -688,7 +703,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
radiusMeters = random().nextDouble() * 333000 + 1.0;
} else {
// So the query can cover at most 50% of the earth's surface:
radiusMeters = random().nextDouble() * GeoProjectionUtils.SEMIMAJOR_AXIS * Math.PI / 2.0 + 1.0;
radiusMeters = random().nextDouble() * GeoUtils.SEMIMAJOR_AXIS * Math.PI / 2.0 + 1.0;
}
// generate a random minimum radius between 1% and 95% the max radius
@ -715,7 +730,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
verifyHits = new VerifyHits() {
@Override
protected Boolean shouldMatch(double pointLat, double pointLon) {
protected boolean shouldMatch(double pointLat, double pointLon) {
if (rangeQuery == false) {
return circleContainsPoint(centerLat, centerLon, radiusMeters, pointLat, pointLon);
} else {
@ -755,7 +770,7 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
verifyHits = new VerifyHits() {
@Override
protected Boolean shouldMatch(double pointLat, double pointLon) {
protected boolean shouldMatch(double pointLat, double pointLon) {
return polyRectContainsPoint(bbox, pointLat, pointLon);
}
@ -834,4 +849,95 @@ public abstract class BaseGeoPointTestCase extends LuceneTestCase {
w.close();
dir.close();
}
/** Run a few iterations with just 10 docs, hopefully easy to debug */
public void testRandomDistance() throws Exception {
for (int iters = 0; iters < 100; iters++) {
doRandomDistanceTest(10, 100);
}
}
/** Runs with thousands of docs */
@Nightly
public void testRandomDistanceHuge() throws Exception {
for (int iters = 0; iters < 10; iters++) {
doRandomDistanceTest(2000, 100);
}
}
private void doRandomDistanceTest(int numDocs, int numQueries) throws IOException {
Directory dir = newDirectory();
IndexWriterConfig iwc = newIndexWriterConfig();
int pointsInLeaf = 2 + random().nextInt(4);
iwc.setCodec(new FilterCodec("Lucene60", TestUtil.getDefaultCodec()) {
@Override
public PointsFormat pointsFormat() {
return new PointsFormat() {
@Override
public PointsWriter fieldsWriter(SegmentWriteState writeState) throws IOException {
return new Lucene60PointsWriter(writeState, pointsInLeaf, BKDWriter.DEFAULT_MAX_MB_SORT_IN_HEAP);
}
@Override
public PointsReader fieldsReader(SegmentReadState readState) throws IOException {
return new Lucene60PointsReader(readState);
}
};
}
});
RandomIndexWriter writer = new RandomIndexWriter(random(), dir, iwc);
for (int i = 0; i < numDocs; i++) {
double latRaw = -90 + 180.0 * random().nextDouble();
double lonRaw = -180 + 360.0 * random().nextDouble();
// pre-normalize up front, so we can just use quantized value for testing and do simple exact comparisons
double lat = quantizeLat(latRaw);
double lon = quantizeLon(lonRaw);
Document doc = new Document();
addPointToDoc("field", doc, lat, lon);
doc.add(new StoredField("lat", lat));
doc.add(new StoredField("lon", lon));
writer.addDocument(doc);
}
IndexReader reader = writer.getReader();
IndexSearcher searcher = newSearcher(reader);
for (int i = 0; i < numQueries; i++) {
double lat = -90 + 180.0 * random().nextDouble();
double lon = -180 + 360.0 * random().nextDouble();
double radius = maxRadius(lat, lon) * random().nextDouble();
BitSet expected = new BitSet();
for (int doc = 0; doc < reader.maxDoc(); doc++) {
double docLatitude = reader.document(doc).getField("lat").numericValue().doubleValue();
double docLongitude = reader.document(doc).getField("lon").numericValue().doubleValue();
double distance = SloppyMath.haversinMeters(lat, lon, docLatitude, docLongitude);
if (distance <= radius) {
expected.set(doc);
}
}
TopDocs topDocs = searcher.search(newDistanceQuery("field", lat, lon, radius), reader.maxDoc(), Sort.INDEXORDER);
BitSet actual = new BitSet();
for (ScoreDoc doc : topDocs.scoreDocs) {
actual.set(doc.doc);
}
try {
assertEquals(expected, actual);
} catch (AssertionError e) {
System.out.println("center: (" + lat + "," + lon + "), radius=" + radius);
for (int doc = 0; doc < reader.maxDoc(); doc++) {
double docLatitude = reader.document(doc).getField("lat").numericValue().doubleValue();
double docLongitude = reader.document(doc).getField("lon").numericValue().doubleValue();
double distance = SloppyMath.haversinMeters(lat, lon, docLatitude, docLongitude);
System.out.println("" + doc + ": (" + docLatitude + "," + docLongitude + "), distance=" + distance);
}
throw e;
}
}
reader.close();
writer.close();
dir.close();
}
}

View File

@ -16,23 +16,11 @@
*/
package org.apache.lucene.spatial.util;
import java.io.PrintWriter;
import java.io.StringWriter;
import java.util.ArrayList;
import java.util.Collection;
import java.util.HashSet;
import java.util.List;
import java.util.Set;
import org.apache.lucene.util.BytesRefBuilder;
import org.apache.lucene.util.LuceneTestCase;
import org.apache.lucene.util.SloppyMath;
import org.junit.BeforeClass;
import com.carrotsearch.randomizedtesting.generators.RandomInts;
import static org.apache.lucene.spatial.util.GeoDistanceUtils.DISTANCE_PCT_ERR;
/**
* Tests class for methods in GeoUtils
*
@ -77,41 +65,6 @@ public class TestGeoUtils extends LuceneTestCase {
assertEquals(0.0, result[1], 0.0);
}
private static class Cell {
static int nextCellID;
final Cell parent;
final int cellID;
final double minLon, maxLon;
final double minLat, maxLat;
final int splitCount;
public Cell(Cell parent,
double minLon, double minLat,
double maxLon, double maxLat,
int splitCount) {
assert maxLon >= minLon;
assert maxLat >= minLat;
this.parent = parent;
this.minLon = minLon;
this.minLat = minLat;
this.maxLon = maxLon;
this.maxLat = maxLat;
this.cellID = nextCellID++;
this.splitCount = splitCount;
}
/** Returns true if the quantized point lies within this cell, inclusive on all bounds. */
public boolean contains(double lon, double lat) {
return lon >= minLon && lon <= maxLon && lat >= minLat && lat <= maxLat;
}
@Override
public String toString() {
return "cell=" + cellID + (parent == null ? "" : " parentCellID=" + parent.cellID) + " lon: " + minLon + " TO " + maxLon + ", lat: " + minLat + " TO " + maxLat + ", splits: " + splitCount;
}
}
public long scaleLon(final double val) {
return (long) ((val-GeoUtils.MIN_LON_INCL) * LON_SCALE);
}
@ -148,263 +101,6 @@ public class TestGeoUtils extends LuceneTestCase {
return result;
}
private void findMatches(Set<Integer> hits, PrintWriter log, Cell root,
double centerLon, double centerLat, double radiusMeters,
double[] docLons, double[] docLats) {
if (VERBOSE) {
log.println(" root cell: " + root);
}
List<Cell> queue = new ArrayList<>();
queue.add(root);
int recurseDepth = RandomInts.randomIntBetween(random(), 5, 15);
while (queue.size() > 0) {
Cell cell = queue.get(queue.size()-1);
queue.remove(queue.size()-1);
if (VERBOSE) {
log.println(" cycle: " + cell + " queue.size()=" + queue.size());
}
if (random().nextInt(10) == 7 || cell.splitCount > recurseDepth) {
if (VERBOSE) {
log.println(" leaf");
}
// Leaf cell: brute force check all docs that fall within this cell:
for(int docID=0;docID<docLons.length;docID++) {
if (cell.contains(docLons[docID], docLats[docID])) {
double distanceMeters = SloppyMath.haversinMeters(centerLat, centerLon, docLats[docID], docLons[docID]);
if (distanceMeters <= radiusMeters) {
if (VERBOSE) {
log.println(" check doc=" + docID + ": match!");
}
hits.add(docID);
} else {
if (VERBOSE) {
log.println(" check doc=" + docID + ": no match");
}
}
}
}
} else {
if (GeoRelationUtils.rectWithinCircle(cell.minLat, cell.maxLat, cell.minLon, cell.maxLon, centerLat, centerLon, radiusMeters, false)) {
// Query circle fully contains this cell, just addAll:
if (VERBOSE) {
log.println(" circle fully contains cell: now addAll");
}
for(int docID=0;docID<docLons.length;docID++) {
if (cell.contains(docLons[docID], docLats[docID])) {
if (VERBOSE) {
log.println(" addAll doc=" + docID);
}
hits.add(docID);
}
}
continue;
} else if (GeoRelationUtils.rectWithin(root.minLat, root.maxLat, root.minLon, root.maxLon,
cell.minLat, cell.maxLat, cell.minLon, cell.maxLon)) {
// Fall through below to "recurse"
if (VERBOSE) {
log.println(" cell fully contains circle: keep splitting");
}
} else if (GeoRelationUtils.rectCrossesCircle(cell.minLat, cell.maxLat, cell.minLon, cell.maxLon,
centerLat, centerLon, radiusMeters, false)) {
// Fall through below to "recurse"
if (VERBOSE) {
log.println(" cell overlaps circle: keep splitting");
}
} else {
if (VERBOSE) {
log.println(" no overlap: drop this cell");
for(int docID=0;docID<docLons.length;docID++) {
if (cell.contains(docLons[docID], docLats[docID])) {
if (VERBOSE) {
log.println(" skip doc=" + docID);
}
}
}
}
continue;
}
// Randomly split:
if (random().nextBoolean()) {
// Split on lon:
double splitValue = cell.minLon + (cell.maxLon - cell.minLon) * random().nextDouble();
if (VERBOSE) {
log.println(" now split on lon=" + splitValue);
}
Cell cell1 = new Cell(cell,
cell.minLon, cell.minLat,
splitValue, cell.maxLat,
cell.splitCount+1);
Cell cell2 = new Cell(cell,
splitValue, cell.minLat,
cell.maxLon, cell.maxLat,
cell.splitCount+1);
if (VERBOSE) {
log.println(" split cell1: " + cell1);
log.println(" split cell2: " + cell2);
}
queue.add(cell1);
queue.add(cell2);
} else {
// Split on lat:
double splitValue = cell.minLat + (cell.maxLat - cell.minLat) * random().nextDouble();
if (VERBOSE) {
log.println(" now split on lat=" + splitValue);
}
Cell cell1 = new Cell(cell,
cell.minLon, cell.minLat,
cell.maxLon, splitValue,
cell.splitCount+1);
Cell cell2 = new Cell(cell,
cell.minLon, splitValue,
cell.maxLon, cell.maxLat,
cell.splitCount+1);
if (VERBOSE) {
log.println(" split cells:\n " + cell1 + "\n " + cell2);
}
queue.add(cell1);
queue.add(cell2);
}
}
}
}
/** Tests consistency of GeoEncodingUtils.rectWithinCircle, .rectCrossesCircle, .rectWithin and SloppyMath.haversine distance check */
public void testGeoRelations() throws Exception {
int numDocs = atLeast(1000);
boolean useSmallRanges = random().nextBoolean();
if (VERBOSE) {
System.out.println("TEST: " + numDocs + " docs useSmallRanges=" + useSmallRanges);
}
double[] docLons = new double[numDocs];
double[] docLats = new double[numDocs];
for(int docID=0;docID<numDocs;docID++) {
docLons[docID] = randomLon(useSmallRanges);
docLats[docID] = randomLat(useSmallRanges);
if (VERBOSE) {
System.out.println(" doc=" + docID + ": lon=" + docLons[docID] + " lat=" + docLats[docID]);
}
}
int iters = atLeast(10);
iters = atLeast(50);
for(int iter=0;iter<iters;iter++) {
Cell.nextCellID = 0;
double centerLon = randomLon(useSmallRanges);
double centerLat = randomLat(useSmallRanges);
// So the circle covers at most 50% of the earth's surface:
double radiusMeters;
// TODO: large exotic rectangles created by BKD may be inaccurate up to 2 times DISTANCE_PCT_ERR.
// restricting size until LUCENE-6994 can be addressed
if (true || useSmallRanges) {
// Approx 3 degrees lon at the equator:
radiusMeters = random().nextDouble() * 333000;
} else {
radiusMeters = random().nextDouble() * GeoProjectionUtils.SEMIMAJOR_AXIS * Math.PI / 2.0;
}
StringWriter sw = new StringWriter();
PrintWriter log = new PrintWriter(sw, true);
if (VERBOSE) {
log.println("\nTEST: iter=" + iter + " radiusMeters=" + radiusMeters + " centerLon=" + centerLon + " centerLat=" + centerLat);
}
GeoRect bbox = GeoUtils.circleToBBox(centerLat, centerLon, radiusMeters);
Set<Integer> hits = new HashSet<>();
if (bbox.maxLon < bbox.minLon) {
// Crosses dateline
log.println(" circle crosses dateline; first left query");
double unwrappedLon = centerLon;
if (unwrappedLon > bbox.maxLon) {
// unwrap left
unwrappedLon += -360.0D;
}
findMatches(hits, log,
new Cell(null,
-180, bbox.minLat,
bbox.maxLon, bbox.maxLat,
0),
unwrappedLon, centerLat, radiusMeters, docLons, docLats);
log.println(" circle crosses dateline; now right query");
if (unwrappedLon < bbox.maxLon) {
// unwrap right
unwrappedLon += 360.0D;
}
findMatches(hits, log,
new Cell(null,
bbox.minLon, bbox.minLat,
180, bbox.maxLat,
0),
unwrappedLon, centerLat, radiusMeters, docLons, docLats);
} else {
// Start with the root cell that fully contains the shape:
findMatches(hits, log,
new Cell(null,
bbox.minLon, bbox.minLat,
bbox.maxLon, bbox.maxLat,
0),
centerLon, centerLat, radiusMeters,
docLons, docLats);
}
if (VERBOSE) {
log.println(" " + hits.size() + " hits");
}
int failCount = 0;
// Done matching, now verify:
for(int docID=0;docID<numDocs;docID++) {
double distanceMeters = SloppyMath.haversinMeters(centerLat, centerLon, docLats[docID], docLons[docID]);
final Boolean expected;
final double percentError = Math.abs(distanceMeters - radiusMeters) / distanceMeters;
if (percentError <= DISTANCE_PCT_ERR) {
expected = null;
} else {
expected = distanceMeters <= radiusMeters;
}
boolean actual = hits.contains(docID);
if (expected != null && actual != expected) {
if (actual) {
log.println("doc=" + docID + " matched but should not with distance error " + percentError + " on iteration " + iter);
} else {
log.println("doc=" + docID + " did not match but should with distance error " + percentError + " on iteration " + iter);
}
log.println(" lon=" + docLons[docID] + " lat=" + docLats[docID] + " distanceMeters=" + distanceMeters + " vs radiusMeters=" + radiusMeters);
failCount++;
}
}
if (failCount != 0) {
System.out.print(sw.toString());
fail(failCount + " incorrect hits (see above)");
}
}
}
/**
* Tests stability of {@link GeoEncodingUtils#geoCodedToPrefixCoded}
*/
@ -613,8 +309,4 @@ public class TestGeoUtils extends LuceneTestCase {
assertFalse(rect.crossesDateline());
}
}
public void testRectCrossesCircle() throws Exception {
assertTrue(GeoRelationUtils.rectCrossesCircle(-90, 0.0, -180, 180, 0.0, 0.667, 88000.0, false));
}
}