1 Introduction

More than one-sixth of the global population rely on glacier and snow melt for their water supply (Barnett et al. 2005). Changes in temperature and precipitation are expected to significantly affect the cryospheric processes and the hydrology of headwater catchments in the Himalayas (Cruz et al. 2007; Immerzeel et al. 2009). Accurate modeling of the hydrological response to climate change in these catchments is complicated due to the large climatic heterogeneity, the lack of data in mountain catchments, the resolution and accuracy of GCM outputs and uncertainty about the response of glacier and snow dynamics (Beniston 2003). Traditionally, glacier dynamics are not linked to other hydrological processes such as evapotranspiration, surface runoff and base flow in a single model (Sharp et al. 1998). Most hydrological impact studies on the contrary deploy simple degree day methods (Hock 2005) and assume hypothetical reduction in future glacier areas (Singh and Bengtsson 2004; Hock 2005; Rees and Collins 2006; Immerzeel et al. 2009). Although these studies provide valuable insights into the possible range of future options, they suffer from large uncertainty about the plausibility of the future evolution of snow and ice. Another issue that is often ignored in modeling melt from glaciers is that both accumulation and ablation of glaciers depend on the glacier area, which does not scale linearly with glacier volume. As larger glaciers have a smaller area to volume ratio they are less sensitive to climate change than smaller glaciers (Van de Wal and Wild 2001).

Recently there has been a strong debate about the melt rate of Himalayan glaciers. The claim that all glaciers in the Himalayas could disappear by 2035 in the AR4 of the IPCC (Cruz et al. 2007) has been the source of great controversy and has been admitted to stem from a wrong quotation of grey literature (Schiermeier 2010). Mass balance studies and regional reviews conclude otherwise, although there are very limited published studies and measurements available (Bolch et al. 2008; Zemp et al. 2009). Obviously, there is large uncertainty and variability in the retreat rates of Himalayan glaciers and to settle the debate in a rational manner there is a strong need for reference studies that reliably model the transient evolution of glaciers under climate change.

In this study we attempt to provide such a reference study by developing a combined cryospheric hydrological model for a glacierized Himalayan catchment in Nepal. The model explicitly simulates glacier movement in combination with major hydrological processes such as evapotranspiration, surface runoff, ablation and groundwater base flow at a high spatial resolution with a daily time step. An innovative two stage calibration procedure is used. First the historical evolution of glaciers is calibrated using the recent location of glacier tongues and secondly the hydrological processes are optimized using discharge observations. The calibrated model is then forced with an ensemble of transient statically downscaled GCM outputs to evaluate the climate change impact on the future evolution of the glaciers and catchment hydrology.

2 Study area

The Langtang river catchment (360 km2) is located approximately 100 km north of Kathmandu (Fig. 1). The elevation ranges from 3800 m.a.s.l. up to the peak of Langtang Lirung at 7234 m.a.s.l with an average altitude of 5169 m.a.s.l.. In total 46% (166 km2) of the catchment is glacierized. The glacier tongues below 5200 m.a.s.l. are generally debris covered (32 km2). The main valley is dissected by the Langtang Khola River and it is typically U-shaped. Several sets of moraines occupy the valley bottom, which have been attributed to the Little Ice Age, Neoglacial and Late glacial (Heuberger et al. 1984). Boulders and scree cover the steep slopes and high plateaus, while the occurrence of forest and grassland in the lower altitudes with less steep slopes along the river is limited to 1.5% of the catchment area. The riverbed consists of boulder and gravel. The area is part of the Main Central Thrust Zone, and its geology consists of granite, gneiss and schist. The characteristic feature of the climate of Nepal is the monsoon circulation with predominant easterly winds in the summer and westerly winds from October to May. In the Langtang Khola catchment, 77% of annual precipitation of 814 mm y−1 falls during monsoon season (June to September; mean over period from 1957 to 2002, (Uppala et al. 2005)). From June to August total precipitation amounts are large and precipitation occurs almost every day. However, the daily amount generally does not exceed 20 mm d−1. In the later part of the monsoon season (September to October) the maximum daily amount is considerably higher, while the number of rainy days decreases gradually compared to previous months. During the dry season (November to May), precipitation (mainly snow) occurs on only a few days, as it is produced by the occasional passage of westerly troughs (B. Wang 2006). In general, the precipitation amounts increase with altitude during both the monsoon and the dry season (Seko 1987). Mean daily air temperature from 1957 to 2002 is 0.5°C from October to June, and 8.4°C during the monsoon season at the meteorological station (3920 ma.s.l.) (Fig. 1) and temperature is strongly correlated with elevation (Shiraiwa et al. 1992; Sakai et al. 2004). The discharge regime can be classified as glacial with maximum discharges in July and August and minimum discharges during the winter season. Winter discharge is characterized by a rather constant base flow with negligible inflows of rainwater or melts water as the air temperature is generally below the melting point.

