A distributed hydrological modeling

Derivation of catchment variables using a DEM -digital elevation Model
Altitude (high-low)
Slope
Aspect
Convergence

Analitical hillshading
Surface specific points
Flow accumulation (upslope area)
Wetness index
Stream power indexSlope length factor
Upstream contributing area for a point
Downstream dispersal area
Solar radiation

Abstract Recent development of remote sensing technology and GIS makes it possible to capture and manage a vast amo unt of data of spatially distributed hydrological parameters and variables. Linking GIS and the distributed hydrolo gical model is of rapidly increasing importance. In this work, an optimal channel routing scheme is developed which keeps downstream channel to be routed after its upstream ones and minimizes the requirement on computer memory. Bas ed on this scheme, a distributed hydrological modelling system including automatic procedures of channel network del ineation, extended Horton-Strahler''s channel ordering, extraction of hydrological attributes is constructed. Impleme ntation of a distributed snowmelt runoff model is demonstrated. By this model, snowmelt analysis of the Uono River u sing remotely sensed snow coverage is carried out.

INTRODUCTION

Recent progresses in remote sensing technology and computer science have been improving the availability of hydrolog ical data and computing resources. Although only few remotely sensed data can be applied to hydrology directly, numerous kinds of hydrologically relevant data, especially spatial information, can and can only be derived from remotely sensed information. This largely stimulates the development of distributed hydrological models which explicitly take into account both the spatial information and the conventional hydrometeorological data. Many distributed hydrological models, for example SHE(European Hydrological System) (Abbott et al., 1986a; Abbott et al., 1986b), IHDM(Institute of Hydrology Distributed Model)(Beven et al., 1987), HYDROTEL(Fortin et al., 1985), WATFLOOD(Kouwen, 1988) and Japanese ones(Takasao et al., 1989; Lu et al., 1989) are now under development.
Utilizing distributed hydrological models and the remotely sensed data require powerful and user friendly data processing hardware and software. Geographic information systems (GIS) have been proved to be very useful in handling point, vector data and raster data. GISs provide the means to geo-referenced data which enables the overlaying, merging and visualization of the data. These are key tasks that simplify distributed hydrological modelling. There is, however, a problem in simulating hydrological processes at the time scale shorter than that of the surface water process. At this time scale, the linkage of GISs and hydrological models becomes difficult because simulation of channel flow depends on the structure of channel network heavily; processing can not be done pixel by pixel. Physically, the channel routing should be done only from upstream to downstream. A channel routing scheme is needed which takes into account the structure of channel network.

The purpose of this study is to develop such a distributed hydrological modelling scheme and construct a modelling system using this scheme. This system will allow closer linkage between GISs and hydrological models.


BASIC MODELING CONCEPTS

In this system, a grid approach is used. The basin topography is expressed by a grid of elevation values (a regular grid DEM). All geographical and meteorological data are stored and processed in this grid system. From field survey s and the result of distributed snowmelt analysis, we found that the hydrological processes can be classified into t he processes within a grid block and processes between grid blocks where the grid spacing is as large as several hun dred meters. Based on this concept, a modelling system shown in Fig.1 is constructed.
The processes within a grid block are simulated grid block by grid block, the interaction between grid blocks is thought to be omittable. GIS facility can be used easily for these processes. And the processes between grid blocks, in this case the channel routing, are simulated by using an optimal channel routing scheme which keeps the downstream channel to be routed after its upstream ones and minimizes the requirement on computer memory.


DISTRIBUTED HYDROLOGICAL MODELING SYSTEM

This system consists of following five parts.
Computation of DEM from contour lines

The first part utilizes an algorithm (Lu et al., 1995) that computes hydrologically sound DEM from contour lines by taking into account the position of valley/ridge lines. A DEM is derived by interpolating the elevation profiles at every concerned co-ordinate value. In computation of profiles and interpolation, linear interpolation is used. However linear interpolation cannot correctly determine the value near the hydrologically important valley/ridge lines. In order to overcome this, the valley/ridge lines are extracted by connecting the central points of flat segments on elevation profiles. Along these lines, the longitudinal elevation profiles are computed, and the elevation values at the grid nodes on these flat segments are re-computed. This algorithm gives more accurate values near these lines and makes the channel delineation and the channel routing more accurate.
Delineation of channel network

