Spectral Loop Functions

Radiative transfer calculations typically loop over many wavenumber bands. Calling the solver repeatedly from Python incurs per-call overhead (GIL acquire/release, type conversions). The spectral loop functions move the entire wavenumber iteration into C++, releasing the GIL once and running all bands natively.

When OpenMP is available, the spectral loops are parallelised across threads, with each thread using its own solver and configuration instance. The number of threads can be controlled via the OMP_NUM_THREADS environment variable.

Per-wavelength overrides 

All four spectral loop functions accept optional keyword arguments that supply per-wavelength values for quantities that commonly depend on wavelength. When an override array is provided, it replaces the corresponding value in the configuration for each spectral point. When omitted (the default), the value from the base configuration is used for all wavelengths.

Available overrides
Keyword	Shape	Description
`delta_tau`	`[n_wl, num_layers]`	Optical depth per layer
`single_scat_albedo`	`[n_wl, num_layers]`	Single-scattering albedo per layer
`phase_function_moments`	`[n_wl, num_layers, nmom+1]`	Phase function Legendre moments per layer
`surface_albedo`	`[n_wl]`	Surface albedo
`emissivity_top`	`[n_wl]`	Top-of-atmosphere emissivity
`direct_beam_flux`	`[n_wl]`	Direct beam flux at the top boundary
`isotropic_flux_top`	`[n_wl]`	Isotropic flux at the top boundary
`isotropic_flux_bottom`	`[n_wl]`	Isotropic flux at the bottom boundary

Here n_wl is the number of wavenumber points (or bands). The outer index always corresponds to the spectral point.

Note

The override arrays are shared read-only across OpenMP threads. Each thread applies the overrides to its own local copy of the configuration, so thread safety is guaranteed.

General solver 

solve_spectral

disortpp.solve_spectral(config, wavenumbers, *, delta_tau=[], single_scat_albedo=[], phase_function_moments=[], surface_albedo=[], emissivity_top=[], direct_beam_flux=[], isotropic_flux_top=[], isotropic_flux_bottom=[])

Solve for a list of single wavenumbers using DisortSolver.

For each wavenumber, wavenumber_low and wavenumber_high are both set to the given value. This is appropriate for monochromatic calculations where thermal emission is not used.

Parameters:

config – DisortConfig instance. A copy is made internally; the original is not modified.
wavenumbers (list[float]) – List of wavenumber values [cm^-1].
delta_tau (list[list[float]]) – Per-wavelength optical depths [n_wl][num_layers]. Empty (default) uses the config value for all wavelengths.
single_scat_albedo (list[list[float]]) – Per-wavelength SSA [n_wl][num_layers].
phase_function_moments (list[list[list[float]]]) – Per-wavelength moments [n_wl][num_layers][nmom+1].
surface_albedo (list[float]) – Per-wavelength surface albedo [n_wl].
emissivity_top (list[float]) – Per-wavelength top emissivity [n_wl].
direct_beam_flux (list[float]) – Per-wavelength beam flux [n_wl].
isotropic_flux_top (list[float]) – Per-wavelength isotropic top flux [n_wl].
isotropic_flux_bottom (list[float]) – Per-wavelength isotropic bottom flux [n_wl].

Returns:

List of DisortResult, one per wavenumber.

Return type:

list[DisortResult]

results = disortpp.solve_spectral(cfg, [5000.0, 10000.0, 15000.0])
flux_up_toa = [r.flux_up[0] for r in results]

solve_spectral_bands

disortpp.solve_spectral_bands(config, wavenumber_bands, *, delta_tau=[], single_scat_albedo=[], phase_function_moments=[], surface_albedo=[], emissivity_top=[], direct_beam_flux=[], isotropic_flux_top=[], isotropic_flux_bottom=[])

Solve for a list of wavenumber bands using DisortSolver.

For each band, wavenumber_low and wavenumber_high are set to the lower and upper edges respectively. This is the appropriate mode for calculations with thermal emission.

Parameters:

config – DisortConfig instance (copied internally).
wavenumber_bands (list[tuple[float, float]]) – List of (low, high) wavenumber pairs [cm^-1].
delta_tau (list[list[float]]) – Per-band optical depths [n_bands][num_layers].
single_scat_albedo (list[list[float]]) – Per-band SSA [n_bands][num_layers].
phase_function_moments (list[list[list[float]]]) – Per-band moments [n_bands][num_layers][nmom+1].
surface_albedo (list[float]) – Per-band surface albedo [n_bands].
emissivity_top (list[float]) – Per-band top emissivity [n_bands].
direct_beam_flux (list[float]) – Per-band beam flux [n_bands].
isotropic_flux_top (list[float]) – Per-band isotropic top flux [n_bands].
isotropic_flux_bottom (list[float]) – Per-band isotropic bottom flux [n_bands].

