Discussion:
[Community] Splitting self-intersecting polygon only returned one polygon in shapely in Python
Yuxiang Wang
2016-01-31 05:28:06 UTC
Permalink
Dear all,

I am using Python 3.5 64 bit in Windows 7 64 bit, shapely version 1.5.13.

I have the following code that returned me a self-intersecting polygon:

import numpy as npfrom shapely.geometry import Polygon,
MultiPolygonimport matplotlib.pyplot as plt

x = np.array([ 0.38517325, 0.40859912, 0.43296919, 0.4583215 , 0.4583215 ,
0.43296919, 0.40859912, 0.38517325, 0.36265506, 0.34100929])
y = np.array([ 62.5 , 56.17977528, 39.39698492, 0. ,
0. , 17.34605377, 39.13341671, 60.4180932 ,
76.02574417, 85.47008547])
polygon = Polygon(np.c_[x, y])
plt.plot(*polygon.exterior.xy)

[image: Self-intersecting polygon] <Loading Image...>

This is correct. Then I tried to obtain the two individual polygons by
using buffer(0):

split_polygon = polygon.buffer(0)
plt.plot(*polygon.exterior.xy)print(type(split_polygon))
plt.fill(*split_polygon.exterior.xy)

Unfortunately, it only returned of the the two polygons:

[image: Only returned one polygon] <Loading Image...>

Could anyone please help? Thanks!


Shawn
Yuxiang Wang
2016-02-01 18:41:41 UTC
Permalink
Dear all,

FYI - The question was perfectly answered by Mike T on StackOverflow. Thank
you so much Mike!

Link:
http://stackoverflow.com/questions/35110632/splitting-self-intersecting-polygon-only-returned-one-polygon-in-shapely-in-pyth/35119152#35119152

Original answer:

The first step is to close the LineString to make a LinearRing, which is
what Polygons are made of.

from shapely.geometry import LineString, MultiPolygonfrom shapely.ops
import polygonize, unary_union
# original data
ls = LineString(np.c_[x, y])# closed, non-simple
lr = LineString(ls.coords[:] + ls.coords[0:1])
lr.is_simple # False

However, note that it is non-simple, since the lines cross to make a
bow-tie. (The widely used buffer(0) trick usually does not work for fixing
bow-ties in my experience). This is unsuitable for a LinearRing, so it
needs further work. Make it simple and MultiLineString with unary_union:

mls = unary_union(lr)
mls.geom_type # MultiLineString'

Then use polygonize to find the Polygons from the linework:

for polygon in polygonize(mls):
print(polygon)

Or if you want one MultiPolygon geometry:

mp = MultiPolygon(list(polygonize(mls)))
Post by Yuxiang Wang
Dear all,
I am using Python 3.5 64 bit in Windows 7 64 bit, shapely version 1.5.13.
import numpy as npfrom shapely.geometry import Polygon, MultiPolygonimport matplotlib.pyplot as plt
x = np.array([ 0.38517325, 0.40859912, 0.43296919, 0.4583215 , 0.4583215 ,
0.43296919, 0.40859912, 0.38517325, 0.36265506, 0.34100929])
y = np.array([ 62.5 , 56.17977528, 39.39698492, 0. ,
0. , 17.34605377, 39.13341671, 60.4180932 ,
76.02574417, 85.47008547])
polygon = Polygon(np.c_[x, y])
plt.plot(*polygon.exterior.xy)
[image: Self-intersecting polygon] <http://i.stack.imgur.com/ieTxz.png>
This is correct. Then I tried to obtain the two individual polygons by
split_polygon = polygon.buffer(0)
plt.plot(*polygon.exterior.xy)print(type(split_polygon))
plt.fill(*split_polygon.exterior.xy)
[image: Only returned one polygon] <http://i.stack.imgur.com/0cRWs.png>
Could anyone please help? Thanks!
Shawn
--
Yuxiang "Shawn" Wang
Gerling Haptics Lab
University of Virginia
***@virginia.edu
+1 (434) 284-0836
https://sites.google.com/a/virginia.edu/yw5aj/
Loading...