⚙️ How-To Urban form & morphology Networks & routing Urban Morphology · step 11 Notebook By · Updated

How to Extract Urban Blocks from a Street Network in Python

On this page

Problem statement

The block — the land enclosed by streets — is the unit at which a city is walkable or not, at which plots are grouped, and at which many density and form measures are reported. No open dataset publishes blocks directly. They have to be derived by polygonising a street network: every closed loop of streets becomes a block.

That is a one-line operation and a data-quality problem. Dual carriageways produce long slivers between their two directions; dead ends and unconnected footpaths fail to close; the choice of which streets count as boundaries decides whether a courtyard with a service road is one block or three.

Quick answer

Polygonise the street lines with momepy, clip to an area, and filter slivers:

import numpy as np, geopandas as gpd, osmnx as ox, momepy
from shapely.geometry import box

ox.settings.http_user_agent = "spatialworkflow-examples/1.0"
ox.settings.overpass_url = "https://maps.mail.ru/osm/tools/overpass/api"   # mirror; overpass-api.de refused
ox.settings.overpass_rate_limit = False
G = ox.graph_from_point((52.0799, 4.3113), dist=1500, network_type="drive")
streets = ox.graph_to_gdfs(ox.convert.to_undirected(G), nodes=False).to_crs(28992)
cx, cy = 81000, 455000
aoi = gpd.GeoSeries([box(cx - 1000, cy - 1000, cx + 1000, cy + 1000)], crs=28992)

blocks = momepy.enclosures(streets, limit=aoi)
blocks["area"] = blocks.area
blocks["elongation"] = momepy.elongation(blocks)
slivers = (blocks.area < 500) | (blocks.elongation < 0.15)
print(f"{len(blocks):,} enclosures from {len(streets):,} street segments in 4 km²")
print(f"slivers (< 500 m² or very elongated): {slivers.sum()} ({slivers.mean():.0%})")
real = blocks[~slivers]
print(f"blocks: median {real.area.median():,.0f} m², 90th percentile "
      f"{real.area.quantile(.9):,.0f} m², {len(real) / 4:.0f} per km²")
819 enclosures from 1,892 street segments in 4 km²
slivers (< 500 m² or very elongated): 345 (42%)
blocks: median 8,315 m², 90th percentile 24,641 m², 118 per km²

Four enclosures in ten fail the filter. Many are the space between the two carriageways of a divided road or the traffic islands at a junction — land that is part of the street, not a block.

Boundary lines, undirected graph, polygonise, remove slivers

1. Decide which lines bound a block

Drivable streets give traffic blocks. Adding footpaths gives pedestrian blocks, smaller and more numerous. Railways and rivers are barriers too and belong in the boundary set.

2. Make the network undirected

Directed graphs contain each two-way street twice; polygonising two identical lines creates zero-area polygons.

3. Polygonise within a limit

momepy.enclosures splits a limit polygon by the street lines. Without a limit, only fully closed loops become polygons and the edge of the area is lost.

4. Remove slivers

Drop polygons that are too small or too thin to be blocks — or, better, collapse dual carriageways to a single centreline before polygonising.

5. Describe the blocks

Area, perimeter, elongation, number of buildings per block, and blocks per km². Block size is one of the strongest single predictors of walkability.

A street network polygonised into blocks, with thin slivers between the two carriageways of a divided road. every closed loop becomes a polygonstreetsdual carriagewayfacessliver: street, not blockFilter slivers by area and elongation, or collapse the carriageways first spatialworkflow.io
Every closed loop becomes a polygon — including the ones between two carriageways.

Drivable against all streets, and buildings per block

The data used on this page

Building footprints with LiDAR heights from the 3D BAG (the Dutch building register with heights from AHN), fetched from its WFS for central Den Haag — 4.285–4.330° E, 52.065–52.095° N — exactly as in How to download building footprints from OSM and Overture, and saved as denhaag_bag.parquet in the Dutch national grid (EPSG:28992). The OpenStreetMap answers are read as the Overpass mirror gave them on 10 October 2026.

import io
import os
import zipfile

import geopandas as gpd
import requests
from pyproj import Transformer

# OpenStreetMap as Overpass returned it on 10 October 2026, saved as osmnx's own cache so the counts
# below do not drift as OSM is edited; delete the cache folder to query live instead
if not os.path.exists("cache"):
    reply = requests.get("https://data.spatialworkflow.io/osm-cache-urban-blocks-python-2026-10-10.zip", timeout=300)
    reply.raise_for_status()
    with zipfile.ZipFile(io.BytesIO(reply.content)) as archive:
        archive.extractall("cache")

if not os.path.exists("denhaag_bag.parquet"):
    W, S, E, N = 4.285, 52.065, 4.330, 52.095
    to_rd = Transformer.from_crs(4326, 28992, always_xy=True)
    (x0, y0), (x1, y1) = to_rd.transform(W, S), to_rd.transform(E, N)
    feats, start = [], 0
    while True:
        fc = requests.get("https://data.3dbag.nl/api/BAG3D/wfs", params={
            "service": "WFS", "version": "2.0.0", "request": "GetFeature", "typeNames": "BAG3D:lod12",
            "count": 5000, "startIndex": start, "outputFormat": "application/json",
            "bbox": f"{x0},{y0},{x1},{y1},EPSG:28992"}, timeout=300).json()
        feats += fc["features"]
        if len(fc["features"]) < 5000:
            break
        start += 5000
    bag = gpd.GeoDataFrame.from_features(feats, crs=28992)
    bag["height"] = bag.b3_h_70p - bag.b3_h_maaiveld
    bag[["identificatie", "b3_h_70p", "b3_h_maaiveld", "height", "b3_opp_grond", "geometry"]] \
        .to_parquet("denhaag_bag.parquet")
