CanopyOptics.jl

CanopyOptics

Canopy radiative-transfer utilities.

Canopy scattering matrices use Z[i_out, j_in]: rows are outgoing streams and columns are incoming streams; Z⁺⁺ is same-sign transmission and Z⁻⁺ is sign-change reflection. Bi-Lambertian kernels follow vSmartMOM's atmospheric phase-matrix convention: the scalar single-scattering albedo is not folded into Z, so conservative bi-Lambertian leaves satisfy

\[\sum_i w_i\,(Z^{++}_{ij} + Z^{-+}_{ij}) \simeq 2 .\]

Specular kernels carry their Fresnel/roughness strength in the returned contribution. Layer solvers that multiply by a separate ϖ need a matching effective albedo when using specular components.

Abstract Type for canopy scattering

AbstractClumping{FT}

Angular clumping model for canopy interception.

Clumping modifies the effective projected area used in gap probability and extinction calculations:

\[G_{eff}(μ) = Ω(μ) G(μ).\]

It does not change the leaf or wood scattering matrix Z; it changes how often canopy elements are encountered along a path.

AbstractHotSpot{FT}

Model for bidirectional hotspot corrections to canopy gap probability.

Hotspot terms are propagation corrections, not leaf scattering matrices: they modify the joint probability that the solar and viewing paths see the same canopy gaps. They should therefore be applied where direct-beam source terms are assembled, not inside compute_Z_matrices.

AbstractLUTWoodReflectance{FT}

Base type for wood reflectance models backed by tabulated spectral data.

Abstract Type for leaf distributions

Abstract Material Type (can be all)

Abstract soil type

Abstract vegetation type (leaves, needles, branches, trunks)

Abstract water type

AbstractWoodReflectance{FT}

Spectral reflectance model for bark or woody canopy elements.

Wood reflectance is kept separate from angular scattering. A reflectance model answers "what is the bark albedo at this wavelength or wavenumber?", while a canopy scattering model answers "how does that woody surface redirect light?".

Model for bi-lambertian canopy leaf scattering

CanopyComponent(scatterer, LAD, area_index; clumping = nothing)
CanopyComponent(; scatterer, LAD, area_index, clumping = nothing)

One canopy scattering population.

scatterer is a canopy scattering model such as BiLambertianCanopyScattering or LambertianWoodCanopyScattering, LAD is that population's element angle distribution, and area_index is its one-sided area index: LAI for leaves, WAI/BAI for woody material. Optional per-component clumping is used by bulk_G, but mixed Z-matrix normalization keeps using the unclumped Ross G convention.

CanopyQuadrature(; n_leaf = 64, n_azimuth = 48)

Integration controls used when projecting canopy scattering models onto Fourier Z matrices.

n_leaf controls the leaf-inclination quadrature used by the bi-Lambertian kernel. n_azimuth controls the outgoing-azimuth quadrature used by the specular kernel and legacy brute-force reference paths.

CompositeCanopyScattering(components...)

Additive canopy scattering model. Components are evaluated independently and their Z matrices are summed, so a leaf can carry both diffuse bi-Lambertian and specular surface terms.

ConstantClumping(; Ω = 1)
ConstantClumping(; Ω₀ = 1)

Direction-independent clumping index.

Values Ω < 1 reduce effective extinction relative to a random leaf distribution. Ω = 1 preserves the unclumped Beer-Lambert limit. Values above one represent more regular spacing and are allowed for completeness.

ConstantWoodReflectance(R)
ConstantWoodReflectance(; R = 0.25)

Wavelength-independent bark/wood reflectance.

EmpiricalDirectionalClumping(; Ω₀ = 0.7, c = 2, e = 2)

Simple empirical angular clumping model:

\[Ω(θ) = \frac{Ω_0}{Ω_0 + (1 - Ω_0)\exp(-c θ^e)}, \qquad θ = \arccos μ.\]

The default shape gives Ω(0) = Ω₀ at nadir and approaches one toward the horizon, where clumping is visually less apparent along the slant path.

This is a compact directional parameterization, not the full Chen-Leblanc row/crown gap-size model. ChenLeblancClumping remains as a compatibility alias for older code that used the previous name.

FourSAILGeometry(LAD; sza_deg, vza_deg, raa_deg, quadrature=CanopyQuadrature())

Geometry and leaf-angle coefficients for a 4SAIL canopy solve.

The expensive geometry/LAD integration is wavelength independent and is done once here. foursail! then evaluates one wavelength per kernel work item from leaf reflectance, leaf transmittance, soil reflectance, and LAI.

FourSAILGeometrySet(LAD; sza_deg, vza_deg, raa_deg, quadrature=CanopyQuadrature())
FourSAILGeometrySet(geometries)

