mitab_coordsys.cpp

在linux环境下

/**********************************************************************
 * $Id: mitab_coordsys.cpp,v 1.30 2005/03/22 23:24:54 dmorissette Exp $
 *
 * Name:     mitab_coordsys.cpp
 * Project:  MapInfo TAB Read/Write library
 * Language: C++
 * Purpose:  Implementation translation between MIF CoordSys format, and
 *           and OGRSpatialRef format.
 * Author:   Frank Warmerdam, [e-mail]
 *
 **********************************************************************
 * Copyright (c) 1999-2001, Frank Warmerdam
 *
 * Permission is hereby granted, free of charge, to any person obtaining a
 * copy of this software and associated documentation files (the "Software"),
 * to deal in the Software without restriction, including without limitation
 * the rights to use, copy, modify, merge, publish, distribute, sublicense,
 * and/or sell copies of the Software, and to permit persons to whom the
 * Software is furnished to do so, subject to the following conditions:
 * 
 * The above copyright notice and this permission notice shall be included
 * in all copies or substantial portions of the Software.
 * 
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.  IN NO EVENT SHALL
 * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER 
 * DEALINGS IN THE SOFTWARE.
 **********************************************************************
 *
 * $Log: mitab_coordsys.cpp,v $
 * Revision 1.30  2005/03/22 23:24:54  dmorissette
 * Added support for datum id in .MAP header (bug 910)
 *
 * Revision 1.29  2004/06/03 19:36:53  fwarmerdam
 * fixed memory leak processing non-earth coordsys
 *
 * Revision 1.28  2003/03/21 14:20:42  warmerda
 * fixed up regional mercator handling, was screwing up transverse mercator
 *
 * Revision 1.27  2003/01/09 17:33:26  warmerda
 * fixed ellipsoid extraction for datum 999/9999
 *
 * Revision 1.26  2002/12/12 20:12:18  warmerda
 * fixed signs of rotational parameters for TOWGS84 in WKT
 *
 * Revision 1.25  2002/10/15 14:33:30  warmerda
 * Added untested support in mitab_spatialref.cpp, and mitab_coordsys.cpp for
 * projections Regional Mercator (26), Polyconic (27), Azimuthal Equidistant -
 * All origin latitudes (28), and Lambert Azimuthal Equal Area - any aspect (29).
 *
 * Revision 1.24  2002/09/23 13:16:04  warmerda
 * fixed leak in MITABExtractCoordSysBounds()
 *
 * Revision 1.23  2002/04/01 19:49:24  warmerda
 * added support for cassini/soldner - proj 30
 *
 * Revision 1.22  2002/03/01 19:00:15  warmerda
 * False Easting/Northing should be in the linear units of measure in MapInfo,
 * but in OGRSpatialReference/WKT they are always in meters.  Convert accordingly.
 *
 * Revision 1.21  2001/04/04 21:43:19  warmerda
 * added code to set WGS84 values
 *
 * Revision 1.20  2001/01/23 21:23:42  daniel
 * Added projection bounds lookup table, called from TABFile::SetProjInfo()
 *
 * Revision 1.19  2001/01/22 16:00:53  warmerda
 * reworked swiss projection support
 *
 * Revision 1.18  2001/01/19 21:56:18  warmerda
 * added untested support for Swiss Oblique Mercator
 *
 * Revision 1.17  2000/12/05 14:55:27  daniel
 * Added some missing unit name aliases in MITABSpatialRef2CoordSys()
 *
 * Revision 1.16  2000/10/16 21:44:50  warmerda
 * added nonearth support
 *
 * Revision 1.15  2000/10/12 22:01:01  daniel
 * Fixed test on number of bound params in MITABExtractCoordSysBounds()
 *
 * Revision 1.14  2000/09/28 16:39:44  warmerda
 * avoid warnings for unused, and unitialized variables
 *
 * Revision 1.13  2000/08/29 22:30:45  daniel
 * Made MITABCoordSys2SpatialRef() return NULL for a NonEarth coordsys since
 * returning an empty SpatialRef created confusion between LatLon and NonEarth
 *
 * Revision 1.12  2000/03/27 03:31:15  daniel
 * Fixed parsing of CoordSys NonEarth broken by previous change
 *
 * Revision 1.11  2000/03/19 23:26:27  daniel
 * Made MITABCoordSys2SpatialRef() remove leading spaces in coordsys string
 * and produce a CPLError() if it fails parsing string.
 *
 * Revision 1.10  2000/02/07 17:43:37  daniel
 * Fixed offset in parsing of custom datum string in SetSpatialRef()
 *
 * Revision 1.9  2000/01/15 22:30:43  daniel
 * Switch to MIT/X-Consortium OpenSource license
 *
 * Revision 1.8  2000/01/14 14:53:32  warmerda
 * fixed generation of MIF %d datum names
 *
 * Revision 1.7  1999/12/19 17:35:16  daniel
 * Fixed memory leaks
 *
 * Revision 1.6  1999/11/20 04:53:41  daniel
 * Fixed warning on return value in MITABExtractCoordSysBounds()
 *
 * Revision 1.5  1999/11/12 05:51:57  daniel
 * Added MITABExtractCoordSysBounds()
 *
 * Revision 1.4  1999/11/11 02:55:25  warmerda
 * fixed problems with stereographic and survey ft
 *
 * Revision 1.3  1999/11/10 20:13:12  warmerda
 * implement spheroid table
 *
 * Revision 1.2  1999/11/10 02:19:05  warmerda
 * fixed up datum support when reading MIF coord sys
 *
 * Revision 1.1  1999/11/09 22:29:38  warmerda
 * New
 *
 **********************************************************************/

