Shapely: Multipolygon.intersection fails
I asked this question on StackOverflow, but I had only two views in a 5-6 hours, so I thought I would try here. Please let me know if it is a better place.
The following code creates different versions of a Multipolygon with shapely, then it calculates their intersection with another Multipolygon. The coordinates are a long list taken by the real program before the crash that bothers me.
Why does multipolygon1.intersection(multipolygon2) fail, but all the other variations of multipolygon1 can successfully intersect with multipolygon2?
And more important, since the polygons were all created by shapely's union and difference operators, how do I avoid this problem? Should I always expect that shapely will create polygons that don't work?
Here are the details
The variations of multipolygon1 are 3:
multipolygon1_bufferedmultipolygon1_without_last_polygonmultipolygon1_with_first_3_and_last
At first I thought that multipolygon1 was an invalid shapely Multipolygon, so I tried to buffer it, and it worked. And also the intersection with the buffered version worked. So multipolygon1 seems to be valid.
Then I removed the last polygon and the intersection worked. I thought I had found the problem: I would need to fix the last polygon.
Then I created another version with the first 3 and the last one, and it worked. So there is nothing to fix on the last one.
Conclusion: multipolygon1 is valid because can be buffered, all the polygons without the last one work, and also the last one with some other polygons work, so I don't know where is the problem.
I noticed that some polygons are very small. For example polygons1[12].area = 1.4210854715202004e-14, but they all were created by a sequence of union and difference. Shapely made them, so they should be valid.
Here is the code
import shapely.geometry
mp1 = [[(22.13354026441686, 0.0),
(0.0, 0.0),
(0.0, 13.242081447963805),
(0.056818181818172775, 13.278846153846153),
(17.25, 13.278846153846153),
(17.25, 11.853846153846154),
(17.375, 11.853846153846154),
(17.375, 13.278846153846153),
(18.8, 13.278846153846153),
(18.8, 13.403846153846153),
(17.375, 13.403846153846153),
(17.375, 24.48472850678734),
(18.551333743585833, 25.24588563498994),
(18.483427333816174, 25.350831904633957),
(17.375, 24.63361421804702),
(17.375, 25.95384615384616),
(17.25, 25.95384615384616),
(17.25, 24.552731865105844),
(0.019903900780485984, 13.403846153846153),
(0.0, 13.403846153846153),
(0.0, 15.285503429011813),
(0.12499999999999992, 15.366385781952989),
(0.125, 14.046153846153844),
(0.25, 14.046153846153844),
(0.24999999999999992, 15.447268134894166),
(17.48009609921951, 26.596153846153857),
(18.8, 26.596153846153857),
(18.8, 26.721153846153857),
(17.67327791740133, 26.721153846153857),
(18.619240153355495, 27.33324705765361),
(18.551333743585836, 27.43819332729763),
(18.25, 27.243212669683263),
(18.25, 28.601114288740316),
(35.480096099219516, 39.75),
(36.8, 39.75),
(36.8, 39.875),
(35.673277917401336, 39.875),
(35.866459735583156, 40.0),
(40.0, 40.0),
(40.0, 28.95022624434389),
(39.94318181818183, 28.913461538461544),
(22.75, 28.913461538461533),
(22.75, 30.338461538461534),
(22.625, 30.338461538461534),
(22.625, 28.913461538461533),
(21.2, 28.913461538461533),
(21.2, 28.788461538461533),
(22.625, 28.788461538461533),
(22.625, 17.70757918552036),
(21.448666256414164, 16.94642205731776),
(21.516572666183823, 16.841475787673744),
(22.625, 17.558693474260682),
(22.625, 16.238461538461536),
(22.75, 16.238461538461536),
(22.75, 17.63957582720186),
(39.980096099219516, 28.788461538461544),
(40.0, 28.788461538461544),
(40.0, 26.906804263295886),
(39.875, 26.825921910354708),
(39.875, 28.146153846153847),
(39.75, 28.146153846153847),
(39.75, 26.74503955741353),
(22.5199039007805, 15.596153846153847),
(21.2, 15.596153846153847),
(21.2, 15.471153846153847),
(22.32672208259868, 15.471153846153847),
(21.38075984664451, 14.859060634654089),
(21.448666256414167, 14.75411436501007),
(21.75, 14.949095022624432),
(21.75, 13.59119340356747),
(21.71308571896216, 13.56730769230769),
(21.2, 13.56730769230769),
(21.2, 13.44230769230769),
(21.519903900780342, 13.44230769230769),
(4.519903900780493, 2.4423076923076925),
(3.1999999999999993, 2.4423076923076925),
(3.1999999999999993, 2.3173076923076925),
(4.3267220825986765, 2.3173076923076925),
(3.3807598466445077, 1.7052144808079315),
(3.4486662564141675, 1.6002682111639128),
(4.55681818181818, 2.3173076923076925),
(21.75, 2.3173076923076925),
(21.75, 1.795248868778278),
(21.448666256414167, 1.600268211163916),
(21.516572666183826, 1.495321941519897),
(21.75, 1.6463631575185977),
(21.75, 0.892307692307693),
(21.875, 0.892307692307693),
(21.875, 1.7272455104597744),
(22.625, 2.2125396281068337),
(22.625, 0.8923076923076909),
(22.75, 0.8923076923076909),
(22.75, 2.2934219810480094),
(22.786914281037692, 2.3173076923076925),
(23.3, 2.3173076923076925),
(23.3, 2.4423076923076925),
(22.98009609921951, 2.4423076923076925),
(39.98009609921951, 13.442307692307693),
(40.0, 13.442307692307693),
(40.0, 11.560650417142032),
(39.875, 11.479768064200856),
(39.875, 12.8),
(39.75, 12.8),
(39.75, 11.39888571125968),
(22.519903900780495, 0.25),
(21.2, 0.25),
(21.2, 0.125),
(22.32672208259868, 0.125),
(22.13354026441686, 0.0)],
[(0.0, 15.434389140271493),
(0.0, 26.596153846153847),
(0.12499999999999921, 26.596153846153847),
(0.1249999999999999, 15.51527149321267),
(0.0, 15.434389140271493)],
[(0.0, 26.721153846153847),
(0.0, 40.0),
(18.125, 40.0),
(18.125, 39.875),
(16.7, 39.875),
(16.7, 39.75),
(18.125, 39.75),
(18.125, 28.66911764705882),
(16.948666256414167, 27.90796051885622),
(17.016572666183826, 27.803014249212204),
(18.125, 28.52023193579914),
(18.125, 27.199999999999996),
(18.183216783216764, 27.199999999999996),
(17.443181818181824, 26.721153846153857),
(0.2499999999999992, 26.721153846153847),
(0.2499999999999991, 28.146153846153847),
(0.12499999999999911, 28.146153846153847),
(0.1249999999999992, 26.721153846153847),
(0.0, 26.721153846153847)],
[(18.25, 40.0),
(35.63636363636365, 40.0),
(35.44318181818183, 39.875),
(18.25, 39.875),
(18.25, 40.0)],
[(21.816783216783378, 14.992307692307794),
(22.556818181818187, 15.471153846153847),
(39.75, 15.471153846153847),
(39.75, 14.046153846153846),
(39.875, 14.046153846153846),
(39.875, 15.471153846153847),
(40.0, 15.471153846153847),
(40.0, 13.604072398190045),
(39.94318181818182, 13.567307692307693),
(22.75, 13.56730769230769),
(22.75, 14.089366515837206),
(23.051333743585833, 14.28434717345157),
(22.983427333816174, 14.389293443095589),
(22.75, 14.238252227096886),
(22.75, 14.99230769230769),
(22.625, 14.99230769230769),
(22.625, 14.15736987415571),
(21.875, 13.672075756508647),
(21.875, 14.992307692307794),
(21.816783216783378, 14.992307692307794)],
[(22.625, 14.008484162896028),
(22.625, 13.56730769230769),
(21.943181818181664, 13.56730769230769),
(22.625, 14.008484162896028)],
[(21.875, 13.44230769230769),
(22.625, 13.44230769230769),
(22.625, 2.4423076923076925),
(21.875, 2.4423076923076925),
(21.875, 13.44230769230769)],
[(22.556818181818187, 2.3173076923076925),
(21.875, 1.8761312217194543),
(21.875, 2.3173076923076925),
(22.556818181818187, 2.3173076923076925)],
[(40.0, 26.7579185520362),
(40.0, 15.596153846153847),
(39.875, 15.596153846153847),
(39.875, 26.677036199095024),
(40.0, 26.7579185520362)],
[(39.75, 0.0),
(22.363636363636367, 0.0),
(22.556818181818183, 0.125),
(39.75, 0.125),
(39.75, 0.0)],
[(40.0, 0.125),
(40.0, 0.0),
(39.875, 0.0),
(39.875, 0.125),
(40.0, 0.125)],
[(40.0, 11.41176470588235),
(40.0, 0.25),
(39.875, 0.25),
(39.875, 11.330882352941174),
(40.0, 11.41176470588235)],
[(39.75, 13.442307692307693),
(22.750000000000004, 2.4423076923076925),
(22.75, 2.4423076923076925),
(39.75, 13.442307692307693)],
[(22.75, 28.788461538461533),
(39.75000000000001, 28.788461538461544),
(34.083333333333336, 25.121794871794876),
(39.75, 28.78846153846154),
(22.75, 28.788461538461533)],
[(4.75, 2.4423076923076925),
(4.749999999999997, 2.4423076923076925),
(11.125, 6.56730769230773),
(4.75, 2.4423076923076925)],
[(35.25, 39.75),
(35.25000000000001, 39.75),
(18.25, 28.749999999999993),
(18.25, 28.749999999999996),
(35.25, 39.75)],
[(17.25, 13.403846153846155),
(17.25, 13.403846153846153),
(0.2499999999999929, 13.403846153846153),
(17.25, 24.403846153846164),
(17.25, 24.40384615384616),
(0.25, 13.403846153846155),
(17.25, 13.403846153846155)],
[(0.25, 15.596153846153847),
(0.2499999999999999, 15.596153846153847),
(0.2499999999999992, 26.596153846153847),
(17.250000000000007, 26.596153846153857),
(8.750000000000002, 21.09615384615385),
(17.25, 26.596153846153854),
(0.25, 26.596153846153847),
(0.25, 15.596153846153847)]]
mp2 = [[(0.25, 0.25),
(39.75, 0.25),
(39.75, 39.75),
(0.25, 39.75),
(0.25, 0.25)]]
polygons1 = [shapely.geometry.Polygon(pg) for pg in mp1]
polygons2 = [shapely.geometry.Polygon(pg) for pg in mp2]
multipolygon1 = shapely.geometry.MultiPolygon(polygons1)
multipolygon1_buffered = multipolygon1.buffer(0.001)
multipolygon1_without_last_polygon = shapely.geometry.MultiPolygon(polygons1[:-1])
multipolygon1_with_first_3_and_last = shapely.geometry.MultiPolygon(polygons1[:3] + polygons1[-1:])
multipolygon2 = shapely.geometry.MultiPolygon(polygons2)
for mp in [multipolygon1,
multipolygon1_buffered,
multipolygon1_without_last_polygon,
multipolygon1_with_first_3_and_last]:
try:
mp.intersection(multipolygon2)
except Exception as e:
print('fail:', e)
else:
print('success')
stenci
>All comments
Closed as not a bug.
At first I thought that multipolygon1 was an invalid shapely Multipolygon, so I tried to buffer it, and it worked. And also the intersection with the buffered version worked. So multipolygon1 seems to be valid.
This logic isn't correct.
We can use the explain_validity function to explore the problem: http://toblerity.org/shapely/manual.html#diagnostics
from shapely.validation import explain_validity
print(multipolygon1.is_valid) # False
print(explain_validity(multipolygon1)) # 'Ring Self-intersection[0.25 26.5961538461538]'
print(multipolygon1[-1].is_valid) # False
print(multipolygon1[0:-1].is_valid) #聽True
So the issue with the polygon is a self intersection in the last ring. The ring is invalid, thus the entire geometry is invalid. Without this last ring the geometry is valid.
Buffering a polygon by zero distance is a common way to fix geometries with issues like this, although you should be aware that it doesn't fix everything and may have unexpected results.
print(multipolygon1.buffer(0.0).is_valid) # True
snorfalorpagus
on 21 Nov 2015
Related issues
LostFan123
路
5Comments
dopplershift
路
3Comments
doctor-ian
路
4Comments
benediktbrandt
路
3Comments
chivasblue
路
3Comments
Most helpful comment
Closed as not a bug.
This logic isn't correct.
We can use the
explain_validityfunction to explore the problem: http://toblerity.org/shapely/manual.html#diagnosticsSo the issue with the polygon is a self intersection in the last ring. The ring is invalid, thus the entire geometry is invalid. Without this last ring the geometry is valid.
Buffering a polygon by zero distance is a common way to fix geometries with issues like this, although you should be aware that it doesn't fix everything and may have unexpected results.