Estimating ice albedo from fine debris cover quantified by a semi-automatic method: the case study of Forni Glacier, Italian Alps

In spite of the quite abundant literature focusing on fine debris deposition over glacier accumulation areas, less attention has been paid to the glacier melting surface. Accordingly, we proposed a novel method based on semiautomatic image analysis to estimate ice albedo from fine debris coverage (d). Our procedure was tested on the surface of a wide Alpine valley glacier (the Forni Glacier, Italy), in summer 2011, 2012 and 2013, acquiring parallel data sets of in situ measurements of ice albedo and high-resolution surface images. Analysis of 51 images yielded d values ranging from 0.01 to 0.63 and albedo was found to vary from 0.06 to 0.32. The estimated d values are in a linear relation with the natural logarithm of measured ice albedo (R =−0.84). The robustness of our approach in evaluating d was analyzed through five sensitivity tests, and we found that it is largely replicable. On the Forni Glacier, we also quantified a mean debris coverage rate (Cr) equal to 6 g m −2 per day during the ablation season of 2013, thus supporting previous studies that describe ongoing darkening phenomena at Alpine debris-free glaciers surface. In addition to debris coverage, we also considered the impact of water (both from melt and rainfall) as a factor that tunes albedo: meltwater occurs during the central hours of the day, decreasing the albedo due to its lower reflectivity; instead, rainfall causes a subsequent mean daily albedo increase slightly higher than 20 %, although it is short-lasting (from 1 to 4 days).


Research motivation and study aims
An understanding of how albedo varies in response to changes in the state of the surface is a crucial component in modeling ice melt and in describing the climate of the ice-covered regions and the climate in general (see Grenfell, 2011).Moreover, in recent climate modeling studies, attention is paid to the ice-albedo feedback and to its action in modulating the changes in the total energy balance of the analyzed area (Grenfell, 2011).One of the most important factors driving albedo changes is the occurrence of debris at the glacier surface, as it influences the features and evolution of glaciers and glacierized areas in numerous ways (Bolch, 2011).Recently, dust and black carbon deposition on glacier accumulation areas (i.e., at the surface of snow and firn) have been of increasing interest to the scientific community due to accelerated snow melting rates affecting glaciers in the highelevation glacierized areas of Asia (Flanner et al., 2009;Yasunari et al., 2010).In addition, Dumont et al. (2014) found that the Greenland springtime darkening since 2009 stems from a widespread increase in the amount of light-absorbing impurities in snow, as well as in the atmosphere.Clarke and Noone (1985) found that the black carbon deposition caused an Arctic snow albedo reduction of 1-3 % in fresh snow and by an additional factor of 3 as the snow ages.Hansen and Nazarenko (2004) modeled this decreased albedo in Arctic snow and sea ice and found this resulted in a hemispheric radiative forcing of +0.3 W m −2 , which may have substantially impacted the Northern Hemisphere climate in recent decades.
Published by Copernicus Publications on behalf of the European Geosciences Union.
In spite of this abundant literature, the effects of fine (mainly dust) debris cover at the melting surface of debrisfree (mountain) glaciers are still poorly debated and sometimes underestimated.
In this contribution, we quantified fine debris coverage at the melting surface of an Alpine debris-free glacier in order to evaluate its seasonal variability and its influence on ice albedo.In particular, to permit comparisons between different glacier zones and among different glaciers, we developed a protocol to standardize field and laboratory analyses.Moreover, we assessed the influence of water (originating from ice melt and liquid precipitation) on the ice albedo variability.Finally, we analyzed the short-term evolution of fine debris occurring at the glacier surface by describing its sedimentological properties and the debris coverage rate during the ice melting season.

