Seasonal and interannual variability of melt-season albedo at Haig Glacier, Canadian Rocky Mountains

In situ observations of summer (June through August, or JJA) albedo are presented for the period 2002–2017 from Haig Glacier in the Canadian Rocky Mountains. The observations provide insight into the seasonal evolution and interannual variability of snow and ice albedo, including the effects of summer snowfall, the decay of snow albedo through the melt season, and the potential short-term impacts of regional wildfire activity on glacier-albedo reductions. Mean JJA albedo (± 1σ ) recorded at an automatic weather station in the upper ablation zone of the glacier was αS = 0.55± 0.07 over this period, with no evidence of longterm trends in surface albedo. Each summer the surface conditions at the weather station undergo a transition from a dry, reflective spring snowpack (αS ∼ 0.8) to a wet, homogeneous midsummer snowpack (αS ∼ 0.5) to exposed, impurity-rich glacier ice, with a measured albedo of 0.21± 0.06 over the study period. The ice albedo drops to ∼ 0.12 during years of intense regional wildfire activity such as 2003 and 2017, but it recovers from this in subsequent years. This seasonal albedo decline is well simulated through a parameterization of snow-albedo decay based on cumulative positive degree days (PDDs), but the parameterization does not capture the impact of summer snowfall events, which cause transient increases in albedo and significantly reduce glacier melt. We introduce this effect through a stochastic parameterization of summer precipitation events within a surface energy balance model. The amount of precipitation and the date of snowfall are randomly selected for each model realization based on a predefined number of summer snow events. This stochastic parameterization provides an improved representation of the mean summer albedo and mass balance at Haig Glacier. We also suggest modifications to conventional degree-day melt factors to better capture the effects of seasonal albedo evolution in temperature-index or positive-degree-day melt models on mountain glaciers. Climate, hydrology, or glacier mass balance models that use these methods typically use a binary rather than continuum approach to prescribing melt factors, with one melt factor for snow and one for ice. As alternatives, monthly melt factors effectively capture the seasonal albedo evolution, or melt factors can be estimated as a function of the albedo where these data are available.

microbial and algal activity require nutrients and meltwater, which increase in association with greater deposition and concentration of impurities, lower albedo, and longer melt seasons. As an example of nutrient delivery, a classical 'spring bloom' of pink algae (Chlamydomonas nivalis) on a glacier can be triggered by walking on a clean, supraglacial snowpack 70 with dirty boots.
In addition, mineral dust deposition on glaciers can increase in association with glacier retreat (Oerlemans et al., 2009), due to exposure of fresh sources of material on the glacier margin as well as melt-concentration effects. Impurities on glacier surfaces are also transported and removed by rainfall and meltwater runoff, as both dissolved and suspended sediment. There 75 is great interest to understand and separate these influences on glacier albedo, to document whether albedo is changing in recent decades, and to quantify the potential impact on glacier mass loss (e.g., Oerlemans et al., 2009;Dumont et al., 2014;Mernild et al., 2015). We examine trends in melt-season and ice albedo at Haig Glacier and discuss the processes governing observed intra-and interannual variations.

80
A final aim of our study is to examine ways in which the seasonal albedo evolution of mountain glaciers can be implicitly included in temperature-index melt models. Models of snow and ice melt frequently employ simplified parameterizations of melt in mountain and polar environments, where essential meteorological data are neither readily available nor easily modelled (Hock 2005;Fausto et al., 2009). Temperature-index or positive-degree-day methods are the most common approach, with melt parameterized as a function of temperature (e.g., Braithwaite 1984) or from a combination of temperature and potential 85 direct shortwave radiation (Cazorzi and Dalla Fontana 1996;Hock 1999). Snow and ice melt are calculated using a melt factor which linearly relates the amount of melt to cumulative positive degree days and potentially other influences, such as incoming shortwave radiation (Hock, 1999). Melt factors are generally taken as constants for snow and for ice, with a higher value for ice due to its lower albedo. This binary treatment of the melt factor does not realistically represent the continuous nature of surface albedo or the systematic seasonal evolution of albedo and melt rates on a glacier. We therefore explore 90 parameterizations of the melt factor that better capture the effects of the seasonal albedo evolution on modelled melt.

Haig Glacier study site
Glaciological and meteorological measurements at Haig Glacier in the Canadian Rocky Mountains were initiated in August 2000 and are ongoing. Surface energy and mass balance characteristics of this site are summarized in Marshall (2014). Haig 95 Glacier is a moderate-sized (2.62 km 2 ) outlet of a small icefield that straddles the North American continental divide between British Columbia and Alberta (Figure 1). It flows southeastwards into Alberta and spans an elevation range from about 2520 m at the terminus to 2750 m at the continental divide. Slopes reaching to 2950 m elevation feed into the upper accumulation area. Local geology is dominated by steeply dipping beds of limestone (CaCO 3 ) and dolostone (MgCO 3 ). https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License.

