Vertical seismic profiling of glaciers: appraising multi-phase mixing models

Abstract

We have investigated the speed of compressional waves in a polythermal glacier by, first, predicting them from a simple three-phase (ice, air, water) model derived from a published ground-penetrating radar study, and then comparing them with field data from four orthogonally orientated walkaway vertical seismic profiles (VSPs) acquired in an 80 m deep borehole drilled in the ablation area of Storglaciären, northern Sweden. The model predicts that the P-wave speed increases gradually with depth from 3700ms^–1 at the surface to 3760ms^–1 at 80m depth, and this change is almost wholly caused by a reduction in air content from 3% at the surface to <0.5% at depth. Changes in P-wave speed due to water content variations are small (<10 ms^–1); the model’s seismic cold–temperate transition surface (CTS) is characterized by a 0.3% decrease downwards in P-wave speed (about ten times smaller than the radar CTS). This lack of sensitivity, and the small contrast at the CTS, makes seismically derived water content estimation very challenging. Nevertheless, for down-going direct-wave first arrivals for zero- and near-offset VSP shots, we find that the model-predicted travel times and field observations agree to within 0.2 ms, i.e. less than the observational uncertainties.

Introduction

Polycrystalline glacier ice near the melting temperature (∼0°C) contains μm-scale water inclusions at the three-grain junctions (Nye and Frank, 1973; Mader, 1992). This water significantly increases the strain rate of temperate glacier ice (Duval, 1977), which has motivated numerous studies targeted at measuring microcrystalline water content in glaciers (e.g. Vallon and others, 1976; Macheret and others, 1993; Murray and others, 2000; Pettersson and others, 2004; Bradford and Harper, 2005; Navarro and others, 2005; Bradford and others, 2009; Endres and others, 2009; Gusmeroli and others, 2010b).

Polythermal glaciers are composed of both cold and temperate ice, and represent a perfect experimental setting to investigate changes in microcrystalline water content. This is because the boundary between cold and temperate ice, the cold-temperate transition surface (CTS; Hutter and others, 1988; Pettersson and others, 2003), is a prominent englacial feature. Across the CTS a jump in water content (from near-zero to single-digit percentages) and hardness (stiffer cold ice to softer temperate ice) is typically observed (Aschwanden and Blatter, 2005; Gusmeroli and others, 2010a). Such an increase in water content on the temperate side of the CTS lowers the propagation speed of both radar and seismic waves. The speed of these signals can therefore theoretically be used to infer microcrystalline water content.

Temperate ice also contains air, mainly in the form of bubbles (Raymond, 1976), and is more appropriately described as a three-phase medium containing solid polycrystals, liquid microcrystalline water and air in bubbles. Such a three-phase approach has recently been demonstrated to be particularly crucial in terms of water content estimates from radar speeds (Bradford and Harper, 2005; Bradford and others, 2009; Gusmeroli and others, 2010a). When the air content of temperate ice is neglected (i.e. when only a two-phase ice/water model is used) water content estimates are typically overestimated (Bradford and Harper, 2005; Gusmeroli and others, 2008). This sensitivity is, however, still unclear from a seismic perspective. Several studies (Navarro and others, 2005; Endres and others, 2009) have interpreted variation in seismic speeds as an indicator of water content using a two-phase (ice/water) approach; however, these studies have not addressed the effects of air content (responsible for bulk density; Kohnen and Bentley 1973) and ice temperature (Kohnen, 1974). Thus the effect of air content on seismic waves propagating in a three-phase mixture of ice, water and air remains to be investigated.

Here we investigate the vertical distribution of compressional-wave (P-wave) speeds using a three-phase framework, by forward modeling the effects that known air content, temperature and water content values have on seismic waves. We simulate the P-wave structure of a vertical section in the ablation area of the polythermal glacier Storglaciaren, Sweden (Fig. 1), using a high-resolution (1 m depth-step) density and water content profile measured from borehole radar (Gusmeroli and others, 2010a). We then compare the modeled speed structure with the results of a series of four orthogonally orientated walkaway vertical seismic profiles (VSPs) acquired in the same area in the summer of 2008.

