Spectroscopy

Cross-section look-up tables σ(ν; T, p): per-layer optical depth for a layered atmosphere (opacity: the OpacityProvider protocol, TableOpacity, and SumOpacity, the sum of several providers on one grid), single-slab interpolation and Beer–Lambert transmittance (tables), v0's table cache and loader (table_cache), the synthetic methane spectra for docs, notebooks and tests (synthetic), the solar irradiance on the model grid from the TSIS-1 HSRS, downloaded once into the table cache, and the Earth–Sun distance (solar), and the offline HITRAN table builder (hitran, the [hitran] extra; python -m radtrax.spectroscopy.hitran CH4 H2O CO2 N2O CO builds v0's tables into the cache). Theory: HITRAN tables.

spectroscopy

Absorption cross sections per layer.

tables interpolates a cross-section look-up table σ(ν; T, p) and applies Beer–Lambert transmittance, including the differential-ratio form used for plume-enhancement retrievals. hitran builds those tables offline from HITRAN through HAPI (the [hitran] extra).

opacity turns them into per-layer optical depth for a layered atmosphere: the OpacityProvider protocol and TableOpacity, bilinear in (ln p, T) (radtrax#33), and SumOpacity, the sum of several providers on one grid (radtrax#184). Line-by-line opacity arrives in v2 of the roadmap. synthetic builds the synthetic methane spectra the docs, notebooks and tests use without HITRAN (radtrax#115). table_cache locates and loads v0's built tables (radtrax#34). solar gives the solar irradiance on the model grid from the TSIS-1 HSRS, fetched once into the same cache, and the Earth–Sun distance (radtrax#107).

opacity

Per-layer optical depth: the OpacityProvider protocol and table opacity.

The optical depth of each layer at each wavenumber, τ (L, G), is radtrax's canonical interface between the atmosphere and a solver (radtrax · 3). v0 computes it from cross-section tables σ_g(ν; p, T), interpolated bilinearly in (ln p, T) at each layer's pressure and temperature:

.. math:: σ_{g,l}(ν) = \operatorname{bilinear}(σ_g;\ \ln p_l, T_l), \qquad τ_l(ν) = \sum_g σ_{g,l}(ν)\, n_l\, x_{g,l}\, Δz_l ,

with n_l the air number density, x_{g,l} the gas's mixing ratio and Δz_l the thickness (:class:~radtrax.atmosphere.AtmosphericState). Sources of optical depth add: :class:SumOpacity holds several providers on one grid and returns Σ_k τ^{(k)} (radtrax#184).

Ported from the radtrax · 3 design notebook's TableOpacity and generalised to a dict of per-gas tables mapped with jax.tree.map, so the gas list lives in the data. The interpolation is :class:jax.scipy.interpolate.RegularGridInterpolator. Shapes: L layers, G wavenumbers, P and T table pressures and temperatures.

OpacityProvider

Bases: Module

Protocol: an atmospheric state in, per-layer optical depth out.

An abstract Equinox module (an eqx.AbstractVar field and an abc-abstract __call__), so filter_jit compiles each concrete provider and a forward model can hold any of them.

Attributes:

Name	Type	Description
nu_cm	AbstractVar[Float[Array, ' G']]	The wavenumber grid [cm⁻¹], (G,).

CrossSectionTable

Bases: Module

One gas's cross sections σ(ν; p, T) on a (ln p, T) grid.

Attributes:

Name	Type	Description
log_p_pa	Float[Array, ' P']	Natural log of the table pressures [ln Pa], strictly increasing, (P,).
t_k	Float[Array, ' T']	Table temperatures [K], strictly increasing, (T,).
sigma_cm2	Float[Array, 'P T G']	Cross sections [cm²/molecule], (P, T, G).

TableOpacity

TableOpacity(nu_cm: Float[Array, ' G'], tables: dict[str, CrossSectionTable])

Bases: OpacityProvider

Optical depth from per-gas cross-section tables.

τ_l(ν) = Σ_g σ_g(ν; p_l, T_l) n_l x_{g,l} Δz_l. The tables are a dict keyed by gas and mapped with jax.tree.map; every table gas must be in the state's vmr (gases only in the state are transparent here). Tables are leaves: to keep them out of a gradient, split them off with eqx.partition. Table opacity needs float64 (jax_enable_x64); a float32 grid warns.

Attributes:

Name	Type	Description
nu_cm	Float[Array, ' G']	Wavenumber grid [cm⁻¹], (G,), shared by every table.
tables	dict[str, CrossSectionTable]	One :class:CrossSectionTable per gas.

Examples:

>>> import jax.numpy as jnp
>>> from radtrax.atmosphere import layers_from_levels
>>> table = CrossSectionTable(
...     log_p_pa=jnp.log(jnp.array([1e4, 1e5])),
...     t_k=jnp.array([200.0, 300.0]),
...     sigma_cm2=jnp.full((2, 2, 3), 1e-21),
... )
>>> opacity = TableOpacity(
...     jnp.array([4300.0, 4300.1, 4300.2]), {"CH4": table}
... )
>>> atm = layers_from_levels(
...     jnp.array([1e5, 5e4]),
...     jnp.array([288.0, 255.0]),
...     {"CH4": jnp.full(2, 1.9e-6)},
... )
>>> opacity(atm).shape
(1, 3)

Source code in src/radtrax/spectroscopy/opacity.py

def __init__(
    self, nu_cm: Float[Array, " G"], tables: dict[str, CrossSectionTable]
) -> None:
    self.nu_cm = jnp.asarray(nu_cm)
    if self.nu_cm.ndim != 1:
        raise ShapeError(
            f"TableOpacity: nu_cm must be (G,), got {self.nu_cm.shape}."
        )
    if self.nu_cm.dtype != jnp.float64:
        warnings.warn(
            f"TableOpacity: nu_cm is {self.nu_cm.dtype}; table opacity needs "
            "float64 (float32 spaces 4300 cm⁻¹ by 4.9e-4 cm⁻¹). Enable "
            "jax_enable_x64 before building the tables.",
            stacklevel=2,
        )
    if not tables:
        raise ParameterError("TableOpacity: need at least one gas table.")
    for gas, table in tables.items():
        n_nu = jnp.shape(table.sigma_cm2)[-1]
        if n_nu != self.nu_cm.shape[0]:
            raise GridError(
                f"TableOpacity: the {gas!r} table has {n_nu} wavenumbers, "
                f"nu_cm has {self.nu_cm.shape[0]}."
            )
    self.tables = dict(tables)

from_datasets classmethod

from_datasets(datasets: dict[str, Dataset], *, var: str | dict[str, str] = 'absorption_cross_section') -> TableOpacity

Build from per-gas look-up tables in radtrax's HITRAN format.

Each dataset is what :func:radtrax.spectroscopy.hitran.build_lut_dataset writes: var on (wavenumber, temperature, pressure) with the pressure in atm, in any order (pressure and temperature are sorted). var may be a dict per gas, e.g. {"CH4": "CH4_cross_section"} for the table of :func:radtrax.spectroscopy.hitran.create_combined_lut. Every gas must share one wavenumber grid.

Raises:

Type	Description
KeyError	If a dataset has no such variable.
ValueError	If the wavenumber grids differ.

Source code in src/radtrax/spectroscopy/opacity.py

@classmethod
def from_datasets(
    cls,
    datasets: dict[str, xr.Dataset],
    *,
    var: str | dict[str, str] = "absorption_cross_section",
) -> TableOpacity:
    """Build from per-gas look-up tables in radtrax's HITRAN format.

    Each dataset is what :func:`radtrax.spectroscopy.hitran.build_lut_dataset`
    writes: ``var`` on ``(wavenumber, temperature, pressure)`` with the
    pressure in atm, in any order (pressure and temperature are sorted).
    ``var`` may be a ``dict`` per gas, e.g. ``{"CH4": "CH4_cross_section"}``
    for the table of :func:`radtrax.spectroscopy.hitran.create_combined_lut`.
    Every gas must share one wavenumber grid.

    Raises:
        KeyError: If a dataset has no such variable.
        ValueError: If the wavenumber grids differ.
    """
    nu: np.ndarray | None = None
    tables: dict[str, CrossSectionTable] = {}
    for gas, ds in datasets.items():
        name = var[gas] if isinstance(var, dict) else var
        if name not in ds:
            raise MissingKeyError(
                f"TableOpacity: {gas!r} has no {name!r} (have {list(ds.data_vars)})."
            )
        da = (
            ds[name]
            .sortby(["pressure", "temperature"])
            .transpose("pressure", "temperature", "wavenumber")
        )
        grid = np.asarray(da["wavenumber"].values, dtype=float)
        if nu is None:
            nu = grid
        elif grid.shape != nu.shape or not np.array_equal(grid, nu):
            raise GridError(
                f"TableOpacity: {gas!r} is on a different wavenumber grid."
            )
        tables[gas] = CrossSectionTable(
            log_p_pa=jnp.log(
                jnp.asarray(da["pressure"].values, dtype=float) * ATM_TO_PA
            ),
            t_k=jnp.asarray(da["temperature"].values, dtype=float),
            sigma_cm2=jnp.asarray(da.values, dtype=float),
        )
    if nu is None:
        raise ParameterError("TableOpacity: need at least one gas table.")
    return cls(jnp.asarray(nu), tables)

SumOpacity

SumOpacity(providers: Sequence[OpacityProvider])

Bases: OpacityProvider

Optical depth from several sources on one wavenumber grid.

τ_l(ν) = Σ_k τ_l^{(k)}(ν): lines, continua and collision-induced absorption are separate sources whose temperature and pressure dependences differ, and they add at the optical-depth level (Korkin et al. 2026, §1, §6). Each provider stays a leaf of the module, so a gradient flows through every one, and each stays testable alone.

Attributes:

Name	Type	Description
nu_cm	Float[Array, ' G']	Wavenumber grid [cm⁻¹], (G,), the first provider's; every provider must share it.
providers	tuple[OpacityProvider, ...]	The opacity providers whose optical depths add.

Raises:

Type	Description
ParameterError	With no providers.
GridError	If a provider's grid differs from the first's in length, or in value when the grids are concrete (under tracing only the lengths are checked).

Examples:

>>> import jax.numpy as jnp
>>> from radtrax.atmosphere import layers_from_levels
>>> def table(sigma):
...     return CrossSectionTable(
...         log_p_pa=jnp.log(jnp.array([1e4, 1e5])),
...         t_k=jnp.array([200.0, 300.0]),
...         sigma_cm2=jnp.full((2, 2, 3), sigma),
...     )
>>> nu = jnp.array([4300.0, 4300.1, 4300.2])
>>> lines = TableOpacity(nu, {"CH4": table(1e-21)})
>>> other = TableOpacity(nu, {"H2O": table(2e-23)})
>>> atm = layers_from_levels(
...     jnp.array([1e5, 5e4]),
...     jnp.array([288.0, 255.0]),
...     {"CH4": jnp.full(2, 1.9e-6), "H2O": jnp.full(2, 1e-2)},
... )
>>> total = SumOpacity((lines, other))
>>> bool(jnp.allclose(total(atm), lines(atm) + other(atm)))
True

Source code in src/radtrax/spectroscopy/opacity.py

def __init__(self, providers: Sequence[OpacityProvider]) -> None:
    providers = tuple(providers)
    if not providers:
        raise ParameterError("SumOpacity: need at least one opacity provider.")
    nu = jnp.asarray(providers[0].nu_cm)
    for k, provider in enumerate(providers[1:], start=1):
        other = jnp.asarray(provider.nu_cm)
        if jnp.shape(other) != jnp.shape(nu):
            raise GridError(
                f"SumOpacity: provider {k} has {jnp.shape(other)} wavenumbers, "
                f"provider 0 has {jnp.shape(nu)}."
            )
        try:
            same = bool(jnp.all(other == nu))
        except (
            jax.errors.TracerBoolConversionError,
            jax.errors.ConcretizationTypeError,
        ):
            continue
        if not same:
            raise GridError(
                f"SumOpacity: provider {k} is on a different wavenumber grid "
                "from provider 0."
            )
    self.nu_cm = nu
    self.providers = providers

solar

The solar irradiance at the top of the atmosphere, on the model grid.

:func:solar_irradiance gives \(E_0(\nu)\) at 1 AU in W m⁻² (cm⁻¹)⁻¹, the solar_irradiance that :class:~radtrax.solvers.BeerLambert takes, on the wavenumber grid before the band integration, so that the Fraunhofer lines are integrated by the instrument rather than averaged into the band beforehand (radtrax#107).

The one source today is the TSIS-1 Hybrid Solar Reference Spectrum, version 2 (Coddington et al. 2021, Geophys. Res. Lett. 48, e2020GL091709; LASP, file hybrid_reference_spectrum_c2022-11-30_with_unc.nc): 202 to 2730 nm in vacuum wavelength at 0.001 nm sampling, in W m⁻² nm⁻¹, at 1 AU. It is 61 MB, so it doesn't ship with radtrax: :func:fetch_solar_spectrum downloads it once into the table cache (:func:~radtrax.spectroscopy.table_cache_dir) and checks its SHA-256, and :func:solar_irradiance fetches it on first use.

Per unit wavenumber, with \(\lambda\) in nm and \(\nu = 10^7 / \lambda\) in cm⁻¹,

\[E_\nu = E_\lambda \left|\frac{d\lambda}{d\nu}\right| = E_\lambda \frac{\lambda^2}{10^7}.\]

At a date other than the mean distance, divide by \(d^2\), with \(d\) from :func:earth_sun_distance.

SolarSource dataclass

SolarSource(name: str, url: str, filename: str, sha256: str, citation: str)

A published solar reference spectrum that radtrax can download.

Attributes:

Name	Type	Description
name	str	The key :func:solar_irradiance takes ("tsis1-hsrs").
url	str	Where the file is downloaded from.
filename	str	The file's name in the cache.
sha256	str	The file's SHA-256, checked after every download.
citation	str	The reference to cite.

solar_spectrum_path

solar_spectrum_path(source: str = 'tsis1-hsrs', *, cache_dir: Path | str | None = None) -> Path

The cache file of a solar spectrum (it may not exist yet).

Parameters:

Name	Type	Description	Default
source	str	A key of :data:SOLAR_SOURCES.	'tsis1-hsrs'
cache_dir	Path | str | None	The cache directory; :func:~radtrax.spectroscopy.table_cache_dir when None.	None

Returns:

Type	Description
Path	<cache_dir>/<filename>.

Raises:

Type	Description
ParameterError	If source is unknown.

Examples:

>>> from radtrax.spectroscopy.solar import solar_spectrum_path
>>> solar_spectrum_path(cache_dir="/cache").name
'hybrid_reference_spectrum_c2022-11-30_with_unc.nc'

Source code in src/radtrax/spectroscopy/solar.py

def solar_spectrum_path(
    source: str = "tsis1-hsrs", *, cache_dir: Path | str | None = None
) -> Path:
    """The cache file of a solar spectrum (it may not exist yet).

    Args:
        source: A key of :data:`SOLAR_SOURCES`.
        cache_dir: The cache directory; :func:`~radtrax.spectroscopy.table_cache_dir`
            when ``None``.

    Returns:
        ``<cache_dir>/<filename>``.

    Raises:
        ParameterError: If ``source`` is unknown.

    Examples:
        >>> from radtrax.spectroscopy.solar import solar_spectrum_path
        >>> solar_spectrum_path(cache_dir="/cache").name
        'hybrid_reference_spectrum_c2022-11-30_with_unc.nc'
    """
    root = Path(cache_dir) if cache_dir is not None else table_cache_dir()
    return root / _source(source).filename

fetch_solar_spectrum

fetch_solar_spectrum(source: str = 'tsis1-hsrs', *, cache_dir: Path | str | None = None) -> Path

Download a solar spectrum into the cache, once, and check its SHA-256.

A file already in the cache is returned without a download.

Parameters:

Name	Type	Description	Default
source	str	A key of :data:SOLAR_SOURCES.	'tsis1-hsrs'
cache_dir	Path | str | None	The cache directory; :func:~radtrax.spectroscopy.table_cache_dir when None.	None

Returns:

Type	Description
Path	The path of the cached file.

Raises:

Type	Description
ParameterError	If source is unknown.
MissingTableError	If the download fails (offline, or the server refuses); the message names the URL and the command to retry.
TableBuildError	If the downloaded file's SHA-256 is not the registered one; nothing is written to the cache.

Source code in src/radtrax/spectroscopy/solar.py

def fetch_solar_spectrum(
    source: str = "tsis1-hsrs", *, cache_dir: Path | str | None = None
) -> Path:
    """Download a solar spectrum into the cache, once, and check its SHA-256.

    A file already in the cache is returned without a download.

    Args:
        source: A key of :data:`SOLAR_SOURCES`.
        cache_dir: The cache directory; :func:`~radtrax.spectroscopy.table_cache_dir`
            when ``None``.

    Returns:
        The path of the cached file.

    Raises:
        ParameterError: If ``source`` is unknown.
        MissingTableError: If the download fails (offline, or the server
            refuses); the message names the URL and the command to retry.
        TableBuildError: If the downloaded file's SHA-256 is not the
            registered one; nothing is written to the cache.
    """
    spec = _source(source)
    path = solar_spectrum_path(source, cache_dir=cache_dir)
    if path.is_file():
        return path
    path.parent.mkdir(parents=True, exist_ok=True)
    fd, tmp_name = tempfile.mkstemp(dir=path.parent, suffix=".part")
    tmp = Path(tmp_name)
    try:
        try:
            with os.fdopen(fd, "wb") as out, urllib.request.urlopen(spec.url) as resp:
                for chunk in iter(lambda: resp.read(1 << 20), b""):
                    out.write(chunk)
        except (urllib.error.URLError, OSError) as err:
            raise MissingTableError(
                f"Could not download the {source} solar spectrum from {spec.url} "
                f"({err}). Retry when online:\n    {_fetch_command(source)}\nor "
                f"download the file by hand to {path}."
            ) from err
        digest = _sha256(tmp)
        if digest != spec.sha256:
            raise TableBuildError(
                f"The {source} solar spectrum from {spec.url} has SHA-256 "
                f"{digest}, expected {spec.sha256}; the file was not cached."
            )
        tmp.replace(path)
    finally:
        tmp.unlink(missing_ok=True)
    return path

solar_irradiance

solar_irradiance(nu_cm: Float[ArrayLike, ' G'], *, source: str = 'tsis1-hsrs', cache_dir: Path | str | None = None, download: bool = True) -> Float[Array, ' G']

Solar irradiance \(E_0(\nu)\) at 1 AU, interpolated onto a wavenumber grid.

Linear in wavenumber between the reference's samples (0.001 nm, i.e. 0.0019 cm⁻¹ at 2.3 µm and 0.0037 cm⁻¹ at 1.6 µm, so finer than v0's 0.01 cm⁻¹ grid: the grid samples the spectrum rather than averaging it). The spectrum is read with NumPy, off the gradient path: call this once, outside jit, and pass the array to :class:~radtrax.solvers.BeerLambert.

Parameters:

Name	Type	Description	Default
nu_cm	Float[ArrayLike, ' G']	The model's wavenumbers [cm⁻¹], (G,), inside the source's range (3663–49505 cm⁻¹, i.e. 202–2730 nm, for TSIS-1 HSRS).	required
source	str	A key of :data:SOLAR_SOURCES.	'tsis1-hsrs'
cache_dir	Path | str | None	The cache directory; :func:~radtrax.spectroscopy.table_cache_dir when None.	None
download	bool	Fetch the spectrum if it is not cached (:func:fetch_solar_spectrum); when False, a missing file raises instead.	True

Returns:

Type	Description
Float[Array, ' G']	\(E_0\) [W m⁻² (cm⁻¹)⁻¹], (G,). Divide by \(d^2\)
Float[Array, ' G']	(:func:earth_sun_distance) for a given date.

Raises:

Type	Description
ParameterError	If source is unknown.
MissingTableError	If the spectrum is not cached and download is False or the download fails; the message names the fetch command.
GridError	If any wavenumber is outside the source's range.

Source code in src/radtrax/spectroscopy/solar.py

def solar_irradiance(
    nu_cm: Float[ArrayLike, " G"],
    *,
    source: str = "tsis1-hsrs",
    cache_dir: Path | str | None = None,
    download: bool = True,
) -> Float[Array, " G"]:
    """Solar irradiance $E_0(\\nu)$ at 1 AU, interpolated onto a wavenumber grid.

    Linear in wavenumber between the reference's samples (0.001 nm, i.e.
    0.0019 cm⁻¹ at 2.3 µm and 0.0037 cm⁻¹ at 1.6 µm, so finer than v0's
    0.01 cm⁻¹ grid: the grid samples the spectrum rather than averaging
    it). The spectrum is read with NumPy, off the gradient path: call this
    once, outside ``jit``, and pass the array to
    :class:`~radtrax.solvers.BeerLambert`.

    Args:
        nu_cm: The model's wavenumbers [cm⁻¹], ``(G,)``, inside the source's
            range (3663–49505 cm⁻¹, i.e. 202–2730 nm, for TSIS-1 HSRS).
        source: A key of :data:`SOLAR_SOURCES`.
        cache_dir: The cache directory; :func:`~radtrax.spectroscopy.table_cache_dir`
            when ``None``.
        download: Fetch the spectrum if it is not cached
            (:func:`fetch_solar_spectrum`); when ``False``, a missing file
            raises instead.

    Returns:
        $E_0$ [W m⁻² (cm⁻¹)⁻¹], ``(G,)``. Divide by $d^2$
        (:func:`earth_sun_distance`) for a given date.

    Raises:
        ParameterError: If ``source`` is unknown.
        MissingTableError: If the spectrum is not cached and ``download`` is
            ``False`` or the download fails; the message names the fetch
            command.
        GridError: If any wavenumber is outside the source's range.
    """
    spec = _source(source)
    path = solar_spectrum_path(source, cache_dir=cache_dir)
    if not path.is_file():
        if not download:
            raise MissingTableError(
                f"No {source} solar spectrum at {path}. Fetch it (61 MB, once):"
                f"\n    {_fetch_command(source)}\n({spec.url})"
            )
        fetch_solar_spectrum(source, cache_dir=cache_dir)
    nu_ref, e_ref = _spectrum_cm(str(path))
    nu = np.asarray(nu_cm, dtype=np.float64)
    if nu.size and (nu.min() < nu_ref[0] or nu.max() > nu_ref[-1]):
        raise GridError(
            f"Wavenumbers {nu.min():.2f}–{nu.max():.2f} cm⁻¹ fall outside the "
            f"{source} solar spectrum ({nu_ref[0]:.2f}–{nu_ref[-1]:.2f} cm⁻¹)."
        )
    return jnp.asarray(np.interp(nu, nu_ref, e_ref))

earth_sun_distance

earth_sun_distance(day_of_year: Real[ArrayLike, '*shape']) -> Float[Array, '*shape']

The Earth–Sun distance [AU] on a day of the year.

\[d = 1 - 0.01673 \cos\big(0.0172\,(t - 4)\big),\]

with \(t\) the day of the year (1 on 1 January; fractional days allowed): the first-order eccentricity term with perihelion on 4 January, the formula georeader's earth_sun_distance_correction_factor uses for EMIT, Sentinel-2 and PRISMA reflectance. The irradiance on that day is \(E_0 / d^2\).

Parameters:

Name	Type	Description	Default
day_of_year	Real[ArrayLike, '*shape']	The day of the year, any shape.	required

Returns:

Type	Description
Float[Array, '*shape']	\(d\) [AU], same shape.

Examples:

>>> from radtrax.spectroscopy.solar import earth_sun_distance
>>> print(f"{float(earth_sun_distance(4.0)):.5f}")
0.98327
>>> print(f"{float(earth_sun_distance(186.0)):.5f}")
1.01673

Source code in src/radtrax/spectroscopy/solar.py

def earth_sun_distance(
    day_of_year: Real[ArrayLike, "*shape"],
) -> Float[Array, "*shape"]:
    """The Earth–Sun distance [AU] on a day of the year.

    $$d = 1 - 0.01673 \\cos\\big(0.0172\\,(t - 4)\\big),$$

    with $t$ the day of the year (1 on 1 January; fractional days allowed):
    the first-order eccentricity term with perihelion on 4 January, the
    formula georeader's ``earth_sun_distance_correction_factor`` uses for
    EMIT, Sentinel-2 and PRISMA reflectance. The irradiance on that day is
    $E_0 / d^2$.

    Args:
        day_of_year: The day of the year, any shape.

    Returns:
        $d$ [AU], same shape.

    Examples:
        >>> from radtrax.spectroscopy.solar import earth_sun_distance
        >>> print(f"{float(earth_sun_distance(4.0)):.5f}")
        0.98327
        >>> print(f"{float(earth_sun_distance(186.0)):.5f}")
        1.01673
    """
    t = jnp.asarray(day_of_year, dtype=float)
    return 1.0 - 0.01673 * jnp.cos(0.0172 * (t - 4.0))

synthetic

Synthetic methane spectra for docs, notebooks and tests, without HITRAN.

Two spectra, each built deterministically from a seed rather than shipped:

• :func:synthetic_table: radtrax · 3's cross-section table, 150 Voigt lines between 4252 and 4348 cm⁻¹ on 4230–4370 cm⁻¹ at 0.01 cm⁻¹, ten pressures and six temperatures, in the format of :func:radtrax.spectroscopy.hitran.build_lut_dataset (about a second to build). Line strengths are log-normal at 296 K, lower-state energies 0–600 cm⁻¹, Lorentz half-widths γ_air (p / 1 atm)(296 / T)^0.75 and Doppler widths of CH₄ at the line's temperature.
• :func:methane_cross_section: radtrax · 1 and · 2's monochromatic cross-section on :data:WL_NM [m²/kg]: Lorentz lines drawn at random in the 1.66 and 2.3 μm bands, scaled so that the averages under two band responses (Sentinel-2's B11 and B12 in those notebooks) are chosen values.

