Lidar snow cover studies on glaciers in the Ötztal Alps ( Austria ) : comparison with snow depths calculated from GPR measurements

The storage of water within the seasonal snow cover is a substantial source of runoff in high mountain catchments. Information about the spatial distribution of snow accumulation is necessary for calibration and validation of hydro-meteorological models. Generally, only a small number of precipitation measurements deliver precipitation input for modelling in mountain areas. The spatial interpolation and extrapolation of measurements of precipitation is still difficult. Multi-temporal application of lidar techniques from aircraft, so-called airborne laser scanning (ALS), provides surface elevations changes even in inaccessible terrain. These ALS surface elevation changes can be used to derive changes in snow depths of the mountain snow cover for seasonal or subseasonal time periods. However, since glacier surfaces are not static over time, ablation, densification of snow, densification of firn and ice flow contribute to surface elevation changes. ALS-derived surface elevation changes were compared to snow depths derived from 35.4 km of ground penetrating radar (GPR) profiles on four glaciers. With this combination of two different data acquisitions, it is possible to evaluate the effect of the summation of these processes on ALS-derived snow depth maps in the high alpine region of the Ötztal Alps (Austria). A Landsat 5 Thematic Mapper image was used to distinguish between snow covered area and bare ice areas of the glaciers at the end of the ablation season. In typical accumulation areas, ALS surface elevation changes differ from snow depths calculated from GPR measurements by −0.4 m on average with a mean standard deviation of 0.34 m. Differences between ALS surface elevation changes and GPR derived snow depths are small along the profiles conducted in areas of bare ice. In these areas, the mean absolute difference of ALS surface elevation changes and GPR snow depths is 0.004 m with a standard deviation of 0.27 m. This study presents a systematic approach to analyze deviations from ALS generated snow depth maps to ground truth measurements on four different glaciers. We could show that ALS can be an important and reliable data source for the spatial distribution of snow depths for most parts of the here inv stigated glaciers. However, within accumulation areas, just utilizing ALS data may lead to systematic underestimation of total snow depth distribution.


Introduction
In alpine catchments, the so-called glacio-nival runoff regime is caused by the storage of water in the seasonal snow cover and in glaciers (Aschwanden et al., 1986;Kuhn and Batlogg, 1998;Kuhn, 2000).High flow rates caused by snow-and ice melt in spring and summer are of interest in terms of water resources management for reservoirs and flood forecasting (e.g.Weingartner et al., 2003;Kirnbauer et al., 2009;Schöber et al., 2010).
For calculating reliable runoff amounts in these mountain catchments, precipitation data in high elevations are needed and the heterogeneous distribution of the snow cover has to Published by Copernicus Publications on behalf of the European Geosciences Union.
be considered in hydro-meteorological models (e.g.Kuhn, 2003;Huss et al., 2008).However, measurements of solid precipitation are known to be affected by considerable errors (e.g.Sevruk, 1985;Lundberg et al., 2010) and only insufficient data are available for model calibration and validation of the snow depth distribution on the catchment scale.Thereby the spatial representativeness of individual measurements of snow depths and snow densities derived from direct measurements in terms of snow probings and snow pits (Fierz et al., 2009) and from automatic measurement systems (Lundberg et al., 2010) has to be questioned (Grünewald and Lehning, 2011).
Whereas various hydro-meteorological models have been tested versus snow cover extent from satellite remote sensing data of MODIS (Moderate Resolution Imaging Spectroradiometer) or Landsat scenes (e.g.Rott and Markl, 1989;Hall et al., 2002;Schöber et al., 2010), to date no snow depth or SWE can be obtained from satellite data in a high spatial resolution of at least 100 m (Lundberg et al., 2010;Nolin, 2011), which is required to capture the spatial variability of the mountain snow cover (Blöschl, 1999).In contrast, laser altimetry was recognized to be useful to derive small scale snow depth distributions in mountain areas (Deems et al., 2013).Several snow cover studies based on the application of lidar from aircraft (airborne laser scanning -ALS) have been performed in mountain areas, mostly in unglacierized terrain (e.g.Deems et al., 2008;Melvold and Skaugen, 2013), but also in glacierized catchments (e.g.Dadic et al., 2010;Helfricht et al., 2012;Sold et al., 2013).Since 2001, ALS surveys of annual and seasonal surface elevation changes ( z ALS ) were performed in glacierized catchment of the Ötztal Alps (Austria, 46 • 50 N, 10 • 50 E) (Geist and Stötter, 2007;Helfricht et al., 2012).In 2010/11 an ALS survey of the whole mountain range of the Ötztal Alps was conducted to estimate seasonal snow depths.
In this study z ALS are investigated for the accumulation season ranging roughly from the beginning of October (t 1 ) to the end of April respective the beginning of May (t 2 ). Figure 1 illustrates the processes causing surface elevation changes on glaciers in this period.The elevation z 1 at the time t 2 corresponds to the primary surface z 1 at the time t 1 .The elevation z 2 corresponds to the actual surface at the time t 2 .z ALS is the result of changes in snow depth ( z HS ), densification of snow and firn ( z DSF ), ablation ( z ABL ) and vertical ice flow ( z ICE ) whereas z HS at the time t 2 is defined as the elevation difference between z 2 and z 1 and the summation of densification, ablation and vertical ice flow is the elevation difference between z 1 and z 1 z DSF + z ABL + z ICE = z 1 − z 1 . (3)

