← Back to branch summary

~ubuntu-branches/ubuntu/vivid/grass/vivid-proposed

« back to all changes in this revision

Viewing changes to vector/v.vol.rst/v.vol.rst.html

  • Committer: Package Import Robot
  • Author(s): Bas Couwenberg
  • Date: 2015-02-20 23:12:08 UTC
  • mfrom: (8.2.6 experimental)
  • Revision ID: package-import@ubuntu.com-20150220231208-1u6qvqm84v430b10
Tags: 7.0.0-1~exp1
* New upstream release.
* Update python-ctypes-ternary.patch to use if/else instead of and/or.
* Drop check4dev patch, rely on upstream check.
* Add build dependency on libpq-dev to grass-dev for libpq-fe.h.
* Drop patches applied upstream, refresh remaining patches.
* Update symlinks for images switched from jpg to png.

Show diffs side-by-side

added added

removed removed

Lines of Context:
vector/v.vol.rst/v.vol.rst.html
 
1
<h2>DESCRIPTION</h2>
 
2
 
 
3
<p><em>v.vol.rst</em> interpolates values to a 3-dimensional raster map from
 
4
3-dimensional point data (e.g. temperature, rainfall data from climatic
 
5
stations, concentrations from drill holes etc.) given in a 3-D vector
 
6
point file named <b>input</b>.&nbsp; The size of the output 
 
7
3D raster map <b>elevation</b> is given by the current 3D region. Sometimes, the
 
8
user
 
9
may want to get a 2-D map showing a modelled phenomenon at a
 
10
crossection surface. In that case, <b>cross_input</b> and <b>cross_output</b>
 
11
options must be specified, with the output 2D raster map <b>cross_output</b>
 
12
containing the crossection of the interpolated volume with a surface 
 
13
defined by <b>cross_input</b>
 
14
2D raster map. As an option, simultaneously with interpolation, 
 
15
geometric parameters of the interpolated
 
16
phenomenon can be computed (magnitude of gradient, direction of
 
17
gradient defined by horizontal and vertical angles), change of gradient,
 
18
Gauss-Kronecker curvature, or mean curvature). These geometric
 
19
parameteres are saved as
 
20
3D raster maps <b>gradient, aspect_horizontal, aspect_vertical, ncurvature, gcurvature, mcurvature</b>,
 
21
respectively. Maps <b>aspect_horizontal</b> and <b>aspect_vertical</b> are in degrees.
 
22
 
 
23
<p>At first, data points are checked for identical positions and points
 
24
that are closer to each other than given <b>dmin</b> are removed.
 
25
Parameters <b>wmult</b> and <b>zmult</b> allow the user to re-scale 
 
26
the w-values and z-coordinates of the point data (useful e.g. for 
 
27
transformation of elevations given in feet to meters, so that the 
 
28
proper values of gradient and curvatures can be computed).
 
29
Rescaling of z-coordinates (<b>zmult</b>) is also needed when the distances
 
30
in vertical direction are much smaller than the horizontal
 
31
distances; if that is the case, the value of <b>zmult</b>
 
32
should be selected so that the vertical and horizontal distances
 
33
have about the same magnitude.
 
34
 
 
35
<p>Regularized spline with tension method is used in the interpolation.
 
36
The <b>tension</b> parameter controls the distance over which 
 
37
each given point influences the resulting volume (with very high tension, 
 
38
each point influences only its close neighborhood and the volume goes
 
39
rapidly to trend between the points).
 
40
Higher values of tension parameter reduce the overshoots that
 
41
can appear in volumes with rapid change of gradient. For noisy data, it
 
42
is possible to define a global smoothing parameter, <b>smooth</b>.
 
43
With the
 
44
smoothing parameter set to zero (<b>smooth=0</b>) the resulting volume
 
45
passes exactly through the data points. 
 
46
When smoothing is used, it is possible to output a vector map <b>deviations</b>
 
