1r.watershed(1) Grass User's Manual r.watershed(1)
2
6 r.watershed - Calculates hydrological parameters and RUSLE factors.
7
9 raster, hydrology, watershed, accumulation, drainage, stream network,
10 stream power index, topographic index
11
13 r.watershed
14 r.watershed --help
15 r.watershed [-s4mab] elevation=name [depression=name] [flow=name]
16 [disturbed_land=name] [blocking=name] [threshold=integer]
17 [max_slope_length=float] [accumulation=name] [tci=name]
18 [spi=name] [drainage=name] [basin=name] [stream=name]
19 [half_basin=name] [length_slope=name] [slope_steepness=name]
20 [convergence=integer] [memory=integer] [--overwrite] [--help]
21 [--verbose] [--quiet] [--ui]
22 Flags:
24 -s
25 SFD (D8) flow (default is MFD)
26 SFD: single flow direction, MFD: multiple flow direction
27 -4
29 Allow only horizontal and vertical flow of water
30 -m
32 Enable disk swap memory option: Operation is slow
33 Only needed if memory requirements exceed available RAM; see manual
34 on how to calculate memory requirements
35 -a
37 Use positive flow accumulation even for likely underestimates
38 See manual for a detailed description of flow accumulation output
39 -b
41 Beautify flat areas
42 Flow direction in flat areas is modified to look prettier
43 --overwrite
45 Allow output files to overwrite existing files
46 --help
48 Print usage summary
49 --verbose
51 Verbose module output
52 --quiet
54 Quiet module output
55 --ui
57 Force launching GUI dialog
58 Parameters:
60 elevation=name [required]
61 Name of input elevation raster map
62 depression=name
64 Name of input depressions raster map
65 All non-NULL and non-zero cells are considered as real depressions
66 flow=name
68 Name of input raster representing amount of overland flow per cell
69 disturbed_land=name
71 Name of input raster map percent of disturbed land
72 For USLE
73 blocking=name
75 Name of input raster map blocking overland surface flow
76 For USLE. All non-NULL and non-zero cells are considered as block‐
77 ing terrain.
78 threshold=integer
80 Minimum size of exterior watershed basin
81 max_slope_length=float
83 Maximum length of surface flow in map units
84 For USLE
85 accumulation=name
87 Name for output accumulation raster map
88 Number of cells that drain through each cell
89 tci=name
91 Name for output topographic index ln(a / tan(b)) map
92 spi=name
94 Stream power index a * tan(b)
95 Name for output raster map
96 drainage=name
98 Name for output drainage direction raster map
99 Directions numbered from 1 to 8
100 basin=name
102 Name for output basins raster map
103 stream=name
105 Name for output stream segments raster map
106 half_basin=name
108 Name for output half basins raster map
109 Each half-basin is given a unique value
110 length_slope=name
112 Name for output slope length raster map
113 Slope length and steepness (LS) factor for USLE
114 slope_steepness=name
116 Name for output slope steepness raster map
117 Slope steepness (S) factor for USLE
118 convergence=integer
120 Convergence factor for MFD (1-10)
121 1 = most diverging flow, 10 = most converging flow. Recommended: 5
122 Default: 5
123 memory=integer
125 Maximum memory to be used with -m flag (in MB)
126 Default: 300
127
129 r.watershed generates a set of maps indicating: 1) flow accumulation,
130 drainage direction, the location of streams and watershed basins, and
131 2) the LS and S factors of the Revised Universal Soil Loss Equation
132 (RUSLE).
133
135 Without flag -m set, the entire analysis is run in memory maintained by
136 the operating system. This can be limiting, but is very fast. Setting
137 this flag causes the program to manage memory on disk which allows very
138 large maps to be processed but is slower.
139 Flag -s force the module to use single flow direction (SFD, D8) instead
141 of multiple flow direction (MFD). MFD is enabled by default.
142 By -4 flag the user allow only horizontal and vertical flow of water.
144 Stream and slope lengths are approximately the same as outputs from
145 default surface flow (allows horizontal, vertical, and diagonal flow of
146 water). This flag will also make the drainage basins look more homoge‐
147 neous.
148 When -a flag is specified the module will use positive flow accumula‐
150 tion even for likely underestimates. When this flag is not set, cells
151 with a flow accumulation value that is likely to be an underestimate
152 are converted to the negative. See below for a detailed description of
153 flow accumulation output.
154 Option convergence specifies convergence factor for MFD. Lower values
156 result in higher divergence, flow is more widely distributed. Higher
157 values result in higher convergence, flow is less widely distributed,
158 becoming more similar to SFD.
159 Option elevation specifies the elevation data on which entire analysis
161 is based. NULL (nodata) cells are ignored, zero and negative values are
162 valid elevation data. Gaps in the elevation map that are located
163 within the area of interest must be filled beforehand, e.g. with
164 r.fillnulls, to avoid distortions. The elevation map need not be
165 sink-filled because the module uses a least-cost algorithm.
166 Option depression specifies the optional map of actual depressions or
168 sinkholes in the landscape that are large enough to slow and store sur‐
169 face runoff from a storm event. All cells that are not NULL and not
170 zero indicate depressions. Water will flow into but not out of depres‐
171 sions.
172 Raster flow map specifies amount of overland flow per cell. This map
174 indicates the amount of overland flow units that each cell will con‐
175 tribute to the watershed basin model. Overland flow units represent the
176 amount of overland flow each cell contributes to surface flow. If omit‐
177 ted, a value of one (1) is assumed.
178 Input Raster map or value containing the percent of disturbed land
180 (i.e., croplands, and construction sites) where the raster or input
181 value of 17 equals 17%. If no map or value is given, r.watershed
182 assumes no disturbed land. This input is used for the RUSLE calcula‐
183 tions.
184 Option blocking specifies terrain that will block overland surface
186 flow. Blocking cells and streams stop the slope length for the RUSLE.
187 All cells that are not NULL and not zero indicate blocking terrain.
188 Option threshold specifies the minimum size of an exterior watershed
190 basin in cells, if no flow map is input, or overland flow units when a
191 flow map is given. Warning: low threshold values will dramactically
192 increase run time and generate difficult to read basin and half_basin
193 results. This parameter also controls the level of detail in the
194 stream segments map.
195 Value given by max_slope_length option indicates the maximum length of
197 overland surface flow in meters. If overland flow travels greater than
198 the maximum length, the program assumes the maximum length (it assumes
199 that landscape characteristics not discernible in the digital elevation
200 model exist that maximize the slope length). This input is used for
201 the RUSLE calculations and is a sensitive parameter.
202 Output accumulation map contains the absolute value of the amount of
204 overland flow that traverses each cell. This value will be the number
205 of upland cells plus one if no overland flow map is given. If the
206 overland flow map is given, the value will be in overland flow units.
207 Negative numbers indicate that those cells possibly have surface runoff
208 from outside of the current geographic region. Thus, any cells with
209 negative values cannot have their surface runoff and sedimentation
210 yields calculated accurately.
211 Output tci raster map contains topographic index TCI, computed as
213 ln(α / tan(β)) where α is the cumulative upslope area
214 draining through a point per unit contour length and tan(β) is the
215 local slope angle. The TCI reflects the tendency of water to accumulate
216 at any point in the catchment and the tendency for gravitational forces
217 to move that water downslope (Quinn et al. 1991). This value will be
218 negative if α / tan(β) < 1.
219 Output spi raster map contains stream power index SPI, computed as
221 α * tan(β) where α is the cumulative upslope area drain‐
222 ing through a point per unit contour length and tan(β) is the
223 local slope angle. The SPI reflects the power of water flow at any
224 point in the catchment and the tendency for gravitational forces to
225 move that water downslope (Moore et al. 1991). This value will be neg‐
226 ative if α < 0, i.e. for cells with possible surface runoff from
227 outside of the current geographic region.
228 Output drainage raster map contains drainage direction. Provides the
230 "aspect" for each cell measured CCW from East. Multiplying positive
231 values by 45 will give the direction in degrees that the surface runoff
232 will travel from that cell. The value 0 (zero) indicates that the cell
233 is a depression area (defined by the depression input map). Negative
234 values indicate that surface runoff is leaving the boundaries of the
235 current geographic region. The absolute value of these negative cells
236 indicates the direction of flow. For MFD, drainage indicates the direc‐
237 tion of maximum flow.
238 The output basin map contains unique label for each watershed basin.
240 Each basin will be given a unique positive even integer. Areas along
241 edges may not be large enough to create an exterior watershed basin.
242 NULL values indicate that the cell is not part of a complete watershed
243 basin in the current geographic region.
244 The output stream contains stream segments. Values correspond to the
246 watershed basin values. Can be vectorized after thinning (r.thin) with
247 r.to.vect.
248 The output half_basin raster map stores each half-basin is given a
250 unique value. Watershed basins are divided into left and right sides.
251 The right-hand side cell of the watershed basin (looking upstream) are
252 given even values corresponding to the values in basin. The left-hand
253 side cells of the watershed basin are given odd values which are one
254 less than the value of the watershed basin.
255 The output length_slope raster map stores slope length and steepness
257 (LS) factor for the Revised Universal Soil Loss Equation (RUSLE).
258 Equations taken from Revised Universal Soil Loss Equation for Western
259 Rangelands (Weltz et al. 1987). Since the LS factor is a small number
260 (usually less than one), the GRASS output map is of type DCELL.
261 The output slope_steepness raster map stores slope steepness (S) factor
263 for the Universal Soil Loss Equation (RUSLE). Equations taken from
264 article entitled Revised Slope Steepness Factor for the Universal Soil
265 Loss Equation (McCool et al. 1987). Since the S factor is a small num‐
266 ber (usually less than one), the GRASS output map is of type DCELL.
267 AT least-cost search algorithm
269 r.watershed uses an AT least-cost search algorithm (see REFERENCES sec‐
270 tion) designed to minimize the impact of DEM data errors. Compared to
271 r.terraflow, this algorithm provides more accurate results in areas of
272 low slope as well as DEMs constructed with techniques that mistake
273 canopy tops as the ground elevation. Kinner et al. (2005), for example,
274 used SRTM and IFSAR DEMs to compare r.watershed against r.terraflow
275 results in Panama. r.terraflow was unable to replicate stream locations
276 in the larger valleys while r.watershed performed much better. Thus, if
277 forest canopy exists in valleys, SRTM, IFSAR, and similar data products
278 will cause major errors in r.terraflow stream output. Under similar
279 conditions, r.watershed will generate better stream and half_basin
280 results. If watershed divides contain flat to low slope, r.watershed
281 will generate better basin results than r.terraflow. (r.terraflow uses
282 the same type of algorithm as ESRI’s ArcGIS watershed software which
283 fails under these conditions.) Also, if watershed divides contain for‐
284 est canopy mixed with uncanopied areas using SRTM, IFSAR, and similar
285 data products, r.watershed will generate better basin results than
286 r.terraflow. The algorithm produces results similar to those obtained
287 when running r.cost and r.drain on every cell on the raster map.
288 Multiple flow direction (MFD)
290 r.watershed offers two methods to calculate surface flow: single flow
291 direction (SFD, D8) and multiple flow direction (MFD). With MFD, water
292 flow is distributed to all neighbouring cells with lower elevation,
293 using slope towards neighbouring cells as a weighing factor for propor‐
294 tional distribution. The AT least-cost path is always included. As a
295 result, depressions and obstacles are traversed with a graceful flow
296 convergence before the overflow. The convergence factor causes flow
297 accumulation to converge more strongly with higher values. The sup‐
298 ported range is 1 to 10, recommended is a convergence factor of 5
299 (Holmgren, 1994). If many small sliver basins are created with MFD,
300 setting the convergence factor to a higher value can reduce the amount
301 of small sliver basins.
302 In-memory mode and disk swap mode
304 There are two versions of this program: ram and seg. ram is used by
305 default, seg can be used by setting the -m flag.
306 The ram version requires a maximum of 31 MB of RAM for 1 million cells.
308 Together with the amount of system memory (RAM) available, this value
309 can be used to estimate whether the current region can be processed
310 with the ram version.
311 The ram version uses virtual memory managed by the operating system to
313 store all the data structures and is faster than the seg version; seg
314 uses the GRASS segmentation library which manages data in disk files.
315 seg uses only as much system memory (RAM) as specified with the memory
316 option, allowing other processes to operate on the same system, even
317 when the current geographic region is huge.
318 Due to memory requirements of both programs, it is quite easy to run
320 out of memory when working with huge map regions. If the ram version
321 runs out of memory and the resolution size of the current geographic
322 region cannot be increased, either more memory needs to be added to the
323 computer, or the swap space size needs to be increased. If seg runs out
324 of memory, additional disk space needs to be freed up for the program
325 to run. The r.terraflow module was specifically designed with huge
326 regions in mind and may be useful here as an alternative, although disk
327 space requirements of r.terraflow are several times higher than of seg.
328 Large regions with many cells
330 The upper limit of the ram version is 2 billion (231 - 1) cells,
331 whereas the upper limit for the seg version is 9 billion-billion (263 -
332 1 = 9.223372e+18) cells.
333 In some situations, the region size (number of cells) may be too large
334 for the amount of time or memory available. Running r.watershed may
335 then require use of a coarser resolution. To make the results more
336 closely resemble the finer terrain data, create a map layer containing
337 the lowest elevation values at the coarser resolution. This is done by:
338 1) Setting the current geographic region equal to the elevation map
339 layer with g.region, and 2) Use the r.neighbors or r.resamp.stats com‐
340 mand to find the lowest value for an area equal in size to the desired
341 resolution. For example, if the resolution of the elevation data is 30
342 meters and the resolution of the geographic region for r.watershed will
343 be 90 meters: use the minimum function for a 3 by 3 neighborhood. After
344 changing to the resolution at which r.watershed will be run, r.water‐
345 shed should be run using the values from the neighborhood output map
346 layer that represents the minimum elevation within the region of the
347 coarser cell.
348 Basin threshold
350 The minimum size of drainage basins, defined by the threshold parame‐
351 ter, is only relevant for those watersheds with a single stream having
352 at least the threshold of cells flowing into it. (These watersheds are
353 called exterior basins.) Interior drainage basins contain stream seg‐
354 ments below multiple tributaries. Interior drainage basins can be of
355 any size because the length of an interior stream segment is determined
356 by the distance between the tributaries flowing into it.
357 MASK and no data
359 The r.watershed program does not require the user to have the current
360 geographic region filled with elevation values. Areas without eleva‐
361 tion data (masked or NULL cells) are ignored. It is NOT necessary to
362 create a raster map (or raster reclassification) named MASK for NULL
363 cells. Areas without elevation data will be treated as if they are off
364 the edge of the region. Such areas will reduce the memory necessary to
365 run the program. Masking out unimportant areas can significantly
366 reduce processing time if the watersheds of interest occupy a small
367 percentage of the overall area.
368 Gaps (NULL cells) in the elevation map that are located within the area
370 of interest will heavily influence the analysis: water will flow into
371 but not out of these gaps. These gaps must be filled beforehand, e.g.
372 with r.fillnulls.
373 Zero (0) and negative values will be treated as elevation data (not
375 no_data).
376 Further processing of output layers
378 Problem areas, i.e. those parts of a basin with a likely underestimate
379 of flow accumulation, can be easily identified with e.g.
380 r.mapcalc "problems = if(flow_acc < 0, basin, null())"
381 If the region of interest contains such problem areas, and this is not
382 desired, the computational region must be expanded until the catchment
383 area for the region of interest is completely included.
384 To isolate an individual river network using the output of this module,
386 a number of approaches may be considered.
387 1 Use a resample of the basins catchment raster map as a MASK.
389 The equivalent vector map method is similar using v.select or
390 v.overlay.
391 2 Use the r.cost module with a point in the river as a starting
393 point.
394 3 Use the v.net.iso module with a node in the river as a starting
396 point.
397 All individual river networks in the stream segments output can be
399 identified through their ultimate outlet points. These points are all
400 cells in the stream segments output with negative drainage direction.
401 These points can be used as start points for r.water.outlet or
402 v.net.iso.
403 To create river mile segmentation from a vectorized streams map, try
405 the v.net.iso or v.lrs.segment modules.
406 The stream segments output can be easily vectorized after thinning with
408 r.thin. Each stream segment in the vector map will have the value of
409 the associated basin. To isolate subbasins and streams for a larger
410 basin, a MASK for the larger basin can be created with r.water.outlet.
411 The stream segments output serves as a guide where to place the outlet
412 point used as input to r.water.outlet. The basin threshold must have
413 been sufficiently small to isolate a stream network and subbasins
414 within the larger basin.
415 Given that the drainage is 8 directions numbered counter-clockwise
417 starting from 1 in north-east direction, multiplying the output by 45
418 (by 45. to get a double precision floating point raster map in r.map‐
419 calc) gives the directions in degrees. For most applications, zeros
420 which indicate depressions specified by depression and negative values
421 which indicate runoff leaving the region should be replaced by NULL
422 (null() in r.mapcalc). The following command performs these replace‐
423 ments:
424 r.mapcalc "drainage_degrees = if(drainage > 0, 45. * drainage, null())"
425 Alternatively, the user can use the -a flag or later the abs() function
426 in r.mapcalc if the runoff is leaving the region.
427
429 These examples use the Spearfish sample dataset.
430 Convert r.watershed streams map output to a vector map
432 If you want a detailed stream network, set the threshold option small
433 to create lots of catchment basins, as only one stream is presented per
434 catchment. The r.to.vect -v flag preserves the catchment ID as the vec‐
435 tor category number.
436 r.watershed elev=elevation.dem stream=rwater.stream
437 r.to.vect -v in=rwater.stream out=rwater_stream
438 Set a different color table for the accumulation map:
440 MAP=rwater.accum
441 r.watershed elev=elevation.dem accum=$MAP
442 eval `r.univar -g "$MAP"`
443 stddev_x_2=`echo $stddev | awk ’{print $1 * 2}’`
444 stddev_div_2=`echo $stddev | awk ’{print $1 / 2}’`
445 r.colors $MAP col=rules << EOF
446 0% red
447 -$stddev_x_2 red
448 -$stddev yellow
449 -$stddev_div_2 cyan
450 -$mean_of_abs blue
451 0 white
452 $mean_of_abs blue
453 $stddev_div_2 cyan
454 $stddev yellow
455 $stddev_x_2 red
456 100% red
457 EOF
458 Create a more detailed stream map using the accumulation map and con‐
460 vert it to a vector output map. The accumulation cut-off, and therefore
461 fractal dimension, is arbitrary; in this example we use the map’s mean
462 number of upstream catchment cells (calculated in the above example by
463 r.univar) as the cut-off value. This only works with SFD, not with MFD.
464 r.watershed elev=elevation.dem accum=rwater.accum
465 r.mapcalc ’MASK = if(!isnull(elevation.dem))’
466 r.mapcalc "rwater.course = \
467 if( abs(rwater.accum) > $mean_of_abs, \
468 abs(rwater.accum), \
469 null() )"
470 r.colors -g rwater.course col=bcyr
471 g.remove -f type=raster name=MASK
472 # Thinning is required before converting raster lines to vector
473 r.thin in=rwater.course out=rwater.course.Thin
474 r.colors -gn rwater.course.Thin color=grey
475 r.to.vect in=rwater.course.Thin out=rwater_course type=line
476 v.db.dropcolumn map=rwater_course column=label
477 Create watershed basins map and convert to a vector polygon map
479 r.watershed elev=elevation.dem basin=rwater.basin thresh=15000
480 r.to.vect -s in=rwater.basin out=rwater_basins type=area
481 v.db.dropcolumn map=rwater_basins column=label
482 v.db.renamecolumn map=rwater_basins column=value,catchment
483 Display output in a nice way
485 r.relief map=elevation.dem
486 d.shade shade=elevation.dem.shade color=rwater.basin bright=40
487 d.vect rwater_course color=orange
488
490 · Ehlschlaeger C. (1989). Using the AT Search Algorithm to
491 Develop Hydrologic Models from Digital Elevation Data, Proceed‐
492 ings of International Geographic Information Systems (IGIS)
493 Symposium ’89, pp 275-281 (Baltimore, MD, 18-19 March 1989).
494 URL: http://chuck.ehlschlaeger.info/older/IGIS/paper.html
495 · Holmgren P. (1994). Multiple flow direction algorithms for
497 runoff modelling in grid based elevation models: An empirical
498 evaluation. Hydrological Processes Vol 8(4), 327-334.
499 DOI: 10.1002/hyp.3360080405
500 · Kinner D., Mitasova H., Harmon R., Toma L., Stallard R. (2005).
502 GIS-based Stream Network Analysis for The Chagres River Basin,
503 Republic of Panama. The Rio Chagres: A Multidisciplinary Pro‐
504 file of a Tropical Watershed, R. Harmon (Ed.), Springer/Kluwer,
505 p.83-95.
506 URL: http://www4.ncsu.edu/~hmitaso/measwork/panama/panama.html
507 · McCool et al. (1987). Revised Slope Steepness Factor for the
509 Universal Soil Loss Equation, Transactions of the ASAE Vol
510 30(5).
511 · Metz M., Mitasova H., Harmon R. (2011). Efficient extraction of
513 drainage networks from massive, radar-based elevation models
514 with least cost path search, Hydrol. Earth Syst. Sci. Vol 15,
515 667-678.
516 DOI: 10.5194/hess-15-667-2011
517 · Moore I.D., Grayson R.B., Ladson A.R. (1991). Digital terrain
519 modelling: a review of hydrogical, geomorphological, and bio‐
520 logical applications, Hydrological Processes, Vol 5(1), 3-30
521 DOI: 10.1002/hyp.3360050103
522 · Quinn P., K. Beven K., Chevallier P., Planchon O. (1991). The
524 prediction of hillslope flow paths for distributed hydrological
525 modelling using Digital Elevation Models, Hydrological Pro‐
526 cesses Vol 5(1), p.59-79.
527 DOI: 10.1002/hyp.3360050106
528 · Weltz M. A., Renard K.G., Simanton J. R. (1987). Revised Uni‐
530 versal Soil Loss Equation for Western Rangelands, U.S.A./Mexico
531 Symposium of Strategies for Classification and Management of
532 Native Vegetation for Food Production In Arid Zones (Tucson,
533 AZ, 12-16 Oct. 1987).
534
536 g.region, r.cost, r.drain, r.fillnulls, r.flow, r.mask, r.neighbors,
537 r.param.scale, r.resamp.interp, r.terraflow, r.topidx, r.water.outlet,
538 r.stream.extract
539
541 Original version: Charles Ehlschlaeger, U.S. Army Construction Engi‐
542 neering Research Laboratory
543 Faster sorting algorithm and MFD support: Markus Metz <markus.metz.gis‐
544 work at gmail.com>
545 Last changed: $Date: 2018-10-18 21:13:18 +0200 (Thu, 18 Oct 2018) $
547
549 Available at: r.watershed source code (history)
550 Main index | Raster index | Topics index | Keywords index | Graphical
552 index | Full index
553 © 2003-2019 GRASS Development Team, GRASS GIS 7.4.4 Reference Manual
555GRASS 7.4.4 r.watershed(1)