I have looked at your data and the book example, the problem is that there are three invalid polygons in data that are processed in the book:
GSHHS_l_L1.shp
ID = 92-W
POLYGON ((-180.0 71.514793999999995,-179.69008299999999 71.577888999999999,-178.648889 71.577416999999997,-178.40644399999999 71.549916999999994,-177.406306 71.244167000000004,-177.877444 71.022889000000006,-179.500111 70.863749999999996,-179.93011100000001 70.979583000000005,-180.0 70.962072000000006))
ID = 486-W
POLYGON ((-180.0 -16.799126,-179.84419399999999 -16.691278,-179.80041700000001 -16.789193999999998,-179.850472 -16.878361000000002,-180.0 -16.959561))
GSHHS_l_L2.shp
ID = 7333-W
POLYGON ((-180.0 65.393473,-179.76583299999999 65.428332999999995,-179.95416700000001 65.385555999999994,-179.90972199999999 65.316389,-180.0 65.321635))
Because this is an example it would be easiest to delete those polygons from the dataset or just add one if statement in you code
if geometry.IsValid():
cursor.execute("INSERT INTO gshhs (level, geom) VALUES (%s, ST_GeomFromText(%s, 4326))", (level, wkt))
You can wrap the polygon in a list
and pass that as an argument to the MultiPolygon
constructor.
Demo:
from shapely.geometry.multipolygon import MultiPolygon
from shapely import wkt
p = wkt.loads(u'POLYGON((0 0,0 1,1 1,0 0))')
m = MultiPolygon([p])
print(m.wkt)
# prints 'MULTIPOLYGON (((0 0, 0 1, 1 1, 0 0)))'
Best Answer
You need to iterate at some level. (Update: I've edited to remove all "for" loops, except for one list comprehension)
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.Here is what
nonoverlap
looks like (via JTS Test Builder):