Fig. 1.
Changes in surface elevation due to snow depth ( z HS ), ablation ( z ABL ) and densification of snow and LSC ( z DSF ) resulting in ALS surface elevation changes z ALS between the two ALS flights (t 1 ) and (t 2 ).Upper part: firn and snow on top of static ice.Lower part: firn and snow on top of a vertically moving ice surface ( z ICE ).
The densification of snow and firn layers existing at t 1 leads to an underestimation of z HS in the case of snow and firn layers on a static ice body at t 2 .Changes of the reference elevation due to vertical ice flow can have both, negative and positive values.Submergence ice flow ( z ICE < 0) leads to an underestimation of z HS by z ALS .Emergent ice flow ( z ICE > 0) causes overestimation of z HS by z ALS .Hence, "ground truth" data are necessary to evaluate z ALS for snow cover and mass balance studies on glacier surfaces (e.g.Huss et al., 2009;Koblet et al., 2010;Fischer, 2011).
The aim of the present study is to show the consequence of differences between z ALS and actual snow depths on the application of multi-temporal ALS to derive seasonal snow depths on the investigated glaciers.The questions is whether and where differences between z ALS and actual snow depths on glacier surfaces are of magnitude, such that they have to be considered in snow-hydrological studies based on ALS surveys.Therefore, we compared z ALS to snow depths calculated from ground penetrating radar (GPR) measurements and analyzed differences larger than the combined uncertainties of ALS and GPR measurements.The presented results will support the application of ALS surveys for snow cover studies on glaciers in this mountain area and in Alpine catchments with comparable glacier behaviour in terms of ice flow velocities and accumulation area extent.Systematic underestimation or overestimation of snow depths by z ALS can be attributed to specific areas for all here investigated glaciers.
After the introduction an overview of methods and data used to process, the GPR and lidar data is given in Sect. 2. Conditions on the glaciers during ALS surveys are treated in more detail in Sect.3. In Sect. 4 the results are presented followed by the discussion in Sect. 5. Finally, conclusions are made about the applicability of surface elevation changes derived from ALS for snow cover studies in glacierized catchments.

Study sites
The study sites are located next to the Alpine main ridge in the Ötztal Alps (46 • 48 N, 10 • 46 E; Fig. 2).Four of the largest glaciers in this mountain range, namely Hintereisferner, Kesselwandferner, Gepatschferner and Vernagtferner were chosen to evaluate z ALS on glacier surface for snow cover studies.All of these glaciers have been subject to scientific research for decades: Hintereisferner (HEF, 46 • 48 N, 10 • 46 E) has one of the world's longest mass balance series starting in 1952 (Hoinkes, 1970;Fischer, 2010;Fischer et al., 2012).Based on these data sets, a variety of models with different complexity were developed to compute glacier surface mass balance (Kuhn et al., 1999;Escher-Vetter et al., 2009).Ice thickness was measured by Span et al. (2005).Since 2001 a series of ALS surveys enabled the comparison of geodetic and direct glaciological measurements of mass balance on Hintereisferner (Geist and Stötter, 2007;Fischer, 2011).However, since the last acceleration of ice flow on Hintereisferner in the 1970s, horizontal velocities have slowed down dramatically to less than 10 m a −1 (Span et al., 1997;Span and Kuhn, 2003).They found that emergence velocities measured on the glacier tongue of Hintereisferner decreased to about 1 m a −1 at the beginning of the 1990s.Haberkorn (2011) observed a depression of the ice surface for a certain area along the tongue of Hintereisferner.Measured horizontal velocities at ablation stakes were 4.6 m a −1 at most.Vertical velocities could not be derived due to small absolute values with changing signs.In 2010/2011 ablation stakes along the flow line on the glacier tongue of Hintereisferner were georeferenced with DGPS (H. Schneider, personal communication, 2013).In combination with the ablation measurements at these positions, vertical velocities were calculated.The mean annual vertical velocities increased from zero in the area of the GPR profile H3 (Fig. 8) to 0.55 m a −1 of annual emergence velocity at about 2 km from the terminus.Annual vertical ice flow measured at an ablation stake in the central flow line near the GPR profile H5 (Fig. 8) was 0.85 m a −1 in 2010/2011.Kesselwandferner (KWF, 46 • 50 N, 10 • 48 E) was subject to scientific studies investigating density profiles and deformation of firn layers (Ambach and Eisner, 1966).Mass balance measurements from Kesselwandferner using the direct glaciological method started in the hydrological year 1952/53 (Fischer et al., 2011) and geodetic mass balances are available since 1964.In Abermann et al. (2007) annual vertical and horizontal velocities at several stake locations along the central flow line are presented.After the glacier's advance in the 1980s, vertical velocities decreased to 1.5 m a −1 at most, both for submergence and emergent flow.Since the measurements of Abermann et al. (2007), annual surveys of the ablation stakes and the accumulation stakes have been continued.For the last 5 yr, mean annual vertical velocities along the flow line remained similar at up to 1 m a −1 (H.Schneider, personal communication, 2013).
Gepatschferner (GF, 46 • 51 N, 10 • 45 E) is the largest glacier in this region.Gepatschferner is the only glacier with an area of more than 10 km 2 in the Ötztal Alps.Its changes in volume and area in the last decades are summarized by Abermann et al. (2009) and ice thickness was measured by Fischer et al. (2007).Glacier velocities have not been investigated on Gepatschferner.Recently a network of ablation stakes was set up to investigate the dynamical behaviour of the glacier tongue.Highest velocities are expected to occur at the crevassed area before the inflow to the glacier tongue.First results show maximum horizontal velocities of up to 50 m a −1 with a vertical ice flow of approx.5 m a −1 measured below the ice fall.Velocities decrease towards the glacier tongue, which itself decays very rapidly.No velocity data exist from the large plateau of the accumulation area.
Since 1964 the mass balance of Vernagtferner (VF, 46 • 52 N, 10 • 49 E) is measured with the direct glaciological method by the Commission for Glaciology, Bavarian Academy of Sciences and Humanities, Munich.In combination with a runoff gauge in front of the glacier tongue, a series of hydro-meteorological models were calibrated for energy balance and hydrological studies (Escher-Vetter et al., 2009).Mayer et al. (2013) showed that the glacier slowed down dramatically.
Apart from the long time series of mass balance and ice flow velocity measurements on these glaciers, no measurements of ice flow in the accumulation season exist.With respect to recently observed slow ice velocities in this region, the influence of ice flow on seasonal surface elevation changes can be assumed to be small.However, even 1 m of submergence would have an enormous effect for the interpretation of z ALS in terms of snow depth without applying corrections for ice dynamics.
For this study, glacier outlines of the four investigated glaciers derived from the Austrian glacier inventory of 2006 ( Abermann et al., 2009) were updated based on ALS data of 2010 according to Abermann et al. (2010).

