CanopyLayers

Structures

CanopyLayers.Canopy4RTType
mutable struct Canopy4RT

A canopy struct for the radiation transfer module

Fields

  • nLayer::Int64

    Number of canopy layers

  • LAI::AbstractFloat

    Leaf Area Index

  • Ω::AbstractFloat

    Clumping factor

  • clump_a::AbstractFloat

    Structure factor a

  • clump_b::AbstractFloat

    Structure factor b

  • leaf_width::AbstractFloat

    Leaf width

  • hc::AbstractFloat

    Vegetation height

  • LIDFa::AbstractFloat

    Leaf Inclination

  • LIDFb::AbstractFloat

    Variation in leaf inclination

  • hot::AbstractFloat

    HotSpot parameter (still need to check!)

  • height::AbstractFloat

    Canopy height [m]

  • z0m::AbstractFloat

    Canopy roughness [m]

  • z0h::AbstractFloat

    Tree roughtnes [m]

  • d::AbstractFloat

    Canopy displacement height [m]

  • Cd::AbstractFloat

    m/sqrt(s) turbulent transfer coefficient

  • litab::Vector{FT} where FT<:AbstractFloat

    List of mean inclination angles [°]

  • litab_bnd::Matrix{FT} where FT<:AbstractFloat

    List of inclination angle boundaries [°]

  • lazitab::Vector{FT} where FT<:AbstractFloat

    List of mean azimuth angles [°]

  • cos_ttlo::Vector{FT} where FT<:AbstractFloat

    Cosine of lazitab

  • cos_philo::Vector{FT} where FT<:AbstractFloat

    Cosine of lazitab - raa (relative azimuth angle), update with time

  • cos_ttli::Vector{FT} where FT<:AbstractFloat

    Cosine of litab

  • sin_ttli::Vector{FT} where FT<:AbstractFloat

    Sine of litab

  • vol_scatt::Vector{FT} where FT<:AbstractFloat

    Cache for volome scatter function

  • lidf::Vector{FT} where FT<:AbstractFloat

    Inclination angles weight distribution

  • xl::Vector{FT} where FT<:AbstractFloat

    List of level location (level = layer + 1)

  • dx::AbstractFloat

    1/nLayer

  • nAzi::Int64

    Number of azimuth angles

  • nIncl::Int64

    Number of inclination angles

source
CanopyLayers.CanopyOpticalsType
mutable struct CanopyOpticals{FT}

A struct for canopy optical properties

Fields

  • nAzi::Int64

    Number of azimuth angles

  • nIncl::Int64

    Number of inclination agles

  • nLayer::Int64

    Number of canopy layers

  • nWL::Int64

    Number of wave lengths

  • sdb::Any

    Solar -> Diffuse backscatter weight

  • sdf::Any

    Solar -> Diffuse forward scatter weight

  • dob::Any

    Diffuse -> Directional backscatter weight

  • dof::Any

    Diffuse -> Directional forward scatter weight

  • ddb::Any

    Diffuse -> Diffuse backscatter weight

  • ddf::Any

    Diffuse -> Diffuse forward scatter weight

  • ks::Any

    Solar beam extinction coefficient weight

  • ko::Any

    Outgoing beam extinction coefficient weight

  • bf::Any

    ?

  • sob::Any

    Weight of specular2directional backscatter coefficient

  • sof::Any

    Weight of specular2directional forward coefficient

  • Ps::Vector

    Probability of directly viewing a leaf in solar direction

  • Po::Vector

    Probability of directly viewing a leaf in viewing direction

  • Pso::Vector

    Bi-directional probability of directly viewing a leaf (solar->canopy->viewing)

  • fs::Matrix

    conversion factor fs to compute irradiance on inclined leaf

  • absfs::Matrix

    abs(fs)

  • absfsfo::Matrix

    abs(fs*fo)

  • fsfo::Matrix

    fs*fo

  • fo::Matrix

    conversion factor fo for angle towards observer (not sun like fs)

  • absfo::Matrix

    abs(fo)

  • cosΘ_l::Matrix

    Cosine of leaf azimuths

  • cos2Θ_l::Matrix

    cos of leaf azimuth sqared

  • sigb::Matrix

    diffuse backscatter scattering coefficient for diffuse incidence

  • sigf::Matrix

    diffuse forward scattering coefficient for diffuse incidence

  • sb::Matrix

    diffuse backscatter scattering coefficient for specular incidence

  • sf::Matrix

    diffuse forward scattering coefficient for specular incidence

  • vb::Matrix

    directional backscatter scattering coefficient for diffuse incidence

  • vf::Matrix

    directional forward scattering coefficient for diffuse incidence

  • w::Matrix

    bidirectional scattering coefficent (directional-directional)

  • a::Matrix

    attenuation

  • Xsd::Matrix

    Effective layer transmittance (direct->diffuse)

  • Xdd::Matrix

    Effective layer transmittance (diffuse->diffuse)

  • R_sd::Matrix

    Effective layer reflectance (direct->diffuse)

  • R_dd::Matrix

    Effective layer reflectance (diffuse->diffuse)

  • Es_::Matrix

    Solar direct radiation per layer)

