Recent accumulation rates of an alpine glacier derived from firn cores and repeated helicopter-borne GPR

Introduction Conclusions References Tables Figures


Introduction
Mountain glaciers are known to be excellent indicators of climatic change (Haeberli and Beniston, 1998;Kaser et al., 2006;IPCC, 2013), but are also important to the hydrology and ecology of alpine regions (Beniston, 2003).In order to better understand their response and adaption, glacier changes are monitored around the globe (Barry, 2006;Zemp et al., 2009) of the two approaches is subject of ongoing research (Huss et al., 2009;Fischer, 2011;Zemp et al., 2013).Field measurements of ablation can be taken accurately at many locations across a glacier.The main drawback is the lack of accumulation measurements that typically involve the time-consuming excavation of snow pits (Østrem and Brugman, 1991).Thus, measurements of accumulation are under-represented and the spatial variability is often not resolved correctly in modern monitoring approaches (e.g., Huss et al., 2014).To compensate for this drawback, an adequate complement to conventional accumulation measurements is necessary, which should preferably be nondestructive and efficient.
In glaciology, ground-penetrating radar (GPR) is used for a variety of purposes (Plewes and Hubbard, 2001).The low conductivity of snow and ice facilitates the deep penetration of the signal, and the use of narrow-bandwidth, high-resolution systems (Ulriksen, 1982).Because the signal is reflected from boundaries with a contrast in dielectric permittivity, a change in the related material properties can be detected, such as from ice to bedrock, water content, impurities, but also changes in density are resolved (Plewes and Hubbard, 2001;Woodward and Burke, 2007).Thus, GPR applications range from ice thickness measurements (Robin et al., 1969;Bauder et al., 2003;Damm, 2004), to the mapping of snow accumulation distribution (Machguth et al., 2006;Sold et al., 2013;Helfricht et al., 2014), to the observation of internal layers within the snow cover (e.g., Heilig et al., 2010).Thereby, airborne GPR surveys can cover large areas in short periods of time.On ice sheets, with rather homogeneous accumulation patterns, internal reflection horizons (IRH) within firn and ice can be tracked over longranging GPR profiles (Arcone et al., 2004;Huybrechts et al., 2009;Miège et al., 2013;Hawley et al., 2014).Estimates of past accumulation rates are typically obtained from the age and density of several pronounced but not necessarily annual IRH, which provides the average mass balance for each sequence.Only few studies have observed and analysed the water equivalent of annually spaced IRH in radar data (e.g., Kohler et al., 1997;Hawley et al., 2006;Kruetzmann et al., 2011;Van Pelt et al., 2014), however, not on temperate mountain glaciers in mid-to low latitudes.Introduction

Conclusions References
Tables Figures

Back Close
Full On temperate glaciers, or glaciers that are constrained by topography and have complex flow fields, the continuous tracking of IRH is often difficult (Karlsson et al., 2009).Therefore, studies have focused on the mapping of facies zones and other general glaciological characteristics (Pälli et al., 2003;Langley et al., 2008;Eisen et al., 2009;Dunse et al., 2009).On the other hand, GPR could not only accompany ongoing monitoring programs to complement the sparse accumulation measurements with high spatial resolution, but could also be used to retrospectively extend newly initiated time series into the past.
To extract past annual accumulation rates from the GPR signal, the IRH must correspond to previous summer surfaces.On mid-latitude glaciers this suggests a high potential of GPR for mass balance studies due to the numerous melt-refreeze cycles that can generate a high-density or ice layer at the snow surface during summer.However, this precondition needs to be verified by independent layer dating information.Additionally, to transform the GPR traveltime signal to the depth-domain and, ultimately, to obtain the water equivalent, an estimate or measurement of the firn density is required (Plewes and Hubbard, 2001).Thus, density introduces an uncertainty, as it also does for all conventional accumulation measurements that are based on thickness determination.
Studies based on GPR that analyse IRH in firn are typically complemented by cores to provide density, layer dating, and potentially, dielectric characteristics (Pälli et al., 2002;Arcone et al., 2004;Eisen et al., 2006).However, firn cores are expensive in terms of cost and effort.Furthermore, they only provide point measurements and little is known about the spatial variability of density (Lundberg et al., 2006).On the other hand, physical modelling of snow or firn density over time requires various input data such as temperature, precipitation, wind, and terrain to obtain the accumulation rate (Ligtenberg et al., 2011), whereas in mountainous terrain, the precipitation data alone can introduce a strong uncertainty (Lopes, 1996).
In this study, we combine 400 MHz helicopter-borne GPR measurements on Findelengletscher, a temperate alpine valley glacier in Switzerland, with a simple transient Introduction

