2

I'm retrieving some geometries from OpenStreetMap with OSMnx library with the following code:

G = ox.graph_from_place('Casco Viejo, Bilbao, Spain', network_type='walk', 
                        retain_all=True, buffer_dist = 50, which_result = 2,
                        infrastructure = 'relation["highway" = "pedestrian"]')

which yields the following graph composed by shapely linestrings:

enter image description here

Then I convert the graph into a geopandas geodataframe:

ped = ox.graph_to_gdfs(G, nodes = False)

I've tried this to convert Linestrings to Points and then Points to Multipolygon

Is there a way to convert this linestrings into a shapely Multipolygon:

from shapely import geometry, ops

combine them into a multi-linestring

multi_line = geometry.MultiLineString(list(ped['geometry']))

merged_line = ops.linemerge(multi_line)

from shapely.geometry import Point, MultiPoint

points = [] for i in range(0, len(merged_line)): points.append((Point(list(merged_line[i].coords[1]))))

coords = [p.coords[:][0] for p in points]

poly = Polygon(coords)

This yields a weird wrong geometry:

shape(poly)

If I try:

MultiPolygon(points)

It gives this Error Message: TypeError: 'Point' object is not subscriptable

Is there a way to transform Linestrings into Multipolygon and this into GeoDataFrame?

Rodrigo Vargas
  • 626
  • 1
  • 6
  • 21

1 Answers1

2

For that you need to use unary_union and polygonize (see Split polygon by MultiLineString - shapely)

....
ped = ox.graph_to_gdfs(G, nodes = False)
from shapely.ops import unary_union, polygonize
multi_line = ped.geometry.values
border_lines = unary_union(multi_line)
result = MultiPolygon(polygonize(border_lines))

The resulting polygons

enter image description here

enter image description hereenter image description here

New

To fix the issue of holes, a solution is to retain only the polygons with holes

polys = []
result = list(ops.polygonize(border_lines))
for poly in result:
    if len(poly.interiors) >= 1: # only polygon with holes
        polys.append(poly)

In one line:

polys = [poly for poly in result if len(poly.interiors) >= 1]

Only two polygons with holes

enter image description here

gene
  • 54,868
  • 3
  • 110
  • 187