source
CanopyLayers.CanopyRadsType
mutable struct CanopyRads{FT}

A struct for canopy radiation information

Fields

  • nAzi::Int64

    Number of azimuth angles

  • nIncl::Int64

    Number of inclination agles

  • nLayer::Int64

    Number of canopy layers

  • nLevel::Int64

    Number of canopy levels

  • nWL::Int64

    Number of wave lengths

  • nWLF::Int64

    Number of wave lengths for SIF

  • intEout::Any

    Integrated TOC outgoing flux [W m⁻²]

  • incomingPAR::Any

    Incident spectrally integrated total PAR [mol m⁻² s⁻¹]

  • incomingPAR_direct::Any

    Incident spectrally integrated direct PAR [mol m⁻² s⁻¹]

  • incomingPAR_diffuse::Any

    Incident spectrally integrated diffuse PAR [mol m⁻² s⁻¹]

  • RnSoil_diffuse::Any

    Net radiation of shaded soil [W m⁻²]

  • RnSoil_direct::Any

    Net Short-wave radiation of sunlit soil [W m⁻²]

  • RnSoil::Any

    Net Short-wave radiation of soil (shaded + sunlit) [W m⁻²]

  • RnSoilLW::Any

    Net long-wave radiation of soil (shaded + sunlit) [W m⁻²]

  • absPAR_shade::Vector

    Net PAR of shaded leaves [mol m⁻² s⁻¹]

  • absPAR_shadeCab::Vector

    Net PAR by Cab+Car of shaded leaves [moles m⁻² s⁻¹]

  • intNetSW_sunlit::Vector

    Spectrally integrated net absorbed direct radiation in each layer per leaf area [W m⁻²]

  • intNetSW_shade::Vector

    Spectrally integrated net absorbed diffuse radiation in each layer per leaf area [W m⁻²]

  • intNetLW_sunlit::Vector

    Spectrally integrated net absorbed direct radiation in each layer per leaf area [W m⁻²]

  • intNetLW_shade::Vector

    Spectrally integrated net absorbed diffuse radiation in each layer per leaf area [W m⁻²]

  • T_sun::Vector

    Leaf temperature (sunlit) [K]

  • T_shade::Vector

    Leaf temperature (shaded) [K]

  • ϕ_shade::Vector

    Fluorescence yield for shaded leaves

  • H_shade::Vector

    Sensible Heat flux H of shaded leaves [W m⁻²]

  • LE_shade::Vector

    Latent Heat flux LE of shaded leaves [W m⁻²]

  • NPQ_shade::Vector

    NPQ of shaded leaves

  • GPP_shade::Vector

    GPP of shaded leaves [μmol m⁻² s⁻¹]

  • gs_shade::Vector

    gs of shaded leaves [mol m⁻² s⁻¹]

  • ψl_shade::Vector

    Leaf water potential of shaded leaves [MPa]

  • Cc_shade::Vector

    Cc of shaded leaves [µmol/mol]

  • Pi_shade::Vector

    internal CO₂ concentration of shaded leaves [µmol/mol]

  • Lo::Vector

    Short-wave TOC outgoing radiance in observation direction [mW m⁻² nm⁻¹ sr⁻¹]

  • Eout::Vector

    Short-wave TOC outgoing radiation [mW m⁻² nm⁻¹]

  • alb_obs::Vector

    Short-wave Albedo in viewing direction

  • alb_direct::Vector

    Short-wave Albedo for direct incoming radiation

  • alb_diffuse::Vector

    Short-wave Albedo for diffuse incoming radiation

  • E_up::Matrix

    Upwelling diffuse short-wave radiation within canopy [mW m⁻² nm⁻¹]

  • E_down::Matrix

    Downwelling diffuse short-wave radiation within canopy [mW m⁻² nm⁻¹]

  • netSW_sunlit::Matrix

    Net absorbed direct radiation in each layer [mW m⁻² nm⁻¹]

  • netSW_shade::Matrix

    net absorbed diffuse radiation in each layer [mW m⁻² nm⁻¹]

  • absPAR_sun::Array{FT, 3} where FT

    net PAR of sunlit leaves [mol m⁻² s⁻¹]

  • absPAR_sunCab::Array{FT, 3} where FT

    net PAR by Cab+Car of sunlit leaves [mol m⁻² s⁻¹]

  • T_sun3D::Array{FT, 3} where FT

    Leaf temperature (sunlit) [K]

  • ϕ_sun::Array{FT, 3} where FT

    Fluorescence yield for sunlit leaves

  • H_sun::Array{FT, 3} where FT

    Sensible Heat flux H of sunlit leaves [W m⁻²]

  • LE_sun::Array{FT, 3} where FT

    Latent Heat flux LE of sunlit leaves [W m⁻²]

  • NPQ_sun::Array{FT, 3} where FT

    NPQ of sunlit leaves

  • GPP_sun::Array{FT, 3} where FT

    GPP of sunlit leaves [μmol m⁻² s⁻¹]

  • gs_sun::Array{FT, 3} where FT

    gs of sunlit leaves [mol m⁻² s⁻¹]

  • ψl_sun::Array{FT, 3} where FT

    Leaf water potential of sunlit leaves [MPa]

  • Cc_sun::Array{FT, 3} where FT

    Cc of sunlit leaves [µmol/mol]

  • Pi_sun::Array{FT, 3} where FT

    Internal CO₂ concentration of sunlit leaves [µmol/mol]

  • SIF_hemi::Vector

    Hemispheric total outgoing SIF flux [mW m⁻² nm⁻¹])

  • SIF_obs::Vector

    Observer-direction outgoing SIF radiance (mW m⁻² nm⁻¹ sr⁻¹))

  • SIF_obs_sunlit::Vector

    Observer-direction outgoing SIF radiance, sunlit leaves (mW m⁻² nm⁻¹ sr⁻¹)

  • SIF_obs_shaded::Vector

    Observer-direction outgoing SIF radiance, shaded leaves (mW m⁻² nm⁻¹ sr⁻¹)

  • SIF_obs_scattered::Vector

    Observer-direction outgoing SIF radiance, scattered (mW m⁻² nm⁻¹ sr⁻¹)

  • SIF_obs_soil::Vector

    Observer-direction outgoing SIF radiance, soil-reflected (mW m⁻² nm⁻¹ sr⁻¹)

  • SIF_sum::Vector

    Total SIF sum of layer sources [mW m⁻² nm⁻¹])

