Introduction

The Paakitsoq area of the Greenland ice sheet has seen a wide variety of focused scientific work for more than half a century, especially following establishment of the Swiss Camp station (69.57^˚ N, 49.31˚W; 1175 ma.s.l.) just below the equilibrium line. Results reported by Reference Van der VeenZwally and others (2002) link glacier acceleration at Swiss Camp to enhanced drainage of meltwater through crevasses and moulins during the summer melt season, suggesting that local hydraulically induced basal sliding can explain seasonal variability in velocity, and implying that thinning in this region can be attributed in part to dynamically induced changes. Reference Thorning, Thomsen and HansenVan de Wal and others (2008) supported this conclusion with a long-term global positioning system (GPS) record of velocities, which showed seasonal velocity increases along the Kangerlussuaq transect (K-transect), a similar ice-marginal area but to the south of Paakitsoq. Das and others (2008) link speed-ups at Swiss Camp with sudden meltwater lake drainage. However, the Paakitsoq region is close to, and may be influenced by, Jakobshavn Isbræ, a fast-moving tidewater outlet glacier whose marine terminus is currently in extended retreat. This implies that the collapse of Jakobshavn may also be influencing the observed thinning at Swiss Camp (Reference Hempel, Thyssen, N., Clausen and MillerJoughin and others, 2004; Reference KitanidisKrabill and others, 2004), a process which may have implications for both the measured velocities and the formation and propagation of crevasses at Swiss Camp, particularly with respect to whether or not crevasses may propagate to the bed. High-resolution information on the basal and surface topography of this region can help to resolve some of these issues. Such data are also important for modelling ice-sheet hydrology and as input to higher-order models of ice dynamics.

The current project was initiated by a feasibility study for a hydroelectricity project intended to feed Ilulissat and requiring detailed analysis of the hydrology and future evolution of the ice-sheet margin. Plans for hydropower development were originally raised during the 1973 oil crisis and recently revisited (Reference AhlstrømAhlstrøm, 2007), but there has not previously been a well-constrained and detailed study of basal topography and ice thickness over the whole region. An early radioglaciological investigation by helicopter using visual navigation in 1985 and 1986 (Reference Thorning and HansenThorning and others, 1986; Reference Thomsen, Thorning and BraithwaiteThorning and Hansen, 1987) produced some useful insights, but was hampered by poor positioning information. In this paper, we present details of basal topography and ice thickness in the region, based on radar measurements, and ice surface elevation based on laser altimeter measurements, both flown in summer 2005. We then apply these in a hydrological model of this part of the ice margin.

Data Acquisition and Processing

Radar data were collected in 2005 by an airborne survey over the region (see Fig. 1) from a Twin Otter aircraft. The region surveyed covers >80km along the ice margin and has a total area of around 2700 km². The map in Figure 1 shows the location of the study area and the density of measurements: the between-track distance was about 1 km near the margin, increasing to around 3 km away from the margin. Aircraft altitude above the ice surface varied from as low as 20m to almost 500m dependent on local conditions.

Fig. 1. Map of the Paakitsoq region illustrating the flight-lines along which radar and lidar data were collected in summer 2005. The coordinates are zone 22 Universal Transverse Mercator (UTM) projection.

The 60MHz ice-sounder developed at the Danish Technical University (DTU) recorded ice surface, ice thickness and subglacial profiles. The ice-sounder uses radio pulses to measure the distance between the ice-sounder antenna and the ice surface and ice bottom (bedrock) respectively, and from these the ice thickness is calculated. The ice-sounder records 50 profiles of 4096 samples per second, and data acquisition consists of both transmitting pulses at a frequency of 5 kHz (sampling in the flight direction) and sampling the returned echo at 75 MHz in range, producing 4096 samples per transmitted pulse. Coordinated universal time (UTC) is recorded by the ice-sounder instrument and is used for positioning the ice-sounder observations with on-board differential GPS units and an inertial navigation system (INS). An example of a radargram obtained by the ice-sounder is displayed in Figure 2a where the horizontal direction represents the time with a spacing of 320 ms per line (i.e. 22.4m spacing at aircraft velocity 70ms^–1). The vertical direction shows propagation time of the radar pulse with a spacing of 80 ns per line of the radar pulse. This represents vertical distances but not by simple scaling because the speed of light within the ice sheet is lower than in free air.

