Radar-Sounding Characterization of the Subglacial Groundwater Table Beneath Hiawatha Glacier, Greenland.

Radar-sounding surveys associated with the discovery of a large impact crater beneath Hiawatha Glacier, Greenland, revealed bright, flat subglacial reflections hypothesized to originate from a subglacial groundwater table. We test this hypothesis using radiometric and hydrologic analysis of those radar data. The dielectric loss between the reflection from the top of the basal layer and subglacial reflection and their reflectivity difference represent dual constraints upon the complex permittivity of the basal material. Either ice-cemented debris or fractured, well-drained bedrock explain the basal layer's radiometric properties. The subglacial reflector's geometry is parallel to isopotential hydraulic head contours, located 7.5-15.3 m below the interface, and 11 ± 7 dB brighter than the ice-basal layer reflection. We conclude that this subglacial reflection is a groundwater table and that its detection was enabled by the wide bandwidth of the radar system and unusual geologic setting, suggesting a path for future direct radar detection of subglacial groundwater elsewhere.

Introduction

Recent airborne radar‐sounding surveys revealed a 31 km‐wide impact crater beneath Hiawatha Glacier, part of the northwestern Greenland Ice Sheet (Kjær et al., 2018). These radar data were collected with a new ultrawideband (UWB) radar sounder that revealed the glacier's bed topography and internal structure in unprecedented detail (Wang et al., 2016; Yan et al., 2017). Within these data, Kjær et al. (2018) also identified a distinct reflection beneath the ice–bed interface that was unusually flat and specular, which they hypothesized to be the groundwater table. However, this observation has yet to be confirmed by radiometric or hydrologic analyses, and the intervening material sandwiched between the glacier and this reflection was not characterized.

The unique geologic setting of a subglacial complex impact crater could be partly responsible for this reflection. While the deglaciated region immediately adjacent to Hiawatha Glacier (Inglefield Land) is composed of highly metamorphosed Paleoproterozoic rock (Kjær et al., 2018), unconsolidated impact breccias are expected to be widespread within the crater floor surrounding the central uplift (Osinski & Pierazzo, 2013). Debris‐rich ice is observed outcropping along the base of ice cliffs along the western margin of Hiawatha Glacier, and Kjær et al. (2018) hypothesized that basal material is being actively entrained within Hiawatha Glacier, based on the UWB radar sounding data. A possible second large subglacial impact crater has also been identified in Greenland (MacGregor et al., 2019), but no comparable subglacial reflection was reported there, although the ice is thicker and that structure is likely older.

Radar sounders are widely deployed to study subglacial and englacial water bodies (e.g., A. Wright & Siegert, 2012; Chu et al., 2018; Jordan et al., 2018; Kendrick et al., 2018; Oswald et al., 2018). This method often succeeds because the complex permittivity contrast between ice and water is large and because such bodies are often specular reflectors (e.g., Schroeder et al., 2015). Subglacial lakes are now regularly found beneath the Antarctica and Greenland ice sheets and inform our understanding of their subglacial hydrology (A. Wright & Siegert, 2012; Livingstone et al., 2019). Groundwater sources are also an important component of glacier hydrology; they can drive water into till, elevate porewater pressures, reduce shear strength and significantly influence ice‐sheet dynamics (Boulton et al., 1995; Gooch et al., 2016; Key & Siegfried, 2017; Siegert et al., 2018). However, relatively few studies have reported detecting the subglacial groundwater table (Christoffersen et al., 2014; Mikucki et al., 2015) despite ground‐penetrating radar surveys being well‐established as a method for identifying groundwater in deglaciated environments (e.g., A. Neal, 2004; Woodward & Burke, 2007).

