Merge branch '3.11' of github.com:libgeos/geos into 3.11

NEWS.md

@@ -10,6 +10,8 @@
   - GEOSConcaveHullOfPolygons, avoid crash on zero-area input (GH-1071, Dan Baston)
   - MinimumClearance, avoid crash on NaN inputs (GH-1079, Dan Baston)
   - Fix TopologyPreservingSimplifier/TaggedLineString to avoid jumping components (JTS-1096, Martin Davis)
+  - Fix BufferOp to increase length of segments removed by heuristic (GH-1200, Martin Davis)
+  - Add ring buffer hole removal heuristic (GH-1233, Martin Davis)
 
 ## Changes in 3.11.4
 2024-06-05

include/geos/operation/buffer/BufferCurveSetBuilder.h

@@ -70,6 +70,8 @@ namespace buffer { // geos.operation.buffer
  *
  */
 class GEOS_DLL BufferCurveSetBuilder {
+    using CoordinateSequence = geos::geom::CoordinateSequence;
+    using Envelope = geos::geom::Envelope;
 
 private:
 
@@ -118,7 +120,7 @@ class GEOS_DLL BufferCurveSetBuilder {
 
     void addPolygon(const geom::Polygon* p);
 
-    void addRingBothSides(const geom::CoordinateSequence* coord, double p_distance);
+    void addLinearRingSides(const geom::CoordinateSequence* coord, double p_distance);
 
     /**
      * Add an offset curve for a polygon ring.
@@ -138,10 +140,14 @@ class GEOS_DLL BufferCurveSetBuilder {
      * @param cwRightLoc the location on the R side of the ring
      *                   (if it is CW)
      */
-    void addRingSide(const geom::CoordinateSequence* coord,
+    void addPolygonRingSide(const geom::CoordinateSequence* coord,
                      double offsetDistance, int side, geom::Location cwLeftLoc,
                      geom::Location cwRightLoc);
 
+    void addRingSide(const geom::CoordinateSequence* coord,
+                     double offsetDistance, int side, geom::Location leftLoc,
+                     geom::Location rightLoc);
+
     /**
      * Tests whether the offset curve for a ring is fully inverted.
      * An inverted ("inside-out") curve occurs in some specific situations
@@ -187,7 +193,10 @@ class GEOS_DLL BufferCurveSetBuilder {
      * @param bufferDistance
      * @return
      */
-    bool isErodedCompletely(const geom::LinearRing* ringCoord,
+    bool isRingFullyEroded(const geom::LinearRing* ring, bool isHole,
+                            double bufferDistance);
+
+    bool isRingFullyEroded(const CoordinateSequence* ringCoord, const Envelope* env, bool isHole,
                             double bufferDistance);
 
     /**

src/operation/buffer/BufferCurveSetBuilder.cpp

@@ -208,16 +208,44 @@ BufferCurveSetBuilder::addLineString(const LineString* line)
      * Singled-sided buffers currently treat rings as if they are lines.
      */
     if (coord->isRing() && ! curveBuilder.getBufferParameters().isSingleSided()) {
-        addRingBothSides(coord.get(), distance);
+        addLinearRingSides(coord.get(), distance);
     }
     else {
         std::vector<CoordinateSequence*> lineList;
         curveBuilder.getLineCurve(coord.get(), distance, lineList);
         addCurves(lineList, Location::EXTERIOR, Location::INTERIOR);
     }
-
 }
 
+/* private */
+void
+BufferCurveSetBuilder::addLinearRingSides(const CoordinateSequence* coord, double p_distance)
+{
+    /*
+     * (f "hole" side will be eroded completely, avoid generating it.
+     * This prevents hole artifacts (e.g. /~https://github.com/libgeos/geos/issues/1223)
+     */
+    //-- distance is assumed positive, due to previous checks
+    Envelope env;
+    coord->expandEnvelope(env);
+    bool isHoleComputed = ! isRingFullyEroded(coord, &env, true, distance);
+
+    bool isCCW = isRingCCW(coord);
+
+    bool isShellLeft = ! isCCW;
+    if (isShellLeft || isHoleComputed) {
+        addRingSide(coord, p_distance,
+                    Position::LEFT,
+                    Location::EXTERIOR, Location::INTERIOR);
+    }
+
+    bool isShellRight = isCCW;
+    if (isShellRight || isHoleComputed) {
+        addRingSide(coord, p_distance,
+                    Position::RIGHT,
+                    Location::INTERIOR, Location::EXTERIOR);
+    }
+}
 
 /*private*/
 void
@@ -235,7 +263,7 @@ BufferCurveSetBuilder::addPolygon(const Polygon* p)
 
     // optimization - don't bother computing buffer
     // if the polygon would be completely eroded
-    if(distance < 0.0 && isErodedCompletely(shell, distance)) {
+    if(distance < 0.0 && isRingFullyEroded(shell, false, distance)) {
 #if GEOS_DEBUG
         std::cerr << __FUNCTION__ << ": polygon is eroded completely " << std::endl;
 #endif
@@ -251,7 +279,7 @@ BufferCurveSetBuilder::addPolygon(const Polygon* p)
         return;
     }
 
-    addRingSide(
+    addPolygonRingSide(
         shellCoord.get(),
         offsetDistance,
         offsetSide,
@@ -264,7 +292,7 @@ BufferCurveSetBuilder::addPolygon(const Polygon* p)
 
         // optimization - don't bother computing buffer for this hole
         // if the hole would be completely covered
-        if(distance > 0.0 && isErodedCompletely(hole, -distance)) {
+        if(distance > 0.0 && isRingFullyEroded(hole, true, distance)) {
             continue;
         }
 
@@ -273,7 +301,7 @@ BufferCurveSetBuilder::addPolygon(const Polygon* p)
         // Holes are topologically labelled opposite to the shell,
         // since the interior of the polygon lies on their opposite
         // side (on the left, if the hole is oriented CCW)
-        addRingSide(
+        addPolygonRingSide(
             holeCoord.get(),
             offsetDistance,
             Position::opposite(offsetSide),
@@ -284,22 +312,7 @@ BufferCurveSetBuilder::addPolygon(const Polygon* p)
 
 /* private */
 void
-BufferCurveSetBuilder::addRingBothSides(const CoordinateSequence* coord, double p_distance)
-{
-    addRingSide(coord, p_distance,
-                Position::LEFT,
-                Location::EXTERIOR, Location::INTERIOR);
-    /* Add the opposite side of the ring
-    */
-    addRingSide(coord, p_distance,
-                Position::RIGHT,
-                Location::INTERIOR, Location::EXTERIOR);
-}
-
-
-/* private */
-void
-BufferCurveSetBuilder::addRingSide(const CoordinateSequence* coord,
+BufferCurveSetBuilder::addPolygonRingSide(const CoordinateSequence* coord,
                                       double offsetDistance, int side, geom::Location cwLeftLoc, geom::Location cwRightLoc)
 {
 
@@ -331,6 +344,14 @@ BufferCurveSetBuilder::addRingSide(const CoordinateSequence* coord,
 #endif
         side = Position::opposite(side);
     }
+    addRingSide(coord, offsetDistance, side, leftLoc, rightLoc);
+}
+
+/* private */
+void
+BufferCurveSetBuilder::addRingSide(const CoordinateSequence* coord,
+                                      double offsetDistance, int side, geom::Location leftLoc, geom::Location rightLoc)
+{
     std::vector<CoordinateSequence*> lineList;
     curveBuilder.getRingCurve(coord, side, offsetDistance, lineList);
     // ASSERT: lineList contains exactly 1 curve (this is teh JTS semantics)
@@ -407,14 +428,22 @@ BufferCurveSetBuilder::maxDistance(const CoordinateSequence*  pts, const Coordin
 
 /*private*/
 bool
-BufferCurveSetBuilder::isErodedCompletely(const LinearRing* ring,
+BufferCurveSetBuilder::isRingFullyEroded(const LinearRing* ring, bool isHole,
         double bufferDistance)
 {
     const CoordinateSequence* ringCoord = ring->getCoordinatesRO();
+    const Envelope* env = ring->getEnvelopeInternal();
+    return isRingFullyEroded(ringCoord, env, isHole, bufferDistance);
+}
 
+/*private*/
+bool
+BufferCurveSetBuilder::isRingFullyEroded(const CoordinateSequence* ringCoord, const Envelope* env, bool isHole,
+        double bufferDistance)
+{
     // degenerate ring has no area
     if(ringCoord->getSize() < 4) {
-        return bufferDistance < 0;
+        return true;
     }
 
     // important test to eliminate inverted triangle bug
@@ -423,10 +452,16 @@ BufferCurveSetBuilder::isErodedCompletely(const LinearRing* ring,
         return isTriangleErodedCompletely(ringCoord, bufferDistance);
     }
 
-    const Envelope* env = ring->getEnvelopeInternal();
-    double envMinDimension = std::min(env->getHeight(), env->getWidth());
-    if(bufferDistance < 0.0 && 2 * std::abs(bufferDistance) > envMinDimension) {
-        return true;
+    bool isErodable = 
+        (  isHole && bufferDistance > 0) ||
+        (! isHole && bufferDistance < 0);
+
+    if (isErodable) {
+      //-- if envelope is narrower than twice the buffer distance, ring is eroded
+        double envMinDimension = std::min(env->getHeight(), env->getWidth());
+        if (2 * std::abs(bufferDistance) > envMinDimension) {
+            return true;
+        }
     }
     return false;
 }

src/operation/buffer/OffsetSegmentGenerator.cpp

@@ -49,7 +49,7 @@ namespace operation { // geos.operation
 namespace buffer { // geos.operation.buffer
 
 /*private data*/
-const double OffsetSegmentGenerator::CURVE_VERTEX_SNAP_DISTANCE_FACTOR = 1.0E-6;
+const double OffsetSegmentGenerator::CURVE_VERTEX_SNAP_DISTANCE_FACTOR = 1.0E-4;
 const double OffsetSegmentGenerator::OFFSET_SEGMENT_SEPARATION_FACTOR = 1.0E-3;
 const double OffsetSegmentGenerator::INSIDE_TURN_VERTEX_SNAP_DISTANCE_FACTOR = 1.0E-3;
 const double OffsetSegmentGenerator::SIMPLIFY_FACTOR = 100.0;