g.proj: consider EPSG code whenever possible

git-svn-id: https://svn.osgeo.org/grass/grass/trunk@72113 15284696-431f-4ddb-bdfa-cd5b030d7da7

general/g.proj/create.c

@@ -6,11 +6,12 @@
 
 #include "local_proto.h"
 
-void create_location(const char *location, const char *epsg)
+void create_location(const char *location)
 {
     int ret;
 
-    ret = G_make_location(location, &cellhd, projinfo, projunits);
+    ret = G_make_location_epsg(location, &cellhd, projinfo, projunits,
+                               projepsg);
     if (ret == 0)
 	G_message(_("Location <%s> created"), location);
     else if (ret == -1)
@@ -23,10 +24,6 @@ void create_location(const char *location, const char *epsg)
 	/* Shouldn't happen */
       G_fatal_error(_("Unable to create location <%s>"), location);
 
-    /* create also PROJ_EPSG */
-    if (epsg)
-        create_epsg(location, epsg);
-        
     G_message(_("You can switch to the new location by\n`%s=%s`"),
 	      "g.mapset mapset=PERMANENT location", location);
 }
@@ -46,7 +43,7 @@ void modify_projinfo()
     
     char path[GPATH_MAX];
 	
-    /* Write out the PROJ_INFO, and PROJ_UNITS if available. */
+    /* Write out the PROJ_INFO, PROJ_UNITS, and PROJ_EPSG if available. */
     if (projinfo != NULL) {
 	G_file_name(path, "", "PROJ_INFO", "PERMANENT");
 	G_write_key_value_file(path, projinfo);
@@ -57,6 +54,11 @@ void modify_projinfo()
 	G_write_key_value_file(path, projunits);
     }
     
+    if (projepsg != NULL) {
+	G_file_name(path, "", "PROJ_EPSG", "PERMANENT");
+	G_write_key_value_file(path, projepsg);
+    }
+
     if ((old_cellhd.zone != cellhd.zone) ||
 	(old_cellhd.proj != cellhd.proj)) {
 	/* Recreate the default, and current window files if projection
@@ -69,29 +71,3 @@ void modify_projinfo()
     }
     G_important_message(_("Projection information updated"));
 }
-
-void create_epsg(const char *location, const char *epsg)
-{
-    FILE *fp;
-    char path[GPATH_MAX];
-    
-    /* if inputs were not clean it should of failed by now */
-    if (location) {
-        G_snprintf(path, sizeof(path), "%s%c%s%c%s%c%s", G_gisdbase(), HOST_DIRSEP, 
-                   location, HOST_DIRSEP,
-                   "PERMANENT", HOST_DIRSEP, "PROJ_EPSG");
-        path[sizeof(path)-1] = '\0';
-    }
-    else {
-        G_file_name(path, "", "PROJ_EPSG", "PERMANENT");
-    }
-    
-    fp = fopen(path, "w");
-    if (!fp)
-        G_fatal_error(_("Unable to create PROJ_EPSG file: %s"), strerror (errno));
-    
-#ifdef HAVE_OGR
-    fprintf(fp, "epsg: %s\n", epsg);
-#endif
-    fclose(fp);
-}

general/g.proj/input.c

@@ -78,6 +78,7 @@ int input_wkt(char *wktfile)
 {
     FILE *infd;
     char buff[8000];
+    OGRSpatialReferenceH hSRS;
     int ret;
 
     if (strcmp(wktfile, "-") == 0)
@@ -100,6 +101,28 @@ int input_wkt(char *wktfile)
     ret = GPJ_wkt_to_grass(&cellhd, &projinfo, &projunits, buff, 0);
     set_default_region();
 
+    hSRS = OSRNewSpatialReference(buff);
+    if (hSRS) {
+	const char *authkey, *authname, *authcode;
+
+	authkey = NULL;
+	if (OSRIsProjected(hSRS))
+	    authkey = "PROJCS";
+	else if (OSRIsGeographic(hSRS))
+	    authkey = "GEOGCS";
+
+	if (authkey) {
+	    authname = OSRGetAuthorityName(hSRS, authkey);
+	    if (authname && *authname && strcmp(authname, "EPSG") == 0) {
+		authcode = OSRGetAuthorityCode(hSRS, authkey);
+		if (authcode && *authcode) {
+		    projepsg = G_create_key_value();
+		    G_set_key_value("epsg", authcode, projepsg);
+		}
+	    }
+	}
+    }
+
     return ret;
 }
 
