Get perimeter out of irregular group of shapes

53 Views Asked by At

I'm working on groups of shapes (think groups of city blocks) that are almost contiguous, like this:

enter image description here

What I'd like to do is to get the shape of the outer perimeter of this group. With more regular forms, this was almost working:

per_df = gpd.GeoSeries(df.geometry.unary_union.convex_hull).boundary

Bringing something like this:

enter image description here

But on an irregular shape it brings this:

enter image description here

Is there some way to "fuse" or join my shapes into one so it is easier to calculate its perimeter/boundary?

Here's a simpler reproducible example:

p1=Polygon([(0,0),(10,0),(10,9.8),(0,9.8)])
p2=Polygon([(10.2,10),(20,10),(20,20),(10.2,20)])
p3=Polygon([(10,10),(9.8,20),(0,10)])

df=gpd.GeoDataFrame(geometry=[p1,p2,p3])
per_df = gpd.GeoSeries(df.geometry.unary_union.convex_hull).boundary
ax = df.plot()

I'd like to get the perimeter of this group of shapes, even though there's a bit of separation between them: enter image description here

Thanks!

3

There are 3 best solutions below

1
Basile On BEST ANSWER

You can go with:

buf = df.geometry.buffer(10) #or a value that suits you better
unary = buf.unary_union
print(unary.length)

You can also buffer with -10 (or the chosen value) the unary union to be closer to the real shape before calculating the length.

0
Sam On

You could try performing a buffer on your geometries so that they overlap, followed by a dissolve, generating a single geometry. You could take the boundary of this geometry. Optionally, a negative buffer can be applied to get the boundary fitting closer to the original data, although some details may be lost

p1=Polygon([(0,0),(10,0),(10,9.8),(0,9.8)])
p2=Polygon([(10.2,10),(20,10),(20,20),(10.2,20)])
p3=Polygon([(10,10),(9.8,20),(0,10)])

df=gpd.GeoDataFrame(geometry=[p1,p2,p3])
df["geometry"] = df.buffer(0.5)
dissolved_df = df.dissolve()
dissolved_df["geometry"] = dissolved_df.buffer(-0.5)
per_df = dissolved_df.boundary
0
Pieter On

If you want to keep the output as close as possible to the original shapes and/or need a solution that scales well for large input files, you could try geofileops.dissolve_within_distance

Sample script:

from pathlib import Path
import geofileops as gfo
import geopandas as gpd
from matplotlib import pyplot as plt
from shapely import Polygon

# Prepare input
p1 = Polygon([(0, 0), (10, 0), (10, 9.8), (0, 9.8)])
p2 = Polygon([(10.2, 10), (20, 10), (20, 20), (10.2, 20)])
p3 = Polygon([(10, 10), (9.8, 20), (0, 10)])
df = gpd.GeoDataFrame(geometry=[p1, p2, p3])
input_path = Path("input.gpkg")
df.to_file("input.gpkg")

# Process
output_path = Path("output.gpkg")
gfo.dissolve_within_distance(input_path, output_path, distance=1, gridsize=0.0)

# Visualize
result_df = gpd.read_file("output.gpkg")
ax = result_df.plot(color="green")
df.plot(ax=ax, edgecolor="blue", facecolor="none")
plt.show()

Result in green, input in blue:

enter image description here

Disclaimer: I'm the developer of geofileops.