#include "mitab.h"
#include "mitab_utils.h"

extern MapInfoDatumInfo asDatumInfoList[200];
extern MapInfoSpheroidInfo asSpheroidInfoList[200];

/************************************************************************/
/*                             GetMIFParm()                             */
/************************************************************************/

static double GetMIFParm( char ** papszFields, int iField, double dfDefault )

{
    if( iField >= CSLCount(papszFields) )
        return dfDefault;
    else
        return atof(papszFields[iField]);
}

/************************************************************************/
/*                      MITABCoordSys2SpatialRef()                      */
/*                                                                      */
/*      Convert a MIF COORDSYS string into a new OGRSpatialReference    */
/*      object.                                                         */
/************************************************************************/

OGRSpatialReference *MITABCoordSys2SpatialRef( const char * pszCoordSys )

{
    char        **papszFields;
    OGRSpatialReference *poSR;

    if( pszCoordSys == NULL )
        return NULL;
    
/* -------------------------------------------------------------------- */
/*      Parse the passed string into words.                             */
/* -------------------------------------------------------------------- */
    while(*pszCoordSys == ' ') pszCoordSys++;  // Eat leading spaces
    if( EQUALN(pszCoordSys,"CoordSys",8) )
        pszCoordSys += 9;
    
    papszFields = CSLTokenizeStringComplex( pszCoordSys, " ,", TRUE, FALSE );

/* -------------------------------------------------------------------- */
/*      Clip off Bounds information.                                    */
/* -------------------------------------------------------------------- */
    int         iBounds = CSLFindString( papszFields, "Bounds" );

    while( iBounds != -1 && papszFields[iBounds] != NULL )
    {
        CPLFree( papszFields[iBounds] );
        papszFields[iBounds] = NULL;
        iBounds++;
    }

/* -------------------------------------------------------------------- */
/*      Create a spatialreference object to operate on.                 */
/* -------------------------------------------------------------------- */
    poSR = new OGRSpatialReference;

/* -------------------------------------------------------------------- */
/*      Fetch the projection.                                           */
/* -------------------------------------------------------------------- */
    char        **papszNextField;
    int nProjection = 0;

    if( CSLCount( papszFields ) >= 3
        && EQUAL(papszFields[0],"Earth")
        && EQUAL(papszFields[1],"Projection") )
    {
        nProjection = atoi(papszFields[2]);
        papszNextField = papszFields + 3;
    }
    else if (CSLCount( papszFields ) >= 2
             && EQUAL(papszFields[0],"NonEarth") )
    {
        // NonEarth Units "..." Bounds (x, y) (x, y)
        nProjection = 0;
        papszNextField = papszFields + 2;

        if( papszNextField[0] != NULL && EQUAL(papszNextField[0],"Units") )
            papszNextField++;
    }
    else
    {
        // Invalid projection string ???
        if (CSLCount(papszFields) > 0)
            CPLError(CE_Warning, CPLE_IllegalArg,
                     "Failed parsing CoordSys: '%s'", pszCoordSys);
        CSLDestroy(papszFields);
        return NULL;
    }

/* -------------------------------------------------------------------- */
/*      Fetch the datum information.                                    */
/* -------------------------------------------------------------------- */
    int         nDatum = 0;
    double      adfDatumParm[8];
    int         nEllipsoid=0;

    if( nProjection != 0 && CSLCount(papszNextField) > 0 )
    {
        nDatum = atoi(papszNextField[0]);
        papszNextField++;
    }

    if( (nDatum == 999 || nDatum == 9999)
        && CSLCount(papszNextField) >= 4 )
    {
        nEllipsoid = atoi(papszNextField[0]);
        adfDatumParm[0] = atof(papszNextField[1]);
        adfDatumParm[1] = atof(papszNextField[2]);
        adfDatumParm[2] = atof(papszNextField[3]);
        papszNextField += 4;
    }

    if( nDatum == 9999
        && CSLCount(papszNextField) >= 5 )
    {
        adfDatumParm[3] = atof(papszNextField[0]);
        adfDatumParm[4] = atof(papszNextField[1]);
        adfDatumParm[5] = atof(papszNextField[2]);
        adfDatumParm[6] = atof(papszNextField[3]);
        adfDatumParm[7] = atof(papszNextField[4]);
        papszNextField += 5;
    }

/* -------------------------------------------------------------------- */
/*      Fetch the units string.                                         */
/* -------------------------------------------------------------------- */
    const char  *pszMIFUnits = NULL;
    const char  *pszUnitsName = NULL;
    double dfUnitsConv = 1.0;
    
    if( CSLCount(papszNextField) > 0 )
    {
        pszMIFUnits = papszNextField[0];
        papszNextField++;
    }

    if( nProjection == 1 || pszMIFUnits == NULL )
        /* do nothing */;
    else if( EQUAL(pszMIFUnits,"km") )
    {
        pszUnitsName = "Kilometer"; 
        dfUnitsConv = 1000.0;
    }
    else if( EQUAL(pszMIFUnits, "in" ) )
    {
        pszUnitsName = "IINCH"; 
        dfUnitsConv = 0.0254; 
    }
    else if( EQUAL(pszMIFUnits, "ft" ) )
    {
        pszUnitsName = SRS_UL_FOOT;
        dfUnitsConv = atof(SRS_UL_FOOT_CONV);
    }
    else if( EQUAL(pszMIFUnits, "yd" ) )
    {
        pszUnitsName = "IYARD";
        dfUnitsConv = 0.9144;
    }
    else if( EQUAL(pszMIFUnits, "mm" ) )
    {
        pszUnitsName = "Millimeter";
        dfUnitsConv = 0.001;
    }
    else if( EQUAL(pszMIFUnits, "cm" ) )
    {
        pszUnitsName = "Centimeter";
        dfUnitsConv = 0.01;
    }
    else if( EQUAL(pszMIFUnits, "m" ) )
    {
        pszUnitsName = SRS_UL_METER;
        dfUnitsConv = 1.0;
    }   
    else if( EQUAL(pszMIFUnits, "survey foot" )
             || EQUAL(pszMIFUnits, "survey ft" ) )
    {
        pszUnitsName = SRS_UL_US_FOOT;
        dfUnitsConv = atof(SRS_UL_US_FOOT_CONV);
    }   
    else if( EQUAL(pszMIFUnits, "nmi" ) )
    {
        pszUnitsName = SRS_UL_NAUTICAL_MILE;
        dfUnitsConv = atof(SRS_UL_NAUTICAL_MILE_CONV);
    }   
    else if( EQUAL(pszMIFUnits, "li" ) )
    {
        pszUnitsName = SRS_UL_LINK;
        dfUnitsConv = atof(SRS_UL_LINK_CONV);
    }
    else if( EQUAL(pszMIFUnits, "ch" ) )
    {
        pszUnitsName = SRS_UL_CHAIN;
        dfUnitsConv = atof(SRS_UL_CHAIN_CONV);
    }   
    else if( EQUAL(pszMIFUnits, "rd" ) )
    {
        pszUnitsName = SRS_UL_ROD;
        dfUnitsConv = atof(SRS_UL_ROD);
    }   
    else if( EQUAL(pszMIFUnits, "mi" ) )
    {
        pszUnitsName = "Mile";
        dfUnitsConv = 1609.344;
    }

/* -------------------------------------------------------------------- */
/*      Handle the PROJCS style projections, but add the datum later.   */
/* -------------------------------------------------------------------- */
    switch( nProjection )
    {
        /*--------------------------------------------------------------
         * NonEarth ... we return with an empty SpatialRef.  Eventually
         * we might want to include the units, but not for now.
         *
         * __TODO__ Changed to return NULL because returning an empty
         * SpatialRef caused confusion between Latlon and NonEarth since
         * empty SpatialRefs do have a GEOGCS set and makes them look like
         * Lat/Lon SpatialRefs.
         *
         * Ideally we would like to return a SpatialRef whith no GEGOCS
         *-------------------------------------------------------------*/
      case 0:
        poSR->SetLocalCS( "Nonearth" );