[PROJ] NTv2 grid file from Bessel ellipsoid to Bessel ellipsoid

[Jochem]

Hi,

I noticed something odd in cs2cs results, which has nothing to do with the
new Proj. version, but could be a old bug. I'm using the OSGeo4W Shell of
QGIS 3.0.1 (Proj.493). 

BACKGROUND
I am testing with different horizontal grid shift files, since we are
defining a new official coordinate tranformation between ETRS89
(ETRF2000(R14)) and the Dutch RD (new) grid and NAP heights. We like to
provide two variants. First, a easy implementation variant with a grid shift
file that handles both the datum transformation (7 parameter Helmert
transformation from Bessel to GRS80) and the datum distortions (correction
grid on Bessel ellipsoid) in a single NTv2 (.gsb) grid file. Since this
method can only be used within the bounding box of the grid file, we also
define a grid without the datum transformation and publish the 7 paramaters
for the datum transformation separately. This 7 parameter transformation
together with the second grid should give the same results as the first grid
within the bounding box of the grid file, but this second method will also
provide useful results outside the bounding box of the grid.

PROBLEM
I do not get the results I expected. This is caused by something odd in
cs2cs. To demonstrate this, I created a NTv2 grid file with only zero
corrections (see ascii grid file below). If I transform coordinates on the
Bessel ellipsoid with this grid (command 1), I do not get unchanged results
(command 2), but I get results which are the same as I get when changing the
ellipsoid (command 3). 

Command 1:
cs2cs +proj=longlat +ellps=bessel +wktext +nadgrids=zeros.gsb +units=m
+no_defs +to +proj=longlat +ellps=bessel +towgs84=0,0,0,0,0,0,0 +units=m
+no_defs -f %.9f
5.000000000     52.000000000 0.0000
result:
5.000000000     51.999419687 698.6070

Command 2:
cs2cs +proj=longlat +ellps=bessel +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
+to +proj=longlat +ellps=bessel +towgs84=0,0,0,0,0,0,0 +units=m +no_defs -f
%.9f
5.000000000     52.000000000 0.0000
result:
5.000000000     52.000000000 0.0000

Command 3:
cs2cs +proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs
+to +proj=longlat +ellps=bessel +towgs84=0,0,0,0,0,0,0 +units=m +no_defs -f
%.9f
5.000000000     52.000000000 0.0000
result:
5.000000000     51.999419686 698.6069

QUESTION
Is there something wrong with my expectations on how cs2cs should work or is
this a bug?


Kind regards, Jochem
J. Lesparre (MSc in geodetic engineering)
Rijksdriehoeksmeting (RD)
Netherlands Cadastre, Land Registry and Mapping Agency (Kadaster)
Netherlands Partnership Geodetic Infrastructure (NSGI)


GRID FILE zeros.gsb
---
NUM_OREC      11
NUM_SREC      11
NUM_FILE       2
GS_TYPE SECONDS 
VERSION NTv2.0  
SYSTEM_FBSSL1841
SYSTEM_TBSSLPSEU
MAJOR_F  6377397.155
MINOR_F  6356078.963
MAJOR_T  6377397.155
MINOR_T  6356078.963
SUB_NAMENL_EEZ  
PARENT  NONE    
CREATED 21-02-19
UPDATED         
S_LAT     180000.000000
N_LAT     201600.000000
E_LONG    -28800.000000
W_LONG     -7200.000000
LAT_INC      360.000000
LONG_INC     360.000000
GS_COUNT    3721
0.000000  0.000000  0.001000  0.001000
0.000000  0.000000  0.001000  0.001000
0.000000  0.000000  0.001000  0.001000
...





--
Sent from: http://osgeo-org.1560.x6.nabble.com/PROJ-4-f3840930.html