El Miércoles, 27 de Mayo de 2009, escribió:
> Jorge:
> > I send a the proposed files.
>
> Hi,
>
> please please file patches as unidiff attachments to a new
> Enhancement ticket in the trac system, diff'd versus SVN trunk.
> otherwise the overhead to deal with everyone's patches gets
> quickly out of hand.
I send the two unidiff file for the ticket
> http://trac.osgeo.org/grass/report
>
>
> thanks,
> Hamish
--
Dr. E. Jorge Tizado
Dpt. Biodiversity and Environmental Management
University of León
--- v.surf.idw/description.html 2009-05-29 16:38:32.000000000 +0200
+++ description.html 2009-05-29 16:37:32.000000000 +0200
@@ -22,11 +22,11 @@
<p>
The amount of memory used by this program is related to the number
-of vector points in the current region. If the vector point map is
+of vector points in the current region. If the vector point map is
very dense (i.e., contains many data points), the program may
not be able to get all the memory it needs from the
-system. The time required to execute is related to the
-resolution of the current region, after an initial delay
+system. The time required to execute is related to the
+resolution of the current region, after an initial delay
determined by the time taken to read the input vector points map.</p>
<p>
@@ -46,9 +46,9 @@
achieve a similar result.</p>
<p>
-If more than <EM>count</EM> points fall into one target raster cell,
+If more than <EM>count</EM> points fall into one target raster cell,
the mean of all the site values will determine the cell value (unless
-the -n flag is specifed, in which case only the <EM>count</EM>
+the -n flag is specifed, in which case only the <EM>count</EM>
points closest to the centre of the cell will be interpolated).</p>
<P>
@@ -56,11 +56,15 @@
Greater values assign greater influence to values closer to the
point to be interpolated. The interpolation function peaks sharply over
the given data points for 0 < <EM>p</EM> < 1 and more smoothly for
-larger values. The default value for the power parameter is 2.
-</p>
+larger values.<br>
+With default mode (Inverse Squared Distance: weight = 1.0/[dist^2^power] );
+the recommended value for the power parameter is 2.<br>
+With flag -g (Gaussian Distance: weight = exp[-dist^2/power^2] );
+no recommended value for the power parameter (it is a distance).
+</P>
<P>
-By setting <EM>npoints=1</EM>, the module can be used
+By setting <EM>npoints=1</EM>, the module can be used
to calculate raster Voronoi diagrams (Thiessen polygons).</p>
@@ -78,8 +82,8 @@
<H2>AUTHOR</H2>
-Michael Shapiro,
-U.S. Army Construction Engineering
+Michael Shapiro,
+U.S. Army Construction Engineering
Research Laboratory
<BR>
Improved algorithm (indexes points according to cell and ignores
--- v.surf.idw/main.c 2009-05-29 16:38:32.000000000 +0200
+++ main.c 2009-05-29 16:42:11.000000000 +0200
@@ -7,7 +7,7 @@
* points outside current region) by Paul Kelly
* further: Radim Blazek <radim.blazek gmail.com>, Huidae Cho <grass4u gmail.com>,
* Glynn Clements <glynn gclements.plus.com>, Markus Neteler <neteler itc.it>
- * PURPOSE:
+ * PURPOSE:
* COPYRIGHT: (C) 2003-2006 by the GRASS Development Team
*
* This program is free software under the GNU General Public
@@ -64,12 +64,13 @@
double dist;
double sum1, sum2, interp_value;
int n, field;
- double p;
+ double d, p, bw;
struct
{
struct Option *input, *npoints, *power, *output, *dfield, *col;
} parm;
- struct cell_list
+ struct Flag *gauss;
+ struct cell_list
{
int row, column;
struct cell_list *next;
@@ -82,9 +83,8 @@
module = G_define_module();
module->keywords = _("vector, interpolation");
- module->description =
- _("Surface interpolation from vector point data by Inverse "
- "Distance Squared Weighting.");
+ module->description = _("Surface interpolation from vector point data by Inverse Distance "
+ "Squared Weighting and Gaussian Radial Basis Functions");
parm.input = G_define_standard_option(G_OPT_V_INPUT);
@@ -98,18 +98,16 @@
parm.npoints->description = _("Number of interpolation points");
parm.npoints->answer = "12";
- parm.power = G_define_option();
- parm.power->key = "power";
- parm.power->type = TYPE_DOUBLE;
- parm.power->answer = "2.0";
- parm.power->description =
- _("Power parameter; greater values assign greater influence to closer points");
+ parm.power = G_define_option();
+ parm.power->key = "power";
+ parm.power->type = TYPE_DOUBLE;
+ parm.power->required = YES;
+ parm.power->description = _("Power parameter");
parm.dfield = G_define_standard_option(G_OPT_V_FIELD);
- parm.dfield->description =
- _("If set to 0, z coordinates are used (3D vector only)");
+ parm.dfield->description = _("If set to 0, z coordinates are used (3D vector only)");
parm.dfield->answer = "1";
- parm.dfield->gisprompt = "old_layer,layer,layer_zero";
+ parm.dfield->gisprompt = "old_layer,layer,layer_zero";
parm.col = G_define_option();
parm.col->key = "column";
@@ -118,6 +116,11 @@
parm.col->label = _("Attribute table column with values to interpolate");
parm.col->description = _("Required if layer > 0");
+ gauss = G_define_flag();
+ gauss->key = 'g';
+ gauss->description = _("Gaussian weights function [ exp(-dist^2/power^2) ]");
+ gauss->answer = 0;
+
noindex = G_define_flag();
noindex->key = 'n';
noindex->label = _("Don't index points by raster cell");
@@ -125,6 +128,7 @@
" less memory and includes points from outside region"
" in the interpolation");
+
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
@@ -144,8 +148,10 @@
list =
(struct list_Point *)G_calloc((size_t) search_points,
sizeof(struct list_Point));
-
- p = atof ( parm.power->answer );
+
+ p = atof( parm.power->answer );
+ if (gauss->answer)
+ p *= p;
/* get the window, dimension arrays */
G_get_window(&window);
@@ -383,19 +389,35 @@
/* interpolate */
sum1 = 0.0;
sum2 = 0.0;
- for (n = 0; n < nsearch; n++) {
- if ((dist = list[n].dist)) {
- sum1 += list[n].z / pow ( dist, p );
- sum2 += 1.0 / pow ( dist, p );
- }
- else {
- /* If one site is dead on the centre of the cell, ignore
- * all the other sites and just use this value.
- * (Unlikely when using floating point numbers?) */
- sum1 = list[n].z;
- sum2 = 1.0;
- break;
- }
+ if (gauss->answer)
+ {
+ for (n = 0; n < nsearch; n++)
+ {
+ d = exp( - list[n].dist / p );
+ sum1 += list[n].z * d;
+ sum2 += d;
+ }
+ }
+ else /* default: inverse squared distance */
+ {
+ for (n = 0; n < nsearch; n++)
+ {
+ if ((dist = list[n].dist))
+ {
+ d = pow ( dist, p );
+ sum1 += list[n].z / d;
+ sum2 += 1.0 / d;
+ }
+ else
+ {
+ /* If one site is dead on the centre of the cell, ignore
+ * all the other sites and just use this value.
+ * (Unlikely when using floating point numbers?) */
+ sum1 = list[n].z;
+ sum2 = 1.0;
+ break;
+ }
+ }
}
interp_value = sum1 / sum2;
}
_______________________________________________
grass-dev mailing list
grass-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-dev