Fig. 1
figure 1

Location of Nepal and the Langtang khola catchment

Full size image

3 Methods

The entire model is developed using the PCRaster environment for numerical modeling in environment science (Karssenberg et al. 2001). The model is setup at a spatial resolution of 90 m (338 × 325 cells). For each cell the following model concepts are simulated at a daily time step.

3.1 Glaciers

We assume that the principal process for glacial movement is basal sliding. There is some evidence that in the north-western part of central Asia in the west Kunlun Shan mountains (Liu et al. 2009), some glaciers are cold based and basal sliding does not occur, but glacial creep is the main process for glacier movement. Hewitt (2007) discusses that glaciers in the Karakoram are both cold and temperate across their length and even height and that shifts in their thermal regime may change their flow rates substantially. He (2003) shows that temperate glaciers are mainly found on the southern-eastern part of the Himalaya’s that is under strong influence of the monsoon. Kumar and Dobhal (1997) show that basal sliding is the main process of glacier movement for the Chhota Shigri glacier, which is also on the southern slopes of the Himalya’s and under strong monsoon influence . Mae et al. (1975) and Ageta and Higuchi (1984) also show that glaciers in the eastern part of Nepal are most likely warm-based. Therefore we assume that on the south-eastern part of the Himalayas, glaciers are mostly temperate given the much wetter climate, glacier morphological characteristics and the synchronous accumulation and ablation during the monsoon season. Basal sliding will be the dominant process of glacier movement in that case and this assumption is evaluated by checking the plausibility of the optimized parameters related to the basal sliding. Glacier bottom motion in that case is modeled as slow, viscous flow using Weertman’s sliding law (Weertman 1957). It is assumed that the basal ice is at the melting point, such that a film of water is conceived to exist between the ice and the underlying bedrock. Two mechanisms are considered. The first is slow viscous deformation and the second is regelation. Regelation occurs because as the ice flows over a bedrock obstacle, the higher upstream pressure causes the ice to melt at the interface, because the melting temperature depends on pressure. The water which is thus formed squirts around the rock and correspondently refreezes on the downstream side. Weertman combines these two processes and derives the Weertman sliding law

$$ {\tau_b} \approx {\nu^2}R{u^{\frac{2}{{n + 1}}}} $$
(1)

here τ b is the basal shear stress (Pa), ν (−) is a measure of the roughness of the bedrock, R is a material roughness coefficient (Pa s1/3), u is the sliding speed (m s-1) and n (−) is the creep constant of Glenn’s flow law (~3 in most cases) Glen (1955). The driving force of glacial movement is gravity and τ b is defined as

$$ {\tau_b} = \rho gH\sin (\beta ) $$
(2)

Where ρ (kg m-3) is the ice density, g is the gravitational acceleration (m s-2), H (m) is the ice thickness and β (°) is the surface slope. By combining Eqs. 1 and 2 under the assumption that glaciers only move when the basal shear stress exceeds the equilibrium shear stress (τ 0 (N m-2)) the sliding speed can be derived

$$ {u^{\frac{2}{{n + 1}}}} = \frac{{\rho gH\sin (\beta ) - {\tau_0}}}{{{\nu^2}R}} $$
(3)

This equation is used to model glacial movement as function of slope, ice thickness and bedrock properties for each daily time step for each cell. Snow accumulates in the upstream parts of the catchment and each time step the sliding speed is calculated for each cell. Based on the sliding speed the ice is transported down the digital elevation model. For numerical stability and given the viscous properties of ice, the ice is not transported to a single downstream cell but is distributed proportionally to all downstream cells according to the slope. As the snow and ice is progressively moving downstream, the temperature increases and snow and ice ablation (Q a ) occurs using a degree day factor (ddf). A fraction (α) of the total ablation leaves the catchment as runoff while, the remainder of the ablation (1- α) is stored in the glacier/snow pack and is released as base flow. The initial ice depth is calculated based on the known extent of the glaciers derived from remote sensing and the depth is estimated by:

