Simulation of Summer Hourly Stream Flow by Applying TOPMODEL and Two Routing Algorithms to the Sparsely Gauged Lhasa River Basin in China

This paper develops a new routing algorithm for improving simulation capacity of physically-based hydrological models applied to sparsely-gauged river basins. The study area is the Lhasa River basin, a large plateau basin with an area of 26,225 km2 in southwest China. In the basin, observations from three hydrological stations are available, and the observed hourly rainfall and summer stream flow data (12 June to 30 September) for the period of 1998–2000 obtained from the three stations were used. TOPMODEL, with its original routing algorithm, which is a distance-related delay function, was applied to the Lhasa River basin. To improve the routing algorithm using a unit hydrograph function and a linear reservoir method, this study proposed a new algorithm; the results revealed that the new algorithm improved the simulation of the variation of hourly stream flow. In addition, to evaluate the influence of rainfall spatial variation on runoff generation, observed rainfall series from the three gauges were used to simulate runoff individually, and it was found that there are significant differences among the three simulated hourly stream flow series.


Introduction
TOPographic Model (TOPMODEL) is a rainfall-runoff model developed using the simple theory of hydrological similarity, with the topographic index designed to represent hydrological similarity [1][2][3][4][5].The original routing algorithm for TOPMODEL simulates the stream flow by transferring the runoff to the basin outlet using a time-delay function.Table 1 presents some example applications of TOPMODEL to different river basins around the world.These applications generally indicate that TOPMODEL has mostly been applied to small river basins or data-rich river basins.In regards to sparsely-gauged large river basins, the exploration of the available routing algorithms for stream flow simulation is particularly lacking.To address this issue, the present study analyses the simulation capacity of TOPMODEL in a sparsely-gauged river basin.The Lhasa River basin in southwestern China is considered as the study area, and summer hourly stream flow in the basin is simulated.A new routing algorithm for TOPMODEL is proposed, and the results from the original and the new routing algorithms in simulating summer hourly stream flow are compared.

Study Area and Related Literature
Located in southwestern China, the Lhasa River (Figure 1) is the longest tributary of the Yarlung Zangbo River, with a mainstream length of approximately 551 km.The river originates at Nyangqentanglha Mountain on the Qinghai-Tibet Plateau [15].The Lhasa River basin has an area of 26,225 km 2 , and its elevation ranges from 3580 to 7000 m, with an average gradient of 2.9%.In the Lhasa River basin, stream flow, and weather data are recorded at three hydrological stations located at Pondo, Tanggya, and Lhasa (see Figure 1 and Table 2).The Lhasa station is located near the river outlet.At this station, the average annual precipitation is 427 mm/year, the average annual temperature is 7.7 °C, and the recorded highest and lowest temperatures are 29.5 °C and −18.0 °C, respectively.Rainfall variability during any year is significant, and its distribution is also uneven, with much of the rainfall concentrating during the summer period of June-September.In this basin, snowmelt mainly occurs during April and May [16,17].Several studies have been conducted to investigate the hydrological processes and water resources of the Qinghai-Tibet Plateau region.Wu et al. [18] applied three methods (the artificial neural network model, the staging calm auto-regressive model and the seasonal auto-regressive model) to a 41-year (1956 to 1968 and 1973 to 2000) average of monthly stream flow data in the Lhasa River and found the artificial neural network model to perform the best at simulating stream flow.In the Niyang River (one of the major tributaries of the Yarlung Zangbo River), a flood forecasting model was developed based on the black-box method using a mixed linear regression model [19].The model was calibrated using the hydrological data of 1997-2002 during the flood season, and verified using the hydrological data observed in 2003.The results showed that the model could be used for flood forecasting.Furthermore, Zhang et al. [20] used CC (cascade correlation) and BP (back propagation) models to forecast the monthly stream flow in the Lhasa River.Such studies have mainly focused on statistical analysis of hydrological and meteorological time series, and there is a general lack of application of physically-based hydrological models to the Yarlung Zangbo River basin, especially the Lhasa River basin.