Conclusions References
Tables Figures

Back Close
Full model for the bulk density of firn layers and the refreezing of meltwater.The model is calibrated at locations where GPR profiles intersect in subsequent years.For each individual GPR trace, this provides the firn layer water equivalents.Comparison with density and impurity profiles from shallow firn cores shows that not all IRH represent former summer surfaces.By taking this into account, we demonstrate that an estimate for the multi-year accumulation of an alpine glacier can be obtained from repeated airborne GPR if complemented with a firn densification model.

Study site and field data
Findelengletscher is a 13 km 2 temperate alpine valley glacier in the vicinity of Zermatt, Switzerland, ranging from 2600 to 3900 m a.s.l.It is situated directly below the main Alpine ridge in north-westerly exposure (Fig. 1).Thus, its accumulation characteristics are strongly determined by the synoptic weather conditions and inflow from the south and south-east.A monitoring program for annual mass balance was started in 2004 (Machguth, 2008) with a network of 13 ablation stakes and two snow pits for density measurements in the accumulation area (Fig. 1).Since 2009, the data record is complemented by measurements of winter mass balance using conventional snow probings and snow pits.The seasonal mass balance measurements are extrapolated to the entire glacier by constraining a distributed accumulation and temperature-index melt model operating on a daily resolution (Huss et al., 2009).The model is calibrated for each year individually using the point measurements of ablation and accumulation and takes into account spatial variations of potential solar radiation and snow accumulation.Meteorological data is available from an automatic weather station at Zermatt (1638 m a.s.l.) and at Stockhorn (3421 m a.s.l.).A comprehensive description of the model is given in Huss et al. (2009Huss et al. ( , 2010)).Since the onset of measurements, the average annual mass balance of Findelengletscher is negative which is accompanied by a retreat of the glacier snout (Glaciological Reports, 2011).Introduction

Conclusions References
Tables Figures

Back Close
Full Helicopter-borne constant-offset 400 MHz GPR measurements were taken for Findelengletscher on 10 May 2012 and 14 April 2013.The surveys were conducted along regular 500 m× 500 m grid lines covering the entire glacier and had a total profile length of 79 km.With a flying speed of approx.10 m s −1 , measurements were taken with a time increment of 0.02 s from 5-10 m above the surface.Together with the position at Differential Global Positioning System (DGPS) accuracy, a time window of 250 ns was recorded for each trace at a waveform sampling interval of 0.24 ns.
To obtain ground reference data, two shallow firn cores were recovered by mechanical drilling in the upper accumulation area of Findelengletscher on the same date as the GPR measurements in May 2012 (Fig. 1).The choice for the two sampling sites on a smooth ridge at 3495 and 3530 m a.s.l., respectively, was based on the availability of GPR profiles of preceding years at these locations and an expected slow glacier flow.The recovered cores were about 10.5 m long and 58 mm in diameter.

GPR data processing
A sequence of processing steps was applied to the GPR data in order to reduce noise and to improve the general visibility of reflectors.The selection of individual processing steps and their parameter settings depend on the field conditions and survey intention and, thus, no universal procedure is available (Annan, 1993;Fisher et al., 1996;Ulriksen, 1982).However, we applied a scheme that was previously used when mapping the snow accumulation distribution on Findelengletscher (Sold et al., 2013), but with some adjustments to the parameter sets.The scheme consists of a spatial interpolation of traces to match their average spacing to correct for variations in the helicopter's velocity; a frequency bandpass filter and background removal to reduce noise; and a gain function to account for signal attenuation with depth.The processing was performed with the software Reflexw (Sandmeier Scientific Software) and was found to Introduction

