Spatial Pipeline EngineeringMRV Automation

Methane Plume Detection from Hyperspectral Imagery

Methane plume detection is the retrieval discipline that turns shortwave-infrared radiance into a per-pixel methane enhancement map, isolates coherent plumes from instrument and surface noise, and attributes each one to a facility — and it is the highest-leverage, highest-risk component of the satellite imagery processing stack. Highest-leverage because a single unlit flare or compressor blowdown can dominate a site’s annual inventory; highest-risk because the same spectral absorption feature that reveals methane is mimicked by bright bare soil, calcite, and paved surfaces, so a careless pipeline reports emissions where none exist. Everything here builds on the same discipline that governs Sentinel-2 and Landsat cloud masking workflows — a mask that is wrong in one direction invents signal, and wrong in the other destroys it.

Methane retrieval chain from shortwave-infrared radiance to an attributed, quantified plume Five sequential stages. Stage one, SWIR radiance from the 2100 to 2400 nanometre window. Stage two, a matched filter against the methane absorption spectrum producing a per-pixel enhancement in parts per million metre. Stage three, plume segmentation that keeps only coherent connected components above the scene noise floor. Stage four, wind-informed quantification converting the integrated mass enhancement to a kilograms per hour flux. Stage five, facility attribution and lineage. Two rejection paths leave the chain: surface-albedo false positives are discarded at segmentation, and plumes without a wind field are held rather than quantified. SWIR radiance 2100–2400 nm window EMIT · PRISMA · S2 B11/B12 Matched filter scene covariance + CH₄ target spectrum Segmentation connected components above scene noise floor Quantification integrated mass ÷ wind → kg h⁻¹ flux Attribution facility polygon join + retrieval lineage Discarded albedo artefacts Held, not quantified no wind field radiance cube ppm·m enhancement plume mask flux estimate reportable record

Role in the MRV Workflow

Methane retrieval sits in the observation layer, immediately after atmospheric correction and immediately before emissions attribution. Its upstream dependency is a radiometrically calibrated, geometrically registered scene: an EMIT or PRISMA hyperspectral cube, an EnMAP acquisition, or — for coarse screening at scale — the Sentinel-2 SWIR band pair. Its downstream consumer is the facility-level inventory, where a plume flux in kilograms per hour becomes an annualised tonnage that flows into corporate reporting through the same channels as any other activity datum, subject to the same emissions data quality validation gates.

What makes methane structurally different from land-carbon work is the time base. A forest-carbon pipeline aggregates observations over months and tolerates missing scenes because the underlying stock changes slowly. A methane plume is an instantaneous snapshot of a process that may have started forty minutes before overpass and stopped an hour after. The retrieval measures a rate, not a stock, and converting a handful of instantaneous rates into an annual mass requires an explicit, documented temporal model — persistence assumptions, detection-limit-aware upper bounds for non-detections, and a duty-cycle estimate for intermittent sources. Pipelines that skip that step produce annual figures whose uncertainty is unbounded, and verifiers reject them.

The second structural difference is the detection limit. Every retrieval has a minimum detectable enhancement set by instrument noise, surface heterogeneity, and plume geometry, and that limit varies pixel-by-pixel across a single scene. A non-detection over dark, spectrally uniform water means something very different from a non-detection over a heterogeneous industrial yard, where the limit may be five times higher. Reporting “no emissions observed” without reporting the scene-specific detection limit alongside it is the single most common way a methane inventory becomes indefensible, and it is why the implementation below carries the per-scene noise floor through to the output record rather than discarding it after thresholding.

Finally, methane retrieval is unusually dependent on geometry that must be exactly right. The plume mask is intersected with facility polygons to assign the emission to an operator, and that join inherits every risk documented under geospatial CRS alignment. A thirty-metre georegistration error is irrelevant when mapping a ten-thousand-hectare forest parcel; it reassigns a plume from one wellpad operator to their neighbour.

Core Failure Modes

Three failure modes account for most bad methane records reaching an inventory. Each has a distinct spectral or geometric root cause and a characteristic signature you can test for.

  1. Albedo-driven false enhancements over bright, mineral-rich surfaces. The matched filter looks for a specific absorption shape in the 2200–2400 nm region, and several common surface materials — calcite, kaolinite, gypsum-rich playa, fresh asphalt, some painted metal roofs — carry absorption features that partially project onto the methane target spectrum. The retrieval reports a positive enhancement that is a property of the ground, not the atmosphere. The signature is diagnostic: a real plume is a coherent, wind-aligned structure anchored to a point source and it moves between overpasses, while an albedo artefact is pixel-locked, repeats in identical shape on every clear acquisition, and correlates with the scene’s broadband reflectance. Pipelines that threshold the enhancement map without a persistence test and an albedo correlation test typically carry 20–60% false positives over arid basins.

  2. Wind-field error dominating the flux uncertainty budget. The integrated mass enhancement (IME) method converts a plume’s total mass to a flux by dividing by an effective length and multiplying by an effective wind speed. The retrieval itself may be good to 15–25%, but a reanalysis wind product at 0.25° resolution, interpolated to an overpass time, routinely carries 40–60% error at the ten-metre level that actually transports the plume. Wind error therefore dominates the total budget, and — critically — it is not reduced by better imagery. Teams that report a flux with a tight uncertainty derived only from retrieval precision are understating their true error by a factor of two or more, which is precisely the misstatement that emission factor uncertainty mapping exists to prevent.

  3. Detection-limit blindness turned into a false zero. Non-detections are data, but only if the detection limit is recorded. A pipeline that emits rows only for detected plumes silently converts “we could not have seen anything below 400 kg h⁻¹ here” into “there was nothing here”, and an annual inventory built from those rows understates emissions by whatever fraction of sources sit below the limit — typically the long tail of small, chronic leaks that dominate site-level totals at many upstream assets. The fix is architectural rather than algorithmic: every scene-facility pair must produce a row, detection or not, carrying the computed limit.

Separating a real plume from an albedo artefact using persistence across three overpasses Two rows of three small scene panels each. The top row, labelled real plume, shows an enhancement blob that changes shape and direction with the wind arrow on each of three dates, and is anchored to a fixed source marker. The bottom row, labelled albedo artefact, shows an identical pixel-locked patch in the same position and shape on all three dates, unrelated to the wind arrows. A verdict column on the right states the test: shape correlation across dates above 0.9 with no wind alignment indicates a surface artefact and is rejected. Persistence test — the cheapest false-positive filter you can run Same facility, three clear overpasses. A plume is transported by wind; a surface feature is not. 2026-03-04 2026-03-19 2026-04-02 Verdict Real plume Albedo artefact Shape varies · aligns with wind Cross-date mask IoU 0.21 · bearing tracks the wind vector → keep, quantify Pixel-locked · identical every date → reject: cross-date IoU 0.94, no wind alignment

Deterministic Implementation Architecture

The detection stage below implements the matched-filter retrieval, the noise-floor threshold, the persistence test, and — most importantly — the rule that every scene-facility pair emits a row whether or not a plume was found. It is written as a Prefect task with structlog telemetry, explicit CRS declaration, and a hard refusal to quantify when the wind field is missing.

from dataclasses import dataclass, asdict

import geopandas as gpd
import numpy as np
import rioxarray
import structlog
import xarray as xr
from prefect import flow, task
from scipy import ndimage

log = structlog.get_logger()

CANONICAL_CRS = "EPSG:6933"      # equal-area: plume areas and IME must be honest
CH4_MOLAR_MASS = 16.04           # g mol-1
NOISE_SIGMA_GATE = 4.0           # keep pixels >= 4 sigma of the scene noise floor
MIN_PLUME_PIXELS = 9             # a coherent plume, not a speckle


@dataclass(frozen=True)
class PlumeRecord:
    """One scene x facility observation. Emitted for NON-detections too — the
    detection limit is the evidence that a zero is a real zero."""
    facility_id: str
    scene_id: str
    acquired: str
    detected: bool
    ime_kg: float | None
    flux_kg_h: float | None
    flux_uncertainty_kg_h: float | None
    detection_limit_kg_h: float
    wind_speed_m_s: float | None
    wind_source: str | None
    retrieval: str


def matched_filter(cube: xr.DataArray, target: np.ndarray) -> xr.DataArray:
    """Column-wise matched filter: enhancement in ppm-m per pixel.

    Covariance is estimated per across-track column, not per scene, because
    push-broom detectors carry column-dependent noise; a scene-wide covariance
    smears that structure into the retrieval and inflates false positives.
    """
    data = cube.values.reshape(cube.sizes["band"], -1)
    out = np.zeros(data.shape[1], dtype="float32")

    for col in range(cube.sizes["x"]):
        idx = np.arange(col, data.shape[1], cube.sizes["x"])
        block = data[:, idx]
        mu = block.mean(axis=1, keepdims=True)
        resid = block - mu
        cov = np.cov(resid) + np.eye(resid.shape[0]) * 1e-6   # ridge: keep it invertible
        inv = np.linalg.inv(cov)
        t = (target * mu.ravel())[:, None]                    # radiance-scaled target
        denom = float(t.T @ inv @ t)
        out[idx] = ((t.T @ inv @ resid) / denom).ravel()

    return xr.DataArray(
        out.reshape(cube.sizes["y"], cube.sizes["x"]),
        coords={"y": cube.y, "x": cube.x}, dims=("y", "x"), name="ch4_ppm_m",
    )


@task(retries=2, retry_delay_seconds=20)
def detect_plumes(
    scene_path: str,
    scene_id: str,
    acquired: str,
    facilities: gpd.GeoDataFrame,
    target_spectrum: np.ndarray,
    wind: dict[str, float] | None,
) -> list[dict]:
    cube = rioxarray.open_rasterio(scene_path, chunks={"x": 512, "y": 512})
    if cube.rio.crs is None:
        raise ValueError(f"{scene_id}: scene carries no CRS; refusing to geolocate a plume")

    enhancement = matched_filter(cube, target_spectrum)
    enhancement = enhancement.rio.write_crs(cube.rio.crs).rio.reproject(CANONICAL_CRS)

    # Scene noise floor from the robust spread of the enhancement field. MAD, not
    # std: a few real plumes must not raise the threshold that finds them.
    mad = float(np.nanmedian(np.abs(enhancement - np.nanmedian(enhancement))))
    sigma = 1.4826 * mad
    threshold = NOISE_SIGMA_GATE * sigma

    mask = (enhancement.values > threshold).astype("uint8")
    labels, n = ndimage.label(mask, structure=np.ones((3, 3)))
    pixel_area_m2 = abs(float(enhancement.rio.resolution()[0])) ** 2

    facilities = facilities.to_crs(CANONICAL_CRS)
    detection_limit = _detection_limit_kg_h(sigma, pixel_area_m2, wind)

    log.info("ch4.scene.retrieved", scene_id=scene_id, sigma_ppm_m=round(sigma, 2),
             threshold_ppm_m=round(threshold, 2), candidates=int(n),
             detection_limit_kg_h=round(detection_limit, 1), crs=CANONICAL_CRS)

    records: list[PlumeRecord] = []
    for _, facility in facilities.iterrows():
        plume = _largest_plume_within(labels, enhancement, facility.geometry, n)

        if plume is None or plume["pixels"] < MIN_PLUME_PIXELS:
            # A non-detection is a row. Without the limit it is an unusable zero.
            records.append(PlumeRecord(
                facility_id=facility.facility_id, scene_id=scene_id, acquired=acquired,
                detected=False, ime_kg=None, flux_kg_h=None, flux_uncertainty_kg_h=None,
                detection_limit_kg_h=detection_limit,
                wind_speed_m_s=(wind or {}).get("speed_m_s"),
                wind_source=(wind or {}).get("source"), retrieval="matched-filter/v2",
            ))
            continue

        ime_kg = plume["sum_ppm_m"] * pixel_area_m2 * CH4_MOLAR_MASS * 1e-9 * 40.87

        if wind is None:
            # Hold, never guess. An unquantified detection is honest; an invented
            # wind speed is a number an auditor cannot trace to an observation.
            log.warning("ch4.flux.withheld", facility_id=facility.facility_id,
                        scene_id=scene_id, reason="no_wind_field", ime_kg=round(ime_kg, 1))
            records.append(PlumeRecord(
                facility_id=facility.facility_id, scene_id=scene_id, acquired=acquired,
                detected=True, ime_kg=round(ime_kg, 1), flux_kg_h=None,
                flux_uncertainty_kg_h=None, detection_limit_kg_h=detection_limit,
                wind_speed_m_s=None, wind_source=None, retrieval="matched-filter/v2",
            ))
            continue

        length_m = np.sqrt(plume["pixels"] * pixel_area_m2)
        u_eff = 0.34 * wind["speed_m_s"] + 0.44          # IME effective-wind scaling
        flux = ime_kg * u_eff / length_m * 3600.0

        # Wind dominates: combine retrieval and wind error in quadrature.
        rel = np.sqrt(0.20 ** 2 + wind.get("rel_error", 0.5) ** 2)

        log.info("ch4.plume.quantified", facility_id=facility.facility_id,
                 scene_id=scene_id, flux_kg_h=round(flux, 1),
                 rel_uncertainty=round(rel, 2), pixels=plume["pixels"])

        records.append(PlumeRecord(
            facility_id=facility.facility_id, scene_id=scene_id, acquired=acquired,
            detected=True, ime_kg=round(ime_kg, 1), flux_kg_h=round(flux, 1),
            flux_uncertainty_kg_h=round(flux * rel, 1),
            detection_limit_kg_h=detection_limit, wind_speed_m_s=wind["speed_m_s"],
            wind_source=wind["source"], retrieval="matched-filter/v2",
        ))

    return [asdict(r) for r in records]


