From d8a173d5aeacfe1ec61058df038e275619d70eab Mon Sep 17 00:00:00 2001
From: Sandro Santilli <strk@kbt.io>
Date: Thu, 4 May 2023 15:09:43 +0200
Subject: [PATCH] Reduce artifacts in single-sided buffer output

Includes unit test

Fixes GH-665
---
 src/operation/buffer/BufferBuilder.cpp       | 87 +++++++++++++++++++-
 tests/unit/operation/buffer/BufferOpTest.cpp | 21 +++++
 2 files changed, 107 insertions(+), 1 deletion(-)

diff --git a/src/operation/buffer/BufferBuilder.cpp b/src/operation/buffer/BufferBuilder.cpp
index 500d729e24..be37585cd0 100644
--- a/src/operation/buffer/BufferBuilder.cpp
+++ b/src/operation/buffer/BufferBuilder.cpp
@@ -411,13 +411,98 @@ BufferBuilder::buffer(const Geometry* g, double distance)
     std::cerr << std::endl << edgeList << std::endl;
 #endif
 
+    std::vector<Edge*> *edges = &edgeList.getEdges();
+    std::vector<Edge*> reducedEdges; // only needed if bufParams.isSingleSided
+
+    if ( bufParams.isSingleSided() )
+    {
+#if GEOS_DEBUG
+        std::cerr << "Single-sided buffer was desired, we'll drop edges not being in the full buffer" << std::endl;
+#endif
+
+        // First, generate the two-sided buffer using a butt-cap.
+        BufferParameters modParams = bufParams;
+        modParams.setEndCapStyle(BufferParameters::CAP_FLAT);
+        modParams.setSingleSided(false);
+        BufferBuilder tmpBB(modParams);
+        std::unique_ptr<Geometry> buf(tmpBB.buffer(g, std::abs(distance)));
+        // Create MultiLineStrings from this polygon.
+        std::unique_ptr<Geometry> bufLineString(buf->getBoundary());
+
+#if GEOS_DEBUG
+        std::cerr << "Boundaries of full buffer:" << *bufLineString << std::endl;
+#endif
+
+        // Get linework of input geom
+        const Geometry *inputLineString = g;
+        std::unique_ptr<Geometry> inputPolygonBoundary;
+        if ( g->getDimension() > 1 )
+        {
+            inputPolygonBoundary = g->getBoundary();
+            inputLineString = inputPolygonBoundary.get();
+        }
+
+#if GEOS_DEBUG
+        std::cerr << "Input linework: " << *inputLineString << std::endl;
+#endif
+
+
+        // Union input line and full buffer line
+        std::unique_ptr<Geometry> usableLines = inputLineString->Union( bufLineString.get() );
+
+#if GEOS_DEBUG
+        std::cerr << "Usable lines: " << *usableLines << std::endl;
+#endif
+
+
+        for ( auto& e : *edges )
+        {
+            std::unique_ptr<Geometry> edgeGeom = geomFact->createLineString(*(e->getCoordinates()));
+
+            // NOTE: we use Snapped overlay because the actual buffer boundary might
+            //       diverge from original offset curves due to the addition of
+            //       intersections with caps and joins curves
+            using geos::operation::overlay::snap::SnapOverlayOp;
+            std::unique_ptr<Geometry> xset = SnapOverlayOp::overlayOp(
+                    *edgeGeom, *usableLines,
+                    overlayng::OverlayNG::INTERSECTION
+            );
+#if GEOS_DEBUG > 1
+            std::cerr << "Intersection: " << *xset << std::endl;
+#endif
+
+            // Drop edges not having linear intersection with
+            // the union of full buffer boundary and input line
+            if ( xset->getDimension() == 1 )
+            {
+#if GEOS_DEBUG > 1
+                std::cerr << "Covered edge: " << *e << std::endl;
+#endif
+                reducedEdges.push_back(e);
+            }
+            else
+            {
+#if GEOS_DEBUG > 1
+                std::cerr << "Non-covered edge: " << *e << std::endl;
+#endif
+            }
+        }
+        if ( reducedEdges.empty() )
+        {
+            // Or we could maybe be tolerant here
+            throw util::GEOSException("Unable to find single-sided Buffer Curves covered by boundary of full Buffer");
+        }
+        edges = &reducedEdges;
+
+    }
+
     std::unique_ptr<Geometry> resultGeom(nullptr);
     std::vector<std::unique_ptr<Geometry>> resultPolyList;
     std::vector<BufferSubgraph*> subgraphList;
 
     try {
         PlanarGraph graph(OverlayNodeFactory::instance());
-        graph.addEdges(edgeList.getEdges());
+        graph.addEdges(*edges);
 
         GEOS_CHECK_FOR_INTERRUPTS();
 
diff --git a/tests/unit/operation/buffer/BufferOpTest.cpp b/tests/unit/operation/buffer/BufferOpTest.cpp
index a40e85a31b..a1ff2e5573 100644
--- a/tests/unit/operation/buffer/BufferOpTest.cpp
+++ b/tests/unit/operation/buffer/BufferOpTest.cpp
@@ -549,4 +549,25 @@ void object::test<20>
     ensure( 0 == dynamic_cast<const geos::geom::Polygon*>(result1.get())->getNumInteriorRing() );
 }
 
+// Test for single-sided buffer
+// See https://github.com/libgeos/geos/issues/665
+template<>
+template<>
+void object::test<21>
+()
+{
+    using geos::operation::buffer::BufferOp;
+    using geos::operation::buffer::BufferParameters;
+
+    std::string wkt("LINESTRING (50 50, 150 150, 150 100, 150 0)");
+    GeomPtr geom(wktreader.read(wkt));
+
+    geos::operation::buffer::BufferParameters bp;
+    bp.setSingleSided(true);
+    geos::operation::buffer::BufferOp op(geom.get(), bp);
+
+    std::unique_ptr<Geometry> result = op.getResultGeometry(-21);
+    ensure_equals(int(result->getArea()), 5055);
+}
+
 } // namespace tut