You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@commons.apache.org by er...@apache.org on 2020/01/11 21:52:30 UTC

[commons-geometry] branch master updated: GEOMETRY-71 Fix spherical barycenter algorithm

This is an automated email from the ASF dual-hosted git repository.

erans pushed a commit to branch master
in repository https://gitbox.apache.org/repos/asf/commons-geometry.git


The following commit(s) were added to refs/heads/master by this push:
     new 6fab539  GEOMETRY-71 Fix spherical barycenter algorithm
     new dd5c19b  Merge branch 'GEOMETRY-71__Baljit'
6fab539 is described below

commit 6fab539c42a561e5d7def8f217e0437963958686
Author: Baljit Singh <ba...@ge.com>
AuthorDate: Tue Jan 7 10:20:26 2020 -0500

    GEOMETRY-71 Fix spherical barycenter algorithm
---
 .../geometry/spherical/twod/ConvexArea2S.java      | 30 ++++++----
 .../geometry/spherical/twod/RegionBSPTree2S.java   |  8 +--
 .../spherical/twod/RegionBSPTree2STest.java        | 66 ++++++++++++++++++++++
 3 files changed, 86 insertions(+), 18 deletions(-)

diff --git a/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/ConvexArea2S.java b/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/ConvexArea2S.java
index d87dceb..c187cbb 100644
--- a/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/ConvexArea2S.java
+++ b/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/ConvexArea2S.java
@@ -128,29 +128,35 @@ public final class ConvexArea2S extends AbstractConvexHyperplaneBoundedRegion<Po
     /** {@inheritDoc} */
     @Override
     public Point2S getBarycenter() {
-        List<GreatArc> arcs = getBoundaries();
-        int numSides = arcs.size();
+        Vector3D weighted = getWeightedBarycenter();
+        return weighted == null ? null : Point2S.from(weighted);
+    }
 