Ground penetrating radar
Ground penetrating radar (GPR) is an active sensing technique to locate targets or internal layers within ground, ice and snow (Daniels, 2004).A transmitted electromagnetic wave is partly reflected at boundaries between two layers with different relative dielectric permittivities ( r ).To convert measured two-way-travel times (TWT) of the GPR signal to depth of the reflecting layer transition, the velocity of propagation of the wave has to be determined.In cryospheric sciences, GPR is used to detect permafrost within rocky ground (e.g.Hinkel et al., 2001;Otto et al., 2012), to determine ice thickness (e.g.Span et al., 2005;Fischer et al., 2007), ice and firn layering (e.g.Spikes et al., 2004) as well as subglacial structures and liquid water content of snow and ice (e.g.Murray et al., 1997;Lundberg et al., 2000).The spatial variability in snow depth and the temporal evolution in snowpack stratigraphy have been analyzed utilizing GPR technology as well (e.g.Machguth et al., 2006;Heilig et al., 2009;Mitterer et al., 2011).

Field campaigns
Three GPR campaigns were conducted close to the respective date of ALS surveys.Dates and snow pit information are displayed in Table 1.
At the end of accumulation season 2010/2011, GPR data were recorded on Vernagtferner by the Commission of Glaciology, Bavarian Academy of Sciences and Humanities (Munich, Germany).A RIS One GPR instrument from IDS (Pisa, Italy) with shielded 600 MHz antennas was used.The Table 1.Attributes of the field measurements in spring.Number of snow depth probings and snow pits shown in brackets indicate measurements not directly located at the GPR profiles.The GPR signal velocities derived from snow depths of the reflecting snow layer in the snow pits (v p ) and GPR signal velocities calculated from measured snow densities (v k ) using Eqs.( 7) and ( 8) with the corresponding coefficient of variation (CV, Eq. 6).GPR was pulled by a snowmobile.At the first and the last point of the sections snow depths were measured by snow depth probing and coordinates were recorded with a handheld GPS.Between these points the GPR data were continuously received while driving the snowmobile at uniform speed, which, however, differ between single sections due to topography.In total more than 15 km of GPR profiles were conducted, divided into four longitudinal profiles, three additional shorter profiles in the accumulation zone and one cross profile (Fig. 6).Snow densities were measured at six locations.

ID
A GPR measurement campaign was conducted on Gepatschferner and Kesselwandferner in 2011.GPR data were recorded with a 3-D system from MALÅ (Malå, Sweden) with shielded 500 MHz antennas along a cross profile in the accumulation zone of Gepatschferner and a longitudinal profile on Kesselwandferner (Fig. 7).In total more than seven kilometres of profiles, divided into 33 sections, were investigated on skis.First and last points of the individual sections were georeferenced by DGPS (Fig. 7 -black dots).Three snow pits were dug on Gepatschferner and one on the tongue of Kesselwandferner.While snow depth probings were performed along the GPR sections on Gepatschferner, snow depth probings and GPR data do not spatially coincide on Kesselwandferner.Sample spacing is larger on Kesselwandferner than on Gepatschferner, due to greater velocities while skiing downhill the glacier (Table 1).
At the end of the accumulation season 2011/2012, a GPR campaign was conducted to measure snow depths along longitudinal sections and cross sections on Hintereisferner and Kesselwandferner (Figs. 8 and 9).Skis were used to pull another IDS RIS One system with shielded 400 MHz antennas.In total over 12 km of GPR data were recorded.Snow depth probings were conducted at every first and last point of the individual sections and seven snow pits were dug.All data points were georeferenced by DGPS.
For reflectors being located at a distance of half the wavelength in the respective medium, constructive and destructive interferences will occur (Daniels, 2004).The applied GPR systems with a centre frequency of 400-600 MHz, hence are capable of distinguishing various reflectors in a vertical distance of more than 0.19 m (600 MHz) to 0.28 m (400 MHz) (Trabant, 1984).All GPR campaigns were conducted in time mode.The trace spacing for each transect was chosen such that 2 consecutive traces are always separated less than 0.5 m.So quasi-continuous information could be tracked along the profiles based on known reflector depths from snow depth probing and from snow pits.

