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

[Ian Erickson]

Nope.

So I'm trying to distill this down...  There's a bunch going on here so 
I thought I should focus on just the one element that's cauing the 
problem - that being msLoadProjectionString().  So with a simple map 
with a projection specifed as epsg:4326, I figured I'd try to alter the 
projection and then save the output map....

Guess what, the projection remains unchanged....  So there appears to be 
a step (or 5) that I'm missing here...

//-------------------------------------------------------------------

    mapObj *map = NULL;
    imageObj *image = NULL;

    msSetup();
    map = msLoadMap( "C:\\original.map", NULL );
    if (!map) {
        cout << "There is an error loading the map..." << endl;
        return;
    }
    //Trying to change the projection to epsg:4269 (and then verify it 
in the output.map file)...
    msLoadProjectionString( &(map->projection), "+init=epsg:4269" );
    msSaveMap( map, "C:\\output.map" );
   
    image = msDrawMap( map );
    msSaveImage( map, image, "C:\\output.png" );

    msFreeImage( image );
    msFreeMap( map );
    msCleanup();
    cout << "Complete." << endl;

//-------------------------------------------------------------------

Ian Erickson
AnalyGIS, LLC
Mancos/Durango, CO  81328
http://www.analygis.com

tel: 877.MAP.DATA (877.627.3282) x 801
fax: 970.565.9268



Stephen Woodbridge wrote:

>
> Ian,
>
> You might try adding a couple of lines I noted below and see if it helps.
>
> -Steve W.
>
> Ian Erickson wrote:
>
>> 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;
>
>
>      msSetup();
>
>>
>>    // This map and every layer in it has the projection of epsg:4326...
>>    map = msLoadMap("C:\\some.map", NULL);
>
>
>      if(!map) writeErrorMessage("Error loading mapfile.");
>
>>
>>    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!
>>
>
>