I am using Shapely in a python script to determine the intersection between LinearRings and symmetric difference between Polygons in 2D. I am running into what I assume are tolerance issues, and I cannot find any information in the Shapely reference regarding the matter. I will now present what I have done so far to try to determine what is going on.
from shapely.geometry import LinearRing, Polygon
ls_blue = LinearRing( [ [25.16,-88.42], [26.24,-87.34], [26.0,-88.42] ] )
ls_red = LinearRing( [ [24.04, -89.54], [27.32, -86.26], [26.0, -89] ] )
inter = ls_blue.intersection(ls_red)
-inter shows no line intersection on the segment of interest (the nearly overlapping segments seen in the image below, [25.16,-88.42], [26.24,-87.34]) and all expected results
p_blue = Polygon( [ [25.16,-88.42], [26.24,-87.34], [26.0,-88.42] ] )
p_red = Polygon( [ [ [24.04, -89.54], [27.32, -86.26], [26.0, -89] ] ] )
p_xor = p_blue.symmetric_difference(p_red)
-segment of interest (the nearly overlapping segments seen in the image below, [25.16,-88.42], [26.24,-87.34]) is included in one of the p_xor polygons and all results are expected
Please note that a LinearRing is the border of the Polygon only.
A diagram of the sample LinearRings/Polygons: https://i.sstatic.net/jcnFa.jpg
In this case, the symmetric difference polygons and intersections are computed as expected.
However, this example is actually a subset of some larger point series. In these larger point series the XOR polygon and line intersections are not computed as expected.
The two point series are as follows:
c_upper = Polygon([[ -38.68 , -116.82 ], [ -37.6 , -116.79 ], [ 45.72 , -116.79 ], [ 45.96 , -116.59 ], [ 46.26 , -115.51 ], [ 46.26 , -74.39 ], [ 45.75 , -73.31 ], [ 45.72 , -73.27 ], [ 44.53 , -74.39 ], [ 43.55 , -75.26 ], [ 42.47 , -76.14 ], [ 41.39 , -76.59 ], [ 40.27 , -77.64 ], [ 39.22 , -78.51 ], [ 37.06 , -79.77 ], [ 35.98 , -80.67 ], [ 33.81 , -81.93 ], [ 32.73 , -82.79 ], [ 31.65 , -83.2 ], [ 29.49 , -85.0 ], [ 28.4 , -85.62 ], [ 27.32 , -86.26 ], [ 24.04 , -89.54 ], [ 23.03 , -90.62 ], [ 22.6 , -91.7 ], [ 21.76 , -92.78 ], [ 20.86 , -93.87 ], [ 20.67 , -94.95 ], [ 19.71 , -96.03 ], [ 19.52 , -97.11 ], [ 18.67 , -98.16 ], [ 18.63 , -98.19 ], [ 18.32 , -99.28 ], [ 17.49 , -100.36 ], [ 17.01 , -101.44 ], [ 16.5 , -102.28 ], [ 16.34 , -102.52 ], [ 15.51 , -103.6 ], [ 15.16 , -104.69 ], [ 14.34 , -105.67 ], [ 14.24 , -105.77 ], [ 13.97 , -106.85 ], [ 13.25 , -107.72 ], [ 13.1 , -107.93 ], [ 12.14 , -109.01 ], [ 11.5 , -110.1 ], [ 11.09 , -110.57 ], [ 10.43 , -111.18 ], [ 8.93 , -112.29 ], [ 7.84 , -112.86 ], [ 6.82 , -113.34 ], [ 6.76 , -113.38 ], [ 5.68 , -113.58 ], [ 4.6 , -113.67 ], [ 3.52 , -113.38 ], [ 2.43 , -112.97 ], [ 1.35 , -112.36 ], [ 1.24 , -112.26 ], [ 0.27 , -111.55 ], [ -0.32 , -111.18 ], [ -0.81 , -110.81 ], [ -1.89 , -109.66 ], [ -2.43 , -109.01 ], [ -3.01 , -107.93 ], [ -4.92 , -105.77 ], [ -5.74 , -104.69 ], [ -6.19 , -103.6 ], [ -7.09 , -102.52 ], [ -7.91 , -101.44 ], [ -8.35 , -100.36 ], [ -9.26 , -99.28 ], [ -9.92 , -98.19 ], [ -10.52 , -97.11 ], [ -11.54 , -96.03 ], [ -12.44 , -94.95 ], [ -13.22 , -93.87 ], [ -13.83 , -92.78 ], [ -14.88 , -91.74 ], [ -15.96 , -91.03 ], [ -17.04 , -90.29 ], [ -17.94 , -89.54 ], [ -19.21 , -88.42 ], [ -20.99 , -87.37 ], [ -21.37 , -87.16 ], [ -22.45 , -86.26 ], [ -23.53 , -85.81 ], [ -24.62 , -85.06 ], [ -25.7 , -84.1 ], [ -26.78 , -83.78 ], [ -27.74 , -83.05 ], [ -28.94 , -82.0 ], [ -30.03 , -81.69 ], [ -31.11 , -80.79 ], [ -32.19 , -79.83 ], [ -33.27 , -79.64 ], [ -34.35 , -78.68 ], [ -35.44 , -77.79 ], [ -36.52 , -77.37 ], [ -37.6 , -76.46 ], [ -38.68 , -75.51 ], [ -39.76 , -75.06 ], [ -40.08 , -75.47 ], [ -40.31 , -76.55 ], [ -40.31 , -77.64 ], [ -40.29 , -78.72 ], [ -40.29 , -114.42 ], [ -40.31 , -115.51 ], [ -39.76 , -116.56 ], [ -39.72 , -116.59 ]])
c_lower = Polygon([[ -38.68 , -116.82 ], [ -37.6 , -116.79 ], [ 45.72 , -116.79 ], [ 45.96 , -116.59 ], [ 46.26 , -115.51 ], [ 46.26 , -73.31 ], [ 45.72 , -72.26 ], [ 44.63 , -73.27 ], [ 43.52 , -74.39 ], [ 42.47 , -75.38 ], [ 41.39 , -76.34 ], [ 40.31 , -77.22 ], [ 39.22 , -77.67 ], [ 38.11 , -78.72 ], [ 37.06 , -79.59 ], [ 36.68 , -79.8 ], [ 34.9 , -80.85 ], [ 33.81 , -81.81 ], [ 32.73 , -82.63 ], [ 31.65 , -83.08 ], [ 30.46 , -84.13 ], [ 29.49 , -84.94 ], [ 29.04 , -85.21 ], [ 27.32 , -86.33 ], [ 26.24 , -87.34 ], [ 25.16 , -88.42 ], [ 24.11 , -89.54 ], [ 23.66 , -90.62 ], [ 22.9 , -91.7 ], [ 21.95 , -92.78 ], [ 21.75 , -93.87 ], [ 20.8 , -94.95 ], [ 20.61 , -96.03 ], [ 19.71 , -97.11 ], [ 19.42 , -98.19 ], [ 18.57 , -99.28 ], [ 18.15 , -100.36 ], [ 17.49 , -101.44 ], [ 16.66 , -102.52 ], [ 16.29 , -103.6 ], [ 15.39 , -104.69 ], [ 15.13 , -105.77 ], [ 14.24 , -106.85 ], [ 13.29 , -107.93 ], [ 12.86 , -109.01 ], [ 11.95 , -110.1 ], [ 11.09 , -111.01 ], [ 10.01 , -111.97 ], [ 8.93 , -112.55 ], [ 7.84 , -113.17 ], [ 6.76 , -113.31 ], [ 5.68 , -113.57 ], [ 4.6 , -113.57 ], [ 3.52 , -113.25 ], [ 2.43 , -112.9 ], [ 0.27 , -111.93 ], [ -0.63 , -111.18 ], [ -0.81 , -111.01 ], [ -2.81 , -109.01 ], [ -2.98 , -108.83 ], [ -3.73 , -107.93 ], [ -4.47 , -106.85 ], [ -5.14 , -105.73 ], [ -6.06 , -104.69 ], [ -6.95 , -103.6 ], [ -7.27 , -102.52 ], [ -8.29 , -101.44 ], [ -9.12 , -100.36 ], [ -9.44 , -99.28 ], [ -10.4 , -98.19 ], [ -11.22 , -97.11 ], [ -11.6 , -96.03 ], [ -12.62 , -94.95 ], [ -13.8 , -93.6 ], [ -14.39 , -92.78 ], [ -14.88 , -92.17 ], [ -15.32 , -91.7 ], [ -15.96 , -91.1 ], [ -16.54 , -90.62 ], [ -18.12 , -89.38 ], [ -19.17 , -88.46 ], [ -20.29 , -87.41 ], [ -21.37 , -87.1 ], [ -22.45 , -86.2 ], [ -23.53 , -85.24 ], [ -24.62 , -84.94 ], [ -25.59 , -84.13 ], [ -26.78 , -83.08 ], [ -27.86 , -82.82 ], [ -30.03 , -80.92 ], [ -30.13 , -80.88 ], [ -31.11 , -80.61 ], [ -32.08 , -79.8 ], [ -33.27 , -78.69 ], [ -34.35 , -78.3 ], [ -35.44 , -77.48 ], [ -36.52 , -76.52 ], [ -37.6 , -75.51 ], [ -38.68 , -74.96 ], [ -39.39 , -74.39 ], [ -39.76 , -74.12 ], [ -40.01 , -74.39 ], [ -40.31 , -75.47 ], [ -40.29 , -77.64 ], [ -40.29 , -114.42 ], [ -40.31 , -115.51 ], [ -39.76 , -116.56 ], [ -39.72 , -116.59 ]])
The segment of interest ([25.16,-88.42], [26.24,-87.34]) is not included in the XOR polygon or as a line intersection when the calculation is done the same way as before, although both of its endpoints are part of the intersection as points.
My algorithm relies on the idea that segments that do not intersect will be part of a polygon produced by the symmetric_difference function, so obviously the above case produces undesirable results and I find this rather disconcerting.
My question is: What is causing the discrepancy between intersection and symmetric difference operations between the cases?
PS: If more images are required for visualization purposes, let me know and I will create them.
This is an issue of precision, since the nodes of red are not part of blue. Interpolating a point on a line is not computationally exact. You can get around this by snapping one geometry to another, within an error of precision. See shapely's snap function to do this.
from shapely.ops import snap
print(ls_blue.intersection(snap(ls_red, ls_blue, 1e-8)))
# LINESTRING (25.16 -88.42, 26.24 -87.34)
# or try
p_blue.symmetric_difference(snap(p_red, p_blue, 1e-8))