5

The following script uses lidR to segment trees and output tree hulls (i.e. polygons representing individual tree canopies). Input data are a directory of laz files and a catalog is used to manage the processing.

library(lidR)
library(rgdal)

my_tree_detection_method <- function(cluster, ws)
{
  las <- readLAS(cluster)
  if (is.empty(las)) return(NULL)

  las_n <- lasnormalize(las, tin())
  chm  <- grid_canopy(las_n, 0.25, pitfree(c(0,2,5,10,15), c(0,1), subcircle = 0.2))
  algo = watershed(chm, th = 4)
  trees <- lastrees(las_n, algo)

  # Remove the buffers
  trees_no_buffer <- lasfilter(trees, buffer == 0)
  return(trees_no_buffer)
}

ws <- "/path/to/laz/directory"
ctg <- readLAScatalog(ws)
lidR:::catalog_laxindex(ctg) 

opt_chunk_buffer(ctg) <- 10
opt_chunk_size(ctg) <- 200 

opt <- list(need_buffer = TRUE) 
output <- catalog_apply(ctg, my_tree_detection_method, ws = 5, .options = opt)
output <- do.call(rbind, output)

hulls  = tree_hulls(output, func = .stdmetrics)
writeOGR(hulls, dsn = "/path/to/dir", layer = "hulls", driver="ESRI Shapefile")

enter image description here

When I write the hulls to shapefile to inspect the results, the hulls appear very distorted. I would expect thousands of small polygons representing tree canopies. What might be causing this distortion?

The following from the documentation is what I would expect to see:

enter image description here

Aaron
  • 51,658
  • 28
  • 154
  • 317

1 Answers1

5

You found why there is no lastrees methods available for a LAScatalog yet in lidR. You did not consider the continuity of the tree IDs. In each chunk the trees are labeled from 1 to n. So in chunk 1 you have a tree labeled 1 and you have another tree labeled 1 in chunks 2,3,4, ...

At the end you merged everything and you have several points tagged 1 spread on the coverage. tree_hull takes all the points tagged 1 and makes the hull of the tree number 1. And this hull is meaningless.

Because the chunks are processed independently and potentially in parallel it is almost impossible to maintain data integrity. You must either find a method to maintain the continuity of the IDs (I'm trying for 2 years) or you must build a more complex function. The following is something I did for a demo. It computes the hulls for each chunk and using the tree tops it removes the hull of the trees that are outside the chunks (i.e in the buffer).

my_process = function(cl) {
  las = readLAS(cl)
  if (is.empty(las)) return(NULL)
  bbox = extent(cl)
  las = lasfilterduplicates(las)
  chm = grid_canopy(las, 0.5, p2r())
  ttops = tree_detection(las, lmf(3, 5))
  las = lastrees(las, dalponte2016(chm, ttops))
  p = tree_metrics(las)
  p = crop(p, bbox)
  m = tree_hulls(las, func = .stdtreemetrics)
  m = m[m$treeID %in% p$treeID,]
  return(m)
}

options = list(automerge = TRUE) # available from v2.2.0
m = catalog_apply(ctg, my_process, .options = options)
JRR
  • 9,389
  • 1
  • 13
  • 28