diff --git a/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Plane.java b/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Plane.java index 6cef2d8528f..5917ce7af3c 100755 --- a/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Plane.java +++ b/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Plane.java @@ -1359,8 +1359,6 @@ public class Plane extends Vector { // m * [- 2*A*ab^2*r + 2*A^2*ab^2*r*q + 2*B^2*ab^2*r*q + 2*C^2*c^2*r*q] + // [ab^2 - 2*A*ab^2*q + A^2*ab^2*q^2 + B^2*ab^2*q^2 + C^2*c^2*q^2] = 0 - //System.err.println(" computing X bound"); - // Useful subexpressions for this bound final double q = A*abSquared*k; final double qSquared = q * q; @@ -1400,29 +1398,33 @@ public class Plane extends Vector { assert Math.abs(a * m1 * m1 + b * m1 + c) < MINIMUM_RESOLUTION; final double m2 = (-b - sqrtResult) * commonDenom; assert Math.abs(a * m2 * m2 + b * m2 + c) < MINIMUM_RESOLUTION; - final double l1 = r * m1 + q; - final double l2 = r * m2 + q; - // x = ((1 - l*A) * ab^2 ) / (2 * m) - // y = (-l*B * ab^2) / ( 2 * m) - // z = (-l*C * c^2)/ (2 * m) - final double denom1 = 0.5 / m1; - final double denom2 = 0.5 / m2; - final GeoPoint thePoint1 = new GeoPoint((1.0-l1*A) * abSquared * denom1, -l1*B * abSquared * denom1, -l1*C * cSquared * denom1); - final GeoPoint thePoint2 = new GeoPoint((1.0-l2*A) * abSquared * denom2, -l2*B * abSquared * denom2, -l2*C * cSquared * denom2); - //Math is not quite accurate enough for this - //assert planetModel.pointOnSurface(thePoint1): "Point1: "+thePoint1+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ - // (thePoint1.x*thePoint1.x*planetModel.inverseAb*planetModel.inverseAb + thePoint1.y*thePoint1.y*planetModel.inverseAb*planetModel.inverseAb + thePoint1.z*thePoint1.z*planetModel.inverseC*planetModel.inverseC); - //assert planetModel.pointOnSurface(thePoint2): "Point1: "+thePoint2+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ - // (thePoint2.x*thePoint2.x*planetModel.inverseAb*planetModel.inverseAb + thePoint2.y*thePoint2.y*planetModel.inverseAb*planetModel.inverseAb + thePoint2.z*thePoint2.z*planetModel.inverseC*planetModel.inverseC); - //assert evaluateIsZero(thePoint1): "Evaluation of point1: "+evaluate(thePoint1); - //assert evaluateIsZero(thePoint2): "Evaluation of point2: "+evaluate(thePoint2); - addPoint(boundsInfo, bounds, thePoint1); - addPoint(boundsInfo, bounds, thePoint2); + if (Math.abs(m1) >= MINIMUM_RESOLUTION || Math.abs(m2) >= MINIMUM_RESOLUTION) { + final double l1 = r * m1 + q; + final double l2 = r * m2 + q; + // x = ((1 - l*A) * ab^2 ) / (2 * m) + // y = (-l*B * ab^2) / ( 2 * m) + // z = (-l*C * c^2)/ (2 * m) + final double denom1 = 0.5 / m1; + final double denom2 = 0.5 / m2; + final GeoPoint thePoint1 = new GeoPoint((1.0-l1*A) * abSquared * denom1, -l1*B * abSquared * denom1, -l1*C * cSquared * denom1); + final GeoPoint thePoint2 = new GeoPoint((1.0-l2*A) * abSquared * denom2, -l2*B * abSquared * denom2, -l2*C * cSquared * denom2); + //Math is not quite accurate enough for this + //assert planetModel.pointOnSurface(thePoint1): "Point1: "+thePoint1+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ + // (thePoint1.x*thePoint1.x*planetModel.inverseAb*planetModel.inverseAb + thePoint1.y*thePoint1.y*planetModel.inverseAb*planetModel.inverseAb + thePoint1.z*thePoint1.z*planetModel.inverseC*planetModel.inverseC); + //assert planetModel.pointOnSurface(thePoint2): "Point1: "+thePoint2+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ + // (thePoint2.x*thePoint2.x*planetModel.inverseAb*planetModel.inverseAb + thePoint2.y*thePoint2.y*planetModel.inverseAb*planetModel.inverseAb + thePoint2.z*thePoint2.z*planetModel.inverseC*planetModel.inverseC); + //assert evaluateIsZero(thePoint1): "Evaluation of point1: "+evaluate(thePoint1); + //assert evaluateIsZero(thePoint2): "Evaluation of point2: "+evaluate(thePoint2); + addPoint(boundsInfo, bounds, thePoint1); + addPoint(boundsInfo, bounds, thePoint2); + } else { + // This is a plane of the form A=n B=0 C=0. We can set a bound only by noting the D value. + boundsInfo.addXValue(-D/A); + } } else { // No solutions } } else if (Math.abs(b) > MINIMUM_RESOLUTION_SQUARED) { - //System.err.println("Not x quadratic"); // a = 0, so m = - c / b final double m = -c / b; final double l = r * m + q; @@ -1569,24 +1571,29 @@ public class Plane extends Vector { assert Math.abs(a * m1 * m1 + b * m1 + c) < MINIMUM_RESOLUTION; final double m2 = (-b - sqrtResult) * commonDenom; assert Math.abs(a * m2 * m2 + b * m2 + c) < MINIMUM_RESOLUTION; - final double l1 = r * m1 + q; - final double l2 = r * m2 + q; - // x = (-l*A * ab^2 ) / (2 * m) - // y = ((1.0-l*B) * ab^2) / ( 2 * m) - // z = (-l*C * c^2)/ (2 * m) - final double denom1 = 0.5 / m1; - final double denom2 = 0.5 / m2; - final GeoPoint thePoint1 = new GeoPoint(-l1*A * abSquared * denom1, (1.0-l1*B) * abSquared * denom1, -l1*C * cSquared * denom1); - final GeoPoint thePoint2 = new GeoPoint(-l2*A * abSquared * denom2, (1.0-l2*B) * abSquared * denom2, -l2*C * cSquared * denom2); - //Math is not quite accurate enough for this - //assert planetModel.pointOnSurface(thePoint1): "Point1: "+thePoint1+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ - // (thePoint1.x*thePoint1.x*planetModel.inverseAb*planetModel.inverseAb + thePoint1.y*thePoint1.y*planetModel.inverseAb*planetModel.inverseAb + thePoint1.z*thePoint1.z*planetModel.inverseC*planetModel.inverseC); - //assert planetModel.pointOnSurface(thePoint2): "Point2: "+thePoint2+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ - // (thePoint2.x*thePoint2.x*planetModel.inverseAb*planetModel.inverseAb + thePoint2.y*thePoint2.y*planetModel.inverseAb*planetModel.inverseAb + thePoint2.z*thePoint2.z*planetModel.inverseC*planetModel.inverseC); - //assert evaluateIsZero(thePoint1): "Evaluation of point1: "+evaluate(thePoint1); - //assert evaluateIsZero(thePoint2): "Evaluation of point2: "+evaluate(thePoint2); - addPoint(boundsInfo, bounds, thePoint1); - addPoint(boundsInfo, bounds, thePoint2); + if (Math.abs(m1) >= MINIMUM_RESOLUTION || Math.abs(m2) >= MINIMUM_RESOLUTION) { + final double l1 = r * m1 + q; + final double l2 = r * m2 + q; + // x = (-l*A * ab^2 ) / (2 * m) + // y = ((1.0-l*B) * ab^2) / ( 2 * m) + // z = (-l*C * c^2)/ (2 * m) + final double denom1 = 0.5 / m1; + final double denom2 = 0.5 / m2; + final GeoPoint thePoint1 = new GeoPoint(-l1*A * abSquared * denom1, (1.0-l1*B) * abSquared * denom1, -l1*C * cSquared * denom1); + final GeoPoint thePoint2 = new GeoPoint(-l2*A * abSquared * denom2, (1.0-l2*B) * abSquared * denom2, -l2*C * cSquared * denom2); + //Math is not quite accurate enough for this + //assert planetModel.pointOnSurface(thePoint1): "Point1: "+thePoint1+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ + // (thePoint1.x*thePoint1.x*planetModel.inverseAb*planetModel.inverseAb + thePoint1.y*thePoint1.y*planetModel.inverseAb*planetModel.inverseAb + thePoint1.z*thePoint1.z*planetModel.inverseC*planetModel.inverseC); + //assert planetModel.pointOnSurface(thePoint2): "Point2: "+thePoint2+"; Planetmodel="+planetModel+"; A="+A+" B="+B+" C="+C+" D="+D+" planetfcn="+ + // (thePoint2.x*thePoint2.x*planetModel.inverseAb*planetModel.inverseAb + thePoint2.y*thePoint2.y*planetModel.inverseAb*planetModel.inverseAb + thePoint2.z*thePoint2.z*planetModel.inverseC*planetModel.inverseC); + //assert evaluateIsZero(thePoint1): "Evaluation of point1: "+evaluate(thePoint1); + //assert evaluateIsZero(thePoint2): "Evaluation of point2: "+evaluate(thePoint2); + addPoint(boundsInfo, bounds, thePoint1); + addPoint(boundsInfo, bounds, thePoint2); + } else { + // This is a plane of the form A=0 B=n C=0. We can set a bound only by noting the D value. + boundsInfo.addYValue(-D/B); + } } else { // No solutions } diff --git a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoCircleTest.java b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoCircleTest.java index 6f9a86b5d23..0190cde7f02 100755 --- a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoCircleTest.java +++ b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoCircleTest.java @@ -405,4 +405,18 @@ public class GeoCircleTest extends LuceneTestCase { assertTrue(solid.isWithin(gp)); } + @Test + public void testBoundsFailureCase2() { + final GeoCircle gc = GeoCircleFactory.makeGeoCircle(PlanetModel.WGS84, -2.7574435614238194E-13, 0.0, 1.5887859182593391); + final GeoPoint gp = new GeoPoint(PlanetModel.WGS84, 0.7980359504429014, 1.5964981068121482); + final XYZBounds bounds = new XYZBounds(); + gc.getBounds(bounds); + System.out.println("Bounds = "+bounds); + System.out.println("Point = "+gp); + final XYZSolid solid = XYZSolidFactory.makeXYZSolid(PlanetModel.WGS84, bounds.getMinimumX(), bounds.getMaximumX(), bounds.getMinimumY(), bounds.getMaximumY(), bounds.getMinimumZ(), bounds.getMaximumZ()); + + assert gc.isWithin(gp)?solid.isWithin(gp):true; + + } + }