In RCBlackburn/lidRmts: Area based lidar attributes

library(rgl)
library(knitr)

A package to be used in conjuction with lidR (https://github.com/Jean-Romain/lidR) that produces area-based variables. Newer functions focus on the use of tree- and voxel-based attributes. These approaches are further outlined in Blackburn et al. 2021 (https://cdnsciencepub.com/doi/10.1139/cjfr-2020-0506). Below are some examples:

Install package

# devtools::install_github("RCBlackburn/lidRattR")

Load in library and example data from lidR package

library(lidR)
library(lidRattR)
library(data.table)
library(terra)
LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
las <- readLAS(LASfile)
las <- filter_poi(las, Classification != LASNOISE)
las <- normalize_height(las, knnidw(k = 8, p = 2))
las <- filter_poi(las, Z < 50 & Z >= 1.37 )
las <- decimate_points(las, random_per_voxel(res = 1, n = 8))

Tree-based approach

The tree based approach summarizes individual tree based information to a given area (e.g., plot, pixel).

First, We segment the trees using a tree segmentation function from lidR. Here, we use the Li et al. 2012 algorithm.

las_tree <- segment_trees(las, li2012())
plot(las_tree, color = "treeID")

Next, we can turn the segmented tree point cloud into a summary statistics of tree attributes to either the extent of the point cloud or within individual raster cells (resolution defined by rast_res).

## entire las summary statistics 
tree_metrics <- tree_summary(las_tree)
kable(tree_metrics)

## rasterized summary statistics 
tree_metrics <- tree_raster(las_tree, rast_res = 20)
plot(tree_metrics)

Voxel-based approach

The voxel-based approach voxelizes the point cloud and summarizes voxel information across a given area. This allows voxel-based variables to be assessed using an area-based approach.

Using these functions, we produce summary statistics of voxel attributes to either the entire las or within individual raster cells (resolution defined by rast_res). These attributes can be computed without defining a specified maximum height for the analysis (max_z = NULL) or setting a maximum height that must be larger than the highest Z value in the study area. Applicability of certain voxel attributes depends on raster and voxel resolution.

## entire las summary statistics 
voxel_metrics <- voxel_summary(las, vox_res = 2, max_z = 50)
kable(voxel_metrics)

## rasterized summary statistics  
voxel_metrics <- voxel_raster(las, vox_res = 2, rast_res = 20, max_z = 50)
plot(voxel_metrics)

Catalog processing

Tree and voxel summaries can also be implemented using lidR's catalog processing engine:

## read in lidar catalog
LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
ctg = readLAScatalog(LASfile)

## set catlog options 
opt_chunk_buffer(ctg) <- 10 
opt_chunk_size(ctg) <- 100    
opt_select(ctg) <- "xyzicr"      
opt_wall_to_wall(ctg) <- TRUE
opt_chunk_alignment(ctg) <- c(0,0)
opt_progress(ctg) <- FALSE
opt_stop_early(ctg) <- TRUE

## check chunks and overlap
plot(ctg, chunk = TRUE)

## define custom function for catalog processing 
all_layers_catalog = function(cluster)
{ 
  las = readLAS(cluster)
  las = las_update(las) 
  if (lidR::is.empty(las)) return(NULL) 
  las = filter_poi(las, Classification != LASNOISE)
  if (lidR::is.empty(las)) return(NULL)
  las = normalize_height(las, knnidw(k = 8, p = 2)) 
  las = filter_poi(las, Z < 49 & Z >= 1.37 ) 
  if (lidR::is.empty(las)) return(NULL)
  las <- decimate_points(las, random_per_voxel(res = 1, n = 8))
  if (lidR::is.empty(las)) return(NULL)

  vmetrics <-voxel_raster(las,vox_res = 2, rast_res = 20, max_z = 40)
  las <- segment_trees(las, li2012())
  tmetrics <-tree_raster(las, rast_res = 20)
  vmetrics <- resample(vmetrics,tmetrics)
  metrics <- c(tmetrics,vmetrics)
  tile0_ext = lidR::ext(las)
  vmetrics <- crop(metrics, tile0_ext)

  return(wrap(vmetrics))
}


layer_list <- catalog_apply(ctg, all_layers_catalog)

raster <- lapply(layer_list, terra::unwrap)
raster <- do.call(terra::merge, raster)

plot(raster)