Spatial Pipeline EngineeringMRV Automation

Handling Anti-Meridian Wrapping in Raster Mosaics

Any MRV pipeline that composites satellite tiles over Fiji, the Aleutians, Chukotka, or the New Zealand EEZ eventually hits the anti-meridian: the ±180° seam where longitude discontinuously flips sign. This guide is the troubleshooting recipe under Temporal Aggregation for Land-Use Change, the compositing discipline within the Satellite Imagery Processing for Emissions Tracking stack. It sits directly alongside monthly temporal aggregation of NDVI for land cover change, because the same tile stack that feeds a monthly composite will silently produce a global-width mosaic the moment it straddles 180°.

The symptom is unmistakable once you know it: a project spanning a few hundred kilometres of Fijian coastline reports a bounding box nearly 360° wide, the mosaic renders as two thin slivers pinned to opposite edges of the map, and the reprojection step either explodes memory or returns an almost-entirely-empty raster. The root cause is not a corrupt tile — it is that EPSG:4326 treats +179.9° and -179.9° as maximally distant when they are physically adjacent. Left unfixed, the defect propagates into area calculations that underpin credit issuance, so it must be caught before the composite reaches the wider Cloud-Optimized Geospatial Formats export path.

Anti-meridian mosaic decision and correction flow Tile footprints enter a detection gate that asks whether the combined extent crosses plus or minus 180 degrees, tested by a bounding-box width greater than 180 or geometries carrying both positive and negative longitudes. If no, the tiles mosaic directly. If yes, one of three corrections applies: shift longitudes to a zero to three-hundred-sixty frame, split geometries at the anti-meridian with shapely, or project to a local UTM or equal-area CRS. All corrected paths converge on a mosaic step, then a validation gate that confirms the bounding-box width is sane before emitting an audited mosaic manifest. Tile footprints 0–360° shift Split at ±180° Project to local CRS Mosaic STAC bboxes longitude re-frame shapely MultiPolygon UTM / EPSG:6933 merge tiles crosses ±180°? width sane? Audited mosaic bbox · CRS · lineage yes no

Root Cause Analysis

The anti-meridian defect is a consequence of representing a spherical coordinate as an unbounded scalar. In EPSG:4326, longitude is a wrapped quantity that discontinuously jumps from +180° to -180°, yet every downstream library — rasterio, shapely, geopandas — treats it as an ordinary real number on a flat plane. Three distinct failure modes follow.

First, bounding boxes explode to near-global width. When you take the union of a tile at +179.7° and a neighbouring tile at -179.6°, the naive extent runs from -179.6° to +179.7° — a reported span of 359.3° for two tiles that physically touch. Any mosaic driver handed that extent allocates a raster wide enough to cover the entire planet at native resolution, so a 10 m Sentinel-2 composite over a Fijian district demands tens of gigabytes of empty pixels and either exhausts memory or writes a mosaic that is 99% nodata with two thin data slivers pinned to the left and right edges.

Second, geometries wrap the wrong way. A project polygon drawn across the seam, if stored with un-normalised longitudes, becomes a shape whose vertices jump from +179° to -179°. shapely connects those vertices along the short numeric path, drawing an edge straight across the entire globe rather than across the narrow strait. Spatial joins against such a geometry match parcels on the wrong side of the planet, and area computations return absurd values — the exact silent corruption that preventing Scope 3 double-counting in spatial joins is designed to catch.

Third, UTM zone discontinuity and CRS ambiguity break area math. The anti-meridian is the boundary between UTM zones 60 and 1, so tiles a few kilometres apart carry different projected CRS codes. A mosaic that assumes a single zone smears the far side by tens of metres, and because EPSG:4326 degrees are not equal-area, any area-weighted carbon figure derived from an un-projected wrapped extent is quietly wrong. This is why disciplined CRS alignment is a precondition, not an afterthought, for seam-crossing regions.

Diagnostic Pipeline / Pre-Flight Validation

Before any merge runs, inspect the tile footprints and detect the crossing explicitly. The pre-flight gate confirms a machine-readable CRS on every input, then flags an anti-meridian crossing using two independent signals: a combined bounding-box width that exceeds 180° in geographic coordinates, or individual footprints that carry both positive and negative longitudes. Reporting both keeps the diagnosis robust when only some tiles straddle the seam.

import geopandas as gpd
import numpy as np
import structlog
from shapely.geometry import box

log = structlog.get_logger()

GEOGRAPHIC_CRS = "EPSG:4326"
WIDTH_GATE_DEG = 180.0  # a real project rarely spans more than a hemisphere


