Handling Datum Shifts Between NAD83 and WGS84 in pyproj

Part of the CRS & Projection Management guide.

Operative rule: EPSG:4269 is not EPSG:4326 — relabelling one as the other is a silent one-to-two metre error, and only a datum transformation removes it.

Why Two Systems of Degrees Disagree

A datum is a model of the Earth’s shape and its orientation in space; a coordinate system built on one datum places the same physical point at slightly different coordinates from one built on another. NAD83 and WGS84 were nearly identical when defined and have since drifted apart, partly because WGS84 is tied to a global reference frame that tracks plate motion and NAD83 is fixed to the North American plate.

The result is a discrepancy that grows over time and varies by location — currently around one to two metres across the contiguous United States, and larger at the continental margins. Because both systems report plain degrees, nothing about the numbers looks wrong. A file labelled EPSG:4269 loaded as EPSG:4326 produces a map that renders perfectly and is quietly displaced.

Two datums, one point, two sets of degrees NAD83 is anchored to the North American plate and moves with it. WGS84 is tied to a global frame that tracks plate motion. A single physical point therefore has two slightly different coordinates, differing by one to two metres today and growing. Same place, two coordinate values, both in degrees EPSG:4269 — NAD83 anchored to the North American plate moves with the continent what most US public data ships in EPSG:4326 — WGS84 tied to a global reference frame tracks plate motion globally what GeoJSON and web maps assume difference today: ≈ 1–2 m across the contiguous US, varying by location and growing invisible on a regional choropleth · decisive for parcels, assets and right-of-way

Production-Ready Implementation

from __future__ import annotations

import geopandas as gpd
import pyproj
from pyproj.transformer import TransformerGroup

MAX_EXPECTED_SHIFT_M = 3.0


def describe_transform(src: str, dst: str, area=None) -> str:
    """Name the transformation PROJ will actually use — not the one you hoped for."""
    group = TransformerGroup(src, dst, area_of_interest=area)
    if not group.transformers:
        raise RuntimeError(f"no transformation available from {src} to {dst}")
    if group.unavailable_operations:
        # A grid is missing; PROJ will silently use a coarser operation instead.
        missing = ", ".join(op.name for op in group.unavailable_operations[:3])
        print(f"warn: better operations unavailable (missing grids): {missing}")
    return group.transformers[0].description


def to_wgs84(gdf: gpd.GeoDataFrame) -> gpd.GeoDataFrame:
    """Transform NAD83 (or anything) to WGS84 THROUGH the datum, and verify."""
    if gdf.crs is None:
        raise ValueError("source CRS undefined — set it before transforming")

    before = gdf.geometry.representative_point()
    out = gdf.to_crs(epsg=4326)
    after = out.geometry.representative_point()

    # Measure the shift in metres with a geodesic, not with degree arithmetic.
    geod = pyproj.Geod(ellps="WGS84")
    _, _, distances = geod.inv(
        before.x.to_numpy(), before.y.to_numpy(),
        after.x.to_numpy(), after.y.to_numpy(),
    )
    worst = float(max(abs(d) for d in distances)) if len(distances) else 0.0
    if worst > MAX_EXPECTED_SHIFT_M:
        raise ValueError(
            f"datum shift of {worst:.2f} m exceeds the expected maximum — "
            "the source CRS is probably mislabelled"
        )
    print(f"transform: {describe_transform(gdf.crs.to_string(), 'EPSG:4326')}")
    print(f"max shift: {worst:.2f} m")
    return out


if __name__ == "__main__":
    parcels = gpd.read_file("parcels_nad83.geojson")   # declared EPSG:4269
    to_wgs84(parcels).to_file("parcels_wgs84.geojson", driver="GeoJSON")

Measuring the shift and asserting a maximum turns an invisible failure into a loud one. A transform that moves features by four hundred kilometres means the source label was wrong; one that moves them by zero means no datum transformation happened at all, which is equally suspicious when the source is NAD83.

Grids, Fallbacks and Reproducibility

PROJ chooses between several possible operations for the same pair of systems, preferring the most accurate one whose grid files are available. When a grid is missing it falls back — silently — to a coarser operation, so two machines with different PROJ data produce different answers for identical input.

