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