Research Data
In the Lhasa River basin, hydrological data measurements are very limited.There are only three hydrological stations, namely Pondo, Tanggya, and Lhasa (see Figure 1 and Table 2).In this study, hourly precipitation and stream flow and daily pan evaporation data observed during the summer (June-September) for the period of 1998-2000 are considered.To investigate the skill of the model in forecasting flood flows in the Lhasa River, this study focuses on the simulation of summer stream flow from 12 June to 30 September, which is the river flood season.Similar to the study of Chen and Kumar [21], digital elevation model (DEM) data with a resolution of 1 km in the Lhasa River basin are extracted from a global DEM dataset (Global Land One-kilometre Base Elevation, National Geophysical Data Center, NOAA), and then the topographic index for TOPMODEL is derived.
It should be noted that the TOPMODEL version used in this study does not include a snowmelt module because the study period is summer (i.e., from 12 June to 30 September), as snow melt occurs only during April and May in this region.To evaluate the influence of snow melt in runoff generation, we investigate the variation of the snow coverage in the study area in the summer.Similar to Chen and Sun [22], the available MODIS (moderate resolution imaging spectroradiometer) snow cover data (approximately 160 remote sensing images obtained from ftp://n4ftl01u.ecs.nasa.gov)from June to September during the period of 2000-2010 are analyzed.It is worth noting that the MODIS data became available in 2000.From the analysis of the MODIS data, we find that the average snow cover is approximately 6% of the study area from June to September during the period of 2000-2010, and the inter-annual variation of the percentage of snow cover is small.Therefore, it can be deduced that rainfall in summer is a dominant factor in generating summer stream flow in the Lhasa River basin, and for simplicity, only rainfall data are used to simulate stream flow.

TOPographic Model (TOPMODEL)
TOPographic Model (TOPMODEL) adopts a topographic index to represent a dynamic saturated area of a basin [1,3,4].For a grid cell of DEM data, i, its topographic index, TIi, is calculated as follows: ln tan where ai is the area of the hillslope per unit contour length, and tanβi is the local surface slope at grid cell i.
In TOPMODEL, the surface runoff at grid i at time t is calculated as follows [2,23]: where uz S is the storage in the unsaturated zone and SD is the local storage deficit.
In the above equation, Sbar is the average storage deficit over the whole basin, zm S is the maximum moisture deficit in the unsaturated zone, * 1 is an average of the topographic indices in the basin, and A is the basin area.Using the topographic index and other TOPMODEL parameters, the base flow, qb, can be calculated as follows [2,23]: where q0 is equal to * λ 0 AT e  and 0 T is the lateral transmissivity when the soil column is saturated.

Original Routing Algorithm
The original routing algorithm for TOPMODEL [1] consists of two parts.The first part is to compute the travel time of runoff from any grid cell to a sub-basin outlet, and the second part is to compute the travel time of runoff in channels from upstream to downstream.For the first part, a flow delay function is normally represented by using a distance-related method.The time taken for runoff to reach a sub-basin outlet from any grid cell is calculated by the equation given below [23]: where j t is the time delay (in units of hours; noted as h hereafter) for runoff from the grid cell j , i x is the flow path length (m), ν is the flow velocity (m/h), and M is the number of the segments between the grid cell j and the sub-basin outlet.For the second part (i.e., the channel routing component), a simple nonlinear algorithm is normally used to compute t c , an average kinematic wave velocity for a river channel: where t Q is the discharge at time t , and CHA and CHB are coefficients, which can be calibrated.The travel time of runoff in a channel from an upstream gauge to a downstream gauge can be calculated using the length of the channel and t c .This approach is an explicit approximation to a kinematic wave channel routing algorithm [23].Most applications of this original routing algorithm have been based on a simple constant wave speed routing algorithm.

New Routing Algorithm
For a sparsely-gauged river basin, a sub-basin for computing the first part of runoff travel time using Equation ( 5) is usually larger, and the assumption of a constant flow velocity v is not reliable; in addition, the channel length from an upstream gauge to a downstream gauge is rather long, and the assumption of a constant wave speed for routing algorithms is not rational.To relax these two assumptions for sparsely-gauged river basins, this paper adopts the UH method to route the overland flow [24]: where ) (k UH is a unit hydrograph, k is the interval index, and m is the number of intervals with non-zero effective rainfall.There is a relationship between UH and Nash's instantaneous unit hydrograph [24], given by where ) (k IUH is Nash's instantaneous unit hydrograph, and t  is the time interval.The mathematical form for Nash's instantaneous unit hydrograph is equivalent to the gamma distribution [25]: where N is the number of linear reservoirs in a series, w is the storage coefficient, and for integer value of N ).According to Nash's instantaneous unit hydrograph and Equation ( 8), a UH can be obtained, and then the overland flow can be computed by substituting UH into Equation (7).Furthermore, the linear reservoir method is used to route the base flow, qb, which is computed by Equation ( 4), and the routed base flow, Qg, is calculated as follows: where g K is the groundwater recession coefficient.Then, the total routed runoff can be computed by the equation below: To use the new routing algorithm, the parameters of N , w and g K require calibration.

Criteria for Evaluating Model Performance
The first criterion used in this study to evaluate the model performance is the Nash-Sutcliffe efficiency index [26], which is defined as follows: where o t Q is the observed stream flow, s t Q is the simulated stream flow, and o Q is the mean value of the observed stream flow.The value of NE is equal to unity if the simulation matches the observation perfectly.
The second efficiency criterion used is the relative error of the volumetric fitness between the observed and simulated stream flow series and is defined by: The value of E R is zero if the total volume of the simulated stream flow is equal to that of the observed ones.
The third criterion is the correlation coefficient, r, to evaluate the match of variations between observations and simulations.The equation for computing r is given as where s Q is the mean value of the simulated stream flow.For model calibration in this study, the Simplex algorithm [27] is used.
Because the simulation is conducted from 12 June to 30 September for the period of data considered (1998-2000), rainfall serves as the driving force for floods, as snowmelt only occurs during April and May.In addition, in the study area, because the flood in the summer of 1998 was the largest flood during the study period, the observed stream flows in 1998 are chosen for verification purposes to better assess the skill of the model in stream flow simulation, especially regarding the flood flows.The floods that occurred during the summer of 1999-2000 are used for calibration.It should be noted that the time interval used in this study is one hour.

Results and Discussion
Because there are only three hydrological stations available in the Lhasa River basin (i.e., the Pondo, Tanggya and Lhasa stations) despite its large area, the basin is separated into three parts in association with these stations.Accordingly, the topography indices for these three sub-basins are calculated.In the new routing algorithm, the values of new three parameters (N, w and Kg) are 16, 3.3 and 0.00033, respectively.
Table 3 lists the stream flow simulation results obtained for the Lhasa station using TOPMODEL and the two routing algorithms mentioned above.Here, Routing 1 and Routing 2 refer to the original and new routing algorithms, respectively.And the average hourly precipitation data at three hydrological stations (namely Pondo, Tanggya, and Lhasa) was used in the simulations.Table 3 shows that the Nash-Sutcliffe efficiencies, NE, obtained using both the routing algorithms are between 0.70 and 0.78, the relative errors of the volumetric fitness, E R , are between −12% and 9%, and the correlation coefficients, r, are between 0.85 and 0.91.Comparison of NE, RE and r reveals that the performances of TOPMODEL with the original and new routing algorithms are similar.Table 4 lists the observation and simulation results of several flood peaks that occurred at the Lhasa station from June to September in the study period.The relative errors between the observed and simulated flood peaks are mostly in the range of −22% to 19%, and the performances of both methods are also similar.Figure 2 shows the observed hourly stream flow and simulated hourly base flows and stream flows for the Lhasa station in the year 1998.Both the original and new routing algorithms can simulate the rising and recessing limbs of the stream flows reasonably well in comparison to the observed flood hydrograph.The new algorithm is insensitive to input, therefore the simulated base flow using the new routing algorithm is smoother than that using the original routing algorithm.Due to sparse data currently in the basin, it is not easy to determine quantity and proportion of base flow in the stream flow.For the further research, isotop experiments [28,29] could be tested in the stream flow, and the actual amount of base flow would remain clear.
The correlation coefficient results obtained for stream flows from 12 June to 30 September for the study period (see Table 3) do not offer a strong indication as to the better routing algorithm (i.e., the new routing algorithm or the original one).However, the use of correlation coefficients can be better presented by examining the stable or rapid variations of stream flow simulated by the two routing algorithms.Table 5 lists the correlation coefficients for flood events in each year, for a total of 11 events over the three-year study period (1998)(1999)(2000).The results show that the correlation coefficients of simulated stream flow from the new routing algorithm are mostly larger than those from the original routing algorithm; specifically, for a flood hydrograph that occurred during 15-17 August 1998, values of 0.59 and 0.71 were found for the original and new algorithms, respectively.Therefore, intuitively, the linear reservoir and the UH method used in the new routing algorithm are more reliable for the Lhasa River basin when compared to the original routing algorithm.Because the observed rainfall at the three hydrological stations (Pondo, Tanggya, and Lhasa) can not sufficiently represent the spatial variation of precipitation over the Lhasa River basin, a sensitivity analysis of the influence of rainfall spatial variation on hourly stream flow simulation is also performed.To this end, four rainfall time series, corresponding to the rainfall recorded at Pondo, Tanggya, Lhasa and their average, are used to represent the rainfall of the whole study area.Figure 3 shows the simulation results of the 1998 summer stream flow at the Lhasa station obtained by using these four rainfall series.The results obtained using the original routing algorithm (Figure 3a) show that a significant difference exists between the observed and the simulated hourly stream flow at the Lhasa station from June to September, 1998.Additionally, the simulated stream flow using the rainfall recorded at the Pondo station significantly overestimates most flood peaks, such as the floods that occurred on 26 July and 7 September.The results further indicate that, among the four rainfall time series considered, the average rainfall of the three stations provide the best stream flow simulation.Compared with the simulations presented in Figure 3a, the simulated stream flows in Figure 3b are more stable.These results are similar to those presented in Figure 2 and Table 5. Figure 3b shows that the simulated stream flows at the outlet of the basin using the rainfall recorded at the Tanggya station based on TOPMODEL are mostly comparable to the observed stream flows.In fact, because the stream flow at the Lhasa station is more closely related to that of the Tanggya station, the Tibetan Bureau of Hydrology (TBH) [30] has developed a correlation diagram method for the Lhasa River on the basis of experience.This diagram has also been used to predict the water level at the Lhasa station, according to the water level relationship between the Lhasa station and the Tanggya station.The simulation results in Figure 3b confirm such a relationship developed by TBH.

Conclusions
Although TOPMODEL has been widely applied to different river basins around the world, such applications have mostly been focused on small basins and data-rich large basins, and there is a lack of application to sparsely-gauged large river basins.In this study, a new routing algorithm was developed for routing the runoff for a relatively longer channel from the upstream station to the downstream one, and both the original and new routing algorithms were used in TOPMODEL in its application to the Lhasa River basin, a sparsely-gauged large basin.Hourly summer stream flow from June to September was simulated.The results for major summer flood peaks at the Lhasa station for the period from 1998 to 2000 were analyzed.In addition, the influence of rainfall spatial variations on stream flow simulation was investigated.The major conclusions drawn from this study are as follows: 1. Application of the original and new routing algorithms to the sparsely gauged Lhasa River basin reveals that the new routing algorithm is suitable for the study area.Furthermore, comparison of the simulated base flows obtained from the original and new routing algorithms showed that the new algorithm can provide a smoother hourly base flow than the original routing algorithm.
2. The simulated stream flows from both the original and new routing algorithms indicate that the three stations available (Pondo, Tanggya, and Lhasa) are not sufficient to represent the rainfall spatial variation in the study area.To evaluate the influence of the spatial variation of rainfall, the mean areal rainfall time series was also used to calibrate the model parameters and the corresponding simulated stream flow from the model were considered as a reference.The rainfall recorded at Pondo, Tanggya, and Lhasa hydrological stations were used as inputs to the model with the mean areal rainfall model parameters, respectively.The simulated stream flow was compared with mean areal rainfall.The results reveal that the hourly stream flow at the outlet of the basin simulated by only using the rainfall recorded at Tanggya station are mostly comparable to the observed stream flow in the same location.
3. In this study, only rainfall data were used to drive TOPMODEL and its related routing algorithms.For a large plateau basin, such as the Lhasa River basin, the snow cover is mostly present throughout the year, and thus, snowmelt can contribute to runoff generation.To investigate the snow cover influence using the available MODIS snow cover data, we found that the average snow cover accounts for approximately 6% of the Lhasa River basin from June to September, and the inter-annual variation of snow cover in the summer is small.Therefore, we assumed that rainfall is a dominant force in generating summer stream flows and that snowmelt is not a particularly important factor.Nevertheless, it would be interesting to study the influence of snowmelt on stream flow simulation throughout the entire year.We will investigate this aspect in our future research.

Figure 1 .
Figure 1.(a) Topographic map; (b) the river network of the Lhasa River basin.

Figure 2 .
Figure 2. Observed and simulated hourly stream flow and base flow at Lhasa in 1998.

Figure
Figure3bpresents the stream flow simulated results obtained using the new routing algorithm.Compared with the simulations presented in Figure3a, the simulated stream flows in Figure3bare more stable.These results are similar to those presented in Figure2and Table5.

Figure 3 .
Figure 3. Observed and simulated hourly stream flow at Lhasa in 1998 by using (a) Routing 1; and (b) Routing 2 based on the different rainfall scenarios.

Table 1 .
Some applications of TOPMODEL to different basins.

Table 2 .
Information of hydrological stations in the Lhasa River basin.

Table 3 .
Simulation results by using TOPMODEL and the original and new routing algorithms.

Table 4 .
Simulation results for major flood peaks at Lhasa station.

Table 5 .
Correlation coefficients for flood events occurred at Lhasa station from 1998 to 2000.