Previous studies and recent literature on fine debris occurring at the glacier surface
One of the most important factors driving glacier albedo (apart from, for instance, the meteorological conditions and light scattering by bubbles and cracks) is light absorption by fine debris and dust (Brock et al., 2000;Brock, 2004;Klok et al., 2003).Fine debris and dust vary across the glacier surface both in space and in time, and consist of mineral and organic fractions with a mean diameter lower than 2 mm.Both components may be autochthonous or allochthonous from both englacial origin and wind transport.The organic elements can originate from bacterial decomposition of organic matter (in situ or outside the glacier), or they can consist of black carbon (so they derive from fossil combustion and fires), as well as living organisms, pollens and other vegetal and organic residuals that remain in the aerosols (Fujita, 2007;Takeuchi et al., 2001;Takeuchi, 2002).The mineral fraction can be locally derived from the weathering of rock outcrops and nunataks or from lateral moraines and debris slopes.In fact, during the summer, when warmer climatic conditions occur, the dry and unconsolidated materials constituting moraines are easily transported by wind gusts and deposited tens to hundreds of meters away and even higher, depending on wind speed and the roughness of the area surrounding the glacier (Oerlemans et al., 2009).In the case of englacial origin, the fine debris can also originate from mechanical disintegration of the bedrock below the glacier or the deformation and weathering processes of the rocks embedded in the ice.Dust and fine debris can also be transported long distances by atmospheric circulation from nonglaciated areas (Ming et al., 2009;Ramanathan, 2007).For instance, the deposition of Saharan dust (Sodemann et al., 2006) or volcanic ash (Conway et al., 1996) on glaciers is a well-known phenomenon that darkens mountain debris-free glaciers (e.g., Paul and Kääb, 2005;Paul et al., 2007;Oerlemans et al., 2009;Diolaiuti and Smiraglia, 2010;Casey, 2012;Painter et al., 2013) thus changing their albedo and affecting melt magnitude and rates.
The occurrence of fine debris at the surface of debris-free glaciers and the role it plays in ice melting rates make correctly determining ice albedo important.The albedo parameterizations used in energy and mass balance models are, however, often inadequate to represent spatial and temporal changes in the surface albedo and are consequently regarded as a major source of errors (e.g., Arnold et al., 1996;Klok and Oerlemans, 2002;Klok et al., 2003).Therefore, studies that combine measurements of fine debris distribution and features with systematic measurements of glacier albedo are needed.
Dust deposition on snowpacks has been well studied (e.g., Qian et al., 2011;Yasunari et al., 2010).A possible snow albedo reduction due to black carbon contamination was revealed by radiation measurements at the snow surface performed at Barrow, Alaska (Aoki et al., 1998(Aoki et al., , 2006)), and in Japanese urban areas (Motoyoshi et al., 2005).In the case of snow and firn, a field procedure was developed, followed by further standardized lab analyses to quantify and describe black carbon presence and features (Yasunari et al., 2010).
However, less attention has been paid to fine debris and dust deposition at the glacier melting surface.A first attempt to parameterize not only snow albedo variability but also the ice albedo on a debris-free glacier was performed by Brock et al. (2000).In spite of the good results they obtained analyzing snow-covered areas, their evaluation of the impact of debris cover on ice albedo was less accurate.In fact, they assessed the debris cover using only a 0.5 m 2 quadrat and basing their investigation on just two criteria (i.e., cumulative melt and number of days, both calculated following exposure of the ice surface).Recently, Pope and Rees (2014) investigated the spectral responses of different ash/debris cover types on the glaciers of Midtre Lovénbreen (Svalbard) and Langjökull (Iceland).These studies suggested the need for further research to standardize the measurements of fine debris and dust at the glacier ice surface, thus avoiding the use of surrogates unable to fully describe debris coverage and its seasonal variability.

Study area
Our experiments were carried out on the ablation tongue of the Forni Glacier (Fig. 1), the widest Italian valley glacier, featuring a surface area of 11.34 km 2 (2007 data, Garavaglia et al., 2012).It is located in the Ortles-Cevedale group, Stelvio National Park, Lombardy Alps.It is widely debrisfree, even if darkening phenomena are ongoing (D'Agata et al., 2014), and some authors have recently pointed out that fine and sparse debris is becoming abundant due to the ongoing glacier shrinkage (Diolaiuti and Smiraglia, 2010;Diolaiuti et al., 2012;Senese et al., 2012a).For this reason, the Forni Glacier can be considered a good laboratory to evaluate fine and sparse debris distribution and seasonal evolution and its influence on ice albedo.The Forni Glacier is facing toward north, is about 3 km long and its altitude ranges from 2600 m to about 3670 m a.s.l.Metamorphic rocks, mostly mica schist rich in quartz, muscovite, chlorite and sericite, constitute the dominant lithology (Montrasio et al., 2008); these rocks emerge from the glacier surface as nunataks (mainly in the accumulation basins) and as rock outcrops (surrounding the glacier tongue).The latter are increasing in size and are becoming very frequent due to the ongoing glacier retreat and thinning (Diolaiuti and Smiraglia, 2010;Diolaiuti et al., 2012).
Studies on short-term changes of the Forni Glacier have been performed through an automatic weather station (named AWS1 Forni) that has been in operation since 2005 at the glacier melting surface.The AWS1 Forni is located on the ablation tongue (∼ 2631 m a.s.l.), about 800 m from the glacier terminus, and it is equipped with sensors for measuring air temperature and humidity, wind speed and direction, atmospheric pressure, liquid precipitation and snow depth, and long-wave and short-wave radiation, both incoming and outgoing (Citterio et al., 2007;Diolaiuti et al., 2009;Senese et al., 2010Senese et al., , 2012aSenese et al., , b, 2014)).

Methods
In the time frame 2011-2013, 51 field measurements in total were obtained on the debris-free ablation tongue of the Forni Glacier (Fig. 1b), both for fine debris quantification (i.e., spatial coverage) and albedo evaluation (Fig. 2).Moreover, we sampled fine debris at the glacier melting surface to assess the debris coverage rate (C r ) and the sedimentological properties (i.e., grain size, humified and total organic carbon and mineralogical properties).The sites for field measurements were chosen considering (i) homogeneity in debris cover, (ii) presence or absence of fine, sparse debris, (iii) diverse debris grain size and (iv) different distances from rock slopes and medial moraines, which are the main debris suppliers, thus assuring that the selected sites are representative of the range of surfaces present at the glacier melting area.At each site, we sampled a 1 m × 1 m parcel with the aim of assuring the effectiveness and the repeatability of the measurements.Larger areas would have required more time for data collection and would have limited the number of sites analyzed, while smaller quadrats would not have captured the spatial heterogeneity of the surface.Each sampling area was selected to be representative of as wide an area as possible.The medial moraines were excluded from this work because we only focused on fine-and sparse-debris-covered ice and not on actual buried ice (i.e., ice covered by a thick and almost continuous debris layer).Figure 1 shows the study area and the positions of the sites analyzed.