Conclusions References
Tables Figures

Back Close
Full sufficiently improve the data quality without the risk of over-processing and generation of artificial reflectors (Fig. 2).Additionally, the profiles were flattened to a horizontal snow surface that is represented by the first and most pronounced reflector.
The processed GPR data revealed several reflectors within the firn column for disjunct parts of the measured profiles in the accumulation area of Findelengletscher (Fig. 2).The penetration depth of the GPR signal depends on the height of the helicopter above ground and the physical properties of the subsurface.In the accumulation area we found IRH down to 20 m below the snow surface.They were digitised manually in the processing software due to noise and lateral variability in layer thickness which impede a reliable automatic tracking of IRH.
For the evaluation of GPR data from consecutive years, we aimed at locating the correspondent of each reflector at its new depth in the subsequent year.A matching procedure was applied on the traveltime-domain at locations where detected reflectors overlap (Fig. 1).Thereby, the uppermost IRH in the last year was omitted, representing the added past year's accumulation.A free scaling of ±20 % of the depth of the uppermost reflector in the earlier year was applied, which experiences the strongest compaction over the year.Then, the one-to-one mapping of reflectors was found with the minimum distance in the time-domain.

Modelling firn density
An assumption of the dielectric permittivity of the firn is required to obtain the radar wave velocity, converting the two-way travel time of the GPR signal from time-to depthdomain.Density strongly affects the dielectric properties of snow and firn and several empirical formulas exist to estimate the permittivity (Kovacs et al., 1993;Frolov and Macheret, 1999).In addition, the density is required to ultimately determine the water equivalent of firn layers.We used the approach to firn layer compaction by Reeh (2008) based on the empirical steady-state model of Herron and Langway (1980), which assumes density Introduction

Conclusions References
Tables Figures

Back Close
Full changes over time t to be linearly related to the pressure of overlying snow or firn: where ρ f is the density of the firn fraction of an annual layer, ρ ice is the density of ice and c is a parameter depending on the mean annual firn temperature T and the average mass balance b that controls the overlying weight.It is split apart by ρ f = 550 kg m −3 to respect different stages during the densification process (Herron and Langway, 1980). For where R is the gas constant and f is an emprirically derived constant controlling the densification rates.This steady-state approach was used by Reeh (2008) but with variable accumulation rates, initial densities and air temperatures.The density ρ f (t, t d ) of an annual firn layer t years after its deposition at t d is then given by solving Eq. ( 1) as an initial value problem with ρ f (0, t d ) being the layer's individual initial density.Thereby, instead of assuming steady-state conditions, an individual c is obtained for each layer from the mean annual mass balances b prev and mean annual temperatures of the t − 1 preceeding years since deposition.
In the context of our study, the local annual mass balances are initially unknown and the model cannot be applied as it is.However, along the GPR profiles, the traveltimethickness of firn layers and the snow cover are given.We used this information to progressively derive the unknown densities ρ f,i and water equivalents ω i of all i = 1, . . ., n layers at the time of the GPR measurements, i.e. t i years after each layer was deposited at t d .According to Reeh (2008), t d is attributed to the end of the respective melting season.Therefore, we set the age of the first annual firn layer (i = 1) as Introduction

Conclusions References
Tables Figures

Back Close
Full t 1 = 0.3 yr.This accounts not only for the GPR measurements taken in the following spring but also for the associated compaction processes which are stronger during the summer months due to higher temperatures and the refreezing of meltwater.For all subsequent layers we set t i +1 = t i + 1 yr.Based on several years of in-situ winter accumulation measurements, we assumed a density of 0.40 g cm −3 for the overlying snow cover.Its water equivalent ω 0 could then be extracted from the GPR traveltime with the radio wave velocity estimated from the density using the empirical relation by Frolov and Macheret (1999).Thus, an estimate for its pressure on the lower firn layers (i ≥ 1) was derived.For the subsequent firn layer the compaction rate c i +1 was then obtained by replacing the average mass balance b in Eq. ( 2) by the mean water equivalents of all overlying layers ω 0 , . . ., ω i .Additionally, we neglected the lower densification stage (ρ f < 550 kg m −3 ) that was suggested by Herron and Langway (1980), because the compaction rate will be calibrated and measured autumn snow densities were not considerably lower (see below).An individual compaction rate c i +1 for each firn layer was then derived from the water equivalent of the overlying layers as where the density of ice was set to ρ ice = 0.918 g cm −3 .The density of the next layer was then derived by solving Eq. ( 1) with the initial density ρ f (0, t d ) as boundary condition.This provided the GPR wave velocity and the water equivalent ω i +1 of that layer.Likewise, the densities and water equivalents of all firn layers were estimated by stepping through the layers i = 1, . . ., n and using the overlying water equivalent to derive its density t i years after deposition, i.e. its age at time of the GPR measurements.Although this approach does not take into account the internal layer densification due to its own weight, spatial variations in a layer's density are obtained through the weight of overlying snow cover and firn layers.Furthermore, our Introduction

