1

I would like to convert this point cloud from US feet to meters. Metadata of the point cloud here.

I firstly though to convert it properly by re-projecting to a metric CRS but this caused me a lot of trouble. See the following posts (in their order of publication) :

But basically, what I want is just metric point clouds! So I came to the idea of simply converting the XYZ coordinates of the point cloud using las_rescale() from lidR R package:

library(lidR)
tmp_las <- readLAS('./USGS_LPC_CA_FEMA_R9_Russian_2017_LAS_2018.laz')
tmp_las_RS <- las_rescale(tmp_las,xscale = 0.3048 , yscale = 0.3048, zscale = 0.3048)
#> X coordinates were moved by 0.0762 on average
#> Y coordinates were moved by 0.0762 on average
#> Z coordinates were moved by 0.0762 on average

But the headers don't take it into account

tmp_las_RS
#> class        : LAS (v1.2 format 3)
#> memory       : 618.7 Mb 
# extent       : 6268117, 6270217, 1998802, 2000902 (xmin, xmax, ymin, ymax)
#> coord. ref.  : +proj=lcc +lat_0=37.6666666666667 +lon_0=-122 +lat_1=39.8333333333333 +lat_2=38.3333333333333 +x_0=2000000.0001016 +y_0=500000.0001016 +ellps=GRS80 +units=us-ft +no_defs 
#> area         : 4.41 kus-ft²
#> points       : 7.05 million points
#> density      : 1.6 points/us-ft²

Is there a better way to do it through lidR?

1
  • Welcome to GIS SE! As a new user be sure to take the Tour to learn about our focussed Q&A format. It emphasizes that there should be only one question asked per question. I recommend reviewing meta.gis.stackexchange.com/a/3353/115 for tips on how to structure a good question.
    – PolyGeo
    Oct 5, 2020 at 4:52

1 Answer 1

2

las_rescale() does not perform the task you think it does. You are not the first one to be confused. See this comprehensive answer: Rescaling and reoffsetting a point cloud with lidR?

Considering how a LAS file stores the data (header, scale factor, offset) it is actually difficult to perform what you want to do because you need to: (1) recompute the coordinates (2) quantize the coordinates (3) update the header (4) update the bounding box. All manually.

But in the version 3.1 things become easier (and better documented). See this thread. So from version 3.1 you can do

ft2m = 1200/3937
las$X = las$X * ft2m
las$Y = las$Y * ft2m
las$Z = las$Z * ft2m

But this does not update the CRS magically. So you must erase the CRS because, according to your previous posts, you don't know what CRS to input.

crs(las) <- sp::CRS()

To install lidR 3.1 which won't be released before Jan 2021.

remotes::install_github("Jean-Romain/lidR", ref = "devel")
2
  • 2
    US survey feet to meters is 1200/3937, not 0.3048, so 0.3048006 (at a minimum). It matters because the state plane coordinates are in the millions.
    – mkennedy
    Oct 4, 2020 at 17:51
  • Thank you for the comment. Definitively true. Using 0.3048 introduces several meters of error.
    – JRR
    Oct 4, 2020 at 18:00

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.