Plume vetting

Plume-candidate vetting on arrays (L1): which of a candidate's pixels are in (the seeds) and which are out (the background pool), and how this scene's clutter scores in the candidate's shape (the null), on radtrax's plume-agnostic absorption-shape test (radtrax.retrieval.absorption: the pairing, the ratio spectrum, the fit and D_norm; see radtrax's absorption-shape test page). Theory: plume vetting; a worked scene: the plume vetting notebook.

Errors are radtrax's public exception classes (ParameterError, DegenerateInputError, ShapeError), each also a ValueError.

vetting

Plume-candidate vetting: is a matched-filter candidate the gas, or clutter?

Candidate-level vetting (Xiang et al. 2024; Růžička et al. 2026) as L1 functions on arrays: no carrier, no geotoolz. A candidate is a boolean mask on a radiance cube and its matched-filter map; this package chooses which of its pixels count as in (the seeds) and which as out (the background pool), and runs radtrax's plume-agnostic absorption-shape test on them (radtrax.retrieval.absorption: the pairing, the ratio spectrum, the fit to the gas's transmittance and D_norm). Theory: the plume-vetting page.

• seeds: select_seeds and background_pool, a candidate's in-plume seeds and the background pool they are paired from (#154).
• null: null_masks and null_pvalue, rotated and translated copies of a candidate's mask off every candidate, and the candidate's p-value among their scores (#155).
• vet: vet_candidate and VettingResult, one candidate vetted in one call: seeds, pool, radtrax's pairs, ratio, fit and D_norm (#156).

vet — one candidate

Plume-candidate vetting: is a matched-filter candidate the gas, or clutter?

Candidate-level vetting (Xiang et al. 2024; Růžička et al. 2026) as L1 functions on arrays: no carrier, no geotoolz. A candidate is a boolean mask on a radiance cube and its matched-filter map; this package chooses which of its pixels count as in (the seeds) and which as out (the background pool), and runs radtrax's plume-agnostic absorption-shape test on them (radtrax.retrieval.absorption: the pairing, the ratio spectrum, the fit to the gas's transmittance and D_norm). Theory: the plume-vetting page.

