Differentiable Atmospheric Scattering Model

• The model is a physically-informed framework that computes smooth gradients through transmittance, optical-depth integrals, and phase functions.
• It generalizes classical haze inversion by incorporating spatially varying, geometry-aware, and multiple-scattering formulations to accurately represent non-uniform atmospheres.
• Numerical realizations include pixelwise inversion, PDE-based solvers, Neural ODEs, and deep unfolding, each enhancing dehazing and atmospheric correction performance.

Searching arXiv for papers on differentiable atmospheric scattering, dehazing, and atmospheric transport. A differentiable atmospheric scattering model is an atmospheric image-formation or radiative-transport model whose outputs are treated as smooth functions of physical parameters, latent geometric variables, or learned fields, so that gradients can be propagated through transmittance, optical-depth integrals, phase functions, inverse operators, and numerical solvers. In the recent literature, the term spans several closely related constructions: pixelwise dehazing models derived from the classical atmospheric scattering model, geometry-aware generalizations for non-uniform media, PDE- and ODE-based differentiable inversion operators, and differentiable or differentiability-oriented formulations of full multiple-scattering transport in the atmosphere (Wang et al., 2021, Colombi et al., 2014, Zheng, 10 Jun 2025, Koch et al., 2024, Fang et al., 15 Jun 2026).

1. Foundational formulations

At the image-formation level, the standard atmospheric scattering model writes the observed hazy image as

$I(x,y) = J(x,y)\, t(x,y) + A\,(1 - t(x,y)),$

with

$t(x,y) = e^{-\beta d(x,y)}.$

Here $I$ is the observed hazy image, $J$ the haze-free scene radiance, $t$ the transmission, $d$ the scene depth, $\beta$ the angular scattering coefficient, and $A$ the atmospheric light factor. In this formulation, differentiability is immediate because the forward map is composed of additions, multiplications, and an exponential (Wang et al., 2021).

At the transport level, the same idea appears as an optical-depth or radiative-transfer problem. A Monte Carlo treatment of atmospheric multiple scattering can be viewed as a path-based stochastic solution of the time-dependent radiative transfer equation,

$$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$

where $\sigma_t$ is total extinction, $t(x,y) = e^{-\beta d(x,y)}.$0 the scattering coefficient, $t(x,y) = e^{-\beta d(x,y)}.$1 the phase function, and $t(x,y) = e^{-\beta d(x,y)}.$2 the source term (Colombi et al., 2014). A related planetary rendering formulation uses the phenomenological radiative transfer equation for radiance $t(x,y) = e^{-\beta d(x,y)}.$3 with extinction $t(x,y) = e^{-\beta d(x,y)}.$4, a phase function $t(x,y) = e^{-\beta d(x,y)}.$5, and transmittance

$t(x,y) = e^{-\beta d(x,y)}.$6

again exposing a smooth dependence on extinction profiles and path geometry (Costa et al., 2020).

These two views are not contradictory. The recent literature instead presents a hierarchy: low-dimensional haze equations for inverse imaging, and higher-dimensional radiative-transfer models when multiple scattering, angular structure, or wavelength dependence must be retained. This suggests that “differentiable atmospheric scattering” is best understood as a modeling class rather than a single canonical equation.

2. From homogeneous haze models to spatially varying and geometry-aware models

The classical atmospheric scattering model assumes a global atmospheric light $t(x,y) = e^{-\beta d(x,y)}.$7 and a global scattering coefficient $t(x,y) = e^{-\beta d(x,y)}.$8, so that only depth varies spatially. This assumption is explicitly challenged in work on real-world dehazing, where non-uniform illumination motivates $t(x,y) = e^{-\beta d(x,y)}.$9 and heterogeneous haze motivates $I$0. The resulting fully non-homogeneous atmospheric scattering model is

$I$1

with inverse

$I$2

For each pixel, $I$3 is therefore a smooth function of $I$4, $I$5, and $I$6, and the model admits analytic partial derivatives such as

$I$7

The same work further shows that image error is linear in ALF bias $I$8, but nonlinear in ASC and depth biases $I$9 and $J$0, which motivates asymmetric, state-dependent losses for differentiable learning (Wang et al., 2021).

A distinct generalization is geometry-aware rather than purely pixelwise. For UAV imagery, the UAV Atmospheric Scattering Model introduces flight altitude $J$1, pitch angle $J$2, and altitude-varying extinction. Along a line of sight parameterized by arc length $J$3,

$J$4

and the optical thickness is

$J$5

For ground-intersecting rays this yields the closed form

$J$6

and transmittance

$J$7

The image-formation equation remains

$J$8

but now the non-uniformity is tied to explicit geometry rather than an unconstrained per-pixel field (Fang et al., 15 Jun 2026).