Processing
GPR raw data were processed applying the ReflexW software from Sandmeier Scientific Software (Karlsruhe, Germany).Coordinates for each trace were calculated at equal distances.The surface signal reflection was set to time zero.Low frequency parts and noise in the spectrum were filtered applying a dewow and bandpass filter.In a next step temporally consistent signals were eliminated utilising a background removal.To compensate for spherical divergence losses, we applied gain functions.No migration was used due to the spatial homogeneity and the low reflector depth relative to the high density of GPR measurements.The boundaries of the snow layers were picked and corrected to the zero phase change.The picks were exported with the attribute of the two-way-travel time.Finally, we merged individual sections to continuous longitudinal and cross profiles.
Snow depth probing and snow pits were used to identify the depth of the seasonal snow layer according to the guidelines for measuring glacier mass balance (Kaser et al., 2003).GPR signal velocities were calculated from measured reflector depths and snow densities at the snow pits only.Assuming that the snow surface and the ice surface are parallel, the slope of the surface has to be considered.The recorded GPR signal corresponds to the snow depths perpendicular to the surface, but snow depths derived from snow probings and snow pits (h m ) are measured vertically.We corrected vertically measured snow depths h m by in accordance to the respective slope angle α.Signal velocities (v) for GPR measurements on Gepatschferner and  Kesselwandferner in 2011 and GPR measurements on Hintereisferner and Kesselwandferner in 2012 were calculated following Daniels (2004): using h cor and the TWT (τ ) at the snow pit locations.A mean signal velocity (v) was calculated for each measurement campaign individually (Table 1).The variation of the signal velocities (CV) is given by as the ratio of the range between the derived maximum velocity v max and minimum velocity v min to the mean velocity v.According to Kovacs et al. (1995), we calculated r of the snow cover where ρ s is the snow density in g cm −1 .Based on r GPR signal velocities v were calculated where c is the velocity of propagation of the electromagnetic wave in a vacuum of approx.0.3 m ns −1 .On Vernagtferner snow pits were not located along the GPR profiles.Therefore GPR signal velocities were calculated from measured snow densities only following Eqs.( 7) and ( 8).Signal velocities derived from measured snow depths in the snow pits (Eq.5) on Hintereisferner and Kesselwandferner were validated with signal velocities calculated from snow densities (Eq.8) (Table 1).The CV of the signal velocities is larger for 2012 data due to higher variations in snow densities measured in the snow pits.Snow depth probings performed along the GPR profiles were used to validate the respective calculated wave speed.All exported TWT τ were converted into snow depths (h v ) using the mean velocities v for each measurement campaign.Again to account for vertically measured z ALS , GPR-determined h v had to be corrected for the prevailing slope α resulting in the final snow depth h GPR Figure 5 shows the magnitude of corrections which have to be applied to h v as a function of snow depth and slope.In the same figure the frequency distribution of snow depths measured with GPR and the frequency distribution of the occurring slopes are shown.The most frequent snow depths were between 1.55 m and 1.65 m.The most frequent slopes were between 5 • and 7 • .Corrections, which therefore had to be applied, were mainly of the order of 0.01 m.

Lidar
Lidar is an active remote sensing technology (Baltsavias, 1999;Wehr and Lohr, 1999;Kraus, 2004).Airborne laser scanning (ALS) can be applied in mountain areas on the catchment scale, because it is hardly affected by topographic shading and travel distance of the signal can be controlled by flight height (e.g.Deems et al., 2013).The DGPS referenced data provide 3-D information of the reflecting surface with a high point density and accuracy.Based on these 3-D point data (the so-called point cloud) digital elevation models (DEMs) can be produced (Lui, 2008), while high point densities allow simple interpolation schemes (Deems et al., 2013).Surface elevation changes can be detected from the difference of the DEMs of multi-temporal laser scan surveys.Accuracy of the 3-D location of each point is affected by the accuracy of the recording of the position of the scanning platform and its orientation (e.g.Joerg et al., 2012;Deems et al., 2013).Accuracy of the resulting DEMs depends on point density, the adjustment of overlapping areas of the stripes and the interpolation algorithm used (e.g.Deems et al., 2013)

Data and processing
In this study four laser scan surveys were used.They were conducted by TopScan using Optech devices mounted on Cesna aircrafts.Technical details and the accuracy measured at a reference surface are listed in Table 2.The last pulses of all ALS points within one pixel of the resulting DEMs were averaged to a grid size of 1 m.Data were analyzed on generally homogeneous glacier surfaces, which, except for very small areas, have a slope of less than 20 • along the investigated profiles (Fig. 5).Surface elevation changes ( z ALS ) were calculated by subtracting a DEM of an almost snowfree ALS survey in fall (t 1 ) from a DEM of an ALS survey at the end of accumulation season in spring (t 2 ), respectively, following Eq.( 1), with z 1 = z t1 and z 2 = z t2 (Fig. 1).

Calculation of differences and uncertainties
Absolute differences ( h abs ) were calculated between ALS surface elevation changes ( z ALS ) and snow depths calculated from GPR (h GPR ) by A moving window approach was used to perform a statistical Mann-Whitney-Wilcoxon test with a number of 300 samples of h abs around each single data point of GPR measurements.According to the mean trace distances shown in Table 1, the number of 300 samples corresponds to a distance of approx.100 m, which is an adequate spatial scale for differences analyzed here.The null hypothesis was tested, if the median of the two samples ( z ALS and h GPR ) is equal with a significance level of p = 0.05.If the null hypothesis is rejected, the medians of the two data sets differ from each other.Additionally the mean value of h abs of the 300 samples was compared to an assessment of absolute uncertainty U abs calculated as follows Uncertainties of gridded lidar data (U ALS = ±0.15m) have to be considered twice, because two elevation models are used to investigate surface elevation changes.Uncertainties caused by the calculated GPR signal velocities (U v ) correspond to half of maximum CV (Eq.6, Table 1), which is ±0.05.U abs are also a function of snow depths (h GPR ).
The results of the statistical analysis show where z ALS and h GPR differ significantly.Furthermore it was calculated where h abs exceed U abs along the profiles.In these areas, z ALS is increasingly influenced by processes more than only snow depth uncertainties (see Sect. 1).

Accumulation areas from Landsat
The densification of firn is supposed to be one source of differences of z ALS from actual snow depths (Sect.1).Also submergence ice flow is very likely in the central parts of typical accumulation areas.Hence, typical accumulation areas on glaciers were expected to show larger differences since both vertical ice flow and densification work vertically in the same direction.Snow covered areas of glaciers can be delineated by application of spaceborne optical data for The Cryosphere, 8, 41-57, 2014

