Delineation of the Strmpfelbach subbasin determination of the

Document Sample
Delineation of the Strmpfelbach subbasin determination of the Powered By Docstoc
					Prof. Dr. D. Schröder
Page 1                         GIS in Hydrology and Water ResourceManagement - ENWAT

                         Tutorial: ArcGIS Tools for Hydrology
                        – Hydrological Analysis of a Watershed
In this tutorial, the tools for hydrologic analysis of ArcGIS for surface water will be
introduced. These tools are part of the extension Spatial Analyst, which is an ArcGIS
extension for raster analysis. Thus most of the analysis will be based on raster maps. The
study area is a small sub basin of the Rems watershed, the “Strümpfelbach” sub basin.

The main tasks of this tutorial are:
    Delineation of the Strümpfelbach sub basin
    Determination of the sub basin characteristics, and
    Calculation of a flow time map and a time area diagram for the outlet point

As input data the following raster and vector layers are available:
    DEM with 25 m resolution
    georeferenced topographic raster map
    digitized streams from topographic map
    Corine Land Cover vector map

Remember: all the rasters we are working with in our analysis should have the same
resolution and the same extent. Here all the analysis is based on the DEM. So it’s a good idea
to use its resolution and extent for all the analysis.

As we will create a number of new maps, it’s a good idea to use a strict file naming
convention scheme!

Step 1: preparation of the DEM
For preparation our DEM for hydrological analysis, we should burn-in the existing stream
network and fill the existing sinks.

Burn in existing channels
In flat areas, it is difficult for the accumulation algorithm to determine the flow channel path.
To improve the result a “burn-in” of the existing stream network can be used, i.e. subtracting
from the DEM a meaningful value for the pixel covered by the streams.
So we need a raster, where all the stream pixels have a value of let’s say 2 m, where as all the
other pixels have the valued 0 (Attention: not noData!)

First we have to “rasterize” the stream ploylines, e.g. convert them to a raster:
Conversion Tools  To Raster  Polyline to Raster
Click the environment button to set the extent of the new raster map:
Environment  General settings  Extent set same as layer DEM
For the cell size, you should set the exact value of the DEM (25.170305836741), best copy it,
e.g. from the cell size setting in the raster analysis settings!

After applying the settings, a new raster layer will be created. The raster cells will have cell
values according the field you have used as input (by default the FID field), all cells not
covered by a stream will just have noData as value.
Prof. Dr. D. Schröder
Page 2                         GIS in Hydrology and Water ResourceManagement - ENWAT

Now we have to reclassify the raster.
All the “stream” cells should have the
value 2, all the others just 0, so that
we can subtract this layer from the
DEM layer.

Spatial Analyst Tools  Reclass 

Now we can burn in our streams
using map Algebra:
Spatial Analyst Tools  Map
Algebra  Single Output Map

The syntax is straight
forward, just subtract
the two layers!

Fill sinks
In the next step we have to correct for errors
in the DEM which will have effects on the
hydrological analysis.
Spatial Analyst Tools  Hydrology  Fill

To check the effect of the fill of sinks,
calculate the difference of the original DEM
and with the filled sinks!

Now we are ready for the hydrological
Prof. Dr. D. Schröder
Page 3                          GIS in Hydrology and Water ResourceManagement - ENWAT

Step 2: Hydrological analysis

Flow direction map and flow accumulation map
Now we will calculate the flow direction map and the flow accumulation. With that we can
extract the (synthetic) stream network and delineate the watershed.

Spatial Analyst Tools  Hydrology  Flow Direction (input: filled_DEM)
Spatial Analyst Tools  Hydrology  Flow Accumulation (input: flow_direction)

Extraction of the stream network
There is no special tool for the extraction of the network. Instead we have to use map algebra
that is creating a new raster layer containing only those pixels, which have an accumulation
value bigger than a fixed threshold.