100
The region has mixed continental and maritime influences, being fed by Pacific air masses that transport moisture to the Canadian Rockies (Sinclair and Marshall, 2009). Snow surveys conducted on the glacier each May indicate a mean winter snowpack of 1.35 m water equivalent (w.e.) on the glacier from 2002-2017, with a standard deviation (σ) of 0.24 m w.e. ( Table   1). The latter provides a measure of the interannual variability. Summer (June through August (JJA)) temperature on the glacier over the period of study averaged 5.3°C from 2002-2017. The site has warm, sunny conditions in the summer months, driving 105 average summer melt totals (± 1σ) of 2.60 ± 0.62 m w.e over this period. Annual mass balance at the site has been negative in every year of the study (Marshall, 2014;Pelto et al., 2019), with the glacier losing all of its winter snow in 9 of 16 years. Mean The forefield AWS has been in place continuously, but sensors can fail, the station was blown down once, and on two other occasions the station was buried by snow from the late spring through early summer. Hence there are occasional data gaps, but there is 92% data coverage for the summer period (June through August, JJA) from 2002-2017. Data collection is ongoing at this site. The glacier AWS is more intermittent, due to the more difficult environment. It was maintained year-round from 120 2001 to 2008, but from 2009 to 2015 the station was set up only in the summer months. It is established at the same site each year. Quality-controlled data represent 79% of JJA days from 2002-2015. Where glacier albedo values are presented in this manuscript, they are restricted to the set of available in situ data. For energy balance and melt modelling, missing glacier meteorological data are estimated from the forefield AWS, with transfer functions (i.e., lapse rates and regression equations) based on the observations that are available from both sites (Marshall, 2014). If data are missing from both stations, gap-filling 125 is based on the mean multi-year value for a given day. This enables a complete estimate of the summer energy and mass balance.

Field Measurements
AWS data are stored at 30-minute intervals, based on the average of 10-second measurements of temperature, humidity, wind 130 speed and direction, incoming and outgoing longwave and shortwave radiation, rainfall, barometric pressure, and snow/ice https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License. surface height. The AWS observations are subject to a manual quality control, and any questionable data are removed from the analysis. The instrumentation and data are described in more detail in Marshall (2014).
This manuscript concentrates on the long-term albedo record at the glacier AWS site, = ↑ / ↓ , for reflected and incoming 135 shortwave radiation ↑ and ↓ . We have 30-minute values for this, but our pragmatic interest in this study is the seasonal and interannual evolution of surface albedo and its influence on mass balance, so we consider just mean daily albedo, calculated from the integrated daily sum of incoming and outgoing shortwave radiation, (1) 140 Note that this gives different values than the average of instantaneous albedo measurements, as it weighs the albedo calculation to the middle of the day, when insolation is highest. This accurately reflects the amount of shortwave energy that is available for glacier melt. It also means that we do not consider zenith angle effects on surface albedo in this study, and measurement uncertainty is reduced through daily averaging. The radiation sensor (Kipp and Zonen CM6B) integrates over the spectral 145 range 0.31 to 2.80 µm, with a manufacturer-reported accuracy of within 5% for mean daily measurements (first class rating from the World Meteorological Organization), although calibration studies indicate mean biases of less than 1% for total daily radiation measurements with this instrument (Myers and Wilcox, 2009). Given that our installation is not maintained on a daily basis and is difficult to maintain an ideal horizontal platform, we take the conservative estimate of an uncertainty of 5% for the mean daily radiation. The uncertainty in mean daily albedo is then 7% (e.g., 0.55 ± 0.04 or 0.20 ± 0.01). 150 Other sources of uncertainty in the albedo measurements include deviations from horizontality for measurements of incoming shortwave radiation, multiple reflections from the undulating glacier surface, reflected radiation from valley walls, and potential covering of upward-looking sensors during snow events, among other effects. The quality control measures identify obvious environmental corruption such as times when fresh snow is covering the sensors or excessive station leaning, and data 155 from these days are omitted from the analysis. The valley walls are steep and are typically free of snow from July through September, and modelling of potential reflected radiation from valleys walls indicates that this is negligible at our AWS site.
Additional data are included in this paper from summer 2017, based on repeat centreline surveys of albedo and chemical analysis of supraglacial snow and ice. These data, described in detail in Miller (2018), provide an additional spatial perspective 160 on Haig Glacier albedo as well as insights about the provenance and concentration of impurities on the glacier surface and their association with albedo. Four surveys were conducted in July and August, 2017, at 33 points on an altitude transect approximating the glacier centreline (see Section 3). Albedo was only measured under clear-sky conditions and within three hours of local solar noon, to minimize the effects of diffuse radiation and high zenith angles. We used a portable pyranometer https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License.
(Jaycar QM1582) for these measurements, taking the average of three upward and three downward shortwave radiation 165 measurements at each point. The sensor was held at a height of ~1.1 m above the glacier surface for all measurements. Based on instrumental fluctuation and repeatability, we assign an uncertainty of 10% to individual radiation measurements. This giving an estimated uncertainty of 8% for the point albedo measurements. Surface snow and ice samples were collected at every third point and were bagged for melting and major ion and organic carbon analyses. The impurity concentrations are referenced here but will be analyzed in detail elsewhere; Miller (2018) provides a complete summary of the chemical analyses. 170

Energy Balance Model
We use a distributed surface energy balance model to examine the influence of seasonal and interannual albedo variations on glacier mass balance and summer runoff for the period 2002-2017 (Ebrahimi and Marshall, 2016). We carry out a survey of the winter snowpack each spring, typically in the second week of May, and use this winter mass balance data as an initial 175 condition for the simulation of summer mass balance. The summer melt model is driven by 30-minute data from the glacier AWS; where these data are lacking we drive the melt model with forefield AWS data, mapped onto the glacier through transfer functions that are well-calibrated from the overlapping data records (Marshall, 2014). Four-component radiation measurements, parameterizations of the turbulent fluxes, and a subsurface model of snow/ice temperature and heat conduction (Ebrahimi and Marshall, 2016) are combined to calculate the net surface energy flux, QN : 180 where Q L ↓ , Q L ↑ , Q H , Q E , and Q C represent incoming and outgoing longwave radiation, sensible and latent heat flux, and subsurface conductive energy flux, respectively. All energy fluxes have units of W m -2 and are defined to be positive when 185 they are sources of energy to the surface. Turbulent fluxes of sensible and latent energy are parameterized from a bulk aerodynamic method (Andreas, 2002;Ebrahimi and Marshall, 2016) and surface temperature and conductive heat flux are internally modelled within a subsurface snow model, which includes calculations of meltwater percolation and refreezing (Samimi and Marshall, 2017). When QN is positive and T s 190 = 0°C, the net energy goes to melting, with melt rate This can be directly related to the classical positive degree day method (Braithwaite 1984), where snow or ice melt m over a period of time τ is calculated from

