Introduction
Over the past several years, momentum has been building in the glaciological and climatological communities to obtain a new deep ice core from central Greenland for research on past climates and atmospheric compositions. Such a project has received very high priority from numerous national and international planning groups (e.g. National Research Council, 1986; National Science Foundation, 1987). An essential recommendation of all these reports is that ancillary studies be done so that a drill site can be carefully selected in order that the scientific results are maximized and interpretation problems are minimized. Paramount among such site-selection activities are radar sounding and altimetry to determine surface and bed topography, internal layer geometry and whether or not there are areas of basal melting. This paper describes the results of a study done in response to these recommendations.
The Committee for Scientific Planning in Greenland, sponsored by the U.S. National Science Foundation, outlined a general plan for ice-radar sounding of the potential drill-site area in central Greenland (Reference Mosely–Thompson, Gow, Herron, Jezek, Kamb and KhahilMosley-Thompson and others, 1985). This plan recommended a large-scale airborne survey to identify the best general area for a core-hole site, followed by a detailed surface-based survey in this area to select the precise drill location.
Our program was directed only at the airborne part of this plan. Nevertheless, we did assume that our results would be used to isolate this smaller area, approximately 25 km by 25 km in size, for subsequent surface-based work. After the field work was done, it became apparent that there would not be enough time to do these additional studies before a site needed to be chosen. Consequently, considerable effort was put into the data analysis so that determination of the surface and bottom topography would be as reliable as possible.
Design Of The Grid
This project was done in close co-operation with a ground-based study of accumulation rates, surface temperature, snow structure, and ice motion done at the same time by Dr J. Bolzan of The Ohio State University. The size and location of the overall region to be studied, called the “site-selection area”, was the same for both parties, except that the radar grid was slightly larger. The location of this area (Fig. 1) was determined by choosing one in which existing data on these surface variables, as well as the known bed topography, appeared to satisfy as best as possible the criteria for a deep drill site (Reference Langway, Oeschger, Dansgaard, Langway, Oeschger and DansgaardLangway and others, 1985). The size of the area, on the other hand, was simply made as large as possible, within the limitations of the available logistics, so that the eventual choice of a drill site would not be artificially predetermined by the initial, somewhat arbitrary, selection of the study area.
Previous suggestions for a core site were generally centered on the old “Crête” location (lat. 71°07.2′N., long. 37°19.0′W.), where a 404 m core was drilled in 1974. However, the site-selection area was shifted further north to take advantage of apparently smoother bed topography, colder surface temperatures, and less chance of surface melting, as well as a relaxing of the minimum acceptable accumulation rate from 250 to about 200 mm/a (personal communication from W. Dansgaard). After several iterations, a square grid, 150 km along each side, and nominally centered on the old “Summit” site of a 31 m core drilled in 1974, was chosen. Our work indicates that the true summit of the Greenland ice sheet is actually about 31 km north, and 9 km east, of this point.
About a month before the radar flights were made, a camp was established on the ice-sheet surface by The Ohio State University group. The coordinates of this camp, referred to here as “Summit/OSU”, were determined by Transit satellite (geoceiver) fixes and used as the precise center of the radar grid. These coordinates were lat. 72°17′ 38.266”N., long. 37°55′ 18.483”W., and altitude 3260.7 m. relative to the WGS-72 ellipsoid (written communication from J. Bolzan). By repeatedly flying over this point, it was possible to make various corrections to the horizontal and vertical coordinates and to transform them into a standard Earth coordinate system. The old Summit site is about 2.0 km west and 0.9 km south of the OSU camp.
A rectangular Cartesian grid, with its origin at the Summit/OSU camp and sides parallel to the true north-south and east-west directions at the origin, was used for the radar flight lines (Fig. I). The +X coordinate is towards the east and the +Y coordinate is towards the north. Based on the previously known surface topography in the region (Quaternary map of Greenland, compiled by A. Weidick and published by the Geological Survey of Greenland), it was thought that the axis of the ice divide was only about 10–15 west of north, so that the north-south lines would be roughly parallel to the ice divide, and the east-west ones normal to it. This would result in one set of radar profiles, the east-west ones, approximately along the ice-flow direction and another set transverse to it. In reality, however, the actual ice divide is almost exactly along a true north-south direction and so the grid lines are fortuitously aligned with the general ice-flow pattern as closely as they could be without having to resort to a curvilinear grid.
A spacing of 12.5 km between flight lines was used, for several reasons. First, this distance corresponded to about three-four ice thicknesses, typical of the topography wavelengths noted around Dye 3, in southern Greenland, by Reference Whillans, Jezek, Drew and GundestrupWhillans and others (1984). Secondly, it was small enough that a 25 km by 25 km sub-region could be selected for further study using data along lines through the interior of the sub-region as well as along its edges. Thirdly, the spacing was large enough that complicated aircraft turns at the end of each line would not be necessary and the total amount of flying required could still be safely accomplished in 6 d, the maximum aircraft time available for the survey. The spacing of 12.5 km gave 13 evenly spaced lines in each coordinate direction over the entire 150 km square grid.
In order to help constrain the topographic interpolation procedures along the borders of the 150 km radar grid, as well as for other reasons discussed below, the flight lines were extended at least 15 km beyond these borders, resulting in continuous data acquisition along straight lines 180 km long. This 180 km by 180 km “extended” grid is the area over which the data analysis was done, and which is depicted in most of the figures.
Fig. 1. The location of the radar grid in central Greenland. The cent er of the grid is al lal. 72°17′ 38”N., long. 37°55′ 18”W., near the old “Summit” site, and is oriented in a true north-south, east-west direction. The expanded view on the right shows the main and extended grids, the d esired flight-line pattern, a sample flight-line loop and the waypoint numbering scheme. The desired flight path, as indicated by the INS readings, is shown dashed; the actual path, after correclion for INS drift, is shown solid. The point M represents the time of passage over the reference paint. NC. North Central; S. Sondestrøm; and CC, Camp Century.
Radar Equipment
The ice radar used for this study was designed and built by the Technical University of Denmark (Reference GudmandsenGudmandsen, 1976; Reference Skou and SondergaardSkou and Sondergaard, 1976) and used to sound the Greenland and Antarctic ice sheets in the mid to late 1970s. It was installed in the same LC-130 aircraft, XD-03, used in this previous work, but unfortunately lost in a crash in Antarctica 6 months after being used for this study. This radar, known informally as the “TUD radar”, was, and still is, the only one in existence which can penetrate, from the air, the thick ice in central Greenland (Reference Bogorodskiy, Bentley and GudmandsenBogorodskiy and others, 1985).
The TUD radar transmits a pulse with a peak power of 10 kW on a carrier frequency of 60 MHz. The pulse length is selectable, but, since experience in Greenland and Antarctica showed that 250 ns almost always gave the best compromise between resolution and sensitivity, this value was used throughout our work as well. With this pulse length, the radio-frequency band width is 4 MHz and the resolution is about 20 m (minimum separation of adjacent layers which can still be individually resolved). Pulses were transmitted continuously at a 12.5 kHz rate, corresponding to one pulse every 80 μs, or about every 10 mm of horizontal travel at typical aircraft speeds of 120 m/s (235 kt). The overall system sensitivity, including antennas and recording, is 218 dB for the 250 ns pulse length.
The received signal was detected, logarithmically amplified, and electronically differentiated to reduce the dynamic range, sharpen edges and enhance features such as internal layering (Reference GudmandsenGudmandsen, 1975). It was then recorded on heat-sensitive, dry-silver paper using a Honeywell “Visicorder”, the same technique used in the last years of the radar’s use in Antarctica and Greenland (Jankowski, unpublished). These records, visible within minutes of acquisition, were used for all the analyses presented here.
In addition to these analog recordings, a new high-speed digital data-acquisition system designed and built by the U.S. Geological Survey for a separate ice-radar program in Antarctica was also used. Even though merging this system with the TUD radar had not been originally envisioned, its development had to be accelerated to meet the Greenland field schedule, and it had never been tested in the field beforehand, over 97% of the ice-radar data were successfully recorded on digital nine-track magnetic tape. These data will be used for future analyses; they were not used here because a large amount of analysis software needs to be developed first. Further details about this system have been given by Reference Wright, Bradley and HodgeWright and others (in press).
Navigation
All navigation was done using the existing instruments on board the LC-130 aircraft. Horizontal position of the aircraft was determined with a Liiton-51 inertial navigation system (INS), altitude above sea-level with a pressure altimeter, and height above terrain with a radar altimeter (model APN-194). Since this study required flying closely spaced grid lines with a series of separate aircraft missions from an airfield over 1000 km away, it was essential to correct for drift of the INS and for changes of atmospheric pressure with time. The Summit-OSU camp was used to provide the necessary horizontal and vertical control, using techniques discussed in later sections.
Navigation data (latitude, longitude, aircraft altitude above sea-level, and aircraft height above the ice surface) were recorded every second on a Digital Acquisition and Display System (DADS) built by the University of Washington for use on this aircraft (Terry and others, unpublished). The clock in the DADS was synchronized with the one used to time-tag the ice-radar data to allow precise merging of the two data sets during post-processing.
A Lambert conformal projection (Reference NewtonNewton, 1985) was used to transform between the rectangular Cartesian coordinates of the radar grid and the geocentric coordinates (latitude, longitude) used by the INS and the aircraft navigator. This projection was centered on the Summit/OSU camp and used a separation of standard parallels of 1°.
It was recognized from the outset that use of an inertial navigation system and a pressure altimeter to position the aircraft would give only barely acceptable results, and then only when combined with periodic resetting over a known reference point. A Global Positioning System (GPS) was therefore installed on the aircraft, and its data output integrated with the DADS and ice-radar data streams. Despite a successful test of the entire GPS data-acquisition hardware and software on a special flight done in California only weeks before the actual Greenland field work, as well as a successful test in central Greenland the year before, satellite-signal strengths were too low to provide any position fixes at any time during the actual radar missions. We do not have any explanation for this unfortunate and frustrating result, but we are reasonably confident that it was not due to equipment or operator malfunction.
Data Collection
The radar flying was done during a period of 6d, from 29 May to 3 June 1987, based out of Thule Air Base in northern Greenland. Since the one-way transit time from Th�le to central Greenland was over 2 h, we were limited to 3–4 h of actual grid flying on each mission. Fortunately, on 4 of the 6 d, the weather at Th�le was sufficiently good that this time was extended somewhat by using up some reserve fuel. Consequently, all 26 lines of the radar grid were successfully completed, with sufficient time left over to repeat eight of them to fill in gaps in the digital radar data. Thus 34 flight lines were flown, each at least 180 km long, for a total of over 6100 km of profiling within the radar grid.
High-quality ice-radar records were obtained over the entire grid, in both digital and analog data streams, with strong, clear, and unambiguous bottom returns everywhere. Internal layering is very prominent, occurring everywhere down to about 60–65% of the total ice depth, and occasionally to within a few hundred meters of the bed. Examples of the ice-radar data records are shown in Figure 2; the particular flight lines used for this figure were chosen because they run through a region of smooth basal topography near the summit of the ice sheet, and thus could be a potential core-hole site. The internal layering is exceptionally strong here, with continuous layers only 600 m above the bed (80% of the ice depth), and traces of layers only 150 m above the bed (95% of the ice depth).
Very few problems were encountered with acquisition of the navigation and analog radar data throughout the flying. Only two gaps, of 25 and 35 km, occurred in the ice-radar records, resulting in over 99% successful coverage of the extended grid with ice-thickness data. Two longer gaps also occurred, one of 66 km in the aircraft height-above-terrain and one of 166 km in the loss of surface return on the radar records, but, as described later, the essential information was effectively recovered by using data from repeated and/or crossing lines.
The precise location of the geoceiver coordinates of the Summit/OSU camp was marked with flag lines and the aircraft was visually flown directly over this point every 1½-2½ h, upon arrival at, and departure from, the radar grid, as well as once, or sometimes twice, during the middle of the grid flying. Ground personnel radioed the exact instant of passage overhead to the aircraft navigator. This allowed the INS to be reset to the correct coordinates, to well within the resolution of the INS (0.1 min of arc). Throughout the 6 d of flying, the weather over central Greenland was generally excellent and very stable. At no time did ground fog obscure the Summit/OSU camp and all control passes over the reference point resulted in good fixes.
Fig. 2. 2. Examples of two radar records, lines N375W and N250E. Line N250E has been photographically inverted so that both lines have west (W) on the left and east (E) on the right. Thus the horizontal scales. D. which are distance along the flight line in the direction actually flown, run in opposite directions. The symbols S and B at the left end indicate the surface and bottom returns, respectively, and the scale at the right end is the ice thickness, H. Both scales are in units of kilometers. The Y-axis of the grid is at approximately D = 90. The end-of-suppression line normally risible on these records just above the ice surface was very faint on these particular records and cannot he seen in this photographic reproduction.
A series of seven INS navigation “waypoints”, four on the outgoing leg and three on the incoming leg, was used to perform the turns at the end of each flight line (Fig. I). The turns were executed by the autopilot, with little or no manual intervention. Two waypoints (1 and 7) were placed at the edges of the main radar grid, 150 km apart, and two (2 and 6) at the edges of the extended grid, 180 km apart; the remaining three waypoints, the outermost of which was 30 km from the edge of the main grid, were used to perform the actual turn. This technique worked well, and all oscillations induced by the autopilot mechanism were usually damped out completely by the time the extended grid was entered. The aircraft was always flown entirely by the autopilot for at least the full 180 km distance across the extended grid (between waypoints 2 and 6), to ensure that each flight line was as straight as possible. Autopilot-induced oscillations did not occur during the resets over the Summit/OSU camp since the aircraft was always flown manually during this time.
Each flight line was done at a constant-pressure altitude, rather than attempt to follow the curvature of the ice surface. The height of the aircraft above the surface was typically about 300–400 m, with the lowest value being about 250 m. Variations in this height had no significant effect on the signal strength of the bottom return, At no time was the aircraft too close to the surface to prevent a surface return from being recorded on the ice-radar records; the line picked always had the expected dome-shaped curvature and was separated from the end-of-supression line by the expected amount.
Data Analysis
Speed of propagation
A constant propagation speed, in ice, of 168 m/¼s, corresponding to a dielectric constant of 3.19, was assumed. This value has been used many times on glaciers and ice sheets by other workers (Reference Robin. deRobin, 1975; Reference PatersonPaterson, 1981; Rose, unpublished). To account for a faster propagation speed in the firn layer at the top of the ice sheet, a constant correction of +10 m was added to all ice thicknesses. This value was computed using the density-velocity relationship of Reference Robin, Evans and BaileyRobin and others (1969) and a measured depth-density curve from “Site A”, about 170 km to the south and almost the same altitude (Reference Alley and KociAlley and Koci, 1988). This value of IOm is identical to that used by Rose (unpublished) in West Antarctica.
Correction for INS drift
After conversion to Lambert coordinates, (X,Y), the horizontal position of the aircraft was corrected for INS drift by assuming a constant drift rate within each flight loop (all flying within two successive reset passes over the Summit/OSU camp). The closure errors, δX and δY, are listed in Table I. These closures were linearly distributed with time to all intermediate navigation data points. The total closure error, δR = (δX 2 + δY 2)½, ranged from 0.5 to 6.4 km, with an average of 3.3 ± 1.6 km. The average drift rate was 1.8 km/h, a typical value for the INS that was used. These average values exclude the data for the first day because they are considerably less reliable than the data for the remaining days.Footnote *
TABLE I. Closure errors. δX and δY are the apparent change in horizontal position at the end of each loop, δZ s is the apparent change in altitude of the ice surface, δR is the total closure error δX 2 + δY 2)½, and 67 is the elapsed time to fly the loop. the asterisk flags an assumed value and the question marks flag data which are considerably less reliable than the rest; all of these values are excluded from the averages (see footnote on p. 20)
The actual flight lines, corrected for FNS drift, are shown on all maps as dashed lines. The dashes are formed by connecting adjacent ice-thickness data points to provide a visualization of the data-sampling density. The flight lines are not drawn where there are gaps in the applicable data, and all flight lines to and from the Summit/OSU camp are omitted.
Correction for atmospheric pressure changes
In a similar manner, the altitude of the aircraft above sea-level, Za, was corrected for changes in atmospheric pressure during the time the loop was flown, again by assuming any such changes were linear with time. The measured pressure altitude and the measured height of the aircraft above the surface, H a, at the start and end of each loop, combined with the known altitude of the Summit/OSU camp, were used to correct all aircraft altitudes to values referenced to the WGS-72 ellipsoid. The altitude closures, δZS, are also listed in Table I. Values range from −5.1 to + 14.6 km, with an average of ±4.2 ± 2.3 m (again excluding the data for day 1).
The second loop on day 3 (not included in Table I) was not closed due to running low on fuel, and so for this loop it was necessary to assume no INS drift and no atmospheric pressure changes with time. † Fortunately, this affected only one flight line, W000SaFootnote †, and this was completed within 50 min of the initial reset pass over the Summit/OSU camp, so any increase in uncertainties for this line should be minimal. On the first loop of day 2, the radar altimeter was not working when the initial reset pass was made, and so for the affected three lines (S125Eb, S250Wb, and S375Ea) it was also necessary to assume no atmospheric pressure changes with time.
The analysis assumed there were no horizontal spatial gradients in atmospheric pressure. This is a potential source of error, but, as discussed in the error analysis, the assumption probably has negligible effect on the final results.
Smoothing the navigation data
The digitized output of the INS, recorded every second, contained a random uncertainty of plus or minus one least significant digit. Since the resolution of the INS was 0.1 min of arc (about 190 m in latitude and 60 m in longitude at the center of the grid), and the aircraft speed was about 120 m/s, this gave rise to an unrealistic aircraft trajectory, with apparent erratic jumping in position, in random directions, of the order of ±200 m between adjacent data points. This caused problems with subsequent calculations, particularly crossing-point adjustments (see below), and so the raw data were smoothed with cubic splines. The amount of smoothing was kept to the minimum necessary to achieve a realistically smooth aircraft trajectory.
The aircraft altitude and height-above-terrain data also contained a digitization uncertainty of the order of ±0.3 m (±1 ft), but they were also subject to a larger oscillation, of the Order of a few meters peak-to-peak, with a period of about 20 s (2–3 km wavelength along the flight line). These were probably induced by the autopilot mechanism since they occurred simultaneously in data from two independent aircraft altimeters as well as the radar altimeter. When combined with the uncertainty in the horizontal position, these oscillations also lead to convergence problems with the crossing-point adjustments, and so they were also smoothed with cubic splines. The smoothing parameters were varied independently of those used for the INS data, so that only the minimum necessary to remove the oscillations was used.
Ice-thickness data
Precise time marks, synchronized with the navigation data clock, were manually placed in real time on the Visicorder ice-radar records at waypoints 6, 7, 1 and 2, as well as approximately every 5 min along the flight line between waypoints 7 and I. From waypoint 6 to waypoint 2, 181 evenly spaced values of ice thickness, H, were measured from these records. Since waypoints 6 and 2 are nominally 180 km apart, this implies the ice thickness was sampled at approximately 1 km intervals along each flight line.
The true location of each time mark was then extracted from the corrected and smoothed INS data set and used to map the sampled ice-thickness points to their true (X,Y) coordinates. Each section between successive time marks, typically 15–35 km long, was treated independently. This process thus makes adjustments for distance scale changes between sections but it assumes that both the recorder-chart speed and the aircraft speed were constant within any given section. However, along any flight line, distance scale changes between sections were always well within the digitizing uncertainty of the INS and so any errors introduced by this assumption are also assumed to be contained within this INS digitizing uncertainty.
A total of 181 ice-thickness values was obtained from each of the 34 flight lines, or 6154 values for the entire grid (less the two short gaps discussed earlier). Each thickness value was resolved to the nearest 0.5 mm (0.25 j�s) on the chart paper, a digitizing error of about ±10 m of ice thickness. This error band of 20 m is approximately the same as the radar resolution.
The ice-radar records were used solely to determine ice thickness, H, whereas the surface topography of the ice sheet, Zs, was determined solely from the aircraft navigation data by subtracting the aircraft height-above-terrain from the aircraft altitude above sea-level: Z s = Z a - H a. The bottom topography of the ice sheet, Z b, was calculated by subtracting the measured ice thickness from the surface altitude of the ice sheet: Z b = Z s - H. Figure 3 shows examples of these profile data, using the same flight lines depicted in Figure 2.
Migration
No migration corrections (Reference HarrisonHarrison, 1970; Reference Brown, Rasmussen and MeierBrown and others, 1986) have been applied to the ice-thickness data. Migration attempts to correct for the fact that the bottom reflection actually comes from the closest point to the aircraft, allowing for refraction effects, rather than from the nadir point directly below the aircraft. Nadir thicknesses are underestimated when this effect is ignored.
This correction is roughly proportional to the bed slope and so in areas of flat bottom topography the effect is small. Figure 4 shows a contour plot, produced with the same techniques described below, of the slopes of the bottom topography over the entire grid. A bed slope of 60 m per kilometer will cause the reflection point to shift about 350 m away from the nadir in 3000 m of ice, and an underestimation of the ice thickness by about 20 m. Thus, only slopes greater than this will always shift the measured ice thickness outside the digitizing error band of ±10 m, and therefore always produce a detectable effect. Areas of bed slope greater than this limit are shown stippled in the figure; they cover a small percentage of the total area (12%) and lie primarily along the mountainous eastern edge of the grid.
The stippled areas in this figure are referred to in the rest of this analysis as the areas of “rough” topography, where, by definition, the bed slope is greater than 60 m/km (3.4°) and migration effects would be detectable. Conversely, the unstippled areas are referred to as the areas of “flat” topography.
Averaging of repeated lines
Ideally, the two lines of a repeated pair should give exactly the same results. However, the uncertainty imposed by the measurement errors, especially the INS drift, produce slight differences. Because the lines are almost, but not exactly coincident, these differences produce many local topographic slopes which are unrealistically large. When the data are contoured, this gives rise to artificial relief along these lines which dominates the overall topography and which cannot be eliminated without also destroying most of the valid information in the rest of the data. Consequently, it was necessary to form a single composite line from the two lines of each repeated pair; this was done using a simple unweighted averaging of each variable (X, Y, H, Za, Ha), point-by-point along the lines. This technique also ensures that the data set is as homogeneous as possible, with the same data density along all flight lines.
Fig. 3. Sample surface- and bed-altitude profiles. The same flight lines used in Figure 2 are depicted here. The navigation data were recorded at approximately 120 m intervals and have been smoothed with cubic splines to remove instrumentai and discretization noise (see text). The thickness data were read from the radar records (Fig. 2) at 1 km intervals. Geoid corrections have not been done Oi litis point and so these data are relative to the WGS-72 ellipsoid.
Adjustment for crossing-point differences
The 13 north-south lines cross the 13 east-west lines at 169 “crossing points”. For the same reasons given for repeated lines, the surface and bottom altitudes are not always exactly the same. This has much the same effect on the final contoured results, with dominant artificial relief at many of the crossing points. To eliminate this, the data values were adjusted as follows.
If ƒ is a data value, for example, the ice thickness, then, at a crossing point of a north-south and east-west line, the crossing-point difference is defined as δ = ƒv - ƒx, where ƒx and ƒy are the data values on the east-west and north-south lines, respectively. For any given line segment between two adjacent crossing points there are thus two values of δ, one at each end of the segment: δ1and δ2, where subscript 1 refers to the crossing point with the smaller spatial coordinate (x for east-west lines and y for north-south lines) and subscript 2 to the other point. The data values on east-west lines were then adjusted with
where
and ƒ is the unadjusted value, and ƒ′ the adjusted value, at an intermediate point x on the segment, x 1 < x < x 2 Data values on north-south lines were adjusted with identical equations, except that the signs of ε1 and e 2 are reversed and x is replaced by y. The adjustment is thus a simple linear one, with the adjustment on any one segment independent of that on any other segment in the grid.
Crossing-point adjustments were applied independently to two variables, the ice thickness, H, and the altitude of the ice surface above sea-level, Z s = Z a - H a. Both of these quantities are required to have the same value at all crossing points. In the latter case, the adjustment was applied to Z a, and H a was left unchanged. Figure 5 shows histograms of the absolute value of the crossing-point differences, |δ|, for these two variables. The average difference in ice thickness is 77.4 ± 80.4 m and the average difference in surface altitude is 4.3 α 3.3 m. These results arc discussed further in the error-analysis section.
Missing data
Complete loss of Visicorder records occurred for 25 km on line N625W and for 35 km on line S750W. Ice-thickness values were interpolated in these gaps from data along the nearest parallel line(s), with a linear scaling so that there was no discontinuity with valid data at the end(s) of the gaps. This data loss affected H and Z b but not Z s.
Fig. 4. The smoothed bottom slopes, in meters of altitude change per kilometer horizontally along the flight lines. The contour interval is 10 m/km. The stippled areas (12% of the total grid) are where the bottom slope exceeds 60 m/km. above which migration corrections would have a detectable effect. The dashed lines on this figure, and on all other maps, indicate the actual flight lines, corrected for INS drift: the dashes are drawn by connecting adjacent ice-thickness data points.
No values of H a were obtained for 66 km at the start of line S125Eb. These were filled in by using the corresponding values from the other line of the repeated pair (S125Ea), again with a constant shift so that there was no discontinuity in Z s at the end of the gap. The only error which is introduced by this process is that caused by differences in the shape of the surface topography along the two lines. These two lines are only 0.8 km apart at the end of the gap; on this scale, such differences are completely insignificant (see Fig. 8). Thus, this adjustment is assumed to have correctly reconstructed the missing values and is not considered to be a true loss of data.
Finally, for 166 km of line N500E there was no surface return on the radar records, due to an incorrect trigger-level setting on the recorder. The offset which resulted when the problem was corrected was used to estimate an initial location for a (constant) surface altitude. Values of H were then determined in the usual manner, but the crossing-point adjustments were modified so that the crossing north-south lines were held fixed along this line and all crossing-point differences were applied to the data on line N500E. This effectively calibrated the values of H with the crossing lines, the only approximation being a linear surface topography between each of the 13 crossing lines.
Interpolation
The final data sets, after averaging of repeated lines and adjustment for crossing-point differences, consisting of 181 × 26 = 4706 values of position (X,Y), ice thickness H, surface altitude Z s, and bottom altitude Z b. To contour these data, they first had to be interpolated on to a square (or at least rectangular) grid. Unfortunately the data, which are very dense along a few discrete lines but completely absent elsewhere, are in a form which is difficult to interpolate realistically, since all available methods assumed a random distribution of points. To help compensate for this problem, both Laplacian and cubic spline interpolation, as well as various relative mixtures of these two methods,* were employed so that the overall effect of the interpolation process on the final results could be judged. The data were interpolated on to a square grid with a 2 km spacing, a total of 912 = 8281 values.
Figure 6 shows the effect of varying the interpolation method. The dimcnsionless parameter ø determines the relative contribution by Laplacian and cubic spline interpolation. An intermediate value, ø = 5, was used for the remaining analysis. This gives a roughly equal mix of the two methods and so their relative advantages and disadvantages should tend to compensate each other.
Smoothing
The final step before contouring was to smooth the interpolated data by passing them several times through a LapIacian smoothing operator.Footnote * If this is not done, the resulting contour maps exhibit a preponderance of topographic details along each flight line, relative to that in the data-free interior zones between the flight lines. Some of these details are caused by measurement noise, but most of them are real. In any case, this effect is a consequence of the interpolation being unable to propagate such features into the interior zones, and it dominates the final results when contoured and viewed at the scale of the full 180 km by 180 km grid.
Fig. 5. Freiiuency distribution of the crossing-point differences. |δ|, for the ice thickness. H, and the surface altitude. Zs. The average difference in ice thickness is 77.4 ± 80.4 m and the average difference in surface altitude is 4.3 ± 3.3 m.
Fig. 6. Effect of varying the interpolation method. The dimensionless parameter ø determines the relative contribution by Laplaciait and cubic spline interpolation, ø = 0 is totally Laplacian and ø = 1000 is totally cubic spline; ø = 5 is an approximately equal mixture of the two methods, and is the value used for all remaining figures. All plots here have been smoothed eight limes (see Fig. 7). The contour interval is 50 m.
Figure 7 shows the effect of varying the amount of smoothing, done by successively doubling the number of passes made with the smoothing function. A smaller contour interval than that justified by the overall error (see next section) was used to ensure that features would not be missed simply because of the arbitrary choice of contour levels. With no smoothing, this clustering of small-scale topography along the flight lines is readily evident. Eight smoothing passes are necessary to eliminate any obvious evidence of this clustering and so this value was used for all contour plots.
Geoid correction
All altitudes were referenced to the GEM-10B geoid model (Reference Lerch, Putney, Wagner and KloskoLerch and others, 1981). This is the same geoid used as the reference for maps of the surface topography of southern Greenland produced from satellite (Seasat) radar altimetry by Reference Bindschadler, Zwally, Major and BrennerBindschadler and others (in press. The correction varies smoothly from a low of +41.9 m in the north-west to a high of +47.4 m in the south-east. All altitudes on the contour maps in this paper are thus altitudes above mean sea-level, as defined by the GEM-10B geoid.
Contouring
Figures 8–10 are contour plots of the surface topography, ice thickness, and bottom topography produced with these techniques using 2 km gridded data sets. The contour intervals are 10, 100, and 100 m, respectively, as determined by the error analysis (next section). Figure 11 is a contour plot of the bottom topography with a 50 m contour interval, presented solely as an aid to visualizing the topographic trends.
Error Analysis
Navigation errors
Uncertainty in the horizontal location of the aircraft, referred to as the “navigation error”, is an indirect source of error in the altitude of both the surface and bed topography. The error in altitude is the product of the error in horizontal position and the slope of the surface or bed, respectively. Because this turns out to be the dominant ice-thickness error term, and thus determines the minimum mapping contour interval allowed, it is discussed in further detail.
The navigation error is caused primarily by the drift inherent in any INS. Errors in the ability of the pilot and navigator to reset precisely over the reference point, as well as errors in the geoceiver coordinates of this point, are at least an order of magnitude less than the drift error, and so are assumed to be negligible. INS drift consists of a linear component and an oscillating component, called the Schuler drift (Rose, unpublished). The period of the Schuler drift is 84.4 min, but its phase and amplitude are variable and difficult, if not impossible, to predict. For the INS used in this study, the peak-to-peak amplitude of the Schuler component is, however, probably of the order of a few kilometers (Rose, unpublished).
Fig. 7. Effect of varying the amount of smoothing, done by making repeated passes with a Laplacian smoothing operator. S is the number of smoothing passes used. The value of S = 8 was used for all other figures. All plots here use an interpolation parameter of ø = 5 (see Fig. 6). The contour interval is 50 m.
Fig. 8. The surface topography, in meters above mean sea-level (GEM-10B geoid). The contour interval is 10 m, and the contours are accurate to about ±6 m. The triangle indicates the location of the n ue summit of the Greenland ice sheet. This has an altitude of 3233 m a.s.l.; its location (lat. 72°34′ N., long. 37°38′ W.) is similarly marked on the remaining maps.
In order to estimate the magnitude of the INS drift error, we assume the linear drift component has been removed successfully with the linear loop-closure corrections described earlier, at least to an accuracy comparable to the INS resolution (approximately ±200 m). We then assume that the standard deviation in the total closure errors (Table I) is, at least in some crude sense, a result of the randomness imposed by performing the INS resets at arbitrary, and unpredictable, points within the 84.4 min Schuler cycle. Combining this value, ±1.6 km, with the resolution uncertainty, ±0.2 km, gives the navigation error, ±1.61 km, used in the remainder of the error analysis.
Surface-altitude errors
The average surface slope over the entire grid, calculated over 1 km intervals from either the raw profile data or the gridded data (Fig. 8), is 0.001 ± 0.0005. The navigation error thus gives rise to an average surface-altitude error of about ±1.6 m. The pressure altimeter is accurate to ±2.1 m (personal communication from H. Terry), and the radar altimeter to ±5 m (personal communication from W. Carver). The latter includes an estimated allowance for local variations in the depth of penetration of the radar-altimeter signal into the surface snow. The geoceiver altitude of the Summit/OSU camp is accurate to about ±2 m (personal communication from J. Bolzan), relative to the GEM-10B geoid. The relative geoid correction applied to the rest of the grid has no error since we define mean sea-level to be this geoid.
The interpolation error is estimated at ±0.4 m, determined by computing the mean difference in the gridded data between the ø = 5 “average” interpolation and the two “extreme” cases of ø = 0 (all Laplacian) and ø = 1000 (all cubic spline). The smoothing error is estimated at ±0.2 m, determined in a similar way by computing the mean difference between the smoothed and unsmoothed gridded data.
The combined effect of all these errors is ±6.0 m. In addition, the pressure-altimeter measurements are subject to four sources of error which are difficult or impossible to estimate: (1) departures of the real atmosphere from a standard atmosphere, (2) spatial variations in pressure, (3) non-linear temporal variations in pressure, and (4) departures of the altimeter’s geopotential reference surface from the WGS-72 ellipsoid.
A “temperature correction” is often done to attempt to correct, at least partially, for non-standard atmospheres (Rose, unpublished). However, since all altitude measurements were made relative to a single known point and the total altitude band covered by the aircraft was only a few hundred meters, this correction was not done, and the resulting errors are assumed to be relatively small. Similarly, a “cross-wind correction” (Rose, unpublished) is often done to allow for spatial variations in pressure. This correction assumes the aircraft was high enough above the terrain where the winds are purely geostrophic. This height, 500–1000 m, is at least twice as high as we typically flew and so the assumptions on which this correction are based probably do not apply well in our situation. Thus, this correction was not done either, and the corresponding errors are also assumed to be relatively small.
Fig. 9. The ice thickness, in meters. The contour interval is 100 m. The contours are accurate to about ±50 m in the flat areas and to about ±125 m in the rough areas (stippled zones in Figure 4). On Figures 9–11 closed contours marked with an “H” indicate local highs and unmarked closed contours indicate local lows
Atmospheric pressure variations at the Summit/OSU camp were less than 3 mbar over the entire 6 d period and less than 1 mbar (11 m at 3000 m altitude) during any 24 h period (personal communication from J. Bolzan). Thus, during a 2–3 h period the errors introduced by non-linear variations in pressure are likely to be very small (the data analysis already allows for linear changes with time).
Even though the errors introduced by (1)-(3) are probably small, they are still largely unknown, and the error resulting from (4) is completely unknown. Fortunately, however, an independent means of estimating the overall error in surface altitude is available from 20 direct measurements of the surface altitude obtained with geoceivers by the OSU group. These geoceiver measurements, accurate to about ±2 m vertically, are distributed reasonably uniformly over the entire grid (written communication from J. Bolzan). The maximum departure of the radar altitudes from the geoceiver altitudes is +14.1 m and the average departure is +4,8 ± 5.4 m. Since both the average and standard deviation are less than the ±6.0 m error already obtained, we consider that any errors introduced by (1)-(4) are contained within this figure and that surface altitudes are indeed accurate to ±6.0 m.
Another check is to note that the crossing-point differences in surface altitude average 4.3 m. Since this is the difference in two values which are effectively averaged in the data analysis, this implies an equivalent altitude error of ±2.15 m, well within the ±6.0 m estimate. A final check is provided by the altitude-closure errors, which averaged +4.2 ± 2.3 m (Table I); since this in effect represents a random sampling of the pressure variations with time, it indicates that any non-linear temporal variations are probably well within the ±6.0 m.
Ice-thickness errors
A similar analysis can be done for the ice thickness. The average bottom slope over the entire grid, calculated over 1 km intervals from either the raw profile data or the gridded data (Fig. 4), is 0.038 ± 0.024. The navigation error thus gives rise to an average ice-thickness error of about ±61.2 m. However, since potential core-hole sites are likely to be chosen from areas of flat bed topography, the error analysis was done separately for both flat and rough areas, as defined earlier. For flat areas, the average bed slope is 0.031 ± 0.014, and the average navigation-induced ice-thickness error is ±49.9 m; for rough areas, the corresponding values are 0.082 ± 0.026 and ±132.0 m.
Since the average bed slupe in the flat areas is half of the 0.060 limit at which migration effects would exceed the radar digitizing error of ±10 m, the average error introduced by not doing a migration correction is assumed to be ±5 m in the flat areas. In the rough areas, the corresponding figure is ±13 m (10 × 0.082/0.060).
Additional errors come from the actual radar measurement. The assumed propagation speed of 168m/μs is probably accurate to about ±0.5 m/μs (Rose, unpublished); in 3000 m of ice this translates to about ±9 m of ice thickness. The firn correction is estimated to be accurate to about ±20%, or about ±2 m. As already noted, digitization errors (identifying the precise time of the start of the surface and bed returns on the radar records, as well as measuring their separation), arc estimated to be ±10 m. Finally, based on crystal-controlled calibration time marks placed by the TUD radar on the records, the assumed recorder sweep rate (0.5 μs/mm of chart paper) is estimated to be accurate to at least ±0.1%, or about ±3 m of ice thickness.
Fig. 10. The bottom topography, in meters above mean sea-level. The contour interval is 100 m. The contours have the same uncertainly as the ice-thickness contours (Fig. 9).
Interpolation and smoothing errors are estimated in the same way done for the surface-altitude error estimate. Interpolation contributes an average error of ±8 m, and smoothing with eight passes, ±16 m.
Combining all these error components, the average error in ice thickness is ±55.0 m for flat areas and ±134.6 m for rough areas. The crossing-point differences in ice thickness average 77.4 m, which implies an equivalent ice-thickness error of ±38.7 m. Like the surface altitudes, this is well within either the ±55.0 or ±134.6 m error estimates.
Bottom-altitude errors
These are simply a combination of the surface-altitude and ice-thickness errors. Flat areas are thus accurate to ±55.3 m, and rough areas to ±134.7 m.
Effect of navigation errors
Over 90% of the error in ice thickness and bottom altitude comes from the navigation error. Even though it has been assumed in this analysis that all error components are statistically random on a scale corresponding to the spatial sample density (about 1 km), this particular component, in fact, is largely non-random on this scale. The reason for this is that the navigation error is dominated by the Schuler cycle. Each 180 km long flight line took about 25 min to fly, or about 30% of the Schuler period. Thus, the navigation error actually imposes a large-scale distortion on each flight line, rather than a small-scale random fluctuation. Only about one-third of a full cycle of distortion is applied to any given flight line.
The overall effect of this is that the ice thicknesses and bottom altitudes are probably significantly more accurate than the above values on a small scale, of the order of a few kilometers. The ice-thickness and bottom-altitude errors are therefore conservatively rounded down to ±50 m for flat areas and ±125 m for rough areas.
Contour intervals
Standard mapping practice requires that 90% of all points tested for altitude are correct to within half the contour interval (Reference WolfWolf, 1983). Since this is a 2σ criterion, and the above errors are Io values, this implies the contour intervals should be 24 m for surface altitude and 200 m for ice thickness and bottom altitude.
In order for the resulting maps to contain some detail, however, we have relaxed this convention to a 1σ value, resulting in a contour interval of 12 m for surface altitude and 100 m for ice thickness and bottom altitude. This is at least partially justified by the fact that not only the navigation error but also many of the other error components are not totally random. Instead, they contain a significant non-random part which does not affect the shape of the topography on scales comparable to the topographic wavelengths (tens of kilometers).
In addition, the contour interval for surface altitude has been decreased slightly to 10 m so that a convenient round number is used. The crossing-point differences and comparison with geoceiver data both support this change.
Summary
The error in surface altitude is ±6 m. This error is controlled largely by the random error in the radar altimeter, and navigation errors are a relatively minor contribution. The error in bottom altitude and ice thickness is ±50 m for flat areas and ±125 m for rough areas (slopes greater than 0.060). This error is dominated by the navigation error, which is mostly non-random on scales less than the topographic wavelengths, and errors from the ice radar itself are a relatively minor contribution. Contour intervals of 10 m for surface altitude and 100 m for bottom altitude and ice thickness are used; the (1σ) error in these maps is approximately one-half the contour interval.
Results
Figure 8 shows a contour plot of the surface topography. The true summit of the Greenland ice sheet is located at about lat. 72°34′N., long. 37°38′W., (X ≃ 9.4, Y ≃ 31.2), and has an altitude of 3233 m a.s.l. This position and altitude were found to be consistent to within about 1 km horizontally and 1 m vertically, regardless of the interpolation method, the amount of smoothing applied, or even whether or not any INS drift corrections were done. Thus we consider these to be very reliable coordinates.
The ice divide south of the true summit is a well-defined ridge which runs almost exactly due south, with a very low slope of about 0.0005, dropping only 55 m in a distance of about 120 km. North of the true summit, however, the ice surface fans out into a broad, north-facing slope. This slope is bounded by a relatively sharp ridge running north-east away from the summit and a much rounder and less well-defined ridge in the north-west direction.
The surface topography surrounding the summit can thus be divided into three sectors, each of which has a reasonably uniform characteristic slope. The north-west-north-east sector has the lowest slope, about 0.0012, the north-east-south sector the steepest, about 0.0020, and the north-west-south sector an intermediate slope of about 0.0015. The west slope dips in a direction about 20° south of west and the east slope about 30° south of east. The lowest surface altitude in the grid is about 3042 m, giving a maximum surface-altitude range of only 191 m over the entire 32 400 km2 of the grid.
The ice thickness (Fig. 9) rangesFootnote * from 1790 to 3375 m, with an average value of 2975 ± 235 m. At the true summit the ice is about 3025 m thick. The thickest ice (3375 m), on the other hand, is near lat. 71°45′N., long. 38°36′W. (X ≃ −24, Y ≃ −60), about 97 km south-southwest of the true summit. Most of the south-west quadrant of the grid, in fact, contains thick ice, typically over 3200 m thick. The ice over most of the rest of the grid is significantly thinner, generally less than 3000 m.
If the original profile data are examined, however, the very thickest, and thinnest, ice occurs in small, localized spots in the north-east corner of the grid, which contains very mountainous bottom relief. These local minimum and maximum values are, approximately: ice thickness, 1550 and 3470 m; bottom altitude, −370 and 1590 m. These values differ considerably from the gridded and smoothed values because, being the extremes, they are subjected to the maximum amount of adjustment. In addition, they are both within a few kilometers of the edge of the grid and so are also subject to large edge effects in the gridding and smoothing.
Because the ice-sheet surface is very flat, the bottom topography (Figs 10 and 11) is essentially a mirror image of the ice thickness. The thick ice in the south-west quadrant is due to a large, flat basin which is almost entirely below sea-level. The lowest point is at the same location as the thickest ice and has an altitude of about −215 m. The relative relief over most of this basin is very flat and smooth, with peak-to-trough amplitudes less than 100 m and wavelengths of tens of kilometers.
The bottom topography over most of the rest of the grid is generally only a few hundred meters above sea-level. The average bed altitude over the entire grid is +180 ± 235 m. The terrain over much of the north-west quadrant is also relatively smooth and drops just below sea-level in a few locations. There is no predominant direction to the basal topography. It appears to be undulating, rolling terrain with no obvious ridge/valley structure.
The true summit of the ice sheet is above the eastern end of a comparatively large plateau. This bench is approximately 10–15 km wide and extends about 50 km to the west. It is about 200–300 m a.s.l. and is reasonably smooth and flat, with relief comparable to that of the south-west basin. The ice thickness above this plateau is 2900–3000 m.
If smoothness of the bottom topography were the sole criterion, the most suitable areas for a core hole would be either this plateau or the south-west basin. If, in addition, the thickest ice was also desired, then a site somewhere in the south-west basin, with X <−15 km and Y <−15km, would be the best choice.
Acknowledgements
A large number of people contributed to the success of this project. Dr P. Gudmandsen participated in the planning and contributed a great amount of invaluable advice and support. Dr D. Drewry gave helpful advice regarding flight planning and aircraft navigation.
We wish to thank Mr H. Terry for the loan of the navigation recording system and for providing us with the navigation data on nine-track tape. We are also grateful to Dr C. Bentley for the loan of the Visicorder used to record the radar data, and to the students at St. Olaf College who measured those records.
The entire crew of the U.S. Navy VXE-6 squadron deserve special mention. All of them went well out of their way to ensure that we obtained the best navigation and greatest amount of data possible in the limited time-available. Lt Cmdr J. Dedon was extremely helpful and in his dual role as lead navigator and project officer-in-charge was instrumental in the overall success. We also wish to express our appreciation to the personnel at Thule Air Base, who provided us with excellent support during our stay in Greenland, and to the Polar Ice Coring Office for their usual exemplary assistance with Greenland logistics.
Finally, we wish to express special thanks to Dr J. Dionne of the U.S. National Science Foundation, without whose tireless and aggressive support the work would probably never have taken place, and to Dr J. Bolzan and the other members of the Summit/OSU camp who struggled out of warm sleeping bags several times during the middle of the night to provide the absolutely essential mark-on-top for our navigation.
This work was supported by the U.S. National Science Foundation agreement number DPP-8520948.
