1. INTRODUCTION
The ongoing Arctic warming results in elevated rates of meltwater production on glaciers and ice sheets (Vaughan and others, Reference Vaughan2013). However, the amount of melt is not necessarily equal to the loss of water from a glacial system, as runoff from accumulation zones of glaciers is reduced by refreezing of meltwater in the snow and firn pore space (e.g. Reijmer and others, Reference Reijmer, van den Broeke, Fettweis, Ettema and Stap2012) or retention in perennial firn aquifers (Forster and others, Reference Forster2014; Christianson and others, Reference Christianson, Kohler, Alley, Nuth and van Pelt2015). Accumulation zones hence act as reservoirs of meltwater in frozen and/or liquid form, effectively reducing and delaying runoff (e.g. Harper and others, Reference Harper, Humphrey, Pfeffer, Brown and Fettweis2012). Accurate description of subsurface processes in snow and firn is important to reduce uncertainties in measurements and simulations of the glacier mass and energy balance.
Snow and firn physical properties are known to exhibit strong lateral variability (e.g. Bell and others, Reference Bell2008; Schweizer and others, Reference Schweizer, Kronholm, Jamieson and Birkeland2008). At the scale of metres the original heterogeneity of a snowpack is created by the irregular pattern of snow deposition influenced by the local topography and the wind field (Sturm and Benson, Reference Sturm and Benson2004). Conditions in snow and firn below the surface continuously evolve, which results in rapid and spatially non-uniform thermal and mechanical changes (Colbeck, Reference Colbeck1982). A major contribution to the local spatial variations in structure of snow and firn comes from the flow and refreezing of water (Jordan and others, Reference Jordan, Albert and Brun2008).
While uncertainties in interpretation and extrapolation of point and profile observations (Machguth and others, Reference Machguth, Eisen, Paul and Hoelzle2006; Bell and others, Reference Bell2008; Trujillo and Lehning, Reference Trujillo and Lehning2015) can be better assessed and reduced by observation of the lateral heterogeneity of snow and firn properties, it is a major challenge to select a reliable field technique. In relatively thin snow and firn packs visual description of pit walls and direct density measurements are conventionally used (Benson, Reference Benson1960; Sturm and Benson, Reference Sturm and Benson2004; Fierz and others, Reference Fierz2009; Sugiyama and others, Reference Sugiyama2012). Alternatively, near infrared imagery (Tape and others, Reference Tape, Rutter, Marshall, Essery and Sturm2010) or a penetrometer, measuring resistance of snow to insertion of a rod (Schneebeli and others, Reference Schneebeli, Pielmeier and Johnson1999; Proksch and others, Reference Proksch, Löwe and Schneebeli2015), can be used for detailed two-dimensional descriptions of the snow properties.
Coring can be used for stratigraphical studies in deeper snow and firn (Brown and others, Reference Brown, Harper, Pfeffer, Humphrey and Bradford2011; Harper and others, Reference Harper, Humphrey, Pfeffer, Brown and Fettweis2012; Machguth and others, Reference Machguth2016). Video observations in boreholes (e.g. Pohjola, Reference Pohjola1994; Malone and others, Reference Malone, Hubbard, Merton-Lyn, Worthington and Zwiggelaar2013) provide a less detailed but also less logistically demanding alternative to coring for description of firn stratigraphy (Langhammer, Reference Langhammer2014), location of summer/autumn erosion surfaces (Hubbard and others, Reference Hubbard, Roberson, Samyn and Merton-Lyn2008) and measurement of firn compaction rate (Hawley and Waddington, Reference Hawley and Waddington2011).
Ground penetrating radar (GPR) measures the variability of dielectric permittivity and electrical conductivity of a medium in space, which serve as a proxy for its physical characteristics (Arcone, Reference Arcone and Jol2009). It offers a relatively affordable and fast alternative to studies in pits and boreholes. GPR is a non-invasive method, which is of particular importance in snow and firn studies as it allows for repeated surveys. However, interpretation of results and validation is often challenging (e.g. Brown and others, Reference Brown, Harper, Pfeffer, Humphrey and Bradford2011). Frequency modulated continuous wave radars having high vertical resolution (Marshall and others, Reference Marshall, Schneebeli and Koh2007) allow to detect seasonal subsurface layers (Richardson and others, Reference Richardson, Aarholt, Hamran, Holmlund and Isaksson1997) and accumulation rates (Forster and others, Reference Forster, Davis, Rand and Moore1991) and may be used to sense flow fingers in a thin seasonal snowpack (Albert and others, Reference Albert, Koh and Perron1999). Pulse radar systems are usually operated at lower central frequencies and can reveal the structure of deeper snow and firn packs. Gridded and profile GPR surveys combined with coring have been used to study layering in firn on the western slope of the Greenland ice sheet (Dunse and others, Reference Dunse2008; Brown and others, Reference Brown, Harper, Pfeffer, Humphrey and Bradford2011), at Kongsvegen, Svalbard (Langley and others, Reference Langley, Lacroix, Hamran and Brandt2009) and in the Canadian Arctic (Gascon and others, Reference Gascon, Sharp, Burgess, Bezeau and Bush2013).
Given the aforementioned limitations in the spatial and temporal extent and resolution of field studies, the need to assess the impact of meltwater flow and retention in snow and firn on glacier mass and energy balance at the scale of glaciers, ice caps or even ice sheets resulted in a number of modelling approaches. Layered models describing water flow through snow and firn (e.g. Bartelt and Lehning, Reference Bartelt and Lehning2002; Reijmer and Hock, Reference Reijmer and Hock2008; Gascon and others, Reference Gascon2014) often utilize the concept of a horizontal wetting front: water advances when the refreezing capacity of the current layer is eliminated and capillary forces cannot accommodate additional water supply. However, both field (e.g. Schneebeli, Reference Schneebeli1995; Bøggild, Reference Bøggild2000; Campbell and others, Reference Campbell, Nienow and Purves2006) and laboratory (e.g. Waldner and others, Reference Waldner, Schneebeli, Schultze-Zimmermann and Flühler2004; Katsushima and others, Reference Katsushima, Yamaguchi, Kumakura and Sato2013) observations show that water percolating through snow and firn pack often follows preferential flow paths originating at stratigraphic boundaries and infiltration along them is more effective than in the adjacent snow and firn volumes. The limited success in multidimensional simulation of water flow in snow (e.g. Illangasekare and others, Reference Illangasekare, Walter, Meier and Pfeffer1990; Katsushima and others, Reference Katsushima, Kumakura and Takeuchi2009; Leroux and Pomeroy, Reference Leroux and Pomeroy2016) and, particularly, thicker firn packs is partly explained by the lack of empirical data on the snow and firn physical properties and their spatial variability.
The purpose of this paper is to present results from a 3 year plot scale (≈10 m × 10 m × 10 m) firn stratigraphy study done at Lomonosovfonna (1200 m a.s.l.), Svalbard, using cores, borehole video and GPR surveys and to report on the strengths, weaknesses and agreement of results from the different methods.
2. STUDY SITE
Lomonosovfonna is a large ice field on central Spitsbergen, Svalbard, nourishing several outlet glaciers including Nordenskiöldbreen (Fig. 1a). Surface energy and mass balance along with firn conditions were studied by Zinger and others (Reference Zinger1966), Zagorodnov and Zotikov (Reference Zagorodnov and Zotikov1980), Pälli and others (Reference Pälli2002), Van Pelt and others (Reference Van Pelt2012, Reference Van Pelt2014). Deep ice cores were drilled at Lomonosovfonna in 1997 (Isaksson and others, Reference Isaksson2001) and in 2009 (Wendl, Reference Wendl2014; Vega and others, Reference Vega2015) and a series of shorter cores was drilled and studied in the recent years (Vega and others, Reference Vega2016).
Fig. 1. (a) location of the study site at the Lomonosovfonna ice field with a map of Svalbard shown in the lower right corner. Also shown is the location of the Soviet expedition in 1965 Zinger1966; (b) the grid of boreholes and GPR profiles used in the study site. Note: the distances between boreholes 2 and 2a, 4 and 4a, 6 and 6a, 8 and 8a are not to the main scale.
The altitude of the equilibrium line at Nordenskiöldbreen was estimated at 719 m a.s.l. (Van Pelt and others, Reference Van Pelt2012). The study site is in the accumulation zone at 78.824°N, 17.432°E, 1200 m a.s.l. on a gently sloped dome. The local ice thickness derived from a GPR survey is 192±5.1 m (R. Pettersson, unpublished data). Ice coring at the same site in 2009 showed a thickness of the firn layer of ~20 m (Wendl, Reference Wendl2014). In accordance with the long-term warming trends observed in Svalbard (Nordli and others, Reference Nordli, Przybylak, Ogilvie and Isaksen2014) and Lomonosovfonna (Van de Wal and others, Reference Van de Wal2002), during 1989–2010 the mean annual air temperature at Nordenskiöldbreen increased by more than 1°C, resulting in a negative glacier-averaged mass balance of −0.39 m w.e. a−1 (Van Pelt and others, Reference Van de Wal2012). Recent field and modelling results in the area showed that at the elevation of the study site the average accumulation rate during 2007–12 is ≈0.75 m w.e. a−1 (Van Pelt and others, Reference Van Pelt2014) and the melt rate during 1989–2010 was on average ≈0.34 m w.e. a−1 (Van Pelt and others, Reference Van Pelt2012). Despite percolation of water below 10 m, runoff is small as most meltwater refreezes in the snow and firn profile exerting a major influence on density, temperature and metamorphism in firn at the field site (Van Pelt and others, Reference Van Pelt2012).
3. METHODS
For presentation of some of the results the glacier surface in April 2014 was defined as the common reference with respect to depth for all 3 years. The vertical offsets of the glacier surfaces in April 2012 and 2013 from that common reference resulting from accumulation are defined as 2.44 and 1.32 m, respectively, based on annual measurements at the nearby stake S11 (Fig. 1a) in April 2012–14. This aligns layers deposited at approximately the same time but recovered in different field campaigns. However, the effect of gravitational densification is not accounted for and additional vertical shifts caused by local-scale lateral variability in accumulation rate are possible. Moreover, we follow the definition of firn as the snow that has survived one or more ablation seasons but has not been transformed to ice (Cogley and others, Reference Cogley2011).
3.1. Firn cores
Firn cores were collected each year using a Kovacs corer with the internal diameter of 10.5 cm. They were used for direct observation of stratigraphy and density. Coring locations in 2012–14 were within 50 m from each other. The length of the retrieved cores was 9.9 m (2012), 11.3 m (2013) and 13.6 m in 2014 (Table 1). The bottom depth of each core section was calculated as the summed lengths of extension rods below the surface plus the length of the coring barrel before extraction to the surface. Coring below ≈10 m was challenging in all three field campaigns as the cold barrel tended to freeze on to the temperate firn at this depth. During extraction and handling of the firn cores 4–7% of the core record was lost in different years.
Table 1. Details on stratigraphy data collected using firn cores, video and GPR surveys
The firn cores were processed in the cold laboratories of the University Centre in Svalbard (Longyearbyen) and of the Norwegian Polar Institute (Tromsø). To enable comparison of firn cores with video observations and GPR measurements the number of distinguished facies was limited to two: porous icy media (snow and firn) and solid ice. Ice lenses thicker than 1 mm were noted. The density of the firn column was derived by measuring the mass and dimensions of core samples. Core pieces with uneven geometry were cut in rectangular samples to reduce uncertainty in calculation of their volume. In 2012 and 2014 density was measured in ≈0.1 m long samples, while in 2013 bulk density of core pieces up to 0.5 m long was measured.
3.2. Borehole video surveys
Optical properties of firn were inspected using a video camera with an internal light source in multiple boreholes with a diameter of 5.5 cm drilled by a Kovacs auger. Brushing the boreholes to remove the drill chips from the walls greatly improved the quality of the video surveys, but at the same time the loose material accumulated at the bottom of the narrow hole and reduced the depth of survey.
The principle behind deriving stratigraphy data from borehole videos is that ice, being less reflective to light than snow and firn, appears darker in the video frames (e.g. Sjögren and others, Reference Sjögren2007; Langhammer, Reference Langhammer2014). Two types of media were distinguished in the videos: firn and ice.
Boreholes surveyed in 2012–14 were positioned on a grid (Fig. 1b): nine holes (labelled 1–9) were placed on a square grid with 3 m spacing, four other holes (labelled 2a, 4a, 6a and 8a) were located 6 m away from the four sides of the grid. In 2012 we filmed all 13 holes during lowering and lifting of the camera, in 2013 five holes (labelled 2–6) and in 2014 four holes (labelled 6–9) within similar grids were drilled and surveyed (Table 1). In 2012 hole 5 was drilled in 0.3 m from the cored borehole and thus the video data grid is centred on it (Fig. 1b). In 2013 and 2014 the video data grids were within 100 m from the coring locations.
Depth marks appearing in video frames were generated by a rotary encoder coaxial with the wheel feeding the camera cable into the borehole. The accuracy of assessing the absolute depth with respect to the surface is estimated to be better than 0.1 m based on the difference between initial and final depth values (i.e. before lowering and after rising of the camera). The resulting videos with the resolution of 720 × 576 pixels were visually analysed and the upper and lower depths of ice layers were recorded with the precision of 5 mm as they appeared in the boreholes.
Multiple stratigraphy profiles derived in each of the three field campaigns in April 2012–14 allowed to derive three laterally averaged datasets revealing the horizontally consistent patterns in firn stratigraphy. For each year the probability of ice layer occurrence for 5 mm depth intervals was defined as the ratio between the number of video surveys that indicated an ice layer there and the number of surveys that reached that particular depth.
3.3. GPR surveys
Grids of GPR profiles were measured using Malå GeoScience ProEx impulse radar and shielded antennas with the central frequency of 800 MHz. In 2012 the GPR grid intersected with the grid of holes surveyed by camera (Fig. 1b), while in other years centres of the radar grids were positioned in 10–50 m from the centres of the boreholes grids (Table 1). At such distances there are no significant variations in the surface conditions, and the general character of the stratigraphy can be expected to be preserved while the horizontal extent of most ice layers in snow and firn is typically smaller and these details are likely to be lost. The radar was triggered every 0.05 m along the profiles by a distance wheel encoder and a stack of four traces consisting of 2024 samples was measured. Spacing of 0.2 m between profiles was controlled by wooden markers plugged in snow at the sides of the grid and a string connecting them. Several tens of radar profiles 6–9 m long each were measured in 2012–14 field campaigns (Table 1). To facilitate the comparison of the GPR data against the video surveys done in 2012, the GPR traces coinciding with the inspected boreholes were marked in the field.
GPR data were processed using custom functions written in MATLAB programming language. First a dewow filter was applied, then the zero time of every trace was adjusted to the first sample exceeding 5% of the trace amplitude range (Yelf, Reference Yelf2004). Gain correction was applied to account for the loss of the radar wave power with depth. Finally, a low-pass filter eliminated the frequencies higher than 1000 MHz from the dataset.
For analysis of the spatial heterogeneity of firn properties we attempt to derive a 3-D GPR dataset. To do so, we first aligned the GPR profiles in the along-profile direction adopting the approach originally suggested by Evangelidis and Psarakis (Reference Evangelidis and Psarakis2008) for alignment of images to radargrams of neighbouring profiles. Then Whittaker–Shannon interpolation was applied to reduce the separation between profiles to 0.05 m preserving the frequency characteristics of the original signal (Oppenheim and others, Reference Oppenheim, Schafer and Buck1999). This multiplied the number of profiles by ca 4. Finally, 3-D Kirchhoff time domain migration (e.g. Hagedoorn, Reference Hagedoorn1954; Schneider, Reference Schneider1978) using an anti-aliasing filter (Gray, Reference Gray1992) was applied to relocate the reflection events to their original position. We then applied automatic gain correction with the window size of 180 samples (≈9% of the total number of samples in a trace) in order to keep amplitude of both noise and signal at the same level throughout the profile.
For the conversion of the vertical scale from two-way travel time to depth the Kovacs and others (Reference Kovacs, Gow and Morey1995) relation was used. The density values measured in the firn cores were averaged over 0.5 m intervals and converted to relative real dielectric constants giving the interval velocity of the GPR signal. While in 2012 and 2014 the chosen sampling frequency allowed sounding to a depth of ≈10.6 m with ≈5 mm spacing between samples, in 2013 the sampling frequency was doubled resulting in equivalent reduction of the depth of sounding and sample spacing. With the mean velocity of the GPR signal of 0.2 m ns−1 and the central frequency of 800 MHz the theoretical wavelength λ is 0.25 m and, consequently, our setup is able to distinguish only between reflectors separated by more than λ/2 = 0.125 m (Hubbard and Glasser, Reference Hubbard and Glasser2005).
A composite GPR signal was calculated to highlight horizontally persistent horizons with increased backscatter for each of the three datasets. For that the sum of the instantaneous amplitudes of several thousand individual traces was calculated using the Hilbert transform (Oppenheim and others, Reference Oppenheim, Schafer and Buck1999) and scaled to have the maximum at 0.5 to facilitate the visual comparison with the data from video inspections.
To facilitate the quantitative analysis of the GPR data and make it comparable with the stratigraphical descriptions in boreholes we use a thresholding technique to derive reflectors in the radar traces in a consistent way compared with subjective visual interpretation of radargrams. First instantaneous amplitudes of all radar traces were calculated, then the samples for which the parameter is larger than the mean plus 2 standard deviations (SDs) of all instantaneous amplitudes were defined as reflecting while others are classified as non-reflecting (Fig. 2). In the example shown in Figure 2 four reflectors (marked by grey bars) were retrieved from the trace, while the rest of the column was assigned to be non-reflecting. The routine effectively cuts off 5% brightest samples on Hilbert filtered radargrams. The choice of the threshold is arbitrary and larger values result in smaller fractions of samples identified as reflecting, consequently, reflectors are fewer and thinner. It has to be noted, though, that the thickness of reflectors rather depends on the strength of the dielectric permittivity and electrical conductivity contrasts between firn layers than on the actual thickness of the layers at the interface with which the backscatter event is initiated.
Fig. 2. Analysis of GPR data. The blue curve is a radar trace from processed and migrated dataset. The green curve shows the corresponding instantaneous amplitudes derived using the Hilbert filter. Similar series were used in Figures 5–7. The violet line shows the threshold defined as the mean plus two standard deviations of all instantaneous amplitudes in traces measured in one field campaign. It was used to decipher reflectors in individual radar traces, which are shown as grey-shaded two-way travel time intervals.
3.4. Firn temperature measurements
During the periods April–October 2012, April–July 2013 and April 2014–April 2015 the evolution of firn temperature was monitored in shallow boreholes using thermistor strings. To minimize preferential percolation of water along them, the holes were backfilled with surface snow. Voltage drop over temperature stable reference resistors connected in series with thermistors was measured by a Campbell Scientific CR10X logger. Accuracy of temperature measurements is estimated as <0.2°C.
4. RESULTS
4.1. Firn temperature and density
The measured evolution of the snow and firn temperature is illustrated by Figure 3. The extent and strength of the cold wave in the upper metres of the firn profile were similar by 10 May in 2012 and 2013. In 2013 the first liquid water infiltration event occurred on the 3 June, which is 17 days earlier than in 2012. In the next few days a warm wave was observed in the upper part of the profile, but after a week it was reduced by the surface cooling. By the end of August 2012 the firn profile down to 12 m was temperate. This was likely also the case in the succeeding year as already by mid-July the upper 5 m of the snow and firn pack reached that state.
Fig. 3. Snow and firn temperature measured during 2012 and 2013.
Figure 4d shows the density distribution measured in the three firn cores. Gravitational densification and freezing of liquid water in the firn matrix can be expected to result in an increase of the firn density between the field campaigns. The most prominent density increase (of ≈130 kg m−3) occurred between April 2013 and April 2014 in the depth range 1.5–5.3 m. The period was also marked by density rise in the lower part of the profile (10–12.7 m) of ≈80 kg m−3. Other density increases over time reach comparable values but have smaller vertical extent (typically <1 m) and, given the local lateral variability in firn conditions, could be a result of cores being drilled at different locations.
Fig. 4. Stratigraphy distribution observed in firn cores drilled in 2012 (a), 2013 (b) and 2014 (c). Different shades of the blue colour mark the depth intervals where snow (light), firn and ice (dark) were observed. Density distribution measured in cores in 2012–2014 (d).
4.2. Firn stratigraphy
4.2.1. Firn cores
The firn core profiles are presented in Figures 4a–c. The visual observation of the firn cores gives the most detailed picture of the firn stratigraphy and allows to detect even very thin ice layers (<1 mm) and to distinguish between layers with different shape and size of the snow and firn grains and pores. Details of the firn structure such as laterally inconsistent ice volumes, ice layers with uneven or inclined boundaries, vertical ice columns above ice layers are apparent from the core face. The entire length of the three cores recovered in 2012–14 consists of snow and firn with several dozen of ice layers (Table 2) from 1 mm to 0.2 m thick. The average thickness of ice layers observed in cores was 1.4–2.3 cm for the different years, but the variability is high with a SD larger than the mean value (Table 2). The mean number of ice layers per metre of core is 3.2–5.8 corresponding to the ice content in firn (i.e. mean ice layer thickness multiplied by the number of ice layers per metre of core) of ≈6–10%.
Table 2. Number and thickness of ice layers/reflectors in individual boreholes/GPR traces
AM, arithmetic mean; SD, standard deviation; RSD, relative SD in % (SD/AM).
Note on SD in thickness: for cores it is for ice layers recovered in a single firn core while for video surveys it is for the mean thicknesses of ice layers from different boreholes.
4.2.2. Borehole video surveys
Interpretation of firn strata from borehole videos brought a less detailed picture of the firn stratigraphy (Fig. 5) than the visual description of the cores: it was only possible to reliably distinguish between two types of facies (firn and ice). In Figure 6 we compare the firn stratigraphy derived in 2012 by video observations in borehole 5 (d) and from the core (c) drilled in 0.3 m from each other. While some ice layers captured by direct observations are matched within 0.15 m by ice layers in the video survey (e.g. at 4.9, 7.1–7.3, 7.5, 8.3 and 9.7 m), most and, primarily, the thinner ice layers found in core were not seen in borehole 5. That could be a result of imperfections in data acquisition techniques as well as a consequence of limited lateral extent of the ice layers.
Fig. 5. Firn stratigraphy observed in 2012 using video surveys. Ice layers found in different boreholes are shown as horizontal black lines. The blue curve to the right illustrates the number of boreholes that indicated an ice layer at a specific depth.
Fig. 6. Vertical profiles of firn stratigraphy derived in 2012 using different methods: firn core stratigraphy (c); ice layers recovered by video surveys in boreholes 2 (a), 5 (d) and 8 (f); reflectors and average instantaneous amplitude calculated for 3 radar traces sounded above boreholes 2 (b), 5 (e) and 8 (g).
The number of ice layers detected by video surveys in multiple boreholes in the three field campaigns was only 25–50% of the number found with core observations (Table 2). The average thickness of ice layers is 3.7–5.9 cm, which is about double the size measured in the firn cores. Thick ice layers are more characteristic for the lower part of the profile (Fig. 5) and in boreholes 5 and 9 at the depth of 12 m, 0.5 m thick ice layers were observed. Ice layers constitute ≈6–9% of the profile (Table 2), a value close to the relative ice content of the cores.
Video observations in multiple boreholes show high lateral variability of the ice volumes in firn. As illustrated by Figure 5, only few layers (for example at 4.5, 5.5 and 7.0 m) can be consistently tracked among neighbouring boreholes, suggesting that horizontal extent of most ice layers is less than the distance between holes (ranging between 3 and 18 m). The high lateral variability of ice layers in firn is further confirmed by the number of boreholes that showed an ice layer at a certain depth (blue curve in Fig. 5), according to which most of the ice layers were observed in less than two boreholes in 2012 and single boreholes in the other 2 years. Furthermore, the SDs in the number and mean thickness of ice layers recovered in the different holes (Table 2) constitute up to 20–40% of the corresponding mean values, which also suggests the low representativeness of a stratigraphy record derived at a single point.
Lateral averaging of the video data from multiple boreholes surveyed in each year allowed to derive more representative descriptions of the ice layers in firn. Up to eight boreholes showed ice layers within the same 5 mm depth interval in 2012 (Fig. 5) and a maximum of four holes had ice layers at the same depth in 2013 and 2014. The maxima in the grid averaged video data reach 0.3–1 and correspond to horizons with high probability of ice layer occurrence. Referenced to the glacier surface in April 2014, in 2012 the horizons were found at the depths of 7.1, 7.8, 8.7, 9.5, 10.7, 12.1, 12.7 and 14.0 m (Fig. 7a), in 2013 at 3.7, 7.1, 8.8, 9.5, 10.8 m (Fig. 7b) and in 2014 at 3.2, 6.7, 7.9, 9.0, 10.6, 11.5 m (Fig. 7c). Some of these depth levels are consistent throughout the 3 years of observations, as will be discussed below. In 2012, when the video inspections were done in the largest number of boreholes (13), a positive trend with respect to depth was revealed in the probability values reached at the horizons of preferential ice formation (Fig. 7a), suggesting more laterally consistent and thick layers in the lower part of the profile.
Fig. 7. Spatially averaged signals from the video and GPR data collected in 2012 (a), 2013 (b) and 2014 (c). Also shown are the summer surfaces 2011–2007 (coloured labels) interpreted from a 500 MHz GPR profile Pelt2014 and corrected for surface accumulation and gravitational firn compaction after April 2012.
4.2.3. GPR surveys
Radargrams of the profiles done in 2012 are presented in Figure 8. Individual radar profiles show regions of increased backscatter that are typically elongated in the horizontal direction (Fig. 8a). However, identification of the preferential reflecting horizons over longer distances (≈0.5 m) is ambiguous due to the limited lateral extent of individual high backscatter regions. The latter particularly applies to the lower part of the profile (below 5 m) where the vertical gradients in backscatter intensity are lower. The latter observation is also confirmed by the averaged radar traces for each of the years (green curves in Fig. 7) having steeper gradients in the upper part of the profile than in the lower. Reflecting horizons become more apparent after calculating the average of more than 100 profiles in the 3-D GPR datasets in the across-profile direction. The averaged profile from 2012 (Fig. 8b) exhibits horizontally consistent reflections at several depths: 1.7, 2.3, 3.2, 3.8, 4.7, 5.6, 6.6, 7.3, 8.6 and 10.1 m, which are less intermittent than the areas with increased backscatter found at same levels in the individual profiles. The ratios between the mean vertical and horizontal gradients in the averaged radar profiles (2012–2014) are 21.2, 37.2 and 23.3, while the corresponding values averaged over all individual profiles collected in each field campaign are 1.7, 3.4 and 1.8. Thus in the averaged radargrams vertical gradients are more dominant than in the radargrams of individual profiles.
Fig. 8. Instantaneous amplitudes of the GPR data collected in April 2012. (a) radar profile crossing locations of boreholes 2, 5 and 8 (shown white, see Fig. 1b); (b) GPR data averaged in across-profile direction, also shown are the summer surfaces 2011–2007 (coloured labels) interpreted from a 500 MHz GPR profile Pelt2014.
The fraction of reflecting samples (5%) resulting from the thresholding routine roughly corresponds to the relative ice content in firn found in cores and videos. However, compared with the ice layers visually observed in cores, reflectors in radar traces composed only 16–24% in number (Table 2). The reflector stratigraphy derived from the GPR traces measured in April 2012 in the closest vicinity of the boreholes 2, 5 and 8, surveyed by camera, is presented in Figure 6. In several cases strong reflectors are found within 0.15 m from the ice layers in the boreholes: 2.8, 4.5, 5.5, 7.2 m in borehole 2 (Figs 6a, b); 4.8 and 9.7 m in borehole 5 (Figs 6d, e); 3.9, 4.9 and 5.5 m in borehole 8 (Figs 6f, g). In some cases amplitude peaks of the GPR signal match with the ice layers found in boreholes (7.2–7.7 m in borehole 8), while the thresholding routine barely detects the reflections. However, the overall similarity between the two data sources is poor: radar reflectors are often detected where no ice layers were observed and most ice layers are not associated with radar reflectors sounded in vicinity of the corresponding boreholes.
The comparison between horizontally averaged radar and video data is presented in Figure 7. The GPR datasets collected in 2012 and 2014 are more useful for the purpose than the 2013 data with the penetration depth of only 4.7 m. The GPR signal (green curves) exhibits multiple spikes, which are expectedly consistent with the bright horizontally elongated areas in the averaged radargrams (Figs 7a and 8b). Black points show the samples in the averaged GPR traces that lie above 0.7 of the amplitude range of each of the data series and at the same time are within 0.15 m of the depth intervals where ice layers were observed in several (three for 2012 and two for 2013 and 2014 data) boreholes using the video surveys. The radar signal from all 3 years contains multiple spikes in the upper part of the profile, which can not be correlated with any maxima in the data from boreholes. However, in 2012 (Fig. 7a) when the largest number of boreholes was inspected using video and the radar survey reached below 10 m such samples are found at six depth intervals. In 2013 the overlap between the two series is limited to the upper 4.7 m and only five boreholes were inspected. The resulting apparent similarity between the radar and video signals is lower (Fig. 7b): only one of the two prominent maxima in video data is repeated by the GPR at ≈3.8 m. In 2014 (Fig. 7c) all but one maxima found in the video data within overlapping depth range are repeated by the averaged GPR signal. Thus the averaged signals from the radar and video surveys show a prominent match in contrast with the data for individual traces and boreholes.
5. DISCUSSION
5.1. Firn temperature and density
Our results show that by the end of the ablation season in 2012 and 2013 the upper 10 m of the firn pack were isothermal at 0°C. Similar measurements done on the 16 August 1965 at 1050 m a.s.l. (Zinger and others, Reference Zinger1966) and on the 1 August 1976 at 1120 m a.s.l. (Zagorodnov and Zotikov, Reference Zagorodnov and Zotikov1980) yielded significantly lower firn temperature at 10 m depth of −3 and 2°C correspondingly. During the last decades the upper accumulation area of Lomonosovfonna has likely seen transformation from the ‘cold infiltration’ zone to the ‘warm infiltration’ zone according to the classification of Shumskii (Reference Shumskii1964) and deep firn water storage and runoff from the domain are now possible.
In 2013 the first refreezing-induced warming in firn was registered as early as the first week of June and the firn profile was isothermal at 0°C 1 month earlier than in 2012 (Fig. 3). This observation may be explained by the change of the wind direction over Svalbard in summer 2013 from the dominant north-westerly to an unusual south-south-westerly, which caused high air temperature (Lang and others, Reference Lang, Fettweis and Erpicum2015).
The gradual increase of firn density measured in cores from 400–450 kg m−3 at the surface to 600 kg m−3 at the depth of 10 m is in correspondence with the previously reported results from the site (Pohjola and others, Reference Pohjola2002; Vega and others, Reference Vega2016). The prominent density increase of ≈130 kg m−3 measured between 1.5 and 5.3 m during April 2013–April 2014 is a result of several processes. Gravitational densification can account for 41 kg m−3 a−1, when calculated following the method suggested by Ligtenberg and others (Reference Ligtenberg, Helsen and van den Broeke2011) using the density measured in 2013 core, temperature evolution measured during April 2014–April 2015 and assuming the annual accumulation rate of 0.75 m a−1 (Van Pelt and others, Reference Van Pelt2014). The amount of meltwater refreezing in spring 2013 can be estimated from the firn refreezing capacity c (Reijmer and others, Reference Reijmer, van den Broeke, Fettweis, Ettema and Stap2012):
$c = \rho \left \vert T \right \vert CL^{ - 1}$
. Assuming that the mean temperature (T) and density (ρ) of firn at the depth 1.5–5.3 m are −9°C and 440 kg m−3, correspondingly (see Figs 3 and 4), relative heat capacity of ice (C) is 2050 J kg−1 K−1, latent heat of water fusion (L) is 334 000 J kg−1 the calculation yields a value of 24 kg m−3. Autumn refreezing of the meltwater suspended in pores of firn with the density of 500 kg m−3 could account for further densification by 25 kg m−3, assuming the liquid water content of 5% (Schneider and Jansson, Reference Schneider and Jansson2004). Apart from the uncertainty in the above estimates, the resulting inconsistency of 40 kg m−3 between the observed density increase and the sum of possible contributions from different densification processes can be explained as an effect of multiple melt-freeze cycles in spring and summer 2013 and/or lateral variability of firn density.
5.2. Firn stratigraphy
Our observations show a high spatial heterogeneity in firn stratigraphy at the scale of several metres (Figs 5 and 8, Table 2). In accordance with earlier studies (e.g. Marsh and Woo, Reference Marsh and Woo1984), we interpret the observation as evidence of the highly heterogeneous nature of water infiltration and refreezing in snow and firn. The latter introduces uncertainty in point observations of firn properties and in results of simulations relying on the hypothesis of a horizontally uniform wetting front in firn (e.g. Gascon and others, Reference Gascon2014).
Spatial averaging of the video signal from multiple boreholes revealed that ice layers tend to group in several preferential horizons (Fig. 7). This finding partly deviates from the conclusions made by Brown and others (Reference Brown, Harper, Pfeffer, Humphrey and Bradford2011) for a thick firn pack with few thin ice lenses studied on the western slope of the Greenland ice sheet (1997 m a.s.l.). In contrast to our results Brown and others (Reference Brown, Harper, Pfeffer, Humphrey and Bradford2011) did not find any horizontally consistent pattern in firn stratigraphy across their eight cores spaced by 1.5–14 m. Although, similar to our findings, their 500 MHz GPR survey on a 20 m× 20 m grid returned several laterally consistent internal reflection horizons (IRHs), which were interpreted as annual features formed by preferential formation of ice layers at the summer surfaces from previous years. Dunse and others (Reference Dunse2008) presented a similar interpretation of horizontally coherent IRHs found in firn at the Greenland ice sheet at a slightly lower site (1940 m a.s.l.).
To test if the maxima in our grid-averaged GPR and also video data likewise originate at the buried summer surfaces, we compare their depths with positions of the summer surfaces 2011–2007 detected by Van Pelt and others (Reference Van Pelt2014) in a 15 km long 500 MHz GPR profile done at Lomonosovfonna and Nordenskiöldbreen in April 2012. The uppermost point of this profile was within 10 m from our plot and is used here for comparison. The summer surfaces are interpreted from IRHs based on their distinct character and lateral continuity at the kilometre scale. To estimate the depths of the summer surfaces after April 2012 we introduce offsets accounting for the gravitational densification of firn following the aforementioned routine (Ligtenberg and others, Reference Ligtenberg, Helsen and van den Broeke2011). Offsets calculated for April 2013 and 2014 vary from 0.13 to 0.53 m a−1 for different depth levels.
In 2012 the summer surfaces found by Van Pelt and others (Reference Van Pelt2014) (coloured labels at 4.0, 5.6, 7.5, 8.9 and 9.9 m in Fig. 7a, also see Fig. 8b) closely match with the maxima in our grid-averaged GPR and video data (green and yellow series in Fig. 7a). Positions of some of the summer surfaces in the following 2 years are also in reasonable correspondence with the maxima in the video and GPR data (7.1, 8.5, 9.4 m in 2013, see Fig. 7b; 3.8, 6.8, 8.1, 9.1 m in 2014, see Fig. 7c), suggesting that the horizons of preferential ice layer formation are persistent in time. The apparent biases can be attributed to the local scale variability in the accumulation rate, effects of firn metamorphism or the inaccuracies in the vertical scales of GPR and video data and the estimated compression rates.
In line with earlier studies (Dunse and others, Reference Dunse2008; Brown and others, Reference Brown, Harper, Pfeffer, Humphrey and Bradford2011; Van Pelt and others, Reference Van Pelt2014; Sold and others, Reference Sold, Huss, Eichler, Schwikowski and Hoelzle2015) we suggest that most of the prominent maxima in the composite data from our radar and video surveys originate at the firn horizons exposed to melt at the end of ablation seasons (summer surfaces). The rate and pattern of water infiltration through snow and firn is highly dependent on the structure of the media (Campbell and others, Reference Campbell, Nienow and Purves2006). The upper decimetres of snow continuously exposed to melt (often referred to as ‘summer surface’) are likely to have an increased density as a result of multiple melt-freeze cycles caused by fluctuations of the surface energy flux around zero at the daily to weekly scale (Dunse and others, Reference Dunse2008). Therefore, it can be expected that formation of ice layers and IRHs associated with them will occur at the buried summer surfaces with pre-existing firn layers that are less permeable to the flow of water.
Some maxima in our data (4.7, 6.2, 7.2 and 8 m in Fig. 7a) fall between the summer surfaces from Van Pelt and others (Reference Van Pelt2014). They are also noticeable on the 500 MHz radargrams (Van Pelt and others, Reference Van Pelt2014) but were not regarded as buried summer surfaces based on the lack of lateral consistency on the kilometre scale. The IRHs falling in between the depths interpreted as buried summer surfaces also appear as less horizontally consistent features in our 800 MHz radargram (Fig. 8b). We interpret them as firn layers that experienced additional densification due to liquid water refreezing events during the accumulation season or were excessively packed by wind (Bøggild, Reference Bøggild2000; Parry and others, Reference Parry2007) before being covered by successive accumulation. Surface melt rates during occasional inter-seasonal melt events are likely to be lower than during the ablation season, causing less prominent perturbations in snow stratigraphy resulting in less laterally consistent IRH. The snow surfaces exposed to strong winds often have a rough topography and the IRHs associated with such buried layers can be expected to be undulating and horizontally intermittent.
5.3. Assessing the methods used to map firn stratigraphy
Direct observation in firn cores and video surveys in multiple boreholes yielded similar values of relative ice content in firn of ≈8%. At the same time video inspections resulted in significantly larger mean thickness and smaller number of ice layers than core observations, suggesting that thin ice layers were omitted by video records or post processing. At the expense of reduced data quality the video surveys offer two main advantages over more conventional ocular observations of firn cores. Firstly, drilling of the narrow boreholes by auger is faster than coring and, secondly, there is no need for transporting the cores to the cold laboratory or arranging the facilities for core analysis in the field. It is suggested that the reproducibility of the coring results can be improved by thorough cleaning of the borehole walls before filming and using video acquisition hardware allowing for faster sampling and possibly manual control on the camera sensitivity.
The point nature and vertical orientation of studies in boreholes bears the risk of over-interpreting recovered ice layers as laterally continuous. However, the observed high spatial variability of firn structure (Table 2) and the poor match between stratigraphy profiles from boreholes separated by 3–18 m and surveyed in one field campaign (Fig. 5) limits the value of stratigraphy studies at single points and calls for spatial averaging over several cores or borehole video surveys to reveal the depth intervals where preferential ice layer formation is observed.
While evidences from studies in boreholes suggest that ice layers tend to be thicker and more laterally consistent with depth (Figs 4a–c, 5 and 7), the regions with increased backscatter of the GPR signal become less distinctive with depth (Figs 6b, e, g, 7 and 8): the vertical gradients in the instantaneous amplitude are less steep and reflectors are less consistent. These observations are expected from the applied 800 MHz GPR system that allow to sounding firn down to ~10.6 m. Due to the gravitational settling and water refreezing in the firn matrix the density gradients between ice layers and adjacent firn decrease with depth and associated differences in electric permittivity are also less pronounced.
Small lateral extent of stratigraphical features in firn, does not justify studies in boreholes for characterization of their shape and size due to associated limitations in horizontal resolution of such surveys. GPR studies can be done with a very high horizontal resolution and combined with effective techniques for post-processing and quantitative analysis, which potentially allows for detailed mapping of the firn structure. At the same time the interpretation of the individual radar reflectors is not straightforward and firn cores and/or video inspections in boreholes remain of great importance for validation of the GPR results and constraining the velocity model used in data post-processing and for conversion of the vertical scale to distance units.
The match observed between the reflectors identified in GPR data and ice layers found in individual boreholes is poor (Fig. 6). However, considering the good match between grid-averaged video and GPR data (Fig. 7), we suggest that radar returns are associated with ice layers in firn but may also come from other targets often related to them. High amplitude of the GPR signal returns is caused by spatial discontinuities of the firn dielectric permittivity and electrical conductivity (Hubbard and Glasser, Reference Hubbard and Glasser2005). In the cold firn the latter can be caused not only by sharp density contrasts (firn/ice), which can be seen in both core and video surveys, but also by less prominent changes in firn density, alternations of firn micro structure and conductive inclusions (Arcone, Reference Arcone and Jol2009; Sold and others, Reference Sold, Huss, Eichler, Schwikowski and Hoelzle2015), detection of which in boreholes would require more detailed studies. The mismatch between visual and radar derived descriptions of firn stratigraphy (Fig. 7) can also be attributed to the uncertainties arising from the lateral variability in firn density not captured by a single firn core. Measured firn density is used in the migration routine and in the time/depth conversion of the vertical scale.
The suggested thresholding routine produces distinct and discrete radar reflectors, which should not be over-interpreted due to several limitations. Firstly the routine is most sensible to the amplitude of the GPR signal, consequently targets with strong contrasts in relative permittivity or electric conductivity with respect to the over- and underlying layers will result in pronounced reflectors independent of the reflecting object's actual thickness. Secondly, a reflection generated by a target in firn is a single wavelet that in our case has the length λ = 0.25 m. Thus even very thin targets in firn may appear as relatively wide reflectors in radargrams and the derived reflectors can overestimate the presence of reflecting objects.
6. CONCLUSIONS
The temporal and spatial dynamics in stratigraphy and density of firn at Lomonosovfonna (1200 m a.s.l.), Svalbard, were studied in April 2012–14 on plots (≈10 m × 10 m × 10 m) using cores, videos surveys in multiple boreholes and 800 MHz GPR profiling. We present results of the field studies and report on the strengths, weaknesses and agreement of the three methods used for stratigraphical descriptions.
Results from three cores and video surveys in multiple boreholes suggest that the upper 12 m of firn contain multiple ice layers 0.1–50 cm thick constituting ≈8% of the profile. Most ice layers are rather thin (≈2 cm) and are laterally discontinuous at the scales >3 m. Ice layers in firn tend to be concentrated in several preferential horizons, which are also consistent with the maxima in the grid-averaged instantaneous amplitude of the GPR signal. We interpret the horizons as buried surfaces that have been previously exposed to melt or wind-driven densification, several of them were consistently recovered in the three field campaigns.
While direct firn observations in cores provided most details about firn stratigraphy, video surveys in boreholes omitted thin layers and allowed to distinguish only between firn and ice. The poor agreement between ice layers visually observed in individual boreholes and reflections found in individual GPR traces is attributed to differences in the physical principles behind the methods and imperfections of the data acquisition routines. The tight correlation between grid-averaged video and GPR data suggests that reflections are associated with the ice layers in firn and that the GPR surveys can be utilized to recover depth intervals with preferential ice layer formation.
The high lateral variability of the firn strata, calls for observations in multiple cores or boreholes to derive more reliable density profiles and distinguish layers with increased probability of ice layer occurrence. Limited lateral extent of ice volumes in firn calls for utilization of high resolution GPR studies for characterization of their shape and size. Such detailed stratigraphic maps along with borehole studies for validation can provide an insight into the processes of water flow and refreezing in firn.
ACKNOWLEDGEMENTS
This publication is contribution N82 of the Nordic Centre of Excellence SVALI, ‘Stability and Variations of Arctic Land Ice’, funded by the Nordic Top-level Research Initiative. This work is also funded by the Vetenskapsrådet grant 621-2014-3735 (V. A. Pohjola). Logistical support was provided by the Swedish Polar Research Secretariat, the Norwegian Polar Institute and the University Centre in Svalbard. The authors are grateful to the Ymer-80 foundation, Arctic Field Grant of the Research Council of Norway, Margit Althins stipend of the Royal Swedish Academy of Sciences for additional funding of the field campaigns. We are also indebted to the two Reviewers who contributed their time and effort to the constructive comments that helped to improve the manuscript.
Open access