205
where f s/i is the degree-day melt factor for snow or ice. This linearly relates the amount of melt to cumulative positive degree days (PDD) over time τ. The integrand can also be modified to include other influences, such as the potential direct incoming shortwave radiation (Hock, 1999).
Eq. (5) is an empirical alternative to the physically-based approach in Eq. (4). It is sometimes helpful because surface energy 210 fluxes are uncertain in the absence of local AWS data, due to poorly-constrained meteorological input variables. Wind, humidity, cloud cover, and radiation fields are difficult to estimate in remote mountain terrain. Eq. (5) requires only temperature, which can be estimated via downscaling or interpolation of regional station data or climate model output. While appealing, it is recognized that this parameterization is over-simplified with respect to its transferability to other locations or times. For instance, there is no direct way to incorporate influences from meteorological variables other than temperature, 215 and melt-albedo feedbacks are not physically represented where fs and f i are taken as constants.
Most temperature-index models approximate this seasonal evolution to first order by using different melt factors for ice and snow; typical values are f i ~ 6-9 mm of water equivalent melt per degree day (mm w.e. °C −1 d −1 ), while f s ~ 3-5 mm w.e. °C −1 d −1 (Braithwaite, 1995;Jóhannesson, 1997;Hock, 2003;Casal et al., 2004;Shea et al., 2009). There is considerable local, 220 regional, and temporal variability in the parameters chosen for different studies, with values sometimes twice as high as this, particularly for glacier ice (see Hock, 2003). Lefebre et al. (2002) also find a large spatial variation in melt factors, through modelling studies of melt patterns in Greenland. This variability is associated with differences in the energy balance and surface conditions that drive melt, much of which may be due to variations in surface albedo.

225
For regions where melting is the dominant process in glacier ablation (cf. Lett et al., 2019), one relatively simple way to improve on temperature-index models is to permit melt factors to vary in space and time, consistent with spatial and temporal https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License. variations in net energy and glacier albedo (Schreider et al., 1997;Arendt and Sharp, 1999). For melting over time τ, one can combine Eqs. (4) and (5) to derive an expression for the melt factor at any location (x, y): Eq. (6) implicitly includes the seasonal evolution of surface albedo, as an important control on the melt energy, but numerous other meteorological influences are embedded in E m , so there is not a direct relation between f (t) and α S (t).
Because absorbed shortwave radiation is the dominant term driving ablation of mountain glaciers (Greuell and Smeets, 2001;235 Klok and Oerlemans, 2002), including Haig Glacier (Marshall, 2014), one can expect that E m ∝ 1/α S . Moreover, melt energy is proportional to PDD, such that the numerator scales with the denominator in Eq. (6). Hence, it should be possible to develop a simple parameterization which includes the lead-order effects of surface albedo on the melt factor. We use Eq. (6) to calculate the seasonal evolution of the melt factor, f (t), at the Haig Glacier AWS site. A compilation of monthly mean values of f and α S then informs a relation f = a − bα S which can better represent the seasonal and spatial evolution of melt 240 factors, where albedo is known or can be estimated. Table 1 gives summary statistics for the observed Haig Glacier mass balance, net energy, and mean summer albedo at the AWS site from 2002-2015. Positive degree day and melt totals are presented for the complete summer melt season, May 245 through September, and temperatures, energy fluxes, and albedo values are given for the core summer months, June through August (JJA), when more than 80% of the melt occurs. Table 2 reports the mean monthly values at the AWS site for the 14year record.