Structure-of-arrays storage for many FourSAILGeometry values. This is the preferred input for phase-space/LUT solves because the kernel can evaluate all (wavelength, angle) pairs in one launch.

FourSAILResult(rddt, rsdt, rdot, rsot, taus, rsos, rsost, rsodt)

Primary 4SAIL reflectance/transmittance factors:

• rddt: bi-hemispherical reflectance over canopy plus soil
• rsdt: directional-hemispherical reflectance for solar illumination
• rdot: hemispherical-directional reflectance for the view direction
• rsot: bidirectional reflectance factor
• taus: total solar-direction canopy transmission (tsd + tss)
• rsos: direct single-scattering canopy contribution
• rsost: direct/direct contribution (rsos plus directly transmitted soil)
• rsodt: diffuse/multiple-scattering contribution, so rsot = rsost + rsodt

KuuskHotSpot(; h = 0.1)
KuuskHotSpot{FT}(; h = FT(0.1))

Kuusk-style hotspot correction with dimensionless size parameter h, roughly leaf size divided by canopy height. Larger h broadens the hotspot. h <= 0 disables the correction and returns the independent-gap probability.

LUTWoodReflectance(grid, R; grid_unit = :nm, extrapolation = :clamp)

Piecewise-linear bark/wood reflectance lookup table.

grid and R must have the same length, with at least two points. The constructor sorts the table by spectral coordinate. grid_unit is :nm for wavelength in nanometers or :cm_inv for wavenumber. extrapolation = :clamp uses the nearest endpoint outside the tabulated range; :error throws.

LambertianWoodCanopyScattering(reflectance)
LambertianWoodCanopyScattering(R)
LambertianWoodCanopyScattering(; R = 0.25, reflectance = nothing)

Opaque Lambertian woody element: reflection only, no transmission.

This is the angular-scattering counterpart to AbstractWoodReflectance. For scalar Z-matrix assembly it delegates to BiLambertianCanopyScattering(R, 0), preserving the existing canopy normalization convention.

struct LeafDistribution{FT<:Real}

A struct that defines the leaf angular distribution in radians (from 0->π/2; here scaled to [0,1])

Fields

• LD: Julia Univariate Distribution from Distributions.js

• scaling: Scaling factor to normalize distribution (here mostly 2/π as Beta distribution is from [0,1])

struct PigmentOpticalProperties{FT}

A struct which stores important absorption cross sections of pigments, water, etc

Fields

• λ: Wavelength [length]

• nᵣ: Refractive index of leaf material

• Kcab: specific absorption coefficient of chlorophyll (a+b) [cm² μg⁻¹]

• Kcar: specific absorption coefficient of carotenoids [cm² μg⁻¹]

• Kant: specific absorption coefficient of Anthocyanins [cm² nmol⁻¹]

• Kb: specific absorption coefficient of brown pigments (arbitrary units)

• Kw: specific absorption coefficient of water [cm⁻¹]

• Km: specific absorption coefficient of dry matter [cm² g⁻¹]

• Kp: specific absorption coefficient of proteins [cm² g⁻¹]

• Kcbc: specific absorption coefficient of carbon based constituents [cm² g⁻¹]

struct LeafProperties{FT}

A struct which stores important variables of leaf chemistry and structure

Fields

• N: Leaf structure parameter [0-3]

• Ccab: Chlorophyll a+b content [µg cm⁻²]

• Ccar: Carotenoid content [µg cm⁻²]

• Canth: Anthocynanin content [nmol cm⁻²]

• Cbrown: Brown pigments content in arbitrary units

• Cw: Equivalent water thickness [cm], typical about 0.002-0.015

• Cm: Dry matter content (dry leaf mass per unit area) [g cm⁻²], typical about 0.003-0.016

• Cprot: protein content [g/cm]

• Ccbc: Carbone-based constituents content in [g/cm⁻²] (cellulose, lignin, sugars...)

LeafUlabyElRayes1987(; M_g=0.5, σ=1.27)

Microwave dielectric model for leaves following Ulaby & El-Rayes (1987). A dual-Debye mixing model combining a non-dispersive residual term, a free-water Debye component, and a bound-water Cole–Cole-like component. Valid for fresh leaves at 0.2–20 GHz with gravimetric moisture M_g ∈ [0, 0.7].

Pass to dielectric together with temperature T [K] (unused by this model — kept for signature consistency) and frequency f [GHz].

Fields

• M_g::FT — gravimetric moisture content ∈ [0, 0.7]
• σ::FT — effective ionic conductivity of free water [S/m] (typ. 1.27)

Reference

