TITLE: Estimating grassy volume from terrestrial LiDAR
DATE: 2021-01-15
AUTHOR: John L. Godlee
====================================================================


As part of my PhD research I have been using terrestrial LiDAR to 
understand woodland tree canopy traits in southern African 
savannas. One of the measurements I wanted to make was an estimate 
of the volume of grassy material below the canopy, so I could make 
inferences about how canopy traits affect the probability of fire. 
In a previous post I've already described how I processed the raw 
point cloud data to produce a .csv with XYZ point coordinates, so 
I'll skip straight to how I used R to estimate grassy volume. I 
adapted this method from Cooper et al. (2017).

  [Cooper et al. (2017)]: https://doi.org/10.3390/rs9060531

After reading in the file with data.table::fread() the first thing 
was to assign each point within a cylinder of grass to 2x2 cm 2D 
bins in the XY plane, then I took the mean height of points within 
each bin and estimated the volume of the column below that mean 
height, assuming that the volume below the mean height was 
completely filled by grass material. That's quite a big assumption 
to make, but for comparison between samples it seems suitable. I 
also have Disc Pasture Meter (DPM) measurements and biomass samples 
on a subset of the sample locations to cross-check the estimates 
from the terrestrial LiDAR.

    # Read file
    dat <- fread(x)

    # Bin into x,y cells
    dat_xy_bin <- dat %>%
    mutate(
      bin_x = cut(.$X, include.lowest = TRUE, labels = FALSE,
        breaks = seq(floor(min(.$X)), ceiling(max(.$X)), by = 
voxel_dim)),
      bin_y = cut(.$Y, include.lowest = TRUE, labels = FALSE,
        breaks = seq(floor(min(.$Y)), ceiling(max(.$Y)), by = 
voxel_dim)))

    # Take mean height of points within a column, then estimate 
volume
    summ <- dat_xy_bin %>%
    group_by(bin_x, bin_y) %>%
    summarise(volume = mean(Z, na.rm = TRUE) * voxel_dim^2)

    # Sum of volumes
    vol <- sum(summ$volume, na.rm = TRUE)