-
Notifications
You must be signed in to change notification settings - Fork 8
Spanner LAI Example
This example script demonstrates the calculation of Leaf Area Index (LAI) using gap fraction profiles, Leaf Area Density (LAD), and related methods in R. It utilizes the lidR
and terra
libraries.
library(lidR)
library(terra)
This function calculates the gap fraction profile of a canopy given a vertical distribution of points.
gap_fraction_profile <- function(z, dz = 1, z0 = 2) {
zrange = range(z)
if (z0 < zrange[1])
z0 = floor((zrange[1] - z0)/dz) * dz + z0
if (z0 >= zrange[2])
return(data.frame(z = numeric(0), gf = numeric(0)))
bk <- seq(z0, ceiling((zrange[2] - z0)/dz) * dz + z0, dz)
histogram <- graphics::hist(z, breaks = c(-Inf, bk), plot = FALSE)
h <- histogram$mids
p <- histogram$counts
cs <- cumsum(p)
i <- cs[1:(length(cs) - 1)]/cs[2:length(cs)]
i[is.na(i)] = 0
z = h[-1]
return(data.frame(z = z, gf = i))
}
This function calculates the Leaf Area Density (LAD) using the gap fraction profile.
LAD <- function(z, dz = 1, k = 0.5, z0 = 2) {
ld <- gap_fraction_profile(z, dz, z0)
if (nrow(ld) == 0)
return(data.frame(z = numeric(0), lad = numeric(0)))
lad <- ld$gf
lad <- -log(lad)/(k * dz)
lad[is.infinite(lad)] = NA
return(data.frame(z = ld$z, lad = lad))
}
The Leaf Area Density (LAD) is calculated using the formula:
Where:
- ( gf ) is the gap fraction,
- ( k ) is the extinction coefficient, and
- ( dz ) is the vertical resolution.
This function computes the Leaf Area Index (LAI) by integrating the LAD over the height.
calculate_lai <- function(z_values, dz = 1, z0 = 2) {
# Compute the Leaf Area Density (LAD)
lad <- LAD(z_values, dz = dz, z0 = z0)
# Integrate the LAD to compute LAI
lai <- sum(lad$lad, na.rm = TRUE) * dz
return(lai)
}
The Leaf Area Index (LAI) is the integral of LAD over height ( z ):
Where:
- ( LAD(z) ) is the leaf area density at height ( z ),
- ( dz ) is the height step used in the numerical integration.
You can use the calculate_lai
function with a set of vertical height values to compute the LAI for a given canopy.
dz = 1
z0 = 2
z_values <- c(1, 3, 5, 7, 9, 11) # Example height values
lai <- calculate_lai(z_values, dz = dz, z0 = z0)
print(lai)
This will print the computed LAI for the provided z_values
.