## Abstract

This paper presents an energy model for determining the overland flow velocity in order to improve the low accuracy problem in flow concentration simulation. It furnishes a novel idea for studying flow concentration in ungauged basins. The model can be widely applied in analysis of spatial velocity field, extraction of instantaneous geomorphic unit hydrograph and development of distributed hydrological model. A distributed flood-forecasting model is constructed for Lianyuan Basin in Hunan Province of China. In the proposed method, gravitational potential energy is transformed into kinetic energy via an analysis of energy distribution of water particles in the basin. Based on the kinetic energy equation, the overland flow velocity simulating the geomorphic unit hydrograph is computed. Rainfall-runoff simulation is then performed by integrating with runoff yield and concentration model. Results indicate that the model based on energy conversion leads to more accurate results. The model has the following advantages: firstly, the spatial distribution of the velocity field is appropriate; secondly, the model has only one parameter, which is easily determined; and finally, flow velocity results can be used for the computation of river network flow concentration.

## HIGHLIGHTS

Propose a geomorphic unit hydrograph based on the principle of energy conversion.

Develop a distributed flood-forecasting model using the proposed geomorphic unit hydrograph.

The process of parameter determination of the developed model is very simple.

The velocity distribution obtained from the analysis is reasonable.

The developed hydrological model has a significant improvement in flood forecasting accuracy.

## INTRODUCTION

Floods of medium and small catchments are among the most costly natural disasters all around the world, posing serious risk of loss of life, physical injury, damage to infrastructure and disruption of economic and social activities (Shi *et al.* 2016). Last decades showed a steady increase of damages due to flooding, which highlighted the need of developing effective measures to reduce flood event impacts (Barbetta *et al.* 2018). Conceptual rainfall-runoff (CRR) models have become a basic tool for flood forecasting and become increasingly important for catchment basin management (Xu *et al.* 2013). However, accurate flood prediction in mountainous watersheds remains a significant challenge for flood-forecasting models. In order to overcome this challenge, a good knowledge of the relationship between rainfall and runoff is entailed. A unit hydrograph is a common method in flood estimation, which is not only applied in peak flow estimation but also in the creation of complicated flood hydrographs (Khaleghi *et al.* 2011) and has been widely used in rainfall-runoff computation since it was proposed by Sherman (1932). Clark (1945) combined the two concepts of the isochrone method and linear reservoir to establish the instantaneous unit hydrograph. Although the Clark unit hydrograph can be a very valuable technique in flood hydrology, the lack of appropriate techniques to estimate the storage coefficient-*R* for ungagged watersheds has diminished the application and utility of this technique (Sabol 1988).

Obviously, there is a close relationship between the hydrological response of a river basin and its geomorphologic characteristics. Since Rodríguez-Iturbe & Valdés (1979) established the geomorphic instantaneous unit hydrograph theory (GIUH) by combining the initial probability of raindrop falling with the period transfer probability, there have been many attempts to propose an Instantaneous Unit Hydrograph that incorporated the geomorphological properties of the watershed (Gupta *et al.* 1980; Rodríguez-Iturbe *et al.* 1982; Agirre *et al.* 2005; Goñi *et al.* 2019). Agnese *et al.* (1988) derived the time scale of a GUIH from an effective stream flow velocity of the highest-order stream and the spatial distribution of velocity throughout the stream network. van der Tak & Bras (1990) incorporated hillslope effects into the gamma GIUH model by assuming that the hillslope travel distance in an area of a given order is approximated by the inverse of twice the local drainage density and introducing a hillslope velocity term. Cudennec *et al.* (2004) investigated the geomorphologic aspect of the unit hydrograph concept and concluded that the use of geomorphologic parameters explained the unit hydrograph and geomorphologic unit hydrograph theories.

Rosso (1984) combined the Nash model with geomorphic parameters based on the Horton Strahler river classification method and analyzed the GIUH by fitting the probability density distribution function of the flow concentration time. This method required computing the average velocity of the basin flow concentration. For areas without sufficient data, the empirical formula proposed by Kirshen & Bras (1983) could be used to estimate the average convergence speed of a basin with the basin area, average rainfall intensity, river gradient, average river width and roughness. In order to reveal the spatially distributed nature of watershed properties, Maidment *et al.* (1996) proposed the concept of distributed spatial velocity field and analyzed the unit hydrograph by computing the spatial distribution of velocity. Tan *et al.* (2018) and Munoth & Goyal (2019) pointed out that when a digital elevation model (DEM) was used to divide sub-basins, DEM resolution and watershed threshold would have a greater impact on geomorphic parameters. Therefore, when using this method to build geomorphic unit hydrograph, it was necessary to analyze and demonstrate the applicability of DEM resolution and watershed division threshold. Kong & Guo (2019) studied the influence of rain intensity on velocity based on the Soil Conservation Service (SCS) velocity formula and proposed a time-varying distributed unit hydrograph computation method. These methods were based on geomorphic parameters.

With the continuous development and improvement, the construction of a geomorphic unit hydrograph based on the spatial velocity field has become a practical technology suitable for the flow concentration computation in ungauged basins. However, it is difficult to get a reasonable velocity distribution if only a single factor is considered. When considering multiple factors to estimate the velocity, it suffers from the limitation that the influence weight and parameter value of each factor are determined. In a flat area, there will be many grids with a slope of 0, which need to be treated with a special method. In fact, the transformation from gravitational potential energy to kinetic energy is the root cause of water flow. In the process of transformation, energy will be lost due to terrain, soil, vegetation, river diversion, water conservancy engineering and other factors.

The objective of this paper is to develop a method of constructing GIUH with clear physical meaning and study a method of determining the value of parameters. Based on this method, the distributed flood-forecasting model of Lianyuan Basin in China is established. In the proposed hydrological model, the spatial velocity field and the instantaneous unit hydrograph can be extracted based on the principle of energy conversion. It has many advantages as mentioned below. Firstly, only one parameter of the model needs to be determined and the process of parameter determination is very simple. Secondly, the error caused by the gradient of 0 in the flat area is avoided and the velocity distribution obtained from the analysis is reasonable. Finally, the velocity results can be used for the simulation of both slope and river network flow concentration, which can ensure the consistency of flow concentration computation to a certain extent and reduce the influence of watershed threshold on model parameters when subwatershed is divided.

