lidR custom grid metric returning infinite values

Question

I am trying to calculate the proportion of points 0.3-2.5m above ground relative to the number of ground points + all points up to 2.5m above ground using lidR grid_metrics(). My custom grid_metrics() calculation works on the lidR demo data provided with the package ("MixedConifer.laz") but fails to work with the data I'm interested in using (New Brunswick, Canada LiDAR tiles).

More specifically, this code successfully produces the metric:

require(lidR)
library(raster)

# read in demo file from lidR package
LASfile <- system.file("extdata", "MixedConifer.laz", package ="lidR")
las <- readLAS(LASfile)

shrub_cover = function(z) {
  n_shrub = sum(z >= 0.3 & z <= 2.5)
  n_veg = sum(z <= 2.5)
  return(list(ratio = n_shrub/n_veg))
}

shrub_cover_metric = grid_metrics(las, ~shrub_cover(Z), res = 10)
writeRaster(shrub_cover_metric, "data/interim/shrub_cover_demo_tile.tif", overwrite=TRUE)

But this script does not... the data I'm working with can be accessed here

las <- readLAS("C:/Users/pende/Downloads/nb_2018_2573000_7416000.laz")

shrub_cover = function(z) {
  n_shrub = sum(z >= 0.3 & z <= 2.5)
  n_veg = sum(z <= 2.5)
  return(list(ratio = n_shrub/n_veg))
}

shrub_cover_metric = grid_metrics(las, ~shrub_cover(Z), res = 10)
writeRaster(shrub_cover_metric, "data/interim/shrub_cover_nb_2018_2573000_7416000.tif", overwrite=TRUE)

I've tried removing duplicate points in the nb_2018_2573000_7416000.laz data, but this doesn't change the output result of a raster full of -Inf data.

What is causing this issue? The grid_metrics() function is the same...

Answer 1

The data I was working with contained raw elevation values and not heights above ground. To solve my problem, I needed to normalize my points. This is the script that is now working:

require(lidR)
library(raster)

las <- readLAS("C:/Users/pende/Downloads/nb_2018_2573000_7416000.laz")
nlas <- normalize_height(las, tin())

shrub_cover = function(z) {
  n_shrub = sum(z >= 0.3 & z <= 2.5)
  n_veg = sum(z <= 2.5)
  return(list(ratio = n_shrub/n_veg))
}

shrub_cover_metric = grid_metrics(nlas, ~shrub_cover(Z), res = 10)
writeRaster(shrub_cover_metric, "data/interim/shrub_cover_nb_2018_2573000_7416000.tif", overwrite=TRUE)