47
containing deviations of the resulting volume from the given data. 
 
48
<p>The user can define a 2D raster map named <b>maskmap</b>, which will
 
49
be used as a mask. The interpolation is skipped for 3-dimensional cells
 
50
whose 2-dimensional projection has a zero value in the mask. Zero values will
 
51
be assigned to these cells in all output 3D raster maps. 
 
52
<p>If the number of given points is greater than 700, segmented
 
53
processing is used. The region is split into 3-dimensional "box"
 
54
segments, each having less than <b>segmax</b> points and interpolation
 
55
is performed on each segment of the region. To ensure the smooth
 
56
connection of segments, the interpolation function for each segment is
 
57
computed using the points in the given segment
 
58
and the points in its neighborhood. The minimum number of points taken
 
59
for interpolation is controlled by <b>npmin</b> , the value of which
 
60
must be larger than <b>segmax</b> and less than 700. This limit of 700 was
 
61
selected to ensure the numerical stability and efficiency of the
 
62
algorithm. 
 
63
 
 
64
 
 
65
<h3>SQL support</h3>
 
66
 
 
67
Using the <b>where</b> parameter, the interpolation can be limited to use
 
68
only a subset of the input vectors.
 
69
 
 
70
<div class="code"><pre>
 
71
# preparation as in above example
 
72
v.vol.rst elevrand_3d wcol=soilrange elevation=soilrange zmult=100 where="soilrange > 3"
 
73
</pre></div>
 
74
 
 
75
 
 
76
<h3>Cross validation procedure</h3>
 
77
 
 
78
Sometimes it can be difficult to figure out the proper values of
 
79
interpolation parameters. In this case, the user can use a
 
80
crossvalidation procedure using <b>-c</b> flag (a.k.a. "jack-knife"
 
81
method) to find optimal parameters for given data. In this method,
 
82
every point in the input point file is temporarily excluded from the
 
83
computation and interpolation error for this point location is
 
84
computed.  During this procedure no output grid files can be
 
85
simultanuously computed.  The procedure for larger datasets may take a
 
86
very long time, so it might be worth to use just a sample data
 
87
representing the whole dataset.
 
88
 
 
89
<p><i>Example 
 
90
 (based on <a href="http://www.grassbook.org/data_menu2nd.php">Slovakia3d dataset</a>):</i>
 
91
<p><div class="code"><pre>
 
92
v.info -c precip3d
 
93
g.region n=5530000 s=5275000 w=4186000 e=4631000 res=500 -p
 
94
v.vol.rst -c input=precip3d wcolumn=precip zmult=50 segmax=700 cvdev=cvdevmap tension=10
 
95
v.db.select cvdevmap
 
96
v.univar cvdevmap col=flt1 type=point
 
97
</pre></div>
 
98
 
 
99
Based on these results, the parameters will have to be optimized. It is
 
100
recommended to plot the CV error as curve while modifying
 
101
the parameters.
 
102
<p>The best approach is to start with <b>tension</b>, <b>smooth</b>
 
103
and <b>zmult</b> with rough steps, or to set <b>zmult</b> to a
 
104
constant somewhere between 30-60. This helps to find minimal RMSE
 
105
values while then finer steps can be used in all parameters. The
 
106
reasonable range is <b>tension</b>=10...100,
 
107
<b>smooth</b>=0.1...1.0, <b>zmult</b>=10...100.
 
108
<p>In <em>v.vol.rst</em> the tension parameter is much more sensitive to
 
109
changes than in <em>v.surf.rst</em>, 
 
110
therefore the user should always check the
 
111
result by visual inspection. Minimizing CV does not always provide the best
 
112
result, especially when the density of data are insufficient. Then
 
113
the optimal result found by CV is an oversmoothed surface.
 
114
 
 
115
<h2>NOTES</h2>
 
116
The vector points map must be a 3D vector map (x, y, z as geometry).
 
