GitHub

4SAIL Hotspot in the Solar Principal Plane

This tutorial walks through a red/NIR principal-plane BRF diagnostic with foursail, using a realistic PROSPECT-PRO leaf and the FourSAILGeometrySet batched entry point. It shows three things in one figure:

  • the full BRF rsot with hotspot off (h = 0)
  • the full BRF rsot with a Kuusk-style hotspot turned on
  • the direct/direct branch rsost from the same hotspot solve

The shape rsot = rsost + rsodt is what makes 4SAIL convenient as a diagnostic tool — the embedded rsost is the path that carries the hotspot, so you can read off the hotspot increment without rerunning the model in a reduced mode.

using CanopyOptics
using CairoMakie

Geometry — signed VZA in the principal plane

The "principal plane" is the vertical plane containing both the sun and the observer. We sweep view zenith from −80° to +80° and split the sign by relative azimuth:

  • signed_vza < 0 → raa = 0°: backscatter side, same azimuth as the sun, the side that carries the hotspot.
  • signed_vza > 0 → raa = 180°: forward-scatter side.

(4SAIL's raa_deg argument is the azimuth difference between sun and view directions, so RAA = 0° puts the observer in the same half-plane as the sun.)

const sza_deg = 30.0
const lai     = 4.0
const hotspot_h = 0.05
const lambda_nm = [685.0, 800.0]

signed_vza = collect(-80.0:2.0:80.0)
vza_abs    = abs.(signed_vza)
raa_deg    = ifelse.(signed_vza .< 0, 0.0, 180.0)

LD    = CanopyOptics.planophile_leaves2(Float64)
geoms = FourSAILGeometrySet(LD;
    sza_deg = sza_deg,
    vza_deg = vza_abs,
    raa_deg = raa_deg,
    quadrature = CanopyOptics.CanopyQuadrature(n_leaf = 64, n_azimuth = 16),
)
FourSAILGeometrySet{Vector{Float64}}([0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386  …  0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386, 0.8660254037844386], [0.17364817766693036, 0.20791169081775934, 0.24192189559966773, 0.27563735581699916, 0.30901699437494745, 0.3420201433256687, 0.374606593415912, 0.4067366430758002, 0.4383711467890774, 0.46947156278589075  …  0.46947156278589075, 0.4383711467890774, 0.4067366430758002, 0.374606593415912, 0.3420201433256687, 0.30901699437494745, 0.27563735581699916, 0.24192189559966773, 0.20791169081775934, 0.17364817766693036], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0  …  3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793, 3.141592653589793], [5.093931550428084, 4.127279840288828, 3.4334306643462185, 2.9100641746512834, 2.5003332679856274, 2.170127150264997, 1.8977365842266702, 1.6686865047145902, 1.4729535723896703, 1.3033761961567063  …  2.458076734535958, 2.627654110768922, 2.823387043093842, 3.052437122605922, 3.3248276886442483, 3.6550338063648793, 4.064764713030535, 4.58813120272547, 5.281980378668081, 6.248632088807336], [0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747  …  0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747, 0.874073911527747], [1.6532675096438016, 1.453346937976536, 1.316625417660598, 1.2188430202059517, 1.1464939769434639, 1.0915895648821745, 1.0491017134803902, 1.015708132401324, 0.9891324578741612, 0.9677676013690241  …  0.9677676013690241, 0.9891324578741612, 1.015708132401324, 1.0491017134803902, 1.0915895648821745, 1.1464939769434639, 1.2188430202059517, 1.316625417660598, 1.453346937976536, 1.6532675096438016], [1.2881039923846997, 1.1824030401838364, 1.109297087214748, 1.056257531792672, 1.0163241863322283, 0.9853776126535339, 0.9608251978738283, 0.9409553530290862, 0.9245964451982709, 0.9109221657472828  …  0.7180820350877849, 0.7184691842983115, 0.7199417755462559, 0.7228928561433312, 0.7278852626260024, 0.7357479567324109, 0.7477448372010822, 0.7659043615917208, 0.7937031819183068, 0.8373557438850313], [0.1514904696290094, 0.10562047057456106, 0.0754603450436938, 0.05481458316162038, 0.040241551475224915, 0.029733101011722513, 0.022040313151641382, 0.016347533770835525, 0.012103525008947135, 0.00892524525015902  …  0.04890042097055637, 0.05978356987328309, 0.07337106019015456, 0.0904992062511663, 0.11235123965346179, 0.1406520569750629, 0.17800925121778177, 0.22856256914842338, 0.2993072169132302, 0.4027907320263704], [0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617  …  0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617, 0.8298315894174617], [0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539  …  0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539, 0.04424232211028539], [1.2194283884754888, 1.1194681026418563, 1.051107342483887, 1.002216143756564, 0.9660416221253201, 0.9385894160946754, 0.9173454903937832, 0.9006486998542501, 0.8873608625906687, 0.8766784343381002  …  0.8766784343381002, 0.8873608625906687, 0.9006486998542501, 0.9173454903937832, 0.9385894160946754, 0.9660416221253201, 1.002216143756564, 1.051107342483887, 1.1194681026418563, 1.2194283884754888], [0.4338391211683127, 0.3338788353346799, 0.2655180751767109, 0.21662687644938772, 0.1804523548181438, 0.15300014878749912, 0.13175622308660695, 0.11505943254707385, 0.10177159528349244, 0.0910891670309239  …  0.0910891670309239, 0.10177159528349244, 0.11505943254707385, 0.13175622308660695, 0.15300014878749912, 0.1804523548181438, 0.21662687644938772, 0.2655180751767109, 0.3338788353346799, 0.4338391211683127], [0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881  …  0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881, 0.8927946336535881], [0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187  …  0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187, 0.10720536634641187])

