2

I'm using the approach described here to change LiDAR point cloud projection in R: Can I re-project an LAS file in LidR But, what about lidR catalog reprojection?

If I used the same reproducible dataset than in the post Can I re-project an LAS file in LidR:

library(lidR)
LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")

# las approach
las = readLAS(LASfile)
las
las2 = spTransform(las, sp::CRS(SRS_string = "EPSG:26918"))
class(las2)
# [1] "LAS"
# attr(,"package")
# [1] "lidR"

# lascatalog approach
ctg <- catalog(LASfile)
ctg2 <- spTransform(ctg, sp::CRS(SRS_string = "EPSG:26918"))
class(ctg2)
# [1] "SpatialPolygonsDataFrame"
# attr(,"package")
# [1] "sp"

The result when using the same command with a lascatalog is a SpatialPolygonsDataFrame and not a lascatalog... Is there a way to reproject lascatalog while remaining in R?

I don't like the las2las solution because you may only work with projection already in pre encoded by Martin within the tool (or maybe I'm missing something with the use or las2las)

1 Answer 1

3

Using spTransform on a LAScatalog should fail. But R being R, it surprisingly worked but not the way you were expecting.

It is not a good idea IMO to re-project a whole point cloud in R. So it is even worse to try to re-project a whole collection of files. The reason is that lidR needs to load the whole point cloud in memory to process it (including all the attributes) and this is done very sub-optimally (see also this vignette). Optimally to re-project you need to load a single point at a time. This is what las2las from LAStool does. It streams from an input file to a output file, processing a single point at a time and is thus way much efficient. Moreover, because lidR uses sp and rgdal and requires at least one copy of the point cloud to handle the format conversion between packages...

This is why re-projection in lidR is not a very developed tools. It exists for convenience on small point cloud but I don't have any interest in making a better tools since it already exists some much better open source tools (LAStools, PDAL). You won't even find spTransform() in lidR's manual.

That being said if you really want to apply spTransform() on a LAScatalog you can design your own function using the catalog processing engine (using the version >= 3.1.2 because I fixed a bug recently. Something like that should work:

reproject_chunk = function(cl, crs, scale, xoffset, yoffset)
{
  las = readLAS(cl)
  if (is.null(las)) return(NULL)
  las = spTransform(las, crs, scale = scale, xoffset = xoffset, yoffset = yoffset)
  return(las)
}

crs = sp::CRS(SRS_string = "EPSG:26918")
opt_chunk_buffer(ctg) <- 0
opt_chunk_size(ctg) <- 0
opt_output_files(ctg) <- "folder/reprojected/{ORIGINALFILENAME}"
reprojected_ctg = catalog_apply(ctg, reproject_chunk, crs = crs, scale = 0.001, xoffset = 12345, yoffset = 12345)

Pay attention to the values of scale, xoffset and yoffset. Depending of your original CRS and target CRS it may be important to specify them correctly. See also this question.

Your Answer

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

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