Phase-function treatment

Legendre moment convention

All phase functions are specified through their Legendre expansion moments \(\chi_\ell\), defined by

(7)\[p(\cos\Theta) = \sum_{\ell=0}^{L-1} (2\ell+1)\,\chi_\ell\, P_\ell(\cos\Theta),\]

where \(P_\ell\) are the Legendre polynomials. The normalisation condition is \(\chi_0 = 1\). The \((2\ell+1)\) weighting is part of the sum; the stored moments \(\chi_\ell\) are the reduced expansion coefficients (so, for example, the asymmetry parameter is \(g = \chi_1\)).

Supported phase functions

The configuration provides helpers for the common cases:

Phase function	Reduced moments
Isotropic	\(\chi_0 = 1\), \(\chi_\ell = 0\) for \(\ell \geq 1\).
Rayleigh	\(\chi_0 = 1\), \(\chi_1 = 0\), \(\chi_2 = 1/10\), \(\chi_\ell = 0\) for \(\ell > 2\). (The familiar \(\tfrac12\) in the \(P_2\) term comes from \(5\chi_2 = \tfrac12\).)
Henyey–Greenstein	\(\chi_\ell = g^\ell\).
Double Henyey–Greenstein	\(\chi_\ell = f_w\, g_1^\ell + (1-f_w)\, g_2^\ell\).

Arbitrary phase functions can be used by filling the Legendre moment arrays in the configuration directly. See Phase functions for the API.

Azimuthally averaged phase matrices

The \(m = 0\) phase matrices are built from the Legendre addition theorem:

(8)\[\begin{split}P^{++}(i,j) &= \frac{1}{2\pi}\sum_{\ell=0}^{L-1} (2\ell+1)\,\chi_\ell\, P_\ell(\mu_i)\,P_\ell(\mu_j), \\ P^{+-}(i,j) &= \frac{1}{2\pi}\sum_{\ell=0}^{L-1} (2\ell+1)\,(-1)^\ell\,\chi_\ell\, P_\ell(\mu_i)\,P_\ell(\mu_j).\end{split}\]

\(P^{++}\) describes scattering within the same hemisphere (forward) and \(P^{+-}\) into the opposite hemisphere (backward); the \((-1)^\ell\) factor implements \(P_\ell(-\mu_j) = (-1)^\ell P_\ell(\mu_j)\).

Hansen normalisation is then applied column-wise: each column \(j\) is rescaled so that

(9)\[\sum_{i=1}^{M} \bigl[P^{++}(i,j) + P^{+-}(i,j)\bigr]\,w_i = \frac{1}{2\pi}.\]

This guarantees exact conservation of scattered energy on the discrete quadrature, which is essential for the doubling recursion to conserve flux. The Legendre polynomials \(P_\ell(\mu_i)\) are precomputed with the standard three-term recurrence and cached.

Solar phase vectors

For a direct beam at \(\mu_0\), the forward- and backward-scattering phase vectors are computed analogously:

(10)\[\begin{split}p^{+}_i &= \frac{1}{2\pi}\sum_{\ell=0}^{L-1} (2\ell+1)\,\chi_\ell\, P_\ell(\mu_i)\,P_\ell(\mu_0), \\ p^{-}_i &= \frac{1}{2\pi}\sum_{\ell=0}^{L-1} (2\ell+1)\,(-1)^\ell\,\chi_\ell\, P_\ell(\mu_i)\,P_\ell(\mu_0),\end{split}\]

with the normalisation \(\sum_i (p^+_i + p^-_i)\,w_i = 1/(2\pi)\).