Automatic datum shift when +nadgrids=@nul is not explicitly given

[proj4-bot]

Reported by janh on 27 Nov 2008 14:27 UTC
An automatic datum shift to WGS84 seems to take place, even when not requested.

Test point 9018/375648 (Dutch RD, EPSG:28992 without WGS84 datum shift)
The source PROJ string is:
+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs.
I tested this string both without and with +nadgrids=@null on four target projections:

1. to latlon Bessel (+proj=longlat +ellps=bessel +no_defs)
   --> (without nadgrids in source) 3.291612 51.352060
   --> (with nadgrids null in source) 3.291612 51.351476

2. to latlon Bessel with nadgrids null (+proj=longlat +ellps=bessel +no_defs +nadgrids=@null)
   --> (without nadgrids in source) 3.291612 51.352643
   --> (with nadgrids null in source) 3.291612 51.352060

3. to latlon Bessel on the WGS84 datum (+proj=longlat +ellps=bessel +no_defs +datum=WGS84)
   --> (without nadgrids in source) 3.291612 51.352060
   --> (with nadgrids null in source) 3.291612 51.351476

4. to latlon on the WGS84 datum (+proj=longlat +no_defs +datum=WGS84)
   --> (without nadgrids in source) 3.291612 51.352643
   --> (with nadgrids null in source) 3.291612 51.352060

5. and 3) give the same results; so do 2) and 4), which seems to inicate that an automatic datum shift to WGS84 takes place, even when not requested. I'm not sure about the effects of +nadgrids=@null. Specifying +ellips=bessel without adding +nadgrids=@null seems to cause an automatic WGS84 datum shift, whether or no datum=WGS84 is specified (1a = 3a). If +nadgrids=@null is specified in the the source projection, but not in the target projection (1b = 3b), the latlon coordinates shift about 65m to the south. When +nadgrids=@null is omitted in the source projection, but added in the target projection, the resulting coordinates are shifted of 65M to the north (2a). Finally, when +nadgrids=@null is specified in both source and target projections, the datum shift is applied nevertheless (2b is equal to 1a and 3a). Finally, specifying a WGS84 datum shift for the target projection without the Bessel ellipsoide gives the same values as 2.

Migrated-From: https://trac.osgeo.org/proj/ticket/22

[proj4-bot]

Comment by warmerdam on 27 Nov 2008 14:36 UTC
Jan,

Could you provide the exact command used for at least a couple of these? Did you use cs2cs? Also, was this 4.6.1?

[proj4-bot]

Comment by janh on 27 Nov 2008 15:23 UTC
Sorry Frank, I have been mixing versions. I use Proj 4.6.1, but the tests were done with an older version of PostGIS, without realizing that it probably uses an older version of the library. I am not sure of all is well with 4.6.1 though. I reran the same test with cs2cs under proj 4.6.1 with the attached shell-file. It gives better results than my first run, but I suspect the projection still uses a WGS84 datum shift by default. Note that all transformations result in the same coordinates as when a datum shift to WGS84 is explicitly requested, except when a) a different datum is requested in the target projection and b) a +nadgrid=@null is specified in the source projection

[proj4-bot]

Attachment added by janh on 27 Nov 2008 15:23 UTC
https://trac.osgeo.org/proj/attachment/ticket/22/testdatum.sh

[kbevers]

Tried to reproduce this with the script below.

echo "1."
echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs +to +proj=longlat +ellps=bessel +no_defs

echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +nadgrids=@null +no_defs +to +proj=longlat +ellps=bessel +no_defs

echo " "
echo "2."
echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs +to +proj=longlat +ellps=bessel  +nadgrids=@null +no_defs

echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +nadgrids=@null +no_defs +to +proj=longlat +ellps=bessel +nadgrids=@null +no_defs

echo " "
echo "3."
echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs +to +proj=longlat +ellps=bessel +datum=WGS84 +no_defs

echo 9018 375648 0 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +nadgrids=@null +no_defs +to +proj=longlat +ellps=bessel +datum=WGS84 +no_defs

echo " "
echo "4."
echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs +to +proj=longlat +datum=WGS84 +no_defs

echo 9018 375648 |./cs2cs -f %.6f +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.999908 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +nadgrids=@null +no_defs +to +proj=longlat +ellps=bessel +no_defs

Which results in the following output:

1.
3.291612        51.352060 0.000000
3.291612        51.352060 0.000000

2.
3.291612        51.352060 0.000000
3.291612        51.352060 0.000000

3.
3.291612        51.352060 0.000000
3.291612        51.351477 699.339164

4.
3.291612        51.352060 0.000000
3.291612        51.352060 0.000000

I conclude that what is discussed in this issue is mostly fixed. The unexpected results in 3) is probably caused by the ellipsoid parameters being adjusted in pj_datum_transform() although the coordinates are not changed. Effectively this results in converting the coordinates from the WGS84 ellipsoid to the Bessel ellipsoid, which is backed up by

echo 3.291612 51.352060 0.000000 |./cct +proj=pipeline +step +proj=cart +ellps=WGS84 +step +proj=cart +ellps=bessel +inv
  3.2916120000   51.3514765406      699.3392           inf