117
The module <a href="v.in.db.html">v.in.db</a> can be used to generate
 
118
a 3D vector map from a table containing x,y,z columns.
 
119
 
 
120
Also, the input data should be in a projected coodinate system, such as
 
121
Univeral Transverse Mercator. The module does not appear to have support for 
 
122
geographic (Lat/Long) coordinates as of May 2009.
 
123
 
 
124
<p><em>v.vol.rst</em> uses regularized spline with tension for
 
125
interpolation from point data (as described in Mitasova and Mitas,
 
126
1993). The implementation has an improved segmentation procedure based
 
127
on Oct-trees which enhances the efficiency for large data sets. 
 
128
 
 
129
<p>Geometric parameters - magnitude of gradient (<b>gradient</b>),
 
130
horizontal (<b>aspect_horizontal</b>) and vertical (<b>aspect_vertical)</b>aspects,
 
131
change of gradient (<b>ncurvature</b>), Gauss-Kronecker (<b>gcurvature</b>) and
 
132
mean curvatures (<b>mcurvature</b>) are computed directly from the
 
133
interpolation function so that the important relationships between
 
134
these parameters are preserved. More information on these parameters
 
135
can be found in Mitasova et al., 1995 or Thorpe, 1979.
 
136
 
 
137
<p>The program gives warning when significant overshoots appear and
 
138
higher tension should be used. However, with tension too high the
 
139
resulting volume will have local maximum in each given point
 
140
and everywhere else the volume goes rapidly to trend. With a smoothing
 
141
parameter greater than zero, the volume will not pass through the data 
 
142
points and the higher the parameter the closer the volume will be to the 
 
143
trend. For theory on smoothing with splines see Talmi and Gilat, 1977 or Wahba, 1990. 
 
144
 
 
145
<p>If a visible connection of segments appears, the program should be
 
146
rerun with higher <b>npmin</b> to get more points from the
 
147
neighborhood of given segment. 
 
148
 
 
149
<p>If the number of points in a vector map is less than 400, <b>segmax</b>
 
150
should be set to 400 so that segmentation is not performed when it is
 
151
not necessary. 
 
152
 
 
153
<p>The program gives a warning when the user wants to interpolate outside the
 
154
"box" given by minimum and maximum coordinates in the input vector map. 
 
155
To remedy this, zoom into the area encompassing the input vector data points.
 
156
 
 
157
<p>For large data sets (thousands of data points), it is suggested to
 
158
zoom into a smaller representative area and test whether the parameters
 
159
chosen (e.g. defaults) are appropriate. 
 
160
 
 
161
<p>The user must run <em>g.region</em> before the program to set the
 
162
3D region for interpolation. 
 
163
 
 
164
 
 
165
<h2>EXAMPLES</h2>
 
166
 
 
167
<!-- TODO: find better data. This example is nonsensical :-) -->
 
168
Spearfish example (we first simulate 3D soil range data):
 
169
 
 
170
<div class="code"><pre>
 
171
g.region -dp
 
172
# define volume
 
173
g.region res=100 tbres=100 res3=100 b=0 t=1500 -ap3
 
174
 
 
175
### First part: generate synthetic 3D data (true 3D soil data preferred)
 
176
# generate random positions from elevation map (2D)
 
177
r.random elevation.10m vector_output=elevrand n=200
 
178
 
 
179
# generate synthetic values
 
180
v.db.addcolumn elevrand col="x double precision, y double precision"
 
181
v.to.db elevrand option=coor col=x,y
 
182
v.db.select elevrand
 
183
 
 
184
# create new 3D map
 
185
v.in.db elevrand out=elevrand_3d x=x y=y z=value key=cat
 
186
v.info -c elevrand_3d
 
187
v.info -t elevrand_3d
 
188
 
 
189
# remove the now superfluous 'x', 'y' and 'value' (z) columns
 