Fig. 1. Map of Storglaciären showing the location of our VSP experiment. (Coordinates Swedish Grid RT90 1616415–7536885.) The inset shows the survey geometry, with the borehole in the middle of the cross. Four walkaway VSP lines were obtained by moving the shot position progressively away from the borehole. The tick-marks on the line indicate shot positions. The approximate position of the equilibrium-line altitude is indicated with a dashed curve.

Storgläciaren: Physical Properties of the Study Area

Storglaciaren (Fig. 1) is a small (3 km long), well-studied valley glacier located in the Kebnekaise massif, Sweden, just north of the Arctic Circle (e.g. Jansson, 1996; Holmlund and others, 2005). Storglaciären is entirely temperate in the accumulation area, whereas it has a polythermal structure in the ablation area. The polythermal structure at our study site is composed of a cold, permanently frozen, surface layer overlying a temperate core (Hooke and others, 1983; Holmlund and Eriksson, 1989; Pettersson and others, 2003). The water content in the temperate ice on Storglaciären varies spatially (Pettersson and others, 2004) but is generally small and homogeneous with depth, typically 0.6–1% (Gusmeroli and others, 2010a). Higher water contents (possibly up to 2%) are suspected for deeper, near-bed ice (Aschwanden and Blatter, 2005, 2009; Gusmeroli and others, 2010a); these higher values result from heat of deformation, which melts some ice and increases the water content.

Gusmeroli and others (2010b) studied the physical properties of a vertical section of the upper ablation area of the glacier (Fig. 1), and obtained vertical profiles of temperature, air content and water content for the upper 80 m of the ice column. A thorough description of the experiments conducted to obtain the englacial measurements is given by Gusmeroli and others (2010b); here we present a brief synopsis of the results.

The temperature profile, measured by thermistors installed within the ice, shows that the CTS is at 21 m depth at this location (Fig. 2a). Ice temperature at the surface is below the melting point and increases with depth. The air content (Fig. 2b) was measured using the near-surface (1m depth) radio-wave speed (Gusmeroli and others, 2010a). This shallow estimate gave an initial value which was subsequently modeled using a pressure/porosity model based on the theory of Bradford and others (2009). The modeled air content profile (Fig. 2b) shows a decrease of density with depth due to air-bubble compression. Gusmeroli and others (2010b) used this air profile to correct the radio-wave speeds. These corrected speeds can be inverted for water content (Fig. 2c). The result shows that water content values measured from cross-borehole radar at Storglaciären are relatively homogeneous with depth in the upper part of the temperate ice, with a tendency towards higher values at greater depths. The mean water content measured within the ice column is 0:6 ± 0 :3%. In this study we simulate the slightly higher water contents at depth with a three-layer water content structure consisting of an upper layer (21-60 m depth) with 0.5% water content and a lower layer (61-85 m depth) with 0.6% water content. The three profiles in Figure 2 represent the starting point of our seismic forward modeling. Although the water content model here is a simplified version of that of Gusmeroli and others (2010b), it sufficiently describes its most important characteristic: a homogeneous layer with only a slight increase in water content at greater depths. These physical properties (Fig. 2) will be used to obtain the P-wave speed structure within the glacier.

Fig. 2. Models of the physical properties of the study area. (a) Temperature distribution with depth. (b) Air content distribution with depth. (c) Water content distribution with depth. Models in this figure are derived from the field measurements described by Gusmeroli and others (2010a).

Modeling the P-Wave Speed Structure

Propagation speed of compressional P-waves within glaciers depends on temperature, air content and water content. Cold temperatures, low air content and low water contents result in high seismic speeds. Based on in situ observations Kohnen (1974) derived the temperature dependence of P-wave speed as

(1)

where T is ice temperature (°C). v _P(T) is often considered a reference value for cold, water-free ice. If glacier ice was not a multi-phase material, then v _P(T) would be the speed structure within the subsurface. However, since temperate glacier ice can also be composed of air (Fig. 2b) and water (Fig. 2c), a three-phase medium should be considered. Gusmeroli and others (2012) reviewed suggestions for temperature correction published since Kohnen’s (1974) work, and concluded that they are all relatively similar. For example, Helgerud and others (2009) suggested a gradient of –2.7. The difference between Helgerud and others (2009) and Kohnen (1974) is almost negligible for the small temperature range used in this study and therefore the Kohnen (1974) relationship is considered adequate.

