- Email: [email protected]
Can the behaviour of different basins be described by the same model’s parameter set? A geomorphologic framework
Physics and Chemistry of the Earth 28 (2003) 289–295 www.elsevier.com/locate/pce
Can the behaviour of different basins be described by the same modelÕs parameter set? A geomorphologic framework F. Giannoni *, G. Roth, R. Rudari CIMA––Centro di ricerca Interuniversitario in Monitoraggio Ambientale and DIAM––Dipartimento di Ingegneria Ambientale, University of Genova, Via Cadorna 7, 17100 Savona, Italy Received 19 October 2002; received in revised form 24 February 2003; accepted 6 March 2003
Abstract Different basin scales and environments usually are a constrain in the calibration phase of hydrological modelsÕ parameters. In the present work this aspect is faced with regard to the calibration and validation of a geomorphologic semi-distributed rainfall– runoff model at different scales. Simulations of the hydrologic response at basin scale are performed here running the Discharge River Forecast model (Phys. Chem. Earth 25 (7–8) (2000) 665). This model is focused on the efficient description of the drainage system in its essential parts. It uses five parameters: two geomorphologic parameters for the drainage network identification; two kinematic parameters to address the time scale of flood formation on hillslopes and in the channel network; the last parameter takes into account the soil antecedent moisture conditions. Parameters calibration and validation have been carried out using intense rainfall events in different basins and sub-basins in a wide range of sizes belonging to different geographic areas: the Liguria region and the upper part of Po basin. Ó 2003 Elsevier Science Ltd. All rights reserved. Keywords: Geomorphology; Rainfall–runoff; Floods; Digital terrain model
1. Introduction The contemporary development of both digital terrain model (DTM) based hydrology––which reflected the necessity of modelling the distribute features of the territory––and of morphological and fractal concepts, led to a synergy between these two fields of study that produced a large amount of scientific and operational works (OÕCallaghan and Mark, 1985; Band, 1986; Tarboton et al., 1988, 1989a,b). In this paper a rainfall– runoff model, which takes advantage of the geomorphological information obtainable from DTMs, is analysed. The model is called Discharge River Forecast (DRiFt) (Giannoni et al., 2000). The importance of a feasible description of the river network in different morphologies is pointed out. This paper is outlined as follows. The global frame of the DRiFt model is given in Section 2; the filtering procedure used to derive the river
* Corresponding author. Address: ARPAL––Agenzia Regionale per lÕAmbiente Ligure, P.zza della Vittoria 1/15, 16100 Genova, Italy. Tel.: +39-019-230271; fax: +39-019-862612. E-mail addresses: [email protected] (F. Giannoni), [email protected] (G. Roth), [email protected] (R. Rudari).
network from data stored in DTMs is explained in detail in Section 3. Section 4 describes the parametersÕ calibration and validation with regard to the Liguria region, that is the morphological region where the DRiFt model has been first applied. In Sections 5 and 6 the extension of the validity of this model to another morphologic region is considered: the parametersÕ estimation in the Tanaro valley together with some results is presented. Finally, the summary and concluding remarks are enlightened in Section 7.
2. Model structure The model frame is quite simple and only the main characteristics of the hydrograph are focused: peak and time to peak. Its parameters have been evaluated by fitting computed hydrographs to observed ones, placing this model in the calibrated parameters models category. The parameters have an intuitive and direct physical meaning; it is therefore possible to control their values after the calibration. This is important because it simplifies the validation procedure, but, moreover, because it shows the goodness of the processesÕ schematisation
1474-7065/03/$ - see front matter Ó 2003 Elsevier Science Ltd. All rights reserved. doi:10.1016/S1474-7065(03)00040-8
290
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
used. We defined it a semi-distributed model, because it takes into account the spatial variation of inputs such as rainfall, morphologic, geological and anthropic characteristics of the basin, but it is lumped in parameters, and outputs can be obtained in a given location wherever in the catchment. From a technical point of view this is the main feature of the model concealing most of the advantages of both distributed and lumped models. The geomorphologic filter, described in detail in Section 3, is used to distinguish hillslopes and channel network starting from the space-filling representation of the drainage system directly obtainable from DTMs. The routing time of each site in the basin is evaluated assigning different typical velocity values in each pixel pertaining to the basin and classified as hillslope or channel. The two velocities, vh and vc , used to describe the flow routing process in each of the two components of the drainage system are assumed here constant; they maintain a physical meaning as the average velocities on hillslopes and in channel network, related to a particular class of extreme events, and therefore to a particular range of discharge, and to the particular class of basin considered. The Instantaneous Unit Hydrograph (IUH) technique (Sherman, 1932), used to determine the basin response is applied with some differences with respect to the classical procedure. Here in fact, the runoff generation is computed in each pixel of the basin at any time step considering the soil infiltration capacity heterogeneity and the rainfall field variability. This leads to a time variant unit response of the basin allowing a better representation of the real basin behaviour excited by a particular rainfall distribution. Schematically, discharge at any location along the drainage network is represented by the expression: Z dh ðxÞ dc ðxÞ QðtÞ ¼ M t ; x dx ð1Þ vh vc B where B is the drainage basin above the specified location, Mðt; xÞ is the runoff rate (mm h1 ) at time t and location x, dh ðxÞ denotes the distance from x to the closest stream channel and dc ðxÞ denotes the distance from x to the outlet of the basin specified by the region B. The runoff rate Mðt; xÞ is computed from the rainfall rate Rðt; xÞ using the Soil Conservation Service method (McCuen, 1982) with the Curve Number parameters locally defined, according to lithology, soil type, land cover and use.
3. Channel identification from DTM Most of the automated procedures available in literature to operationally derive the drainage network starting from the information retrieved by DTMs data lead, as final result, to a space-filling description. As it
is well known, this type of representation of the river pattern does not interpret properly the transport processes on the hillslopes and in the channelled network, especially when this kind of approach is targeted to rainfall–runoff modelling. In this case it is necessary to set a physical limit to the upstream development of the channelled network. From a physical point of view, this means to set a minimum hillslope length for the channel initiation. The filtering procedure used in this work to obtain a more suitable description of the river network from the DTM, is a further development of the area-threshold concept introduced by Tarboton et al., 1992. Montgomery and Dietrich (1992) investigated with a field study this topic pointing out the existence of a channelization area and slope threshold able to set a finite scale to the landscape. This constrains hillslopes to assume the length required by overland flow to initiate surface erosion and, consequently, to form and maintain a channel. The threshold criterion involves the interaction between climate and soil throughout erosion and sediment processes. Further studies related to this topic (Rodriguez-Iturbe et al., 1992; La Barbera and Roth, 1994) agree in considering the constancy of the quantity AS k for points belonging to the channel-network, where A is the drained area, S is the local slope and k is a geomorphologic exponent. In some previous works (Roth et al., 1996) this filter has been applied in order to describe the topological characteristics of the network constraining the drainage density to be equal to a particular value established a priori in order to guarantee physical control; this is due also to the fact that no literature background, regarding the threshold value to use, is available. In the present work it has been possible to apply the filter procedure without imposing any external constrains because of its direct meaning in the rainfall–runoff model context. The threshold value in fact influences the calibration of the model on the basis of recorded hydrographs. The filtering procedure can be synthesized in the following steps: (i) pits removal, (ii) evaluation of a quasi-local slope (Giannoni et al., 2000)––by taking into account the slope of up and downstream links of the cell studied in place of a local slope, problems due to the sampling DTMs errors are partially avoided; (iii) computation of the quantity AS k in each pixel of the basin grid and application of the slope-area filter by comparing the threshold and the local AS k values, pixels characterized by threshold value greater than the chosen threshold are assumed to represent the effective channelled structure. Actually, in order to eliminate unrealistic channel path interruptions due mainly to DTM errors, the matrix is scanned from mountain to valley zones and when a channel-pixel is found, the channel is traced to the outlet following the maximum slope direction.
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
291
Table 1 Liguria region and Tanaro basin parameters
k AS k (m2 ) vc (m s1 ) vh (m s1 )
Liguria region
Tanaro basin Mountain
Pede mountain
Valley
1.7 105 2.5 0.16
1.7 105 2.5 0.16
1.7 105 2.2 0.16
1.7 105 1.5 0.12
4. Models parameters: Liguria environment To stress the correct leading processes and their representation is the core of a robust semi-distributed model. The major advantages of this class of models is their ability in describing the whole catchment system with a limited set of conceptual, but physically understandable, parameters. One of the most interesting qualities of this model is its simplicity in inputs needed and in the limited number of parameters used. To be semi-distributed here means to be lumped in parameters: this is without any doubt an important characteristic of the model itself. In fact the capability of using only some lumped parameters to describe the average behaviour of the basin is a way to make easier the hydrologic simulation reducing the uncertainty due to parametersÕ calibration and validation. On the other hand, DRiFt is distributed in considering the actual morphology of the basin from a topological point of view, and in considering the spatial variability of the soil characteristics and spatial and temporal variability of the precipitation fields used to stress the model. This model has five parameters. Two morphologic parameters are used to illustrate the morphology of the environment where the model is applied. In particular they are the threshold value AS k , called morphologic threshold, and the geomorphologic exponent, k. Both of them have a long literature background with regard to morphologic studies (Roth et al., 1996), but only some recent studies applied this methodology within the rainfall–runoff modelling context (Giannoni et al., 2000). The other two parameters used by the model are the two constant velocities applied to interpret the flow routing on the hillslope elements and afterwards in the channel. The use of two velocities allowed us to introduce the scale of the Geomorphologic Instantaneous Unit Hydrograph (GIUH), (Rodriguez-Iturbe and Valdes, 1979), considered the weak point of this kind of approach (Shamseldin and Nash, 1998). This scale is here defined within a physically based framework, directly linked to pluviometric and hydrometric records. The model uses also another parameter to set up the correct degree of soil moisture to start the simulation. As a matter of fact this last is not a real parameter but must be more reasonably interpreted as a degree of freedom left to the user. The parameters estimation, performed on several historical
observations in different Ligurian basins has led to the set of parameters schematised in Table 1.
5. Liguria and Tanaro: comparison between two different morphologies Liguria region spreads in an arch from the mouth of Roja to that of the Magra Rivers, embracing the south side of the Ligurian Alps and Apennines as well as a large part of the Po Valley flanks. The Ligurian lithology is quite heterogeneous: metamorphic rocks and green stones prevail in the Alps and in the Apennines, while calcareous and dolomiteous rocks, sand stone, clayey rocks and several parallel crests are present to form longitudinal valleys in the east part. Most of the territory is mountainous or hilly with narrow strips of fairly low terrain along tracts of the coast or near several low alluvial valleys. Numerous valleys penetrate the mountains: those to the south cut mainly across the lie of the mountains, and their rivers are generally fastflowing torrents. To the north of the watershed, the mountains are broken high in the valleys by tributaries of the Po, the major is the Tanaro River. The Tanaro source is in the Ligurian Alps (Fig. 1), very close to the Marguareis mountains; in the flat area when it joints the Po River it drains about 8000 km2 , more than 500 km2
Fig. 1. Localization of Liguria region and Tanaro river basin.
292
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
are in the mountain part while the others are in the plain. The Tanaro basin dimensions (length and area) make this basin very heterogeneous, showing characteristics quite far from the Alpine basins and from the Apennine basins too. The major tributaries from the mountain towards the valley region are Stura di Demonte (1430 km2 ), Alto Tanaro (1840 km2 ), Belbo (480 km2 ), Bormida (1640 km2 ) and Orba (840 km2 ). The morphologic variability shown allows the identification of three main zones having a peculiar behaviour. A mountain part that shows a 6% average slope, very steep catchments and very deep river bed, the mild part with 1% average slope, with a shallow river bed and steep catchments, and finally the alluvial part, characterized by very low slope values.
6. Models parameters calibration and validation: Tanaro region 6.1. Morphologic parameters estimation In this study, as far as in the Liguria region, the morphologic exponent k has been assumed a priori equal to 1.7 according to literature background. In many previous works (Mandelbrot, 1983; Tarboton et al., 1990; Rosso et al., 1991) the fractal measure of the whole drainage system was assumed equal to the product D d, where D represents the river branching fractal dimension, while d is the fractal dimension of single rivers. The evidence that the network tends to drain the whole river basin constrains its fractal dimension to approach 2. Being usually d equal 1.136 from field experiments (Tarboton et al., 1990), D tends to be equal to 1.7. Roth et al. (1996) derived in a theoretical framework the link between k and D. Moreover, in the study by (Flint, 1974) the exponent introduced in this work took an average value of 1.7. The threshold value, AS k , calibration has consisted in some subsequent phases, following the same procedure already applied in the Ligurian environment. First the goodness of the network identified has been checked by comparing it with the cartography blue lines, then other important features, directly depending from the network such as hillslope length and drainage density distribution, IUH shape and travelling time distribution, have been considered. In order to make clear the morphologic threshold estimation, we present the comparison of all the features we have focused on showing the results obtained using only three threshold values: 104 , 105 and 106 m2 . In fact we are convinced that these three values, different for one order of magnitude, are sufficient to sum up the results of the work, describing the different behaviours we met. First of all we want to verify the river network consistency: to be in accord with a natural one, the
identified drainage network must show a non uniform drainage density, higher in the upper part and lower in the plain zone, where the gentle slope does not support any more the channel initiation. In Fig. 2 the drainage networks obtained using the three different threshold values mentioned above are compared. The lowest threshold value (104 m2 ) does not allow the identification of the river network in the mountain part, and also the plain part shows a non-reasonable drainage density value. The value 105 m2 instead allows a good identification of the network; the drainage density is clearly decreasing from mountain towards valley (Fig. 3) and all the morphologic characteristics concerning the topology of the network are preserved. The network obtained with the highest threshold (106 m2 ) value maintains the major characteristics of the river network, but besides this, it loses a lot of information about the topology of the network itself; especially in the north-western part of the basin. Here the drainage network is not any more recognizable: the Borbore subbasin (8E–45N) is completely lost. A quantitative analysis concerning the different drainage network obtained is summarized in Fig. 4 that shows, respectively, channel and hillslope percentages. For values that are less than 105 m2 the percentages are very sensitive to the threshold value changes, while once this value is exceeded a substantial constancy is obtained. This first piece of information leads to choose the threshold value in the field of low variability. In fact the theory linked to this type of study states the constancy of the AS k quantity in the channelled network. This theory predicts that the percentage of channel reaches asymptotically a quasi-stable value by increasing the threshold and then drops dramatically when the threshold value is exceeded. Even if this kind of behaviour is not strictly respected in nature, actually it should be possible to observe an asymptotic trend in the neighbourhoods of the threshold value. The second test that has to be performed once identified some networks is to analyse the distribution of concentration time to the outlet through the inspection of the IUH. In Fig. 5 it is quite clear how a different threshold value might affect the hydrologic response. In a distributed physically based model, the coherence between the shape of the basin itself and the profile of the instantaneous hydrologic response must be respected; in particular the complexity of the basin must be reflected by the shape of the IUH: for the Tanaro the presence of two distinct, same importance major sub-basins must be observable in the IUH that in this case shows a clear bimodality. Using the higher threshold value, 106 m2 , this bimodal behaviour is lost and the shape of the IUH is unnaturally skewed to the right. To put some emphasis on this concept, in Fig. 6 the two IUH corresponding to the main tributaries of the Tanaro River are presented. It is evident that the global IUH is the result
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
293
Fig. 2. River network filtered using different morphologic threshold values: AS k ¼ 104 m2 (a), AS k ¼ 105 m2 (b) and AS k ¼ 106 m2 (c). In the middle panel some points are highlighted where the drainage density is computed along the Stura sub-catchment.
Fig. 3. Dimensionless drainage density, Stura sub-catchment.
of the linear combination of the tributariesÕ IUH, highlighted in Fig. 7. 6.2. Kinematic parameters estimation The estimation of kinematic parameters has been performed using different recorded hydrographs in several hydrometric stations. The best fit calibration technique identified three set of kinematic parameters representing the variability of the Tanaro region morphology: the presence of a mountain zone characterized by very steep catchments, a Pede-mountain part and a
Fig. 4. Hillslope and channel percentages as functions of the threshold value.
valley alluvial part. The average velocities decrease from steep mountain sub-catchments to valley basins both on hillslopes and channels. Their values maintain always a physical meaning with regard to the morphology considered (Mancini et al., 1994) supporting the conviction of a good description capacity of the catchment system offered by the model. In Fig. 8 a validation run is presented. As it is possible to observe the model performance, evaluated on the basis of the peak and the time to peak capturing, is completely satisfying. The runoff volumes are consistent and only
294
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
Fig. 5. IUH correspondent to the network obtained using different threshold values, AS k ¼ 104 m2 , (continuous thin line), AS k ¼ 105 m2 , (continuous thick line), AS k ¼ 106 m2 (dotted thick line).
Fig. 8. Validation event: Tanaro at Alba.
7. Conclusions
Fig. 6. IUH correspondent to different sub-catchments: Orba and Stura (dotted line), total IUH (continuous line).
Fig. 7. TributariesÕ sub-catchments: Stura (light grey), Orba (dark grey).
in the initial and final part of the hydrograph the simulated differs from the observed one consistently. This single case was here chosen to represent the validation phase being the other results comparable with this one. In Table 1 parameters estimation results for the Tanaro River are summarised.
The semi-distributed approach to rainfall–runoff modelling here presented allows the comprehensive use of commonly available distributed information, i.e. digital elevation models, rainfall data, soil characteristics and land use. On the basis of this knowledge, a lumped parameter set allows drainage network characterization and surface runoff routing to the site of interest for hydrologic modelling purposes. To this aim, four lumped parameters are defined. Two parameters are used to identify the channelled portion of the drainage network. The drainage density redistribution in the mountain and valley parts of the basin is related to k; its value was derived, in a theoretical framework, from drainage network fractal properties. The quantity AS k is calibrated on the basis of catchment drainage density and hydrologic response characteristics. Surface runoff routing is then performed with reference to calibrated kinematic parameters that maintain a physical meaning as average velocities, vh and vc , on the two drainage components: hillslopes and channels. The modelÕs properties outlined above are of foremost importance when the model is used for flood discharge prediction in non-gauged sites, where parameters direct estimation from recorded events is not possible. In particular, the capacity to treat different, measured or predicted, distributed rainfall inputs, the invariance of the parameters used to identify the channelled portion of the drainage and the strong physical meaning of the kinematic parameters, allows the use of the model for flood prediction on the basis of available quantitative precipitation forecasts. Results obtained from the calibration and validation phases show that the parameters used to characterize the drainage network can be assumed constant for the entire basin set here analysed, i.e. from 100 to 1000 km2 for basins in the Liguria region to 8000 km2 for the Tanaro basin. On the other side, values of the kinematic parameters appear to be sensitive to the average drainage slope. With respect to this, basins pertaining to the Liguria environment can be assumed to be homoge-
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
neous and characterized by behaviour similar to that of sub-basins pertaining to the mountain part of the Tanaro basin. The Tanaro basin is not homogeneous; its sub-basins are placed in environments––mountain, mild and alluvial––characterized, from a hydrologic modelling point of view, by different kinematic parameters. In this case, the physical meaning assigned to the parameters is of paramount importance: the slow decreasing of the two average velocities moving from mountain towards alluvial catchments is a strong indication of the coherence of the calibration phase. Acknowledgements The research was partially supported by the Italian National Research Council, through GNDCI, the Italian National Group for the Protection from Hydrological Hazards. References Band, L.E., 1986. Topographic partition of watershed with digital elevation models. Water Resour. Res. 22 (1), 15–24. Flint, J.J., 1974. Stream gradient as a function of order, magnitude and discharge. Water Resour. Res. 10 (5), 969–973. Giannoni, F., Roth, G., Rudari, R., 2000. A semi-distributed rainfall– runoff model based on a geomorphologic approach. Phys. Chem. Earth 25 (7–8), 665–671. La Barbera, P., Roth, G., 1994. Invariance and scaling properties in the distributions of contributing area and energy in drainage basins. Hydrol. Proc. 8, 125–135. Mancini, M., Troch, P.A., Wood, E.F., 1994. Overland flow routing over a range of catchment scales. In: Rosso, R., Peano, A., Becchi,
295
I., Bemporad, G.A. (Eds.), Advances in Distributed Hydrology. Water Resource Publication, Highlands Ranch, CO, pp. 225–244. Mandelbrot, B.B., 1983. The Fractal Geometry of Nature. Freeman, New York. McCuen, R.H., 1982. A Guide to Hydrologic Analysis using SCS Methods. Prentice-Hall, Englewood Cliffs, NY. Montgomery, D.R., Dietrich, W.E., 1992. Channel initiation and the problem of landscape scale. Science 225, 826–830. OÕCallaghan, J.F., Mark, D.M., 1985. The extraction of drainage networks from digital elevation data. Comput. Vision Graph. Image Process 28, 323–344. Rodriguez-Iturbe, I., Valdes, J.B., 1979. The geomorphologic structure of hydrologic response. Water Resour. Res. 15 (6), 1409–1420. Rodriguez-Iturbe, I., Rinaldo, A., Rigon, R., Bras, R.L., Marani, A., Ijjasz-Vasquez, E., 1992. Energy dissipation, runoff production, and the three-dimensional structure of river basins. Water Resour. Res. 28 (4), 1095–1103. Rosso, R., Bacchi, B., Barbera, P., 1991. Fractal relation of mainstream length to catchment area in river networks. Water Resour. Res. 27 (3), 381–387. Roth, G., Barbera, P., Greco, M., 1996. On the description of the basin drainage structure. J. Hydrol. 187, 119–135. Shamseldin, A.Y., Nash, J.E., 1998. The geomorphological unit hydrograph––a critical review. Hess 2 (1), 1–8. Sherman, L.K., 1932. Streamflow from rainfall by the unit hydrograph method. Eng. News Rec. 108, 501–505. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1988. The fractal nature of river networks. Water Resour. Res. 24 (8), 1317–1322. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1989a. Scaling and elevation in river networks. Water Resour. Res. 25 (9), 2037–2051. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1989b. The analysis of river basins and channel network using digital terrain data. Technical Report 326, Ralph. M. Parsons Laboratory, MIT, Cambridge, MA. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1990. Comment on: On the fractal dimension of stream networks by P. La Barbera and R. Rosso. Water Resour Res. 6 (8), 2243–2244. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1992. A physical basis for drainage density. Geomorphology 5, 59–76.
Can the behaviour of different basins be described by the same modelÕs parameter set? A geomorphologic framework F. Giannoni *, G. Roth, R. Rudari CIMA––Centro di ricerca Interuniversitario in Monitoraggio Ambientale and DIAM––Dipartimento di Ingegneria Ambientale, University of Genova, Via Cadorna 7, 17100 Savona, Italy Received 19 October 2002; received in revised form 24 February 2003; accepted 6 March 2003
Abstract Different basin scales and environments usually are a constrain in the calibration phase of hydrological modelsÕ parameters. In the present work this aspect is faced with regard to the calibration and validation of a geomorphologic semi-distributed rainfall– runoff model at different scales. Simulations of the hydrologic response at basin scale are performed here running the Discharge River Forecast model (Phys. Chem. Earth 25 (7–8) (2000) 665). This model is focused on the efficient description of the drainage system in its essential parts. It uses five parameters: two geomorphologic parameters for the drainage network identification; two kinematic parameters to address the time scale of flood formation on hillslopes and in the channel network; the last parameter takes into account the soil antecedent moisture conditions. Parameters calibration and validation have been carried out using intense rainfall events in different basins and sub-basins in a wide range of sizes belonging to different geographic areas: the Liguria region and the upper part of Po basin. Ó 2003 Elsevier Science Ltd. All rights reserved. Keywords: Geomorphology; Rainfall–runoff; Floods; Digital terrain model
1. Introduction The contemporary development of both digital terrain model (DTM) based hydrology––which reflected the necessity of modelling the distribute features of the territory––and of morphological and fractal concepts, led to a synergy between these two fields of study that produced a large amount of scientific and operational works (OÕCallaghan and Mark, 1985; Band, 1986; Tarboton et al., 1988, 1989a,b). In this paper a rainfall– runoff model, which takes advantage of the geomorphological information obtainable from DTMs, is analysed. The model is called Discharge River Forecast (DRiFt) (Giannoni et al., 2000). The importance of a feasible description of the river network in different morphologies is pointed out. This paper is outlined as follows. The global frame of the DRiFt model is given in Section 2; the filtering procedure used to derive the river
* Corresponding author. Address: ARPAL––Agenzia Regionale per lÕAmbiente Ligure, P.zza della Vittoria 1/15, 16100 Genova, Italy. Tel.: +39-019-230271; fax: +39-019-862612. E-mail addresses: [email protected] (F. Giannoni), [email protected] (G. Roth), [email protected] (R. Rudari).
network from data stored in DTMs is explained in detail in Section 3. Section 4 describes the parametersÕ calibration and validation with regard to the Liguria region, that is the morphological region where the DRiFt model has been first applied. In Sections 5 and 6 the extension of the validity of this model to another morphologic region is considered: the parametersÕ estimation in the Tanaro valley together with some results is presented. Finally, the summary and concluding remarks are enlightened in Section 7.
2. Model structure The model frame is quite simple and only the main characteristics of the hydrograph are focused: peak and time to peak. Its parameters have been evaluated by fitting computed hydrographs to observed ones, placing this model in the calibrated parameters models category. The parameters have an intuitive and direct physical meaning; it is therefore possible to control their values after the calibration. This is important because it simplifies the validation procedure, but, moreover, because it shows the goodness of the processesÕ schematisation
1474-7065/03/$ - see front matter Ó 2003 Elsevier Science Ltd. All rights reserved. doi:10.1016/S1474-7065(03)00040-8
290
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
used. We defined it a semi-distributed model, because it takes into account the spatial variation of inputs such as rainfall, morphologic, geological and anthropic characteristics of the basin, but it is lumped in parameters, and outputs can be obtained in a given location wherever in the catchment. From a technical point of view this is the main feature of the model concealing most of the advantages of both distributed and lumped models. The geomorphologic filter, described in detail in Section 3, is used to distinguish hillslopes and channel network starting from the space-filling representation of the drainage system directly obtainable from DTMs. The routing time of each site in the basin is evaluated assigning different typical velocity values in each pixel pertaining to the basin and classified as hillslope or channel. The two velocities, vh and vc , used to describe the flow routing process in each of the two components of the drainage system are assumed here constant; they maintain a physical meaning as the average velocities on hillslopes and in channel network, related to a particular class of extreme events, and therefore to a particular range of discharge, and to the particular class of basin considered. The Instantaneous Unit Hydrograph (IUH) technique (Sherman, 1932), used to determine the basin response is applied with some differences with respect to the classical procedure. Here in fact, the runoff generation is computed in each pixel of the basin at any time step considering the soil infiltration capacity heterogeneity and the rainfall field variability. This leads to a time variant unit response of the basin allowing a better representation of the real basin behaviour excited by a particular rainfall distribution. Schematically, discharge at any location along the drainage network is represented by the expression: Z dh ðxÞ dc ðxÞ QðtÞ ¼ M t ; x dx ð1Þ vh vc B where B is the drainage basin above the specified location, Mðt; xÞ is the runoff rate (mm h1 ) at time t and location x, dh ðxÞ denotes the distance from x to the closest stream channel and dc ðxÞ denotes the distance from x to the outlet of the basin specified by the region B. The runoff rate Mðt; xÞ is computed from the rainfall rate Rðt; xÞ using the Soil Conservation Service method (McCuen, 1982) with the Curve Number parameters locally defined, according to lithology, soil type, land cover and use.
3. Channel identification from DTM Most of the automated procedures available in literature to operationally derive the drainage network starting from the information retrieved by DTMs data lead, as final result, to a space-filling description. As it
is well known, this type of representation of the river pattern does not interpret properly the transport processes on the hillslopes and in the channelled network, especially when this kind of approach is targeted to rainfall–runoff modelling. In this case it is necessary to set a physical limit to the upstream development of the channelled network. From a physical point of view, this means to set a minimum hillslope length for the channel initiation. The filtering procedure used in this work to obtain a more suitable description of the river network from the DTM, is a further development of the area-threshold concept introduced by Tarboton et al., 1992. Montgomery and Dietrich (1992) investigated with a field study this topic pointing out the existence of a channelization area and slope threshold able to set a finite scale to the landscape. This constrains hillslopes to assume the length required by overland flow to initiate surface erosion and, consequently, to form and maintain a channel. The threshold criterion involves the interaction between climate and soil throughout erosion and sediment processes. Further studies related to this topic (Rodriguez-Iturbe et al., 1992; La Barbera and Roth, 1994) agree in considering the constancy of the quantity AS k for points belonging to the channel-network, where A is the drained area, S is the local slope and k is a geomorphologic exponent. In some previous works (Roth et al., 1996) this filter has been applied in order to describe the topological characteristics of the network constraining the drainage density to be equal to a particular value established a priori in order to guarantee physical control; this is due also to the fact that no literature background, regarding the threshold value to use, is available. In the present work it has been possible to apply the filter procedure without imposing any external constrains because of its direct meaning in the rainfall–runoff model context. The threshold value in fact influences the calibration of the model on the basis of recorded hydrographs. The filtering procedure can be synthesized in the following steps: (i) pits removal, (ii) evaluation of a quasi-local slope (Giannoni et al., 2000)––by taking into account the slope of up and downstream links of the cell studied in place of a local slope, problems due to the sampling DTMs errors are partially avoided; (iii) computation of the quantity AS k in each pixel of the basin grid and application of the slope-area filter by comparing the threshold and the local AS k values, pixels characterized by threshold value greater than the chosen threshold are assumed to represent the effective channelled structure. Actually, in order to eliminate unrealistic channel path interruptions due mainly to DTM errors, the matrix is scanned from mountain to valley zones and when a channel-pixel is found, the channel is traced to the outlet following the maximum slope direction.
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
291
Table 1 Liguria region and Tanaro basin parameters
k AS k (m2 ) vc (m s1 ) vh (m s1 )
Liguria region
Tanaro basin Mountain
Pede mountain
Valley
1.7 105 2.5 0.16
1.7 105 2.5 0.16
1.7 105 2.2 0.16
1.7 105 1.5 0.12
4. Models parameters: Liguria environment To stress the correct leading processes and their representation is the core of a robust semi-distributed model. The major advantages of this class of models is their ability in describing the whole catchment system with a limited set of conceptual, but physically understandable, parameters. One of the most interesting qualities of this model is its simplicity in inputs needed and in the limited number of parameters used. To be semi-distributed here means to be lumped in parameters: this is without any doubt an important characteristic of the model itself. In fact the capability of using only some lumped parameters to describe the average behaviour of the basin is a way to make easier the hydrologic simulation reducing the uncertainty due to parametersÕ calibration and validation. On the other hand, DRiFt is distributed in considering the actual morphology of the basin from a topological point of view, and in considering the spatial variability of the soil characteristics and spatial and temporal variability of the precipitation fields used to stress the model. This model has five parameters. Two morphologic parameters are used to illustrate the morphology of the environment where the model is applied. In particular they are the threshold value AS k , called morphologic threshold, and the geomorphologic exponent, k. Both of them have a long literature background with regard to morphologic studies (Roth et al., 1996), but only some recent studies applied this methodology within the rainfall–runoff modelling context (Giannoni et al., 2000). The other two parameters used by the model are the two constant velocities applied to interpret the flow routing on the hillslope elements and afterwards in the channel. The use of two velocities allowed us to introduce the scale of the Geomorphologic Instantaneous Unit Hydrograph (GIUH), (Rodriguez-Iturbe and Valdes, 1979), considered the weak point of this kind of approach (Shamseldin and Nash, 1998). This scale is here defined within a physically based framework, directly linked to pluviometric and hydrometric records. The model uses also another parameter to set up the correct degree of soil moisture to start the simulation. As a matter of fact this last is not a real parameter but must be more reasonably interpreted as a degree of freedom left to the user. The parameters estimation, performed on several historical
observations in different Ligurian basins has led to the set of parameters schematised in Table 1.
5. Liguria and Tanaro: comparison between two different morphologies Liguria region spreads in an arch from the mouth of Roja to that of the Magra Rivers, embracing the south side of the Ligurian Alps and Apennines as well as a large part of the Po Valley flanks. The Ligurian lithology is quite heterogeneous: metamorphic rocks and green stones prevail in the Alps and in the Apennines, while calcareous and dolomiteous rocks, sand stone, clayey rocks and several parallel crests are present to form longitudinal valleys in the east part. Most of the territory is mountainous or hilly with narrow strips of fairly low terrain along tracts of the coast or near several low alluvial valleys. Numerous valleys penetrate the mountains: those to the south cut mainly across the lie of the mountains, and their rivers are generally fastflowing torrents. To the north of the watershed, the mountains are broken high in the valleys by tributaries of the Po, the major is the Tanaro River. The Tanaro source is in the Ligurian Alps (Fig. 1), very close to the Marguareis mountains; in the flat area when it joints the Po River it drains about 8000 km2 , more than 500 km2
Fig. 1. Localization of Liguria region and Tanaro river basin.
292
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
are in the mountain part while the others are in the plain. The Tanaro basin dimensions (length and area) make this basin very heterogeneous, showing characteristics quite far from the Alpine basins and from the Apennine basins too. The major tributaries from the mountain towards the valley region are Stura di Demonte (1430 km2 ), Alto Tanaro (1840 km2 ), Belbo (480 km2 ), Bormida (1640 km2 ) and Orba (840 km2 ). The morphologic variability shown allows the identification of three main zones having a peculiar behaviour. A mountain part that shows a 6% average slope, very steep catchments and very deep river bed, the mild part with 1% average slope, with a shallow river bed and steep catchments, and finally the alluvial part, characterized by very low slope values.
6. Models parameters calibration and validation: Tanaro region 6.1. Morphologic parameters estimation In this study, as far as in the Liguria region, the morphologic exponent k has been assumed a priori equal to 1.7 according to literature background. In many previous works (Mandelbrot, 1983; Tarboton et al., 1990; Rosso et al., 1991) the fractal measure of the whole drainage system was assumed equal to the product D d, where D represents the river branching fractal dimension, while d is the fractal dimension of single rivers. The evidence that the network tends to drain the whole river basin constrains its fractal dimension to approach 2. Being usually d equal 1.136 from field experiments (Tarboton et al., 1990), D tends to be equal to 1.7. Roth et al. (1996) derived in a theoretical framework the link between k and D. Moreover, in the study by (Flint, 1974) the exponent introduced in this work took an average value of 1.7. The threshold value, AS k , calibration has consisted in some subsequent phases, following the same procedure already applied in the Ligurian environment. First the goodness of the network identified has been checked by comparing it with the cartography blue lines, then other important features, directly depending from the network such as hillslope length and drainage density distribution, IUH shape and travelling time distribution, have been considered. In order to make clear the morphologic threshold estimation, we present the comparison of all the features we have focused on showing the results obtained using only three threshold values: 104 , 105 and 106 m2 . In fact we are convinced that these three values, different for one order of magnitude, are sufficient to sum up the results of the work, describing the different behaviours we met. First of all we want to verify the river network consistency: to be in accord with a natural one, the
identified drainage network must show a non uniform drainage density, higher in the upper part and lower in the plain zone, where the gentle slope does not support any more the channel initiation. In Fig. 2 the drainage networks obtained using the three different threshold values mentioned above are compared. The lowest threshold value (104 m2 ) does not allow the identification of the river network in the mountain part, and also the plain part shows a non-reasonable drainage density value. The value 105 m2 instead allows a good identification of the network; the drainage density is clearly decreasing from mountain towards valley (Fig. 3) and all the morphologic characteristics concerning the topology of the network are preserved. The network obtained with the highest threshold (106 m2 ) value maintains the major characteristics of the river network, but besides this, it loses a lot of information about the topology of the network itself; especially in the north-western part of the basin. Here the drainage network is not any more recognizable: the Borbore subbasin (8E–45N) is completely lost. A quantitative analysis concerning the different drainage network obtained is summarized in Fig. 4 that shows, respectively, channel and hillslope percentages. For values that are less than 105 m2 the percentages are very sensitive to the threshold value changes, while once this value is exceeded a substantial constancy is obtained. This first piece of information leads to choose the threshold value in the field of low variability. In fact the theory linked to this type of study states the constancy of the AS k quantity in the channelled network. This theory predicts that the percentage of channel reaches asymptotically a quasi-stable value by increasing the threshold and then drops dramatically when the threshold value is exceeded. Even if this kind of behaviour is not strictly respected in nature, actually it should be possible to observe an asymptotic trend in the neighbourhoods of the threshold value. The second test that has to be performed once identified some networks is to analyse the distribution of concentration time to the outlet through the inspection of the IUH. In Fig. 5 it is quite clear how a different threshold value might affect the hydrologic response. In a distributed physically based model, the coherence between the shape of the basin itself and the profile of the instantaneous hydrologic response must be respected; in particular the complexity of the basin must be reflected by the shape of the IUH: for the Tanaro the presence of two distinct, same importance major sub-basins must be observable in the IUH that in this case shows a clear bimodality. Using the higher threshold value, 106 m2 , this bimodal behaviour is lost and the shape of the IUH is unnaturally skewed to the right. To put some emphasis on this concept, in Fig. 6 the two IUH corresponding to the main tributaries of the Tanaro River are presented. It is evident that the global IUH is the result
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
293
Fig. 2. River network filtered using different morphologic threshold values: AS k ¼ 104 m2 (a), AS k ¼ 105 m2 (b) and AS k ¼ 106 m2 (c). In the middle panel some points are highlighted where the drainage density is computed along the Stura sub-catchment.
Fig. 3. Dimensionless drainage density, Stura sub-catchment.
of the linear combination of the tributariesÕ IUH, highlighted in Fig. 7. 6.2. Kinematic parameters estimation The estimation of kinematic parameters has been performed using different recorded hydrographs in several hydrometric stations. The best fit calibration technique identified three set of kinematic parameters representing the variability of the Tanaro region morphology: the presence of a mountain zone characterized by very steep catchments, a Pede-mountain part and a
Fig. 4. Hillslope and channel percentages as functions of the threshold value.
valley alluvial part. The average velocities decrease from steep mountain sub-catchments to valley basins both on hillslopes and channels. Their values maintain always a physical meaning with regard to the morphology considered (Mancini et al., 1994) supporting the conviction of a good description capacity of the catchment system offered by the model. In Fig. 8 a validation run is presented. As it is possible to observe the model performance, evaluated on the basis of the peak and the time to peak capturing, is completely satisfying. The runoff volumes are consistent and only
294
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
Fig. 5. IUH correspondent to the network obtained using different threshold values, AS k ¼ 104 m2 , (continuous thin line), AS k ¼ 105 m2 , (continuous thick line), AS k ¼ 106 m2 (dotted thick line).
Fig. 8. Validation event: Tanaro at Alba.
7. Conclusions
Fig. 6. IUH correspondent to different sub-catchments: Orba and Stura (dotted line), total IUH (continuous line).
Fig. 7. TributariesÕ sub-catchments: Stura (light grey), Orba (dark grey).
in the initial and final part of the hydrograph the simulated differs from the observed one consistently. This single case was here chosen to represent the validation phase being the other results comparable with this one. In Table 1 parameters estimation results for the Tanaro River are summarised.
The semi-distributed approach to rainfall–runoff modelling here presented allows the comprehensive use of commonly available distributed information, i.e. digital elevation models, rainfall data, soil characteristics and land use. On the basis of this knowledge, a lumped parameter set allows drainage network characterization and surface runoff routing to the site of interest for hydrologic modelling purposes. To this aim, four lumped parameters are defined. Two parameters are used to identify the channelled portion of the drainage network. The drainage density redistribution in the mountain and valley parts of the basin is related to k; its value was derived, in a theoretical framework, from drainage network fractal properties. The quantity AS k is calibrated on the basis of catchment drainage density and hydrologic response characteristics. Surface runoff routing is then performed with reference to calibrated kinematic parameters that maintain a physical meaning as average velocities, vh and vc , on the two drainage components: hillslopes and channels. The modelÕs properties outlined above are of foremost importance when the model is used for flood discharge prediction in non-gauged sites, where parameters direct estimation from recorded events is not possible. In particular, the capacity to treat different, measured or predicted, distributed rainfall inputs, the invariance of the parameters used to identify the channelled portion of the drainage and the strong physical meaning of the kinematic parameters, allows the use of the model for flood prediction on the basis of available quantitative precipitation forecasts. Results obtained from the calibration and validation phases show that the parameters used to characterize the drainage network can be assumed constant for the entire basin set here analysed, i.e. from 100 to 1000 km2 for basins in the Liguria region to 8000 km2 for the Tanaro basin. On the other side, values of the kinematic parameters appear to be sensitive to the average drainage slope. With respect to this, basins pertaining to the Liguria environment can be assumed to be homoge-
F. Giannoni et al. / Physics and Chemistry of the Earth 28 (2003) 289–295
neous and characterized by behaviour similar to that of sub-basins pertaining to the mountain part of the Tanaro basin. The Tanaro basin is not homogeneous; its sub-basins are placed in environments––mountain, mild and alluvial––characterized, from a hydrologic modelling point of view, by different kinematic parameters. In this case, the physical meaning assigned to the parameters is of paramount importance: the slow decreasing of the two average velocities moving from mountain towards alluvial catchments is a strong indication of the coherence of the calibration phase. Acknowledgements The research was partially supported by the Italian National Research Council, through GNDCI, the Italian National Group for the Protection from Hydrological Hazards. References Band, L.E., 1986. Topographic partition of watershed with digital elevation models. Water Resour. Res. 22 (1), 15–24. Flint, J.J., 1974. Stream gradient as a function of order, magnitude and discharge. Water Resour. Res. 10 (5), 969–973. Giannoni, F., Roth, G., Rudari, R., 2000. A semi-distributed rainfall– runoff model based on a geomorphologic approach. Phys. Chem. Earth 25 (7–8), 665–671. La Barbera, P., Roth, G., 1994. Invariance and scaling properties in the distributions of contributing area and energy in drainage basins. Hydrol. Proc. 8, 125–135. Mancini, M., Troch, P.A., Wood, E.F., 1994. Overland flow routing over a range of catchment scales. In: Rosso, R., Peano, A., Becchi,
295
I., Bemporad, G.A. (Eds.), Advances in Distributed Hydrology. Water Resource Publication, Highlands Ranch, CO, pp. 225–244. Mandelbrot, B.B., 1983. The Fractal Geometry of Nature. Freeman, New York. McCuen, R.H., 1982. A Guide to Hydrologic Analysis using SCS Methods. Prentice-Hall, Englewood Cliffs, NY. Montgomery, D.R., Dietrich, W.E., 1992. Channel initiation and the problem of landscape scale. Science 225, 826–830. OÕCallaghan, J.F., Mark, D.M., 1985. The extraction of drainage networks from digital elevation data. Comput. Vision Graph. Image Process 28, 323–344. Rodriguez-Iturbe, I., Valdes, J.B., 1979. The geomorphologic structure of hydrologic response. Water Resour. Res. 15 (6), 1409–1420. Rodriguez-Iturbe, I., Rinaldo, A., Rigon, R., Bras, R.L., Marani, A., Ijjasz-Vasquez, E., 1992. Energy dissipation, runoff production, and the three-dimensional structure of river basins. Water Resour. Res. 28 (4), 1095–1103. Rosso, R., Bacchi, B., Barbera, P., 1991. Fractal relation of mainstream length to catchment area in river networks. Water Resour. Res. 27 (3), 381–387. Roth, G., Barbera, P., Greco, M., 1996. On the description of the basin drainage structure. J. Hydrol. 187, 119–135. Shamseldin, A.Y., Nash, J.E., 1998. The geomorphological unit hydrograph––a critical review. Hess 2 (1), 1–8. Sherman, L.K., 1932. Streamflow from rainfall by the unit hydrograph method. Eng. News Rec. 108, 501–505. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1988. The fractal nature of river networks. Water Resour. Res. 24 (8), 1317–1322. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1989a. Scaling and elevation in river networks. Water Resour. Res. 25 (9), 2037–2051. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1989b. The analysis of river basins and channel network using digital terrain data. Technical Report 326, Ralph. M. Parsons Laboratory, MIT, Cambridge, MA. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1990. Comment on: On the fractal dimension of stream networks by P. La Barbera and R. Rosso. Water Resour Res. 6 (8), 2243–2244. Tarboton, D.G., Bras, R.L., Rodriguez-Iturbe, I., 1992. A physical basis for drainage density. Geomorphology 5, 59–76.