Using the "C API" and modifying the projection...

[Ian Erickson]

Gentlemen,

I know that the C API was never really intended for use directly as a 
developer's platform, although when evaluating all of my options it was 
clear that what I needed to do was simple and the C API (map.h, 
mapproject.h) were ideal for my situation....however....

I've put together a simple console application to test the API and 
(whatever I try) I cannot get the map to change projections and specify 
the extent of the new map in projected coordinates.

Here's an extract...
//----------------------------------------------------------------------------
void simpleMap( double minx = -180.0, double miny = -90.0, double maxx = 
180.0, double maxy = 90.0 )
{
    mapObj *map = NULL;
    imageObj *image = NULL;

    // This map and every layer in it has the projection of epsg:4326...
    map = msLoadMap("C:\\some.map", NULL);

    pointObj *pt_ll = pointObj_new();
    pointObj *pt_ur = pointObj_new();

    // *** Assign the coordinates to be reprojected.
    pt_ll->x = minx;
    pt_ll->y = miny;
    pt_ur->x = maxx;
    pt_ur->y = maxy;

    cout << "BEFORE:" << endl;
    cout << "pt_ll: (" << pt_ll->x << ", " << pt_ll->y << ")" << endl;
    cout << "pt_ur: (" << pt_ur->x << ", " << pt_ur->y << ")" << endl;
    cout << endl;

    // *** Just to confirm, that the mapObj is working as expected, I 
test this statement and indeed,
    // *** the map is recentered correctly....
    msMapSetExtent( map, pt_ll->x, pt_ll->y, pt_ur->x, pt_ur->y );


    // *** Create projectionObj pointers to convert between geographic 
coordinates
    // *** into the Transverse Mercator projection required for USGS quads.
    // *** projectionObj_new is defined elsewhere and functions 
properly...as the
    // *** transformed coordinates are correct.
    projectionObj *projIn = projectionObj_new( "+init=epsg:4326" );

    string s_projOut = "+proj=tmerc +lat_0=0.0";
    s_projOut.append( " +lon_0=" );
    s_projOut.append( doubleToString( (double)((pt_ll->x + pt_ur->x) / 
2.0) ) );
    s_projOut.append( " +k=0.9996 +x_0=0.0 +y_0=0.0 +ellps=GRS80" );
    projectionObj *projOut = projectionObj_new( _strdup( 
s_projOut.c_str() ) );

    cout << "Output Projection: " << msGetProjectionString( projOut ) << 
endl;


    // *** Reproject the input points.  One point for the lower-left 
coordinate
    // *** and one for the upper-right corner.
    msProjectPoint( projIn, projOut, pt_ll );
    msProjectPoint( projIn, projOut, pt_ur );

    // *** These points are correct! No problem here....
    cout << "AFTER: " << endl;
    cout << "p_min: (" << pt_ll->x << ", " << pt_ll->y << ")" << endl;
    cout << "p_max: (" << pt_ur->x << ", " << pt_ur->y << ")" << endl;
    cout << endl;

    // *** This is where things fall apart.  Taking a page from the 
php_mapscript.c and
    // *** mapscript_i.c sources, I init the projectionObj pointers, run 
msLoadProjectionString,
    // *** and recompure the extent of the map....msLoadProjectionString 
by itself doesn't work either.
    int nStatus = 0;
    int nUnits = MS_METERS;
    projectionObj in;
    projectionObj out;
    rectObj sRect;
    int bSetNewExtents = 0;
    int bSetUnitsAndExtents = 1;

    msInitProjection( &in );
    in = map->projection;
    msInitProjection( &out );

    msLoadProjectionString( &(out), msGetProjectionString( projOut ) );
    sRect = map->extent;

    if (in.proj != NULL && out.proj != NULL)
    {
        if (msProjectionsDiffer( &in, &out ) )
        {
            if (msProjectRect( &in, &out, &sRect) == 0)
                bSetNewExtents = 1;
        }
    }

    msLoadProjectionString( &map->projection, msGetProjectionString( 
projOut ) );

    nUnits = GetMapserverUnitUsingProj( &map->projection );
    map->units = (MS_UNITS)nUnits;
    if (bSetNewExtents == 1)
    {
        map->extent = sRect;
        map->cellsize = msAdjustExtent( &map->extent, map->width, 
map->height );
        msCalculateScale( map->extent, map->units, map->width, 
map->height, map->resolution, &map->scale );
    }


    // *** Compute the size of the output map in meters.  This is later 
used to
    // *** determine the output size of the map in pixels.
    double width_m = pt_ur->x - pt_ll->x;
    double height_m = pt_ur->x - pt_ll->y;

    double width = 500.0;
    double height = ((height_m / width_m) * width);
    cout << "EARTH SIZE [Width, Height] m: " << ((double)width_m) << ", 
" << ((double)height_m) << endl;
    cout << endl;


    // *** Set the size of the output map.
    msMapSetSize( map, (int)width, (int)height );
    cout << "MAP SIZE [Width, Height]: " << map->width << ", " << 
map->height << endl;


    // *** Set the extent of the map as per the newly projected coordinates
    // *** computed above.

    // *** This call is where things go wrong.  Apparently the 
msLoadProjectionString is not
    // *** working properly, because if I use the original coordinates 
supplied as parameters
    // *** the map moves to the correct location -- not the projected 
location.
    msMapSetExtent( map, pt_ll->x, pt_ll->y, pt_ur->x, pt_ur->y );


    image = msDrawMap( map );
    msSaveImage( map, image, "C:\\simple.png");

    pointObj_destroy( pt_ll );
    pointObj_destroy( pt_ur );
    msFreeImage( image );
    msFreeMap( map );
    msCleanup();
    cout << "Complete!" << endl;
}
// ------------------------------------------------------------------

Please! I'm looking for some guidance here as it doesn't seem too 
difficult - and the same order and function calls work in PHP/MapScript!

-- 
Ian Erickson