Radar‐sounder designs range from ground‐based impulse and frequency‐modulated continuous‐wave systems to multi‐channel chirped airborne systems (Li et al., 2019). The UWB chirped radar system developed by the Center for Remote Sensing Ice Sheets and deployed by Kjær et al. (2018) is a relatively new version of the Multi‐channel Coherent Radar Sounder (MCoRDS v5) characterized by a much larger bandwidth than previous versions, weaker sidelobes and a higher signal‐to‐noise ratio (SNR) (Rodríguez‐Morales et al., 2013; Wang et al., 2016). The technical advances of MCoRDS v5 raise the possibility that its capabilities alone are what enabled the detection of the hypothesized subglacial groundwater table.

A reflection from a subglacial groundwater table ought to posses a radiometric signature distinct from off‐nadir bed reflections, because the dielectric contrast that induces the reflection should be due to the contrast between unsaturated and saturated sediment, rather than between ice and more typical subglacial interfaces (marine sediment or bedrock). This reflection's subglacial depth should also be consistent with the predicted hydrology of groundwater flow through such systems. If the interface is indeed a water table, then it should also be conformal to isopotential contours of hydraulic head (Flowers, 2015; Rutishauser et al., 2018; Wright et al., 2008). Here we test the hypothesis that the subglacial reflection at Hiawatha Glacier is indeed that of a groundwater table using both radiometric and hydrologic analyses. These tests inform an assessment of the cause of this reflection, prospects of its detectability elsewhere and provide a framework for future investigations of subglacial groundwater systems beneath glaciers and ice sheets.

Data and Methods

The MCoRDS v5 data used in this study were collected in May 2016 (Kjær et al., 2018). The radar system is described in detail by Wang et al. (2016) and consists of three eight‐element arrays, operating between 150–520 MHz at a 10 kHz pulse repetition frequency. These arrays were mounted on the Alfred Wegener Institute's Polar 6, a Basler BT‐67 aircraft. After pulse compression and synthetic aperture processing, the data have a vertical (range) resolution of 0.5 m and an along‐track (azimuth) resolution of 15 m. Figure 1 shows the four flight tracks and radargrams from this survey where the putative groundwater table was detected. The peak power of the putative groundwater‐table reflection was extracted using a local depth window that was selected manually to bound this reflection. This depth window was also used to re‐track some regions where the original ice–basal layer peak power picks corresponded with the hypothesized groundwater table.

Figure 1

Radargrams from the 2016 Hiawatha Glacier survey and location of the hypothesized groundwater table. Panels a–d show full radargrams across the crater (20160517_03_008, 20160512_02_009, 20160516_02_006 and 20160512_02_007, respectively). Panels e–h are portions of a–d zoomed in on the hypothesized groundwater table. The vertical scale bar in e–h corresponds to the depth range in ice (assuming ϵ′ = 3.15), not sediment or rock. Panel i shows the bed topography (Morlighem et al., 2017) overlain by all the 2016 survey flights over the crater (white), with black segments (panels a–d) representing those with a potential subglacial groundwater table beneath a portion thereof (colors, panels e–h).

Radiometric Analysis

Our radiometric analysis aims to test whether the received power from the ice–basal layer and hypothesized basal layer–groundwater interface, along with the dielectric loss within the intervening basal layer, are consistent with subglacial groundwater table. We assume that this system can be represented by a three‐layer dielectric model, where the relative complex permittivity, is uniform within each layer, where j 2 = −1 and tan δ = ϵ″/ϵ′ is the loss tangent (dielectric loss). Two different models for the basal layer are explored in parallel, and in each case the basal layer is described by a three‐component mixture. We show below that the combination of reflectivity difference between the ice–basal layer and subglacial reflections and tan δ within the basal layer provide dual constraints upon .

For our first model, the top layer is polar ice with an assumed permittivity expressed in the form (Fujita et al., 2000; Peters et al., 2005). The middle layer sandwiched between the ice–basal layer interface and the putative groundwater table has an unknown permittivity and is assumed to be a mixture of granitic sand, groundwater‐saturated till and ice (with initially unknown volume fractions). The bottom layer (groundwater table) is modeled using the dielectric properties of thawed, groundwater‐saturated granitic till (Christianson et al., 2016). These assumed layer compositions were selected a posteriori by testing a range of plausible subglacial permittivity values to maximize overlap between the dual constraints. We note that the dielectric contrast between the middle layer and the putative groundwater table is so great that our conclusions are not significantly affected by the plausible range of permittivity values for the middle layer.

