ERS tandem study of glacier dynamics in NE-Greenland

Abstract

This paper describes a technique for decomposition of the interferometric phase from multiple ERS Tandem acquisitions into terrain height and surface displacement. Descending and ascending data are combined, and the 2-D surface parallel flow of the glacier extracted. Due to the utilization of precision orbit data only four tie-points are needed, and those need not be identified in the image. Mean SAR derived velocities, extracted in small patches on the surrounding bed rocks, are generally less than 5 m/y, with a standard deviation on the order of 1.4-2.0 m/y in the N-S direction and 0.7-1.0 m/y in the E-W direction. The spatially dense sampled height and velocity data provide an outstanding opportunity to model the dynamics of the glacier and thus enhance our understanding of buildup and surge.
Keywords: SAR Interferometry, Glacier, ERS

Introduction

The presented work is performed within the frame of the announcement of opportunity study "An Investigation of the Utility of ERS SAR Data for Studies of Glacier dynamics" conducted by the Danish Center for Remote Sensing (DCRS). The study concerns the use of ERS-1/2 Tandem SAR data to provide spatial measurements of instantaneous ice velocity and also microtopography using a multibaseline technique.

The primary test field is Storstrømmen, a major outlet glacier in North East Greenland. In the period 1978-1984 the glacier went through a large surge, with a frontal advance of approximately 12 km, (Reeh, Bøggild, Oerter, 1994). During six field seasons in the period 1989-1995 field measurements of ice velocities and meteorological conditions have been acquired by the Alfred Wegener Institute and the Danish Polar Center. Polarimetric and interferometric SAR imagery have also been acquired over an area near the terminus with the airborne EMISAR operated by DCRS.

A unified approach to processing the airborne EMISAR data and the satellite ERS-1/2 has been developed. For airborne repeat track interferometry the important issue is motion compensation - an integrated part of the focusing. Our ERS study uses raw unfocussed SAR data, even though ESA provides processed images (SLC data) as a standard product, because the optimal choice of processing parameters (e.g. Doppler centroid and Doppler bandwidth) require knowledge of which data sets is to be combined in the later processing steps. Full control of the entire processing chain and access to raw unfocussed SAR data, is also important for calibration and for assessment of system artifacts such as missing lines, phase drift etc.

Repeat track interferometry

The phase of a pixel in an interferogram formed from two images denoted and are (to first order)

(1)

where denotes wavelength, parallel baseline, perpendicular baseline, is the angular deviation from the reference surface (the topography), and is the unknown displacement in the line of sight direction (the target motion). Random noise as well as systematic errors e.g. phase unwrapping errors or path length differences induced by atmospheric inhomogenities is included in .

In order to extract the displacement term from (1) the topography phase must either be supplied by an external digital elevation model (DEM), or multiple interferograms can be used to decompose into elevation and displacement, provided that the baselines (spatial and temporal) for two or more interferograms give a linearly independent equation system, (Massonnet, et.al., 1993), (Gabriel, Goldstein, Zebker, 1989), (Zebker, et.al., 1994), (Kwok, Fahnestock, 1996).

The interferometric technique only provides the displacement in the line of sight direction, but the flow vector consist of three components. Thus a direct measurement of the flow require measurements from three different directions. This can, in principle, be accomplished by an airborne SAR system, but not with the ERS-1/2 satellites only capable of viewing from two directions (ascending and descending orbits.) However, the full three dimensional flow can be derived, if the flow can be confined. Assuming that the flow is parallel to the surface the three equations for the determination of the flow vector (of each pixel) become

(2)

where and are the displacements measured from the ascending and descending orbits, and are the corresponding line of sight vectors and and are the observation intervals. The surface normal can, if multiple interferograms are utilized, be derived from the radar derived DEM.

The application of (2) is straightforward, provided that the measurements from each orbit are geocoded to a common reference system, such as geographical latitude and longitude or UTM. With the advent of precision orbit (PRC) data, this geocoding can be performed without identifying tie-points in the images, provided that the PRC data can be synchronized with the radar data, (Massmann, et.al., 1997). The raw SAR data from ESA is timetagged using a 3.9 ms counter, corresponding to an along track distance of approximately 30 m. Identification of tie-points in imagery focused to zero Doppler compared to a map indicates that the time tagging used for the PRCs and the raw data is consistent to within 50 m. However, this is an upper bound (i.e. synchronization might be better) as identification of tie-points in regions without man-made targets is extremely difficult. Similarly the range-offset provided in the CEOS leader files has been assessed. Again, the analysis can only provide a crude upper bound on the absolute accuracy, in this case on the order of 75 m. Thus, presently the range-offset provided by ESA is used.

