Working with EPSG and WKT2 CRS Definitions in Python

Anchoring a cadastral coordinate to an unambiguous reference frame is the first obligation of every workflow in Core Transformation Fundamentals & Standards: before a single datum shift or projection runs, the pipeline must state, in a form a reviewer can reproduce, exactly which coordinate reference system a number belongs to. A bare authority code such as EPSG:4326 names a system but leaves three properties open — the axis order, the linear or angular unit, and which member of a datum ensemble the coordinate realizes — and any of the three can silently invert or displace a boundary corner. This page narrows the standards framework to one concrete discipline: turning an EPSG authority code into a fully-specified pyproj.CRS object, inspecting its datum, ellipsoid, coordinate system, axis order and area of use, and emitting a WKT2:2019 string that freezes all of it for the permanent record. The complementary task of resolving those codes into a running transformer with the correct operation selected is covered in setting up high-precision coordinate reference systems; here the focus is the identity of the CRS itself, not the operation between two of them.

Figure — an EPSG code is a lookup key into the same datum/ellipsoid/CS/axis tree that WKT2:2019 writes out in full.

Why an Ambiguous CRS Fails a Survey Deliverable #

An authority code is a promise that a registry entry exists, not a self-contained description. The reason this distinction is not pedantic is that the two most common shorthands a surveyor inherits — the string EPSG:4326 and a bare PROJ string like +proj=longlat +datum=WGS84 +no_defs — each leave a load-bearing property undefined at the point of use.

The authoritative definition of EPSG:4326 places geodetic latitude first, longitude second, both in degrees. A great deal of software, however, ingests GeoJSON, shapefiles and GML in longitude-latitude order and tags them EPSG:4326 regardless. The code is therefore consistent with two mutually exclusive physical orderings, and the disagreement does not raise an error — it silently transposes northings and eastings, placing a parcel corner on the wrong side of a diagonal. A bare PROJ string is worse: it discards the axis metadata entirely, cannot distinguish the several realizations of a datum ensemble, and carries no authority identifier to trace back to a published record. For legally defensible work under the ISO 19111 information model, neither is acceptable as the recorded CRS of a deliverable.

WKT2:2019 closes the gap. It is a text encoding of the complete CRS object — datum, ellipsoid with its defining parameters, coordinate system with an explicit axis list and direction, unit, and an embedded authority ID — so a reader recovers the exact frame without a registry lookup. The correct pattern is to treat the EPSG code as an ingress key, resolve it once into a pyproj.CRS, assert the axis order and datum you expected, and then serialize WKT2 as the frozen record that travels with the coordinates. That principle — resolve codes to explicit WKT2 to remove implicit axis and operation ambiguity — is the same one that governs the whole reference and reappears whenever CRS mismatches must be reconciled in cadastral datasets.

The CRS Object Model #

A pyproj.CRS mirrors the ISO 19111 object graph. Understanding the four nodes the code must inspect prevents the class of error where two coordinates that “have the same EPSG code” are nonetheless incompatible:

• Datum — the physical realization that ties the coordinate system to the Earth. NAD83 and WGS 84 differ by up to roughly two metres horizontally and the offset is not uniform; the difference between two realizations of NAD83 (for example the CORS96 and 2011 adjustments) is smaller but still exceeds cadastral tolerance in active-deformation zones.
• Ellipsoid — owned by the datum, defined by a semi-major axis $a$ and either the inverse flattening $1/f$ or the semi-minor axis $b$. The semi-minor axis follows from $b = a\,(1 - f)$, and the first eccentricity squared, which projection math consumes directly, is $e^2 = 2f - f^2$.
• Coordinate system (CS) — the abstract axes: ellipsoidal (latitude/longitude, angular) for a geographic CRS, Cartesian (easting/northing, linear) for a projected CRS.
• Axis order and direction — the ordered list of axes with their directions (north, east) and units. This is the property a bare code cannot pin down and the one most likely to corrupt a deliverable.

Two CRS objects are equivalent for cadastral purposes only when all four agree. pyproj exposes a strict equality and a looser equals(..., ignore_axis_order=True), and choosing between them deliberately is itself a compliance decision, examined in the verification section below.

Step-by-Step Implementation #

The routines below use only pyproj, dataclasses, logging and the standard library. Each step is a runnable, type-hinted unit, and each carries a comment naming the ISO 19111 or EPSG concept it enforces.

Step 1 — Build a CRS from an authority code, defensively #

from __future__ import annotations

import logging

from pyproj import CRS
from pyproj.exceptions import CRSError

logger = logging.getLogger("crs")


