ExamplesInversions

Inversion Examples

Spectral retrieval (field spectrometer)

from biosnicar.emulator import Emulator
from biosnicar.inverse import retrieve
 
emu = Emulator.load("data/emulators/glacier_ice_8_param_default.npz")
 
# measured_albedo = your 480-band spectral albedo measurement
result = retrieve(
    observed=measured_albedo,
    parameters=["ssa", "black_carbon", "glacier_algae", "dust"],
    emulator=emu,
    fixed_params={"direct": 1, "solzen": 50, "snow_algae": 0},
)
print(result.summary())
print(f"SSA = {result.best_fit['ssa']:.2f} ± {result.uncertainty['ssa']:.2f} m²/kg")

True vs retrieved spectrum

Satellite band retrieval

import numpy as np
 
result = retrieve(
    observed=np.array([0.82, 0.78, 0.75, 0.45, 0.03]),
    parameters=["ssa", "black_carbon", "glacier_algae"],
    emulator=emu,
    platform="sentinel2",
    observed_band_names=["B2", "B3", "B4", "B8", "B11"],
    obs_uncertainty=np.array([0.02, 0.02, 0.02, 0.03, 0.05]),
    fixed_params={"direct": 1, "solzen": 50, "dust": 1000, "snow_algae": 0},
)
print(result.summary())

Retrieval with known density

When density is known from field measurements, retrieve rds directly:

result = retrieve(
    observed=measured_albedo,
    parameters=["rds", "black_carbon", "glacier_algae"],
    emulator=emu,
    fixed_params={"direct": 1, "solzen": 50, "rho": 750, "dust": 500, "snow_algae": 0},
)

Comparing optimisation methods

for method in ["L-BFGS-B", "Nelder-Mead", "differential_evolution"]:
    result = retrieve(
        observed=measured_albedo,
        parameters=["ssa", "glacier_algae"],
        emulator=emu,
        method=method,
        fixed_params={"direct": 1, "solzen": 50, "black_carbon": 0,
                      "dust": 0, "snow_algae": 0},
    )
    print(f"{method:25s}  SSA={result.best_fit['ssa']:.4f}  "
          f"cost={result.cost:.6f}  evals={result.n_function_evals}")

MCMC for publication-quality uncertainty

result = retrieve(
    observed=measured_albedo,
    parameters=["ssa", "black_carbon", "glacier_algae"],
    emulator=emu,
    method="mcmc",
    mcmc_walkers=32,
    mcmc_steps=5000,
    mcmc_burn=1000,
    fixed_params={"direct": 1, "solzen": 50, "snow_algae": 0},
)
print(result.summary())
 
# Full chains for corner plots
flat_chain = result.chains.reshape(-1, 3)
# import corner; corner.corner(flat_chain, labels=["SSA", "BC", "algae"])

With regularization

Add prior information when available:

result = retrieve(
    observed=measured_albedo,
    parameters=["ssa", "black_carbon", "glacier_algae", "dust"],
    emulator=emu,
    regularization={"dust": (500, 200)},  # prior: 500 ± 200 ppb
    fixed_params={"direct": 1, "solzen": 50, "snow_algae": 0},
)