source
CanopyLayers.IncomingRadiationType
mutable struct IncomingRadiation{FT}

Incoming radiation information

Fields

  • E_direct::Vector

    Direct incoming radiation [mW m⁻² nm⁻¹]

  • E_diffuse::Vector

    Diffuse incoming radiation [mW m⁻² nm⁻¹]

source
CanopyLayers.LeafBiosType
mutable struct LeafBios{FT}

A struct of leaf biological parameters

Fields

  • nWL::Int64

    Number of wave length

  • nWLE::Int64

    Number of wave length for excitation

  • nWLF::Int64

    Number of wave length for SIF

  • N::Any

    Leaf structure parameter

  • Cab::Any

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

  • Car::Any

    Carotenoid content [µg cm⁻²]

  • Ant::Any

    Anthocynanin content [µg cm⁻²]

  • Cs::Any

    Senescent material fraction

  • Cw::Any

    Equivalent water thickness [cm]

  • Cm::Any

    Dry matter content (dry leaf mass per unit area) [g cm⁻²]

  • Cx::Any

    Fractionation between Zeaxanthin and Violaxanthin in Car (1=all Zeaxanthin) (-)

  • fqe::Any

    Leaf fluorescence efficiency (Fo standard)

  • ρ_LW::Any

    Broadband thermal reflectance (-)

  • τ_LW::Any

    Broadband thermal transmission (-)

  • ρ_SW::Vector

    Shortwave leaf reflectance

  • τ_SW::Vector

    Shortwave leaf transmission

  • α_SW::Vector

    Shortwave absorption

  • kChlrel::Vector

    Relative absorbtion by Chlorophyll+Car

  • kChlrel_old::Vector

    Relative absorbtion by Chlorophyll

  • Mb::Matrix

    Fluorescence excitation matrix backwards

  • Mf::Matrix

    Fluorescence excitation matrix forwards

  • ndub::Int64

    Doubling adding layers