Fig. 2. Radargrams illustrating bed and surface topographic traces taken from two different locations within the region. (a) The initial radargram, obtained from the sounder. (b) The processed radargram, with the digitized traces of the surface layer shown in blue and the bottom returns shown in red.

Reflection and attenuation within the ice sheet reduces the strength of the returned echo. Substantial processing is therefore required to produce a radargram that makes detection of the ice-sheet bottom echo possible and is done both online during acquisition and offline using software developed at the Microwave and Remote Sensing section, DTU-Space. Processing consisted of coherent filtering to suppress unwanted signal components and noise, followed by incoherent averaging to reduce speckle noise. Figure 2b shows the radargram produced after processing, with the bottom and surface layers clearly represented. Based on known instrumental tolerances, the accuracy with which elevations can be retrieved from the radar data is conservatively estimated to be ±80 m.

Following processing, the basal and surface layers were identified and digitized using ReflexW# software. The two-way travel time of the radar pulse from the plane to the basal and surface layers and the elevation of the plane above the World Geodetic System 1984 (WGS84) ellipsoid was used to calculate the surface elevation (Equation (1)) from which the ice thickness was derived (Equation (2)).

The surface elevation, H _s, is determined using the measured elevation of the plane, H _p, the two-way travel time to the surface, t _s, and the velocity of light through air, v _air . The factor 0.5 is applied to give the single travel time. In this study we used the lidar-obtained data to calculate the surface elevation.

The ice thickness, H _i, is derived from t _b, the two-way travel time to the bed, and v _ice, the velocity of light through ice. The bed elevation, H _b, is then simply calculated by subtracting the ice thickness from the surface elevation:

All elevations are given relative to the standard WGS84 ellipsoid, as the local correction to sea level is not well known. The speed of the radar pulse through air was assumed to be 0.3mns^–1 and through ice 0.169mns^–1 (Reference Bogorodsky, Bentley and GudmandsenBogorodsky and others, 1985; Reference Hagen, Etzelmüller and NuttallHempel and others, 2000). As the surface layer was not easy to identify at many points for which basal returns were clear, a further correction was applied using the laser surface altimetry measurements to calculate aircraft altitude above the ice. During analysis of the radar data, returns from ice were distinguished and separated out from the returns from water and bedrock, in order to create digital terrain models (DTMs) of both glaciated and unglaciated terrain.

Surface elevation data were derived from laser altimeter measurements made concurrently with the radar measurements and using the same INS and GPS set-up to derive positioning information. A Riegl scanning laser was used to make the elevation measurements, which provides crosstrack scans with a range accuracy better than 5 cm. The laser operates in the near-infrared wavelength band, and has a scan angle of 60^˚, giving a swath width similar to the flight elevation above the ground. Post-processing to combine positional and elevation information in vector format files was carried out at the Geodynamics Department of DTU-Space.

