Introduction

The ice-sheet topography plays an important role in ice-sheet-flow studies. On the one hand, large-scale topography controls flow directions and its accurate mapping may allow us to derive balance velocity of ice sheets. On the other hand, the knowledge of the topographic slope and of the ice thickness allows an estimate of the basal shear stress which is related to the deformational velocity of the ice.

The best way to obtain a global and precise topography of ice sheets is by satellite altimetry. Its coverage and precision allow us to improve our knowledge of ice-sheet flow; for example, a few attempts have been made to estimate rheological parameters, the Glen flow-law exponent n and the activation energy Q (Young and others, 1989; Rémy and Minster, 1993). The latter demonstrates the role of the topography in such studies. Using a precise altimetric topography deduced from Seasat, they found n = 1 and Q = 70 kJ mol⁻¹ for temperatures lower than −10°C, but using earlier maps (Drewry, 1983) they found n = 1.8 and Q = 30 kj mol⁻¹ for the same conditions. For a typical value of stress of 0.5 bar, this leads to a difference of 40% between estimates of deformation velocity. Moreover, satellite altimetry provides time series that could yield estimates for long-term or seasonal variations in ice-sheet elevation. Zwally and others (1989) measured a mean increase in the southern Greenland ice sheet of 23 cm year⁻¹, using Geosat and Seasat altimeter data, while Andersen (1994) extracted an annual signal of mean amplitude of 40 cm over Greenland from the analysis of ERS-1 preliminary data.

Satellite altimetry also allows us to map short-wavelength features on ice sheets such as undulations (λ = 20 km). The presence of these short-wavelength features has been attributed to the process of ice flowing over an irregular bedrock (Budd, 1970; Hutter and others, 1981) and it appears that in places they are a direct reflection of the basal topography. Their occurrence is also linked with basal conditions; whereas these features are enhanced in areas of high friction, the presence of subglacial lakes suppresses them (Ridley and others, 1993). Amplitude and wavelength of such features can be directly estimated by using satellite altimetry.

Nevertheless, altimetric mapping for ice-sheet dynamics presents some difficulties. First, it is important to correct the data for instrumental, orbit and surface slope-induced errors. Among other things, we will see in the next section that, with insufficient care in the altimetric treatment, errors can induce a global underestimation of the surface slope. Secondly, the distribution of altimetric data is not homogeneous (350 m separation along track, several kilometres between tracks), leading to distortions in the map of undulations, To prevent this, we must extrapolate carefully height measurements in the across-track directions: the inverse method, that will be described later, partially corrects this effect. Finally, the height measurement is derived from a beam with a footprint of several kilometres and is thus distorted by the same scale topographic features (Wingbam, 1995). A direct analysis of the radar-altimeter wave-form echoes will provide an independent estimation.

Altimetric Processing

We used the ERS-1 35 d repeat cycle waveform altimeter product (WAP), delivered by the UK-Processing and Archiving Facility (UK-PAF).

Fig. 1. Global coverage of the ERS-1 ice-mode altimeter data (cycle 88). Each dot represents 100 measurements of height. The average density varies from 0.2 data km⁻² at 68° S to 6 data km⁻² at 81.5° S. The worst-covered zones are located near the coast, where there are rock outcrops.

Among other instruments, the ERS-1 payload comprises an altimeter working with two operational modes, designed to operate well over both the ocean and the continental ice. It covers all latitudes up to ± 81.4°, whereas Geosat and Seasat, its predecessors, were limited to ± 72.1°. ERS-I therefore allows a global and synoptic control of most of Greenland and 80% of Antarctica.

In this paper, we have processed all data collected during the 88th 35 d repeat cycle, an ice-mode cycle in October 1992, and the whole 89th repeat cycle which was an ocean-mode cycle in November 1992. The ice-mode tracking is devoted to the tracking of steep topography and allows complete coverage of Antarctica; while ocean-mode tracking, having a 4 times higher band width, is more accurate but loses track close to the margins. The mixing of both data sets provides a good compromise: ice mode for large-scale and global analysis and ocean mode for precise and short-scale analysis. In order to ensure maximum coverage, by minimizing seasonal effects, some missing tracks of the 89th cycle were replaced by the same tracks from neighbouring ocean cycles (83rd, 85th and 87th). The instrument delivers 20 wave forms per second, or every 350 m along track, each wave form being the time series of samples of integrated echo power received by the altimeter. This provides about 3 × 10⁶ height measurements on each complete cycle (Fig. 1). These raw data must be processed and corrected for many instrumental and physical effects. Only geophysical effects will be discussed here.