A common misconception is that differentiability requires abandoning physical parameterization. The current literature instead tends to preserve an explicit optical-depth structure—global $J$9 in the classical model, $t$0 in FNH-ASM, or altitude–pitch factorization in UASM—and then differentiates through that structure.

3. Phase functions, vertical profiles, and multiple scattering

In full atmospheric transport, differentiability depends not only on transmittance but also on how densities and phase functions are parameterized. A representative multiple-scattering study models a vertically inhomogeneous, laterally homogeneous atmosphere with two populations—molecules and aerosols—using height-dependent attenuation lengths $t$1 and $t$2 of exponential form. Typical molecular parameters are $t$3 and $t$4, while aerosols use $t$5 and variable $t$6. Because these profiles are exponential in altitude and scale parameters, they are explicitly smooth with respect to $t$7, $t$8, and height (Colombi et al., 2014).

The same study uses the Rayleigh phase function for molecules,

$t$9

and a Double Henyey–Greenstein phase function for aerosols,

$d$0

with asymmetry parameter $d$1 and backward-scattering correction $d$2. In the reported simulations, $d$3 is fixed and $d$4 is varied across regimes loosely associated with aerosols, haze, mist, fog, and rain. The extracted analysis emphasizes that $d$5 is smooth in $d$6, $d$7, and $d$8 except near the formal limit $d$9, making it suitable for gradient-based estimation of aerosol size via $\beta$0 (Colombi et al., 2014).

Differentiability becomes more intricate when multiple scattering and detector observables are included. The atmospheric point spread function is defined through the angular deviation $\beta$1 between direct and indirect photons and a cumulative fraction

$\beta$2

The direct beam obeys Beer–Lambert attenuation,

$\beta$3

which is smooth in the vertical-profile parameters. By contrast, stochastic branching, detector-hit classification, and arrival-time thresholds are not differentiable in the ordinary sense, and the same analysis identifies reparameterization gradients, likelihood-ratio estimators, and adjoint or path-space differentiation as the relevant remedies for differentiable rendering (Colombi et al., 2014).

Related atmospheric transport work extends differentiable structure beyond haze imaging. Finite plane-parallel atmospheres with isotropic scattering and thermal emission are written in optical-depth form, introducing modified reflection and transmission operators $\beta$4 and $\beta$5, together with emission-generalized response functions $\beta$6 and $\beta$7. Their integral equations are smooth in total optical thickness $\beta$8, single-scattering albedo $\beta$9, and boundary emission parameters $A$0 and $A$1, and the extracted analysis explicitly notes that they are suitable for automatic differentiation after angular discretization (Sengupta, 2022).

4. Numerical realizations and differentiable solvers

Recent work implements differentiable atmospheric scattering through several numerical templates. The principal families can be summarized as follows.

Family	Core mechanism	Representative papers
Pixelwise inversion layer	Closed-form haze inversion with elementwise exponentials	(Wang et al., 2021)
PDE-based inversion	Elliptic PDE with differentiable diffusion and Gaussian regularization	(Zheng, 10 Jun 2025)
Neural ODE transmission	Forward and inverse atmospheric transmission operator $A$2 via ODE solve	(Koch et al., 2024)
Deep unfolding	Unrolled proximal gradient descent with physics-based data term	(Fang et al., 15 Jun 2026)

In the CNN-based dehazing setting, FNHD-Net uses a shared U-Net style encoder, three decoders for $A$3, $A$4, and $A$5, and an Estimating Haze-Free Image Module that implements

$A$6

The module is followed by a convolution layer and a BReLU activation, while the losses include FWB-Loss for $A$7, the exponential $A$8-Loss,

$A$9

and the analogous D-Loss for depth. Training uses a two-stage alternating strategy in PyTorch with SGD, learning rate $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$0, weight decay $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$1, and momentum $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$2 (Wang et al., 2021).

A PDE-based alternative defines the dehazed image $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$3 as the weak solution of

$$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$4

with

$$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$5

and

$$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$6

The paper proves existence and uniqueness of weak solutions in $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$7 using the Lax–Milgram theorem, and implements a fixed-point iteration scheme in PyTorch using Scharr filters, Gaussian convolution, and unrolled updates. The extracted discussion emphasizes that the resulting map $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$8 is differentiable almost everywhere and can be backpropagated through as an unrolled optimization layer (Zheng, 10 Jun 2025).

A spectrally resolved approach is DINSAT, which replaces explicit radiative transfer with a differentiable atmospheric transmission operator $$\frac{1}{c}\frac{\partial I(\mathbf{x}, \mathbf{\Omega}, t)}{\partial t} + \mathbf{\Omega}\cdot\nabla I(\mathbf{x}, \mathbf{\Omega}, t) = -\sigma_t(\mathbf{x}) I(\mathbf{x}, \mathbf{\Omega}, t) + \int_{4\pi} \sigma_s(\mathbf{x}) p(\mathbf{\Omega}'\to\mathbf{\Omega}) I(\mathbf{x}, \mathbf{\Omega}', t)\,\mathrm{d}\mathbf{\Omega}' + S(\mathbf{x},\mathbf{\Omega},t),$$9 defined by an ODE,