In polar ice sheets, P-wave speed increases with depth, mirroring the increase of density, p (Kohnen and Bentley, 1973). This behavior reflects the reduction of the relative proportion of air in the two-phase firn. A simple, well-established method to estimate the bulk speed in a multiphase medium is Wyllie’s time-average equation (Wyllie and others, 1958). The bulk propagation speed in the medium will be proportional to the relative volumetric concentrations of the components. Thus, by defining vP_a and v _Pw as the P-wave speed in air (330 ms–¹; Sheriff and Geldart, 1999) and water (1450 ms–¹; Sheriff and Geldart, 1999), respectively, we can obtain v _P(T, a, w), the bulk speed in the three-phase glacier ice,

(2)

where Φ_a and Φ _w are volumetric air and water fractions, respectively. Both vP_a and vP_w may vary with depth in the borehole, due to temperature and pressure change (Mackenzie, 1981; Sheriff and Geldart, 1999); however, these variations are ∼2ms–¹, which is insignificant for our study. We thus assume constant v _Pa and vP_w. Equation (2) includes the strong dependence of seismic speeds on density (Kohnen and Bentley, 1973).

Vertical Seismic Profiling (Vsp)

Borehole VSP is an established technique in exploration seismology which allows quantification of subsurface seismic speed (Sheriff and Geldart, 1999). In typical geological applications heterogeneity in the stratified subsurface results in high P-wave speed contrasts, which allow differentiation of material properties. Speed contrasts in glaciological applications are smaller than in rock applications, so an additional aim of this study is to examine the usefulness of VSP experiments in alpine glaciers. In VSP a series of receivers, located in a borehole, record seismic energy generated at the surface at some variable distance from the borehole (Fig. 3).

Fig. 3. Schematic representation of a VSP experiment using the negligible-ray-bending assumption. The shot position, at offset x = is indicated by a star located at A. Receivers, indicated by circles, are located in the borehole at depths z ₁, z ₂ and z ₃. A, B, C and E are the trigonometric points used in the calculations.

Travel-time modeling of VSP

In a hypothetical VSP experiment the first prominent arrival is the direct wave, which propagates directly from the source to the in-hole receiver. If the source is located at some offset, x, from the borehole (Fig. 3), the path of the propagating direct wave is oblique; otherwise, when the source is at the position of the borehole and the borehole is straight, the path is vertical and in this case the VSP may be called zero-offset VSP (Sheriff and Geldart, 1999). In layered media the propagation of seismic waves obeys Snell’s law: ray paths refract at the interface between two different layers, according to the incidence angle and the magnitude of the speed contrast. Such behavior causes rays to bend at interfaces.

In common geological situations seismic speeds in the subsurface vary over a very large range (<300 to >2000m s–¹; Miller and Xia, 1998; Sheriff and Geldart, 1999). Since speed contrasts are thus very high, and ray curvature is not negligible, a bent-ray approach is required to model VSP travel times (Moret and others, 2004). A bent-ray approach should also be applied in polar ice sheets or in the accumulation area of mountain glaciers, where the presence of a strong density gradient leads to P-wave speeds <1500ms–¹ in snow, rising to ∼3000ms–¹ in firn and reaching ∼3800 ms–¹ in ice (e.g. King and Jarvis, 2007). Such a significant speed gradient is not present in the ablation area of glaciers, where snow and firn are absent. Consequently for the ablation area of Storglaciären we can model ray paths using the negligible-ray-bending assumption (Schuster and others, 1988) and ray paths can be modeled as straight.

The travel-time model is structured using the physical properties illustrated in Figure 2. An 85 m VSP is modeled with 1 m depth-step sampling of the direct wave (the first-arriving compressional wave that, following a near-straight ray path, propagates between the source and the receiver). In our model the shot position is located at distance x from the borehole top. Receivers are located at depths z_i, where subscript i indicates the depth of the layer. By considering a VSP experiment of a simple two-layer case (Fig. 3) the travel times, t, of the direct wave recorded at positions z ₁ and z ₂ are given by

(3)

and

(4)

respectively. The path lengths, AE = EC, are easily found by trigonometrical relationships on the triangle ABC (Fig. 3).

The propagation angle, a, changes according to receiver depth. Thus for each z

(5)

Similarly, AE, the travel path in each layer, will diminish as a increases:

