Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

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

share|improve this question
    
Why don't you have the towgs84 parameter in the cs2cs version? –  mkennedy Apr 19 '13 at 17:22
    
Hm, seems I forgot it there. But that doesn't change anything. See the solving answer below. –  t1mg Apr 22 '13 at 8:25

1 Answer 1

up vote 1 down vote accepted

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
share|improve this answer
    
Thanks! Totally makes sense, since RD doesn't use angles. –  t1mg Apr 22 '13 at 8:24

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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