Dry tropospheric error

The interaction of the travelling wave with the dry troposphere creates a propagation delay. The resulting error is corrected using a theoretical formula derived from Saastromoinen (1972):

where h is the local ice-sheet altitude. This correction, cor _tropo, is added to the altitude estimate.

Slope error

This error is due to the shift of the reflecting point of the radar wave on the ice surface, in the upslope direction. Before proposing a correction for this error, let us briefly recall its characteristics, since it is the most critical for our study. First, it is dominant, except near the domes where the slopes are less than 1 m km⁻¹ . Secondly, it is difficult to correct with a single track, because we have no information on the slope in the across-track direction. Thirdly, this error not only depends on the surface slope of the reflecting surface but also on the curvature (Rémy and others, 1989): if f(x) is the true topography along the upslope direction and g(x) is the apparent one. “seen” by the satellite, then we can write, to second order,

(1)

where H is the satellite height, almost 800 km. fʹ(x) and fʺ(x), respectively, what we call surface slope and surface curvature, represent the first and second derivatives of the topography with distance in the upslope direction. fʹʺ(x) is supposed to be negligible.

The residual error is always positive and is about 40 cm for a slope of 1 m km⁻¹ but it reaches 1.6 m for a slope of 2 in km⁻¹, as found typically within the continent and 10 m for a slope of 5 m km⁻¹.

Using Equation (m1), we can see th a t the apparent slope g'(x) IS, to second order,

(2)

If one does not lake the curvature fʺ(x) into account, the residual error on the slope gʹ(x) is about O.1 fʹ(x) for a typical curvature value of 10⁻⁷ m⁻¹. Moreover, because of the near-parabolic shape of the Antarctic ice sheet, slope and curvature are geographically linked, a poor slope correction will cause a general underestimation of the derived slope, and thus an underestimation of the derived shear stress.

Our method consists of computing the exact surface by fining the two-dimensional apparent topography with a biquadratic form. Slope and curvature parameters are then determined with nine neighbouring grid points. The derivation of Equation (2) leads to

(3)

so that

(4)

Equations (1), (2) and (4) then lead to

(5)

This technique has the advantage of being easy to implement. One can note that when fʺ(x) = 1/H. Equation (1) is unbounded. Indeed, in such a case, the curvature is too important and the first return of the radar wave has not been reflected by the bottom of the undulation. The determination of topography is actually limited to this curvature value. Fortunately, in Antarctica, the curvature of the topography does not reach this critical value, except near the coast.

Mapping the Topography and Surface Slope Using a Total Inverse Method

At this stage of processing, the altimetric data have been interpolated with a precise orbit and corrected for instrumental, atmospheric and tracking errors. Any residual errors still remaining are within the 80 cm r.m.s. cross-over error. In order to reduce this remainder error, we compute the topography with a total inverse technique, developed for ocean mapping with altimetric data (Mazzega and Houry, 1989) and adapted for ice-sheet mapping (Rémy and others, 1989). Each term of the error budget is taken into account (except slope error which is corrected at last), by describing its statistics in terms of the covariance function and the decorrelation length. The height signal we are seeking is modelled with a Gaussian function, defined by its maximum value that we call variance, and by its decorrelation length. This technique has several advantages. First, the knowledge of the statistics of the various components of the error budget help us to separate the error from the topographic signal. Secondly, the a posteriori error can be simultaneously computed. Finally, the full benefit of the inversion technique lies in the possibility of separating large-scale features from short-scale undulations. The spectrum of the huge-scale signal is very different from that of the remainder errors and topographic undulations and we recover the topography in two steps: a large-scale inversion of the ice-mode cycle provides a mean surface covering most of Antarctica and a second iteration allows us to determine a very precise topography with the ocean-mode cycle.

Fig. 2. Region of 500 km × 800 km, indicated in Figure 1, showing the large-scale topography, contoured at 10 m intervals. The 100 km wavelength signal is removed.

Fig. 3. Residuals obtained by the second-step inversion, superimposed on satellite tracks, for three small regions defined in Figure 2. We can see at 128° E (b) a strong kilometre-wavelength signal which is also detected on the back-scattered energy variations map (Fig. 8). Contour interval is 2 m.

