locationtech / jts

The JTS Topology Suite is a Java library for creating and manipulating vector geometry.
Other
1.99k stars 443 forks source link

Problem with Intersection between two polygons #1046

Closed roarbra closed 7 months ago

roarbra commented 7 months ago

I'm having two polygons, where one is totally within the other. When I try to make a call outer.intersection(inner) I would assume to get the inner as a result. Instead I get a TopologyException. In the FAQ I saw something about using the PresicionModel of the GeometryFactory to minimize the scale, and so I implemented a tiny little program that would test that. The essential part is here:

      double scale = 1.0;
      GeometryFactory gf = new GeometryFactory(new PrecisionModel(scale));
      WKTReader reader = new WKTReader(gf);
      Geometry clip = reader.read(STR_CLIP);
      Geometry geom = reader.read(STR_GEOM);
      System.out.println("Intersects = " + geom.intersects(clip));
      Geometry intersection = geom.intersection(clip);
      System.out.println("Equals = " + intersection.equals(clip));

I've experienced with the scale from 1.0 - 0.0001 and here are the results: scale = 1.0

Intersects = true
org.locationtech.jts.geom.TopologyException: found non-noded intersection between LINESTRING ( 307695 6604259, 307725 6604250 ) and LINESTRING ( 307763 6604236, 307725 6604250 ) [ (307725.0, 6604250.0, NaN) ]
        at org.locationtech.jts.noding.FastNodingValidator.checkValid(FastNodingValidator.java:139)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:80)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:45)
        at org.locationtech.jts.operation.overlay.OverlayOp.computeOverlay(OverlayOp.java:229)
        at org.locationtech.jts.operation.overlay.OverlayOp.getResultGeometry(OverlayOp.java:181)
        at org.locationtech.jts.operation.overlay.OverlayOp.overlayOp(OverlayOp.java:84)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.getResultGeometry(SnapIfNeededOverlayOp.java:75)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.overlayOp(SnapIfNeededOverlayOp.java:37)
        at org.locationtech.jts.geom.GeometryOverlay.overlay(GeometryOverlay.java:76)
        at org.locationtech.jts.geom.GeometryOverlay.intersection(GeometryOverlay.java:119)
        at org.locationtech.jts.geom.Geometry.intersection(Geometry.java:1348)
        at app.JTSInterseciton.main(JTSInterseciton.java:59)

scale = 0.1

Intersects = true
org.locationtech.jts.geom.TopologyException: found non-noded intersection between LINESTRING ( 307690 6604260, 307730 6604250 ) and LINESTRING ( 307760 6604240, 307730 6604250 ) [ (307730.0, 6604250.0, NaN) ]
        at org.locationtech.jts.noding.FastNodingValidator.checkValid(FastNodingValidator.java:139)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:80)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:45)
        at org.locationtech.jts.operation.overlay.OverlayOp.computeOverlay(OverlayOp.java:229)
        at org.locationtech.jts.operation.overlay.OverlayOp.getResultGeometry(OverlayOp.java:181)
        at org.locationtech.jts.operation.overlay.OverlayOp.overlayOp(OverlayOp.java:84)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.getResultGeometry(SnapIfNeededOverlayOp.java:75)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.overlayOp(SnapIfNeededOverlayOp.java:37)
        at org.locationtech.jts.geom.GeometryOverlay.overlay(GeometryOverlay.java:76)
        at org.locationtech.jts.geom.GeometryOverlay.intersection(GeometryOverlay.java:119)
        at org.locationtech.jts.geom.Geometry.intersection(Geometry.java:1348)
        at app.JTSInterseciton.main(JTSInterseciton.java:59)

scale = 0.01

Intersects = true
Equals = false

scale = 0.001

Intersects = true
org.locationtech.jts.geom.TopologyException: found non-noded intersection between LINESTRING ( 304000 6601000, 304000 6602000 ) and LINESTRING ( 305000 6602000, 304000 6602000 ) [ (304000.0, 6602000.0, NaN) ]
        at org.locationtech.jts.noding.FastNodingValidator.checkValid(FastNodingValidator.java:139)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:80)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:45)
        at org.locationtech.jts.operation.overlay.OverlayOp.computeOverlay(OverlayOp.java:229)
        at org.locationtech.jts.operation.overlay.OverlayOp.getResultGeometry(OverlayOp.java:181)
        at org.locationtech.jts.operation.overlay.OverlayOp.overlayOp(OverlayOp.java:84)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.getResultGeometry(SnapIfNeededOverlayOp.java:75)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.overlayOp(SnapIfNeededOverlayOp.java:37)
        at org.locationtech.jts.geom.GeometryOverlay.overlay(GeometryOverlay.java:76)
        at org.locationtech.jts.geom.GeometryOverlay.intersection(GeometryOverlay.java:119)
        at org.locationtech.jts.geom.Geometry.intersection(Geometry.java:1348)
        at app.JTSInterseciton.main(JTSInterseciton.java:59)

scale = 0.0001

Intersects = true
Equals = true

This isn't what I would call robust. The whole program with the polygons is available as a gist.

dr-jts commented 7 months ago

The more recent OverlayNG algorithm is more robust. See the doc for how to enable it via a command-line switch, or use the OverlayNGRobust entry point directly.

dr-jts commented 7 months ago

Also, are you sure you understand how the PrecisionModel scale factor works? It is the reciprocal of the grid size. So a scale factor of 1 rounds to integers, and a scale factor of 0.1 rounds to the nearest 10 (not to the first decimal place). Using scale factors < 1 is unusual.

roarbra commented 7 months ago

Also, are you sure you understand how the PrecisionModel scale factor works? It is the reciprocal of the grid size. So a scale factor of 1 rounds to integers, and a scale factor of 0.1 rounds to the nearest 10 (not to the first decimal place). Using scale factors < 1 is unusual.

Yes, that was what I understood from the docs. And you can see it from the exception message as well.

roarbra commented 7 months ago

The more recent OverlayNG algorithm is more robust. See the doc for how to enable it via a command-line switch, or use the OverlayNGRobust entry point directly.

Yes, this seems to work a lot better. I will try to set that command-line switch in the startup of GeoServer. In which version was this added?

dr-jts commented 7 months ago

In which version was this added?

1.18

roarbra commented 7 months ago

Solved by setting -Djts.overlay=ng in Geoserver startup.