Leaf optics from PROSPECT-PRO

Evaluate PROSPECT once on a 1 nm grid in the 400–2500 nm range, then pull the two requested wavelengths out by nearest-neighbor index.

leaf = CanopyOptics.LeafProspectProProperties{Float64}(
    N = 1.5, Ccab = 40.0, Ccar = 8.0,
    Canth = 0.0, Cbrown = 0.0,
    Cw = 0.012, Cm = 0.009,
    Cprot = 0.0, Ccbc = 0.0,
)

opti = CanopyOptics.createLeafOpticalStruct(400.0:1.0:2500.0)
T_grid, R_grid = CanopyOptics.prospect(leaf, opti)
λ_grid   = [Float64(v.val) for v in opti.λ]
grid_idx = [argmin(abs.(λ_grid .- λ)) for λ in lambda_nm]
R_leaf   = R_grid[grid_idx]
T_leaf   = T_grid[grid_idx]

println("PROSPECT leaf optics at the two probe wavelengths:")
for (λ, ρ, τ) in zip(lambda_nm, R_leaf, T_leaf)
    println("  λ = ", λ, " nm   ρ = ", round(ρ, digits = 4),
            "   τ = ", round(τ, digits = 4))
end
┌ Warning: Assuming unit of [nm] here for optical struct
└ @ CanopyOptics ~/work/CanopyOptics.jl/CanopyOptics.jl/src/initialization/loadProspect.jl:64
PROSPECT leaf optics at the two probe wavelengths:
  λ = 685.0 nm   ρ = 0.0386   τ = 0.0121
  λ = 800.0 nm   ρ = 0.4425   τ = 0.4746

Two 4SAIL solves — hotspot off, hotspot on

foursail returns a FourSAILResult. For the batched FourSAILGeometrySet path each field is a (nwavelength × nangle) matrix.

soil_albedo = 0.0

sail_nohot = foursail(R_leaf, T_leaf, soil_albedo, geoms, lai;
                      hotspot = 0.0, mode = :full)
sail_hot   = foursail(R_leaf, T_leaf, soil_albedo, geoms, lai;
                      hotspot = hotspot_h, mode = :full)

size(sail_hot.rsot)
(2, 81)

Quick numeric readout at three reference angles: the hotspot (signed_vza = -sza), nadir, and the symmetric forward-scatter angle.

hotspot_idx = findfirst(==(-sza_deg), signed_vza)
nadir_idx   = findfirst(==(0.0), signed_vza)
forward_idx = findfirst(==(sza_deg), signed_vza)

for (iλ, λ) in pairs(lambda_nm)
    println("λ = ", λ, " nm")
    for (label, idx) in (("hotspot (-sza)", hotspot_idx),
                         ("nadir",          nadir_idx),
                         ("forward (+sza)", forward_idx))
        idx === nothing && continue
        println("  ", rpad(label, 16),
                "  rsot(h=0) = ",      round(sail_nohot.rsot[iλ, idx], digits = 4),
                "  rsot(h=$(hotspot_h)) = ", round(sail_hot.rsot[iλ, idx], digits = 4),
                "  rsost = ",          round(sail_hot.rsost[iλ, idx], digits = 4))
    end
