[pdal] Regarding Reprojection Between Vertical Datums

Tue Jul 12 14:43:08 PDT 2022

The following have worked for me in the past with PDAL 2.3.0, GDAL 3.4, Proj 8.2.0 (although I have edited it since I was going from EPSG:6344+5703 to EPSG:4326 at the time):

--filter filters.reprojection --filters.reprojection.out_srs="EPSG:32610"

or

--filter filters.reprojection --filters.reprojection.in_srs="+init=epsg:6350 +geoidgrids=us_noaa_g2018u0.tif" --filters.reprojection.out_srs="+init=epsg:32610"

or

--filter filters.reprojection --filters.reprojection.in_srs="+init=epsg:6350 +geoidgrids=us_noaa_g2012bu0.tif" --filters.reprojection.out_srs="+init=epsg:32610"

All three shifted the data about 28m "lower" in Z (in Minnesota).

On 7/12/22 16:12, Nicholas Stanley wrote:
>
> Thanks,
>
> Yes, that was an oversight in respect to the +4326 - thanks.
>
> To reiterate:
>
> I am trying to get the /NAD83(2011) / Conus Albers + NAVD88 height/ data *out* of its native elevation model, converted *into* a wgs84 EPSG 32610 ellipsoidal projection.
>
> When I run the conversion command /  //pdal translate filename.laz filename_utm.laz reprojection --filters.reprojection.out_srs="EPSG:32610"//
>
> /The resulting file lines up in Easting/Northing precisely, but seems to have retained the old elevation data (which is ~30M+ above our UTM10 WGS84 "Ground Truth")
>
> On 7/12/22 14:03, Jim Klassen wrote:
>>
>>
>> On 7/12/22 15:41, Nicholas Stanley wrote:
>>>
>>> Greetings,
>>>
>>> My goal is to reproject/align a third party GRS80 point cloud data source to an existing UTM10 WGS84 (EPSG 32610) georeferenced point cloud file/"ground truth".
>>>
>>> Running Lasinfo on the third party data source (to be corrected/reprojected/aligned to Ground Truth) shows the following:
>>>
>>>
>>> /// 'OGC Transformation Record'//
>>> //    WKT OGC COORDINATE SYSTEM://
>>> //    COMPD_CS["NAD83(2011) / Conus Albers + NAVD88 height",PROJCS["NAD83(2011) / Conus Albers",GEOGCS["NAD83(2011)",DATUM["NAD83 (National Spatial Reference System 2011)",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","1116"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","6318"]],PROJECTION["Albers_Conic_Equal_Area"],PARAMETER["standard_parallel_1",29.5],PARAMETER["standard_parallel_2",45.5],PARAMETER["latitude_of_center",23],PARAMETER["longitude_of_center",-96],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["meter",1,AUTHORITY["EPSG","9001"]],AXIS["X",EAST],AXIS["Y",NORTH],AUTHORITY["EPSG","6350"]],VERT_CS["NAVD88 height",VERT_DATUM["North American Vertical Datum 1988",2005,AUTHORITY["EPSG","5103"]],UNIT["meter",1,AUTHORITY["EPSG","9001"]],AXIS["Up",UP],AUTHORITY["EPSG","5703"]]] ///
>>>
>>>
>>> I am able to successfully reproject the file in PDAL so that the relative 2D alignment in Easting/Northing is corrected - using a command like this:
>>>
>>> ///pdal translate filename.laz filename_utm.laz reprojection --filters.reprojection.out_srs="EPSG:32610"///
>>>
>>>
>>> *Reviewing the output, this shows proper 2D Grid alignment between ground truth and third party data, but a delta of ~30M in elevation between the ground planes. **This needs to be rectified.
>>> *
>>>
>>> *I've recently found this site describing a process to reproject a .las file in respect to various vertical datums.*
>>>
>>>
>>> http://www.gadom.ski/2017/05/31/vdatum-with-pdal.html
>>>
>>> ///In order to use//PDAL <https://www.pdal.io/>//to convert vertical coordinates to a new vertical datum, we’ll need to fetch the appropriate geoid. PDAL uses//proj.4 <http://proj4.org/>//under the hood, and proj.4 accepts geoid’s in NOAA VDatum’s//|gtx|//file format. We can fetch the 12B geoid from the//VDatum download website <https://vdatum.noaa.gov/download.php>//. Inside the 12B zipfile there is a file called//|g2012b_conus.gtx|//; this is the geoid for the Continental United States (CONUS). We’ll use this file to convert our data.///
>>>
>>>
>>> I've downloaded the VDatum collection as stated - That said, I've been unable to get the process to work properly.
>>>
>>> ///pdal translate filename.laz filename_utm.laz reprojection --filters.reprojection.out_srs="+init=EPSG:32610 +geoidgrids=~/Downloads/vdatum/g2012b_conus.gtx" --writers.las.a_srs=EPSG:32610+4326//
>>> //(pdal translate Error) GDAL failure (1) SetCompoundCS() fails, vertical component is not VERT_CS.//
>>> //(pdal translate writers.las Error) GDAL failure (1) PROJ: proj_crs_get_coordinate_system: Object is not a SingleCRS///
>>>
>>>
>>> Does anyone have any Tips/Pointers as to how to resolve/debug this issue?
>>>
>>> Thanks.
>>>
>>>
>>> _______________________________________________
>>> pdal mailing list
>>> pdal at lists.osgeo.org
>>> https://lists.osgeo.org/mailman/listinfo/pdal
>>
>>
>> --writers.las.a_srs=EPSG:32610+4326
>>
>> There may be more going on here, but "+4326" is likely not correct.  You need a vert_cs, something like 5703, there.  4326 is a horizontal CS.  I'm not sure if there is a way to explicitly state you want ellipsoid heights.  In my experience that's the default if a separate vertical component isn't given./
>>
>> /Also if you have recent PDAL built against recent PROJ (6+) you probably want to pull the TIFF files from cdn.proj.org (instead of gtx) for the vertical datum conversion./
>>
>> /
>>
>> _______________________________________________
>> pdal mailing list
>> pdal at lists.osgeo.org
>> https://lists.osgeo.org/mailman/listinfo/pdal
>
> _______________________________________________
> pdal mailing list
> pdal at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/pdal
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/pdal/attachments/20220712/6d499ead/attachment-0001.htm>