lasdatum issues

110 views
Skip to first unread message

Michael Perdue

unread,
Oct 3, 2022, 2:58:39 PM10/3/22
to LAStools - efficient command line tools for LIDAR processing
I'm trying to isolate an odd issue that I'm having while trying to apply a seven parameter transformation to data that is in longlat coordinates. When I try to transform a file that is in -longlat to utm, the program generates incorrect answers. But if I first reproject the file to ECEF, then the answers are correct.

What really baffles me is that when I manually extract a single point from a text file as geographic coordinates and run it through the flow, then it generates the correct answer. So, I'm a little stumped. Below is a series of commands, their output and comments that illustrate what I'm seeing. The laz file being used can be retrieved from here:

#  Convert the laz file to txt
MichaelPerdue@Hactar /cygdrive/e/Projects/6322267a/Output/testing
> las2txt -parse txyz -i geo_ortho.laz -o geo_ortho.txt

 MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> head geo_ortho.txt
348064598.000008 -84.687553064 32.848354023 377.080
348064598.000008 -84.687555197 32.848351052 373.866
348064598.000008 -84.687557562 32.848347757 370.303
348064598.000008 -84.687559393 32.848345206 367.543
348064598.000008 -84.687366817 32.848514384 387.439
348064598.000008 -84.687368447 32.848511366 382.943
348064598.000014 -84.687551102 32.848355418 376.728
348064598.000014 -84.687364777 32.848515872 387.157
348064598.000020 -84.687549351 32.848356741 376.613
348064598.000020 -84.687551130 32.848354251 373.904

#  Transform the longlat file to utm with a 7param Helmert Transformation 
 MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> lasdatum -pvr -seven 1.0026,-1.9083,-0.5164,-0.0265849,-0.0018569,-0.0110891,-0.0012048 -longlat -target_epsg 6345 -i geo_ortho.laz -o utm_ortho.laz
shifted 401149 points in 'geo_ortho.laz'. took 0.992 sec.

#  Convert to txt to check the answers 
 MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> las2txt -parse txyz -i utm_ortho.laz -o utm_ortho.txt

# These are wrong! X & Y are wrong, but elevation looks right
MichaelPerdue@Hactar  /cygdrive/e/Projects/6322267a/Output/testing
> head utm_ortho.txt
348064598.000008 716477.294692753 3636887.548518884 378.469
348064598.000008 716477.102246283 3636887.214676205 375.255
348064598.000008 716476.888870044 3636886.844427837 371.692
348064598.000008 716476.723672583 3636886.557780305 368.932
348064598.000008 716477.159274430 3636888.533745605 388.828
348064598.000008 716477.014027458 3636888.195721587 384.332
348064598.000014 716477.474961051 3636887.707239874 378.117
348064598.000014 716477.346617516 3636884.407972819 388.546
348064598.000020 716477.635653146 3636887.857543868 378.002
348064598.000020 716477.475175109 3636887.577767536 375.293

# Extract the first long lat point from txt file to a laz file  
MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> echo '348064598.000008 -84.687553064 32.848354023 377.080' | txt2las -parse txyz -itxt -stdin -o foo.laz -set_scale 0.000000001 0.000000001 0.001
done with 'foo.laz'. total time 0.015 sec.

# Perform the same 7p Helmert transformation to the extracted point 
MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> lasdatum -pvr -seven 1.0026,-1.9083,-0.5164,-0.0265849,-0.0018569,-0.0110891,-0.0012048 -longlat -target_epsg 6345 -i foo.laz -o bar.laz
shifted 1 point in 'foo.laz'. took 0.015 sec.

# Not sure why, but now it works. These are the right coordinates!
MichaelPerdue@Hactar  r /cygdrive/e/Projects/6322267a/Output/testing
> las2txt -parse txyz -stdout -i bar.laz
348064598.000008 716412.870183312 3636844.598845925 378.469

# Now we convert the longlat file to ECEF  
MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> las2las -longlat -target_ecef -i geo_ortho.laz -o ecef_ortho.laz -target_precision 0.001 -target_elevation_precision 0.001

