You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@sis.apache.org by js...@apache.org on 2026/10/14 15:00:14 UTC
(sis) 07/12: feat(Geometry): implement default polygon area computation
This is an automated email from the ASF dual-hosted git repository.
asf-gitbox-commits pushed a commit to branch geoapi-4.0
in repository https://gitbox.apache.org/repos/asf/sis.git
commit 42152ec38249a43068359f3df73250c16ef67f8a
Author: jsorel <jo...@geomatys.com>
AuthorDate: Wed Oct 14 12:31:43 2026 +0200
feat(Geometry): implement default polygon area computation
---
.../sis/geometries/internal/shared/Curves.java | 69 +++++++++++++++++++++-
.../geometries/internal/shared/DefaultPolygon.java | 15 ++++-
.../apache/sis/geometries/surface/Triangle.java | 18 +++---
3 files changed, 88 insertions(+), 14 deletions(-)
diff --git a/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/Curves.java b/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/Curves.java
index cda7562b28..47e6c59ea8 100644
--- a/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/Curves.java
+++ b/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/Curves.java
@@ -21,11 +21,13 @@ import javax.measure.Quantity;
import javax.measure.Unit;
import org.apache.sis.geometries.DataPoints;
import org.apache.sis.geometries.Point;
+import org.apache.sis.maths.Vector;
import org.apache.sis.measure.Units;
/**
- * Helper methods for locating positions at a distance along a curve.
+ * Helper methods for computations on curves, such as locating positions at a distance along a curve
+ * or measuring the area enclosed by a ring.
*
* @author Johann Sorel (Geomatys)
*/
@@ -34,6 +36,71 @@ public final class Curves {
private Curves() {
}
+ /**
+ * Returns the area enclosed by the given ring of positions.
+ *
+ * <ul>
+ * <li>In 2 dimensions, the area is computed with the shoelace formula.</li>
+ * <li>In 3 dimensions or more, the area is computed on the first three coordinates with Newell's
+ * method: half the norm of the sum of the cross products of consecutive positions. This is
+ * exact for a ring lying in a plane, in any orientation of that plane, and gives the area
+ * projected on the best-fitting plane otherwise.</li>
+ * </ul>
+ *
+ * <p>The ring may repeat its first position at the end or not. The positions are taken relative
+ * to the first one, which avoids losing precision with large coordinate values.</p>
+ *
+ * @param ring positions of the ring, in order.
+ * @return the area enclosed by the ring, in the square of the unit of the coordinates.
+ */
+ public static double ringArea(final DataPoints ring) {
+ final int n = ring.size();
+ final int dim = ring.getDimension();
+ if (n < 3 || dim < 2) {
+ return 0;
+ }
+ final Vector<?> origin = ring.getPosition(0);
+ if (dim == 2) {
+ final double ox = origin.get(0);
+ final double oy = origin.get(1);
+ Vector<?> p = ring.getPosition(n - 1);
+ double ax = p.get(0) - ox;
+ double ay = p.get(1) - oy;
+ double sum = 0;
+ for (int i = 0; i < n; i++) {
+ p = ring.getPosition(i);
+ final double bx = p.get(0) - ox;
+ final double by = p.get(1) - oy;
+ sum += ax*by - ay*bx;
+ ax = bx;
+ ay = by;
+ }
+ return Math.abs(sum) / 2;
+ } else {
+ final double ox = origin.get(0);
+ final double oy = origin.get(1);
+ final double oz = origin.get(2);
+ Vector<?> p = ring.getPosition(n - 1);
+ double ax = p.get(0) - ox;
+ double ay = p.get(1) - oy;
+ double az = p.get(2) - oz;
+ double nx = 0, ny = 0, nz = 0;
+ for (int i = 0; i < n; i++) {
+ p = ring.getPosition(i);
+ final double bx = p.get(0) - ox;
+ final double by = p.get(1) - oy;
+ final double bz = p.get(2) - oz;
+ nx += ay*bz - az*by;
+ ny += az*bx - ax*bz;
+ nz += ax*by - ay*bx;
+ ax = bx;
+ ay = by;
+ az = bz;
+ }
+ return Math.sqrt(nx*nx + ny*ny + nz*nz) / 2;
+ }
+ }
+
/**
* Returns whether the given distance is expressed in percent of the curve length.
*
diff --git a/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/DefaultPolygon.java b/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/DefaultPolygon.java
index 42ab75bc80..5df2377c5c 100644
--- a/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/DefaultPolygon.java
+++ b/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/internal/shared/DefaultPolygon.java
@@ -22,6 +22,7 @@ import javax.measure.Quantity;
import org.opengis.geometry.Envelope;
import org.opengis.referencing.crs.CoordinateReferenceSystem;
import org.apache.sis.geometries.Curve;
+import org.apache.sis.geometries.DataPoints;
import org.apache.sis.geometries.Geometries;
import org.apache.sis.geometries.Geometry;
import org.apache.sis.geometries.GeometryFactory;
@@ -87,10 +88,20 @@ public non-sealed class DefaultPolygon extends AbstractGeometry implements Polyg
return getExteriorRing().getEnvelope();
}
+ /**
+ * Returns the area enclosed by the exterior ring minus the areas enclosed by the interior rings.
+ * The area is measured in the plane of the polygon, which may have any orientation in a
+ * 3-dimensional space.
+ *
+ * @see Curves#ringArea(DataPoints)
+ */
@Override
public Quantity<?> getArea() {
- //TODO : fallback on JTS until implemented
- return Quantities.create(asJTS().getArea(), Geometries.getAreaUnit(this));
+ double area = Curves.ringArea(exterior.getDataPoints());
+ for (final LinearRing interior : interiors) {
+ area -= Curves.ringArea(interior.getDataPoints());
+ }
+ return Quantities.create(area, Geometries.getAreaUnit(this));
}
@Override
diff --git a/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/surface/Triangle.java b/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/surface/Triangle.java
index 8fb960cc00..30baca299e 100644
--- a/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/surface/Triangle.java
+++ b/incubator/src/org.apache.sis.geometry/main/org/apache/sis/geometries/surface/Triangle.java
@@ -25,10 +25,12 @@ import org.opengis.geometry.Envelope;
import org.opengis.referencing.crs.CoordinateReferenceSystem;
import org.apache.sis.geometries.BBox;
import org.apache.sis.geometries.DataPoints;
+import org.apache.sis.geometries.Geometries;
import org.apache.sis.geometries.GeometryType;
import org.apache.sis.geometries.Point;
import org.apache.sis.geometries.DataPointsType;
import org.apache.sis.geometries.curve.LinearRing;
+import org.apache.sis.geometries.internal.shared.Curves;
import org.apache.sis.geometries.internal.shared.DefaultTriangle;
import org.apache.sis.geometries.operation.simplify.greedyinsert.WTriangle;
import org.apache.sis.geometries.operation.triangulate.delaunay.OrientedTriangle;
@@ -37,7 +39,6 @@ import org.apache.sis.maths.Vector;
import org.apache.sis.maths.Vector2D;
import org.apache.sis.maths.Vector3D;
import org.apache.sis.measure.Quantities;
-import org.apache.sis.measure.Units;
/**
@@ -124,18 +125,13 @@ public sealed interface Triangle extends Polygon
return env;
}
+ /**
+ * Returns the area of this triangle, measured in its plane, which may have any orientation
+ * in a 3-dimensional space.
+ */
@Override
default Quantity<?> getArea() {
- final DataPoints points = getExteriorRing().getDataPoints();
- final Vector<?> a = points.getPosition(0);
- final Vector<?> b = points.getPosition(1);
- final Vector<?> c = points.getPosition(2);
- final double area = (
- a.get(0) * (b.get(1) - c.get(1))
- + b.get(0) * (c.get(1) - a.get(1))
- + c.get(0) * (a.get(1) - b.get(1))
- ) / 2.0;
- return Quantities.create(Math.abs(area), Units.SQUARE_METRE);
+ return Quantities.create(Curves.ringArea(getExteriorRing().getDataPoints()), Geometries.getAreaUnit(this));
}
/**