(6)

Thus the travel time of the received P-wave at depth z, t_z , is simply given by the sum of the travel times in each layer of different speed, v(z):

(7)

where n is the total number of layers considered in the model.

Field experiment

We collected four orthogonally orientated walkaway VSPs in an 80m deep, water-filled borehole located in the upper ablation area of Storglaciären in summer 2008 (Fig. 1). The borehole was drilled using the Stockholm University hot-water drill (Jansson and Näslund, 2009). A walkaway is simply obtained by moving the shot position progressively away from the borehole. Seismic energy was generated with a sledgehammer and received by a 12-channel, 1 m spacing, 10 Hz hydrophone string, which was lowered down the borehole to different depths. The borehole was entirely water-filled. The array was lowered down-hole so the travel times of the direct wave could be estimated at a 1 m depth spacing. The shot position was located using a tape measure. The area investigated by the VSP experiment was generally flat or had very little topography. A flat area was primarily chosen because it facilitated drilling operations. Some topography (up to ∼2 m height difference in 60 m distance) was, however, present in the along-flow direction. In other words, the shot located 30 m west of the borehole was ∼2 m higher than the shot located 30 m east of the borehole (Fig. 1). The effect of this topographic gradient in the along-flow direction is discussed below.

Sample data for a VSP with source located 3 m away from the borehole are shown in Figure 4. The collected seismic profiles are of good quality. They show clear P-wave direct arrival first breaks, and are rich in other arrivals (i.e. converted-wave direct arrivals; englacial and subglacial P-wave reflections; tube waves, a ubiquitous noise source in VSPs) that have varying dominant frequencies as revealed by bandpass filtering (Fig. 4). Here we use only the direct arrivals first breaks to compare seismic travel times with our a priori model. The other seismic signals observed in the VSP records and the full velocity inversion will be the subject of future investigations.

Fig. 4. Sample VSP data collected at Storglaciären. The source (hammer and plate) was located 30 m from the borehole. (a) Filtered data (cut off at 150 and 800 Hz, cosine ramps between 150–300 and 400–800 Hz). (b) Sketch of the main events: P and converted P arrivals (dashed curves), tube waves (dotted curves) and internal reflections (solid curves).

Results

Modeled P-wave speed structure in the ablation area

The modeled P-wave speed through the 80 m ice column imaged (Eqns (1) and (2)) varies by ∼100ms–¹ (Fig. 5, asterisks). The lower speed (∼3700ms^–1) is found in the upper part of the ice column, due to the presence of air within the cold ice. The transition from cold ice to temperate ice is marked by a ∼10 m s–¹ decrease in speed (Fig. 5). The seismic CTS at Storglaciären is thus characterized by a 0.3% decrease in speed, ∼10 times lower than the decrease in speed observed at the radar CTS (Gusmeroli and others, 2010a). The air content and its decrease with depth dominates the overall modeled seismic speed structure (Figs 5 and 2b). The temperature-induced decrease of speed from the cold to the temperate side of the CTS is very small (5 m s–¹; Fig. 5, curve T). In the bulk part of the ice column, speed gradually increases with depth, reaching the highest value of 3760 m s–¹ at 80 m depth.

Fig. 5. P-wave speed structure of the ablation area of Storglaciären inferred from the physical properties illustrated in Figure 2. The final speed model used to predict VSP arrivals (asterisks) was calculated from Eqn (2) by considering temperature, T, air and water. The speed structure calculated by temperature only (Eqn (1)) is shown with crosses. Combinations of T, water and T, air are also indicated with circles and dots, respectively.

VSP travel times

The arrival times of the propagating P-wave at different receiver depths were obtained from four orthogonally orientated VSP lines (Fig. 1). In this way, potential azimuthal variation in P-wave speed (e.g. ice anisotropy) can be detected, since four different propagation directions were sampled (up-flow, down-flow and across-flow in two opposite directions). We varied the source-offset distance, x, in 5 m increments from 5 to 30 m, and collected six VSPs for each line. Here we only show the differences in modeled and measured travel time for a representative 30 m offset, since it sufficiently summarizes the main findings of our experiment.

