Automating Coordinate Validation with Python and Shapely

This page builds a single, runnable Python function that screens every incoming store coordinate for WGS84 legality, swapped axes, default-geocoder artifacts, and operational-region containment before the record is allowed into an analytics pipeline.

When location intelligence teams onboard hundreds of lease candidates from fragmented CRM exports or third-party broker feeds, manual coordinate checking becomes the bottleneck that silently corrupts trade area modeling, demographic overlays, and cannibalization forecasts. A deterministic validation layer — built with Shapely’s geometry engine and an explicit coordinate reference system — catches those defects at the front door instead of after they have propagated into financial models. This task is one concrete implementation of the broader Data Validation Rules for Store Coordinates, and it sits inside the wider Location Intelligence Architecture & Data Foundations framework as the gate every record passes through before reaching the spatial database.

Prerequisites #

Before running the validator you need a Python 3.10+ environment with the geometry and array libraries installed, plus an optional reference polygon describing where your stores can legitimately exist. The pipeline writes nothing to disk; it consumes a DataFrame (or single coordinate pair) and emits a structured result you can route downstream.

pip install "shapely>=2.0" "pyproj>=3.4" "numpy>=1.24" "pandas>=2.0"

The operational boundary is the only domain-specific input. For a Midwest retailer it might be the dissolved outline of the states you operate in; for a national chain it is the CONUS landmass. Storing that geometry alongside curated site data — for example as GeoParquet in your data lake — keeps the validation contract reproducible across pipeline runs.

Configuration and execution parameters #

The validator exposes a small set of flags that control how strict each gate is. All comparisons assume the input is EPSG:4326 (WGS84) — latitude in degrees on [-90, 90], longitude on [-180, 180] — and that convention is asserted, not assumed, so a later reprojection cannot misread the axis order.

The axis-swap heuristic exploits the geometry of the continental US: stores sit at negative longitudes (roughly -66 to -125) and positive latitudes (24 to 49), so a swapped CONUS point usually pushes lat above 90 — an unambiguous signal that the pair was written (lon, lat). Outside CONUS you should narrow the heuristic with your operational_boundary rather than trusting the swap alone.

Annotated implementation #

The function performs five gates in order: type coercion, numeric safety, WGS84 bounds (with optional swap repair), Null Island filtering, and containment. It returns a dataclass carrying a machine-readable error_code, a human message, and cleaned coordinates — ready to branch on in an ingestion service.

import logging
from dataclasses import dataclass
from typing import Optional

import numpy as np
from pyproj import CRS
from shapely.geometry import Point, Polygon

logging.basicConfig(
    level=logging.INFO,
    format="%(asctime)s | %(levelname)s | %(module)s | %(message)s",
)
logger = logging.getLogger(__name__)

# Assert the CRS once at module load so every comparison below is unambiguous.
WGS84 = CRS.from_epsg(4326)
assert WGS84.is_geographic, "Validator assumes geographic lat/lon degrees"


@dataclass
class ValidationResult:
    is_valid: bool
    error_code: Optional[str] = None
    message: Optional[str] = None
    cleaned_lat: Optional[float] = None
    cleaned_lon: Optional[float] = None


def validate_coordinate(
    lat: float,
    lon: float,
    operational_boundary: Optional[Polygon] = None,
    allow_null_island: bool = False,
    repair_axis_swap: bool = True,
    precision: int = 6,
) -> ValidationResult:
    """Validate one EPSG:4326 lat/lon pair. Returns a structured result."""

    # 1. Type coercion & numeric safety -----------------------------------
    try:
        lat_val = float(lat)
        lon_val = float(lon)
    except (ValueError, TypeError):
        return ValidationResult(False, "TYPE_MISMATCH", "Non-numeric coordinate values.")

    if np.isnan(lat_val) or np.isnan(lon_val) or np.isinf(lat_val) or np.isinf(lon_val):
        return ValidationResult(False, "NUMERIC_INVALID", "Coordinate is NaN or infinite.")

    # 2. WGS84 bounds enforcement (with optional swap repair) -------------
    in_bounds = -90.0 <= lat_val <= 90.0 and -180.0 <= lon_val <= 180.0
    if not in_bounds:
        looks_swapped = -180.0 <= lat_val <= 180.0 and -90.0 <= lon_val <= 90.0
        if looks_swapped and repair_axis_swap:
            logger.warning("Repairing swapped axes: (%s, %s)", lat_val, lon_val)
            lat_val, lon_val = lon_val, lat_val  # promote to (lat, lon)
        elif looks_swapped:
            return ValidationResult(
                False, "AXIS_SWAP_DETECTED",
                f"Axes appear swapped: ({lat_val}, {lon_val}). Expected (lat, lon).",
            )
        else:
            return ValidationResult(
                False, "OUT_OF_BOUNDS",
                f"Coordinates exceed WGS84 limits: ({lat_val}, {lon_val}).",
            )

    # 3. Null Island / default geocoder filter ---------------------------
    if not allow_null_island and lat_val == 0.0 and lon_val == 0.0:
        return ValidationResult(
            False, "NULL_ISLAND",
            "Coordinate defaults to (0.0, 0.0). Verify the upstream geocoder.",
        )

    # 4. Operational-boundary containment (optional) ---------------------
    # Shapely uses (x, y) = (longitude, latitude) per GIS convention.
    if operational_boundary is not None:
        if not operational_boundary.contains(Point(lon_val, lat_val)):
            return ValidationResult(
                False, "OUTSIDE_BOUNDARY",
                "Coordinate falls outside the designated operational polygon.",
            )

    # 5. Precision normalization -----------------------------------------
    return ValidationResult(
        is_valid=True,
        cleaned_lat=round(lat_val, precision),
        cleaned_lon=round(lon_val, precision),
    )