It is worth noting that the viewing geometry of the ERS-1/2 system is not optimal for extracting the essentially horizontal glacier flow vectors. The ratio of the vertical errors compared to the horizontal is approximately 0.4 (tangent to the incidence angle), see fig. 1. Similarly the E-W errors vs. the N-S errors have a ratio of , where is the angular separation of the ascending and descending orbit tracks, see fig. 2. For the Storstrømmen glacier, at a latitude of 76° N, the angular separation 56° is corresponding to a ratio of 0.5. For areas more towards the equator the geometry deteriorates, for areas more towards the poles the geometry become more favorable.

Fig. 1. Vertical vs. horizontal errors. The incidence angle for ERS-1/2 is 23° at mid-swath.

Fig. 2. E-W vs. N-S errors. The angular separation between descending and ascending orbits is 56° at 76° N.

Dependent on the approach, and the temporal decorrelation, a different number of images is required. It is our experience that winter data usually exhibit good correlation during the 1-day repeat period of the ERS-1/2 tandem data. Interferograms formed from 1-day summer data are generally difficult to unwrap due to temporal decorrelation (melting is one important decorrelation mechanism) and 35-day data typically show no correlation at all over glaciers.

Decomposing multibaseline data into the topography and displacement terms, for descending as well as ascending orbit data requires in total 8 images of each frame. This approach is presently used. The advantage is that no external DEM is required, the disadvantage that 8 high quality images must be available. The other extreme, requiring only 4 good images, is to supply the topography phase from an external DEM. Although more robust, this method has the disadvantage that a DEM must be available. For Greenland this is an important constraint, as the horizontal resolution of DEMs in many areas is several kilometers, and in the poorest mapped areas no map data exist at all, even at a basic 1:250.000 scale, (Ekholm, 1996). The ideal solution is obviously to implement both options. In that case a the DEM extracted from one look direction can be used to remove the topographic phase from the other viewing direction. This only requires 6 good images and does not require an external DEM.

Note that the DEMs is necessary not only for removing the topography phase, but also for the geocoding. This is critical as an vertical height error, translates into a horizontal positioning error a factor 1/0.4=2.5 greater.

Processing Techniques

The processing chain has been designed for airborne as well as satellite data, see fig. 3. For the ERS-1/2 data the difference from the conventional processing of satellite data is that the utilization of precision orbit data is integrated in the basic SAR focusing via a motion compensation scheme. This procedure automatically aligns the images before interferogram formation, and additionally remove the flat earth phase. The other major difference from conventional interferometric processing schemes is that the output from the topography/motion decomposition for each pixel is i) the 3-D position, ii) the radar line of sight vector, iii) the displacement in the line of sight direction and iv) radar brightness, coherence etc. This is important when measurements from different viewing geometry’s is combined (ascending descending orbits).

Fig. 3. Block diagram of processing chain. Each box corresponds to one program unit. Two arrows indicates that multiple (two) images is combined. Each program writes a log-file with parameters needed for the subsequent processing steps, but basic setup (of for example generic file names or multi-looking factors) is done separately for each unit.

Initially the raw SAR data is cleaned to add dummy data for missing lines, and to align data in the case of sampling window start time (SWST) changes. Multiple frames can be merged, easing the processing of strips.

Precision orbit data is interpolated and a local track is chosen for each image pair. An additional global reference line is chosen to facilitate a combination of interferograms from different 35-day cycles.

In the SAR focusing each image is resampled to the local reference line, similar to motion compensation in the airborne case. A spherical earth reference surface is used. However, this is not critical as the motion compensation phase (the flat earth phase) is removed prior to the inversion of the interferometric phase to elevation/displacement. A common Doppler is used for processing and the image is band-pass filtered (and decimated) in azimuth prior to range compression. Processing of several ascending and descending scenes have revealed a consistent internal 5 m slant range delay difference between ERS-1 and ERS-2. Due to the coarse timetagging of the raw data, a small image patch is initially processed and the differential azimuth timing error estimated (typical less than 1 ms).