The second model again includes a top layer of polar ice and an unknown middle layer, but differs in assuming the bottom layer is a groundwater aquifer in porous or fractured granitic rock (rather than an unsorted till) with higher water content and thus higher permittivity (). In this case, the middle basal layer is assumed to be a mixture of fractured bedrock, water and air (rather than ice).

These two distinct models effectively test whether detection of the subglacial reflection is due to a thermal transition from frozen to thawed material (model 1), or due to a hydraulic transition from drained to saturated bedrock (model 2). Only their differences in the assumed composition of the bottom layer affect the resulting mixture ratios for the sandwiched basal layer. We consider these two cases to be the most plausible, with the primary goal of testing for the existence of the groundwater table rather than robustly identifying nature of the middle layer. For simplicity, the following sections use equations with subscripts and descriptions for model one only, but the analysis for model two is the same except using groundwater gw instead of groundwater‐saturated till gwt.

The difference in received power between the basal layer–groundwater till (b–gwt) and the ice–basal layer (ice–b) reflectors is given by where Δ[P] = [P ] − [P ], Δ[R] = [R ] − [R ] is the reflectivity difference between the basal layer–groundwater and the ice–basal layer reflections, [L ] is the dielectric attenuation within the basal layer material, and the notation [X] = 10 log 10 X is used for power in decibels. For this relation, birefringence loss and the radar system performance cancel out (Fujita, 2006; Haynes, 2020; Matsuoka et al., 2012). VHF birefringence within subglacial materials has not been reported so we neglect this possible confounding factor (Jordan, et al., 2020). Given the small traveltime differences and plausible range of permittivities, suggesting a basal layer ∼10 m thick, the difference in geometric spreading loss between the ice–basal layer interface and groundwater interface is negligible and also ignored. To assess the potential effect of interface roughness upon reflection scattering loss and how it might impact interpretation of Δ[R], we compared the spread of the reflectivity distributions for the ice‐basal layer interface and the basal layer‐groundwater interface (see supplement) (Grima et al., 2019; Jordan et al., 2016).

Spatial variation in the thickness of the basal layer between the ice–basal layer interface and the putative groundwater table can be further exploited to estimate the mean dielectric attenuation rate within this layer (Campbell et al., 2008). This regression assumes that (linear) power decays exponentially with travel time t (equivalent to layer thickness for uniform ), resulting in a linear relationship between Δ[P] and t. This method also assumes that tan δ is uniform, that volume scattering within the layer is negligible, and that roughness‐induced scattering losses and [R ] are uncorrelated with t. By neglecting volume scattering within this middle layer, our estimate of water content in the basal layer is a conservative upper bound. [L ] is obtained for each along‐track sample using the regression slope (−Δ[P]/Δt) and Δ[R] is obtained using Equation 1. The loss tangent of the sandwiched basal layer is thus where λ is the radar wavelength in the vacuum (0.9 m), c is its speed in the vacuum and Δt is the two‐way travel time (Campbell et al., 2008). This approach differs from the typical procedure to estimate englacial attenuation rates, where power is regressed against ice thickness for an assumed value of (Jacobel et al., 2009). The rationale for our approach, which was originally applied to the subsurface of Mars (Campbell et al., 2008), is that we cannot assume a value for a priori.