-        if (numSides == 0) {
+    /**
+     * Returns the weighted vector for barycenter.
+     *
+     * @return Weighted barycenter
+     */
+    Vector3D getWeightedBarycenter() {
+        List<GreatArc> arcs = getBoundaries();
+        switch (arcs.size()) {
+        case 0:
             // full space; no barycenter
             return null;
-        } else if (numSides == 1) {
+        case 1:
             // hemisphere; barycenter is the pole of the hemisphere
-            return arcs.get(0).getCircle().getPolePoint();
-        } else {
+            GreatArc singleArc = arcs.get(0);
+            return singleArc.getCircle().getPole().withNorm(singleArc.getSize());
+        default:
             // 2 or more sides; use an extension of the approach outlined here:
             // https://archive.org/details/centroidinertiat00broc
             // In short, the barycenter is the sum of the pole vectors of each side
             // multiplied by their arc lengths.
             Vector3D barycenter = Vector3D.ZERO;
-
             for (GreatArc arc : getBoundaries()) {
-                barycenter = Vector3D.linearCombination(
-                        1, barycenter,
-                        arc.getSize(), arc.getCircle().getPole());
+                barycenter = barycenter.add(arc.getCircle().getPole().withNorm(arc.getSize()));
             }
-
-            return Point2S.from(barycenter);
+            return barycenter;
         }
     }
 
diff --git a/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2S.java b/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2S.java
index 8c721aa..17ddd45 100644
--- a/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2S.java
+++ b/commons-geometry-spherical/src/main/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2S.java
@@ -169,12 +169,8 @@ public class RegionBSPTree2S extends AbstractRegionBSPTree<Point2S, RegionBSPTre
 
         double size;
         for (ConvexArea2S area : areas) {
-            size = area.getSize();
-
-            sizeSum += size;
-            barycenterVector = Vector3D.linearCombination(
-                    1, barycenterVector,
-                    size, area.getBarycenter().getVector());
+            sizeSum += area.getSize();
+            barycenterVector = barycenterVector.add(area.getWeightedBarycenter());
         }
 
         Point2S barycenter = barycenterVector.eq(Vector3D.ZERO, precision) ?
diff --git a/commons-geometry-spherical/src/test/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2STest.java b/commons-geometry-spherical/src/test/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2STest.java
index 7a7943d..3e0e290 100644
--- a/commons-geometry-spherical/src/test/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2STest.java
+++ b/commons-geometry-spherical/src/test/java/org/apache/commons/geometry/spherical/twod/RegionBSPTree2STest.java
@@ -706,6 +706,54 @@ public class RegionBSPTree2STest {
         checkBarycenterConsistency(france);
     }
 
+    @Test
+    public void testCircleToPolygonBarycenter() {
+        double radius = 0.0001;
+        Point2S center = Point2S.of(1.0, 1.0);
+        int numPts = 200;
+
+        // counterclockwise
+        RegionBSPTree2S ccw = circleToPolygon(center, radius, numPts, false);
+        SphericalTestUtils.assertPointsEq(center, ccw.getBarycenter(), BARYCENTER_EPS);
+
+        // clockwise; barycenter should just be antipodal for the circle center
+        RegionBSPTree2S cw = circleToPolygon(center, radius, numPts, true);
+        SphericalTestUtils.assertPointsEq(center.antipodal(), cw.getBarycenter(), BARYCENTER_EPS);
+    }
+
+    @Test
+    public void testCircleToPolygonSize() {
+        double radius = 0.0001;
+        Point2S center = Point2S.of(1.0, 1.0);
+        int numPts = 200;
+
+        // https://en.wikipedia.org/wiki/Spherical_cap
+        double ccwArea = 4.0 * PlaneAngleRadians.PI * Math.pow(Math.sin(radius / 2.0), 2.0);
+        double cwArea = 4.0 * PlaneAngleRadians.PI - ccwArea;
+
+        RegionBSPTree2S ccw = circleToPolygon(center, radius, numPts, false);
+        Assert.assertEquals("Counterclockwise size", ccwArea, ccw.getSize(), TEST_EPS);
+
+        RegionBSPTree2S cw = circleToPolygon(center, radius, numPts, true);
+        Assert.assertEquals("Clockwise size", cwArea, cw.getSize(), TEST_EPS);
+    }
+
+    @Test
+    public void testCircleToPolygonBoundarySize() {
+        double radius = 0.0001;
+        Point2S center = Point2S.of(1.0, 1.0);
+        int numPts = 200;
+
+        // boundary size is independent from winding
+        double boundary = PlaneAngleRadians.TWO_PI * Math.sin(radius);
+
+        RegionBSPTree2S ccw = circleToPolygon(center, radius, numPts, false);
+        Assert.assertEquals("Counterclockwise boundary size", boundary, ccw.getBoundarySize(), 1.0e-7);
+
+        RegionBSPTree2S cw = circleToPolygon(center, radius, numPts, true);
+        Assert.assertEquals("Clockwise boundary size", boundary, cw.getBoundarySize(), 1.0e-7);
+    }
+
     /**
      * Insert convex subhyperplanes defining the positive quadrant area.
      * @param tree
@@ -791,4 +839,22 @@ public class RegionBSPTree2STest {
             SphericalTestUtils.assertPointsEq(barycenter, computedBarycenter, BARYCENTER_EPS);
         }
     }
+
+    private static RegionBSPTree2S circleToPolygon(Point2S center, double radius, int numPts, boolean clockwise) {
+        List<Point2S> pts = new ArrayList<>(numPts);
+
+        // get an arbitrary point on the circle boundary
+        pts.add(Transform2S.createRotation(center.getVector().orthogonal(), radius).apply(center));
+
+        // create the list of boundary points by rotating the previous point around the circle center
+        double span = PlaneAngleRadians.TWO_PI / numPts;
+
+        // negate the span for clockwise winding
+        Transform2S rotate = Transform2S.createRotation(center, clockwise ? -span : span);
+        for (int i = 1; i < numPts; ++i) {
+            pts.add(rotate.apply(pts.get(i - 1)));
+        }
+
+        return GreatArcPath.fromVertexLoop(pts, TEST_PRECISION).toTree();
+    }
 }