Positioning information and the calculated elevations were merged and interpolated across a grid to create a DTM using the kriging technique. A block kriging algorithm (Reference CarrCarr, 1995; Reference Joughin, Abdalati and FahnestockKitanidis, 1997), based on a rational quadratic model and incorporating an adjustment for the nugget effect (i.e. the measurement error in positioning, ±80m for radar data and ±1m for lidar data), was developed using plotted semivariograms in eight directions to take into account a priori knowledge of bed topographic variability (Reference AhlstrømAhlstrøm and others, 2002, Reference Ahlstrøm, Mohr, Reeh, Christensen and Hooke2005, Reference Ahlstrøm, Mottram, Nielsen, Reeh and Andersen2008). The model was then used to create a grid of basal topographic elevations over the whole area. The DTM also incorporated an additional topographic dataset of the unglaciated foreland, constructed from aerial photographs taken in 1985 by the Geological Survey of Denmark and Greenland (GEUS), (personal communication from H. Jepsen, 2008), as well as bathymetric data from two ice-marginal lakes supplied by The Greenland Survey (ASIAQ) (personal communication from D. Petersen, 2008). The interpolated bed elevation grid is 1 km by 1 km, a compromise between the high resolution along track and the relatively low resolution across track, particularly in areas further away from the margin.

Surface data from laser altimeter measurements made concurrently with the radar measurements were used to construct a surface terrain model. In this case, a linear model was found to be the best fit when applying the kriging interpolation method to the surface data. The ice-surface DTM was then also merged with the DTM of the unglaciated terrain in order to develop a hydrological model of the surface run-off.

The hydrological model, developed using the Rivertools# software, predicts the location of the watershed and river channels within each drainage basin, based on surface slope, curvature and relief of the modelled surface. Further modelling within the package was used to delineate drainage basins and predict meltwater discharge into the lakes proposed as sources for the hydropower project. Delineation of hydrological ice-sheet drainage basins is complicated because the watershed is not simply defined by the surface drainage, but must include a formulation of the englacial water routing as shown by Reference Price, Payne, Catania and NeumannShreve (1972). This formulation is simplified by the assumption of Reference BjörnssonBjörnsson (1982) that all meltwater reaches the bedrock through moulins and crevasses and drains along the base of the ice sheet, which is assumed to be impermeable. This is a reasonable assumption for large-scale flow in regions with basal ice near the pressure-melting point and has been applied by, for example, Reference ShreveThomsen and others (1988), Reference Dass.B.Hagen and others (2000) and Reference AhlstrømAhlstrøm and others (2002). Recent observations at Swiss Camp by Das and others (2008) have proven that surface meltwater can penetrate 1 km of cold ice (i.e. ice at subfreezing temperature). Das and others (2008) observed how a large surface meltwater lake drained catastrophically to the basal drainage system through cracks formed by filling existing crevasses with meltwater from the lake. The weight and melting capacity of the meltwater from the lake caused the crack to penetrate all the way to the bed, a process that has previously been proven to be theoretically possible (Reference Van de WalVan der Veen, 2007). The Paakitsoq region includes the study area of Das and others (2008) and is known to be drained by a well-developed network of surface streams feeding into local moulins, a map of which was produced, based on 1985 data, by Reference ShreveThomsen and others (1988). The application of a simplified model assuming instant meltwater transport to the bed through thick, cold ice is therefore justified by observations. The simplified model implies that the direction of the water flow at the base of the glacier is determined by a water-pressure potential Φ_b, given by

where ρ _w and ρ _i are the densities of water and ice, respectively, Z _b is the bedrock elevation, Z _s the elevation of the ice-sheet surface and g the gravitational acceleration. The last term in the equation is the subglacial water pressure which is proportional to the pressure of the overlying ice, with the k factor ranging from k = 0, corresponding to atmospheric pressure in subglacial channels, to k = 1 for the situation where the subglacial water pressure equals the overburden pressure exerted by the ice. The water at the base of the glacier will flow in the direction of the maximum gradient of the water-pressure potential. Basin delineation and a drainage pattern was then calculated from the potential surface by the RiverTools# package. Thus, this model predicts, based on surface and subglacial topography, the direction of the flow of water beneath the ice.

Results