Debris cover quantification
The quantification of sparse and fine debris at the glacier melting surface was performed by acquiring high-resolution digital images at each site analyzed (Fig. 2c) and processing them with image analysis software ImageJ following Irvine-Fynn et al. (2010) (Fig. 3).Digital RGB (red-green-blue, in color composite) photographs of the 1 m × 1 m parcel were www.the-cryosphere.net/10/665/2016/The Cryosphere, 10, 665-679, 2016 taken using a digital camera (Nikon D40, 6.1 megapixels).
The images affected by shadows, deformations and photographic imperfections (e.g., poor exposure, incorrect focus) were excluded from the analysis, and for each measurement site (total 51), we selected the image that best captured sharp differences between bare ice and fine-debris-covered ice.The selected images were first cropped delimiting the 1 m × 1 m parcel (Fig. 3a).Second, we converted them to 8 bit greyscale in order to highlight the contrast between glacier ice and debris/dust.Third, as a darker grey pixel denotes the presence of debris or water or shadow (the latter due to surface roughness), we assumed that debris granules could be isolated by thresholding for those pixels with brightness values that fall below a specified grayscale threshold level (T GS , Fig. 3b), specified by a supervised classification.In particular, the threshold was iteratively adjusted until the isolated image pixels best coincided with the debris/dust (Fig. 3c).An 8 bit image is composed of 256 grey tones ranging from 0 (black) to 255 (white), and ice surfaces can be isolated by selecting the pixels with brightness values higher than a specified T GS .For instance, if the T GS value is fixed at 100, pixels with a grey tone from 0 to 100 represent debris and pixels with a grey tone from 101 to 255 represent ice.For each image, the pixels with a value lower than T GS were changed to black and the other ones were changed to white (Fig. 3d).Finally, the ratio of the surface covered by debris (d) was obtained as d = number of black pixels total number of pixels , where the total number of pixels is 6.1 × 10 6 .As only the greyscale threshold choice is manual and the other steps are automatic ones, the proposed method for quantifying the ratio of the glacier surface covered by debris (d) can be defined as semi-automatic.The reliability of this method was evaluated through five tests.More precisely, we selected 10 images from the 51 photos used in this work and we compared the d values obtained from the chosen T GS data (i.e., applying our semi-automatic procedure) to d values derived from changed T GS or from the application of different methods.These d values are (i) d 10PI data derived from the point intercept approach (i.e., another largely applied method based on a visual estimation by placing a grid on the investigated area; for a detailed discussion of this approach see Elzinga et al., 2001), (ii) d 10IJ data derived from the application of our procedure by several untrained users (thus showing the sensitivity of our method to changes in the user), (iii) d +10 % and d −10 % values obtained varying the selected T GS up to ±10 % of its initial value (T GS+10 % and T GS−10 % , respectively), (iv) d MOD data obtained selecting the modal greyscale value as the threshold (T GS-MOD ) and (v) d AVE data derived averaging all 51 T GS values to obtain an average threshold (T GS-AVE ).