During the interferogram formation images are resampled to the global reference line (no phase correction) in effect aligning all interferograms. Examples on an ascending orbit (track 244) and a descending orbit (track 382) interferogram and correlation image covering a 100 km by 200 km area around Storstrømmen is attached, see plates. In particular note the phase artifacts and the corresponding correlation streaks near the terminus of the glacier in the descending orbit interferogram.

Fig. 4. Example: Ascending orbit arc raw interferogram and correlation image. ERS-1: April 10th, 1996. ERS-2: April 11th, 1997. Track=244.

Fig. 5. Example: Descending orbit arc raw interferogram and correlation image. ERS-1: October 28th, 1995. ERS-2: October 29th, 1995. Track=382.

Phase unwrapping is performed using a Goldstein cutline line approach, with optional use of neutrons (guidance to formation of cutlines).

The decomposition of multiple interferograms is performed on a pixel by pixel basis, as the baseline is varying in the azimuth direction.

The 3-D positions, displacements etc. are now (optionally) resampled to a ground coordinate system and (optionally) converted to geographical coordinates. The resampling to a regular grid is done by a local mean, plane or parabolic fit using the weights

(3)

where d is the distance from the input pixel to the output position, and (k,s,m) are constants; for ERS-1/2 data averaged to 100 m by 100 m typical values used for (k,s,m) are k=2, s=100 and m=3. The final derivation of flow vectors is done by a direct application of (2).

Results

A initial assessment of the accuracy can be performed by analyzing the measured velocities in small 5 by 5 pixel patches on the bed-rock. The mean value of the velocities is typically less than 5 m/y with a standard deviation on the order of 1.4-2.0 m/y in the N-S direction and 0.7-1.0 m/y in the E-W direction. A more detailed assessment and interpretation of the results is given in (Mohr, Reeh, Madsen, 1997).

Conclusion

A processing system for extraction of 3-D glacier velocities and surface elevations from ascending and descending orbit ERS-1/2 data has been developed. Due to the utilization of precision orbit data from the D-PAF, images can automatically be aligned prior to interferogram formation. Also the geocoding benefits from the use of PRCs as only 4 tie-points, which need not be identified in the images is needed. This is important for mapping in arctic regions, where identification is extremely difficult. It is anticipated, that the direct measurements of surface velocities and elevations, will provide important basis for our understanding of glacier and ice sheet dynamics.

References

Ekholm, 1996:
A Full Coverage, High-Resolution, Topographic Model of Greenland Computed from a Variety of Digital Elevation Data. J. Geophysical Research, 101, no. B10, pp. 21961-21972.

Gabriel, Goldstein, Zebker, 1989:
Mapping Small Elevation Changes over Large Areas. J. Geophysical Research, 94, pp. 9183-9191.

Kwok, Fahnestock, 1996:
Ice Sheet Motion and Topography from Radar Interferometry. IEEE trans. Geosci. Remote Sensing, 34, pp. 189-200.

Massmann, et.al., 1997:
Quality of the D-PAF ERS Orbits Before and After the Inclusion of Prare Data. These proceedings.

Massonnet, et.al., 1993:
The Displacement Field of the Landers Earthquake Mapped by Radar Interferometry. Nature, 364, pp. 138-142.

Mohr, Reeh, Madsen, 1997:
Three Dimensional Glacial Flow and Surface Elevation Measured with Radar Interferometry. Submitted to Nature.

Reeh, Bøggild, Oerter, 1994:
Surge of Storstrømmen, a Large Outlet Glacier from the Inland Ice of North-East Greenland. Rapp. Grønlands Geol. Unders., 162, pp. 201-209.

Zebker, et.al., 1994:
On the Derivation of Coseismic Displacement Fields using Differential Radar Interferometry; the Landers Earthquake. J. Geophysical Research, 99, pp. 19617-19634.

Keywords: ESA European Space Agency - Agence spatiale europeenne, observation de la terre, earth observation, satellite remote sensing, teledetection, geophysique, altimetrie, radar, chimique atmospherique, geophysics, altimetry, radar, atmospheric chemistry