def diagnose_antimeridian(footprints: gpd.GeoDataFrame) -> dict:
    """Inspect tile footprints for anti-meridian crossing before mosaicking.

    Detects: missing/mismatched CRS, a combined extent wider than 180 degrees,
    and individual footprints spanning both longitude signs.
    """
    issues: list[str] = []

    if footprints.crs is None:
        raise ValueError("untagged footprints; refusing to guess a datum.")
    if str(footprints.crs) != GEOGRAPHIC_CRS:
        # Detection logic below assumes geographic degrees.
        footprints = footprints.to_crs(GEOGRAPHIC_CRS)

    minx, miny, maxx, maxy = footprints.total_bounds
    combined_width = float(maxx - minx)

    # Signal 1: the union spans an implausible fraction of the globe.
    width_flag = combined_width > WIDTH_GATE_DEG

    # Signal 2: at least one footprint carries both signs of longitude.
    per_tile_span = footprints.geometry.bounds
    sign_flag = bool(
        ((per_tile_span["minx"] < -170) & (per_tile_span["maxx"] > 170)).any()
    )

    crosses = width_flag or sign_flag
    if width_flag:
        issues.append(f"bbox_width_exceeds_gate:{combined_width:.2f}deg")
    if sign_flag:
        issues.append("footprint_spans_both_longitude_signs")

    report = {
        "n_tiles": int(len(footprints)),
        "combined_bbox_width_deg": round(combined_width, 3),
        "crosses_antimeridian": crosses,
        "issues": issues,
    }
    log.warning("mosaic.antimeridian_preflight", **report) if crosses else \
        log.info("mosaic.antimeridian_preflight", **report)
    return report

A report with crosses_antimeridian=True does not abort the run — it routes the tiles to the correction path below and records why, so the seam handling is visible to a verifier rather than being an undocumented reprojection buried in the mosaic driver.

Deterministic Transformation Logic

Two corrections resolve the wrap, and the right one depends on the downstream product. For a raster mosaic you normalise longitudes into a continuous 0–360° frame so the tiles become numerically contiguous; for vector geometries that must stay valid in EPSG:4326 you split the shape at ±180° into a MultiPolygon. In both cases the final compositing and any area math happen in a projected, area-honest CRS — a local UTM zone for a compact project, or the global equal-area EPSG:6933 for a wide EEZ — never in raw degrees.

The routine below normalises footprint longitudes into the 0–360° frame, splits any seam-crossing geometry with shapely, reprojects to the chosen projected CRS, and re-validates the extent through a hard gate before returning the corrected frame.

import geopandas as gpd
import numpy as np
import pyproj
import structlog
from datetime import datetime, timezone
from shapely.geometry import MultiPolygon, Polygon, box
from shapely.ops import split, transform

log = structlog.get_logger()

GEOGRAPHIC_CRS = "EPSG:4326"
EQUAL_AREA_CRS = "EPSG:6933"       # global equal-area fallback for wide extents
SANE_WIDTH_GATE_M = 5_000_000.0    # projected extents wider than this are rejected


def _shift_to_0_360(x: float, y: float, z=None):
    """Map longitudes from [-180, 180) into [0, 360) so the seam is contiguous."""
    return (x % 360.0), y


def _split_at_antimeridian(geom: Polygon) -> MultiPolygon:
    """Cut a seam-crossing polygon into east/west parts valid in EPSG:4326."""
    seam = box(180.0, -90.0, 180.0, 90.0)  # degenerate line at +180 in 0-360 frame
    shifted = transform(_shift_to_0_360, geom)
    pieces = split(shifted, seam.exterior)
    parts = []
    for part in pieces.geoms:
        # Return east-of-seam pieces to negative longitudes.
        parts.append(transform(lambda x, y, z=None:
                               (x - 360.0 if x > 180.0 else x, y), part))
    return MultiPolygon(parts) if len(parts) > 1 else MultiPolygon([parts[0]])