print(len(gpd.read_parquet("denhaag_bag.parquet")), "buildings")
23929 buildings

Example 1 — drivable streets against all streets

G_all = ox.graph_from_point((52.0799, 4.3113), dist=1500, network_type="all")
streets_all = ox.graph_to_gdfs(ox.convert.to_undirected(G_all), nodes=False).to_crs(28992)
blocks_all = momepy.enclosures(streets_all, limit=aoi)
blocks_all = blocks_all[(blocks_all.area >= 500) & (momepy.elongation(blocks_all) >= 0.15)]
for name, g in (("drive network", real), ("all streets and paths", blocks_all)):
    print(f"{name:<22} {len(g):>5} blocks, median {g.area.median():>7,.0f} m², "
          f"{len(g) / 4:>4.0f} per km²")
drive network            474 blocks, median   8,315 m²,  118 per km²
all streets and paths   1304 blocks, median   2,746 m²,  326 per km²

Footpaths through courtyards and across parks cut drive-network blocks into much smaller pedestrian blocks. Which one counts as "the block" depends on who is moving through it.

Example 2 — buildings per block

b = gpd.read_parquet("denhaag_bag.parquet")
pts = gpd.GeoDataFrame(geometry=b.representative_point(), crs=28992)
counts = gpd.sjoin(pts, real[["geometry"]], predicate="within").groupby("index_right").size()
real = real.assign(buildings=counts.reindex(real.index).fillna(0).astype(int))
print(f"blocks with no building: {(real.buildings == 0).sum()} "
      f"(parks, squares, water, rail yards)")
print(f"buildings per built block: median {real[real.buildings > 0].buildings.median():.0f}, "
      f"max {real.buildings.max()}")
blocks with no building: 106 (parks, squares, water, rail yards)
buildings per built block: median 32, max 757

A block holding 757 buildings is the other extreme: on the drive network, a quarter whose inner lanes cannot be driven is one block, however many streets a pedestrian would count. Empty blocks are real too: they are the parks, squares and yards of the city. Leaving them in changes averages of block density; flagging them keeps the information.

Example 3 — intersection density alongside blocks

nodes = ox.graph_to_gdfs(G, edges=False).to_crs(28992)
Gp = ox.project_graph(ox.convert.to_undirected(G))
spn = ox.stats.count_streets_per_node(Gp)
junctions = nodes[nodes.index.map(spn) >= 3]
inside = junctions[junctions.within(aoi.iloc[0])]
print(f"intersections (3+ streets) per km²: {len(inside) / 4:.0f}")
print(f"mean block perimeter: {real.length.mean():,.0f} m")
intersections (3+ streets) per km²: 126
mean block perimeter: 469 m

Block size and intersection density are two views of the same grain: small blocks mean many junctions. Walkability indices typically include one or the other.

Median block area for blocks from the drive network and from all streets and paths in 4 km² of Den Haag. median block area, m²drive network (474)8,315 m² · 118 per km²all streets (1,304)2,746 m² · 326 per km²045009000Footpaths through courtyards and parks cut blocks to a third of the size spatialworkflow.io
The pedestrian grain is finer than the traffic grain.

Blocks come from polygonising

Why blocks come from polygonising

A street network is a set of lines; a block is the area they enclose. Polygonising builds every face of the planar graph the lines form. That is exact when the lines are correct and meet at junctions, and fragile when they do not — which is why cleaning the network is most of the work.

Why dual carriageways produce slivers

OSM maps divided roads as two one-way lines a few metres apart. Polygonising puts a long thin face between them. It is street, not block, and it distorts every statistic of block size and count. momepy offers tools to consolidate such lines before polygonising; a sliver filter afterwards is the blunt alternative.

Why the boundary set is a choice

A service road through a housing estate either splits it into blocks or does not, depending on whether it is included. Different studies make different choices; block counts are only comparable with the same boundary set.

From enclosures polygonised from the drive network to blocks and to blocks with buildings in 4 km² of Den Haag. enclosures to blocks, 4 km²enclosures819−345 sliversafter the sliver filter474−106 parks, squares, waterblocks with buildings368Median 32 buildings per built block; the largest holds 757 spatialworkflow.io
Four enclosures in ten are slivers; a fifth of the rest hold no building.

Dangling lines, bridges and roundabouts

  • Unconnected lines (dangling ends) do not split a block; snap them to the network or extend them with momepy.extend_lines.
  • Bridges and tunnels cross without junctions; polygonising 2D lines treats them as crossings.
  • Roundabouts create small central polygons; filter them like slivers.
  • Edge blocks clipped by the limit are partial; flag them in any statistics.
  • Railways and water as additional barriers: pass them to momepy.enclosures(additional_barriers=...).

FAQ

How do I get city blocks from OpenStreetMap?

Download the street network, make it undirected, and polygonise it within a limit with momepy.enclosures.

Why do I get thin sliver blocks?

Dual carriageways and junction islands. Filter by area and elongation or consolidate the lines first.

Should footpaths count as block boundaries?

For pedestrian blocks, yes; for traffic blocks, no. Choose and report.

What is a typical block size?

It varies enormously: dense historic centres have blocks of a few thousand m²; car-oriented suburbs of tens of thousands.

Why are some blocks empty?

They are parks, squares, water or rail land enclosed by streets. Flag rather than delete them.

Smaller blocks mean more route choices and shorter detours, which is why block size and intersection density appear in walkability indices.