Fixing OpenStreetMap Lakes and Parks That Show as Lines Instead of Polygons

On this page

Problem statement

You extract lakes, parks or buildings from OpenStreetMap and they come out as outlines that cannot be filled, have no area and cannot be used in a spatial join. Worse, the largest features, a wide river, a big park, a lake with islands, are missing altogether. Both symptoms have the same cause: the reader was not asked to assemble areas. Without area assembly, closed ways come back as rings of line, and multipolygon relations, which is how large and complex areas are mapped, are not built at all.

This page reproduces both symptoms with the water bodies of Washington, D.C. from the OpenStreetMap extract of 1 January 2026, measures what is lost, and fixes it with pyosmium's area assembly.

Quick answer

Enable area assembly and keep only the areas:

fp = (osmium.FileProcessor(str(PBF)).with_areas()
      .with_filter(osmium.filter.EntityFilter(osmium.osm.AREA))
      .with_filter(osmium.filter.TagFilter(("natural", "water")))
      .with_filter(osmium.filter.GeoInterfaceFilter(tags=["natural", "water", "name"])))
good = gpd.GeoDataFrame.from_features(fp, crs=4326)
print(f"water areas: {len(good):,}; {good.geom_type.value_counts().to_dict()}; {good.to_crs(32618).area.sum() / 1e6:.2f} km²")
water areas: 325; {'MultiPolygon': 325}; 23.88 km²

With area assembly, the extract yields 325 water polygons covering 23.88 km², including the parts of the extract outside the District.

Two maps of water around the Potomac and Anacostia: read without areas, 303 outlines and gaps where the rivers should be; read with areas, 325 filled polygons including the full river surfaces. Two maps of water around the Potomac and Anacostia: read without areas, 303 outlines and gaps where the rivers should be; read with areas, 325 filled polygons including the full river surfaces.
The same water, read two ways.

Read without areas, see what went missing, why rings are not enough

1. Reproduce: read without areas

fp = (osmium.FileProcessor(str(PBF)).with_locations()
      .with_filter(osmium.filter.TagFilter(("natural", "water")))
      .with_filter(osmium.filter.GeoInterfaceFilter(tags=["natural", "water", "name"])))
bad = gpd.GeoDataFrame.from_features(fp, crs=4326)
print(f"natural=water without areas: {len(bad):,} features; {bad.geom_type.value_counts().to_dict()}")
print(f"area of those geometries: {bad.to_crs(32618).area.sum():.0f} m²")
natural=water without areas: 303 features; {'LineString': 303}
area of those geometries: 0 m²

Without it, the same filter returns 303 lines with a total area of zero square metres.

2. What was missing entirely

Multipolygon relations are not built without area assembly. Count the water relations in the file:

con = duckdb.connect(); con.sql("INSTALL spatial; LOAD spatial")
rel = con.sql(f"""SELECT id, tags['name'] AS name, tags['water'] AS water, len(refs) AS members FROM ST_ReadOSM('{PBF}')
                  WHERE kind = 'relation' AND tags['natural'] = 'water' ORDER BY members DESC""").df()
print(f"natural=water multipolygon relations: {len(rel)}; largest by members:")
print(rel.head(5).to_string(index=False))
natural=water multipolygon relations: 22; largest by members:
     id               name     water  members
1857635                NaN     river       15
1858408                NaN     river       13
9370920                NaN     river       12
1870168       Kingman Lake reservoir        7
5991610 Washington Channel     river        6

The file has 22 water multipolygon relations, the largest of them unnamed sections of the Potomac and Anacostia surfaces. Kingman Lake and the Washington Channel are among the named ones. None of them is built without area assembly.

Water area recovered from the D.C. extract by three methods. water area by methodlines, no areas0 km²rings made into polygons1.74 km²with_areas()23.88 km²08 km²16 km²24 km²22 water relations are only built by area assembly spatialworkflow.io
Water area by method.

3. Turning rings into polygons yourself is not enough

Converting each closed line to a polygon recovers simple ponds, but not the relations, and it fills islands:

from shapely.geometry import Polygon
rings = bad[bad.geometry.is_closed].geometry.apply(lambda g: Polygon(g.coords))
print(f"closed rings converted: {len(rings):,}, {gpd.GeoSeries(rings, crs=4326).to_crs(32618).area.sum() / 1e6:.2f} km², against {good.to_crs(32618).area.sum() / 1e6:.2f} km² with area assembly")
closed rings converted: 303, 1.74 km², against 23.88 km² with area assembly

Converting the 303 rings to polygons recovers only 1.74 km², less than a tenth of the 23.88 km² that area assembly finds, because the rivers are relations.

Parks, the same problem