def correct_mosaic_extent(
    footprints: gpd.GeoDataFrame,
    strategy: str = "shift",
    projected_crs: str | None = None,
) -> tuple[gpd.GeoDataFrame, dict]:
    """Normalise seam-crossing footprints and reproject before mosaicking.

    strategy:
      - "shift":  re-frame longitudes to 0-360 (preferred for raster mosaics)
      - "split":  cut geometries at +/-180 into MultiPolygons (vector-safe)
    Returns: (corrected_footprints, audit_manifest)
    """
    if str(footprints.crs) != GEOGRAPHIC_CRS:
        footprints = footprints.to_crs(GEOGRAPHIC_CRS)

    if strategy == "shift":
        corrected = footprints.copy()
        corrected["geometry"] = corrected.geometry.apply(
            lambda g: transform(_shift_to_0_360, g))
    elif strategy == "split":
        corrected = footprints.copy()
        corrected["geometry"] = corrected.geometry.apply(_split_at_antimeridian)
    else:
        raise ValueError(f"Unknown strategy: {strategy}. Use 'shift' or 'split'.")

    # Choose an area-honest projected CRS before any width validation or merge.
    target = projected_crs or EQUAL_AREA_CRS
    corrected = corrected.set_crs(GEOGRAPHIC_CRS, allow_override=True).to_crs(target)

    minx, miny, maxx, maxy = corrected.total_bounds
    projected_width_m = float(maxx - minx)

    # Validation gate: a corrected extent must be physically plausible.
    if projected_width_m > SANE_WIDTH_GATE_M or projected_width_m <= 0:
        raise RuntimeError(
            f"corrected extent width {projected_width_m:.0f} m is implausible; "
            "anti-meridian correction did not resolve the wrap.")

    manifest = {
        "pipeline_version": "1.3.0-mrv",
        "strategy": strategy,
        "source_crs": GEOGRAPHIC_CRS,
        "projected_crs": target,
        "n_tiles": int(len(corrected)),
        "corrected_width_m": round(projected_width_m, 1),
        "crossed_antimeridian": True,
        "generated_utc": datetime.now(timezone.utc).isoformat(),
        "compliance_standard": "GHG Protocol LULUCF Activity Data",
    }
    log.info("mosaic.extent_corrected", **manifest)
    return corrected, manifest

The shift strategy is the default for raster compositing: once longitudes live in a contiguous 0–360° frame, the mosaic driver sees one narrow extent and allocates a sensibly sized output. The split strategy is reserved for vector products that must round-trip through EPSG:4326 — a GeoJSON or GeoParquet footprint that a verifier will open in a standards-compliant reader, which under RFC 7946 requires anti-meridian-crossing polygons to be split into a MultiPolygon rather than stored with wrapped coordinates.

Compliance Gating & Audit Trail Generation

The correction routine embeds a machine-readable manifest recording the strategy applied, the source and projected CRS, and the corrected extent width. That record is what turns a reprojection from an invisible driver default into an auditable decision. The gates it enforces are:

  1. Extent-sanity gate. A corrected projected width exceeding SANE_WIDTH_GATE_M, or collapsing to zero, raises rather than proceeding — a wrap that survives correction can never reach the composite.
  2. CRS declaration. Every input is refused unless it carries a machine-readable datum, and the area-honest projected CRS is stamped into the manifest so an auditor can confirm areas were computed on an equal-area or UTM grid, not in degrees.
  3. Strategy transparency. The strategy and crossed_antimeridian fields make the seam handling explicit, so a verifier can trace exactly which tiles were re-framed and how, satisfying the reproducibility expectation of ISO 14064-3.

These attributes flow downstream into MRV data lineage and provenance tracking, where the corrected extent and its CRS become part of the queryable record an auditor traces. For the authoritative convention on encoding seam-crossing geometries, consult RFC 7946 §3.1.9 (Antimeridian Cutting).

Production Integration

Deploy the correction within the tile-processing framework following a fixed ingest → diagnose → transform → validate → export → submit sequence:

  1. Ingest. Query the STAC API for tile footprints over the reporting extent, keeping each item’s declared CRS. For regions near the seam, request footprints in EPSG:4326 so the diagnostic works in geographic degrees.
  2. Diagnose. Run diagnose_antimeridian to test the combined bbox width and per-tile longitude signs; route any crossing to the correction path and log the report.
  3. Transform. Call correct_mosaic_extent with strategy="shift" for a raster mosaic (or "split" for a vector product), passing a local UTM zone for compact projects or letting it fall back to EPSG:6933.
  4. Validate. Let the extent-sanity gate assert the corrected width is plausible; reject any composite whose extent did not resolve.
  5. Export. Merge the corrected tiles and serialize to a Cloud-Optimized GeoTIFF with the manifest attached, following the Cloud-Optimized Geospatial Formats conventions so the audited extent survives downstream reads.
  6. Submit. Forward the composite into the temporal aggregation stage — the same monthly reducers described in monthly temporal aggregation of NDVI for land cover change — and carry the lineage through to registry submission.

By detecting the crossing on two independent signals, choosing the correction that matches the product, and gating the corrected extent before any merge, anti-meridian handling turns a class of silent, area-corrupting mosaic failures into a documented, reproducible step. The Fijian district that once reported a 359° bounding box now composites into a compact, area-honest mosaic whose every hectare an auditor can recompute — the deterministic foundation automated MRV compliance depends on.