# Do the 7p helmert transformation on the ECEF file
MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> lasdatum -pvr -seven 1.0026,-1.9083,-0.5164,-0.0265849,-0.0018569,-0.0110891,-0.0012048 -ecef -target_epsg 6345 -i ecef_ortho.laz -o utm_ortho.laz
shifted 401149 points in 'ecef_ortho.laz'. took 1.106 sec.

# Convert the file to txt for inspection
MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> las2txt -parse txyz -i utm_ortho.laz -o utm_ortho.txt

# Now the answers are right!
MichaelPerdue@Hactar   /cygdrive/e/Projects/6322267a/Output/testing
> head utm_ortho.txt
348064598.000008 716412.871 3636844.599 378.469
348064598.000008 716412.677 3636844.265 375.255
348064598.000008 716412.464 3636843.895 371.692
348064598.000008 716412.299 3636843.608 368.932
348064598.000008 716429.914 3636862.764 388.828
348064598.000008 716429.770 3636862.426 384.332
348064598.000014 716413.050 3636844.758 378.117
348064598.000014 716430.102 3636862.933 388.547
348064598.000020 716413.211 3636844.908 378.001
348064598.000020 716413.051 3636844.628 375.293

What am I missing???

Cheers,

Mike

Dave Stoll

unread,
Oct 7, 2022, 1:01:23 PM10/7/22
to LAStools - efficient tools for LiDAR processing

Hi Mike,
Having a hard time understanding your issue. This is near Manchester, Georgia?
Where did you get this data?
7-parameter transforms are usually meant for Datum to Datum. I've never heard of using it to transform from Lat/Long to Northing/Easting.(GCS to PCS)

LASInfo reveals:
reporting minimum and maximum for all LAS point record entries ...
  X           -87955709  -85772434
  Y            47795519   49629756
  Z              358566     392250
scale factor x y z: 0.000000001 0.000000001 0.001 (Huh?)
offset x y z: -84.599999999999994 32.799999999999997 0 (Huh?)
offset to point data: 43476 (Huh?)
The X & Y scale factors look wrong, and I totally don't understand the offsets.
Dave

Michael Perdue

unread,
Oct 11, 2022, 3:18:59 PM10/11/22
to last...@googlegroups.com
Hi Dave,

This is a datum to datum transformation. More specifically, it is an ITRF2000 epoch 1997.0 to NAD83(2011) datum transformation. The full datum transformation is a 14 parameter datum transformation, but since both epoch's are fixed (delta T is fixed, so all rotational/translational rates become constants), I have manually reduced the 14 parameters to 7 parameters so that they will work with lasdatum.

This data is produced from raw sensor observations using Cloudpro (Leica's point cloud processing package for their legacy ALS line of instruments). The scale factors and offsets are determined by Cloudpro on export. Nine decimal places in decimal degrees is probably a little excessive and it certainly didn't produce nice clean offsets, but I have no control over that. I would be the first to agree that they are ugly, but they are not wrong or illegal (violate the LAS specification) and they should not affect the datum transformation. The large offset to point data is because Leica likes to write a bunch of VLRs in the files. At Least one of those VLR's stores a histogram of intensity distribution and it takes up some room.

The helmert transformation is not being used to go from lat/long to Northing/Easting. The datum transformation is supposed to be applied to ECEF coordinates, and the proper workflow for this transformation would be:
  1. Unproject the data to ECEF,
  2. apply the helmert transformation
  3. reproject the coordinates to the destination projection
  4. write the points to a new file.
My understanding is that this is what lasdatum does (or is supposed to do). If I manually complete the first step with las2las (conversion to ECEF) then everything works properly. But if I feed the lat/long file directly, things go sideways and it is possible that the unusual offsets/scale factors are triggering the issue and potentially explain why it works when I manually extract one of the points. But if these scale factors and offsets are not illegal or wrong, then wouldn't this constitute a bug?

Cheers,

Mike

--
Download LAStools at
https://rapidlasso.de
Manage your settings at
https://groups.google.com/g/lastools/membership
---
You received this message because you are subscribed to the Google Groups "LAStools - efficient tools for LiDAR processing" group.
To unsubscribe from this group and stop receiving emails from it, send an email to lastools+u...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/lastools/5215ffee-1e4a-4379-af40-9d04e3fdd725n%40googlegroups.com.
Reply all
Reply to author
Forward
0 new messages