[GRASS-SVN] r53297 - grass/trunk/general/g.proj
svn_grass at osgeo.org
svn_grass at osgeo.org
Mon Oct 1 14:37:51 PDT 2012
Author: pkelly
Date: 2012-10-01 14:37:50 -0700 (Mon, 01 Oct 2012)
New Revision: 53297
Modified:
grass/trunk/general/g.proj/datumtrans.c
grass/trunk/general/g.proj/g.proj.html
grass/trunk/general/g.proj/local_proto.h
grass/trunk/general/g.proj/main.c
Log:
Add new datum= option allowing the datum in the input co-ordinate
system to be overridden (or added if missing). Should be used very
sparingly, as it doesn't check that the correct ellipsoid is used
which could potentially result in ellipsoid/datum mis-matches.
Modified: grass/trunk/general/g.proj/datumtrans.c
===================================================================
--- grass/trunk/general/g.proj/datumtrans.c 2012-10-01 20:53:01 UTC (rev 53296)
+++ grass/trunk/general/g.proj/datumtrans.c 2012-10-01 21:37:50 UTC (rev 53297)
@@ -18,6 +18,7 @@
#include <stdlib.h>
#include <string.h>
#include <grass/gis.h>
+#include <grass/glocale.h>
#include <grass/gprojects.h>
#include "local_proto.h"
@@ -27,6 +28,9 @@
* \brief Add or replace datum transformation parameters to the current
* co-ordinate system definition
*
+ * \param datum Use this datum (overrides any datum found in the
+ * current co-ordinate system definition).
+ *
* \param datumtrans Index number of parameter set to use, 0 to leave
* unspecifed (or remove specific parameters, leaving just
* the datum name), -1 to list the available parameter
@@ -38,18 +42,25 @@
* \return 1 if a change was made, 0 if not.
**/
-int set_datumtrans(int datumtrans, int force)
+int set_datumtrans(char *datum, int datumtrans, int force)
{
- char *params, *datum = NULL;
int paramsets, status;
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)
+ /* datum name provided; set status to 1 to indicate no parameters available */
+ status = 1;
+ else {
+ /* extract datum name and parameter status from input co-ordinate system */
+ char *params;
+ 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 */
@@ -136,7 +147,8 @@
/* Copy old PROJ_INFO, skipping out any keys related
* to datum parameters */
for (i = 0; i < projinfo->nitems; i++) {
- if (strcmp(projinfo->key[i], "dx") == 0
+ 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
@@ -148,7 +160,9 @@
temp_projinfo);
}
- /* Finally add new parameters (if we have them) */
+ /* Finally (re-)add datum name and add new parameters (if we have them) */
+ G_set_key_value("datum", datum, temp_projinfo);
+ G_message(_("Datum set to <%s>"), datum);
if (chosenparams != NULL) {
/* Now split 'chosenparams' into key/value format */
paramkey = strtok(chosenparams, "=");
Modified: grass/trunk/general/g.proj/g.proj.html
===================================================================
--- grass/trunk/general/g.proj/g.proj.html 2012-10-01 20:53:01 UTC (rev 53296)
+++ grass/trunk/general/g.proj/g.proj.html 2012-10-01 21:37:50 UTC (rev 53297)
@@ -7,8 +7,8 @@
<ul>
<li>Reporting the projection information for the current location,
either in conventional GRASS (-p flag) or PROJ.4 (-j flag) format</li>
-<li>Reporting and modifying the datum transformation parameters for
-the current location</li>
+<li>Changing the datum, or reporting and modifying the datum transformation
+parameters, for the current location</li>
</ul>
<p>When compiled with OGR, functionality is increased and allows output of
@@ -33,7 +33,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>
@@ -47,6 +47,18 @@
to support future revisions of the EPSG database.</dd>
</dl>
+<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. Note that the -t flag is
+implicit when a value for <em>datum</em> is specified. 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
@@ -167,6 +179,14 @@
g.proj -c wkt=irish_grid.prj location=irish_grid
</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>
Modified: grass/trunk/general/g.proj/local_proto.h
===================================================================
--- grass/trunk/general/g.proj/local_proto.h 2012-10-01 20:53:01 UTC (rev 53296)
+++ grass/trunk/general/g.proj/local_proto.h 2012-10-01 21:37:50 UTC (rev 53297)
@@ -21,7 +21,7 @@
#endif
/* datumtrans.c */
-int set_datumtrans(int, int);
+int set_datumtrans(char *, int, int);
/* create.c */
void create_location(char *);
Modified: grass/trunk/general/g.proj/main.c
===================================================================
--- grass/trunk/general/g.proj/main.c 2012-10-01 20:53:01 UTC (rev 53296)
+++ grass/trunk/general/g.proj/main.c 2012-10-01 21:37:50 UTC (rev 53297)
@@ -48,6 +48,7 @@
*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;
@@ -156,6 +157,17 @@
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;
@@ -204,6 +216,20 @@
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 */
@@ -233,7 +259,7 @@
G_fatal_error(_("Projection files missing"));
/* Set Datum Parameters if necessary or requested */
- set_datumtrans(atoi(dtrans->answer), forcedatumtrans->answer);
+ set_datumtrans(datum->answer, atoi(dtrans->answer), forcedatumtrans->answer);
/* Output */
More information about the grass-commit
mailing list