• seeds: select_seeds and background_pool, a candidate's in-plume seeds and the background pool they are paired from (#154).
• null: null_masks and null_pvalue, rotated and translated copies of a candidate's mask off every candidate, and the candidate's p-value among their scores (#155).
• vet: vet_candidate and VettingResult, one candidate vetted in one call: seeds, pool, radtrax's pairs, ratio, fit and D_norm (#156).

VettingResult

Bases: Module

One candidate's scores and the spectra behind them.

Attributes:

Name	Type	Description
delta	Float[Array, '']	ΔX̂, the fitted column enhancement [mol/m²], ().
d_norm	Float[Array, '']	D_norm of the fit (dimensionless), ().
n_pairs	int	The number of pairs the ratio averages (static).
ratio	Float[Array, ' Bf']	R, the ratio spectrum on the window's fitting bands, (Bf,); NaN in a band with no usable pair.
fit	TransmittanceFit	radtrax's TransmittanceFit of ratio (radtrax.retrieval.absorption): its continuum P, transmittance E and residual fit.residual_ss at ΔX̂.
residual_ss_no_gas	Float[Array, '']	The residual sum of squares of the continuum alone (E = 1, ΔX = 0) on the same bands and degree, (); with fit.residual_ss, the likelihood ratio SSE(0) / SSE(ΔX̂) of #158's third option.

vet_candidate

vet_candidate(cube: Float[ArrayLike, 'H W B'], mf: Float[ArrayLike, 'H W'], mask: Bool[ArrayLike, 'H W'], wavelength_nm: Float[ArrayLike, ' B'], table: Float[ArrayLike, 'D Bf'], delta_grid: Float[ArrayLike, ' D'], window: AbsorptionWindow, *, exclude: Bool[ArrayLike, ' C H W'] | None = None, k: int = 40, keep: float = 0.5, degree: int = 2, box: int = 200) -> VettingResult

Vet one candidate mask: seeds, pool, pairs, ratio, fit and D_norm.

With F the window's fitting bands and O its pairing bands,

I = select_seeds(MF, M), P = background_pool(MF, I) below the window's mf_background_mol_m2;

(i_n, o_n) = pair_pixels(L, I, P, O), the best keep of a one-to-one assignment of the two sets as flat pixel indices;

R(λ) = ∑_n L_in,n(λ) / ∑_n L_out,n(λ) for λ ∈ F;

(ΔX̂, ĉ) = argmin ∑_λ (R(λ) − P(λ; c) E(λ; ΔX))² over the grid, and D_norm of that fit.

Runs eagerly, not under jit (pairing is NumPy and SciPy), and is not differentiable as a whole (seeds, pairs and ΔX̂ are discrete); the pieces that are differentiable say so.

Parameters:

Name	Type	Description	Default
cube	Float[ArrayLike, 'H W B']	Radiance, (H, W, B).	required
mf	Float[ArrayLike, 'H W']	Matched-filter column ΔX [mol/m²], (H, W), in the units of window.mf_background_mol_m2.	required
mask	Bool[ArrayLike, 'H W']	The candidate M, (H, W).	required
wavelength_nm	Float[ArrayLike, ' B']	The cube's band centres [nm], (B,).	required
table	Float[ArrayLike, 'D Bf']	E(λ; ΔX_d) on the window's fitting bands, (D, Bf) with Bf = window.fit_bands(wavelength_nm).sum(), from radtrax's transmittance_table.	required
delta_grid	Float[ArrayLike, ' D']	The columns ΔX_d of the table's rows [mol/m²], (D,), starting at 0.	required
window	AbsorptionWindow	The gas's AbsorptionWindow (radtrax.retrieval.absorption).	required
exclude	Bool[ArrayLike, ' C H W'] | None	The scene's C candidate masks, (C, H, W), kept out of the background pool; pass every candidate of the scene.	None
k	int	Number of seeds before growth; the reference code's 40.	40
keep	float	Fraction of the pairs kept, in (0, 1].	0.5
degree	int	The continuum's Legendre degree.	2
box	int	Margin, in pixels, around the seeds' bounding box that the pool is drawn from.	200

Returns:

Type	Description
VettingResult	The candidate's VettingResult.

Raises:

Type	Description
ShapeError	If wavelength_nm is not (B,) with the cube's B, or table does not have Bf columns; and whatever the steps raise (select_seeds, background_pool, pair_pixels, fit_transmittance) on their own inputs.

Examples:

>>> import numpy as np
>>> from radtrax.retrieval.absorption import AbsorptionWindow
>>> from plumax.vetting import vet_candidate
>>> wl = np.array([1000.0, 1100.0, 1200.0, 2300.0, 2310.0, 2320.0])
>>> eps = np.array([2.0, 5.0, 3.0])  # m²/mol on the fitting bands
>>> grid = np.linspace(0.0, 0.2, 21)  # mol/m²
>>> table = np.exp(-grid[:, None] * eps)  # (21, 3)
>>> mask = np.zeros((6, 6), bool)
>>> mask[2:4, 2:4] = True
>>> mf = np.where(mask, 0.1, 0.0)
>>> e = np.where(mask[..., None], np.exp(-0.1 * eps), 1.0)  # (6, 6, 3)
>>> cube = np.concatenate([np.ones((6, 6, 3)), e], axis=-1)
>>> window = AbsorptionWindow(
...     "CH4", ((2290.0, 2330.0),), ((2290.0, 2330.0),), 0.01
... )
>>> r = vet_candidate(cube, mf, mask, wl, table, grid, window, degree=0)
>>> print(f"{r.delta:.2f} {r.d_norm:.2f} {r.n_pairs}")
0.10 0.00 2

Source code in src/plumax/vetting/vet.py

def vet_candidate(
    cube: Float[ArrayLike, "H W B"],
    mf: Float[ArrayLike, "H W"],
    mask: Bool[ArrayLike, "H W"],
    wavelength_nm: Float[ArrayLike, " B"],
    table: Float[ArrayLike, "D Bf"],
    delta_grid: Float[ArrayLike, " D"],
    window: AbsorptionWindow,
    *,
    exclude: Bool[ArrayLike, " C H W"] | None = None,
    k: int = 40,
    keep: float = 0.5,
    degree: int = 2,
    box: int = 200,
) -> VettingResult:
    """Vet one candidate mask: seeds, pool, pairs, ratio, fit and ``D_norm``.

    With ``F`` the window's fitting bands and ``O`` its pairing bands,

    ``I = select_seeds(MF, M)``, ``P = background_pool(MF, I)`` below the
    window's ``mf_background_mol_m2``;

    ``(i_n, o_n) = pair_pixels(L, I, P, O)``, the best ``keep`` of a
    one-to-one assignment of the two sets as flat pixel indices;

    ``R(λ) = ∑_n L_in,n(λ) / ∑_n L_out,n(λ)`` for ``λ ∈ F``;

    ``(ΔX̂, ĉ) = argmin ∑_λ (R(λ) − P(λ; c) E(λ; ΔX))²`` over the grid, and
    ``D_norm`` of that fit.

    Runs eagerly, not under ``jit`` (pairing is NumPy and SciPy), and is not
    differentiable as a whole (seeds, pairs and ``ΔX̂`` are discrete); the
    pieces that are differentiable say so.

    Args:
        cube: Radiance, ``(H, W, B)``.
        mf: Matched-filter column ``ΔX`` [mol/m²], ``(H, W)``, in the units
            of ``window.mf_background_mol_m2``.
        mask: The candidate ``M``, ``(H, W)``.
        wavelength_nm: The cube's band centres [nm], ``(B,)``.
        table: ``E(λ; ΔX_d)`` on the window's fitting bands, ``(D, Bf)``
            with ``Bf = window.fit_bands(wavelength_nm).sum()``, from
            radtrax's ``transmittance_table``.
        delta_grid: The columns ``ΔX_d`` of the table's rows [mol/m²],
            ``(D,)``, starting at 0.
        window: The gas's ``AbsorptionWindow``
            (``radtrax.retrieval.absorption``).
        exclude: The scene's ``C`` candidate masks, ``(C, H, W)``, kept out
            of the background pool; pass every candidate of the scene.
        k: Number of seeds before growth; the reference code's 40.
        keep: Fraction of the pairs kept, in (0, 1].
        degree: The continuum's Legendre degree.
        box: Margin, in pixels, around the seeds' bounding box that the pool
            is drawn from.

    Returns:
        The candidate's [`VettingResult`][plumax.vetting.VettingResult].

    Raises:
        ShapeError: If ``wavelength_nm`` is not ``(B,)`` with the cube's
            ``B``, or ``table`` does not have ``Bf`` columns; and whatever
            the steps raise (`select_seeds`, `background_pool`,
            `pair_pixels`, `fit_transmittance`) on their own inputs.

    Examples:
        >>> import numpy as np
        >>> from radtrax.retrieval.absorption import AbsorptionWindow
        >>> from plumax.vetting import vet_candidate
        >>> wl = np.array([1000.0, 1100.0, 1200.0, 2300.0, 2310.0, 2320.0])
        >>> eps = np.array([2.0, 5.0, 3.0])  # m²/mol on the fitting bands
        >>> grid = np.linspace(0.0, 0.2, 21)  # mol/m²
        >>> table = np.exp(-grid[:, None] * eps)  # (21, 3)
        >>> mask = np.zeros((6, 6), bool)
        >>> mask[2:4, 2:4] = True
        >>> mf = np.where(mask, 0.1, 0.0)
        >>> e = np.where(mask[..., None], np.exp(-0.1 * eps), 1.0)  # (6, 6, 3)
        >>> cube = np.concatenate([np.ones((6, 6, 3)), e], axis=-1)
        >>> window = AbsorptionWindow(
        ...     "CH4", ((2290.0, 2330.0),), ((2290.0, 2330.0),), 0.01
        ... )
        >>> r = vet_candidate(cube, mf, mask, wl, table, grid, window, degree=0)
        >>> print(f"{r.delta:.2f} {r.d_norm:.2f} {r.n_pairs}")
        0.10 0.00 2
    """
    wavelength_nm = np.asarray(wavelength_nm)
    cube = jnp.asarray(cube)
    table = jnp.asarray(table)
    if cube.ndim != 3 or wavelength_nm.shape != cube.shape[2:]:
        raise ShapeError(
            f"vet_candidate: wavelength_nm has shape {wavelength_nm.shape}, the "
            f"cube {cube.shape}; need (B,) with the cube's (H, W, B)"
        )
    fit_b = window.fit_bands(wavelength_nm)  # (B,)
    n_fit = int(fit_b.sum())
    if table.ndim != 2 or table.shape[1] != n_fit:
        raise ShapeError(
            f"vet_candidate: table must be (D, {n_fit}), one column per band of "
            f"the window's {n_fit} fitting bands, got {table.shape}"
        )

    seeds = select_seeds(mf, mask, k=k)  # (H, W)
    pool = background_pool(
        mf, seeds, threshold=window.mf_background_mol_m2, exclude=exclude, box=box
    )  # (H, W)
    i_in, i_out = pair_pixels(
        cube,
        np.flatnonzero(seeds),
        np.flatnonzero(pool),
        window.pairing_bands(wavelength_nm),
        keep=keep,
    )  # (N,), (N,)
    ratio = ratio_spectrum(cube, i_in, i_out)[fit_b]  # (Bf,)
    wl_fit = wavelength_nm[fit_b]  # (Bf,)
    fit = fit_transmittance(ratio, table, delta_grid, wl_fit, degree=degree)
    no_gas = fit_transmittance(
        ratio,
        jnp.ones_like(table[:1]),
        jnp.zeros(1, table.dtype),
        wl_fit,
        degree=degree,
    )
    return VettingResult(
        delta=fit.delta,
        d_norm=d_norm(ratio, fit),
        n_pairs=int(i_in.size),
        ratio=ratio,
        fit=fit,
        residual_ss_no_gas=no_gas.residual_ss,
    )

seeds — the seeds and the background pool

Plume-candidate vetting: is a matched-filter candidate the gas, or clutter?

Candidate-level vetting (Xiang et al. 2024; Růžička et al. 2026) as L1 functions on arrays: no carrier, no geotoolz. A candidate is a boolean mask on a radiance cube and its matched-filter map; this package chooses which of its pixels count as in (the seeds) and which as out (the background pool), and runs radtrax's plume-agnostic absorption-shape test on them (radtrax.retrieval.absorption: the pairing, the ratio spectrum, the fit to the gas's transmittance and D_norm). Theory: the plume-vetting page.