AWS Albedo Measurements, 2002-2015
The mean JJA surface albedo at the AWS site is 0.55, with a marked decrease through the summer months and a minimum of 250 0.38 in August ( Figure 2a). The AWS site is in the upper ablation zone of the glacier, and seasonal snow gives way to bare glacier ice at some point in late summer. This is attended by a sharp drop in albedo, to the bare-ice value of 0.21 ± 0.06. Mean August albedo values represent an average of aged, wet snow and bare ice, with year-to-year variability associated with the timing of seasonal snow depletion. The average date of seasonal snow depletion at the AWS site is August 3, ranging from July 21 to August 20 over the study period. New snow accumulation at the AWS site ('winter snow') begins in September in 255 most years, accounting for the albedo increase this month (Table 2). Persistent snow began to accumulate at the AWS site between August 30 and September 25 during our study period. On average, there are 25 ± 10 days with glacier ice exposed at the AWS site during the melt season. The site was selected because it is near the equilibrium mass balance point of the glacier: https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License. the elevation at which annual snow accumulation is equal to summer melt, for the glacier to be in balance. The observations of bare ice exposure are consistent with the persistent negative mass balance of the glacier over the period of observations. 260 Intermittent summer snow events also temporarily refresh the glacier surface (e.g., Figure 2b), reducing the number of snowfree days on the glacier surface. The average melt-season albedo evolution at the site in Figure 2a averages out the impact of episodic summer snowfall events that refresh the snow or ice surface (cf. Figure 4d). As seen in Figure 2b, these cause an immediate increase in albedo to a fresh-snow value of ~0.9, followed by a decay back to the albedo of the underlying surface 265 over the course of hours to a few days. Figure 3 provides a more detailed illustration of summer snowfall events over exposed glacier ice. This plot covers the period August 3-28, 2015, during which there were three distinct summer snow events, each of which increased the surface albedo for two to three days. The events can be predicted somewhat from the meteorological conditions, where temperature drops below 0°C and relative humidity reaches 100% (Figures 3a,b), but they are most clearly evident in the albedo record ( Figure 3c). The accumulation of new snow is also apparent in the glacier surface height (SR50) 270 data ( Figure 3d), attended by an interruption in surface melting.
Even a modest amount of fresh snow has a strong albedo impact; there were roughly 4 cm of accumulation in the first two snow events and 8 cm for the third event in Figure 3. The latter event had a longer impact, roughly three days before the surface albedo returned to values typical of bare ice. Total surface ablation at the AWS site was 1.05 m over this 25-day period ( Figure  275 3d), equivalent to about 0.95 m w.e. Based on the observed bare-ice vs. actual average albedo values over the 25-day period, 0.15 vs. 0.27, we calculate that the snow events reduced the average net energy by 24 W m −2 , equivalent to 0.16 m w.e. or a 17% reduction in melting over this period. Hence, the direct impact of summer snowfall on glacier mass balance is generally minor (estimated at ~0.03 m w.e. for the events in Figure 3), but the indirect impact through increased albedo and reduced melting is important. 280 Based on analysis of the SR50 and albedo data over the full study period, an average of 9.3 ± 2.6 ephemeral snowfall events per year occurred at the AWS site from May to September. This included 6.3 ± 2.2 summer (JJA) events. Our main criteria to identify summer snowfall events is a mean daily albedo jump of at least 0.15. This may not capture trace precipitation events that are too minor to be seen in either of the SR50 or albedo measurements (i.e., too ephemeral or not enough snow to mask 285 the underlying surface).

Relation Between Summer Albedo and Mass Balance
The broader relations between glacier mass balance and albedo, summer snow events, and temperature at the Haig Glacier AWS site are summarized in Table 3 Table 3 are correlated, with numerous interactions, but the importance of albedo is clear. Monthly mean albedo is highly correlated with monthly melt (r = −0.89) and net energy (r = −0.84), in addition to strong negative correlations with other melt indicators such as mean monthly temperature and PDD (r = −0.74). Monthly albedo values are also significantly correlated with the optimal monthly degree-day melt factor calculated from Eq. (6) (r = −0.66), which we discuss further in Section 4.2. 295 Correlation coefficients for the mean summer conditions are generally weaker, though there are fewer samples so these statistics are less robust. Mean melt-season albedo remains strongly correlated with summer mass balance (melt), net energy, and temperature, but is only weakly associated with total melt-season PDD. The influence of winter mass balance is also evident through positive correlations with albedo and net mass balance; deeper snowpacks take longer to melt out, delaying 300 the transition to the low-albedo summer surface. Mean summer albedo is also positively correlated with the number of summer snow events (r = 0.66). Due to these compounding influences, mean summer albedo is highly correlated with net annual mass balance (r = 0.76). Summer snow events have a significant overall influence on the summer and net mass balance (r = −0.73 and r = 0.70, respectively).

