[PROJ] PROJ 6.1.0RC1

[Even Rouault]

Bas,

> IIRC epsg:4258 is interpreted as a Geographic 2D CRS which likely causes
> the 3D component to be ignored.

Yes, you should use EPSG:4937 (ETRS89 Geographic 3D CRS), but that will not 
solve the issue with the Z value. The issue here is that in the rdnap:rdnap 
string, the +geoidgrids=naptrans2008.gtx implies a transformation between 
WGS84 3D and the vertical CRS implied by naptrans2008.gtx. But you want to 
transform between ETRS89 3D and that vertical CRS. As there is no 
transformation registered in the EPSG dataset between WGS84 3D and ETRS89 3D 
(there is one only on the 2D part), PROJ cannot use the vertical grid

What works however is completely dropping the use of rdnap:rdnap and using 
EPSG:28992+5709 instead, that is "Amersfoort / RD New + NAP height", as there 
is a transformation registered between ETRS89 and NAP height, which uses the 
naptrans2008.gtx grid

$ echo "53.160753042 4.824761912 42.8614" | \
      PROJ_LIB=$PWD:data src/cs2cs -f "%.4f" EPSG:4937 EPSG:28992+5709
117380.1203	575040.3399 0.9985

Which is close to what you get with PROJ 5.2.0 with
cs2cs -r  -f %.4f +init=epsg:4258 +to +init=rdnap:rdnap:
117380.1203	575040.3398 1.0000

> 
> The projinfo output is not something I can make sense of, it is as
> follows:
> 
> $ PROJ_LIB=.:/usr/share/proj projinfo -s EPSG:4258 -t rdnap:rdnap
> target CRS: parsing of user string failed: crs not found

Note: you can use +init=rdnap:rdnap to resolved in the rdnap file.

Even

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com