Ulaby, F. T. & El-Rayes, M. A. (1987). Microwave Dielectric Spectrum of Vegetation – Part II: Dual-Dispersion Model. IEEE TGRS 25(5), 550–557. Reproduced in Ulaby & Long (2014), §11-9.

Pure liquid water

Salty liquid water

MixedCanopy(components...)
MixedCanopy((components...,))

Area-weighted canopy made from multiple scattering populations.

Each component may have its own scatterer, angle distribution, area index, and clumping model. Z matrices are mixed column-by-column using the incoming-stream projected area weights

\[w_c(μ_j) = \frac{AI_c G_c(μ_j)} {\sum_d AI_d G_d(μ_j)}.\]

The optional component clumping models are intentionally not included in this Z normalization; use bulk_G with clumped = true for propagation.

NoClumping()
NoClumping{FT}()

Default clumping model with Ω(μ) = 1.

NoHotSpot()
NoHotSpot{FT}()

Default model with independent sun and view gaps. The joint gap probability is exp(-(k_s + k_o)L).

PolynomialWoodReflectance(coeffs; grid_unit = :nm, x_offset = 0, x_scale = 1)

Polynomial bark/wood reflectance model.

Coefficients use ascending order:

\[R(x) = c_0 + c_1 \hat{x} + c_2 \hat{x}^2 + \cdots, \qquad \hat{x} = (x - x_{offset}) / x_{scale}.\]

Use x_offset and x_scale to fit on a normalized wavelength or wavenumber coordinate instead of raw nanometers.

Pure Ice

Soil MW properties

Model for specular canopy leaf scattering

dirVector{FT}

Struct for spherical coordinate directions in θ (elevation angle) and ϕ (azimuth angle)

Fields

• θ

• ϕ

dirVector_μ{FT}

Struct for spherical coordinate directions in θ (elevation angle) and ϕ (azimuth angle)"

Fields

• μ

• ϕ

BranchComponent(; BAI, R = 0.30, scatterer = nothing, LAD = nothing,
                clumping = nothing) -> CanopyComponent

Convenience constructor for the branch scattering population.

Similar to StemComponent but with a default plagiophile (oblique, on-average ~45°-tilted) angle distribution that represents branches hanging off the main stem at intermediate inclinations. The default R = 0.30 is a broadband typical younger-bark value, slightly higher than the mature-stem default to reflect smoother branch bark; override with any AbstractWoodReflectance for band-resolved spectra.

When scatterer and LAD are left at their defaults, they are constructed in the float type of BAI — so BranchComponent(BAI = 0.1f0) returns a CanopyComponent{Float32, ...} end-to-end.

BAI is the one-sided branch area index [m² branch / m² ground]. Must be non-negative (throws DomainError otherwise).

G(μ::Array{FT}, LD::AbstractLeafDistribution; nLeg=20)

Returns the integrated projection of leaf area in the direction of μ, assumes azimuthally uniform distribution and a LD distribution for leaf polar angle θ. This function is often referred to as the function O(B) (Goudriaan 1977) or G(Ζ) (Ross 1975,1981), see Bonan modeling book, eqs. 14.21-14.26.

Arguments

• μ an array of cos(θ) (directions [0,1])
• LD an AbstractLeafDistribution type struct, includes a leaf distribution function
• nLeg an optional parameter for the number of legendre polynomials to integrate over the leaf distribution (default=20)

Examples

julia> μ,w = CanopyOptics.gauleg(10,0.0,1.0);       # Create 10 quadrature points in μ      
julia> LD  = CanopyOptics.spherical_leaves()        # Create a default spherical leaf distribution
julia> G   = CanopyOptics.G(μ, LD)                  # Compute G(μ)
10-element Vector{Float64}:
 0.5002522783000879
 0.5002715115149204
 0.5003537989277846
 0.5004432798701134
 0.5005134448870893
 0.5003026448466977
 0.4999186257540982
 0.4994511190721635
 0.49907252201082375
 0.49936166823681594

G2(μ::AbstractArray{FT}, LD::AbstractLeafDistribution; nLeg=40)

Reference implementation of G that integrates in leaf-normal cosine μ_L instead of leaf inclination angle θ_L.

This path is retained for cross-checking the production θ_L quadrature in G. New code should normally call G.

LeafComponent(; LAI, scatterer = nothing, LAD = nothing,
              clumping = nothing) -> CanopyComponent

Convenience constructor for the leaf scattering population of a tree canopy. Defaults to a BiLambertianCanopyScattering (reflection

• transmission, the canonical leaf scatterer) and spherical_leaves

