9

I have a shapefile and I would like the exploded view with GeoPandas in several polygons, I cannot use ArcGIS and keep the polygons that have a given.

import csv 
import glob
import time
import sys
import datetime
import numpy as np
import pandas as pd
import geopandas as gpd
from shapely.geometry import Polygon

gdf = gpd.read_file("dossier/essais.shp")

time.sleep(10)
xmin,ymin,xmax,ymax = gdf.total_bounds

lenght = 2065
wide = 3069

cols = list(range(int(np.floor(xmin)), int(np.ceil(xmax)), wide))
rows = list(range(int(np.floor(ymin)), int(np.ceil(ymax)), lenght))
rows.reverse()

polygons = []
for x in cols:
    for y in rows:
        polygons.append( Polygon([(x,y), (x+wide, y), (x+wide, y+lenght), (x, y+lenght)]) )

grid = gpd.GeoDataFrame({'geometry':polygons}, crs = gdf.crs.to_string())
grid.to_file("dossier/grid.shp")

#df4.head(1000).to_file("dossier/"+file)
gdf.to_file("dossier/grille.shp")

enter image description here

4
  • What do you mean by "explode"? Do you mean "extract" the empty polygons? Commented Jan 23, 2021 at 8:42
  • Hello thank you for your answer indeed the purpose and to have the information if the polygon and empty or filled. Commented Jan 23, 2021 at 8:46
  • is explode means separate the individual areas ? Commented Jan 23, 2021 at 19:37
  • Please, do not forget about "What should I do when someone answers my question?"
    – Taras
    Commented Mar 21, 2022 at 19:49

1 Answer 1

10
gdf = gpd.read_file("dossier/essais.shp")
gdf_merged_geom = gdf.unary_union

#
# other lines
#

grid = gpd.GeoDataFrame({'geometry':polygons}, crs=gdf.crs.to_string())

# add new CASE column
grid["CASE"] = "FILLED"

for index, polygon in grid.iterrows():
    if polygon.geometry.disjoint(gdf_merged_geom):
        grid.loc[index, "CASE"] = "EMPTY"

print(grid.head())

# OUTPUT
#           geometry      CASE
# 0 POLYGON ((282...    EMPTY
# 1 POLYGON ((282...    FILLED
# 2 POLYGON ((282...    FILLED
# 3 POLYGON ((282...    FILLED
# 4 POLYGON ((282...    FILLED

Sample data:

enter image description here

Result:

custom_draw(grid[grid["CASE"]=="EMPTY"].geometry)
# OR 
# grid[grid["CASE"]=="EMPTY"].plot()

enter image description here

2
  • 3
    Please use unary_union instead of cascaded_union as recommended in the Shapely manual (The unary union function can operate on different geometry types, not only polygons as is the case for the older cascaded union.)
    – gene
    Commented Jan 23, 2021 at 12:18
  • @gene Done. Thank you. I misunderstood the note in the shapely manual. Commented Jan 23, 2021 at 15:50

Your Answer

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

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