def load_crs(identifier: int | str) -> CRS:
    """Resolve an EPSG code (or authority string) into a fully-specified CRS.

    ISO 19111: a CRS must carry an authority identifier. We reject deprecated
    entries loudly rather than let a superseded code reach a deliverable.
    """
    crs = CRS.from_user_input(identifier)  # accepts 4326, "EPSG:4326", WKT, PROJ
    auth = crs.to_authority()              # (authority, code) or None
    if auth is None:
        logger.warning("CRS %r has no authority ID; it cannot be traced to a register", identifier)
    if crs.is_deprecated:
        # EPSG marks superseded CRSs deprecated; using one silently is a defect.
        logger.warning("CRS %s is deprecated; resolve its replacement before use", auth)
    return crs

Step 2 — Inspect datum, ellipsoid, CS and axis order #

from dataclasses import dataclass


@dataclass(frozen=True)
class CRSFingerprint:
    authority: str | None      # e.g. "EPSG:6318", or None if untraceable
    name: str
    kind: str                  # "geographic", "projected", "compound", ...
    datum_name: str
    ellipsoid_name: str
    semi_major_m: float        # a, metres
    inverse_flattening: float  # 1/f (0.0 for a sphere)
    axis_abbrevs: tuple[str, ...]     # ("Lat", "Lon") or ("E", "N")
    axis_directions: tuple[str, ...]  # ("north", "east") or ("east", "north")
    unit_names: tuple[str, ...]       # ("degree", "degree") or ("metre", "metre")


def fingerprint(crs: CRS) -> CRSFingerprint:
    """Extract the identity-bearing components of a CRS for assertion/logging."""
    ell = crs.ellipsoid
    axes = crs.axis_info  # ordered exactly as the CS declares them
    kind = (
        "geographic" if crs.is_geographic
        else "projected" if crs.is_projected
        else "compound" if crs.is_compound
        else "other"
    )
    auth = crs.to_authority()
    return CRSFingerprint(
        authority=f"{auth[0]}:{auth[1]}" if auth else None,
        name=crs.name,
        kind=kind,
        datum_name=crs.datum.name if crs.datum else "<undefined>",
        ellipsoid_name=ell.name,
        semi_major_m=float(ell.semi_major_metre),
        inverse_flattening=float(ell.inverse_flattening or 0.0),
        axis_abbrevs=tuple(ax.abbrev for ax in axes),
        axis_directions=tuple(ax.direction for ax in axes),
        unit_names=tuple(ax.unit_name for ax in axes),
    )

Step 3 — Read the area of use #

@dataclass(frozen=True)
class AreaOfUse:
    name: str
    west: float   # degrees, positive east
    south: float
    east: float
    north: float


def area_of_use(crs: CRS) -> AreaOfUse | None:
    """Return the CRS's declared domain of validity as a lon/lat bounding box.

    Applying a CRS outside its area of use is an out-of-domain error: the
    datum realization is only defined where its adjustment was computed.
    """
    aou = crs.area_of_use
    if aou is None:
        return None
    return AreaOfUse(
        name=aou.name,
        west=float(aou.west), south=float(aou.south),
        east=float(aou.east), north=float(aou.north),
    )


def contains_point(aou: AreaOfUse, lon_deg: float, lat_deg: float) -> bool:
    """True if a lon/lat point lies inside the CRS domain of validity."""
    return (aou.west <= lon_deg <= aou.east) and (aou.south <= lat_deg <= aou.north)

Step 4 — Compare two CRS for equivalence #

def crs_match(a: CRS, b: CRS, *, strict_axis_order: bool = True) -> dict[str, bool]:
    """Report component-wise agreement between two CRS objects.

    strict_axis_order=True is the survey default: a lat/lon vs lon/lat CRS
    pair is treated as NON-equivalent, because feeding one to the other's
    consumer transposes the coordinates.
    """
    return {
        "exact_equal": a == b,  # pyproj strict equality (axis order included)
        "equal_ignoring_axis_order": a.equals(b, ignore_axis_order=True),
        "same_datum": (a.datum == b.datum) if (a.datum and b.datum) else False,
        "same_ellipsoid": a.ellipsoid == b.ellipsoid,
        "usable": (a == b) if strict_axis_order
                  else a.equals(b, ignore_axis_order=True),
    }

Step 5 — Emit WKT2:2019 for the record #

from pyproj.enums import WktVersion


def to_wkt2(crs: CRS, *, pretty: bool = True) -> str:
    """Serialize a CRS to a WKT2:2019 string for the permanent deliverable.

    WKT2_2019 is requested explicitly so the output carries the modern
    GEOGCRS/PROJCRS keywords, axis directions and the embedded authority ID,
    rather than the ambiguous legacy WKT1 GEOGCS form.
    """
    wkt = crs.to_wkt(version=WktVersion.WKT2_2019, pretty=pretty)
    if not wkt.lstrip().startswith(("GEOGCRS", "PROJCRS", "GEODCRS",
                                    "BOUNDCRS", "COMPOUNDCRS", "VERTCRS",
                                    "PROJCS")):
        raise ValueError("Serialized WKT does not start with a WKT2 CRS keyword")
    return wkt