The derived values of Δ[R] and tan δ provide two independent constraints upon and hence its composition. To relate these constraints to , we consider a three‐component mixture of granite (considering a range of granitic sand to rock permittivities), ice, and groundwater‐saturated till using a power‐law mixing model of the form where is the complex permittivity of granite and ϕ , ϕ , ϕ are the fractional volumes of granite, ice and groundwater‐saturated till, respectively (Wilhelms, 2005; Nerozzi & Holt, 2019). We assume γ = 3 following Looyenga (1965) and a range of values for between 5(1 − j6.8 × 10−5) and 9(1 − j0.068), with a mean value of 7(1 − j0.034). These values were determined by converting electrical conductivity σ values for granite from between 10−5 − 10−2 S m−1 (Bogorodsky et al., 1985) and considering a real permittivity range of 5–9 (Martinez & Barnes, 2001; Nerozzi & Holt, 2019), i.e., where f is the radar center frequency (335 MHz), ϵ 0 is the permittivity of the vacuum and to 9. Expected values of Δ[R] and tan δ were modeled using Equation 3, assuming a specular Fresnel reflection for Δ[R], and evaluated for all possible fractional combinations of ϕ , ϕ and ϕ to produce a ternary diagram for each constraint.

Hydraulic Analysis

Hydraulic head is a measure of liquid potential and its spatial pattern determines where groundwater flows (Freeze & Cherry, 1979). Here, we apply a traditional Darcian approach of evaluating the aquifers. We assume flow in one direction and a homogeneous aquifer in both models. For model 1, the total potential at the potentiometric surface (the elevation of the water table) is evaluated assuming the groundwater till aquifer is confined and partially driven by pressure from the top and middle layers (frozen basal layer and the overlying ice sheet). Model 2 considers flow through fractures dominated by gravity. While flow through such fractures is not well described by a simple Darcian model, for our purposes these simple hydraulic models enable the examination of the range of depths and shapes of the groundwater potential compared to the groundwater echo, rather than robustly modeling or supporting either of the two hydraulic scenarios.

To determine isopotential contours of hydraulic equilibrium, we calculate the hydraulic head as a function of elevation head and pressure head (neglecting velocity head) as where z is elevation, ρ is density and g is the acceleration due to Earth's gravity (Rutishauser et al., 2018; Shreve, 1972). Because the basal layer is thin (∼10 m), uncertainty in the density of the basal layer does not significantly affect the hydraulic head calculation and is neglected. We assume the basal layer density is comparable to granite (ρ  = 2700 kg m−3) used in prior studies of northwest Greenland and note that granitic sand at lower densities would produce similar results (Corbett et al., 2015; Vermassen et al., 2019). The density of groundwater‐saturated till (ρ ) and groundwater (ρ ) are both assumed to be 997 kg m−3 assuming the till will not flow (see supplement). To bound these two end members, both possibilities are shown in Figure 2.

Figure 2

Radiometric data analysis for profiles shown in Figure. 1e–h. Panels a–d show power loss in the basal layer material versus two‐way travel time. The loss tangents are obtained from the regression slopes. Reflectivity difference between the basal layer–groundwater and ice–basal layer reflectors are shown in e‐h and isopotential hydraulic head contours for the putative groundwater table in i‐l.

Radar System Analysis

We evaluate the performance parameters of MCoRDS v5 against those of other commonly deployed radar sounders to address whether MCoRDS v5 itself is responsible for the detection of a subglacial reflection. The characteristic bandwidths, center frequency, pulse length and windowing techniques are incorporated to generate and compare the sidelobe patterns of HiCARS (Peters et al., 2007), MCoRDS v3 (Shi et al., 2010) and MCoRDS v5 (see Table S1 for radar‐system parameters used). This comparison tests whether the observed subglacial reflections are likely to be “visible”, or stronger than the sidelobes from basal layer echoes.

Results

For our radiometric analysis, the four profiles where the subglacial reflection was detected were initially analyzed separately (Figure 1). Best‐fit loss tangents range between 0.0102 and 0.0128 (Figure 2). The reflectivity difference between the ice–basal layer and basal layer–groundwater reflection, Δ[R], and hydraulic head are shown for each track section in Figure 2. The mean reflectivity difference (Δ[R]) for all four profiles is 11.1 ± 6.8 dB, accounting for both the standard deviation of each segment and the propagated error in [L ] from the regression slope. The four reflectivity distributions all satisfy Lilliefors and Jarque‐Bera tests for normality, and their mean spread (one standard deviation about Δ[R]) is 6.0 dB.