$$ H = \frac{{{\tau_0}}}{{\rho g\sin (\beta )}} $$
(4)

The ddf is not uniform in space but depends on both exposure and in case of glaciers on debris cover. Debris covered glaciers are insulated and melt at a much lower rate than clean ice glaciers and on south faces the incident shortwave radiation is higher and thus the ddf is larger than on north faces (Konz et al. 2007). The exposure dependence of the ddf on the aspect is calculated by

$$ dd{f_c} = ddf \cdot ({1} - {R_{exp}} \cdot {\text{cos(}}aspect{)}) $$
(5)

Where ddf c is the corrected ddf and R exp is a factor quantifying the aspect dependence of the ddf. For the debris covered glaciers a multiplicative reduction factor has been used and for the Lirung glacier a separate multiplicative reduction factor for the ddf has been used. The relative low position of the glacial snout of the Lirung glacier can only be explained by a relative low melt rate of the glacier tongue due to a thick insulating debris cover. The model parameters related to glacier modeling are shown in Table 1.

Table 1 Model parameters
Full size table

3.2 Hydrological processes

The model is forced by daily precipitation and temperature data at Kyangjing from 2000 to 2007. Temperature is spatially differentiated using a vertical lapse (λ t ) rate, which is a calibration parameter. For precipitation both a positive vertical (λ p,v ) and a negative horizontal gradient (λ p,h ) from west to east is applied similar to (Konz et al. 2007). Precipitation is partitioned in either snow (S) or rain (P) using these daily fields, the lapsed temperature fields and a threshold temperature.

Reference evapotranspiration (ET 0 ) is calculated based on minimum, maximum and average temperature data at Kyangjing from 2000 to 2007 according to the Hargreaves equation (Hargreaves et al. 1985). ET 0 is spatially differentiated per cell by applying a temperature dependent correction factor (cET 0 ). This factor is derived by using the temperature dependence of the Hargreaves equation. Actual ET (ET a ) is derived by first deriving potential evapotranspiration by multiplying ET 0 by a crop reduction factor (K e ). By limiting potential evapotranspiration with the actual soil water content actual evapotranspiration is calculated.

On snow and ice free cells surface runoff (Q s ) is calculated according to the curve number method (SCS USDA 1972). The curve number (CN) parameter is a calibration parameter. The sum of P-ET a -Q s is then added to the soil water storage and if the maximum soil water storage (θ m ), which is calibrated, is exceeded recharge to the groundwater occurs. Groundwater base flow (Q b ) is modeled similar to the SWAT model (Neitsch et al. 2005) that is based on the work of (Smedema and Rycroft 1983) to quantify the non-steady response of groundwater flow to recharge. Finally, for each cell the total runoff (Q) is calculated as the sum of Q a   + Q s   + Q b and routed to the catchment outlet following a recession equation similar to the Snow Melt Runoff Model (Martinec 1975):

$$ {Q_{out,t}} = (1 - k) \cdot {Q_t} + k \cdot {Q_{out,t - 1}} $$
(6)

Where k is a recession coefficient that is calibrated and Q out,t (m3 s-1) is the river discharge at the catchment outlet on day t. The model parameters related to the hydrological processes are shown in Table 1.

3.3 Calibration

The model is calibrated using the Parameter ESTimation (PEST) software. PEST is able to run a model as many times as it needs to while adjusting its parameters until the discrepancies between selected model outputs and a set of observations is reduced to a minimum. The PEST algorithm method is based on non-linear parameter estimation theory (Doherty 2005).

The calibration procedure follows a two step approach. First the parameters related to the glacier modeling are manually calibrated as glacier evolution is a much slower process than the hydrological processes such as evapotranspiration, surface runoff, ablation and base flow. Four parameters that influence glacier evolution are calibrated (τ 0 , R, λ t , ddf). A bias corrected time series of precipitation and temperature from 1957–2002, based on the ERA40 dataset (Uppala et al. 2005), is used to force the model with the aim to reproduce the location of glaciers and permanent snow in the year 2000 based on remote sensing (Konz et al. 2007). The bias correction is done similar to the corrections applied for the GCMs Eq. 7/Eq. 9. In the second step, PEST is used and four hydrological parameters are calibrated (θ m , CN, α and k) by simulating the period 2000–2006 for which observations of precipitation and temperature at Kyangjing are available. The parameters are calibrated against observed daily discharges from 2000–2006 at the outlet of the catchment.