Results from VSP arrival times for the four lines (Fig. 6 for a 30 m offset VSP) show that the modeled P-wave structure generally matches the field experiment. This is especially evident for the two transverse profiles (LINE1 and LINE3; Fig. 6) as the misfit between modeled and measured travel time for these two lines is within 0.2 ms (Fig. 6). This misfit is within the travel-time uncertainties for the P-wave direct arrival first picks, which is approximated as S _t = 0:25 ms (twice our sampling rate of 0.125 ms). _t can be approximated to twice the sampling rate, because the sampled first break might actually be arriving one or two samples earlier. The uncertainty in travel time causes speed uncertainties which vary with the travel path. Our resultant P-wave velocity uncertainty, _v, is a function of the uncertainties in our observed travel time, t, and the ray-path length, d, for a given source-to-receiver arrangement,

(8)

Fig. 6. Example of a misfit plot (difference between modeled and measured arrival time) for a 30m offset VSP acquired at Storglaciären.

and the error in speed varies inversely with distance travelled, d. v is very high (between ±400 and ±100ms–¹ ) in the upper 10 m of the ice column, but it reduces to ±50 m s–¹ in the bulk of the ice column. We can thus be confident that the measured travel times at depths >10m in the across-flow direction agree with the speed structure indicated by asterisks in Figure 5, within a margin of ±50 ms^-1.

In the along-flow direction, the observed direct arrival travel times show a deviation of up to 1 ms from our speed model, the up-glacier travel times being consistently slower than the model predictions (LINE4; Figs 6 and 7). The misfit between our model and the seismic observations is indeed higher for these two profiles and well above the 0.25 uncertainty in our travel-time picks, with a deviation of up to ±1.0ms below 40 m receiver depth (Fig. 7). These results suggest apparent anisotropy in the flow direction; P-wave speed appears to be lower when propagating along-flow, whereas it is higher when propagating in the opposite direction. Apparent anisotropy is also confirmed by looking at the mean misfit measured in each survey (Fig. 7).

Fig. 7. Mean misfit for the 24 VSP surveys acquired. The error bars represent the standard deviation of the misfit for each survey. The gray area is 0 ± _t where _t is our mean standard deviation (0.18 ms) which is about twice the sampling rate (_δt = 0:25 ms); a reasonable error for our travel-time estimates. LINE1 and LINE3 are stable on the value of 0 ± _t , whereas those waves sampled in LINE2 and LINE4 are slower and faster, respectively. Theout-of-trend data point in LINE1 is believed to be caused by poor data quality for that particular VSP. Each line has six data points, which correspond to the six walkaway positions covered (5, 10, 15, 20, 25 and 30 m).

Discussion

For down-going direct-wave first arrivals, for zero- and near-offset VSP shots, we find that the travel times predicted by the model agree with field observations to within 0.2 ms, i.e. less than the observational uncertainties (Fig. 7). This is particularly true for the VSP surveys acquired perpendicular to the glacier flow, where we found negligible differences between modeled and measured P-wave arrival times. Potential along-flow seismic anisotropy was also observed because waves propagating down-glacier arrived ∼0.8–1.0 ms after those propagating up-glacier. In this section we discuss the two main issues that arose in this study: (1) the relationship between seismic speeds and physical properties in the ablation area of a polythermal glacier and (2) the presence of apparent along-flow seismic anisotropy.

Is the straight-ray approach safe and appropriate?

The ice of the ablation area is homogeneous, old metamorphic glacier ice, already densified. We are thus correct to assume an a priori absence of a significant density gradient that leads to an equally important velocity gradient. However, our choice of a straight-ray approach is simplistic and violates Snell’s law (Fig. 3). Figure 5 demonstrates that seismic velocities can vary by ∼-100ms–¹ in the upper 20-30 m of the ice column; a potentially important velocity gradient thus exists and whether or not this impacts our analysis needs to be determined.

We thus undertook further analysis; based on the model observed in Figure 5, we compared modeled travel times following the straight-ray assumption and utilizing a framework that bends rays for each layer and does not violate Snell’s law. The results show that there is no practical difference between the two (Fig. 8) and that the straight-ray approach described above is applicable for the ablation area of a glacier where the near-surface speed gradient is minimal.

Fig. 8. Comparison between the straight-ray and the bent-ray travel-time modeling of a VSP. Dots and solid lines indicate straight-ray and bent-ray model, respectively. There is no practical difference between the two for the small velocity gradients observed in our study.