The paper is organized as follows: Section ‘Methodology’ describes the methodology proposed for the area rainfall calculation and runoff generation model and presents the flow concentration model in this study. Section ‘Case study’ gives details on the selected case study and model evaluation criteria. Section ‘Results’ presents and discusses the results obtained through the proposed methodology. Finally, we provide the conclusions of this work in the last section.

## METHODOLOGY

### Areal rainfall computation

The areal rainfall is computed by an inverse distance square (IDS) model. It can be approximately considered that the rainfall distribution in the basin is uniform because the area of sub-basin is very small. The result of the interpolation of the central point of the subwatershed is taken as the areal rainfall of the subwatershed.

*Z*is the estimated value of the interpolation point,

*Z*is the measured sample value,

_{i}*n*is the number of measured samples involved in the computation,

*D*is the distance between the interpolation point and the

_{i}*i*th station,

*p*is the power of the distance, which significantly affects the interpolation results, and its selection standard is the minimum average absolute error. When

*p*is taken as 0, the weights of all sample points are equal, and the method degenerates to an arithmetic mean model. When

*p*is infinite, the result is approximately equal to the value of the nearest sample point, and the method is equivalent to the Voronoi model. In practical application,

*p*is usually taken as 2 under the IDS model, which is the method adopted in this paper.

### Runoff generation model

The Xinanjiang model is adopted for runoff generation computation. It is a CRR model with distributed parameters (Wang *et al.* 2012), which was developed by Zhao *et al.* (1980) and has been successfully and widely applied to the humid and semi-humid regions of China for flood forecasting (Zhao 1992). The model computes the basin evaporation and emission according to a three-layer evaporation and emission model. The total runoff produced by rainfall is computed according to the concept of full storage and runoff yield, and the influence of uneven underlying surface on the runoff yield area is considered by the water storage curve of the basin. In the aspect of runoff component division, according to the runoff production theory of ‘hillside hydrology’, the total runoff is divided into saturated surface runoff, soil water runoff and groundwater runoff by a free water storage reservoir with limited volume, measuring hole and bottom hole. For detailed descriptions and explanations of the Xinanjiang model, refer to Zhao (1992). The model structure is shown in Figure 1. It has 13 parameters and their physical descriptions are listed in Table 1. The value of each parameter is usually within certain ranges according to physical and mathematical constraints, information about watershed characteristics and from modelling experiences (Wang *et al.* 2012).

. | Parameter . | Physical description . | Unit . |
---|---|---|---|

Runoff generating parameter | |||

1 | K | Ratio of potential evapotranspiration to pan evaporation | – |

2 | U_{m} | Averaged soil moisture storage capacity of the upper layer | mm |

3 | L_{m} | Averaged soil moisture storage capacity of the lower layer | mm |

4 | D_{m} | Averaged soil moisture storage capacity of the deep layer | mm |

5 | C | Coefficient of the deep layer that depends on the proportion of the basin area covered by vegetation with deep roots | – |

6 | B | Exponential parameter with a single parabolic curve, which represents the nonuniformity of the spatial distribution of the soil moisture storage capacity over the catchment | – |

7 | I_{m} | Percentage of impervious and saturated areas in the catchment | – |

Runoff routing parameter | |||

8 | S_{m} | Areal mean free water capacity of the surface soil layer, which represents the maximum possible deficit of free water storage | mm |

9 | E _{x} | Exponent of the free water capacity curve influencing the development of the saturated area | – |

10 | K_{g} | Outflow coefficients of the free water storage to groundwater relationships | – |

11 | K_{i} | Outflow coefficients of the free water storage to interflow relationships | – |

12 | C_{i} | Recession constants of the lower interflow storage | – |

13 | C_{g} | Recession constants of the groundwater storage | – |

. | Parameter . | Physical description . | Unit . |
---|---|---|---|

Runoff generating parameter | |||

1 | K | Ratio of potential evapotranspiration to pan evaporation | – |

2 | U_{m} | Averaged soil moisture storage capacity of the upper layer | mm |

3 | L_{m} | Averaged soil moisture storage capacity of the lower layer | mm |

4 | D_{m} | Averaged soil moisture storage capacity of the deep layer | mm |

5 | C | Coefficient of the deep layer that depends on the proportion of the basin area covered by vegetation with deep roots | – |

6 | B | Exponential parameter with a single parabolic curve, which represents the nonuniformity of the spatial distribution of the soil moisture storage capacity over the catchment | – |

7 | I_{m} | Percentage of impervious and saturated areas in the catchment | – |

Runoff routing parameter | |||

8 | S_{m} | Areal mean free water capacity of the surface soil layer, which represents the maximum possible deficit of free water storage | mm |

9 | E _{x} | Exponent of the free water capacity curve influencing the development of the saturated area | – |

10 | K_{g} | Outflow coefficients of the free water storage to groundwater relationships | – |

11 | K_{i} | Outflow coefficients of the free water storage to interflow relationships | – |

12 | C_{i} | Recession constants of the lower interflow storage | – |

13 | C_{g} | Recession constants of the groundwater storage | – |

### Flow concentration model

The concentration model is composed of two parts: the overland flow concentration and the river network flow concentration. The overland flow concentration denotes the process of net rain gathering along the slope to the river channel, and the river network concentration means that the flow converges to the downstream along all levels of river channels. The flow concentration time mentioned below denotes the time required for a water particle at any position in the basin to flow to the outlet section of the basin.

#### Overland flow concentration

According to D8 algorithm (Fairfield & Leymarie 1991), a water particle in a grid is assumed to flow to the adjacent grid in the direction of the maximum gradient in DEM. Figure 2(a) shows the elevation of each grid. Figure 2(b) displays the flow direction of each grid, in which water particles in any grid flow to its lowest adjacent grid, and thus, the flow concentration routes of each grid can be determined. The flow routes from a grid along the flow direction to the outlet of the basin are illustrated in Figure 2(c).

*L*is the length of the edge of the mesh converted to the actual distance and

*v*is the flow velocity in the grid. Along the flow concentration route, the flow concentration time from each grid to the outlet of the drainage basin can be computed by the following formula:where

*m*is the number of grids on the flow concentration route. Since the area of each grid is known, the drainage area corresponding to different flow concentration times can be computed, so as to obtain the flow concentration time area relationship and compute the unit hydrograph according to the method in Kong

*et al.*(2007).

*et al.*1996; Kang

*et al.*2006; Song

*et al.*2010):where

*S*is the slope of two adjacent grid points,

*A*is the grid catchment area,

*h*is the net rainfall of the unit grid,