@@ -171,6 +194,7 @@ int input_proj4(char *proj4params)
 int input_epsg(int epsg_num)
 {
     OGRSpatialReferenceH hSRS;
+    char epsgstr[100];
     int ret = 0;
 
     /* Set finder function for locating OGR csv co-ordinate system tables */
@@ -182,6 +206,10 @@ int input_epsg(int epsg_num)
 
     ret = GPJ_osr_to_grass(&cellhd, &projinfo, &projunits, hSRS, 0);
 
+    sprintf(epsgstr, "%d", epsg_num);
+    projepsg = G_create_key_value();
+    G_set_key_value("epsg", epsgstr, projepsg);
+
     OSRDestroySpatialReference(hSRS);
 
     set_default_region();
@@ -211,6 +239,7 @@ int input_epsg(int epsg_num)
 int input_georef(char *geofile)
 {
     OGRDataSourceH ogr_ds;
+    OGRSpatialReferenceH hSRS;
     int ret = 0;
 
     /* Try opening file with OGR first because it doesn't output a
@@ -218,16 +247,16 @@ int input_georef(char *geofile)
     G_message(_("Trying to open with OGR..."));
     OGRRegisterAll();
 
+    hSRS = NULL;
     if ((ogr_ds = OGROpen(geofile, FALSE, NULL))
 	&& (OGR_DS_GetLayerCount(ogr_ds) > 0)) {
 	OGRLayerH ogr_layer;
-	OGRSpatialReferenceH ogr_srs;
 
 	G_message(_("...succeeded."));
 	/* Get the first layer */
 	ogr_layer = OGR_DS_GetLayer(ogr_ds, 0);
-	ogr_srs = OGR_L_GetSpatialRef(ogr_layer);
-	ret = GPJ_osr_to_grass(&cellhd, &projinfo, &projunits, ogr_srs, 0);
+	hSRS = OGR_L_GetSpatialRef(ogr_layer);
+	ret = GPJ_osr_to_grass(&cellhd, &projinfo, &projunits, hSRS, 0);
 	set_default_region();
 
 	OGR_DS_Destroy(ogr_ds);
@@ -249,6 +278,7 @@ int input_georef(char *geofile)
 				 0);
 
 	    set_gdal_region(gdal_ds);
+	    hSRS = OSRNewSpatialReference(wktstring);
 	}
 	else
 	    G_fatal_error(_("Could not read georeferenced file %s using "
@@ -260,6 +290,27 @@ int input_georef(char *geofile)
 		    "projection information. 'XY (unprojected)' will be used"),
 		  geofile);
 
+    if (hSRS) {
+	const char *authkey, *authname, *authcode;
+
+	authkey = NULL;
+	if (OSRIsProjected(hSRS))
+	    authkey = "PROJCS";
+	else if (OSRIsGeographic(hSRS))
+	    authkey = "GEOGCS";
+
+	if (authkey) {
+	    authname = OSRGetAuthorityName(hSRS, authkey);
+	    if (authname && *authname && strcmp(authname, "EPSG") == 0) {
+		authcode = OSRGetAuthorityCode(hSRS, authkey);
+		if (authcode && *authcode) {
+		    projepsg = G_create_key_value();
+		    G_set_key_value("epsg", authcode, projepsg);
+		}
+	    }
+	}
+    }
+
     return ret;
 
 }

general/g.proj/local_proto.h

@@ -13,7 +13,7 @@ int input_georef(char *);
 #endif
 
 /* output.c */
-void print_projinfo(int, const char *);
+void print_projinfo(int);
 void print_datuminfo(void);
 void print_proj4(int);
 #ifdef HAVE_OGR
@@ -25,6 +25,5 @@ int set_datum(char *);
 int set_datumtrans(int, int);
 
 /* create.c */
-void create_location(const char *, const char *);
+void create_location(const char *);
 void modify_projinfo();
-void create_epsg(const char *, const char *);

general/g.proj/main.c

@@ -233,6 +233,7 @@ int main(int argc, char *argv[])
     }
 
     epsg = inepsg->answer;