@flow(name="ch4_plume_detection")
def run(scenes: list[dict], facility_path: str, target_spectrum_path: str) -> list[dict]:
    facilities = gpd.read_file(facility_path)
    if facilities.crs is None:
        raise ValueError("facility layer has no CRS; attribution would be unverifiable")
    target = np.load(target_spectrum_path)

    rows: list[dict] = []
    for scene in scenes:
        rows.extend(detect_plumes(
            scene["path"], scene["scene_id"], scene["acquired"],
            facilities, target, scene.get("wind"),
        ))
    log.info("ch4.run.complete", scenes=len(scenes), rows=len(rows),
             detections=sum(1 for r in rows if r["detected"]))
    return rows

Four decisions in that code are the ones worth defending in a design review. Column-wise covariance rather than scene-wide, because push-broom instruments carry per-column noise structure that a scene covariance smears into apparent enhancement. MAD-based noise estimation rather than standard deviation, so a scene containing several genuine plumes does not raise the threshold that would find them. Withholding rather than guessing when the wind field is absent — an unquantified detection is a defensible record, an invented wind speed is not. And a row for every scene-facility pair, which is what turns non-detections from silence into evidence.

Validation, Debugging & Compliance Mapping

Methane records enter the inventory as activity data and are scrutinised on three axes: whether the detection is real, whether the flux is defensible, and whether the annualisation from instantaneous rates is documented.

  • Detection validity → controlled-release calibration. The standard evidence is performance against a blind controlled-release experiment, in which known volumes are released and the pipeline’s detection rate and quantification bias are measured without the operator knowing the true values. Record the resulting detection-limit curve and quantification bias per instrument and surface class, and cite it in the methodology annex. A retrieval with no controlled-release provenance is, to a verifier, an unvalidated model.
  • Flux uncertainty → honest error propagation. The uncertainty attached to each flux must include wind error, not just retrieval precision. The quadrature combination above is the minimum; where a site drives a material fraction of a reported inventory, propagate the wind distribution through a Monte Carlo pass in the manner set out under Monte Carlo uncertainty propagation for emission factors.
  • Annualisation → an explicit temporal model. Converting N instantaneous observations into an annual mass requires a stated model: persistent-source assumption, duty-cycle estimate, or a survival-analysis treatment of intermittent sources. Whichever you choose, the assumption belongs in the record’s lineage, wired through MRV data lineage and provenance tracking, so an auditor can re-derive the annual figure from the same observations.
  • Attribution → geometry you can defend. The facility join must record the buffer distance used, the CRS, and the tie-break rule when two facilities are within the plume’s uncertainty ellipse. Undocumented attribution is where methane inventories most often become contested rather than merely uncertain.

The detection limit is not one number per instrument — it is a curve that varies with surface heterogeneity and wind. Reporting a single headline sensitivity flatters the pipeline and misleads the reader; reporting the curve, and the surface class each observation fell into, is what makes a non-detection quantitative.

