Diagnostic calibration of a hydrological model in a mountain 1 area by hydrograph partitioning 2

Hydrological modeling can exploit informative signatures extracted from long time sequences of observed streamflow for parameter calibration and model diagnosis. In this study we explore the diagnostic potential of hydrograph partitioning for model calibration in mountain areas, where meltwater from snow and glaciers is an important source for river runoff (in addition to rainwater). We propose an index-based method to partition the hydrograph according to dominant runoff water sources, and a diagnostic approach to calibrate a mountain hydrological model. First, by accounting for the seasonal variability of precipitation and the altitudinal variability of temperature and snow/glacier coverage, we develop a set of indices to indicate the daily status of runoff generation from each type of water source (i.e., glacier meltwater, snow meltwater, rainwater, and groundwater). Second, these indices are used to partition a hydrograph into four parts associated with four different combinations of dominant water sources (i.e., groundwater, groundwater + snow meltwater, groundwater + snow meltwater + glacier meltwater, and groundwater + snow meltwater + glacier meltwater + rainwater). Third, the hydrological model parameters are grouped by the associated runoff sources, and each group is calibrated to match the corresponding hydrograph partition in a stepwise and iterative manner. Similar to use of the regime curve to diagnose seasonality of streamflow, the hydrograph partitioning curve based on a dominant runoff water source (more briefly called the partitioning curve, not necessarily continuous) can serve as a diagnostic signature that helps relate model performance to model components. The proposed methods are demonstrated via application of a semi-distributed hydrological model (THREW, Tsinghua Representative Elementary Watershed) to the Tailan River basin (TRB) (1324 km 2 ) in the Tianshan Mountains of China. Results show that the proposed calibration approach performed reasonably well. Cross-validation and comparison to an automatic calibration method indicated its robustness.


