Merging Adjacent OSM Polygons for Low Zoom Jump to heading
Turn ten thousand individual residential parcels into forty built-up areas a reader can actually see, without merging a park into a car park or inventing an area that spans a river.
Prerequisites Jump to heading
Conceptual minimum Jump to heading
Aggregation replaces several source features with one output feature that never existed in the data. That is a bigger step than simplification, and it needs three decisions made explicitly.
What may merge. Adjacency alone is not a reason to merge: a park touching an industrial estate is two things, and merging them produces an area that is neither. Grouping must be by class, with the permitted classes named rather than inferred.
How close counts as adjacent. Polygons in OSM rarely share exact boundaries; a road, a hedge or a sliver of unclassified land usually sits between them. A gap tolerance — buffer outwards, union, buffer back — bridges those gaps. The tolerance must come from the target zoom: at zoom 8 a twenty-metre gap is invisible and should be bridged; at zoom 14 it is a street and must not be.
What the result is worth keeping. A merged area smaller than a few pixels at the target zoom contributes nothing but bytes. A minimum-area threshold derived from the same zoom removes them.
Runnable solution Jump to heading
from __future__ import annotations
import logging
import math
from dataclasses import dataclass
import geopandas as gpd
from shapely.geometry import MultiPolygon, Polygon
from shapely.ops import unary_union
logging.basicConfig(level=logging.INFO, format="%(levelname)s: %(message)s")
logger = logging.getLogger("osm.tiles.aggregate")
EARTH_CIRCUMFERENCE = 40_075_016.686
EXTENT = 4096
TILE_PIXELS = 512
# Which classes may merge with which. Anything absent merges only with itself.
MERGE_GROUPS: dict[str, str] = {
"residential": "built_up", "retail": "built_up", "commercial": "built_up",
"industrial": "industrial",
"forest": "green", "wood": "green", "meadow": "green", "grass": "green",
"park": "park", # deliberately its own group: parks are not "green"
}
@dataclass(frozen=True)
class ZoomBudget:
zoom: int
def metres_per_pixel(self, latitude_deg: float = 50.0) -> float:
tile_m = EARTH_CIRCUMFERENCE / (2 ** self.zoom)
return tile_m / TILE_PIXELS * math.cos(math.radians(latitude_deg))
def gap_tolerance_m(self) -> float:
"""Gaps narrower than about three pixels are invisible; bridge them."""
return 3.0 * self.metres_per_pixel()
def min_area_m2(self) -> float:
"""Areas smaller than about four pixels square are not worth a feature."""
return (4.0 * self.metres_per_pixel()) ** 2
def aggregate(frame: gpd.GeoDataFrame, zoom: int,
class_column: str = "class") -> gpd.GeoDataFrame:
"""Merge touching same-group polygons into low-zoom areas.
`frame` must be in a metric projection: buffering in degrees is meaningless
because a degree of longitude varies with latitude.
"""
if frame.crs is None or frame.crs.is_geographic:
raise ValueError("reproject to a metric CRS before aggregating")
budget = ZoomBudget(zoom)
gap = budget.gap_tolerance_m()
floor = budget.min_area_m2()
logger.info("zoom %d: bridging gaps under %.0f m, dropping areas under %.0f m2",
zoom, gap, floor)
frame = frame.copy()
frame["_group"] = frame[class_column].map(MERGE_GROUPS).fillna(
frame[class_column])
rows: list[dict] = []
for group, members in frame.groupby("_group"):
# Grow, dissolve, shrink: half the tolerance each way keeps the outline honest.
grown = [g.buffer(gap / 2.0, join_style=2) for g in members.geometry]
dissolved = unary_union(grown).buffer(-gap / 2.0, join_style=2)
parts = (list(dissolved.geoms)
if isinstance(dissolved, MultiPolygon) else [dissolved])
kept = 0
for part in parts:
if part.is_empty or part.area < floor:
continue
# Attributes are RECOMPUTED from the members, never inherited from one.
inside = members[members.geometry.intersects(part)]
rows.append({
"class": group,
"member_count": int(len(inside)),
"source_area_m2": float(inside.geometry.area.sum()),
"merged_area_m2": float(part.area),
"geometry": part,
})
kept += 1
logger.info("group %-12s %5d member(s) -> %3d merged area(s)",
group, len(members), kept)
out = gpd.GeoDataFrame(rows, geometry="geometry", crs=frame.crs)
if not out.empty:
inflation = out["merged_area_m2"].sum() / out["source_area_m2"].sum()
logger.info("merged area is %.2fx the source area", inflation)
if inflation > 1.35:
logger.warning("aggregation inflated area by more than 35%% — the gap "
"tolerance is probably too large for this zoom")
return out
if __name__ == "__main__":
logger.info("reproject to a metric CRS, then call aggregate(frame, zoom=8)")
Step-by-step walkthrough Jump to heading
- Refuse geographic coordinates. Buffering by a distance requires a metric projection; buffering degrees produces a shape that is wrong by the cosine of the latitude and silently so.
- Map classes to merge groups. The table is explicit, and a class not in it merges only with itself. Putting parks in their own group rather than lumping them with generic green space is the kind of decision that has to be visible.
- Derive both thresholds from the zoom. Gap tolerance is a few screen pixels’ worth of ground distance; minimum area is a small number of pixels squared. Both then scale correctly across the zoom range without a second table.
- Buffer by half in each direction. Growing by half the tolerance and shrinking by the same amount closes gaps up to the full tolerance while leaving the outer boundary approximately where it started.
- Use a mitre join. A round join adds vertices at every corner; a mitre join keeps the merged outline closer to the source shapes and produces far less geometry.
- Recompute attributes from members. The merged feature carries a member count and the summed source area, which are meaningful. Inheriting the name or identifier of an arbitrary member is not.
- Drop results below the floor. An area under a few pixels square at the target zoom costs bytes and contributes nothing.
- Watch the inflation ratio. Merged area meaningfully larger than the sum of source areas means the buffer is bridging gaps that are real features, which is the characteristic failure of a tolerance set too high.
Verification Jump to heading
- No cross-class merges. Group the output by class and confirm every merged area’s members share one group.
- Inflation is modest. Merged area should be within roughly a third of the summed source area; more means the tolerance is bridging real gaps.
- Rivers and motorways still separate areas. Pick a built-up area split by a river and confirm the merge did not jump it.
- Member counts are plausible. A merged area with one member is not an aggregation; a very large count at high zoom suggests the tolerance is too coarse.
- The output is materially smaller. Vertex and feature counts should fall by an order of magnitude at low zoom; if they do not, aggregation is not earning its cost.
Common errors and fixes Jump to heading
| Symptom | Root cause | One-line fix |
|---|---|---|
| A park merged into an industrial area | Grouped on adjacency, not class | Map classes to explicit merge groups |
| Merged areas span rivers | Gap tolerance larger than the river | Derive the tolerance from metres per pixel at the zoom |
| Outlines visibly inflated | Buffer applied outwards only | Buffer out and back by the same distance |
| Enormous vertex counts | Round join style on every corner | Use a mitre join for the buffer operations |
| Results wrong by a latitude factor | Buffered in geographic degrees | Reproject to a metric CRS before buffering |
| Merged feature carries a random name | Attributes inherited from one member | Recompute attributes from all members |
| Tiny slivers in the output | No minimum-area threshold | Drop results below a few pixels squared |
Specification reference Jump to heading
A buffer operation offsets a geometry by a fixed distance, and a positive buffer followed by a negative buffer of the same magnitude — a morphological closing — joins shapes separated by less than twice that distance while approximately restoring the original outline. Shapely implements both through
buffer, with the join style controlling corner treatment. See the Shapely documentation forbuffer,unary_unionand the join-style options used here.
Frequently Asked Questions Jump to heading
Why merge by class rather than by adjacency alone?
Because an area that is part park and part industrial estate describes nothing. Adjacency tells you two polygons touch; it says nothing about whether the thing they would become is a thing. Grouping by an explicit class map keeps merged areas meaningful and, just as importantly, makes the cartographic decision visible in one table a colleague can review rather than implicit in the geometry code.
How do I choose the gap tolerance?
Derive it from metres per pixel at the target zoom. A gap of a few pixels is invisible to a reader and should be bridged; a gap of many pixels is a real feature — a river, a motorway, a railway corridor — and must not be. Because metres per pixel changes by a factor of two per zoom level, one derived formula covers the whole range while a hand-picked distance is right at exactly one zoom.
Why buffer outwards and then inwards instead of just buffering out?
Because buffering outwards alone inflates every outline by the tolerance, so a merged built-up area spills across its real boundary by the same distance it used to bridge gaps. Applying the inverse buffer afterwards restores the outline to approximately where it began while leaving the members joined, which is the whole point of the closing operation.
What attributes should a merged feature carry?
Ones computed from the members: the class, how many members were merged, and the summed source area alongside the merged area. Those are meaningful and they support the verification. What it must not carry is a name or identifier taken from one arbitrary member, which implies the merged area is that feature when it is a new object representing many.
Related Jump to heading
- Cartographic Generalization of OSM Data — the parent topic and where aggregation sits among the operators.
- Simplifying OSM Geometry per Zoom Level — the step that usually precedes aggregation.
- Snapping Near-Duplicate OSM Nodes with a Tolerance — the same tolerance-choosing discipline at vertex level.
- Measuring Area Accurately on OSM Polygons — why the metric projection matters before buffering.
- Mapping OSM Tags to a Fixed Schema with YAML — producing the class attribute the merge groups key on.
Up one level: Cartographic Generalization of OSM Data.