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.

The buffer, union, unbuffer sequence and what each step does Four steps. The group step partitions features by class so only compatible areas can merge. The buffer step grows each polygon outwards by half the gap tolerance, closing the small gaps between neighbours that a road or hedge creates. The union step dissolves the overlapping buffered shapes into connected components. The unbuffer step shrinks the result back by the same distance, restoring approximately the original outline while keeping the members joined. Grow, dissolve, shrink — grouped by class group partition by class never merge across buffer out half the gap tolerance closes small gaps union dissolve overlaps connected components buffer in same distance back outline restored Buffering half the tolerance each way leaves the outline almost unchanged while gaps up to the full tolerance still close.
The symmetry of the two buffers is what keeps the merged outline honest rather than inflated.

Runnable solution Jump to heading

python
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

  1. 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.
  2. 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.
  3. 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.
  4. 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.
  5. 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.
  6. 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.
  7. Drop results below the floor. An area under a few pixels square at the target zoom costs bytes and contributes nothing.
  8. 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.
How the gap tolerance and minimum area change across zoom levels A grid of four zoom levels against three derived quantities. At zoom six the ground size of one screen pixel is about a kilometre and a half, the gap tolerance is several kilometres and the minimum area covers tens of square kilometres. At zoom eight those fall to about four hundred metres, a kilometre and several square kilometres. At zoom ten they fall to about a hundred metres, three hundred metres and under a square kilometre. At zoom twelve they fall to about twenty five metres, eighty metres and a few hectares. Both thresholds derive from one number: metres per pixel m per pixel Gap bridged Minimum area zoom 6 about 1,500 m about 4,500 m about 36 km2 zoom 8 about 390 m about 1,170 m about 2.4 km2 zoom 10 about 98 m about 294 m about 0.15 km2 zoom 12 about 24 m about 73 m about 1 hectare At zoom 12 the bridged gap is about the width of a wide street, which is the point at which aggregation should stop being applied at all.
Reading down the middle column shows exactly which real features each zoom is willing to ignore.
Feature and vertex reduction from aggregating a country land-use layer at zoom eight Four measurements before and after aggregation for a country-sized land-use layer at zoom eight. The source layer holds several hundred thousand polygons. After grouping and merging, a few thousand areas remain. Source vertex count is in the tens of millions. Merged vertex count is in the hundreds of thousands, a reduction of roughly two orders of magnitude, because shared interior boundaries between merged members disappear entirely. Two orders of magnitude, mostly from vanished interior edges Source polygons about 420,000 Merged areas about 5,000 Source vertices about 18 million Merged vertices about 220,000 Most of the saving is interior boundaries between merged members, which existed only to separate parcels a reader cannot distinguish.
Aggregation beats selection here because the removed geometry was never conveying anything at this zoom.

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 for buffer, unary_union and 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.

Up one level: Cartographic Generalization of OSM Data.