I have some point clouds from mobile terrestrial laser scanning of circular sampling plots in woodland environments. I am interested in studying the structural complexity of the woodland understorey. I have used the clip_circle() function to crop the point clouds to the 15 metre radius of the sampling plot.

Now, I would like to horizontally slice the point cloud so that I can obtain the first two metres of points from the woodland floor. Therefore, I would like to extract a 15-metre radius cylinder of two metres depth, positioned on the floor.

I managed to do this in CloudCompare and reimport the cloud back into R like this:

See image 2

However I would like to do the whole pipeline in R to save time when I have more clouds to process.


There are several options to achieve this task (assuming the point cloud is height normalized)

The simplest options is to read only the point below 2 meters and then clip

las = readLAS(file, filter = "-drop_z_above 2")
plot = clip_circle(las, x, y, radius)

On other option is to do the same from a LAScatalog. This is more suitable for multiple queries an more memory optimized if the original file is particularly big

ctg = readLAScatalog(file, filter = "-drop_z_above 2")
plot = clip_circle(ctg, x, y, radius)

The last option and the most trivial but the least efficient is to filter the points of interest with filter_poi()

las = readLAS(file)
las = filter_poi(las, Z <= 2)
plot = clip_circle(las, x, y, radius)

It's not memory-safe like catalog_apply, but you could process scans in parallel as shown below.

 #install packages and libs w/ pacman ----------------------------------------
if (!require("pacman")) install.packages("pacman"); library(pacman)
p_load(lidR, TreeLS, doParallel, parallel, here)

# read files into a list of data frames ------------------------------------
dat = as.data.frame(list.files(here("path/to/normalized/data"), pattern="laz$"))
seq_id_all = seq_along(1:nrow(dat))

# define function -----------------------------------------------------------
my_fun <- function(...) {
 i <- (...)
# load all libs into function
library(TreeLS); library(lidR); library(here); library(tidyverse)  
zeb_data <- readTLS(here("path/to/normalized/data", dat[i,]))
(a <- zeb_data %>% filter_poi(Z >= 0 & Z <1)) # filter by Z (0 to 1m)
(b <- zeb_data %>% filter_poi(Z >= 1 & Z <2))

# write output
writeLAS(a, file= here("output/path", paste0(dat[i,], "_0to1m.laz")))
writeLAS(b, file= here("output/path", paste0(dat[i,], "_1to2m.laz")))

# get number of available cores for parallel processing ------------
n_cores = detectCores(logical = TRUE)
cl = makeCluster(n_cores-1)  # leave 1 core availible
registerDoParallel(cl)  # register the cluster

# export function and 'dat' to the cluster and run ----------------
clusterExport(cl, list('my_fun','dat'))
# create a output directory
system.time(results <- c(parLapply(cl, seq_id_all, fun=my_fun)))

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.