Conclusions References
Tables Figures

Back Close
Full procedure neglects the variability of the loading history, affecting the rate of compaction with depth (Reeh, 2008).
The inputs to define are f controlling the compaction rate, the firn temperature T , and the initial density ρ f (0, t d ).In order to account for the temperate nature of Findelengletscher, temperature was set to T = 0 • C = 273.15K (constant).Although k is dependent on the firn temperature (Eq.2), any change to the constant value will be compensated by the calibration of f .On the other hand, a dynamic implementation of the heat transfer from air into the firn column would not only require a fully coupled model for refreezing, but also for melting at the surface and the percolation of meltwater.For the initial density we used ρ f (0, t d ) = 0.49 g cm −3 as the mean density measured within the conventional glaciological accumulation measurements at the end of the summer in 2010-2013.We obtained a conservative uncertainty estimate of ±0.061 g cm −3 from the standard deviation of (1) multiple density measurements within the same snowpits (measurement error), (2) the spatial variability within single years and (3) the temporal variability.The initial density is assumed to be spatially and temporarily constant.However, as the uppermost firn layer is of age t = 0.3 yr, it experienced compaction and thus, this assumption does not imply a spatially homogeneous layer density, but rather defines the initial (i.e. the preceding autumn) model condition.

Refreezing of meltwater
Due to cooling from the surface during winter, firn can periodically have temperatures below 0 • C even on temperate glaciers.Thus, the refreezing of meltwater percolating into the firn column has to be taken into account by the firn densification scheme, because the relocation and refreezing of water does not only alter the water equivalent of firn layers but also affects their density.We numerically modelled the one-dimensional firn temperature profile at each GPR measurement location using heat conduction from the surface.This is a reasonable approximation as no melting occurs in the accumulation area of Findelengletscher over the winter season.The thermal diffusitivity of snow was estimated from an empirical relationship with snow density (Calonne et al., 2011).Introduction

Conclusions References
Tables Figures

Back Close
Full To obtain a conservative estimate for the end-of-winter temperature profile, the respective modelled end-of-winter density profile was used.The winter snow cover was distributed linearly between 1 October and 1 May and its density was assumed to increase linearly from 0.1 to 0.4 g cm −3 .Due to the temperate conditions, we could assume that all negative temperatures are compensated by the refreezing of meltwater.We applied the model to every annual layer and all years within its individual history using the air temperature record from the automatic weather station at Stockhorn, corrected for elevation with a moist adiabatic lapse rate of −0.006K m −1 .In order to account for differences between the air and snow surface temperatures we applied a constant offset of −4.9 K as the observed long-term mean difference estimated from outgoing longwave radiation at the Stockhorn weather station during the winter months.The amount of refrozen meltwater was annually added to the density and water equivalent of the affected layers.However, the two models for temperature and density were not fully coupled because the densification model primarily steps through the layers from top to bottom, whereas the heat conduction was solved over time.Thus, each year's density profile remains static in the temperature model, whereas the layer water equivalents stay constant within the densification model but are adjusted separately.This is in line with the model itself neglecting variations in the loading history of each firn layer.