*i*is the net rainfall intensity and

*a*,

*b*,

*c*,

*d*are some parameters, where parameter

*a*summarizes the combined effect of all other factors not described in the equations.

*E*

_{ideal}is the total kinetic energy of flow in the grid under the ideal conditions;

*n*is the total number of grids in the basin upstream of the target grid (including the target grid);

*E*

_{p}is the gravitational potential energy of each upstream grid relative to the target grid;

*m*is the mass of unit net rain in the grid;

*g*is the gravity acceleration and

*H*

_{Avg}is the average elevation difference between all upstream grids and the target grid.

*μ*. According to the kinetic energy formula, the following formula is derived:where

*E*

_{k}is the kinetic energy of the flow from the target grid. The following equation is obtained by combining Equations (9) and (10):

*μ′*is the energy residual coefficient,

*θ*is the slope angle of the grid outflow direction and Δ

*h*is the elevation difference between the target grid and the outflow grid. For grids that are not on the watershed, the kinetic energy of water flow comes from two parts, one is the kinetic energy of the inflow grids, and the other is the gravitational potential energy difference between the target grid and the outflow grid. The energy equation is as follows:where

*N*is the number of inflow grids of the target grid. The following direct equation of

*v*can be obtained by combining Equations (10) and (13):

#### River network flow concentration

*I*

_{t}_{−1}and

*Q*

_{t}_{−1}are the inflow and outflow of the channel at the beginning of the simulation, respectively.

*I*and

_{t}*Q*are the inflow and outflow of the channel at the end of the simulation, respectively.

_{t}*C*

_{0},

*C*

_{1}and

*C*

_{2}are the parameters of the simulation and they can be determined by the following formula:

*dt*is the time step of the simulation. In order to determine the parameter

*K*, the wave velocity must be computed first. The spatial velocity field is used to compute the average velocity

*v*, and a flow equation (Todini 2007) is used to compute the wave velocity of flood

*C*.where

*λ*is the conversion coefficient. According to Tewolde & Smithers (2006),

*λ*for rectangular

*n*, parabola and triangle sections are taken as 5/3, 13/9 and 4/3, respectively.

*x*is given by Song et al

*.*(2011) where

*L*is the length of the reach,

*n*is the roughness,

*Q*is the average of the maximum and minimum flow values in the inflow process and

*J*is the gradient of the reach.

## CASE STUDY

### Study area and data

As shown in Figure 3, the study area is located in Hunan province of southern China and is the source of Lianshui, which is one of the tributary rivers in the Xiangjiang River. It originates in the south of Guanyin Mountain and the river flows through Lianyuan City from southwest to northeast. The catchment area of the basin is 160 km^{2} with an average slope of 25.2%. The longest flow concentration route is 25.7 km. The average elevation of the basin is 268 m, and the maximum drop is 658 m. The soil in the basin is mainly sand clay on the hillside and clay loam on the riverside. The vegetation distribution consists of 60% of farmland and 40% of forest. Lianyuan Basin is located in the mid subtropical continental monsoon humid climate area, with the annual average precipitation of 1377.1 mm and 45% of the total rainfall falls between April and June.

The study area is divided into 11 sub-basin units by DEM data of 25 m × 25 m resolution, including 256,000 grids. There are five rain gauging stations and one hydrological station in the study area. The rainfall weights of each rain gauging station computed by IDS for each sub-basin are listed in Table 2. A total of 29 historical floods between 1979 and 2006 are employed to calibrate the model parameters while 10 floods between 2010 and 2017 are utilized to verify these parameters. The time step in the model simulation is half an hour.

Sub-basin . | Weight of rainfall station . | |||||
---|---|---|---|---|---|---|

. | LianYuan . | XiaShanQiao . | XieJia . | ChangRong . | TuanJia . | TuZhu . |

LianYuan2_R120W120 | 0.881 | 0.005 | 0.013 | 0.041 | 0.051 | 0.009 |

LianYuan2_R140W140 | 0.039 | 0.048 | 0.058 | 0.182 | 0.055 | 0.618 |

LianYuan2_R150W150 | 0.349 | 0.018 | 0.052 | 0.299 | 0.248 | 0.034 |

LianYuan2_R160W160 | 0.015 | 0.007 | 0.018 | 0.913 | 0.030 | 0.019 |

LianYuan2_R170W170 | 0.003 | 0.001 | 0.004 | 0.007 | 0.984 | 0.001 |

LianYuan2_R180W180 | 0.037 | 0.036 | 0.091 | 0.642 | 0.090 | 0.104 |

LianYuan2_R190W190 | 0.007 | 0.014 | 0.907 | 0.029 | 0.031 | 0.011 |

LianYuan2_R290W290 | 0.034 | 0.136 | 0.132 | 0.164 | 0.068 | 0.466 |

LianYuan2_R300W300 | 0.031 | 0.270 | 0.368 | 0.123 | 0.081 | 0.128 |

LianYuan2_R320W320 | 0.019 | 0.748 | 0.097 | 0.042 | 0.037 | 0.058 |

LianYuan2_R330W330 | 0.025 | 0.572 | 0.217 | 0.061 | 0.058 | 0.066 |

Sub-basin . | Weight of rainfall station . | |||||
---|---|---|---|---|---|---|

. | LianYuan . | XiaShanQiao . | XieJia . | ChangRong . | TuanJia . | TuZhu . |

LianYuan2_R120W120 | 0.881 | 0.005 | 0.013 | 0.041 | 0.051 | 0.009 |

LianYuan2_R140W140 | 0.039 | 0.048 | 0.058 | 0.182 | 0.055 | 0.618 |

LianYuan2_R150W150 | 0.349 | 0.018 | 0.052 | 0.299 | 0.248 | 0.034 |

LianYuan2_R160W160 | 0.015 | 0.007 | 0.018 | 0.913 | 0.030 | 0.019 |

LianYuan2_R170W170 | 0.003 | 0.001 | 0.004 | 0.007 | 0.984 | 0.001 |

LianYuan2_R180W180 | 0.037 | 0.036 | 0.091 | 0.642 | 0.090 | 0.104 |

LianYuan2_R190W190 | 0.007 | 0.014 | 0.907 | 0.029 | 0.031 | 0.011 |

LianYuan2_R290W290 | 0.034 | 0.136 | 0.132 | 0.164 | 0.068 | 0.466 |

LianYuan2_R300W300 | 0.031 | 0.270 | 0.368 | 0.123 | 0.081 | 0.128 |

