0

I have a Geopandas data frame with five multi-polygons geometries. Now, I want to intersect all the geometries at once and get the intersected polygon out of it. I tried using unary union and polygonize but this is giving me a list of polygons but I want only one polygon that has the intersection of the set of multi-polygon polygons. How can we intersect a set of multi-polygon or polygons together and get the intersected polygon?

df=
location    geometry    
1          MULTIPOLYGON (((-0.304766 51.425882, -0.304904...    
2          MULTIPOLYGON (((-0.305968 51.427425, -0.30608 ...    
3          MULTIPOLYGON (((-0.358358 51.423471, -0.3581 5...    
4          MULTIPOLYGON (((-0.357654 51.413925, -0.357604...    
    rows=[]
    listpoly = [a.intersection(b) for a, b in combinations(df['geometry'], 2)]
    rings = []
    for poly in listpoly:
      if type(poly) == MultiPolygon:
         exterior_lines = LineString()
         for p in poly:
            exterior_lines = exterior_lines.union(p.exterior)
         rings.append(exterior_lines)
      elif type(poly) == Polygon:
         rings.append(LineString(list(poly.exterior.coords)))
    union = unary_union(rings)
    result = [geom for geom in polygonize(union)]
    print(result)
MULTILINESTRING((-0.0345 54.900...))
MULTILINESTRING((-0.045 54.200...))
MULTILINESTRING((-0.05 54.650...))
MULTILINESTRING((-0.04 54.750...))
1
  • Can you add a screenshot showing your data and your expected result
    – BERA
    May 11 at 5:48

1 Answer 1

1

This should work for a set of polygons:

import os
os.environ['USE_PYGEOS'] = '0'  # this is required by shapely when using python 3.10
import geopandas


file = geopandas.read_file("C:/Users/user/Desktop/all_multi.shp")

# Set initials
intersect_list = []
uniq = []
any_intersect = []

# Loop unitil only one polygon left
while True:
    intersect_list = []
    for i, geom in enumerate(file.geometry):
        for j, next_geom in enumerate(file.geometry):
            if i != j:
                if geom.intersects(next_geom):
                    intersect_list.append(geom.intersection(next_geom))
                                      
    # Set last unique list
    last = uniq

    # Remove duplicate geometries
    uniq = []
    for poly in intersect_list:
        if not any(p.equals(poly) for p in uniq):
            uniq.append(poly)    
    
    # Check if there are anymore intersections
    if len(uniq) == 0:
        break

    # Update GeoDataFrame
    file = geopandas.GeoDataFrame(geometry=uniq)

    
geopandas.GeoDataFrame(geometry=last).to_file("C:/Users/user/Desktop/intersection.shp")

enter image description here

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.