1. Introduction
Snow covered and glacierised catchments play an important role for controlling seasonal discharge in river systems (van Tiel and others, Reference van Tiel2023). Winter snowfall is stored in the catchment area and generally released the following spring and early summer when temperatures are high and dry conditions prevail. In glacierised catchments, some fraction of the winter snowfall can be transformed to ice and stored for longer periods before it is released to the water cycle again. Both constituents regulate streamflow on different time scales and have the potential of filling precipitation gaps during hot and dry periods (Fountain and Tangborn, Reference Fountain and Tangborn1985; Singh and Singh, Reference Singh and Singh2001; Kaser and others, Reference Kaser, Großhauser and Marzeion2010).
Non-glacierised catchments are governed by precipitation events, while runoff from glaciers depends on temperature and radiation variations (Chen and Ohmura, Reference Chen and Ohmura1990). Snow-dominated catchments show a maximum in streamflow during the spring snowmelt season, while glaciers provide an additional source of water, especially in summer after the snow cover in the nonglacierised part has melted (Moore and others, Reference Moore2009). In addition, strong diurnal streamflow variations occur, especially during the high melt season, due to the daily temperature and radiation variability (Lane and Nienow, Reference Lane and Nienow2019). Increasingly negative glacier mass balances in many high mountain regions, including the Alps, will lead to a higher melt contribution from glaciers, which will reach a maximum and then decrease due to the shrinking glacier area (Immerzeel and others, Reference Immerzeel, Pellicciotti and Bierkens2013; Huss and Hock, Reference Huss and Hock2018; Pritchard, Reference Pritchard2019). At the same time, the gradual removal of the firn layer results in an acceleration of the water flow through the glacier system (de Woul and others, Reference de Woul2006), which affects the seasonality of the discharge (high spring and lower autumn discharge).
Due to the important influence of snow and ice ablation on the spring and summer discharge in particular, runoff from high alpine catchments has been investigated in many studies (e.g. Hanzer and others, Reference Hanzer, Helfricht, Marke and Strasser2016; Wortmann and others, Reference Wortmann, Bolch, Su and Krysanova2019; Eidhammer and others, Reference Eidhammer2021). Several model approaches have been developed to simulate the discharge from such catchments (e.g. Hanzer and others, Reference Hanzer, Helfricht, Marke and Strasser2016; Jobst and others, Reference Jobst, KingstonD, Cullen and Schmid2018; Biemans and others, Reference Biemans2019; Shahgedanova and others, Reference Shahgedanova2020). However, the model performance is often affected by high uncertainties in input datasets and parameters. Especially precipitation input is affected by considerable errors, and there exists the risk of incorrect parameterisations in accumulation and melt processes due to error compensation in the model (van Tiel and others, Reference van Tiel, Stahl, Freudiger and Seibert2020a). Unfortunately, observations of melt rates and streamflow, which are essential for constraining models, often do not exist in high alpine catchments. The lack of meteorological observations in these regions frequently requires driving snow-hydrological models with information from remote stations, mainly at low altitudes (van Tiel and others, Reference van Tiel, Stahl, Freudiger and Seibert2020a), and relying on estimated lapse rates for e.g. air temperature and precipitation. Therefore, the performance of models does not necessarily improve with the degree of sophistication, especially in data scarce regions (Muñoz and others, Reference Muñoz, Huggel, Drenkhan, Vis and Viviroli2021). In addition, there exist general uncertainties regarding the quality of forcing data (e.g. precipitation), as well as the evaluation data (e.g. measured discharge) (Shannon and others, Reference Shannon2023). On the other hand, well-performing physical models can be used for a process-based analysis of cause-effect relationships and thus can provide improved results for e.g. trend attribution, compared to knowledge derived from measurements alone (Duethmann and others, Reference Duethmann2015).
In the light of these results, it can be expected that high quality in situ observations have the potential to considerably improve the model performance, and thus allow a better insight into the significance of glaciers to the runoff characteristics in glacierised catchments. Especially the compensating effect of ice melt during dry and hot periods, as well as the timing and magnitude of the snowmelt contribution to streamflow depend also on local characteristics, which might not be captured well by generalised model approaches (van Tiel and others, Reference van Tiel, Kohn, Van Loon and Stahl2020b, Reference van Tiel, Van Loon A, Seibert and Stahl2021).
The quantification of the individual runoff components is highly relevant for characterising the runoff dynamics and seasonality of high alpine basins (He and others, Reference He, Duethmann and Tian2021). The component separation can be done by stable isotope investigations (e.g. Crespo and others, Reference Crespo, Fernandoy, Cara, Klarian and Lavergne2020), using constrained hydrological models (e.g. Fatima and others, Reference Fatima, Hassan, Hasson, Ahmad and Ali2020), or quantifying observations of the water balance terms (e.g. Kumar and others, Reference Kumar, Singh and Singh2007). We take advantage of the rather unique setting at Vernagtferner in the Ötztal Alps to investigate the contribution of runoff components to streamflow in detail. The glacierised catchment provides a high density of in situ observations of meteorological and hydrological parameters, but also information about the evolution of the snow and ice cover from the long-term monitoring programme. We established additional automatic weather stations (AWS) on the glacier, in order to obtain direct observations of important mass and energy balance parameters on the ice, but also from the accumulation region.
The study aimed to investigate runoff at high temporal resolution and the dependence of streamflow variability on the different melt contributions, by using high-quality local observations as model input and for validation. One focus was to test the robustness of simulating hourly to daily discharge amounts in comparison to measured values, while predicting the contribution of the runoff components. For this purpose, we performed a detailed model study on the evolution of the snow pack across the glacier for the glacier mass balance years 2019/20 through 2021/22. Due to the very different nature of these three years (moderately negative glacier mass balance, close to long term mean balance and extremely negative mass balance), we could reconstruct the range of variability of the discharge from this glacierised basin, not only by magnitude, but also in relation to the water sources. We implemented the SNOWPACK and Alpine3D software packages (Lehning and others, Reference Lehning, Bartelt, Brown and Fierz2002b; Lehning and others, Reference Lehning2006) to simulate the distributed snow pack evolution and firn/ice ablation. In addition, we use a distributed and site-calibrated version of the HBV-ETH model (Mayr and others, Reference Mayr, Hagg, Mayer and Braun2013) to determine the basin discharge at the gauging station Vernagtbach. Snow thickness measurements were used to calibrate the SNOWPACK model, while meltwater production was validated against the continuous discharge measurements at the basin scale.
2. Study area
Vernagtferner in the Ötztal Alps has a long history of glaciological observations. It was visited in the 17th century, due to its episodic surges (Walcher, Reference Walcher1773). Seasonal mass balance has been measured since 1964 (Reinwarth and Escher-Vetter, Reference Reinwarth and Escher-Vetter1999), while geodetic investigations reach back to 1889, when the first detailed map of the glacier was derived from terrestrial photogrammetry (Finsterwalder, Reference Finsterwalder1897). Monitoring of discharge from the Vernagtferner basin was started in 1974 by the construction of a river gauge, together with meteorological observations (Oerter and others, Reference Oerter, Baker, Moser and Reinwarth1981) at Pegelstation Vernagtbach (PS in Fig. 1). The drainage basin has an area of 11.7 km2, of which the glacier covered about 6.9 km2 in 2022. The altitude in the drainage basin ranges from 2640 m a.s.l. to 3610 m a.s.l.. The region is an ideal test site for investigating the glacial snow and firn cover dynamics and the associated meltwater production in detail, due to the wealth of direct observations of the glacier evolution and basin discharge.
To better constrain the accumulation and melt conditions, we performed additional observations on the glacier. Snow thickness and ice ablation is monitored continuously close to the eastern terminus since 2003 (Ablatometer in Fig. 1). An automatic weather station collects meteorological data, snow thickness and ice ablation in the central eastern part of the glacier (AWS-A in Fig. 1) at 3070 m a.s.l.. Another station was installed at the Hochvernagtplateau, the highest plateau of the glacier at 3450 m a.s.l. in 2018 (HVP in Fig. 1), where meteorological data, including snow thickness are continuously observed. In addition, an upward looking radar (upGPR) at HVP provides information about the internal structure of the snow and firn pack, which allows also the retrieval of the snow water equivalent (SWE) variability over time (Heilig and others, Reference Heilig, Eisen and Schneebeli2010; Schmid and others, Reference Schmid2015). The glacier has had negative mass balances since the early 1990s and accumulation basins are currently restricted to less than 20% of the glacier area. Discharge from the glacier varies strongly during the year, but also on a daily scale, due to the changing melt and accumulation conditions. While total discharge in winter is the order of 0.01–0.1 m3 s−1, peak discharge can reach more than 10 m3 s−1 during strong melt periods. Rainfall events (afternoon thunderstorms during hot spells) can lead to an even higher total discharge.
3. Data and methods
3.1 Data
3.1.1 Meteorological data
We use data from the different AWS for this study. The stations record meteorological data and snow cover properties every five minutes (Table 1). Data recording started in different years but since 2018 the data has been collected in parallel. The availability of meteorological observations at different altitudes allows us to adjust altitude dependent lapse rates used in the model to the local conditions. Furthermore, we use selected data sets from the stations for validating the model results.
TA, air temperature; RH, relative humidity; HS, snow thickness; VW, wind velocity; DW, wind direction; ISWR, incoming short wave radiation; ILWR, incoming longwave radiation; RSWR, reflected short wave radiation; RLWR, reflected longwave radiation; PSUM, precipitation.
The abbreviations of the meteorological parameters are explained below the table.
To fill for data gaps, which occurred due to sensor problems or power failure (periods are given in Table S1), data from operating stations are used for the interpolation of parameters such as air temperature, relative humidity, wind direction and wind velocity, as well as snow thickness. The gaps are filled by using the lapse rates or characteristic differences calculated from periods when the stations were operating in parallel. Data gaps only occurred at HVP and AWS-A for some periods.
3.1.2 Glaciological data
Seasonal mass balance data exist since 1964 for Vernagtferner. Annual balance is determined at the end of September, while winter balance is measured at the end of April (fixed date system) (Moser, Reference Moser1980; Reinwarth and Escher-Vetter, Reference Reinwarth and Escher-Vetter1999). Stake measurements of ablation, snow density from shallow snow pits in the accumulation region and the geometry of the glacier tongues and the equilibrium line from terrestrial GNSS measurements and remote sensing data are used for determining the annual mass balance. The spatial distribution is interpolated from these point and line data. The winter balance is based on distributed manual snow thickness probing and several snow pits across the basin for determining the mean snow density. The summer balance results from the difference of annual and winter balance. Density profiles are measured in snow pits at the AWS locations in order to relate the recorded snow thickness to SWE. The overall mean uncertainties of the glaciological mass balance measurements are in the order of 300 mm w.e. a−1 (Zemp and others, Reference Zemp2013)
3.1.3 Runoff measurements
The discharge from the catchment area is recorded at Pegelstation Vernagtbach (PS) since 1974. For this purpose, the water level in the gauge is observed continuously, while the rating curve is manually determined every summer to account for changing water flow conditions (e.g. due to the relocation of large boulders upstream of the gauge). Winter runoff is usually below the sensitivity of the gauge and contributes only very minor amounts to the annual runoff. Reliable measurements therefore start with the onset of snowmelt during May and usually end in October or November, when a continuous snow cover is established and air temperatures are below the freezing point in general. Early manual measurements and isotope analysis in the 1970s show that winter runoff is very low (~0.09 m3 s−1 in November and ~0.01–0.03 m3 s−1in April) and the contribution from groundwater is in the order of 0.1 m3 s−1 in summer and as little as 0.01 m3 s−1 in winter (Behrens and others, Reference Behrens, Moser, Oerter, Rauert and Stichler1979; Oerter and others, Reference Oerter, Baker, Moser and Reinwarth1981). This relates to a relative contribution of flow through the groundwater storage to total discharge of 13% in summer and up to 100% in winter. The mean estimated error of the discharge observations is about 10%.
3.1.4 Geometry of the drainage basin
We use a digital elevation model (DEM) based on aerial photogrammetry recorded in summer 2018 by 3D RealityMaps, with a resampled ground resolution of 1 m as input for the modelling study. Horizontal location errors are less than the original pixel size (0.2 m) and vertical errors are in the order of a few decimetres (Geissler and others, Reference Geissler, Mayer, Jubanski, Münzer and Siegert2021). The accompanying orthophoto was used to derive the glacier extent in 2018 and to determine the surface classification for the model. The resampled land surface class map and DEM (pixel size: 20 m) form the basis for the spatially distributed model simulations.
3.2 Methods
3.2.1 SNOWPACK model
SNOWPACK is a physical, 1-D snow and land surface model, originally developed by the Swiss Federal Institute for Snow and Avalanche Research for avalanche warning purposes. The model considers the exchange of mass and energy between the snow and the atmosphere, as well as soil and vegetation. SNOWPACK treats snow as a three-component porous medium (ice – water – air), simulates its microstructural properties in relation to its thermal and mechanical properties (Bartelt and Lehning, Reference Bartelt and Lehning2002; Steger and others, Reference Steger2017) and considers phase-change processes as well as transport of heat and water in the snow cover (Lehning and others, Reference Lehning1999; Fierz and Lehning, Reference Fierz and Lehning2001; Lehning and others, Reference Lehning, Bartelt, Brown, Fierz and Satyawali2002a, Reference Lehning, Bartelt, Brown and Fierz2002b). The original water bucket model for vertical water transport in snow and soil has been further developed to extend its applicability, introducing different water transport schemes and explicitly solving Richard's equation (Wever and others, Reference Wever, Fierz, Mitterer, Hirashima and Lehning2014, Reference Wever2015). Also, the simulation of ice layer formation has been introduced (Wever and others, Reference Wever, Comola, Bavay and Lehning2017), as well as including soil-specific characteristics. The mass and energy transport and phase-change processes within the soil are treated in the same way as in the snow layers for investigation of e.g. permafrost sensitivity (Luetschg and others, Reference Luetschg, Lehning and Haeberli2008). Because the model is one-dimensional, the input consists of meteorological data for a specific location (e.g. an automatic weather station), as well as information about the initial vertical snow and soil profiles. The most important input parameters are air temperature, relative humidity, wind speed and direction, short- and long-wave incoming radiation, snow thickness and precipitation. We use the bucket model for vertical water transport and the routines for phase-change processes like melting, refreezing and ice layer formation. Detailed analysis has shown that SNOWPACK simulates the temporal evolution of the snow thickness and SWE realistically (Michel and others, Reference Michel2022).
3.2.2 Alpine3D model
Alpine3D is a physically based, 3-component alpine surface process model, with SNOWPACK at its core, constructed to offer a high-resolution description of mountain surface processes, many of which show a high variability on a small spatial scale (Lehning and others, Reference Lehning2006). The physical processes of mass and energy transport, as well as resulting phase changes are identical to SNOWPACK. Atmospheric processes, such as radiation, shading and reflection (Helbig and Löwe, Reference Helbig and Löwe2014) and preferential deposition and snow transport (Lehning and others, Reference Lehning, Löwe, Ryser and Raderschall2008) are simulated using 3-dimensional model modules, whereas below surface processes in snow and soil are described by a 1-dimensional module without lateral exchange. Depending on the goals of the simulation and data availability, Alpine3D can be driven with data from a single AWS within or close to the model domain, multiple AWS data or gridded data from a regional climate model. Additionally, a DEM, a land-use classification and, as for SNOWPACK, an initial vertical snow and soil point profile must be provided for the initialisation of the model domain. Ice area and thickness need to be provided for glacierised regions as well. For our model runs, we used an ice layer of constant thickness and thicker than the maximum ablation during the model period. Ice dynamics are not considered, because the typical maximum velocities of less than 2 m yr−1 have no influence on the results for the chosen simulation period. The standard parameters for the individual physical processes included in Alpine3D were used as described in the model documentation. In addition, water transport on, in and below the glacier in this small catchment is fast (maximum routing length is less than 2 km) and the mean residence times of a few hours for the glacierised part and in the order of days in the rest (Oerter and others, Reference Oerter, Baker, Moser and Reinwarth1981) are treated in the routing of the HBV model. Subsequently, the data are statistically interpolated over the DEM of the domain, by the module METEOIO (Bavay and others, Reference Bavay and Egger2014). Details about the initial parameter values can be found in the supplementary material (S2). Grid resolution can have some impact on the determined SWE (Schlögl and others, Reference Schlögl, Marty, Bavay and Lehning2016), which indicates the need for a careful simulation design. However, Alpine3D has been successfully applied by a number of investigations (e.g. Michlmayr and others, Reference Michlmayr2008; Mott and others, Reference Mott2008; Haberkorn and others, Reference Haberkorn2017), also demonstrating its suitability for assessing the sensitivity of catchments to flood risks from snowmelt and rainfall events (Wever and others, Reference Wever, Comola, Bavay and Lehning2017).
3.2.3 Distributed HBV model
We calculate discharge by using the modelled meltwater production on the 20 m model grid as input to the routing routine of a distributed version of the HBV-ETH model (Mayr and others, Reference Mayr, Hagg, Mayer and Braun2013). The original HBV-ETH model (Braun and Aellen, Reference Braun and Aellen1990) is a conceptual runoff model for glacierised basins, based on a lumped temperature-index model, which builds on altitude, aspect and land cover classification of the drainage basin. It contains routing through several storage units, representing the different sensitivity of sediment, soil and groundwater bodies. The model was refined by Mayr and others, Reference Mayr, Hagg, Mayer and Braun2013 to run on a distributed grid, in order to improve the runoff dynamics according to the real basin geometry. We use the runoff routine of this modified HBV-ETH model instead of the routing scheme of Alpine3D, because it is already calibrated to the conditions of the Vernagtferner basin (Mayr and others, Reference Mayr, Hagg, Mayer and Braun2013). This runoff routine contains an upper and lower storage compartment and a soil module. Outflow from the upper compartment occurs by surface runoff and interflow and depends on the maximum fill level and the storage level of the lower compartment. The lower compartment drains linearly like a groundwater storage, with a calibrated drainage rate. Outflow from the different compartments and surface runoff is lumped into a catchment runoff, which we use for the discharge analysis.
3.2.4 Snow, ice and runoff modelling
The snowpack evolution at Vernagtferner is modelled at the locations of AWS-A and HVP, using SNOWPACK with a time step of 15 min for the period 09/2019–09/2022. The measured snow thicknesses are used as input in addition to recorded (or reconstructed) meteorological parameters. The resulting precipitation output, the parameter PSUM (accumulated total precipitation), is then used as input in Alpine3D. The Alpine3D model was calibrated against field measurements (snow thickness, SWE and point mass balance) for the entire period of three years, in order to obtain an optimised parameter set for the best overall performance (i.e. reaching a maximum correlation coefficient).
In addition to using the METEOIO module within Alpine3D for the creation of filtered and homogenised time series of the meteorological variables (Fig. 2), the model is used to simulate the distributed snowpack over Vernagtferner glacier, at an hourly time step on a 20 × 20 m spatial resolution. The input consists of a DEM (Geissler and others, Reference Geissler, Mayer, Jubanski, Münzer and Siegert2021), a land-use model derived from the DEM and orthophotos and meteorological input from the three AWSs (AWS-A, HVP, PS) in the domain.
The modelled distributed grid-cell runoff from Alpine3D is used as input for the runoff routine of the HBV-ETH model, in order to produce the total discharge at the location of the Pegelstation Vernagtbach gauge measurements.
4. Results
We focus on the evolution of the snow and glacier resources in the Vernagtferner basin during the three mass balance years from autumn 2019 until autumn 2022. Conditions during these years are rather different and cover a large part of the recent variability of mass balances of Vernagtferner (Table 2). While the annual mass balance of 2019/20 is almost equal to the 30-year mean value (1992–2022), the mass balance of 2020/21 is considerably less negative and the mass balance of 2021/22 is the most negative mass balance recorded since the beginning of the measurements in 1964. The winter balances for 2019/20 and 2020/21 are close to the mean winter accumulation, while 2021/22 shows 60% winter accumulation only, compared to the 30-year mean value.
4.1 Evaluation of the simulated snowpack (calibration/validation)
4.1.1 Snowpack evolution at point scale
A critically important input variable for the snowpack evolution is total precipitation (i.e. the water equivalent of precipitation, either liquid or solid). Because total precipitation (parameter PSUM in the models) is not observed at the glacier stations HVP and AWS-A (only snow thickness is measured), it needs to be derived separately. PSUM is modelled with SNOWPACK at the locations of HVP and AWS-A, forced by the locally observed snow thickness and the additional meteorological parameters, as described in 3.2.4. For the summer, we additionally used the directly measured precipitation at PS, which in general accounts for the liquid phase of precipitation. The resulting PSUM time series acts as input for a validation run to evaluate the agreement between resulting modelled snow thicknesses and the measurements at the Ablatometer site (ABL). The main tuning parameter during this process is the snow/rain transition factor T_s, while the other meteorological parameters are constrained by the observations. T_s is manually varied until an optimal correlation is found for simulated and observed snow thicknesses (better than 0.99 at both locations, with a mean bias of 14 cm at HVP and 3 cm at AWS-A, final value T_s = 1.8°C). The error of the snow thickness measurements with ultrasonic rangers mainly originates from the necessary readjustment of the instrument level and to a minor degree from the measurement itself. The combined error is estimated to about 5 cm. The basin-wide Alpine3D simulation is then based on this optimised PSUM time series. The comparison of modelled and observed snow thicknesses, as well as ice melt in the summer season at the ABL station serves as an independent validation of the simulation. The very high correlation between modelled and observed snow thicknesses/ice surface level at the calibration locations HVP and AWS-A (Fig. 3) provides a measure of model quality.
As can be expected, the modelled snow thicknesses are very similar to the observations at the calibration locations. However, during summer 2021 the simulated snow thickness decrease at HVP is stronger than the observations, while the opposite is observed in summer 2022.
SWE/ice melt values at point scale are available from the annual mass balance measurements. Involved errors are in the order of 10%. The modelled SWE/ice melt values agree rather well with the local observations (in general better than 10% and thus within the error range, Table 3), also for the HVP station in 2021. This indicates that the snow thickness offset at HVP station is not related to discrepancies in the modelled and observed mass balance, but rather to differences in snow densities.
4.1.2 Spatial variability
In a next step, the spatial distribution of the snow pack evolution and ice ablation based on Alpine3D simulations is evaluated against available measurements. The snow thicknesses from the regular winter mass balance measurements at Vernagtferner, as well as the accumulation area extent determined at the end of the ablation season, are used for this purpose. The agreement between the manual end-of-winter snow soundings and the modelled local snow thickness is rather good in general (Fig. 4). Across most of the ablation area, the comparison shows only rather minor differences of less than 18 cm on average. However, there are some characteristic larger discrepancies: a few regions close to the main ridge (e.g. Sexenjoch and Taschachhochjoch; see Fig. 1 for location names) show considerably higher measured snow thicknesses compared to the model result (mean difference of 35 cm). Another discrepancy occurs in the upper part of the Schwarzwand region, where observations are decidedly lower than the simulated results.
The differences in simulated and observed snow thicknesses at the end of winter directly translate to the spatial distribution of the accumulation area at the end of the ablation season. The region with too low modelled snow thicknesses result in a locally smaller accumulation area, while the region with larger model snow thicknesses leads to a larger extent of the accumulation area, compared to the observations. In general, however, the spatial pattern of snow thickness and ice ablation distribution is captured rather well.
The glacier-wide, annual mass balances (Table 2) are modelled reasonably well, with an overestimation of the magnitude by 24% in 2019/20, 20% in 2020/21 and 4% in 2021/22. The absolute offsets of 202 mm and less are smaller than the general error of the mass balance measurements. The winter balances are captured better, with an underestimation of 13% in 2020/21, but only 5% in 2019/20 and 4% in 2021/22.
4.2 Runoff simulations and comparison with discharge
Runoff on a grid-cell scale from Alpine3D is evaluated by combining Alpine3D and the routing routine of the gridded version of HBV-ETH (Mayr and others, Reference Mayr, Hagg, Mayer and Braun2013) and comparing the resulting discharge with the observations at the Pegelstation gauge (Fig. 5). In general, the simulated discharge correlates very well with the observations from the gauging station at Pegelstation (orange and red lines in Fig. 5), with a Nash-Sutcliffe efficiency coefficient of 0.924. The comparison shows, however, that snowmelt in the early melt season (May and June) does not immediately lead to an increase in river discharge.
A larger discrepancy is observed in early summer 2021 (mid-June to mid-July), when the associated discharge from modelled snow and ice melt is considerably higher than the observations.
Modelled and measured discharge can only be compared for measuring periods at the gauging station PS and therefore not during winter (Table 4). The comparison of the seasonal sums shows that there is a rather good agreement between the simulation of the discharge in the river and the measurements at PS, even though offsets can be up to 14%. However, the accuracy of the discharge measurements is estimated to be at best ± 10% and, due to COVID, calibration of the rating curve was not possible in 2021. The ratio between the modelled and measured runoff shows an excess of available water during the early melt period. The total annual basin precipitation is calculated from the reconstructed PSUM values at the monitoring sites and results in 21.92 (2020), 22.99 (2021) and 15.70 million m3 (2022).
Also, the total modelled annual discharge (Q_HBV_t) in the balance year (01 October – 30 September) and the summer precipitation at PS is shown.
4.3 Melt characteristics in the glacierised basin
The individual runoff components show different temporal characteristics. Therefore, we analysed the composition of the total runoff in terms of these components, especially during the ablation season.
There exists a typical pattern of runoff contributions over the balance year (Fig. 6). Starting with the onset of melt, usually in May, the snow cover begins to reduce. During the first weeks, large amounts of meltwater will be likely retained in the depleted groundwater storage. Direct discharge of almost the entire melt into the river is reached after about a month. Ice melt starts when the snow line reaches the lowermost glacier margin, usually in July. The snowmelt volume decreases rapidly during this time and within about a month, ice melt is the dominating contribution to discharge. Melt from the firn area contributes to the discharge when the snow line migrates up-glacier across the firn line. High discharge values occur during the peak of the snowmelt season, when there is still snow cover in the entire basin. Major peaks occur as well when the largest part of the glacier is snow free in July and August. The contribution of rain generally does not play a major role on a daily basis during the ablation season, even though high, short-term peaks can occur during summer thunderstorms or heavy precipitation events (e.g. 29.08.2020, where rain contributes more than 50% to the direct discharge). Melt water from snow, firn and ice usually strongly reduces around mid to end of September, when the ablation period comes to an end.
4.4 Differences of interannual storage (snow, firn, ice, rain)
The storage magnitude of individual melt units depends on the seasonal weather patterns and can be rather different from year to year. The timing and magnitude of snowfall and melt characterize the runoff dynamics and seasonality of high alpine basins and are thus relevant for water management planning (He and others, Reference He, Duethmann and Tian2021). The storage size of different units and their interannual evolution is of high interest, especially with respect to the glacier reaction on climatic changes.
The role of ice melt for the total discharge from the glacier basin increases rapidly once the lower parts of the glacier are exposed to the atmosphere. The relative contribution of firn to the total discharge is less than 2% in the investigated period. However the total firn melt is about twice in 2022, compared to 2020 and 2021 (0.58 × 106 m3 compared to 0.33 × 106 m3 and 0.24 × 106 m3, respectively). Given that the recent firn area (2021) of Vernagtferner only covers a size of 1.5 km2, the volume loss equates to a mean height loss of 0.38 m w.e. in 2022, 0.22 m w.e. in 2020 and 0.16 m w.e. 2021. The cumulative loss of 0.76 m w.e. (or 1.38 m firn) is considerable, considering that the remaining firn body at Vernagtferner is probably less than ten metres thick on average. At the same time, it is clear that the high discharge amounts in 2022 are based on a massive loss of glacier ice and the ice contribution to the total runoff exceeds 75% (Fig. 7).
Today, snow on the glacier accounts for 63% of the total snowmelt in the basin (see Fig. 8 for the temporal distribution), even though only 58% of the basin is covered by glacier, because of the lower snow thicknesses in the low and nonglacierised areas. However, the snowmelt contribution from the entire basin to the total discharge varies strongly from year to year and shows a similar pattern as the snowmelt contribution from the glacierised area, contributing 42% (2020), 47% (2021) and 24% in 2020 (2022).
5. Discussion
The model simulations provide a detailed insight into the dynamics of the runoff components of the glacierised catchment of Vernagtferner, but also reveal valuable details about the suitability of the chosen model setup, as well as the consistency of the observations.
Snow thickness and ice surface level at the observation sites HVP and AWS-A were used for calibrating the total precipitation for the model period, because continuous observations of SWE do not exist. The validation at the Ablatometer site reveals a very good agreement between simulation and observations (Fig. 3). The net difference after three years is less than one metre or about 7% at the ABL site. However, the results for HVP show some systematic offset during the ablation periods in 2021 and 2022. These discrepancies are likely due to differences in the compaction of multi-year firn layers in the model compared to reality in summer 2021. In summer 2022, the entire firn layer down to the reference level is removed at HVP by melting, including the firn pack with a potentially different degree of compaction, and thus the difference in snow thickness disappears again. This possible source of deviation is supported by the annual point mass balance observations, which agree very well with the modelled SWE values (Table 3). Therefore, we conclude that Alpine3D simulates the SWE correctly, while the density distribution in the snow and firn pack does not agree completely with the observations. A possible reason could be the initialisation of the model, where deeper firn layers are not included due to a lack of information for these layers. The quality of the model results is also confirmed by comparing the seasonal and annual mass balance values, where the deviations are smaller than the error of the balance calculations.
The correct spatial distribution of snow accumulation, as well as snow, firn and ice melt is important for the temporal evolution of meltwater production. The annual accumulation measurements at Vernagtferner confirm that the model simulations represent the real conditions rather well, even though some regions show clear deviations. At the regions of underestimated snow thicknesses, the differences are very likely to originate from poor representation of snowdrift across high ridges in the Alpine3D model due to a lack of a detailed wind field. However, this offset is not observed at the HVP station, because the station data were used for model calibration at this high-elevated plateau. The likely reason for the underestimation of snow accumulation in the upper Schwarzwand region is high wind erosion in this section of the glacier, which results in very low snow thicknesses compared to other regions at similar elevations. These local differences in the accumulation pattern lead to a slightly different extent of the accumulation regions, but the influence on the general mass balance and thus the meltwater production are small.
The comparison of the modelled discharge at the gridcell scale (melt and rain) with the measured runoff at the PS gauge allows us to investigate the basin characteristics. In general, the simulations agree rather well with the observations and the differences of the annual sums are within the measurement errors. This indicates that the combination of Alpine3D and the site-calibrated, distributed HBV-ETH model enables a realistic representation of dynamic processes in a highly glacierised catchment.
Still, there exists a characteristic offset during the early melt season, where the simulated discharge values are generally higher than the measurements. This is also indicated by the ratio of the glacier runoff and the measured gauge discharge, which shows a clear signal of water excess in the first weeks of the melt period. We interpret this spring surplus as a likely source to replenish the groundwater storage, which drained during the preceding winter. There is a rising awareness about the role of groundwater flow for streamflow compensation during drought periods (e.g. Somers and others, Reference Somers2016; Saberi and others Reference Saberi2019; Hayashi, Reference Hayashi2020; Somers and McKenzie, Reference Somers and McKenzie2020). However, it is difficult to quantify the groundwater contribution to streamflow, as well as the renewal rates of groundwater during the course of the year, especially in high alpine catchments with strong flow variations (Hayashi, Reference Hayashi2020). There exist some dedicated investigations for the Vernagtferner basin, which tried to disentangle the composition of the flow already in the 1970s (Behrens and others, Reference Behrens, Moser, Oerter, Rauert and Stichler1979; Oerter and others, Reference Oerter, Baker, Moser and Reinwarth1981). The glacier covered about 84% of the basin at that time and the analysis of isotope concentrations in the water indicated that the contribution of groundwater to streamflow was on the order of 0.01–0.1 m3 s−1, depending on the season (Behrens and others, Reference Behrens, Moser, Oerter, Rauert and Stichler1979). The winter streamflow was assumed to consist of 100% groundwater, while the maximum contribution to the mean summer discharge of 1–2 m3 s−1 was in the order to 5–10%. The reason for such a rather small contribution is the very limited storage volume in the proglacial sediment deposits (Moser and others, Reference Moser, Escher-Vetter, Oerter, Reinwarth and Zunke1986). The glacier cover in the basin has since reduced from 84% to 58% and additional sediment deposits are now uncovered. These deposits might already have contributed to the groundwater storage in the 1970s, because Vernagtferner is a temperate glacier with abundant basal water flow. The water excess between simulated runoff and measurements in spring amounts to about 2–3 × 106 m3 in each year. Unfortunately, low flow and thus winter discharge cannot be determined precisely at the gauge, but estimates result to about 1.6 × 106 m3 between October and May. This is similar in magnitude to the estimated spring melt excess, indicating that refilling the groundwater storage in spring explains the majority of this discrepancy. The role of the groundwater storage for compensating low flow periods is not as large as has been estimated for other catchments (Hayashi, Reference Hayashi2020) due to the rather limited storage volume compared to the total discharge volumes.
The time of complete filling of the groundwater storage is obviously related to the melt intensity. The early summer conditions in 2021 were considerably warmer than in 2022, even though 2022 was the most negative balance year so far. The positive degree day sum from the beginning of the snowmelt until the modelled melt equals the measured discharge is almost 40% higher in 2021, compared to 2022 (207°C vs 149°C). Also, the mean incoming short wave radiation is about 28% higher in 2021 compared to 2022 for the respective period. This results in considerably higher simulated melt rates during the early summer 2021. However, it seems that the HBV-ETH model (with the calibration of Mayr and others, Reference Mayr, Hagg, Mayer and Braun2013) cannot account for such high water fluxes into the storage system and produces a too high discharge, which is not observed at the PS gauge in early summer 2021.
Remarkably, the match between meltwater generation and runoff (blue and red curves in Fig. 5) and thus the assumed end of the refilling phase occurs for rather similar snow cover conditions on the glacier (Fig. 9), independent of the date.
The three ablation periods show rather different characteristics: the central glacier tongue started to become snow free around the 6 July in 2020 and the 8 July in 2021, but already at the 12 June in 2022. The snow line roughly reached the long-term equilibrium line on 14 of August in 2021 and the 16 of August in 2020, but in 2022 already on the 2nd of July. This demonstrates that the ice ablation period starts almost a month earlier in 2022 compared to the two previous years, and ice melt affects the major part of the glacier even more than a month earlier.
Even though the glacier mass balance was almost 40% more negative in 2020 compared to 2021 (Table 3), the modelled discharge in the two years is very similar. The lower contribution from glacier melt to discharge in 2021 is very likely compensated by considerably higher summer precipitation compared to 2020 (Table 4).
It is also noteworthy that the excess discharge (i.e. the difference between annual basin precipitation and annual basin discharge) and thus the contribution of glacier imbalance varies greatly in the three years. Excess runoff shows large differences from 7.4% in 2021 to 57.1% in 2022, while 10.6% of runoff in 2020 is due to uncompensated glacier loss. Previous work showed that the long-term excess discharge from the Vent catchment (glacier cover 35%), which includes the Vernagt basin, is 9% for the period 1969–1997 (Lambrecht and Mayer, Reference Lambrecht and Mayer2009). Even though this result cannot be directly compared to the Vernagtferner basin due to the different glacier cover (58% compared to 35%), it supports the finding that 2022 was an outstanding melt year for the Vernagtferner basin, where the majority of basin runoff was due to a net loss of glacier ice.
The general seasonality of the discharge dynamics of a high-alpine drainage basin is characterised by the total amount of winter snow accumulation, the snowmelt in spring, the glacier ice melt during the summer and episodic strong rain events. In addition, the groundwater share of total discharge is highly variable (Hayashi, Reference Hayashi2020). The contribution of the different cryospheric storage units to melt is of particular interest, with respect to glacier mass balance. In years with balanced conditions, the remaining amount of snow on the glacier at the end of the ablation season should equal the melted glacier ice. For negative mass balance years, not only snow but also firn is removed from the upper regions of the glacier. Quantifying the individual contributions to the discharge by direct measurements is difficult (e.g. Oerter and others, Reference Oerter, Baker, Moser and Reinwarth1981; Baraer and others, Reference Baraer2015), whereas a well-calibrated snow and ice model provides such information.
Our analysis of the discharge characteristics focuses on the ablation seasons 2020 and 2022 as representatives of a mean and an extreme melt year. The basin wide snowmelt (green line in Fig. 6) and the contribution from melt on the glacier to discharge (blue line, defined as Q 0 in Fig. 6) show a rather different seasonality. Melt starts in the middle of May for both years, but it is restricted to lower elevations and short periods in 2020, while there is already a strong initial melt event in 2022 (Fig. 6). The major snowmelt period is reached at the end of June/beginning of July in 2020 for about two and a half weeks, followed by a phase of reduced melt due snowfall between 15 and 17 July. In contrast, snowmelt is almost complete by the end of June in 2022, due to considerably less winter accumulation (539 mm w.e. vs 882 mm w.e.) and a persistent, intensive melt period. The model results show clearly that ice melt plays the dominant role in runoff generation for the major part of the ablation season in 2022, from the beginning of July until mid-September. In 2020, the dominance of ice melt starts around mid of August, but at a low rate, while the ablation period is about one week longer into autumn. A strong precipitation event occurred from the 28 August until the 30 August in 2020, with a total rainfall of more than 70 mm. This event led to the highest discharge values of the entire ablation season, but ended with snowfall and a strong reduction of discharge during the following days. The year 2022 shows almost twice as high mean discharge values in the main ablation season (almost 3 m3 s−1 on average from mid-July until the end of August), compared to 2020. Rain, however, plays a negligible role, apart from some single events at the end of June. The continuously high Q 0 values during August 2022 are in line with the stable sunny weather conditions during this period. However, the shape of Q 0 in 2022 also demonstrates that the melt potential in August is already decidedly lower than in June and July, when the highest melt peaks are reached, due to less radiative energy. In contrast, melt is considerably reduced by unfavourable weather conditions during August 2020 with snowfall. The highest discharge peak in the end of August is dominated by rain, rather than melt. This analysis is in line with earlier investigations to the compensating effect of glacier melt on streamflow (Pohl and others, Reference Pohl, Gloaguen, Andermann and Knoche2017; van Tiel and others, Reference van Tiel, Van Loon A, Seibert and Stahl2021). If 2020 represents a mean post-2000 mass balance year, the contrast to 2022 is obvious. Ice melt clearly overcompensates the lack of snowmelt due to reduced winter precipitation, with about 80% higher August discharge values in 2022, compared to 2020, while the accumulation area ratio (AAR) is reduced to zero at the end of the balance year.
A closer look on the meltwater composition reveals some insight into the cryospheric reservoirs on the glacier. Firn contributes to the discharge only in negative mass balance years and its role is rather minor (Fig. 7). However, the extent of the firn is important for long-term mass balance simulations, as firn has a different albedo and density than snow or ice (Naegeli and Huss, Reference Naegeli and Huss2017). In this context, the extreme loss of firn cover in 2022 represents a strong break in the mass balance evolution, because a considerably larger dark ice area is exposed after snowmelt in the future and thus influences the melt conditions, especially for a further changing climate.
It is worthwhile to note that, even for years with higher snowmelt amounts compared to the total ice melt (e.g. in 2021, Fig. 7), the glacier experiences a negative mass balance. This indicates that the relative snowmelt contribution to the total mass loss must be well above 60% under neutral-balance glacier conditions. Such conditions could be met either by high winter snowfall, or less summer melt. The simulations confirm that melting ice today accounts for the largest share of runoff in the glacierised part of the catchment, while snow contributes less than 50% in average years and as little as 25% in years with extreme melt. Only for winters with clearly above-average snow accumulation, does snowmelt dominate the runoff during the ablation season. This result is in line with observations from past periods with positive mass balance conditions. In 1978, ice melt contributed only about 7.5% to the summer runoff, while the snow melt contribution was 80% (Behrens and others, Reference Behrens, Moser, Oerter, Rauert and Stichler1979), even though the glacier cover was 84% of the catchment at that time, compared to 58% today. The mass balance in 1977/78 was +288 mm w.e. and the annual runoff was about 10 × 106 m3, compared to almost 30 × 106 m3 during the last years. Similar conditions existed until the middle of the 1980s, where both, the snow and firn contribution to discharge was always higher, compared to ice melt, except for the slightly negative balance years 1982 and 1983 (Escher-Vetter and Oerter, Reference Escher-Vetter and Oerter2013), when ice melt was slightly higher than the two other components. It seems evident that there exists a general relationship between the snow fraction of the total melt and the glacier mass balance, because the snow layer thickness in spring and the available melt energy during summer both influence the speed and degree of snow cover reduction and thus the final magnitude of ice loss. This is also in line with a rather well established correlation between total discharge and the AAR (Fig. S3), which indicates that for strong imbalances (AAR < 0.5) the discharge is considerably higher than the mean annual basin precipitation and thus is dominated by ice loss. It would be interesting to see similar relationships computed for different glacierised drainage basins.
We also observe that the contribution of firn to the melt does not increase steadily during the ablation season, but reaches its maximum shortly after the start of firn melt at the beginning of September, even though additional firn areas might be exposed by further snow ablation (Fig. 8). This is probably due to the already low sun elevation and thus generally reduced melt rates at the end of the ablation season. Compared to 2020 and 2021, the firn loss in 2022 is extensive, which leads to a larger area of glacier ice which is exposed and thus a future decrease of glacier albedo.
We did not test the relative impact of individual parameter uncertainties on the quality of the model results. However, we recognise that local observations considerably improve the reliability of the model output. The possibility to reconstruct the total precipitation at different elevations and locations of the domain is especially important to obtain a rather realistic distribution of the SWE evolution over time. This is also confirmed by a detailed study on the relevance of input parameter errors on model output quality for Alpine3D (Schlögl and others, Reference Schlögl, Marty, Bavay and Lehning2016), which emphasises that station coverage in the basin is important for obtaining reliable results. This study also demonstrates that typical input uncertainties have a similar effect on Alpine3D results as uncertainties of SWE measurements. Our high-resolution basin geometry and the use of several monitoring stations with high-accuracy sensors enable the physical model to produce robust results. This is also confirmed by the validation of the results against measured snow thickness, SWE and discharge, which typically show uncertainties in the range of the observation error. A critical parameter is liquid precipitation, because of known problems of undercatch and site-representativity. However, liquid precipitation plays only a small role for the general discharge and the mass budget of the basin due to the infrequent occurrence of rain conditions during summer. Even though the observed summer precipitation amounts seem to fit well into the seasonal and annual water budget, we concentrate our study on the comparison of cryospheric components of the discharge. The measured snow thicknesses and SWE values have a considerably smaller error than rain gauge observations and show a much better spatial distribution. Based on our analysis and the findings of Schlögl and others (Reference Schlögl, Marty, Bavay and Lehning2016), we conclude that the presented results are robust and representative for the conditions in highly glacierised Alpine basins.
6. Conclusions
Our results demonstrate that melt and discharge in and from glacierised basins can be simulated rather well with the SNOWPACK/Alpine3D model and an appropriate calibration using local meteorological variables. We relied on field observations for model forcing, while total precipitation (rain and snow) was reconstructed from observed snow thickness and additional meteorological data in a SNOWPACK model run. Subsequent evaluation of local snow pack evolution, as well as the seasonal, glacier wide mass balance distribution showed a good agreement with the observations. Some minor deviations in the spatial distribution are likely connected to wind redistribution of snow at certain critical areas and indicate that wind induced deposition patterns require special attention. This is in agreement with other investigations (e.g. Freudiger and others, Reference Freudiger, Kohn, Seibert, Stahl and Weiler2017; Mott and others, Reference Mott2023), who emphasise the importance of snow redistribution.
The observed mismatch between modelled discharge and the measured runoff during the early melt season is assumed to primarily replenish the seasonal to multi-seasonal water storage units (groundwater recharge). This mechanism reduces the discharge at the gauge during the early snowmelt period for six to eight weeks, after which the glacier has lost 30–40% of its snow cover. The subsequent discharge is resolved very well by the model runs. Winter runoff can only be estimated, but is in the order of 0.01 to 0.1 m3 s−1 and is mainly attributed to the depletion of the groundwater storage in the absence of other water sources during this season.
Snowmelt dominates discharge for the first weeks of the melt period. However, even the highest early-season, basin-wide snowmelt peaks reach hardly the magnitude of midsummer glacier melt peaks, despite the much larger melt area. The annual contribution of snowmelt to discharge during recent years is less than 50% in general. The mean annual glacier melt is decidedly higher than the mean snowfall in the basin, which is in line with the consistently negative mass balances during the last decades. The mean melt share from the glacier area of 73% in 2020–2022 is comparable to earlier studies, where the contribution ranges from 6–15% for basins with 5% glacier cover (Jost and others, Reference Jost, Moore, Menounos and Wheate2012; Soruco and others, Reference Soruco2015) to 55% for a glacier cover of 44% (Young and others, Reference Young2021). The melting of the still extensive glacier surface clearly overcompensates the lack of precipitation as has been observed in other highly glacierised regions (e.g. Pohl and others, Reference Pohl, Gloaguen, Andermann and Knoche2017; van Tiel and others, Reference van Tiel, Kohn, Van Loon and Stahl2020b). The groundwater contribution to summer discharge is, however, small (less than 5%) even in dry years, in contrast to other, less glacierised catchments (Hayashi, Reference Hayashi2020; Somers and others, Reference Somers and McKenzie2020). Even though the snow melt magnitude varies strongly during June and July, depending on winter snow amounts and weather conditions, ice melt controls the runoff during August and September even for years with high winter precipitation. More than 80% of the total melt or about 62% of the total discharge originated from glacier ice in 2022, which emphasises the dominating role of glacier melt for discharge. This is in strong contrast to observations a few decades ago with slightly positive mass balance conditions, when ice melt played only a small role in the total discharge (Behrens and others, Reference Behrens, Moser, Oerter, Rauert and Stichler1979). The hydrological character of Vernagtferner basin has therefore changed fundamentally from a snowmelt dominated regime to a regime dominated by ice melt intensity on the glacier today, with persistently high river levels in July and August, sometimes even into September.
The evolution of the firn body might be crucial for a realistic simulation of future glacier response to climate variations, even though the contribution of the firn area to discharge is negligible. The firn body is a dominant part of the accumulation area and a loss of firn extent during extreme melt years (e.g. 2022) leads to the exposure of larger areas of glacier ice with a lower albedo.
Our detailed analysis shows that the chosen approach is able to simulate the hydrological balance of a glacierised basin rather well and with a high spatial (20 m) and temporal (1 h) resolution. Local meteorological observations at the glacier are essential as reliable input for the calculation of snow pack formation and snow and ice melt.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/jog.2024.48
Acknowledgements
Adrien Michel was highly supportive, even when asking basic questions about SNOWPACK/Alpine3D. We are also very grateful for support from the SNOWPACK/Alpine3D team. Financial support from the Bavarian State Ministry for the Environment and Consumer Protection (grant TUS01UFS-77317) is gratefully acknowledged. We sincerely thank Evan Miles for improving the English language. We thank Lizz Ultee and two anonymous reviewers, as well as the scientific editor Shad O'Neel for providing highly valuable recommendations, which considerably improved the quality of the manuscript.
Data and code availability
The source codes of MeteoIO, SNOWPACK and Alpine3D are available at https://gitlabext.wsl.ch/public. The code of the distributed HBV-ETH version is available on request only by Elisabeth Mayr (Mayr and others. 2013). The HBV-light version is available at https://www.geo.uzh.ch/en/units/h2k/Services/HBV-Model/HBV-Download.html. Real-time data of the monitoring stations at Vernagtferner are found at https://geo.badw.de/vernagtferner-digital.html. The meteorological and discharge data are available at https://doi.pangaea.de/10.1594/PANGAEA.775113 until 2013. The more recent data will be added to this repository in due time.