For DataFrame-scale ingestion, wrap the row-level validator so it appends structured columns you can filter and quarantine on:

import pandas as pd


def batch_validate_coordinates(
    df: pd.DataFrame,
    lat_col: str,
    lon_col: str,
    boundary: Optional[Polygon] = None,
) -> pd.DataFrame:
    """Validate every row and append result columns. See the performance note
    below for the vectorized path on datasets over 100k rows."""
    logger.info("Validating %d coordinate records.", len(df))

    results = df.apply(
        lambda r: validate_coordinate(r[lat_col], r[lon_col], boundary), axis=1
    )

    out = df.copy()
    out["validation_status"] = results.apply(lambda x: x.is_valid)
    out["error_code"] = results.apply(lambda x: x.error_code)
    out["error_message"] = results.apply(lambda x: x.message)
    out["cleaned_lat"] = results.apply(lambda x: x.cleaned_lat)
    out["cleaned_lon"] = results.apply(lambda x: x.cleaned_lon)

    passed = int(out["validation_status"].sum())
    logger.info("Done: %d/%d records passed all gates.", passed, len(df))
    return out

Failure modes and debugging #

The gates fail in predictable ways; each maps to one error_code so you can quarantine and triage records without reading every message by hand.

The most expensive recurring bug is a CRS mismatch on the boundary itself: Polygon.contains is purely numeric and will happily compare degrees against meters, rejecting every store. Always assert the boundary is geographic — assert CRS.from_user_input(boundary_crs).is_geographic — before the containment loop. The same CRS discipline governs the downstream PostGIS retail database these cleaned records flow into.

Performance #

For datasets above ~100,000 rows the DataFrame.apply pattern becomes the bottleneck because it materializes one Python ValidationResult object per row. Shapely 2.0 ships array-based geometry operations, so once the schema is stable you can replace the boundary loop with a single vectorized call:

import numpy as np
import shapely

def fast_boundary_mask(lat: np.ndarray, lon: np.ndarray, boundary: Polygon) -> np.ndarray:
    """Vectorized containment for large batches. Bounds/swap checks stay as
    cheap NumPy comparisons; only containment uses the Shapely fast path."""
    in_bounds = (lat >= -90) & (lat <= 90) & (lon >= -180) & (lon <= 180)
    inside = shapely.contains_xy(boundary, lon, lat)  # (x, y) = (lon, lat)
    return in_bounds & inside

This runs containment at near-C speed and keeps the bounds and swap logic as masked NumPy comparisons. Reserve the per-row dataclass path for sources small enough that auditable per-record messages matter more than throughput.

Verification #

Confirm the validator behaves before trusting it on production feeds. Three cheap checks catch almost every regression:

• Pass-rate sanity: out["validation_status"].mean() on a known-clean source should sit near 1.0; a sudden drop signals an upstream format change, not bad stores.
• Error histogram: out["error_code"].value_counts() should be dominated by the defects you expect (mostly AXIS_SWAP_DETECTED from a known-swapped broker, near-zero TYPE_MISMATCH from a typed source).
• Bounding-box reality check: after validation, out["cleaned_lat"].agg(["min", "max"]) and the same for longitude must fall inside your operational region — a max latitude of 80 in a US-only dataset means a swap slipped through.

For a stricter gate, route the cleaned coordinates through a held-out reference set and confirm they land in the right polygons, the same approach described in validating spatial join accuracy with ground truth. Coordinates that clear all five gates are then safe inputs for catchment work such as point-in-polygon joins for store catchments.

Related #

• Performing point-in-polygon joins for store catchments
• Validating spatial join accuracy with ground truth
• Setting up PostGIS for retail analytics

← Back to Data Validation Rules for Store Coordinates