3.4 Downscaling GCM output

Future climate GCM data for the A1B SRES scenario (Meehl et al. 2007) were used. Monthly temperature and precipitation from 2001 to 2099 were downloaded from the IPCC data distribution centre (http://www.ipcc-data.org) of five different GCMs: (i) the CGCM3.1(T47) model of the Canadian Centre for Climate Modelling and Analysis (CCMA-CGCM3), the CM2.0 model of the Geophysical Fluid Dynamics Laboratory in the USA (GFDL-CM2), the ECHAM5 model of the Max Planck Institute for Meteorology in Germany (MPIM-ECHAM5), the high resolution version of the MIROC3.2 model of the National Institute for Environmental Studies in Japan (NIES-MIROC3) and the HADGEM1 model of the Hadley Centre for Climate Prediction and Research in the United Kingdom (UKMO-HADGEM1). More information and details on these models are given in (Randall et al. 2007). Given that climate is highly variable in mountain areas and the resolution of the GCMs is coarse, the monthly GCM data were statistically downscaled using the 2001–2007 precipitation and temperature reference dataset of Kyangjing. The GCM data were corrected such that the monthly averages and variability during 2001–2007 matches the observations at Kyangjing similar to the procedure described in (Aerts and Droogers 2004). First the monthly averages and standard deviations are calculated for the reference dataset and for each GCM then a correction factor is determined for the average Eq. 7 and the variation Eq. 8 in each climate parameter.

$$ {a_{adj,M}} = \frac{{\overline {{a_{obs,M}}} }}{{\overline {{a_{gcm,M}}} }} $$
(7)

and

$$ {\sigma_{adj,M}} = \frac{{{\sigma_{obs,M}}}}{{{\sigma_{gcm,M}}}} $$
(8)

Where aadj is the correction factor for averages, \( \overline {{a_{obs}}} \) the average observed climate parameter (either average temperature or precipitation), \( \overline {{a_{gcm}}} \) the average simulated climate parameter, σadj is the correction factor for standard deviations, σobs the standard deviation of the observed climate parameter, σgcm the standard deviation of the simulated climate parameter and the subscript M indicates the month. From these two adjustment factors the monthly GCM values from 2001–2099 are calculated according to

$$ a_{gcm,M}^\prime = \left( {{a_{gcm,M}} - \overline {{a_{gcm,M}}} } \right) \cdot {\sigma_{adj,M}} + (\overline {{a_{gcm,M}}} \cdot {a_{adj,M}}) $$
(9)

Where a’gcm, M is the corrected climate parameter.

The monthly downscaled GCM time series from 2001 to 2099 were subsequently disaggregated into daily values by using monthly distributions of the daily 45 year ERA40 dataset from 1957–2002. The total monthly sums are equal to the downscaled GCM time series but the partitioning into daily values is based on the statistical distribution of the ERA40 dataset. This procedure resulted in a daily statistically downscaled precipitation and temperature time series from 2001 to 2099 for the five different GCMs for the Kyangjing station, which were then used to force the calibrated model. It can be expected that in the transient GCM runs both the effects as climate change (trends) as well as multi-year variability are present, such that they will be accounted for in our projections. Finally using ERA40 to simulate daily variation will assure that within month variation is properly represented. Note that further distribution of precipitation and temperature of the downscaled Kyangjing time series within the catchment is achieved by applying vertical temperature lapse rate and horizontal and vertical precipitation lapse rates (Table 1) similar to the reference run from 2000–2007.

4 Results and discussion

Figure 2 compares remotely sensed glacier tongue positions (panel A) with simulation results (panel B). The model is able to simulate glacier evolution with high accuracy. After a 45 year simulation the glacier tongues are nearly in the same position as is observed by remote sensing. The green dot in the left figure also shows the snout position of the Lirung glacier as derived from GPS ground observations in August 2009. The difference between the remotely sensed position and the field observations can be explained by the 9 year time difference between the satellite image (year 2000) and 2009. However, the simulation of the glacier snout position corresponds reasonably well with the ground data and the slight difference could be due to the 7 years difference between the end of the 45 year simulation and the field work in 2009. The calibrated parameters related to glacier movement are all within realistic and plausible ranges and this confirms our assumption that glacier movement occurs mainly by basal sliding (Weertman 1957; Zhang et al. 2006; Konz et al. 2007). We have assumed that basal sliding is the major process of glacier movement and that the resistance to glacier slip is dominated by drag over bed irregularities. Experiments have revealed however that in case of debris covered glaciers additional debris-bed resistance may further slowdown the flow velocity of debris covered glaciers (Iverson et al. 2003). We use an overall resistance with the bedrock, but further research is recommended as to how debris-bed resistance could be parameterized in data-scarce environments.

Fig. 2
figure 2

Glacier extent in the year 2000 based on visual interpolation of Landsat image a and simulated ice depth in 2002 based on a 1957–2002 simulation with the calibrated model b

Full size image

Figure 3 shows the average annual ice balance of the catchment from the 45 year simulation. Snow falls on the high elevations and is progressively transported downwards using the glacier movement equation Eq. 2. Most melt occurs on the glacier tongues and the lower melt rate on the Lirung glacier also shows in the ablation map. The net budget after 45 years of simulation is a slight thickening on the glacier tongues and a very small retreat of the snouts of most of the glacier tongues compared with the initial ice conditions that are based on the equilibrium shear stress. From these results we conclude that the cryospheric processes have been represented sufficiently realistic to allow for a climate change impact analysis.

Fig. 3
figure 3

Average annual snow fall A, ablation B, net mass balance by sliding C and net ice balance D based on a 1956–2002 simulation

Full size image

Figure 4 shows the daily simulated and observed discharges at Kyangjing station. Based on these time series a number of model performance indicators were determined; the Nash-Sutcliffe model efficiency coefficient equals 0.76, the Pearson correlation coefficient is 0.87 and the bias is 3%. In addition the calibrated parameters (Table 1) are all in a plausible range. Hence, the model is able to simulate the discharge at the outlet of the catchment realistically and with relatively high accuracy.

Fig. 4
figure 4

Simulated and observed river discharge at the Kyangjing discharge station from 2001 to 2006

Full size image

Figure 5 shows that between 2000 and 2100 an annual increase in temperature of 0.06°C y-1 is projected. There is considerably variability between the different GCMs but the average trend is distinctly positive and in line with global projections for the northern hemisphere. The precipitation prediction exhibits much more capricious behavior and the variability between the 5 GCMs is large. However, on average a positive precipitation trend is predicted of 1.9 mm y-1 consistent with the acceleration of the hydrological cycle due to the increased atmospheric temperature (Fig. 5). The positive temperature and precipitation trend will influence the hydrological cycle in several ways. Firstly, the higher temperatures will result in a higher reference evapotranspiration and increased melt of snow and ice. Secondly in addition to an increase in total precipitation the fraction of that falls in the form of rain instead of snow will increase. The net effect on the glacier area, total discharge and its composition, and temporal shifts in the hydrograph shifts are further quantified by our hydrological model analysis.

Fig. 5
figure 5

Downscaled annual total precipitation and average temperature for Kyangjing from 2001 to 2006. Bold line indicates the multi-model average (MMA) of 5 GCMs. Error bars indicate the 1σ of the 5 GCMs

Full size image

Figures 6 and 7 show the future decrease in glacier area and glacier volume. A rapid and steady decline of the glaciers can be observed and from 2075 onwards permanent ice can only be found at the highest elevations in the catchment. By 2035 the glacier area has decreased by 32%, by 50% in 2055 and by 75% in 2088. The total ice volume shows a similar behavior but the decrease is slightly faster than the area. This catchment is representative of high altitude glacierized catchments in the central and eastern Himalayas and the findings confirm recent statements that the Himalayan glaciers are likely not to have disappeared by 2035 (Schiermeier 2010). However it should be stressed that it is not possible to give a single date for the entire Himalayas. Across the Himalayas there are west to east and north to south transitions towards a warmer and wetter climate and this shows in the pattern and complexity of glacier changes being observed. Overall the Himalayan glaciers are melting, but regional anomalies exist and the rate of change is highly variable because glacier dynamics depend on the precipitation and temperature regime, topography, volume-area ratios and debris cover.

Fig. 6
figure 6

Simulated glacier extent in 2025 A, 2050 B, 2075 C and 2100 D based on forcing with MMA time series from 2001–2100

Full size image
Fig. 7
figure 7

Anomalies in glacier area (top) and glacier ice volume (bottom) (% deviation from 2000). Error bars indicate 1σ of the 5 GCMs. Glaciers are defined as all pixels with an ice thickness of more than 5 m

Full size image

The eventual impact of climate change on the hydrology does not merely depend on the melting of the glaciers but other factors are equally important. Our results reveal surprisingly that the total discharge is increasing between 2000 and 2100 (Fig. 8). The mean ensemble discharge significantly increases annually with 0.05 m3 s-1 (~4 mm) and the increase slows down towards the end of the century. By 2050 the discharge has increased by 32%. Part of this increase can be attributed to an annual increase in net precipitation (0.03 m3 s-1). The remainder is explained by the gradual change from a melt water river towards a rain river that causes the runoff coefficient to increase from 50% in 2000 to 65% in 2100.

Fig. 8
figure 8

Total MMA simulated discharge partitioned into rain runoff, glacier runoff, snow runoff and baseflow

Full size image

Figures 8 and 9 also show the partitioning of stream flow in rain runoff, base flow, snow and glacier runoff. Over a period of 100 years an increase in rain runoff and baseflow is observed, the snow runoff remains constant and the glacier runoff decreases. It is interesting to note that the glacier runoff remains relatively constant until 2040 as the reduction in glacier area is compensated by the increased melt per unit area. However, after 2040 the decrease in glacier area becomes dominant leading to runoff decline. The increase in rain runoff and base flow is related to an overall increase in precipitation, a shift from rain to snow and a decrease in glacier area. It should be noted that an increase in the relative contribution of rain runoff will also results in a more direct response of river runoff to extreme dry or wet years. Figure 9 also shows that there is a large variation between the different GCMs in the partitioning of streamflow, in particular in base flow and glacier runoff. These differences can be explained by the large variation in temperature and precipitation projections between the GCMs.

Fig. 9
figure 9

Streamflow partitioning in 2005, 2025, 2050 and 2075 of the MMA. Error bars indicate the 1σ of the 5 GCMs

Full size image

Figure 10 shows the monthly distribution of stream flow in 2000, 2025, 2050 and 2075 revealing a similar trend in partitioning of runoff components as in Figs. 8 and 9. A steady base flow is simulated in winter and the stream flow is highest in July during the peak of the melt and monsoon season. Shifts from glacier to rain runoff and base flow are most pronounced during the summer months. It is however interesting to note that the monthly distribution of total runoff remains relatively unchanged because the melt season coincides with the monsoon season so that no shift in runoff regime is point is lacking after observed.

Fig. 10
figure 10

Monthly partitioning in 2005, 2025, 2050 and 2075 of the MMA simulated discharges

Full size image

5 Conclusions

In this study a combined cryospheric hydrological model is developed that explicitly simulates glacier movement in combination with major hydrological processes. The model is used to model the transient evolution of glaciers and hydrology under climate change. From an application of the model in the Langtang catchment in Nepal the following conclusions are drawn:

  • The combination of hydrological processes and glacier movement in a high resolution raster based model enables the accurate simulation of both glacier evolution and river flow in a high altitude glacierized catchment in the Himalayas.

  • Climate change analysis using downscaled data from 5 different GCMS shows that temperatures are projected to increase by 0.06°C y-1 and precipitation by 1.9 mm y-1. The analysis also reveals a large variability among the different GCMs in particular for precipitation.

  • In the catchment the glaciers are retreating steadily under climate change and it is estimated that in 2035 the glacier area has been reduced by 32%. This catchment is representative for the southern slopes of central and eastern Himalayas where glacier systems are dynamic, moderate in size and often characterized by debris covered tongues.

  • The positive temperature and precipitation trends will increase evapotranspiration and snow and ice melt while more precipitation will fall as rain instead of snow. The net result is an increase in stream flow by 4 mm y-1 that can be attributed to the increase in precipitation and the change from melt river to rain river. The partitioning of stream flow is indeed showing strong changes. Rain runoff and base flow are increasing, snow runoff remains more or less constant and glacier runoff is eventually decreasing.

This study shows an extensive analysis of a glacierized catchment in the central Himalayas but the results are not representative for the entire Himalayas. To arrive at a comprehensive assessment on how climate change is affecting the hydrology of the Himalayas it is recommendable to perform this analysis in reference catchments covering the east–west and north–south gradient in climatology and glacier and hydrological dynamics.