Multi-Spatial Resolution Rainfall-Runoff Modelling—A Case Study of Sabari River Basin, India

: One of the challenges in rainfall-runoff modeling is the identiﬁcation of an appropriate model spatial resolution that allows streamﬂow estimation at customized locations of the river basin. In lumped modeling, spatial resolution is not an issue as spatial variability is not accounted for, whereas in distributed modeling grid or cell resolution can be related to spatial resolution but its application is limited because of its large data requirements. Streamﬂow estimation at the data-poor customized locations is not possible in lumped modeling, whereas it is challenging in distributed modeling. In this context, semi-distributed modeling offers a solution including model resolution and estimation of streamﬂow at customized locations of a river basins with less data requirements. In this study, the Hydrologic Engineering Center-Hydrologic Modeling System (HEC-HMS) model is employed in semi-distribution mode on river basins of six different spatial resolutions. The model was calibrated and validated for ﬁfteen and three selected ﬂood events, respectively, of three types, i.e., single peak (SP), double peak (DP)- and multiple peaks (MP) at six different spatial resolution of the Sabari River Basin (SRB), a sub-basin of the Godavari basin, India. Calibrated parameters were analyzed to understand hydrologic parameter variability in the context of spatial resolution and ﬂood event aspects. Streamﬂow hydrographs were developed, and various veriﬁcation metrics and model scores were calculated for reference- and calibration- scenarios. During the calibration phase, the median of correlation coefﬁcient and NSE for all 15 events of all six conﬁgurations was 0.90 and 0.69, respectively. The estimated streamﬂow hydrographs from six conﬁgurations suggest the model’s ability to simulate the processes efﬁciently. Parameters obtained from the calibration phase were used to generate an ensemble of streamﬂow at multiple locations including basin outlet as part of the validation. The estimated ensemble of streamﬂows appeared to be realistic, and both single-valued and ensemble veriﬁcation metrics indicated the model’s good performance. The results suggested better performance of lumped modeling followed by the semi-distributed modeling with a ﬁner spatial resolution. Thus, the study demonstrates a method that can be applied for real-time streamﬂow forecast at interior locations of a basin, which are not necessarily data rich.


Introduction
Hydrologic science, including rainfall-runoff modeling, has been an important discipline for several decades because of intriguing science questions and associated socioeconomic issues [1][2][3]. Increased data access over a wider area for longer periods, and availability of high computing power in combination with various advanced techniques has allowed hydrologic sciences to evolve over different paradigms, i.e., empiricism, theory, and computational simulation, and now data-intensive sciences [4][5][6]. It enhanced the understanding of rainfall-runoff processes as well as assisted in developing solutions of various aspects ranging from engineering to policy. However, heterogeneity of hydrologic processes in space and time, and unavailability of data either due to a lack of monitoring networks or absence of measurements during extremes make rainfall-runoff modeling The SRB receives an average seasonal rainfall (June to September) of 1077 mm (calculated for 1966-2017). The SRB has four gauge-discharge (GD) stations and two-level flood forecasting (FF) stations. While GD stations record streamflow and gauge values, FF stations provide the flood warning issued as water level. Konta is the most downstream station among all the four GD stations, and it is in close proximity to the outlet of the basin. Therefore, Konta is treated as an outlet of the basin. The other three GD stations, Injaram, Potteru and Saradaput are the three stations in the interior of the basin. Average daily rainfall and average daily streamflow at Konta for the south-west monsoon period (June to September) are 8 mm/day and 944 m 3 /s, respectively, and these values are relatively high as compared to the annual averages (3.8 mm/day and 460 m 3 /s).

Data
Daily streamflow data at Konta, Injaram, Saradaput, and Potteru GD stations are available for the period of record 1966 to 2017, 1966 to 2006, 1997 to 2005 and 1970 to 2015, respectively. The streamflow data was obtained from the Krishna Godavari Board Organization (KGBO), CWC, India. Among all GD stations, Konta has the longest period of the record from 1966 to 2017, and consideration of the common period of the record for all stations decreased the length of available data, and consequently, the number of major flood events available for calibration were reduced. Therefore, streamflow data from Konta was considered to define the flood events for model calibration. However, data of all four GD stations was considered to define flood events for model validation.
Daily gridded rainfall data of 0.25 • × 0.25 • spatial resolution was obtained from the India Meteorological Department (IMD) [47]. The gridded data was used to extract the weighted average rainfall for each sub-basin from six configurations for all the flood events. Sub-basin specific rainfall amounts for the selected events were given as input to the Hydrologic Engineering Center-Hydrologic Modeling System (HEC-HMS) model. The Shuttle Radar Topography Mission (SRTM) based digital elevation model (DEM) of 30 m × 30 m [48] spatial resolution was used to extract the river streams and basin's physical characteristics. Land use and land cover (LULC) data of 100 m spatial resolution for the years 1985, 1995, and 2005 is available for the Indian region [45] from which the values corresponding to the SRB were obtained. Soil data for the SRB of approximately 2.5 km spatial resolution were obtained from Harmonized World Soil Database Version 3.0 [46]. The LULC and soil data were used in obtaining the excess rainfall related parameters. Similar to rainfall, sub-basin specific LULC and soil data were extracted from the corresponding gridded data.

