5

I try to convert between dutch RD New (EPSG:28992) and WGS84 (EPSG:4326) via PROJ4 C++ API. On the console this works with:

cs2cs +proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +units=m +no_defs +to +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
138494.92605 467792.640021
5d8'46.283"E    52d11'56.38"N 0.000

However, doing the same thing in C++ with:

#include <iostream>
#include <proj_api.h>

using namespace std;

int main()
{
  double x = 138494.92605;
  double y = 467792.640021;

  char *srid28992 = "+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +towgs84=565.04,49.91,465.84,-1.9848,1.7439,-9.0587,4.0772 +units=m +no_defs";
  char *srid4326 = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs";

  projPJ source = pj_init_plus(srid28992);
  projPJ target = pj_init_plus(srid4326);

  if(source==NULL || target==NULL)
    return false;

  x *= DEG_TO_RAD;
  y *= DEG_TO_RAD;

  int success = pj_transform(source, target, 1, 1, &x, &y, NULL );

  x *= RAD_TO_DEG;
  y *= RAD_TO_DEG;

  cout << success << endl << x << endl << y << endl;
}

results in wrong numbers:

./main 
0
3.34135
48.0483

I wasn't able to find a working example for this. What am I missing here?

PROJ4 Version is 4.6.1

3
  • Why don't you have the towgs84 parameter in the cs2cs version?
    – mkennedy
    Apr 19, 2013 at 17:22
  • Hm, seems I forgot it there. But that doesn't change anything. See the solving answer below.
    – t1mg
    Apr 22, 2013 at 8:25
  • @t1mg: Hi! I had problems building your program. After I finally did it, I summarized my solution here: github.com/thomasfermi/proj4_cmake_example Please let me know if you are not ok with me using your sample program in that github repo. The repo is mainly about the CMake and vcpkg troubles and not the actual program. Nov 6, 2019 at 11:46

1 Answer 1

3

The source (inputs) are not in degrees. You shouldn't convert them before the pj_transform() call:

#include <iostream>
#include <proj_api.h>

using namespace std;

int main()
{
  double x = 138494.92605;
  double y = 467792.640021;

  char *srid28992 = "+proj=sterea +lat_0=52.15616055555555 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel +towgs84=565.04,49.91,465.84,-1.9848,1.7439,-9.0587,4.0772 +units=m +no_defs";
  char *srid4326 = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs";

  projPJ source = pj_init_plus(srid28992);
  projPJ target = pj_init_plus(srid4326);

  if(source==NULL || target==NULL)
    return false;

//  x *= DEG_TO_RAD;
//  y *= DEG_TO_RAD;

  int success = pj_transform(source, target, 1, 1, &x, &y, NULL );

  x *= RAD_TO_DEG;
  y *= RAD_TO_DEG;

  cout << success << endl << x << endl << y << endl;
}

Then building and running it:

$ g++ -Wall -W test.cc -lproj
test.cc: In function ‘int main()’:
test.cc:11:21: warning: deprecated conversion from string constant to ‘char*’ [-Wwrite-strings]
test.cc:12:20: warning: deprecated conversion from string constant to ‘char*’ [-Wwrite-strings]
bradh@ornata:~/proj4$ ./a.out 
0
5.14431
52.1976
1
  • Thanks! Totally makes sense, since RD doesn't use angles.
    – t1mg
    Apr 22, 2013 at 8:24

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.