LianYuan2_R320W320 | 0.019 | 0.748 | 0.097 | 0.042 | 0.037 | 0.058 |

LianYuan2_R330W330 | 0.025 | 0.572 | 0.217 | 0.061 | 0.058 | 0.066 |

### Evaluation criteria

To evaluate the suitability of the proposed model for the studied basin, three criteria are chosen to analyze the degree of goodness of fit. These criteria are Root Mean Square Error (RMSE), Nash–Sutcliffe Coefficient of Efficiency (NSCE) and three statistical ratios of acceptable criteria relative to the peak discharge, peak time and total runoff volume, which are introduced as follows:

*et al.*2000; Cooper

*et al.*2007; Wang

*et al.*2012; Xu

*et al.*2013):where

*Q*

_{0}(

*i*) and

*Q*

_{s}(

*i*) are the measured and simulated runoff or corresponding runoff, respectively, and

*N*is the number of data points involved.

*n*represents the total number of floods that satisfy the acceptable criteria relative to the peak discharge, peak time and total runoff volume, respectively; and

*m*is the total number of the calibrated floods or validated ones. Three statistics (QR

_{peak_discharge}, QR

_{peak_tm}and QR

_{runoff}) are used to evaluate the effect of parameter calibration of a rainfall-runoff model (Wang

*et al.*2012; Xu

*et al.*2013).

The scheme is excellent when QR reaches 85%. The scheme is good when QR is greater than 75% and less than 85%. Otherwise, the results of the performances of parameter calibration are unsatisfactory for online flood forecasting.

## RESULTS

29 floods are selected for model parameter calibration from 1979 to 2006. Table 3 gives the results of parameter calibration. According to three statistical ratios of acceptable national criteria relative to the peak discharge, peak time and total runoff volume among the calibrated and validated historical flood events for flood forecasting in China, Table 4 presents the detailed results of statistics performance using calibrated parameters during the calibration for flood simulation. It indicates that the qualified quantity is 24 and the ratio of qualifying simulation is 82.76% relative to the error of peak discharge; the qualified quantity is 27 and the ratio of qualifying simulation is 93.1% relative to the error of peak time; and the qualified quantity is 25 and the ratio of qualifying simulation is 86.21% relative to the error of total runoff volume. Ten floods are used for validation from 2010 to 2017. Table 5 gives the detailed results of statistics performance using the calibrated parameters during the validation for flood forecasting. It shows that the qualified quantity is 9 and the ratio of qualifying simulation is 90.0% relative to the error of peak discharge; the qualified quantity is 10 and the ratio of qualifying simulation is 100% relative to the error of peak time; and the qualified quantity is 9 and the ratio of qualifying simulation is 90.0% relative to the error of total runoff volume.

Parameter | U_{m} (mm) | L_{m} (mm) | D_{m} (mm) | B | I_{m} | K | C | S_{m} (mm) |

Value | 10 | 40 | 15 | 0.15 | 0.01 | 0.9 | 0.15 | 30 |

Parameter | E _{x} | K_{g} | K_{i} | C_{g} | C_{i} | μ′ | ||

Value | 1.1 | 0.5 | 0.4 | 0.9 | 0.3 | 0.005 |

Parameter | U_{m} (mm) | L_{m} (mm) | D_{m} (mm) | B | I_{m} | K | C | S_{m} (mm) |

Value | 10 | 40 | 15 | 0.15 | 0.01 | 0.9 | 0.15 | 30 |

Parameter | E _{x} | K_{g} | K_{i} | C_{g} | C_{i} | μ′ | ||

Value | 1.1 | 0.5 | 0.4 | 0.9 | 0.3 | 0.005 |

Floods . | Observed (m^{3}/s)
. | Simulated (m^{3}/s)
. | Percentage error (%) . | Observed peak time (yyyy-mm-dd hh) . | Simulated peak time (yyyy-mm-dd hh) . | Error (number) . | Total volume error (%) . |
---|---|---|---|---|---|---|---|

19790627 | 115 | 91.4 | −20.5 | 1979-06-27 15:00 | 1979-06-27 16:30 | −1.5 | 15.7 |

19800812 | 108 | 103 | −8.3 | 1980-08-12 08:00 | 1980-08-12 17:00 | −9 | 12.8 |

19810407 | 114 | 109.3 | −4.1 | 1981-04-07 07:30 | 1981-04-07 09:00 | −1.5 | −4.5 |

19820616 | 170 | 151.2 | −11.1 | 1982-06-16 19:42 | 1982-06-16 17:30 | 2.2 | 7.1 |

19870512 | 109 | 149.2 | 36.9 | 1987-05-13 02:30 | 1987-05-13 03:00 | −0.5 | 25.6 |

19880829 | 99.9 | 95.3 | −4.6 | 1988-09-03 15:36 | 1988-09-03 15:30 | 0.1 | 18.1 |

19900607 | 108 | 128.5 | 19.0 | 1990-06-07 11:00 | 1990-06-07 09:00 | 2 | 16.3 |

19900615 | 185 | 154.9 | −16.3 | 1990-06-15 15:00 | 1990-06-12 23:00 | 64 | −5.6 |

19920322 | 130 | 133.1 | 2.4 | 1992-03-20 06:00 | 1992-03-20 06:00 | 0 | 1.4 |

19920516 | 172 | 158.4 | 11.3 | 1992-05-17 00:12 | 1992-05-16 23:30 | 0.7 | 11.5 |

19920615 | 109 | 102.2 | −6.2 | 1992-06-15 12:00 | 1992-06-15 13:00 | −1 | −5.2 |

19920622 | 259 | 250.9 | −3.1 | 1992-06-22 09:30 | 1992-06-22 09:30 | 0 | 1.9 |

19930704 | 162 | 134.6 | −17.5 | 1993-07-04 20:00 | 1993-07-04 20:30 | −0.5 | −28.1 |

19930720 | 170 | 263.6 | 55.1 | 1993-07-20 10:00 | 1993-07-20 08:00 | 2 | 49.2 |

19940425 | 98 | 82.7 | −15.6 | 1994-04-25 08:00 | 1994-04-25 08:00 | 0 | −6.0 |

19940718 | 206 | 184.3 | −10.5 | 1994-07-18 17:12 | 1994-07-18 17:30 | −0.3 | −4.4 |

