ogr_srs_proj4.cpp

mitab,读取MapInfo的地图文件

/******************************************************************************
 * $Id: ogr_srs_proj4.cpp 12537 2007-10-26 03:16:26Z warmerdam $
 *
 * Project:  OpenGIS Simple Features Reference Implementation
 * Purpose:  OGRSpatialReference interface to PROJ.4.
 * Author:   Frank Warmerdam <[e-mail]>
 *
 ******************************************************************************
 * Copyright (c) 1999,  Les Technologies SoftMap Inc. 
 *
 * 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.
 ****************************************************************************/

#include "ogr_spatialref.h"
#include "ogr_p.h"
#include "cpl_conv.h"

extern int EPSGGetWGS84Transform( int nGeogCS, double *padfTransform );

CPL_CVSID("$Id: ogr_srs_proj4.cpp 12537 2007-10-26 03:16:26Z warmerdam $");

/* -------------------------------------------------------------------- */
/*      The following list comes from osrs/proj/src/pj_ellps.c          */
/*      ... please update from time to time.                            */
/* -------------------------------------------------------------------- */
static const char *ogr_pj_ellps[] = {
"MERIT",        "a=6378137.0", "rf=298.257", "MERIT 1983",
"SGS85",        "a=6378136.0", "rf=298.257",  "Soviet Geodetic System 85",
"GRS80",        "a=6378137.0", "rf=298.257222101", "GRS 1980(IUGG, 1980)",
"IAU76",        "a=6378140.0", "rf=298.257", "IAU 1976",
"airy",         "a=6377563.396", "b=6356256.910", "Airy 1830",
"APL4.9",       "a=6378137.0.",  "rf=298.25", "Appl. Physics. 1965",
"NWL9D",        "a=6378145.0.",  "rf=298.25", "Naval Weapons Lab., 1965",
"mod_airy",     "a=6377340.189", "b=6356034.446", "Modified Airy",
"andrae",       "a=6377104.43",  "rf=300.0",    "Andrae 1876 (Den., Iclnd.)",
"aust_SA",      "a=6378160.0", "rf=298.25", "Australian Natl & S. Amer. 1969",
"GRS67",        "a=6378160.0", "rf=298.2471674270", "GRS 67(IUGG 1967)",
"bessel",       "a=6377397.155", "rf=299.1528128", "Bessel 1841",
"bess_nam",     "a=6377483.865", "rf=299.1528128", "Bessel 1841 (Namibia)",
"clrk66",       "a=6378206.4", "b=6356583.8", "Clarke 1866",
"clrk80",       "a=6378249.145", "rf=293.4663", "Clarke 1880 mod.",
"CPM",          "a=6375738.7", "rf=334.29", "Comm. des Poids et Mesures 1799",
"delmbr",       "a=6376428.",  "rf=311.5", "Delambre 1810 (Belgium)",
"engelis",      "a=6378136.05", "rf=298.2566", "Engelis 1985",
"evrst30",  "a=6377276.345", "rf=300.8017",  "Everest 1830",
"evrst48",  "a=6377304.063", "rf=300.8017",  "Everest 1948",
"evrst56",  "a=6377301.243", "rf=300.8017",  "Everest 1956",
"evrst69",  "a=6377295.664", "rf=300.8017",  "Everest 1969",
"evrstSS",  "a=6377298.556", "rf=300.8017",  "Everest (Sabah & Sarawak)",
"fschr60",  "a=6378166.",   "rf=298.3", "Fischer (Mercury Datum) 1960",
"fschr60m", "a=6378155.",   "rf=298.3", "Modified Fischer 1960",
"fschr68",  "a=6378150.",   "rf=298.3", "Fischer 1968",
"helmert",  "a=6378200.",   "rf=298.3", "Helmert 1906",
"hough",        "a=6378270.0", "rf=297.", "Hough",
"intl",         "a=6378388.0", "rf=297.", "International 1909 (Hayford)",
"krass",        "a=6378245.0", "rf=298.3", "Krassovsky, 1942",
"kaula",        "a=6378163.",  "rf=298.24", "Kaula 1961",
"lerch",        "a=6378139.",  "rf=298.257", "Lerch 1979",
"mprts",        "a=6397300.",  "rf=191.", "Maupertius 1738",
"new_intl",     "a=6378157.5", "b=6356772.2", "New International 1967",
"plessis",      "a=6376523.",  "b=6355863.", "Plessis 1817 (France)",
"SEasia",       "a=6378155.0", "b=6356773.3205", "Southeast Asia",
"walbeck",      "a=6376896.0", "b=6355834.8467", "Walbeck",
"WGS60",    "a=6378165.0",  "rf=298.3", "WGS 60",
"WGS66",        "a=6378145.0", "rf=298.25", "WGS 66",
"WGS72",        "a=6378135.0", "rf=298.26", "WGS 72",
"WGS84",    "a=6378137.0",  "rf=298.257223563", "WGS 84",
"sphere",   "a=6370997.0",  "b=6370997.0", "Normal Sphere (r=6370997)",
0, 0, 0, 0,
};

