las2las coodinate transformation

936 views
Skip to first unread message

Drexel Peter

unread,
Dec 5, 2019, 9:27:30 AM12/5/19
to LAStools - efficient tools for LiDAR processing
Hi Martin,
I tried to use

las2las -i SynthetischePunkte_neuer_rhein.laz -o 31254-SynthetischePunkte_neuer_rhein.laz -olaz -epsg 25832 -target_epsg 31254

to convert my EPSG 25832-Pointcloud to EPSG 31254.
The results are off by ~500m...
Thanks
Peter

Martin Isenburg

unread,
Dec 5, 2019, 10:56:58 AM12/5/19
to LAStools - efficient command line tools for LIDAR processing
Hello Peter,

because going from EPSG 25832 to EPSG 31254 means not only a reprojection but also a datum transform from ETRS99 (which uses the GRS1980 ellipsoid) to MGI (which uses the Bessel 1841 ellipsoid) things are a bit more complicated and need to be done in a number of steps. Below we walk you through it step by step with the example data you provided me.

 DATUM["European_Terrestrial_Reference_System_1989",
            SPHEROID["GRS 1980",6378137,298.257222101,
                AUTHORITY["EPSG","7019"]],
            TOWGS84[0,0,0,0,0,0,0],
            AUTHORITY["EPSG","6258"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4258"]],
 DATUM["Militar_Geographische_Institute",
            SPHEROID["Bessel 1841",6377397.155,299.1528128,
                AUTHORITY["EPSG","7004"]],
            TOWGS84[577.326,90.129,463.919,5.137,1.474,5.297,2.4232],
            AUTHORITY["EPSG","6312"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4312"]],

We do this in THREE steps:

(1) convert the data in EPSG 25832 to ECEF (Earth Centered Earth Fixed) but still with the ETRS99 datum.

las2las -i NeuerRhein.laz -epsg 25832 -target_ecef -odix _ecef -olaz

(2) Then apply the Helmert Transform that goes from ETSR99 to WGS84. This is a null transform because TOWGS84[0,0,0,0,0,0,0] so we skip it.  apply the ***inverse*** of the Helmert Transform that goes from MGI to WGS84 to go from WGS84 to MGI. This should be the inverse of TOWGS84[577.326,90.129,463.919,5.137,1.474,5.297,2.4232] aka the Helmert Transform -577.326,-90.129,-463.919,-5.137,-1.474,-5.297,-2.4232 (at least I hope that inverse is simply negated).

las2las -i NeuerRhein_ecef.laz -transform_helmert -577.326,-90.129,-463.919,-5.137,-1.474,-5.297,-2.4232 -odix _mgi -olaz

(3) Now that we are on the correct MGI datum with ECEF coordinates we project them into the EPSG 31254 projection with:

las2las -i NeuerRhein_ecef_mgi.laz -ecef -target_epsg 31254 -set_ogc_wkt -ocut 9 -odix _31254 -olaz

The use of seven parameter Helmert Transforms is of course not the best way to do a Datum transform as it is not so precise. The best way for doing this would be high-resolution datum shift grids (if available) but this is not yet implemented in LAStools.

==================================================

E:\software\LAStools\bin>lasinfo -i NeuerRhein.laz
lasinfo (191127) report for 'NeuerRhein.laz'
reporting all LAS header entries:
  file signature:             'LASF'
  file source ID:             0
  global_encoding:            16
  project ID GUID data 1-4:   00000000-0000-0000-0000-000000000000
  version major.minor:        1.4
  system identifier:          'LAStools (c) by rapidlasso GmbH'
  generating software:        'las2las (version 190623)'
  file creation day/year:     60/2017
  header size:                375
  offset to point data:       1154
  number var. length records: 2
  point data format:          6
  point data record length:   30
  number of point records:    0
  number of points by return: 0 0 0 0 0
  scale factor x y z:         0.01 0.01 0.01
  offset x y z:               500000 5200000 0
  min x y z:                  548728.04 5262303.85 439.77
  max x y z:                  548972.16 5262625.59 446.01
  start of waveform data packet record: 0
  start of first extended variable length record: 0
  number of extended_variable length records: 0
  extended number of point records: 155469
  extended number of points by return: 155321 148 0 0 0 0 0 0 0 0 0 0 0 0 0
variable length header record 1 of 2:
  reserved             0
  user ID              'LASF_Projection'
  record ID            34735
  length after header  40
  description          'by LAStools of rapidlasso GmbH'
    GeoKeyDirectoryTag version 1.1.0 number of keys 4
      key 1024 tiff_tag_location 0 count 1 value_offset 1 - GTModelTypeGeoKey: ModelTypeProjected
      key 3072 tiff_tag_location 0 count 1 value_offset 25832 - ProjectedCSTypeGeoKey: ETRS89 / UTM 32N
      key 3076 tiff_tag_location 0 count 1 value_offset 9001 - ProjLinearUnitsGeoKey: Linear_Meter
      key 4099 tiff_tag_location 0 count 1 value_offset 9001 - VerticalUnitsGeoKey: Linear_Meter
variable length header record 2 of 2:
  reserved             0
  user ID              'LASF_Projection'
  record ID            2112
  length after header  629
  description          'by LAStools of rapidlasso GmbH'
    WKT OGC COORDINATE SYSTEM:
    PROJCS["ETRS89 / UTM 32N",GEOGCS["ETRS89",DATUM["European_Terrestrial_Reference_System_1989",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","701
9"]],AUTHORITY["EPSG","6258"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4258
"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",9],PARAMETER["scale_factor",0.9996],PARAMETER["false_eastin
g",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","25832"]]
the header is followed by 2 user-defined bytes
LASzip compression (version 3.4r1 c3 50000): POINT14 3
reporting minimum and maximum for all LAS point record entries ...
  X             4872804    4897216
  Y             6230385    6262559
  Z               43977      44601
  intensity           5       6487
  return_number       1          2
  number_of_returns   1          2
  edge_of_flight_line 0          0
  scan_direction_flag 1          1
  classification      0          9
  scan_angle_rank   -28         21
  user_data          29        142
  point_source_ID   185        187
  gps_time 121970.520615 122644.581305
  extended_return_number          1      2
  extended_number_of_returns      1      2
  extended_classification         0      9
  extended_scan_angle         -4667   3500
  extended_scanner_channel        0      0
number of first returns:        155321
number of intermediate returns: 0
number of last returns:         155301
number of single returns:       155153
overview over extended number of returns of given pulse: 155153 316 0 0 0 0 0 0 0 0 0 0 0 0 0
histogram of classification of points:
            1095  never classified (0)
             107  unclassified (1)
          103916  ground (2)
             419  low vegetation (3)
           49932  water (9)

==================

E:\software\LAStools\bin>lasinfo -i NeuerRhein_31254.laz
lasinfo (191127) report for 'NeuerRhein_31254.laz'
reporting all LAS header entries:
  file signature:             'LASF'
  file source ID:             0
  global_encoding:            16
  project ID GUID data 1-4:   00000000-0000-0000-0000-000000000000
  version major.minor:        1.4
  system identifier:          'LAStools (c) by rapidlasso GmbH'
  generating software:        'las2las (version 191127)'
  file creation day/year:     60/2017
  header size:                375
  offset to point data:       1161
  number var. length records: 2
  point data format:          6
  point data record length:   30
  number of point records:    0
  number of points by return: 0 0 0 0 0
  scale factor x y z:         0.01 0.01 0.01
  offset x y z:               0 200000 0
  min x y z:                  -51658.35 263964.79 391.41
  max x y z:                  -51417.13 264288.64 397.65
  start of waveform data packet record: 0
  start of first extended variable length record: 0
  number of extended_variable length records: 0
  extended number of point records: 155469
  extended number of points by return: 155321 148 0 0 0 0 0 0 0 0 0 0 0 0 0
variable length header record 1 of 2:
  reserved             0
  user ID              'LASF_Projection'
  record ID            34735
  length after header  40
  description          'by LAStools of rapidlasso GmbH'
    GeoKeyDirectoryTag version 1.1.0 number of keys 4
      key 1024 tiff_tag_location 0 count 1 value_offset 1 - GTModelTypeGeoKey: ModelTypeProjected
      key 3072 tiff_tag_location 0 count 1 value_offset 31254 - ProjectedCSTypeGeoKey: MGI / Austria GK West
      key 3076 tiff_tag_location 0 count 1 value_offset 9001 - ProjLinearUnitsGeoKey: Linear_Meter
      key 4099 tiff_tag_location 0 count 1 value_offset 9001 - VerticalUnitsGeoKey: Linear_Meter
variable length header record 2 of 2:
  reserved             0
  user ID              'LASF_Projection'
  record ID            2112
  length after header  636
  description          'by LAStools of rapidlasso GmbH'
    WKT OGC COORDINATE SYSTEM:
    PROJCS["MGI / Austria GK West",GEOGCS["MGI",DATUM["Militar-Geographische Institut",SPHEROID["Bessel 1841",6377397.155,299.1528128,AUTHORITY["EPSG","7004"]],
AUTHORITY["EPSG","6312"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.01745329251994328,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4312"]],P
ROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",10.3333333333333],PARAMETER["scale_factor",1],PARAMETER["false_e
asting",0],PARAMETER["false_northing",-5000000],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","31254"]]
the header is followed by 2 user-defined bytes
LASzip compression (version 3.4r3 c3 50000): POINT14 3
reporting minimum and maximum for all LAS point record entries ...
  X            -5165835   -5141713
  Y             6396479    6428864
  Z               39141      39765
  intensity           5       6487
  return_number       1          2
  number_of_returns   1          2
  edge_of_flight_line 0          0
  scan_direction_flag 1          1
  classification      0          9
  scan_angle_rank   -28         21
  user_data          29        142
  point_source_ID   185        187
  gps_time 121970.520615 122644.581305
  extended_return_number          1      2
  extended_number_of_returns      1      2
  extended_classification         0      9
  extended_scan_angle         -4667   3500
  extended_scanner_channel        0      0
number of first returns:        155321
number of intermediate returns: 0
number of last returns:         155301
number of single returns:       155153
overview over extended number of returns of given pulse: 155153 316 0 0 0 0 0 0 0 0 0 0 0 0 0
histogram of classification of points:
            1095  never classified (0)
             107  unclassified (1)
          103916  ground (2)
             419  low vegetation (3)
           49932  water (9)

--
Download LAStools at
http://lastools.org
http://rapidlasso.com
Be social with LAStools at
http://facebook.com/LAStools
http://twitter.com/LAStools
http://linkedin.com/groups/LAStools-4408378
Manage your settings at
http://groups.google.com/group/lastools/subscribe
---
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/a85cdfc5-93f9-407c-a8c2-9045c92b272b%40googlegroups.com.

Drexel Peter

unread,
Dec 9, 2019, 11:33:40 AM12/9/19
to LAStools - efficient tools for LiDAR processing
Hello Martin,
works fine so far, thanks.
2D is - as expected - shifted by ~1m,
but what happens with the elevations? Are these Elevations above Bessel 1841 Elipsoid?
On Thu, Dec 5, 2019 at 3:27 PM Drexel Peter <peter...@vorarlberg.at> wrote:
Hi Martin,
I tried to use

las2las -i SynthetischePunkte_neuer_rhein.laz -o 31254-SynthetischePunkte_neuer_rhein.laz -olaz -epsg 25832 -target_epsg 31254

to convert my EPSG 25832-Pointcloud to EPSG 31254.
The results are off by ~500m...
Thanks
Peter

--
Download LAStools at
http://lastools.org
http://rapidlasso.com
Be social with LAStools at
http://facebook.com/LAStools
http://twitter.com/LAStools
http://linkedin.com/groups/LAStools-4408378
Manage your settings at
http://groups.google.com/group/lastools/subscribe
---
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 last...@googlegroups.com.

Martin Isenburg

unread,
Dec 10, 2019, 2:52:31 PM12/10/19
to LAStools - efficient command line tools for LIDAR processing
Hello,

Honestly I do not know that details as I have just recently learned about those seven parameter transforms myself. Can anyone else weigh in here? What happens with the elevation values when you do a DATUM change in an Earth-Centered Earth-Fixed representation using a Helmert Transform with all 7 parameters as described in my answer earlier?

Regards,

Martin

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/b48ab34f-9291-4622-ba5f-e952732521da%40googlegroups.com.

Kirk Waters - NOAA Federal

unread,
Dec 10, 2019, 3:51:37 PM12/10/19
to LAStools - efficient command line tools for LIDAR processing
In general, you're using the Helmert transform to go between different reference frames that are based on ellipsoids. You go through ECEF when doing that. The changes are 3D as a rule, so if you started in ellipsoid heights in reference frame A and transformed to reference frame B, you are now in ellipsoid heights for reference frame B (whatever its ellipsoid might be). Even if the ellipsoids are the same between the two frames (e.g. WGS84 and NAD83 both use GRS80), your heights will change (for example, see this blog post). If you started in orthometric heights and applied the Helmert transform, you probably have a mess vertically. You would want to reduce to ellipsoid heights first, Helmert transform, and then apply the appropriate geoid model for your new frame to get back to some orthometric datum. 

Kirk Waters, PhD, BJCP          | NOAA Office for Coastal Management
Applied Sciences Program      | 2234 South Hobson Ave
843-740-1227                          | Charleston, SC 29405    



Sverre Froyen

unread,
Dec 10, 2019, 9:19:53 PM12/10/19
to last...@googlegroups.com
Interesting…

The Helmert transform parameters are time dependent so you usually use a 14 parameter transform (the original 7 plus their time derivatives). In addition you need the epoch date.

Two questions:

1) Are the geoid parameters time dependent too (or is the geoid resolution low enough that time dependence can be ignored)?

2) Is there a standard reference frame used by GPS devices when reporting lat/lon and height?

Sverre Froyen
> To view this discussion on the web visit https://groups.google.com/d/msgid/lastools/CADm%3DQrQazeDGuZKifkiFPgvucL6_hPKrHA0RRNB%3DFtBhFq1d%2BA%40mail.gmail.com.

Kirk Waters - NOAA Federal

unread,
Dec 11, 2019, 7:57:20 AM12/11/19
to LAStools - efficient command line tools for LIDAR processing
Sverre,
You'll find both 7 parameter (no time dependence) and 14 parameter (time dependent) Helmert transform parameters. To my knowledge (I'm not a geodesist), geoid models do not have a time dependency, but they are often associated with a particular realization of a reference frame. I'd guess that we're generally not a point where we have enough confidence in the models that we would be able to assign velocities. Many geoid models are grid systems, so you'd probably need to come up with a velocity for each cell. I think your point about the resolution is valid too.

Regarding the GPS or GNSS systems, I'm doing some guessing and I hope someone will correct me if I get this wrong.  I'm pretty sure they are working in an ECEF reference frame at their most fundamental level. Basically doing something like a triangulation based on distance from each satellite and our understanding of their orbits. This is why you might wait a couple weeks until we have more precise information about what the orbits really were at the time of the collection before processing. The reference frame most aligned with GPS is probably ITRF14 (see http://itrf.ensg.ign.fr/general.php for more on ITRF). I expect most GPS/GNSS collection systems have the equations built in to give you coordinates in other systems.

I think one of the things to think about is why there might be time components. It's the same earth being modeled for every reference frame, so why would the relationship between frames change over time? I think at least part of this has to do with plate tectonics. A frame like NAD83 for the U.S. is tied to the North American plate and moves with it, so the coordinates for a piece of ground don't change (too much) over time. A frame like WGS84 is not tied to any plate, it is fixed in space. So, those two frames have a changing relationship as the North American plate moves. A different frame that is also fixed in space might have all zeros for the time dependence because the relationship of the axes between it and WGS84 aren't changing.

There are probably some geodesists lurking out there that can give better answers. I hope some will chime in as I'd like to understand some of this better too.


Kirk Waters, PhD, BJCP          | NOAA Office for Coastal Management
Applied Sciences Program      | 2234 South Hobson Ave
843-740-1227                          | Charleston, SC 29405    


Drexel Peter

unread,
Dec 16, 2019, 9:39:56 AM12/16/19
to LAStools - efficient tools for LiDAR processing
Hello,
do you have any plans to include NTV2-Grids as an alternative to helmert?

Also a great improvement would be to use the helmert only in 2D but leave the z as it is. So I could do the 2. step (applying the geoid) with lasheight.exe

Regards
Peter


Am Dienstag, 10. Dezember 2019 20:52:31 UTC+1 schrieb Martin Isenburg:
Hello,

Honestly I do not know that details as I have just recently learned about those seven parameter transforms myself. Can anyone else weigh in here? What happens with the elevation values when you do a DATUM change in an Earth-Centered Earth-Fixed representation using a Helmert Transform with all 7 parameters as described in my answer earlier?

Regards,

Martin

Martin Isenburg

unread,
Dec 17, 2019, 10:15:21 AM12/17/19
to LAStools - efficient command line tools for LIDAR processing
Hello,

as Kirk points out you need to use ellipsoidal heights before converting to ECEF, so if you have orthometric heights now you will need first to use lasheight to go from orthometric to  ellipsoidal height, then do the entire transform I described before, and at the very end you need to use lasheight again - with a different set of GEOID correction grids - to go from the new ellipsoid that you are on then back to orthometric.

Are the NTV2-Grids for your area available and described somewhere?

Regards,

Martin

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/568e7031-add4-45ca-89e6-c3acb9ae9ff9%40googlegroups.com.

michaelperdue99

unread,
Dec 17, 2019, 4:05:21 PM12/17/19
to LAStools - efficient tools for LiDAR processing
Hello,

I stumbled upon the NTv2 Users Guide a number of years ago. It describes the file format, but I no longer have a local copy and can no longer find the pdf on the web. It appears that NRCAN has closed the barn door again. You can get a copy from NRCAN if you sign an agreement not to redistribute it, but I that agreement bars it's use in open source projects (laslib). Unfortunately, NRCAN is less that friendly when it comes to sharing the NTv2 file format specifications (among other things) with the open source community.

GDAL supports reading NTv2 file format and you could probably build your own reader using the information in it's source code but I'm not sure you want to dive down that rabbit hole.

Cheers,

Mike


On Tuesday, December 17, 2019 at 8:15:21 AM UTC-7, Martin Isenburg wrote:
Hello,

as Kirk points out you need to use ellipsoidal heights before converting to ECEF, so if you have orthometric heights now you will need first to use lasheight to go from orthometric to  ellipsoidal height, then do the entire transform I described before, and at the very end you need to use lasheight again - with a different set of GEOID correction grids - to go from the new ellipsoid that you are on then back to orthometric.

Are the NTV2-Grids for your area available and described somewhere?

Regards,

Martin

Michael Stimson

unread,
Dec 17, 2019, 6:27:54 PM12/17/19
to last...@googlegroups.com

There’s nothing wrong with using GDAL to read an NTv2 transformation file as a raster, the API is fairly easy to use, however doing so would taint LasTools with open source so it’s a big decision for commercially licensed software. I know that the author of LasTools has GeoTIFF reading/writing ability so perhaps it would be better for users to use GDAL_Translate or similar to export the NTv2 GSB to GeoTIFF 3 band (xShift, yShift, zOffset), read about it https://gis.stackexchange.com/questions/272693/analyze-the-extent-of-a-ntv2-grid , which would then remove LasTools from the open source operation and protect the well deserved commercial status of the LasTools suite.

 

Michael Stimson

Senior GIS Analyst
RPS | Australia Asia Pacific
+61 7 3539 9694
michael...@rpsgroup.com.au

 

From: last...@googlegroups.com <last...@googlegroups.com> On Behalf Of michaelperdue99
Sent: Wednesday, 18 December 2019 2:58 AM
To: LAStools - efficient tools for LiDAR processing <last...@googlegroups.com>
Subject: Re: [LAStools] las2las coodinate transformation

 

CAUTION: This email originated from outside of RPS.

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/13baa430-5de5-4dbd-ada2-9f054443d537%40googlegroups.com.

This email and its attachments may contain confidential and/or privileged information and is for the sole use of the intended recipient(s). The contents of this email must not be disclosed to or used by or copied in any way by anyone other than the intended recipient(s). If you are not the intended recipient, any use, distribution or copying of the information contained in this email and its attachments is strictly prohibited. Confidentiality and/or privilege in the content of this email is not waived. If you have received this email in error, please email the sender by replying to this message and immediately delete and destroy any copies of this email and any attachments. Please note that neither RPS Consultants Pty Ltd, any subsidiary, related entity ('RPS') nor the sender accepts any responsibility for viruses and it is your responsibility to scan or otherwise check this email and any attachments. The views or opinions expressed are the author's own and may not reflect the views or opinions of RPS

Newcomb, Doug

unread,
Dec 18, 2019, 8:19:51 AM12/18/19
to last...@googlegroups.com
As open source licenses go, gdal's is very  permissive
https://gdal.org/license.html  

Many of the companies on this list, https://trac.osgeo.org/gdal/wiki/SoftwareUsingGdal , do not share your concerns about gdal's license "tainting" their software.




--
Doug Newcomb - Cartographer
USFWS
551F Pylon Dr
Raleigh, NC
---------------------------------------------------------------------------------------------------------

NOTE: This email correspondence and any attachments to and from this sender is subject to the Freedom of Information Act (FOIA) and may be disclosed to third parties.

Nicolas Cadieux

unread,
Dec 18, 2019, 10:06:06 AM12/18/19
to last...@googlegroups.com, Michael Stimson

Hi,

I'm not sure if Mike is Michael here but I am sending you the description of byn Grid Format.  Last year, I collaborated with NRCAN and a programmer to add the .byn geoid files to gdal to permit reading the file format.  NRCAN   You can also download a version at https://www.nrcan.gc.ca/sites/www.nrcan.gc.ca/files/earthsciences/pdf/gpshgrid_e.pdf.

.byn Files are covered by Canada's Open License now.

Nicolas

DescriptionGridFormats_e.pdf

Michael Perdue

unread,
Dec 18, 2019, 12:02:45 PM12/18/19
to LAStools - efficient command line tools for LIDAR processing
Hi Michael,

Only select tools of the LAStools package are closed source. The underlying library being used by LAStools (laslib and laszip) to read and write files is released under a copyleft LGPL license.

GDAL and PROJ are both released under the X/MIT license which, as Doug has already mentioned, is very permissive. As a result, there are no licensing issues preventing LAStools from linking into these libraries. Martin has expressed his reasoning for not wanting to use these external license in the past and I understand and respect his decision.

The NTv2 file format is not simply a flat binary file. Yes, GDAL can read and convert this file to other formats and for small nations or limited extents, there may be a single grid stored within the *.GSB file. But that is definitely not true for other locals. For instance, look at the gdalinfo output for Canada's NAD27 <-> NAD83 NTv2 conversion grid.
 gdalinfo.exe NTV2_0.GSB
Driver: NTv2/NTv2 Datum Grid Shift
Files: NTV2_0.GSB
Size is 529, 241
Coordinate System is:
GEOGCRS["WGS 84",
    DATUM["World Geodetic System 1984",
        ELLIPSOID["WGS 84",6378137,298.257223563,
            LENGTHUNIT["metre",1]]],
    PRIMEM["Greenwich",0,
        ANGLEUNIT["degree",0.0174532925199433]],
    CS[ellipsoidal,2],
        AXIS["geodetic latitude (Lat)",north,
            ORDER[1],
            ANGLEUNIT["degree",0.0174532925199433]],
        AXIS["geodetic longitude (Lon)",east,
            ORDER[2],
            ANGLEUNIT["degree",0.0174532925199433]],
    ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
Origin = (-88.041666666666671,60.041666666666664)
Pixel Size = (0.083333333333333,-0.083333333333333)
Metadata:
  CREATED=95-06-30
  DATUM_F=NAD27
  DATUM_T=NAD83
  GS_TYPE=SECONDS
  MAJOR_F=6378206.4
  MAJOR_T=6378137
  MINOR_F=6356583.8
  MINOR_T=6356752.314
  PARENT=NONE
  SUB_NAME=CAeast
  UPDATED=95-07-05
  VERSION=NTv2.0
Subdatasets:
  SUBDATASET_0_NAME=NTv2:0:NTV2_0.GSB
  SUBDATASET_0_DESC=CAeast
  SUBDATASET_1_NAME=NTv2:1:NTV2_0.GSB
  SUBDATASET_1_DESC=CAwest
  SUBDATASET_2_NAME=NTv2:2:NTV2_0.GSB
  SUBDATASET_2_DESC=CAnorth
  SUBDATASET_3_NAME=NTv2:3:NTV2_0.GSB
  SUBDATASET_3_DESC=CAarctic
  SUBDATASET_4_NAME=NTv2:4:NTV2_0.GSB
  SUBDATASET_4_DESC=ALbanff
  SUBDATASET_5_NAME=NTv2:5:NTV2_0.GSB
  SUBDATASET_5_DESC=ALbarhed
  SUBDATASET_6_NAME=NTv2:6:NTV2_0.GSB
  SUBDATASET_6_DESC=ALbonvil
  SUBDATASET_7_NAME=NTv2:7:NTV2_0.GSB
  SUBDATASET_7_DESC=ALbowisl
  SUBDATASET_8_NAME=NTv2:8:NTV2_0.GSB
  SUBDATASET_8_DESC=ALbrooks
  SUBDATASET_9_NAME=NTv2:9:NTV2_0.GSB
  SUBDATASET_9_DESC=ALcalgry
  SUBDATASET_10_NAME=NTv2:10:NTV2_0.GSB
  SUBDATASET_10_DESC=ALcamros
  SUBDATASET_11_NAME=NTv2:11:NTV2_0.GSB
  SUBDATASET_11_DESC=ALcanmor
  SUBDATASET_12_NAME=NTv2:12:NTV2_0.GSB
  SUBDATASET_12_DESC=ALcardst
  SUBDATASET_13_NAME=NTv2:13:NTV2_0.GSB
  SUBDATASET_13_DESC=ALcarsta
  SUBDATASET_14_NAME=NTv2:14:NTV2_0.GSB
  SUBDATASET_14_DESC=ALclarho
  SUBDATASET_15_NAME=NTv2:15:NTV2_0.GSB
  SUBDATASET_15_DESC=ALcoldlk
  SUBDATASET_16_NAME=NTv2:16:NTV2_0.GSB
  SUBDATASET_16_DESC=ALcrowps
  SUBDATASET_17_NAME=NTv2:17:NTV2_0.GSB
  SUBDATASET_17_DESC=ALdraytn
  SUBDATASET_18_NAME=NTv2:18:NTV2_0.GSB
  SUBDATASET_18_DESC=ALdrumhl
  SUBDATASET_19_NAME=NTv2:19:NTV2_0.GSB
  SUBDATASET_19_DESC=ALedmntn
  SUBDATASET_20_NAME=NTv2:20:NTV2_0.GSB
  SUBDATASET_20_DESC=ALedson
  SUBDATASET_21_NAME=NTv2:21:NTV2_0.GSB
  SUBDATASET_21_DESC=ALfairvw
  SUBDATASET_22_NAME=NTv2:22:NTV2_0.GSB
  SUBDATASET_22_DESC=ALftmacl
  SUBDATASET_23_NAME=NTv2:23:NTV2_0.GSB
  SUBDATASET_23_DESC=ALftmcmr
  SUBDATASET_24_NAME=NTv2:24:NTV2_0.GSB
  SUBDATASET_24_DESC=ALgrcach
  SUBDATASET_25_NAME=NTv2:25:NTV2_0.GSB
  SUBDATASET_25_DESC=ALgrimsh
  SUBDATASET_26_NAME=NTv2:26:NTV2_0.GSB
  SUBDATASET_26_DESC=ALgrprar
  SUBDATASET_27_NAME=NTv2:27:NTV2_0.GSB
  SUBDATASET_27_DESC=ALhanna
  SUBDATASET_28_NAME=NTv2:28:NTV2_0.GSB
  SUBDATASET_28_DESC=ALhilevl
  SUBDATASET_29_NAME=NTv2:29:NTV2_0.GSB
  SUBDATASET_29_DESC=ALhinton
  SUBDATASET_30_NAME=NTv2:30:NTV2_0.GSB
  SUBDATASET_30_DESC=ALhiprai
  SUBDATASET_31_NAME=NTv2:31:NTV2_0.GSB
  SUBDATASET_31_DESC=ALinnsfl
  SUBDATASET_32_NAME=NTv2:32:NTV2_0.GSB
  SUBDATASET_32_DESC=ALjasper
  SUBDATASET_33_NAME=NTv2:33:NTV2_0.GSB
  SUBDATASET_33_DESC=ALlacbic
  SUBDATASET_34_NAME=NTv2:34:NTV2_0.GSB
  SUBDATASET_34_DESC=ALlacomb
  SUBDATASET_35_NAME=NTv2:35:NTV2_0.GSB
  SUBDATASET_35_DESC=ALletbrg
  SUBDATASET_36_NAME=NTv2:36:NTV2_0.GSB
  SUBDATASET_36_DESC=ALlkloui
  SUBDATASET_37_NAME=NTv2:37:NTV2_0.GSB
  SUBDATASET_37_DESC=ALlydmin
  SUBDATASET_38_NAME=NTv2:38:NTV2_0.GSB
  SUBDATASET_38_DESC=ALmedhat
  SUBDATASET_39_NAME=NTv2:39:NTV2_0.GSB
  SUBDATASET_39_DESC=ALolds
  SUBDATASET_40_NAME=NTv2:40:NTV2_0.GSB
  SUBDATASET_40_DESC=ALoyen
  SUBDATASET_41_NAME=NTv2:41:NTV2_0.GSB
  SUBDATASET_41_DESC=ALpeacer
  SUBDATASET_42_NAME=NTv2:42:NTV2_0.GSB
  SUBDATASET_42_DESC=ALpinchr
  SUBDATASET_43_NAME=NTv2:43:NTV2_0.GSB
  SUBDATASET_43_DESC=ALponoka
  SUBDATASET_44_NAME=NTv2:44:NTV2_0.GSB
  SUBDATASET_44_DESC=ALraymnd
  SUBDATASET_45_NAME=NTv2:45:NTV2_0.GSB
  SUBDATASET_45_DESC=ALredeer
  SUBDATASET_46_NAME=NTv2:46:NTV2_0.GSB
  SUBDATASET_46_DESC=ALrockmt
  SUBDATASET_47_NAME=NTv2:47:NTV2_0.GSB
  SUBDATASET_47_DESC=ALslavlk
  SUBDATASET_48_NAME=NTv2:48:NTV2_0.GSB
  SUBDATASET_48_DESC=ALstetlr
  SUBDATASET_49_NAME=NTv2:49:NTV2_0.GSB
  SUBDATASET_49_DESC=ALstpaul
  SUBDATASET_50_NAME=NTv2:50:NTV2_0.GSB
  SUBDATASET_50_DESC=ALstramr
  SUBDATASET_51_NAME=NTv2:51:NTV2_0.GSB
  SUBDATASET_51_DESC=ALswanhi
  SUBDATASET_52_NAME=NTv2:52:NTV2_0.GSB
  SUBDATASET_52_DESC=ALtaber
  SUBDATASET_53_NAME=NTv2:53:NTV2_0.GSB
  SUBDATASET_53_DESC=ALtrehil
  SUBDATASET_54_NAME=NTv2:54:NTV2_0.GSB
  SUBDATASET_54_DESC=ALvegvil
  SUBDATASET_55_NAME=NTv2:55:NTV2_0.GSB
  SUBDATASET_55_DESC=ALvermil
  SUBDATASET_56_NAME=NTv2:56:NTV2_0.GSB
  SUBDATASET_56_DESC=ALwanwgt
  SUBDATASET_57_NAME=NTv2:57:NTV2_0.GSB
  SUBDATASET_57_DESC=ALweslok
  SUBDATASET_58_NAME=NTv2:58:NTV2_0.GSB
  SUBDATASET_58_DESC=ALwetask
  SUBDATASET_59_NAME=NTv2:59:NTV2_0.GSB
  SUBDATASET_59_DESC=ALwhitec
  SUBDATASET_60_NAME=NTv2:60:NTV2_0.GSB
  SUBDATASET_60_DESC=BCcambel
  SUBDATASET_61_NAME=NTv2:61:NTV2_0.GSB
  SUBDATASET_61_DESC=BCcranbk
  SUBDATASET_62_NAME=NTv2:62:NTV2_0.GSB
  SUBDATASET_62_DESC=BCdawson
  SUBDATASET_63_NAME=NTv2:63:NTV2_0.GSB
  SUBDATASET_63_DESC=BCelkfrd
  SUBDATASET_64_NAME=NTv2:64:NTV2_0.GSB
  SUBDATASET_64_DESC=BCfield
  SUBDATASET_65_NAME=NTv2:65:NTV2_0.GSB
  SUBDATASET_65_DESC=BCgranil
  SUBDATASET_66_NAME=NTv2:66:NTV2_0.GSB
  SUBDATASET_66_DESC=BCkamlop
  SUBDATASET_67_NAME=NTv2:67:NTV2_0.GSB
  SUBDATASET_67_DESC=BCkelwna
  SUBDATASET_68_NAME=NTv2:68:NTV2_0.GSB
  SUBDATASET_68_DESC=BClogan
  SUBDATASET_69_NAME=NTv2:69:NTV2_0.GSB
  SUBDATASET_69_DESC=BCmacknz
  SUBDATASET_70_NAME=NTv2:70:NTV2_0.GSB
  SUBDATASET_70_DESC=BCnanimo
  SUBDATASET_71_NAME=NTv2:71:NTV2_0.GSB
  SUBDATASET_71_DESC=BCnelson
  SUBDATASET_72_NAME=NTv2:72:NTV2_0.GSB
  SUBDATASET_72_DESC=BCparkvl
  SUBDATASET_73_NAME=NTv2:73:NTV2_0.GSB
  SUBDATASET_73_DESC=BCpentic
  SUBDATASET_74_NAME=NTv2:74:NTV2_0.GSB
  SUBDATASET_74_DESC=BCportal
  SUBDATASET_75_NAME=NTv2:75:NTV2_0.GSB
  SUBDATASET_75_DESC=BCpowell
  SUBDATASET_76_NAME=NTv2:76:NTV2_0.GSB
  SUBDATASET_76_DESC=BCprigeo
  SUBDATASET_77_NAME=NTv2:77:NTV2_0.GSB
  SUBDATASET_77_DESC=BCroslnd
  SUBDATASET_78_NAME=NTv2:78:NTV2_0.GSB
  SUBDATASET_78_DESC=BCtrail
  SUBDATASET_79_NAME=NTv2:79:NTV2_0.GSB
  SUBDATASET_79_DESC=BCtumblr
  SUBDATASET_80_NAME=NTv2:80:NTV2_0.GSB
  SUBDATASET_80_DESC=BCvancvr
  SUBDATASET_81_NAME=NTv2:81:NTV2_0.GSB
  SUBDATASET_81_DESC=BCvernon
  SUBDATASET_82_NAME=NTv2:82:NTV2_0.GSB
  SUBDATASET_82_DESC=BCvictor
  SUBDATASET_83_NAME=NTv2:83:NTV2_0.GSB
  SUBDATASET_83_DESC=NFstjohn
  SUBDATASET_84_NAME=NTv2:84:NTV2_0.GSB
  SUBDATASET_84_DESC=NWclyder
  SUBDATASET_85_NAME=NTv2:85:NTV2_0.GSB
  SUBDATASET_85_DESC=NWftgood
  SUBDATASET_86_NAME=NTv2:86:NTV2_0.GSB
  SUBDATASET_86_DESC=NWhayriv
  SUBDATASET_87_NAME=NTv2:87:NTV2_0.GSB
  SUBDATASET_87_DESC=NWinuvik
  SUBDATASET_88_NAME=NTv2:88:NTV2_0.GSB
  SUBDATASET_88_DESC=NWiqulit
  SUBDATASET_89_NAME=NTv2:89:NTV2_0.GSB
  SUBDATASET_89_DESC=NWpondin
  SUBDATASET_90_NAME=NTv2:90:NTV2_0.GSB
  SUBDATASET_90_DESC=NWrankin
  SUBDATASET_91_NAME=NTv2:91:NTV2_0.GSB
  SUBDATASET_91_DESC=NWyellow
  SUBDATASET_92_NAME=NTv2:92:NTV2_0.GSB
  SUBDATASET_92_DESC=ONkinstn
  SUBDATASET_93_NAME=NTv2:93:NTV2_0.GSB
  SUBDATASET_93_DESC=ONottawa
  SUBDATASET_94_NAME=NTv2:94:NTV2_0.GSB
  SUBDATASET_94_DESC=ONsarnia
  SUBDATASET_95_NAME=NTv2:95:NTV2_0.GSB
  SUBDATASET_95_DESC=ONsault
  SUBDATASET_96_NAME=NTv2:96:NTV2_0.GSB
  SUBDATASET_96_DESC=ONthundr
  SUBDATASET_97_NAME=NTv2:97:NTV2_0.GSB
  SUBDATASET_97_DESC=ONtimins
  SUBDATASET_98_NAME=NTv2:98:NTV2_0.GSB
  SUBDATASET_98_DESC=ONtronto
  SUBDATASET_99_NAME=NTv2:99:NTV2_0.GSB
  SUBDATASET_99_DESC=ONwinsor
  SUBDATASET_100_NAME=NTv2:100:NTV2_0.GSB
  SUBDATASET_100_DESC=SAestvan
  SUBDATASET_101_NAME=NTv2:101:NTV2_0.GSB
  SUBDATASET_101_DESC=SAmelfrt
  SUBDATASET_102_NAME=NTv2:102:NTV2_0.GSB
  SUBDATASET_102_DESC=SAmelvil
  SUBDATASET_103_NAME=NTv2:103:NTV2_0.GSB
  SUBDATASET_103_DESC=SAmosjaw
  SUBDATASET_104_NAME=NTv2:104:NTV2_0.GSB
  SUBDATASET_104_DESC=SAnbatle
  SUBDATASET_105_NAME=NTv2:105:NTV2_0.GSB
  SUBDATASET_105_DESC=SApralbt
  SUBDATASET_106_NAME=NTv2:106:NTV2_0.GSB
  SUBDATASET_106_DESC=SAregina
  SUBDATASET_107_NAME=NTv2:107:NTV2_0.GSB
  SUBDATASET_107_DESC=SAsatoon
  SUBDATASET_108_NAME=NTv2:108:NTV2_0.GSB
  SUBDATASET_108_DESC=SAswiftc
  SUBDATASET_109_NAME=NTv2:109:NTV2_0.GSB
  SUBDATASET_109_DESC=SAweybrn
  SUBDATASET_110_NAME=NTv2:110:NTV2_0.GSB
  SUBDATASET_110_DESC=SAyorktn
  SUBDATASET_111_NAME=NTv2:111:NTV2_0.GSB
  SUBDATASET_111_DESC=YUdawson
  SUBDATASET_112_NAME=NTv2:112:NTV2_0.GSB
  SUBDATASET_112_DESC=YUrossri
  SUBDATASET_113_NAME=NTv2:113:NTV2_0.GSB
  SUBDATASET_113_DESC=YUwhiteh
Corner Coordinates:
Upper Left  ( -88.0416667,  60.0416667) ( 88d 2'30.00"W, 60d 2'30.00"N)
Lower Left  ( -88.0416667,  39.9583333) ( 88d 2'30.00"W, 39d57'30.00"N)
Upper Right ( -43.9583333,  60.0416667) ( 43d57'30.00"W, 60d 2'30.00"N)
Lower Right ( -43.9583333,  39.9583333) ( 43d57'30.00"W, 39d57'30.00"N)
Center      ( -66.0000000,  50.0000000) ( 66d 0' 0.00"W, 50d 0' 0.00"N)
Band 1 Block=529x1 Type=Float32, ColorInterp=Undefined
  Description = Latitude Offset (arc seconds)
Band 2 Block=529x1 Type=Float32, ColorInterp=Undefined
  Description = Longitude Offset (arc seconds)
Band 3 Block=529x1 Type=Float32, ColorInterp=Undefined
  Description = Latitude Error
Band 4 Block=529x1 Type=Float32, ColorInterp=Undefined
  Description = Longitude Error


It contains 113 sub datasets that differ in their spatial extents, resolution and accuracy. These datasets are arranged in a hierarchical structure where the values in one dataset will supersede the values in a subordinate grid and it is the softwares responsibility to determine which is the best grid to derive values from. If a user simply does;
gdal_translate NTv2.gsb NTv2.tif
gdal_translate will simply pull the first grid, which will probably be wrong. What you are describing could still be done but getting the correct answers would take some data investigation/messaging that is beyond the abilities of most users.

Cheers,

Mike

Michael Perdue

unread,
Dec 18, 2019, 12:11:17 PM12/18/19
to LAStools - efficient command line tools for LIDAR processing
Hi Nicolas,

This document describes the binary formats typically used to store geoids. Maybe you can use your contacts to enquire about the equivalent document for the NTv2 file format. It is a very popular format for storing grid shift information. I think it is only the US which uses a different file format.

Cheers,

Mike

Nicolas Cadieux

unread,
Dec 18, 2019, 1:13:27 PM12/18/19
to last...@googlegroups.com, Michael Perdue
Reply all
Reply to author
Forward
0 new messages