[gdal-dev] Upgrading to GDAL 3.4.2 from 1.10.0 shows different outputs from gdaltransform and corresponding C++ API

Mon Aug 26 09:26:47 PDT 2024

Shawn,

The results depends on which PROJ grids are available on your system (or 
if you allowed network access).

1) The result you get with GDAL 1.10 can be reproduced with modern PROJ 
versions by allowing (and restricting) to using the CONUS / 
us_noaa_conus.tif grid which corresponds to the "NAD27 to WGS 84 (79)" 
EPSG transformation, with an advertized accuracy of 5.0 m:

$ echo 32.804191 -117.258911 0 | PROJ_NETWORK=ON cct -d 8 +proj=pipeline 
+step +proj=axisswap +order=2,1 +step +proj=unitconvert +xy_in=deg 
+xy_out=rad +step +proj=hgridshift +grids=us_noaa_conus.tif +step 
+proj=unitconvert +xy_in=rad +xy_out=deg +step +proj=axisswap +order=2,1
32.80423938   -117.25978145    0.00000000           inf

2) The result you get with GDAL 3.4.2 / PROJ 8.1.1 can be actually 
reproduced by using only the "NAD27 to WGS 84 (6)" EPSG Helmert 
transformation, with an advertized accuracy of 7.0 m, which is the most 
precise one in the absence of any grid:

$ echo 32.804191 -117.258911 | cs2cs -d 8 NAD27 WGS84
32.80423884    -117.25976448 0.00000000

or explicitly with the pipeline show by "projinfo NAD27 WGS84 --bbox 
-118,32,-117,33 --spatial-test intersects --single-line":

$ echo 32.804191 -117.258911 0 |  cct -d 8 +proj=pipeline +step 
+proj=axisswap +order=2,1 +step +proj=unitconvert +xy_in=deg +xy_out=rad 
+step +proj=push +v_3 +step +proj=cart +ellps=clrk66 +step +proj=helmert 
+x=-8 +y=159 +z=175 +step +inv +proj=cart +ellps=WGS84 +step +proj=pop 
+v_3 +step +proj=unitconvert +xy_in=rad +xy_out=deg +step +proj=axisswap 
+order=2,1
32.80423884   -117.25976448    0.00000000           inf

3) You would get the """"most precise result"""" (with lots of quotes, 
because this is highly subjective, and assumes that NAD83(HARN) is 
equivalent to WGS 84, according to the comment in the EPSG database) 
with the "NAD27 to WGS 84 (41)" EPSG transformation, with an advertized 
accuracy of 1.5m, which combines the CONUS grid and us_noaa_cshpgn:

$ echo 32.804191 -117.258911 | PROJ_NETWORK=ON cs2cs -d 8 NAD27 WGS84
32.80424065    -117.25978105 0.00000000

or explicitly with the pipeline show by "projinfo NAD27 WGS84 --bbox 
-118,32,-117,33 --spatial-test intersects --single-line":

$ echo 32.804191 -117.258911 0 | PROJ_NETWORK=ON cct -d 8 +proj=pipeline 
+step +proj=axisswap +order=2,1 +step +proj=unitconvert +xy_in=deg 
+xy_out=rad +step +proj=hgridshift +grids=us_noaa_conus.tif +step 
+proj=hgridshift +grids=us_noaa_cshpgn.tif +step +proj=unitconvert 
+xy_in=rad +xy_out=deg +step +proj=axisswap +order=2,1
32.80424065   -117.25978105    0.00000000           inf

==> All those results are valid. Some are more precise than others, 
under certain assumptions... But basically speaking about WGS 84 
introduces an intrinsic inaccuracy of about 2 meters (cf 
https://www.youtube.com/watch?v=M2ck3cAGvhg). I don't think there are 
any published time-dependent transformations for NAD27 to WGS84 or any 
of its realization.

Even

Le 24/08/2024 à 03:14, Fox, Shawn D (US) via gdal-dev a écrit :
>
> Hello GDAL users.
>
> I have noticed that from GDAL 1.10 – GDAL 3.1.0 that gdaltransform is 
> producing different results for NAD27 GCS to WGS84 GCS 
> transformations.  The differences in longitude in the output is more 
> significant.  Having tested the differences in feet between the points 
> output from each version I can see a difference of ~5 feet primarily 
> due to the change in the longitude value.
>
> Would anyone look at the following results and share your thoughts 
> about why the result of this transformation is different in a newer 
> version?  Could this be related to the concept of WGS 84 realization 
> being different in newer versions of GDAL?  I am also linking GDAL 
> 3.4.2 with PROJ 8.1.1 at this time.  These results were produced with 
> GDAL 3.4.2, but the same results occur with GDAL 3.1.0 as well.  GDAL 
> 3.4.2 is what I am currently trying to upgrade to.  The results from 
> the C++ programming API are similar in difference.  I just used 
> gdaltransform for the purposes of this message because it was a more 
> convenient way of producing something that others might be able to 
> reproduce.
>
> *GDAL 1.10.0 results*
>
> ./gdaltransform -s_srs NAD27 -t_srs WGS84
>
> -117.258911 32.804191
>
> -117.25978144845 32.8042393819884 1.32247805595398e-07
>
> *GDAL 3.4.2 results*
>
> ./gdaltransform -s_srs NAD27 -t_srs WGS84
>
> -117.258911 32.804191
>
> -117.259764476223 32.8042388409449 0
>
> I have read that the support for specifying epochs related to 
> realizations is not available until GDAL 3.8.0.
>
> Here are some references I’ve read so far.
>
> Transforming between WGS84 Realizations | Journal of Surveying 
> Engineering | Vol 148, No 2 (ascelibrary.org) 
> <https://ascelibrary.org/doi/10.1061/%28ASCE%29SU.1943-5428.0000389>
>
> RFC 81: Support for coordinate epochs in geospatial formats — GDAL 
> documentation 
> <https://gdal.org/development/rfc/rfc81_coordinate_epoch.html>
>
> Thanks,
>
> Shawn Fox
>
> San Diego, CA
>
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev

-- 
http://www.spatialys.com
My software is free, but my time generally not.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20240826/6677cc17/attachment-0001.htm>