A procedure (Lu et al., 1989) is developed to automatically delineate a continuous channel network out of both a DEM and channel data digitized from a contour map. For each grid node, the central point of a grid block, its eight neighbours are examined, and the runoff generated in this grid block is considered to flow to the grid node forming the steepest slope. This direction is referred as flow direction, and the link between these two grid nodes is modelled as a channel. Though this process seems to be reasonable and to work successfully, there are cases it does fail. Depressions where there are no positive slopes in all eight directions, intersecting of two channels and closed loop of channels make this procedure fail to compute a continuous channel network. Channel routing is also prohibited. This procedure determines the flow directions at these grid nodes by using primary and secondary processing iteratively. In the primary processing cycle, the elevation of the current grid node is removed by the average value of its neighbours and the flow directions for this grid node and its neighbours are re-computed. In the secondary processing cycle, the neighbour with shortest path to main channels is selected as the exit. This procedure produces a channel network in which all grid nodes are hydrologically connected to the basin outlet

From this channel network and the DEM, many hydrological attributes can be calculated. The slope and aspect of each grid block are computed from the normal vector at the grid node; the channel slope is computed from elevation data of two grid nodes it links; drainage area at each grid node is calculated from the number of grid nodes flowing to it. The Horton-Strahler''s stream order is computed by using an extended Horton-Strahler stream ordering procedure which allows the junctions of more than two streams. The extended ordering procedure can be expressed by

streams which originate at a source are defined to be streams of the first order;

where denotes the stream order of downstream channel, expresses the stream order of the i-th upstream channel which is sorted in descending order. Fig. 2 shows a simple example of this ordering procedure. The number in circle shows the Horton-Strahler order resulted from this procedure. In conjunction with the flow direction, the stream order is used to represent the structure of channel network.

Determination of optimal channel routing order

The routing within a channel network should be processed from upstream to downstream. For a channel network with several thousands of channels, it is necessary to develop an automatic procedure to determine the routing order based on the structure of the channel network. The theorem about optimal routing order is proved theoretically by following discussion.
[ Theorem. ] n buffers are necessary and sufficient for routing a channel network of order n .

[ Proof. ] We prove this theorem by induction on , the order of a channel network.

When , this theorem is an axiom. Here, we will show that it is valid for ( ) whenever it is valid for ( ).

Suppose the channel network of order n is resulted by the junction of m sub-networks. These sub-networks and their orders are referred to as , and , ( i=1, 2, , m ), respectively. For convenience, suppose . Considering Horton-Strahler ordering system, there are two cases:

( ) for all i;
( ) and [ ( ), ( ) ].
In the first case, each sub-network can be routed by using buffer No. 1 through No. n-1, and adding the computed output to buffer No. n. Thus the upstream inflow can be loaded from buffer No. n in the routing of the most downstream channel of this channel network.
In the second case, we consider first. For , there are also two cases. Suppose case 1 is achieved after tracing back l times from . Here the time traced back is represented by superscripts. From previous result, can be routed by using n buffers. Then can be routed by routing and other sub-networks , ( ) successively. This will make the sequence stored in buffer No. n free from conflicting usage. By repeating this procedure, the studied network can be routed by using n buffers.

Conversely, suppose that it is possible to route a channel network of order n with n-1 buffers. To show this is not always possible, consider channel network of order n = 1 consisting more than two channels. It is obviously impossible to route such a channel network without any buffers. consequently, n buffers are needed for channel networks of order n . [ end ]

Obviously, the optimal routing order is not unique. There will be many ways to get such an order. In this study, following computer algorithm that satisfies the requirements in the above proof is used to determine the optimal routing order.

Place current point (CP) at the outlet of the basin;
Search channels flowing to CP, and find the channel with the highest Horton-Strahler order. If channel net work vanishes, finish this procedure. If no inflow channel exists, jump to step 4;
Push CP into a stack, a one dimensional array large enough, and place the new CP at the upstream end of the channel specified in step 2. Then return to step 2;
Write out the CP and its input/output buffer number to a sequential file and cut off that channel from the tree structure. Then pop up previous CP from the stack and return to step 2.
The CP also represents a channel starting from it except for the first one located at the basin outlet. The numbers of input and output buffer of each channel are equal to the stream order of itself and the stream order of the channel it flows to. The optimal routing order (h, i, m, j, e, d, f, b, k, l, g, c, a) can be derived by applying this algorithm to the channel network shown in Fig. 2.

Implementation of channel routing -- buffering technique