Albedo
The bihemispherical reflectance, generally called albedo (α), is defined as the ratio of the radiant flux reflected from a unit surface area into the whole hemisphere to the incident radiant flux of hemispherical angular extent (Schaepman-Strub et al., 2006) in the approximate spectral range 350-3000 nm (Grenfell, 2011).The albedo is an apparent optical property.This means that it depends on the angular distribution and spectral composition of the ambient radiation field as well as on the inherent optical properties, which depend only on the structural and optical properties of the medium (Grenfell, 2011).Thus, it is important to consider solar elevation, cloudiness, presence of liquid water, crystal structure, ice surface conditions and the presence or absence of materials at the surface (rock debris, dust, organic matter, etc.).It is estimated as the ratio of measured outgoing shortwave (SW out ) to measured incoming shortwave (SW in ): For this study, the albedo was calculated from radiation data measured using two pyranometers (the ones installed in the net radiometer CNR1, Kipp & Zonen; see Fig. 2).The sensor features an accuracy of ±5 %.The net radiometer was equipped with a waterproof box containing a data logger, a 5 Ah battery and a 10 W solar panel on the lateral face.Moreover, a tripod was used to raise the net radiometer for short periods (∼ 20-30 min for each measurement) above the ice surface.Tests regarding the influence of the height of the sensor above the surface on albedo values were performed, installing the sensor at various distances from the surface.Since we chose sites featuring homogeneous surfaces, we did not find appreciable differences between albedo values measured at the same site varying the distance between the instrument and the ice surface, thus suggesting a negligible influence of the sensor height.The CNR1 net radiometer was chosen for its accuracy and resolution in measuring short-wave radiation, and it is also the same type as the one running at the AWS1 Forni (Citterio et al., 2007;Senese et al., 2012a), thus assuring the comparability between the two data sets (accomplishing the recommendation described by Grenfell, 2011).Then the radiation data collected using the portable instrument were crosschecked and analyzed against the data acquired by the AWS1 Forni.The measurements were carried out following the guidelines of the WMO (2008).
The radiation data were acquired every second, and every minute the minimum, average, maximum and standard deviation values were calculated.Albedo measurements were taken in the central hours of clear-sky days (i.e., from 11:00 to 15:00 LT (local time), when the solar incidence angles are smaller), thus ensuring the greatest possible accuracy and reliability of the albedo calculations (Brock et al., 2000;Brock, 2004;Oerlemans, 2010).The mean geographic coordinates (WGS84 datum) for each measurement site were recorded by a GPS receiver and the features characterizing the local ice surface were also noted.Fifty-one measurements were carried out from the beginning of the ice ablation period (when snow coverage at the melting tongue disappeared, exposing ice to solar radiation and dust/debris deposition) to the end of the ice melting season (before the occurrence of the first snowfall event covering the glacier ice and preventing dust/debris deposition): on 30 June and 25 August 2011, 4 July, 7 August and 9 September 2012 and 31 July and 6 September 2013.
In addition to debris, water plays a significant role in changing ice albedo; water washes out the finer sediments on the glacier (Oerlemans et al., 2009) and makes the ice surface smooth.Thus, the effect of water (derived from both melting processes and rainfall) on glacier albedo variability was assessed during each ice ablation season from 2011 to 2013.The length of the ice ablation period was investigated, coupling albedo and melting data (i.e., considering the time window featuring melt and with an albedo lower than 0.40, more details in Senese et al., 2012a).The occurrence of melting was investigated by applying the energy balance model from meteorological data and energy fluxes measured by the AWS1 Forni (for more details regarding the melting model, see Senese et al., 2012aSenese et al., , b, 2014)).Finally, the temporal length (i.e., number of rainy days) and amount (i.e., mm of rain) of liquid precipitation were measured by an unheated pluviometer installed at the AWS1 Forni (DQA035, LSI-Lastem).The effect of liquid precipitation was quantified by comparing albedo values before, during and after the occurrence of liquid precipitation.Any event featuring an hourly liquid precipitation higher than 0.2 mm (i.e., the threshold to activate the toggle switch of the rain gauge) was considered to be rainfall.

Sedimentological analyses and debris coverage rate evaluation
Several bulk samples of sediment were collected from the glacier surface (Fig. 2a) and divided into subsamples for physical and chemical analyses.In 2011, eight samples were collected, choosing surfaces with diverse debris grain size and different distances from rock slopes and medial moraines, and these samples were used to characterize the spatial variability of debris at the glacier melting surface.
Then in 2012 (4 July, 7 August and 9 September), the temporal evolution of debris features was studied by sampling three sites (identified by ablation stakes) with different conditions of debris cover: (i) samples 9a, 9b and 9c fine and sparse sediment, (ii) samples 10a, 10b and 10c widespread debris cover, and (iii) samples 11a, 11b and 11c coarse debris.Finally, in 2013 we assessed the debris coverage rate (C r , the fine debris amount reaching the surface over a defined time frame).The samples were collected four times (samples 12a, 12b, 12c, 12d): 11 July, 31 July, 6 September and 4 October 2013.
Each sample was collected by scraping the glacier surface with a cleaned chisel, completely removing the surface layer (from 2 to 5 cm deep, depending on the surface roughness); the collected material was preserved in appropriate holders.A cold chain (ice boxes) was used to preserve sediment samples at cold temperature conditions (lower than +4 • C) during transport to the laboratory, where further analyses were carried out.
For evaluating the debris coverage rate, debris samples were periodically collected from the same sites.First, it was necessary to clean the 1 m × 1 m parcel, completely removing surface debris (i.e., scraping at least 2 cm of surface ice).Second, about 1 month later, the sampling of the surface sediments was repeated on the same glacier parcel, which was marked on the field.Then in the lab, debris samples were dried and weighed.The ratio between the weight of the debris deposited at the ice surface (in g m −2 ) and the time frame (days) permitted evaluation of the debris coverage rate (C r in g m −2 per day): The samples collected in 2011 and 2012 for evaluating the spatial and temporal variability were subjected to the analytical procedures, summarized as follows.Grain size analyses (Gale and Hoare, 1991) were performed after removing organics using hydrogen peroxide (130 vol) treatment; sediments were wet-sieved (diameter from 1000 to 63 µm), then the finer fraction (63 µm) was determined by aerometer on the basis of Stokes's law.Humified organic carbon was iden- tified by means of the Walkley and Black (1934) method, using chromic acid to measure the oxidizable organic carbon (titration).Total organic carbon (TOC) was estimated by loss on ignition (Heiri et al., 2001), with an uncertainty margin of ±0.1 %; samples were air-dried and organic matter was oxidized at 500-550 • C to carbon dioxide and ash; the samples were weighed before and after heating to calculate the quantity lost during the reaction.Additionally, we performed several XRD (X-ray diffraction) and SEM (scanning electron microscope) analyses on randomly oriented powder from the bulk debris samples to investigate the mineralogical properties of the fine debris and dust and the occurrence of micro features (e.g., pollen, spores, micro-and mesofauna, algae).