$\sigma_t$0

solved on $\sigma_t$1. In the linear case,

$\sigma_t$2

while the nonlinear case uses

$\sigma_t$3

with encoder and decoder MLPs. The forward model

$\sigma_t$4

and its inverse

$\sigma_t$5

are both evaluated with a differentiable ODE solver, so gradients can be propagated through forward transmission, inverse transmission, and spectral losses (Koch et al., 2024).

A deep-unfolding formulation appears in GP-DUN for UAV dehazing. Its Latent Geometry Estimator predicts

$\sigma_t$6

where $\sigma_t$7 and $\sigma_t$8. The Geometry-aware Gradient Descent Module then applies the closed-form data-term gradient

$\sigma_t$9

through the stage update

$t(x,y) = e^{-\beta d(x,y)}.$00

The learned proximal step is implemented by the Pooling-Expert Proximal Mapping Module, and the $t(x,y) = e^{-\beta d(x,y)}.$01-stage architecture is trained end-to-end with $t(x,y) = e^{-\beta d(x,y)}.$02 (Fang et al., 15 Jun 2026).

5. Observables, tasks, and representative systems

Differentiable atmospheric scattering models are used to predict or invert different observables depending on domain. In dehazing, the target is the haze-free image $t(x,y) = e^{-\beta d(x,y)}.$03, together with intermediate fields such as $t(x,y) = e^{-\beta d(x,y)}.$04, $t(x,y) = e^{-\beta d(x,y)}.$05, $t(x,y) = e^{-\beta d(x,y)}.$06, transmittance, or geometry proxies. In atmospheric optics and astronomical imaging, the observables include transmission spectra, atmospheric point spread functions, photometry, background, and differential astrometry.

For dehazing, the fully non-homogeneous model reports that FNH-ASM strictly generalizes both ASM and NH-ASM, and that PSNR gains from replacing ASM/NH-ASM by FNH-ASM are approximately $t(x,y) = e^{-\beta d(x,y)}.$07, $t(x,y) = e^{-\beta d(x,y)}.$08, and $t(x,y) = e^{-\beta d(x,y)}.$09 dB for base $t(x,y) = e^{-\beta d(x,y)}.$10 levels $t(x,y) = e^{-\beta d(x,y)}.$11, $t(x,y) = e^{-\beta d(x,y)}.$12, and $t(x,y) = e^{-\beta d(x,y)}.$13, respectively; $t(x,y) = e^{-\beta d(x,y)}.$14-Loss adds approximately $t(x,y) = e^{-\beta d(x,y)}.$15–$t(x,y) = e^{-\beta d(x,y)}.$16 dB and D-Loss approximately $t(x,y) = e^{-\beta d(x,y)}.$17–$t(x,y) = e^{-\beta d(x,y)}.$18 dB. The same work states that performance improvement becomes more obvious in dense and heterogeneous haze scenes (Wang et al., 2021).

For UAV imagery, UASM-HazeSet provides $t(x,y) = e^{-\beta d(x,y)}.$19 paired hazy/clear synthetic images and $t(x,y) = e^{-\beta d(x,y)}.$20 real UAV haze images for testing. The corresponding GP-DUN formulation is explicitly built so that the synthetic data generator and the unfolded data-fidelity term share the same UASM, which the paper argues improves generalization to real UAV haze (Fang et al., 15 Jun 2026).

For hyperspectral atmospheric correction, DINSAT is designed to perform atmospheric correction, recast spectral data between various modalities such as radiance and reflectance at the surface and at the sensor, and infer atmospheric transmission profiles including absorbing bands and their relative magnitudes. On the HyMap Cooke City scene, the reported supervised linear model gives surface reflectance MSE of approximately $t(x,y) = e^{-\beta d(x,y)}.$21 and at-sensor radiance synthesis MSE of approximately $t(x,y) = e^{-\beta d(x,y)}.$22 for a red fabric target, while the supervised nonlinear model gives surface reflectance MSE of approximately $t(x,y) = e^{-\beta d(x,y)}.$23 and at-sensor radiance synthesis MSE of approximately $t(x,y) = e^{-\beta d(x,y)}.$24 (Koch et al., 2024).

