[GRASS-SVN] r61520 - in grass/trunk: raster/r.resamp.bspline scripts/r.fillnulls vector/v.surf.bspline
svn_grass at osgeo.org
svn_grass at osgeo.org
Mon Aug 4 15:01:28 PDT 2014
Author: neteler
Date: 2014-08-04 15:01:28 -0700 (Mon, 04 Aug 2014)
New Revision: 61520
Modified:
grass/trunk/raster/r.resamp.bspline/crosscorr.c
grass/trunk/raster/r.resamp.bspline/main.c
grass/trunk/raster/r.resamp.bspline/r.resamp.bspline.html
grass/trunk/scripts/r.fillnulls/r.fillnulls.py
grass/trunk/vector/v.surf.bspline/crosscorr.c
grass/trunk/vector/v.surf.bspline/main.c
grass/trunk/vector/v.surf.bspline/v.surf.bspline.html
Log:
r.resamp.bspline/v.surf.bspline: standardized se/sn and sie/sin parameters to ew_step/ns_step (trac #2299); manual: r.resamp.bspline example added
Modified: grass/trunk/raster/r.resamp.bspline/crosscorr.c
===================================================================
--- grass/trunk/raster/r.resamp.bspline/crosscorr.c 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/raster/r.resamp.bspline/crosscorr.c 2014-08-04 22:01:28 UTC (rev 61520)
@@ -195,7 +195,7 @@
if (nparam_spl > 22900)
G_fatal_error(_("Too many splines (%d x %d). "
- "Consider changing spline steps \"sn=\" \"se=\"."),
+ "Consider changing spline steps \"ew_step=\" \"ns_step=\"."),
nsplx, nsply);
BW = P_get_BandWidth(bilin, nsply);
Modified: grass/trunk/raster/r.resamp.bspline/main.c
===================================================================
--- grass/trunk/raster/r.resamp.bspline/main.c 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/raster/r.resamp.bspline/main.c 2014-08-04 22:01:28 UTC (rev 61520)
@@ -100,7 +100,7 @@
mask_opt->required = NO;
stepE_opt = G_define_option();
- stepE_opt->key = "se";
+ stepE_opt->key = "ew_step";
stepE_opt->type = TYPE_DOUBLE;
stepE_opt->required = NO;
stepE_opt->description =
@@ -108,7 +108,7 @@
stepE_opt->guisection = _("Settings");
stepN_opt = G_define_option();
- stepN_opt->key = "sn";
+ stepN_opt->key = "ns_step";
stepN_opt->type = TYPE_DOUBLE;
stepN_opt->required = NO;
stepN_opt->description =
@@ -183,7 +183,7 @@
if (stepE_opt->answer) {
stepE = atof(stepE_opt->answer);
if (stepE <= .0)
- G_fatal_error(_("se must be positive"));
+ G_fatal_error(_("ew_step must be positive"));
}
else
stepE = src_reg.ew_res * 1.5;
@@ -191,7 +191,7 @@
if (stepN_opt->answer) {
stepN = atof(stepN_opt->answer);
if (stepN <= .0)
- G_fatal_error(_("sn must be positive"));
+ G_fatal_error(_("ns_step must be positive"));
}
else
stepN = src_reg.ns_res * 1.5;
@@ -363,7 +363,7 @@
else {
G_debug(1, "Cross validation finished correctly");
- G_done_msg(_("Cross validation finished for se = %f and sn = %f"), stepE, stepN);
+ G_done_msg(_("Cross validation finished for ew_step = %f and ns_step = %f"), stepE, stepN);
segment_release(&in_seg); /* release memory */
close(in_fd);
Modified: grass/trunk/raster/r.resamp.bspline/r.resamp.bspline.html
===================================================================
--- grass/trunk/raster/r.resamp.bspline/r.resamp.bspline.html 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/raster/r.resamp.bspline/r.resamp.bspline.html 2014-08-04 22:01:28 UTC (rev 61520)
@@ -1,35 +1,39 @@
<h2>DESCRIPTION</h2>
+
<em>r.resamp.bspline</em> performs a bilinear/bicubic spline interpolation with
Tykhonov regularization. The input is a raster surface map, e.g. elevation,
temperature, precipitation etc. Output is a raster map. Optionally, only
input NULL cells are interpolated, useful to fill NULL cells, an alternative
to <em><a href="r.fillnulls.html">r.fillnulls</a></em>. Using the <b>-n</b> flag to only
interpolate NULL cells will considerably speed up the module.
-<p>The input raster map is read at its native resolution, the output raster
+<p>
+The input raster map is read at its native resolution, the output raster
map will be produced for the current computational region set with
<em><a href="g.region.html">g.region</a></em>. Any MASK will be respected, masked
values will be treated as NULL cells in both the input and the output map.
-<p>Spline step values <b>se</b> for the east-west direction and
-<b>sn</b> for the north-south direction should not be smaller than
+<p>Spline step values <b>ew_step</b> for the east-west direction and
+<b>ns_step</b> for the north-south direction should not be smaller than
the east-west and north-south resolutions of the input map. For a raster
map without NULL cells, 1 * resolution can be used, but check for
undershoots and overshoots. For very large areas with missing values
(NULL cells), larger spline step values may be required, but most of the
time the defaults (1.5 x resolution) should be fine.
-<p>The Tykhonov regularization parameter (<b>lambda</b>) acts to
+<p>
+The Tykhonov regularization parameter (<b>lambda</b>) acts to
smooth the interpolation. With a small <b>lambda</b>, the
interpolated surface closely follows observation points; a larger value
will produce a smoother interpolation. Reasonable values are 0.0001,
0.001, 0.005, 0.01, 0.02, 0.05, 0.1 (needs more testing). For seamless
NULL cell interpolation, a small value is required and default is set to 0.005.
-<p>From a theoretical perspective, the interpolating procedure takes place in two
+<p>
+From a theoretical perspective, the interpolating procedure takes place in two
parts: the first is an estimate of the linear coefficients of a spline function;
these are derived from the observation points using a least squares regression; the
second is the computation of the interpolated surface (or interpolated vector
points). As used here, the splines are 2D piece-wise non-zero polynomial
functions calculated within a limited 2D area. The length of each spline step
-is defined by <b>se</b> for the east-west direction and
-<b>sn</b> for the north-south direction. For optimal performance, the
+is defined by <b>ew_step</b> for the east-west direction and
+<b>ns_step</b> for the north-south direction. For optimal performance, the
spline step values should be no less than the east-west and north-south
resolutions of the input map. Each non-NULL cell observation is modeled as a
linear function of the non-zero splines in the area around the observation.
@@ -61,6 +65,7 @@
<h3>Interpolation of NULL cells and patching</h3>
+General procedure:
<div class="code"><pre>
# set region to area with NULL cells, align region to input map
g.region n=north s=south e=east w=west align=input -p
@@ -72,6 +77,36 @@
r.patch input=input_raster,interpolated_nulls output=input_raster_gapfilled
</pre></div>
+<h3>Interpolation of NULL cells and patching (NC data)</h3>
+
+In this example, the SRTM elevation map in the
+North Carolina sample dataset location is filtered for outlier
+elevation values; missing pixels are then re-interpolated to obtain
+a complete elevation map:
+
+<div class="code"><pre>
+g.region rast=elev_srtm_30m -p
+d.mon wx0
+d.histogram elev_srtm_30m
+
+r.univar -e elev_srtm_30m
+
+# remove too low elevations (esp. lakes)
+# Threshold: thresh = Q1 - 1.5 * (Q3 - Q1)
+r.mapcalc "elev_srtm_30m_filt = if(elev_srtm_30m < 50.0, null(), elev_srtm_30m)"
+
+# verify
+d.histogram elev_srtm_30m_filt
+d.erase
+d.rast elev_srtm_30m_filt
+
+r.resamp.bspline -n input=elev_srtm_30m_filt output=elev_srtm_30m_complete \
+ method=bicubic
+
+d.histogram elev_srtm_30m_complete
+d.rast elev_srtm_30m_complete
+</pre></div>
+
<h3>Estimation of <b>lambda</b> parameter with a cross validation proccess</h3>
A random sample of points should be generated first with
Modified: grass/trunk/scripts/r.fillnulls/r.fillnulls.py
===================================================================
--- grass/trunk/scripts/r.fillnulls/r.fillnulls.py 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/scripts/r.fillnulls/r.fillnulls.py 2014-08-04 22:01:28 UTC (rev 61520)
@@ -8,13 +8,13 @@
# Updated to GRASS 6.0 by Markus Neteler
# Ring and zoom improvements by Hamish Bowman
# Converted to Python by Glynn Clements
-# Add support to v.surf.bspline by Luca Delucchi
+# Added support for r.resamp.bspline by Luca Delucchi
# Per hole filling with RST by Maris Nartiss
# Speedup for per hole filling with RST by Stefan Blumentrath
# PURPOSE: fills NULL (no data areas) in raster maps
# The script respects a user mask (MASK) if present.
#
-# COPYRIGHT: (C) 2001-2012 by the GRASS Development Team
+# COPYRIGHT: (C) 2001-2014 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
@@ -349,12 +349,12 @@
if usermask:
grass.run_command('r.resamp.bspline', input = input, mask = usermask,
output = prefix + 'filled', method = method,
- se = 3 * reg['ewres'], sn = 3 * reg['nsres'],
+ ew_step = 3 * reg['ewres'], ns_step = 3 * reg['nsres'],
_lambda = 0.01, flags = 'n')
else:
grass.run_command('r.resamp.bspline', input = input,
output = prefix + 'filled', method = method,
- se = 3 * reg['ewres'], sn = 3 * reg['nsres'],
+ ew_step = 3 * reg['ewres'], ns_step = 3 * reg['nsres'],
_lambda = 0.01, flags = 'n')
# restoring user's mask, if present:
Modified: grass/trunk/vector/v.surf.bspline/crosscorr.c
===================================================================
--- grass/trunk/vector/v.surf.bspline/crosscorr.c 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/vector/v.surf.bspline/crosscorr.c 2014-08-04 22:01:28 UTC (rev 61520)
@@ -149,7 +149,7 @@
if (nparam_spl > 22900)
G_fatal_error(_("Too many splines (%d x %d). "
- "Consider changing spline steps \"sie=\" \"sin=\"."),
+ "Consider changing spline steps \"ew_step=\" \"ns_step=\"."),
nsplx, nsply);
BW = P_get_BandWidth(bilin, nsply);
Modified: grass/trunk/vector/v.surf.bspline/main.c
===================================================================
--- grass/trunk/vector/v.surf.bspline/main.c 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/vector/v.surf.bspline/main.c 2014-08-04 22:01:28 UTC (rev 61520)
@@ -131,7 +131,7 @@
mask_opt->required = NO;
stepE_opt = G_define_option();
- stepE_opt->key = "sie";
+ stepE_opt->key = "ew_step";
stepE_opt->type = TYPE_DOUBLE;
stepE_opt->required = NO;
stepE_opt->answer = "4";
@@ -140,7 +140,7 @@
stepE_opt->guisection = _("Settings");
stepN_opt = G_define_option();
- stepN_opt->key = "sin";
+ stepN_opt->key = "ns_step";
stepN_opt->type = TYPE_DOUBLE;
stepN_opt->required = NO;
stepN_opt->answer = "4";
@@ -301,7 +301,7 @@
Vect_close(&In);
- G_done_msg(_("Cross validation finished for sie = %f and sin = %f"), stepE, stepN);
+ G_done_msg(_("Cross validation finished for ew_step = %f and ns_step = %f"), stepE, stepN);
exit(EXIT_SUCCESS);
}
}
Modified: grass/trunk/vector/v.surf.bspline/v.surf.bspline.html
===================================================================
--- grass/trunk/vector/v.surf.bspline/v.surf.bspline.html 2014-08-04 21:54:45 UTC (rev 61519)
+++ grass/trunk/vector/v.surf.bspline/v.surf.bspline.html 2014-08-04 22:01:28 UTC (rev 61520)
@@ -18,14 +18,14 @@
of the interpolated surface (or interpolated vector points). As used
here, the splines are 2D piece-wise non-zero polynomial functions
calculated within a limited, 2D area. The length of each spline step
-is defined by <b>sie</b> for the east-west direction and
-<b>sin</b> for the north-south direction. For optimum performance, the
+is defined by <b>ew_step</b> for the east-west direction and
+<b>ns_step</b> for the north-south direction. For optimal performance, the
length of spline step should be no less than the distance between observation
points. Each vector point observation is modeled as a linear function of the
non-zero splines in the area around the observation. The least squares
regression predicts the the coefficients of these linear functions.
-Regularization, avoids the need to have one one observation and one coefficient
-for each spline (in order to avoid instability).
+Regularization, avoids the need to have one observation and one
+coefficient for each spline (in order to avoid instability).
<p>With regularly distributed data points, a spline step corresponding
to the maximum distance between two points in both the east and north
@@ -40,7 +40,7 @@
account for some degree of anisotropy in the distribution of
observation points. Short spline step lengths - especially spline step
lengths that are less than the distance between observation points -
-can greatly increase processing time.
+can greatly increase the processing time.
<p>Moreover, the maximum number of splines for each direction at each
time is fixed, regardless of the spline step length. As the total
@@ -102,7 +102,7 @@
<h3>Basic interpolation and raster output with a longer spline step</h3>
<div class="code"><pre>
-v.surf.bspline input=point_vector raster=interpolate_surface sie=25 sin=25
+v.surf.bspline input=point_vector raster=interpolate_surface ew_step=25 ns_step=25
</pre></div>
A bilinear spline interpolation will be done with a spline step length
@@ -127,7 +127,8 @@
<h3>Using attribute values instead z-coordinates</h3>
<div class="code"><pre>
-v.surf.bspline input=point_vector raster=interpolate_surface layer=1 column=attrib_column
+v.surf.bspline input=point_vector raster=interpolate_surface layer=1 \
+ column=attrib_column
</pre></div>
The interpolation will be done using the values
@@ -136,7 +137,9 @@
<h3>North carolina location example using z-coordinates for interpolation</h3>
<div class="code"><pre>
-v.surf.bspline input=elev_lid792_bepts raster=elev_lid792_rast sie=5 sin=5 method=bicubic lambda_i=0.1
+g.region region=rural_1m res=2 -p
+v.surf.bspline input=elev_lid792_bepts raster=elev_lid792_rast \
+ ew_step=5 ns_step=5 method=bicubic lambda_i=0.1
</pre></div>
More information about the grass-commit
mailing list