Model calibration
Aside from the refreezing, the change in density over time is controlled by the density contrast to the density of ice and the weight of overlying layers (Eq.4).The model sensitivity to those variables is controlled by the change rate c that is scaled by the factor f (Eq.3).Changes in f propagate through the entire model as they do not only affect the densities of individual layers, but also their water equivalent, their weight and, thus, the modelled densities of the following layers.Because the initial density ρ f (0, t d ) can be constrained with conventional measurements in autumn, we used the parameter f to calibrate the model to our site characteristics.Introduction

Conclusions References
Tables Figures

Back Close
Full At locations were GPR profiles intersected in the two consecutive years, the density change of each layer was modelled over the time of one year.In order to obtain the new overlying weight on each layer, the traveltime thickness of the snow cover and the new annual accumulation layer were considered.The resulting layer densities for the second year were used to estimate the respective hypothetical GPR signal traveltimes.
By comparing the modelled and measured traveltime thicknesses of 35 layers at 12 locations (Fig. 1), the optimal scaling factor f = 2900 within the change rate formulation c (Eq. 2) was found.However, the temporal density change also depends on the density contrast that is linked to the initial density.Therefore, the optimal value of f is related to the choice of ρ f (0, t d ) (Fig. 3).

Firn sampling, analysis and dating
The two 10.5 m firn cores recovered during the campaign in spring 2012 using the electromechanical drill FELICS small (Ginot et al., 2002) were stored in polyethylene tubes in sections of 0.5-0.7 m length at −20 • C.After a visual inspection for melt features or dust layers, the outer part was removed and the sections were cut with a resolution of ca.6 cm using a band-saw with stainless steel blades and Teflon-covered tabletops and saw guides to derive decontaminated samples of approx.60 mm × 17 mm × 17 mm dimension.They were melted in individual pre-cleaned 50 mL plastic tubes under inert gas (N 2 ) before the analysis.The density of each sample was determined gravimetrically from its measured volume.The concentration of the major soluble inorganic cations and anions including ammonium (NH + 4 ) and sulfate (SO 2− 4 ) was determined by ion chromatography (Metrohm 850 Professional IC combined with a 872 Extension Module and a 858 Professional Sample Processor autosampler).Additionally, the water stable isotope ratio (δD) was measured using Wavelength-Scanned Cavity Ring Down Spectrometry (Picarro L2130-i).
Dating of the two firn cores was performed by annual layer counting, using the pronounced seasonality of ammonium concentrations and water stable isotope ratios δD, both peaking in summer (Dansgaard, 1964;Eichler et al., 2000) (Fig. 4).A distinct dust Introduction

Conclusions References
Tables Figures

Back Close
Full

Local results at firn core locations
We observed a variable density profile for both firn cores that increased with depth (Fig. 4).Maximum densities corresponded to melt features of bubble free or poor ice.They are more transparent and form when surface snow melts and the meltwater percolates into deeper layers where it fills the pores and refreezes.In many cases melt features involve enhanced concentrations of NH + 4 and δD.Due to the temperate conditions at the drilling site, refreezing of meltwater contributes to the densification of firn and alters the water equivalents of firn layers, but also strongly influences the chemical and stable isotope composition (Ginot et al., 2010).The drilling took place in spring 2012 and thus, the snow accumulation of winter 2011/12 did not experience major melting.In contrast, the firn layers below experienced one or more cycles of melting and refreezing, leading to an enhanced smoothing of seasonal variations and potentially, a partial relocation of chemical species.Nevertheless, in agreement with studies of the nearby polythermal Grenzgletscher (Eichler et al., 2001), NH + 4 and δD seasonal cycles were still detectable, even under temperate conditions.
The measured density profile allowed for a direct estimation of the GPR signal velocity using the empirical relation by Frolov and Macheret (1999).On the traveltimedomain, we observed agreement between the GPR signal and pronounced changes in the density profile (Fig. 4), suggesting that the reflectors were composed of dense layers.This underlines the high potential of GPR for the detection of summer layers.However, not all density maxima and melt layers correspond to a former glacier summer Introduction

Conclusions References
Tables Figures

