I'm trying to capture all the non-overlapping polygons indicated below using Shapely (given polygons A,B & C). Moreover, I'm hoping to do so without iteration, testing for intersect etc. The accepted answer to this question expresses the PostGIS method but it would seem that 'union' means different things to different people.
2 Answers
You need to iterate at some level. (Update: I've edited to remove all "for" loops, except for one list comprehension)
# imports used throughout this example
from shapely.geometry import Point
from shapely.ops import cascaded_union
from itertools import combinations
# Here are your input shapes (circles A, B, C)
A = Point(3, 6).buffer(4)
B = Point(6, 2).buffer(4)
C = Point(1, 2).buffer(4)
# list the shapes so they are iterable
shapes = [A, B, C]
First you need the union of all intersections (use a cascaded union), using the combination pair of each shape. Then you remove (via difference
) the intersections from the union of all shapes.
# All intersections
inter = cascaded_union([pair[0].intersection(pair[1]) for pair in combinations(shapes, 2)])
# Remove from union of all shapes
nonoverlap = cascaded_union(shapes).difference(inter)
Here is what nonoverlap
looks like (via JTS Test Builder):
After a few years there seems to be a better solution through shapely
:
# imports used throughout this example
from shapely.geometry import Point
from shapely.ops import polygonize, unary_union
# Here are your input shapes (circles A, B, C)
A = Point(3, 6).buffer(4)
B = Point(6, 2).buffer(4)
C = Point(1, 2).buffer(4)
...
# list the shapes so they are iterable
shapes = [A, B, C, ...]
# generate the overlay
list(polygonize(unary_union(list(x.exterior for x in shapes))))
It supports any lenght of geometry, the only problem it's about the calculation time, and not supporting polygon with holes.
-
Out of curiosity, why do you think your solution is better than the one by @MikeT? I can only read about a problem in terms of calculation times.– mgriCommented Jul 3, 2017 at 8:08