Seismic speeds and physical properties of glacier ice

Seismic speeds in temperate ice are important because they can be used to infer the microcrystalline water content and fabric, which are fundamental parameters in modeling ice viscosity. Navarro and others (2005) and Endres and others (2009) inferred water contents as high as 2% by inverting low seismic speed for water content in two-phase (ice/water) models. Our model and its agreement with the VSP travel times suggest that the P-wave speed in the cold surface layer is relatively low, ∼3700ms–¹. This relatively low value is also supported by the speed analysis of the refraction surveys collected in the same area (Gusmeroli and others, 2010b). Gusmeroli (2010) provides further agreement, measuring an average P-wave speed of 3617 ±57 m s–¹ in the uppermost ice. Here we suggest that this relatively low seismic speed is primarily caused by the air content of the ice. Such low seismic speeds have traditionally been interpreted as being typical of temperate ice with some volumetric percentages of microcrystalline water (Benjumea and others, 2003; Endres and others, 2009). A simple (two-phase) interpretation of speeds as low as 3600 m s–¹ could lead to the wrong conclusion that the upper 20 m of the ice column are composed of temperate ice with a volumetric water content of ∼3%. This value is calculated using Eqn (2) for a water/ice mixture. Other two-phase models, such as the Riznichencko model (Benjumea and others, 2003), do not differ significantly.

Instead seismic speeds measured and modeled in our study area suggest that relatively low P-wave speeds (e.g. 3600-3700 ms–¹) are to be expected for cold ice in the ablation area. These low speeds are not due to the presence of water but, instead, to the presence of small percentages of air. As observed in polar ice masses, v _P increases because of increasing pressure which decreases bubble size and air content. Although dominant, the magnitude of the air-content-induced change is, however, very small (only 100ms–¹ between the lowest and the highest speed in Fig. 5). An even smaller, perhaps undetectable, change (only 10ms–¹) is caused by the low water contents measured at Storglaciären. We believe that the change in vP due to water content might be hard to detect given the uncertainty often observed in many glaciological seismic speed estimates (tens of m s–¹).

Equation (2) indirectly incorporates the densification with depth using air content. We thus consider the bulk density of the polycrystalline glacier ice to be primarily conditioned by the volumetric air content. Kohnen and Bentley (1973) pointed out that a universal density/speed function applicable to different glaciological settings might be difficult to obtain. Most of our knowledge of density/speed relationships comes from observations in polar firn (Kohnen and Bentley, 1973; King and Jarvis, 2007), where density changes up to 0.5 gcm^{- 3} lead to P-wave speed variations of ∼1500 m s–¹ (King and Jarvis, 2007). In the ablation area of a mountain glacier (e.g. this study) there is no firn and the density gradient is significantly lower.

The necessity of including air content when inferring water content from radio-wave speed estimates appears to have been established in recent radioglaciological literature (Gusmeroli and others, 2008, 2010b; Bradford and others, 2009). The reason is that the presence of air increases the bulk radio-wave speed and thus causes underestimation of water content. From a seismic perspective, this behavior is quite the opposite; the presence of air within the ice matrix decreases the bulk speed, meaning that water content can be overestimated if air content is ignored. The glaciological consequences of this are important; by ignoring the presence of air, bulk speeds in the water-free cold ice of the ablation area of Storglaciaren suggest water content up to 3%.

We stress that the P-wave structure modeled in this study is valid only for the ablation area of the glacier. In the accumulation area the presence of snow and firn might lead to speeds considerably lower than those reported in Figure 5. Additionally, we have very little knowledge about the water content of the ice in the accumulation area. Schneider and Jansson (2004) found values as high as 3% and these values are likely to be significant. Acquiring borehole geophysical data in the ablation area of glaciers is easier than working in the accumulation areas, where the occurrence of snow-covered crevasses and the presence of porous firn result in logistical challenges that must be overcome to improve our ability to use ground-based geophysics to detect ice properties.

Azimuthal variation in seismic results