• seeds: select_seeds and background_pool, a candidate's in-plume seeds and the background pool they are paired from (#154).
• null: null_masks and null_pvalue, rotated and translated copies of a candidate's mask off every candidate, and the candidate's p-value among their scores (#155).
• vet: vet_candidate and VettingResult, one candidate vetted in one call: seeds, pool, radtrax's pairs, ratio, fit and D_norm (#156).

select_seeds

select_seeds(mf: Float[ArrayLike, 'H W'], mask: Bool[ArrayLike, 'H W'], *, k: int = 40, cap_percentile: Float[ArrayLike, ''] = 99.0, grow: int = 1) -> Bool[Array, 'H W']

The k strongest matched-filter pixels of a candidate, grown inside it.

With M the mask's pixels of finite score and q the cap_percentile-th percentile of the score over M,

I = grow(top_k{p ∈ M : MF(p) ≤ q}) ∩ M,

where each growth step is one 3 × 3 dilation intersected with M (the reference code's Gaussian filter with σ = 1 and threshold 0.05 is the same dilation). With grow=0, |I| = min(k, |{p ∈ M : MF(p) ≤ q}|) exactly: ties at the k-th score are broken by the lower flat index (jax.lax.top_k).

Parameters:

Name	Type	Description	Default
mf	Float[ArrayLike, 'H W']	Matched-filter score, e.g. ΔX in mol/m², (H, W); a non-finite pixel is never a seed.	required
mask	Bool[ArrayLike, 'H W']	The candidate, (H, W).	required
k	int	Number of seeds before growth (static); the reference code's 40.	40
cap_percentile	Float[ArrayLike, '']	Scores above this percentile of the mask's are not seeds, in (0, 100]; the reference code's 99.	99.0
grow	int	Growth steps inside the mask (static); 0 keeps the top k.	1

Returns:

Type	Description
Bool[Array, 'H W']	The seeds, (H, W), a subset of mask.

Raises:

Type	Description
ParameterError	If k < 1, cap_percentile is outside (0, 100] or grow < 0. k and grow are static under jit; cap_percentile may be traced, and then its range check runs with eqx.error_if, as the empty-mask check below.
ShapeError	If mf is not 2-D or mask has another shape.
DegenerateInputError	If the mask holds no pixel with a finite score. Under jit this check runs with eqx.error_if and raises a RuntimeError instead (equinox.EquinoxRuntimeError under eqx.filter_jit, jax.errors.JaxRuntimeError under jax.jit).

Examples:

>>> import jax.numpy as jnp
>>> from plumax.vetting import select_seeds
>>> mf = jnp.arange(25.0).reshape(5, 5)
>>> mask = jnp.zeros((5, 5), bool).at[1:4, 1:4].set(True)
>>> print(
...     select_seeds(mf, mask, k=2, cap_percentile=90.0, grow=0).astype(int)
... )
[[0 0 0 0 0]
 [0 0 0 0 0]
 [0 0 0 0 0]
 [0 1 1 0 0]
 [0 0 0 0 0]]
>>> print(select_seeds(mf, mask, k=2, cap_percentile=90.0).astype(int))
[[0 0 0 0 0]
 [0 0 0 0 0]
 [0 1 1 1 0]
 [0 1 1 1 0]
 [0 0 0 0 0]]

Source code in src/plumax/vetting/seeds.py

def select_seeds(
    mf: Float[ArrayLike, "H W"],
    mask: Bool[ArrayLike, "H W"],
    *,
    k: int = 40,
    cap_percentile: Float[ArrayLike, ""] = 99.0,
    grow: int = 1,
) -> Bool[Array, "H W"]:
    """The ``k`` strongest matched-filter pixels of a candidate, grown inside it.

    With ``M`` the mask's pixels of finite score and ``q`` the
    ``cap_percentile``-th percentile of the score over ``M``,

    ``I = grow(top_k{p ∈ M : MF(p) ≤ q}) ∩ M``,

    where each growth step is one 3 × 3 dilation intersected with ``M`` (the
    reference code's Gaussian filter with σ = 1 and threshold 0.05 is the
    same dilation). With ``grow=0``, ``|I| = min(k, |{p ∈ M : MF(p) ≤ q}|)``
    exactly: ties at the ``k``-th score are broken by the lower flat index
    (``jax.lax.top_k``).

    Args:
        mf: Matched-filter score, e.g. ``ΔX`` in mol/m², ``(H, W)``; a
            non-finite pixel is never a seed.
        mask: The candidate, ``(H, W)``.
        k: Number of seeds before growth (static); the reference code's 40.
        cap_percentile: Scores above this percentile of the mask's are not
            seeds, in (0, 100]; the reference code's 99.
        grow: Growth steps inside the mask (static); ``0`` keeps the top
            ``k``.

    Returns:
        The seeds, ``(H, W)``, a subset of ``mask``.

    Raises:
        ParameterError: If ``k < 1``, ``cap_percentile`` is outside
            (0, 100] or ``grow < 0``. ``k`` and ``grow`` are static under
            ``jit``; ``cap_percentile`` may be traced, and then its range
            check runs with ``eqx.error_if``, as the empty-mask check below.
        ShapeError: If ``mf`` is not 2-D or ``mask`` has another shape.
        DegenerateInputError: If the mask holds no pixel with a finite
            score. Under ``jit`` this check runs with ``eqx.error_if`` and
            raises a ``RuntimeError`` instead (``equinox.EquinoxRuntimeError``
            under ``eqx.filter_jit``, ``jax.errors.JaxRuntimeError`` under
            ``jax.jit``).

    Examples:
        >>> import jax.numpy as jnp
        >>> from plumax.vetting import select_seeds
        >>> mf = jnp.arange(25.0).reshape(5, 5)
        >>> mask = jnp.zeros((5, 5), bool).at[1:4, 1:4].set(True)
        >>> print(
        ...     select_seeds(mf, mask, k=2, cap_percentile=90.0, grow=0).astype(int)
        ... )
        [[0 0 0 0 0]
         [0 0 0 0 0]
         [0 0 0 0 0]
         [0 1 1 0 0]
         [0 0 0 0 0]]
        >>> print(select_seeds(mf, mask, k=2, cap_percentile=90.0).astype(int))
        [[0 0 0 0 0]
         [0 0 0 0 0]
         [0 1 1 1 0]
         [0 1 1 1 0]
         [0 0 0 0 0]]
    """
    if k < 1:
        raise ParameterError(f"select_seeds: k must be >= 1, got {k}")
    cap_percentile = _check_percentile(cap_percentile)
    if grow < 0:
        raise ParameterError(f"select_seeds: grow must be >= 0, got {grow}")
    mf = jnp.asarray(mf)
    mask = jnp.asarray(mask, dtype=bool)
    _check_image("select_seeds", mf, mask, "mask")
    inside = _require_any(
        mask & jnp.isfinite(mf),
        "select_seeds: the mask holds no pixel with a finite matched-filter score",
    )
    cap = jnp.nanpercentile(jnp.where(inside, mf, jnp.nan), cap_percentile)
    score = jnp.where(inside & (mf <= cap), mf, -jnp.inf).ravel()
    _, top = jax.lax.top_k(score, min(k, score.size))
    seeds = jnp.zeros(score.shape, dtype=bool).at[top].set(True)
    seeds = (seeds & jnp.isfinite(score)).reshape(mf.shape)
    for _ in range(grow):
        seeds = dilate(seeds) & inside
    return seeds

background_pool

background_pool(mf: Float[ArrayLike, 'H W'], seeds: Bool[ArrayLike, 'H W'], *, threshold: Float[ArrayLike, ''], exclude: Bool[ArrayLike, ' C H W'] | None = None, box: int = 200) -> Bool[Array, 'H W']

Background pixels near the seeds, outside every candidate, with small scores.

O = box(I) ∖ I ∖ ⋃_c M_c ∩ {p : |MF(p)| < τ}, where box(I) is the seeds' bounding box widened by box pixels on every side (clipped to the image) and M_c the C candidate masks in exclude. The box keeps each pair under the same sun, elevation and water vapour; the threshold keeps plume pixels out of the pool. Pass every candidate of the scene, this one included, as exclude: without it only the seeds themselves are kept out.

Parameters:

Name	Type	Description	Default
mf	Float[ArrayLike, 'H W']	Matched-filter score, (H, W), in the units of threshold; a non-finite pixel is never in the pool.	required
seeds	Bool[ArrayLike, 'H W']	The in-plume set, (H, W), from select_seeds.	required
threshold	Float[ArrayLike, '']	τ: pixels with |MF| strictly below it join the pool, in mol/m² for a ΔX score (radtrax's AbsorptionWindow.mf_background_mol_m2; the reference code's 30 ppm·m).	required
exclude	Bool[ArrayLike, ' C H W'] | None	The scene's C candidate masks, (C, H, W); none of their pixels joins the pool.	None
box	int	Margin, in pixels, around the seeds' bounding box (static); the reference code's 200.	200

Returns:

Type	Description
Bool[Array, 'H W']	The pool, (H, W), disjoint from seeds and from every mask in
Bool[Array, 'H W']	exclude.

Raises:

Type	Description
ParameterError	If box < 0.
ShapeError	If mf is not 2-D, seeds has another shape, or exclude is not (C, H, W).
DegenerateInputError	If the pool is empty (no seeds, or no pixel near them passes the cuts). Under jit this check runs with eqx.error_if and raises a RuntimeError instead (equinox.EquinoxRuntimeError under eqx.filter_jit, jax.errors.JaxRuntimeError under jax.jit).

Examples:

>>> import jax.numpy as jnp
>>> from plumax.vetting import background_pool
>>> mf = jnp.zeros((6, 6)).at[0, :].set(5.0)
>>> seeds = jnp.zeros((6, 6), bool).at[2, 2].set(True)
>>> candidate = jnp.zeros((1, 6, 6), bool).at[0, 2:4, 2:4].set(True)
>>> pool = background_pool(
...     mf, seeds, threshold=1.0, exclude=candidate, box=1
... )
>>> print(pool.astype(int))
[[0 0 0 0 0 0]
 [0 1 1 1 0 0]
 [0 1 0 0 0 0]
 [0 1 0 0 0 0]
 [0 0 0 0 0 0]
 [0 0 0 0 0 0]]

Source code in src/plumax/vetting/seeds.py

def background_pool(
    mf: Float[ArrayLike, "H W"],
    seeds: Bool[ArrayLike, "H W"],
    *,
    threshold: Float[ArrayLike, ""],
    exclude: Bool[ArrayLike, " C H W"] | None = None,
    box: int = 200,
) -> Bool[Array, "H W"]:
    """Background pixels near the seeds, outside every candidate, with small scores.

    ``O = box(I) ∖ I ∖ ⋃_c M_c ∩ {p : |MF(p)| < τ}``, where ``box(I)`` is the
    seeds' bounding box widened by ``box`` pixels on every side (clipped to
    the image) and ``M_c`` the ``C`` candidate masks in ``exclude``. The box
    keeps each pair under the same sun, elevation and water vapour; the
    threshold keeps plume pixels out of the pool. Pass every candidate of
    the scene, this one included, as ``exclude``: without it only the seeds
    themselves are kept out.

    Args:
        mf: Matched-filter score, ``(H, W)``, in the units of ``threshold``;
            a non-finite pixel is never in the pool.
        seeds: The in-plume set, ``(H, W)``, from
            [`select_seeds`][plumax.vetting.select_seeds].
        threshold: ``τ``: pixels with ``|MF|`` strictly below it join the
            pool, in mol/m² for a ΔX score
            (radtrax's ``AbsorptionWindow.mf_background_mol_m2``; the
            reference code's 30 ppm·m).
        exclude: The scene's ``C`` candidate masks, ``(C, H, W)``; none of
            their pixels joins the pool.
        box: Margin, in pixels, around the seeds' bounding box (static); the
            reference code's 200.

    Returns:
        The pool, ``(H, W)``, disjoint from ``seeds`` and from every mask in
        ``exclude``.

    Raises:
        ParameterError: If ``box < 0``.
        ShapeError: If ``mf`` is not 2-D, ``seeds`` has another shape, or
            ``exclude`` is not ``(C, H, W)``.
        DegenerateInputError: If the pool is empty (no seeds, or no pixel
            near them passes the cuts). Under ``jit`` this check runs with
            ``eqx.error_if`` and raises a ``RuntimeError`` instead
            (``equinox.EquinoxRuntimeError`` under ``eqx.filter_jit``,
            ``jax.errors.JaxRuntimeError`` under ``jax.jit``).

    Examples:
        >>> import jax.numpy as jnp
        >>> from plumax.vetting import background_pool
        >>> mf = jnp.zeros((6, 6)).at[0, :].set(5.0)
        >>> seeds = jnp.zeros((6, 6), bool).at[2, 2].set(True)
        >>> candidate = jnp.zeros((1, 6, 6), bool).at[0, 2:4, 2:4].set(True)
        >>> pool = background_pool(
        ...     mf, seeds, threshold=1.0, exclude=candidate, box=1
        ... )
        >>> print(pool.astype(int))
        [[0 0 0 0 0 0]
         [0 1 1 1 0 0]
         [0 1 0 0 0 0]
         [0 1 0 0 0 0]
         [0 0 0 0 0 0]
         [0 0 0 0 0 0]]
    """
    if box < 0:
        raise ParameterError(f"background_pool: box must be >= 0, got {box}")
    mf = jnp.asarray(mf)
    seeds = jnp.asarray(seeds, dtype=bool)
    _check_image("background_pool", mf, seeds, "seeds")
    height, width = mf.shape
    rows = jnp.any(seeds, axis=1)  # (H,)
    cols = jnp.any(seeds, axis=0)  # (W,)
    r0, r1 = jnp.argmax(rows), height - 1 - jnp.argmax(rows[::-1])
    c0, c1 = jnp.argmax(cols), width - 1 - jnp.argmax(cols[::-1])
    i = jnp.arange(height)[:, None]
    j = jnp.arange(width)[None, :]
    near = (i >= r0 - box) & (i <= r1 + box) & (j >= c0 - box) & (j <= c1 + box)
    pool = near & jnp.any(rows) & ~seeds & (jnp.abs(mf) < threshold)
    if exclude is not None:
        exclude = jnp.asarray(exclude, dtype=bool)
        if exclude.ndim != 3 or exclude.shape[1:] != mf.shape:
            raise ShapeError(
                f"background_pool: exclude must be (C, {height}, {width}), "
                f"got {exclude.shape}"
            )
        pool &= ~jnp.any(exclude, axis=0)
    return _require_any(
        pool,
        "background_pool: the pool is empty; no pixel near the seeds is "
        "outside every candidate with |MF| below the threshold",
    )

null — a null from the same scene

Plume-candidate vetting: is a matched-filter candidate the gas, or clutter?

Candidate-level vetting (Xiang et al. 2024; Růžička et al. 2026) as L1 functions on arrays: no carrier, no geotoolz. A candidate is a boolean mask on a radiance cube and its matched-filter map; this package chooses which of its pixels count as in (the seeds) and which as out (the background pool), and runs radtrax's plume-agnostic absorption-shape test on them (radtrax.retrieval.absorption: the pairing, the ratio spectrum, the fit to the gas's transmittance and D_norm). Theory: the plume-vetting page.

• seeds: select_seeds and background_pool, a candidate's in-plume seeds and the background pool they are paired from (#154).
• null: null_masks and null_pvalue, rotated and translated copies of a candidate's mask off every candidate, and the candidate's p-value among their scores (#155).
• vet: vet_candidate and VettingResult, one candidate vetted in one call: seeds, pool, radtrax's pairs, ratio, fit and D_norm (#156).

null_masks

null_masks(mask: Bool[ArrayLike, 'H W'], n: int, *, key: PRNGKeyArray, exclude: Bool[ArrayLike, ' C H W'] | None = None, min_fraction: float = 0.9) -> Bool[np.ndarray, 'S H W']

n rotated and translated copies of a mask, off every candidate.

Each attempt draws an angle θ ~ U[0°, 360°) and rotates mask about the image centre by θ (nearest neighbour, the image's shape kept); a rotation that keeps fewer than min_fraction · |M| pixels is discarded. Among the integer shifts that keep the rotated mask inside the image, the free ones are those where

M_s ∩ (M ∪ ⋃_c M_c) = ∅,

with M_c the C candidate masks in exclude; one free shift is drawn uniformly, and an angle with none is discarded. So a copy never touches a candidate, this one included, and given its angle its place is uniform over the free places (found exactly, by correlating the rotated mask with the candidates, rather than by rejecting random shifts, which a mask covering a large part of the scene almost never passes). Copies may overlap one another. They are returned in the order they were drawn, so the first m of n copies are the m copies the same key gives.

Parameters:

Name	Type	Description	Default
mask	Bool[ArrayLike, 'H W']	The candidate M, (H, W).	required
n	int	Number of copies S; the reference package's 100.	required
key	PRNGKeyArray	A jax.random key, the only source of randomness: the same key gives the same copies.	required
exclude	Bool[ArrayLike, ' C H W'] | None	The scene's C candidate masks, (C, H, W); pass every candidate (mask itself is always avoided).	None
min_fraction	float	A copy keeps at least this fraction of M's pixels (rotation can push some off the image or merge them), in (0, 1].	0.9

Returns:

Type	Description
Bool[ndarray, 'S H W']	The copies, (S, H, W), a NumPy array.

Raises:

Type	Description
ParameterError	If n < 1 or min_fraction is outside (0, 1].
ShapeError	If mask is not 2-D or exclude is not (C, H, W).
DegenerateInputError	If mask is empty, or 100 · n angles place fewer than n copies (the scene is too small, or too crowded with candidates, for copies of this mask).

Examples:

>>> import jax
>>> import numpy as np
>>> from plumax.vetting import null_masks
>>> mask = np.zeros((32, 32), bool)
>>> mask[4:8, 4:12] = True
>>> copies = null_masks(mask, 5, key=jax.random.key(0), exclude=mask[None])
>>> copies.shape
(5, 32, 32)
>>> bool((copies & mask).any())
False
>>> bool((copies.sum(axis=(1, 2)) >= 0.9 * mask.sum()).all())
True

Source code in src/plumax/vetting/null.py

def null_masks(
    mask: Bool[ArrayLike, "H W"],
    n: int,
    *,
    key: PRNGKeyArray,
    exclude: Bool[ArrayLike, " C H W"] | None = None,
    min_fraction: float = 0.9,
) -> Bool[np.ndarray, "S H W"]:
    """``n`` rotated and translated copies of a mask, off every candidate.

    Each attempt draws an angle ``θ ~ U[0°, 360°)`` and rotates ``mask``
    about the image centre by ``θ`` (nearest neighbour, the image's shape
    kept); a rotation that keeps fewer than ``min_fraction · |M|`` pixels is
    discarded. Among the integer shifts that keep the rotated mask inside
    the image, the free ones are those where

    ``M_s ∩ (M ∪ ⋃_c M_c) = ∅``,

    with ``M_c`` the ``C`` candidate masks in ``exclude``; one free shift is
    drawn uniformly, and an angle with none is discarded. So a copy never
    touches a candidate, this one included, and given its angle its place is
    uniform over the free places (found exactly, by correlating the rotated
    mask with the candidates, rather than by rejecting random shifts, which
    a mask covering a large part of the scene almost never passes). Copies
    may overlap one another. They are returned in the order they were drawn,
    so the first ``m`` of ``n`` copies are the ``m`` copies the same key
    gives.

    Args:
        mask: The candidate ``M``, ``(H, W)``.
        n: Number of copies ``S``; the reference package's 100.
        key: A ``jax.random`` key, the only source of randomness: the same
            key gives the same copies.
        exclude: The scene's ``C`` candidate masks, ``(C, H, W)``; pass
            every candidate (``mask`` itself is always avoided).
        min_fraction: A copy keeps at least this fraction of ``M``'s pixels
            (rotation can push some off the image or merge them), in (0, 1].

    Returns:
        The copies, ``(S, H, W)``, a NumPy array.

    Raises:
        ParameterError: If ``n < 1`` or ``min_fraction`` is outside (0, 1].
        ShapeError: If ``mask`` is not 2-D or ``exclude`` is not
            ``(C, H, W)``.
        DegenerateInputError: If ``mask`` is empty, or ``100 · n`` angles
            place fewer than ``n`` copies (the scene is too small, or too
            crowded with candidates, for copies of this mask).

    Examples:
        >>> import jax
        >>> import numpy as np
        >>> from plumax.vetting import null_masks
        >>> mask = np.zeros((32, 32), bool)
        >>> mask[4:8, 4:12] = True
        >>> copies = null_masks(mask, 5, key=jax.random.key(0), exclude=mask[None])
        >>> copies.shape
        (5, 32, 32)
        >>> bool((copies & mask).any())
        False
        >>> bool((copies.sum(axis=(1, 2)) >= 0.9 * mask.sum()).all())
        True
    """
    if n < 1:
        raise ParameterError(f"null_masks: n must be >= 1, got {n}")
    if not 0.0 < min_fraction <= 1.0:
        raise ParameterError(
            f"null_masks: min_fraction must be in (0, 1], got {min_fraction}"
        )
    mask = np.asarray(mask, dtype=bool)
    if mask.ndim != 2:
        raise ShapeError(f"null_masks: mask must be (H, W), got {mask.shape}")
    height, width = mask.shape
    banned = mask.copy()
    if exclude is not None:
        exclude = np.asarray(exclude, dtype=bool)
        if exclude.ndim != 3 or exclude.shape[1:] != mask.shape:
            raise ShapeError(
                f"null_masks: exclude must be (C, {height}, {width}), got "
                f"{exclude.shape}"
            )
        banned |= np.any(exclude, axis=0)
    size = int(mask.sum())
    if size == 0:
        raise DegenerateInputError("null_masks: the mask is empty")

    from scipy.ndimage import rotate
    from scipy.signal import correlate

    needed = min_fraction * size
    image = mask.astype(np.float64)
    forbidden = banned.astype(np.float64)
    copies: list[np.ndarray] = []
    for _ in range(100 * n):
        key, k_angle, k_place = jax.random.split(key, 3)
        angle = float(jax.random.uniform(k_angle, maxval=360.0))
        rotated = rotate(image, angle, reshape=False, order=0) > 0.5  # (H, W)
        rows, cols = np.nonzero(rotated)
        if rows.size < needed:
            continue
        crop = rotated[rows.min() : rows.max() + 1, cols.min() : cols.max() + 1]
        # Candidate pixels under the crop at each top-left corner that keeps
        # it inside the scene: (H - h + 1, W - w + 1).
        overlap = correlate(forbidden, crop.astype(np.float64), mode="valid")
        free = np.flatnonzero(overlap < 0.5)
        if free.size == 0:
            continue
        pick = free[int(jax.random.randint(k_place, (), 0, free.size))]
        top, left = np.unravel_index(pick, overlap.shape)
        copy = np.zeros_like(mask)
        copy[top : top + crop.shape[0], left : left + crop.shape[1]] = crop
        copies.append(copy)
        if len(copies) == n:
            return np.stack(copies)
    raise DegenerateInputError(
        f"null_masks: {100 * n} angles placed {len(copies)} of {n} copies of a "
        f"{size}-pixel mask in a {height} × {width} scene off every candidate"
    )

null_pvalue

null_pvalue(score: float, null_scores: Float[ArrayLike, ' S']) -> float

The fraction of null copies that score at least as well as the candidate.

With D the candidate's score and D_s the S' finite scores of its null copies (non-finite ones, copies that could not be scored, are left out),

p̂ = (1/S') ∑_s 1[D_s ≤ D],

ties counting as at least as good. Lower scores are better, as for D_norm; for a score where higher is better (a likelihood ratio), pass both negated. p̂ moves in steps of 1/S' and p̂ = 0 says only that no copy did as well; the Monte Carlo p-value, never zero and valid when the candidate and the copies are exchangeable under the null (Phipson and Smyth 2010), is (1 + S'·p̂) / (1 + S').

Parameters:

Name	Type	Description	Default
score	float	The candidate's score D, e.g. its D_norm.	required
null_scores	Float[ArrayLike, ' S']	The scores of its null copies, (S,), from null_masks vetted as the candidate was.	required

Returns:

Type	Description
float	p̂, in [0, 1].

Raises:

Type	Description
ShapeError	If null_scores is not 1-D.
DegenerateInputError	If score is not finite or no null score is.

Examples:

>>> import numpy as np
>>> from plumax.vetting import null_pvalue
>>> null_pvalue(0.3, np.array([0.2, 0.3, 0.5, 0.9]))
0.5
>>> null_pvalue(0.1, np.array([0.2, np.nan, 0.5]))
0.0

Source code in src/plumax/vetting/null.py

def null_pvalue(score: float, null_scores: Float[ArrayLike, " S"]) -> float:
    """The fraction of null copies that score at least as well as the candidate.

    With ``D`` the candidate's score and ``D_s`` the ``S'`` finite scores of
    its null copies (non-finite ones, copies that could not be scored, are
    left out),

    ``p̂ = (1/S') ∑_s 1[D_s ≤ D]``,

    ties counting as at least as good. Lower scores are better, as for
    ``D_norm``; for a score where higher is better (a likelihood ratio),
    pass both negated. ``p̂`` moves in steps of ``1/S'`` and ``p̂ = 0`` says
    only that no copy did as well; the Monte Carlo p-value, never zero and
    valid when the candidate and the copies are exchangeable under the null
    (Phipson and Smyth 2010), is ``(1 + S'·p̂) / (1 + S')``.

    Args:
        score: The candidate's score ``D``, e.g. its ``D_norm``.
        null_scores: The scores of its null copies, ``(S,)``, from
            [`null_masks`][plumax.vetting.null_masks] vetted as
            the candidate was.

    Returns:
        ``p̂``, in [0, 1].

    Raises:
        ShapeError: If ``null_scores`` is not 1-D.
        DegenerateInputError: If ``score`` is not finite or no null score
            is.

    Examples:
        >>> import numpy as np
        >>> from plumax.vetting import null_pvalue
        >>> null_pvalue(0.3, np.array([0.2, 0.3, 0.5, 0.9]))
        0.5
        >>> null_pvalue(0.1, np.array([0.2, np.nan, 0.5]))
        0.0
    """
    null_scores = jnp.asarray(null_scores)
    if null_scores.ndim != 1:
        raise ShapeError(
            f"null_pvalue: null_scores must be (S,), got {null_scores.shape}"
        )
    score = float(score)
    if not np.isfinite(score):
        raise DegenerateInputError(f"null_pvalue: the score is {score}, not finite")
    finite = jnp.isfinite(null_scores)  # (S,)
    n_finite = int(jnp.sum(finite))
    if n_finite == 0:
        raise DegenerateInputError("null_pvalue: no null score is finite")
    at_least_as_good = jnp.sum(finite & (null_scores <= score))
    return float(at_least_as_good) / n_finite