19940905 | 129 | 178.2 | 38.1 | 1994-09-06 11:30 | 1994-09-06 11:30 | 0 | 13.7 |

19950701 | 423 | 401.5 | −5.1 | 1995-06-30 12:12 | 1995-06-30 12:00 | 0.2 | 2.6 |

19960601 | 108 | 108.8 | 4.5 | 1996-06-01 13:42 | 1996-06-01 15:30 | −1.8 | 0.0 |

19960717 | 140 | 128.1 | −17.4 | 1996-07-18 08:00 | 1996-07-18 08:30 | −0.5 | −17.9 |

19980522 | 293 | 309.6 | 5.7 | 1998-05-22 09:00 | 1998-05-22 09:30 | −0.5 | −0.3 |

19980624 | 169 | 144.8 | −14.3 | 1998-06-25 02:30 | 1998-06-25 03:30 | −1 | −9.8 |

19990716 | 97.6 | 94.3 | −3.4 | 1999-07-17 12:00 | 1999-07-17 10:00 | 2 | −23.0 |

20020818 | 140 | 115 | −17.9 | 2002-08-19 18:18 | 2002-08-19 16:00 | 2.3 | −17.0 |

20030516 | 109 | 92.9 | −14.8 | 2003-05-16 12:00 | 2003-05-16 12:30 | −0.5 | −14.3 |

20030605 | 81.5 | 56.8 | −30.3 | 2003-06-05 20:42 | 2003-06-05 22:00 | −1.3 | −19.7 |

20040515 | 206 | 186.7 | −9.4 | 2004-05-15 16:12 | 2004-05-15 16:30 | −0.3 | −4.6 |

20050601 | 183 | 161.2 | −11.9 | 2005-06-01 10:30 | 2005-06-01 11:30 | −1 | −13.2 |

20060411 | 108 | 118.8 | 10.0 | 2006-04-12 12:48 | 2006-04 − 12 13:30 | −0.7 | −2.3 |

Floods . | Observed (m^{3}/s)
. | Simulated (m^{3}/s)
. | Percentage error (%) . | Observed peak time (yyyy-mm-dd hh) . | Simulated peak time (yyyy-mm-dd hh) . | Error (number) . | Total volume error (%) . |
---|---|---|---|---|---|---|---|

19790627 | 115 | 91.4 | −20.5 | 1979-06-27 15:00 | 1979-06-27 16:30 | −1.5 | 15.7 |

19800812 | 108 | 103 | −8.3 | 1980-08-12 08:00 | 1980-08-12 17:00 | −9 | 12.8 |

19810407 | 114 | 109.3 | −4.1 | 1981-04-07 07:30 | 1981-04-07 09:00 | −1.5 | −4.5 |

19820616 | 170 | 151.2 | −11.1 | 1982-06-16 19:42 | 1982-06-16 17:30 | 2.2 | 7.1 |

19870512 | 109 | 149.2 | 36.9 | 1987-05-13 02:30 | 1987-05-13 03:00 | −0.5 | 25.6 |

19880829 | 99.9 | 95.3 | −4.6 | 1988-09-03 15:36 | 1988-09-03 15:30 | 0.1 | 18.1 |

19900607 | 108 | 128.5 | 19.0 | 1990-06-07 11:00 | 1990-06-07 09:00 | 2 | 16.3 |

19900615 | 185 | 154.9 | −16.3 | 1990-06-15 15:00 | 1990-06-12 23:00 | 64 | −5.6 |

19920322 | 130 | 133.1 | 2.4 | 1992-03-20 06:00 | 1992-03-20 06:00 | 0 | 1.4 |

19920516 | 172 | 158.4 | 11.3 | 1992-05-17 00:12 | 1992-05-16 23:30 | 0.7 | 11.5 |

19920615 | 109 | 102.2 | −6.2 | 1992-06-15 12:00 | 1992-06-15 13:00 | −1 | −5.2 |

19920622 | 259 | 250.9 | −3.1 | 1992-06-22 09:30 | 1992-06-22 09:30 | 0 | 1.9 |

19930704 | 162 | 134.6 | −17.5 | 1993-07-04 20:00 | 1993-07-04 20:30 | −0.5 | −28.1 |

19930720 | 170 | 263.6 | 55.1 | 1993-07-20 10:00 | 1993-07-20 08:00 | 2 | 49.2 |

19940425 | 98 | 82.7 | −15.6 | 1994-04-25 08:00 | 1994-04-25 08:00 | 0 | −6.0 |

19940718 | 206 | 184.3 | −10.5 | 1994-07-18 17:12 | 1994-07-18 17:30 | −0.3 | −4.4 |

19940905 | 129 | 178.2 | 38.1 | 1994-09-06 11:30 | 1994-09-06 11:30 | 0 | 13.7 |

19950701 | 423 | 401.5 | −5.1 | 1995-06-30 12:12 | 1995-06-30 12:00 | 0.2 | 2.6 |

19960601 | 108 | 108.8 | 4.5 | 1996-06-01 13:42 | 1996-06-01 15:30 | −1.8 | 0.0 |

19960717 | 140 | 128.1 | −17.4 | 1996-07-18 08:00 | 1996-07-18 08:30 | −0.5 | −17.9 |

19980522 | 293 | 309.6 | 5.7 | 1998-05-22 09:00 | 1998-05-22 09:30 | −0.5 | −0.3 |

19980624 | 169 | 144.8 | −14.3 | 1998-06-25 02:30 | 1998-06-25 03:30 | −1 | −9.8 |

19990716 | 97.6 | 94.3 | −3.4 | 1999-07-17 12:00 | 1999-07-17 10:00 | 2 | −23.0 |

20020818 | 140 | 115 | −17.9 | 2002-08-19 18:18 | 2002-08-19 16:00 | 2.3 | −17.0 |

20030516 | 109 | 92.9 | −14.8 | 2003-05-16 12:00 | 2003-05-16 12:30 | −0.5 | −14.3 |

20030605 | 81.5 | 56.8 | −30.3 | 2003-06-05 20:42 | 2003-06-05 22:00 | −1.3 | −19.7 |

20040515 | 206 | 186.7 | −9.4 | 2004-05-15 16:12 | 2004-05-15 16:30 | −0.3 | −4.6 |

20050601 | 183 | 161.2 | −11.9 | 2005-06-01 10:30 | 2005-06-01 11:30 | −1 | −13.2 |

20060411 | 108 | 118.8 | 10.0 | 2006-04-12 12:48 | 2006-04 − 12 13:30 | −0.7 | −2.3 |

