Watershed analysis

At this point of the tutorial it may be useful to experiment with some more advanced GRASS functionality for hydrologic analysis.

First of all reset the region settings with g.region making them match the DTM map elevation (from the window g.region set elevation in the field Set region to match this raster map, making sure that the resolution is 30 meters, in the Resolution tab).
For this exercise the map elevation will be used as DTM, aspect as aspect and slope as slope.

Sub basins extraction

Display the map elevation (DTM). The area can be divided into a set of basins and sub basins. This operation can be performed automatically using the r.watershed (Raster -> Hydrologic modeling -> Watershed analysis) module to identify all the basins with a minimum area of 1 square km. The r.watershed module can perform further analysis which will be skipped now: some of the parameters can be left to their default values and some of the fields can be left blank (the corresponding maps will not be created).
Conversely, the following parameters must be provided:

r.watershed module interface

r.watershed can be also run in background; you can check the r.watershed module options using g.manual (note that in some cases r.watershed performs better on DEMs which have been smoothed using a neighborhood filter using the mean operator. This can be done using the r.neighbors module; this operation is not necessary in our case.).

basin raster map

The output map from r.watershed identifies each basin with a different integer number, assigned to all the cells within the basin; this value can be detected using the already described procedure for querying raster maps.

r.watershed gives to the map a "random" color table. Using the r.colors module (Raster -> Manage colors -> Color tables), it is possible to change the color table, for example using the "rainbow" color table. Select the map basins in Tab Required and "rainbow" in Colors . Remember to refresh the display by clicking the icon on the Map Display, the color sequence hints at how the basins are numerated starting from the gauge point of each main drainage.

r.colors interface Tab Colors

basins raster map and its sequential sorting of the watershed basins

It can be useful, to better understand the results, to "merge" the output map with the aspect map while displaying. To do this, in the GIS Layer Manager select the button to add a raster map choose Add shaded relief map layer: in the current case select the options as shown in figure.

Shadowing effect

Merge of the aspect map and basins

River network extraction

The accumulation map evaluates the number of cells draining in each cell. Negative values indicate that a cell receives flux from an area outside the current region.

The accumulation map

It is possible to detect the main channels, but not much more. Using r.mapcalc it is possible to compute the logarithm of the absolute value of the accumulation, an useful parameter for hydrologic calculations. On the console, write:

> r.mapcalc 'log_accumulation=log(abs(accumulation)+1)'

Display, then, the new map (log_accumulation):

Accumulation logarithm a parameter of hydrologic model

Use again r.mapcalc to extract from this map the river network: a good threshold value is 6. The value 1 is assigned to the cells with a value whose logarithm is larger than 6, 0 to the other ones.

r.mapcalc 'inf_rivers=if(log_accumulation>6)'

The inf_rivers raster map

Finally, to create a vector map of the river network, use r.thin from the menu, selecting Raster -> Transform features -> Thin to thin out to a single pixel width the lines on the binary raster map (inf_rivers), creating the new raster map thin_rivers.

Using the module r.to.vect (selecting File -> Map type conversions -> Raster to vector on the menu) it is possible to create a new vector map of the river network, named here rivers.

r.to.vect module interface

Display the vector map rivers using the icon in the Display Manager (as seen before when transforming contour lines from raster to vector).

rivers vector map

The rivers map contains, near the upper and left borders, three small lines which are not part of the river network: they belong to basins only for a small part inside the current region (less than 1 square km, in fact they do not have a category in the basins map). To avoid this pathologic situation, only the basins falling completely inside the region should be evaluated. In the menu located in the upper right of the Map Display select the option Digitize: a toolbar containing some buttons to modify vector maps appears. The map to be modified must be selected in the menu in the left. Click on rivers.

The digitizing module in Map Display

Select the icon and click on the small lines with the left mouse button to select them, click on the small lines with the right button to delete them. At the end of the operation, exit from the digitizing module by clicking on the icon , saving all changes.

Overlay this new vector map of the river network to the basins (raster) map: it is possible now to detect the watersheds. The result is shown below:

Results of the processing

It may be interesting to create a vector map showing the basins boundaries.
Select on the menu File -> Map type conversions -> Raster to vector map and set the parameters as shown in the next figure:

r.to.vect interface

Display the new vector map basins: try to hide the "centroids" elements and to assign random colors to the areas (in the next figure the map rivers have been overlaid to give a more complete view).

Display of the basins and rivers vector maps


The suggested resolution for this exercise is 30 meters: it may be interesting to try the same analysis with different resolutions, for example 10 or 50 meters.