[Proj] Creating a nadgrid shift file

Andre Joost andre+joost at nurfuerspam.de
Sun Jan 20 10:28:17 PST 2013


Hi,

I think I found a solution, if I have understood your test data correct:

What you have and what you want:

source:
  id |  lon   |   lat   | deg E |  deg N |   dms E  | dms N
----+--------+---------+-------+--------+----------+----------+
  ll | 5.3876 | 52.1561 | 19395 | 187762 | 5d23'15" | 52d9'22"
  lr | 5.4022 | 52.1561 | 19448 | 187762 | 5d24'08" | 52d9'22"
  ur | 5.4022 | 52.1651 | 19448 | 187794 | 5d24'08" | 52d9'54"
  ul | 5.3876 | 52.1651 | 19395 | 187794 | 5d23'15" | 52d9'54"

destination:
  id |  lon   |   lat   | deg E |  deg N |   dms E  | dms N
----+--------+---------+-------+--------+----------+----------+
  ll | 5.37299| 52.14711| 19343 | 187730 | 5d22'23" | 52d8'50"
  lr | 5.41681| 52.14711| 19501 | 187730 | 5d25'1"  | 52d8'50"
  ur | 5.41682| 52.17409| 19501 | 187826 | 5d25'1"  | 52d10'27"
  ul | 5.37298| 52.17409| 19343 | 187826 | 5d22'23" | 52d10'27"

I expanded the data a bit:
S 52.13811°N =187697"N
N 52.18309°N =187859"N
E -5.43142°E =-19553"E
W -5.35838°E =-19290"E
with zero values for the border cells

My ASCII ntv2 file has the following header:

NUM_OREC 11
NUM_SREC 11
NUM_FILE  1
GS_TYPE SECONDS
VERSION
SYSTEM_F
SYSTEM_T
MAJOR_F        0.000
MINOR_F        0.000
MAJOR_T        0.000
MINOR_T        0.000
SUB_NAME
PARENT  NONE
CREATED
UPDATED
S_LAT     187697.000000
N_LAT     187859.000000
E_LONG    -19553.000000
W_LONG    -19290.000000
LAT_INC       32.400000
LONG_INC      52.600000
GS_COUNT    36

dx and dy are swapped in the file, upper right first, and dx values 
inverted:
-32.370100 -52.59180  0.000000  0.000000
-32.363600 52.602400  0.000000  0.000000
   32.35040 -52.64480  0.000000  0.000000
   32.35610 52.634210  0.000000  0.000000
(and a lot of null value lines between)

Converted with GDAy, because the old GDAit won't run on Windows 7 
anymore :-(

For cs2cs, I made an input file with your WGS84 values, and used 
Amersfoort lat/lon instead of the projected CRS.

Put the following in a batch file:

echo epsg4326 RDgrid.txt >out.txt
cs2cs +init=epsg:4326 +to +proj=longlat +ellps=bessel 
+nadgrids=D:\Karten\gdal\ntv2\RD.gsb +no_defs RD-grid.txt >>out.txt

echo Amersf ll nowgs84 - grid >>out.txt
cs2cs +proj=longlat +ellps=bessel +towgs84=0,0,0,0,0,0,0 +no_defs +to 
+proj=longlat +ellps=bessel +nadgrids=D:\Karten\gdal\ntv2\RD.gsb 
+no_defs RD-grid.txt >>out.txt

echo Amersf ll - grid >>out.txt
cs2cs +proj=longlat +ellps=bessel 
+towgs84=565.417,50.3319,465.552,-0.398957,0.343988,-1.8774,4.0725 
+no_defs +to +proj=longlat +ellps=bessel 
+nadgrids=D:\Karten\gdal\ntv2\RD.gsb +no_defs RD-grid.txt >>out.txt

echo grid -   Amersf ll nowgs84 >>out.txt
cs2cs +proj=longlat +ellps=bessel +nadgrids=D:\Karten\gdal\ntv2\RD.gsb 
+no_defs +to +proj=longlat +ellps=bessel +towgs84=0,0,0,0,0,0,0 +no_defs 
RD-grid.txt >>out.txt

echo grid 2 Amersf ll >>out.txt
cs2cs +proj=longlat +ellps=bessel +nadgrids=D:\Karten\gdal\ntv2\RD.gsb 
+no_defs +to +proj=longlat +ellps=bessel 
+towgs84=565.417,50.3319,465.552,-0.398957,0.343988,-1.8774,4.0725 
+no_defs RD-grid.txt >>out.txt

gives this result:
epsg4326 RDgrid.txt
5d23'15.36"E	52d9'21.96"N 0.000
5d24'7.92"E	52d9'21.96"N 0.000
5d24'7.92"E	52d9'54.36"N 0.000
5d23'15.36"E	52d9'54.36"N 0.000
Amersf ll nowgs84 - grid
5d23'15.36"E	52d9'24.046"N -698.431
5d24'7.92"E	52d9'24.046"N -698.431
5d24'7.92"E	52d9'56.446"N -698.420
5d23'15.36"E	52d9'56.446"N -698.420
Amersf ll - grid
5d23'13.793"E	52d9'18.402"N 43.347
5d24'6.345"E	52d9'18.402"N 43.346
5d24'6.345"E	52d9'50.798"N 43.331
5d23'13.793"E	52d9'50.799"N 43.331
grid -   Amersf ll nowgs84
5d22'23.077"E	52d8'47.83"N 698.441
5d25'0.512"E	52d8'47.823"N 698.441
5d25'0.565"E	52d10'24.624"N 698.411
5d22'23.046"E	52d10'24.63"N 698.411
grid 2 Amersf ll
5d22'24.637"E	52d8'53.47"N -43.360
5d25'2.095"E	52d8'53.464"N -43.359
5d25'2.148"E	52d10'30.276"N -43.313
5d22'24.606"E	52d10'30.282"N -43.314

Compared to the start, the "grid 2 Amersfoort" is the way to get the 
desired coordinates. Don't ask me why...

Maybe a case of "If you give me rubbish, I do nothing"-behaviour of proj.

HTH,
André Joost






More information about the Proj mailing list