Back Close
Full surface, but can also originate from exceptional melting and refreezing in other seasons, or from percolation of melt water into deeper parts of the firn.In core 1 at a depth of 9.2 m (7.2 m in core 2) the accumulation layer of 2009 contained a pronounced ice layer of 0.12 m (0.15 m) thickness, corresponding to maximum measured densities of 0.86 g cm −3 (0.85 g cm −3 ) (Fig. 4).The low δD ratio implies low temperatures during the formation of the precipitation that was incorporated into the ice layer by in-situ densification or refreezing.We suppose that this layer was formed by surface meltwater that was hindered from further percolation by a blocking layer, such as a boundary with a change in the grain structure or a thin ice layer (Pfeffer and Humphrey, 1996).This is supported by a high concentration of most of the measured ions right below the ice layer (e.g.SO 2− 4 ) that typically just extends over one or two samples.By preventing meltwater from percolating and refreezing underneath the ice layer, the chemical composition is preserved and does not experience enhanced smoothing.After a period of major snow falls, high wind-speeds and mild temperatures were registered on 20 December 2008.It was followed by high radiation weather until mid January (Etter et al., 2011).We suggest that this allowed the formation of a dense layer at the surface of a fresh snow pack, creating a hydraulic barrier for meltwater formed later in the season.Because no defined clear ice layer was built up in the following summer of 2009 (Fig. 4), meltwater of several seasons could have penetrated down to, and been retained above the blocking layer, resulting in the thick ice layer with a considerable water equivalent within the firn layer of 2009 and a strong IRH in the GPR signal.

Spatial extraction of recent accumulation rates
The firn core profiles of density, NH + 4 concentrations, and δD ratios suggested that, aside from the ice lens inherent to the 2009 accumulation layer, the IRH can be attributed to dense summer layers.In fact, the related reflection in the GPR signal was observed across the entire accumulation area and was more pronounced at higher altitudes.The spatial continuity of the ice layer is further supported by its evidence in both firn cores.By disregarding the respective IRH, the water equivalents of annual Introduction

Conclusions References
Tables Figures

Back Close
Full accumulation layers were extracted.As described above, the firn densification model was applied to each GPR trace individually to estimate local layer densities based on the weight from overlying layers.Layer water equivalents were then derived from the IRH traveltime-thicknesses and the density-based GPR wave velocity estimate.
The water equivalents of annual firn layers derived from GPR exceeded those obtained by extrapolation of conventional mass balance measurements, albeit reproducing their temporal pattern (Fig. 5a).This can be attributed to the limited number of direct observations in the mass balance monitoring programme and the extrapolation scheme applied.Layer water equivalents derived from the firn cores revealed similar discrepancies with extrapolated mass balance.These findings highlight the need for a further investigation and better representation of accumulation processes in current glacier mass balance monitoring.
The spatial distribution of each annual layer's water equivalent showed a considerable small-and medium-scale variability, following the general topography to some degree (Fig. 6).The large-scale pattern of areas with increased annual accumulation rates recurs in all four layers.However, with increasing depth, older IRHs were more difficult to locate in the radargrams and appear in shorter profiles sections, thus hampering a direct temporal comparison of the accumulation patterns.Modelled layer densities increased with age as expected (Fig. 5b) reaching a mean value of 0.74 g cm −3 four years after deposition.The cumulative effect of spatial variations in layer water equivalent, together with the fixed initial density in the densification model, results in a higher variability of the density of deeper, i.e. older, firn layers.Refreezing is a considerable contributor to firn densification on alpine glaciers.It accounted for approx.9 % of the density increase over time and depth (Fig. 5b) and corresponds to a mean of 28 kg m −2 yr −1 over the entire firn column in the accumulation area of Findelengletscher.The conservative uncertainty estimate for the initial model density that covers measurement uncertainty, spatial and temporal variability, and the related uncertainty introduced within the model calibration affects the modelled layer densities and thus, the layer water equivalents by 7 % (Fig. 5a).A conventional glaciological point accumu-

TCD Introduction Conclusions References
Tables Figures

Back Close
Full lation measurement is only slightly better (±0.021 g cm −3 density measurement error), but lacks any information about the spatial distribution.By combining GPR with the model for firn density a considerable spatial coverage for four annual accumulation layers can be obtained with a small trade-off in terms of uncertainty.

