TransWikia.com

Reprojecting a lascatalog

Geographic Information Systems Asked on July 21, 2021

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)

One Answer

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.

Correct answer by JRR on July 21, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP