[GRASS-SVN] r40990 - grass/trunk/general/g.transform
svn_grass at osgeo.org
svn_grass at osgeo.org
Sun Feb 14 02:05:41 EST 2010
Author: hamish
Date: 2010-02-14 02:05:40 -0500 (Sun, 14 Feb 2010)
New Revision: 40990
Modified:
grass/trunk/general/g.transform/g.transform.html
grass/trunk/general/g.transform/main.c
Log:
add new options to dump coeffs and transform arbitrary points, both fwd and reverse (trac #929; merge from devbr6)
Modified: grass/trunk/general/g.transform/g.transform.html
===================================================================
--- grass/trunk/general/g.transform/g.transform.html 2010-02-14 07:01:54 UTC (rev 40989)
+++ grass/trunk/general/g.transform/g.transform.html 2010-02-14 07:05:40 UTC (rev 40990)
@@ -3,13 +3,65 @@
<EM>g.transform</EM> is an utility to compute transformation
based upon GCPs and output error measurements.
+
+<H2>NOTES</H2>
+
+For coordnaites given with the <b>coords</b> file option or fed from
+<tt>stdin</tt>, the input format is "x y" with one coordinate pair per
+line.
+
+<p>
+The transformations are:
+<p>
+order=1:
+<div class="code"><pre>
+ e = [E0 E1][1].[1]
+ [E2 0][e] [n]
+
+ n = [N0 N1][1].[1]
+ [N2 0][e] [n]
+</pre></div>
+
+order=2:
+<div class="code"><pre>
+ e = [E0 E1 E3][1 ] [1 ]
+ [E2 E4 0][e ].[n ]
+ [E5 0 0][e²] [n²]
+
+ n = [N0 N1 N3][1 ] [1 ]
+ [N2 N4 0][e ].[n ]
+ [N5 0 0][e²] [n²]
+</pre></div>
+
+order=3:
+<div class="code"><pre>
+ e = [E0 E1 E3 E6][1 ] [1 ]
+ [E2 E4 E7 0][e ].[n ]
+ [E5 E8 0 0][e²] [n²]
+ [E9 0 0 0][e³] [n³]
+
+ n = [N0 N1 N3 N6][1 ] [1 ]
+ [N2 N4 N7 0][e ].[n ]
+ [N5 N8 0 0][e²] [n²]
+ [N9 0 0 0][e³] [n³]
+</pre></div>
+
+["." = dot-product, (AE).N = N'EA.]
+<p>
+In other words, order=1 and order=2 are equivalent to order=3 with
+the higher coefficients equal to zero.
+
+
<H2>SEE ALSO</H2>
<EM><A HREF="i.rectify.html">i.rectify</A></EM>
+
<H2>AUTHORS</H2>
Brian J. Buckley<br>
-Glynn Clements
+Glynn Clements<br>
+Hamish Bowman
-<p><i>Last changed: $Date$</i></p>
+<p>
+<i>Last changed: $Date$</i></p>
Modified: grass/trunk/general/g.transform/main.c
===================================================================
--- grass/trunk/general/g.transform/main.c 2010-02-14 07:01:54 UTC (rev 40989)
+++ grass/trunk/general/g.transform/main.c 2010-02-14 07:05:40 UTC (rev 40990)
@@ -2,11 +2,12 @@
/****************************************************************************
*
* MODULE: g.transform
- * AUTHOR(S): Brian J. Buckley<br>
+ * AUTHOR(S): Brian J. Buckley
* Glynn Clements
+ * Hamish Bowman
* PURPOSE: Utility to compute transformation based upon GCPs and
* output error measurements
- * COPYRIGHT: (C) 2006-2008 by the GRASS Development Team
+ * COPYRIGHT: (C) 2006-2010 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
@@ -44,12 +45,16 @@
static char *name;
static int order;
static int summary;
+static int forward;
static char **columns;
static int need_fwd;
static int need_rev;
static int need_fd;
static int need_rd;
+static char *coord_file;
+static double E12[10], N12[10], E21[10], N21[10];
+
static struct Control_Points points;
static int equation_stat;
@@ -83,7 +88,6 @@
static void compute_transformation(void)
{
static const int order_pnts[3] = { 3, 6, 10 };
- double E12[10], N12[10], E21[10], N21[10];
int n, i;
equation_stat =
@@ -94,7 +98,7 @@
order_pnts[order - 1]);
if (equation_stat <= 0)
- return;
+ G_fatal_error(_("Error conducting transform (%d)"), equation_stat);
count = 0;
@@ -134,7 +138,7 @@
update_stats(&rev, n, rx, ry, rd, rd2);
}
- if (!columns)
+ if (!columns[0])
continue;
for (i = 0;; i++) {
@@ -236,10 +240,74 @@
}
}
+static void dump_cooefs(void)
+{
+ int i;
+ static const int order_pnts[3] = { 3, 6, 10 };
+
+ for (i = 0; i < order_pnts[order - 1]; i++)
+ fprintf(stdout, "E%d=%.15g\n", i, forward ? E12[i] : E21[i]);
+
+ for (i = 0; i < order_pnts[order - 1]; i++)
+ fprintf(stdout, "N%d=%.15g\n", i, forward ? N12[i] : N21[i]);
+}
+
+static void xform_value(double east, double north)
+{
+ double xe, xn;
+
+ if(forward)
+ CRS_georef(east, north, &xe, &xn, E12, N12, order);
+ else
+ CRS_georef(east, north, &xe, &xn, E21, N21, order);
+
+ fprintf(stdout, "%.15g %.15g\n", xe, xn);
+}
+
+static void do_pt_xforms(void)
+{
+ double easting, northing;
+ int ret;
+ FILE *fp;
+
+ if (strcmp(coord_file, "-") == 0)
+ fp = stdin;
+ else {
+ fp = fopen(coord_file, "r");
+ if (!fp)
+ G_fatal_error(_("Unable to open file <%s>"), coord_file);
+ }
+
+ for (;;) {
+ char buf[64];
+
+ if (!G_getl2(buf, sizeof(buf), fp))
+ break;
+
+ if ((buf[0] == '#') || (buf[0] == '\0'))
+ continue;
+
+ /* ? sscanf(buf, "%s %s", &east_str, &north_str)
+ ? G_scan_easting(,,-1)
+ ? G_scan_northing(,,-1) */
+ /* ? muliple delims with sscanf(buf, "%[ ,|\t]", &dummy) ? */
+
+ ret = sscanf(buf, "%lf %lf", &easting, &northing);
+ if (ret != 2)
+ G_fatal_error(_("Invalid coordinates: [%s]"), buf);
+
+ xform_value(easting, northing);
+ }
+
+ if (fp != stdin)
+ fclose(fp);
+}
+
+
int main(int argc, char **argv)
{
- struct Option *grp, *val, *fmt;
- struct Flag *sum;
+ struct Option *grp, *val, *fmt, *xfm_pts;
+ struct Flag *sum, *rev_flag, *dump_flag;
struct GModule *module;
G_gisinit(argv[0]);
@@ -283,13 +351,32 @@
sum->key = 's';
sum->description = _("Display summary information");
+ xfm_pts = G_define_standard_option(G_OPT_F_INPUT);
+ xfm_pts->key = "coords";
+ xfm_pts->required = NO;
+ xfm_pts->label =
+ _("File containing coordinates to transform (\"-\" to read from stdin)");
+ xfm_pts->description = _("Local x,y coordinates to target east,north");
+
+ rev_flag = G_define_flag();
+ rev_flag->key = 'r';
+ rev_flag->label = _("Reverse transform of coords file or coeff. dump");
+ rev_flag->description = _("Target east,north coordinates to local x,y");
+
+ dump_flag = G_define_flag();
+ dump_flag->key = 'x';
+ dump_flag->description = _("Display transform matrix coefficients");
+
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
+
name = grp->answer;
order = atoi(val->answer);
summary = !!sum->answer;
columns = fmt->answers;
+ forward = !rev_flag->answer;
+ coord_file = xfm_pts->answer;
I_get_control_points(name, &points);
@@ -301,5 +388,11 @@
analyze();
+ if(dump_flag->answer)
+ dump_cooefs();
+
+ if(coord_file)
+ do_pt_xforms();
+
return 0;
}
More information about the grass-commit
mailing list