(Bonan/Norman default — projection function G(μ) ≈ 0.5). Override either via kwarg: e.g. scatterer = some_prospect_driven_scatterer or LAD = planophile_leaves(FT).

When scatterer and LAD are left at their defaults, they are constructed in the float type of LAI — so LeafComponent(LAI = 3.5f0) returns a CanopyComponent{Float32, ...} end-to-end, with no silent Float64 promotion.

LAI is the one-sided leaf area index [m² leaf / m² ground]. Must be non-negative (throws DomainError otherwise).

The resulting CanopyComponent plugs directly into MixedCanopy and the existing compute_Z_matrices, bulk_G, effective_G paths — no other call sites change.

StemComponent(; SAI, R = 0.25, scatterer = nothing, LAD = nothing,
              clumping = nothing) -> CanopyComponent

Convenience constructor for the stem scattering population.

Stems are modeled as opaque Lambertian (reflection only, no transmission — the "no T, just R" point) with a default erectophile (mostly vertical) angle distribution. This captures the fact that trunks intercept little overhead sun (small G at high μ) but a lot of low-angle sun (large G at small μ).

R may be a scalar (broadband bark albedo) or any AbstractWoodReflectance (e.g. LUTWoodReflectance for band-resolved bark spectra). The default R = 0.25 is a broadband typical mature-bark value (Bonan 2019, Climate Change and Terrestrial Ecosystem Modeling, Table 14.1). Pass scatterer = ... to override the entire scattering model (rarely needed).

When scatterer and LAD are left at their defaults, they are constructed in the float type of SAI — so StemComponent(SAI = 0.9f0) returns a CanopyComponent{Float32, ...} end-to-end.

SAI is the one-sided stem area index [m² stem / m² ground]. Must be non-negative (throws DomainError otherwise).

TreeCanopy(; LAI, SAI = 0, BAI = 0,
           leaf_kwargs = (;), stem_kwargs = (;), branch_kwargs = (;))
    -> MixedCanopy

One-line constructor for "leaves + stems + branches" tree canopies. Returns a MixedCanopy with up to three CanopyComponents.

A SAI or BAI of zero omits the corresponding component entirely — MixedCanopy doesn't carry spurious zero-area populations in its projected-area weighting. When SAI = BAI = 0 the result is a single-component MixedCanopy whose bulk_G matches the leaf-only G(μ, ::LeafDistribution) exactly.

Each population's defaults can be overridden via the *_kwargs NamedTuples, which are spliced into the corresponding LeafComponent / StemComponent / BranchComponent call. For example, to use band-resolved bark reflectance:

using CanopyOptics
canopy = TreeCanopy(;
    LAI = 4.0,
    SAI = 0.9,
    BAI = 0.1,
    stem_kwargs = (; R = LUTWoodReflectance(
        grid = [400, 800, 2500], R = [0.10, 0.20, 0.45], grid_unit = :nm)),
)

The resulting MixedCanopy plugs into all existing paths (bulk_G, compute_Z_matrices, etc.) and is the recommended way to construct tree canopies in upstream consumers (RRTMGP, ClimaLand multi-layer canopies, CanopyColumn).

beta_leaves(α::FT, β::FT) where FT

Creates a LeafDistribution from an arbitrary Beta(α, β) distribution over [0, π/2]. The parameters α and β correspond to ν and μ in the RAMI benchmark convention. Use planophile_leaves, spherical_leaves, etc. for standard named distributions.

bfG(μ::Array{FT}, LD::AbstractLeafDistribution; nLeg=20)

Brute-force two-angle reference for the Ross projection factor G(μ).

This directly integrates abs(Ω ⋅ Ω_L) over leaf inclination and leaf azimuth. It is useful for tests or debugging quadrature changes, but it is intentionally slower than G.

bulk_G(canopy, μ; clumped = true)

Area-index-weighted canopy projected area:

\[G_{bulk}(μ) = \sum_c AI_c G_c(μ),\]

or \sum_c AI_c Ω_c(μ)G_c(μ) when clumped = true.

canopy_extinction(G_value, μ)

Convert projected area G(μ) into the directional extinction coefficient k = G(μ) / μ used by Beer-Lambert canopy gap probabilities.

clumping_index(model, μ)

Evaluate Ω(μ) for one direction cosine or an array of direction cosines.

component_G(component, μ; clumped = false)

Projected area for one CanopyComponent. Set clumped = true to apply the component's clumping model.

compute_Z_matrices(canopy::MixedCanopy, μ, m; quadrature, npol = 1)
compute_Z_matrices(canopy::MixedCanopy, μ, ignored_LD, m; quadrature, npol = 1)