Background
Parameter calibration has been singled out as one of the major issues in the application of hydrological models (Johnston and Pilgrim, 1976;Gupta and Sorooshian, 1983;Beven and Binley, 1992;Boyle et al., 2000).Commonly, one or more objective functions are selected as criteria to evaluate the similarity between observed and simulated hydrographs (Nash and Sutcliffe, 1970;Brazil, 1989;Gupta et al., 1998;van Griensven and Bauwens, 2003).As model complexity increases, parameter dimensionality also increases significantly, which makes it much more difficult to calibrate model parameters manually.For this reason, automatic calibration procedures have been developed to identify the optimal parameter set (Gupta and Sorooshian, 1985;Gan and Biftu, 1996;Vrugt et al, 2003a, b).However, due to limitations in process understanding and measurement technologies, one can find different parameter sets within a chosen space that may acceptably reproduce the observed aspects of the catchment system (Sorooshian and Gupta, 1983;Beven and Freer, 2001).This phenomenon, which has been called equifinality, causes uncertainty in simulation and prediction (Duan et Published by Copernicus Publications on behalf of the European Geosciences Union.al., 1992;Beven, 1993Beven, , 1996)), and highlights the need for methods that are powerful enough to diagnostically evaluate and correct models, i.e., that are capable of indicating to what degree a realistic representation of the real world has been achieved and pointing towards how the model should be improved (Spear and Hornberger, 1980;Gupta et al., 1998Gupta et al., , 2008)).
Traditional regression-based model evaluation strategies (e.g., based on the use of mean squared error or Nash-Sutcliffe efficiency (NSE) as performance criteria) are demonstrably poor in their ability to identify the roles of various model components or parameters in the model output (Van Straten and Keesman, 1991;Zhang et al., 2008;Gupta et al., 2008;Yilmaz et al., 2008;Hingray et al., 2010), which is due in part to the loss of meaningful information when projecting from the high dimension of the data set (like hydrograph) down to the low (often one) dimension of the measurement (Yilmaz et al., 2008;Gupta et al., 2009).A diagnostic evaluation method should match the number of unknowns (parameters) with the number of pieces of information by making use of multiple measures of model performance (Gupta et al., 1998(Gupta et al., , 2008(Gupta et al., , 2009;;Yilmaz et al., 2008).One way to exploit hydrological information is to analyze the spatiotemporal characteristics of hydrological variables that can be related to specific hydrological processes in the form of signature indices (Richter et al., 1996;Sivapalan et al., 2003;Gupta et al., 2008;Yilmaz et al., 2008).Ideally, a signature should represent some invariant property of the system, be readily identifiable from available data, directly reflect some system function, and be maximally related to some structure or parameter in the model.
Attention to hydrological signatures, therefore, constitutes the natural basis for model diagnosis (Gupta et al., 2008).Placed in this context, the body of literature on the topic is indeed large.Jothityangkoon et al. (2001) proposed a downward approach to evaluate the model's performance against appropriate signatures at progressively refined timescale.Signatures that govern the evaluation of model complexity are the inter-annual variability, mean monthly variation in runoff (called regime curve), and the flow duration curve (FDC).Farmer et al. (2003) evaluated the climate, soil, and vegetation controls on the variability of water balance through four signatures: gradient of the annual yield frequency graph, average yield over many years for each month, FDC, and magnitude and shape of the hydrograph.Shamir et al. (2005a) described a parameter estimation method based on hydrograph descriptors (total flow, range between the extreme values, monthly rising limb density (RLD) of the hydrograph, monthly maximum flow, and negative/positive change) that characterize dominant streamflow patterns at three timescales (monthly, yearly, and record extent).Detenbeck et al. (2005) calculated several hydrologic indices including daily flow indices (mean, median, coefficient of variation, and skewness), overall flood indices (flood frequency, magnitude, duration, and flood timing of various levels), low flow variables (mean annual daily minimum), and ranges of flow percentiles to study the relationship of the streamflow regime to watershed characteristics.Shamir et al. (2005b) presented two streamflow indices to describe the shape of the hydrograph (RLD and declining limb density -DLD) for parameter estimation in 19 basins of United States.Yadav et al. (2007) used similarity indices and hydrological signatures (runoff ratio and slope of the FDC) to classify catchments.Westerberg et al. (2011) selected several evaluation points on the FDC to calibrate models, and compared two selection methods to evaluate their effects on parameter calibration.
Generally, the reported signatures have the following two characteristics: (1) they concentrate on the extraction of hydrologically meaningful information contained in hydrographs, and (2) they focus on either an entire study period or a special continuous section of the entire period.They have occasionally considered temporal variability of runoff components and dominance of different runoff sources during different periods (e.g., the seasonal switching of runoff sources discussed in Tian et al., 2012).However, a hydrograph could be dominated by various components or water sources at different response times (Haberlandt et al., 2001;Eder et al., 2005).With this in mind, a few studies have explored the use of hydrological information in time dimension for stepwise calibration.For example, Schaefli et al. (2005) presented a stepwise calibration method for seven parameters in a high mountainous area: snow and ice melt degree-day factors were conditioned by mass balance, slow reservoir parameters were determined by base flow, reservoir coefficients were calibrated by summer runoff, and the direct runoff coefficient was used to control discharge during precipitation events.Another notable example is Hingray et al. (2010), in which the authors estimated the value of the snowmelt degree-day factor in a mountain basin by progressively minimizing the differences between observed and simulated values of different magnitude hydrographs.There are also many other follow up studies.
In mountain areas, streamflow is composed of both snow/glacier meltwater and rainwater.The energy-based and temperature-index models are two principal approaches to simulate snow and glacier melt (Rango and Martinec, 1979;Howard, 1996;Kane et al., 1997;Singh et al., 2000;Fierz et al., 2003).To describe significant heterogeneity of temperature, precipitation, snow, and glacier, distributed hydrological models are generally used for precipitation-runoff modeling in mountain regions (Daly et al., 2000;Klok et al., 2001).Also, the utilization of remote-sensing products of precipitation and snow cover data in mountain runoff modeling has become more popular in recent years (Swamy and Brivio, 1997;Akyurek et al., 2011;Liu et al., 2012).Most of these studies report sound simulation results.However, the need to develop an appropriate calibration strategy for precipitation-runoff modeling in mountain areas remains a key issue for two reasons: first, hydrological processes are usually more complex (with snow/glacier melt and possibly soil freezing/thawing) than those in warmer areas, which implies a larger dimension of parameter (R P ) in the corresponding hydrological model; second, the measured data set useful for model identification is usually limited due to a sparse gauge network.To address this problem, related studies are putting effort into two particular directions of research.One is to reduce the calibrated R P by estimating some of the parameters based on basin characteristics a priori.For example, Gurtz et al. (1999) proposed a parameterization method based on elevation, slope, and shading derived from basin terrain.Gomez-Landesa and Rango (2002) obtained model parameters of ungauged basins from gauged basins by basin size, proximity of location, and shape similarities.Eder et al. (2005) estimated most of the parameters a priori from basin physiography before an automatic calibration is applied.The parameterization method may involve some uncertainties but is useful for determining insensitive parameters.
The second direction is to exploit hydrological information from implicit measurement data.For instance, Dunn and Colohan (1999) used baseflow data as an additional criteria for model evaluation.Mendoza et al. (2003) exploited recession-flow data to estimate hydraulic parameters.Stahl et al. (2008) used glacier mass balance information combined with stream hydrographs to constrain melt factors.Huss et al. (2008) used annual ice volume change data for optimizing melt and radiation factors, and glacier equilibrium line altitude for precipitation correction factors.Schaefli and Huss (2011) integrated the seasonal information of point glacier mass balance for model calibration by modifying the GSM-SOCONT model.Jost et al. (2012) introduced glacier volume loss calculated by high-resolution digital elevation models (DEMs) to calibrate a hydrologic model.Knowledge acquired from the aforementioned research indicates that the use of additional information (e.g., baseflow, recession flow, and glacier mass balance) can effectively help reduce parameter uncertainty.
However, glacier mass data and baseflow data are usually not available in some mountain basins.In these cases, hydrograph partitioning is another possible way to exploit information from available data.Information about dominant hydrological processes contained in a hydrograph can be extracted by hydrograph partitioning or separation: this has long been a topic of interest in hydrology.Several different kinds of methods have been proposed (Pinder and Jones, 1969;Mc-Cuen, 1989;Nathan and McMahon, 1990;Arnold et al., 1995;Arnold and Allen, 1999;Vivoni et al., 2007), which can generally be classified into graphical methods, analytical methods, empirical methods, geochemical methods, and automated program techniques (Nejadhashemi et al., 2009).Most of them primarily focus on the partitioning of baseflow and are not capable of identifying more than two components.With the advent of isotope methods, multi-component hydrograph separation models have been developed.However, these models need to run for an extended period of time (usually a minimum of one hydrologic year) for the assumption that the isotopes of components are conserved to hold (Hooper and Shoemaker, 1986) and call for volumes of field data that are seldom available in poorly gauged and difficult to access mountain basins.

Objectives and scope
This paper explores the benefits of partitioning the hydrograph into several parts, each related to one combination of dominant water sources for runoff generation.The parameter group controlling each type of runoff generation is then calibrated using the corresponding partitioning hydrographic curves via a stepwise approach, and model deficiencies are diagnosed by evaluating the model simulations associated with each partitioning curve (as a diagnostic signature).We demonstrate the potential of this approach in a mountain area where streamflow is the result of complex runoff generation processes arising from combinations of storm events and snow/glacier melt.The influence of each type of water source (groundwater, snow meltwater, glacier meltwater, or rainwater) varies in time and can be determined by an analysis of the dynamic spatiotemporal information in the available data series.
The paper is organized as follows.Section 2 contains a description of the geographic and hydrological characteristics of the study basin, including the main data sources and data preprocessing.Section 3 details the proposed method of hydrograph partitioning and parameter calibration based on a semi-distributed model coupled with the temperature-index method.Section 4 presents the results and discusses the possible sources of uncertainty.Section 5 provides a summary of this study and discusses further applications of the partitioning strategy.

Overview of the study area
The study mountain area (Tailan River basin, TRB) is on the south slope of the Tianshan Mountains (one of the highest mountain areas in China) in the Xinjiang Uygur autonomous region of China and extends from 41 • 35 to 42 • 05 N and 80 • 04 to 80 • 35 E, covering a drainage area of 1324 km 2 .Elevation ranges from 1600 m to 7100 m a.s.l. with an average height of 4100 m a.s.l.Precipitation occurs mainly in summer and rarely in winter, and winter precipitation always comes in the form of snowfall.Snow coverage accumulates in winter and ablates from spring into late summer when it melts away completely; the snow coverage dynamics can be obtained from MODIS data (see Fig. 4).The basin is highly glacierized with approximately 33 % of the basin area covered by glacier ice (see Fig. 1).The glacier coverage stretches from approximately 3000 to 7100 m a.s.l. and exists mainly at an altitude range of 4000 to 5000 m a.s.l.Glacier melt and snowmelt form runoff as long as the temperature rises above a certain threshold and provide the primary source of downstream discharge.
TRB is a heavily studied mountain watershed in northwestern China.The relevant literature (Kang et al., 1980;Shen et al., 2003;Xie et al., 2004;Gao et al., 2011;Sun et al., 2012) are reviewed below, and the main conclusions about the hydrometeorological characteristics are summarized as follows: 1.The climate presents strong altitudinal variability.The mean annual precipitation in higher mountain areas is approximately 1200 mm (Kang et al., 1980), while it is approximately only 180 mm in the outlet plain area (Xie et al., 2004).The mean annual temperature ranges from below 0 • C in mountain areas to approximately 9 • C at the basin outlet (Sun et al., 2012).
2. Meltwater is the principal source of streamflow.Snow and glacier meltwater account for approximately 63 % of the annual runoff (Kang et al., 1980).The contribution of rainwater is relatively lower and occurs mainly in the storm rain period (SRP) (May to September) (Xie et al., 2004).Groundwater baseflow is smaller but dominates the streamflow in the winter (January, February and December), during which either rainfall or melt rarely occur (Kang et al., 1980).
3. The TRB river network is a simple fan system.Given large topographic drop and moderate drainage area, the runoff concentration time is no longer than 1 day (Xie et al., 2004).Melting and falling water can quickly flow into the main channel and reach the basin outlet.

Data & preprocessing
The Tailan hydrologic station (THS; 1602 m a.s.l.) is located in the outlet of the watershed, where runoff, precipitation, and temperature have been measured since 1957.which indicates that TRB glacier coverage is relatively stable during the study period.The DEM, river system, gauging stations and glacier distribution are shown in Fig. 1.

Temperature lapse rate
Altitudinal distribution of temperature can be estimated through the lapse rate (Rango and Martinec, 1979;Tabony, 1985).According to Aizen et al. (2000), rates of temperature decrease with increasing elevation are quite different in various months, and ignoring this difference may lead to significant errors in the simulation of snow accumulation and melt.The lapse rate was therefore estimated for each month.Temperature variations with altitude can be estimated by the following equation where T o is the temperature value at low altitude (THS in this study), T p is the temperature lapse rate (usually nega- where i indicates the ith day in the analyzed month, T i is the observed temperature in AWS, which is the mean value of the TG AWS and XT AWS in this study. The temperature series data from 1 July 2011 to 31 December 2012 at THS, TG AWS, and XT AWS were used to estimate the temperature lapse rate.The results (Table 1) indicate significant month-to-month variation ranging from −0.30 • C 100 m −1 in December to −0.86 • C 100 m −1 in August.To validate the temperature lapse rates, the estimated and observed temperature data at BT AWS were compared (Fig. 2).We also compared the estimated temperature by an annual constant lapse rate (−0.62 • C 100 m −1 , a similar value to previous studies, e.g., Tabony, 1985;Tahir et al., 2011).This constant value is optimized by the same method in Eq. ( 2) but using all daily temperature measurements.Figure 2 indicates that the monthly lapse rate method performs better than the annual constant rate method at the BT station for all months throughout the year.Further, the temperature curves estimated by monthly lapse rates for April to August match the observed ones rather well.Note that the estimated temperatures tend to underestimate observed ones for the rest of the months, which, however, will not affect the melt runoff significantly due to the general freezing condition during this period.

Precipitation lapse rate
Based on the precipitation series measured at THS, the monthly precipitation to annual precipitation ratio (Fig. 3) for the study period (2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012) indicates that precipitation  occurs mainly in May to September.The lapse rate of precipitation was also estimated monthly, and a similar procedure as temperature was applied.The difference is that the precipitation analysis was conducted at a weekly rather than daily time step, and the maximum measured precipitation of the two installed AWS was used instead of the mean value.The analyzed period is limited to the SRP (May to September).Other months are not included due to the relatively small amount of precipitation.The weekly precipitation lapse rates are listed in Table 2. Daily precipitation differences between higher and lower altitudes can be estimated as the weekly precipitation lapse multiplied by the ratio of daily precipitation to the corresponding weekly amount in THS.The precipitation lapse rate was not validated against BT AWS because of significant differences in precipitation distribution between the two basins (i.e., Tailan and Kumalak).

Filtering of MODIS snow cover area data
Snow cover extent was obtained from MODIS products.The MOD10A2 and MYD10A2 products were downloaded from http://reverb.echo.nasa.gov.In total, we obtained 460 8-day images (two tiles, h23v04 and h24v04) from 2003 to 2012 for  2005).The term mod is the snow cover area from MOD10A2 products, myd is MYD10A2 products, combined is the combined result from step1, spatial-comb from step2 and temporal-comb from step3.See Sect.2.2.3 for details.
each product.Given that the accuracy of the MODIS SCA product is affected by cloud coverage to a significant degree, the remotely sensed images should be filtered to avoid the noise from clouds before using it for hydrological modeling (Ackerman et al., 1998).The following three successive steps are adopted to filter the products based on previous reports (Gafurov and Bardossy, 2009;Wang et al., 2009;Lopez-Burgos et al., 2013): 1. Satellite combination: the snow cover products of two satellites, Terra (MOD10A2) and Aqua (MYD10A2), were combined.As long as the value of a pixel is marked as snow in either satellite, the pixel value is marked as snow.
2. Spatial combination: inspecting the values of the nearest four pixels around one center pixel marked as cloud, if at least three of the four surrounding pixels are marked as snow, the center pixel is modified as snow.
3. Temporal combination: if one pixel is marked as cloud, its values in the previous and following observations are investigated.If both of the two observed values are snow, then the present value of the same pixel is snow.
As an example, the filtered results from 2004 to 2005 shown in Fig. 4   when the snow has been melt away completely.The filtered results indicate a relatively stable coverage of glacier in TRB.

Altitudinal cumulative melt curve
The daily temperature of each cell in MODIS SCA images can be estimated by a temperature lapse rate based on its elevation and daily temperature measured at THS.As long as the temperature exceeds a specific threshold value for melt (assumed to be 0 • C in this study), a given cell was labeled as an active cell in terms of melt.The land cover type for each cell was classified into glacier, snow, and other land cover according to the CGI and MODIS SCA product.To obtain the area covered by snow only, we subtracted the glacier area in CGI from the SCA (a similar procedure can be found in Luo et al., 2013).When a glacier or snow cover cell is active, it is labeled as a melt cell, and the melt area is computed as the number of active cells multiplied by the area of a cell.
Organizing the melt area by elevation from low to high and summing the melt area at each elevation, we can get the altitudinal cumulative melt curve, which can be used to describe the spatiotemporal distribution of melt area.The altitudinal cumulative melt curves calculated from 2003 to 2012 for all months (Fig. 5) show that melt mainly occurs from May to September, which coincides with the precipitation period.Snowmelt starts at an elevation of approximately 1650 m a.s.l., while glacier melt starts at an elevation of approximately 2950 m a.s.l., which has an important implication for hydrograph partitioning.

Methodology
Theoretically, every drop of water in the streamflow comes ultimately from precipitation.Practically, we can consider water sources for runoff generation in mountain areas as mainly consisting of meltwater from snow and glacier, rainwater, and groundwater.Groundwater at the basin scale is recharged by direct infiltration and run-on infiltration of meltwater or rainwater, and it is mainly discharged as baseflow via a subsurface flow path (especially in mountain areas where the large elevation gradient favors baseflow discharge).For the purpose of hydrograph partitioning, we can consider recharge to be a separate water source for streamflow, independent of meltwater and rainwater, which principally forms the baseflow part of a hydrograph.The remaining part of a hydrograph is principally formed by meltwater and rainwater via surface flow path (Blöschl et al., 2013).We develop three indices to indicate the water sources for runoff generation at the daily timescale.The hydrograph is further partitioned into several sub-parts based on the indices' values.Each sub-part is dominated by one or more water sources for runoff generation.With the partitioning hydrographic curves, the parameters of hydrological models are correspondingly grouped by runoff sources and calibrated in a stepwise fashion.We use the Tsinghua Representative Elementary Watershed (THREW) model coupled with a temperature-index module as an exploratory tool.To better demonstrate usefulness of the proposed methods, only the runoff generation related parameters, which are also significantly sensitive parameters (see Sect. 4.6), are calibrated.Other insensitive parameters are fixed at their initial values, specified a priori from the literature or by expert knowledge.

An index-based method for hydrograph partitioning
In mountain areas, the relative contribution of different runoff water sources to the total streamflow varies throughout the year (Martinec et al., 1982;Dunn and Colohan, 1999;Yang et al., 2007).For rainwater sources, Fig. 3 shows that precipitation in TRB presents strong seasonality and is primarily concentrated (more than 76 %) in the SRP from May to September.During the relatively dry period from October to April, mean precipitation gauged at the THS is just 43 mm, while precipitation in higher mountainous regions is mainly snowfall.Therefore, surface runoff induced by rainwater can rarely occur during relative dry period.It is reasonable to assume that the rainwater source can only contribute to the surface runoff part of a hydrograph on the same day during the SRP (May to September) except for the baseflow occurring much later.
For the meltwater sources, the altitudinal cumulative melt curves (Fig. 5 months, i.e., glacier melt can only occur in the areas higher than 2950 m (the lower elevation limit of glacier coverage), while snowmelt can occur in areas higher than 1650 m.It can be deduced that snowmelt generally occurs at lower elevations than glacier melt.Remember that temperature decreases with increase in altitude.There should exist a period of time during which temperature at 1650 m is higher than the snowmelt threshold, while temperature above 2950 m is lower than the glacier threshold and thus snowmelt occurs but glacier melt does not.
The groundwater source should be a dominant source for the baseflow part of a hydrograph and, of course, it dominates the recession limb of a hydrograph (part of a baseflow partition) when no rainfall or melting occurs.
Based on the above physical understanding, we can partition the hydrograph using the following three indices: 1. Date index (D i ): D i is used to distinguish the dates on which rainfall and thus possible rainwater direct runoff processes occur.For simplicity, in this study we use D i to distinguish the dry period and the SRP and assume no rainfall runoff in the dry period, i.e., 1, for days in SRP from May to September 0, for days in the relative dry period from October to April.
(3) 2. Snowmelt index (S i ): S i indicates whether snowmelt possibly occurs on a given day: 1, for days when temperature at altitude 1650 m is higher than 0 • C 0, for other days.
(4) 3. Glacier melt index (G i ): G i is used to identify days when glacier melt possibly occurs: 1, for days when temperature at altitude 2950 m is higher than 0 • C 0, for other days. (5) The hydrograph is then partitioned according to the three indices by using the following rules: and where Q is the overall streamflow series, Q SB stands for the baseflow generated by groundwater source, Q SM for snow meltwater runoff, Q GM for glacier meltwater runoff, and Q R for rainwater direct runoff.The partitioning principles are described as follows: when both meltwater and rainwater direct runoff do not occur.This condition requires S i = 0, G i = 0, and D i = 0.
2. Snow meltwater and groundwater are the dominant components (Q = Q SB +Q SM )when the temperature is higher than 0 • C at 1650 m a.s.l. and lower than 0 • C at 2950 m a.s.l.(requires S i = 1, G i = 0, and D i = 0).
3. Snow meltwater and glacier meltwater coupled with groundwater dominate (Q = Q SB + Q SM + Q GM ) on days when the temperature at 2950 m a.s.l.exceeds 0 • C in October to April.This means G i = 1, D i = 0, and S i = 1, noting that S i must be equal to 1 when G i = 1 for the decreasing nature of temperature along altitude.
4. Finally, all sources are mixed (Q for other days in the SRP (May to September, D i = 1).Each category contains days that could be continuous or discontinuous in time and could lie within different weeks due to temporal variability of precipitation and temperature.

Tsinghua representative elementary watershed hydrological model
The THREW model used for the hydrological simulation in this study, has been successfully applied in many watersheds in both China and the United States (see Tian et al., 2008Tian et al., , 2012;;Li et al., 2012;Liu et al., 2012), including application to a high mountainous catchment of Urumqi River basin by Mou et al. (2008).The THREW model adopts the Representative Elementary Watershed (REW) approach to conceptualize a watershed, where REW is the sub-catchment unit for hydrological modeling.The study basin was divided into several units (REW) based on a DEM.Sub-catchment units were further divided into a surface and subsurface layer, each layer containing several sub-zones.The subsurface layer is composed of two zones: saturated zone and unsaturated zone, and the surface layer consists of six zones: vegetated zone, bare soil zone, snow covered zone, glacier covered zone, substream-network zone, and main channel reach; see Tian et al. (2006) for further details.
The main runoff generation processes simulated by the THREW model include rainfall surface runoff, groundwater baseflow, snowmelt, and glacier melt.Rainfall surface runoff is simulated by a Xinanjiang module, which adopts a water storage capacity curve to describe non-uniform distribution of water storage capacity of a sub-catchment (Zhao, 1992).The storage capacity curve is determined by two parameters (spatial averaged storage capacity W M and shape coefficient B).Rainfall surface runoff forms on areas where storage is replete.Replete areas are calculated by the antecedent storage and current rainfall.The saturation excess runoff is computed based on water balance.(Jiang, 1987;Yang et al., 1987;Mu and Jiang, 2009), which indicates its applicability in our study area.
For the simulation of melt processes in this study, the THREW model was modified to couple with the temperature-index method, given the easy accessibility of air temperature data and generally good model performance of the temperature-index model (Hock, 2003;Singh et al., 2000).Snowmelt and glacier melt are simulated using separate degree-day factors (snowmelt degree-day factor D s and glacier melt degree-day factor D g ).Glacier melt only occurs in glacier areas according to CGI, which remains stable during the study period (2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012); see discussion in Sect.2.2.3).Precipitation in the snow and glacier zone is divided into rainfall and snowfall according to two threshold temperature values (0 and 2.5 • C are adopted in this study according to Wu and Li, 2007), i.e., when temperature is higher than 2.5 • C, all precipitation is rainfall, when temperature is lower than 0 • C, all precipitation is snowfall, and when temperature falls between the two thresholds, precipitation is divided into half rainfall and half snowfall (a simple division scheme adopted here).Rainfall on glacier areas forms runoff and flows into the stream network directly without infiltration into soil.The snow water equivalent (SWE) on glacier areas is updated by combining snowfall and snowmelt, and for simplicity, snow is assumed to cover all glacier areas when the corresponding SWE is not zero.Snowmelt in glacier areas is simulated using snow degree-day factor D s until it melts away completely.The snow cover area in a non-glacier area is updated using MODIS data.To be noted, snowfall in each sub-catchment is calculated according to the daily precipitation and temperature, and snowmelt is simulated using the degree-day method.However, the SWE in the snow cover zone (non-glacier area) is not computed.The existing snow cover in each sub-catchment is only determined by MODIS snow image.When the MODIS image indicates the existence of snow cover and at the same time the daily temperature is higher than 0 • C, then snowmelt will occur; otherwise, snowmelt will not occur.The identification of snow cover by the MODIS image is in accordance with the fact that the partitioning of a snowmelt dominant hydrograph is based on MODIS snow products.If the existence of snow cover is determined by the SWE, the temperature parameters to calculate snowfall can have significant effects on the estimation of the degree-day factor for snowmelt.To partly reduce this effect, we calibrate the degree-day factor for snowmelt on the basis of MODIS snow cover products.Although in this way, the water balance of snow cover is not taken into account in the snow cover zone: it should not impact the calibration of the degree-day factor for snowmelt.Since MODIS SCA products (i.e., MYD10A2) have been available since 2003, the model simulation period is from 2003 to 2012, of which 2003-2007 is for calibration and 2008-2012 is for evaluation.The time step for simulation is daily.

Stepwise calibration of grouped parameters on partitioning curves
Model parameters are grouped a priori according to their connection with causal physical mechanisms (see Table 3).
According to Xie et al. (2004) and Kang et al. (1980), parameters that control groundwater baseflow, snowmelt, glacier melt, and rainwater surface runoff should be the most sensitive parameters for the runoff simulation (also see our sensitivity analysis in Sect.4.6).These parameters are subjected to calibration in this study.They are related to the corresponding hydrograph parts and then calibrated in a stepwise manner.First, groundwater baseflow parameters (K A and K D ) are estimated based on the Q SB part of the hydrograph.Second, the snowmelt degree-day factor (D s ) is calibrated on the Q SB + Q SM part.Third, glacier melt degree-day factor (D g ) is determined according to the Q SB + Q SM + Q GM part.Finally, rainfall surface runoff parameters (B, W M ) are calibrated on days when D i equals to 1, i.e., the In each step, only the specific parameter group is subjected to calibration.The parameters determined in the previous steps are kept constant, and all other parameters that will be calibrated in the next steps adopt their initial values.As the simulation in each step can, to some degree, be affected by the initial conditions produced in the preceding step, an iterative procedure is implemented to progressively minimize this influence.The parameter groups are first calibrated based on the corresponding hydrograph parts, and then the stepwise sequence is repeated until the calibrated parameters converge, i.e., the difference in parameter values between two contiguous iterations is less than 10 %.In each calibration step, we use RMSEln (Eq. 7, emphasizing low flow) or RMSE (Eq.8, emphasizing high flow) as an objective function for parameter optimization.The remaining, insensitive, parameters are determined a priori according to previous modeling experience (mainly from Sun et al., 2012) and listed in Table 3.The initial values of the calibrated parameters are also determined a priori according to Sun et al. (2012) and Tian et al. (2012).
The overall streamflow can be simulated with all calibrated parameters, which is evaluated with the NSE and NSEln (logarithm Nash criterion) values.Given that it is relatively easier to obtain high-evaluation merit values in snowmeltdriven basins due to strong seasonality of streamflow, we further adopt a simple benchmark model (the inter-annual    9) for detail.
4 Results and discussion

Partitioning hydrographic curves
The hydrograph from 2003 to 2012 was partitioned based on Eq. ( 6).In total, we obtained four kinds of partitioning The numbers of non-melt days (i.e., the Q SB part, due to glacier melt, generally occurs in the Q SB +Q SM +Q GM +Q R part) in each year of this 5 year period is 114, 80, 89, 96, and 68; correspondingly, the mean temperatures in those years gauged at the THS are 8.9, 10.1, 9.9, 10.4, and 11.3 • C. A lower mean annual temperature causes a longer non-melt period in that year and vice versa.Note that the partitioning curves can be discontinuous in time due to the spatiotemporal variability of temperature.

Model calibration by the stepwise method
The six key parameters (K A , K D , D s , D g , W M , and B) were initially calibrated by the proposed stepwise and iterative method.To focus on baseflow generated by the groundwater source during the Q SB period, the RMSEln metric that emphasizes low flow is chosen as the evaluation criterion for the calibration of parameters K A and K D .Conversely, high flow is our focus for the remaining periods ( and the RMSE metric is chosen as the evaluation criterion for calibration of parameters D s , D g , and W M and B. To deal with interaction between steps, an iterative calibration approach was adopted.A total of five iterations were implemented until the parameter estimates became stable; the simulation of each kind of partitioning curve in each step of the last iteration is presented in Fig. 7.The calibrated parameters are shown in Table 4 and the evaluation merits are listed in Table 5.
Figure 7a shows that the magnitude of baseflow in the Q SB part was captured well at most of the times.The RMSEln merit is 0.302 m 3 s −1 , and the parameters K A and K D are 1.1 and 0.002, respectively.Streamflow in the Q SB + Q SM part is dominated by both snow meltwater and groundwater.Figure 7b shows that melt peak flow events have also been captured well by a calibrated D s as 2.5 mm •−1 C day −1 after the determination of K A and K D in the first step.For the Q SB + Q SM + Q GM part, glacier meltwater began to control the streamflow in combination with snow meltwater and groundwater.Snowmelt and baseflow were determined a priori by previously calibrated parameters.The remaining residual between the simulated and observed discharge can be attributed to glacier melt alone, which was thus used for the calibration of glacier melt factor D g .The RMSE value for this hydrograph partition was optimized as 4.784 m 3 s −1 and we obtained a sound simulation by a calibrated D g as 7.2 mm • C −1 day −1 as shown in Fig. 7c.During the SRPs ( , rainwater direct runoff is an additional important component of river runoff.Similarly, parameters W M and B can be calibrated separately after a priori determination of melt runoff and groundwater baseflow.The simulated RMSE value in this period is 12.650 m 3 s −1 , with calibrated W M = 10.50 cm and B = 0.80.The overall daily streamflow simulation is obtained by combining the four partitions together (see Fig. 8a).The corresponding NSE index is 0.881 and NSEln is 0.929.Generally, the results suggest a sound simulation compared to the observation.
To be noted, the calibrated values of melt degree-day factors D s (2.5 mm • C −1 day −1 ) and D g (7.2 mm • C −1 day −1 ) are similar to the values obtained in other studies in the Tainshan area, e.g., D s is calibrated as 2.5 mm • C −1 day −1 by Liu et al. (2012), and D s and D g are estimated as 3.1 and 7.3 mm • C −1 day −1 , respectively, based on observed mass balance data by Liu et al. (1999), which indicates the robustness of our calibration method.

Comparison to automatic calibration method
For comparison, we also carry out an automatic calibration with the help of the ε-NSGAII algorithm, an optimization method developed by Deb et al. (2002) and Kollat and Reed (2006).The six parameters were calibrated together and evaluated by the NSE value of the overall hydrograph.The run time of the automatic algorithm is about 5 weeks (840 h on a desktop equipped with an Intel Core i7 CPU with 2.8 GHz).The NSE value for the final optimized parameters is 0.868, and the NSEln value is 0.846 (Fig. 8b), both of which are lower than the values obtained by the proposed stepwise method.The parameters calibrated by ε-NSGAII are listed in Table 4, and are different from those calibrated by the stepwise method.Specifically, the snowmelt degreeday factor (D s ) and groundwater baseflow parameters (K A and K D ) obtained by ε-NSGAII are 2.03 mm • C −1 day −1 and 5.6 and 99.1, respectively.The evaluation merits of RMSE and RMSEln for each partitioning curve are also shown in Table 5.In general, the simulation by the automatic algorithm is not as good as that by the stepwise method, especially for the low and middle flow partitions (Q SB + Q SM and Q SB +Q SM +Q GM ).This may be due to the tendency of NSE-based automatic calibration to emphasize high flows.To make a further evaluation, a benchmark model suggested by Schaefli and Gupta (2007) is used for the comparison, which simply simulates daily runoff as the interannual daily mean value.Simulation results by the benchmark model are shown in the Fig. 8c, which shows the NSE value as 0.815 and NSEln value as 0.923.The high NSE and NSEln values can be attributed to the strong seasonality of stream discharge in the study basin (Schaefli and Gupta, 2007).The BE index (Eq.9, see Table 5) is used to measure the improvement of simulations by the calibration methods compared to the benchmark model.A positive value for BE means that the evaluated method outperforms the benchmark model.Figure 8 shows the simulations of daily streamflow by the three methods (Fig. 8a by stepwise calibration method, Fig. 8b by automatic calibration method, and Fig. 8c by benchmark model), which shows better simulation by the two calibration runs with the THREW model than the benchmark model (BE values are both positive).The stepwise calibration run obtained a BE value of 0.355, while the BE of the automatic calibration run is 0.271.The benchmark model describes the mean value of daily discharge on each calendar day.The higher the BE value is, the better the seasonal variability of the hydrograph is captured by the evaluation method.The higher BE value in the stepwise calibration method can be attributed to the better simulation of middle and low flows which are dominated by groundwater and meltwater (Fig. 8a).However, BE values simulated by two calibrated parameter sets are both relatively low, which is attributed to the poor mimic of the (rapidly rising and falling) peaks.
Note that the automatic calibration method based on the NSE value of the overall hydrograph adopts 1-D measurement information to optimize four parameter groups.Benefitting from the partitioning curves, however, the stepwise calibration method increases the dimension of hydrological signature to four.The signature dimension is now equal to the number of parameter groups, and the grouped parameters can be optimized according to their corresponding runoff components separately.A sound simulation of the overall hydrograph is obtained by the reasonable reproduction of the separate partitioning curves.Therefore, parameters calibrated by the stepwise method are inclined to have a more explicit physical basis.
In regards to computation efficiency, the stepwise calibration required 385 runs of the model to be completed, with each model run taking about 1.5 min and the total computation time being about 10 h.In contrast, the state-of-the-art automatic calibration algorithm required about 5 weeks of CPU time consumption on a desktop equipped with an Intel Core i7 CPU and 2.8 GHz.The comparison indicates that the stepwise calibration method is both more physically based as well as more computationally efficient.
It is worth noting, the performance of the automatic calibration algorithm can increase if the algorithm keeps on running, and even be higher than that of the stepwise calibration method.The comparison here is intending to show that the stepwise calibration method based on hydrograph partition can achieve considerable performance more effectively.The automatic algorithm here treats all the parameters equally during the calibration period.Each parameter should be optimized when searching for the optimal parameter set.This searching algorithm hampers the efficiency of the calibration procedure without identifying the dominant sub-periods for different parameters.In the stepwise calibration method, only parameters that are responsible for the simulation of corresponding hydrograph partition are optimized in each step.Furthermore, the calibration of the parameter by this method reflects the role of each parameter for the basin runoff generation.

Evaluation for the stepwise calibration method
The parameter set calibrated by the stepwise method is applied to the evaluation period (2008)(2009)(2010)(2011)(2012), and the daily discharge simulation is shown in Fig. 9a.The evaluation merits are listed in Table 5.The NSE, NSEln, and RMSE values for the whole period indicate sound evaluation results but generally lower performance compared to calibration pe- riod.However, the evaluation results by the stepwise method are still significantly better than the benchmark model, which obtained a NSE value as low as 0.577 (Fig. 9b and Table 5).The BE value in evaluation period by the stepwise calibration method is 0.413.Furthermore, from the partition perspective, the RMSEln and RMSE values for four partitions in Table 5 show that the low flow simulations (Q SB , Q SB + Q SM , and Q SB + Q SM + Q GM parts) are pretty good and even outperform the calibration simulations.The high flow simulation (Q SB +Q SM +Q GM +Q R part) is, however, insufficient, with RMSE 16.727 m 3 s −1 (compared to 12.65 m 3 s −1 in calibration period).The lower performance of overall evaluation should be attributed to the insufficiency in storm rain days, especially for some extreme storm events in the summer of 2010 (see Fig. 9a).The underestimation of these events is likely due to inadequate observations of rainfall, which are principally due to the strong spatial variability of rainfall in mountainous areas.It is widely acknowledged that the extreme runoff events are difficult to capture in mountain areas, where gauged stations are scarce, on the daily scale (Aizen et al., 2000;Jasper et al., 2002).However, the accuracy of our results is similar to Li and Williams (2008) (used the SRM model) and Liu et al. (2012) (who used the MIKE-SHE model) who performed similar work in a basin that is close to TRB in the Tianshan Mountains.Their Nash values for daily discharge varied from 0.51 to 0.78, and also failed to simu-late the peak flows in summer.They also attributed the low efficiency to the heavy precipitation.
To further evaluate the robustness of the stepwise calibration method based on partitioning curves, cross-validation was implemented.The hydrograph in the evaluation period was partitioned based on dominant runoff sources, as was done in the calibration years 2003-2007.

Sensitivity analysis on the index-based partitioning method
The stepwise calibration method relies heavily on the hydrograph partition for different runoff components.The in-dices defined in Sect.3.1 are key to identifying the dominant days for meltwater and rainwater.The definitions for elevation bands for the 0 • C isotherm and for storm rain days in the year producing rainwater runoff should have significant influence on the parameter calibration.In this study, the elevation band of 0 • C isotherm for snowmelt is fixed and defined as 1650 m.This value should have minimal effect on the snowmelt simulation, as the occurrence of snowmelt is actually determined by the MODIS snow cover data.The glacier cover area is assumed as constant, which is very rough as we have only one CGI data.In this section, we define different elevation bands of 0 • C isotherm for glacier to analyze the effect of glacier area variation on the model calibration.We also select different seasons as the SRP to analyze its sensitive effect.
According to the CGI data, the glacier area extends from the altitude of 2950 m in 2002.Considering the possible variability, we define four different lowest elevation bands for the glacier area (LEGs), i.e., −500 m (2450 m), −200 m (2750 m), +200 m (3150 m), and +500 m (3450 m).As an example, various hydrograph partition patterns in year 2003 are shown in Fig. 10.For the SRP, new seasons are defined as April to October, April to September, May to October, and June to August compared to the benchmark period May to September.A new hydrograph partition pattern in year 2003 is also shown in Fig. 10.The left column in Fig. 10 shows that the Q SB + Q SM + Q GM partition becomes longer while the Q SB + Q SM partition becomes shorter when the LEG is lower.Therefore, glacier melt starts earlier and ends later in the years with lower LEGs.In the right column, the Q SB +Q SM +Q GM partition becomes longer with the shorter SRP, while the variation of the Q SB + Q SM partition can be negligible.Parameters were re-calibrated according to the new partition curves, and the results are shown in Table 6, indicating the increase of degree-day factor for glacier melt (D g ) with the increase of the LEG.The value of D g is also found to become higher when the SRP falls in the warmer months.The variation of LEG imposes significant impacts on the calibration of D g , with a result ranging from 5.8 to 8.0 mm • C −1 day −1 , while the variation of the SRP principally impacts the calibration of parameter W M , with a result ranging from 8.2 to 10.5 cm.However, the NSE values (see Table 6) for different settings show minimal differences.This can be attributed to the fact that parameters are optimized on separate partitioning curves in the stepwise calibration method.Each hydrograph partition can be well simulated by adjusting the parameter values.The partition patterns can influence the value of parameters significantly but only slightly influence the discharge simulation.Among various LEGs, the setting of 2950 m leads to the highest NSE value.The glacier melt degree-day factor (D g ) calibrated with this LEG is 7.2 mm • C −1 day −1 , which is very close to the value estimated as 7.3 mm • C −1 day −1 by Liu et al. (1999), in which the D g is estimated according to the observed glacier mass balance data in the Tianshan area.This can further demon-strate the reasonability of the assumption in Sect.3.2 that the glacier area is stable and its lowest elevation is fixed at 2950 m during the study period.For the various SRPs, when the May to October period is adopted, the discharge simulation is slightly better than the benchmark setting of the SRP, i.e., May to September.This phenomenon seems to indicate the importance of precipitation measurement as discussed in Sect.4.4.Given that the hydrograph partition in Fig. 6 is on the basis of setting the SRP as May to September, some small rain events in April are not taken into account.Sensitive analysis in Table 6 indicates that taking these events into account (i.e., defining SRP as April to October and April to September), the calibrated value of parameter W M can be significantly different.With the help of more advanced precipitation measurement, the SRP can be determined more precisely to improve the model simulation.
To evaluate the relative dominance of multiple runoff components on the total runoff, we compute their contributions to total runoff by various LEGs and SRPs in Fig. 11.The mean contributions of every runoff component are as follows: groundwater contributes 17 %, snow meltwater contributes 16.5 %, glacier meltwater contributes 40 %, and rainwater direct runoff contributes 26.5 %.Total meltwater (snowmelt and glacier melt) occupies approximately 56.5 % and is close to the ratio of 63 % suggested by Kang et al. (1980).

Sensitivity analysis on parameters
The number of parameters to be calibrated is determined by the parameter sensitivity and a priori analysis.To evaluate the effect of different parameters on the simulation of different hydrograph partitions, we implemented a simple parameter sensitivity procedure that is carried out by a one-at-atime approach.Parameters from different groups in Table 3 are selected for sensitivity analysis, including saturated hydraulic conductivity for unsaturated zone K u s , saturated hydraulic conductivity for saturated zone K s s , subsurface flow coefficient K A and K D , manning roughness coefficient for hillslope n t , spatial heterogeneous coefficient for infiltration capacity α IFL , ground surface depression storage capacity Fmax b , shape coefficient to calculate the saturation excess runoff area from the Xinanjiang model B, spatially averaged tension water storage capacity in the Xinanjiang model W M , glacier degree-day factor D g , and snowmelt degree per day factor D s .Parameters varied from −50 to +50 % of the calibrated values using the stepwise method in Table 4.The relative change (R MS ) of simulated measurement merits (RM-SEln or RMSE) for different hydrograph partitions are used to evaluate the sensitivity (Eq.10), where MS is the value of measurement merits by the calibrated parameter, MS + is the merits value obtained by the parameter +50 % of the calibrated one, and MS − is the merits value obtained by the parameter −50 % of the calibrated one.The sensitivity simulation results are shown in Table 7, which demonstrates the dominant control of parameter K A , K D , W M , B, D s , and D g .Some parameters have significant effects on simulation of multi-hydrograph partitions.For example, parameters controlling the Q SB + Q SM + Q GM + Q R period can also have significant effect on the other periods.To minimize this interaction, iterative calibration was implemented in the calibration procedure.The number of calibrated parameters is determined as six, which control the main runoff components (i.e., groundwater baseflow, snowmelt, glacier melt, and rainwater direct runoff).Note that the low dimension of parameter calibration should not account for the low efficiency of peak flow simulation, referring to the similar study in the Tianshan Mountains by Li and Williams (2008) and Liu et al. (2012), in which the models have a higher parameter dimension (higher than six), and the peak flow simulations are still inadequate.

Summary and conclusion
This study proposes a diagnostic calibration approach to extract hydrological signatures from available data series in a mountain area, which can be further used to partition the hydrograph into dominant runoff sources.The parameters  of a hydrological model were grouped according to runoff sources and then related to the corresponding hydrologic partitioning curve.Each parameter group was calibrated to improve the simulation of the corresponding partitioning curve in a stepwise way.In this way, the dimension of the hydrological signature is expanded to equal the number of parameter groups.The parameter uncertainty due to interaction of parameters is reduced via an iterative calibration procedure.Application to a mountain watershed in the Tianshan Mountains in northwestern China showed that the approach performed reasonably well.Cross-validation and comparison to an automatic calibration method indicated its applicability.Note that a semi-distributed hydrological model was utilized to illustrate the proposed diagnostic calibration approach in the high mountainous TRB.Glacier mass balance is not simulated in the model and the glacier coverage was kept fixed during the study period, which can be subject to significant change in the context of global warming.According to existing studies (Stahl et al., 2008;Schaefli and Huss, 2011;Jost et al., 2012), glacier mass balance data were useful to constrain the parameter uncertainty for hydrological modeling in a glaciered basin.While arguing that our assumption of unchanged glacier coverage will not weaken the importance of the proposed approach, we acknowledge that an improved model coupled with glacier mass balance equations will improve the accuracy of hydrological simulation aided by glacier mass balance observations.This is left for future research.
A prerequisite for the proposed approach is hydrograph partitioning based on dominant runoff sources.The key to the partition procedure is to identify the functional domain of each runoff source from signature information extracted from easily available data.A partition can be achieved in which the relative roles of different runoff components in the basin runoff vary significantly with time.The mountain watershed is an area in which the runoff source can be separated by the combination of topography, ground-gauged temperature and precipitation, and remotely sensed snow and glacier coverage.Other areas with strong temporal variability of catchment wetness along with precipitation (e.g., monsoon zones) could also be suitable for the proposed approach.The Dunne runoff is prone to dominate the hydrograph when the catchment is wet and it could switch to Hortonian runoff rapidly under the combination of high evaporative demand and less precipitation, as shown by Tian et al. (2012) in the Blue River basin of Oklahoma.This is, however, also left for future research.

Z
. H.He et al.: Diagnostic calibration of a hydrological model

Figure 1 .
Figure 1.Location of the Tailan River basin in the Xinjiang Uygur autonomous region, China.Two automatic weather stations (TG is the Tagelake automatic weather station (AWS) at 2381 m a.s.l. and XT is the Xiaotailan AWS at 2116 m a.s.l.) were set up in an upstream mountain area in July 2011.Additionally, the BT is the Bingtan AWS (3950 m a.s.l.) located in the adjacent Kumalak River basin was used to validate the estimated temperature lapse rates.The Tailan hydrologic station (THS) has gauged streamflow data at the catchment outlet since 1957 (a).Glacier occupies approximately 33 % of the total basin area (b).

Figure 2 .
Figure2.Evaluation of the estimated temperature lapse rate at the BT station.The black solid line is the observed temperature series at BT (Obs.tem); the red solid line is the estimated temperature by monthly lapse rate (Mrate.tem).The red dotted line indicates the estimated temperature based on annual constant rate (Yrate.tem).The goodness of fit between the observed and estimated temperature is measured by RMSE M for monthly lapse rate and RMSE Y for annual constant rate, respectively.The temperature series in September and October are absent at BT.

Figure 3 .
Figure 3. Proportion of monthly precipitation to annual amount (2003-2012).The red line in each box represents the median value for each month from 2003 to 2012.Red crosses indicate abnormal values that exceed 1.5 times the inter-quartile range.

Figure 5 .
Figure 5. Altitudinal Cumulative Melt Curve: (a) cumulative monthly snowmelt area distribution by elevation (2003-2012); (b) cumulative monthly glacier melt area distribution by elevation (2003-2012).The snowmelt areas in December and January and the glacier melt areas in November, December, January, and February are zero and are not shown in this figure.
) show that the areas experiencing glacier melt and snowmelt change significantly with elevation.Melting of glacier and snow begins at different elevations in different www.hydrol-earth-syst-sci.net/19/1807/2015/ Hydrol.Earth Syst.Sci., 19, 1807-1826, 2015 • C −1 day −1 Glacier melt degree-day factor Calibrated D s mm • C −1 day −1 Snowmelt degree-day factor Calibrated mean value for every calendar day) to evaluate the performance of the proposed method by subtracting streamflow seasonality.This benchmark model is proposed by Schaefli and Gupta (2007) for basins having a relatively constant seasonality.The improvement of a model compared to the benchmark model is quantified by the BE index; see Eq. (

Figure 6 .
Figure 6.Hydrograph partition in 2003.Q SB stands for subsurface baseflow generated by groundwater, Q SM and Q GM for snow meltwater and glacier meltwater, respectively, and Q R for rainwater direct runoff.

Figure 7 .
Figure 7. Stepwise calibration of grouped parameters on partitioning curves.(a) Partitioning curves after calibrating K A and K D on Q SB .(b) Partitioning curves after calibrating D s on Q SB + Q SM .(c) Partitioning curves after calibrating D g on Q SB + Q SM + Q GM .(d) Partitioning curves after calibrating W M and B on Q SB +Q SM +Q GM +Q R .The goodness of fit between observed and simulated discharge is measured by RMSEln (for the Q SB part) or RMSE (for other parts).

Figure 8 .
Figure 8. Simulation of daily streamflow by different methods from 2003 to 2007.(a) by the proposed stepwise method, (b) by the automatic calibration method, and (c) by the benchmark model.The performance of the simulations is measured in NSE, NSEln, and RMSE.

Figure 9 .
Figure 9. Evaluation of the stepwise calibration method.(a) discharge simulation in evaluation period 2008-2012 using the stepwise calibrated parameters in calibration period 2003-2007.(b) discharge simulation in evaluation period 2008-2012 by the benchmark model.(c) Cross-validation simulation of daily discharge in 2003-2007.x coordinate presents the simulated daily discharges by parameters calibrated in the period 2003-2007.y coordinate presents the simulated daily discharges by parameters calibrated in the period 2008-2012.(d) Crossvalidation simulation of daily discharge in 2008-2012.x coordinate presents the simulated daily discharges by parameters calibrated in the period 2008-2012.y coordinate presents the simulated daily discharges by parameters calibrated in the period 2003-2007.
We calibrated the model to 2008-2012 and evaluated it for 2003-2007.The new calibrated parameter values are K A = 0.9, K D = 0.003, D s = 2.2 mm • C −1 day −1 , D g = 7.4 mm • C −1 day −1 , W M = 10.2 cm, and B = 0.77, which are similar to the values calibrated in 2003-2007 and listed in Table 4.The NSE, NSEln, and RMSE values for the calibration period 2008-2012 and the evaluation period 2003-2007 are 0.757, 0.900, 10.892 m 3 s −1 and 0.883, 0.910, 8.589 m 3 s −1 , respectively, using this new calibrated parameter set.The simulations of the two periods by cross-validation are presented in Fig. 9cd, which show similar performance by two calibrated parameter sets and further demonstrates the robustness of the proposed stepwise calibration method.

Figure 10 .
Figure 10.Sensitivity analysis for the hydrograph partition.The first column is the hydrograph partition pattern using different lowest elevation bands of the glacier area (LEGs).The second column is the hydrograph partition pattern using different storm rain periods (SRPs).

Figure 11 .
Figure 11.Sensitivity analysis on the contributions of different runoff sources to total runoff.(a) is the contribution pattern under different lowest elevation bands of glacier area (LEGs), where the storm rain period (SRP) is fixed as May to September.(b) is the contribution pattern under different SRPs, where the LEG is fixed as 2950 m.The red line stands for the mean contribution for each runoff source, and the top/bottom end of each plot presents the highest/lowest contribution ratio.SB is groundwater baseflow, SM is snowmelt, GM is glacier melt, and R is rainwater direct runoff.

Table 1 .
Estimated monthly temperature lapse rate in the TRB.
• C day −1 /100 m) tive), and H and h are the elevation values at high and low positions, i.e., the mean elevation of two AWS and the elevation of THS, respectively.The values of T p in different months are obtained by minimizing the error function:

Table 2 .
Estimated weekly precipitation lapse rate in storm rain months.
Precipitation lapseMonth rate (mm week −1 /100 m) can infiltrate into soil and become additional contributions to groundwater.Groundwater forms baseflow that is separately calculated by two coefficients: K A and K D .K A and K D are outflow coefficients of groundwater storage.Their sum determines the flow rate of groundwater baseflow and their ratio (K D /K A ) dominate the proportion of free groundwater storage.Infiltration and storage should have effects on the calibration of the two parameters.The Xinanjiang module has been successfully applied to the Qiedeke, Kaidu, Manasi, and Kahai basins in the Tianshan Mountains by different authors The remainder of rainfall Hydrol.Earth Syst.Sci., 19, 1807-1826, 2015 www.hydrol-earth-syst-sci.net/19/1807/2015/

Table 3 .
Grouped parameters in the THREW model.Parameters subjected to calibration are highlighted in bold.

Table 4 .
Calibrated parameters by the stepwise and automatic methods.January, and February and are denoted by black dots.The rainwater surface runoff occurs in the SRP only (May to September, denoted by blue dots).The total number of days of the Q SB + Q SM part in the period from 2003 to 2007 is 365, and that of the Q SB +Q SM +Q GM part is 249, while the

Table 5 .
Evaluation merits for the stepwise and automatic calibration methods.

Table 7 .
R MS (%) for parameter sensitivity (R MS values indicating the most sensitive parameters are labeled in bold).SB + Q SM + Q GM ) SB + Q SM + Q GM + Q R ) 0.17