They were tests/_synthetic_table.py here and notebooks/geomodels/synthetic_ch4.py in geostack; both are bit for bit the same here (radtrax#115). SciPy (the Voigt profile) is imported inside :func:synthetic_table.

Shapes: G wavenumbers or wavelengths; T temperatures; P pressures.

synthetic_table

synthetic_table(seed: int = 0) -> Dataset

radtrax · 3's σ(ν; T, p) [cm²/molecule] on (wavenumber, temperature, pressure [atm]).

Parameters:

Name	Type	Description	Default
seed	int	Seed of the line draw (numpy.random.default_rng).	0

Returns:

Type	Description
Dataset	absorption_cross_section on :data:NU_CM, :data:T_K and
Dataset	data:P_ATM.

Examples:

>>> table = synthetic_table()
>>> dict(table.sizes)
{'wavenumber': 14000, 'temperature': 6, 'pressure': 10}

Source code in src/radtrax/spectroscopy/synthetic.py

def synthetic_table(seed: int = 0) -> xr.Dataset:
    """*radtrax · 3*'s σ(ν; T, p) [cm²/molecule] on (wavenumber, temperature, pressure [atm]).

    Args:
        seed: Seed of the line draw (``numpy.random.default_rng``).

    Returns:
        ``absorption_cross_section`` on :data:`NU_CM`, :data:`T_K` and
        :data:`P_ATM`.

    Examples:
        >>> table = synthetic_table()
        >>> dict(table.sizes)
        {'wavenumber': 14000, 'temperature': 6, 'pressure': 10}
    """
    from scipy.special import voigt_profile

    rng = np.random.default_rng(seed)
    line_nu = rng.uniform(4252.0, 4348.0, _N_LINES)
    line_s = rng.lognormal(0.0, 1.2, _N_LINES) * 2e-22  # cm/molecule at 296 K
    line_e = rng.uniform(0.0, 600.0, _N_LINES)  # lower-state energy, cm⁻¹
    sigma = np.zeros((NU_CM.size, T_K.size, P_ATM.size))
    for j, t in enumerate(T_K):
        # S(T) = S₀ (296/T)^1.5 exp(−c₂ E (1/T − 1/296)); σ_D = ν √(kT/m) / c
        s_t = line_s * (296.0 / t) ** 1.5 * np.exp(-_C2 * line_e * (1 / t - 1 / 296.0))
        sig_d = (
            line_nu * np.sqrt(BOLTZMANN_J_PER_K * t / (CH4_MASS_U * _AMU_KG)) / _C_LIGHT
        )
        for i, p in enumerate(P_ATM):
            gamma = _GAMMA_AIR * p * (296.0 / t) ** _N_AIR
            for k in range(_N_LINES):
                lo, hi = np.searchsorted(NU_CM, [line_nu[k] - 5, line_nu[k] + 5])
                sigma[lo:hi, j, i] += s_t[k] * voigt_profile(
                    NU_CM[lo:hi] - line_nu[k], sig_d[k], gamma
                )
    return xr.Dataset(
        {
            "absorption_cross_section": (
                ("wavenumber", "temperature", "pressure"),
                sigma,
                {"units": "cm^2 / molecule"},
            )
        },
        coords={"wavenumber": NU_CM, "temperature": T_K, "pressure": P_ATM},
        attrs={
            "description": "radtrax · 3 synthetic table: 150 Voigt lines near 4300 cm⁻¹"
        },
    )

band_mean

band_mean(x: Float[ndarray, '*batch G'], f: Float[ndarray, ' G']) -> float | Float[ndarray, '*batch']

⟨x⟩_b = ∫ x f dλ / ∫ f dλ on a shared uniform grid.

Parameters:

Name	Type	Description	Default
x	Float[ndarray, '*batch G']	Values on the grid, (..., G).	required
f	Float[ndarray, ' G']	A band response on the same grid, (G,).	required

Returns:

Type	Description
float | Float[ndarray, '*batch']	The response-weighted mean over the last axis (a float for a 1-D x).

Source code in src/radtrax/spectroscopy/synthetic.py

def band_mean(
    x: Float[np.ndarray, "*batch G"], f: Float[np.ndarray, " G"]
) -> float | Float[np.ndarray, "*batch"]:
    """``⟨x⟩_b = ∫ x f dλ / ∫ f dλ`` on a shared uniform grid.

    Args:
        x: Values on the grid, ``(..., G)``.
        f: A band response on the same grid, ``(G,)``.

    Returns:
        The response-weighted mean over the last axis (a float for a 1-D ``x``).
    """
    return (x * f).sum(-1) / f.sum()

methane_cross_section

methane_cross_section(rng: Generator, f11: Float[ndarray, ' G'], f12: Float[ndarray, ' G'], sig12: float = 0.4, ratio11: float = 0.1) -> Float[ndarray, ' G']

radtrax · 1 and · 2's synthetic CH₄ cross-section σ(λ) on :data:WL_NM [m²/kg].

Draws 3,000 lines in the 2.3 μm band, then 600 in the 1.66 μm band, from rng (in that order, so a caller's later draws are unchanged), and scales them so that ⟨σ⟩₁₂ = sig12 and ⟨σ⟩₁₁ = ratio11 · sig12 under the band responses f12 and f11.

Parameters:

Name	Type	Description	Default
rng	Generator	The caller's generator, advanced by the draws.	required
f11	Float[ndarray, ' G']	The 1.66 μm band's response on :data:WL_NM (Sentinel-2 B11 in the notebooks), (G,).	required
f12	Float[ndarray, ' G']	The 2.3 μm band's response on :data:WL_NM (B12), (G,).	required
sig12	float	The 2.3 μm band's mean cross-section [m²/kg].	0.4
ratio11	float	The 1.66 μm band's mean over the 2.3 μm band's.	0.1

Returns:

Type	Description
Float[ndarray, ' G']	σ(λ) [m²/kg], (G,).

Examples:

>>> f11 = ((WL_NM >= 1569) & (WL_NM <= 1658)).astype(float)
>>> f12 = ((WL_NM >= 2113) & (WL_NM <= 2286)).astype(float)
>>> sigma = methane_cross_section(np.random.default_rng(0), f11, f12)
>>> (
...     round(float(band_mean(sigma, f12)), 12),
...     round(float(band_mean(sigma, f11)), 12),
... )
(0.4, 0.04)

Source code in src/radtrax/spectroscopy/synthetic.py

def methane_cross_section(
    rng: np.random.Generator,
    f11: Float[np.ndarray, " G"],
    f12: Float[np.ndarray, " G"],
    sig12: float = 0.40,
    ratio11: float = 0.10,
) -> Float[np.ndarray, " G"]:
    """*radtrax · 1* and *· 2*'s synthetic CH₄ cross-section σ(λ) on :data:`WL_NM` [m²/kg].

    Draws 3,000 lines in the 2.3 μm band, then 600 in the 1.66 μm band, from
    ``rng`` (in that order, so a caller's later draws are unchanged), and
    scales them so that ``⟨σ⟩₁₂ = sig12`` and ``⟨σ⟩₁₁ = ratio11 · sig12``
    under the band responses ``f12`` and ``f11``.

    Args:
        rng: The caller's generator, advanced by the draws.
        f11: The 1.66 μm band's response on :data:`WL_NM` (Sentinel-2 B11
            in the notebooks), ``(G,)``.
        f12: The 2.3 μm band's response on :data:`WL_NM` (B12), ``(G,)``.
        sig12: The 2.3 μm band's mean cross-section [m²/kg].
        ratio11: The 1.66 μm band's mean over the 2.3 μm band's.

    Returns:
        σ(λ) [m²/kg], ``(G,)``.

    Examples:
        >>> f11 = ((WL_NM >= 1569) & (WL_NM <= 1658)).astype(float)
        >>> f12 = ((WL_NM >= 2113) & (WL_NM <= 2286)).astype(float)
        >>> sigma = methane_cross_section(np.random.default_rng(0), f11, f12)
        >>> (
        ...     round(float(band_mean(sigma, f12)), 12),
        ...     round(float(band_mean(sigma, f11)), 12),
        ... )
        (0.4, 0.04)
    """
    band_23 = _lorentz_sum(*_lines(rng, 3000, 2165.0, 2420.0, 2310.0, 70.0))
    band_166 = _lorentz_sum(*_lines(rng, 600, 1620.0, 1710.0, 1666.0, 15.0))
    sigma = band_23 * sig12 / band_mean(band_23, f12)
    return sigma + band_166 * (ratio11 * sig12 - band_mean(sigma, f11)) / band_mean(
        band_166, f11
    )

table_cache

The cross-section tables of v0: where they live and how to load them.

The tables are hundreds of MB, so they don't ship with radtrax (decision D6). Each is built once with the [hitran] extra,

.. code-block:: console

python -m radtrax.spectroscopy.hitran CH4 H2O CO2 N2O CO

and cached under the user cache directory, keyed by gas, grid, HITRAN version and table version; :func:load_table and :func:load_opacity read them back, and a missing table raises :class:~radtrax.MissingTableError naming that command. Loading needs no extra (xarray and netCDF4 are base dependencies).

V0_GRID is the grid v0 builds on: 4000–6500 cm⁻¹ (Sentinel-2's bands 11 and 12 and the 1.6 and 2.3 µm windows of imaging spectrometers) at 0.01 cm⁻¹, 16 pressures from 0.01 to 1.1 atm (the surface to ~31 km, log spaced, ratio 1.37) and 8 temperatures from 180 to 320 K: 128 MB per gas as float32. It is a proposal that gate V2 checks: band transmittance within 0.1 % of HAPI through the AFGL US Standard atmosphere (tests/spectroscopy/test_v2_hapi.py, which also reports the monochromatic difference). A single-line study puts the monochromatic error of this grid near 4.5e-3 and that of 32 pressures × 15 temperatures near 1e-3 (480 MB per gas); refine the grid if the band gate fails.

TableGrid dataclass

TableGrid(name: str, nu_min_cm: float, nu_max_cm: float, nu_step_cm: float, p_atm: tuple[float, ...], t_k: tuple[float, ...])

The (ν, p, T) grid a cross-section table is built on.

Attributes:

Name	Type	Description
name	str	A short name for messages and the CLI ("v0").
nu_min_cm	float	First wavenumber [cm⁻¹].
nu_max_cm	float	Last wavenumber [cm⁻¹].
nu_step_cm	float	Wavenumber step [cm⁻¹].
p_atm	tuple[float, ...]	Table pressures [atm], increasing.
t_k	tuple[float, ...]	Table temperatures [K], increasing.

key

key() -> str

A short hash of the grid's values (not its name), for file names.

Source code in src/radtrax/spectroscopy/table_cache.py

def key(self) -> str:
    """A short hash of the grid's values (not its name), for file names."""
    values = {k: v for k, v in asdict(self).items() if k != "name"}
    digest = hashlib.sha1(json.dumps(values, sort_keys=True).encode()).hexdigest()
    return digest[:10]

table_cache_dir

table_cache_dir() -> Path

Where tables are cached.

$RADTRAX_TABLE_DIR if set, else $XDG_CACHE_HOME/radtrax/tables, else ~/.cache/radtrax/tables.

Source code in src/radtrax/spectroscopy/table_cache.py

def table_cache_dir() -> Path:
    """Where tables are cached.

    ``$RADTRAX_TABLE_DIR`` if set, else ``$XDG_CACHE_HOME/radtrax/tables``,
    else ``~/.cache/radtrax/tables``.
    """
    env = os.environ.get(TABLE_DIR_ENV)
    if env:
        return Path(env).expanduser()
    base = os.environ.get("XDG_CACHE_HOME")
    root = Path(base).expanduser() if base else Path.home() / ".cache"
    return root / "radtrax" / "tables"

table_path

table_path(gas: str, grid: TableGrid = V0_GRID, *, hitran_version: str = HITRAN_VERSION, cache_dir: Path | str | None = None) -> Path

The cache file of gas on grid: <gas>_<grid key>_<HITRAN>_v<table>.nc.

Source code in src/radtrax/spectroscopy/table_cache.py

def table_path(
    gas: str,
    grid: TableGrid = V0_GRID,
    *,
    hitran_version: str = HITRAN_VERSION,
    cache_dir: Path | str | None = None,
) -> Path:
    """The cache file of ``gas`` on ``grid``: ``<gas>_<grid key>_<HITRAN>_v<table>.nc``."""
    root = Path(cache_dir) if cache_dir is not None else table_cache_dir()
    return root / f"{gas}_{grid.key()}_{hitran_version}_v{TABLE_VERSION}.nc"

build_command

build_command(gases: Sequence[str], grid: TableGrid = V0_GRID) -> str

The shell command that builds gases on grid.

Source code in src/radtrax/spectroscopy/table_cache.py

def build_command(gases: Sequence[str], grid: TableGrid = V0_GRID) -> str:
    """The shell command that builds ``gases`` on ``grid``."""
    grid_arg = "" if grid == V0_GRID else f" --grid {grid.name}"
    return f"python -m radtrax.spectroscopy.hitran {' '.join(gases)}{grid_arg}"

load_table

load_table(gas: str, grid: TableGrid = V0_GRID, *, hitran_version: str = HITRAN_VERSION, cache_dir: Path | str | None = None) -> Dataset

Open gas's cached table (lazily, with xarray).

Raises:

Type	Description
MissingTableError	If it isn't built; the message names the command.

Source code in src/radtrax/spectroscopy/table_cache.py

def load_table(
    gas: str,
    grid: TableGrid = V0_GRID,
    *,
    hitran_version: str = HITRAN_VERSION,
    cache_dir: Path | str | None = None,
) -> xr.Dataset:
    """Open ``gas``'s cached table (lazily, with xarray).

    Raises:
        MissingTableError: If it isn't built; the message names the command.
    """
    import xarray as xr

    path = table_path(gas, grid, hitran_version=hitran_version, cache_dir=cache_dir)
    if not path.is_file():
        raise MissingTableError(
            f"No {gas} cross-section table for grid {grid.name!r} ({hitran_version}) "
            f"at {path}. Build it (needs the [hitran] extra and network access to "
            f"HITRAN):\n    pip install 'radtrax[hitran]'\n    "
            f"{build_command([gas], grid)}"
        )
    return xr.open_dataset(path)

load_opacity

load_opacity(gases: Sequence[str], grid: TableGrid = V0_GRID, *, hitran_version: str = HITRAN_VERSION, cache_dir: Path | str | None = None) -> TableOpacity

A :class:~radtrax.spectroscopy.TableOpacity from cached tables.

Raises:

Type	Description
MissingTableError	Naming every gas that isn't built, with one command that builds them all.

Source code in src/radtrax/spectroscopy/table_cache.py

def load_opacity(
    gases: Sequence[str],
    grid: TableGrid = V0_GRID,
    *,
    hitran_version: str = HITRAN_VERSION,
    cache_dir: Path | str | None = None,
) -> TableOpacity:
    """A :class:`~radtrax.spectroscopy.TableOpacity` from cached tables.

    Raises:
        MissingTableError: Naming every gas that isn't built, with one
            command that builds them all.
    """
    from radtrax.spectroscopy.opacity import TableOpacity

    missing = [
        gas
        for gas in gases
        if not table_path(
            gas, grid, hitran_version=hitran_version, cache_dir=cache_dir
        ).is_file()
    ]
    if missing:
        raise MissingTableError(
            f"No cross-section tables for {missing} on grid {grid.name!r} "
            f"({hitran_version}). Build them (needs the [hitran] extra):\n    "
            f"pip install 'radtrax[hitran]'\n    {build_command(missing, grid)}"
        )
    return TableOpacity.from_datasets(
        {
            gas: load_table(
                gas, grid, hitran_version=hitran_version, cache_dir=cache_dir
            )
            for gas in gases
        }
    )

tables

LUT interpolation + Beer-Lambert transmittance for one-shot forward models.

This module closes the loop between an offline-generated σ(ν, T, P) LUT and the naïve single-layer radiative-transfer model used in the tutorial notebooks. The physics is:

N_total = p / (k_B T)              [molecules/cm³]  (ideal gas)
N_gas   = N_total · VMR_gas        [molecules/cm³]
α(ν)    = σ(ν; T, p) · N_gas       [1/cm]
AMF     = 1/cos(SZA) + 1/cos(VZA)  [plane-parallel, two-way]
τ(ν)    = exp[ -α(ν) · L_vert · AMF ]   [Beer-Lambert]

For plume-enhancement retrievals the differential form is

τ_ratio(ν) = τ_total(ν) / τ_background(ν)
           = exp[ -σ(ν; T, p) · N_total · L_vert · AMF · (VMR_total - VMR_bg) ]

— any factor common to the "plume" and "background" pixel (solar spectrum, broadband albedo, aerosol) cancels, leaving the narrow-band gas signature.

Units: the LUT is in HAPI's HITRAN convention (σ in cm², p in atm). We work in cgs throughout — pressures are converted to Pa only to evaluate number density via the ideal-gas law with k_B = 1.380649e-23 J/K.

Ported from plumax 0.1.2 at 28f257e (src/plumax/hapi_lut/beers.py); only import paths changed.

number_density

number_density(p_atm: Real[ArrayLike, '*#shape'], T_K: Real[ArrayLike, '*#shape']) -> Float[ArrayLike, '*#shape']

Total number density from the ideal-gas law [molecules / cm³].

\[ n = \frac{p}{k_B T} \]

with \(p\) converted from atm to Pa and the result from m⁻³ to cm⁻³. Element-wise on arrays; no validation (radtrax.slab.number_density_cm3 adds it for scalar inputs).

Parameters:

Name	Type	Description	Default
p_atm	Real[ArrayLike, '*#shape']	Pressure [atm].	required
T_K	Real[ArrayLike, '*#shape']	Temperature [K].	required

Returns:

Type	Description
Float[ArrayLike, '*#shape']	Number density [molecules / cm³], the broadcast shape of the inputs.

Examples:

>>> f"{number_density(1.0, 296.0):.4e}"
'2.4794e+19'

Source code in src/radtrax/atmosphere/gas.py

def number_density(
    p_atm: Real[ArrayLike, "*#shape"], T_K: Real[ArrayLike, "*#shape"]
) -> Float[ArrayLike, "*#shape"]:
    """Total number density from the ideal-gas law [molecules / cm³].

    $$
    n = \\frac{p}{k_B T}
    $$

    with $p$ converted from atm to Pa and the result from m⁻³ to cm⁻³.
    Element-wise on arrays; no validation (``radtrax.slab.number_density_cm3``
    adds it for scalar inputs).

    Args:
        p_atm: Pressure [atm].
        T_K: Temperature [K].

    Returns:
        Number density [molecules / cm³], the broadcast shape of the inputs.

    Examples:
        >>> f"{number_density(1.0, 296.0):.4e}"
        '2.4794e+19'
    """
    p_pa = p_atm * ATM_TO_PA
    n_per_m3 = p_pa / (BOLTZMANN_J_PER_K * T_K)
    return n_per_m3 / CM3_PER_M3

air_mass_factor

air_mass_factor(sza_deg: Real[ArrayLike, '*#shape'], vza_deg: Real[ArrayLike, '*#shape']) -> Float[ArrayLike, '*#shape']

Plane-parallel two-way air-mass factor, sun → surface → sensor.

\[ m = \frac{1}{\cos\theta_0} + \frac{1}{\cos\theta} \]

with \(\theta_0\) the solar and \(\theta\) the viewing zenith angle. The plane-parallel form breaks down past about 75°, where a spherical (Chapman) correction is needed; that arrives with v1 of the roadmap.

Parameters:

Name	Type	Description	Default
sza_deg	Real[ArrayLike, '*#shape']	Solar zenith angle [degrees].	required
vza_deg	Real[ArrayLike, '*#shape']	Viewing zenith angle [degrees].	required

Returns:

Type	Description
Float[ArrayLike, '*#shape']	The two-way air-mass factor (dimensionless).

Examples:

>>> air_mass_factor(0.0, 0.0)
np.float64(2.0)
>>> round(float(air_mass_factor(60.0, 0.0)), 12)
3.0

Source code in src/radtrax/atmosphere/geometry.py

def air_mass_factor(
    sza_deg: Real[ArrayLike, "*#shape"], vza_deg: Real[ArrayLike, "*#shape"]
) -> Float[ArrayLike, "*#shape"]:
    """Plane-parallel two-way air-mass factor, sun → surface → sensor.

    $$
    m = \\frac{1}{\\cos\\theta_0} + \\frac{1}{\\cos\\theta}
    $$

    with $\\theta_0$ the solar and $\\theta$ the viewing zenith angle. The
    plane-parallel form breaks down past about 75°, where a spherical
    (Chapman) correction is needed; that arrives with v1 of the roadmap.

    Args:
        sza_deg: Solar zenith angle [degrees].
        vza_deg: Viewing zenith angle [degrees].

    Returns:
        The two-way air-mass factor (dimensionless).

    Examples:
        >>> air_mass_factor(0.0, 0.0)
        np.float64(2.0)
        >>> round(float(air_mass_factor(60.0, 0.0)), 12)
        3.0
    """
    sza = np.deg2rad(sza_deg)
    vza = np.deg2rad(vza_deg)
    return 1.0 / np.cos(sza) + 1.0 / np.cos(vza)

interpolate_cross_section

interpolate_cross_section(ds: Dataset, T_K: float, p_atm: float, *, var: str = 'absorption_cross_section') -> DataArray

Return σ(ν) at a given (T, p) via bilinear LUT interpolation.

The LUT dataset must carry temperature [K] and pressure [atm] coords as in build_lut_dataset. No extrapolation — callers should clip to the LUT range before interpolating.

Source code in src/radtrax/spectroscopy/tables.py

def interpolate_cross_section(
    ds: xr.Dataset,
    T_K: float,
    p_atm: float,
    *,
    var: str = "absorption_cross_section",
) -> xr.DataArray:
    """Return σ(ν) at a given (T, p) via bilinear LUT interpolation.

    The LUT dataset must carry ``temperature`` [K] and ``pressure`` [atm]
    coords as in ``build_lut_dataset``. No extrapolation — callers should
    clip to the LUT range before interpolating.
    """
    if var not in ds:
        raise MissingKeyError(
            f"variable {var!r} not in dataset (have {list(ds.data_vars)})"
        )
    return ds[var].interp(temperature=T_K, pressure=p_atm)

absorption_coefficient

absorption_coefficient(sigma: DataArray | Float[ndarray, '*shape'], *, vmr: float, p_atm: float, T_K: float) -> DataArray | Float[ndarray, '*shape']

α(ν) = σ(ν) · N_gas [1/cm].

Source code in src/radtrax/spectroscopy/tables.py

def absorption_coefficient(
    sigma: xr.DataArray | Float[np.ndarray, "*shape"],
    *,
    vmr: float,
    p_atm: float,
    T_K: float,
) -> xr.DataArray | Float[np.ndarray, "*shape"]:
    """α(ν) = σ(ν) · N_gas [1/cm]."""
    n_total = number_density(p_atm, T_K)
    return sigma * (n_total * vmr)

transmittance

transmittance(alpha: DataArray | Float[ndarray, '*shape'], *, l_vert_cm: float, sza_deg: float, vza_deg: float) -> DataArray | Float[ndarray, '*shape']

τ(ν) = exp(-α · L_vert · AMF).

Source code in src/radtrax/spectroscopy/tables.py

def transmittance(
    alpha: xr.DataArray | Float[np.ndarray, "*shape"],
    *,
    l_vert_cm: float,
    sza_deg: float,
    vza_deg: float,
) -> xr.DataArray | Float[np.ndarray, "*shape"]:
    """τ(ν) = exp(-α · L_vert · AMF)."""
    amf = air_mass_factor(sza_deg, vza_deg)
    return np.exp(-alpha * l_vert_cm * amf)

beers_law_from_lut

beers_law_from_lut(ds: Dataset, *, vmr: float, T_K: float, p_atm: float, l_vert_cm: float, sza_deg: float, vza_deg: float, var: str = 'absorption_cross_section') -> DataArray

End-to-end: σ-from-LUT → α → τ(ν) for a single homogeneous layer.

Returns an xarray.DataArray carrying the wavenumber + wavelength coords so callers can plot τ(λ) without re-introducing units.

Source code in src/radtrax/spectroscopy/tables.py

def beers_law_from_lut(
    ds: xr.Dataset,
    *,
    vmr: float,
    T_K: float,
    p_atm: float,
    l_vert_cm: float,
    sza_deg: float,
    vza_deg: float,
    var: str = "absorption_cross_section",
) -> xr.DataArray:
    """End-to-end: σ-from-LUT → α → τ(ν) for a single homogeneous layer.

    Returns an ``xarray.DataArray`` carrying the ``wavenumber`` + ``wavelength``
    coords so callers can plot τ(λ) without re-introducing units.
    """
    sigma = interpolate_cross_section(ds, T_K=T_K, p_atm=p_atm, var=var)
    alpha = absorption_coefficient(sigma, vmr=vmr, p_atm=p_atm, T_K=T_K)
    tau = transmittance(alpha, l_vert_cm=l_vert_cm, sza_deg=sza_deg, vza_deg=vza_deg)
    return tau.assign_attrs(
        description="Single-layer Beer-Lambert transmittance from LUT",
        vmr=vmr,
        T_K=T_K,
        p_atm=p_atm,
        l_vert_cm=l_vert_cm,
        sza_deg=sza_deg,
        vza_deg=vza_deg,
    )

plume_ratio_spectrum

plume_ratio_spectrum(ds: Dataset, *, vmr_background: float, vmr_total: float, T_K: float, p_atm: float, l_vert_cm: float, sza_deg: float, vza_deg: float, var: str = 'absorption_cross_section') -> DataArray

Differential-Beer-Lambert ratio τ_total(ν) / τ_bg(ν).

Callers typically pass vmr_total = vmr_background * (1 + enhancement) where enhancement is the fitting parameter in the plume retrieval.

Source code in src/radtrax/spectroscopy/tables.py

def plume_ratio_spectrum(
    ds: xr.Dataset,
    *,
    vmr_background: float,
    vmr_total: float,
    T_K: float,
    p_atm: float,
    l_vert_cm: float,
    sza_deg: float,
    vza_deg: float,
    var: str = "absorption_cross_section",
) -> xr.DataArray:
    """Differential-Beer-Lambert ratio τ_total(ν) / τ_bg(ν).

    Callers typically pass ``vmr_total = vmr_background * (1 + enhancement)``
    where ``enhancement`` is the fitting parameter in the plume retrieval.
    """
    tau_bg = beers_law_from_lut(
        ds,
        vmr=vmr_background,
        T_K=T_K,
        p_atm=p_atm,
        l_vert_cm=l_vert_cm,
        sza_deg=sza_deg,
        vza_deg=vza_deg,
        var=var,
    )
    tau_total = beers_law_from_lut(
        ds,
        vmr=vmr_total,
        T_K=T_K,
        p_atm=p_atm,
        l_vert_cm=l_vert_cm,
        sza_deg=sza_deg,
        vza_deg=vza_deg,
        var=var,
    )
    ratio = tau_total / tau_bg
    return ratio.assign_attrs(
        description="Plume/background transmittance ratio (differential Beer-Lambert)",
        vmr_background=vmr_background,
        vmr_total=vmr_total,
    )

HITRAN table building

Needs the [hitran] extra; HAPI is imported inside the builders.

hitran

Offline HITRAN cross-section table building through HAPI.

Submodules

• config — GasConfig / LUTGridConfig / ATMOSPHERIC_GASES.
• generator — single-gas fetch / compute / wrap / save pipeline.
• multi — multi-gas orchestration (separate-per-gas + combined).
• build — v0's tables into the cache, and its command line, python -m radtrax.spectroscopy.hitran GAS [GAS ...] (radtrax#34).

HAPI (hitran-api, the [hitran] extra) is imported lazily inside the generator and multi pipelines: importing this subpackage does not need HAPI, only building or fetching a table does. Interpolating a built table is radtrax.spectroscopy.tables.

Ported from plumax 0.1.2 at 28f257e (plumax.hapi_lut).

config

Gas + grid configuration for HAPI absorption-cross-section look-up tables.

A LUT is an N-dimensional array of absorption cross-sections

σ(ν, T, P) = Σ_l S_l(T) f_V(ν − ν_l; Γ_D, Γ_L) · 1[|ν − ν_l| ≤ w]    [cm² / molecule]

computed line-by-line via HAPI's Voigt routine, over the lines l of every isotopologue of the gas, each cut at w = ±25 cm⁻¹ from its centre, the convention continua are defined against (Korkin et al. 2025, eqs. 16–18). HITRAN's intensities S_l are already weighted by each isotopologue's natural abundance, so the sum needs no rescaling. GasConfig identifies the HITRAN line set (molecule + isotopologues + spectral range) and LUTGridConfig describes the common (T, P, ν-step, cut-off) grid used across gases.

The ATMOSPHERIC_GASES registry covers the six species most commonly needed for SWIR/NIR remote sensing (CH₄, CO₂, H₂O, O₂, N₂O, CO). Add entries for other HITRAN species by constructing a GasConfig directly.

Ported from plumax 0.1.2 at 28f257e (src/plumax/hapi_lut/config.py); every isotopologue and the ±25 cm⁻¹ cut-off since radtrax#172.

GasConfig dataclass

GasConfig(name: str, molecule_id: int, nu_min: float, nu_max: float, description: str = '', *, isotopologue_ids: tuple[int, ...] | None = None, isotopologue_id: int | None = None)

Configuration for a single HITRAN gas species.

Attributes:

Name	Type	Description
name	str	Short name used as HAPI SourceTable id (e.g. 'CH4').
molecule_id	int	HITRAN molecule id (see hitran.org).
nu_min	float	Lower wavenumber bound [cm⁻¹].
nu_max	float	Upper wavenumber bound [cm⁻¹].
description	str	Free-text description for dataset metadata.
isotopologue_ids	tuple[int, ...] | None	HITRAN global isotopologue ids (the id of hapi.ISO, e.g. 32–35 for CH₄); None (the default) takes every isotopologue of the molecule. Restrict it only on purpose, e.g. for an isotopologue retrieval.
isotopologue_id	int | None	Deprecated: one local isotopologue number (1 = most abundant); use isotopologue_ids.

Raises:

Type	Description
ParameterError	If both isotopologue_ids and isotopologue_id are given, or description is not a string (before radtrax#172 the positional order was name, molecule_id, isotopologue_id, nu_min, nu_max).

resolve_isotopologue_ids

resolve_isotopologue_ids() -> tuple[int, ...]

The HITRAN global isotopologue ids to fetch (needs [hitran]).

isotopologue_ids as given; else the deprecated local isotopologue_id mapped through hapi.ISO; else every global id whose hapi.ISO entry has this molecule_id, in HITRAN's isotopologue order.

Raises:

Type	Description
MissingExtraError	If the ids must come from hapi.ISO and the [hitran] extra is not installed.
KeyError	If hapi.ISO has no such molecule or isotopologue.

Source code in src/radtrax/spectroscopy/hitran/config.py

def resolve_isotopologue_ids(self) -> tuple[int, ...]:
    """The HITRAN global isotopologue ids to fetch (needs ``[hitran]``).

    ``isotopologue_ids`` as given; else the deprecated local
    ``isotopologue_id`` mapped through ``hapi.ISO``; else every global id
    whose ``hapi.ISO`` entry has this ``molecule_id``, in HITRAN's
    isotopologue order.

    Raises:
        MissingExtraError: If the ids must come from ``hapi.ISO`` and the
            ``[hitran]`` extra is not installed.
        KeyError: If ``hapi.ISO`` has no such molecule or isotopologue.
    """
    if self.isotopologue_ids is not None:
        return self.isotopologue_ids
    try:
        from hapi import ISO
    except ModuleNotFoundError as err:
        if err.name != "hapi":
            raise
        raise MissingExtraError(
            f"GasConfig({self.name!r}).resolve_isotopologue_ids() needs the "
            "[hitran] extra: pip install 'radtrax[hitran]'"
        ) from err

    global_ids = {
        local: int(entry[0])
        for (molecule, local), entry in ISO.items()
        if molecule == self.molecule_id and self.isotopologue_id in (None, local)
    }
    if not global_ids:
        number = "" if self.isotopologue_id is None else f" {self.isotopologue_id}"
        raise MissingKeyError(
            f"GasConfig({self.name!r}): hapi.ISO has no isotopologue{number} "
            f"of molecule {self.molecule_id}."
        )
    return tuple(global_ids[local] for local in sorted(global_ids))

LUTGridConfig dataclass

LUTGridConfig(T_grid: ndarray = (lambda: array([200, 220, 240, 260, 280, 300, 320], dtype=float))(), P_grid: ndarray = (lambda: array([0.1, 0.2, 0.4, 0.6, 0.8, 1.0], dtype=float))(), nu_step: float = 0.01, diluent_composition: dict[str, float] | None = None, wing_cm: float = 25.0, wing_halfwidths: float = 0.0)

Grid configuration shared across gases for a single LUT build.

Attributes:

Name	Type	Description
T_grid	ndarray	Temperature knots [K]. Should bracket the atmospheric range of interest (upper-troposphere ~200 K to surface ~320 K).
P_grid	ndarray	Pressure knots [atm]. Should bracket the altitude range of interest (upper-tropo ~0.1 atm to surface ~1.0 atm).
nu_step	float	Spectral step [cm⁻¹] on the output wavenumber grid. HITRAN lines have HWHMs down to ~0.005 cm⁻¹ at low pressure — step ≤ 0.01 cm⁻¹ for hyperspectral work, 0.05 cm⁻¹ for tutorials.
diluent_composition	dict[str, float] | None	Optional override for the air diluent fractions. If None, each gas gets {'air': 1-VMR, 'self': VMR}.
wing_cm	float	Absolute line cut-off w [cm⁻¹]: each line is evaluated within ±w of its centre (HAPI's WavenumberWing). The default 25 cm⁻¹ is the convention HITRAN line parameters are used with and continua (MT_CKD) are defined against.
wing_halfwidths	float	Relative cut-off in half-widths max(Γ_L, Γ_D) (HAPI's WavenumberWingHW); the larger of the two applies. 0 by default, so the absolute cut-off alone applies.

get_diluent_for_gas

get_diluent_for_gas(gas_name: str, vmr_nominal: float | None = None) -> dict[str, float]

Return the HAPI Diluent dict for a given gas.

If an explicit diluent_composition override was supplied at construction time it is returned verbatim; otherwise a standard {'air': 1 - VMR, 'self': VMR} split is built from DEFAULT_VMR_NOMINAL[gas_name] (falling back to 1e-6).

Source code in src/radtrax/spectroscopy/hitran/config.py

def get_diluent_for_gas(
    self,
    gas_name: str,
    vmr_nominal: float | None = None,
) -> dict[str, float]:
    """Return the HAPI ``Diluent`` dict for a given gas.

    If an explicit ``diluent_composition`` override was supplied at
    construction time it is returned verbatim; otherwise a standard
    ``{'air': 1 - VMR, 'self': VMR}`` split is built from
    ``DEFAULT_VMR_NOMINAL[gas_name]`` (falling back to 1e-6).
    """
    if self.diluent_composition is not None:
        return dict(self.diluent_composition)
    if vmr_nominal is None:
        vmr_nominal = DEFAULT_VMR_NOMINAL.get(gas_name, 1e-6)
    return {"air": 1.0 - vmr_nominal, "self": vmr_nominal}

generator

Single-gas absorption-cross-section LUT builder.

Given a GasConfig and a LUTGridConfig this module fetches HITRAN line data (cached under a user-provided directory), evaluates a Voigt-profile line-by-line sum on every (T, P) grid cell, and wraps the result as a CF-conformant xarray.Dataset ready for NetCDF persistence.

HAPI is imported lazily inside the functions that need it so import radtrax.spectroscopy.hitran does not require hitran-api to be installed.

Ported from plumax 0.1.2 at 28f257e (src/plumax/hapi_lut/generator.py); every isotopologue (fetch_by_ids) and the ±25 cm⁻¹ cut-off since radtrax#172.

default_cache_dir

default_cache_dir() -> Path

Return the HITRAN line-data cache directory.

Resolution order

1. $HAPI_CACHE_DIR if set.
2. ~/.cache/radtrax/hitran — a per-user cache, so that invoking the generator from anywhere on disk does not leave unignored artefacts in the repository tree. The projects/plumax/data/hitran_cache path used by the notebooks is always passed explicitly.

Source code in src/radtrax/spectroscopy/hitran/generator.py

def default_cache_dir() -> Path:
    """Return the HITRAN line-data cache directory.

    Resolution order:
      1. ``$HAPI_CACHE_DIR`` if set.
      2. ``~/.cache/radtrax/hitran`` — a per-user cache, so that
         invoking the generator from anywhere on disk does not leave
         unignored artefacts in the repository tree. The
         ``projects/plumax/data/hitran_cache`` path used by the
         notebooks is always passed explicitly.
    """
    env = os.environ.get(HAPI_CACHE_ENV)
    if env:
        return Path(env).expanduser().resolve()
    return Path.home() / ".cache" / "radtrax" / "hitran"

fetch_hitran_data

fetch_hitran_data(gas_config: GasConfig, cache_dir: Path | str | None = None) -> Path

Fetch HITRAN line parameters for gas_config into cache_dir.

Fetches the lines of :meth:GasConfig.resolve_isotopologue_ids <radtrax.spectroscopy.hitran.GasConfig.resolve_isotopologue_ids>, every isotopologue of the molecule by default, with HAPI's fetch_by_ids. HAPI writes a <GasConfig.name>.data + .header pair under cache_dir, overwriting any existing pair (HAPI does not no-op on a cache hit). A successful fetch also records the fetched ν range and isotopologue ids in a <name>.range.json sidecar.

When the fetch fails (no network, API error) the previously-cached pair is only reused if its recorded range covers the requested [nu_min, nu_max] and it holds the same isotopologues — otherwise a narrowed cache (e.g. left by :func:radtrax.spectroscopy.hitran.multi.create_combined_lut, which fetches over the gas intersection range) or one isotopologue's lines would silently serve a truncated line list.

Raises:

Type	Description
RuntimeError	if the fetch raises and either no cached pair is present or its recorded range or isotopologues do not match the request. The underlying HAPI exception is chained so the network / API error is visible.

Returns:

Type	Description
Path	The resolved cache directory.

Source code in src/radtrax/spectroscopy/hitran/generator.py

def fetch_hitran_data(
    gas_config: GasConfig,
    cache_dir: Path | str | None = None,
) -> Path:
    """Fetch HITRAN line parameters for ``gas_config`` into ``cache_dir``.

    Fetches the lines of :meth:`GasConfig.resolve_isotopologue_ids
    <radtrax.spectroscopy.hitran.GasConfig.resolve_isotopologue_ids>`, every
    isotopologue of the molecule by default, with HAPI's ``fetch_by_ids``.
    HAPI writes a ``<GasConfig.name>.data`` + ``.header`` pair under
    ``cache_dir``, **overwriting** any existing pair (HAPI does not no-op on
    a cache hit). A successful fetch also records the fetched ν range and
    isotopologue ids in a ``<name>.range.json`` sidecar.

    When the fetch fails (no network, API error) the previously-cached pair is
    only reused if its recorded range **covers** the requested
    ``[nu_min, nu_max]`` and it holds the same isotopologues — otherwise a
    narrowed cache (e.g. left by
    :func:`radtrax.spectroscopy.hitran.multi.create_combined_lut`, which
    fetches over the gas intersection range) or one isotopologue's lines would
    silently serve a truncated line list.

    Raises:
        RuntimeError: if the fetch raises and either no cached pair is
            present or its recorded range or isotopologues do not match the
            request. The underlying HAPI exception is chained so the network /
            API error is visible.

    Returns:
        The resolved cache directory.
    """
    from hapi import fetch_by_ids

    ids = gas_config.resolve_isotopologue_ids()
    cache = Path(cache_dir) if cache_dir is not None else default_cache_dir()
    _init_hapi_cache(cache)

    data_path = cache / f"{gas_config.name}.data"
    header_path = cache / f"{gas_config.name}.header"
    already_cached = data_path.exists() and header_path.exists()

    logger.info(
        "Fetching HITRAN %s (M%d, isotopologue ids %s, %.1f-%.1f cm^-1)%s",
        gas_config.name,
        gas_config.molecule_id,
        ids,
        gas_config.nu_min,
        gas_config.nu_max,
        " [cache present]" if already_cached else "",
    )

    fetched = True
    try:
        fetch_by_ids(gas_config.name, list(ids), gas_config.nu_min, gas_config.nu_max)
    except Exception as exc:
        # HAPI raises plain Exceptions on network failures and API errors. A
        # cached pair is only trustworthy if its recorded range covers the
        # request and it holds the same isotopologues — otherwise proceeding
        # would feed a truncated line list to the Voigt sum.
        fetched = False
        cached = _read_cached_range(cache, gas_config.name)
        if not already_cached:
            detail = "no cached data was found"
        elif cached is None:
            detail = "the cached line list has no recorded range"
        elif not _range_covers(cached, gas_config):
            detail = (
                f"the cached line list covers only {cached[0]:.1f}-{cached[1]:.1f} "
                "cm^-1"
            )
        elif cached[2] != ids:
            detail = (
                f"the cached line list holds isotopologue ids {cached[2]}, not {ids}"
            )
        else:
            detail = ""
        if detail:
            raise TableBuildError(
                f"HITRAN fetch failed for {gas_config.name} "
                f"(M{gas_config.molecule_id}, isotopologue ids {ids}, "
                f"{gas_config.nu_min:.1f}-{gas_config.nu_max:.1f} cm^-1) and "
                f"{detail} at {cache}. Re-fetch online or point `cache_dir` at "
                "a cache whose range spans the request."
            ) from exc
        logger.warning(
            "fetch_by_ids() raised %s — proceeding with cached %s (its range "
            "and isotopologues match the request).",
            type(exc).__name__,
            data_path,
        )

    if not (data_path.exists() and header_path.exists()):
        raise TableBuildError(
            f"HITRAN fetch for {gas_config.name} returned without error but "
            f"no {gas_config.name}.data/.header pair was produced under {cache}."
        )
    if fetched:
        _write_cached_range(
            cache, gas_config.name, gas_config.nu_min, gas_config.nu_max, ids
        )
    return cache

compute_absorption_lut

compute_absorption_lut(gas_config: GasConfig, grid_config: LUTGridConfig, *, cache_dir: Path | str | None = None, progress: bool = True) -> tuple[Float[ndarray, ' G'], Float[ndarray, 'G T P'], Float[ndarray, ' G']]

Evaluate σ(ν, T, P) on the full LUT grid.

σ(ν) = Σ_l S_l(T) f_V(ν − ν_l) · 1[|ν − ν_l| ≤ w] over every line in the fetched table, each cut at w = max(wing_cm, wing_halfwidths · max(Γ_L, Γ_D)) from its centre (HAPI's rule; ±25 cm⁻¹ by default). Assumes fetch_hitran_data has already populated the cache, over the grid plus w so lines just outside it add their wings. Fails fast — a HAPI failure at any (T, P) knot, or any non-finite cross-section in the result, raises rather than leaving silent NaN holes in the LUT.

Parameters:

Name	Type	Description	Default
gas_config	GasConfig	Target gas (molecule + isotopologues + ν range).	required
grid_config	LUTGridConfig	Common T/P/ν-step grid and line cut-off.	required
cache_dir	Path | str | None	Where the fetched HITRAN .data/.header pair lives.	None
progress	bool	Emit per-(T, P) log lines at INFO.	True

Raises:

Type	Description
RuntimeError	if absorptionCoefficient_Voigt raises at any grid knot, or if the returned cross-section contains any non-finite values.

Returns:

Type	Description
Float[ndarray, ' G']	(nu_grid, sigma_lut, wavelength_nm):
Float[ndarray, 'G T P']	nu_grid: wavenumbers [cm⁻¹], shape (n_nu,)
Float[ndarray, ' G']	sigma_lut: cross-sections [cm²/molecule], shape (n_nu, n_T, n_P).
tuple[Float[ndarray, ' G'], Float[ndarray, 'G T P'], Float[ndarray, ' G']]	wavelength_nm: wavelengths [nm], shape (n_nu,).

Source code in src/radtrax/spectroscopy/hitran/generator.py

def compute_absorption_lut(
    gas_config: GasConfig,
    grid_config: LUTGridConfig,
    *,
    cache_dir: Path | str | None = None,
    progress: bool = True,
) -> tuple[
    Float[np.ndarray, " G"], Float[np.ndarray, "G T P"], Float[np.ndarray, " G"]
]:
    """Evaluate σ(ν, T, P) on the full LUT grid.

    ``σ(ν) = Σ_l S_l(T) f_V(ν − ν_l) · 1[|ν − ν_l| ≤ w]`` over every line in
    the fetched table, each cut at ``w = max(wing_cm, wing_halfwidths ·
    max(Γ_L, Γ_D))`` from its centre (HAPI's rule; ±25 cm⁻¹ by default).
    Assumes ``fetch_hitran_data`` has already populated the cache, over the
    grid plus ``w`` so lines just outside it add their wings. Fails
    fast — a HAPI failure at any ``(T, P)`` knot, or any non-finite
    cross-section in the result, raises rather than leaving silent NaN
    holes in the LUT.

    Args:
        gas_config:   Target gas (molecule + isotopologues + ν range).
        grid_config:  Common T/P/ν-step grid and line cut-off.
        cache_dir:    Where the fetched HITRAN ``.data/.header`` pair lives.
        progress:     Emit per-(T, P) log lines at INFO.

    Raises:
        RuntimeError: if ``absorptionCoefficient_Voigt`` raises at any
            grid knot, or if the returned cross-section contains any
            non-finite values.

    Returns:
        ``(nu_grid, sigma_lut, wavelength_nm)``:
        - ``nu_grid``: wavenumbers [cm⁻¹], shape ``(n_nu,)``
        - ``sigma_lut``: cross-sections [cm²/molecule], shape ``(n_nu, n_T, n_P)``.
        - ``wavelength_nm``: wavelengths [nm], shape ``(n_nu,)``.
    """
    from hapi import absorptionCoefficient_Voigt

    if cache_dir is not None:
        _init_hapi_cache(Path(cache_dir))

    nu_grid = np.arange(gas_config.nu_min, gas_config.nu_max, grid_config.nu_step)
    wavelength_nm = 1e7 / nu_grid

    n_T = len(grid_config.T_grid)
    n_P = len(grid_config.P_grid)
    sigma_lut = np.zeros((nu_grid.size, n_T, n_P), dtype=np.float64)

    vmr_nominal = DEFAULT_VMR_NOMINAL.get(gas_config.name, 1e-6)
    diluent = grid_config.get_diluent_for_gas(gas_config.name, vmr_nominal)

    total = n_T * n_P
    for i_T, T in enumerate(grid_config.T_grid):
        for i_P, P in enumerate(grid_config.P_grid):
            if progress:
                logger.info(
                    "  [%d/%d] %s T=%.0fK P=%.2fatm",
                    i_T * n_P + i_P + 1,
                    total,
                    gas_config.name,
                    T,
                    P,
                )
            try:
                _, coef = absorptionCoefficient_Voigt(
                    SourceTables=gas_config.name,
                    WavenumberGrid=nu_grid,
                    Environment={"T": float(T), "p": float(P)},
                    Diluent=diluent,
                    HITRAN_units=True,  # cm²/molecule
                    WavenumberWing=float(grid_config.wing_cm),
                    WavenumberWingHW=float(grid_config.wing_halfwidths),
                )
            except Exception as exc:
                raise TableBuildError(
                    f"absorptionCoefficient_Voigt failed for {gas_config.name} "
                    f"at T={T:.1f} K, p={P:.3f} atm."
                ) from exc
            sigma_lut[:, i_T, i_P] = coef

    if not np.all(np.isfinite(sigma_lut)):
        n_bad = int(np.sum(~np.isfinite(sigma_lut)))
        raise TableBuildError(
            f"compute_absorption_lut produced {n_bad} non-finite entries for "
            f"{gas_config.name}; refusing to return a partial LUT."
        )

    return nu_grid, sigma_lut, wavelength_nm

build_lut_dataset

build_lut_dataset(gas_config: GasConfig, grid_config: LUTGridConfig, nu_grid: Float[ndarray, ' G'], sigma_lut: Float[ndarray, 'G T P'], wavelength_nm: Float[ndarray, ' G']) -> Dataset

Wrap the raw LUT arrays in a CF-conformant xarray.Dataset.

The attributes record the line set and its cut-off: isotopologue_ids (HITRAN global ids, from :meth:~radtrax.spectroscopy.hitran.GasConfig.resolve_isotopologue_ids, so it needs [hitran]), line_wing_cm and line_wing_halfwidths.

Source code in src/radtrax/spectroscopy/hitran/generator.py

def build_lut_dataset(
    gas_config: GasConfig,
    grid_config: LUTGridConfig,
    nu_grid: Float[np.ndarray, " G"],
    sigma_lut: Float[np.ndarray, "G T P"],
    wavelength_nm: Float[np.ndarray, " G"],
) -> xr.Dataset:
    """Wrap the raw LUT arrays in a CF-conformant ``xarray.Dataset``.

    The attributes record the line set and its cut-off: ``isotopologue_ids``
    (HITRAN global ids, from
    :meth:`~radtrax.spectroscopy.hitran.GasConfig.resolve_isotopologue_ids`,
    so it needs ``[hitran]``), ``line_wing_cm`` and ``line_wing_halfwidths``.
    """
    now = datetime.now(UTC).isoformat()
    ds = xr.Dataset(
        data_vars={
            "absorption_cross_section": (
                ["wavenumber", "temperature", "pressure"],
                sigma_lut,
                {
                    "units": "cm^2 / molecule",
                    "long_name": f"{gas_config.name} absorption cross-section",
                    "description": (
                        "Voigt-profile line-by-line absorption cross-section "
                        "from HITRAN, computed with HAPI."
                    ),
                    "standard_name": (
                        f"absorption_cross_section_of_{gas_config.name.lower()}_in_air"
                    ),
                },
            ),
        },
        coords={
            "wavenumber": (
                ["wavenumber"],
                nu_grid,
                {
                    "units": "cm^-1",
                    "long_name": "Wavenumber",
                    "standard_name": "wavenumber",
                },
            ),
            "wavelength": (
                ["wavenumber"],
                wavelength_nm,
                {
                    "units": "nm",
                    "long_name": "Wavelength",
                    "standard_name": "radiation_wavelength",
                },
            ),
            "temperature": (
                ["temperature"],
                np.asarray(grid_config.T_grid, dtype=float),
                {
                    "units": "K",
                    "long_name": "Temperature",
                    "standard_name": "air_temperature",
                },
            ),
            "pressure": (
                ["pressure"],
                np.asarray(grid_config.P_grid, dtype=float),
                {
                    "units": "atm",
                    "long_name": "Pressure",
                    "standard_name": "air_pressure",
                },
            ),
        },
        attrs={
            "title": f"{gas_config.name} Absorption Cross-Section Look-Up Table",
            "institution": "Generated using HITRAN API (HAPI)",
            "source": "HITRAN molecular spectroscopic database",
            "molecule": gas_config.name,
            "molecule_id": gas_config.molecule_id,
            "isotopologue_ids": list(gas_config.resolve_isotopologue_ids()),
            "line_wing_cm": float(grid_config.wing_cm),
            "line_wing_halfwidths": float(grid_config.wing_halfwidths),
            "spectral_range_cm_1": f"{nu_grid.min():.1f}-{nu_grid.max():.1f}",
            "spectral_range_nm": f"{wavelength_nm[-1]:.1f}-{wavelength_nm[0]:.1f}",
            "temperature_range_K": f"{float(grid_config.T_grid.min()):.0f}-{float(grid_config.T_grid.max()):.0f}",
            "pressure_range_atm": f"{float(grid_config.P_grid.min()):.2f}-{float(grid_config.P_grid.max()):.2f}",
            "line_shape": "Voigt profile",
            "description": gas_config.description,
            "conventions": "CF-1.8",
            "creation_date": now,
        },
    )
    return ds

save_lut

save_lut(ds: Dataset, output_dir: Path | str, gas_name: str, *, filename: str | None = None) -> Path

Persist ds to <output_dir>/<filename> (zlib, σ as float32).

filename defaults to <gas_name>_absorption_lut.nc. The table is written to a temporary file in output_dir and renamed into place, so a crashed or concurrent write never leaves a half table under the name.

Source code in src/radtrax/spectroscopy/hitran/generator.py

def save_lut(
    ds: xr.Dataset,
    output_dir: Path | str,
    gas_name: str,
    *,
    filename: str | None = None,
) -> Path:
    """Persist ``ds`` to ``<output_dir>/<filename>`` (zlib, σ as float32).

    ``filename`` defaults to ``<gas_name>_absorption_lut.nc``. The table is
    written to a temporary file in ``output_dir`` and renamed into place, so
    a crashed or concurrent write never leaves a half table under the name.
    """
    out = Path(output_dir)
    out.mkdir(parents=True, exist_ok=True)
    path = out / (filename or f"{gas_name.lower()}_absorption_lut.nc")
    with tempfile.NamedTemporaryFile(
        dir=out, prefix=f".{path.name}.", suffix=".partial", delete=False
    ) as tmp:
        partial = Path(tmp.name)
    try:
        ds.to_netcdf(
            partial,
            encoding={
                "absorption_cross_section": {
                    "zlib": True,
                    "complevel": 4,
                    "dtype": "float32",
                }
            },
        )
        partial.replace(path)
    finally:
        partial.unlink(missing_ok=True)
    logger.info(
        "Saved %s LUT to %s (%.1f MB)", gas_name, path, path.stat().st_size / 1e6
    )
    return path

generate_single_gas_lut

generate_single_gas_lut(gas_config: GasConfig, grid_config: LUTGridConfig | None = None, *, cache_dir: Path | str | None = None, output_dir: Path | str | None = None, save: bool = False) -> Dataset

End-to-end single-gas pipeline: fetch → compute → wrap → (optionally) save.

Source code in src/radtrax/spectroscopy/hitran/generator.py

def generate_single_gas_lut(
    gas_config: GasConfig,
    grid_config: LUTGridConfig | None = None,
    *,
    cache_dir: Path | str | None = None,
    output_dir: Path | str | None = None,
    save: bool = False,
) -> xr.Dataset:
    """End-to-end single-gas pipeline: fetch → compute → wrap → (optionally) save."""
    grid_config = grid_config if grid_config is not None else LUTGridConfig()
    cache = Path(cache_dir) if cache_dir is not None else default_cache_dir()

    fetch_hitran_data(gas_config, cache_dir=cache)
    nu_grid, sigma_lut, wavelength_nm = compute_absorption_lut(
        gas_config, grid_config, cache_dir=cache
    )
    ds = build_lut_dataset(gas_config, grid_config, nu_grid, sigma_lut, wavelength_nm)
    if save:
        if output_dir is None:
            raise ParameterError("save=True requires output_dir.")
        save_lut(ds, output_dir, gas_config.name)
    return ds

build

Build v0's cross-section tables into the cache (the [hitran] extra).

:func:build_table runs the existing pipeline (:func:~radtrax.spectroscopy.hitran.generate_single_gas_lut: fetch the lines from HITRAN, compute σ with HAPI, wrap them) on a :class:~radtrax.spectroscopy.table_cache.TableGrid, over the grid's wavenumber range for every gas so their tables share one grid, and writes the result where :func:~radtrax.spectroscopy.table_cache.load_table looks. python -m radtrax.spectroscopy.hitran GAS [GAS ...] is its command line.

gas_config_for

gas_config_for(gas: str, grid: TableGrid = V0_GRID) -> GasConfig

gas's HITRAN configuration over grid's wavenumber range.

Every gas is built on the grid's range, not its own, so the tables share one wavenumber grid (:class:~radtrax.spectroscopy.TableOpacity needs it).

Raises:

Type	Description
KeyError	If gas has no HITRAN configuration.

Source code in src/radtrax/spectroscopy/hitran/build.py

def gas_config_for(gas: str, grid: TableGrid = V0_GRID) -> GasConfig:
    """``gas``'s HITRAN configuration over ``grid``'s wavenumber range.

    Every gas is built on the grid's range, not its own, so the tables share
    one wavenumber grid (:class:`~radtrax.spectroscopy.TableOpacity` needs it).

    Raises:
        KeyError: If ``gas`` has no HITRAN configuration.
    """
    if gas not in ATMOSPHERIC_GASES:
        raise MissingKeyError(
            f"no HITRAN configuration for {gas!r} (have {sorted(ATMOSPHERIC_GASES)})."
        )
    return dataclasses.replace(
        ATMOSPHERIC_GASES[gas], nu_min=grid.nu_min_cm, nu_max=grid.nu_max_cm
    )

fetch_lines

fetch_lines(gas: str, grid: TableGrid = V0_GRID, *, line_cache_dir: Path | str | None = None) -> None

Download gas's HITRAN lines over the grid plus a margin ([hitran]).

Source code in src/radtrax/spectroscopy/hitran/build.py

def fetch_lines(
    gas: str, grid: TableGrid = V0_GRID, *, line_cache_dir: Path | str | None = None
) -> None:
    """Download ``gas``'s HITRAN lines over the grid plus a margin (``[hitran]``)."""
    config = gas_config_for(gas, grid)
    padded = dataclasses.replace(
        config,
        nu_min=config.nu_min - LINE_MARGIN_CM,
        nu_max=config.nu_max + LINE_MARGIN_CM,
    )
    _with_extra(fetch_hitran_data, padded, cache_dir=line_cache_dir)

build_table

build_table(gas: str, grid: TableGrid = V0_GRID, *, hitran_version: str = HITRAN_VERSION, cache_dir: Path | str | None = None, line_cache_dir: Path | str | None = None, overwrite: bool = False) -> Path

Build gas's table on grid and write it to the table cache.

Parameters:

Name	Type	Description	Default
gas	str	A gas of :data:~radtrax.spectroscopy.hitran.ATMOSPHERIC_GASES.	required
grid	TableGrid	The table grid.	V0_GRID
hitran_version	str	Recorded in the file and its cache key.	HITRAN_VERSION
cache_dir	Path | str | None	Table cache directory (default :func:~radtrax.spectroscopy.table_cache.table_cache_dir).	None
line_cache_dir	Path | str | None	Where HAPI keeps the downloaded line lists (default :func:~radtrax.spectroscopy.hitran.default_cache_dir).	None
overwrite	bool	Rebuild a table that already exists.	False

Returns:

Type	Description
Path	The table's path.

Raises:

Type	Description
KeyError	If gas has no HITRAN configuration.
ValueError	If the line cut-off exceeds :data:LINE_MARGIN_CM, so lines outside the fetched range would be silently dropped.
ImportError	Without the [hitran] extra, naming it.

Source code in src/radtrax/spectroscopy/hitran/build.py

def build_table(
    gas: str,
    grid: TableGrid = V0_GRID,
    *,
    hitran_version: str = HITRAN_VERSION,
    cache_dir: Path | str | None = None,
    line_cache_dir: Path | str | None = None,
    overwrite: bool = False,
) -> Path:
    """Build ``gas``'s table on ``grid`` and write it to the table cache.

    Args:
        gas: A gas of :data:`~radtrax.spectroscopy.hitran.ATMOSPHERIC_GASES`.
        grid: The table grid.
        hitran_version: Recorded in the file and its cache key.
        cache_dir: Table cache directory (default
            :func:`~radtrax.spectroscopy.table_cache.table_cache_dir`).
        line_cache_dir: Where HAPI keeps the downloaded line lists (default
            :func:`~radtrax.spectroscopy.hitran.default_cache_dir`).
        overwrite: Rebuild a table that already exists.

    Returns:
        The table's path.

    Raises:
        KeyError: If ``gas`` has no HITRAN configuration.
        ValueError: If the line cut-off exceeds :data:`LINE_MARGIN_CM`, so
            lines outside the fetched range would be silently dropped.
        ImportError: Without the ``[hitran]`` extra, naming it.
    """
    gas_config = gas_config_for(gas, grid)
    lut_grid = LUTGridConfig(
        T_grid=np.asarray(grid.t_k, dtype=float),
        P_grid=np.asarray(grid.p_atm, dtype=float),
        nu_step=grid.nu_step_cm,
    )
    if lut_grid.wing_cm > LINE_MARGIN_CM:
        raise ParameterError(
            f"line cut-off {lut_grid.wing_cm} cm⁻¹ exceeds the fetch margin "
            f"LINE_MARGIN_CM = {LINE_MARGIN_CM} cm⁻¹."
        )
    path = table_path(gas, grid, hitran_version=hitran_version, cache_dir=cache_dir)
    if path.is_file() and not overwrite:
        logger.info("%s table already built at %s", gas, path)
        return path
    fetch_lines(gas, grid, line_cache_dir=line_cache_dir)
    ds = _with_extra(
        generate_single_gas_lut, gas_config, lut_grid, cache_dir=line_cache_dir
    )
    ds.attrs.update(
        radtrax_grid=grid.name,
        radtrax_grid_key=grid.key(),
        hitran_version=hitran_version,
        radtrax_table_version=TABLE_VERSION,
    )
    save_lut(ds, path.parent, gas, filename=path.name)
    logger.info("Built %s table at %s (%.1f MB)", gas, path, path.stat().st_size / 1e6)
    return path

main

main(argv: Sequence[str] | None = None) -> int

python -m radtrax.spectroscopy.hitran GAS [GAS ...] [--grid v0].

Source code in src/radtrax/spectroscopy/hitran/build.py

def main(argv: Sequence[str] | None = None) -> int:
    """``python -m radtrax.spectroscopy.hitran GAS [GAS ...] [--grid v0]``."""
    parser = argparse.ArgumentParser(
        prog="python -m radtrax.spectroscopy.hitran",
        description="Build radtrax cross-section tables from HITRAN into the cache.",
    )
    parser.add_argument("gases", nargs="+", choices=sorted(ATMOSPHERIC_GASES))
    parser.add_argument("--grid", default=V0_GRID.name, choices=sorted(GRIDS))
    parser.add_argument("--cache-dir", type=Path, default=None)
    parser.add_argument("--line-cache-dir", type=Path, default=None)
    parser.add_argument("--hitran-version", default=HITRAN_VERSION)
    parser.add_argument("--overwrite", action="store_true")
    args = parser.parse_args(argv)
    logging.basicConfig(level=logging.INFO, format="%(message)s")
    for gas in args.gases:
        path = build_table(
            gas,
            GRIDS[args.grid],
            hitran_version=args.hitran_version,
            cache_dir=args.cache_dir,
            line_cache_dir=args.line_cache_dir,
            overwrite=args.overwrite,
        )
        print(path)
    return 0

multi

Multi-gas LUT orchestration — separate-per-gas and combined builds.

Ported from plumax 0.1.2 at 28f257e (src/plumax/hapi_lut/multi.py); only import paths changed.

create_multi_gas_luts

create_multi_gas_luts(gases: list[str], grid_config: LUTGridConfig | None = None, *, output_dir: Path | str, cache_dir: Path | str | None = None, force_recompute: bool = False) -> dict[str, Path]

Build one separate <gas>_absorption_lut.nc per gas.

Gases are looked up in ATMOSPHERIC_GASES. Unknown names are skipped with a warning rather than raising — lets users pass a superset list.

Parameters:

Name	Type	Description	Default
gases	list[str]	List of gas names (e.g. ['CH4', 'CO2', 'H2O']).	required
grid_config	LUTGridConfig | None	Common grid; defaults to LUTGridConfig().	None
output_dir	Path | str	Required. Where to write NetCDF LUTs. No default is provided to avoid dropping unignored artefacts into the current working directory; callers should point this at the project-local git-ignored projects/plumax/data/hapi_lut/ path (or an out-of-tree scratch directory).	required
cache_dir	Path | str | None	HITRAN cache dir; defaults to default_cache_dir().	None
force_recompute	bool	If False, skip gases whose NetCDF already exists.	False

Returns:

Type	Description
dict[str, Path]	Mapping {gas_name: Path} for successfully-built LUTs.

Source code in src/radtrax/spectroscopy/hitran/multi.py

def create_multi_gas_luts(
    gases: list[str],
    grid_config: LUTGridConfig | None = None,
    *,
    output_dir: Path | str,
    cache_dir: Path | str | None = None,
    force_recompute: bool = False,
) -> dict[str, Path]:
    """Build one separate ``<gas>_absorption_lut.nc`` per gas.

    Gases are looked up in ``ATMOSPHERIC_GASES``. Unknown names are skipped
    with a warning rather than raising — lets users pass a superset list.

    Args:
        gases:           List of gas names (e.g. ``['CH4', 'CO2', 'H2O']``).
        grid_config:     Common grid; defaults to ``LUTGridConfig()``.
        output_dir:      Required. Where to write NetCDF LUTs. No default
                         is provided to avoid dropping unignored artefacts
                         into the current working directory; callers should
                         point this at the project-local git-ignored
                         ``projects/plumax/data/hapi_lut/`` path
                         (or an out-of-tree scratch directory).
        cache_dir:       HITRAN cache dir; defaults to ``default_cache_dir()``.
        force_recompute: If False, skip gases whose NetCDF already exists.

    Returns:
        Mapping ``{gas_name: Path}`` for successfully-built LUTs.
    """
    grid_config = grid_config if grid_config is not None else LUTGridConfig()
    cache = Path(cache_dir) if cache_dir is not None else default_cache_dir()
    out = Path(output_dir)
    out.mkdir(parents=True, exist_ok=True)

    lut_files: dict[str, Path] = {}
    for gas_name in gases:
        if gas_name not in ATMOSPHERIC_GASES:
            logger.warning("Unknown gas %r — skipping.", gas_name)
            continue
        gas_config = ATMOSPHERIC_GASES[gas_name]

        path = out / f"{gas_name.lower()}_absorption_lut.nc"
        if path.exists() and not force_recompute:
            logger.info("LUT already exists for %s (%s); skipping.", gas_name, path)
            lut_files[gas_name] = path
            continue

        fetch_hitran_data(gas_config, cache_dir=cache)
        nu_grid, sigma_lut, wavelength_nm = compute_absorption_lut(
            gas_config, grid_config, cache_dir=cache
        )
        ds = build_lut_dataset(
            gas_config, grid_config, nu_grid, sigma_lut, wavelength_nm
        )
        lut_files[gas_name] = save_lut(ds, out, gas_name)

    logger.info("Built %d of %d requested LUTs.", len(lut_files), len(gases))
    return lut_files

create_combined_lut

create_combined_lut(gases: list[str], grid_config: LUTGridConfig | None = None, *, output_dir: Path | str, cache_dir: Path | str | None = None, filename: str = 'combined_atmospheric_lut.nc') -> Path

Build a single NetCDF with one cross-section variable per gas.

The shared wavenumber axis is the intersection of the per-gas ranges: [max(nu_min), min(nu_max)]. All variables therefore live on one (wavenumber, temperature, pressure) coordinate system — convenient for matched-filter retrievals where CH₄ and CO₂ share the same SWIR window.

output_dir has no default for the same reason as :func:create_multi_gas_luts — forcing an explicit path keeps generated artefacts out of the current working directory.

Source code in src/radtrax/spectroscopy/hitran/multi.py

def create_combined_lut(
    gases: list[str],
    grid_config: LUTGridConfig | None = None,
    *,
    output_dir: Path | str,
    cache_dir: Path | str | None = None,
    filename: str = "combined_atmospheric_lut.nc",
) -> Path:
    """Build a single NetCDF with one cross-section variable per gas.

    The shared wavenumber axis is the intersection of the per-gas ranges:
    ``[max(nu_min), min(nu_max)]``. All variables therefore live on one
    ``(wavenumber, temperature, pressure)`` coordinate system — convenient
    for matched-filter retrievals where CH₄ and CO₂ share the same SWIR window.

    ``output_dir`` has no default for the same reason as
    :func:`create_multi_gas_luts` — forcing an explicit path keeps
    generated artefacts out of the current working directory.
    """
    grid_config = grid_config if grid_config is not None else LUTGridConfig()
    cache = Path(cache_dir) if cache_dir is not None else default_cache_dir()
    out = Path(output_dir)
    out.mkdir(parents=True, exist_ok=True)

    missing = [g for g in gases if g not in ATMOSPHERIC_GASES]
    if missing:
        raise MissingKeyError(f"Unknown gases: {missing}")

    nu_min = max(ATMOSPHERIC_GASES[g].nu_min for g in gases)
    nu_max = min(ATMOSPHERIC_GASES[g].nu_max for g in gases)
    if nu_min >= nu_max:
        raise GridError(
            f"Gases {gases} have no common spectral range ({nu_min:.1f} >= {nu_max:.1f})."
        )
    nu_grid = np.arange(nu_min, nu_max, grid_config.nu_step)
    wavelength_nm = 1e7 / nu_grid

    logger.info(
        "Building combined LUT for %s on %.1f-%.1f cm^-1 (%d pts).",
        gases,
        nu_min,
        nu_max,
        nu_grid.size,
    )

    data_vars: dict[str, tuple] = {}
    for gas_name in gases:
        base = ATMOSPHERIC_GASES[gas_name]
        sub = dataclasses.replace(base, nu_min=nu_min, nu_max=nu_max)
        fetch_hitran_data(sub, cache_dir=cache)
        _, sigma_lut, _ = compute_absorption_lut(sub, grid_config, cache_dir=cache)
        data_vars[f"{gas_name}_cross_section"] = (
            ["wavenumber", "temperature", "pressure"],
            sigma_lut,
            {
                "units": "cm^2 / molecule",
                "long_name": f"{gas_name} absorption cross-section",
                "molecule_id": base.molecule_id,
                "isotopologue_ids": list(sub.resolve_isotopologue_ids()),
                "line_wing_cm": float(grid_config.wing_cm),
                "line_wing_halfwidths": float(grid_config.wing_halfwidths),
            },
        )

    ds = xr.Dataset(
        data_vars=data_vars,
        coords={
            "wavenumber": (["wavenumber"], nu_grid, {"units": "cm^-1"}),
            "wavelength": (["wavenumber"], wavelength_nm, {"units": "nm"}),
            "temperature": (
                ["temperature"],
                np.asarray(grid_config.T_grid, dtype=float),
                {"units": "K"},
            ),
            "pressure": (
                ["pressure"],
                np.asarray(grid_config.P_grid, dtype=float),
                {"units": "atm"},
            ),
        },
        attrs={
            "title": "Combined Atmospheric Absorption LUT",
            "gases": ", ".join(gases),
            "line_shape": "Voigt profile",
            "conventions": "CF-1.8",
        },
    )

    path = out / filename
    ds.to_netcdf(
        path,
        encoding={
            var: {"zlib": True, "complevel": 4, "dtype": "float32"} for var in data_vars
        },
    )
    logger.info("Saved combined LUT to %s (%.1f MB)", path, path.stat().st_size / 1e6)
    return path