rlidr

Attempting to Calculate Leaf Area Index in lidR package results in "Duplicate elements found error"


I am attempt to write a function to calculate leaf area index using discrete return aerial lidar data using the lidR package in R. My function follows the Beer-Lambert Law with the assumption of randomly oriented leaves so that the Ross-Nilson gap function value = 0.5 in the calculation of the extinction coefficient. With discrete return data, I am calculating the ratio of intercepted to transmitted photosynthetically active radiation as the ratio of canopy returns (defined as above 3 meters) to ground returns (defined as below 2 meters). I have created a "low-level api" function called Calc_LAI() that does the following:

  1. Counts the number of ground returns
  2. Counts the number of vegetation returns
  3. Calculates leaf area index within a call to pixel_metrics()
  4. Returns a raster of leaf area index values

I am stumped by two issues:

First, I have a helper function called lai_func() that does the arithmetical calculation of leaf area index. Currently I have it nested in the conditional, for if(is(las,"LAS")){}. However, when I run the whole Calc_LAI() function, I get an error saying that the lai_func() isn't found. If I move the helper function out of Calc_LAI() this error goes away, but I have had helper functions in other "low-level api" calls to catalog_apply() where they were nested within the if(is(las,"LAS")){} conditional that worked without error.

Second, the call to pixel_metrics() results in the following error Error: Duplicated elements found. At least one of the metrics was not a number. Each metric should be a single number. which perplexes me because I ensure that I only get real values from the calculation of leaf area index by adding 1 to both the numerator and denominator in the ratio of canopy returns to ground returns to ensure that I never divide by zero or take the natural log of 0.

Below is my reproducible example:

##Loading Necessary Package##
library(lidR)

##Reading in Example Lidar Dataset##
lasfile<- system.file("extdata", "MixedConifer.laz", package="lidR")
MixedConifer<-readLAS(lasfile)


##Low-level API function to calculate leaf area index##
Calc_LAI <- function(las)
{
  if (is(las, "LAScatalog"))  {
    options <- list(automerge = TRUE, need_buffer = TRUE)
    LAI <- catalog_apply(las, Calc_LAI, .options = options)
    return(LAI)
  }
  else if (is(las, "LAScluster")) {
    bbox <- st_bbox(las)
    las <- readLAS(las)
    if (is.empty(las)) return(NULL) 
    LAI <- Calc_LAI(las, param1, param2)
    LAI <- sf::st_crop(LAI, bbox)
    return(LAI)
  }
  else if (is(las, "LAS")) {
    las<-add_attribute(las, paste(las$gpstime, las$ReturnNumber, sep="_"), "UID")
    dem<-rasterize_terrain(las, res=1, algorithm = knnidw())
    norm<-las-dem
    
    lai_func<-function(ScanAngle, Z, UID){
      nground<-length(unique(UID[Z <= 6.6]))
      nveg<-length(unique(UID[Z >= 9.9]))
      out<-log(nveg+1/(nground+1))*cos((ScanAngle*pi/180)/0.5)
      return(list(LAI=out))
    }
    
    LAI<-pixel_metrics(las=norm, func=lai_func(ScanAngle=ScanAngleRank, Z=Z, UID=UID), res=9.9)
    return(LAI)
  }
  else {
    stop("This type is not supported.")
  }
}

##First Attempt with lai_func() nested within Calc_LAI()##
Calc_LAI(MixedConifer)#throws error that lai_func() not found

##Second Attempt with lai_func() outside of Calc_LAI##
lai_func<-function(ScanAngle, Z, UID){
  nground<-length(unique(UID[Z <= 6.6]))
  nveg<-length(unique(UID[Z >= 9.9]))
  out<-log(nveg+1/(nground+1))*cos((ScanAngle*pi/180)/0.5)
  return(list(LAI=out))
}

Calc_LAI(MixedConifer)#throws error that Duplicate elements are found and at least one metric was not a number

##Running the function elements directly on example data set
las<-MixedConifer
las<-add_attribute(las, paste(las$gpstime, las$ReturnNumber, sep="_"), "UID")
dem<-rasterize_terrain(las, res=1, algorithm = knnidw())
norm<-las-dem

lai_func<-function(ScanAngle, Z, UID){
  nground<-length(unique(UID[Z <= 6.6]))
  nveg<-length(unique(UID[Z >= 9.9]))
  out<-log(nveg+1/(nground+1))*cos((ScanAngle*pi/180)/0.5)
  return(list(LAI=out))
}

LAI<-pixel_metrics(las=norm, func=lai_func(ScanAngle=ScanAngleRank, Z=Z, UID=UID), res=9.9) #throw the duplicate elements error again


Solution

  • throws error that Duplicate elements are found and at least one metric was not a number

      out<-log(nveg+1/(nground+1))*cos((ScanAngle*pi/180)/0.5)
      return(list(LAI=out))
    

    Here, clearly LAI is not a metric (i.e a number) but is a vector of size length(ScanAngle)