You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@lucene.apache.org by kw...@apache.org on 2018/01/25 12:06:05 UTC
lucene-solr:branch_7x: LUCENE-8136: Adopt iterative convergence for
construction of planes from two vectors. Thanks, Ignacio Vera!
Repository: lucene-solr
Updated Branches:
refs/heads/branch_7x eb6467391 -> a69f34084
LUCENE-8136: Adopt iterative convergence for construction of planes from two vectors. Thanks, Ignacio Vera!
Project: http://git-wip-us.apache.org/repos/asf/lucene-solr/repo
Commit: http://git-wip-us.apache.org/repos/asf/lucene-solr/commit/a69f3408
Tree: http://git-wip-us.apache.org/repos/asf/lucene-solr/tree/a69f3408
Diff: http://git-wip-us.apache.org/repos/asf/lucene-solr/diff/a69f3408
Branch: refs/heads/branch_7x
Commit: a69f34084f8852f04c325ec0ae3817725d94f772
Parents: eb64673
Author: Karl Wright <Da...@gmail.com>
Authored: Thu Jan 25 07:03:45 2018 -0500
Committer: Karl Wright <Da...@gmail.com>
Committed: Thu Jan 25 07:05:58 2018 -0500
----------------------------------------------------------------------
.../apache/lucene/spatial3d/geom/GeoPoint.java | 8 ++
.../apache/lucene/spatial3d/geom/Vector.java | 85 +++++++++++++++-----
.../lucene/spatial3d/geom/GeoPolygonTest.java | 24 +++---
.../lucene/spatial3d/geom/RandomPlaneTest.java | 74 +++++++++++++++++
4 files changed, 160 insertions(+), 31 deletions(-)
----------------------------------------------------------------------
http://git-wip-us.apache.org/repos/asf/lucene-solr/blob/a69f3408/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/GeoPoint.java
----------------------------------------------------------------------
diff --git a/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/GeoPoint.java b/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/GeoPoint.java
index de51a2b..66c6226 100755
--- a/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/GeoPoint.java
+++ b/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/GeoPoint.java
@@ -190,6 +190,14 @@ public class GeoPoint extends Vector implements SerializableObject {
}
/** Compute whether point matches another.
+ *@param p is the other point.
+ *@return true if the same.
+ */
+ public boolean isIdentical(final GeoPoint p) {
+ return isIdentical(p.x, p.y, p.z);
+ }
+
+ /** Compute whether point matches another.
*@param x is the x value
*@param y is the y value
*@param z is the z value
http://git-wip-us.apache.org/repos/asf/lucene-solr/blob/a69f3408/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Vector.java
----------------------------------------------------------------------
diff --git a/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Vector.java b/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Vector.java
index 8fad521..454439b 100755
--- a/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Vector.java
+++ b/lucene/spatial3d/src/java/org/apache/lucene/spatial3d/geom/Vector.java
@@ -27,7 +27,7 @@ public class Vector {
* Values that are all considered to be essentially zero have a magnitude
* less than this.
*/
- public static final double MINIMUM_RESOLUTION = 1.5e-12;
+ public static final double MINIMUM_RESOLUTION = 1.0e-12;
/**
* Angular version of minimum resolution.
*/
@@ -72,20 +72,78 @@ public class Vector {
* @param BZ is the Z value of the second
*/
public Vector(final Vector A, final double BX, final double BY, final double BZ) {
+ // We're really looking at two vectors and computing a perpendicular one from that.
+ // We can think of this as having three points -- the origin, and two points that aren't the origin.
+ // Normally, we can compute the perpendicular vector this way:
// x = u2v3 - u3v2
// y = u3v1 - u1v3
// z = u1v2 - u2v1
+ // Sometimes that produces a plane that does not contain the original three points, however, due to
+ // numerical precision issues. Then we continue making the answer more precise using the
+ // Gram-Schmidt process: https://en.wikipedia.org/wiki/Gram%E2%80%93Schmidt_process
+
+ // Compute the naive perpendicular
final double thisX = A.y * BZ - A.z * BY;
final double thisY = A.z * BX - A.x * BZ;
final double thisZ = A.x * BY - A.y * BX;
+
final double magnitude = magnitude(thisX, thisY, thisZ);
- if (Math.abs(magnitude) < MINIMUM_RESOLUTION) {
+ if (magnitude < MINIMUM_RESOLUTION) {
throw new IllegalArgumentException("Degenerate/parallel vector constructed");
}
- final double inverseMagnitude = 1.0 / magnitude;
- this.x = thisX * inverseMagnitude;
- this.y = thisY * inverseMagnitude;
- this.z = thisZ * inverseMagnitude;
+
+ final double inverseMagnitude = 1.0/magnitude;
+
+ double normalizeX = thisX * inverseMagnitude;
+ double normalizeY = thisY * inverseMagnitude;
+ double normalizeZ = thisZ * inverseMagnitude;
+ // For a plane to work, the dot product between the normal vector
+ // and the points needs to be less than the minimum resolution.
+ // This is sometimes not true for points that are very close. Therefore
+ // we need to adjust
+ int i = 0;
+ while (true) {
+ final double currentDotProdA = A.x * normalizeX + A.y * normalizeY + A.z * normalizeZ;
+ final double currentDotProdB = BX * normalizeX + BY * normalizeY + BZ * normalizeZ;
+ if (Math.abs(currentDotProdA) < MINIMUM_RESOLUTION && Math.abs(currentDotProdB) < MINIMUM_RESOLUTION) {
+ break;
+ }
+ // Converge on the one that has largest dot product
+ final double currentVectorX;
+ final double currentVectorY;
+ final double currentVectorZ;
+ final double currentDotProd;
+ if (Math.abs(currentDotProdA) > Math.abs(currentDotProdB)) {
+ currentVectorX = A.x;
+ currentVectorY = A.y;
+ currentVectorZ = A.z;
+ currentDotProd = currentDotProdA;
+ } else {
+ currentVectorX = BX;
+ currentVectorY = BY;
+ currentVectorZ = BZ;
+ currentDotProd = currentDotProdB;
+ }
+
+ // Adjust
+ normalizeX = normalizeX - currentDotProd * currentVectorX;
+ normalizeY = normalizeY - currentDotProd * currentVectorY;
+ normalizeZ = normalizeZ - currentDotProd * currentVectorZ;
+ // Normalize
+ final double correctedMagnitude = magnitude(normalizeX, normalizeY, normalizeZ);
+ final double inverseCorrectedMagnitude = 1.0 / correctedMagnitude;
+ normalizeX = normalizeX * inverseCorrectedMagnitude;
+ normalizeY = normalizeY * inverseCorrectedMagnitude;
+ normalizeZ = normalizeZ * inverseCorrectedMagnitude;
+ //This is probably not needed as the method seems to converge
+ //quite quickly. But it is safer to have a way out.
+ if (i++ > 10) {
+ throw new IllegalArgumentException("Plane could not be constructed! Could not find a normal vector.");
+ }
+ }
+ this.x = normalizeX;
+ this.y = normalizeY;
+ this.z = normalizeZ;
}
/**
@@ -98,20 +156,7 @@ public class Vector {
* @param B is the second
*/
public Vector(final Vector A, final Vector B) {
- // x = u2v3 - u3v2
- // y = u3v1 - u1v3
- // z = u1v2 - u2v1
- final double thisX = A.y * B.z - A.z * B.y;
- final double thisY = A.z * B.x - A.x * B.z;
- final double thisZ = A.x * B.y - A.y * B.x;
- final double magnitude = magnitude(thisX, thisY, thisZ);
- if (Math.abs(magnitude) < MINIMUM_RESOLUTION) {
- throw new IllegalArgumentException("Degenerate/parallel vector constructed");
- }
- final double inverseMagnitude = 1.0 / magnitude;
- this.x = thisX * inverseMagnitude;
- this.y = thisY * inverseMagnitude;
- this.z = thisZ * inverseMagnitude;
+ this(A, B.x, B.y, B.z);
}
/** Compute a magnitude of an x,y,z value.
http://git-wip-us.apache.org/repos/asf/lucene-solr/blob/a69f3408/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoPolygonTest.java
----------------------------------------------------------------------
diff --git a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoPolygonTest.java b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoPolygonTest.java
index 8892111..6291057 100755
--- a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoPolygonTest.java
+++ b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/GeoPolygonTest.java
@@ -1035,17 +1035,19 @@ shape:
GeoPoint point1 = new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434456), Geo3DUtil.fromDegrees(14.459204));
GeoPoint point2 = new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.43394), Geo3DUtil.fromDegrees(14.459206));
GeoPoint check = new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434067), Geo3DUtil.fromDegrees(14.458927));
- SidedPlane plane = new SidedPlane(check, point1, point2);
- assertTrue(plane.isWithin(check));
- assertTrue(plane.isWithin(point1));
- assertTrue(plane.isWithin(point2));
- //POLYGON((14.459204 -23.434456, 14.459206 -23.43394,14.458647 -23.434196, 14.458646 -23.434452,14.459204 -23.434456))
- List<GeoPoint> points = new ArrayList<>();
- points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434456), Geo3DUtil.fromDegrees(14.459204)));
- points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees( -23.43394), Geo3DUtil.fromDegrees(14.459206)));
- points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434196), Geo3DUtil.fromDegrees(14.458647)));
- points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434452), Geo3DUtil.fromDegrees(14.458646)));
- GeoPolygonFactory.makeGeoPolygon(PlanetModel.SPHERE, points);
+ if (!point1.isIdentical(point2) && !check.isIdentical(point1) && !check.isIdentical(point2)) {
+ SidedPlane plane = new SidedPlane(check, point1, point2);
+ assertTrue(plane.isWithin(check));
+ assertTrue(plane.isWithin(point1));
+ assertTrue(plane.isWithin(point2));
+ //POLYGON((14.459204 -23.434456, 14.459206 -23.43394,14.458647 -23.434196, 14.458646 -23.434452,14.459204 -23.434456))
+ List<GeoPoint> points = new ArrayList<>();
+ points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434456), Geo3DUtil.fromDegrees(14.459204)));
+ points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees( -23.43394), Geo3DUtil.fromDegrees(14.459206)));
+ points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434196), Geo3DUtil.fromDegrees(14.458647)));
+ points.add(new GeoPoint(PlanetModel.SPHERE, Geo3DUtil.fromDegrees(-23.434452), Geo3DUtil.fromDegrees(14.458646)));
+ GeoPolygonFactory.makeGeoPolygon(PlanetModel.SPHERE, points);
+ }
}
http://git-wip-us.apache.org/repos/asf/lucene-solr/blob/a69f3408/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/RandomPlaneTest.java
----------------------------------------------------------------------
diff --git a/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/RandomPlaneTest.java b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/RandomPlaneTest.java
new file mode 100644
index 0000000..1e19194
--- /dev/null
+++ b/lucene/spatial3d/src/test/org/apache/lucene/spatial3d/geom/RandomPlaneTest.java
@@ -0,0 +1,74 @@
+/*
+ * 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.spatial3d.geom;
+
+import java.util.ArrayList;
+import java.util.List;
+
+import com.carrotsearch.randomizedtesting.annotations.Repeat;
+import org.junit.Test;
+
+/**
+ * Random test for planes.
+ */
+public class RandomPlaneTest extends RandomGeo3dShapeGenerator {
+
+ @Test
+ @Repeat(iterations = 10)
+ public void testPlaneAccuracy() {
+ PlanetModel planetModel = randomPlanetModel();
+ GeoPoint point1 = randomGeoPoint(planetModel);
+ for (int i= 0; i < 1000; i++) {
+ double dist = random().nextDouble() * Vector.MINIMUM_ANGULAR_RESOLUTION + Vector.MINIMUM_ANGULAR_RESOLUTION;
+ double bearing = random().nextDouble() * 2 * Math.PI;
+ GeoPoint point2 = planetModel.surfacePointOnBearing(point1, dist, bearing );
+ GeoPoint check = randomGeoPoint(planetModel);
+ if (!point1.isNumericallyIdentical(point2)) {
+ SidedPlane plane = new SidedPlane(check, point1, point2);
+ String msg = dist + " point 1: " + point1 + ", point 2: " + point2 + " , check: " + check;
+ assertTrue(msg, plane.isWithin(check));
+ assertTrue(msg, plane.isWithin(point2));
+ assertTrue(msg, plane.isWithin(point1));
+ }
+ else {
+ assertFalse("numerically identical", true);
+ }
+ }
+ }
+
+ /*
+ @Test
+ @Repeat(iterations = 10)
+ public void testPolygonAccuracy() {
+ PlanetModel planetModel = randomPlanetModel();
+ GeoPoint point1 = randomGeoPoint(planetModel);
+ for (int i= 0; i < 1000; i++) {
+ double dist = random().nextDouble() * 1e-6 + Vector.MINIMUM_ANGULAR_RESOLUTION;
+ GeoPoint point2 = planetModel.surfacePointOnBearing(point1, dist, 0);
+ GeoPoint point3 = planetModel.surfacePointOnBearing(point1, dist, 0.5 * Math.PI);
+
+ List<GeoPoint> points = new ArrayList<>();
+ points.add(point1);
+ points.add(point2);
+ points.add(point3);
+ GeoPolygonFactory.makeGeoPolygon(planetModel, points);
+
+ }
+ }
+ */
+}