Transport backends

The Tier II and Tier III transports behind CoupledForward. A TransportBackend gives each instrument's column per unit emission through the shared ColumnObservation: LagrangianSource from met-ensemble backward footprints, EulerianSource from the les_fvm finite-volume run at the overpass time. PlumeSource is the Tier I backend. Tier IV.

backends

Tier II and Tier III transport backends (implementation in plumax.coupled._src.backends).

EulerianSource dataclass

EulerianSource(forward: EulerianForward4DVar, background: float = 0.0)

Tier III transport: the les_fvm finite-volume run through the shared column.

r_i = A_i Σ_k Δz c_k(x_i, t_inst; Q = 1),

the field of a unit constant emission at the overpass time t_inst, read in the T-cell that holds each receptor. The forward is linear in the emission with a zero initial field, so Q · r is the run at constant emission Q.

Attributes:

Name	Type	Description
forward	EulerianForward4DVar	The differentiable transport (les_fvm.build_forward); its save times must include each instrument's overpass time.
background	float	Additive background column c_bg.

unit_history cached property

unit_history: Float[Array, 'time level y x']

Interior field history of a unit constant emission, one solve.

plume_grid property

plume_grid: PlumeGrid3D

The analysis grid the transport runs on.

time_index

time_index(instrument: Instrument) -> int

Index of the save time that equals the instrument's overpass time.

Raises:

Type	Description
GridMismatchError	If no save time matches it.

Source code in src/plumax/coupled/_src/backends.py

def time_index(self, instrument: Instrument) -> int:
    """Index of the save time that equals the instrument's overpass time.

    Raises:
        GridMismatchError: If no save time matches it.
    """
    save_times = np.asarray(self.forward.save_times, dtype=float)
    matches = np.flatnonzero(np.isclose(save_times, instrument.observation_time))
    if matches.size == 0:
        raise GridMismatchError(
            f"EulerianSource: overpass time {instrument.observation_time} s of "
            f"{instrument.name!r} is not one of the forward's save times; "
            "build the forward with it in `save_times`."
        )
    return int(matches[0])

column_response

column_response(instrument: Instrument) -> Float[Array, ' obs']

Column per unit emission at the instrument's receptors.

Parameters:

Name	Type	Description	Default
instrument	Instrument	The observing instrument.	required

Returns:

Type	Description
Float[Array, ' obs']	The response r_inst.

Raises:

Type	Description
GridMismatchError	If the overpass time is not a save time.
DomainError	If a receptor lies outside the grid.

Source code in src/plumax/coupled/_src/backends.py

def column_response(self, instrument: Instrument) -> Float[Array, " obs"]:
    """Column per unit emission at the instrument's receptors.

    Args:
        instrument: The observing instrument.

    Returns:
        The response r_inst.

    Raises:
        GridMismatchError: If the overpass time is not a save time.
        DomainError: If a receptor lies outside the grid.
    """
    index = self.time_index(instrument)
    observation = ColumnObservation.for_grid(self.plume_grid, instrument=instrument)
    return observation.from_field(self.unit_history[index])

LagrangianSource dataclass

LagrangianSource(emission_pattern: Float[Array, ' cell'] | Float[Array, 'cell param'], winds: tuple[GriddedWind | Callable[[Array], Array], ...], turbulence: TurbulenceModel, domain_x: tuple[float, float, int], domain_y: tuple[float, float, int], column_z: Float[Array, ' level'], n_particles: int = 1000, t_back: float = 600.0, dt: float = 1.0, pbl_height: float = 1000.0, pbl_fraction: float = 0.5, air_density: float = 1.2, seed: int = 0, background: float = 0.0)

Tier II transport: met-ensemble backward footprints through the shared column.

c_k(x_i) = ρ_air F̄(x_i, z_k) q, F̄ = mean_e F_e, r_i = A_i Σ_k h_k c_k(x_i),

with F_e the footprint of member e at receptor (x_i, z_k), q the surface flux of a unit source and h_k the trapezoid weights over column_z.

Attributes:

Name	Type	Description
emission_pattern	Float[Array, ' cell'] | Float[Array, 'cell param']	Surface flux per unit emission q [m⁻²] on the footprint grid, cells flattened row-major over (x, y). A second axis holds one pattern per emission parameter.
winds	tuple[GriddedWind | Callable[[Array], Array], ...]	One mean wind per met-ensemble member (a GriddedWind or a t -> (u, v, w) callable). One member is a deterministic run.
turbulence	TurbulenceModel	Turbulence model shared by the members.
domain_x	tuple[float, float, int]	(start, stop, n_cells) of the footprint grid along x.
domain_y	tuple[float, float, int]	(start, stop, n_cells) of the footprint grid along y.
column_z	Float[Array, ' level']	Receptor levels z_k [m] of the column, strictly increasing.
n_particles	int	Particles per receptor level and member.
t_back	float	Backward integration time [s].
dt	float	Particle time step [s].
pbl_height	float	Boundary-layer height [m].
pbl_fraction	float	Fraction of the boundary layer that sees the surface flux.
air_density	float	Air density ρ_air [kg m⁻³].
seed	int	Base PRNG seed of the particle draws.
background	float	Additive background column c_bg.