Transformation selection with and without grids With grids installed PROJ uses a grid-based operation accurate to a few centimetres. Without them it falls back to a coarse operation accurate to about a metre, or in the worst case to a null transformation that changes nothing at all. No warning is raised by geopandas in either case. The same call, three possible operations grids installed → grid-based operation accurate to a few centimetres what you want, and what CI must guarantee grids missing → coarse fallback accurate to about a metre no exception, no warning from geopandas — only TransformerGroup reveals it worst case → null transformation coordinates unchanged, error preserved a measured shift of exactly zero from NAD83 is the signature of this case

For a pipeline that anyone will ever need to reproduce, pin the PROJ data version in the environment image, log the transformation description on every run into the run record, and fail the build when a required grid is unavailable rather than accepting the fallback.

Does It Matter for Your Map?

When a two-metre difference matters For regional choropleths the difference is far below one pixel and can be ignored. For street-level operational maps it is roughly one pixel at maximum zoom and rarely matters. For asset, parcel and right-of-way work it is decisive. For survey and engineering use it is unacceptable and the transformation is mandatory. Two metres, four contexts regional choropleth far below one pixel — safe to ignore, but label the CRS honestly anyway street-level operations about one pixel at z18 — rarely decisive, occasionally embarrassing assets, parcels, right-of-way decisive — two metres crosses boundaries and moves assets off-lot survey and engineering unacceptable — the transformation is mandatory, and so is recording it

Verification Steps

  • Print the transformation description on every run and confirm it names a grid-based operation, not a fallback.
  • Measure the maximum shift geodesically and assert it is within the expected band — neither zero nor kilometres.
  • Run the same transform in CI and locally and compare outputs to sub-centimetre precision.
  • Check a known control point against its published coordinates in both systems.
  • Confirm the source CRS was declared by the provider rather than assumed by your loader.

Common Errors & Fixes

The shift is exactly zero

No datum transformation was applied — either the source was already WGS84, or PROJ used a null transformation because grids were unavailable. Inspect TransformerGroup to see which.

Features move by hundreds of kilometres

The source label was wrong, not the datum. This is the UTM zone class of error, not a datum issue.

to_crs is slow on large frames

The transformer is being rebuilt per call inside a loop. Transform the whole GeoDataFrame in one pass; geopandas constructs the transformer once and vectorises.

Results differ between a container and a laptop

Different PROJ data. Pin the version in the image, and treat an unavailable operation as a build failure rather than a warning.

Recording the Decision, Not Just the Result

Datum handling is the part of a pipeline most likely to be revisited by someone who was not there when it was written. A year later, a colleague looking at a transformation step cannot tell from the code whether the two-metre shift was applied deliberately, applied accidentally, or skipped because someone decided it did not matter. All three produce code that looks identical.

The remedy is to write the decision down where the code lives. Three facts are enough: which datum the source is in and how that was established, whether a transformation is applied, and what accuracy the dashboard requires. A short comment block above the transform, repeated in the pipeline’s own documentation, turns an invisible assumption into a reviewable statement.

The second half of the practice is to record what actually happened at run time. The transformation description PROJ selected, the maximum shift measured, and the version of the PROJ data in use all belong in the run record. That combination answers the question that arrives during an incident — “did this change?” — without anyone needing to reconstruct the environment.

There is also a practical argument for asserting rather than trusting. Source metadata is frequently wrong: files labelled WGS84 that are plainly NAD83, shapefiles whose projection file was written by hand, exports whose CRS was inherited from a template. An assertion on the measured shift catches all three, because each produces a distance that does not match the transformation that was supposed to occur. A file that claims to be NAD83 and moves by zero metres is telling you something, and the pipeline should say so out loud rather than continuing.

Finally, be explicit about what “good enough” means for this particular dashboard. Two metres is nothing on a national choropleth and decisive on a parcel map. Writing that judgement down — with the use case, not just the number — is what lets the next person extend the dashboard without silently inheriting a tolerance that was chosen for something else entirely.

Gotchas & Edge Cases

  • NAD83 has several realisations, and they differ from one another by centimetres to decimetres; for most dashboards that is noise, but for survey work the realisation matters as much as the datum name.
  • A projection file written by hand — or copied from another dataset — is a claim, not evidence. Treat the declared CRS as a hypothesis and check it against a known control point before trusting it.
  • Transforming a large frame twice in a pipeline (once to a projected system, once back) accumulates rounding; do it once, from the source CRS to the target, in a single call.
  • Grid files are versioned and occasionally updated. Pinning the PROJ data version means a rebuild produces the same numbers next year, which is what makes an archived artifact reproducible.
  • A dataset assembled from several sources may carry several datums. Transform each source before concatenating, never after — a merged frame has one CRS attribute and no memory of where each row came from.