To estimate the subglacial material composition, we first used the arithmetic mean to combined the individual profile‐mean estimates of tan δ and Δ[R], yielding tan δ  = 0.0115 ± 0.0013 and Δ[R] = 11.1 ± 6.8 dB. The regions of the ternary diagrams consistent with these estimates are shown in Figure 3. The upper and lower bounds for the volume fractions consider the intersection of the outlined regions in Figures 3b and 3c and Figures 3e and 3f, which account for the full ranges of uncertainty in tan δ , ΔR and possible complex permittivity values of granite and groundwater till.

Figure 3

(a) Three‐layer dielectric model of ice, a basal layer (constrained by the radiometric analysis), and groundwater‐saturated till. (b) Ternary diagrams for tan δ and (c) ΔR. The second row shows (d) a second hypothesized three‐layer dielectric model of ice, bedrock, and groundwater. (e) Ternary diagrams for tan δ and (f) ΔR with respect to the second model. Black outlined regions show most likely basal layer volume fractions, and dotted white outline showing the overlapping area of both the loss‐tangent‐ and reflectivity‐analysis probability regions. These ternary diagrams assume , and (Christianson et al., 2016), but the outlined regions encompass permittivity and conductivity ranges (Bogorodsky et al., 1985; Martinez & Barnes, 2001; Nerozzi & Holt, 2019) and (Christianson et al., 2016).

The first hypothesized model (Figure 3a) results in material volume fractions of ϕ  = 16 ± 9%, ϕ  = 74 ± 10%, ϕ  = 10 ± 7%. Substituting these volume fractions into Equation 3, we derive . Using the estimates for ϵ ′ and the mean travel times (Figures 2a–2d), the mean thickness of the basal layer between the ice and groundwater table, averaged over the four sections, is 13 ± 4.7 m. Both this loss and reflectivity analysis indicate the presence of a debris‐laden basal ice layer above a groundwater table. The second hypothesized model (Figure 3d) resulted in material volume fractions of ϕ  = 1 ± 1%, ϕ  = 39 ± 23%, ϕ  = 60 ± 23%, bedrock permittivity of , and indicate the presence of drained, fractured bedrock 9.8 ± 2.8 m thick. While the two models differ, both are consistent with a groundwater table located ∼10 m below overlying material of either (1) ice‐cemented debris or (2) drained, fractured bedrock.

Equipotential hydraulic head lines were compared against the groundwater interfaces for all segments, and the interfaces often followed isocontours (Figure. 2i‐l). Thus, the interfaces are qualitatively consistent with a groundwater table in hydrologic equilibrium. The deviations from these isopotential lines appear to be related with deviations of the flight tracks from the local ice‐flow direction, especially toward the northwestern margin of the ice sheet. The lateral extent of the groundwater system is ∼15 km2.

Analysis of sidelobe patterns shows the potential of MCoRDS v5 and other systems to detect similar subglacial groundwater tables (Figure 4). MCoRDS v3 could plausibly detect nearly all Hiawatha subglacial groundwater reflections, but many would be on the edge of detectability for HiCARS due to its narrower bandwidth (15 MHz). This is because the subglacial groundwater table reflections only slightly exceed the sidelobes generated by basal layer echoes from these two systems. However, for MCoRDS v5 these subglacial reflections are consistently tens of decibels higher than the sidelobes. Therefore, the combination of high SNR and wide bandwidth–resulting in faster sidelobe fall‐off–is likely a significant factor in explaining why the subglacial groundwater table was detected beneath a portion of Hiawatha Glacier. The lack of detection of the subglacial groundwater table in other regions of the crater could be because: 1. No groundwater table is present there; 2. The basal layer is insufficiently frozen or drained to permit substantial radio‐wave penetration; 3. The groundwater table is not sufficiently contiguous to identify in the radargrams; or 4. The interface is too deep to be detected.

