Cortando horizontalmente la nube de puntos en un cilindro

Question

Tengo algunas nubes de puntos obtenidas de escaneo láser terrestre móvil de parcelas circulares de muestreo en entornos boscosos. Estoy interesado en estudiar la complejidad estructural del sotobosque del bosque. He utilizado la función clip_circle() para recortar las nubes de puntos al radio de 15 metros de la parcela de muestreo.

Ahora, me gustaría cortar horizontalmente la nube de puntos para obtener los primeros dos metros de puntos desde el suelo del bosque. Por lo tanto, me gustaría extraer un cilindro de 15 metros de radio y dos metros de profundidad, posicionado en el suelo.

Logré hacer esto en CloudCompare y volver a importar la nube a R de esta manera:

Sin embargo, me gustaría hacer todo el proceso en R para ahorrar tiempo cuando tenga más nubes que procesar.

Answer 1

Existen varias opciones para lograr esta tarea (asumiendo que la nube de puntos está normalizada en altura)

La opción más simple es leer solo el punto por debajo de 2 metros y luego recortar

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

Otra opción es hacer lo mismo desde un LAScatalog. Esto es más adecuado para múltiples consultas y más optimizado en memoria si el archivo original es particularmente grande

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

La última opción y la más trivial pero la menos eficiente es filtrar los puntos de interés con filter_poi()

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

Answer 2

No es seguro en memoria como catalog_apply, pero puedes procesar escaneos en paralelo como se muestra a continuación.

 #instalar paquetes y librerías con pacman---------------------------------------
if (!require("pacman")) install.packages("pacman"); library(pacman)
p_load(lidR, TreeLS, doParallel, parallel, here)

# leer archivos en una lista de data frames-------------------------------------
dat = as.data.frame(list.files(here("ruta/hacia/los/datos/normalizados"), pattern="laz$"))
seq_id_all = seq_along(1:nrow(dat))

# definir función ------------------------------------------------------------
my_fun <- function(...) {
 i <- (...)
# cargar todas las librerías dentro de la función
library(TreeLS); library(lidR); library(here); library(tidyverse)  
zeb_data <- readTLS(here("ruta/hacia/los/datos/normalizados", dat[i,]))
(a <- zeb_data %>% filter_poi(Z >= 0 & Z <1)) # filtrar por Z (0 a 1m)
(b <- zeb_data %>% filter_poi(Z >= 1 & Z <2))

# escribir salida
writeLAS(a, file= here("ruta/de/salida", paste0(dat[i,], "_0to1m.laz")))
writeLAS(b, file= here("ruta/de/salida", paste0(dat[i,], "_1to2m.laz")))
}

# obtener el número de núcleos disponibles para el procesamiento paralelo--------
n_cores = detectCores(logical = TRUE)
cl = makeCluster(n_cores-1)  # dejar disponible 1 núcleo
registerDoParallel(cl)  # registrar el cluster

# exportar la función y 'dat' al cluster y ejecutarlo-----------------------------
clusterExport(cl, list('my_fun','dat'))
# crear un directorio de salida
dir.create(here('ruta/de/salida'))
system.time(results <- c(parLapply(cl, seq_id_all, fun=my_fun)))
stopCluster(cl)