GRASS Programmer's Manual  6.4.2(2012)
convert.c
Go to the documentation of this file.
00001 
00016 #include <grass/config.h>
00017 
00018 #ifdef HAVE_OGR
00019 #include <stdio.h>
00020 #include <stdlib.h>
00021 #include <string.h>
00022 #include <math.h>
00023 #include <grass/gis.h>
00024 #include <grass/gprojects.h>
00025 #include <grass/glocale.h>
00026 #include <cpl_csv.h>
00027 #include "local_proto.h"
00028 
00029 /* GRASS relative location of OGR co-ordinate system lookup tables */
00030 #define CSVDIR "/etc/ogr_csv"
00031 
00032 static void DatumNameMassage(char **);
00033 
00052 char *GPJ_grass_to_wkt(struct Key_Value *proj_info,
00053                        struct Key_Value *proj_units,
00054                        int esri_style, int prettify)
00055 {
00056     OGRSpatialReferenceH hSRS;
00057     char *wkt, *local_wkt;
00058 
00059     hSRS = GPJ_grass_to_osr(proj_info, proj_units);
00060 
00061     if (hSRS == NULL)
00062         return NULL;
00063 
00064     if (esri_style)
00065         OSRMorphToESRI(hSRS);
00066 
00067     if (prettify)
00068         OSRExportToPrettyWkt(hSRS, &wkt, 0);
00069     else
00070         OSRExportToWkt(hSRS, &wkt);
00071 
00072     local_wkt = G_store(wkt);
00073     CPLFree(wkt);
00074     OSRDestroySpatialReference(hSRS);
00075     return local_wkt;
00076 }
00077 
00088 OGRSpatialReferenceH GPJ_grass_to_osr(struct Key_Value * proj_info,
00089                                       struct Key_Value * proj_units)
00090 {
00091     struct pj_info pjinfo;
00092     char *proj4, *proj4mod, *wkt, *modwkt, *startmod, *lastpart;
00093     OGRSpatialReferenceH hSRS, hSRS2;
00094     OGRErr errcode;
00095     struct gpj_datum dstruct;
00096     struct gpj_ellps estruct;
00097     size_t len;
00098     char *ellps, *ellpslong, *datum, *params, *towgs84, *datumlongname,
00099         *start, *end, *unit, *unfact, *buff;
00100     const char *sysname, *osrunit, *osrunfact;
00101     double a, es, rf;
00102     int haveparams = 0;
00103 
00104     if ((proj_info == NULL) || (proj_units == NULL))
00105         return NULL;
00106 
00107     hSRS = OSRNewSpatialReference(NULL);
00108 
00109     if (pj_get_kv(&pjinfo, proj_info, proj_units) < 0) {
00110         G_warning(_("Unable parse GRASS PROJ_INFO file"));
00111         return NULL;
00112     }
00113 
00114     if ((proj4 = pj_get_def(pjinfo.pj, 0)) == NULL) {
00115         G_warning(_("Unable get PROJ.4-style parameter string"));
00116         return NULL;
00117     }
00118 
00119     unit = G_find_key_value("unit", proj_units);
00120     unfact = G_find_key_value("meters", proj_units);
00121     if (unfact != NULL && (strcmp(pjinfo.proj, "ll") != 0))
00122         G_asprintf(&proj4mod, "%s +to_meter=%s", proj4, unfact);
00123     else
00124         proj4mod = proj4;
00125 
00126     if ((errcode = OSRImportFromProj4(hSRS, proj4mod)) != OGRERR_NONE) {
00127         G_warning(_("OGR can't parse PROJ.4-style parameter string: "
00128                     "%s (OGR Error code was %d)"), proj4mod, errcode);
00129         return NULL;
00130     }
00131 
00132     if ((errcode = OSRExportToWkt(hSRS, &wkt)) != OGRERR_NONE) {
00133         G_warning(_("OGR can't get WKT-style parameter string "
00134                     "(OGR Error code was %d)"), errcode);
00135         return NULL;
00136     }
00137 
00138     ellps = G_find_key_value("ellps", proj_info);
00139     GPJ__get_ellipsoid_params(proj_info, &a, &es, &rf);
00140     haveparams = GPJ__get_datum_params(proj_info, &datum, &params);
00141 
00142     if ((datum == NULL) || (GPJ_get_datum_by_name(datum, &dstruct) < 0)) {
00143         G_asprintf(&datumlongname, "unknown");
00144         if (ellps == NULL)
00145             G_asprintf(&ellps, "unnamed");
00146     }
00147     else {
00148         datumlongname = G_store(dstruct.longname);
00149         if (ellps == NULL)
00150             ellps = G_store(dstruct.ellps);
00151         GPJ_free_datum(&dstruct);
00152     }
00153     if (GPJ_get_ellipsoid_by_name(ellps, &estruct) > 0) {
00154         ellpslong = G_store(estruct.longname);
00155         DatumNameMassage(&ellpslong);
00156         GPJ_free_ellps(&estruct);
00157     }
00158     else
00159         ellpslong = G_store(ellps);
00160 
00161     startmod = G_strstr(wkt, "GEOGCS");
00162     lastpart = G_strstr(wkt, "PRIMEM");
00163     len = strlen(wkt) - strlen(startmod);
00164     wkt[len] = '\0';
00165     if (haveparams == 2) {
00166         /* Only put datum params into the WKT if they were specifically
00167          * specified in PROJ_INFO */
00168         char *paramkey, *paramvalue;
00169 
00170         paramkey = strtok(params, "=");
00171         paramvalue = params + strlen(paramkey) + 1;
00172         if (G_strcasecmp(paramkey, "towgs84") == 0)
00173             G_asprintf(&towgs84, ",TOWGS84[%s]", paramvalue);
00174         else
00175             towgs84 = "";
00176     }
00177     else
00178         towgs84 = "";
00179 
00180     sysname = OSRGetAttrValue(hSRS, "PROJCS", 0);
00181     if (sysname == NULL) {
00182         /* Not a projected co-ordinate system */
00183         start = "";
00184         end = "";
00185     }
00186     else {
00187         if ((strcmp(sysname, "unnamed") == 0) &&
00188             (G_find_key_value("name", proj_info) != NULL))
00189             G_asprintf(&start, "PROJCS[\"%s\",",
00190                        G_find_key_value("name", proj_info));
00191         else
00192             start = G_store(wkt);
00193 
00194         osrunit = OSRGetAttrValue(hSRS, "UNIT", 0);
00195         osrunfact = OSRGetAttrValue(hSRS, "UNIT", 1);
00196 
00197         if ((unfact == NULL) || (G_strcasecmp(osrunit, "unknown") != 0))
00198             end = "";
00199         else {
00200             double unfactf = atof(unfact);
00201 
00202             G_asprintf(&buff, ",UNIT[\"%s\",", osrunit);
00203 
00204             startmod = G_strstr(lastpart, buff);
00205             len = strlen(lastpart) - strlen(startmod);
00206             lastpart[len] = '\0';
00207 
00208             if (unit == NULL)
00209                 G_asprintf(&unit, "unknown");
00210             G_asprintf(&end, ",UNIT[\"%s\",%.16g]]", unit, unfactf);
00211         }
00212 
00213     }
00214 
00215     G_asprintf(&modwkt,
00216                "%sGEOGCS[\"%s\",DATUM[\"%s\",SPHEROID[\"%s\",%.16g,%.16g]%s],%s%s",
00217                start, ellps, datumlongname, ellpslong, a, rf, towgs84,
00218                lastpart, end);
00219 
00220     hSRS2 = OSRNewSpatialReference(modwkt);
00221 
00222     OSRDestroySpatialReference(hSRS);
00223     G_free(modwkt);
00224     CPLFree(wkt);
00225     pj_dalloc(proj4);
00226     if (proj4 != proj4mod)
00227         G_free(proj4mod);
00228     G_free(datum);
00229     G_free(params);
00230     G_free(datumlongname);
00231     pj_free(pjinfo.pj);
00232     G_free(ellpslong);
00233     /* Other string pointers may or may not need to be freed depending
00234      * on sequence of execution so just leave them. */
00235 
00236     return hSRS2;
00237 }
00238 
00258 int GPJ_osr_to_grass(struct Cell_head *cellhd, struct Key_Value **projinfo,
00259                      struct Key_Value **projunits, OGRSpatialReferenceH hSRS,
00260                      int datumtrans)
00261 {
00262     struct Key_Value *temp_projinfo;
00263     char *pszProj4 = NULL, *pszRemaining;
00264     char *pszProj = NULL;
00265     char *datum = NULL;
00266     struct gpj_datum dstruct;
00267 
00268     if (hSRS == NULL)
00269         goto default_to_xy;
00270 
00271     /* Set finder function for locating OGR csv co-ordinate system tables */
00272     SetCSVFilenameHook(GPJ_set_csv_loc);
00273 
00274     /* Hopefully this doesn't do any harm if it wasn't in ESRI format
00275      * to start with... */
00276     OSRMorphFromESRI(hSRS);
00277 
00278     /* -------------------------------------------------------------------- */
00279     /*      Set cellhd for well known coordinate systems.                   */
00280     /* -------------------------------------------------------------------- */
00281     if (!OSRIsGeographic(hSRS) && !OSRIsProjected(hSRS))
00282         goto default_to_xy;
00283 
00284     if (cellhd) {
00285         int bNorth;
00286 
00287         if (OSRIsGeographic(hSRS)) {
00288             cellhd->proj = PROJECTION_LL;
00289             cellhd->zone = 0;
00290         }
00291         else if (OSRGetUTMZone(hSRS, &bNorth) != 0) {
00292             cellhd->proj = PROJECTION_UTM;
00293             cellhd->zone = OSRGetUTMZone(hSRS, &bNorth);
00294             if (!bNorth)
00295                 cellhd->zone *= -1;
00296         }
00297         else {
00298             cellhd->proj = PROJECTION_OTHER;
00299             cellhd->zone = 0;
00300         }
00301     }
00302 
00303     /* -------------------------------------------------------------------- */
00304     /*      Get the coordinate system definition in PROJ.4 format.          */
00305     /* -------------------------------------------------------------------- */
00306     if (OSRExportToProj4(hSRS, &pszProj4) != OGRERR_NONE)
00307         goto default_to_xy;
00308 
00309     /* -------------------------------------------------------------------- */
00310     /*      Parse the PROJ.4 string into key/value pairs.  Do a bit of      */
00311     /*      extra work to "GRASSify" the result.                            */
00312     /* -------------------------------------------------------------------- */
00313     temp_projinfo = G_create_key_value();
00314 
00315     /* Create "local" copy of proj4 string so we can modify and free it
00316      * using GRASS functions */
00317     pszRemaining = G_store(pszProj4);
00318     CPLFree(pszProj4);
00319     pszProj4 = pszRemaining;
00320     while ((pszRemaining = strstr(pszRemaining, "+")) != NULL) {
00321         char *pszToken, *pszValue;
00322 
00323         pszRemaining++;
00324 
00325         /* Advance pszRemaining to end of this token[=value] pair */
00326         pszToken = pszRemaining;
00327         while (*pszRemaining != ' ' && *pszRemaining != '\0')
00328             pszRemaining++;
00329 
00330         if (*pszRemaining == ' ') {
00331             *pszRemaining = '\0';
00332             pszRemaining++;
00333         }
00334 
00335         /* parse token, and value */
00336         if (strstr(pszToken, "=") != NULL) {
00337             pszValue = strstr(pszToken, "=");
00338             *pszValue = '\0';
00339             pszValue++;
00340         }
00341         else
00342             pszValue = "defined";
00343 
00344 
00345         if (G_strcasecmp(pszToken, "proj") == 0) {
00346             /* The ll projection is known as longlat in PROJ.4 */
00347             if (G_strcasecmp(pszValue, "longlat") == 0)
00348                 pszValue = "ll";
00349 
00350             pszProj = pszValue;
00351             continue;
00352         }
00353 
00354         /* Ellipsoid and datum handled separately below */
00355         if (G_strcasecmp(pszToken, "ellps") == 0
00356             || G_strcasecmp(pszToken, "a") == 0
00357             || G_strcasecmp(pszToken, "b") == 0
00358             || G_strcasecmp(pszToken, "es") == 0
00359             || G_strcasecmp(pszToken, "rf") == 0
00360             || G_strcasecmp(pszToken, "datum") == 0)
00361             continue;
00362 
00363         /* We will handle units separately */
00364         if (G_strcasecmp(pszToken, "to_meter") == 0
00365             || G_strcasecmp(pszToken, "units") == 0)
00366             continue;
00367 
00368         G_set_key_value(pszToken, pszValue, temp_projinfo);
00369     }
00370 
00371     *projinfo = G_create_key_value();
00372 
00373     /* -------------------------------------------------------------------- */
00374     /*      Derive the user name for the projection.                        */
00375     /* -------------------------------------------------------------------- */
00376     if (pszProj) {
00377         char path[4095];
00378         char name[80];
00379 
00380         sprintf(path, "%s/etc/projections", G_gisbase());
00381         if (G_lookup_key_value_from_file(path, pszProj, name, sizeof(name)) >
00382             0)
00383             G_set_key_value("name", name, *projinfo);
00384         else
00385             G_set_key_value("name", pszProj, *projinfo);
00386 
00387         G_set_key_value("proj", pszProj, *projinfo);
00388     }
00389     else
00390         G_warning(_("No projection name! Projection parameters likely to be meaningless."));
00391 
00392 
00393     /* -------------------------------------------------------------------- */
00394     /*      Find the GRASS datum name and choose parameters either          */
00395     /*      interactively or not.                                           */
00396     /* -------------------------------------------------------------------- */
00397 
00398     {
00399         const char *pszDatumNameConst = OSRGetAttrValue(hSRS, "DATUM", 0);
00400         struct datum_list *list, *listhead;
00401         char *dum1, *dum2, *pszDatumName;
00402         int paramspresent =
00403             GPJ__get_datum_params(temp_projinfo, &dum1, &dum2);
00404 
00405         if (pszDatumNameConst) {
00406             /* Need to make a new copy of the string so we don't mess
00407              * around with the memory inside the OGRSpatialReferenceH? */
00408 
00409             pszDatumName = G_store(pszDatumNameConst);
00410             DatumNameMassage(&pszDatumName);
00411 
00412             list = listhead = read_datum_table();
00413 
00414             while (list != NULL) {
00415                 if (G_strcasecmp(pszDatumName, list->longname) == 0) {
00416                     datum = G_store(list->name);
00417                     break;
00418                 }
00419                 list = list->next;
00420             }
00421             free_datum_list(listhead);
00422 
00423             if (datum == NULL) {
00424                 if (paramspresent < 2)
00425                     /* Only give warning if no parameters present */
00426                     G_warning(_("Datum <%s> not recognised by GRASS and no parameters found"),
00427                               pszDatumName);
00428             }
00429             else {
00430                 G_set_key_value("datum", datum, *projinfo);
00431 
00432                 if (paramspresent < 2) {
00433                     /* If no datum parameters were imported from the OSR
00434                      * object then we should use the set specified by datumtrans */
00435                     char *params, *chosenparams = NULL;
00436                     int paramsets;
00437 
00438                     paramsets =
00439                         GPJ_get_default_datum_params_by_name(datum, &params);
00440 
00441                     if (paramsets < 0)
00442                         G_warning(_("Datum <%s> apparently recognised by GRASS but no parameters found. "
00443                                    "You may want to look into this."), datum);
00444                     else if (datumtrans > paramsets) {
00445 
00446                         G_warning(_("Invalid transformation number %d; valid range is 1 to %d. "
00447                                    "Leaving datum transform parameters unspecified."),
00448                                   datumtrans, paramsets);
00449                         datumtrans = 0;
00450                     }
00451 
00452                     if (paramsets > 0) {
00453                         struct gpj_datum_transform_list *list, *old;
00454 
00455                         list = GPJ_get_datum_transform_by_name(datum);
00456 
00457                         if (list != NULL) {
00458                             do {
00459                                 if (list->count == datumtrans) {
00460                                     chosenparams = G_store(list->params);
00461                                     break;
00462                                 }
00463                                 old = list;
00464                                 list = list->next;
00465                                 G_free(old);
00466                             } while (list != NULL);
00467                         }
00468                     }
00469 
00470                     if (chosenparams != NULL) {
00471                         char *paramkey, *paramvalue;
00472 
00473                         paramkey = strtok(chosenparams, "=");
00474                         paramvalue = chosenparams + strlen(paramkey) + 1;
00475                         G_set_key_value(paramkey, paramvalue, *projinfo);
00476                         G_free(chosenparams);
00477                     }
00478 
00479                     if (paramsets > 0)
00480                         G_free(params);
00481                 }
00482 
00483             }
00484         }
00485     }
00486 
00487     /* -------------------------------------------------------------------- */
00488     /*   Determine an appropriate GRASS ellipsoid name if possible, or      */
00489     /*   else just put a and es values into PROJ_INFO                       */
00490     /* -------------------------------------------------------------------- */
00491 
00492     if ((datum != NULL) && (GPJ_get_datum_by_name(datum, &dstruct) > 0)) {
00493         /* Use ellps name associated with datum */
00494         G_set_key_value("ellps", dstruct.ellps, *projinfo);
00495         GPJ_free_datum(&dstruct);
00496         G_free(datum);
00497     }
00498     else {
00499         /* If we can't determine the ellipsoid from the datum, derive it
00500          * directly from "SPHEROID" parameters in WKT */
00501         const char *pszSemiMajor = OSRGetAttrValue(hSRS, "SPHEROID", 1);
00502         const char *pszInvFlat = OSRGetAttrValue(hSRS, "SPHEROID", 2);
00503 
00504         if (pszSemiMajor != NULL && pszInvFlat != NULL) {
00505             char *ellps = NULL;
00506             struct ellps_list *list, *listhead;
00507             double a = atof(pszSemiMajor), invflat = atof(pszInvFlat), flat;
00508             double es;
00509 
00510             /* Allow for incorrect WKT describing a sphere where InvFlat 
00511              * is given as 0 rather than inf */
00512             if (invflat > 0)
00513                 flat = 1 / invflat;
00514             else
00515                 flat = 0;
00516 
00517             es = flat * (2.0 - flat);
00518 
00519             list = listhead = read_ellipsoid_table(0);
00520 
00521             while (list != NULL) {
00522                 /* Try and match a and es against GRASS defined ellipsoids;
00523                  * accept first one that matches. These numbers were found
00524                  * by trial and error and could be fine-tuned, or possibly
00525                  * a direct comparison of IEEE floating point values used. */
00526                 if ((a == list->a || fabs(a - list->a) < 0.1 || fabs(1 - a / list->a) < 0.0000001) && ((es == 0 && list->es == 0) ||    /* Special case for sphere */
00527                                                                                                        (invflat == list->rf || fabs(invflat - list->rf) < 0.0000001))) {
00528                     ellps = G_store(list->name);
00529                     break;
00530                 }
00531                 list = list->next;
00532             }
00533             if (listhead != NULL)
00534                 free_ellps_list(listhead);
00535 
00536             if (ellps == NULL) {
00537                 /* If we weren't able to find a matching ellps name, set
00538                  * a and es values directly from WKT-derived data */
00539                 char es_str[100];
00540 
00541                 G_set_key_value("a", (char *)pszSemiMajor, *projinfo);
00542 
00543                 sprintf(es_str, "%.16g", es);
00544                 G_set_key_value("es", es_str, *projinfo);
00545             }
00546             else {
00547                 /* else specify the GRASS ellps name for readability */
00548                 G_set_key_value("ellps", ellps, *projinfo);
00549                 G_free(ellps);
00550             }
00551 
00552         }
00553 
00554     }
00555 
00556     /* -------------------------------------------------------------------- */
00557     /*      Finally append the detailed projection parameters to the end    */
00558     /* -------------------------------------------------------------------- */
00559 
00560     {
00561         int i;
00562 
00563         for (i = 0; i < temp_projinfo->nitems; i++)
00564             G_set_key_value(temp_projinfo->key[i],
00565                             temp_projinfo->value[i], *projinfo);
00566 
00567         G_free_key_value(temp_projinfo);
00568     }
00569 
00570     G_free(pszProj4);
00571 
00572     /* -------------------------------------------------------------------- */
00573     /*      Set the linear units.                                           */
00574     /* -------------------------------------------------------------------- */
00575     *projunits = G_create_key_value();
00576 
00577     if (OSRIsGeographic(hSRS)) {
00578         /* We assume degrees ... someday we will be wrong! */
00579         G_set_key_value("unit", "degree", *projunits);
00580         G_set_key_value("units", "degrees", *projunits);
00581         G_set_key_value("meters", "1.0", *projunits);
00582     }
00583     else {
00584         char szFormatBuf[256];
00585         char *pszUnitsName = NULL;
00586         double dfToMeters;
00587         char *pszUnitsPlural, *pszStringEnd;
00588 
00589         dfToMeters = OSRGetLinearUnits(hSRS, &pszUnitsName);
00590 
00591         /* Workaround for the most obvious case when unit name is unknown */
00592         if ((G_strcasecmp(pszUnitsName, "unknown") == 0) &&
00593             (dfToMeters == 1.))
00594             G_asprintf(&pszUnitsName, "meter");
00595 
00596         G_set_key_value("unit", pszUnitsName, *projunits);
00597 
00598         /* Attempt at plural formation (WKT format doesn't store plural
00599          * form of unit name) */
00600         pszUnitsPlural = G_malloc(strlen(pszUnitsName) + 3);
00601         strcpy(pszUnitsPlural, pszUnitsName);
00602         pszStringEnd = pszUnitsPlural + strlen(pszUnitsPlural) - 4;
00603         if (G_strcasecmp(pszStringEnd, "foot") == 0) {
00604             /* Special case for foot - change two o's to e's */
00605             pszStringEnd[1] = 'e';
00606             pszStringEnd[2] = 'e';
00607         }
00608         else if (G_strcasecmp(pszStringEnd, "inch") == 0) {
00609             /* Special case for inch - add es */
00610             pszStringEnd[4] = 'e';
00611             pszStringEnd[5] = 's';
00612             pszStringEnd[6] = '\0';
00613         }
00614         else {
00615             /* For anything else add an s at the end */
00616             pszStringEnd[4] = 's';
00617             pszStringEnd[5] = '\0';
00618         }
00619 
00620         G_set_key_value("units", pszUnitsPlural, *projunits);
00621         G_free(pszUnitsPlural);
00622 
00623         sprintf(szFormatBuf, "%.16g", dfToMeters);
00624         G_set_key_value("meters", szFormatBuf, *projunits);
00625 
00626     }
00627 
00628     return 2;
00629 
00630     /* -------------------------------------------------------------------- */
00631     /*      Fallback to returning an ungeoreferenced definition.            */
00632     /* -------------------------------------------------------------------- */
00633   default_to_xy:
00634     if (cellhd != NULL) {
00635         cellhd->proj = PROJECTION_XY;
00636         cellhd->zone = 0;
00637     }
00638 
00639     *projinfo = NULL;
00640     *projunits = NULL;
00641 
00642     return 1;
00643 }
00644 
00645 
00665 int GPJ_wkt_to_grass(struct Cell_head *cellhd, struct Key_Value **projinfo,
00666                      struct Key_Value **projunits, const char *wkt,
00667                      int datumtrans)
00668 {
00669     int retval;
00670 
00671     if (wkt == NULL)
00672         retval =
00673             GPJ_osr_to_grass(cellhd, projinfo, projunits, NULL, datumtrans);
00674     else {
00675         OGRSpatialReferenceH hSRS;
00676 
00677         /* Set finder function for locating OGR csv co-ordinate system tables */
00678         SetCSVFilenameHook(GPJ_set_csv_loc);
00679 
00680         hSRS = OSRNewSpatialReference(wkt);
00681         retval =
00682             GPJ_osr_to_grass(cellhd, projinfo, projunits, hSRS, datumtrans);
00683         OSRDestroySpatialReference(hSRS);
00684     }
00685 
00686     return retval;
00687 }
00688 
00689 
00690 /* GPJ_set_csv_loc()
00691  * 'finder function' for use with OGR SetCSVFilenameHook() function */
00692 
00693 const char *GPJ_set_csv_loc(const char *name)
00694 {
00695     const char *gisbase = G_gisbase();
00696     static char *buf = NULL;
00697 
00698     if (buf != NULL)
00699         G_free(buf);
00700 
00701     G_asprintf(&buf, "%s%s/%s", gisbase, CSVDIR, name);
00702 
00703     return buf;
00704 }
00705 
00706 
00707 /* The list below is only for files that use a non-standard name for a 
00708  * datum that is already supported in GRASS. The number of entries must be even;
00709  * they are all in pairs. The first one in the pair is the non-standard name;
00710  * the second is the GRASS name. If a name appears more than once (as for
00711  * European_Terrestrial_Reference_System_1989) then it means there was more
00712  * than one non-standard name for it that needs to be accounted for. 
00713  *
00714  * N.B. The order of these pairs is different from that in 
00715  * ogr/ogrfromepsg.cpp in the GDAL source tree! GRASS uses the EPSG
00716  * names in its WKT representation except WGS_1984 and WGS_1972 as
00717  * these shortened versions seem to be standard
00718  */
00719 
00720 static const char *papszDatumEquiv[] = {
00721     "Militar_Geographische_Institute",
00722     "Militar_Geographische_Institut",
00723     "World_Geodetic_System_1984",
00724     "WGS_1984",
00725     "World_Geodetic_System_1972",
00726     "WGS_1972",
00727     "European_Terrestrial_Reference_System_89",
00728     "European_Terrestrial_Reference_System_1989",
00729     "European_Reference_System_1989",
00730     "European_Terrestrial_Reference_System_1989",
00731     "ETRS_1989",
00732     "European_Terrestrial_Reference_System_1989",
00733     "ETRS89",
00734     "European_Terrestrial_Reference_System_1989",
00735     "ETRF_1989",
00736     "European_Terrestrial_Reference_System_1989",
00737     "NZGD_2000",
00738     "New_Zealand_Geodetic_Datum_2000",
00739     "Monte_Mario_Rome",
00740     "Monte_Mario",
00741     "MONTROME",
00742     "Monte_Mario",
00743     "Campo_Inchauspe_1969",
00744     "Campo_Inchauspe",
00745     "S_JTSK_Ferro",
00746     "Militar_Geographische_Institut",
00747     "Potsdam_Datum_83",
00748     "Deutsches_Hauptdreiecksnetz",
00749     "South_American_1969",
00750     "South_American_Datum_1969",
00751     NULL
00752 };
00753 
00754 /************************************************************************/
00755 /*                      OGREPSGDatumNameMassage()                       */
00756 /*                                                                      */
00757 /*      Massage an EPSG datum name into WMT format.  Also transform     */
00758 /*      specific exception cases into WKT versions.                     */
00759 
00760 /************************************************************************/
00761 
00762 static void DatumNameMassage(char **ppszDatum)
00763 {
00764     int i, j;
00765     char *pszDatum = *ppszDatum;
00766 
00767     /* -------------------------------------------------------------------- */
00768     /*      Translate non-alphanumeric values to underscores.               */
00769     /* -------------------------------------------------------------------- */
00770     for (i = 0; pszDatum[i] != '\0'; i++) {
00771         if (!(pszDatum[i] >= 'A' && pszDatum[i] <= 'Z')
00772             && !(pszDatum[i] >= 'a' && pszDatum[i] <= 'z')
00773             && !(pszDatum[i] >= '0' && pszDatum[i] <= '9')) {
00774             pszDatum[i] = '_';
00775         }
00776     }
00777 
00778     /* -------------------------------------------------------------------- */
00779     /*      Remove repeated and trailing underscores.                       */
00780     /* -------------------------------------------------------------------- */
00781     for (i = 1, j = 0; pszDatum[i] != '\0'; i++) {
00782         if (pszDatum[j] == '_' && pszDatum[i] == '_')
00783             continue;
00784 
00785         pszDatum[++j] = pszDatum[i];
00786     }
00787     if (pszDatum[j] == '_')
00788         pszDatum[j] = '\0';
00789     else
00790         pszDatum[j + 1] = '\0';
00791 
00792     /* -------------------------------------------------------------------- */
00793     /*      Search for datum equivalences.  Specific massaged names get     */
00794     /*      mapped to OpenGIS specified names.                              */
00795     /* -------------------------------------------------------------------- */
00796     for (i = 0; papszDatumEquiv[i] != NULL; i += 2) {
00797         if (EQUAL(*ppszDatum, papszDatumEquiv[i])) {
00798             G_free(*ppszDatum);
00799             *ppszDatum = G_store(papszDatumEquiv[i + 1]);
00800             break;
00801         }
00802     }
00803 }
00804 
00805 #endif /* HAVE_OGR */
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines