Create a function that can process a LAScatalog

The following demonstrates how to write your own functions that are fully applicable on a wide catalog of point clouds and based on the available lidR tools. We will create a simple filter_noise function. This example should not be considered as the reference method for filtering noise, but rather as a demonstration to help understand the logic behind the design of lidR, and as a full example of how to create a user-defined function that is fully operational. The code design is the one used internally in lidR and relies on S3 method dispatch. For more details, we recommend reading the chapter about S3 method dispatch from the Advanced R book.

Create a generic filter_noise function

First we create a generic function filter_noise that will be usable on different classes.

filter_noise = function(las, ...)
  UseMethod("filter_noise", las)

Create a filter_noise for LAS objects

A simple (perhaps too simplistic) way to detect outliers is to measure the 95th percentile of height in 10 x 10-m pixels (area-based approach) and then remove the points that are above the 95th percentile in each pixel plus, for example, 20%. This can easily be built in lidR using grid_metrics, merge_spatial and filter_poi, and should work either on a normalized or a raw point cloud. Let’s create a function method filter_noise for LAS objects:

filter_noise.LAS = function(las, sensitivity)
  p95 <- grid_metrics(las, ~quantile(Z, probs = 0.95), 10)
  las <- merge_spatial(las, p95, "p95")
  las <- filter_poi(las, Z < p95*sensitivity)
  las$p95 <- NULL

This function is fully functional on a point cloud loaded in memory

las <- readLAS("file.las")
las <- filter_noise(las, sensitivity = 1.2)
writeLAS(las, "denoised-file.las")

Extend the filter_noise function to a LAScatalog

Users can access the catalog processing engine with the function catalog_apply i.e. the engine used internally. It can be applied to any function over an entire catalog. Here we will apply our custom filter_noise function. To use our function filter_noise on a LAScatalog we must create a compatible function (see documentation of catalog_apply). In the lidR package we usually create an intermediate method. Here filter_noise for LAScluster objects (see also the documentation for catalog_apply):

filter_noise.LAScluster = function(las, sensitivity)
  # The function is automatically fed with LAScluster objects
  # Here the input 'las' will a LAScluster

  las <- readLAS(las)                          # Read the LAScluster
  if (is.empty(las)) return(NULL)              # Exit early (see documentation)
  las <- filter_noise(las, sensitivity)      # Filter the noise
  las <- filter_poi(las, buffer == 0)           # Don't forget to remove the buffer
  return(las)                                  # Return the filtered point cloud

This function can be used in catalog_apply. We can then create a method filter_noise for a LAScatalog:

filter_noise.LAScatalog = function(las, sensitivity)
   catalog_apply(las, filter_noise, sensitivity = sensitivity)

And it just works. This function filter_noise is now fully compatible with the catalog processing engine and supports all the options of the engine.

myproject <- catalog("folder/to/lidar/data/")

opt_filter(myproject)       <- "-drop_z_below 0"
opt_chunk_buffer(myproject) <- 10
opt_chunk_size(myproject)   <- 0
opt_cores(myproject)        <- 2
opt_output_files(myproject) <- "folder/to/lidar/data/denoised/{ORIGINALFILENAME}_denoised"

output <- filter_noise(myproject, tolerance = 1.2)

Finalize the functions

As is, the function filter_noise.LAScatalog is not actually complete. Indeed:

  1. The processing options were not checked. For example, this function should not allow the output to be returned into R otherwise the whole point cloud will be returned.
  2. The output is a list of written files that can be simplified into a LAScatalog.

In lidR the functions usually look like this:

filter_noise.LAScatalog = function(las, sensitivity)
   opt_select(las) <-  "*"         # Do not respect the select argument by overwriting it

   options <- list(
       need_output_file = TRUE,    # Throw an error if no output template is provided
       need_buffer = TRUE,         # Throw an error if buffer is 0
       automerge = TRUE)           # Automatically merge the output list (here into a LAScatalog)
   output  <- catalog_apply(las, filter_noise, sensitivity = sensitivity, .options = options)

Now you know how to build your custom functions that work either on a LAS or a LAScatalog object.