+    projinfo = projunits = projepsg = NULL;
 
     /* Input */
     /* We can only have one input source, hence if..else construct */
@@ -292,7 +293,7 @@ int main(int argc, char *argv[])
 		      create->key);
 
     if (printinfo->answer || shellinfo->answer)
-        print_projinfo(shellinfo->answer, epsg);
+        print_projinfo(shellinfo->answer);
     else if (datuminfo->answer)
 	print_datuminfo();
     else if (printproj4->answer)
@@ -302,7 +303,7 @@ int main(int argc, char *argv[])
 	print_wkt(esristyle->answer, dontprettify->answer);
 #endif
     else if (location->answer)
-	create_location(location->answer, epsg);
+	create_location(location->answer);
     else if (create->answer)
 	modify_projinfo();
     else
@@ -316,21 +317,13 @@ int main(int argc, char *argv[])
 		      printinfo->key, shellinfo->key, printproj4->key);
 #endif
 
-#ifdef HAVE_OGR
-    if (create->answer && epsg) {
-#else
-    if (create->answer){ 
-#endif
-	/* preserve epsg code for user records only (not used by grass's pj routines) */
-        create_epsg(location->answer, epsg);
-    }
-
-
     /* Tidy Up */
     if (projinfo != NULL)
 	G_free_key_value(projinfo);
     if (projunits != NULL)
 	G_free_key_value(projunits);
+    if (projepsg != NULL)
+	G_free_key_value(projepsg);
 
     exit(EXIT_SUCCESS);
 

general/g.proj/output.c

@@ -28,10 +28,10 @@
 
 static int check_xy(int shell);
 
-void print_projinfo(int shell, const char *force_epsg)
+/* print projection information gathered from one of the possible inputs */
+void print_projinfo(int shell)
 {
     int i;
-    char path[GPATH_MAX];
 
     if (check_xy(shell))
 	return;
@@ -46,22 +46,12 @@ void print_projinfo(int shell, const char *force_epsg)
 	    fprintf(stdout, "%-11s: %s\n", projinfo->key[i], projinfo->value[i]);
     }
 
-    /* EPSG code is preserved for historical metadata interest only:
-	the contents of this file are not used by pj_*() routines at all */
-    G_file_name(path, "", "PROJ_EPSG", "PERMANENT");
-    if (access(path, F_OK) == 0) {
+    if (projepsg) {
         const char *epsg_value, *epsg_key;
-        struct Key_Value *in_epsg_key;
-        
-        if (force_epsg) {
-            epsg_key = "epsg";
-            epsg_value = force_epsg;
-        }
-        else {
-            in_epsg_key = G_read_key_value_file(path);
-            epsg_key = in_epsg_key->key[0];
-            epsg_value = in_epsg_key->value[0];
-        }
+
+	epsg_key = projepsg->key[0];
+	epsg_value = projepsg->value[0];
+
 	if (!shell) {
 	    fprintf(stdout,
 		"-PROJ_EPSG-------------------------------------------------\n");
@@ -69,21 +59,20 @@ void print_projinfo(int shell, const char *force_epsg)
 	}
 	else
 	    fprintf(stdout, "%s=%s\n", epsg_key, epsg_value);
-
-        if (!force_epsg)
-            G_free_key_value(in_epsg_key);
     }
  
-    if (!shell)
-	fprintf(stdout,
-		"-PROJ_UNITS------------------------------------------------\n");
-    for (i = 0; i < projunits->nitems; i++) {
-	if (shell)
-	    fprintf(stdout, "%s=%s\n",
-		    projunits->key[i], projunits->value[i]);
-	else
-	    fprintf(stdout, "%-11s: %s\n",
-		    projunits->key[i], projunits->value[i]);
+    if (projunits) {
+	if (!shell)
+	    fprintf(stdout,
+		    "-PROJ_UNITS------------------------------------------------\n");
+	for (i = 0; i < projunits->nitems; i++) {
+	    if (shell)
+		fprintf(stdout, "%s=%s\n",
+			projunits->key[i], projunits->value[i]);
+	    else
+		fprintf(stdout, "%-11s: %s\n",
+			projunits->key[i], projunits->value[i]);
+	}
     }
     
     return;