The parsing side of this round trip — reading such a WKT2 string back into a validated CRS and confirming it is WKT2 rather than legacy WKT1 — is a self-contained routine in parsing WKT2:2019 CRS strings in Python.

Parameter and Return-Value Reference #

Worked Example: NAD83(2011) Geographic and its UTM Projection #

Take two systems a U.S. cadastral job routinely pairs: EPSG:6318, NAD83(2011) geographic 2D, and EPSG:6339, NAD83(2011) / UTM zone 10N in metres. Resolving both, inspecting them, and serializing WKT2 makes their differences explicit rather than assumed.

geographic = load_crs("EPSG:6318")   # NAD83(2011) lat/lon, degrees
projected = load_crs("EPSG:6339")    # NAD83(2011) / UTM zone 10N (metre)

fp_geo = fingerprint(geographic)
fp_proj = fingerprint(projected)

print(fp_geo.kind, fp_geo.axis_abbrevs, fp_geo.unit_names)
# geographic ('Lat', 'Lon') ('degree', 'degree')
print(fp_proj.kind, fp_proj.axis_abbrevs, fp_proj.unit_names)
# projected ('E', 'N') ('metre', 'metre')

print(round(fp_geo.semi_major_m, 3), round(fp_geo.inverse_flattening, 9))
# 6378137.0 298.257222101   (GRS 1980, shared by both)

aou = area_of_use(projected)
print(contains_point(aou, lon_deg=-121.5, lat_deg=37.3))
# True   -> the point lies inside the zone's domain of validity

match = crs_match(geographic, projected)
print(match["same_datum"], match["usable"])
# True False   -> same datum, but NOT interchangeable (geographic vs projected)

The pair shares the GRS 1980 ellipsoid and the NAD83(2011) datum, yet usable is False: one is angular latitude/longitude, the other linear easting/northing, and treating them as interchangeable because “they are both NAD83(2011)” is exactly the confusion the fingerprint is designed to surface. Note too that the geographic frame’s own domain of validity spans the Aleutian antimeridian, so a naive longitude-bounds test must be reserved for a single-zone extent like the projected CRS above. Emitting the WKT2 record makes the contrast permanent:

wkt = to_wkt2(geographic)
print(wkt.splitlines()[0])
# GEOGCRS["NAD83(2011)",
assert 'ID["EPSG",6318]' in wkt  # authority survives into the record

Verification and CRS Equality Audit #

CRS handling is verified not against a control monument but against the intent recorded at ingress. The check that matters is a two-part assertion: that the resolved axis order matches what the pipeline expects to feed its transformer, and that a CRS declared equivalent to a reference truly is. Emit a structured record so the decision is auditable.

import json


def audit_crs(crs: CRS, expected_authority: str,
              expected_axis: tuple[str, ...]) -> dict[str, object]:
    """Assert a resolved CRS matches ingress expectations; emit an audit record.

    Raises when the axis order or authority disagrees, because both faults are
    silent at transform time and only surface as a displaced coordinate.
    """
    fp = fingerprint(crs)
    record = {
        "authority": fp.authority,
        "expected_authority": expected_authority,
        "axis": list(fp.axis_abbrevs),
        "expected_axis": list(expected_axis),
        "authority_ok": fp.authority == expected_authority,
        "axis_ok": fp.axis_abbrevs == expected_axis,
        "is_deprecated": crs.is_deprecated,
    }
    logger.info("crs_audit %s", json.dumps(record))
    if not (record["authority_ok"] and record["axis_ok"] and not record["is_deprecated"]):
        raise ValueError(f"CRS audit failed: {record}")
    return record


rec = audit_crs(load_crs("EPSG:6318"), "EPSG:6318", ("Lat", "Lon"))
assert rec["axis_ok"]

The subtle failure this guards against is CRS equality itself. Strict == in pyproj treats a latitude/longitude CRS and its longitude/latitude twin as unequal, which is correct for a deliverable, whereas equals(ignore_axis_order=True) treats them as equal, which is convenient only when you have already normalized the coordinate order elsewhere. Recording which comparison you used, and why, is part of the audit trail the parent reference requires alongside every transformation. The residual-level counterpart — checking transformed coordinates against surveyed monuments once the CRS is pinned — is developed in validating datum alignment with control points.

Troubleshooting and Gotchas #

Frequently Asked Questions #

Related References #

• Parsing WKT2:2019 CRS strings in Python — the parsing and validation half of the WKT2 round trip.
• Setting up high-precision coordinate reference systems — turning resolved CRSs into a transformer with the right operation selected.
• Handling CRS mismatches in cadastral datasets — reconciling mixed or mislabelled CRS tags across a dataset.
• Validating datum alignment with control points — the residual-level check once the CRS identity is pinned.
• Core Transformation Fundamentals & Standards — the parent reference on CRS resolution, grid selection, and ISO 19111 compliance.