[PROJ] National CRS to ITRS

[Even Rouault]

Le 21/02/2023 à 11:00, Lesparre, Jochem a écrit :
>
> Dear Even,
>
> Sorry, I was using PROJ 8.2.1 as I didn’t expect the PROJ behaviour 
> had changed since that version (see updated command results below).
>
> We use the ETRF2000 realisation like most other NMAs, as this is 
> recommended by EUREF for georeferencing, not ETRF2014! Therefore, 
> ETRS89 should not be transformed to ETRF2014 with a null transformation.
>
The null transformation between ETRS89 and ETRF2014 comes with the 
definition of the datum ensemble ETRS89 (cf change 
https://github.com/OSGeo/PROJ/issues/3263)


> The correct transformation routes are:
>
>   * RD -> Amersfoort -> ETRF2000 -> ITRF2000 -> ITRF2014
>   * NAP -> ETRF2000 -> ITRF2000 -> ITRF2014
>
> I expected this could be achieved by changing the Amersfoort-ETRS89 
> transformation to Amersfoort-ETRF2000 “instead” in EPSG. I don’t 
> understand what advantage a transformation to the ETRS89 datum 
> ensemble “in addition” to a transformation to the specific ETRF2000 
> realisation in EPSG would give.
>
Changing from ETRS89 to ETRF2000 could (potentially) negatively impact 
workflows RDNAP to ETRS89 that behave properly currently. Anyway you can 
simulate that by hacking your proj.db and modifying the existing records 
or adding ones to see the effect of your changes
>
> *In PROJ 9.1.1. RDNAP to ITRF2014:*
>
> echo 128410.0957 445806.4960 -0.4754 2023.00 | cs2cs epsg:7415 
> epsg:7912 -d 9
>
> Gives:
>
> 52.000004918 5.000008844 43.0018 2023.00
>
> This should be:
>
> 52.000005512 5.000008256 43.0197 2023.00
>
> Like RDNAP via ETRF2000 to ITRF2014:
>
> echo 128410.0957 445806.4960 -0.4754 2023.00 | cs2cs epsg:7415 
> epsg:7931 -d 9 | cs2cs epsg:7931 epsg:7912 -d 9
>
> *In PROJ 9.1.1. RDNAP to WGS84:*
>
> echo 128410.0957 445806.4960 -0.4754 2023.00 | cs2cs epsg:7415 
> epsg:4326 -d 9
>
> Gives:
>
> 51.999999889 5.000000153 -0.4754 2023.00
>
> This should be:
>
> 52.000000000 5.000000000 43.0000 2023.00
>
> Like RDNAP via ETRF2000 to WGS84:
>
> echo 128410.0957 445806.4960 -0.4754 2023.00 | cs2cs epsg:7415 
> epsg:7931 -d 9 | cs2cs epsg:7931 epsg:4326 -d 9
>
If you only specify EPSG:4326, the vertical transformation is not used, 
so you get the following pipeline:

unknown id, Inverse of RD New + Amersfoort to WGS 84 (4), 1 m, 
Netherlands - onshore, including Waddenzee, Dutch Wadden Islands and 
12-mile offshore coastal zone.

PROJ string:
+proj=pipeline
   +step +inv +proj=sterea +lat_0=52.1561605555556 +lon_0=5.38763888888889
         +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel
   +step +proj=push +v_3
   +step +proj=cart +ellps=bessel
   +step +proj=helmert +x=565.4171 +y=50.3319 +z=465.5524 
+rx=0.398957388243134
         +ry=-0.343987817378283 +rz=1.87740163998045 +s=4.0725
         +convention=coordinate_frame
   +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


To get a 3D transformation, you need to specify EPSG:4979 as target:

unknown id, Inverse of RD New + Amersfoort to ETRS89 (9) + Inverse of 
ETRS89 to NAP height (2) + ETRS89 to WGS 84 (1), 1.002 m, Netherlands - 
onshore, including Waddenzee, Dutch Wadden Islands and 12-mile offshore 
coastal zone., at least one grid missing

PROJ string:
+proj=pipeline
   +step +inv +proj=sterea +lat_0=52.1561605555556 +lon_0=5.38763888888889
         +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel
   +step +proj=hgridshift +grids=nl_nsgi_rdtrans2018.tif
   +step +proj=vgridshift +grids=nl_nsgi_nlgeo2018.tif +multiplier=1
   +step +proj=unitconvert +xy_in=rad +xy_out=deg
   +step +proj=axisswap +order=2,1

Even

-- 
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/proj/attachments/20230221/a81691de/attachment-0001.htm>