*Notes*: The total number of floods is **29**, which are qualificatory relative to the error of peak discharge, is **24** and the ratio of qualifying simulation is **82.76%**, which are qualificatory relative to the error of peak time, is **27** and the ratio of qualifying simulation is **93.1%**, which are qualificatory relative to the error of total runoff volume, is **25** and the ratio of qualifying simulation is **86.21%**.

Floods . | Observed (m^{3}/s)
. | Simulated (m^{3}/s)
. | Percentage error (%) . | Observed peak time (yyyy-mm-dd hh) . | Simulated peak time (yyyy-mm-dd hh) . | Error (number) . | Total volume error (%) . |
---|---|---|---|---|---|---|---|

20100510 | 227 | 266 | 17.2 | 2010-05-13 15:35 | 2010-05-13 14:00 | 1.6 | 17.2 |

20100608 | 167 | 217.5 | 30.3 | 2010-06-08 16:00 | 2010-06-08 16:30 | −0.5 | 20.0 |

20100622 | 205 | 213 | 3.9 | 2010-06-24 10:15 | 2010-06-24 10:30 | −0.3 | 10.4 |

20110612 | 184 | 199.3 | 8.3 | 2011-06-15 06:35 | 2011-06-15 04:00 | 2.6 | 24.0 |

20120610 | 125 | 138.1 | 10.5 | 2012-06-11 05:15 | 2012-06-11 03:30 | 1.8 | 4.2 |

20140704 | 117 | 111.8 | −4.4 | 2014-07-05 05:20 | 2014-07-05 05:00 | 0.3 | −3.7 |

20150605 | 206 | 242.1 | 17.5 | 2015-06-08 16:10 | 2015-06-08 16:30 | −0.3 | 17.0 |

20150619 | 241 | 275.3 | 14.2 | 2015-06-21 14:30 | 2015-06-21 15:00 | −0.5 | 13.3 |

20160704 | 170 | 194.6 | 14.5 | 2016-07-04 00:30 | 2016-07-04 00:30 | 0 | 1.8 |

20170627 | 254 | 207.4 | −18.4 | 2017-07-01 11:30 | 2017-07-01 12:30 | −1 | −3.0 |

Floods . | Observed (m^{3}/s)
. | Simulated (m^{3}/s)
. | Percentage error (%) . | Observed peak time (yyyy-mm-dd hh) . | Simulated peak time (yyyy-mm-dd hh) . | Error (number) . | Total volume error (%) . |
---|---|---|---|---|---|---|---|

20100510 | 227 | 266 | 17.2 | 2010-05-13 15:35 | 2010-05-13 14:00 | 1.6 | 17.2 |

20100608 | 167 | 217.5 | 30.3 | 2010-06-08 16:00 | 2010-06-08 16:30 | −0.5 | 20.0 |

20100622 | 205 | 213 | 3.9 | 2010-06-24 10:15 | 2010-06-24 10:30 | −0.3 | 10.4 |

20110612 | 184 | 199.3 | 8.3 | 2011-06-15 06:35 | 2011-06-15 04:00 | 2.6 | 24.0 |

20120610 | 125 | 138.1 | 10.5 | 2012-06-11 05:15 | 2012-06-11 03:30 | 1.8 | 4.2 |

20140704 | 117 | 111.8 | −4.4 | 2014-07-05 05:20 | 2014-07-05 05:00 | 0.3 | −3.7 |

20150605 | 206 | 242.1 | 17.5 | 2015-06-08 16:10 | 2015-06-08 16:30 | −0.3 | 17.0 |

20150619 | 241 | 275.3 | 14.2 | 2015-06-21 14:30 | 2015-06-21 15:00 | −0.5 | 13.3 |

20160704 | 170 | 194.6 | 14.5 | 2016-07-04 00:30 | 2016-07-04 00:30 | 0 | 1.8 |

20170627 | 254 | 207.4 | −18.4 | 2017-07-01 11:30 | 2017-07-01 12:30 | −1 | −3.0 |

*Notes*: The total number of floods is **10**, which are qualificatory relative to the error of peak discharge, is **9** and the ratio of qualifying simulation is **90.0%**, which are qualificatory relative to the error of peak time, is **10** and the ratio of qualifying simulation is **100%**, which are qualificatory relative to the error of total runoff volume, is **9** and the ratio of qualifying simulation is **90.0%**.

In order to demonstrate the comparative performance of the proposed Distributed Model based on Spatially Distributed Velocity Field using Equation (14) (DMSDVF) in this study, the Distributed Model based on Slope and Net Rainfall Intensity using Equation (7) (DMSNR) and the lumped Xinanjiang model (LXM) are employed for comparison. The calibrating and validating datasets of the two models are the same as those of the proposed model. The parameters of the two employed models are set in accordance with the best performance based on the calibrating programme. The performance statistics of the three models are summarized in Table 6 according to the evaluation criteria. It can be observed from Table 6 that the proposed DMSDVF attains the best RMSE, NSCE, QR_{peak_discharge}, QR_{peak_time} and QR_{runoff} statistics of 12.6, 0.858, 82.76%, 93.1% and 86.21%, respectively, during calibration and attains the best RMSE, NSCE, QR_{peak_discharge} and QR_{peak_time} statistics of 11.3, 0.918, 90% and 100%, respectively, during validation. The results of this analysis indicate that the proposed model is able to obtain the best result in terms of different evaluation measures during the calibration phase; and except QR_{runoff} value, the proposed model is also able to obtain the best results in terms of other evaluation measures during the validation phase. The performances of the three models developed in this paper during the calibration and validation periods in the studied site are shown in Figures 4 and 5.

. | Evaluation result during calibration . | Evaluation result during validation . | ||||||||
---|---|---|---|---|---|---|---|---|---|---|

Model . | RMSE (m^{3}/s)
. | NSCE . | QR_{peak_discharge} (%)
. | QR_{peak_time} (%)
. | QR_{runoff} (%)
. | RMSE (m^{3}/s)
. | NSCE . | QR_{peak_discharge} (%)
. | QR_{peak_time} (%)
. | QR_{runoff} (%)
. |

DMSDVF | 12.6 | 0.858 | 82.76 | 93.1 | 86.21 | 11.3 | 0.918 | 90.0 | 100.0 | 90.0 |

DMSNR | 14.3 | 0.804 | 48.3 | 93.1 | 55.2 | 11.7 | 0.910 | 90.0 | 80.0 | 100.0 |

