Validating Shared Boundaries Between Adjacent Polygons

9 min read

Two parcels that render as neighbours with no visible gap are, in most datasets, not sharing a boundary at all. They share a region of space to within a few micrometres, with different vertex counts along the common edge and coordinates that differ in the last few decimal places. This guide sits within topology rule enforcement and covers testing that an edge is genuinely shared rather than merely close.

The distinction matters because near-coincident boundaries are stable under rendering and unstable under everything else. Dissolve them and slivers appear. Simplify them and the two sides diverge. Reproject them and the divergence grows. A gap-and-overlap check will pass throughout, because the discrepancy is far below any area threshold worth setting.

Root cause: adjacency is asserted by area, not by geometry

The standard no-gaps check unions the layer and measures the holes. A hole of 10−1110^{-11} square metres is below every sane tolerance, so a boundary that disagrees at the eleventh decimal place passes. But the failure is not about area — it is that the two polygons have different vertex sequences along an edge they both claim.

Three boundary relationships that all pass an area-based gap check Three ways two neighbouring polygons can meet are compared, all of which satisfy an area-based no-gaps check. In the first, the boundary is genuinely shared: both polygons use exactly the same vertices along the common edge, so dissolving, simplifying and reprojecting all behave predictably. In the second, the boundary is near-coincident: the two sides differ in the final decimal places, so the enclosed sliver has negligible area and the gap check passes, but simplification pulls the sides apart and dissolving leaves visible slivers. In the third, both sides trace the same line but one carries additional vertices, so the geometries agree in space yet cannot be compared vertex by vertex and simplification moves them differently. truly shared identical vertices on both sides dissolve, simplify, reproject: safe gap check passes — correctly near-coincident sides differ in the last decimals simplify pulls them apart; slivers gap check passes — wrongly differing vertex density same line, one side has extra points no vertex-by-vertex comparison possible gap check passes — ambiguously

The three cases need different tests. The first is what you want. The second is a defect. The third is usually acceptable but must be a deliberate decision, because it rules out the cheapest form of the check.

Step-by-step implementation

Step 1 — Find the adjacent pairs efficiently

Comparing every pair is O(n2)O(n^2) and unnecessary. A spatial index reduces it to the pairs whose envelopes interact:

import geopandas as gpd
from shapely import STRtree


def adjacent_pairs(gdf: gpd.GeoDataFrame, buffer: float = 1e-6):
    """Yield (i, j) index pairs whose geometries touch or nearly touch, i < j."""
    geoms = list(gdf.geometry)
    tree = STRtree(geoms)

    for i, geom in enumerate(geoms):
        # Query on a slightly grown envelope so near-misses are caught too.
        for j in tree.query(geom.buffer(buffer)):
            j = int(j)
            if j <= i:
                continue
            if geom.buffer(buffer).intersects(geoms[j]):
                yield i, j

Growing the query geometry by a small buffer is what makes the near-coincident case visible: two polygons separated by 10−910^{-9} metres do not touch(), so an unbuffered query misses exactly the pairs this check exists to find. The trade-off between index strategies is covered in R-tree vs GiST index performance.

Step 2 — Extract the shared edge and measure the disagreement

For each pair, the shared portion is the intersection of the two boundaries. What matters is how far the two sides wander from it:

from shapely.geometry.base import BaseGeometry


def boundary_disagreement(a: BaseGeometry, b: BaseGeometry) -> float:
    """Maximum separation between the two polygons' versions of their common edge.

    Returns 0.0 for a truly shared boundary and the sliver width otherwise.
    """
    shared = a.boundary.intersection(b.boundary)
    if shared.is_empty:
        return float("inf")            # not adjacent at all

    # The symmetric difference of the two polygons along the edge is the sliver.
    sliver = a.symmetric_difference(b).difference(a.union(b).buffer(-1e-9))
    if sliver.is_empty:
        return 0.0

    # Hausdorff distance from the sliver back to the shared edge is its width.
    return shared.hausdorff_distance(sliver)

Hausdorff distance is the right measure here rather than area, because it is scale-free with respect to edge length. A one-micrometre disagreement along a 10 km boundary and along a 10 m boundary are the same defect and produce the same number, whereas the sliver areas differ by three orders of magnitude and would need different thresholds.

Formally, for the two boundary representations EaE_a and EbE_b the quantity being bounded is

dH(Ea,Eb)=max⁡{sup⁡p∈Eainf⁡q∈Eb∥p−q∥, sup⁡q∈Ebinf⁡p∈Ea∥p−q∥}d_H(E_a, E_b) = \max\left\{\sup_{p \in E_a} \inf_{q \in E_b} \lVert p - q \rVert,\ \sup_{q \in E_b} \inf_{p \in E_a} \lVert p - q \rVert\right\}

which is zero exactly when the two sides describe the same point set.

Step 3 — Choose a threshold from the data’s precision, not from taste

The threshold should come from the source data’s coordinate precision. A dataset captured at millimetre precision in a projected CRS has a natural floor of 10−310^{-3} m; anything below that is noise and anything above it is a real disagreement.

Source precision CRS Threshold Rationale
Survey-grade, projected EPSG:27700 1e-3 m Matches capture precision
Digitised from imagery EPSG:3857 1e-2 m Below one pixel at max zoom
Geographic degrees EPSG:4326 1e-9 deg ~0.1 mm at the equator
Generalised for display any not applicable Boundaries are not expected to match

Deriving the threshold rather than guessing it is what keeps the gate from being loosened every time it fires. Note the degree row: a tolerance expressed in degrees is not a distance, and the same numeric value means a different physical distance at different latitudes — a reason to run this check in a projected CRS wherever possible.

