[GRASS-SVN] r53453 - in grass/branches/develbranch_6: general/g.proj include lib/proj
svn_grass at osgeo.org
svn_grass at osgeo.org
Wed Oct 17 13:20:25 PDT 2012
Author: neteler
Date: 2012-10-17 13:20:24 -0700 (Wed, 17 Oct 2012)
New Revision: 53453
Added:
grass/branches/develbranch_6/general/g.proj/create.c
Modified:
grass/branches/develbranch_6/general/g.proj/datumtrans.c
grass/branches/develbranch_6/general/g.proj/description.html
grass/branches/develbranch_6/general/g.proj/input.c
grass/branches/develbranch_6/general/g.proj/local_proto.h
grass/branches/develbranch_6/general/g.proj/main.c
grass/branches/develbranch_6/general/g.proj/output.c
grass/branches/develbranch_6/include/gprojects.h
grass/branches/develbranch_6/lib/proj/datum.c
Log:
backport with GRASS 6 adaptations of datum handling patches by Paul Kelly: r53297, r53305, r53310 from trunk
Added: grass/branches/develbranch_6/general/g.proj/create.c
===================================================================
--- grass/branches/develbranch_6/general/g.proj/create.c (rev 0)
+++ grass/branches/develbranch_6/general/g.proj/create.c 2012-10-17 20:20:24 UTC (rev 53453)
@@ -0,0 +1,68 @@
+#include <errno.h>
+#include <string.h>
+
+#include <grass/gis.h>
+#include <grass/glocale.h>
+
+#include "local_proto.h"
+
+void create_location(char *location)
+{
+ int ret;
+
+ ret = G__make_location(location, &cellhd, projinfo, projunits, NULL);
+ if (ret == 0)
+ G_message(_("Location <%s> created"), location);
+ else if (ret == -1)
+ G_fatal_error(_("Unable to create location <%s>: %s"),
+ location, strerror(errno));
+ else if (ret == -2)
+ G_fatal_error(_("Unable to create projection files: %s"),
+ strerror(errno));
+ else
+ /* Shouldn't happen */
+ G_fatal_error(_("Unspecified error while creating new location"));
+
+ G_message(_("You can switch to the new location by\n`%s=%s`"),
+ "g.mapset mapset=PERMANENT location", location);
+}
+
+void modify_projinfo()
+{
+ const char *mapset = G_mapset();
+ struct Cell_head old_cellhd;
+
+ if (strcmp(mapset, "PERMANENT") != 0)
+ G_fatal_error(_("You must select the PERMANENT mapset before updating the "
+ "current location's projection (current mapset is <%s>)."),
+ mapset);
+
+ /* Read projection information from current location first */
+ G_get_default_window(&old_cellhd);
+
+ char path[GPATH_MAX];
+ int stat;
+
+ /* Write out the PROJ_INFO, and PROJ_UNITS if available. */
+ if (projinfo != NULL) {
+ G__file_name(path, "", "PROJ_INFO", "PERMANENT");
+ G_write_key_value_file(path, projinfo, &stat);
+ }
+
+ if (projunits != NULL) {
+ G__file_name(path, "", "PROJ_UNITS", "PERMANENT");
+ G_write_key_value_file(path, projunits, &stat);
+ }
+
+ if ((old_cellhd.zone != cellhd.zone) ||
+ (old_cellhd.proj != cellhd.proj)) {
+ /* Recreate the default, and current window files if projection
+ * number or zone have changed */
+ G__put_window(&cellhd, "", "DEFAULT_WIND");
+ G__put_window(&cellhd, "", "WIND");
+ G_message(_("Default region was updated to the new projection, but if you have "
+ "multiple mapsets `g.region -d` should be run in each to update the "
+ "region from the default"));
+ }
+ G_important_message(_("Projection information updated"));
+}
Modified: grass/branches/develbranch_6/general/g.proj/datumtrans.c
===================================================================
--- grass/branches/develbranch_6/general/g.proj/datumtrans.c 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/general/g.proj/datumtrans.c 2012-10-17 20:20:24 UTC (rev 53453)
@@ -16,13 +16,77 @@
#include <stdio.h>
#include <stdlib.h>
+#include <string.h>
#include <grass/gis.h>
+#include <grass/glocale.h>
#include <grass/gprojects.h>
#include "local_proto.h"
/**
*
+ * \brief Add or replace datum to the current co-ordinate system definition
+ *
+ * \param datum Use this datum (overrides any datum found in the
+ * current co-ordinate system definition).
+ *
+ * \return 1 if a change was made, 0 if not.
+ **/
+
+int set_datum(char *datum)
+{
+ struct gpj_datum dstruct;
+ struct Key_Value *temp_projinfo;
+ int i;
+
+ if (cellhd.proj == PROJECTION_XY)
+ return 0;
+
+ if (!datum || GPJ_get_datum_by_name(datum, &dstruct) < 0)
+ {
+ G_fatal_error(_("Invalid datum code <%s>"), datum);
+ return 0;
+ }
+
+ temp_projinfo = G_create_key_value();
+
+ /* Copy old PROJ_INFO, skipping out any keys related
+ * to datum or ellipsoid parameters */
+ for (i = 0; i < projinfo->nitems; i++) {
+ if (strcmp(projinfo->key[i], "datum") == 0
+ || strcmp(projinfo->key[i], "dx") == 0
+ || strcmp(projinfo->key[i], "dy") == 0
+ || strcmp(projinfo->key[i], "dz") == 0
+ || strcmp(projinfo->key[i], "datumparams") == 0
+ || strcmp(projinfo->key[i], "nadgrids") == 0
+ || strcmp(projinfo->key[i], "towgs84") == 0
+ || strcmp(projinfo->key[i], "ellps") == 0
+ || strcmp(projinfo->key[i], "a") == 0
+ || strcmp(projinfo->key[i], "b") == 0
+ || strcmp(projinfo->key[i], "es") == 0
+ || strcmp(projinfo->key[i], "f") == 0
+ || strcmp(projinfo->key[i], "rf") == 0)
+ continue;
+
+ G_set_key_value(projinfo->key[i], projinfo->value[i],
+ temp_projinfo);
+ }
+
+ /* Finally add datum and ellipsoid names */
+ G_set_key_value("datum", dstruct.name, temp_projinfo);
+ G_message(_("Datum set to <%s>"), dstruct.name);
+ G_set_key_value("ellps", dstruct.ellps, temp_projinfo);
+ G_message(_("Ellipsoid set to <%s>"), dstruct.ellps);
+
+ /* Destroy original key/value structure and replace with new one */
+ G_free_key_value(projinfo);
+ projinfo = temp_projinfo;
+
+ return 1;
+}
+
+/**
+ *
* \brief Add or replace datum transformation parameters to the current
* co-ordinate system definition
*
@@ -34,27 +98,25 @@
* \param force Force editing of parameters even if current co-ordinate
* system already contains fully specified parameters.
*
- * \param interactive Interactively prompt user for parameters, ignoring
- * value of datumtrans variable.
- *
* \return 1 if a change was made, 0 if not.
**/
-int set_datumtrans(int datumtrans, int force, int interactive)
+int set_datumtrans(int datumtrans, int force)
{
char *params, *datum = NULL;
int paramsets, status;
- struct gpj_datum dstruct;
if (cellhd.proj == PROJECTION_XY)
return 0;
status = GPJ__get_datum_params(projinfo, &datum, ¶ms);
G_debug(3, "set_datumtrans(): GPJ__get_datum_params() status=%d", status);
+ G_free(params);
if (datum) {
/* A datum name is specified; need to determine if
* there are parameters to choose from for this datum */
+ struct gpj_datum dstruct;
if (GPJ_get_datum_by_name(datum, &dstruct) > 0) {
char *defparams;
@@ -62,6 +124,8 @@
paramsets =
GPJ_get_default_datum_params_by_name(dstruct.name,
&defparams);
+ G_free(defparams);
+ GPJ_free_datum(&dstruct);
G_debug(3, "set_datumtrans(): datum transform terms found "
"with %d options", paramsets);
@@ -89,15 +153,9 @@
struct Key_Value *temp_projinfo;
int i;
- /* First of all obtain the new parameters either interactively
- * or through the supplied transform number index */
- if (interactive) {
- /* Print to stderr here as the function call will also */
- fprintf(stderr, "Datum name: %s (%s)", dstruct.name,
- dstruct.longname);
- GPJ_ask_datum_params(datum, &chosenparams);
- }
- else {
+ /* First of all obtain the new parameters
+ * through the supplied transform number index */
+ {
struct gpj_datum_transform_list *list;
if (datumtrans > paramsets)
@@ -113,21 +171,23 @@
if (list != NULL) {
if (datumtrans == -1) {
do {
+ struct gpj_datum_transform_list *old = list;
fprintf(stdout, "---\n%d\nUsed in %s\n%s\n%s\n",
list->count, list->where_used,
list->params, list->comment);
list = list->next;
+ GPJ_free_datum_transform(old);
} while (list != NULL);
exit(EXIT_SUCCESS);
}
else {
do {
- if (list->count == datumtrans) {
+ struct gpj_datum_transform_list *old = list;
+ if (list->count == datumtrans)
chosenparams = G_store(list->params);
- break;
- }
list = list->next;
+ GPJ_free_datum_transform(old);
} while (list != NULL);
}
}
@@ -165,6 +225,7 @@
projinfo = temp_projinfo;
}
+ G_free(datum);
return force;
}
Modified: grass/branches/develbranch_6/general/g.proj/description.html
===================================================================
--- grass/branches/develbranch_6/general/g.proj/description.html 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/general/g.proj/description.html 2012-10-17 20:20:24 UTC (rev 53453)
@@ -1,18 +1,24 @@
<h2>DESCRIPTION</h2>
-<p><em>g.proj</em> provides a means of converting a co-ordinate system
+<em>g.proj</em> provides a means of converting a co-ordinate system
description (i.e. projection information) between various formats.
-It requires <a href="http://www.gdal.org/ogr/">OGR</a> to compile. The basic
-functionality of the module is to report the projection information for the
-current location, either in conventional GRASS (-p flag) or PROJ.4 (-j flag)
-format.
+If compiled without <a href="http://www.gdal.org/ogr/">OGR</a> present, the
+functionality is limited to:
+<ul>
+<li>Reporting the projection information for the current location,
+either in conventional GRASS (-p flag) or PROJ.4 (-j flag) format</li>
+<li>Changing the datum, or reporting and modifying the datum transformation
+parameters, for the current location</li>
+</ul>
-<p>Projection information may also be output in the Well-Known Text (WKT)
-format popularised
+<p>
+When compiled with OGR, functionality is increased and allows output of
+the projection information in the Well-Known Text (WKT) format popularised
by proprietary GIS. In addition, if one of the parameters <em>georef</em>,
<em>wkt</em>, <em>proj4</em> or <em>epsg</em> is specified, rather than the
projection information being read from the current location it is imported
from an external source as follows:
+
<dl>
<dt>georef=<em>filename</em></dt>
<dd><em>g.proj</em> attempts to invoke GDAL and OGR in turn to read a
@@ -28,7 +34,7 @@
<dt>proj4=<em>description</em> or <em>-</em></dt>
<dd><em>description</em> should be a projection description in
-<a href="http://remotesensing.org/proj/">PROJ.4</a> format, enclosed in
+<a href="http://proj.osgeo.org/">PROJ.4</a> format, enclosed in
quotation marks if there are any spaces. If <em>-</em> is given for
<em>description</em>, the PROJ.4 description will be read from stdin rather
than as a directly-supplied command-line parameter.</dd>
@@ -38,36 +44,51 @@
co-ordinate system in the <a href="http://www.epsg.org/CurrentDB.html">EPSG
database</a>. EPSG code support is based upon a local copy of the GDAL CSV
co-ordinate system and datum information files, stored in the directory
-${GISBASE}/etc/ogr_csv. These can be updated if necessary to support future
-revisions of the EPSG database.</dd>
+<tt>$GISBASE/etc/proj/ogr_csv</tt>. These can be updated if necessary
+to support future revisions of the EPSG database.</dd>
</dl>
-<p>The -p, -j, -w, etc. flags are all functional when importing projection
+<p>
+If datum information is incorrect or missing in the input
+co-ordinate system definition (e.g. PROJ.4 descriptions have very limited
+support for specifying datum names), a GRASS datum abbreviation can instead be
+supplied using the <em>datum</em> parameter. This will override any
+datum contained in the input co-ordinate system, and discard
+any datum transformation parameters. Enter datum=<em>list</em> to return a
+list of all the datums supported by GRASS. Since any
+existing datum transformation parameters will have been discarded, the
+<em>datumtrans</em> parameter should in general always be used in
+conjunction with <em>datum</em>.
+
+<p>
+The -p, -j, -w, etc. flags are all functional when importing projection
information from an external source, meaning that <em>g.proj</em> can be
used to convert between representations of the information. It is
<strong>not</strong> required that either the input or output be in GRASS
format.
-<p>In addition however, if the -c flag is specified, <em>g.proj</em> will
+<p>
+In addition however, if the -c flag is specified, <em>g.proj</em> will
create new GRASS projection files (PROJ_INFO, PROJ_UNITS, WIND and
DEFAULT_WIND) based on the imported information. If the <em>location</em>
parameter is specified in addition to -c, then a new location will be created.
Otherwise the projection information files in the current location will be
-overwritten. The program will warn before doing this only if command-line
-interactive mode (<em>-i</em> flag) is selected.
+overwritten. The program will <strong>not</strong> warn before doing this.
-<p>The final mode of operation of g.proj is to report on the datum
+<p>
+The final mode of operation of g.proj is to report on the datum
information and datum transformation parameters associated with the
co-ordinate system. The -d flag will report a human-readable summary of
this.
<h2>NOTES</h2>
-<p>If the input co-ordinate system contains a datum name but no
+<p>
+If the input co-ordinate system contains a datum name but no
transformation parameters, and there is more than one suitable parameter set
available (according to the files datum.table and datumtransform.table in
-${GISBASE}/etc), g.proj will check the value of the <em>datumtrans</em>
-option and act according to the following:<br>
+<tt>$GISBASE/etc/proj</tt>), <em>g.proj</em> will check the value of
+the <em>datumtrans</em> option and act according to the following:<br>
<strong>-1:</strong> List available parameter sets in a GUI-parsable (but also
human-readable) format and exit.<br>
<strong>0 (default):</strong> Continue without specifying parameters - if
@@ -76,24 +97,24 @@
<strong>Any other number less than or equal to the number of parameter sets
available for this datum:</strong> Choose this parameter set and add it to the
co-ordinate system description.<br>
-If the module is being used from the command-line through an interactive
-terminal, the <em>-i</em> flag can be specified to enable interactive
-selection of the parameter set, and the value of <em>datumtrans</em> (if
-specified) is ignored.<br>
If the <em>-t</em> flag is specified, the module will attempt to change the
datum transformation parameters using one of the above two methods
<strong>even if</strong> a valid parameter set is already specified in the
-input co-ordinate system.
+input co-ordinate system. This can be useful to change the datum information
+for an existing location.
-<p>Output is simply based on the input projection information. g.proj does
+<p>
+Output is simply based on the input projection information. g.proj does
<strong>not</strong> attempt to verify that the co-ordinate system thus
described matches an existing system in use in the world. In particular,
this means there are no EPSG Authority codes in the WKT output.
-<p>WKT format shows the false eastings and northings in the projected unit
+<p>
+WKT format shows the false eastings and northings in the projected unit
(e.g. meters, feet) but in PROJ format it should always be given in meters.
-<p>The maximum size of input WKT or PROJ.4 projection descriptions is
+<p>
+The maximum size of input WKT or PROJ.4 projection descriptions is
limited to 8000 bytes.
<h2>EXAMPLES</h2>
@@ -105,6 +126,7 @@
</pre></div>
<p>
+
Create a '.prj' file in ESRI format corresponding to the current location:<br>
<div class="code"><pre>
@@ -154,18 +176,17 @@
<p>
Create a new location with the same co-ordinate system as the current
-location:<br>
+location, but forcing a change to datum transformation parameter set no. 1:<br>
<div class="code"><pre>
-g.proj -c location=newloc
+g.proj -c location=newloc -t datumtrans=1
</pre></div>
<p>
-Interactively change/update the datum transformation parameters for the
-current location:<br>
+List the possible datum transformation parameters for the current location:<br>
<div class="code"><pre>
-g.proj -itc
+g.proj -t datumtrans=-1
</pre></div>
<p>
@@ -177,6 +198,15 @@
</pre></div>
<p>
+Create a new location from a PROJ.4 description, explicitly
+specifying a datum and using the default datum transformation
+parameters:<br>
+
+<div class="code"><pre>
+g.proj -c location=spain proj4="+proj=utm +zone=30 +ellps=intl" datum=eur50 datumtrans=0
+</pre></div>
+
+<p>
Reproject external raster map to current GRASS projection (does not always make sense!)
using the GDAL 'gdalwarp' tool. We recommend to use the ERDAS/Img format and not
to use the ESRI style of WKT:<br>
@@ -217,4 +247,5 @@
Paul Kelly
-<p><i>Last changed: $Date$</i>
+<p>
+<i>Last changed: $Date$</i>
Modified: grass/branches/develbranch_6/general/g.proj/input.c
===================================================================
--- grass/branches/develbranch_6/general/g.proj/input.c 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/general/g.proj/input.c 2012-10-17 20:20:24 UTC (rev 53453)
@@ -22,15 +22,22 @@
#include <grass/gis.h>
#include <grass/gprojects.h>
#include <grass/glocale.h>
-#include <gdal.h>
-#include <ogr_api.h>
-#include <cpl_csv.h>
+#include <grass/config.h>
+#ifdef HAVE_OGR
+# include <gdal.h>
+# include <ogr_api.h>
+# include <cpl_csv.h>
+#endif
+
#include "local_proto.h"
static void set_default_region(void);
+
+#ifdef HAVE_OGR
static void set_gdal_region(GDALDatasetH);
static void set_ogr_region(OGRLayerH);
+#endif
/**
* \brief Read projection and region information from current location
@@ -50,6 +57,8 @@
return;
}
+#ifdef HAVE_OGR
+
/**
* \brief Read projection information in WKT format from stdin or a file
*
@@ -86,7 +95,7 @@
G_squeeze(buff);
}
else
- G_fatal_error(_("Unable to open file [%s] for reading"), wktfile);
+ G_fatal_error(_("Unable to open file '%s' for reading"), wktfile);
ret = GPJ_wkt_to_grass(&cellhd, &projinfo, &projunits, buff, 0);
set_default_region();
@@ -254,6 +263,7 @@
return ret;
}
+#endif /* HAVE_OGR */
/**
* \brief Populates global cellhd with "default" region settings
@@ -283,6 +293,8 @@
return;
}
+#ifdef HAVE_OGR
+
/**
* \brief Populates global cellhd with region settings based on
* georeferencing information in a GDAL dataset
@@ -365,3 +377,4 @@
return;
}
+#endif /* HAVE_OGR */
Modified: grass/branches/develbranch_6/general/g.proj/local_proto.h
===================================================================
--- grass/branches/develbranch_6/general/g.proj/local_proto.h 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/general/g.proj/local_proto.h 2012-10-17 20:20:24 UTC (rev 53453)
@@ -1,25 +1,29 @@
-#ifdef G_PROJ_MAIN
-# define G_PROJ_GLOBAL
-#else
-# define G_PROJ_GLOBAL extern
-#endif
+#include <grass/config.h>
-G_PROJ_GLOBAL struct Key_Value *projinfo, *projunits;
-G_PROJ_GLOBAL struct Cell_head cellhd;
+extern struct Key_Value *projinfo, *projunits;
+extern struct Cell_head cellhd;
/* input.c */
void input_currloc(void);
+#ifdef HAVE_OGR
int input_wkt(char *);
int input_proj4(char *);
int input_epsg(int);
int input_georef(char *);
+#endif
/* output.c */
-void print_projinfo(void);
+void print_projinfo(int);
void print_datuminfo(void);
void print_proj4(int);
+#ifdef HAVE_OGR
void print_wkt(int, int);
-void create_location(char *, int);
+#endif
/* datumtrans.c */
-int set_datumtrans(int, int, int);
+int set_datum(char *);
+int set_datumtrans(int, int);
+
+/* create.c */
+void create_location(char *);
+void modify_projinfo();
Modified: grass/branches/develbranch_6/general/g.proj/main.c
===================================================================
--- grass/branches/develbranch_6/general/g.proj/main.c 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/general/g.proj/main.c 2012-10-17 20:20:24 UTC (rev 53453)
@@ -3,14 +3,15 @@
*
* MODULE: g.proj
* AUTHOR(S): Paul Kelly - paul-grass at stjohnspoint.co.uk
+ * Shell script style by Martin Landa <landa.martin gmail.com>
* PURPOSE: Provides a means of reporting the contents of GRASS
* projection information files and creating
* new projection information files.
- * COPYRIGHT: (C) 2003-2007 by the GRASS Development Team
+ * COPYRIGHT: (C) 2003-2007, 2011 by the GRASS Development Team
*
- * This program is free software under the GNU General Public
- * License (>=v2). Read the file COPYING that comes with GRASS
- * for details.
+ * This program is free software under the GNU General
+ * Public License (>=v2). Read the file COPYING that
+ * comes with GRASS for details.
*
*****************************************************************************/
@@ -18,33 +19,41 @@
#include <string.h>
#include <grass/gis.h>
#include <grass/glocale.h>
+#include <grass/config.h>
-#define G_PROJ_MAIN
#include "local_proto.h"
+struct Key_Value *projinfo, *projunits;
+struct Cell_head cellhd;
+
int main(int argc, char *argv[])
{
struct Flag *printinfo, /* Print contents of PROJ_INFO & PROJ_UNITS */
- *printproj4, /* Print projection in PROJ.4 format */
- *datuminfo, /* Check if datum information is present */
- *create, /* Create new projection files */
- *printwkt, /* Print projection in WKT format */
- *esristyle, /* Use ESRI-style WKT format */
- *dontprettify, /* Print 'flat' output (no linebreaks) */
- *forcedatumtrans, /* Force override of datumtrans parameters */
- *interprompt; /* Enable interactive prompting */
-
+ *shellinfo, /* Print in shell script style */
+ *printproj4, /* Print projection in PROJ.4 format */
+ *datuminfo, /* Check if datum information is present */
+ *create, /* Create new projection files */
+#ifdef HAVE_OGR
+ *printwkt, /* Print projection in WKT format */
+ *esristyle, /* Use ESRI-style WKT format */
+#endif
+ *dontprettify, /* Print 'flat' output (no linebreaks) */
+ *forcedatumtrans; /* Force override of datumtrans parameters */
+
struct Option *location, /* Name of new location to create */
- *inepsg, /* EPSG projection code */
- *dtrans, /* index to datum transform option */
- *inwkt, /* Input file with projection in WKT format */
- *inproj4, /* Projection in PROJ.4 format */
- *ingeo; /* Input geo-referenced file readable by
+#ifdef HAVE_OGR
+ *inepsg, /* EPSG projection code */
+ *inwkt, /* Input file with projection in WKT format */
+ *inproj4, /* Projection in PROJ.4 format */
+ *ingeo, /* Input geo-referenced file readable by
* GDAL or OGR */
+#endif
+ *datum, /* datum to add (or replace existing datum) */
+ *dtrans; /* index to datum transform option */
struct GModule *module;
+
+ int formats;
- int importformats;
-
G_set_program_name(argv[0]);
G_no_gisinit(); /* We don't call G_gisinit() here because it validates the
* mapset, whereas this module may legitmately be used
@@ -52,18 +61,29 @@
module = G_define_module();
module->keywords = _("general, projection");
+#ifdef HAVE_OGR
module->label =
- _("Converts co-ordinate system descriptions (i.e. projection "
- "information) between various formats (including GRASS format).");
+ _("Prints and manipulates GRASS projection information files "
+ "(in various co-ordinate system descriptions).");
module->description =
- _("Can also be used to create GRASS locations.");
+ _("Can also be used to create new GRASS locations.");
+#else
+ module->description =
+ _("Prints and manipulates GRASS projection information files.");
+#endif
printinfo = G_define_flag();
printinfo->key = 'p';
printinfo->guisection = _("Print");
printinfo->description =
- _("Print projection information (in conventional GRASS format)");
+ _("Print projection information in conventional GRASS format");
+ shellinfo = G_define_flag();
+ shellinfo->key = 'g';
+ shellinfo->guisection = _("Print");
+ shellinfo->description =
+ _("Print projection information in shell script style");
+
datuminfo = G_define_flag();
datuminfo->key = 'd';
datuminfo->guisection = _("Print");
@@ -76,6 +96,17 @@
printproj4->description =
_("Print projection information in PROJ.4 format");
+ dontprettify = G_define_flag();
+ dontprettify->key = 'f';
+ dontprettify->guisection = _("Print");
+ dontprettify->description =
+ _("Print 'flat' output with no linebreaks (applies to "
+#ifdef HAVE_OGR
+ "WKT and "
+#endif
+ "PROJ.4 output)");
+
+#ifdef HAVE_OGR
printwkt = G_define_flag();
printwkt->key = 'w';
printwkt->guisection = _("Print");
@@ -86,21 +117,14 @@
esristyle->guisection = _("Print");
esristyle->description =
_("Use ESRI-style format (applies to WKT output only)");
-
- dontprettify = G_define_flag();
- dontprettify->key = 'f';
- dontprettify->guisection = _("Print");
- dontprettify->description =
- _("Print 'flat' output with no linebreaks (applies to WKT and "
- "PROJ.4 output)");
-
+
ingeo = G_define_option();
ingeo->key = "georef";
ingeo->type = TYPE_STRING;
ingeo->key_desc = "file";
ingeo->required = NO;
- ingeo->guisection = _("Input");
- ingeo->description = _("Georeferenced data file to read projection "
+ ingeo->guisection = _("Specification");
+ ingeo->description = _("Name of georeferenced data file to read projection "
"information from");
inwkt = G_define_option();
@@ -108,29 +132,45 @@
inwkt->type = TYPE_STRING;
inwkt->key_desc = "file";
inwkt->required = NO;
- inwkt->guisection = _("Input");
- inwkt->description = _("ASCII file containing a WKT projection "
- "description (- for stdin)");
+ inwkt->guisection = _("Specification");
+ inwkt->label = _("Name of ASCII file containing a WKT projection "
+ "description");
+ inwkt->description = _("'-' for standard input");
inproj4 = G_define_option();
inproj4->key = "proj4";
inproj4->type = TYPE_STRING;
inproj4->key_desc = "params";
inproj4->required = NO;
- inproj4->guisection = _("Input");
- inproj4->description = _("PROJ.4 projection description (- for stdin)");
+ inproj4->guisection = _("Specification");
+ inproj4->label = _("PROJ.4 projection description");
+ inproj4->description = _("'-' for standard input");
inepsg = G_define_option();
inepsg->key = "epsg";
inepsg->type = TYPE_INTEGER;
+ inepsg->key_desc = "code";
inepsg->required = NO;
inepsg->options = "1-1000000";
- inepsg->guisection = _("Input");
+ inepsg->guisection = _("Specification");
inepsg->description = _("EPSG projection code");
+#endif
+ datum = G_define_option();
+ datum->key = "datum";
+ datum->type = TYPE_STRING;
+ datum->key_desc = "name";
+ datum->required = NO;
+ datum->guisection = _("Datum");
+ datum->label =
+ _("Datum (overrides any datum specified in input co-ordinate system)");
+ datum->description =
+ _("Accepts standard GRASS datum codes, or \"list\" to list and exit");
+
dtrans = G_define_option();
dtrans->key = "datumtrans";
dtrans->type = TYPE_INTEGER;
+ dtrans->key_desc = "index";
dtrans->required = NO;
dtrans->options = "-1-100";
dtrans->answer = "0";
@@ -147,48 +187,57 @@
create = G_define_flag();
create->key = 'c';
- create->guisection = _("Create/Edit");
+ create->guisection = _("Modify");
create->description = _("Create new projection files (modifies current "
- "location unless 'location' option specified)");
+ "location)");
location = G_define_option();
location->key = "location";
location->type = TYPE_STRING;
location->key_desc = "name";
location->required = NO;
- location->guisection = _("Create/Edit");
+ location->guisection = _("Create");
location->description = _("Name of new location to create");
- interprompt = G_define_flag();
- interprompt->key = 'i';
- interprompt->description =
- _("Enable interactive prompting (for command-line use only)");
-
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
/* Initialisation & Validation */
+#ifdef HAVE_OGR
/* -e implies -w */
if (esristyle->answer && !printwkt->answer)
printwkt->answer = 1;
- projinfo = NULL;
- projunits = NULL;
-
- importformats = ((ingeo->answer ? 1 : 0) + (inwkt->answer ? 1 : 0) +
- (inproj4->answer ? 1 : 0) + (inepsg->answer ? 1 : 0));
- if (importformats > 1)
+ formats = ((ingeo->answer ? 1 : 0) + (inwkt->answer ? 1 : 0) +
+ (inproj4->answer ? 1 : 0) + (inepsg->answer ? 1 : 0));
+ if (formats > 1)
G_fatal_error(_("Only one of '%s', '%s', '%s' or '%s' options may be specified"),
ingeo->key, inwkt->key, inproj4->key, inepsg->key);
+ /* List supported datums if requested; code originally
+ * from G_ask_datum_name() (formerly in libgis) */
+ if (datum->answer && strcmp(datum->answer, "list") == 0) {
+ const char *dat;
+ int i;
+
+ for (i = 0; (dat = G_datum_name(i)); i++) {
+ fprintf(stdout, "---\n%d\n%s\n%s\n%s ellipsoid\n",
+ i, dat, G_datum_description(i), G_datum_ellipsoid(i));
+ }
+
+ exit(EXIT_SUCCESS);
+ }
+
/* Input */
/* We can only have one input source, hence if..else construct */
- if (importformats == 0)
+ if (formats == 0)
+#endif
/* Input is projection of current location */
input_currloc();
+#ifdef HAVE_OGR
else if (inwkt->answer)
/* Input in WKT format */
input_wkt(inwkt->answer);
@@ -201,41 +250,93 @@
else
/* Input from georeferenced file */
input_georef(ingeo->answer);
+#endif
-
/* Consistency Check */
if ((cellhd.proj != PROJECTION_XY)
&& (projinfo == NULL || projunits == NULL))
G_fatal_error(_("Projection files missing"));
+ /* Override input datum if requested */
+ if(datum->answer)
+ set_datum(datum->answer);
+
/* Set Datum Parameters if necessary or requested */
- set_datumtrans(atoi(dtrans->answer), forcedatumtrans->answer,
- interprompt->answer);
+ set_datumtrans(atoi(dtrans->answer), forcedatumtrans->answer);
/* Output */
- /* We can output the same information in multiple formats if
- * necessary, hence multiple if statements */
+ /* Only allow one output format at a time, to reduce confusion */
+ formats = ((printinfo->answer ? 1 : 0) + (shellinfo->answer ? 1 : 0) +
+ (datuminfo->answer ? 1 : 0) +
+ (printproj4->answer ? 1 : 0) +
+#ifdef HAVE_OGR
+ (printwkt->answer ? 1 : 0) +
+#endif
+ (create->answer ? 1 : 0));
+ if (formats > 1)
+ G_fatal_error(_("Only one of -%c, -%c, -%c, -%c"
+#ifdef HAVE_OGR
+ ", -%c"
+#endif
+ " or -%c flags may be specified"),
+ printinfo->key, shellinfo->key, datuminfo->key, printproj4->key,
+#ifdef HAVE_OGR
+ printwkt->key,
+#endif
+ create->key);
- if (printinfo->answer)
- print_projinfo();
-
- if (datuminfo->answer)
+ if (printinfo->answer || shellinfo->answer)
+ print_projinfo(shellinfo->answer);
+ else if (datuminfo->answer)
print_datuminfo();
-
- if (printproj4->answer)
+ else if (printproj4->answer)
print_proj4(dontprettify->answer);
-
- if (printwkt->answer)
+#ifdef HAVE_OGR
+ else if (printwkt->answer)
print_wkt(esristyle->answer, dontprettify->answer);
+#endif
+ else if (location->answer)
+ create_location(location->answer);
+ else if (create->answer)
+ modify_projinfo();
+ else
+#ifdef HAVE_OGR
+ G_fatal_error(_("No output format specified, define one "
+ "of flags -%c, -%c, -%c, or -%c"),
+ printinfo->key, shellinfo->key, printproj4->key, printwkt->key);
+#else
+ G_fatal_error(_("No output format specified, define one "
+ "of flags -%c, -%c, or -%c"),
+ printinfo->key, shellinfo->key, printproj4->key);
+#endif
- if (create->answer)
- create_location(location->answer, interprompt->answer);
+#ifdef HAVE_OGR
+ if (create->answer && inepsg->answer) {
+#else
+ if (create->answer){
+#endif
+ /* preserve epsg code for user records only (not used by grass's pj routines) */
+ FILE *fp;
+ char path[GPATH_MAX];
+ /* if inputs were not clean it should of failed by now */
+ if (location->answer) {
+ snprintf(path, sizeof(path), "%s/%s/%s/%s", G_gisdbase(),
+ location->answer, "PERMANENT", "PROJ_EPSG");
+ path[sizeof(path)-1] = '\0';
+ }
+ else
+ G__file_name(path, "", "PROJ_EPSG", "PERMANENT");
+ fp = fopen(path, "w");
+#ifdef HAVE_OGR
+ fprintf(fp, "epsg: %s\n", inepsg->answer);
+#endif
+ fclose(fp);
+ }
/* Tidy Up */
-
if (projinfo != NULL)
G_free_key_value(projinfo);
if (projunits != NULL)
Modified: grass/branches/develbranch_6/general/g.proj/output.c
===================================================================
--- grass/branches/develbranch_6/general/g.proj/output.c 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/general/g.proj/output.c 2012-10-17 20:20:24 UTC (rev 53453)
@@ -15,34 +15,70 @@
*****************************************************************************/
#include <stdio.h>
+#include <unistd.h>
#include <string.h>
-#include <errno.h>
+#include <proj_api.h>
+
#include <grass/gis.h>
#include <grass/gprojects.h>
#include <grass/glocale.h>
+#include <grass/config.h>
#include "local_proto.h"
-static int check_xy(void);
+static int check_xy(int shell);
-void print_projinfo(void)
+void print_projinfo(int shell)
{
int i;
+ char path[GPATH_MAX];
- if (check_xy())
+ if (check_xy(shell))
return;
- fprintf(stdout,
- "-PROJ_INFO-------------------------------------------------\n");
- for (i = 0; i < projinfo->nitems; i++)
- fprintf(stdout, "%-11s: %s\n", projinfo->key[i], projinfo->value[i]);
+ if (!shell)
+ fprintf(stdout,
+ "-PROJ_INFO-------------------------------------------------\n");
+ for (i = 0; i < projinfo->nitems; i++) {
+ if (shell)
+ fprintf(stdout, "%s=%s\n", projinfo->key[i], projinfo->value[i]);
+ else
+ fprintf(stdout, "%-11s: %s\n", projinfo->key[i], projinfo->value[i]);
+ }
- fprintf(stdout,
- "-PROJ_UNITS------------------------------------------------\n");
- for (i = 0; i < projunits->nitems; i++)
- fprintf(stdout, "%-11s: %s\n",
- projunits->key[i], projunits->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) {
+ struct Key_Value *in_epsg_key;
+ int stat;
+ in_epsg_key = G_read_key_value_file(path, &stat);
+ if (!shell) {
+ fprintf(stdout,
+ "-PROJ_EPSG-------------------------------------------------\n");
+ fprintf(stdout, "%-11s: %s\n", in_epsg_key->key[0],
+ in_epsg_key->value[0]);
+ }
+ else
+ fprintf(stdout, "%s=%s\n", in_epsg_key->key[0],
+ in_epsg_key->value[0]);
+ if (in_epsg_key != NULL)
+ 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]);
+ }
+
return;
}
@@ -52,7 +88,7 @@
struct gpj_datum dstruct;
int validdatum = 0;
- if (check_xy())
+ if (check_xy(FALSE))
return;
GPJ__get_datum_params(projinfo, &datum, ¶ms);
@@ -93,13 +129,15 @@
void print_proj4(int dontprettify)
{
struct pj_info pjinfo;
- char *proj4, *proj4mod, *i, *unfact;
+ char *proj4, *proj4mod, *i;
+ const char *unfact;
- if (check_xy())
+ if (check_xy(FALSE))
return;
pj_get_kv(&pjinfo, projinfo, projunits);
proj4 = pj_get_def(pjinfo.pj, 0);
+ pj_free(pjinfo.pj);
/* GRASS-style PROJ.4 strings don't include a unit factor as this is
* handled separately in GRASS - must include it here though */
@@ -107,7 +145,8 @@
if (unfact != NULL && (strcmp(pjinfo.proj, "ll") != 0))
G_asprintf(&proj4mod, "%s +to_meter=%s", proj4, unfact);
else
- proj4mod = proj4;
+ proj4mod = G_store(proj4);
+ pj_dalloc(proj4);
for (i = proj4mod; *i; i++) {
/* Don't print the first space */
@@ -120,15 +159,17 @@
fputc(*i, stdout);
}
fputc('\n', stdout);
+ G_free(proj4mod);
return;
}
+#ifdef HAVE_OGR
void print_wkt(int esristyle, int dontprettify)
{
char *outwkt;
- if (check_xy())
+ if (check_xy(FALSE))
return;
outwkt = GPJ_grass_to_wkt(projinfo, projunits, esristyle,
@@ -142,121 +183,13 @@
return;
}
+#endif
-void create_location(char *location, int interactive)
+static int check_xy(int shell)
{
- int ret;
-
- if (location) {
- ret = G__make_location(location, &cellhd, projinfo, projunits, NULL);
- if (ret == 0)
- G_message(_("Location %s created!"), location);
- else if (ret == -1)
- G_fatal_error(_("Unable to create location: %s"),
- strerror(errno));
- else if (ret == -2)
- G_fatal_error(_("Unable to create projection files: %s"),
- strerror(errno));
- else
- /* Shouldn't happen */
- G_fatal_error(_("Unspecified error while creating new location"));
- }
- else {
- /* Create flag given but no location specified; overwrite
- * projection files for current location */
-
- int go_ahead = 0;
- char *mapset = G_mapset();
- struct Key_Value *old_projinfo = NULL, *old_projunits = NULL;
- struct Cell_head old_cellhd;
-
- if (strcmp(mapset, "PERMANENT") != 0)
- G_fatal_error(_("You must select the PERMANENT mapset before updating the "
- "current location's projection. (Current mapset is %s)"),
- mapset);
-
- /* Read projection information from current location first */
- G_get_default_window(&old_cellhd);
-
- if (old_cellhd.proj != PROJECTION_XY) {
- old_projinfo = G_get_projinfo();
- old_projunits = G_get_projunits();
- }
-
- if (old_projinfo && old_projunits) {
- if (interactive) {
- /* Warn as in g.setproj before overwriting current location */
- fprintf(stderr,
- _("\n\nWARNING! A projection file already exists for this location\n"));
- fprintf(stderr,
- "\nThis file contains all the parameters for the\nlocation's projection: %s\n",
- G_find_key_value("proj", old_projinfo));
- fprintf(stderr,
- "\n Overriding this information implies that the old projection parameters\n");
- fprintf(stderr,
- " were incorrect. If you change the parameters, all existing data will be\n");
- fprintf(stderr,
- " interpreted differently by the projection software.\n%c%c%c",
- 7, 7, 7);
- fprintf(stderr,
- " GRASS will not re-project your data automatically\n\n");
-
- if (G_yes
- (_("Would you still like to overwrite the current projection information "),
- 0))
- go_ahead = 1;
- }
- else
- go_ahead = 1;
- }
- else {
- /* Projection files missing for some reason;
- * go ahead and update */
- go_ahead = 1;
- }
-
- if (go_ahead) {
- int out_stat;
- char path[4096];
-
- /* Write out the PROJ_INFO, and PROJ_UNITS if available. */
- if (projinfo != NULL) {
- G__file_name(path, "", "PROJ_INFO", "PERMANENT");
- G_write_key_value_file(path, projinfo, &out_stat);
- if (out_stat != 0)
- G_fatal_error(_("Error writing PROJ_INFO"));
- }
-
- if (projunits != NULL) {
- G__file_name(path, "", "PROJ_UNITS", "PERMANENT");
- G_write_key_value_file(path, projunits, &out_stat);
- if (out_stat != 0)
- G_fatal_error(_("Error writing PROJ_UNITS"));
- }
-
- if ((old_cellhd.zone != cellhd.zone) ||
- (old_cellhd.proj != cellhd.proj)) {
- /* Recreate the default, and current window files if projection
- * number or zone have changed */
- G__put_window(&cellhd, "", "DEFAULT_WIND");
- G__put_window(&cellhd, "", "WIND");
- G_message(_("N.B. The default region was updated to the new projection, but if you have "
- "multiple mapsets g.region -d should be run in each to update the region from "
- "the default."));
- }
- G_message(_("Projection information updated!"));
- }
- else
- G_message(_("The projection information will not be updated."));
-
- }
-
- return;
-}
-
-static int check_xy(void)
-{
if (cellhd.proj == PROJECTION_XY) {
+ if (shell)
+ fprintf(stdout, "name=");
fprintf(stdout, "XY location (unprojected)\n");
return 1;
}
Modified: grass/branches/develbranch_6/include/gprojects.h
===================================================================
--- grass/branches/develbranch_6/include/gprojects.h 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/include/gprojects.h 2012-10-17 20:20:24 UTC (rev 53453)
@@ -95,6 +95,7 @@
int GPJ_get_datum_params(char **, char **);
int GPJ__get_datum_params(struct Key_Value *, char **, char **);
void GPJ_free_datum(struct gpj_datum *);
+void GPJ_free_datum_transform(struct gpj_datum_transform_list *);
int GPJ_ask_datum_params(const char *, char **);
struct gpj_datum_transform_list *GPJ_get_datum_transform_by_name(const char
*inputname);
Modified: grass/branches/develbranch_6/lib/proj/datum.c
===================================================================
--- grass/branches/develbranch_6/lib/proj/datum.c 2012-10-17 20:15:48 UTC (rev 53452)
+++ grass/branches/develbranch_6/lib/proj/datum.c 2012-10-17 20:20:24 UTC (rev 53453)
@@ -85,7 +85,7 @@
int GPJ_get_default_datum_params_by_name(const char *name, char **params)
{
struct gpj_datum_transform_list *list, *old;
- int count = 1;
+ int count = 0;
list = GPJ_get_datum_transform_by_name(name);
@@ -98,14 +98,13 @@
* (will normally be a 3-parameter transformation) */
*params = G_store(list->params);
- while (list->next != NULL) {
+ while (list != NULL) {
count++;
old = list;
list = list->next;
- G_free(old);
+ GPJ_free_datum_transform(old);
}
- G_free(list);
return count;
}
@@ -459,6 +458,21 @@
}
/**
+ * \brief Free the memory used by a gpj_datum_transform_list struct
+ *
+ * \param item gpj_datum_transform_list struct to be freed
+ **/
+
+void GPJ_free_datum_transform(struct gpj_datum_transform_list *item)
+{
+ G_free(item->params);
+ G_free(item->where_used);
+ G_free(item->comment);
+ G_free(item);
+ return;
+}
+
+/**
* \brief Read the current GRASS datum.table from disk and store in
* memory
*
More information about the grass-commit
mailing list