Selection of Flood Events
An event is defined as when the daily streamflow is greater than the 90th percentile threshold at least for five consecutive time steps and continues until the streamflow is less than the given threshold for four consecutive time steps [49]. Note that the events defined by the algorithm were plotted, analyzed and then revised with respect to start and endpoints. The final set of 70 events consist of 21 single peak (SP), 13 double peak (DP) and 36 multiple peak (MP) events. The top five events of higher magnitude in each category were selected and used to calibrate the model. Table 1 details the characteristics of the events including peak flow, duration, time to peak and accumulated rainfall of the event. Three events, viz., 2003SP, 1997DP and 2004MP of duration 10, 15 and 30 days were selected for model validation. Event characteristics for three interior locations, i.e., Injaram, Saradaput, Potteru, along with the basin outlet Konta are mentioned in Table 2. Table 1. Flood aspects at Konta GD station, i.e., peak flows (Qp), flood duration (t), time to peak (Tp) and accumulated rainfall (Rsum) for the selected 15 events used in the model's calibration. The letters in the event name, i.e., SP, DP and MP, correspond to single-, double-and multiple-peak, respectively.

Event
Peak Flows (Q p , m 3

Model Setup
The HEC-HMS model is a physically-based, semi-distributed hydrological model developed by the U.S. Army Corps of Engineers [50]. The model can be run in lumped and semi-distributed modes at different time intervals ranging from minute to daily time steps, and parameters can be estimated in both event and continuous modes [51]. The HEC-HMS grouped key rainfall-runoff processes into four components and wrapped multiple methods as one method to model each of these components. In the present study, the HEC-HMS model was applied on a daily time scale in the event mode. The Soil Conservation Service-Curve Number (SCS-CN) method [52] was used to compute the excess rainfall. The Clark's unit hydrograph (UH) method [53] was adopted to compute the direct runoff. An exponential recession model was employed to estimate baseflow, and the Muskingum method [54] was used for the channel flow routing.
In this study, the SRB was delineated according to the stream network with five different area thresholds. The resulting set of sub-basins represents the spatial resolution, and each spatial resolution was treated as a configuration in the modeling context as per the number of sub-basins. The model was run in the semi-distributed mode for these configurations. Moreover, the entire basin was considered as a single unit, and the model was run in the lumped mode for the configuration. Thus, a total of six configurations, i.e., a lumped model and five semi-distributed models were considered in this study. Basin delineation and calculation of physical characteristics of river and sub-basin such as length, width, area, basin lag and slope for all configurations were performed using the HEC-GeoHMS tool [55]. The above physical characteristics were used to estimate the initial values of parameters using empirical equations [50]. All six configurations were exported to the HEC-HMS model environment to perform the event-based rainfall-runoff modeling. The model calibration procedure is mentioned below.

Model Calibration
The model calibration consists of the following two components, i.e., (1) estimation of initial model parameter values using well-established empirical relationships, and (2) estimation of optimal values of model parameters using only the HEC-HMS model's Neldar-Mead simplex algorithm and in combination with trial-and-error.

Model Initial Parameter Estimation
The set of parameters consist of initial abstractions (I a ), curve number (CN) and percentage of impervious area (imp%) of SCS-CN method, time of concentration (T c ) and storage coefficient (R) of Clark's UH method, initial baseflow (Q 0 ), recession constant (k), and ratio to peak (R p ) of recession baseflow method, and travel time (K) for flood wave and dimensionless weight (X) of the Muskingum method. The percentage of impervious area was calculated from satellite maps. Therefore, this parameter was not included in the model calibration. The equations to estimate the initial values of parameters are shown in Table 3. Additional information related to rainfall-runoff methods of the HEC-HMS model can be found in the Supplementary Material (Section SI1).
Ratio to peak (R p ) 0.1 0.1-1 9 Travel time (K) in hr Weighting factor (X) [57], c [50] Where CN composite is the composite curve number, i is the index of sub-basins of uniform land use and soil type, A i is the area of i th sub-basin, CN i is the curve number for i th sub-basin, S is the potential maximum retention that represents mainly the infiltration occurring after runoff has started, A Built−up is the built-up area of a basin or a sub-basin, L is the longest flow path of a sub-basin, s is the slope of the channel, l is the channel length, V W is the wave celerity, B is the top width of the water surface, dQ dy is the slope of the rating curve at a known cross-section, Q 0 is the reference flow from the inflow hydrograph, an average flow, midway between the baseflow and peak flow, Q t is the baseflow at any time t, S 0 is the frictional slope, and ∆x is the length of the reach.

Model Parameter Optimization
The parameters I a , CN, T c , R, Q 0 , k and R p were optimized simultaneously using the simplex algorithm for the lumped model. For the remaining five configurations, i.e., for five semi-distributed models, only I a and CN were optimized using the simplex algorithm and the remaining seven parameters were estimated in the trial-and-error mode one-by-one in the following order, i.e., T c , R, Q 0 , k, R p , K and X. The parameter R was calibrated first in case of either high volumetric error or no error in time to peak flow corresponding to optimal values of I a and CN. The range of values for which parameters can be varied during its calibration suggested by the HEC-HMS user manual is mentioned in Table 3.
Optimal parameter estimation in the trial-and-error mode was stopped when the Nash Sutcliffe Efficiency (NSE) of simulations was greater than or equivalent to 0.4. The chosen NSE threshold is subjective and can be modified.
The HEC-HMS's optimization routine Nelder and Mead Simplex algorithm was chosen to minimize the objective function, Peak-Weighted RMSE (PWRMSE) [50]. The PWRMSE is a measure that implicitly compares magnitudes of flow peak, volume and time to peak between observed and simulated hydrographs. As it is seen from the equation (Table 4), the measure is a weighted sum of squared differences of simulated flows and corresponding observed flows. The weight of each ordinate is proportional to its observed streamflow magnitude.
Streamflow hydrographs that were derived based on initial values of parameters referred to as the reference scenario, whereas streamflow hydrographs that were derived based on optimal values of parameters treated as the calibration scenario. Streamflow from both reference and calibration scenarios were evaluated by calculating multiple verification metrics, i.e., the correlation coefficient (r), Nash-Sutcliffe Efficiency (NSE), Root Mean Square Error (RMSE), Mean Absolute Error (MAE) [58], Percentage Bias (PBIAS) [59], Percentage Error in Peak flows (PEP) [50] and Error in Time to Peak (ETP). Table 4 shows the mathematical formulae of the objective function and verification metrics. Table 4. Mathematical formulae of the objective function and verification metrics used for the model calibration.

Metric Equation Range
Peak-Weighted Root Mean Square Error (PWRMSE) (m 3 /s) Root Mean Squared Error, RMSE (m 3 /s) Error in Time to Peak, ETP (number of days) Where Z is the objective function, NQ is the number of hydrograph ordinates, N indicates the number of observations or the duration in days of the event; O and P are the observed and simulated values; O and P correspond to the average of observed and simulated values; Q PP is the peak flow of simulated hydrograph and Q PO is the peak flow of the observed hydrograph.

Model Validation
As part of the model validation, streamflow for the three flood events was estimated and then verified at all four selected locations, i.e., the basin's outlet and three interior locations of the basin, i.e., Injaram, Saradaput and Potteru. The model parameter set consists of CN, I a , T c , R, k, Q 0 , R p , K and X, and calibration of the model for each event yields an optimal value for each parameter. Thereby, 15 events together result in a set of 15 values for each of the model parameters. In this way, first, an ensemble set of parameters was created, and then, a corresponding ensemble of streamflow was produced for the three selected events at all four selected locations. However, noting that the CN and I a values as event specific as they depend on antecedent moisture conditions (AMC), the values of these parameters were not taken from calibration, consequently, the parameters were not part of the ensemble set of parameters. Similar to the calibration phase, the CN values were estimated for each event using rainfall five days prior to the event and LULC [60], and I a was calculated as 30% of the potential maximum retention (S) ( Table 1). The ensemble of streamflows was evaluated using the verification metrics shown in Table 4 as well as the Continuous Ranked Probability Score (CRPS), which is an ensemble verification metric [61]. The metrics in Table 4 are single-valued measures, and the ensemble mean of all ensemble members was used here. Thus, the model output in the validation mode was evaluated in both single-valued and ensemble modes.

Results and Discussions
5.1. The Six Configurations of Sabari River Basin Figure 2 displays sub-basins that resulted from six configurations. The configuration that has no sub-basins resulted from treating the entire basin as a single unit and is named as S1. The configurations that were generated with five different area thresholds produced as small as eight to 132 sub-basins (Table 5). Accordingly, the basin configuration was named S8, S20, S48, S52 and S132. Sub-basins in all configurations except S8 were algorithm-based and automatic. The S8 was resulted from merging of a few sub-basins with a small stream order so that basins with only the highest stream order were retained. The logic behind the S8 is to perform the rainfall-runoff calculations by considering the tributaries that are directly joining the Sabari river. As the spatial resolution becomes finer, a few sub-basins are further divided into smaller sub-basins, and a few remained intact. Note that a sub-basin of S8 contains one or multiple reservoirs, whereas a subbasin of a finer configuration typically contains only a part of the reservoir. Thus, basins with finer configuration require multiple sub-basins to encompass completely or a part of the reservoir.  Several topographical parameters, i.e., area, longest flow path, elevation, average channel slope, average basin slope and basin lag, were calculated for all the sub-basins for the six spatial resolutions. Table 5 displays basin/sub-basin average values of topographical parameters. As the spatial resolution became finer, the sub-basin average values of parameters were decreased. The values for the six configurations differed significantly suggesting basins of a different order. The relation between area and physical characteristics was found to be a power-law; thus, the characteristics exhibited scaling with respect to the basin or sub-basin area. Table 5. Average values of topographical parameters of Sabari River Basin. The values for six configurations were calculated by taking the average value of all sub-basins of the corresponding configuration. Here, At is the threshold area used to delineate the configuration, A̅ is the mean area of the basin, L̅ is the mean longest flow path, and Lag is the basin lag.

Threshold Area
Mean Longest Flows Mean Basin Lag S8 S1 Konta Konta Konta S132 S58 S48 Several topographical parameters, i.e., area, longest flow path, elevation, average channel slope, average basin slope and basin lag, were calculated for all the sub-basins for the six spatial resolutions. Table 5 displays basin/sub-basin average values of topographical parameters. As the spatial resolution became finer, the sub-basin average values of parameters were decreased. The values for the six configurations differed significantly suggesting basins of a different order. The relation between area and physical characteristics was found to be a power-law; thus, the characteristics exhibited scaling with respect to the basin or sub-basin area. Table 5. Average values of topographical parameters of Sabari River Basin. The values for six configurations were calculated by taking the average value of all sub-basins of the corresponding configuration. Here, A t is the threshold area used to delineate the configuration, A is the mean area of the basin, L is the mean longest flow path, and Lag is the basin lag.

Configuration
Threshold

Results from the Model Calibration
As mentioned in Section 4.3.2, the set of seven parameters were optimized simultaneously for S1. For the remaining five configurations, the set of nine parameters were optimized. The two parameters, I a and CN were optimized simultaneously, and the remaining seven parameters were optimized one-by-one through trial-and-error. The calibration was continued until the NSE of the simulations was greater than or equal to 0.4. It was observed that both baseflow and routing related parameters (k, Q 0 , R p , K and X) have less influence in improving the NSE of simulations, therefore these parameters were always assigned with their initial values. This leaves excess rainfall and direct runoff related parameters (CN, I a , T c and R) as active in terms of calibration; hence optimal values of these parameters were analyzed for six spatial resolutions with respect to 15 flood events and event characteristics. The following sub-sections detail how these optimized parameters vary in multiple aspects.

Spatial Resolution
The set of calibrated parameters was analyzed for 15 events for six spatial resolutions in the context of basin response as well as basin area. Figure 3 has four panels, and each panel corresponds to a parameter. Each panel has six boxplots, and each boxplot corresponds to a configuration. A filled circle in each boxplot correspond to a sub-basin's parameter value for an event of a configuration. Thus, filled circles represent all values of a parameter for 15 events and all sub-basins of the configuration. Importantly, a boxplot with filled circles highlights whether values take only a certain portion of the range or vary continuously throughout the range, thus provides information on their distribution. Increased sample size (number of filled circles) was observed for finer resolution configurations because of an increase in the number of sub-basins.
A skewed and a horizontal band of values were observed for all the configurations and four parameters (I a , CN, T c , R). A skewed distribution implies a unique response of a sub-basin, where a long right or left tail of values was observed. Whereas the horizontal band of values suggests either a similar response of a sub-basin irrespective of the event or a similar response of all sub-basins to an event, hence the same set of values was observed. In this regard, the values of the parameters were analyzed with respect to the sub-basin area (Supplementary Materials Figure SI1-SI5) and an event. Although the sample size of the parameters for each configuration is different, the following was observed. The median of the values for parameters I a , CN, T c and R is high (low) for the S20 (S132, S48 and S1), S48 (S58 and S1), S1 (S58 and S8) and S1 (S132, S58, S4, S20 and S8) configurations, respectively. The results indicate the dominance of no single configuration across all parameters either for high or low values.
Noting that each parameter corresponds to a different physical process, the results suggest a different response of sub-basins of a configuration for each process and an association among parameters of a configuration. For example, relatively large I a values of the S20 configuration indicate dry antecedent moisture conditions and are in accordance with the corresponding moderate CN values. Moreover, relatively low T c values are in accordance with observed high R values for the S58 configuration. Other notable observations are decreased R values for configurations of finer spatial resolution (S8 to S132) and exponentially increasing R values with respect to the basin area, e.g., for S58 configuration ( Figure SI6). The higher the R, the longer it takes for attenuation of peak flow and this was observed for larger basins. Thus, the results are consistent. vary in multiple aspects.

Spatial Resolution
The set of calibrated parameters was analyzed for 15 events for six spatial resolutions in the context of basin response as well as basin area. Figure 3 has four panels, and each panel corresponds to a parameter. Each panel has six boxplots, and each boxplot corresponds to a configuration. A filled circle in each boxplot correspond to a sub-basin's parameter value for an event of a configuration. Thus, filled circles represent all values of a parameter for 15 events and all sub-basins of the configuration. Importantly, a boxplot with filled circles highlights whether values take only a certain portion of the range or vary continuously throughout the range, thus provides information on their distribution. Increased sample size (number of filled circles) was observed for finer resolution configurations because of an increase in the number of sub-basins.  Figure SI1-SI5) and an event. Although the sample size of the parameters for each configuration is different, the following was observed. The median of the values for parameters Ia, CN, Tc and R is high (low) for the S20 (S132, S48 and S1), S48 (S58 and S1), S1 (S58 and S8) and S1 (S132, S58, S4, S20 and S8) configurations, Apart from the above analysis, the values of parameters for sub-basins that are common in both S8 and S58 configurations were analyzed ( Figure 4). A total of five sub-basins were chosen and were named differently in two different configurations. The common basins for both configurations were named in S58 (S8) as W1110 (W4), W1140 (W5), W1150 (W6), W1160 (W7) and W1170 (W8). The boxplot of values for the same basin of two configurations indicates that the parameters are not similar in their magnitude and distribution. The above implies that the response of a basin is different when it is modelled with a different configuration. However, this is not valid because the basin's response will be the same irrespective of the configuration. All parameters, except I a , exhibited an association with the sub-basin area (Supplementary Materials Figure SI1-SI5). However, a medium or large spread of parameter values implies a significantly different basin response and raises questions on the transferability of parameters for a basin for all events. In this regard, the values of the parameters were analyzed for each event separately. tion. The above implies that the response of a basin is different when it is modelled with a different configuration. However, this is not valid because the basin's response will be the same irrespective of the configuration. All parameters, except Ia, exhibited an association with the sub-basin area (Supplementary Materials Figure SI1-SI5). However, a medium or large spread of parameter values implies a significantly different basin response and raises questions on the transferability of parameters for a basin for all events. In this regard, the values of the parameters were analyzed for each event separately.
10 -1 10 0 10 1 10 2 10 -1 10 0 10 1 10 2    Figure 5 has four panels, and each panel corresponds to a parameter. Calibrated parameter values of all sub-basins of a configuration specific to an event are displayed in a boxplot. A set of six boxplots, of which each boxplot is in a different color and corresponds to a configuration, were developed for each event, and are arranged in chronological order of events. A filled circle in a boxplot corresponds to a value of the parameter for a sub-basin. Thus, the number of filled circles in each boxplot is equivalent to the number of sub-basins of the configuration. A notable observation among all plots is the pattern exhibited by the CN, i.e., high and low CN values for events of the early-and recent-part of the record. High CN values correspond to less infiltration, which may potentially result in higher runoff. However, it is a counter-intuitive observation to see high CN values in the early part of the record before any significant development in the basin. Apart from CN, no other parameters exhibited a pattern except that the range of values for a few parameters is high for certain events.

Selected Flood Events
Except for I a and CN values of S58 configuration, the boxplots broadly suggest no large variations in parameter values among sub-basins of the same configuration (indicated by boxplot of small size) and all configurations of an event (indicated by a set of boxplots corresponding to an event). However, the boxplots of values indicate a large event-to-event variability, particularly for I a values (Figure 5a). Moreover, patterns with small variability in CN values (Figure 5b) and moderate to small variability in T c and R values (Figure 5c,d) are mainly responsible for the boxplot of values to be different for six configurations were observed. Nevertheless, excluding 1966SP and 1986SP, no event with a combination of a parameter set with either high or low values was observed. In this regard, parameters were analyzed with respect to flood event characteristics. boxplots corresponding to an event). However, the boxplots of values indicate a large event-to-event variability, particularly for Ia values (Figure 5a). Moreover, patterns with small variability in CN values (Figure 5b) and moderate to small variability in Tc and R values (Figure 5c,d) are mainly responsible for the boxplot of values to be different for six configurations were observed. Nevertheless, excluding 1966SP and 1986SP, no event with a combination of a parameter set with either high or low values was observed. In this regard, parameters were analyzed with respect to flood event characteristics. Boxplots of red, green, orange, blue, purple and yellow colors correspond to S132, S48, S8, S58, S20 and S1, respectively. The letters in the event name, i.e., SP, DP and MP, correspond to single peak, double peak and multiple peak flood events. Panels in the clockwise direction from the top left, (a-d) correspond to I a , CN, T c and R, respectively.