Compute canopy Z matrices for a mixture of components with their own angle distributions and area indices. The ignored_LD compatibility form exists so callers that already carry (scatterer, μ, LAD, m) can pass a MixedCanopy as the scatterer without changing their call shape.

compute_Z_matrices(mod::AbstractCanopyScatteringType,
                   μ, LD, m_range::AbstractUnitRange; quadrature, npol = 1)

Compute a stack of canopy Z matrices for Fourier moments in m_range.

The scalar default returns arrays with shape (length(μ), length(μ), length(m_range)). Passing npol = 3 or 4 returns Stokes-expanded arrays with the first two dimensions multiplied by npol. In both cases the package convention is Z[i_out, j_in, k], where k = m - first(m_range) + 1.

compute_Z_matrices(mod::BiLambertianCanopyScattering,
                   μ::AbstractVector, LD::AbstractLeafDistribution,
                   m::Integer; quadrature = CanopyQuadrature(), npol = 1) -> (Z⁺⁺, Z⁻⁺)

Compute one cosine Fourier moment of the bi-Lambertian canopy scattering matrices. Rows are outgoing streams and columns are incoming streams: Z[i_out, j_in].

npol = 1 returns scalar Stokes-I matrices. npol = 3 and npol = 4 expand the scalar diffuse result into vSmartMOM-style Stokes blocks with only I→I populated.

This method uses the closed-form Fourier moments from compute_Z_matrices_aniso_analytic. For all moments in one call, pass a range such as 0:m_max.

compute_Z_matrices(mod::CompositeCanopyScattering, μ, LD, m; quadrature, npol = 1)

Compute one Fourier moment for an additive canopy-scattering model.

Each component is evaluated independently and the resulting (Z⁺⁺, Z⁻⁺) matrices are summed. Diffuse components use the vSmartMOM convention where scalar single-scattering albedo is divided out of Z; specular components carry their Fresnel/roughness strength in their kernel. Consumers that multiply by a separate canopy single-scattering albedo must account for that when using specular components in a full RT layer.

compute_Z_matrices(mod::LambertianWoodCanopyScattering,
                   μ, LD, m; spectral_coordinate = nothing, grid_unit = :nm,
                   npol = 1)

Compute one Fourier moment for opaque Lambertian wood.

The wood reflectance is evaluated, then the angular Z construction delegates to BiLambertianCanopyScattering(R, 0). For non-constant wood reflectance models, pass spectral_coordinate in the unit specified by grid_unit.

compute_Z_matrices(mod::SpecularCanopyScattering,
                   μ::AbstractVector, LD::AbstractLeafDistribution,
                   m::Integer; quadrature = CanopyQuadrature())

Convenience method for non-Array vector inputs. It materializes μ with collect and delegates to the specular azimuth-integration method.

compute_Z_matrices_aniso_analytic(mod::BiLambertianCanopyScattering,
                                  μ::AbstractVector{FT},
                                  LD::AbstractLeafDistribution,
                                  m_max::Int;
                                  quadrature = CanopyQuadrature(),
                                  npol = 1) -> (Z⁺⁺, Z⁻⁺)

Compute all scalar Fourier moments m = 0:m_max of the bi-Lambertian canopy phase matrices without azimuthal quadrature.

The scalar default returns arrays of shape (length(μ), length(μ), m_max + 1). Passing npol = 3 or npol = 4 expands those scalar moments into Stokes-block matrices of shape (npol*length(μ), npol*length(μ), m_max + 1). Because bi-Lambertian scattering is unpolarized here, only the I→I entry of each stream block is populated.

All cases use the vSmartMOM convention

\[Z[i_{out}, j_{in}, m+1],\]

where Z⁺⁺ is same-sign transmission and Z⁻⁺ is sign-change reflection. The single-scattering albedo ϖ = R + T is divided out, so for conservative leaves the m = 0 column integral of Z⁺⁺ + Z⁻⁺ targets ≈ 2 under the same quadrature and G(μ) convention used by vSmartMOM.

Origin trace

Shultis and Myneni (1988), Eq. 45 writes the azimuthally averaged Lambertian canopy kernels as products of one-direction projected-area functions. For each leaf inclination this implementation generalizes that factorization from m = 0 to arbitrary cosine Fourier order by using the closed-form moments returned by _one_sided_projection_moments.

For signed incoming/outgoing directions, let (Pᵢ, Nᵢ) and (Pₒ, Nₒ) be the positive-face and negative-face projection moments. The circular-convolution theorem gives the per-leaf Fourier coefficients

\[Ψ^+_m = P_{i,m}P_{o,m} + N_{i,m}N_{o,m}, \qquad Ψ^-_m = P_{i,m}N_{o,m} + N_{i,m}P_{o,m}.\]