Figures 3–5 present maps of the ice surface, the basal topography and ice thickness for the Paakitsoq region, based predominantly on the flight-lines shown in Figure 1 as well as the other datasets described above. The red line indicates the ice margin mapped from aerial photographs taken in 1985 (Reference ShreveThomsen and others, 1988). The Jakobshavn outlet glacier, which has since retreated from the terminus position shown, was not covered by this study; nevertheless the surface contours to the north of the basin in Figure 3 suggest that it has a substantial influence on ice flow in the Paakitsoq area. The position of the Swiss Camp research base is shown on both maps. Figure 6 shows the modelled drainage of the area based on the mapped topography, while Figure 7 shows the modelled drainage basins for the Paakitsoq area of the ice-sheet margin, assuming a range of different values for the factor k.

Fig. 3. Map of ice surface topography and unglaciated foreland in the Paakitsoq region. The red line indicates the ice margin mapped from aerial photographs taken in 1985. SA indicates the position of the calving front of Sermeq Akuliarusersua glacier. The area of white space indicates missing data. The map is in zone 22 UTM coordinates; note the black and red dot showing the position of Swiss Camp. Colour scale refers only to the ice surface; for elevations on the unglaciated foreland, see Figure 4.

Fig. 4. Map of basal topography and ice-free terrain in the Paakitsoq region. The red line indicates the mapped ice margin in 1985, based on aerial photographs. SA indicates the position of the calving front of Sermeq Akuliarusersua glacier. White areas indicate missing or poor-quality data. The map uses zone 22 UTM coordinates. The red and white dot marks the position of Swiss Camp.

Fig. 5. Contour map of ice thickness in the Paakitsoq region. The red line indicates the mapped ice margin in 1985, based on aerial photographs. The white dot indicates the position of Swiss Camp. Areas not contoured have insufficient data. The image uses zone 22 UTM coordinates.

Fig. 6. Hydrological model of the surface drainage in the Paakitsoq area. The red line indicates the 1985 ice margin. The model does not incorporate a correction for sea level or for englacial or basal drainage, hence the appearance of rivers following the fjords. Note the two different trends, with surface channels running from north to south in the eastern part and east to west in the western part of the area.

Fig. 7. Delineated drainage basins in the Paakitsoq region, incorporating surface, basal and englacial flow and using a number of different values for the k factor. Borehole measurements suggest a value of 0.70–1.05 is appropriate for the Paakitsoq area. The surface and bed topography measured in 2005 were used as input to the hydrological model. The three lakes marked are currently under consideration as reservoirs for a hydropower project. Note that Swiss Camp is significantly to the east of these delineated drainage basins of the area and therefore does not appear on this map.

Topography and Ice Thickness

Surface topographic contours show the large-scale trend of the ice surface which slopes down towards the margin, and also south towards the Jakobshavn basin. The resolution of the grid does not allow small-scale features such as melt ponds and run-off channels to be discerned, although these are certainly present (Reference ShreveThomsen and others, 1988). There is one obvious artefact in this dataset, an east–west-trending ridge which does not feature in any of the aerial photographs or satellite imagery of the area. This is the product of an erroneous surface dataset, and the data were excluded in the DEMs used to calculate surface hydrology. Other than this feature, the shape of the surface conforms well to that indicated in existing optical imagery of the region.

The basal topographic map in Figure 4 shows that the hilly ice-free terrain continues, at least initially, underneath the current ice cover, with a large plateau-like feature reaching elevations of around 500m in the north of the area. Further inland, however, the topography becomes much lower-lying, with large parts of the area at or below the ellipsoidal datum, by as much as 200 m. Total relief in the area is thus around 700 m. A large valley system with a significant overdeepening extending several kilometres back is shown behind the current calving front of the outlet glacier, Sermeq Akuliarusersua, marked SA in Figure 4. A large trough is also identifiable close to Swiss Camp, and this implies that the location, which is on a bedrock slope, may be experiencing substantial tensile stresses which could significantly increase the chances of deep crevasses forming (Reference KrabillPrice and others, 2008). Both these large features also show on the subglacial topography given by Reference ShreveThomsen and others (1988) and more recently shown in Reference KrabillPrice and others (2008), giving additional confidence in the results. As with the outlet glacier trough, the Swiss Camp feature also trends broadly northeast to southwest; however, note that there is poor data coverage in the area to the north of Swiss Camp. Although the kriging routine has taken into account the resolution of the radar data (±80 m), some of the small irregular features may be artefacts of the kriging and gridding processes. However, comparison with the earlier survey data (Reference ShreveThomsen and others, 1988) also shows similar features in the same location, so it is likely that, at least in part, these features are real.