Event Characteristics
Four characteristics of a flood event, i.e., peak flow (Q p ), event duration (t) and time to peak (T p ), and accumulated rainfall (R sum ), corresponding to the event were considered as part of the analysis (see Table 1). While the time-to-peak values are uniformly distributed, flood duration values are relatively less right-skewed than the flood peak values ( Figure 6). Rainfall amount increased with an increase in the flood duration of an event. With an event exclusion, two aspects, i.e., flood duration and flood peak exhibited a moderate linear association. However, flood magnitude exhibited no apparent association. As no clear pattern is observed among different aspects of flood events, parameters of all configurations were analyzed with event characteristics (Figure 7). as part of the analysis (see Table 1). While the time-to-peak values are uniformly distributed, flood duration values are relatively less right-skewed than the flood peak values ( Figure 6). Rainfall amount increased with an increase in the flood duration of an event.
With an event exclusion, two aspects, i.e., flood duration and flood peak exhibited a moderate linear association. However, flood magnitude exhibited no apparent association. As no clear pattern is observed among different aspects of flood events, parameters of all configurations were analyzed with event characteristics (Figure 7). In Figure 7, each row corresponds to an event characteristic, and each column corresponds to a parameter. As parameter Ia was not exhibiting any specific pattern with any event characteristic, the corresponding association plot is not shown in Figure 7. Relatively high CN values were observed for high peak flow events and value greater than 70 were observed for approximately eight events for most of the configurations. High CN values with a small spread (i.e., small boxplot) were observed for the two highest peak flow events (Figure 7a). The large spread in CN values for three specific flood durations (16, 28 and 37 days) suggests that the CN values are different for flood events with the In Figure 7, each row corresponds to an event characteristic, and each column corresponds to a parameter. As parameter I a was not exhibiting any specific pattern with any event characteristic, the corresponding association plot is not shown in Figure 7. Relatively high CN values were observed for high peak flow events and value greater than 70 were observed for approximately eight events for most of the configurations. High CN values with a small spread (i.e., small boxplot) were observed for the two highest peak flow events (Figure 7a). The large spread in CN values for three specific flood durations (16, 28 and 37 days) suggests that the CN values are different for flood events with the same flood duration (Table 1; Figure 7b). The difference in CN values reflects the different infiltration conditions, which can be attributed to different antecedent rainfall and soil moisture conditions of the events. Unlike for flood duration, the large spread in CN values for three different time to peak values (6, 9 and 15 days) was seen in all configurations and it indicates different responses from sub-basins of a configuration (Figure 7c).
Most T c values for many configurations fall within a range of values from 3 to 10 h. High T c values (~10 h) were observed for two flood events of high peak flow for all configurations (Figure 7d). Typically, high T c values can be related to flood events of large magnitude for which most of the basin contributes as well as large flood events that result from high rainfall events in mid-and down-stream sub-basins. Thus, spatial variability of rainfall plays a key role in T c values, and it is dominant for all configurations except S1. This also explains the relatively high T c values of the S1 configuration. T c values of the S1 were increased with respect to flood duration (Figure 7e), and it can be attributed to high rainfall amount ignoring spatial variability; typically, in these conditions, increased rainfall is associated with increased flood duration. Similar to CN values, the large spread in T c values for the same flood duration and time to peak was observed and it is due to combining values from different flood events of the same flood duration or same time to peak (Table 1;  The rows from top to bottom correspond to flood peak, flood duration and time to peak, respectively. Boxplots of red, blue, green, purple, orange and yellow colors correspond to S132-, S58-S48-, S20-, S8-and S1-configurations, respectively.

Performance Metrics
Model performance was analyzed by calculating various verification metrics using estimated and corresponding observed streamflow values for six configurations, for the 15 events. The verification metrics consist of r, NSE, RMSE, MAE, PBIAS, PEP and ETP, and were calculated separately for streamflow of reference and calibrated scenarios (Table  A1 and   to T c and (g-i) correspond to R, respectively. The rows from top to bottom correspond to flood peak, flood duration and time to peak, respectively. Boxplots of red, blue, green, purple, orange and yellow colors correspond to S132-, S58-S48-, S20-, S8-and S1-configurations, respectively.
Storage coefficient (R) values of S1 increased with peak flow magnitude; however, approximately similar values were observed for S8 (Figure 7g). Relatively, smaller R values were observed for S58 and S132, and they exhibited event-to-event variation only for a few events. The R values indicate the number of hours that a basin can temporarily store the excess rainfall, and it explains relatively small values for the basins of finer resolution as compared to the basins of coarser resolution. The three large sub-basins of S8 (see Figure 2) have one or more hydraulic structures which can store the excess rainfall during the event; therefore, similar R values were observed for S8. Similar to CN and T c values, the R values for the flood events of the same flood duration or time to peak differ significantly (Table 1; Figure 7h,i). The two events with the same flood duration and time to peak resulted in a large spread (large boxplot) of R values. Note that there are also a few flood events for which R values are the same. Therefore, no large spread was observed, e.g., 1972SP, 2005DP,1988DP and 1969MP.

Performance Metrics
Model performance was analyzed by calculating various verification metrics using estimated and corresponding observed streamflow values for six configurations, for the 15 events. The verification metrics consist of r, NSE, RMSE, MAE, PBIAS, PEP and ETP, and were calculated separately for streamflow of reference and calibrated scenarios (Tables A1 and A2). As seen in the table, each of these scenarios has 90 values for each verification metric and result from a combination of six configurations and 15 events.
Positive correlation coefficient (r) values for both scenarios were observed. For the calibration scenarios, approximately 92% of the values are greater than 0.7 and 50% of the values are greater than 0.9. High values of r,~0.85 or greater values were observed for at least one of the six configurations when 1988DP and 2014MP are excluded. However, the r values for the reference scenario suggest its poor performance, yet 33% of values are greater than 0.9 and nine of 15 events exhibited a value of~0.85 at least for one of the configurations. In terms of NSE, approximately 85% of the values were greater than 0.4 for the calibration scenario. However, both positive and negative values were observed for the reference scenario, i.e., 50% of the NSE values are smaller than zero, and 20% of the values are greater than zero and less than 0.4. For the calibration scenario, the highest NSE values were greater than 0.4 for at least one of the six configurations for all the events. Whereas, for the reference scenario, the highest values of greater than 0.4 were observed only for seven of the 15 events. These suggest the better performance of the calibration scenario. While the linear association among estimated and corresponding observed streamflow and performance of estimated streamflow with respect to climatology are essential, the difference in streamflow magnitude should not be large. Therefore, metrics based on error terms were estimated.
The error term corresponds to a flow magnitude, and it is challenging to set a threshold to define the best or acceptable value. In this regard, first, the configuration with the lowest error was identified. Second, the score "1" was assigned for configurations with error values under 25% of the lowest error; otherwise "0" was assigned. Thus, if the error is within 25% of the lowest error for all configurations, then all configurations get a score of "1" and are accepted as equally performing. Finally, the above criterion was repeated for all 15 events and configurations. A model score can range from 0 to 15, and a high score indicates a relatively better performing configuration. The scheme was extended for both r and NSE metrics, but configurations with the metric values within 10% of the best metric were considered. Similarly, the model scores were calculated for error in time to peak (ETP), and configurations with non-zero values were penalized with a score of "0". The model scores were tabulated for all configurations of both reference and calibration scenarios, and the label of a column without the asterisk (*) corresponds to the reference scenario (Table 6). Table 6. Model scores were calculated for six configurations of both reference and calibration scenarios. A model score can range from 0 to 15, and a high score indicates a relatively better performing configuration. The label of a column without and with asterisk (*) corresponds to streamflow simulations of reference and calibration scenarios, respectively. Verification Metric S132 S132 * S58 S58 * S48 S48 * S20 S20 * S8 S8 * S1 S1 * Among all configurations of calibrated scenario, S1 stands out with very high model score values for four verification metrics (r, NSE, RMSE and MAE) followed by S132. Except for S8, all configurations exhibited similar scores for r and NSE. Whereas S1, S132 and S58 exhibited slightly higher values than S20 and S48 configurations in the case of RMSE and MAE. For PBIAS and PEP, all configurations exhibited relatively smaller model score, and among them, S132 and S58 have slightly higher values followed by S1. Low model score indicates no consensus on the best configuration, i.e., a few configurations performed better for certain events. The model score values for ETP suggest the relatively better performance of S20, S132, S58 and S48 compared to S1 and S8. When the threshold values of 25% and 10% were modified to 15% and 5%, respectively, no significant changes in the model scores were observed.
In addition to verification metrics, streamflow hydrographs were developed and used in the qualitative assessment of streamflow simulations. Figure 8 shows hydrographs of observed and estimated streamflow of reference and calibrated scenarios for six configurations, for 15 events. The plots are arranged in descending order based on magnitude of observed peak flow. Each plot has 13 hydrographs, of which one corresponds to observed streamflow, and two sets each with six hydrographs correspond to streamflow of reference and calibration scenarios. Observed flow S132 S132Ref S58 S58Ref S48 S48Ref S20 S20Ref S8 S8Ref S1 S1Ref Hydrographs of all events suggest that estimated streamflow from all configurations are realistic and exhibit day-to-day variations like observed streamflow. However, the plots also highlighted underestimated peak flow and volume for the two highest flood events, i.e., 1986SP and1966SP, and for a multi-peak flood event, i.e., 2014MP, for all six configurations, for both the scenarios. The large error for the two highest flood events indicates the model's inability to simulate excess rainfall and its translation into direct runoff. For 1986SP, underestimated streamflow is attributed the backwater effect in the Sabari river due to the high streamflow in the Godavari river (main tributary). The backwater effect increased the streamflow at the Konta GD station, but it is not accounted for by the HEC-HMS's regular direct runoff related calculations. Therefore, underestimated streamflow was produced. For the 1966SP, dry AMC resulted in incorrect CN values, which in combination with other parameters produced underestimated streamflow. The time to peak (Tp) was advanced by a day for S1, S8 and S48 (Table A1). Streamflow of the Hydrographs of all events suggest that estimated streamflow from all configurations are realistic and exhibit day-to-day variations like observed streamflow. However, the plots also highlighted underestimated peak flow and volume for the two highest flood events, i.e., 1986SP and1966SP, and for a multi-peak flood event, i.e., 2014MP, for all six configurations, for both the scenarios. The large error for the two highest flood events indicates the model's inability to simulate excess rainfall and its translation into direct runoff. For 1986SP, underestimated streamflow is attributed the backwater effect in the Sabari river due to the high streamflow in the Godavari river (main tributary). The backwater effect increased the streamflow at the Konta GD station, but it is not accounted for by the HEC-HMS's regular direct runoff related calculations. Therefore, underestimated streamflow was produced. For the 1966SP, dry AMC resulted in incorrect CN values, which in combination with other parameters produced underestimated streamflow. The time to peak (T p ) was advanced by a day for S1, S8 and S48 (Table A1). Streamflow of the reference scenario was underestimated for six events, overestimated for another six events, and showed less bias for three events. In a similar proportion, time to peak was either advanced or delayed; however, no association was observed with either under-or over-estimated streamflow. Whereas streamflow of the calibration scenario was relatively better, i.e., exhibited less bias and correctly estimated time to peak, which was reflected in the model score as well. For the two events (1988DP and 2014MP), the time to peak was largely advanced or delayed for most of the configurations. The same was observed for 2006MP for the S8. The error can be attributed to the presence of two or more peaks of the same magnitude at longer days apart and overestimated peak flows for other than the actual peak. However, there are other events of two or more peaks for which most of the configurations have zero error in time to peak e.g., 1966DP, 1975DP and 1984MP.

Results from the Model Validation
The model calibration results suggested a better performance of S1 followed by S132. Given that S1 estimates streamflow only at the basin's outlet and S132 estimates streamflow at both basins' outlet and interior locations, S132 was considered for the model validation. The HEC-HMS model was run for three events at four locations (Konta, Injaram, Saradaput and Potteru) using the ensemble of parameter values that were obtained during the calibration. The parameter set consists of nine parameters, of which a total of five parameters from baseflow (Q 0 , k and R p ) and routing (K and X) methods exhibited no significant variation irrespective of an event. Therefore, the same value of the parameter was used in the validation, i.e., no ensemble of values for these parameters. In addition, the excess rainfall method parameters are events specific, therefore, the parameters are not part of the ensemble set of parameter values. Thus, it results in an ensemble set of parameter values for direct runoff method, i.e., T c and R.
The excess rainfall method parameters, i.e., I a and CN, vary with each event, therefore the parameters were estimated separately for three events. The Clark's UH parameters (T c and R) of all events excluding two extreme events, i.e., 1966SP and 1986SP, were used. Thus, it potentially yields 13 sets of T c and R. However, values for a set of events are the same. In this regard, T c and R values of events, 1968SP, 2015SP, 2005DP, 1988DP, 1975DP and 1969MP were grouped as set 1, and events 1969DP and 1966DP were grouped as set 2. The values for the remaining five events (1972SP, 2006MP, 2014MP, 1984MP and 1978MP) are different and treated as a separate ensemble member. Thus, a total of seven sets of parameters were created.
The ensemble of seven sets of parameters was used to generate an ensemble of streamflow for four locations, for three events (Figure 9). Both single-valued and mean continuous ranked probability score (CRPS) [61] were calculated using ensemble median and all ensemble members, respectively ( Table 7). The mean CRPS considers all members of an ensemble forecast, and it reduces to the mean absolute error (MAE) of the single-valued forecasts. Thus, both mean CRPS and MAE can be compared. Performance of estimated streamflow varied with each event; however, relatively better performance was observed for 2003SP in terms of linear association and NSE (   As the need for the basin-wide real-time forecast is imminent, reliable and accurate runoff estimates at customized locations of a basin including ungauged basins are necessary. The estimation of streamflow for the interior part of a basin requires a scientific un-  In addition to the verification metrics, hydrographs for an ensemble of estimated streamflow were developed (Figure 9). Each ensemble member and all members together showed realistic estimated streamflows. The estimated streamflow exhibited day-to-day variability and produced single-, double and multiple peaks encompassing of observed streamflow for three events for all locations excluding Potteru. Note that the estimated streamflow corresponds to set 1 and set 2 were of the same values; therefore, they are not distinguished in the plot. The estimated ensemble of streamflows exhibited poor performance at the Potteru for events 2003SP and 2004MP; however, they appeared realistic for the event 1997DP with a systematic bias. The poor model performance at the Potteru for events 2003SP and 2004MP could be attributed to the controlled streamflows from upstream at the Jalaput reservoir, which has been operational since 2000. While no systematic pattern was seen for events for the remaining three locations, a few verification metrics suggest good performance with error metrics emphasizing post-processing of streamflows [62][63][64][65].
As the need for the basin-wide real-time forecast is imminent, reliable and accurate runoff estimates at customized locations of a basin including ungauged basins are necessary. The estimation of streamflow for the interior part of a basin requires a scientific understanding of hydrologic processes with respect to resolution and integration of information in the modeling context so that a robust forecast system can be developed. Insights gained from the analysis of parameter values and model performance with respect to sub-basin area, model spatial resolution and events assist in developing the system. The results suggest that the runoff estimates are chiefly driven by excess rainfall and direct runoff components. Equal performance of all configurations in Percentage Error in Peak (PEP, see Tables A1 and A2) agree with [66], who observed smaller differences in the model performance metrics with an increase in spatial resolution for two flood events. However, no consensus was observed among configurations with respect to other verification metrics, and it requires further analysis.
As part of the future work to develop an operational real-time flood forecast system, the following needs to be considered. Sub-basin specific rainfall amounts with respect to events highlight the influence of spatial variability of rainfall. However, it warrants a detailed analysis (e.g., [36,37,51]). Channel routing is an essential component of modeling, and apparently, the Muskingum method did not exhibit any significant role. In this regard, other methods of channel routing need to be explored. In addition, a detailed analysis of rainfall-runoff modeling with finer and coarser resolutions to understand the anthropogenic effects (e.g., land-use change and the role of reservoirs) and integration of hydrologic information as per its availability, including for a short period of the record, is much needed.

Conclusions
The Sabari River Basin was delineated into six configurations, and event-based rainfallrunoff modeling was carried out using the HEC-HMS model. The model calibration and validation were performed for 15 and three flood events of single, double and multiple peaks. The calibrated parameters were analyzed with respect to spatial resolution (subbasin area), flood events and event characteristics, and then performance metrics were calculated for the reference and the calibration scenarios. It was observed that the parameters vary among the sub-basins and all configurations for a given event. Values of the parameters suggest the key role of initial abstractions (I a ), curve number (CN), unit hydrograph related parameters, and less influence of baseflow and channel routing parameters. While I a and CN varied event-to-event without any pattern, T c and R exhibited similar values for most of the events for finer resolution configurations. Being unit hydrograph parameters, T c and R represent the basin's geomorphological characteristics. Therefore, approximate T c and R values for most of the events are agreeable. The performance metrics highlight the better performance of S1 and S132 during the calibration. The hydrographs suggest realistic estimates of streamflow for most of the configurations and events. Noting that S1 estimates streamflow only at the basin outlet, S132 was used to estimate streamflow at four locations as part of the validation. S132 was implemented with a set of seven unique values of calibrated parameters and produced an ensemble of streamflow for three events at four locations. Three of the four locations are interior locations of the basin; thus, the model validation mimics basin-wide streamflow estimates. While performance metrics for a few configurations suggested an agreeable performance, the estimated streamflow was not of similar accuracy as in calibration because of input forecast uncertainties [36,67,68]. As in calibration, streamflow hydrographs suggested realistic flows. However, streamflow estimation may be improved with post-processing including simple techniques for bias correction. Thus, this study shows the estimation of an ensemble of streamflows that can be further enhanced at the basin s interior locations at multiple spatial resolutions. In addition, these types of studies need to be extended in the context of other aspects of modeling, for example, input-including rainfall and model structure-related uncertainties. To address the above limitations, a systematic analysis of rainfall-runoff modeling using other well-established hydrologic models and techniques is necessary, and this study can be treated as one of the steps in this direction. Table A1. Verification metrics for streamflow simulations of the reference scenario; verification metrics from left to right column are r, NSE, RMSE, MAE, PBIAS and PEP for six configurations, and each row corresponds to an event.