For setting the threshold we should look at our topographic map, checking where here the
streams start. On the other hand, we have to consider the resolution of our DEM of about 25
m, too. For the extracted streams later gully flow will be assumed, for all the other pixels
sheet flow. But our Strümpfelbach is only a small brook, with a top with of about 5 m at its
junction to its headwater. Thus we have to use some compromise for fixing a meaningful

For Map Algebra, we have to use the conditional function con, with the general syntax:
       Con(<condition>, <true_expression>, <false_expression>)

Here e.g.
          con(flowacc_flow2 > 1500,1,0)
with the threshold set to 1500 accumulated pixels corresponding to a drainage area of about
94ha. As the true expression is 1, all pixel values with an accumulation bigger than 1500 will
be assigned 1, all others 0.

With the stream network we now are able to calculate:
- different segments for the network (Spatial Analyst Tools  Hydrology  Stream Link)
- stream order according to Strahler or Shreve (Spatial Analyst Tools  Hydrology  Stream
- convert the (raster) stream network to vector polylines (Spatial Analyst Tools  Hydrology
 Stream to Feature)

Watershed Delineation
Now we are ready to delineate the watershed. In the Spatial Analyst, there are two tools: basin
and watershed. The first one will create all sub basins for the whole study area, where as with
the second you can define watershed upstream to defined outlet points (pour points). So we
will use watershed, as we want to delineate the watershed upstream to the outlet point, i.e. the
junction of the Strümpfelbach to its headwater.

The pour point hast to be located exactly on the pixel of the highest accumulation to
determine the correct watershed. For that, the tool Snap Pour Point has to be used.

Here the steps in detail:
- use ArcCatalog to create a new shape file (geometry type Point)
- digitize the pout point
Prof. Dr. D. Schröder
Page 4                          GIS in Hydrology and Water ResourceManagement - ENWAT

- use the snap to Pour point tool to drag the point exactly on the pixel with highest
accumulation at the outlet.
- delineate the watershed using the watershed tool
- convert the raster to a (vector) polygon (Conversions Tool  Raster to Polygon)

Step 3: Calculate watershed/stream network characteristics
After extracting the streams and delineating the watershed, we can calculate some general
parameters of the watershed:
- Channel length: (see also tutorial 1 on vector data!) add a new field to the table of your
vector streams and use Calculate Geometry (only in ArcGIS 9.2! look in the help of the field
calculator for older versions how to calculate the length). Select the main channel segments
and use Statistics to get the sum of the segments (about 5 km)
- Channel slope: use the identify tool to get the height values for the outlet point and the
source of your main stream (2.7%)
- Drainage area: calculate Geometry for the polygon or number of pixels x area of a pixel
(949ha polygon, pixel 14981 x 625= 936ha)
- Drainage density: Select all segments of your streams inside the polygon using Selection by
location. Divide the length by the drainage area (0.6 17 1/km)
- Hydrological length: Extract by mask from the (filled) DEM the pixels inside the watershed,
recalculate flow direction map, use the tool Flow length, the maximum value is the
hydrological length (6.3km).
- Watershed slope: From the pixels inside the watershed take the difference of the highest and
the lowest value and divide it by the hydrological length (480-247)/6300 = 3.6%
- Length of the centroid: use Data Management Tools  Features  Feature to Point to
calculate the centroid. Use the measurement tool to measure the approximate distance to the
perpendicular toe point on the stream. (There are tools to make it better, but remember, a GIS
is not a CAD, so it is rather cumbersome to de it exact!) (2.9km)
- Shape Factor Ll = (LLca)0.3  2.39 km2/3
- Circularity ratio Fc = P/(4πA)0.5 (P=14392m; Fc=1.32)
- Circularity ration Rc = A/Ao (U=2π r  r = 2274 A0 = pi r2  A0 = 162 Rc =5.85)
- Ruggedness Number = basin relief times drainage density.(m/km) 143

Step 4: Flow velocities
Based on the flow velocities, the travel time from each pixel to the outlet can be calculated.
Calculation of the velocity map is the most crucial step in the workflow. The flow velocity for
each pixel depends on the slope and the roughness of the surface. To calculate the velocities,
for instance the Manning-Strickler equation can be used. Here we have carefully to
distinguish between sheet and gully flow, i.e. different parameters for the roughness as well as
for the hydraulic radius should be used.
The roughness itself is determined by the land use. In the literature table for the Manning
value n can be found (Kst is the inverse of the Manning value n), for instance Haestad
Methods: Floodplain modelling using HEC-RAS, 2003.

Sheet flow
Surface description                              Kst in m1/3s -1
Smooth surface (concrete, asphalt, gravel, or    90
bare soil)
Fallow                                           20
Cultivated soils                                 10
Grass                                            5
Wood                                             17
Prof. Dr. D. Schröder
Page 5                          GIS in Hydrology and Water ResourceManagement - ENWAT

Channel flow
Surface description                               Kst in m1/3s -1
Minor streams, some grass and weeds               30

The Corine Land Cover layer can be used for assigning meaningful averaged roughness
values by adding a new field (Kst) to the attribute table and edit it manually.
After adding the values, we are ready to rasterize the map using the Kst values. Now we have
to do some map algebra to calculate the flow velocities. In complex expression, It is a good
idea to split up the calculation into different steps.
Use map Algebra to calculate first the Kst x r 2/3 map, where r=0.03, so r2/3=0.0965

For modelling the channel flow, we have to keep in mind that the map resolution is 25m, so
that only a part of a pixel is really covered by the channel. For the rest we will still have sheet
flow. So again, some averaged value should be used as roughness and hydraulic radius.
For the first order streams, let’s assume a top width of 1 m and a depth of 0.3, so that for a
rectangular channel shape we have a wetted parameter of 1.6m and the hydraulic radius
(A/Pw) becomes 0.1875 m. For second order streams assume a top width of 2 m and a depth
of 0.5 m, i.e. a wetted parameter of 3m and a hydraulic radius 0.333m.
As not the whole pixel will be covered by the channel, we use a “corrected” Kst of 20.

Thus we have as kst r2/3 for our channels: first order 7 and second order 10.

Now we can Reclassify the Strahler order map assigning the values above to the stream pixels,
all other pixels should have the value 0.

Now we have to combine the two maps with sheet and gully flow. One possibility is to first
“burn out” the streams in the sheet flow map by using a binary stream map (1 for no stream, 0
for stream) and multiply this with the sheet flow map. Secondly just sum the two maps to get
the final kst r2/3 map for combined sheet and gully flow.

To get the final velocities map, we have to multiply this map with the square root slope map.
Creating the slope map, you should take care of the units. Using percentage we have to divide
the numbers first by 100, than we can multiply to get the final velocity map.

Step 5 Calculating flow time and the simple IUH
As the flow time is the quotient of path and velocity, we can use Flow length to calculate it.
We have only to introduce the inverse of the velocity as weights. As a weight raster has to be
an integer raster, we have first to convert it to integer.

Dividing by 3600 we will get the travel time in h as an isochron map. The time area diagram
we get from the histogram of this map. For an overview, use the histogram (accumulated)
from the symbology dialog. If you need more details, export (copy and paste) the VAT of the
raster layer to Excel. As the ordinates of the IUH are the differential changes of the drainage
area with respect to time, the ordinate of the IUH are the gradients of the time area diagram
and can be easily calculated.
Prof. Dr. D. Schröder
Page 6                                 GIS in Hydrology and Water ResourceManagement - ENWAT

                                                                   Flow time map

 Time area diagram (histogram from ArcGIS)



    pixel count

                  800                                                  IUH in minutes
                                                                       (Excel based)



                         0   10   20       30       40   50   60
                                       time (min)

Shared By:
dominic.cecilia dominic.cecilia http://