Returns:

List of DisortResult, one per band.

Return type:

list[DisortResult]

bands = [(100.0, 200.0), (200.0, 300.0), (300.0, 400.0)]
results = disortpp.solve_spectral_bands(cfg, bands)

Flux solver 

solve_flux_spectral

disortpp.solve_flux_spectral(config, wavenumbers, *, delta_tau=[], single_scat_albedo=[], phase_function_moments=[], surface_albedo=[], emissivity_top=[], direct_beam_flux=[], isotropic_flux_top=[], isotropic_flux_bottom=[])

Solve for fluxes at a list of single wavenumbers using the flux solver. Dispatches to the appropriate DisortFluxSolver<NStr> based on config.num_streams.

Parameters:

config – DisortFluxConfig instance.
wavenumbers (list[float]) – List of wavenumber values [cm^-1].
delta_tau (list[list[float]]) – Per-wavelength optical depths [n_wl][num_layers].
single_scat_albedo (list[list[float]]) – Per-wavelength SSA [n_wl][num_layers].
phase_function_moments (list[list[list[float]]]) – Per-wavelength moments [n_wl][num_layers][nmom+1].
surface_albedo (list[float]) – Per-wavelength surface albedo [n_wl].
emissivity_top (list[float]) – Per-wavelength top emissivity [n_wl].
direct_beam_flux (list[float]) – Per-wavelength beam flux [n_wl].
isotropic_flux_top (list[float]) – Per-wavelength isotropic top flux [n_wl].
isotropic_flux_bottom (list[float]) – Per-wavelength isotropic bottom flux [n_wl].

Returns:

List of FluxResult, one per wavenumber.

Return type:

list[FluxResult]

solve_flux_spectral_bands

disortpp.solve_flux_spectral_bands(config, wavenumber_bands, *, delta_tau=[], single_scat_albedo=[], phase_function_moments=[], surface_albedo=[], emissivity_top=[], direct_beam_flux=[], isotropic_flux_top=[], isotropic_flux_bottom=[])

Solve for fluxes at a list of wavenumber bands using the flux solver.

Parameters:

config – DisortFluxConfig instance.
wavenumber_bands (list[tuple[float, float]]) – List of (low, high) wavenumber pairs [cm^-1].
delta_tau (list[list[float]]) – Per-band optical depths [n_bands][num_layers].
single_scat_albedo (list[list[float]]) – Per-band SSA [n_bands][num_layers].
phase_function_moments (list[list[list[float]]]) – Per-band moments [n_bands][num_layers][nmom+1].
surface_albedo (list[float]) – Per-band surface albedo [n_bands].
emissivity_top (list[float]) – Per-band top emissivity [n_bands].
direct_beam_flux (list[float]) – Per-band beam flux [n_bands].
isotropic_flux_top (list[float]) – Per-band isotropic top flux [n_bands].
isotropic_flux_bottom (list[float]) – Per-band isotropic bottom flux [n_bands].

Returns:

List of FluxResult, one per band.

Return type:

list[FluxResult]

Performance considerations 

The C++ spectral loop avoids Python overhead for each wavenumber band. The speedup depends on the number of bands and the per-band computation time:

For small per-band cost (few layers, few streams), the overhead reduction is significant – typical speedups of 2–5x over an equivalent Python loop.
For large per-band cost (many layers, many streams), the overhead is a smaller fraction and the speedup is more modest.
With OpenMP enabled, the spectral loops parallelise across CPU cores, providing near-linear scaling up to the number of available cores.

Example measuring the speedup:

import time

# C++ spectral loop
t0 = time.monotonic()
results_cpp = disortpp.solve_spectral_bands(cfg, bands)
t_cpp = time.monotonic() - t0

# Equivalent Python loop
solver = disortpp.DisortSolver()
t0 = time.monotonic()
results_py = []
for wn_lo, wn_hi in bands:
    cfg.wavenumber_low = wn_lo
    cfg.wavenumber_high = wn_hi
    results_py.append(solver.solve(cfg))
t_py = time.monotonic() - t0

print(f"C++ loop: {t_cpp:.3f} s, Python loop: {t_py:.3f} s")
print(f"Speedup: {t_py / t_cpp:.1f}x")

Spectral Loop Functions

Per-wavelength overrides

General solver

solve_spectral

solve_spectral_bands

Flux solver

solve_flux_spectral

solve_flux_spectral_bands

Performance considerations

Per-wavelength overrides 

General solver 

Flux solver 

Performance considerations 