190
v.db.dropcolumn elevrand_3d col=x
 
191
v.db.dropcolumn elevrand_3d col=y
 
192
v.db.dropcolumn elevrand_3d col=value
 
193
 
 
194
# add attribute to have data available for 3D interpolation
 
195
# (Soil range types taken from the USDA Soil Survey)
 
196
d.mon wx0
 
197
d.rast soils.range
 
198
d.vect elevrand_3d
 
199
v.db.addcolumn elevrand_3d col="soilrange integer"
 
200
v.what.rast elevrand_3d col=soilrange rast=soils.range
 
201
 
 
202
# fix 0 (no data in raster map) to NULL:
 
203
v.db.update elevrand_3d col=soilrange value=NULL where="soilrange=0"
 
204
v.db.select elevrand_3d
 
205
 
 
206
# optionally: check 3D points in Paraview
 
207
v.out.vtk input=elevrand_3d output=elevrand_3d.vtk type=point dp=2
 
208
paraview --data=elevrand_3d.vtk
 
209
 
 
210
### Second part: 3D interpolation from 3D point data
 
211
# interpolate volume to "soilrange" voxel map
 
212
v.vol.rst input=elevrand_3d wcol=soilrange elevation=soilrange zmult=100
 
213
 
 
214
# visualize I: in GRASS GIS wxGUI
 
215
g.gui
 
216
# load: 2D raster map: elevation.10m
 
217
#       3D raster map: soilrange
 
218
 
 
219
# visualize II: export to Paraview
 
220
r.mapcalc "bottom = 0.0"
 
221
r3.out.vtk -s input=soilrange top=elevation.10m bottom=bottom dp=2 output=volume.vtk
 
222
paraview --data=volume.vtk
 
223
</pre></div>
 
224
 
 
225
<h2>BUGS</h2>
 
226
<b>deviations</b> file is written as 2D and deviations are not written as attributes.
 
227
 
 
228
<h2>REFERENCES</h2>
 
229
<p>Hofierka J., Parajka J., Mitasova H., Mitas L., 2002, Multivariate
 
230
Interpolation of Precipitation Using Regularized Spline with Tension.
 
231
Transactions in
 
232
GIS&nbsp; 6, pp. 135-150. 
 
233
<p><a href="http://www4.ncsu.edu/~hmitaso/gmslab/">Mitas, L.,
 
234
Mitasova, H.</a>, 1999, Spatial Interpolation. In: P.Longley, M.F.
 
235
Goodchild, D.J. Maguire, D.W.Rhind (Eds.), Geographical Information
 
236
Systems: Principles, Techniques, Management and Applications, Wiley,
 
237
pp.481-492 
 
238
<p>Mitas L., Brown W. M., Mitasova H., 1997,
 
239
<a href="http://www4.ncsu.edu/~hmitaso/gmslab/lcgfin/cg-mitas.html">Role
 
240
of dynamic cartography in simulations of landscape processes based on
 
241
multi-variate fields.</a> Computers and Geosciences, Vol. 23, No. 4,
 
242
pp. 437-446 (includes CDROM and WWW: www.elsevier.nl/locate/cgvis) 
 
243
<p>Mitasova H., Mitas L.,&nbsp; Brown W.M.,&nbsp; D.P. Gerdes, I.
 
244
Kosinovsky, Baker, T.1995, Modeling spatially and temporally
 
245
distributed phenomena:
 
246
New methods and tools for GRASS GIS. International Journal of GIS, 9
 
247
(4),
 
248
special issue on Integrating GIS and Environmental modeling, 433-446. 
 
249
<p> Mitasova, H., Mitas, L., Brown, B., Kosinovsky, I., Baker, T.,
 
250
Gerdes, D. (1994):
 
251
<a href="http://www4.ncsu.edu/~hmitaso/gmslab/viz/ches.html">Multidimensional
 
252
interpolation and visualization in GRASS GIS</a> 
 