From layer to reported boundary disagreement The check is drawn as a pipeline. A spatial index is built over the layer. Each geometry is queried against the index using a slightly buffered envelope so that near-misses, which do not strictly touch, are still returned as candidates. Each candidate pair has the Hausdorff distance computed between the two polygons' versions of their common edge. That distance is compared against a threshold derived from the source data's capture precision rather than chosen by preference. Pairs measuring exactly zero are truly shared, pairs within the threshold are accepted, and pairs above it are reported together with both feature identifiers and the measured separation. STRtree over the layer buffered query near-misses included Hausdorff distance per candidate pair compare to threshold derived from capture precision distance = 0 truly shared — the vertices are identical on both sides 0 < distance ≤ threshold accepted — within the precision the source data was captured at distance > threshold reported — both feature ids plus the measured separation

Step 4 — Write the gate

import pytest

THRESHOLD_M = 1e-3          # survey-grade capture, EPSG:27700


def test_adjacent_boundaries_agree(parcels):
    assert parcels.crs.is_projected, "run this check in a projected CRS"

    failures = []
    for i, j in adjacent_pairs(parcels):
        d = boundary_disagreement(parcels.geometry.iloc[i],
                                  parcels.geometry.iloc[j])
        if d == float("inf"):
            continue                    # envelopes interact, boundaries do not
        if d > THRESHOLD_M:
            failures.append((parcels.id.iloc[i], parcels.id.iloc[j], d))

    failures.sort(key=lambda f: -f[2])
    assert not failures, (
        f"{len(failures)} adjacent pairs disagree beyond {THRESHOLD_M} m; "
        f"worst: {failures[:5]}"
    )

Sorting by descending separation before the assertion is a small thing that changes how the failure is used: the message names the five worst offenders, which are almost always one upstream problem rather than five, so the first investigation resolves the batch.

Repairing a disagreement

When the check fires, the repair is a snap to a common precision — but snapping the whole layer is a blunt instrument that moves vertices which were correct. The targeted version snaps only the pairs that failed:

from shapely import set_precision


def snap_pair(a, b, grid: float = 1e-3):
    """Force both polygons onto the same coordinate grid so the edge becomes shared."""
    return set_precision(a, grid), set_precision(b, grid)

set_precision collapses coordinates onto a fixed grid and rebuilds the geometry, which makes two near-coincident edges become identical when the disagreement is smaller than the grid. It can also collapse genuinely small features, so run the validity and area checks afterwards rather than trusting the result.

What a dissolve does to an exact boundary and to a near-coincident one Two dissolve outcomes are compared. Where two neighbouring polygons share their boundary exactly, the union removes the internal edge cleanly and produces a single polygon with no interior rings. Where the two sides disagree by a micrometre, the union preserves a hairline interior ring at every point of disagreement, producing a polygon carrying many invisible holes that no visual inspection will find and that every later area calculation, containment test and rendering pass will carry. A closing note observes that running the boundary check before the dissolve converts that outcome into a named failure with both feature identifiers attached. exact shared boundary dissolve one polygon, no interior rings — the internal edge is gone This is the outcome the boundary check exists to guarantee before a dissolve is ever run. near-coincident boundary → the union keeps a hairline interior ring at every disagreement Dozens of holes smaller than a grain of sand: invisible to inspection, carried by every later area and containment result. Check first and it is one named failure with two feature ids; check afterwards and it is a polygon with sixty rings and no cause.

Failure modes and edge cases

Corner touches are adjacency of a different kind. Two polygons meeting at a single point produce a Point intersection of their boundaries, not a LineString. boundary_disagreement will report zero because there is no sliver, which is correct but means corner-only neighbours contribute nothing to the check. Filter to pairs whose boundary intersection has non-zero length if you want line-adjacency specifically.

Multi-part polygons produce several shared edges. A MultiPolygon neighbour may share edges in two disconnected places; the intersection is then a MultiLineString and the Hausdorff distance is taken over all of it, which is the behaviour you want but makes the reported figure a maximum over several edges rather than one.

The check is quadratic in dense clusters. An index keeps the average case near linear, but a layer where thousands of small polygons share a bounding box degenerates. Cap candidate pairs per feature and report the cap rather than letting the suite run for an hour — the runtime budget pattern applies directly.

Coincident boundaries across a layer join are a separate problem. Two layers from different sources will almost never agree at the vertex level, and requiring them to is a specification error, not a data defect. Run this check within a layer; between layers, assert containment or overlap rather than boundary identity.

Reprojection breaks previously shared boundaries. Vertices that were identical in the source CRS transform to values differing in the final bits, so a layer that passed before reprojection can fail after. Run the check in the CRS the data is delivered in, and if the pipeline reprojects, run it again afterwards.

Why this check belongs upstream of the dissolve

The clearest argument for running a boundary check at all is what happens when it is skipped and a dissolve runs anyway. Dissolving a layer on an attribute merges polygons sharing that value, and the merge is a union: where the two sides agreed exactly, the internal edge disappears cleanly, and where they disagreed by a micrometre, the union keeps a hairline interior ring. The result is a polygon with dozens of invisible holes, each smaller than a grain of sand, which no visual inspection will find and which every subsequent area calculation, containment test and rendering pass will carry.

Running the boundary check before the dissolve turns that into a named failure with two feature identifiers attached. Running it afterwards means investigating a polygon with sixty interior rings and no record of which neighbours produced them. The ordering costs nothing and is the difference between a defect that names its cause and one that does not.

Conclusion

An area-based gap check answers “is there a hole”, which is not the question a shared boundary poses. Index the layer, query with a small buffer so near-misses are candidates, measure the Hausdorff distance between the two sides of each common edge, and compare it against a threshold derived from the source data’s capture precision. Repair by snapping only the failing pairs to a common grid, and re-run validity checks afterwards because snapping is itself a transformation.