For physically based atmospheric blur, the multiple-scattering study emphasizes that aerosols can have a non-negligible effect on the atmospheric point spread function even though they are usually in lower density than molecules, especially for ground-based detectors when large aerosols are present. The reported sensitivities indicate that the indirect-light fraction near the direct path increases with larger $t(x,y) = e^{-\beta d(x,y)}.$25, longer integration time $t(x,y) = e^{-\beta d(x,y)}.$26, and smaller inclination angle $t(x,y) = e^{-\beta d(x,y)}.$27, making PSF measurements informative for aerosol asymmetry and density estimation (Colombi et al., 2014).

Astronomical transport simulations extend the same paradigm to broader observables. A self-consistent atmosphere representation coupled to photon Monte Carlo yields differential astrometry and elliptical point spread functions predicted completely to the horizon, arcminute-scale spatial-dependent photometry variations at $t(x,y) = e^{-\beta d(x,y)}.$28 mmag for short exposures, excess background spatial variations at $t(x,y) = e^{-\beta d(x,y)}.$29 due the atmosphere, and a point spread function wing due to water droplets. The same study reports correlations of altitude with absolute photometry and relative photometry, anti-correlations of altitude with differential astrometry and background levels, and a significant correlation of PSF wing with cloud depth (Peterson et al., 2024).

These examples show that the “output” of a differentiable atmospheric scattering model is not restricted to a restored RGB image. It may instead be a PSF, a transmission curve, a dehazed solution of a PDE, or a hyperspectral inverse estimate, depending on which part of atmospheric transport is treated as primary.

6. Assumptions, misconceptions, and limitations

A recurrent limitation is that many differentiable dehazing models are built on single-scattering or Koschmieder-style assumptions. The PDE-based dehazing framework explicitly assumes a homogeneous participating medium, a single-scattering approximation, and spatially constant atmospheric light $t(x,y) = e^{-\beta d(x,y)}.$30, with transmission map $t(x,y) = e^{-\beta d(x,y)}.$31 estimated by the dark channel prior (Zheng, 10 Jun 2025). Likewise, the classical ASM and its immediate variants do not model explicit phase functions, vertical layering, or multi-angle geometry (Wang et al., 2021).

At the opposite end, full multiple-scattering formulations remain computationally and algorithmically harder to differentiate. The Monte Carlo PSF study identifies the main difficulties as random draws, discrete branching between scatter and no-scatter, aerosol versus molecular interactions, and discontinuous detector classification or arrival-time thresholds; its extracted analysis therefore points to reparameterization, score-function estimators, and adjoint or path-space differentiation as necessary tools (Colombi et al., 2014). A hybrid Monte Carlo/deterministic transport method for a scattering and absorbing atmosphere similarly suggests that deterministic adjoint surrogates can reduce variance, and the extracted discussion explicitly frames this as compatible with differentiable programming, although the original work is not itself a differentiable-rendering paper (Bal et al., 2011).

Even when the forward model is analytic, strict smoothness is not universal. The PDE framework uses $t(x,y) = e^{-\beta d(x,y)}.$32 and clamping, which are only piecewise smooth; the planetary visualization model relies on texture interpolation, precomputed LUTs, and conditional occlusion logic; and the Cerro Paranal sky model includes piecewise spectral fits and discrete seasonal and night-time bins (Zheng, 10 Jun 2025, Costa et al., 2020, Noll et al., 2012). These are not fatal to gradient-based methods, but they constrain where exact derivatives exist and where subgradients or smoothed surrogates may be needed.

Another misconception is that a differentiable model must always be physically complete. The surveyed literature shows the opposite design pattern: models are often deliberately low-dimensional. FNH-ASM keeps only $t(x,y) = e^{-\beta d(x,y)}.$33, $t(x,y) = e^{-\beta d(x,y)}.$34, and $t(x,y) = e^{-\beta d(x,y)}.$35; UASM factorizes transmittance as $t(x,y) = e^{-\beta d(x,y)}.$36; DINSAT collapses transport into a single Beer–Lambert-like ODE layer; and finite-atmosphere scattering-with-emission can be compressed into boundary response functions $t(x,y) = e^{-\beta d(x,y)}.$37, $t(x,y) = e^{-\beta d(x,y)}.$38, $t(x,y) = e^{-\beta d(x,y)}.$39, and $t(x,y) = e^{-\beta d(x,y)}.$40 (Wang et al., 2021, Fang et al., 15 Jun 2026, Koch et al., 2024, Sengupta, 2022). This suggests that differentiability and physical fidelity are typically balanced through carefully chosen state reduction rather than maximal mechanistic detail.

The current literature therefore presents differentiable atmospheric scattering models as a spectrum of physics-informed operators. At one end lie explicit haze inversion layers embedded in CNNs; at the other lie transport solvers with angular, spectral, or vertical structure. What unifies them is not a particular equation, but a shared commitment to parameterized atmospheric physics that remains amenable to gradient-based estimation, optimization, or learning.