Debris coverage ratio (d) and ice albedo (α)
Image analysis yielded 51 d values ranging from 0.01 to 0.63 (Fig. 4).The ice albedo acquired by the portable net radiometer varied from 0.06 to 0.32.The two data records (i.e., d and α) are highly correlated.A plot showing the natural logarithm of ice albedo (ln α, y axis) vs. d values (x axis) is reported in Fig. 5.The regression line is given by The correlation is −0.84 (the 95 and 99 % confidence intervals ranging from −0.91 to −0.74 and from −0.92 to −0.69, respectively); the p value of the correlation coefficient is lower than 10 −9 .For low values of debris cover ratio, the correlation appears less accurate.This can be due to the occurrence of other influencing parameters that become dominant whenever debris is limited or absent (i.e., d < 0.10).Among the most important factors, bubbles and other air inclusions modulate the volume scattering and then albedo (see Mullen and Warren, 1988).Moreover, Grenfell (2011) reported ice inhomogeneities (also at a microscale) to be significant in determining albedo.Nevertheless, these other factors become negligible whenever the debris cover ratio is higher than 0.10.The most frequently occurring d value was 0.03, but the wide variability of surface features indicates that a large number of samples is required to describe a d pattern and albedo distribution.Points featuring a high debris cover ratio are less numerous than the ones showing low values.In fact, in our study we considered glacier areas featuring fine and sparse debris coverage and not zones characterized by an almost continuous debris coverage (d > 0.60), such as medial moraines or actual buried ice sectors.
In addition to debris occurrence, ice albedo depends on water presence.To evaluate the impact of water on albedo, we considered both melting processes and rainfall occurrence during each ice ablation season from 2011 to 2013.The beginning and the end of each ice melting period are shown in Table 1.The meteorological data from 3 to 13 July 2013 are lacking; however, on 3 July, the albedo was equal to 0.55, indicating snow cover, while on 13 of the same month the albedo was 0.18.The latter value is characteristic of bare ice, thus indicating that the ablation season started between these dates.
By coupling melt and albedo data, higher ablation rates were found to correspond to decreases in ice albedo (Fig. 6).When melting was less intense, the albedo decrease was smaller (e.g., on 20 June 2012, see Fig. 6).In general, we can deduce that meltwater occurs during the central hours of the day (when the solar radiation input is higher, the melting processes are more intense), decreasing the albedo due to its lower reflectivity (i.e., equal to 0.05-0.10;Hartmann, 1994).This lower albedo implies a more intense absorption of incoming solar radiation, which leads to more energy being available for melting.As a consequence, these factors (i.e., solar radiation input, melt and water) have a positive feedback in influencing albedo.This trend is also found when analyzing the liquid precipitation.In particular, we report in Table 1 the mean daily albedo values before, during and after rainfall events.First, the days before the rainfall featured a mean daily albedo equal to 0.22.Second, whenever precipitation occurred the mean daily reflectivity was reduced to 0.20, probably due to water albedo being lower than ice.This phenomenon occurred in 18 of 30 events.Of the remaining 12 cases, 6 featured an albedo increase and 6 steady-state albedo conditions.This variable trend can be attributed to the rain amount; in fact, a misty rain decreases the surface albedo less than a heavy liquid precipitation.Third, once the rain event has washed out the dust and smoothed the surface, the mean daily albedo resulted to be 0.26.Almost all rain events (28 over a total number of 30) showed a mean daily albedo increase slightly higher than 20 %.In contrast, when albedo before the rainfall was higher than 0.30, the water effect was not so appreciable.In fact, this reflectivity value is typical of bare ice with no fine debris coverage.
The occurrence of the rainfall washing out effect, and the consequent reflectivity increase, was found to be shortlasting.The mean time period to restore the previous albedo value was 1.8 days (ranging from 1 to 4 days), occurring over 10 events out of a total of 30.The rain effect is less evident and effective whenever the time interval between two rainfalls is very restricted (i.e., 1 day).

