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))) On Sun, Jan 31, 2016 at 12:28 AM, Yuxiang Wang <yw...@virginia.edu> wrote: > 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] <http://i.stack.imgur.com/ieTxz.png> > > 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] <http://i.stack.imgur.com/0cRWs.png> > > Could anyone please help? Thanks! > > > Shawn > -- Yuxiang "Shawn" Wang Gerling Haptics Lab University of Virginia yw...@virginia.edu +1 (434) 284-0836 https://sites.google.com/a/virginia.edu/yw5aj/
_______________________________________________ Community mailing list Community@lists.gispython.org http://lists.gispython.org/mailman/listinfo/community