Results from our VSP surveys (Figs 6 and 7) might suggest an up-/down-glacier oriented seismic anisotropy. P-waves propagating from west to east (down-flow) propagate more slowly than those waves propagating in the opposite direction (Fig. 6). This behavior is consistent with P-wave arrival times, which are higher than modeled in the down-glacier direction but are lower than modeled in the up-glacier direction (LINE2 and LINE4 in Figs 6 and 7). In this subsection we briefly evaluate whether or not intrinsic anisotropy is observable in our data. Seismic P-wave speed is strongly anisotropic in single-crystal ice (see Gusmeroli and others, 2012, for a recent review). When the wave propagates parallel to the c-axis it propagates ∼5% faster than when it propagates parallel to the basal plane (Bentley, 1972; Rothlisberger, 1972; Blankenship and Bentley, 1987; Anandakrishnan and others, 1994; Gusmeroli and others, 2012). In a polycrystalline aggregate of glacier ice this seismic anisotropy is possible when the majority of crystals are aligned along certain preferential directions. A response such as the one observed in our dataset (slower wave from up-glacier) can be explained by a fabric style which is not typically observed in glacier ice. Another factor that could potentially contribute to along-flow differences is variations in the percentage of water. A slower wave from up-glacier implies higher water content up-glacier from the center of the VSP experiment. Pettersson and others (2004) found distinct patterns of water content at the CTS, with lower and higher water content on either side of the glacier center line. We thus may expect across-glacier rather than along-glacier variability. Additionally the horizontal and vertical variability documented by both Pettersson and others (2004) and Gusmeroli and others (2010b) still lies within (/>_w<1%. For these reasons we do not believe that the observed azimuthal variation can be explained with different cj) _w.

Inclinometry data from the borehole show very little deviation from the vertical. The maximum deviation is <1 m and the main deviation is towards the northeast. The upper 30 m of the borehole exhibits a deviation of <0.5 m; a maximum deviation of 1 m is observed at the base of the borehole (80 m depth), resulting in a 0.01 ms change in travel times from a truly vertical borehole. Deviation of ∼1 m has a negligible effect (∼0.01 ms) on the arrival times. We used differential GPS (dGPS) to measure differences in topography of ∼2 m for the greatest offset distance covered in our study (30 m). It thus follows that a receiver located at depth z is in reality located at depth z + 1 m for the wave propagating down-glacier and at depth z – 1 m for the wave propagating up-glacier. For z = 40 m we obtain travel-time differences of ∼0.5 ms between the two directions. Topography might thus explain most of the observed apparent travel-time anisotropy, because singular profiles (Fig. 6) and mean misfits are <1 ms (Fig. 7). It is therefore likely that the observed misfits are due to topography instead of real physical anisotropy.

Conclusions and Outlook

We have studied the seismic P-wave structure of the upper ablation area of Storglaciaren, and found that the seismic velocity structure of the ice is largely dominated by the presence of air. P-wave speed typically increases with depth from 3700 m s–¹ at the surface to 3760 m s–¹ at 80 m depth; this change is almost entirely explained by a reduction in air content from 3% at the surface to <0.5% at depth. Changes in P-wave speed due to water content are very small (∼10ms^–1). The seismic CTS is characterized by a 0.3% decrease in speed, about ten times lower than the decrease in radar-wave speed observed for the radar CTS (Gusmeroli and others, 2010a). These results highlight the importance of including air content when inferring physical properties of temperate glacier ice; if air is neglected, a two-phase (ice/water) interpretation of relatively low P-wave speed in cold, water-free ice could lead to the conclusion that the ice was temperate with some volumetric percentage of water.

Our modeled arrival times match the observed arrival times to within 0.25 ms in the transverse direction, corresponding to a seismic speed uncertainty of 50ms–¹. This uncertainty in constraining the seismic speed is high for detecting water content. Future studies should focus on reducing uncertainties in seismic speed determination as well as applying a joint inversion of radar-seismic datasets, as in Endres and others (2009). We believe that future application of VSP to glaciological problems, perhaps combined with the use of explosives to increase the signal-to-noise ratio and increase sampling frequency, might enhance important advances in glacier geophysics. These advances include a detailed analysis of the shear waves collected here (as they are more sensitive to water than P-waves), and ongoing efforts to seismically infer englacial temperature (Gusmeroli and others, 2010a; Peters and others, 2012), model ice rheology and, ultimately, obtain reliable datasets that can be useful for modelers when predicting the future dynamics of mountain glaciers and ice sheets.