The ice-thickness map, computed from the surface and bed topographic models, shows that some of the thickest ice in the area is directly to the east of Swiss Camp, largely due to the deep trough identified there, while ice at the camp location itself is around 1200m thick. The rapid thinning close to the margin reflects the relatively steep surface slope of the ice sheet.

Glacier hydrology

The surface run-off in the Paakitsoq area predicted by this initial analysis (Fig. 6) shows a strong bidirectional trend. To the east, the influence of the Jakobshavn basin is evident, with the predicted surface meltwater streams running from north to south, clearly reflecting the surface contours shown in Figure 3. The rest of the area shows an east–west surface drainage pattern, consistent with ice flow from the ice divide westward to the margin. Aerial photographs were not taken during the field campaign of this project, but the predicted drainage pattern based on surface contours is supported by analysis of old aerial photographs carried out by Reference ShreveThomsen and others (1988). Interestingly, this model suggests that meltwater from Swiss Camp will drain out towards the Jakobshavn basin, rather than the ice-sheet margin. This in turn suggests that the measured thinning and seasonal speed-ups at Swiss Camp may, at least in part, be a complex response to the recent collapse and rapid retreat of Jakobshavn Isbræ.

Drainage basins delineated using the full surface, basal and englacial hydrology model are shown in Figure 7. In order to delineate the basins, a value for the k factor, the ratio between ice overburden pressure and basal meltwater pressure, must be specified. Measurements in boreholes drilled through the ice near the ice-sheet margin within the basin late in the melt season point at basal water pressures ranging between 79% and 105% of the ice overburden pressure, corresponding to a value of the k factor between 0.79 and 1.05 (Reference AhlstrømAhlstrøm, 2007). Figure 7 presents results for three values of the k factor, 0.5, 0.7 and 0.9, as a combined drainage area for the three lakes selected for the hydropower project. As the project was mainly concerned with the use of hydropower development, only the drainage basins relevant for the ice-marginal lakes at Paakitsoq were computed. Figure 7 shows that Swiss Camp is outside the catchment of the three Paakitsoq lakes, and this supports, though not conclusively, the interpretation of Figure 6 suggesting Swiss Camp is hydrologically linked to the Jakobshavn drainage basin.

Reference Van der VeenZwally and others (2002) and Das and others (2008) link the velocity variability at Swiss Camp with enhanced meltwater drainage during the summer melt season. Recent modelling work by Reference KrabillPrice and others (2008) suggests that the observed velocity variability at Swiss Camp may also be a response to enhanced meltwater drainage up to 12 km lower down on the ice-sheet margin, which influences the up-glacier velocity through longitudinal coupling. The findings presented here are not in conflict with the hypothesis that meltwater draining to the bed causes such accelerations, but they do suggest that the rates of thinning measured at Swiss Camp may also be a dynamic response to the rapid terminus retreat of Jakobshavn Isbræ.

Conclusions and Implications

The basal and ice-surface topographic datasets presented here are of interest to the glaciological community and will be especially useful as input to higher-order models of ice dynamics as well as to research on the local hydrology of the ice sheet around Swiss Camp.

Furthermore, the results presented here indicate that the Swiss Camp region of the ice sheet has a complex relationship with the Jakobshavn Isbræ outlet glacier. Our results imply that the meltwater from this region is directed through the Jakobshavn drainage basin. Further work modelling the hydrology of the Swiss Camp area in detail is required to firmly resolve this problem and is a focus of future work in this project.