LXM | 15.0 | 0.792 | 65.5 | 79.3 | 58.6 | 13.5 | 0.887 | 40.0 | 90.0 | 80.0 |

. | Evaluation result during calibration . | Evaluation result during validation . | ||||||||
---|---|---|---|---|---|---|---|---|---|---|

Model . | RMSE (m^{3}/s)
. | NSCE . | QR_{peak_discharge} (%)
. | QR_{peak_time} (%)
. | QR_{runoff} (%)
. | RMSE (m^{3}/s)
. | NSCE . | QR_{peak_discharge} (%)
. | QR_{peak_time} (%)
. | QR_{runoff} (%)
. |

DMSDVF | 12.6 | 0.858 | 82.76 | 93.1 | 86.21 | 11.3 | 0.918 | 90.0 | 100.0 | 90.0 |

DMSNR | 14.3 | 0.804 | 48.3 | 93.1 | 55.2 | 11.7 | 0.910 | 90.0 | 80.0 | 100.0 |

LXM | 15.0 | 0.792 | 65.5 | 79.3 | 58.6 | 13.5 | 0.887 | 40.0 | 90.0 | 80.0 |

## DISCUSSION

In the process of model parameter calibration, the most important parameter of the Xinanjiang model is the areal mean free water capacity of the surface soil layer *S*_{m}, which controls the redistribution of the total runoff under the assumption of the surface soil layer as a reservoir. As such, it is divided into surface water, soil-free water and groundwater. This parameter can be determined by simulation of the secondary flood peak flow. Secondly, the parameters that should be concerned are *K _{g}* and

*K*,

_{i}*C*and

_{g}*C*, which can be estimated by analyzing the duration of flood recession and calibrated by simulating the recession curve. The other parameters of the Xinanjiang model are rather insensitive, which only have a weak impact on the local details of the flood simulation process.

_{i}There are few parameters in the concentration model, but they significantly affect the simulation results of peak discharge and peak time. The geomorphic unit hydrograph parameters in Equations (5)–(8) are mostly empirical parameters, which have no clear physical meaning and can only be determined by calibration. The improved geomorphic unit hydrograph model has only one parameter *μ*′, which can be determined both by calibration and measurement.

In order to analyze the effect of different velocity formulas, we use the same runoff generation model and parameters, and use Equations (7) and (14) to generate different geomorphic unit hydrographs to simulate rainfall-runoff. As can be seen from the previous results, the simulation and forecasting results of the distributed model use Equation (14) based on the energy conversion method are obviously better than the other models. When Equation (7) is used to compute the flow velocity in the grid, the velocity is often considered to be positively correlated with slope and rainfall intensity, that is, under the same conditions, the greater the slope and rainfall intensity, the faster the velocity will be. In this paper, parameter *b* of velocity formula is taken as 0.5 according to Maidment *et al.* (1996), parameter *d* is taken as 0.4 according to Li and Wen (1989) and the parameter *a* is calibrated as 0.5 based on the hydrological data. Although this method can also get good results, but its velocity distribution results are very different from the energy conversion method.

Figure 6(a) shows the velocity distribution computed by Equation (14) based on energy conversion. The velocity is generally 1–3 m/s on the steep slope, no more than 1 m/s on the flat land and mostly 1–2 m/s in the river channel. The velocity distribution considering slope and net rainfall intensity (in mm/h) is adopted by Equation (7) in Figure 6(b). It can be seen from Figure 6 that when the rain intensity is 1 mm/h, the velocity is generally 5–10 m/s on the steep slope, generally 2–5 m/s on the flat land and generally not more than 0.3 m/s in the river channel (many grid values are affected by a slope of 0). It is because the velocity under other rain intensities is in multiple ratios with the result under 1 mm/h with a ratio coefficient is *i ^{d}*. When the rainfall intensity is 10 mm/h, the ratios are 2.51, 3.31 for 20 mm/h and 3.90 for 30 mm/h, respectively. From the numerical results, it can be seen that the velocity result is high on the slope and reasonable in the river by considering the slope and net rainfall intensity.

If *a* = 0.5, *b* = 0.5 and *d* = 0.4 are selected for the velocity estimation method of Equation (8) considering slope, rainfall intensity and catchment area, the velocity result close to Equation (7) can be obtained when parameter *c* tends to 0. When *c* is greater than 0, it can improve the situation of the velocity difference between slope land and river channel in Equation (7). The larger the catchment area is, the larger the velocity gain of the grid point is, but it will cause a very complicated calibration process.

It is easy to know from Equation (15) that velocity *v* is directly proportional to the square root of parameter *μ*′ for any point in space, which is with a value of 15.38 m/s when *μ*′ = 1. The corresponding relationship between discharge and velocity is analyzed according to both the section measurement data and the stage–discharge relation of Lianyuan station. The average velocity of the section is computed to be 1.085 m/s by analyzing the total water volume and total kinetic energy of 39 floods. Table 7 presents the average velocity of each flood.

Floods . | Flow peak (m^{3}/s)
. | Total volume (10^{6} m^{3})
. | Average velocity (m/s) . | Floods . | Flow peak (m^{3}/s)
. | Total volume (10^{6} m^{3})
. | Average velocity (m/s) . |
---|---|---|---|---|---|---|---|

19790627 | 115 | 7.523 | 0.97 | 19980522 | 293 | 15.561 | 1.29 |

19800812 | 108 | 7.973 | 0.97 | 19980624 | 169 | 29.19 | 1.13 |

19810407 | 114 | 12.316 | 1.01 | 19990716 | 97.6 | 12.481 | 1.01 |

19820616 | 170 | 12.765 | 1.09 | 20020818 | 140 | 19.351 | 1.05 |

19870512 | 109 | 6.458 | 0.93 | 20030516 | 109 | 11.952 | 1.01 |

19880829 | 99.9 | 27.268 | 0.9 | 20030605 | 81.5 | 5.049 | 0.9 |

19900607 | 108 | 8.697 | 0.95 | 20040515 | 206 | 11.734 | 1.18 |

19900615 | 185 | 26.319 | 1.08 | 20050601 | 183 | 15.716 | 1.15 |

19920322 | 130 | 24.745 | 0.89 | 20060411 | 108 | 10.373 | 0.92 |

19920516 | 172 | 10.263 | 1.13 | 20100510 | 227 | 26.876 | 1.1 |