/************************************************************************/
/*                          OSRProj4Tokenize()                          */
/*                                                                      */
/*      Custom tokenizing function for PROJ.4 strings.  The main        */
/*      reason we can't just use CSLTokenizeString is to handle         */
/*      strings with a + sign in the exponents of parameter values.     */
/************************************************************************/

char **OSRProj4Tokenize( const char *pszFull )

{
    char *pszStart = NULL;
    char *pszFullWrk;
    char **papszTokens = NULL;
    int  i;

    if( pszFull == NULL )
        return NULL;

    pszFullWrk = CPLStrdup( pszFull );

    for( i=0; pszFullWrk[i] != '\0'; i++ )
    {
        switch( pszFullWrk[i] )
        {
          case '+':
            if( i == 0 || pszFullWrk[i-1] == '\0' )
            {
                if( pszStart != NULL )
                {
                    if( strstr(pszStart,"=") != NULL )
                        papszTokens = CSLAddString( papszTokens, pszStart );
                    else
                    {
                        CPLString osAsBoolean = pszStart;
                        osAsBoolean += "=yes";
                        papszTokens = CSLAddString( papszTokens, osAsBoolean );
                    }
                }
                pszStart = pszFullWrk + i + 1;
            }
            break;

          case ' ':
          case '\t':
          case '\n':
            pszFullWrk[i] = '\0';
            break;

          default:
            break;
        }
    }

    if( pszStart != NULL && strlen(pszStart) > 0 )
        papszTokens = CSLAddString( papszTokens, pszStart );

    CPLFree( pszFullWrk );

    return papszTokens;
}


/************************************************************************/
/*                         OSRImportFromProj4()                         */
/************************************************************************/

OGRErr OSRImportFromProj4( OGRSpatialReferenceH hSRS, const char *pszProj4 )

{
    return ((OGRSpatialReference *) hSRS)->importFromProj4( pszProj4 );
}

/************************************************************************/
/*                              OSR_GDV()                               */
/*                                                                      */
/*      Fetch a particular parameter out of the parameter list, or      */
/*      the indicated default if it isn't available.  This is a         */
/*      helper function for importFromProj4().                          */
/************************************************************************/

static double OSR_GDV( char **papszNV, const char * pszField, 
                       double dfDefaultValue )

{
    const char * pszValue;

    pszValue = CSLFetchNameValue( papszNV, pszField );

    // special hack to use k_0 if available.
    if( pszValue == NULL && EQUAL(pszField,"k") )
        pszValue = CSLFetchNameValue( papszNV, "k_0" );

    if( pszValue == NULL )
        return dfDefaultValue;
    else
        return CPLDMSToDec(pszValue);
}

/************************************************************************/
/*                          importFromProj4()                           */
/************************************************************************/

/**
 * Import PROJ.4 coordinate string.
 *
 * The OGRSpatialReference is initialized from the passed PROJ.4 style
 * coordinate system string.  In addition to many +proj formulations which
 * have OGC equivelents, it is also possible to import "+init=epsg:n" style
 * definitions.  These are passed to importFromEPSG().  Other init strings
 * (such as the state plane zones) are not currently supported.   
 *
 * Example:
 *   pszProj4 = "+proj=utm +zone=11 +datum=WGS84" 
 *
 * This method is the equivelent of the C function OSRImportFromProj4().
 *
 * @param pszProj4 the PROJ.4 style string. 
 *
 * @return OGRERR_NONE on success or OGRERR_CORRUPT_DATA on failure.
 */

OGRErr OGRSpatialReference::importFromProj4( const char * pszProj4 )

