2

following up on this question: Capture and parse output of lascheck from lidR... (while lidR 3 is in development, I thought I'd continue with a workaround.)

I'm using catalog_apply to perform lascheck on individual las files in the catalog. I'm getting the right output for warnings/errors, but I can't identify the name of the las file (chunk) because although the warnings/errors are there in the output of lascheck, the name isn't.

However, a list of names is available in the LASCatalog. I would like to extract the corresponding name, but I can't increment a variable to count (for reasons of scope I think). I.e. capture.output(ctg@data$filename[x]) where x is position of the name in the list.

E.g. the following code does return the warnings/errors but x gets stuck at 1.

library(lidR)

ctg <- readLAScatalog("D:/Projects2020/Sherbrooke/") #, progress = TRUE, list.files())

x <- 0

errwrn <- function(chunk)
{
  las = readLAS(chunk)
  if (is.empty(las)) return(NULL)

  las@header@PHB$`Point Data Format ID` <- 12L # Introduce errors.

  ew <- capture.output(lascheck(las))
  err_pos  = grep("+2717", ew)  
  errs = ew[err_pos]
  warn_pos = grep("+26A0", ew)
  warns = ew[warn_pos]

  x <- x + 1
  lasName <- capture.output(ctg@data$filename[x])
  ews = c(x, lasName, errs, warns)

  return (ews)
}

lascheckList <- catalog_apply(ctg, errwrn)

Can anyone see how I could increment x? Or maybe I should use a loop to go through each file in the directory, bypassing the use of a catalog?

3

There is no documented way for what you are trying to do simply because it does not always makes sense. In your case you want to process by file with no buffer so it actually makes sense because what you are actually doing is a for loop. Below several options. I'm assuming you are processing with:

opt_chunk_size(ctg) <- 0
opt_chunk_buffer(ctg) <- 0

Otherwise nothing makes sense in your example. By default the processing options include a buffer and thus you are not checking each file but an object built out of several files.

Attribute the names after the process

Processing order is preserved so the following is valid

lascheckList <- catalog_apply(ctg, errwrn)
names(laschekList) <- basename(ctg$filename)

Extract data from the LAScluster

This returns all the files that will be read. In this specific case a single one

errwrn <- function(chunk)
{
  las = readLAS(chunk)
  if (is.empty(las)) return(NULL)
  file <- basename(chunk@files)
  ...
}

Regular loop

You are processing by file without buffer so actually you don't need the engine you can use a regular loop to go through each file.

errwrn <- function(f)
{
  las = readLAS(f)
  ...
}

lascheckList <- lapply(ctg$filename, errwrn)
names(laschekList) <- basename(ctg$filename)
| improve this answer | |
  • Okay I see; this is very helpful. For this case then, a regular loop is easier (although we don't get a visual with the plot). I modified your suggestion slightly so that results could be stored in a file. – Ray J Mar 3 at 18:32
  • Yes, in a regular loop you loose the monitoring (yet you can add a progress bar) and the multi-core capabilities but in this case it is not critical. – JRR Mar 3 at 18:46
0

Following JRR suggestion of using a regular loop, I wrote the code below which produces a text file listing errors and warnings and what type (i.e. just error, just warning, both, or none). It is inefficient in it's present state (but does what I need)... a shapefile of the catalog index was exported outside of this script; I used a semicolon as a column separator so column headers are added in Excel and I save the file as .xls. With the shapefile and the .xls joined in ArcMap you can see which tiles give problems.

library(lidR)
library(pbapply)

options(max.print=1000000)
LASlist <- list.files(path = "D:/Projects2020/Sherbrooke/", pattern = "*.las", full.names = TRUE, recursive = FALSE)
file.create("D:/Projects2020/Sherbrooke/ews.txt")

pblapply(LASlist, function(lasName)
{
  lasF <- readLAS(lasName, select = "*", filter = "-drop_withheld -drop_overlap -keep_class 2 3 4 5 9")
  if (is.empty(lasF)) return(NULL)
  #lasF@header@PHB$`Point Data Format ID` <- 12L # Introduce errors for testing.
  #lasF$ReturnNumber[1] <- 0L # Introduce warnings for testing.
  ew <- capture.output(lascheck(lasF))
  err_pos <- grep("+2717", ew)  
  errs <- ew[err_pos]
  warn_pos <- grep("+26A0", ew)
  warns <- ew[warn_pos]

  ewType <- "n"
  if (length(err_pos) > 0 & length(warn_pos) == 0)
  {
    ewType <- "e"
  }
  if (length(warn_pos) > 0 & length(err_pos) == 0)
  {
    ewType <- "w"
  }
  if (length(warn_pos) > 0 & length(err_pos) > 0)
  {
    ewType <- "ew"
  }

  ews <- (paste(c(lasName, ewType, errs, warns), collapse = ";"))
  cat(ews, file = "D:/Projects2020/Sherbrooke/ews.txt", sep = " ", fill = TRUE, labels = NULL, append = TRUE)
  return(ews)
})
| improve this answer | |

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.