[Proj] please help C++ idiot with datum conversion

[Frank Warmerdam]

Mike Zulauf wrote:
> Hi all,
> 
> I have what I believe to be a simple question.  I've got a very simple 
> bit of code C++ that converts a vector containing NAD27 long-lat to 
> NAD27 State Plane coordinates.  It's written so that I can call it from 
> Fortran.  It appears to be working, and giving correct results.
> 
> I need to modify it, however, so that it will also do a datum 
> conversion, from WGS84 lat-lon to NAD27 State Plane coordinates.  I'm 
> sure this is simple, but I'm a newbie with PROJ.  To make matters worse, 
> I'm not a C++ programmer.  I generally work with Fortran.
> 
> Below is the code as it stands now.  What's the best way to modify the 
> code to include the datum conversion?  It's ok to just have it 
> hard-wired - I don't need it to be flexible (ie, it doesn't need to be 
> able to handle NAD27 _or_ WGS84, just WGS84).  But I still need the 
> output to be NAD27.
> 
> Even though it works as it is now (without datum conversion), it's 
> entirely likely I have some stupid error or nonstandard code, since I'm 
> not a C++ programmer.  If there's anything of that sort which might 
> cause problems down the line, I wouldn't mind hearing about that either.

Mike,

I have included an alternate implementation using pj_transform().  The
code compiled but I haven't tested it.  Note:

  o I use proj_api.h, the public include file instead of projects.h which is
    really only for internal use (or backward compatibility).

  o I use pj_transform() instead of pj_fwd() as it includes support for two
    arbitrary coordinate systems and extra stuff like datum conversion.

#include <stdio.h>
#include <proj_api.h>

void conv_coords(double coord_1[], double coord_2[], int *ni,
                  int *conv_stat)
{
     static const char *src_crs = "+proj=latlong +datum=WGS84";
     static const char *dst_crs = "+init=nad27:3601 +datum=NAD27";

     projPJ src_PJ, dst_PJ;
     projUV data;
     int i;

     src_PJ = pj_init_plus( src_crs );
     dst_PJ = pj_init_plus( dst_crs );

     *conv_stat = 1;

     if ( src_PJ == NULL || dst_PJ == NULL )
     {
         printf("pj_init error\n");
         *conv_stat = 0;
     }

     for (i = 0; i <= *ni-1; i++)
     {
         double dummyz = 0.0;

         if( pj_transform( src_PJ, dst_PJ,  1, 0,
                           coord_1 + i, coord_2 + i, &dummyz ) != 0 )
         {
             printf("pj_fwd error\n");
             *conv_stat = 0;
         }

         if (coord_1[i] != HUGE_VAL)
         {
             coord_1[i] = data.u;
             coord_2[i] = data.v;
         }
         else
         {
             printf("pj_fwd error\n");
             *conv_stat = 0;
         }
     }

     pj_free(src_PJ);
     pj_free(dst_PJ);
}


Good luck,

-- 
---------------------------------------+--------------------------------------
I set the clouds in motion - turn up   | Frank Warmerdam, warmerdam at pobox.com
light and sound - activate the windows | http://pobox.com/~warmerdam
and watch the world go round - Rush    | Geospatial Programmer for Rent