In the first step, variance of the signal has been set at 700 m and the decorrelation length is 150 km. This distance is larger than the distance between tracks, so the inversion algorithm is not sensitive to the poor distribution of data. The resulting smooth topography is restored over a regular grid with steps of 0.5° in longitude and 0.1° in latitude (about 15 × 10 km²). This map (Fig. 2) is used as a reference surface for the second step.

After subtracting the surface derived in the first step from the measured altitude of the ocean-cycle data, the resulting residuals were processed through the inversion algorithm. This residual signal, mostly a small-scale signal, shows the undulations typically of 10 m amplitude and 20 km wavelength. In the second-step inversion, the variance of the signal is set to 7 m, when the variance of the orbit error is set to 1 m and the noise is set to 1 m. The decorrelation length is very different for each component (20 km for signal, 1 cycle/revolution for orbit error, 350 m for noise) and both signal and errors are well separated. In fact, the space time structure of the various data errors is difficult to characterize, so that we prefer slightly to overestimate the a priori error variances. This choice will tend to slightly oversmooth the topographic signal. Figure 3 shows three small maps of residual topography from the area indicated on Figure 2. One can see on figure 3a that coverage is dense, except at 140° E, where strong undulations are visible. In such areas, the lack of data is due to failure of the onboard tracker over rough surfaces, and the same parts of the track are missing in other ocean cycles. At high latitudes, where the inter-track distance is less than the decorrelation length, the inversion algorithm is able to estimate the geolocation and value of a local extreme between two tracks (Fig. 3a and b). In that case, the a posteriori error, which is the difference between the a priori error and the information gained through the inversion (Fig. 4a and b), is homogeneous and less than 1 m except where there are no data. In figure 3c, undulations are found under the tracks; this is a consequence of the particular geographic sampling. Indeed, for latitudes lower than 75° S, the distance between two tracks is too large to restore the small-scale features and so the algorithm fits topography on the reference only. This effect also explains the shape of the formal error shown in figure 4c, where the error is maximum between two tracks. In fact, over such zones, undulations are not well recovered, as we do not have sufficient information on the across-track direction slope.

Fig. 4. The a posteriori error maps calculated within the inversion, for the areas defined in Figure 2. The satellite tracks are also drawn. Breaks in these tracks indicate that no data were collected. In (a) and (b), this error is almost homogeneous and sub-metric, whereas in (c) it is minimum only under tracks.

We will see in the next section that the estimation of the undulation characteristics from the along-track profiles induces a global underestimation of these features.

Finally, the second-step residuals are added to the large-scale topography in order to estimate the local surface slope and curvature and to compute the corresponding correction to be applied to the data (Equation (5)). Figure 5 shows this topography for our three regions. The same process yielded the final topography of the whole of Antarctica, which is presented on a grid of about 5 km (Fig. 6). We believe that the accuracy of this map is better than previous maps. There are certainly significant changes, for instance, Dome C, a precise description of which was needed for the Italian–French ice-core project Concordia, appears at 75.1° S, 123.3° W, 60 km from its previous position. This has been confirmed by ground-based measurements. Other topographic features such as ice divides, domes and drainage basins are clearly identifiable.

Surface slopes derived from Figure 6 are shown in Figure 7. Slopes are less than 0.5 m km⁻¹ nearer the dome or divide and increase in the downslope direction. Some anomalies are also apparent: for instance, flat regions such as Vostok lake (76–78° S, 105° E), already mapped by Ridley and others (1993), Astrolabe Lake (70° S, 135° E), or the great slope region along longitude 135° E which corresponds to an elongated bedrock feature (Drewry, 1983).

Mapping Small-Scale Features

As explained in the introduction, small-scale features that correspond to kilometre-scale roughness directly reflect bedrock irregularities and basal conditions. Moreover, they are related to the small-scale shear-stress variations, that affect ice-sheet-flow modelling (Rémy and others, 1996b).

The mapping of these structures with satellite-altimeter data is critically dependent on the inter-track distance. As shown in figure 3c when the inter-track distance is large, the maximum values of the short-scale features are found under satellite tracks. This effect will be Stronger at lower latitudes and as such undulations are strong. Moreover, the across-track widths of these structures depends strongly on the covariance length used in the mapping.

One can show that, if undulations are isotrapically dome-shaped, as suggested by McIntyre and Drewry (1984) (i.e. if h = Aexp(−2((x ² + y ²)/λ²)), where λ is the wavelength, A is the amplitude and y is the across-track distance), then the maximum amplitude as seen by the satellite will be Aexp(−2(d ²/λ²)), where d is the nearest distance of the satellite tract to the dome centre. If we do not consider the shift of the reflecting point due to the regional slope, d is statistically uniformly distributed between 0 to λ, yielding to an average amplitude of Aexp(−1/2) or about 60% of the real amplitude.