Figure 4

Antenna patterns for three radar systems (Table S1) compared against the putative groundwater echoes as a function of traveltime through the basal layer. The groundwater echoes are shown as in Figures 2a–2d. For each radar system, the potential detectability of the any echo increases with the difference in power between the echo and the radar system's antenna pattern, for example at a traveltime of 0.1 μs, the putative groundwater echoes are ∼70 dB above the noise floor of MCoRDS v5, but < 15 dB above that of HiCARS.

Discussion and Conclusions

Both our radiometric and hydrologic analysis are consistent with the anomalous subglacial reflection originating from a groundwater table beneath either a well‐drained or partially frozen basal layer within the Hiawatha impact crater floor. Our radiometric analysis shows the groundwater‐table reflection is typically over 10 dB stronger than the overlying ice–basal layer reflection, strongly indicative of the presence of water‐saturated material, i.e., a groundwater aquifer.

Our first hypothesized model (Figure 3a) indicates that an ice‐cemented debris layer lies above thawed, saturated groundwater till, consistent of a mixture of groundwater till, granite, and ice, with ratios of approximately 16%, 10%, and 74%, respectively. In this model, water can exist both above and below the aquitard of the frozen basal layer, and the low attenuation rate of the basal layer is the result of its thermal state, i.e., the pores are filled with ice rather than water. The underlying thawed layer might also be trapped by frozen layers above it, a feature observed in firn hydrology (Koenig et al., 2014; Chu et al., 2018). This layer could be liquid because it is confined and pressurized (Steinbrügge et al., 2020), due to refreezing, heat advected by subglacial water flow or higher salinity (Rutishauser et al., 2018).

Our second hypothesized model (Figure 3d) indicates a basal layer of porous, well‐drained rock above the groundwater table consisting of 1% groundwater, 60% granite, and 39% ice. In this case, the low‐loss basal layer is the result of efficient vertical drainage rather than freezing. The impact should have produced a thick layer of impact breccia, which would be permeable and conducive to rapid subglacial drainage from the overlying glacier into a groundwater system, and this second model indirectly assumes that this layer is still present. Thus, the unique detection of this subglacial groundwater table could be in part due to the uniqueness of its geologic setting.

The radar profiles where the groundwater table was detected are found close to each other, within the northwestern section of the crater (Kjær et al., 2018). A map of the ice–basal layer reflectivity from the 2016 survey indicates increasing relative reflectivity from the southeastern corner of the crater toward its northwestern corner (Fig. s3). Our hydrologic analyses are also consistent with this drainage pattern, in that deviations of the observed groundwater table depth from equipotential hydraulic contours show a pattern of decreasing pressure gradient toward the northwest. This pattern suggests that groundwater is indeed flowing through the crater toward the ice‐sheet margin in the same direction indicated by Kjær et al. (2018). Further investigation into character and flow of the subglacial and groundwater hydrology of the Hiawatha Glacier region of Greenland will require more sophisticated modeling, such as considering variations in hydraulic conductivity and intrinsic permeability.

As ​unique as the geologic setting of Hiawatha Glacier may be, observation of its groundwater was also partly enabled by the large bandwidth and SNR of the MCoRDS v5 system. Additional surveys by similar wideband sounders over other sites with known or hypothesized groundwater, or surveys of Hiawatha Glacier region by other radar sounders could validate the potential for wider applications of this work. This conclusion raises the possibility that other subglacial groundwater systems could be mapped using wideband radar sounders, providing new insights into the poorly understood role of groundwater in the subglacial hydrology of Greenland, Antarctica and other glaciated regions (Key & Siegfried, 2017; Siegert et al., 2018; Williams et al., 2020).

Supporting Information S1

Click here for additional data file.