Prevailing local conditions during lidar surveys
Surface elevation changes between field campaigns and the dates of ALS surveys (t 1 , t 2 ) contribute to uncertainties in h abs .Snow accumulation before the ALS surveys in the higher elevated glacier areas in fall has to be considered when identifying the snow layers in the GPR analysis.Photographs and measurements of snow depths and snow densities taken during annual glacier mass balance field campaigns on Hintereisferner, Kesselwandferner and Vernagtferner were included for data analysis.Information on the weather conditions before the ALS surveys, between the ALS surveys and the field campaigns as well as on the weather conditions during the accumulation seasons were obtained from two automatic weather stations located next to Hintereisferner (46 • 47 55 N,10 • 45 37 E; 3027 m a.s.l.) and next to Vernagtferner (46 • 51 23 N, 10 • 49 43 E; 2640 m a.s.l., Fig. 2).The course of daily mean temperatures and cumulative precipitation at the two automatic weather stations is presented in Figs. 3 and 4.
On 12 October 2010 snow depths and snow densities of a snow layer covering the 2010 firn layer on Gepatschferner and Kesselwandferner was measured by 43 snow depth prob-ings and three snow pits.This snow layer was accumulated in September 2010.In a warm period before the ALS survey (Fig. 3), this layer was transformed into a snowpack with a high bulk density and a melt-freeze crust at the surface.This snow layer could be identified in the stratigraphies of spring snow pits and in the snow layering derived from GPR data in 2011.Snow depths of the same snow layer were recorded by 34 snow probings within the annual mass balance measurements at Vernagtferner. Between the ALS survey in spring 2011 and the corresponding GPR field campaigns a shallow snow cover of about 0.1 m with a low density was accumulated.
In September 2011 glacier tongues and steep southerly exposed parts on Hintereisferner and Kesselwandferner appeared to be snow free on photographs taken alongside the annual mass balance measurements on 1 October 2011.Only a shallow snow layer was covering the typical accumulation areas (A.Fischer, personal communication, 2012).Shortly after the ALS survey in October 2011 a snowfall event preserved the surface recorded by ALS from melt.

Results
The spatial distribution of z ALS and snow depths calculated from GPR measurements (h GPR ) along the profiles are illustrated in Figs.6-9.Locations of snow pits and snow depth probings of the spring field campaigns are included.Data corresponding to the individual profiles are listed in Table 3.
Figure 10 shows the total number of analyzed data points for LSC areas and bare ice areas (Sect.2.5) from the profiles and the statistical distribution of all absolute differences ( h abs , Eq. 11) in terms of the median (red line), the scatter (25 % and 75 % percentiles, blue box), the range of extreme values (1.5 times the interquartile range added to the 25 % percentile and to the 75 % percentile, black whiskers) and the outliers (red crosses).
The profiles H1 (Fig. 11) and K4 (Fig. 12) along the central flow lines of the glaciers show distinct negative h abs (Eq.11) in the LSC areas.h abs are constantly small along most of the glacier tongue of Hintereisferner.At the very front of the glacier tongue, h abs increases towards positive values.Along the profile K4 on KWF distinctly less h GPR match with the corresponding z ALS .Towards the tongue of the glacier, h abs is constantly positive for about 500 m in www.the-cryosphere.net/8/41/2014/The Cryosphere, 8, 41-57, 2014 distance of the profile.For most parts of the profiles from Vernagtferner, h abs do not exceed the uncertainty range (Sect.2.4).
The profiles H2 on Hintereisferner, K1, K2, K3 on Kesselwandferner, V8 on Vernagtferner and G1 on Gepatschferner are mostly located in the accumulation area of the glaciers.Only very few parts of the respective profiles cover bare ice areas.For each of those profiles, the median of h abs is distinctly negative.The cross profiles H3, H4, H5, K5 and V7, located in areas with bare ice at the surface, show up with a median of h abs , which is close to 0. The interquartile ranges of these profiles are comparably small in contrast to most LSC areas (Fig. 10).
Figure 13 shows the comparison of h abs to calculated uncertainties U abs (Eq.12) in combination with the results of the Mann-Whitney-Wilcoxon test (Sect.2.4).h abs calculated for LSC areas of several profiles exceed U abs , which is ±0.24 m for a mean snow depth of 2 m (Table 3).Along most parts of the profiles located in ice areas on the investigated glaciers, differences between z ALS and h GPR are significant (p = 0.05) but within the range of the calculated uncertainties (Eq.12).Distinct positive h abs were found on the very front of the glacier tongue of Hintereisferner and on the tongue of Kesselwandferner.Small areas of positive h abs were found in crevassed zones located in the higher elevated parts of the glaciers.
Table 3 highlights the increased negative values of h abs calculated for ALS measurements in LSC areas compared to h abs for ice areas.Figure 10 shows that only a small number of h abs values in ice areas are as negative as the majority of h abs determined in LSC areas.In these LSC areas mean h abs is −0.40 m with a standard deviation of 0.34 m, which is approximately 20 % of the mean value of h GPR with a standard deviation of 17 %.
Very locally, absolute differences are in the magnitude of about 1 m, which corresponds to relative differences of up to −40 % or more.For the ice areas of the profiles, the mean absolute difference is 0.004 m with a standard deviation of 0.27 m, in relative values −0.002 % with a standard deviation of 14.4 %.Depending on the location of the profiles, negative h abs were also found for ice areas (e.g.K2, G1; Table 3) or increased positive h abs of up to 1 m (K0, K4; Fig. 12) were calculated.

