libgeos / geos

Geometry Engine, Open Source
https://libgeos.org
GNU Lesser General Public License v2.1
1.13k stars 346 forks source link

BUG: polygon intersection sensitive to winding order since GEOS 3.10 #827

Open brendan-ward opened 1 year ago

brendan-ward commented 1 year ago

First reported at Shapely #1767.

Given the attached WKB geometries: wkbs.zip

Intersection should produce a polygon output.

The first polygon has clockwise winding order, the second is counterclockwise. Intersection used to work as expected for GEOS 3.9, and only works properly for GEOS >= 3.10 if both polygons are oriented the same way.

# GEOS 3.12.0dev
geosop intersection -a /tmp/p1.wkb -b /tmp/p2.wkb -f wkt
# MULTIPOINT (-0.0232907037082979 43.207097646087036, -0.0570036 43.0674921)

# GEOS 3.10
geosop intersection -a /tmp/p1.wkb -b /tmp/p2.wkb -f wkt
# MULTIPOINT (-0.0232907037082979 43.207097646087036, -0.0570036 43.0674921)

# GEOS 3.9 (expected result)
geosop intersection -a /tmp/p1.wkb -b /tmp/p2.wkb -f wkt
# POLYGON ((0.009930363585773504 43.20709764608704, 0.009930363585773504 43.05895741058383, -0.0590476798926261 43.05895741058383, -0.0570036 43.0674921, -0.02329070370829787 43.20709764608704, 0.009930363585773504 43.20709764608704))

Reversing the first so they have the same winding order produces valid results:

# GEOS 3.12.0dev
geosop reverse -a /tmp/p1.wkb -f wkt | geosop intersection -a - -b /tmp/p2.wkb -f wkt
# POLYGON ((0.0099303635857735 43.207097646087036, 0.0099303635857735 43.05895741058383, -0.0590476798926261 43.05895741058383, -0.0570036 43.0674921, -0.0232907037082979 43.207097646087036, 0.0099303635857735 43.207097646087036))

Note that testing with WKT does not work, because it slightly modifies the precision and produces expected outputs rather than the errors above.

dr-jts commented 1 year ago

Confirmed this is also a bug in JTS. I'll investigate there.

dbaston commented 1 year ago

Here is an XML test to help with diagnosing the issue. geos-gh-827.xml.txt](https://github.com/libgeos/geos/files/10845933/issue-geos-gh-827.xml.txt) (added .txt so GitHub will allow it)

I looked into it a bit and found the breaking commit in GEOS is 3320d96a88a1d8492016c37f23983de164470670. In JTS, changing the corresponding object from a HashMap to a TreeMap also makes the test pass. I am suspicious that there is a bug in how the overlay labels are being propagated in OverlayLabeller, such that the order in which the edges are processed affects the final labeling. In particular, two edges are labeled interior-right or exterior-right depending on the order in which they are processed.

skygering commented 1 year ago

Is there any outlook on this being fixed? Thank you!

dr-jts commented 1 year ago

Possibly related to #741 ?