The routing computation of each channel is done in the above optimal order using following buffering technique:
reading upstream inflow from the input buffer and filling it with zero;
adding local runoff to upstream inflow to derive the input to the channel;
routing the input to the downstream end of the channel and
adding the resulted sequence to the output buffer.
For the routing model, there may be many alternatives available. At current point, kinematic wave approximation is implemented.

IMPLEMENTATION OF DISTRIBUTED HYDROLOGICAL MODEL

Using this modelling system, the implementation of a distributed hydrological model becomes applying a hydrological model for each grid block and creating a data record for each grid block which includes hydrological attributes used in this model. A sequential file containing these records in the optimal order controls the execution of the distributed model in the following way:
reading a record;
computing runoff from this grid block using the hydrological model; and
routing resulted runoff and inflow from upstream channels using buffering technique.
This paper shows the schematic of the implementation of a distributed snowmelt runoff model.
This model accepts snow cover derived from Landsat TM data and other meteorological data. For each grid block, snowmelt is computed by using a snowmelt model (Koike et al., 1985) in which snowmelt due to radiation, degree-hour and rainfall are computed. The melted water is separated into direct runoff and subsurface runoff by final infiltration capacity . The baseflow is then computed by using the following simple storage drainage model expressing the recession of un-confined groundwater. F71

Here , , are storage, baseflow and subsurface runoff at time t respectively, and is the recession constant. The baseflow is then inputted to the channel network together with the direct runoff

A snowmelt analysis is carried out by using this model on Uono River. The study period is from April 23 to May 30,1993. Fig. 3 shows the hydrological attributes extracted by using this system. In this analysis, information about snow covered area is extracted from 3 consecutive Landsat TM images on April 26, May 12 and May 28 (Fig. 4). Fig. 5 shows the hourly observed and calculated hydrographs.


CONCLUSIONS

An optimal channel routing scheme including the procedure of routing runoff through a channel, the optimal routing order and a buffering technique is developed considering the structure of a channel network. This scheme keeps downstream channel to be routed after its upstream ones and minimizes the requirement on computer memory. Even a personal computer using MS-DOS can be used to run such a distributed model. Based on this scheme, a distributed hydrological modelling system is constructed which makes it very easy to implement a distributed model incorporating remotely sensed data and digitized geographic information as well as conventional hydrological inputs. A distributed snowmelt runoff model is implemented to accept remote sensed hydrological data, and shows very good performance.

REFERENCES

Abbott,M.B., Bathurst,J.C., Cunge,J.A., O''Connell,P.E. and Rasmussen,J.(1986a): An introduction to the European Hydrological System-- Systeme Hydrologique European,''SHE'', 1:history and philosophy of a physically-based distributed modelling system, J. Hydrol., 87,45-59.

Abbott,M.B., Bathurst,J.C., Cunge,J.A., O''Connell,P.E. and Rasmussen,J.(1986b): An introduction to the European Hydrological System-- Systeme Hydrologique European,''SHE'', 2:structure of a physically-based distributed modelling system, J. Hydrol., 87,61-77.

Beven,K., Calver, A. and Morris, E. M.(1987): Institute of hydrology distributed model, Internal Report, Institute of Hydrology, Wallingford.

Fortin,J.P., Villeneuve,J.P., Guilbot,A. and Seguin,B.(1986): Development of a modular hydrological forecasting model based on remotely sensed data, for interactive utilization on a microcomputer, IAHS Publ., 160, 307-319.

Koike,T., Takahasi,Y. and Yosino,S.(1985): A study on modeling of snowmelt distribution(in Japanese), Proc. JSCE, 363, 165-174.

Kouwen,N.(1988): WATFLOOD: a micro-computer based flood forecasting system based on real-time weather radar, Can. Wat. Res. J. , 13-1, 62-77.

Lu,M., Koike,T. and Hayakawa,N.(1989): A model using distributed data of radar rain and altitude(in Japanese), Proc. JSCE, 411/II-12, 135-142.

Lu,M., Koike,T. and Hayakawa,N.(1995): Computation of hydrologically sound grid-based DEM from contour lines considering the position of valleys and ridges(in Japanese), Annual J. Hydraul. Engng, JSCE, 39, 127-132.

Takasao, T., Shiiba, M. and Tachikawa, Y.(1989): Quasi-three-dimensional slope runoff model taking account of Topography of a natural watershed and automatic generation of a basin(in Japanese), Proc. Japanese Conf. on Hydraul., 33, 139-144.

A distributed hydrological modeling system linking GIS and hydrological models