{
    char **papszNV = NULL;
    char **papszTokens;
    int  i;
    char *pszCleanCopy;

/* -------------------------------------------------------------------- */
/*      Strip any newlines or other "funny" stuff that might occur      */
/*      if this string just came from reading a file.                   */
/* -------------------------------------------------------------------- */
    pszCleanCopy = CPLStrdup( pszProj4 );
    for( i = 0; pszCleanCopy[i] != '\0'; i++ )
    {
        if( pszCleanCopy[i] == 10 
            || pszCleanCopy[i] == 13 
            || pszCleanCopy[i] == 9 )
            pszCleanCopy[i] = ' ';
    }

/* -------------------------------------------------------------------- */
/*      Try to normalize the definition.  This should expand +init=     */
/*      clauses and so forth.                                           */
/* -------------------------------------------------------------------- */
    char *pszNormalized;

    pszNormalized = OCTProj4Normalize( pszCleanCopy );
    CPLFree( pszCleanCopy );
    
/* -------------------------------------------------------------------- */
/*      If we have an EPSG based init string, try to process it         */
/*      directly to get the fullest possible EPSG definition.           */
/* -------------------------------------------------------------------- */
    if( strstr(pszNormalized,"init=epsg:") != NULL )
    {
        OGRErr eErr;
        const char *pszNumber = strstr(pszNormalized,"init=epsg:") + 10;

        eErr = importFromEPSG( atoi(pszNumber) );
        if( eErr == OGRERR_NONE )
        {
            CPLFree( pszNormalized );
            return eErr;
        }
    }

/* -------------------------------------------------------------------- */
/*      Parse the PROJ.4 string into a cpl_string.h style name/value    */
/*      list.                                                           */
/* -------------------------------------------------------------------- */
    papszTokens = OSRProj4Tokenize( pszNormalized );
    CPLFree( pszNormalized );
    
    for( i = 0; papszTokens != NULL && papszTokens[i] != NULL; i++ )
    {
        char *pszEqual = strstr(papszTokens[i],"=");

        if( pszEqual == NULL )
            papszNV = CSLAddNameValue(papszNV, papszTokens[i], "" );
        else
        {
            pszEqual[0] = '\0';
            papszNV = CSLAddNameValue( papszNV, papszTokens[i], pszEqual+1 );
        }
    }

    CSLDestroy( papszTokens );

/* -------------------------------------------------------------------- */
/*      Extract the prime meridian, if there is one set.                */
/* -------------------------------------------------------------------- */
    const char *pszPM = CSLFetchNameValue( papszNV, "pm" );
    double dfFromGreenwich = 0.0;
    int    nPMCode = -1;

    if( pszPM != NULL )
    {
        if( EQUAL(pszPM,"lisbon") )
        {
            dfFromGreenwich = CPLDMSToDec( "9d07'54.862\"W" );
            nPMCode = 8902;
        }
        else if( EQUAL(pszPM,"paris") )
        {
            dfFromGreenwich = CPLDMSToDec( "2d20'14.025\"E" );
            nPMCode = 8903;
        }
        else if( EQUAL(pszPM,"bogota") )
        {
            dfFromGreenwich = CPLDMSToDec( "74d04'51.3\"W" );
            nPMCode = 8904;
        }
        else if( EQUAL(pszPM,"madrid") )
        {
            dfFromGreenwich = CPLDMSToDec( "3d41'16.48\"W" );
            nPMCode = 8905;
        }
        else if( EQUAL(pszPM,"rome") )
        {
            dfFromGreenwich = CPLDMSToDec( "12d27'8.4\"E" );
            nPMCode = 8906;
        }
        else if( EQUAL(pszPM,"bern") )
        {
            dfFromGreenwich = CPLDMSToDec( "7d26'22.5\"E" );
            nPMCode = 8907;
        }
        else if( EQUAL(pszPM,"jakarta") )
        {
            dfFromGreenwich = CPLDMSToDec( "106d48'27.79\"E" );
            nPMCode = 8908;
        }
        else if( EQUAL(pszPM,"ferro") )
        {
            dfFromGreenwich = CPLDMSToDec( "17d40'W" );
            nPMCode = 8909;
        }
        else if( EQUAL(pszPM,"brussels") )
        {
            dfFromGreenwich = CPLDMSToDec( "4d22'4.71\"E" );