2

I have a big point cloud (3.5 GB) and want to create a big DTM from it using R. In the "old" lidR version this was possible by turning on automerge = TRUE. If the tiles could not be merged, a virtual raster was returned. If they could, a raster object was returned. I might be wrong, because I am new to the terra package, but this did not work for me anymore. I tried making my DTM with the following code:

# load data
ctg <- readTLSLAScatalog(ctg_path)

# set options
opt_output_files(ctg) <- paste0(dtm_path, "/tile_{ID}")
opt_chunk_size(ctg) <- 20
opt_chunk_buffer(ctg) <- 5
opt_merge(ctg) <- TRUE

# my function
make_dtm <- function(chunk){ 
  las <- readLAS(chunk)
  if (is.empty(las)) return(NULL)
  las <- classify_ground(las, csf(class_threshold = 0.3, cloth_resolution = 0.1))
  dtm <- rasterize_terrain(las, res = 0.01, algorithm = tin())
  dtm <- crop(dtm, ext(chunk)) # remove buffer
  return(dtm)
}

# execution
dtm <- catalog_apply(ctg, make_dtm)
writeRaster(dtm, file.path(dtm_path, "final.tif"), overwrite = T)

But the return is always a list with the filenames. I guess I am just not up-to-date anymore on lidR, but how can I use catalog_apply() to obtain a SpatRaster which I can save as a single tif?

Edit: The new working solution:

# load data
ctg <- readTLSLAScatalog(ctg_path, chunk_size = 20, chunk_buffer = 5)

# my functions
make_dtm <- function(chunk){ 
  las <- readLAS(chunk)
  if (is.empty(las)) return(NULL)
  las <- classify_ground(las, csf(class_threshold = 0.3, cloth_resolution = 0.3))
  dtm <- rasterize_terrain(las, res = 0.01, algorithm = tin())
  dtm <- crop(dtm, ext(chunk)) # remove buffer
  return(dtm)
}
make_dtm_ctg <- function(las) {
  options <- list(
    need_buffer = TRUE, # buffer necessary
    automerge = TRUE) # combine outputs
  output <- catalog_apply(las, make_dtm, .options = options)
  return(output)
}

# execution
dtm <- make_dtm_ctg(ctg)
writeRaster(dtm, file.path(dtm_path, "DTM.tif"))

1 Answer 1

1

Nothing changed on this aspect. The option auto-merge is still supported but you did not use it here.

options = list(automerge = TRUE)
dtm <- catalog_apply(ctg, make_dtm, .options = options)

It is different from opt_merge. The idea is that catalog_apply is not intended to be used bare bone but is expected to be used inside a wrapper function. For example:

make_dtm <- function(las){ 
  if (is(las, "LAS")) {
    las <- classify_ground(las, csf(class_threshold = 0.3, cloth_resolution = 0.1))
    dtm <- rasterize_terrain(las, res = 0.01, algorithm = tin())
    return(dtm)
  } else if (is(las, "LAScluster")) {
    las <- readLAS(chunk)
    if (is.empty(las)) return(NULL)
    dtm <- make_dtm(las)
    dtm <- crop(dtm, ext(chunk)) # remove buffer
    return(dtm)
  } else if (is(las, "LAScatalog")) {
    options = list(automerge = TRUE)
    dtm <- catalog_apply(ctg, make_dtm, .options = options)
    return(dtm)
  }
}

There are options for the developers who write the function and options for the users who use the function. Here you are the developer and with the options automerge you state I want catalog_apply to automatically merge the chunks in a standard way. If you do not add this options, the function will never do it, no mater what is in opt_merge() because if you do not use automerge we assume that your output is complex and cannot be automatically be merged with simple standard code so you must perform this task yourself. With opt_merge() you are on the user side and you state I want do disable the automerge option to get a raw list. opt_merge() is more for debugging purposes.

Also from lidR v4.0.0 you have catalog_map() that simplifies catalog_apply() for simple cases like your. Something like that should work.

make_dtm <- function(las){ 
  if (is(las, "LAS")) {
    las <- classify_ground(las, csf(class_threshold = 0.3, cloth_resolution = 0.1))
    dtm <- rasterize_terrain(las, res = 0.01, algorithm = tin())
    return(dtm)
  } else if (is(las, "LAScatalog")) {
    dtm <- catalog_map(ctg, make_dtm)
    return(dtm)
  }
}
2
  • Thanks a lot! I wasn't aware that the opt_...() options only work if explicitly included in the function. I always get confused whether to use .options or these opt_...()-functions...
    – Zoe
    Commented Jun 22, 2022 at 10:58
  • 1
    Think this way: as the developer of the function how I want the function to behave (e.g. fail if buffer = 0, automerge, etc.), and as a user of the function what can I control (chunk size, buffer, parallelization, etc.)
    – JRR
    Commented Jun 22, 2022 at 12:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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