Debris composition and debris coverage rate (C r )
The spatial variability of fine debris cover is highlighted from 2011 data.The sediment analysis performed in the laboratory indicates significant variability in total organic carbon (TOC, from 0.6 to 5.9 %; see Table 2).The highest content of organic matter was found in samples 5, 6 and 8; in particular, sample 5 was wholly cryoconite, where generally the development of algae and bacteria communities is extremely favored (Takeuchi et al., 2000(Takeuchi et al., , 2005)).The lowest value of total organic carbon was found in sample 2, which was collected on a glacier area located close to the flank of the nesting rock walls, a site which receives a high amount of debris originating from rock weathering processes such as the macrogelivation, which in this area has been reported by Guglielmin and Notarpietro (1997).Rock debris deposits in the area are active and unstable; they are continuously suffering renewal of the surface, and are therefore poorly colonized by supraglacial organisms.Moreover, the grain-size analysis shows that samples collected at these sites are characterized by coarser sediments, in keeping with their origin, mostly due to mechanical weathering.
Regarding the lithology of the debris, X-ray diffraction indicates that the samples are enriched with quartz, muscovite, chlorite, sericite and albite.This reflects the local geological bedrock, which corresponds to mica schist (Montrasio et , 2008;Chiesa et al., 2011).According to regional geological maps (Montrasio et al., 2008;Chiesa et al., 2011), we may exclude a contribution in the formation of debris from nearby localities.In fact, outside the Forni Glacier basin at about 8 km northward, a lithological and tectonic discontinu-ity (namely the Zebrù Line) is located, where the carbonatebearing sedimentary rocks (i.e., dolomite) outcrop.Moreover, at about 18 km southward, the intrusive rocks of the Adamello pluton are present.SEM analyses (Fig. 7) revealed algae, spores, pollen and mesofauna.Moreover, we also observed spherical structures (Fig. 7d) characterized (from energy dispersive spectroscopy (EDS) analysis) by an abundance of FeO (45.3 %) and Al 2 O 3 (19.8%).This composition confirms that these structures are cenospheres (i.e., a residual product of carbon combustion; see Kolay and Singh, 2001).These cenospheres may have been transported by wind and are probably derived from diesel fuel combustion or also from siderurgic factories, suggesting allochthonous inputs and human impacts at the glacier surface, even if limited.However, the assessment of the source and origin of fine debris was not of particular interest in this study and therefore remains at least partly open.
To assess the evolution of the fine debris cover, the debris coverage rate (C r ) was evaluated (see Eq. 3) from 11 July to 4 October 2013 (i.e., sample 12, Table 2) and it was found to be equal to 6 g m −2 per day.Immediately after each sampling, the cleaned 1 m 2 parcel could be clearly distinguished from the glacier areas nearby (see also the portion of cleaned ice in Fig. 2a).However, in the following survey the sampled parcel became completely covered by fine debris and it was identified only thanks to the signals we put on the field.This suggests that the development of debris coverage occurs at a fast rate.This evolution is also highlighted from the sedimentological analyses performed on the 2012 samples.During the ablation season, the grain size remained almost equivalent with a slight increase of finer sedimentological classes (i.e., silt and clay).However, a rapid increase in total organic carbon (TOC) along the season was observed.At the beginning of July the TOC ranged from 1.6 g kg −1 (at sample 9a) to 26.3 g kg −1 (at sample 10a); at the end of the ablation period the organic carbon increased up to 41.9 g kg −1 .The higher values correspond to finer debris (i.e., samples 9a, 9b, 9c, 11a and 11c, enriched in silt and clay); conversely, in coarser samples 10a, 10b and 10c, the TOC is lower.
The debris evolution is also analyzed through SEM observations.At the beginning of the melting time frame, the sediment was characterized by sharp and angular clasts, suggesting a supraglacial mass transport; the samples collected in September 2012 featured more rounded shapes, suggesting an englacial mass transport.

Accuracy assessment of semi-automatic debris cover quantification
We performed several tests to evaluate the robustness of our semi-automatic method to quantify α from the estimation of d.Firstly, we asked a representative sample of users (10 geologists, but non-glaciologists) to apply different approaches to 10 images randomly selected from the whole sample (see Fig. 4).In particular, they estimated the debris coverage ratio applying both the point intercept method obtaining d 10PI (see an example in Fig. 8b) and the ImageJ procedure we proposed obtaining d 10IJ from T GS-10IJ (see Fig. 8c).The average d 10PI and d 10IJ for each of 51 images (d 10PI-AVE and d 10IJ-AVE , light-blue and green circles in Fig. 9) were compared with the debris coverage ratio values obtained by the well-trained operator and used in Eq. ( 4) (d, red diamonds in Fig. 9).The d 10PI values were affected by a very high standard deviation (from ±0.09 to ±0.34, see the light-blue errors bars in Fig. 9, displayed only for positive errors), sug- Secondly, as almost all T GS-10IJ values fell within ±10 % of the corresponding T GS values, we also selected the extremes of this interval (T GS−10 % and T GS+10 % , respectively) to investigate the sensitivity of our method to changes in the chosen T GS .For example, whenever the applied T GS value was 100, we recalculated d with 90 (T GS−10 % , obtaining d −10 % , black dashes in Fig. 9) and with 110 (T GS+10 % , obtaining d +10 % , black crosses in Fig. 9).We applied this test also to all the measurements taken (51 field data), thus obtaining 51 d −10 % values and 51 d +10 % values.Subsequently, we evaluated the departures of d −10 % from d (i.e., d − d −10 % ) which resulted up to −0.07 (with a mean value of −0.02).Moreover, we calculated the departures of d +10 % from d (i.e., d − d +10 % ), which were found to be lower than +0.09 (with a mean value of +0.02).We discovered that whenever d was higher than 0. In addition, we used d −10 % and d +10 % values to obtain two more relations with α (reported in Table 3 together with Eq. 4).Applying these three equations to the d data set, it was possible to estimate glacier albedo; the albedo values modeled with the three relations were compared to the albedo data obtained from the field radiation measurements and the departures between the modeled and the observed records proved to be very small, with a mean value lower than ±0.01 and a standard deviation of 0.04 (Fig. 10).This suggests that besides giving evidence of a good performance in estimating albedo, our method is also robust to changes in the applied threshold, supporting the replicability of the results.
Finally, we tested the utility of a user-defined threshold, thus different for each image.For this experiment, we performed two other sensitivity tests: we considered the most frequent grey tone for each image (the top of the curve in Fig. 11c), obtaining 51 T GS-MOD values, and we used a unique value, averaging all the 51 T GS values (T GS-AVE ).This latter threshold was found to be equal to 92.
Applying T GS-MOD gave an incorrect selection of the pixels; in particular, some pixels with clean ice are selected as debris-covered ones, thus overestimating the d value (an example is shown in Fig. 11 where d MOD is 0.50 instead of the actual d value of 0.10).In fact, the most frequent grey tone could correspond to pixels of ice featuring a thin film of water and not covered by debris.This is also evident from Fig. 9 (red stars).
As far as T GS-AVE is concerned, Fig. 12 shows the comparison between d (obtained from T GS ) and d AVE (obtained from T GS-AVE ) values.The relation between the two data sets was not negligible, thus suggesting that a unique threshold value could be sufficient to describe debris distribution on different images (see also yellow triangles in Fig. 9); then we applied the obtained 51 d AVE values to calculate a new relation with α data:  (5) Equation ( 5) features a R value of −0.68 (p value < 10 −6 , Table 3), meaningful but lower than the one given by Eq. ( 4), thus suggesting that the different T GS values we found for each image (even if they required spending more time in the image analysis) permit a better and more detailed determination of debris distribution, and consequently a more accurate d evaluation.