The code multiplies these products by 2 before applying the existing f_0 = 2, f_{m>0} = 4 factors. That leading 2 is not leaf physics; it converts the 1/(2π) Fourier coefficients above to the half-range cosine moment used by the historical compute_Z_matrices_aniso implementation: (2/π)∫_0^π Γ(Δϕ)cos(mΔϕ)dΔϕ. This is the normalization expected by vSmartMOM's elemental kernels, where ϖ is applied outside Z.

compute_Z_matrices_aniso_analytic(mod::CompositeCanopyScattering,
                                  μ, LD, m_max; quadrature)

Compute all moments 0:m_max for a composite canopy-scattering model by delegating to the stacked compute_Z_matrices API.

compute_reflection(mod::SpecularCanopyScattering, Ωⁱⁿ::dirVector_μ{FT}, Ωᵒᵘᵗ::dirVector_μ{FT},
                   LD::AbstractLeafDistribution) where FT

Returns the specular bidirectional scattering coefficient for direction pair (Ωⁱⁿ, Ωᵒᵘᵗ), following Knyazikhin & Marshak, Discrete Ordinates Method for Photon Transport in Leaf Canopies, Eq. 2.39:

$f_s(Ω' \to Ω) = \frac{1}{8}\,g_L(θ^*)\,K(κ, α^*)\,F_r(n_r, α^*)$

where:

• θ* = polar angle of the specular leaf normal (from getSpecularΩ)
• α* = incidence half-angle = arccos(Ωⁱⁿ ⋅ Ωᵒᵘᵗ) / 2
• K(κ, α*) = Nilson–Kuusk roughness factor (from K)
• Fᵣ(nᵣ, α*) = unpolarized Fresnel reflectance (from Fᵣ)

For Stokes-vector propagation, use compute_reflection_mueller, which constructs the Fresnel Mueller matrix and rotates it into the scattering-plane basis used by the canopy Fourier Z matrices.

Note: only reflection is currently modelled; specular transmission is not yet implemented.

compute_reflection(mod::SpecularCanopyScattering,
                   Ωⁱⁿ::dirVector{FT}, Ωᵒᵘᵗ::dirVector{FT}, LD) where FT

Angle-based variant of the specular leaf bidirectional scattering coefficient.

This uses the same Knyazikhin-Marshak specular-leaf model as the dirVector_μ method below, but accepts polar-angle direction vectors. The dirVector_μ method is the one used by Z-matrix assembly; this method is retained for interactive angle-grid checks.

compute_reflection_mueller(mod, Ωin, Ωout, LD, npol)

Specular leaf bidirectional scattering Mueller matrix for one direction pair.

The scalar I→I element reduces to compute_reflection. For npol = 3 or 4, the Fresnel Mueller matrix is rotated between the scattering plane and the local leaf incidence/reflection planes before the Nilson-Kuusk roughness and leaf-normal PDF factors are applied.

createLeafOpticalStruct(λ_bnds)

Loads in the PROSPECT-PRO database of pigments (and other) absorption cross section in leaves, returns a [`LeafOpticalProperties`](@ref) type struct with spectral units attached.

Arguments

- `λ_bnds` an array (with or without a spectral grid unit) that defines the upper and lower limits over which to average the absorption cross sections

Examples

julia> using Unitful                                                               # Need to include Unitful package 
julia> opti = createLeafOpticalStruct((400.0:5:2400)*u"nm");                       # in nm
julia> opti = createLeafOpticalStruct((0.4:0.1:2.4)*u"μm");                        # in μm
julia> opti = CanopyOptics.createLeafOpticalStruct((10000.0:100:25000.0)u"1/cm");  # in wavenumber (cm⁻¹)

createLeafOpticalStruct()
As in createLeafOpticalStruct(λ_bnds) but reads in the in Prospect-PRO database at original resolution (400-2500nm in 1nm steps)

dielectric(mod::LeafUlabyElRayes1987, T, f)

Complex relative permittivity ε(f) of a fresh leaf with gravimetric moisture mod.M_g, at frequency f [GHz]. Implements the dual-Debye mixing model of Ulaby & El-Rayes (1987):

\[\varepsilon = \varepsilon_r + v_{fw}\,\varepsilon_{fw} + v_b\,\varepsilon_b\]

where the residual, free-water and bound-water volume fractions are fitted from M_g, and $\varepsilon_{fw}$, $\varepsilon_b$ are the free- and bound-water dispersions. The temperature T [K] is accepted for signature consistency with the other dielectric methods but is not used by this model.

Sign convention: +i ε'' (loss positive in imaginary part), matching LiquidPureWater and SoilMW in this package.

Examples

julia> leaf = LeafUlabyElRayes1987(M_g = 0.5);

julia> dielectric(leaf, 295.0, 5.0)         # C-band, mid-moisture
14.400774… + 4.689454…im

Reference

Ulaby & El-Rayes 1987, IEEE TGRS 25(5), 550–557 (Eqs. 11–14); Ulaby & Long 2014, §11-9.

dielectric(mod::LiquidPureWater, T::FT,f::FT)

Computes the complex dieletric of liquid pure walter (Ulaby & Long book)

Arguments

• mod an LiquidSaltWater type struct
• T Temperature in [⁰K]
• f Frequency in [GHz]

Examples

julia> w = CanopyOptics.LiquidPureWater()     # Create struct for salty seawater
julia> CanopyOptics.dielectric(w,283.0,10.0)
53.45306674215052 + 38.068642430090044im

dielectric(mod::LiquidSaltWater, T::FT,f::FT)

Computes the complex dieletric of liquid salty walter (Ulaby & Long book)

Arguments

• mod an LiquidSaltWater type struct, provide Salinity in PSU
• T Temperature in [⁰K]
• f Frequency in [GHz]
• S Salinity in [PSU] comes from the mod LiquidSaltWater struct

Examples

julia> w = CanopyOptics.LiquidSaltWater(S=10.0)     # Create struct for salty seawater
julia> CanopyOptics.dielectric(w,283.0,10.0)
51.79254073931596 + 38.32304044382495im

dielectric(mod::PureIce, T::FT,f::FT)

Computes the complex dielectric constant of pure ice (Ulaby & Long book, Eqs. 4.22–4.25). Valid temperature range: 233–273.15 K.

Arguments

• mod a PureIce type struct
• T Temperature in [K], must be ∈ [233, 273.15]
• f Frequency in [GHz]

Examples

julia> w = CanopyOptics.PureIce()
julia> CanopyOptics.dielectric(w, 253.0, 10.0)   # T = -20°C, f = 10 GHz
# Real part ≈ 3.17 (ice has low imaginary dielectric compared to liquid water)

dielectric(mod::SoilMW, T::FT,f::FT)

Computes the complex dielectric constant of moist soil (Ulaby & Long book, Eqs. 4.66–4.70).

Arguments

• mod a SoilMW type struct (fields: sand_frac, clay_frac, mᵥ, ρ)
• T Temperature in [K]
• f Frequency in [GHz]

Examples

julia> w = SoilMW(sand_frac=0.2, clay_frac=0.2, mᵥ=0.3, ρ=1.7)
julia> dielectric(w, 283.0, 10.0)

effective_G(clumping, G_value, μ)

Apply clumping to a Ross projection factor: G_eff = Ω(μ) G.

This helper is intended for propagation/extinction code. Z-matrix construction should generally keep using the unclumped leaf-angle G normalization.

erectophile_leaves(FT=Float64)

Creates a LeafDistribution following a erectophile (mostly vertical) distribution using a Beta distribution

Standard values for θₗ (63.24) and s (18.5) from Bonan https://doi.org/10.1017/9781107339217.003, Table 2.1

flat_leaves(FT=Float64)

Creates a LeafDistribution representing nearly horizontal (flat) leaves using a strongly planophile Beta(1, 100) distribution, which peaks sharply near θ_L ≈ 0.

foursail!(rddt, rsdt, rdot, rsot, taus, leaf_ρ, leaf_τ, soil_ρ,
          geom, lai; hotspot=0, mode=:full)
foursail!(result::FourSAILResult, leaf_ρ, leaf_τ, soil_ρ,
          geom, lai; hotspot=0, mode=:full)

In-place 4SAIL evaluation: writes canopy + soil reflectance and transmittance factors into the supplied output arrays (or FourSAILResult) for one geometry or a batched FourSAILGeometrySet. Outputs are dispatched on the shape of the inputs:

• vector outputs + a FourSAILGeometry → one wavelength per vector element, one geometry.
• matrix outputs + a FourSAILGeometrySet → wavelengths along rows, geometries along columns; the GPU-friendly batched path.

The full 4SAIL bidirectional reflectance factor rsot is decomposed as rsot = rsost + rsodt, where rsost is the direct/direct branch (single scattering plus directly transmitted soil reflectance) and rsodt is the diffuse / multiple-scattering remainder — see FourSAILResult.

Allocating wrappers: see foursail.

Reference

Verhoef, W. (1984). Light scattering by leaf layers with application to canopy reflectance modeling: the SAIL model. Remote Sens. Env. 16, 125–141. Extended in Verhoef et al. (2007), 4SAIL.

foursail(leaf_reflectance, leaf_transmittance, soil_reflectance,
         geometry, LAI; hotspot=0)
foursail(leaf_reflectance, leaf_transmittance, soil_reflectance,
         geometry_set, LAI; hotspot=0, mode=:full)

Allocate outputs and evaluate 4SAIL over all wavelengths. For array inputs the batched kernel path is used; for scalar inputs this returns a scalar FourSAILResult.

hotspot_correction(model, k_s, k_o, μ_s, μ_o, dϕ, L)

Return the multiplicative correction C_hs applied to the independent joint gap probability. k_s and k_o are directional extinction coefficients, μ_s and μ_o are positive cosines for sun and view paths, dϕ is relative azimuth in radians, and L is cumulative LAI or area index.

hotspot_separation(μ_s, μ_o, dϕ, h)

Angular separation parameter used by the Kuusk hotspot correction:

\[α = \frac{1}{h}\sqrt{\tan^2 θ_s + \tan^2 θ_o - 2\tan θ_s\tan θ_o\cos Δϕ}.\]

joint_gap_probability(model, k_s, k_o, μ_s, μ_o, dϕ, L)

Return the full correlated bidirectional gap probability

\[P_{so}(L) = \exp[-(k_s + k_o)L] C_{hs}(L).\]

Use NoHotSpot for the independent-path Beer-Lambert result.

Hotspot corrections belong in the direct-beam source assembly of the canopy RT consumer, not in compute_Z_matrices. A typical per-layer skeleton is:

G_s = G([μ_s], LAD)[1]
G_o = G([μ_o], LAD)[1]
k_s = canopy_extinction(G_s, μ_s)
k_o = canopy_extinction(G_o, μ_o)

P_so = joint_gap_probability(hotspot, k_s, k_o, μ_s, μ_o, dϕ, L_cumulative)
P_o  = exp(-k_o * L_cumulative)
direct_source_scale = P_so / P_o

The consumer then applies direct_source_scale to the direct solar source term before multiplying by the canopy scattering kernel.

plagiophile_leaves(FT=Float64)

Creates a LeafDistribution following a plagiophile (mostly 45degrees) distribution using a Beta distribution

Standard values for θₗ (45.0) and s (16.27) from Bonan https://doi.org/10.1017/9781107339217.003, Table 2.1

planophile_leaves(FT=Float64)

Creates a LeafDistribution following a planophile (mostly horizontal) distribution using a Beta distribution

Standard values for θₗ (26.76) and s (18.51) from Bonan https://doi.org/10.1017/9781107339217.003, Table 2.1

planophile_leaves2(FT=Float64)

Creates a LeafDistribution following a planophile distribution using Beta(1, 2.768), an analytically convenient approximation to Bunnik's planophile distribution. See Shultis & Myneni (1988) Eq. (59): gL(θL) = (2/π)(1 + cos 2θ_L).

prospect(leaf::LeafProspectProProperties{FT},
                optis) where {FT<:AbstractFloat}

Computes leaf optical properties (reflectance and transittance) based on PROSPECT-PRO

Arguments

• leaf LeafProspectProProperties type struct which provides leaf composition
• optis LeafOpticalProperties type struct, which provides absorption cross sections and spectral grid

Examples

julia> opti = createLeafOpticalStruct((400.0:5:2400)*u"nm");
julia> leaf = LeafProspectProProperties{Float64}(Ccab=30.0);
julia> T,R = prospect(leaf,opti);

spherical_leaves(FT=Float64)

Creates a LeafDistribution following a spherical (distributed as facets on a sphere) distribution using a Beta distribution

Standard values for θₗ (57.3) and s (21.55) from Bonan https://doi.org/10.1017/9781107339217.003, Table 2.1

uniform_leaves(FT=Float64)

Creates a LeafDistribution following a uniform distribution (all angles equally likely)

wood_reflectance(model, x; grid_unit = :nm)
wood_reflectance(model, grid; grid_unit = :nm)

Evaluate a wood reflectance model at one spectral coordinate or a vector of coordinates. grid_unit describes the input coordinate unit.

βparameters(θₗ, s)

Returns Beta distribution parameters (α, β) from the mean leaf inclination angle θₗ and its standard deviation s (both in radians), following Bonan Eqs. 2.9–2.10:

$α = \frac{1 - (s^2 + θₗ^2)/(θₗ π/2)}{(s^2 + θₗ^2)/θₗ^2 - 1}$

$β = \left(\frac{π/2}{θₗ} - 1\right) α$

The resulting Beta(α, β) distribution is defined on [0, 1], which maps to leaf zenith angles [0, π/2] via θ_L = (π/2) u for u ~ Beta(α, β).