end
λ = 685.0 nm
  hotspot (-sza)    rsot(h=0) = 0.0184  rsot(h=0.05) = 0.0355  rsost = 0.0352
  nadir             rsot(h=0) = 0.0176  rsot(h=0.05) = 0.0209  rsost = 0.0207
  forward (+sza)    rsot(h=0) = 0.0169  rsot(h=0.05) = 0.0187  rsost = 0.0184
λ = 800.0 nm
  hotspot (-sza)    rsot(h=0) = 0.4824  rsot(h=0.05) = 0.6779  rsost = 0.4032
  nadir             rsot(h=0) = 0.4738  rsot(h=0.05) = 0.5117  rsost = 0.2371
  forward (+sza)    rsot(h=0) = 0.4655  rsot(h=0.05) = 0.4863  rsost = 0.2116

Plot — full BRF and the hotspot increment

Top row: BRF curves themselves (with and without hotspot, plus the direct/direct branch rsost). Bottom row: the hotspot increment rsot(h) − rsot(0) is concentrated near vza = −sza, exactly where the camera looks back along the sun direction.

fig = Figure(size = (1100, 700))
colors = (sail_nohot = :dodgerblue3,
          sail_hot   = :seagreen4,
          sail_ss    = :darkorange3,
          increment  = :purple4)

for (iλ, λ) in pairs(lambda_nm)
    ax = Axis(fig[iλ, 1];
        xlabel = iλ == length(lambda_nm) ?
                 "signed VZA in solar principal plane (deg)" : "",
        ylabel = "BRF",
        title  = "$(Int(round(λ))) nm — principal-plane BRF")
    lines!(ax, signed_vza, sail_nohot.rsot[iλ, :];
           color = colors.sail_nohot, linewidth = 2.2, linestyle = :dash,
           label = "rsot, h = 0")
    lines!(ax, signed_vza, sail_hot.rsot[iλ, :];
           color = colors.sail_hot,   linewidth = 2.5,
           label = "rsot, h = $(hotspot_h)")
    lines!(ax, signed_vza, sail_hot.rsost[iλ, :];
           color = colors.sail_ss,    linewidth = 2.0, linestyle = :dot,
           label = "rsost (direct/direct), h = $(hotspot_h)")
    vlines!(ax, [-sza_deg]; color = :gray30, linestyle = :dashdot, linewidth = 1.2,
            label = "hotspot direction")
    vlines!(ax, [0.0];      color = :gray60, linestyle = :dot,    linewidth = 1.0)
    axislegend(ax; position = :lt, framevisible = false)

    axd = Axis(fig[iλ, 2];
        xlabel = iλ == length(lambda_nm) ?
                 "signed VZA in solar principal plane (deg)" : "",
        ylabel = "Δ BRF",
        title  = "$(Int(round(λ))) nm — hotspot increment")
    lines!(axd, signed_vza, sail_hot.rsot[iλ, :] .- sail_nohot.rsot[iλ, :];
           color = colors.increment, linewidth = 2.4,
           label = "rsot(h) − rsot(0)")
    hlines!(axd, [0.0]; color = :gray70, linewidth = 1.0)
    vlines!(axd, [-sza_deg]; color = :gray30, linestyle = :dashdot, linewidth = 1.2)
    axislegend(axd; position = :rt, framevisible = false)
end

Label(fig[0, 1:2],
      "4SAIL principal-plane BRF and hotspot increment " *
      "(SZA = $(Int(sza_deg))°, LAI = $(lai), planophile LAD, soil albedo = $(soil_albedo))";
      fontsize = 16)
fig

Reading the figure

  • Without the hotspot (h = 0) the BRF is smooth across the principal plane. The slight back–forward asymmetry is just the LAD-weighted single-scattering shape.
  • Turning the hotspot on lifts the curve sharply near vza = -sza. The increment is essentially zero on the forward-scatter side, so the multiple-scattering remainder rsodt is hardly perturbed — the hotspot lives in the rsost branch.
  • rsost itself sits well below the full BRF in the NIR because most of the 800 nm signal arrives through diffuse multiple scattering, not single scattering. In the red band leaves absorb most of what they scatter once, so rsost and the full rsot are much closer.

This page was generated using Literate.jl.