Discussion
The main task of this study was to evaluate the role of fine debris in modulating ice albedo, thus in driving ice ablation.Accordingly, we measured and analyzed the present debris coverage at the melting surface of a wide and representative Alpine glacier, finding an actual relation with measured ice albedo.An advantage of our approach is that a thorough analysis of debris origin (autochthonous or allochthonous) and history (from englacial origin or transported by wind) is not required.Nevertheless, research focused on the human impact, if any, on glacier ice (mainly black carbon deposition) should also consider these issues.In this case, distinguishing between local debris and particles transported by wind, and also considering the paleo implications of emerging englacial debris would be desirable.Moreover, our findings support the recent literature describing darkening phenomena occurring at the Alpine debris-free glaciers (Paul and Kääb, 2005;Paul et al., 2007;Oerlemans et al., 2009;Painter et al., 2013;Diolaiuti et al., 2012;Gabbi et al., 2014;Brun et al., 2015;Naegeli et al., 2015) and in particular at the Forni Glacier tongue (Diolaiuti and Smiraglia, 2010;D'Agata et al., 2014).Indeed, we found that the albedo range (i.e., 0.06-0.32)agrees with the values characteristic for debris-rich ice albedo found by Cuffey and Paterson (2010).They reported a range of 0.06-0.30,whereas clean ice is described as varying from 0.30 to 0.46.Moreover, on average from our results the ice albedo decreases along the ablation tongue, becoming more absorptive.This can be due to high melt rates and long exposure times (Klok and Oerlemans, 2002).Effectively, debris can be concentrated not only over a single melting season, but cumulated over many years (Cuffey and Paterson, 2010).In addition to spatial variability, we also observed that albedo decreases over time in the melt season (see also Senese et al., 2012a).This entails that the period of most effective energy absorption occurs later than the peak of insolation.Maximum of exoatmospheric insolation is at the solstice (except in the tropics), but much of the winter snow still covers glaciers at this time.As the melt season progresses, darkening of the surface increases the absorbed short-wave radiation available for melt by a factor of 3 to 4 (Cuffey and Paterson, 2010).This is one possible reason (in addition to warmer air in midsummer) to explain why the peak of melt rates occurs on midlatitude glaciers 1 to 2 months after the solstice (Cuffey and Paterson, 2010).

Conclusions
In spite of the quite abundant literature focusing on fine debris deposition over glacier accumulation areas, less attention has been paid to the glacier melting surface.Accordingly, we proposed a novel method based on semi-automatic image analysis to quantify fine debris coverage (d) on glacier ice.We tested this procedure over the widest Italian valley glacier (Forni, Stelvio National Park).Analysis of 51 images yielded d values ranging from 0.01 to 0.63.Together with image acquisition, we also measured ice albedo, which varied from 0.06 to 0.32.The estimated d values were found to be in a linear relation with the natural logarithm of measured ice albedo (R = −0.84).We performed five sensitivity tests to investigate the robustness of our approach in evaluating d. (i) 10 users estimated d 10PI using the point intercept method and their results were affected by a high variability, thus suggesting strong subjectivity.(ii) The same 10 users quantified d 10IJ using our approach and these values featured a variability lower than 10 %, thus supporting a wide replicability by untrained operators.(iii) We computed d −10 % and d +10 % , varying the selected threshold up to ±10 % for each image (T GS−10 % and T GS+10 % ), and we found that more than 70 % of departures of d −10 % from d were up to −0.02 and of d +10 % from d were lower than +0.03.(iv) We derived d MOD by applying as T GS the most frequent greyscale value (T GS-MOD ), but in several cases this caused an incorrect selection of the debris pixels.(v) We computed d AVE , applying a unique averaged threshold (T GS-AVE ) which led to a slight loss of accuracy.In summary, these tests suggest that the method we propose gives results featuring a lower variability thus suggesting a lower subjectivity than the point intercept method.Moreover, the sensitivity of our method to changes in the applied threshold is not so high as to affect the reliability of the results also supporting the replicability of the approach.In addition, even if T GS analysis is more time-consuming, the specific T GS value found for each image permits a better and more detailed determination of debris distribution and therefore a more accurate d evaluation.Then, to look for the most reliable relation between d and α, the best and most suitable solution is to apply a specific T GS value for each image.
On the Forni Glacier, we also evaluated the debris coverage rate (C r ) during ablation season 2013.Our data show a mean C r equal to 6 g m −2 per day, thus supporting the recent literature which describes the ongoing darkening phenomena at the Alpine debris-free glaciers (Paul and Kääb, 2005;Paul et al., 2007;Oerlemans et al., 2009;Diolaiuti and Smiraglia, 2010;Painter et al., 2013;D'Agata et al., 2014).In addition to debris coverage, we considered both melting processes and rainfall occurrence during each ice ablation season from 2011 to 2013.Our analyses indicate that meltwater occurs during the central hours of the day, decreasing the albedo due to its lower reflectivity.This lower albedo implies a more intense absorption of incoming solar radiation, which leads to more energy being available for melting.Consequently, these factors (i.e., solar radiation input, melt and water) play a positive feedback in influencing albedo.In addition, we found that almost all rain events caused a mean daily albedo increase slightly higher than 20 %, although it was short-lasting (from 1 to 4 days).
In conclusion, the semi-automatic image analysis method we are proposing could potentially be applied to images acquired by an unmanned aerial vehicle (UAV), covering the entire ablation tongue of the Forni Glacier.The UAV could provide an excellent basis for capturing high-resolution images of a much larger expanse of the glacier surface than would have been possible using ground-based surveillance (Hodson et al., 2007;Fugazza et al., 2015).Once the extension of the surface covered by fine debris is quantified, the relation between d and α could be applied, thus deriving an albedo map.The latter will feature a higher resolution (pixel size of a few centimeters) than the ones derived from other approaches such as Landsat images (pixel size of 15 m × 15 m, Fugazza et al., 2016).

Figure 1 .
Figure 1.(a) The position of Forni Glacier in the Italian Alps and of AWS1 Forni (black star, in both panels).(b) Enlarged view of the Forni Glacier (image credit: GoogleEarth TM ), showing the location of field measurements (black dots).

Figure 2 .
Figure 2. Series of pictures illustrating (a) sampling supraglacial debris, (b) measuring albedo and (c) acquiring high-resolution digital images of the glacier surface.

Figure 3 .
Figure 3. Example of the procedure followed in image analysis: (a) original cut frame; (b) 8 bit conversion and discrimination between the debris-covered (in blue) and debris-free ice surface; (c) definition of the threshold; (d) calculated debris cover ratio.

Figure 4 .
Figure 4. Examples of Forni Glacier surfaces; debris cover ratio (d) and measured albedo (α) values are shown.The first and last images show the minimum and maximum d values.Images 1-10 are used for further sensitivity tests (see Sect. 6).

Figure 6 .
Figure 6.Comparison between albedo measured by the AWS1 Forni (grey line) and melting amount estimated by the energy balance model (black line).

Figure 7 .
Figure 7. SEM investigation on bulk samples from Forni Glacier evidenced the presence of (a) organisms of collembolan order, (b) spores, (c) diatoms and (d) cenospheres (a residual product of carbon combustion).

Figure 8 .
Figure 8.An example of the application of different approaches to the same image.The original image is shown in (a).In (b) and (c), the application of the point intercept method and our approach, respectively.

Figure 9 .
Figure 9. Values of debris cover ratio from analyzing 10 randomly selected images (shown in Fig. 4).

Figure 10 .
Figure 10.Comparison between measured and modeled albedo values.The applied equations (Table3) are shown in the legend.
25, the d − d −10 % and d − d +10 % values reached their maxima.In fact, d − d −10 % was −0.07 with a d value equal to 0.28 and d − d +10 % was +0.09 with a d value of 0.60.Considering the sample as a whole, slightly more than 70 % of d featured d − d −10 % values of up to −0.02 and d − d +10 % values lower than +0.03.

Figure 11 .
Figure 11.Comparison between the application of the threshold chosen correctly by the user, T GS (b), and the one corresponding to the most frequent grey tone, T GS-MOD (d), deducted by the frequency distribution curve (c).The relative d and d MOD values are shown in (b) and (d).In (a), the analyzed surface is presented.

Figure 12 .
Figure 12.Scatter plot reporting d (obtained from T GS ) vs. d AVE (obtained from T GS-AVE ) values.

Table 1 .
Influence of rainfall on ice albedo (α) measured from the AWS1 Forni.In the table the 30 rain events (and relative rain amount) that occurred during the 2011, 2012 and 2013 ablation seasons, and the albedo values before, during and after every rainfall are reported.

Table 2 .
(Wentworth, 1922)es sampled for sedimentological analyses and results.The grain-size classes refer to Krumbein's scale(Wentworth, 1922).An asterisk indicates a debris sample insufficient for grain-size analysis.

Table 3 .
Depending on the four different ratio data sets (considering d +10 % , d, d −10 % and d AVE ), the differences between measured (α M ) and calculated (α C ) albedo data are shown.