253
<p><a href="http://www4.ncsu.edu/~hmitaso/gmslab/papers/lmg.rev1.ps">Mitasova
 
254
H. and Mitas L. 1993</a>: Interpolation by Regularized Spline with
 
255
Tension: I. Theory and Implementation, <i>Mathematical Geology</i> 25,
 
256
641-655. 
 
257
<p><a href="http://www4.ncsu.edu/~hmitaso/gmslab/papers/hmg.rev1.ps">Mitasova
 
258
H. and Hofierka J. 1993</a>: Interpolation by Regularized Spline with
 
259
Tension: II. Application to Terrain Modeling and Surface Geometry
 
260
Analysis, <i>Mathematical Geology</i> 25, 657-667. 
 
261
<p>Mitasova, H., 1992 : New capabilities for interpolation and
 
262
topographic analysis in GRASS, GRASSclippings 6, No.2 (summer), p.13. 
 
263
<p>Wahba, G., 1990 : Spline Models for Observational Data, CNMS-NSF
 
264
Regional Conference series in applied mathematics, 59, SIAM,
 
265
Philadelphia, Pennsylvania. 
 
266
<p>Mitas, L., Mitasova H., 1988 : General variational approach to the
 
267
interpolation problem, Computers and Mathematics with Applications 16,
 
268
p. 983 
 
269
<p>Talmi, A. and Gilat, G., 1977 : Method for Smooth Approximation of
 
270
Data, Journal of Computational Physics, 23, p.93-123. 
 
271
<p>Thorpe, J. A. (1979): Elementary Topics in Differential Geometry.
 
272
Springer-Verlag, New York, pp. 6-94.
 
273
 
 
274
<h2>SEE ALSO</h2>
 
275
 
 
276
<em>
 
277
<a href="g.region.html">g.region</a>,
 
278
<a href="v.in.ascii.html">v.in.ascii</a>,
 
279
<a href="r3.mask.html">r3.mask</a>,
 
280
<a href="v.in.db.html">v.in.db</a>,
 
281
<a href="v.surf.rst.html">v.surf.rst</a>,
 
282
<a href="v.univar.html">v.univar</a>
 
283
</em>
 
284
 
 
285
<h2>AUTHOR</h2>
 
286
 
 
287
Original version of program (in FORTRAN) and GRASS enhancements: <br>
 
288
Lubos Mitas, NCSA, University of Illinois at Urbana-Champaign,
 
289
Illinois, USA, since 2000 at Department of Physics, 
 
290
North Carolina State University, Raleigh, USA
 
291
<a href="mailto:lubos_mitas@ncsu.edu">lubos_mitas@ncsu.edu</a>
 
292
<br>
 
293
Helena Mitasova, Department of Marine, Earth and Atmospheric Sciences,
 
294
North Carolina State University, Raleigh, USA,
 
295
<a href="mailto:hmitaso@unity.ncsu.edu">hmitaso@unity.ncsu.edu</a>
 
296
<p>
 
297
Modified program (translated to C, adapted for GRASS, new
 
298
segmentation procedure): <br>
 
299
Irina Kosinovsky, US Army CERL, Champaign, Illinois, USA <br>
 
300
Dave Gerdes, US Army CERL, Champaign, Illinois, USA 
 
301
<p>
 
302
Modifications for g3d library, geometric parameters,
 
303
cross-validation, deviations: <br>
 
304
Jaro Hofierka, Department of Geography and Regional Development,
 
305
University of Presov, Presov, Slovakia,
 
306
<a href="mailto:hofierka@fhpv.unipo.sk">hofierka@fhpv.unipo.sk</a>,
 
307
<a href="http://www.geomodel.sk">http://www.geomodel.sk</a>
 
308
 
 
309
<p><i>Last changed: $Date: 2014-11-28 17:25:40 +0100 (Fri, 28 Nov 2014) $</i>