Discussion
Generally, the formation of ice layers is not limited to surface processes.At boundaries with changes in the grain structure or on top of a blocking layer, e.g. a thin ice layer, meltwater can be retained (e.g., Pfeffer and Humphrey, 1996).Thus, an outstanding weather event during winter is not a necessary criterion to generate an ice lens that can potentially be misinterpreted as a former summer surface in the GPR signal.In the case of Findelengletscher, such a layer could be identified within the 2009 accumulation with the help of reference data from firn cores.Ultimately, the presented method of deriving past accumulation rates from GPR is based on the interpretation of GPR profiles by an observer, with the IRH relating to high-density layers (Fig. 4).Therefore, the interpretation of GPR data can be ranked similar to a visual and density-based interpretation of a firn profile.Conventional glaciological measurements, comprising snow probings and snow pits for density measurements, are likewise prone to misinterpretations of the stratigraphy.This can be either visually in a snow pit, by hitting an internal ice layer, or by penetrating the designated summer surface during probing.In contrast to the presented approach, conventional accumulation measurements are conducted on an annual basis, cover individual years and thus, do not require interpretation of multiple layers.On the other hand, a misinterpretation of an IRH affects not only the estimated water equivalent of the respective accumulation layer, but in fact, of all following layers, by shifting the chronology of the layers.However, this source of error is inherent in all approaches that rely on annual layer counting.An independent verification of the annual character of IRH requires either in-situ measurements from firn cores or a model that can simulate

Conclusions References
Tables Figures

Back Close
Full the relevant small-scale processes.This comprises solving the surface energy balance, heat transfer, water percolation, and formation of internal hydraulic barriers such as ice layers (Mitterer et al., 2011).Modelled or independently measured mass balance data could be used for a plausibility check (e.g., Van Pelt et al., 2014).We showed that the weak representation of the spatial distribution of accumulation is a drawback in current mass balance extrapolation schemes (Fig. 5a) and thus, such a comparison requires a cautious interpretation.In our case, the GPR-derived layer water equivalents exceed the extrapolated mass balance measurements while reproducing the temporal pattern.Thus, relative deviations could still be used for a first comparison.By relying only on the firn cores, we avoid an unstable verification approach and remain independent from measurements taken in years before the GPR survey and the extrapolation of sparse measurements.
The firn densification model was calibrated with the observed change in traveltime from spring 2012 to 2013 at 12 locations and for 35 firn layers.The fitting procedure is not exceptionally stable because an adjustment of the change rate c with the scale factor f (Eq. 3) affects the densities of layers, their water equivalent and thus, their pressure on the underlying firn mass (Eq.4).On the other hand, the density change is linearly dependent on the contrast to the density of ice (Eq.2) and the two effects counteract each other in the calibration scheme.In order to obtain a reliable change rate in the model we did not run the calibration on individual locations to account for spatial variations in the densification regime.Rather, the derived parameter set was assumed to represent the mean temporal and spatial condition for the investigated site.Comparison with the densities and layer water equivalents derived from the firn cores showed that the local deviations are within the uncertainty range provided by the initial density (Fig. 7a).When hypothetical GPR traveltimes of the layers are calculated from the density and water equivalent of the firn core, the model revealed an adequate density for the topmost layer.For the following layers of 2010 and older, densities were too high (Fig. 7b), accounting for 9 % of their water equivalent.This indicates that for this individual location the given initial density is reasonable, but the change rate is Introduction

Conclusions References
Tables Figures