The data used on this page

The Geofabrik OpenStreetMap extract for the District of Columbia dated 1 January 2026 (© OpenStreetMap contributors, ODbL), read with pyosmium 4 and DuckDB's spatial extension.

from pathlib import Path

import duckdb
import geopandas as gpd
import osmium
import requests

Path("data").mkdir(exist_ok=True)
PBF = Path("data/district-of-columbia-260101.osm.pbf")
if not PBF.exists():
    r = requests.get("https://download.geofabrik.de/north-america/us/district-of-columbia-260101.osm.pbf", timeout=600)
    r.raise_for_status()
    PBF.write_bytes(r.content)
print(f"{PBF.name} ready")
district-of-columbia-260101.osm.pbf ready

The dated D.C. extract is ready.

Example 1: parks, the same problem

Parks are mapped the same way. Compare leisure=park with and without area assembly:

def read(areas):
    fp = osmium.FileProcessor(str(PBF))
    fp = fp.with_areas().with_filter(osmium.filter.EntityFilter(osmium.osm.NODE | osmium.osm.AREA)) if areas else fp.with_locations()
    fp = fp.with_filter(osmium.filter.TagFilter(("leisure", "park"))).with_filter(osmium.filter.GeoInterfaceFilter(tags=["name"]))
    return gpd.GeoDataFrame.from_features(fp, crs=4326)
for label, g in (("without areas", read(False)), ("with areas", read(True))):
    print(f"{label:>13}: {len(g):,} parks; {g.geom_type.value_counts().to_dict()}; polygon area {g.to_crs(32618).area.sum() / 1e6:.2f} km²")
without areas: 730 parks; {'LineString': 691, 'Point': 39}; polygon area 0.00 km²
   with areas: 775 parks; {'MultiPolygon': 736, 'Point': 39}; polygon area 24.73 km²

Parks behave the same: without areas there are 730 features, 691 of them lines, with no area at all. With areas there are 775, of which 736 are polygons covering 24.73 km².

The two ways OpenStreetMap stores an area. how OSM stores areasmultipolygon relationouter and inner member waysrivers, big parks, holesclosed waya ring of nodesponds, small parksnodescoordinatesshared by bothA way has at most 2,000 nodes and no holes spatialworkflow.io
How OSM stores areas.

Closed ways and multipolygon relations

How pyosmium assembles areas from ways and relations. area assemblyfirst passcollect relationmemberssecond passnode locationsringsouter and innerpolygonsholes cut outareasMultiPolygonfeatureswith_areas() in pyosmium 4 spatialworkflow.io
Area assembly.

OpenStreetMap stores areas in two ways. A simple area is a closed way: a ring of nodes with tags such as natural=water or leisure=park. A complex area is a multipolygon relation: a relation whose members are ways with the roles outer and inner, assembled into polygons with holes. Large features, such as a river surface hundreds of kilometres long, a lake with islands, or a park made of many parcels, are almost always relations, because a single way has a limit of 2,000 nodes and cannot have holes.

A reader that only processes ways sees the closed way as a sequence of points, so it returns a line. It never assembles relations, so multipolygons disappear, and with them the biggest features. Converting rings to polygons afterwards fixes the simple cases but cannot recover relations and ignores inner rings. The only correct fix is to let the reader assemble areas, with_areas() in pyosmium, and to keep the areas rather than the ways. Then every lake and park arrives as a polygon, holes and all.

Incomplete relations, coastlines and invalid areas

Checks when extracting areas from OpenStreetMap. before trusting OSM polygonswith_areas() onor relations vanishkeep AREA, drop WAYno line copiescount the relations22 water relations herecompare areas0 vs 23.88 km²extract edgesincomplete relations dropD.C., OSM 1 January 2026 spatialworkflow.io
Before trusting OSM polygons.
  • Incomplete relations: a relation whose members run outside the extract cannot be assembled and is silently dropped; see polygons missing at extract edges.
  • Coastlines: seas are not natural=water; coastline polygons are built separately from natural=coastline lines.
  • Invalid areas: broken rings in the data fail assembly; pyosmium skips them, so compare counts with the raw relations.
  • Duplicates: with areas enabled, also drop the way copies; see features counted twice.

FAQ

Why are OSM lakes lines instead of polygons?

The reader returned closed ways as lines because area assembly was off; enable it (with_areas() in pyosmium) and keep the areas.

Why are large lakes and parks missing from my OSM data?

They are multipolygon relations, which are only built when areas are assembled.

Can I just convert closed lines to polygons?

Only for simple areas; relations stay missing and inner holes are filled.

What is a multipolygon relation?

A relation of outer and inner member ways that together describe one area, possibly with holes.