19920615 | 109 | 6.748 | 0.97 | 20100608 | 167 | 9.536 | 1.13 |

19920622 | 259 | 15.39 | 1.23 | 20100622 | 205 | 26.258 | 1.1 |

19930704 | 162 | 18.345 | 1.16 | 20110612 | 184 | 20.241 | 1.07 |

19930720 | 170 | 10.528 | 1.16 | 20120610 | 125 | 10.353 | 1.04 |

19940425 | 98 | 10.005 | 0.98 | 20140704 | 117 | 8.606 | 0.99 |

19940718 | 206 | 11.44 | 1.2 | 20150605 | 206 | 20.019 | 1.08 |

19940905 | 129 | 8.244 | 1.04 | 20150619 | 241 | 21.294 | 1.14 |

19950701 | 423 | 29.805 | 1.27 | 20160704 | 170 | 23.326 | 1.11 |

19960601 | 108 | 8.571 | 1.03 | 20170627 | 254 | 37.894 | 1.1 |

19960717 | 140 | 16.216 | 1.04 |

Floods . | Flow peak (m^{3}/s)
. | Total volume (10^{6} m^{3})
. | Average velocity (m/s) . | Floods . | Flow peak (m^{3}/s)
. | Total volume (10^{6} m^{3})
. | Average velocity (m/s) . |
---|---|---|---|---|---|---|---|

19790627 | 115 | 7.523 | 0.97 | 19980522 | 293 | 15.561 | 1.29 |

19800812 | 108 | 7.973 | 0.97 | 19980624 | 169 | 29.19 | 1.13 |

19810407 | 114 | 12.316 | 1.01 | 19990716 | 97.6 | 12.481 | 1.01 |

19820616 | 170 | 12.765 | 1.09 | 20020818 | 140 | 19.351 | 1.05 |

19870512 | 109 | 6.458 | 0.93 | 20030516 | 109 | 11.952 | 1.01 |

19880829 | 99.9 | 27.268 | 0.9 | 20030605 | 81.5 | 5.049 | 0.9 |

19900607 | 108 | 8.697 | 0.95 | 20040515 | 206 | 11.734 | 1.18 |

19900615 | 185 | 26.319 | 1.08 | 20050601 | 183 | 15.716 | 1.15 |

19920322 | 130 | 24.745 | 0.89 | 20060411 | 108 | 10.373 | 0.92 |

19920516 | 172 | 10.263 | 1.13 | 20100510 | 227 | 26.876 | 1.1 |

19920615 | 109 | 6.748 | 0.97 | 20100608 | 167 | 9.536 | 1.13 |

19920622 | 259 | 15.39 | 1.23 | 20100622 | 205 | 26.258 | 1.1 |

19930704 | 162 | 18.345 | 1.16 | 20110612 | 184 | 20.241 | 1.07 |

19930720 | 170 | 10.528 | 1.16 | 20120610 | 125 | 10.353 | 1.04 |

19940425 | 98 | 10.005 | 0.98 | 20140704 | 117 | 8.606 | 0.99 |

19940718 | 206 | 11.44 | 1.2 | 20150605 | 206 | 20.019 | 1.08 |

19940905 | 129 | 8.244 | 1.04 | 20150619 | 241 | 21.294 | 1.14 |

19950701 | 423 | 29.805 | 1.27 | 20160704 | 170 | 23.326 | 1.11 |

19960601 | 108 | 8.571 | 1.03 | 20170627 | 254 | 37.894 | 1.1 |

19960717 | 140 | 16.216 | 1.04 |

*μ*′ using Equation (24) is 0.00498, which is basically consistent with the calibration result

*μ*′ = 0.005.The average elevation is 268 m and the outlet elevation is 137 m in the studied area. Assuming that the gravity potential energy is all converted into kinetic energy under an ideal condition, the average velocity of the outlet section of study area

*v*

_{ideal}= is about 51 m/s. The energy conversion

*μ*

_{efficiency}can be computed by Equation (25):

When *μ*′ = 0.005, the energy conversion efficiency is 0.45‰, and it can be seen that a large amount of energy is lost in the process of water flow movement, denoting an inefficient conversion rate.

## CONCLUSIONS

Hydrologic models are widely used tools for flood forecasting in the world, as the demand for accurate and reliable flood forecasts has increased in the management of water resources. This paper presents a new method for watershed concentration computation, which can be used to construct the hydrological model in areas without enough data. In order to simplify the computation, the residual energy coefficient is assumed to be constant in this paper and its value can be determined by calibration or measurement of the velocity. From the application point of view, the discretization of residual energy coefficient in spatial distribution is realized by measuring the velocity of different control nodes. As a matter of fact, it translates the flow concentration problem into a measurement problem (we can measure the flow velocity directly, or compute the flow velocity according to the hydraulic formula by using the measured channel section), and more precise flow concentration parameters can be obtained by more measurement data.

To demonstrate the applicability of the proposed method, a real-world catchment, Lianyuan River Basin is selected. For a 0.5 h time step, 29 historical floods in 28 years (1979–2006) are used for calibration, whereas 10 historical floods in 8 years (2010–2017) are utilized for verification. Three criteria are employed to evaluate the performances of the developed model. The result comparison indicates a significant improvement in the flood forecasting accuracy of the established distribution model, which uses energy conversion analysis of the geomorphic unit hydrograph. Hence, the velocity distribution computed by the energy conversion equation is more reasonable. The results of flow velocity can be used for the computation of both slope and river network flow concentration, which can ensure the consistency of flow concentration computation to a certain extent and reduce the influence of watershed threshold on model parameters in subwatershed division.

In terms of energy conversion efficiency, the ratio of gravitational potential energy to kinetic energy is very low in the process of flood movement, and most of it is converted into other forms of energy. There are many factors that affect the energy conversion. Future work should focus on further exploring the quantitative relationship between the energy coefficient and relevant influencing factors.

## ACKNOWLEDGMENTS

The authors are grateful to the Project of Key Science and Technology of the Henan Province (No. 202102310588), Henan Province University Scientific and Technological Innovation Team (No. 18IRTSTHN009) and the supports of Hunan Hydrological and Water Resources Survey Center and Loudi Hydrological and Water Resources Survey Center. We gratefully acknowledge the critical comments and corrections of the anonymous reviewers and editors, which improved the presentation of the paper considerably.

## DATA AVAILABILITY STATEMENT

Data cannot be made publicly available; readers should contact the corresponding author for details.