Back Close
Full too high.We suggest that this discrepancy is due to spatial variations in the density regime, related to processes not incorporated in the model, such as the snow grain structure, the temperature-dependence of densification processes, or the particular weather conditions in general.
Throughout the entire study we refer to the water equivalent of firn layers in the sense of annual accumulation rates.We modelled refreezing melt water within the firn each year by using a daily time series of temperature data to calculate the temperature profile of the firn column that is entirely compensated by refreezing.However, this approach only takes into account the effect on bulk layer density and water equivalent but does not shift mass between layers.This is in line with conventional glaciological measurements that only cover the uppermost accumulation layer and thus, similarly miss the additional accumulation that is generated by refreezing in firn layers of previous years deeper in the firn column.As our GPR measurements were conducted in spring 2012, however, the winter accumulation of 2011/12 did not experience major melting and therefore, no external refreezing is expected to have occurred in the topmost firn layer representing the accumulation of 2011.In order to consistently quantify vertical mass fluxes in the firn throughout the year, a coupled energy-balance -snow model can be complemented with GPR data.Van Pelt et al. (2014) calibrated the accumulation routine of such a model by matching modelled and observed IRH traveltimes.This provides the adjusted accumulation distribution, but requires numerous input data in order to solve the surface energy balance.
The water equivalent stored as firn between IRH also differs from the accumulation part within the mass balance term regarding the temporal breaks between years.Mass balance monitoring approaches typically use the hydrologic year, with 30 September marking the end of the ablation period on the Northern Hemisphere.In contrast, the analysis of the GPR data is based on the firn stratigraphy where the build-up of a highdensity layer due to particular weather conditions provides the chronology of firn layers.Thus, a comparison or combination of such data sets requires conversions between the Introduction

Conclusions References
Tables Figures

Back Close
Full stratigraphic and the fixed-date mass balance (Mayo et al., 1972;Huss et al., 2009), but does not generally conflict with the mass balance concept.

Conclusions
We established a new approach to derive past accumulation rates on temperate alpine valley glaciers from repeated GPR surveys.Our method is based on the interpretation of IRH as previous summer surfaces.By linking the measured GPR traveltime of each layer to a simple model for firn densification, an estimate of the bulk density was derived for each layer.Thereby, we take into account the refreezing of meltwater by modelling the end-of-winter temperature profile that is entirely compensated by refreezing under temperate conditions.The model was calibrated by comparing the modelled and measured change in traveltime thickness of 35 layers at locations where GPR profiles intersect in two subsequent years.Thus, our approach is independent from external information such as firn cores.It was applied to each GPR measurement location independently to obtain the distribution of the mass of annual firn layers along the GPR profiles.
The strongest limitation to the methodology is the non-annual character of IRH for our study site as it was found from the analysis of two shallow firn cores.A potential misinterpretation of IRH affects the results for all subsequent layers by changing their ages.We suggest that, if no firn cores are available, measured mass balance data or modelled accumulation rates based on meteorological information could serve for a plausibility check.However, the formation of high-density and ice layers is linked to smallscale physical characteristics of snow and firn and is driven by particular weather conditions.Thus, a general conclusion about the validity of a dating-by-counting scheme for different study sites cannot be drawn.
We showed that GPR measurements can be used to estimate multi-year accumulation rates on a temperate valley glacier.Thereby, measurements can be conducted from helicopter, facilitating efficient data acquisition even in inaccessible areas.Be-Introduction

Conclusions References
Tables Figures

Back Close
Full cause the GPR penetrates several layers of accumulation, it could be used to complement existing monitoring programs, but also for a retrospective extension of newly initiated time series.Introduction

Conclusions References
Tables Figures

TCD Introduction
Full . Conventional mass balance measurements are prone to large uncertainties when compared to geodetically derived volume change.The reconciliation Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Screen / Esc Printer-friendly Version Interactive Discussion Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | layer from the Saharan dust fall on 28 May 2008(Hostettler and Bader, 2008) found in 10.2 m (Core 1) and 8.5 m (Core 2) depth was used to corroborate the dating.Core 1 and 2 cover the time periods summer 2008 to spring 2012 and summer 2007 to spring 2012, respectively.

Figure 1 .Figure 2 .
Figure1.Findelengletscher with the setup of annual mass balance measurements, the mean equilibrium line altitude (ELA), the firn core locations and all analysed GPR profiles measured in 2012, as well the locations used for the model calibration.The overview map shows its location within the European Alps.

Figure 7 .
Figure 7. Layer water equivalents (a) and densities (b) at firn core 1, derived from GPR, directly from the firn core, and modelled using the hypothetical traveltime obtained from the firn core layer water equivalent and density.Error bars show the uncertainty derived from variations in the initial density model parameter ρ f (0, t d ) and the linked scaling for the compaction rate f .