Fig. 5. The total topography in the three regions defined in Figure 2 by adding the second-step residuals to the large-scale topography. Contour interval is 10 m.

Fig. 6. Global topography of the Antarctic ice sheet as determined over significant regions of the ERS-1 cycle 89 of the altimetric data set. Thin lines correspond to 100 m interval contour lines, medium-width lines to 1000 m interval contour lines and the thick one is 3000 m.

This fact suggests we should search for another way to estimate the amplitude of kilometre-scale undulations. Currently, we tend to ignore the back-scattered energy but the time series recorded by the altimeter can be thought of as a histogram of the received energy vs the arrival time and corresponds to a footprint of a few kilometres radius, depending on the local surface. This energy results from reflection both by surface micro-roughness, largely created by katabatic wind (Rémy and others, 1990) and by volume-scattering inside the snow-pack (Ridley and Partington, 1988), while the short-scale spatial variations of this parameter are mostly due to the variations of the illuminated surface and to the kilometre-scale curvature of the topography. Note that other mechanisms such as surface temperature and drifting snow also play a role but these meteorological phenomena only play a role on the 100 km spatial scale. Legresy and Rémy (paper in preparation) show that the amplitude of the variations of the energy δσ, calculated over a 50 km scale, are related to undulations by

Fig. 7. Global surface slope derived from the topographic map. This map shines the pattern of the principal surface features. Fig. 7. Global surface slope derived from the topographic map. This map shows the pattern of the principal surface features.

Fig. 8. Amplitude of the variation of the returned energy (in dB), averaged over 50 × 50 km² squares. This parameter is directly related to the undulation curvatures.

(6)

where c is the undulation curvature, H is the satellite altitude and δσ is expressed in dB. δσ yields the true characteristics of the undulations at the footprint scale. It is totally independent of the altimetric height, which is distorted near the coast or in the presence of strong undulations. This kilometre-scale curvature is directly linked to the small-scale variations of the slope and therefore is related to the accuracy of the restored surface slope.

The large-scale geographical distribution of this parameter, averaged over areas of 50 × 50 km² (Fig. 8), is strongly correlated with the surface slope (Fig. 7) (correlation of 0.7 for the set of 12 000 grid points), as theoretically suggested by Hutter and others (1981). The region indicated in Figure 8, with a strong gradient, is the region mapped in Figure 2. This gradient is also visible in figure 3b, where the lefthand side has far more undulations than the righthand part, and in figure 3c, where a flat region is seen around 133° E, corresponding to the minimum of δσ.

West Antarctica shows large values of δσ (Fig. 8), as does East Antarctica near the coast, on the great convergence area, and where the bedrock is important (see 135° E). A few low values can be observed, which correspond to the known subglacial lakes, as for example the Vostok lake (76–78° S. 105° E; Ridley and others, 1993), the Dome C region or the Astrolabe Basin (Oswald and Robin, 1973).

Conclusions

Satellite altimetric observations above the ice sheet provide surface topographic and surface-slope maps, and contribute to glaciological understanding, by producing a small-scale curvature map which may be related to basal conditions. However, to be useful, these parameters must be carefully corrected for the principal altimetric errors. Error due to surface slope is critical, because it contains a systematic component not only related to surface slope but also to surface curvature, which generally correlate over the ice sheet. The use of the two-stage formulation presented here reduces the residual error and prevents systematic bias. The inverse technique used allows the restitution of the a posteriori error, which is around 1 m at high latitude but reaches several metres near 70° S. Because the linear sampling of altimeter data does not allow the correct retrieval of the small-scale topography, we prefer to use the variations in back-scattered energy which are directly linked with the topographic curvature at kilometre scale. We show that this parameter is strongly related to surface slope and to the presence of several subglacial lakes, and so constitutes a valuable complement to height measurements.

Finally, the most important problem preventing sub-metre precision altimetry and high-precision slope mapping of the whole continent is the inadequate coverage. The map is, however, better than the earlier ones. The 168 d repeat cycle will give better resolution than the 35 d repeat cycle and it promises to provide a coverage five times as dense at the 35 d repeat one, and so it will allow good estimation of kilometre-scale roughness and its spatial distribution at low latitudes.