Ice Albedo Values 305
The progressive decline in glacier surface albedo through the melt season has been reported in many previous studies (e.g., Brock et al., 2000, Klok and Oerlemans, 2002, 2004. Within the seasonal snow, this can generally be related to cumulative melting, with its associated effects on snow grain size, liquid water content, and increasing concentration of impurities (Warren and Wiscombe, 1980). We discuss modelling of this seasonal evolution in Section 5. At Haig Glacier, the ripened and saturated July snowpack typically asymptotes at an albedo values of about 0.5, before surface albedo drops sharply to a value of ~0.2 310 once bare ice is exposed. Figure 4a visually captures this transition. The mean glacier albedo value from a composite of 224 bare-ice days in summer (JJA) is 0.21 ± 0.06. Including the month of September, the number of bare-ice days increases to 272 and the mean ice albedo is 0.22 ± 0.07.
Our values for glacier ice albedo are in line with other mid-latitude glacier observations (e.g., Brock et al., 2000;Gerbaux et 315 al., 2005;Naegeli et al., 2019) and the value of 0.2 recommended by Cuffey and Paterson (2010) for impurity-rich ice.
Particulate concentrations are high in the old snow and glacier ice on Haig Glacier, and include a combination of mineral dust, black carbon, and organic material (see Section 4.4). Ice albedo values of 0.07 have been measured on the lower glacier in multiple years, in association with high impurity loads (Figure 4b,c). Indeed, during spatial albedo surveys, measured albedo is generally higher on the proglacial limestone than in the lower ablation zone (e.g., Figure 1c and Figure 4b). No part of Haig 320 Glacier is considered to be debris-covered, where material covering the glacier is thick enough to insulate the ice surface from ablation; rather, supraglacial particulate matter takes the form of a thin film, with considerable spatial heterogeneity and https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License. temporal variability in impurity concentrations. The heterogeneity is presumably associated with variable patterns of atmospheric deposition, flushing (cleansing of the glacier surface through rain events or meltwater runoff), and microbial/algal activity. Temporal changes in these processes may underlie the range of ice albedo values measured at the AWS, from 0.11 to 325 0.34 ( Figure 5).
We sorted all bare-ice days into subsets of clear-sky and overcast conditions, based on incoming shortwave radiation measurements at the AWS (specifically, the ratio of the total daily and potential direct incoming solar radiation). The spectral reflectance of snow is dependent on the solar incidence angle (Hubley, 1955;Wiscombe and Warren, 1980), hence differs for 330 direct vs. diffuse radiation. As a result, mean daily albedo values can be expected to be higher on cloudy days, when diffuse radiation is dominant (Cutler and Munro, 1996;Brock, 2004;Abermann et al., 2014). However, we found no difference between the mean ice albedo values for clear-sky and overcast days in our dataset (a mean value of 0.21 for each subset).
Glacier ice albedo is less sensitive to the zenith angle than snow (Cutler and Munro, 1996), and this appears to be particularly true for impurity-rich glacier ice, possibly due to isotropic absorption by impurities and liquid water on the glacier surface and 335 in the intergranular interstices.
However, another phenomenon may contribute to some of the higher ice-albedo values recorded at our AWS site. Values above 0.3 are most common in September, after ephemeral snowfall events have melted away; these serve to increase albedo by 0.1 to 0.2 above the minimum seasonal values attained in August. We interpret this to be due to either residual, refrozen 340 (i.e. superimposed) ice that is more reflective or because of meltwater flushing of some of the local impurities. The albedo record in Figure 3c illustrates this temporary increase in albedo in the days following fresh snowmelt, particularly for the third snow event on August 22. The storm snow was melted away by August 24 (Figure 3d), but the exposed ice albedo values remained above their early-August 'baseline' value of ~0.13 until August 28. Albedo values after this returned to the baseline, indicating that a potential crust of superimposed or flushed ice had been melted away. This pattern is typical of the albedo 345 evolution following summer snow events. There is no trend of ice albedo decrease through the summer melt season; hence, no evidence of increasing impurities that cause progressive darkening once the glacier surface ice is exposed. In contrast, ice albedo increases in late 350 August and September, perhaps associated with the brightening influence of superimposed ice or meltwater flushing, as hypothesized above. The summer 2003 was an interesting exception, plotted in red in Figure 6a. Ice albedo declined through July and the first week of August in 2003, reaching a minimum of 0.11 (the lowest mean daily value on record at the Haig AWS) and remaining at ~0.13 until strong melting caused the AWS to lean beyond a condition that ensures reliable data after August 22. Severe wildfire conditions in southwestern Canada that summer resulted in a an evacuation order for the region in mid-August, so we were forced to leave the site and could not maintain the AWS. These same wildfire conditions may have resulted in deposition of soot and black carbon that produced the extremely low albedo values that summer.

Albedo Transects and Snow/Ice Impurity Data 365
To supplement the AWS albedo record from a single point on the glacier, we conducted spatial albedo surveys across the glacier during different seasons. In summer 2017 we completed four centrelines albedo surveys through July and August, in conjunction with collection of snow and ice samples to analyze the chemistry and concentration of particulate matter on the glacier surface. Figure 7 plots the location of the survey sites and the centreline albedo data from these four surveys, with summary data provided in Table 4. 370 The characteristic decrease in surface albedo on Haig Glacier over the summer melt season is evident in Figure 7b. For the initial survey on July 13, the glacier surface was still completely snow-covered, with a relatively uniform albedo typical of old, wet snow. The average albedo value (± 1 standard deviation) on the July 13 survey was 0.48 ± 0.04, and albedo declined through each ~two-week period (Table 4). Glacier ice was exposed as the seasonal snowline moved upglacier in the following 375 weeks, with albedo values dropping to 0.16 ± 0.11 on August 22. The toe of the glacier is a high-accumulation area due to wind-blown snow deposition in the lee of a convexity (Adhikari and Marshall, 2013); it retained seasonal snow through mid-August, but was snow-free by August 22. On this final day of sampling, only the uppermost sampling site retained seasonal snow cover, possibly refreshed by a snow event on August 13-14.

380
The seasonal albedo decline is partly associated with the transition from snow to bare ice and partly because the glacier ice albedo systematically decreased over the course of the melt season, from an average value of 0.21 ± 0.07 on July 25 to 0.13 ± 0.05 on August 22 (Table 4). The snow albedo at the glacier toe also decreased from July 13 to August 9 (Figure 7b), but a decline in snow albedo was not apparent on the upper glacier. Much of the glacier surface had an albedo of ~0.1 by the end of the melt season, and the lower glacier had values of 0.07. As reported in previous studies (Brock et al., 2000;Klok and 385 Oerlemans, 2002), ice albedo generally increases with altitude on Haig Glacier. Impurity measurements from snow and ice samples collected during each glacier visit indicate a ~four-fold increase in total carbon on the glacier surface from July 26 to August 9, from average concentrations of 5.6 to 22.7 mg/L, respectively (Miller, 2018). These data will be described in detail elsewhere (Miller and Marshall, in preparation). The increase in impurities was 400 evident in both inorganic and organic carbon, which are present in roughly equal proportions, and are consistent with the potential impacts of wildfire fallout on surface albedo. Particulate matter from local terrigenous dust also increased over this period, but by a factor of about two. The mineral dust load is dominated by calcium and magnesium carbonate, with an average summer carbon concentration associated with carbonaceous dust, [Cdust], of 2.2 mg/L. This compares with an average total carbon concentration of 14.5 mg/L, indicating that up to 85% of the carbon on the glacier has a source other than local, 405 terrigenous dust, although we recognize that Ca and Mg are highly soluble and may have been preferentially removed by meltwater. We are not able to partition the non-dust carbon between algal, wildfire, or other potential sources such as British Columbia industrial activity.

Albedo Modelling
Given the strong variation of surface albedo through the summer melt season -a typical decline from ~0.9 to ~0.2 from May to August -it is important to capture the seasonal albedo evolution in glacier hydrological and mass balance models. Numerous other researchers have tackled this problem, and physically-based snow albedo models have been proposed (e.g., Marshall and Oglesby, 1994;Flanner and Zender, 2006;Gardner and Sharp, 2010;Aoki et al., 2011). 415 In glacier modelling, the decrease in supraglacial snow albedo through the melt season can be approximated by a proxy for snow age or cumulative melting, to capture the systematic decline in albedo due to rounding and growth of snow grains, the effects of liquid water content in the snowpack, and increasing concentration of impurities (Brock et al., 2000). Following https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License. Hirose and Marshall (2013), a parameterization based on cumulative PDD gives a reasonable fit to observed albedo at Haig 420 Glacier, where 0 is the fresh-snow albedo (~0.85) and k is an albedo decay coefficient. The free parameter, k, can be tuned to fit a 425 given observational record such as those plotted in Figure 2b. The influence of snow grain size or impurity concentration could also be incorporated in k in this type of model. Once the seasonal snow has melted, the albedo drops to that of firn or ice.
Additional details can be added to the snow albedo model for thin snowpacks (e.g., less than ~10 cm), to capture the influence of the underlying firn or ice albedo when it begins to show through (Oerlemans and Knap, 1998).

430
This simple parameterization works reasonably well, with minimal inputs, but fails to capture the albedo impact of fresh snow events, which temporarily brighten the glacier surface. The albedo impact of summer snowfall is ephemeral, but these events significantly increase the mean summer albedo and reduce the total summer melt, as discussed in Section 4.2. This is a difficult thing to model remotely or in future projections, as precipitation events can be extremely local in the mountains and the phase of precipitation (rain vs. snow) is difficult to predict. Rain and snow events are both are common on Haig Glacier in the 435 summer months, often mixed on the glacier as a function of elevation.
As a simple approach to address this, we recommend the introduction of summer snow events as a stochastic process. A normal distribution characterized by the mean and standard deviation of the expected number of summer precipitation events (Table   1) can be randomly sampled, with the phase of precipitation determined by the current local temperature (e.g., T < 2°C for 440 snowfall). Total event precipitation can also be treated as a random variable. Each summer snow event then resets the surface albedo to the fresh-snow value, 0 , and the albedo decay begins anew with this fresh summer snow, until it is ablated and the underlying, darker surface re-emerges (old snow, firn, or ice). The albedo of the underlying glacier firn or ice can be held constant or can be parameterized to decay with the amount of time exposed or as a function of impurity concentration. Lacking a good understanding or independent model of these processes, we assign the ice albedo to be equal to the observed longterm 445 mean at the Haig Glacier AWS site, 0.21, once again including some stochastic variability based on the observed distribution of ice albedo values (Table 1). in a glacier energy balance model (Ebrahimi and Marshall, 2016) that calculates PDD and melting, forced by the observed 450 AWS data. The model is seeded with the observed winter snowpack, measured on April 13 of that year, and is run from May different model realizations, illustrating the differences in timing of summer snow events. These are not completely random, as a stochastic precipitation event will only register as a snowfall when temperatures are cold enough. Because of this, snow 455 events are more common in early and late summer and also correlated in different model realizations. This temperature control also helps the model to capture the end of the summer melt season (beginning of winter snow accumulation), although not always. For instance, the end of summer, which occurs around September 17 this year, is well captured in Figure 8a but is ~one week late in Figure 8b.

Implications for Glacier Mass Balance 460
The seasonal albedo decline, summer snow events, and realistic values of firn and ice albedo are all important to resolve in models of glacier energy and mass balance. As an example, the energy balance model driven by the observed AWS meteorological data and surface albedo gives a total melt of 2.97 m w.e. for the 2007 melt season (May to September), the case study in Figure 8, corresponding to a mean surface albedo of 0.59. With the random snow events as illustrated in Figure   8, the mean of five realizations gives an estimated summer melt of 2.82 ± 0.25 m w.e., corresponding to a mean melt-season 465 albedo of 0.60. This slightly underestimates but is within the uncertainty of the observation-based estimate; the energy balance model is well-calibrated to this site. Without the summer snow events, however, modelled melt equals 3.61 m w.e., with a mean melt-season albedo of 0.48. This overestimates summer ablation and runoff by 22%.
There are also important implications for modelling of glacier mass balance at remote sites or in future projections. As 470 temperatures warm, summer precipitation events can be expected to shift to rain rather than snowfall, with positive feedbacks on glacier melting. Without an explicit treatment of summer snow events and their impact on albedo, models calibrated to present-day conditions will not capture this feedback. Similar caveats can be raised about assignment of a constant, observationally-based ice albedo in mass balance models; conditions vary between glaciers and may change in the future as a function of changing particulate loads, and possibly other factors. Physically-based models of impurity deposition and washout 475 and the relation to ice albedo are needed for more reliable regional models and future projections. Similar efforts are underway to improve mass balance models for debris-covered glaciers (e.g., Reid and Brock, 2010;Rowan et al., 2015), although the processes differ for transport and dispersal of coarse debris such as rockfall.
Most studies to date involving regional or future models of glacier mass balance employ degree-day or temperature-index melt 480 models, since in situ meteorological data are unavailable (e.g., Marzeion et al., 2014;Clarke et al., 2015). A full surface energy balance requires a large suite of variables such as wind speed, humidity, and cloud conditions. These are difficult to downscale from climate models with fidelity, relative to temperature fields. In this case, albedo does not appear directly in most formulations of the melt parameterization, but is implicit in the assignment of different degree-day melt factors for snow and ice, fsnow and f ice .
Observations of surface albedo evolution on mountain glaciers make it clear that a continuum approach is more appropriate, with changing degree-day melt factors that track the seasonal albedo evolution (e.g., Arendt and Sharp, 1999). The monthly melt factor f was calculated from Eq. (6) using mean monthly values of melt energy and temperature at the Haig Glacier AWS from 2002-2015, giving a range of values that can be compared directly with mean monthly albedo ( Figure 9). As expected, 490 there is a relatively strong inverse relation, with a linear correlation coefficient of −0.66 and the regression equation f = 7.98 -6.16 α s . This relation could be applied if one has independent estimates of the surface albedo, e.g., from remote sensing or UAV surveys. Alternatively, Table 2 includes mean monthly values of the melt factor for the full dataset, which would be preferable to using single values for snow and for ice. Equivalent monthly factors could also be calculated for the radiation melt coefficient in enhanced temperature-index melt models which use potential direct solar radiation as an input (e.g., Hock, 495 1999;Clarke et al., 2015;Carenzo et al., 2016). Table 3 summarizes the broader relations between albedo, temperature, and mass balance conditions on Haig Glacier. The monthly net energy and melt are strongly correlated with temperature and PDD (r ~ 0.9), implying that temperature-index melt models could give good estimates of monthly melt at this site, given a judicious choice of melt factor, 500

The correlation matrix in
f. This relationship is weaker but still significant for the total summer melt (r ~ 0.7). Annual mass balance is highly correlated the summer balance (r = 0.94), emphasizing the importance of the summer melt season, which in turn is highly sensitive to surface albedo. Mean summer albedo is highly correlated with net annual mass balance (r = 0.76). This is stronger than the correlation of net balance with mean summer temperature or PDD totals. Closely related to this, summer snow events have a significant association with the summer and net mass balance (r = −0.73 and r = 0.70, respectively). The amount of mass added 505 to the glacier is small in summer relative to the winter snowpack (less than 5%), but net mass balance is more strongly correlated with the number of summer snow events than the winter balance, due to the albedo impact.

Temporal Variability and Trends in Ice Albedo 510
Glacier ice albedo is low at this site, in association with high concentrations of supraglacial impurities. The impurities are a combination of mineral dust, primarily calcium and magnesium carbonate, other sources of inorganic carbon, and organic carbon, including active algal populations. The mean value of summer ice albedo at the AWS site is 0.21, but this dips to 0.11 in some years (2003,2015,2017), possibly in association with regional wildfire activity. The summers of 2003 and 2017 were particularly active wildfire seasons in southern British Columbia, upwind of our study site, and ice albedo declined through 515 the summer melt season in these two summers. This was unusual, however; overall, there is no evidence of decreases in ice albedo over the melt season (Figure 6a). In contrast, bare-ice albedo increases slightly in the late summer and early autumn, perhaps in association with superimposed ice formation and/or meltwater rinsing of the glacier surface after transient summer snow events.

520
Similarly, there is no multi-year trend for ice albedo at this site (Figure 6b), although this record is limited to a relatively short period (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015) at just one location on the glacier. Based on the available data, however, there is no evidence of glacier 'darkening' over the period of study, despite years such as 2003 which experienced heavy deposition and accumulation of particulate matter. This stands in contrast to reported glacier albedo reductions over the last two decades in other regions (e.g., Mernild et al., 2015;Naegeli et al., 2019). These results imply that there is some degree of effective cleansing and refreshing 525 of the glacier surface through rainfall and meltwater runoff, although the baseline albedo remains low, in connection with a multiyear supraglacial accumulation of impurities.
The transect data from 2017 indicate lower albedo values and greater impurity loads near the glacier terminus (Figure 7b; Miller and Marshall, in preparation). This is consistent with increased concentration of residual particulate matter due to 530 cumulative melting, within a given summer or over many years, as well as the possibility of greater mineral dust loading on the lower glacier (and associated nutrients to support algal activity), as reported by Oerlemans et al. (2009). We do not have the data to assess multiyear albedo or impurity trends on the lower glacier, to test whether the terminus zone is darkening as a feedback to negative mass balance trends, as has been reported elsewhere (Oerlemans et al., 2009;Naegeli et al., 2019). which should lead to leaching and removal of some dissolved ions in the meltwater runoff, but concentrations increased 2-to 540 6-fold. Subsurface snow samples were essentially clean (concentrations below the detection limit of 0.1 mg/L), so the increased particulate matter must have been due to deposition on the glacier. Further work is needed to quantify deposition and rinsing (leaching, washout) of particulates, to develop models for these processes, and to characterize their influence on temporal and spatial variations in albedo.

Conclusions 545
Albedo measurements from the upper ablation area of Haig Glacier over the period 2002-2017 indicate significant interannual variability in mean melt-season and glacier ice albedo, but there is no temporal trend in surface albedo over this period. This runs counter to documented albedo reductions elsewhere (Oerlemans et al., 2009;Mernild et al., 2015;Williamson et al., 2019;di Mauro et al., 2020), and to anecdotal evidence of darkening glaciers in the Canadian Rockies. The result may just be specific to the Haig Glacier AWS site; we do not have data to constrain albedo trends on the lower glacier, where changes and melt 550 feedbacks have been strongest over the observation period. Moreover, the record is short for trend detection and Haig Glacier https://doi.org/10.5194/tc-2020-87 Preprint. Discussion started: 17 April 2020 c Author(s) 2020. CC BY 4.0 License. mass balance has been negative through this whole period. It is plausible that there have not been significant changes in albedo in the 2000s, but there may be substantial darkening over a multi-decadal time frame (e.g., since the 1970s).
The baseline summer ice albedo at the Haig Glacier AWS site is 0.21 ± 0.06, so it has been relatively low for the whole period 555 of study. It drops to values as low as 0.11 in certain years though (2003,2017), in association with strong wildfire seasons in southern British Columbia, upwind of our study site. Large increases in total and organic carbon concentrations measured on the glacier in August 2017 support this association. The glacier ice appears to recover from these low-albedo summers, however, returning to albedo values of ~0.2 in subsequent years. This is evidence of cleansing of the glacier surface by rainfall and meltwater runoff. Impurities at Haig Glacier are dominated by fine particulate matter, mineral dust in particular and much 560 of this may be effectively leached as dissolved sediment load. The mass balance of supraglacial particulate matter is not well understood, and requires further study.
Other processes controlling the variability in ice albedo also require further study. We find no relation between mean daily ice albedo and cloud conditions in our data, as reported elsewhere (e.g., Brock, 2004;Abermann et al., 2014) and as theoretically 565 expected. This may be because the albedo is relatively low, with a high concentration of impurities; particulate matter and liquid water content act as isotropic absorbers, reducing the sensitivity of specular reflection to zenith angle (hence, diffuse vs. direct radiation). We also see no evidence of ice albedo reductions through the melt season, unlike in seasonal snow, although the major wildfire years provide an exception to this. This further argues for effective rinsing of the glacier surface in most summers, save when dry deposition of particulate matter is unusually high. 570 We do see evidence of temporary increases in bare-ice albedo to values of ~0.3 following melting and runoff of fresh summer snow. The post-snowfall glacier ice albedo is commonly about 0.15 higher than before the snow event. This may be due to a reflective, superimposed ice crust that temporarily forms after snow events, or it could be a result of effective washing of the glacier surface from the melting of clean snow. The increase in ice albedo is transient, but the effect persists for two or more 575 days after the new snow has melted away. This effect is subtle, but overall, summer snow events at Haig Glacier have a large impact on mean summer albedo and glacier mass balance. An average of 9.3 ± 2.6 such events were recorded each summer, resulting in a mean melt-season albedo increase of about 0.1 (e.g., from 0.48 to 0.60 in 2007). Such events are particularly significant when they occur late in the summer, 580 temporarily brightening the low-albedo ice. Based on both energy balance modelling and direct AWS observations, we estimate that summer snow events reduce summer melting and runoff by about 20% at Haig Glacier. This is an important potential feedback and sensitivity to climate change, as warming is likely to cause more of these summer precipitation events to shift to rainfall rather than snow in the coming decades.
A stochastic model of summer snow events helps to capture the typical melt-season albedo evolution at Haig Glacier. This is necessary for realistic mass balance modelling at this site, and the approach could be adapted for use elsewhere, with some knowledge of precipitation frequency during the melt season. The seasonal albedo evolution on glaciers governs how effectively incoming shortwave radiation is converted to melt, and it is important to capture this influence in simplified melt and mass balance modelling applications where local meteorological data are not available. We suggest ways in which the 590 melt factor in temperature-index melt models can be parameterized as a function of albedo, to better capture the conversion of positive degree days to melt.
Haig Glacier albedo values and summer snow conditions may not be broadly applicable, particularly given regional differences in the provenance and concentration of impurities. Particulate loading is highly variable in space, even within a given glacier. 595 The processes discussed in this contribution and the general pattern of melt-season albedo evolution are relevant to most mountain glaciers, however. These observations can help to inform regional models of glacier mass balance and assessments of glacier response to climate change. We emphasize the need for process studies of particulate mass balance (deposition, accumulation, transport, and removal) in supraglacial environments, including the potential effects of forest fire fallout on glacier albedo and mass balance. 600

Acknowledgements
We are grateful to the Natural Science and Engineering Research Council of Canada (NSERC) and the Canada Research Chairs program for sustained, long-term support of the Haig Glacier project. Rick Smith of the University of Calgary Weather Research Station has been instrumental in helping to maintain and calibrate our sensors. Numerous graduate and undergraduate 605 students assisted with the Haig Glacier fieldwork since 2000, and we particularly thank Patrick Coulas for his assistance with the summer 2017 albedo and supraglacial snow/ice sampling.

Code/Data Availability
The automatic weather station data from Haig Glacier and MATLAB code for the surface energy balance model used in this study are available from the authors on request. 610

Author Contributions
SM initiated the Haig Glacier field study, led the field effort and data collection, wrote the MATLAB code for the surface energy balance modelling, was responsible for the data analysis and wrote the manuscript. KM collected the summer 2017 surface albedo and supraglacial chemistry data as part of her Masters research at the University of Calgary. The authors declare no competing interests and no conflict of interest with this research or its conclusions. Figures   Table 1. Mean summer albedo and mass balance conditions at Haig Glacier, 2002-2017, based on glacier-wide simulations driven by local AWS data. Bw, Bs, and Bn are the winter, summer and annual specific mass balance, Nm and Nss are the number of melt days and summer 760 snowfall days from May through September, PDD are the positive degree days over the summer melt season, T and QN are the mean JJA air temperature and net energy flux, and Em is the total summer melt energy. Albedo values are as measured at the glacier AWS, where α S is the mean JJA surface albedo and αi is the measured ice albedo for the composite of snow-free days (N = 224).   (17) August 9 0.23 ± 0.11 (33) 0.41 ± 0.03 (9) 0.17 ± 0.05 (24) August 22 0.16 ± 0.11 (33) 0.47 ± 0.08 (3) 0.13 ± 0.05 (30)