Port dr-jts/jts@2a385fe Change line output behaviour to preserve noding

include/geos/index/kdtree/KdNode.h

@@ -29,20 +29,20 @@ class GEOS_DLL KdNode {
 private:
 
     geom::Coordinate p;
-    const void* data;
+    void* data;
     KdNode* left;
     KdNode* right;
     std::size_t count;
 
 public:
 
-    KdNode(double p_x, double p_y, const void* p_data);
-    KdNode(const geom::Coordinate& p_p, const void* p_data);
+    KdNode(double p_x, double p_y, void* p_data);
+    KdNode(const geom::Coordinate& p_p, void* p_data);
 
     double getX() { return p.x; }
     double getY() { return p.y; }
     const geom::Coordinate& getCoordinate() { return p; }
-    const void* getData() { return data; }
+    void* getData() { return data; }
     KdNode* getLeft() { return left; }
     KdNode* getRight() { return right; }
     void increment() { count++; }

include/geos/index/kdtree/KdTree.h

@@ -62,7 +62,7 @@ class GEOS_DLL KdTree {
     double tolerance;
 
     KdNode* findBestMatchNode(const geom::Coordinate& p);
-    KdNode* insertExact(const geom::Coordinate& p, const void* data);
+    KdNode* insertExact(const geom::Coordinate& p, void* data);
 
     void queryNode(KdNode* currentNode, const geom::Envelope& queryEnv, bool odd, KdNodeVisitor& visitor);
     KdNode* queryNodePoint(KdNode* currentNode, const geom::Coordinate& queryPt, bool odd);
@@ -71,7 +71,7 @@ class GEOS_DLL KdTree {
     * Create a node on a locally managed deque to allow easy
     * disposal and hopefully faster allocation as well.
     */
-    KdNode* createNode(const geom::Coordinate& p, const void* data);
+    KdNode* createNode(const geom::Coordinate& p, void* data);
 
 
     /**
@@ -158,7 +158,7 @@ class GEOS_DLL KdTree {
     * Inserts a new point in the kd-tree.
     */
     KdNode* insert(const geom::Coordinate& p);
-    KdNode* insert(const geom::Coordinate& p, const void* data);
+    KdNode* insert(const geom::Coordinate& p, void* data);
 
     /**
     * Performs a range search of the points in the index and visits all nodes found.

include/geos/noding/snapround/HotPixel.h

@@ -76,6 +76,7 @@ class GEOS_DLL HotPixel {
     geom::Coordinate ptHot;
     geom::Coordinate originalPt;
 
+    bool hpIsNode;
     double scaleFactor;
 
     double minx;
@@ -134,6 +135,16 @@ class GEOS_DLL HotPixel {
     bool intersects(const geom::Coordinate& p0,
                     const geom::Coordinate& p1) const;
 
+    /**
+    * Tests whether a coordinate lies in (intersects) this hot pixel.
+    *
+    * @param p the coordinate to test
+    * @return true if the coordinate intersects this hot pixel
+    */
+    bool intersects(const geom::Coordinate& p) const;
+
+    bool isNode() const { return hpIsNode; };
+    void setToNode() { hpIsNode = true; };
 
     std::ostream& operator<< (std::ostream& os);
 };

include/geos/noding/snapround/HotPixelIndex.h

@@ -68,15 +68,24 @@ class GEOS_DLL HotPixelIndex {
 
     /* methods */
     geom::Coordinate round(const geom::Coordinate& c);
-    const HotPixel* find(const geom::Coordinate& pixelPt);
+    HotPixel* find(const geom::Coordinate& pixelPt);
 
 public:
 
     HotPixelIndex(const geom::PrecisionModel* p_pm);
-    const HotPixel* add(const geom::Coordinate& pt);
+    HotPixel* add(const geom::Coordinate& pt);
     void add(const geom::CoordinateSequence* pts);
     void add(const std::vector<geom::Coordinate>& pts);
-    void query(const geom::Coordinate& p0, const geom::Coordinate& p1, index::kdtree::KdNodeVisitor& visitor);
+    void addNodes(const geom::CoordinateSequence* pts);
+    void addNodes(const std::vector<geom::Coordinate>& pts);
+
+    /**
+    * Visits all the hot pixels which may intersect a segment (p0-p1).
+    * The visitor must determine whether each hot pixel actually intersects
+    * the segment.
+    */
+    void query(const geom::Coordinate& p0, const geom::Coordinate& p1,
+               index::kdtree::KdNodeVisitor& visitor);
 
 };
 

include/geos/noding/snapround/SnapRoundingIntersectionAdder.h

@@ -46,10 +46,29 @@ namespace geos {
 namespace noding { // geos::noding
 namespace snapround { // geos::noding::snapround
 
+/**
+ * Finds intersections between line segments which will be snap-rounded,
+ * and adds them as nodes to the segments.
+ *
+ * Intersections are detected and computed using full precision.
+ * Snapping takes place in a subsequent phase.
+ *
+ * The intersection points are recorded, so that HotPixels can be created for them.
+ *
+ * To avoid robustness issues with vertices which lie very close to line segments
+ * a heuristic is used:
+ * nodes are created if a vertex lies within a tolerance distance
+ * of the interior of a segment.
+ * The tolerance distance is chosen to be significantly below the snap-rounding grid size.
+ * This has empirically proven to eliminate noding failures.
+ */
 class GEOS_DLL SnapRoundingIntersectionAdder: public SegmentIntersector { // implements SegmentIntersector
 
 private:
-
+    /**
+    * The division factor used to determine
+    * nearness distance tolerance for interior intersection detection.
+    */
     static constexpr int NEARNESS_FACTOR = 100;
 
     algorithm::LineIntersector li;
@@ -59,10 +78,10 @@ class GEOS_DLL SnapRoundingIntersectionAdder: public SegmentIntersector { // imp
 
     /**
     * If an endpoint of one segment is near
-    * the <i>interior</i> of the other segment, add it as an intersection.
+    * the interior of the other segment, add it as an intersection.
     * EXCEPT if the endpoint is also close to a segment endpoint
     * (since this can introduce "zigs" in the linework).
-    * <p>
+    *
     * This resolves situations where
     * a segment A endpoint is extremely close to another segment B,
     * but is not quite crossing.  Due to robustness issues
@@ -84,7 +103,7 @@ class GEOS_DLL SnapRoundingIntersectionAdder: public SegmentIntersector { // imp
     * This method is called by clients
     * of the {@link SegmentIntersector} class to process
     * intersections for two segments of the {@link SegmentString}s being intersected.
-    * Note that some clients (such as <code>MonotoneChain</code>s) may optimize away
+    * Note that some clients (such as MonotoneChains) may optimize away
     * this call for segment pairs which they have determined do not intersect
     * (e.g. by an disjoint envelope test).
     */

include/geos/noding/snapround/SnapRoundingNoder.h

@@ -41,6 +41,33 @@ namespace geos {
 namespace noding { // geos::noding
 namespace snapround { // geos::noding::snapround
 
+/**
+ * Uses Snap Rounding to compute a rounded,
+ * fully noded arrangement from a set of {@link SegmentString}s,
+ * in a performant way, and avoiding unnecessary noding.
+ *
+ * Implements the Snap Rounding technique described in
+ * the papers by Hobby, Guibas & Marimont, and Goodrich et al.
+ * Snap Rounding enforces that all output vertices lie on a uniform grid,
+ * which is determined by the provided {@link PrecisionModel}.
+ *
+ * Input vertices do not have to be rounded to the grid beforehand;
+ * this is done during the snap-rounding process.
+ * In fact, rounding cannot be done a priori,
+ * since rounding vertices by themselves can distort the rounded topology
+ * of the arrangement (i.e. by moving segments away from hot pixels
+ * that would otherwise intersect them, or by moving vertices
+ * across segments).
+ *
+ * To minimize the number of introduced nodes,
+ * the Snap-Rounding Noder avoids creating nodes
+ * at edge vertices if there is no intersection or snap at that location.
+ * However, if two different input edges contain identical segments,
+ * each of the segment vertices will be noded.
+ * This still provides fully-noded output.
+ * This is the same behaviour provided by other noders,
+ * such as {@link MCIndexNoder} and {@link SnappingNoder}.
+ */
 class GEOS_DLL SnapRoundingNoder : public Noder {
 
 private:
@@ -51,12 +78,23 @@ class GEOS_DLL SnapRoundingNoder : public Noder {
     std::vector<SegmentString*> snappedResult;
 
     // Methods
-    void snapRound(const std::vector<SegmentString*>& inputSegStrings, std::vector<SegmentString*>& resultNodedSegments);
+    void snapRound(std::vector<SegmentString*>& inputSegStrings, std::vector<SegmentString*>& resultNodedSegments);
 
-    static void createNodedStrings(const std::vector<SegmentString*>& segStrings,
-        std::vector<SegmentString*>& nodedStrings);
+    /**
+    * Creates HotPixels for each vertex in the input segStrings.
+    * The HotPixels are not marked as nodes, since they will
+    * only be nodes in the final line arrangement
+    * if they interact with other segments (or they are already
+    * created as intersection nodes).
+    */
+    void addVertexPixels(std::vector<SegmentString*>& segStrings);
 
-    void addVertexPixels(const std::vector<SegmentString*>& segStrings);
+    /**
+    * Detects interior intersections in the collection of {@link SegmentString}s,
+    * and adds nodes for them to the segment strings.
+    * Also creates HotPixel nodes for the intersection points.
+    */
+    void addIntersectionPixels(std::vector<SegmentString*>& segStrings);
 
     void round(const geom::Coordinate& pt, geom::Coordinate& ptOut);
 
@@ -69,25 +107,15 @@ class GEOS_DLL SnapRoundingNoder : public Noder {
     */
     std::unique_ptr<std::vector<geom::Coordinate>> round(const std::vector<geom::Coordinate>& pts);
 
-    /**
-    * Computes all interior intersections in the collection of {@link SegmentString}s,
-    * and returns their {@link Coordinate}s.
-    *
-    * Also adds the intersection nodes to the segments.
-    *
-    * @return a list of Coordinates for the intersections
-    */
-    std::unique_ptr<std::vector<geom::Coordinate>> findInteriorIntersections(std::vector<SegmentString*>& inputSS);
-
     /**
     * Computes new segment strings which are rounded and contain
-    * any intersections added as a result of snapping segments to snap points (hot pixels).
+    * intersections added as a result of snapping segments to snap points (hot pixels).
     *
     * @param segStrings segments to snap
     * @return the snapped segment strings
     */
     void computeSnaps(const std::vector<SegmentString*>& segStrings, std::vector<SegmentString*>& snapped);
-    NodedSegmentString* computeSnaps(NodedSegmentString* ss);
+    NodedSegmentString* computeSegmentSnaps(NodedSegmentString* ss);
 
     /**
     * Snaps a segment in a segmentString to HotPixels that it intersects.
@@ -99,6 +127,13 @@ class GEOS_DLL SnapRoundingNoder : public Noder {
     */
     void snapSegment(geom::Coordinate& p0, geom::Coordinate& p1, NodedSegmentString* ss, size_t segIndex);
 
+    /**
+    * Add nodes for any vertices in hot pixels that were
+    * added as nodes during segment noding.
+    */
+    void addVertexNodeSnaps(NodedSegmentString* ss);
+
+    void snapVertexNode(const geom::Coordinate& p0, NodedSegmentString* ss, size_t segIndex);
 
 public:
 
@@ -112,6 +147,10 @@ class GEOS_DLL SnapRoundingNoder : public Noder {
     */
     std::vector<SegmentString*>* getNodedSubstrings() const override;
 
+    /**
+    * Computes the nodes in the snap-rounding line arrangement.
+    * The nodes are added to the {@link NodedSegmentString}s provided as the input.
+    */
     void computeNodes(std::vector<SegmentString*>* inputSegStrings) override; //override
 
 };

include/geos/operation/overlayng/LineBuilder.h

@@ -109,6 +109,9 @@ class GEOS_DLL LineBuilder {
     geom::Location effectiveLocation(int geomIndex, const OverlayLabel* lbl) const;
 
     void addResultLines();
+    void addResultLinesMerged();
+
+    std::unique_ptr<geom::LineString> toLine(OverlayEdge* edge);
 
     void addResultLinesForNodes();
 

include/geos/operation/overlayng/OverlayLabel.h

@@ -235,7 +235,6 @@ class GEOS_DLL OverlayLabel {
     Location getLocation(int index, int position, bool isForward) const;
     bool hasSides(int index) const;
 
-    OverlayLabel copyFlip() const;
     OverlayLabel copy() const;
 
 

src/index/kdtree/KdNode.cpp

@@ -20,14 +20,14 @@ namespace geos {
 namespace index { // geos.index
 namespace kdtree { // geos.index.kdtree
 
-KdNode::KdNode(double p_x, double p_y, const void* p_data) :
+KdNode::KdNode(double p_x, double p_y, void* p_data) :
     p(p_x, p_y),
     data(p_data),
     left(nullptr),
     right(nullptr),
     count(1) {}
 
-KdNode::KdNode(const Coordinate& p_p, const void* p_data) :
+KdNode::KdNode(const Coordinate& p_p, void* p_data) :
     p(p_p),
     data(p_data),
     left(nullptr),

src/index/kdtree/KdTree.cpp

@@ -51,7 +51,7 @@ KdTree::toCoordinates(std::vector<KdNode*>& kdnodes, bool includeRepeated)
 
 /*private*/
 KdNode*
-KdTree::createNode(const Coordinate& p, const void* data)
+KdTree::createNode(const Coordinate& p, void* data)
 {
     auto it = nodeQue.emplace(nodeQue.end(), p, data);
     return &(*it);
@@ -66,7 +66,7 @@ KdTree::insert(const Coordinate& p)
 
 /*public*/
 KdNode*
-KdTree::insert(const Coordinate& p, const void* data)
+KdTree::insert(const Coordinate& p, void* data)
 {
     if (root == nullptr) {
         root = createNode(p, data);
@@ -98,7 +98,7 @@ KdTree::findBestMatchNode(const Coordinate& p) {
 }
 
 KdNode*
-KdTree::insertExact(const geom::Coordinate& p, const void* data)
+KdTree::insertExact(const geom::Coordinate& p, void* data)
 {
     KdNode* currentNode = root;
     KdNode* leafNode = root;

src/noding/snapround/HotPixel.cpp

@@ -16,8 +16,9 @@
  *
  **********************************************************************/
 
-#include <geos/algorithm/CGAlgorithmsDD.h>
 #include <geos/noding/snapround/HotPixel.h>
+
+#include <geos/algorithm/CGAlgorithmsDD.h>
 #include <geos/noding/NodedSegmentString.h>
 #include <geos/algorithm/LineIntersector.h>
 #include <geos/geom/Coordinate.h>
@@ -39,10 +40,10 @@ namespace noding { // geos.noding
 namespace snapround { // geos.noding.snapround
 
 HotPixel::HotPixel(const Coordinate& newPt, double newScaleFactor)
-    :
-    ptHot(newPt),
-    originalPt(newPt),
-    scaleFactor(newScaleFactor)
+    : ptHot(newPt)
+    , originalPt(newPt)
+    , hpIsNode(false)
+    , scaleFactor(newScaleFactor)
 {
     if(scaleFactor <= 0.0) {
         throw util::IllegalArgumentException("Scale factor must be non-zero");
@@ -59,10 +60,27 @@ HotPixel::HotPixel(const Coordinate& newPt, double newScaleFactor)
 
 /*public*/
 const geom::Coordinate&
-HotPixel::getCoordinate() const {
+HotPixel::getCoordinate() const
+{
     return originalPt;
 }
 
+/* public */
+bool
+HotPixel::intersects(const Coordinate& p) const
+{
+    double x = scale(p.x);
+    double y = scale(p.y);
+    if (x >= maxx) return false;
+    // check Left side
+    if (x < minx) return false;
+    // check Top side
+    if (y >= maxy) return false;
+    // check Bottom side
+    if (y < miny) return false;
+    // finally
+    return true;
+}
 
 /*public*/
 bool