point classmethod

point(location: tuple[float, float], *, winds: Sequence[GriddedWind | Callable[[Array], Array]], turbulence: TurbulenceModel, domain_x: tuple[float, float, int], domain_y: tuple[float, float, int], column_z: Float[Array, ' level'], **settings: float | int) -> LagrangianSource

A point source: unit emission spread over the cell that holds it.

q_s = 1 / A_cell in the containing cell s, zero elsewhere.

Parameters:

Name	Type	Description	Default
location	tuple[float, float]	Source (x, y) [m].	required
winds	Sequence[GriddedWind | Callable[[Array], Array]]	One mean wind per met member.	required
turbulence	TurbulenceModel	Turbulence model.	required
domain_x	tuple[float, float, int]	(start, stop, n_cells) along x.	required
domain_y	tuple[float, float, int]	(start, stop, n_cells) along y.	required
column_z	Float[Array, ' level']	Receptor levels z_k [m].	required
**settings	float | int	Any other field of :class:LagrangianSource.	{}

Returns:

Type	Description
LagrangianSource	The source, with a one-cell emission pattern.

Raises:

Type	Description
DomainError	If the location lies outside the footprint grid.

Source code in src/plumax/coupled/_src/backends.py

@classmethod
def point(
    cls,
    location: tuple[float, float],
    *,
    winds: Sequence[GriddedWind | Callable[[jax.Array], jax.Array]],
    turbulence: TurbulenceModel,
    domain_x: tuple[float, float, int],
    domain_y: tuple[float, float, int],
    column_z: Float[Array, " level"],
    **settings: float | int,
) -> LagrangianSource:
    """A point source: unit emission spread over the cell that holds it.

    q_s = 1 / A_cell in the containing cell s, zero elsewhere.

    Args:
        location: Source ``(x, y)`` [m].
        winds: One mean wind per met member.
        turbulence: Turbulence model.
        domain_x: ``(start, stop, n_cells)`` along x.
        domain_y: ``(start, stop, n_cells)`` along y.
        column_z: Receptor levels z_k [m].
        **settings: Any other field of :class:`LagrangianSource`.

    Returns:
        The source, with a one-cell emission pattern.

    Raises:
        DomainError: If the location lies outside the footprint grid.
    """
    x_edges = np.linspace(domain_x[0], domain_x[1], int(domain_x[2]) + 1)
    y_edges = np.linspace(domain_y[0], domain_y[1], int(domain_y[2]) + 1)
    ix = int(np.searchsorted(x_edges, location[0], side="right")) - 1
    iy = int(np.searchsorted(y_edges, location[1], side="right")) - 1
    if not (0 <= ix < x_edges.size - 1 and 0 <= iy < y_edges.size - 1):
        raise DomainError(
            f"LagrangianSource.point: {location} lies outside the footprint grid."
        )
    # A_cell = Δx Δy
    cell_area = (x_edges[1] - x_edges[0]) * (y_edges[1] - y_edges[0])
    pattern = np.zeros((x_edges.size - 1, y_edges.size - 1))
    # q_s = 1 / A_cell in the containing cell
    pattern[ix, iy] = 1.0 / cell_area
    return cls(
        emission_pattern=jnp.asarray(pattern.reshape(-1)),
        winds=tuple(winds),
        turbulence=turbulence,
        domain_x=domain_x,
        domain_y=domain_y,
        column_z=jnp.asarray(column_z),
        **settings,
    )

footprints

footprints(instrument: Instrument) -> Float[Array, 'member point cell']

Member footprints F_e at every receptor level, level fastest.

Parameters:

Name	Type	Description	Default
instrument	Instrument	The observing instrument; its overpass time is the receptor time of the backward runs.	required

Returns:

Type	Description
Float[Array, 'member point cell']	One row per (receptor, level) pair, units s m² kg⁻¹.

Source code in src/plumax/coupled/_src/backends.py