source
CanopyLayers.LeafOpticalsType
mutable struct LeafOpticals{FT}

Struct for leaf optical properties

Fields

  • nr::Vector

  • Km::Vector

  • Kab::Vector

  • Kant::Vector

  • Kcar::Vector

  • Kw::Vector

  • KBrown::Vector

  • phi::Vector

  • KcaV::Vector

  • KcaZ::Vector

  • lambda::Vector

    Wave length [nm], same as WL in WaveLengths`

source
CanopyLayers.RTDimensionsType
mutable struct RTDimensions

Struct that stores matrix dimension information

Fields

  • nAzi::Int64

    Number of azimuth angles

  • nIncl::Int64

    Number of inclination agles

  • nLayer::Int64

    Number of canopy layers

  • nLevel::Int64

    Number of canopy layer boundaries nLayer+1

  • nPAR::Int64

    Number of PAR wave lengths

  • nWL::Int64

    Number of wave lengths

  • nWLE::Int64

    Number of wave length for excitation

  • nWLF::Int64

    Number of wave lengths for SIF

source
CanopyLayers.SoilOpticalsType
mutable struct SoilOpticals{FT}

A struct of soil optical parameters

Fields

  • T::Any

    Soil surface temperature

  • color::Int64

    Soil color class

  • ρ_NIR::Any

    Shortwave albedo for NIR

  • ρ_PAR::Any

    Shortwave albedo for PAR

  • ρ_SW::Vector

    Shortwave albedo that matches WL from Wavelengths

  • ρ_SW_SIF::Vector

    Shortwave albedo that matches WLF from Wavelengths

  • ε_SW::Vector

    Shortwave absorption that equals 1 - ρ_SW

  • SW_mat_raw_4::Matrix

    Shortwave albedo matrix from 4 bands, wavelengths are 400:10:2500 nm

  • SW_mat_raw_2::Matrix

    Shortwave albedo matrix from 2 bands, wavelengths are 400:10:2500 nm

  • SW_mat_4::Matrix

    Shortwave albedo matrix from 4 bands with WL from Wavelengths

  • SW_mat_2::Matrix

    Shortwave albedo matrix from 2 bands with WL from Wavelengths

  • SW_vec_4::Vector

    Shortwave albedo weight from 4 bands

  • SW_vec_2::Vector

    Shortwave albedo weight from 2 bands

  • ρ_SW_raw::Vector

    Shortwave albedo, wavelengths are 400:10:2500 nm

  • ρ_LW::Vector

    Longtwave albedo

  • dry_NIR::Any

    Mean value for day band 1 in NIR region

  • dry_PAR::Any

    Mean value for day band 1 in PAR region

  • wet_NIR::Any

    Mean value for day band 1 in NIR region

  • wet_PAR::Any

    Mean value for day band 1 in PAR region

source
CanopyLayers.SolarAnglesType
struct SolarAngles{FT}

Struct for observation and solar angles

Fields

  • hza::Any

    Hill zenith angle [°], hill slope angle

  • haa::Any

    Hill azimuth angle [°], 0 for North, 180 for south

  • saa::Any

    Solar azimuth angle [°], a function of time

  • sza::Any

    Solar zenith angle [°], a function of lat and time

  • vaa::Any

    Viewing azimuth angle [°]

  • vza::Any

    Viewing zenith angle [°]

  • raa::Any

    Relative azimuth angle [°], difference between saa and vaa

source
CanopyLayers.WaveLengthsType
mutable struct WaveLengths{FT}

Struct for pre-set wave length parameters

Fields

  • minwlPAR::Any

    Minimal WL for PAR [nm]

  • maxwlPAR::Any

    Maximal WL for PAR [nm]

  • maxwlNIR::Any

    Maximal WL for NIR [nm]

  • minwle::Any

    Minimal WL for SIF excitation [nm]

  • maxwle::Any

    Maximal WL for SIF excitation [nm]

  • minwlf::Any

    Minimal WL for SIF emission/fluorescence [nm]

  • maxwlf::Any

    Maximal WL for SIF emission/fluorescence [nm]

  • sWL::Vector

    Standard wave length [nm]

  • dWL::Vector

    Differential wavelength

  • optis::LeafOpticals

    Leaf optical parameter set

  • WL::Vector

    Wave length [nm]

  • iWLE::Vector{Int64}

    Index of WLE in WL

  • iWLF::Vector{Int64}

    Index of WLF in WL

  • iPAR::Vector{Int64}

    index of wlPAR in WL

  • iNIR::Vector{Int64}

    index of wlNIR in WL

  • WLE::Vector

    excitation wave length [nm]

  • WLF::Vector

    Fluorescence wave length [nm]

  • WL_iPAR::Vector

    Wave length for PAR

  • dWL_iPAR::Vector

    Differential wave length for PAR

  • dWL_iWLE::Vector

    Differential wave length for iWLE

  • nPAR::Int64

    Length of WL_iPAR

  • nWL::Int64

    Length of WL

  • nWLE::Int64

    length of WLE

  • nWLF::Int64

    length of WLF

source

Caches

CanopyLayers.CFCacheType
mutable struct CFCache{FT}

Cache to speed canopy_fluxes! by pre-allocating arrays

Fields

  • abs_wave::Vector

    absorbed energy from wave lengths

  • absfs_lidf::Vector

    absfs' * lidf [nAzi]

  • E_all::Vector

    wave length energy [same as dWL]

  • E_iPAR::Vector

    wave length energy [same as iPAR]

  • lPs::Vector

    lPs [nLayer]

  • kChlrel::Vector

    kChlrel [same as iPAR]

  • PAR_diff::Vector

    diffusive PAR [same as iPAR]

  • PAR_diffCab::Vector

    diffusive PAR for photosynthesis [same as iPAR]

  • PAR_dir::Vector

    direct PAR [same as iPAR]

  • PAR_dirCab::Vector

    diffusive PAR for photosynthesis [same as iPAR]

source
CanopyLayers.CGCacheType
mutable struct CGCache{FT}

Cache to speed canopy_geometry! by pre-allocating arrays

Fields

  • _Co::Vector

    cos_ttli .* cos(vza) dim: nIncl

  • _Cs::Vector

    cos_ttli .* cos(sza) dim: nIncl

  • _So::Vector

    sin_ttli .* sin(vza) dim: nIncl

  • _Ss::Vector

    sin_ttli .* sin(sza) dim: nIncl

  • _1s::Matrix

    maxtrix filled with 1 dim: (1, nAzi)

  • _2d::Matrix

    2D array to speed up _cds and _cdo dim: (nIncl, nAzi)

  • _cdo::Matrix

    Co * _1s .+ _So * cosphilo' dim: (nIncl, nAzi)

  • _cds::Matrix

    Cs * _1s .+ _Ss * costtlo' dim: (nIncl, nAzi)

source
CanopyLayers.SFCacheType
mutable struct SFCache{FT}

Cache to speed SIF_fluxes! by pre-allocating arrays

Fields

  • M⁻_sun::Vector

  • M⁺_sun::Vector

  • wfEs::Vector

  • sfEs::Vector

  • sbEs::Vector

  • M⁺⁻::Vector

  • M⁺⁺::Vector

  • M⁻⁺::Vector

  • M⁻⁻::Vector

  • sun_dwl_iWlE::Vector

  • tmp_dwl_iWlE::Vector

  • ϕ_cosΘ_lidf::Vector

  • vfEplu_shade::Vector

  • vbEmin_shade::Vector

  • vfEplu_sun::Vector

  • vbEmin_sun::Vector

  • sigfEmin_shade::Vector

  • sigbEmin_shade::Vector

  • sigfEmin_sun::Vector

  • sigbEmin_sun::Vector

  • sigfEplu_shade::Vector

  • sigbEplu_shade::Vector

  • sigfEplu_sun::Vector

  • sigbEplu_sun::Vector

  • zeroB::Vector

  • tmp_1d_nWlF::Vector

  • tmp_1d_nLayer::Vector

  • dnorm::Vector

  • τ_dd::Matrix

    transmission of diffusive light?

  • ρ_dd::Matrix

    extinction of diffuse light?

  • Xdd::Matrix

  • Rdd::Matrix

  • Y::Matrix

  • U::Matrix

  • S⁻::Matrix

  • S⁺::Matrix

  • piLs::Matrix

  • piLd::Matrix

  • Fsmin::Matrix

  • Fsplu::Matrix

  • Fdmin::Matrix

  • Fdplu::Matrix

  • Femo::Matrix

  • M⁺::Matrix

  • M⁻::Matrix

  • ϕ_cosΘ::Matrix

  • F⁻::Matrix

  • F⁺::Matrix

  • net_diffuse::Matrix

  • tmp_2d_nWlF_nLayer::Matrix

  • tmp_2d_nWlF_nLayer_2::Matrix

source
CanopyLayers.SWCacheType
mutable struct SWCache{FT}

Cache to speed short_wave! by pre-allocating arrays

Fields

  • dnorm::Vector

    dnorm?

  • piLo::Vector

    pi * Lo

  • piLoc::Vector

    pi * Lo from canopy

  • piLos::Vector

    pi * Lo from soil

  • piLoc2::Matrix

    pi * Lo from canopy 2D matrix

  • ρ_dd::Matrix

    extinction of diffuse light?

  • ρ_sd::Matrix

    extinction of direct light?

  • τ_dd::Matrix

    transmission of diffusive light?

  • τ_sd::Matrix

    transmission of direct light?

source
CanopyLayers.RTCacheType
mutable struct RTCache{FT}

Collection of caches to speed up RT module

Fields

  • cf_con::CanopyLayers.CFCache

    CFCache type cache

  • cg_con::CanopyLayers.CGCache

    CGCache type cache

  • sf_con::CanopyLayers.SFCache

    SFCache type cache

  • sw_con::CanopyLayers.SWCache

    SWCache type cache

source

Initialization of Structures

CanopyLayers.create_incoming_radiationFunction
create_incoming_radiation(
            wls::WaveLengths{FT},
            wlfn::String = FILE_SUN
) where {FT<:AbstractFloat}

Create an AbstractIncomingRadiation struct, given

  • wls WaveLengths type struct
  • wlfn File that saves incoming wave information
source
CanopyLayers.create_leaf_opticalsFunction
create_leaf_opticals(
            sWL::Array{FT,1},
            file::String = OPTI_2021
) where {FT<:AbstractFloat}

Create an AbstractLeafOptiPara struct, given

  • sWL Standard wave length
  • opfn File that saves optical parameters
source

Big Leaf Model

CanopyLayers.big_leaf_partitionFunction
big_leaf_partition(
            lai::FT,
            zenith::FT,
            r_all::FT,
            r_dir::FT = FT(0.8)
) where {FT<:AbstractFloat}

Partition the big-leaf canopy into sunlit and shaded layers, given

  • lai Leaf area index
  • zenith Zenith angle in degree
  • r_all Total radiation in [W m⁻²]
  • r_dir Direct radiation partition in r_all

The function returns

  • ratio ratio of sunlit leaves out of all leaves
  • q_slm Mean sunlit layer PAR
  • q_shm Mean shaded layer PAR
  • e_sl Mean sunlit layer absorbed total energy
  • e_sh Mean shaded layer absorbed total energy
source

SCOPE Model

CanopyLayers.canopy_fluxes!Function
canopy_fluxes!(
            can::Canopy4RT{FT},
            can_opt::CanopyOpticals{FT},
            can_rad::CanopyRads{FT},
            in_rad::IncomingRadiation{FT},
            soil::SoilOpticals{FT},
            leaves::Array{LeafBios{FT},1},
            wls::WaveLengths{FT},
            rt_con::RTCache{FT}
) where {FT<:AbstractFloat}

Computes a variety of integrated fluxes from the spectrally resolved computations in the short-wave Canopy RT (e.g. absorbed soil radiation, absorbed direct and diffuse PAR by layer (and angles for direct), net direct and diffuse energy balance per layer), given

source
CanopyLayers.canopy_geometry!Function
canopy_geometry!(
            can::Canopy4RT{FT},
            angles::SolarAngles{FT},
            can_opt::CanopyOpticals{FT},
            rt_con::RTCache{FT}
) where {FT<:AbstractFloat}

Computes canopy optical properties (extinction coefficients for direct and diffuse light) based on the SAIL model. Most important input parameters are leaf inclination and azimuth distribution functions and sun-sensor geometry. Canopy clumping Ω is implemented as in Pinty et al (2015), given

source
CanopyLayers.canopy_matrices!Function
canopy_matrices!(
            leaves::Array{LeafBios{FT},1},
            can_opt::CanopyOpticals{FT}
) where {FT<:AbstractFloat}

Compute scattering coefficient matrices for direct and diffuse light given geometry dependent overall extinction coefficients and pigment dependent leaf reflectance and transmission (computed via fluspect). This function has to be called before short_wave! can be used.

source
CanopyLayers.diffusive_SFunction
diffusive_S(τ_dd::Array{FT},
            ρ_dd::Array{FT},
            S⁻::Array{FT},
            S⁺::Array{FT},
            boundary_top::Array{FT},
            boundary_bottom::Array{FT},
            rsoil::Array{FT}
) where {FT<:AbstractFloat}

Computes 2-stream diffusive radiation transport (used for thermal and SIF) given:

  • τ_dd A 2D Array with layer reflectances
  • ρ_dd A 2D Array with layer transmissions
  • S⁻ A 2D Array with layer source terms in the downwelling direction
  • S⁺ A 2D Array with layer source terms in the upwelling direction
  • boundary_top A 1D array with downwelling radiation at the top (top of canopy)
  • boundary_bottom A 1D array with upwnwelling radiation at the bottom (soil)
  • rsoil A 1D array with soil reflectance
source
CanopyLayers.diffusive_S!Function
diffusive_S!(
            sf_con::SFCache{FT},
            soil::SoilOpticals{FT},
            rt_dim::RTDimensions
) where {FT<:AbstractFloat}

Computes 2-stream diffusive radiation transport (used for thermal and SIF), given

source
CanopyLayers.fluspect!Function
fluspect!(leaf::LeafBios{FT},
          wls::WaveLengths{FT};
          APAR_car::Bool = true
) where {FT<:AbstractFloat}

Computes leaf optical properties (reflectance and transittance) based on pigment concentrations. Also computes Fluorescence excitation matrices. Mostly based on PROSPECT-D for leaf reflectance/transmission and FluSpec for fluorescence.

  • leaf LeafBios type struct
  • wls WaveLengths type struct
  • APAR_car If true, include Car absorption in APAR for photosynthesis
source
CanopyLayers.SIF_fluxes!Function
SIF_fluxes!(leaves::Array{LeafBios{FT},1},
            can_opt::CanopyOpticals{FT},
            can_rad::CanopyRads{FT},
            can::Canopy4RT{FT},
            soil::SoilOpticals{FT},
            wls::WaveLengths{FT},
            rt_con::RTCache{FT},
            rt_dim::RTDimensions;
            photon::Bool = true
) where {FT<:AbstractFloat}

Computes 2-stream diffusive radiation transport for SIF radiation (calls [diffusive_S!] internally). Layer reflectance and transmission is computed from SW optical properties, layer sources from absorbed light and SIF efficiencies. Boundary conditions are zero SIF incoming from atmosphere or soil.

source
SIF_fluxes!(leaf::LeafBios{FT},
            in_rad::IncomingRadiation{FT},
            wls::WaveLengths{FT},
            rt_con::RTCache{FT},
            fqe::FT = FT(0.01);
            photon::Bool = true
) where {FT<:AbstractFloat}

Leaf level SIF, given

  • leaf LeafBio type struct
  • in_rad IncomingRadiation type struct
  • wls WaveLengths type struct
  • rt_con RTCache type cache
  • fqe Fluorescence quantum yield (default at 1%)
  • photon If true, use photon unit in the matrix conversion

Note that in_rad assumes direct light with zenith angle of 0, and a zenith angle correction needs to be made before passing it to this function. The up- and down-ward SIF are stored in sf_con as M⁻_sun and M⁺_sun.

source
CanopyLayers.thermal_fluxes!Function
thermal_fluxes!(
            leaves::Array{LeafBios{FT},1},
            can_opt::CanopyOpticals{FT},
            can_rad::CanopyRads{FT},
            can::Canopy4RT{FT},
            soil::SoilOpticals{FT},
            incLW::Array{FT},
            wls::WaveLengths{FT}
) where {FT<:AbstractFloat}

Computes 2-stream diffusive radiation transport for thermal radiation (calls [diffusive_S] internally). Layer reflectance and transmission is computed from LW optical properties, layer sources from temperature and Planck law, boundary conditions from the atmosphere and soil emissivity and temperature. Currently only uses Stefan Boltzmann law to compute spectrally integrated LW but can be easily adjusted to be spectrally resolved.

source

Indicies

CanopyLayers.REF_WLFunction
REF_WL(wls::WaveLengths{FT},
       can_rad::CanopyRads{FT}
       wls::WaveLengths{FT},
       twl::FT
) where {FT<:AbstractFloat}

Return the Reflectance, given

source
CanopyLayers.SIF_757Function
SIF_757(can_rad::CanopyRads{FT},
        wls::WaveLengths{FT};
        oco::Int = 2
) where {FT<:AbstractFloat}

Return the SIF @ 758.68 (OCO2) or 758.77 (OCO3) nm, given

source
CanopyLayers.SIF_771Function
SIF_771(can_rad::CanopyRads{FT},
        wls::WaveLengths{FT};
        oco::Int = 2
) where {FT<:AbstractFloat}

Return the SIF @ 769.94 (OCO2) or 770.005 (OCO3) nm, given

source
CanopyLayers.SIF_WLFunction
SIF_WL(wls::WaveLengths{FT},
       can_rad::CanopyRads{FT}
       wls::WaveLengths{FT},
       twl::FT
) where {FT<:AbstractFloat}

Return the SIF, given

source
CanopyLayers.SPECTRUM_WLFunction
SPECTRUM_WL(wls::WaveLengths{FT},
       can_rad::CanopyRads{FT}
       wls::WaveLengths{FT},
       twl::FT
) where {FT<:AbstractFloat}

Return the Reflectance, given

source

Utils

CanopyLayers.calctavFunction
calctav(α::FT, nr::FT) where {FT<:AbstractFloat}

Computes transmission of isotropic radiation across an interface between two dielectrics (Stern F., 1964; Allen W.A., 1973)). From calctav.m in PROSPECT-D

  • α angle of incidence
  • nr Index of refraction
source
CanopyLayers.dcumFunction
dcum(a::FT, b::FT, t::FT) where {FT<:AbstractFloat}

TODO Add function description

source
CanopyLayers.dladgenFunction
dladgen(a::FT, b::FT, litab_bnd::Array{FT,2}) where {FT<:AbstractFloat}

TODO Calculate the freqency of WHAT?

source
CanopyLayers.e2photFunction
e2phot(λ::Array{FT}, E::Array{FT}) where {FT<:AbstractFloat}

Calculates the number of moles of photons, given

  • λ An array of wave length in [nm], converted to [m] by _FAC
  • E Joules of energy
source
CanopyLayers.e2phot!Function
e2phot!(
            λ::Array{FT,1},
            E::Array{FT,1},
            cache::Array{FT,1}
) where {FT<:AbstractFloat}

Calculates the number of moles of photons, given

  • λ An array of wave length in [nm], converted to [m] by _FAC
  • E Joules of energy
  • cache Cache to avoid memory allocations
source
CanopyLayers.psofunctionFunction
psofunction(K::FT,
            k::FT,
            Ω::FT,
            LAI::FT,
            q::FT,
            dso::FT,
            xl::FT
) where {FT<:AbstractFloat}

TODO explain the variables

Return the probability of observing a sunlit leaf at depth xl (pso, see eq 31 in vdT 2009), given

  • xl Leaf depth in the canopy
source
CanopyLayers.volscatt!Function
volscatt!(cache::Array{FT,1},
          sza::FT,
          vza::FT,
          raa::FT,
          ttl::FT
) where {FT<:AbstractFloat}

Calculate interception parameters (chi_s and chi_s) and leaf reflectance multiplier (frho) and transmittance multiplier (ftau), given

  • cache Array cache for results
  • sza Solar zenith angle
  • vza Viewing zenith angle
  • raa Relative azimuth angle
  • ttl Leaf inclination angle
source