1. Introduction
Modelling melt rates across an Alpine glacier is essential for studying the intra- and interannual variability of snow and ice cover and for the estimation of the total runoff from glacierized areas. Melt can be computed by means of physically based energy-balance models or by more empirical models, usually known as temperature-index models (see Reference HockHock, 2005; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005, for a definition). Temperature-index models have been used in a wide range of applications (e.g. Reference Huybrechts, Letréguilly and ReehHuybrechts and others, 1991; Reference JóhannessonJóhannesson, 1997; Reference OerlemansOerlemans and others, 1998; Reference Braithwaite and ZhangBraithwaite and Zhang, 2000; Reference Zappa, Pos, Strasser, Warmerdam and GurtzZappa and others, 2003), and a number of enhanced versions of this approach have been proposed (Reference Cazorzi and Dalla FontanaCazorzi and Dalla Fontana, 1996; Reference HockHock, 1999; Reference Daly, Davis, Ochs and PangburnDaly and others, 2000; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). These provided important improvements over the simple standard version, and their performance can be compared to those of energy-balance models when applied to one season or a single catchment (Reference HockHock, 1999; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). However, a criticism often directed at this type of approach is that such models depend on empirical parameters that are calibrated for a single glacier or ablation season and are therefore not transferable. For this reason, their use in the simulation of the impact of climate change on glaciers has been questioned. Investigating the parameters’ transferability and model generalization helps to address this criticism.
The aim of this study is to test the transferability and robustness of an enhanced temperature-index (ETI) model of glacier melt that was developed for Haut Glacier d’Arolla, Switzerland, using data measured on the glacier at five automatic weather stations (AWSs) during the 2001 ablation season (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). The ETI model includes the shortwave radiation balance in addition to air temperature, and is therefore more physically based than the classical method which uses air temperature alone. It could also be regarded as a simplified energy-balance model in which the temperature-dependent components of the energy balance (longwave radiation and sensible heat fluxes) are parameterized in one term, whereas the shortwave radiation flux, which is the dominant source of energy for Alpine glaciers, is accounted for explicitly (Reference Greuell, Genthon, Bamber and PayneGreuell and Genthon, 2004). It still includes, however, a threshold temperature that remains an empirical element of the model.
In the original paper (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005), the model parameters (temperature factor, TF, and shortwave radiation factor, SRF) were optimized for one particular AWS location by calibrating the ETI model output against output from a physically based energy-balance model, which had itself been tested against measured melt rates. The ETI model was then applied to other AWS locations and tested against the output from the physically based energy-balance model run at those locations.
This paper investigates the transferability of the original model (calibrated against energy-balance model output at one location on one glacier over one melt season) in terms of its success at reproducing energy-balance model output: (1) at other locations on the same glacier in the same year; (2) at the original location on the same glacier but when the year is divided into distinct time periods containing different surface characteristics (e.g. snow vs ice) or different weather types (e.g. sunny vs cloudy); (3) at locations on the same glacier in different years; and (4) at locations on other glaciers in different years.
The model transferability is assessed in two ways. First, the original model with the original parameter values (i.e. the parameters of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005) is run for the other locations, time periods and glaciers and the outputs are compared with those from the energy-balance model. The Reference Nash and SutcliffeNash and Sutcliffe (1970) efficiency criterion, R 2, is used to compare the results. Second, the model’s parameters are recalibrated the for the new locations, time periods and glaciers, and the model performance obtained in this way is compared with that using the original parameters. We also compare the new parameter values with the original ones and explain the observed changes by relating the values of the parameters to the predominant meteorological conditions and partition of the energy fluxes at the glacier–atmosphere interface. For this purpose, we analyse in detail the components of the energy balance at the locations on Haut Glacier d’Arolla.
A second aim of this work is to test the transferability of the ETI model version in which both albedo and incoming shortwave radiation are modelled (model E of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). This uses only air temperature as input data and is therefore comparable to the standard degree-day approach and to the model of Reference HockHock (1999) in terms of data requirements.
The ETI model is calibrated against reference hourly melt rates computed by a physically based energy-balance model (Reference Brock and ArnoldBrock and Arnold, 2000), which in turn is tested against ablation readings obtained from measurements of surface lowering recorded at ultrasonic depth gauges (UDGs). The energy-balance model is used for calibration and validation of the ETI model because ablation readings at stakes are affected by significant errors (Reference Müller and KeelerMüller and Keeler, 1969, Reference MunroMunro, 1990) and are not available at hourly resolution, while UDG measurements, although continuous in time, are often noisy and provide a good test for the cumulative ablation but not for high temporal resolution (i.e. hourly) melt rates. Both types of data are also affected by time lags between melt and its expression as a lowering of the surface (Reference MunroMunro, 1990). This follows the strategy adopted by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005). This study is conducted at the point scale, in order to avoid inaccuracies introduced by the extrapolation or parameterization of the meteorological input data. The ETI and energy-balance models are run using input data measured at the location of the AWS, and simulations at hourly resolution are used. The dataset used is the same as that of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) for the 2001 ablation season, plus meteorological data recorded on Haut Glacier d’Arolla in the 2005 and 2006 seasons at the location of the lowest automatic weather station (Fig. 2), and data collected on both Gornergletscher and Tsa de la Tsan glacier (Fig. 1).
Fig. 1. Location of the three glaciers. The thick lines indicate the glacier borders.
Fig. 2. Map of Haut Glacier d’Arolla showing the position of the five AWSs installed in 2001 and of the fixed station in the proglacial valley. Background image is a relief-shaded digital elevation model of 10 m resolution showing the glacier border (solid curve) and the watershed (dashed curve). Northing and easting are in the Swiss coordinate system. The station names correspond to station numbers as follows: central station (17), uppermost station (15), north-central station (13), south-central station (14) and lowest station (16) (see Table 2).
2. Study Sites
Haut Glacier d’Arolla is a valley glacier in the southern Alps of Switzerland, with an elevation range of ∼2560–3500 m a.s.l. and an area of ∼6.3 km2. The glacier consists of an upper basin with northwesterly aspect that feeds a glacier tongue flowing to the north (Fig. 2). The glacier surface meteorology, melt and hydrology has been extensively investigated in the last decade (e.g. Reference Arnold, Richards, Willis and SharpArnold and others, 1998; Reference Brock and ArnoldBrock and others, 2000, Reference Brock, Willis and Sharp2006; Reference Willis, Arnold and BrockWillis and others, 2002; Reference Strasser, Corripio, Pellicciotti, Burlando, Brock and FunkStrasser and others, 2004; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005, Reference Pellicciotti, Carenzo, Helbing, Rimkus and Burlando2009). A glaciological program currently being conducted on the glacier has meant that meteorological and ablation data are available over several seasons at the same location, which is optimal for the aim of this work.
Gornergletscher, also in the Swiss Alps, is a much larger valley glacier with a total area of ∼57.5 km2 and length of ∼14 km, consisting of several tributaries flowing down from the highest peaks of Switzerland. The glacier is polythermal and is characterized by an extensive debris cover on both sides of the tongue. It is located ∼20 km east of Haut Glacier d’Arolla (Fig. 1), but meteorological conditions differ, especially in terms of precipitation, which is lower (Table 2). Field campaigns were conducted on Gornergletscher in the 2005 and 2006 melt seasons.
Table 1. Characteristics and period of functioning of the AWSs at the three study sites: Haut Glacier d’Arolla (2001, 2005 and 2006), Gornergletscher (2005 and 2006) and Tsa de la Tsan glacier (2006). HGdA denotes Haut Glacier d’Arolla. The coordinates are in the Swiss coordinate system (easting and northing of Tsa de la Tsan glacier in the Italian UTM ED50 32 N system are 388785 and 5093025, respectively)
Table 2. Main features (climatic conditions and surface and topography characteristics) of the AWS locations at the three study sites. HGdA denotes Haut Glacier d’Arolla. T is air temperature, h is height of precipitation and I is incoming shortwave radiation
Tsa de la Tsan glacier, in the upper Valpelline valley in the Italian Pennine Alps, is also located close to Haut Glacier d’Arolla (Fig. 1), on the other side of the Col du Mont Brulé, Pointe de la Grande Arête and Bouquetins ridge across the Italian border. The glacier has undergone a significant recession and only its upper basin remains. The glacier has an area of ∼2.9 km2, an elevation range of 3040–790 m a.s.l. and prevailing southwest exposition. A field campaign was carried out on the glacier in 2006 (Table 1). The three glaciers differ substantially in terms of area, surface characteristics, elevation range and meteorological conditions, thus offering a perfect case study for the analysis of the model transferability and application to other sites.
3. Data
The meteorological and ablation datasets collected in the 2001 ablation season on Haut Glacier d’Arolla are described in detail elsewhere (Reference PellicciottiPellicciotti, 2004; Reference Strasser, Corripio, Pellicciotti, Burlando, Brock and FunkStrasser and others, 2004; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). Five AWSs set up in the upper basin and on the glacier tongue for the duration of the melt season (Fig. 2; Table 1) provided 5 min averages of air temperature (°C), relative humidity (%), incoming and reflected shortwave radiation (W m−2), wind speed (m s−1) and direction (°). All data were recorded at a nominal height of 2 m, and incoming and reflected shortwave radiation were measured parallel to the surface (Reference Greuell, Genthon, Bamber and PayneGreuell and Genthon, 2004; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). An ultrasonic depth gauge was installed at the location of the central station, and recorded measurements of surface lowering every 5 min. Snow density measurements were also carried out at the central station (see Reference PellicciottiPellicciotti, 2004, for details). In the 2005 and 2006 seasons, only the lowest station operated (Table 1). In both years, the AWS was equipped with a UDG, which, combined with density measurements, provided the data necessary for validation of the energy-balance model. An AWS installed in the proglacial valley provided 5 min measurements of precipitation (Fig. 2).
One AWS with the same characteristics as those on Haut Glacier d’Arolla was set up on Gornergletscher in both the 2005 and 2006 ablation seasons on the tributary Grenzgletscher at ∼2600 m a.s.l. (Table 1). The station recorded the same variables as at Haut Glacier d’Arolla at the same temporal resolution. The set-up also included a UDG. An AWS sitting on the slopes facing the glacier recorded hourly precipitation.
One AWS was installed in 2006 on Tsa de la Tsan glacier at 3250 m a.s.l. along the glacier flowline. The station functioned from 26 July to 10 October (Table 1). Air-temperature measurements were not ventilated. A UDG was installed together with the AWS and recorded 5 min observations of surface lowering for the same period of record as the AWS. Precipitation was also measured at the glacier AWS.
Data from all the AWSs were treated in the same way: 5 min data were aggregated into hourly data which were used as input to the melt models. UDG observations of surface lowering were smoothed and then converted into continuous time series of snow or ice water equivalent (w.e.) (section 4.3).
4. Methods
4.1. The enhanced temperature-index model
The ETI model is described in detail by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005), so only the main features of the model are recalled here. Melt rate, M (mm w.e. h−1), is computed as
where T is air temperature (°C), α is albedo, I is incoming shortwave radiation (Wm−2), and TF and SRF are empirical coefficients, called temperature factor and short-wave radiation factor, expressed in mm h−1°C−1 and mmh−1 W−1 m2, respectively. T T is the threshold temperature above which melt is assumed to occur (1°C). The two empirical parameters were optimized at the central station on Haut Glacier d’Arolla (Fig. 2) during the 2001 ablation season, and their values are TF = 0.05 mm h−1°C−1 and SRF = 0.0094 mm h−1 W−1 m2 (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). Daily albedo is used in Equation (1) (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). Both I and α can be measured or modelled values (see below), depending on the availability of the input data. If both variables are modelled, the ETI model requires no more data than the standard degree-day approach.
The main characteristic of the ETI model, as compared to the standard temperature-index (TI) method, is the inclusion in the melt calculations of the shortwave radiation balance (Equation (1)), which accounts for the fact that incoming shortwave radiation is the main source of energy for melt on Alpine glaciers (e.g. Reference Arnold, Willis, Sharp, Richards and LawsonArnold and others, 1996; Reference Greuell and SmeetsGreuell and Smeets, 2001; Reference Willis, Arnold and BrockWillis and others, 2002). Separation of temperature-dependent and temperature-independent terms, which follows from consideration of the energy balance (Reference Greuell, Genthon, Bamber and PayneGreuell and Genthon, 2004), provides an increase in model performance over the standard TI method and other approaches proposed in the literature (e.g. Reference HockHock, 1999) by reducing the model over-sensitivity to temperature fluctuations (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005).
Two versions of the ETI model are applied, following Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005): their model D, in which measured air temperature, incoming shortwave radiation and albedo are used as input to the model; and their model E, in which only measurements of air temperature are needed, since incoming shortwave radiation and albedo are parameterized. Incoming solar radiation is modelled using a parametric model for clear-sky conditions (Reference IqbalIqbal, 1983; Reference CorripioCorripio, 2003b), together with a cloud-factor parameterization, in which the effect of clouds is computed as a function of daily temperature variations (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). Daily albedo is parameterized as a function of accumulated daily maximum temperature since snowfall (Reference Brock, Willis and SharpBrock and others 2000). The detailed explanation of the albedo parameterization and solar radiation model can be found in Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005). We follow the naming convention of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005). Model D is used to optimize the model empirical parameters, since these depend only on the meteorological and surface conditions typical of the location and not on the parameterizations used for the input variables. Model E, however, is more operationally applicable because of the lower data requirement, and therefore we also tested the performance of this model.
4.2. The energy-balance model
The physically based energy-balance model is described in detail by Reference Brock and ArnoldBrock and Arnold (2000) and was used by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005). The model calculates the energy balance and surface melt rate at an hourly time resolution for a point on snow or ice, and is based on the assumption that there is no conduction of heat into the snowpack or glacier (zero-degree assumption). The model is forced by hourly measurements of incoming shortwave radiation (W m−2), reflected shortwave radiation (W m−2), air temperature (°C), air vapour pressure (Pa) and wind speed (m s−1).
Surface melt rate, Q M, is computed as residual in the surface energy-balance equation:
where Q I is the net shortwave radiation flux, L is the net longwave radiation flux, Q H is the turbulent sensible heat flux and Q L is the latent heat flux (Reference Brock and ArnoldBrock and Arnold, 2000). Fluxes directed toward the surface are assumed to be positive (e.g. Reference Röthlisberger, Lang, Gurnell and ClarkRöthlisberger and Lang, 1987). In this study, Q I is computed from measurements of incoming and reflected shortwave radiation. L is the sum of incoming (L↓) and outgoing (L↑) longwave radiation, both modelled. Assuming that the surface is at 0°C and radiates as a black body (emissivity equal to 1), L↑ is a constant flux-, time- and space-independent, and equal to −316 W m−2 (Reference OkeOke, 1987). L↓ is calculated from the Stefan–Boltzmann relationship, in which the emissivity is a function of temperature, cloud type and cloud amount (Reference Brock and ArnoldBrock and Arnold, 2000; Reference PellicciottiPellicciotti and others, 2008).
Q H and Q L are computed using the bulk aerodynamic method, which requires wind speed, air temperature and humidity to be measured at only one height above the surface (usually 2 m). The approach is discussed in detail by Reference MunroMunro (1989), Reference Braithwaite, Konzelmann, Marty and OlesenBraithwaite and others (1998), Reference Brock and ArnoldBrock and Arnold (2000) and Reference Denby and GreuellDenby and Greuell (2000). The two fluxes depend additionally on the stability-correction factors for momentum, heat and humidity, the Monin–Obukhov length scale (Reference ObukhovObukhov, 1971) and the scaling lengths for aerodynamic roughness (z 0), temperature (z t) and humidity (z e). In the model, z t and z e are computed as functions of z 0 using the roughness Reynolds number, Re*, following Reference AndreasAndreas (1987). The aerodynamic roughness, z 0, has to be provided to the model, and its evaluation follows the scheme employed by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005), which assigns different values of z 0 to fresh snow, snow after snowfall when melting has taken place, and ice (z 0 = 0.1, 1.0 and 2.0 mm, respectively). These values are well in agreement with values reported in the literature (Reference Brock, Willis and SharpBrock and others, 2006).
The analysis in this paper is conducted using both the original Reference Brock and ArnoldBrock and Arnold (2000) model and an updated version that includes minor changes in the computation of potential incoming shortwave radiation (K in Reference Brock and ArnoldBrock and Arnold, 2000), which is used to derive the cloud amount for the longwave flux calculation. The two model versions do not show any substantial difference in the simulated melt rates (correlation coefficients between simulations with the two model versions are 0.97 at the Haut Glacier d’Arolla central station and above 0.99 at all other Haut Glacier d’Arolla stations). Optimization of the 2001 model parameters with the updated version resulted in a TF of 0.04 and an SRF of 0.0094 (as compared to the original parameter values of TF = 0.05 and SRF = 0.0094; see section 4.1). There is no difference between the two values of SRF, and TF differs very little. In this work, we use the updated model, except for the analysis of model transferability in the 2001 ablation season, which is conducted with the original version of the model in order to ensure a fair comparison with the original work of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005).
4.3. Validation of the energy-balance model
Validation of the energy-balance model is important since hourly melt rates simulated by the model are used as reference to calibrate the ETI model parameters and therefore confidence on their accuracy is fundamental. The energy-balance model has already been tested on Haut Glacier d’Arolla, during the 1993 ablation season (Reference Brock and ArnoldBrock and Arnold, 2000) and for the 2001 ablation season at the central station (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). Here we test the model calculations for other seasons and glaciers when reliable ablation observations are available. Validation is conducted for Haut Glacier d’Arolla in 2006 and Gornergletscher in 2005 and in 2006. No reliable ablation observations were available for Tsa de la Tsan glacier.
Melt simulated by the energy balance is compared with continuous readings of the UDG. The surface lowering observations at the UDG are smoothed to remove noise using a Hamming window (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005) and are converted into mm w.e. using snow or ice density (assuming ρ = 900 kg m−3 for ice, ρ = 100 kg m−3 for fresh snow and ρ = 500 kg m−3 for old snow, from measurements on Haut Glacier d’Arolla; see Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005).
Figure 3 compares energy-balance (EB) simulations with UDG measurements for Haut Glacier d’Arolla in the 2006 season and for Gornergletscher in the 2005 and 2006 seasons. Agreement between energy-balance simulations and cumulative melt from the UDG readings is very good in all three cases (see also Table 3), especially considering the uncertainties in the UDG record and conversion to w.e. These are mainly due to errors in the observations of surface lowering associated with vibration of the sensor, particularly strong for high wind velocities, and to errors associated with the value of snow density (and fresh snow in particular) used to convert the surface lowering into w.e. Additional errors can be introduced by snowdrift and snow compaction.
Fig. 3. Validation of the energy-balance simulations on (a) Haut Glacier d’Arolla in the 2006 ablation season, (b) Gornergletscher in the 2005 ablation season and (c) Gornergletscher in the 2006 ablation season. The three plots show comparisons between cumulative melt simulated by the energy-balance model and measured at the UDG (from measurements of surface lowering converted into water equivalent through density values).
Table 3. Total ablation computed by the energy-balance model (EB) and obtained from UDG measurements for the three sites where validation is conducted. Total ablation is in mm w.e., and is the value at the end of the period for which observations are available
On Haut Glacier d’Arolla, agreement between UDG cumulative melt and energy-balance simulations is very good for the entire season, with a slight overestimation by the model at the beginning of the period. It is lower during snowfalls and days with low temperature (in August and at the beginning of September), when the energy-balance model underestimates melt (Fig. 3a). We observed snowfalls on 4, 12–13 and 27–29 August. Both the UDG record and energy-balance simulations indicate a reduced melt rate in these periods, and the disagreement might be explained by either errors in the conversion of the UDG record into ablation or errors in the energy-balance calculations. The overall agreement between model and measurements, however, is very good, with ∼2% difference in total melt (Table 3).
Comparison of modelled and observed ablation on Gornergletscher from 22 May to 22 July 2006 also indicates a very good agreement between computed and measured cumulative melt (Fig. 3c). Discrepancies are evident on 9 June, 23–24 June and 9 July, when the observed melt rate is smaller than the modelled one, and lower than the average observed melt rate. The energy-balance model predicts higher melt rates than observed, but the discrepancies are minor, and lead to a difference in total ablation at the end of the season of <1% (Table 3). Validation of the 2005 Gornergletscher simulations also shows that agreement between cumulative melt computed with the energy-balance model and calculated from the UDG readings is very good over a period of more than 2 months (Fig. 3b). Agreement between the total melt values at the end of the season is also very good, with a difference of <1% (Table 3), confirming that the energy-balance model can correctly simulate melt over the periods of investigation at different locations on Alpine glaciers. No substantial differences in the comparison of modelled and measured melt at the three sites could be detected, pointing to the general validity of the assumption that the energy-balance simulations can be used as reference melt.
4.4. ETI model recalibration and validation
Four modelling experiments were performed to test the transferability of the model parameters:
-
1. Model transferability across Haut Glacier d’Arolla. The model parameters are recalibrated at the five AWSs on Haut Glacier d’Arolla for the 2001 ablation season. The simulation period covers more than 3 months, from 30 May to 11 September. The experiment allows us to assess the spatial variability of the parameters across a single glacier, and to establish their dependence on site characteristics, such as elevation, location (accumulation or ablation area), topographic factors such as slope and aspect, and surface characteristics such as albedo. The parameter values obtained in this way are mean values typical of the average meteorological conditions of the melt season.
-
2. Robustness of the parameters within one season at one location on Haut Glacier d’Arolla. The optimal values of the parameters are identified for subperiods of the Haut Glacier d’Arolla 2001 ablation season. This analysis is conducted only at the central station.
Three sets of different subperiods are considered. First, we divided the ablation season into two subperiods, depending on whether the surface was snow (30 May–20 August) or ice (21–30 August).
Second, we identified three short subperiods characterized by distinct meteorological conditions. This experiment allows evaluation of the temporal variations in the parameter values associated with changes in meteorological conditions more than with surface and topographic characteristics. It has to be noted, however, that this set of subperiods is not strictly independent of the first one (snow/ice). The subperiods are: (a) a period with both sunny and cloudy days (25 July–2 August); (b) a period of clear-sky conditions with cold spells and temperatures below the melt threshold value of 1°C (10–15 August); and (c) a period of clear-sky conditions with temperatures always above the threshold value (21–29 August).
Third, we looked at the effect of clear-sky and overcast conditions independently from other meteorological characteristics (such as air temperature) and over the entire season. An algorithm based on incoming short-wave radiation measurements (explained in detail below) was developed to discriminate between sunny and cloudy days. Then the model parameters were recalibrated separately for these two types of meteorological conditions.
-
3. Model transferability in time. The ETI model empirical parameters are recomputed for the 2005 and 2006 ablation seasons at the lowest station on Haut Glacier d’Arolla.
-
4. Model transferability in space to other glaciers. The ETI model is applied to both the 2005 and 2006 ablation seasons on Gornergletscher and to the 2006 ablation season on Tsa de la Tsan glacier, and for each season and glacier the model parameters are optimized.
For all the experiments described above, the analysis procedure is the same as that used by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005): hourly surface melt rates are computed using the ETI model and the results are validated against the simulations of the physically based energy-balance model. The model performance is evaluated for each simulation using the Reference Nash and SutcliffeNash and Sutcliffe (1970) efficiency criterion, R 2, defined as
where M is the melt rate, and the subscripts ‘r’ and ‘s’ denote reference (energy-balance model) and simulated (ETI model) values, respectively.
In all modelling experiments the ETI model is also run with the original parameter values obtained from the original calibration procedure at the Haut Glacier d’Arolla central station in 2001 (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005), and the results of these simulations are compared with those of the recalibration. The calibration procedure is an optimal one, in the sense that all parameter values are allowed in the search for the optimal combination (no trial and error or manual adjustments). We defined a range and step size for each parameter and used parameter ranges based on the values adopted by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005).
The algorithm employed in experiment 2 for discrimination of sunny and cloudy days is a simple criterion based on the number of times, within a day (defined from sunrise to sunset), in which the solar radiation curve changes its slope. On an ideal sunny day, this number is equal to one, whereas it is higher on a cloudy day. In this way, we identified the cloudy days in the season from the record of hourly incoming shortwave radiation. From visual inspection of the record, however, we observed that this simple algorithm cannot select all the cloudy days, and we therefore imposed a threshold of 5 W m−2 below which values of incoming shortwave radiation are not considered, in order to remove erroneous incoming shortwave radiation values in the night hours; and we increased the counter used to discriminate between sunny and cloudy days from one to three (allowing more than one peak in a sunny day). In this way, the presence of a single hour of overcast conditions on a clear-sky day does not mean the the day is classified as overcast.
5. Results
5.1. Analysis of the energy fluxes: Haut Glacier d’Arolla 2001
Table 4 shows the daily energy fluxes averaged over the period of record for the five stations on Haut Glacier d’Arolla in 2001. The shortwave radiation flux is the dominant source of energy for melt at all locations, which agrees with several other energy-balance studies on Alpine glaciers (e.g. Reference OerlemansOerlemans, 2000; Reference Greuell and SmeetsGreuell and Smeets, 2001; Reference Willis, Arnold and BrockWillis and others, 2002). The longwave radiation flux and the turbulent flux of latent heat represent a sink of energy (negative contribution to total melt), whereas the turbulent flux of sensible heat generally contributes positively to the energy balance of the glacier surface (Table 4; Fig. 4). The sign of the contribution of the longwave radiation and turbulent sensible heat fluxes also agrees with most energy-balance studies on Alpine glaciers (e.g. Reference OerlemansOerlemans, 2000), whereas the results for the latent heat are less conclusive, since both positive and negative contributions have been found (see Reference HockHock, 2005, for a review of studies). Negative latent heat fluxes were observed, for example, by Reference Van de Wal, Oerlemans and Van der HageVan de Waal and others (1992), Reference Wagnon, Ribstein, Francou and PouyaudWagnon and others (1999) and Reference PellicciottiPellicciotti and others (2008).
Table 4. Average daily energy fluxes (W m−2) measured or computed over the period of record using the energy-balance model at the locations of the AWSs on Haut Glacier d’Arolla in 2001, and for the separate cases of overcast and clear-sky conditions at the central station: net shortwave radiation, Q I; net longwave radiation, sensible and latent heat fluxes (L, Q H and Q L, respectively). Q M is the resulting energy available for melt (Equation (2)). Also shown is the average daily ablation, M D (mm w.e. d−1). All mean values are computed over the period 30 May–11 September
Fig. 4. Hourly means of the shortwave radiation (a) and sensible heat flux (b) at the five AWSs on Haut Glacier d’Arolla in the 2001 ablation season. The values are hourly means from 0 to 23 h over the entire season, both in W m−2 and converted to the corresponding melt value in mm w.e. h−1.
The shortwave radiation flux exhibits strong diurnal variations, and differences can be observed between the stations (Fig. 4); at the sites in the ablation area (lowest, central and north-central stations) its contribution to total melt is higher than at those in the accumulation area (uppermost and south-central stations). This is the result of different factors, such as surface conditions (albedo is higher in the upper part, so more incoming shortwave radiation is reflected at the surface), shading by the surrounding peaks (reducing the incoming shortwave radiation reaching the ground, which is particularly evident at the south-central station under La Vierge peak) and exposition. The main factor, however, is the higher albedo typical of the accumulation area. Incoming solar radiation peaks are also shifted in time (Fig. 4), reflecting the different aspects and locations of the five stations with respect to the solar path and geometry (Fig. 2). The lowest mean shortwave radiation flux is observed at the south-central station, due to a combination of high albedo and shading.
The longwave radiation flux is very similar at all stations and always negative (Table 4). This is also in good agreement with most energy-balance studies (Reference OerlemansOerlemans, 2000; Reference Greuell and SmeetsGreuell and Smeets, 2001; Reference OhmuraOhmura, 2001).
The turbulent flux of latent heat has little or no diurnal variation, and is always negative and relatively small in comparison with the other fluxes. Its highest negative values are found at the lowest station (Table 4). The turbulent sensible heat flux exhibits a pronounced diurnal cycle, mostly positive and of the same magnitude as the turbulent latent heat flux (Fig. 4). The higher values in the afternoon hours are associated with higher temperature gradients in the boundary layer and the development of a glacier wind, which develops late in the day because of the transport of cold, heavier air from the upper to the lower sections of the glacier due to gravity (e.g. Reference Denby and GreuellDenby and Greuell, 2000; Reference Greuell and SmeetsGreuell and Smeets, 2001; Reference PellicciottiPellicciotti and others, 2008). The presence of glacier wind explains the differences in turbulent fluxes observed between stations, since fluxes are higher at the lowest and north-central station, where the glacier wind is stronger (Fig. 5). The sensible flux is lowest at . the uppermost station, where the lack of an adequate fetch hinders the development of the glacier wind (Fig. 5; Reference Strasser, Corripio, Pellicciotti, Burlando, Brock and FunkStrasser and others, 2004). The second lowest mean value of the sensible heat flux is observed at the south-central station, which is marginal to the glacier flowline and therefore less affected by glacier wind (Fig. 5). Analysis of the relative contribution of the energy fluxes at the five sites is important in the interpretation of the values assumed by the two empirical factors in the analysis of spatial transferability across Haut Glacier d’Arolla, discussed below.
Fig. 5. Distribution of 2 m hourly measurements of wind direction at the proglacial station and at the five AWSs on Haut Glacier d’Arolla, 2001, as a percentage over the entire season (indicated on the diagonal axes). Wind direction indicates the direction where the wind comes from. Down-glacier direction is ∼90°C at the proglacial station, 150°C at the lowest and 120°C at the central, north-central and uppermost stations, while the south-central station is outside the glacier flowline (see Fig. 2).
5.2. Model runs with original parameters from Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005)
The efficiency criteria, R 2, obtained when running the model with the parameter values of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) are listed in Table 5, for all experiments. The efficiency criterion reported in that paper, optimizing the model parameters against the energy-balance outputs at Haut Glacier d’Arolla central station, was R 2 = 0.911. Applying the model at all other sites with these (original) parameters we observe both higher and lower model performances (Table 5). The mean of all experiments is R 2 = 0.925 and the standard deviation is σ = 0.044, with a coefficient of variation of 4.8%. Values of R 2 higher than at the original site are typical of ablation-area sites and of clear-sky conditions, while a reduction in model performance is evident for snow-covered sites in the accumulation area and overcast conditions (Table 5). The maximum R 2 is for the period of ice exposure and clear-sky conditions at Haut Glacier d’Arolla central station (R 2 = 0.986). The minimum efficiency criterion corresponds to overcast conditions at Haut Glacier d’Arolla central station (R 2 = 0.78) and is significantly lower than the average value (Table 5) and lower than the values at the snow-covered sites, where a reduction in model performance is also evident: south-central station at Haut Glacier d’Arolla 2001 (R 2 = 0.895), uppermost station at Haut Glacier d’Arolla 2001 (R 2 = 0.904), snow period at Haut Glacier d’Arolla central station 2001 (R 2 = 0.893) and Tsa de la Tsan glacier 2006 (R 2 = 0.899).
Table 5. Efficiency criterion, R 2, obtained using the original parameters of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) and following their recalibration, for the experiments described in section 4.4. HGdA denotes Haut Glacier d’Arolla
Analysis of the model performance across seasons at Haut Glacier d’Arolla lowest station shows that R 2 is highest for the 2001 season (R 2 = 0.951), the season in which the original parameters were calibrated (for central station), and it is lower for both 2005 and 2006, the decrease being stronger for the 2005 ablation season (with R 2 of 0.925 as compared to 0.935 in 2006; Table 5). For this season, however, only data from the end of July and August were available, a period which was characterized by an initial short very warm spell followed by a cold spell with frequent snowfalls (see Tables 1 and 2). Simulations with the original model parameters underestimate high melt rates, as can be seen in Figure 6, which depicts hourly melt rates simulated by model D and reference melt rates computed by the energy-balance model. This translates to an underestimation of total melt at the end of the season for both years. The decrease in model performance, however, is not large and R 2 stays above 0.92.
Fig. 6. Hourly melt rate simulated by the ETI model (model D) using the empirical parameters calibrated at central station in the 2001 ablation season vs the reference melt rate computed by the energy-balance model. The scatter plots refer to: (a) Haut Glacier d’Arolla lowest station 2001; (b) Haut Glacier d’Arolla lowest station 2005; (c) Haut Glacier d’Arolla lowest station 2006; (d) Gornergletscher 2005; (e) Gornergletscher 2006; and (f) Tsa de la Tsan glacier 2006. The points on the x axis are due to the assumption that melt occurs only if the temperature is above the threshold (1°C), and correspond to those hours in which temperature is below the threshold value for melt onset, where the ETI model computes zero melt.
Analysis of the model transferability to other locations on other glaciers shows that the model performance stays high for the Gornergletscher location in both 2005 and 2006 (Table 5). The model, however, underestimates high melt rates, as for the simulations on Haut Glacier d’Arolla (Fig. 6).
5.3. Model runs with recalibrated parameters
Model performance
Recalibration of model parameters results in higher efficiency criteria than those obtained without recalibration, with a mean R 2 = 0.932, standard deviation σ = 0.037, and therefore a coefficient of variation of 4% (Table 5). A higher mean and smaller coefficient of variation are the result of the recalibration, which leads to an increase in model performance at all sites and thus to a smaller dispersion of the values. This increase, however, is on average not high, with a mean difference of 0.007 (with differences between experiments that are discussed below), and indicates that the model is robust and transferable on the whole, although it performs better in some cases than in others.
To establish whether the difference between the R 2 obtained with the original and recalibrated parameters is statistically significant we performed a paired t test on the two groups of R 2 in Table 5. Results show that there is no statistically significant difference between the mean performance of the model with the parameters of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) and that with recalibrated parameters at the 5% significance level. Moreover, the pattern of efficiency criteria remains the same as in the model runs with the parameter values of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005): very high R 2 values are evident at ablation-area locations (Haut Glacier d’Arolla north-central station in 2001, Haut Glacier d’Arolla lowest station in 2001, ice conditions on Haut Glacier d’Arolla central station in 2001, Haut Glacier d’Arolla lowest station in 2005 and 2006) and for clear-sky conditions, whereas lower values are typical of snow-covered sites (Haut Glacier d’Arolla south-central station in 2001, Haut Glacier d’Arolla uppermost station in 2001, snow conditions on Haut Glacier d’Arolla central station in 2001, Tsa de la Tsan station in 2006), and a distinctly low model performance is evident for overcast days (Table 5). This is the same pattern of results as obtained for the model application with the original parameters. It also indicates that the lower model performance obtained in some experiments, in particular for overcast conditions, cannot be overcome by recalibration, but is an intrinsic characteristic of the model, which works less well for overcast conditions (this point is discussed in more detail below, in relation to the values of the model parameters).
Recalibration gives the highest increase in model performance for those glaciers, locations and periods where R 2 were lowest, in particular for overcast conditions, when R 2 increases from 0.78 to 0.822. Increases in R 2 above the mean increase are observed for the two accumulation sites at Haut Glacier d’Arolla 2001 (south-central and uppermost station), for overcast conditions at Haut Glacier d’Arolla central station 2001, for Haut Glacier d’Arolla lowest station in 2005 and 2006, and for Gornergletscher station in 2006.
Recalibration of the parameters at the four Haut Glacier d’Arolla 2001 locations (experiment 1) leads to an increase in model performance at the two accumulation sites, but not at the locations in the ablation area, where the differences in R 2 between simulations with the optimized parameters and those of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) are minimal. Differences are also very small for all the Haut Glacier d’Arolla 2001 subperiods (experiment 2), except in the case of overcast conditions, discussed above. This is a strong indication of the model’s robustness with respect to the varying meteorological and surface conditions typical of one ablation season. As far as the influence of meteorological conditions is concerned (experiment 2, three subperiods), the model performance (both with recalibrated and original parameters) is highest in the period with clear-sky days and high temperatures (21–29 August, R 2 = 0.988), followed by the period with both sunny and cloudy days (25 July–2 August, R 2 = 0.956) and finally by the period of clear-sky conditions with a cold spell (10–15 August), in which the model performance is the lowest (R 2 = 0.932). This clearly indicates that a problem is related to the days on which temperature is below the threshold or close to it. On clear-sky days with high air temperatures (21–29 August), total melt is largely due to shortwave radiation (SRF is high and TF tends to 0), and the model simulates surface melt with very high accuracy (R 2 = 0.988), nearing the performance of the energy-balance model. The difference in model performance that follows from the recalibration of the parameters for the clear-sky period is very small.
The model performance over the three years investigated at Haut Glacier d’Arolla lowest station (experiment 3) is very high, although the R 2 values obtained with the original parameters were also high. Recalibration leads to the same model performance for the two seasons 2005 and 2006, but the increase in R 2 resulting from recalibration is stronger for the 2005 ablation season (R 2 increases from 0.925 to 0.949), because conditions in 2005 were not representative of the entire ablation season (the 2005 dataset is shorter and includes only the ‘anomalous’ August month; see above and Table 2). On Gornergletscher, a larger increase in model performance is evident for the 2006 season (R 2 increases from 0.916 to 0.93), where simulations with the original parameters underestimated total melt (Fig. 6). In 2005, the efficiency criterion was very high with both original and recalibrated parameters (0.941 and 0.947, respectively), probably due to the higher mean temperatures typical of the 2005 season and more clear-sky days (Table 2), so recalibration does not bring any major improvement over more than 3 months.
Model parameters
The recalibrated model parameters are reported in Table 6 for all experiments (new sites, subperiods, new years and new glaciers), together with the original values of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005). TF values have larger variations around the value found by calibration at the 2001 Haut Glacier d’Arolla central station than SRF values. TF values are rather scattered, with a mean value of 0.055 and a coefficient of variation of 56.3%, while SRFs exhibit lower variability, with a mean of 0.0093 (compared to the optimal value at Haut Glacier d’Arolla central station of 0.0094) and a coefficient of variation of 5.8%.
Table 6. Model parameters TF and SRF obtained by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) (original parameters) and recalibrated in this work for all the experiments discussed in section 4.4
The lowest TF value (0.1) is obtained in three cases: (1) for the ice condition at Haut Glacier d’Arolla 2001 central station; (2) for the 21–29 August period at Haut Glacier d’Arolla 2001 central station and (3) for clear-sky days at Haut Glacier d’Arolla 2001 central station (Table 6). The highest values correspond to overcast conditions at Haut Glacier d’Arolla 2001 central station (TF =0.11) and to Haut Glacier d’Arolla lowest station in 2005 (TF = 0.12) (Table 6). Higher SRFs are typical of clear-sky conditions, whereas the lowest SRF is found for overcast conditions, with an SRF of 0.0079, a difference of 16% from the mean value. The second lowest SRF is at Tsa de la Tsan glacier (SRF of 0.0087), where snowfalls and overcast conditions were more frequent than at all other sites and counted for more than half the recording period (Table 2). A clear pattern is thus evident for both parameters: low TFs and high SRFs correspond to clear-sky conditions, in which incoming shortwave radiation is predominant, while high TFs and low SRFs are typical of overcast days.
High TFs are also observed for ablation-area sites, compared to accumulation-area locations, given the same meteorological forcing (Haut Glacier d’Arolla ablation season 2001): TF is higher than the original at the two ablation sites, lowest and north-central (TF =0.07 at both sites), and lower at the stations in the accumulation area, south-central and uppermost station (TF = 0.03 and 0.04, respectively). Although the parameters vary between 20% and 40% of the mean value, the increase in model performance associated with these values is small, and negligible for lowest and north-central stations (R 2 increases from 0.955 to 0.956 and from 0.940 to 0.941, respectively; Table 5). This pattern of higher TF for ablation-area sites is confirmed by the recalibration of TF values at Haut Glacier d’Arolla lowest station in 2005 and 2006 and at Gornergletscher 2005 and 2006, where all recalibrated values are higher than the original TF at Haut Glacier d’Arolla central station. These four sites are located in the ablation area of the two glaciers, with the number of days of exposed ice >50% (Table 2).
The low TF value obtained for ice conditions at Haut Glacier d’Arolla central station in 2001 is due to the prevalence of clear-sky conditions in the period, and not to the surface characteristics. Although, in this case, recalibration leads to rather different parameter values, these changes do not correspond to significant increases in model performance (Table 6). This points to a possible case of equifinality that we discuss below.
The TF found at the ablation site of Haut Glacier d’Arolla lowest station in 2005 is higher than at all other ablation sites (TF = 0.12), and comparable to that of overcast conditions, despite the relatively high number of clear-sky days (Table 2). This might seem to contradict the finding discussed above. It may be due, however, to the high air temperatures typical of the first half of that season, which result in a higher contribution from the turbulent fluxes and thus of the temperature-dependent term of Equation (1). The first part of the 2005 recording period (and in particular the period 25–30 July) was characterized by very high air temperatures (mean value of 7.5°C; 4°C above the mean for the entire period), and melt at night caused by high turbulent fluxes and a less negative net longwave flux. There followed a period of lower temperatures and snowfalls, which is reflected in the large number of days with T < 0°C (Table 2) which, however, cannot compensate for the earlier very high temperatures, which result in an overall mean temperature over the entire period which is still high compared to the 2001 and 2006 seasons (Table 2). Since the observation period is short compared to the other seasons, the specific meteorological conditions of the first half of the period count more than over a long season, and lead to the rather high TF. We recalibrated the model parameters for the period from 31 July to the end, excluding the very warm spell of 20–30 July, and obtained TF = 0.06 and SRF = 0.0094, in agreement with the values in the other two seasons and at the other ablation sites.
The high values of TF at ablation-area sites are explained by the higher temperatures typical of lower sites, combined with stronger glacier winds that cause higher turbulent fluxes, while lower TF in the accumulation area is associated with lower temperature and lower turbulent fluxes. At both north-central and lowest station in the ablation area, glacier wind was well developed (Fig. 5) and the turbulent fluxes are the highest, while they are lowest at the two accumulation-area sites (Fig. 4; Table 4), due to less glacier wind (Fig. 5). At these two sites, lower temperatures and less sensible heat fluxes together decrease the contribution to surface melt of the temperature-dependent fluxes.
A clear result of this analysis is that remarkably distinct values are assumed by the two empirical factors for overcast and clear-sky conditions. A distinct model performance is also associated with these two types of conditions. In order to explain these values, we looked at the components of the energy balance on clear-sky days (when the model can simulate melt rate with high accuracy, independently of recalibration) and for overcast conditions (when the model performance quickly deteriorates, with both the original and recalibrated parameters). Under these two conditions, different processes take place, reflected in a distinctive partition of the energy fluxes available for melt. Table 4 shows the mean daily value of the four fluxes over the ablation season, for overcast and clear-sky days separately. The major difference between sunny and cloudy days is in the shortwave radiation flux, which is much smaller on overcast days than on clear-sky days, although it remains the main source of energy available for melt. Conversely, contribution to the energy balance by the longwave radiation and the sensible and latent heat fluxes is higher on overcast days, either as higher positive values (sensible fluxes) or as a smaller negative contribution (longwave radiation and latent heat fluxes), thus favouring the melt process in both cases. The recalibrated model parameters reflect this changed contribution of the individual fluxes, as on clear-sky days the shortwave radiation factor, SRF, tends to become higher, nearing the conversion factor from units of W m−2 to mmw.e. h−1 of 0.01078 (Table 6) and reflecting the higher contribution of solar radiation to the total melt. Conversely, the temperature factor, TF, decreases on clear-sky days (from 0.04 to 0.01 for clear-sky conditions at Haut Glacier d’Arolla central station), but it increases for overcast days (from 0.04 to 0.11) as a result of a much higher contribution from the temperature-dependent terms of the energy budget (Table 4). The optimal SRF for cloudy conditions decreases (from 0.0094 to 0.0079) and assumes its lowest value among those for the different subperiods and conditions considered, as a consequence of the much reduced energy flux that is due to incoming solar radiation on overcast days (Table 4).
Variations in the recalibrated model parameters did not result in significant changes in model performance in any of our experiments, as confirmed by the results of the t test. The largest variation in R 2 obtained by recalibration is for overcast conditions at Haut Glacier d’Arolla central station, where recalibration leads to an increase of 4.5% in model performance. The variation in the associated model parameters is 128% in TF and 16% in SRF. For several experiments, such as for clear-sky days, variations in TF as large as 55% lead to a minimal change in R 2, from 0.943 to 0.945 (Table 5). We therefore looked at the model sensitivity to both TF and SRF. Figure 7 shows R 2 computed for each pair of parameters used in the optimization procedure for six of the cases analysed in this work, representative of different parameter behaviours: (a) Haut Glacier d’Arolla central station 2001, where the original parameters were calibrated by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005), and representative of the central tendency; (b) clear-sky days at Haut Glacier d’Arolla central station 2001, characterized by high SRF and low TF, where the model performance is very high; (c) overcast days at Haut Glacier d’Arolla central station 2001, characterized by high TF and low SRF and low model performance; (d) Haut Glacier d’Arolla lowest station 2001, representative of an ablation site, with very high model performance, relatively high TF values and average SRF values; (e) Haut Glacier d’Arolla lowest station 2005 and (f) Haut Glacier d’Arolla lowest station 2006, the same two ablation-area sites but in different years, with equal model performance but rather different factors because of the short and anomalous period in 2005. Similar matrices were obtained at all sites but are not reported here for reasons of space.
Fig. 7. R 2 values corresponding to all possible combinations in the space of the two model parameters for: (a) Haut Glacier d’Arolla central station 2001; (b) clear-sky conditions at Haut Glacier d’Arolla central station 2001; (c) overcast conditions at Haut Glacier d’Arolla central station 2001; (d) Haut Glacier d’Arolla lowest station 2001; (e) Haut Glacier d’Arolla lowest station 2005 and (f) Haut Glacier d’Arolla lowest station 2006.
The surfaces in Figure 7 do not show a clear peak corresponding to a specific pair of parameter values in the space of the two parameters, since R 2 assumes very similar values for a rather broad range of TF and SRF. This points to a problem of equifinality, where several pairs of parameters produce very similar model performances. Independently of the values of the two parameters, which vary as discussed above, the same shape of the surface is obtained for all cases, including the two extreme cases of clear-sky days and overcast conditions (Fig. 7b and c). Increases in TF correspond to a decrease in SRF, and the slope of this relationship is fairly constant across all cases, except for overcast conditions at Haut Glacier d’Arolla central station (Fig. 7c), where the slope of the ridge of equally performing parameter pairs is less steep than in the other cases. Thus, in this case the model is less sensitive to the SRF, probably because the radiation term of Equation (1) counts less under overcast conditions. The model seems otherwise to be equally and little sensitive to changes in SRF and in TF, confirming the results obtained in recalibration and indicating that the model is fairly robust.
5.4. Model E analysis
Table 5 lists the R 2 obtained using the parameterizations of albedo and incoming shortwave radiation and the parameters calibrated at central station during the 2001 ablation season by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) (model E). The model performance decreases for all cases investigated in this work (transferability in space across one glacier, within one season, to other seasons and to other glaciers), as was expected given the high degree of parameterization introduced in this way. The lower input data requirement translates into performance loss (see Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005). The mean R 2 is 0.779, and variability of model performance is higher than for model D, with a coefficient of variation of 13%. However, the same pattern of model performance as observed for model D is evident: the model performs better for clear-sky conditions and at ablation-area sites, and less well for snow conditions and particularly for overcast days (with an R 2 of 0.633 for overcast conditions at Haut Glacier d’Arolla central station). We obtain a significantly smaller R 2 for the period 10–15 August 2001 at Haut Glacier d’Arolla central station (Table 5), but one has to consider that these are only five days, the first of which was cold. The model performance on 10 August is very low, and this effect is felt over the entire short period because of the way the R 2 is computed. Except for this case, the ETI model performance is also quite good when using parameterizations of albedo and incoming shortwave radiation, and comparable with that of the best available model with the same data requirements (see Reference HockHock, 1999; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005).
The decrease in model performance from model D (measured input data) to model E (modelled albedo and incoming shortwave radiation) was also observed by Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005), and a discussion of the inaccuracies introduced by the parameterization of the input variables can be found there. Here we quantify the percentage of errors attributable to the parameterization of albedo and to the solar radiation model, respectively, in order to indicate directions for future research.
Table 7 shows results of simulations with the original model parameters with different degrees of parameterization of the input variables albedo and incoming shortwave radiation. Four model runs were conducted: (1) the standard model D with all measured input data; (2) a run with measured temperature and albedo but parameterized incoming shortwave radiation (indicated as E* in Table 7); (3) a run with measured temperature and incoming shortwave radiation but parameterized albedo (indicated as E** in Table 7); and (4) model E with measured temperature but both albedo and incoming shortwave radiation parameterized. Results show that the albedo parameterization affects the performance of the model only marginally, and its influence on R 2 is less than that of the solar radiation model (Table 7). The albedo model explains between ∼3% and 34% of the difference in performance between models D and E, whereas the radiation model explains between ∼66% and 97% (Table 7). Given the high accuracy of the parametric clear-sky solar radiation model (Reference CorripioCorripio, 2003a; Reference Strasser, Corripio, Pellicciotti, Burlando, Brock and FunkStrasser and others, 2004; Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005), these results suggest that future research should address the modelling of the impact of cloud cover on the radiation receipt, starting by testing and recalibration of the cloud-factor parameterization used in this study (Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005), which was applied as such to other seasons and glaciers.
Table 7. Efficiency criterion, R 2, obtained computing surface melt using the ETI model with the original model parameters (TF = 0.04 and SRF = 0.0094). The performance is calculated using: model D (all input measured data); model E including both the albedo and incoming shortwave radiation parameterizations (and measured temperature); model E* using only the parameterization for incoming shortwave radiation (and measured temperature and albedo); and model E** with only the albedo parameterization (and measured temperature and incoming shortwave radiation). The relative influence of the modelling components is expressed in per cent
6. Concluding Remarks
We have analysed a broad range of conditions across a single glacier, within a single ablation season, over several seasons and at locations on different glaciers in order to test the robustness of the empirical parameters of an ETI model. For this purpose, the same type and quality of input data were used (those provided by AWSs with the same set-up across the three glaciers and several seasons investigated). A thorough validation of the energy-balance simulations used to test the ETI model was also conducted.
The model, on average, works very well when applied to other sites, years and glaciers, and it can be applied without recalibration with only a small drop in model performance in most of the cases considered. The biggest drop in model performance (corresponding to the least model transferability) was obtained for overcast conditions, for which the model performance decays rapidly compared to all other cases (R 2 = 0.78, compared to an average R 2 for all cases considered of 0.925). This result was obtained considering only the cloudy days of the 2001 Haut Glacier d’Arolla season. Overcast conditions, however, are rarely encountered over an entire melt season, so the R 2 of 0.78 may be regarded as a case limit. For clear-sky conditions the model performance (transferability) is extremely high. In general, higher model performance is typical of ablation-area sites and clear-sky conditions, whereas the model performance is slightly lower for accumulation-area sites. Application of the model to other seasons and locations without recalibrating the empirical parameters (i.e. using the original values of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005) has the effect of underestimating melt rates, in particular high melt rates, and therefore the total melt over the entire season. The largest underestimate was obtained for HGdA 2005, equal to ∼14% of the total melt, which is a reasonably low value, and well within the margin of accuracy of ablation measurements. In all other cases underestimation was lower. The 2005 season on Haut Glacier d’Arolla was short and anomalous, with a few very warm days which were followed by a period of frequent snowfalls and cold days, not entirely representative of the average meteorological conditions of a melt season (the efficiency criterion was still >0.92). We therefore conclude that the model is robust, i.e. it performs well for conditions other than those for which it was calibrated, except for the case of overcast conditions, and it is therefore generally transferable with only a small loss in model performance. This also suggests that the stronger physical basis of the ETI model makes its parameters more independent of meteorological variables than the classical temperature-index approach.
Recalibration of model parameters does not result in a great increase in model performance. When recalibrating the two parameters, we obtained an increase in model performance which is not statistically significant for the group of experiments of this study. In the case of overcast conditions, discussed above, for which model performance was also low with the original parameters, recalibration leads to the highest increase in R 2. This, however, stays below the average model performance and remains much lower than in all other cases both with original and recalibrated parameters, indicating that the low model performance is not related to recalibration or the parameter values but is an intrinsic characteristic of the model, which seems not to be particularly suitable for these conditions. The same pattern of model performance is obtained regardless of whether the parameters are recalibrated or the coefficients of Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others (2005) are used.
The values of the optimal model parameters were shown to vary little (Table 6). Larger variations in parameter values were observed when considering different climatic characteristics (clear-sky, overcast), whereas variability within a single glacier was smaller. The recalibrated parameters also exhibit a clear pattern: low TF and high SRF are typical of clear-sky conditions, in which incoming shortwave radiation is dominant, whereas high TF and low SRF correspond to overcast days (Fig. 8).
Fig. 8. Optimal values of TF and SRF in the parameters space for the study sites and seasons examined in this work. TF and SRF are expressed in mmh−1°C−1 and mm h−1 W−1 m2, respectively. The bars denote the range of parameter values for which R 2 is higher than a certain value (equal to R 2 − 0.2% R 2). For the station names corresponding to the numbers see Figure 2.
This work has provided a range of values of model parameters that could be used in the application of the model to glaciers with limited datasets. They do not reflect the conditions typical of one single glacier or meteorological period, since their robustness has been proven over a large range of surface and meteorological conditions. An advantage of this work is that recalibration was conducted over the entire ablation season for almost all cases considered (except where differently stated for the aim of this study), whereas most studies extend only over timescales of days or weeks (Reference OhmuraOhmura, 2001; Reference HockHock, 2005).
Most important, we have related the values of model parameters to specific conditions of the energy balance at the glacier–atmosphere interface, thus providing the physical basis for the interpretation of their values. Higher TF values were found to be associated with higher turbulent fluxes due to the well-developed presence of glacier wind, whereas sites where the glacier wind is less important (i.e. accumulation-area sites where the lack of an adequate fetch prevents the development of glacier wind) exhibited much lower TF. High SRF values, very similar to the physically based value, were obtained in conditions where clear-sky days were prevalent. In this way, we have clearly established a range of parameters that can be used for modelling studies of the melt regime across glaciers where the lack of adequate meteorological or ablation data prevents accurate calibration of the model parameters.
Analysis of model sensitivity to variations in the two empirical factors also confirms the model’s robustness, since its performance does not deteriorate significantly for changes over a wide range of both model parameters. It has also indicated, however, an equifinality problem, in that several pairs of model parameters result in a similar model performance. An analogous result was obtained by Reference HockHock (1999).
We suggest that, given the robustness of model parameters, the ETI model can be used for studies of climate-change impact on glacier melt, when energy-balance models cannot be used because of their complexity or higher data requirements. The ETI model will be more transferable than the standard degree-day approach, which uses only air temperature as an integrated index of the complex processes of energy exchange at the glacier surface. In the standard degree-day approach, the empirical calibration factor is more site-dependent and its physical meaning is more difficult to identify and therefore to relate to the prevailing meteorological conditions. Indeed, values of the degree-day factor have been shown to vary over a broad range and no regional patterns could be identified in the degree-day factors reported in the literature, which are also affected by the different integration periods used in the studies and by the way in which they are computed (Reference HockHock, 2003). Whereas standard or classical temperature-index or degree-day models suffer from limitations related to the need for recalibrating empirical parameters that have very little physical basis (Reference HockHock, 2003), the model presented in this paper uses parameters with a clearer physical meaning.
Most studies of melt regime across glaciers, especially those interested in high temporal- or spatial-resolution simulations, focus on the application of models to one or two ablation seasons (e.g. Reference Kustas, Rango and UijlenhoetKustas and others, 1994; Reference Hamlin, Pietroniro, Prowse, Soulis and KouwenHamlin and others, 1998; Reference HockHock, 1999). Other modelling approaches, including the standard temperature-index method, have been used for simulations of runoff from mountainous catchments over several seasons (Reference Klok, Jasper, Roelofsma, Gurtz and BadouxKlok and others, 2001; Reference Zappa, Pos, Strasser, Warmerdam and GurtzZappa and others, 2003; Reference Schaefli, Hingray, Niggli and MusySchaefli and others, 2005). In these studies, however, calibration involves so many parameters that identifying meaningful values for those related to melt simulation is difficult. Although we have focused on a melt model that is suitable for high-temporal-resolution modelling, we have explored its robustness over several seasons and meteorological conditions, partly bridging the gap between the two types of strategies.
We have also tested the transferability of a version of the model that uses only air temperature as measured input data and incorporates parameterizations of albedo, incoming clear-sky solar radiation and cloud factor. This model, despite the high degree of parameterization introduced (see Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005), has a relatively high performance, comparable to that of the best alternative temperature-dependent model (Reference HockHock, 1999). Furthermore, in these tests none of the above-mentioned parameterizations was recalibrated, and an increase in model performance could be achieved if data were available for optimization of these parameterizations. Computation of the cloud factor is the main cause of uncertainty in the model, and further research could be devoted to estimate the effect of cloud cover on shortwave radiation (see also Reference Pellicciotti, Brock, Strasser, Burlando, Funk and CorripioPellicciotti and others, 2005).
This study has focused on the model transferability over seasons and to different glaciers in the same large climatic setting of the European Alps. The main findings summarized above indicate that the model can be used to investigate the spatial and temporal pattern of melt rate for a range of elevations and meteorological, surface and topographic conditions typical of the Alpine region. A next step in the direction of the research presented here is to test the model transferability to different latitudes and climatic settings, where the energy balance at the glacier–atmosphere interface and, therefore, the factors driving the melt process are markedly different to those in the European Alps. One such area is the dry Andes of central Chile, where parameter recalibration has pointed to a substantial difference in the TF factor (Reference PellicciottiPellicciotti and others, 2008).
Acknowledgements
We thank A. Kretz for joining the fieldwork on Gornergletscher in 2006, E. Cremonese and U. Morra di Cella for collecting the Tsa de la Tsan data in 2006 and R. Dadić for providing the 2005 and 2006 Haut Glacier d’Arolla dataset. Reviewers I. Willis and T. Jóhannesson provided useful comments that improved the original manuscript. Thanks also to an anonymous reviewer.