Minimum detectable methane flux against wind speed, for three surface classes A line chart with wind speed from 1 to 9 metres per second on the horizontal axis and minimum detectable flux from 0 to 900 kilograms per hour on the vertical axis. Three rising curves: uniform dark surfaces such as water and dense canopy stay lowest, from about 60 to 300 kilograms per hour; moderate surfaces such as grassland and cropland run from about 120 to 520; heterogeneous bright surfaces such as an industrial yard or arid playa run highest, from about 240 to 860. A shaded band marks the region below 200 kilograms per hour where most chronic leaks sit, showing that only dark uniform surfaces resolve them at moderate wind. Minimum detectable flux is a curve, not a constant Representative hyperspectral retrieval at 30 m ground sampling; higher wind disperses the plume and raises the floor. 0 200 400 600 800 1 3 5 7 9 wind speed (m s⁻¹) min. detectable flux (kg h⁻¹) chronic-leak band (< 200 kg h⁻¹) dark, uniform water · dense canopy moderate grassland · cropland bright, heterogeneous industrial yard · playa

For debugging, three diagnostics resolve most incidents. Plot the enhancement field against broadband SWIR reflectance for the scene: a correlation above roughly 0.4 signals albedo contamination rather than atmosphere. Compute cross-date mask intersection-over-union for repeat detections at a fixed location, as in the persistence figure above; values near 1.0 with no wind alignment mean a surface feature. And re-run one scene with the wind speed perturbed by ±50% — if the reported flux moves less than the stated uncertainty, the uncertainty is understated.

Frequently Asked Questions

Can Sentinel-2 detect methane, or do I need a hyperspectral instrument?

Sentinel-2’s B11/B12 SWIR band pair can detect large plumes — typically above roughly 1–3 tonnes per hour depending on surface and wind — using multi-band or multi-temporal ratio methods. That is enough to catch major blowouts and unlit flares over spectrally uniform bright surfaces, and its five-day revisit is a real advantage for persistence testing. It is not enough for routine site-level MRV of chronic small leaks, where you need EMIT, PRISMA, EnMAP, or a commercial targeted-observation instrument. The practical architecture uses Sentinel-2 as a wide-area screening tier that triggers a tasked hyperspectral acquisition.

Why does the wind speed dominate the uncertainty, and can I fix it with better imagery?

No — and that is the point. The IME method converts mass to flux by multiplying by an effective wind speed, so relative wind error passes into relative flux error essentially one-for-one. A retrieval precise to 15% divided by a wind field uncertain to 50% yields a flux uncertain to about 52%. Better imagery improves the numerator only. The fixes are meteorological: use the highest-resolution wind product available for the overpass time, prefer on-site anemometry where a facility provides it, and record which source was used so the uncertainty can be re-derived if a better reanalysis is published later.

How should non-detections be represented in the inventory?

As rows, with the scene-specific detection limit attached. A non-detection means “no source above L kg h⁻¹ was present at overpass”, and L varies by an order of magnitude across surfaces within a single scene. Storing only detections converts a bounded statement into an unbounded one and biases the annual total low. When the inventory is aggregated, the non-detection rows let you compute a defensible upper bound on the unobserved tail instead of implicitly assuming it is zero.

What causes a plume to appear at a facility that had no activity that day?

Three causes, in rough order of frequency. First, an albedo artefact — run the persistence and reflectance-correlation tests. Second, misattribution from georegistration error or an over-generous facility buffer, which reassigns a neighbour’s plume; check the join distance and the CRS on both layers. Third, genuine transport from an off-frame source upwind, which the wind-alignment test will reveal because the plume’s upwind extent leaves the facility polygon. Only the third is a real emission, and it belongs to someone else.

Does a methane plume record need different lineage from a land-carbon record?

It needs the same lineage plus three fields that land-carbon work does not have: the retrieval algorithm and version, the wind product and its timestamp, and the detection limit. Those three are what make the record reproducible, because the same radiance cube processed with a different target spectrum or a different wind source yields a materially different flux. Everything else — source scene identifiers, CRS, transformation chain, code version — follows the standard provenance contract.

Conclusion

Methane plume detection rewards engineering discipline more than algorithmic novelty. The matched filter is well understood and widely implemented; what separates a defensible inventory from a contested one is whether the pipeline tests for albedo artefacts before believing an enhancement, whether it carries wind error into the uncertainty it reports, whether it withholds a flux instead of inventing a wind speed, and whether it emits non-detections with their limits attached. Build those four behaviours into the stage and the retrieval becomes an evidence generator rather than a source of numbers that cannot survive review. For the quantification mathematics in full, work through quantifying methane plume emission rates in Python; for the false-positive problem in depth, see troubleshooting false methane detections over bright surfaces.