def footprints(self, instrument: Instrument) -> Float[Array, "member point cell"]:
    """Member footprints F_e at every receptor level, level fastest.

    Args:
        instrument: The observing instrument; its overpass time is the
            receptor time of the backward runs.

    Returns:
        One row per (receptor, level) pair, units s m² kg⁻¹.
    """
    levels = np.asarray(self.column_z, dtype=float)
    receptors = np.asarray(instrument.receptors, dtype=float)
    # every receptor (x_i, y_i) at every level z_k, level fastest
    horizontal = np.repeat(receptors, levels.size, axis=0)
    heights = np.tile(levels, receptors.shape[0])
    points = np.column_stack([horizontal, heights])
    return ensemble_footprints(
        list(self.winds),
        jnp.asarray(points),
        self.turbulence,
        self.domain_x,
        self.domain_y,
        n_particles=self.n_particles,
        t_back=self.t_back,
        dt=self.dt,
        pbl_height=self.pbl_height,
        pbl_fraction=self.pbl_fraction,
        air_density=self.air_density,
        receptor_time=instrument.observation_time,
        seed=self.seed,
    )

column_response

column_response(instrument: Instrument) -> Float[Array, ' obs'] | Float[Array, 'obs param']

Column per unit emission at the instrument's receptors.

Algorithm

1. Backward footprints from every receptor level, per member.
2. Average over the members: F̄ = mean_e F_e.
3. Concentration of a unit source: c = ρ_air F̄ q.
4. Column and averaging kernel: r_i = A_i Σ_k h_k c_k(x_i).

Parameters:

Name	Type	Description	Default
instrument	Instrument	The observing instrument.	required

Returns:

Type	Description
Float[Array, ' obs'] | Float[Array, 'obs param']	The response r_inst, one column per emission pattern when the
Float[Array, ' obs'] | Float[Array, 'obs param']	pattern has a parameter axis.

Source code in src/plumax/coupled/_src/backends.py

def column_response(
    self, instrument: Instrument
) -> Float[Array, " obs"] | Float[Array, "obs param"]:
    """Column per unit emission at the instrument's receptors.

    Algorithm:
        1. Backward footprints from every receptor level, per member.
        2. Average over the members: F̄ = mean_e F_e.
        3. Concentration of a unit source: c = ρ_air F̄ q.
        4. Column and averaging kernel: r_i = A_i Σ_k h_k c_k(x_i).

    Args:
        instrument: The observing instrument.

    Returns:
        The response r_inst, one column per emission pattern when the
        pattern has a parameter axis.
    """
    n_levels = int(np.shape(self.column_z)[0])
    footprints = self.footprints(instrument)
    # F̄ = mean_e F_e
    mean_footprint = einx.mean("member point cell -> point cell", footprints)
    pattern = jnp.asarray(self.emission_pattern, dtype=mean_footprint.dtype)
    observation = ColumnObservation.trapezoid(
        jnp.asarray(self.column_z, dtype=mean_footprint.dtype),
        averaging_kernel=instrument.ak,
    )
    if pattern.ndim == 1:
        # F̄ q
        mixing_ratio = einx.dot(
            "point cell, cell -> point", mean_footprint, pattern
        )
        # c = ρ_air F̄ q
        concentration = self.air_density * mixing_ratio
        profiles = einx.rearrange(
            "(obs level) -> obs level", concentration, level=n_levels
        )
        return observation.column(profiles)
    # F̄ q_p, one per pattern
    mixing_ratio = einx.dot(
        "point cell, cell param -> point param", mean_footprint, pattern
    )
    # c_p = ρ_air F̄ q_p
    concentration = self.air_density * mixing_ratio
    profiles = einx.rearrange(
        "(obs level) param -> param obs level", concentration, level=n_levels
    )
    columns = jax.vmap(observation.column)(profiles)
    return einx.rearrange("param obs -> obs param", columns)

TransportBackend

Bases: Protocol

A transport tier as the coupled forward sees it.

r_inst = A_inst col_z(c(·; Q = 1)), the column per unit emission at each of the instrument's receptors.

Attributes:

Name	Type	Description
background	float	Additive background column c_bg, in the units of the observation.

column_response

column_response(instrument: Instrument) -> Float[Array, ' obs'] | Float[Array, 'obs param']

Column per unit emission at the instrument's receptors.

Parameters:

Name	Type	Description	Default
instrument	Instrument	The observing instrument (receptors and kernel).	required

Returns:

Type	Description
Float[Array, ' obs'] | Float[Array, 'obs param']	The response r_inst: one value per receptor for a scalar emission
Float[Array, ' obs'] | Float[Array, 'obs param']	rate, or one column per emission parameter.

Source code in src/plumax/coupled/_src/backends.py

def column_response(
    self, instrument: Instrument
) -> Float[Array, " obs"] | Float[Array, "obs param"]:
    """Column per unit emission at the instrument's receptors.

    Args:
        instrument: The observing instrument (receptors and kernel).

    Returns:
        The response r_inst: one value per receptor for a scalar emission
        rate, or one column per emission parameter.
    """
    ...