[pdal] Regarding Reprojection Between Vertical Datums
Paul Harwood
runette at gmail.com
Tue Jul 12 23:34:22 PDT 2022
Pdalinfo reports the bbox in the native CRS AND 4326 - so that 4326 tag is
created by pdal and is not in the original. Look at the "native" tag under
bbox.
On Wed, 13 Jul 2022, 06:33 Nicholas Stanley, <
nicholas.stanley at luminartech.com> wrote:
>
>
>
>
>
>
> * For the fun of it, here's a link to the Ground Truth and Data to be
> reprojected in alignment to it. If Someone thinks they can crack this,
> please give it a shot!
> https://drive.google.com/drive/folders/1KdUc1xnOtNl0RkWfr1-C_gh7jzU7q9U1?usp=sharing
> <https://drive.google.com/drive/folders/1KdUc1xnOtNl0RkWfr1-C_gh7jzU7q9U1?usp=sharing>
> *
>
> *The Source Materials (to be corrected) are described here:*
> https://www.fisheries.noaa.gov/inport/item/59131
>
>
> Coordinate Reference System
> CRS Type: Geographic 3D
> *EPSG Code: EPSG:6319*
> EPSG Name: NAD83(2011)
> Datum Name: NAD83 (National Spatial Reference System 2011)
> Ellipsoid: GRS 1980
> Semimajor Axis: 6378137.0
> Flattening Ratio: 298.257222101
> *"ldrgeoid" : "GEOID 12B",*
>
> *That said, when i run pdal info on the file, i see the following:*
>
> "file_size": 10024109,
> "filename": "USGS_LPC_CA_NoCAL_Wildfires_B5b_2018_w2274n1959.laz",
> "now": "2022-07-12T21:37:52-0700",
> "pdal_version": "2.2.0 (git-version: Release)",
> "reader": "readers.las",
> "stats":
> {
> "bbox":
> {
> * "EPSG:4326":*
> {
> "bbox":
> {
> "maxx": -122.4022284,
> "maxy": 37.81572701,
> "maxz": 99.83,
> "minx": -122.4164013,
> "miny": 37.80469945,
> "minz": -70.28
> },
> "boundary": { "type": "Polygon", "coordinates": [ [ [
> -122.413254958928221, 37.804699453399095, -70.28 ], [ -122.416401346373235,
> 37.813279215545926, -70.28 ], [ -122.405373596515062, 37.815727007213006,
> 99.83 ], [ -122.402228384472423, 37.807146965417473, 99.83 ], [
> -122.413254958928221, 37.804699453399095, -70.28 ] ] ] }
> },
> "native":
> {
> "bbox":
> {
> "maxx": -2273000,
> "maxy": 1959999.99,
> "maxz": 99.83,
> "minx": -2273999.99,
> "miny": 1959000,
> "minz": -70.28
> },
> "boundary": { "type": "Polygon", "coordinates": [ [ [
> -2273999.990000000223517, 1959000.0, -70.28 ], [ -2273999.990000000223517,
> 1959999.99, -70.28 ], [ -2273000.0, 1959999.99, 99.83 ], [ -2273000.0,
> 1959000.0, 99.83 ], [ -2273999.990000000223517, 1959000.0, -70.28 ] ] ] }
> }
> },
>
>
>
>
> *EPSG 4326 is a WGS84 code, but this is merely defining the bounding box
> container it seems. *
>
>
>
>
>
>
> *Lasinfo shows: description '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"]]]
> variable length header record 4 of 4: *
> ---------
>
> *//What I am also suggesting is that you obtain a point manually (some
> lat/lon and the height) from your dataset, at one of the ground truth
> dataset points, and then take that and do the various conversions
> manually using NGS and other tools, such as//*
>
> Having done this, I see coordinates that more closely align in Elevation
> to our ground Truth (although they still appear to be off by ~2-3 meters,
> not 30)
>
> Output from GEOID18
>
> latitude longitude N error (95% confidence interval)
> Station Name ddd mm ss.sssss ddd mm ss.sssss meters meters
> USER LOCATION 37 48 25.00920 237 35 38.63040 -32.486 0.031
>
> Our "Ground Truth" in this area shows an Elevation at -29.17m.
>
>
>
> On 7/12/22 16:35, Greg Troxel wrote:
>
>
> Nicholas Stanley <nicholas.stanley at luminartech.com> <nicholas.stanley at luminartech.com> writes:
>
>
> Our "ground truth" is using WGS84 Ellipsoidal Height. The survey was
> created in ~2017, but unfortunately I don't have much additional info
> about the Epoch.
>
>
> One would assume that the coordinate epoch is the date of the
> observations.
>
> I am very curious about this survey and how it was done. I do not know
> of any mechanisms to obtain precise coordinates in WGS84 (assuming you
> aren't working within NGA). If it's uncorrected data I don't see how
> it would be that useful. If it's RTK or static carrier phase, I don't
> see how it could be in WGS84.
>
>
> Technically, yes, I am in California - is it fair to say this
> complicates things somewhat? The data for conversion is a USGS survey
> of the San Francisco Bay Area.
>
>
> As you know California is unstable :-) Most stations in the US are
> located on the North American plate and their NAD83 velocities are low,
> on the order of a few mm/year. Parts of CA are on the Pacific plate and
> have much higher velocities. And, NAD83 has a velocity relative to
> WGS84/ITRF. So if you are trying to get things right at the 10 cm
> level, you'll have to be careful about this. NGS has lots of
> information about seismic/plate-motion coordinate changes. I don't
> really understand the details as I'm in Massachusetts, far from a plate
> boundary. But, if you are going from NAD83(2011) epoch 2010.0 to ITRF,
> I think the Pacific plate issues won't be that much trouble.
>
>
> Following the logic of your process:
>
> *Step 1: Convert the **/NAD83(2011) / Conus Albers + NAVD88 height
> /**to NAD83 Ellipsoidal using GEOID18 *
>
>
> E.G.: pdal translate
> USGS_LPC_CA_NoCAL_Wildfires_B5b_2018_w2274n1959.laz Step1_NAD83.laz
> reprojection --filters.reprojection.out_srs="+init=EPSG:xxxx
> +geoidgrids=~/Downloads/vdatum/core/geoid18/v2018prvi.gtx"
> --writers.las.a_srs=EPSG:xxxx
>
>
> v2018prvi looks like a grid for the Puerto Rico region. Is there a
> v2018conus.gtx? Why are you using the prvi one?
>
> If you run pdalinfo on the .laz file, does it output a CRS, and does
> that CRS match the documented CRS for the LIDAR data?
>
> What I am also suggesting is that you obtain a point manually (some
> lat/lon and the height) from your dataset, at one of the ground truth
> dataset points, and then take that and do the various conversions
> manually using NGS and other tools, such as
> https://www.ngs.noaa.gov/NCAT/
> https://www.ngs.noaa.gov/GEOID/GEOID18/
> https://www.ngs.noaa.gov/GEOID/GEOID18/computation.html
> https://geographiclib.sourceforge.io/cgi-bin/GeoidEval
>
> (Not sure which EPSG codes to use here - any guess?)
>
>
> I don't think you can assume that there is an EPSG code for a compound
> CRS that's e.g. a combination of a UTM zone and a height.
>
>
> *Step 2. Convert resulting Step 1. file into destination projection
> (WGS84 UTM 10 Ellipsoidal https://epsg.io/32610, or alternatively,https://epsg.io/7789 ITRF 2014)
>
> E.G.: pdal translate
> USGS_LPC_CA_NoCAL_Wildfires_B5b_2018_w2274n1959.laz Step1_NAD83.laz
> reprojection --filters.reprojection.out_srs="+init=EPSG:xxxx
> +geoidgrids=~/Downloads/vdatum/core/geoid18/v2018prvi.gtx"
> --writers.las.a_srs=EPSG:xxxx
>
>
> I think you need to define a compound CRS that includes height and avoid
> handling 3D data with a 2D CRS.
>
> I am not really familiar with pdal translate options; your arguments
> look like old-style proj rather than modern proj which lets you query
> for transformations and see them.
>
> _______________________________________________
> 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/20220713/523ebd60/attachment.htm>
More information about the pdal
mailing list