<div dir="ltr">Nicholas,<div>Based on your link to the metadata, it looks like you have at least found a reference to the data on the NOAA system. That link describes the data as it exists in the backend of the <a href="https://coast.noaa.gov/dataviewer/#/lidar/search/where:ID=9067">Data Access Viewer</a> (the link will zoom to that data set). But, if your original data is in Albers, it sounds like you picked up the data from USGS. The quick cheat is to use the Data Access Viewer and request the data in WGS84/ITRF ellipsoid heights for the area you care about.</div><div><br></div><div>I think your primary issue is that the geoid model isn't getting removed. As mentioned, you have to specify the grid to use in PDAL. The CONUS area is split into 8 grids as distributed, though it's possible that someone has combined them. You'll need to make sure you specify the right grid. Looking at the metadata you linked, I see that we removed GEOID12B, so you need to make sure that's the grid set you use (not GEOID18). The gird name will likely look like g2012buX.gtx where X is a number from 1 to 8. The u is what is indicating that it is CONUS. PDAL may allow you to specify multiple grids at a time so you don't have to figure out which of the 8 it is. I suspect it is either grid 1 or 2 for California. I noticed that VDatum also has a g2012bu0.gtx which appears to be all 8 grids combined, so that would be an easy pick.</div><div><br></div><div>Once your sure you've got the geoid taken care of, I think there are still some other issues. As was mentioned by others, you need to specify that you're going to transform to ITRF2008 or ITRF2014 instead of simply WGS84. That's because PROJ is using a ballpark offset to convert to WGS84 (you can see that with projinfo -s EPSG:6318 -t EPSG:4326) with an accuracy around 2 meters. I suspect it's using the NAD83(86) to WGS84(Transit), but doesn't have a way to convert NAD83(2011) to NAD83(86) because the only published way is a multi-grid hop. Something like NAD83(2011)->NAD83(2007)->NAD83(HARN)->NAD83(86), but I might have missed a hop in there. However, there are published Helmert transforms from NAD83(2011) to ITRF(2008) and ITRF(2014). The accuracy of those transforms will be much better.</div><div><br></div><div>One last issue. It can be hard to get the georeferencing right for ellipsoid heights. There is no vertical code in EPSG for it. There isn't, AFAIK, either a convenient compound CRS EPSG code or a combination of horizontal and vertical codes to give for what you want as your output. You can get to geographic coordinates and ellipsoid heights using the 3D codes (e.g. for ITRF2014 it's EPSG:7912), but not projected coordinates. You may need to specify your transform in pdal to go to 7912 and then transform to UTM while avoiding a vertical transformation. The file georeferencing on your output may be missing the vertical, though I haven't tried it. </div><div><br></div><div><div><div dir="ltr" class="gmail_signature" data-smartmail="gmail_signature"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div dir="ltr"><div><font face="arial, helvetica, sans-serif">Kirk Waters, PhD </font></div><div><font face="arial, helvetica, sans-serif">NOAA Office for Coastal Management<br></font></div><div><font face="arial, helvetica, sans-serif">Applied Sciences Program </font></div><div>Phone: 843-284-6962 (New as of 2022) <font face="arial, helvetica, sans-serif"> </font><div><font face="arial, helvetica, sans-serif"><a href="http://coast.noaa.gov/digitalcoast" target="_blank">coast.noaa.gov/digitalcoast</a></font><br></div></div><div><br></div><div>Thunder is good, thunder is impressive; but it is lightning that does the work - Mark Twain</div><div><br><br></div></div></div></div></div></div></div></div></div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Wed, Jul 13, 2022 at 1:33 AM Nicholas Stanley <<a href="mailto:nicholas.stanley@luminartech.com">nicholas.stanley@luminartech.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div>
<p><b><br>
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!<br>
<br>
<a href="https://drive.google.com/drive/folders/1KdUc1xnOtNl0RkWfr1-C_gh7jzU7q9U1?usp=sharing" target="_blank">https://drive.google.com/drive/folders/1KdUc1xnOtNl0RkWfr1-C_gh7jzU7q9U1?usp=sharing</a><br>
<br>
<br>
</b></p>
<p><b>The Source Materials (to be corrected) are described here:</b>
<a href="https://www.fisheries.noaa.gov/inport/item/59131" target="_blank">https://www.fisheries.noaa.gov/inport/item/59131<br>
</a><br>
<br>
Coordinate Reference System<br>
CRS Type: Geographic 3D<br>
<b>EPSG Code: EPSG:6319</b><br>
EPSG Name: NAD83(2011)<br>
Datum Name: NAD83 (National Spatial Reference System 2011)<br>
Ellipsoid: GRS 1980<br>
Semimajor Axis: 6378137.0<br>
Flattening Ratio: 298.257222101<br>
<b>"ldrgeoid" : "GEOID 12B",</b><br>
<br>
<b>That said, when i run pdal info on the file, i see the
following:</b><br>
<br>
"file_size": 10024109,<br>
"filename":
"USGS_LPC_CA_NoCAL_Wildfires_B5b_2018_w2274n1959.laz",<br>
"now": "2022-07-12T21:37:52-0700",<br>
"pdal_version": "2.2.0 (git-version: Release)",<br>
"reader": "readers.las",<br>
"stats":<br>
{<br>
"bbox":<br>
{<br>
<b> "EPSG:4326":</b><br>
{<br>
"bbox":<br>
{<br>
"maxx": -122.4022284,<br>
"maxy": 37.81572701,<br>
"maxz": 99.83,<br>
"minx": -122.4164013,<br>
"miny": 37.80469945,<br>
"minz": -70.28<br>
},<br>
"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 ] ] ] }<br>
},<br>
"native":<br>
{<br>
"bbox":<br>
{<br>
"maxx": -2273000,<br>
"maxy": 1959999.99,<br>
"maxz": 99.83,<br>
"minx": -2273999.99,<br>
"miny": 1959000,<br>
"minz": -70.28<br>
},<br>
"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 ] ] ] }<br>
}<br>
},<br>
<br>
<br>
<b>EPSG 4326 is a WGS84 code, but this is merely defining the
bounding box container it seems.<br>
<br>
</b><br>
<b>Lasinfo shows:<br>
<br>
description 'OGC Transformation Record'<br>
WKT OGC COORDINATE SYSTEM:<br>
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"]]]<br>
variable length header record 4 of 4:<br>
</b><br>
---------<br>
</p>
<pre><i>//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//</i></pre>
<i>
</i>
<p>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)<br>
<br>
</p>
<h2>Output from GEOID18</h2>
<pre style="color:rgb(0,0,0);font-style:normal;font-variant-ligatures:normal;font-variant-caps:normal;font-weight:400;letter-spacing:normal;text-align:start;text-indent:0px;text-transform:none;word-spacing:0px;text-decoration-style:initial;text-decoration-color:initial"> 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</pre>
<p>Our "Ground Truth" in this area shows an Elevation at -29.17m.<br>
<br>
<br>
</p>
<p><br>
</p>
<div>On 7/12/22 16:35, Greg Troxel wrote:<br>
</div>
<blockquote type="cite">
<pre>Nicholas Stanley <a href="mailto:nicholas.stanley@luminartech.com" target="_blank"><nicholas.stanley@luminartech.com></a> writes:
</pre>
<blockquote type="cite">
<pre>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.
</pre>
</blockquote>
<pre>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.
</pre>
<blockquote type="cite">
<pre>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.
</pre>
</blockquote>
<pre>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.
</pre>
<blockquote type="cite">
<pre>Following the logic of your process:
*Step 1: Convert the **/NAD83(2011) / Conus Albers + NAVD88 height
/**to NAD83 Ellipsoidal using GEOID18 *
</pre>
</blockquote>
<pre></pre>
<blockquote type="cite">
<pre>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
</pre>
</blockquote>
<pre>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
<a href="https://www.ngs.noaa.gov/NCAT/" target="_blank">https://www.ngs.noaa.gov/NCAT/</a>
<a href="https://www.ngs.noaa.gov/GEOID/GEOID18/" target="_blank">https://www.ngs.noaa.gov/GEOID/GEOID18/</a>
<a href="https://www.ngs.noaa.gov/GEOID/GEOID18/computation.html" target="_blank">https://www.ngs.noaa.gov/GEOID/GEOID18/computation.html</a>
<a href="https://geographiclib.sourceforge.io/cgi-bin/GeoidEval" target="_blank">https://geographiclib.sourceforge.io/cgi-bin/GeoidEval</a>
</pre>
<blockquote type="cite">
<pre>(Not sure which EPSG codes to use here - any guess?)
</pre>
</blockquote>
<pre>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.
</pre>
<blockquote type="cite">
<pre>*Step 2. Convert resulting Step 1. file into destination projection
(WGS84 UTM 10 Ellipsoidal <a href="https://epsg.io/32610" target="_blank">https://epsg.io/32610</a>, or alternatively,
<a href="https://epsg.io/7789" target="_blank">https://epsg.io/7789</a> 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
</pre>
</blockquote>
<pre>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.
</pre>
</blockquote>
</div>
_______________________________________________<br>
pdal mailing list<br>
<a href="mailto:pdal@lists.osgeo.org" target="_blank">pdal@lists.osgeo.org</a><br>
<a href="https://lists.osgeo.org/mailman/listinfo/pdal" rel="noreferrer" target="_blank">https://lists.osgeo.org/mailman/listinfo/pdal</a><br>
</blockquote></div>