Differences in accumulation areas
The Landsat image of snow covered areas from 31 August 2009 cannot give direct proof for the general distinction between areas, where z ALS is considerably influenced by snow densification and vertical ice flow, and areas with only small h abs .In addition, only the areal extent of the snow cover (LSC) can be derived from Landsat imagery (30 m pixel resolution).No information on snow depths and, hence, spatial distribution of densification is provided by this data source.It does, however, give a broad picture of the distribution of typical firn areas and ice in the period of investigation.Areal snow cover extent on 31 August 2009 is generally in good agreement with the small firn areas observed for the last decade, which caused low accumulation area ratios (AAR) on these Alpine glaciers (e.g.http://www.glaziologie.de/)LSC area probably do not display the total zone of submergence ice flow, because flow dynamics are controlled by a variety of properties and usually adapt to mass changes within decades.Superimposed ice zones, which are part of the accumulation zone, are likely not visible to Landsat either.Furthermore, the applicability of Landsat scenes to distinguish typical accumulation areas on glaciers has to be proved carefully in terms of the state of seasonal ablation and the presence of new snow at the date of survey.However, LSC areas at the end of ablation season used here show consistently higher h abs compared to ice areas (Table 3,Fig.13), since, as mentioned in Sect.2.5, both vertical ice flow and densification of firn work in the same direction in typical accumulation areas of the glaciers.On Vernagtferner firn areas detected from Landsat data are very small.Even for these small firn areas the distinctly negative h abs are evident.The total profile on Gepatschferner shows negative differences for subtracting h GPR from z ALS .This can be attributed to the location of the profile in the accumulation zone of the glacier.However, even in this accumulation zone different magnitudes of h abs calculated for ice areas and h abs calculated for LSC areas could be detected (Table 3).Ambach and Eisner (1966) measured snow densities in a firn pit of 20 m depth on Kesselwandferner at an elevation of 3200 m a.s.l.Firn densities increased roughly linear from 630 kg m −3 in the uppermost layer to 836 kg m −3 in a layer at the age of 10 yr.The firn compression rate of 3 % per year in this 10 yr firn cover lead to a vertical compression of 16 m firn depth by approx.0.5 m per year.Similar rates of densification are presented in a density profile of upper Seward Glacier (Yukon, Canada) in Cuffey and Paterson (2010) (page 17, Fig. 2.3).However, percolating meltwater may contribute to annual densification of firn on Alpine glaciers.Because h abs were investigated for the comparatively shorter accumulation season, the contribution of densification of firn to h abs is assumed to be less than 0.5 m.

Comparison to measured ice flow velocities
Our data comparison does not allow to assess the quantitative proportion of the particular processes contributing to differences of z ALS from actual snow depths (Fig. 1).However, the four investigated glaciers show different dynamical behaviour, which is also the case for annual ice flow velocities presented in Sect.2.1.
In the uppermost parts of Hintereisferner negative h abs (Eq. 11) occur in between the crevassed areas (Fig. 11).These h abs exceed the uncertainties of the applied methods (Eq.12).However, no measurements of vertical ice flow exist for this area.Positive h abs could be detected only at the very front of the glacier tongue of Hintereisferner.These results coincide with the findings that mean annual vertical ice flow velocities are around 0.5 m s −1 and less apart from the front of the glacier tongue (Sect.2.1) .
Distinct positive and negative h abs are evident on Kesselwandferner (Fig. 12).Maximum h abs in the LSC areas are within the magnitude of the measured annual vertical ice flow velocity (Sect.2.1).h abs of the profile K3 of 0.6 m (Table 3) is slightly more than half of the measured annual vertical ice flow velocities.Locally, derived h abs are larger than half of the observed annual emergence of the ice surface which is about 1 m at the glacier tongue of Kesselwandferner (Sect.2.1 and Fig. 9).For the accumulation zone of Gepatschferner no measurements of ice flow velocities exists.On this glacier h abs are comparable to h abs of the uppermost profiles on Hintereisferner (H1) and Kesselwandferner (K3).Differences between z ALS and h GPR are within the uncertainty range of the measurements for most parts of Vernagtferner.This coincides with the decreased ice flow of the glacier observed by Mayer et al. (2013).

Influence of prevailing conditions
In glacierized mountain catchments, it is generally difficult to coordinate the dates of maximum ablation in fall and maximum accumulation in spring with the ALS surveys.In this study, both the temporal offsets between GPR and ALS measurements and the dates of ALS surveys in relation to the weather conditions are in good accordance to the required accuracy based on the previous assumptions for this study (Sect. 3,Figs. 3 and 4).
Due to the short delay between ALS surveys and GPR campaigns, simple assumptions were made for processes causing surface elevation changes within this time period.The thin layer of new snow in spring 2011 was neglected, because this layer was moved aside when pulling the GPR antenna over the snow surface on Gepatschferner and Kesselwandferner.On Vernagtferner this new-snow layer was assumed to be depleted by the snowmobile.In 2012 the GPR campaign took place before the ALS survey.The snow depth reduction due to melt conditions is assumed to be in the magnitude of the immersion of the GPR antenna into the uppermost snow layer during the measurements.
To account for the densification of snow accumulated before the ALS survey in fall 2010, measurements of this snow layer can give a rough estimate.From snow depth and density measurements in fall 2010 and spring 2011, a compaction of approx.0.07 m of a corresponding snow depth at ALS (t 1 ) of 1 m was calculated.However, measured snow depths at time of ALS(t 1 ) in 2010 was considerably less with mean values of 0.64 m on Gepatschferner and Kesselwandferner and 0.38 m on Vernagtferner.

Uncertainties
U abs was calculated accounting for stochastic uncertainties of the combination of the applied measurement methods (Eq.12).Further sources for systematic errors and uncertainties have to be mentioned: data were analyzed in the spatial resolution of 1 m for DEM data and for each GPR trace separated by the trace distance (Table 1).Both, in ALS and in GPR data, surface heterogeneities (e.g.crevasses) in the scale of 1 m are mapped.Horizontal ice flow shifts surface heterogeneities along the flow direction within one accumulation season.Hence, location of crevasses in fall are included in the data of z ALS , whereas GPR data show the location of the crevasse in spring (e.g.Fig. 12 at a profile distance of 700 m).In the analysis of h abs the horizontal shift is not considered, because ALS data are validated with respect to all processes contributing to surface elevation change.Resulting maximum differences are displayed as outliers in the Fig. 10.As shown by Joerg et al. (2012), stochastic uncertainties of ALS-derived surface elevations are strongly reduced on less inclined, homogeneous glacier surfaces.Mean systematic offset of each ALS survey was calculated based on a reference surface (Table 2).However, this systematic offset is not identical with a systematic offset of the ALS data at the investigated glaciers.Due to the fact that the systematic offsets are smaller than the assumed vertical uncertainty of the DEMs of ±0.15 m, no corrections have been applied to ALS data.
The quality of the GPR data analysis depends strongly on the amount of additional information.While stochastic uncertainties are due to the vertical resolution of the used GPR frequency, systematic errors and stochastic uncertainties result from velocity assumptions and the non-automated analysis of the GPR data.In a snowpack, r is a function of the mixing ratio between ice, air and liquid water (e.g.Lundberg et al., 2006).However, the proportion of liquid water in a wet snowpack can alter r and hence the signal velocity in a significant way (Frolov and Macheret, 1999;Sundström et al., 2012).Snow depths calculated from GPR data (h v , Eq. 9) and slope corrected snow depths (h cor , Eq. 4) from snow probings were compared to validate the used GPR velocities.On Vernagtferner 95 snow probings show a mean offset of +0.08 m with a standard deviation of 0.13 m, accounting for lower snow depths calculated from GPR data.This can be related to the usage of a snowmobile for GPR measurements.While GPR data were recorded in the track of the snowmobile, snow depth probings were performed beside the snowmobile in the undisturbed snow.On Hintereisferner and Kesselwandferner in 2012, the comparison of snow depths calculated from GPR measurements and snow depth probings show a mean difference of 0.01 m with a standard deviation of 0.13 m.

Consequence for ALS snow cover studies on glaciers
For most of the ice areas of the four investigated glaciers differences between z ALS and observed "ground truth" snow depths are less than the assumed uncertainties of the applied methods.The corresponding relative differences are less than 15 % for a mean snow depth of 1.5 m, decreasing with increasing snow depths.This is considerably less than the relative errors of measurements of solid precipitation, which may be up to 50 % when using rain gauges in high mountain catchments (Sevruk, 1985).In contrast, rain gauge measurements have the advantage of a high temporal resolution in comparison to ALS surveys.Hence the combination of both, ALS data for total volume and spatial distribution of the snow cover and rain gauge data for the timing of precipitation events, ensure a high potential for simulating the seasonal hydrology of mountain catchments more realistically.
Typical accumulation areas have to be considered when relating z ALS to snow depth.For these areas, z ALS has to be considered as a minimum snow depth value, as deviations to h GPR appear to be systematic.However, at the end of ablation seasons, recent firn areas cover considerably smaller parts of the investigated glaciers than bare ice areas.Accumulationarea ratios (AAR) derived from LSC are 0.14 for Vernagtferner, 0.15 for Hintereisferner and 0.41 for Kesselwandferner.Mean h abs of −0.4 m on average corresponds to a Fig. 13.Results of similarity between z ALS and h GPR applying a Mann-Whitney-Wilcoxon test (level of significance p = 0.05) and the comparison of h abs (Eq.11) to uncertainties of the two methods: ALS and GPR (Eq.12) are given.Each midpoint of a window including n = 300 samples is coloured indicating not significant h abs (black), significant h abs (blue), significant and negative h abs larger than uncertainty (red) and significant and positive h abs larger than uncertainty (green).The outlines of the investigated glaciers (black) and the snow cover derived from Landsat image of 31 August 2009 (LSC, striped) are shown.mean relative underestimation of approx.20 % for the investigated LSC areas.Regarding the whole investigated glacierized areas, this systematic underestimation results in a mean relative error of 3 to 7 % for the respective glaciers.Hence, the information on snow depths from z ALS is supposed to be a more reliable measurement not only of the spatial distribution of snow depths, but also for total volume of the snow cover in the investigated glacierized region at the time of the snow-on ALS survey (t 2 ).
In addition ablation in fall or already prevailing seasonal accumulation before the ALS survey at t 1 influence the accuracy of the derived ALS-derived snow depths on glacier surfaces.Furthermore, in contrast to non-glaciated areas, the ALS survey for t 1 has to be repeated annually, because glacier surface elevation varies annually.To calculate SWE from z ALS , additional assumptions on the relation of snow density on snow depth (e.g.Jonas et al., 2009) or a spatially uniform snow density (e.g.Lehning et al., 2011) have to be made.This results in additional uncertainties of the derived accumulation based on the conversion used and the prevailing conditions of the snowpack.
The data presented in this study were analyzed for differences between z ALS and actual snow depths, because recently an increased number of research projects using ALS data to investigate the spatial and temporal variability of the cryosphere are conducted in the mountain range of Ötztal Alps.The results may be transferable in the study region and to Alpine glaciers of similar size and with similar mass turnover.Nevertheless, various process parameters contributing to z ALS (Sect. 1) have to be considered, for example, on larger or steeper glaciers, on glaciers with faster ice flow, on glaciers gaining more accumulation or on glaciers with a less distinctive seasonal cycle of ablation and accumulation period.
In difficult mountain terrain, the acquisition of actual snow depths by GPR along profiles delivers less information on spatial snow distribution compared to ALS surveys covering total mountain catchments.However, systematic differences between z ALS and actual snow depths have to be analyzed separately for individual glaciers when using ALS for recording snow depths even in glacierized catchments.This study demonstrated that the usage of ALS data and GPR measurements appears to be a promising combination of methods to identify differences of z ALS from actual snow depths and corresponding uncertainties not only for single points, but along profiles on glacier surfaces worldwide.

Conclusions
The study presented here was accomplished to evaluate surface elevation changes derived from airborne laser scanning ( z ALS ) on glacier surfaces for snow cover studies in the Ötztal Alps.Therefore, differences of z ALS from snow depths calculated from more than 35 km of GPR profiles (h GPR ) on the four glaciers Hintereisferner, Kesselwandferner, Vernagtferner and Gepatschferner were analyzed.It is shown that GPR can be applied to validate z ALS regarding snow depths and be used to analyze systematic differences.The sum of all processes contributing to z ALS in addition to seasonal snow depths on glacier surfaces could be detected by combining ALS and GPR data.A more extensive measurement set-up would be necessary for detailed analysis of the partitioning and magnitude of the individual processes causing these differences.The delineation of snow covered areas at the end of ablation season from Landsat imagery highlights areas of increased differences between z ALS and "ground truth" snow depth derived from h GPR , snow pits and snow probing.This study shows that on four different glaciers systematic differences between z ALS and h GPR greater than calculated uncertainties of the applied techniques were observable almost exclusively within the respective accumulation areas.Most of ice areas on the investigated glaciers do not show differences between z ALS and actual snow depths larger than the uncertainty of the used combination of measurements.Distinct overestimation of actual snow depths by z ALS could be detected at the very front of the glacier tongues.In summary, relative differences of z ALS from actual snow depths are small for analysis of the mountain snow cover in the study area and comparable glacierized catchments.Hence, spatial distributed snow depths derived from z ALS at glacier surface provide an information on minimum snow depths which must be simulated by hydro-meteorological models, particularly in the higher elevated parts of the glaciers in this region.

Fig. 5 .
Fig. 5. Contour lines of the correction of GPR snow depth as a function of slope and snow depth.The relative frequency distribution of GPR snow depths h v and the frequency distribution of the occurring slopes are presented in horizontal and vertical bars, respectively.

Fig. 6 .
Fig. 6.Spatial distribution of ALS-derived surface elevation changes for for 10 October 2010 to 23 April 2011.Snow depths calculated from GPR measurements along the profiles at Vernagtferner plotted in the same colour-scale.Black dots indicate locations of snow depth probing.Locations of snow pits are shown as white squares.Snow covered areas derived from Landsat image of 31 August 2009 (LSC) are striped.

Fig. 7 .
Fig. 7. Spatial distribution of ALS-derived surface elevation changes for 10 October 2010 to 23 April 2011.Snow depths calculated from GPR measurements along the profiles at Gepatschferner and Kesselwandferner shown in the same colour-scale.Black dots indicate locations of snow depth probing.Locations of snow pits are shown as white squares.Snow covered areas derived from Landsat image of 31 August 2009 (LSC) are striped.

Fig. 8 .
Fig. 8. Spatial distribution of ALS-derived surface elevation changes from 4 October 2011 to 11 May 2012.Snow depths calculated from GPR measurements along the profiles at Hintereisferner shown in the same colour-scale.Black dots indicate locations of snow depth probing.Locations of snow pits are shown as white squares.Snow covered areas derived from Landsat image of 31 August 2009 (LSC) are striped.

Fig. 9 .
Fig. 9. Spatial distribution of ALS-derived surface elevation changes from 4 October 2011 to 11 May 2012.Snow depths calculated from GPR measurements along the profiles at Kesselwandferner shown in the same colour-scale.Black dots indicate locations of snow depth probing.Locations of snow pits are shown as white squares.Snow covered areas derived from Landsat image of 31 August 2009 (LSC) are striped.

Fig. 10 .
Fig.10.Box plots of the absolute differences h abs at the profiles for snow covered areas derived from Landsat (LSC, upper plot) and ice areas (lower plot) separately.Statistical distribution is shown in terms of the median (red line), the scatter (25 % and 75 % percentiles, blue box), the range of values within 1.5 times the interquartile range added to the 25 % percentile and to the 75 % percentile (black whiskers) and the outliers (red crosses).The number of points included in the statistics are plotted as grey bars (right scale).

Fig. 11 .
Fig. 11.Elevation of the surface in fall 2011 (blue curve), ALS surface elevation change z ALS (red curve) and GPR snow depth h GPR (green curve) along the longitudinal profile H1 at Hintereisferner.Snow covered areas derived from Landsat image of 31 August 2009 (LSC) are shaded in blue.

Fig. 12 .
Fig. 12. Elevation of the surface in fall 2011 (blue curve), ALS surface elevation change z ALS (red curve) and GPR snow depth h GPR (green curve) along the longitudinal profile K4 at Kesselwandferner.Snow covered areas derived from Landsat image of 31 August 2009 (LSC) are shaded in blue.
. For DEMs derived from ALS data of the Hintereisferner region, Bollmann et al.(2011)showed that for slopes below 35 • a vertical accuracy of ±0.15 m can be assumed.

Table 2 .
Hall et al., 1987)S flight campaigns investigated in this study.The vertical accuracy is obtained from differences between ALS surface elevations and elevations of a known reference surface.The vertical accuracy is shown in terms of a mean offset (mean) and the standard deviation (σ ).Hall et al., 1987).The application of Landsat imagery is limited by the discrepancy of the date of surveys to the date of interest and the occurrence of clouds masking the surface information.The conditions of this scene fulfill all requirements in terms of no clouds, recent new-snow precipitation on the glaciers and an advanced state of snow ablation in Cogley et al., 2011)accumulation area ratios (AAR,Cogley et al., 2011)derived using the NDSI for separation of snow and bare ice are comparable to those observed within the annual mass balance measurements on the glaciers in recent years.Hence the delineated snow covered areas were assumed to indicate recent accumulation areas on the investigated glaciers and are referred to as Landsat snow cover (LSC) hereafter.

Table 3 .
Properties of the GPR profiles in terms of ID, profile length, number of sections, the mean section length (l S ), minimum elevation and maximum elevation.Mean value of the surface elevation changes derived from ALS data z ALS (Eq. 1) and the mean absolute difference of this value from mean GPR snow depth h abs (Eq.11) are shown for LSC areas and ice areas of each profile.