Mike Seymour and Ian Cumming
University of British Columbia
Dept. of Electrical and Computer Engineering
2356 Main Mall, Vancouver, B.C.,Canada V6T 1Z4
tel: +1-604-822-4988     fax: +1-604-822-5949
http://www.ee.ubc.ca/sar/
Presented at the International Geoscience and Remote Sensing Symposium, IGARSS'97, Singapore, August 3-8, 1997.

1. Introduction
2. Overview of Algorithm
3. Algorithm Performance
4. Simulations
5. ERS Example
6. Conclusions

We have developed a method for employing DEMs to estimate the geometry of SAR images used in interferometric SAR processing. The geometry of the interferometer is important because it controls the mapping of topography to interferometric phase. We give an overview of the method and then describe the performance of the algorithm. Simulations are presented as examples of the algorithm performance. Finally, comparisons of geometry derived using medium to high quality DEMs for an ERS FRINGE test site with precision orbit data are also made. We found the prime limitation on the performance of the algorithm is the influence of DEM errors on estimating the error due to improper baseline estimates.

Synthetic aperture radar interferometry (InSAR) is a technique for estimating topography using the phase difference between coincident SAR images. One important part of the InSAR algorithm is calibrating the interferometer by estimating the geometry (e.g. the baseline magnitude and orientation) of the coincident SAR images. In addition, InSAR's emphasis on coherently processing overlapping data sets has spawned a number of other applications which usually require relative geometry of SAR images.

There are a number of different methods for estimating the baseline magnitude and orientation at various accuracies [1,2,3]. To obtain accurate baseline estimates, existing height control points [1] must be used along with the unwrapped phase data. For ERS data, precision orbit data provide baselines which are sufficiently accurate to perform differential interferometry to detect centimeter magnitude topographic changes [4] but not accurate enough for topographic estimation.

There is actually a great range of quality and resolution in existing digital elevation model (DEM) data See Bruce Gitting's DEM catalog at http://www.geo.ed.ac.uk/home/ded.html . At one end of the spectrum there is the GTOPO30 (See http://edcwww.cr.usgs.gov/landdaac/gtopo30/gtopo30.html ) data set which has coverage over most of the Earth's land surface and is currently available publically from the United States Geological Survey (USGS). GTOPO30 has low resolution with nominal posting of 1 km and nominal 90 significance level of 160 meters. This DEM is ideally suited to large scale modelling of the Earth's topography. At the other end of the spectrum are DEMs of high resolution and very high accuracy (RMS error <1 meter) which generally cover small areas and are not publically available. In between these two extremes there is a range of DEM data usually with regional coverage and sampling spacing of 30 to 200 meters and accuracy from 5 to 30 meters. For example, in the province of British Columbia in Canada, TRIM trim digitial maps provide data sampled at nominally 100 m posting with accuracy of 5m. Note that if the DEM quality or resolution is poorer than the nominal 30 m posting, 10 m RMS error, performance of SAR interferometry, the InSAR DEM can be used to enhance the resolution and quality of the existing DEM data.

Existing DEM data can be applied in both the flattening and calibration stages of processing [6,7]. In particular, the calibration operation with an existing DEM can proceed without phase unwrapping, thereby avoiding a potentially difficult processing step. In the following, we investigate the performance of the baseline estimation algorithm. We start by giving a brief overview of the algorithm in Section 2 followed by an error analysis in Section 3. Simulations results are described in Section 4. An example of processing with an ERS data set is given in Section 5. Finally, some conclusions are given in Section 6.

An interferogram is formed by complex conjugate multiplication of two registered SAR images. If the images are sufficiently correlated, the interferogram phase is a function of the slant range difference between the positions of the SAR sensors used to form the interferogram and the common ground patch:

  …   (1)

Assuming a constant travel time for the microwave energy, the reference slant range (r) can be calculated from the radar timing information while the second slant range (r2) can be accurately calculated via phase unwrapping. However, can also be estimated as a function of the interferometer geometry if the topography is known. We exploit the geometry and the existence of a DEM by solving the minimization problem:

  …   (2)

as a function of the baseline length (B) and orientation (Q). Estimates of r2 are made from the interferogram and processed data as the algorithm progresses.

The algorithm proceeds in two stages. In the first stage, a coarse estimate of r2 is derived from the registration relation between the two images used to form the interferogram. Using this value in (2) provides reasonable initial estimates of B and Q . In the second stage, a residual interferogram is analyzed to provide a more accurate estimate of r2 to increase the accuracy of the B and Q estimates.

After removing the preliminary estimate of the interferogram phase from the interferogram a residual interferogram having phase Y _res is formed. The residual interferogram consists of three phase terms: a term due to the estimate of the geometry error (Y _geom), a term due to uncorrelated noise in the interferogram (Y _noise), and a term due to the errors in the DEM (Y _dem):

  …   (3)

The geometry error phase is approximately:

  …   (4)
where:

The elements of Y_geom can be estimated using standard frequency domain techniques. Then, the estimate of r2 can be corrected to provide a high accuracy second slant range term to fit against in (2). The final estimate of the geometry is made by solving the minimization problem with the estimated slant range from the coarse estimate plus the correction from the residual interferogram. In practice this procedure must be iterated, but the number of iterations required for convergence is usually low (between 2 and 7).

A critical step in the algorithm is estimating the parameters of (4) to update the r2 estimate used for the fitting process. Both the phase noise (Y _noise) and DEM errors (Y _dem) can possibly interfere with the estimate of Y _res. Usually, the effect of the noise is minimized by the large number of samples in the residual interferogram used to estimate w_r,a , w_r , w_a.

However, DEM errors can cause errors in two ways: correlation between DEM errors and the geometry error phase terms of (4); and smearing of the spectrum estimate of the residual interferogram. The DEM errors and the residual phase due to geometry errors must be uncorrelated for the residual interferogram phase to reflect geometry errors alone. In practice, this means that the DEM errors cannot have linear or bi-linear trends. This condition should be satisified by processing a sufficiently large-sized interferogram and having DEMs of good quality to process.

Because the DEM error term and the geometry error term are multiplied, any spectral method actually processes the convolution of the spectrum of and .Y _dem is a function of the interferometer geometry and the DEM errors (Dh):

  …   (5)

A necessary condition for the residual interferogram peak to correspond to geometric phase error is that the frequency spectrum of have a maximum at DC. Statistiscal characterization of the spectrum of is difficult because of the influence of the unknown topography and the impact of interpolation from noisy DEM data to the radar data grid. In general, as the accuracy of the DEM decreases, it becomes less and less likely that the spectrum of will have an absolute maximum at 0 frequency. Note that the spectrum of the residual interferogram can serve as a diagnostic about the performance of the algorithm. If there are a large number of discrete peaks near 0 frequency, then it is likely that an error will have ocurred in the geometrical phase estimate.

There are other sources of error in the geometry estimate. The two most serious are the altitude of the satellite and the estimate of the mean second slant range from the interferogram fit procedure. For ERS data, the precision orbit data product [4] provides more than accurate enough altitude values. Adequate estimation of the mean r2 value for the preliminary estimates of B and Q hinges on using enough samples during the initial registration of the SAR images used to calculate the interferogram.

To demonstrate the algorithm, a synthetic interferogram was generated using terrain heights from TRIM [5] data maps covering an approximately 40 by 10 kilometer continuous patch of mapsheet 92-O, Taseko Lakes. For the purpose of testing, the TRIM data maps were also degraded by deleting data points from the DEM data. The resulting DEMs have significantly less data points than the orginal DEM with a manufactured error [8]. Data sets were generated from which one could reconstruct ground range topography with nominal maximum errors of 10, 25, 50, and 100 meter errors respectively. In addition, overlapping coverage of 30 Arc-second Digital Chart of the World (DCW) data [9] was also available for the scene. The statistics of the errors for the various data sets are summarized in Table 1. Note that the accuracy requirement is defined assuming that the mean height error of the DEM is 0 and that the net residual slope due to baseline error is less than 1 in 10000.

Since the error performance of the algorithm depends both on the baseline magnitude and orientation, simulations were performed with baseline magnitudes of 70, 170, 270 and 370 meters respectively. The orientation of the baseline was chosen to be Q » 1.17 radians as a representative case. To investigate the dependence of the algorithm on the baseline orientation, the algorithm performance at other values of was also investigated and similar trends in the results to those reported here were found. Results are summarized in Table 2.

From Table 2, it is clear that the algorithm works well with the higher quality input DEMs: the 10m and 25m data sets. With these DEMS, the algorithm has no problem in estimating the baseline magnitude and orientation within the requirements of having less than 10 meters slope in the InSAR DEM error from near to far range. The 50m data set nearly always meets the criterion. For the 100 m data set, the algorithm begins to fail when the baseline magnitude changes from 70 meter to 170 meters. The far coarser and higher error DCW data set did not ever provide the correct geometry. In further experiments we found that the normal baseline had to be less than 20 m for the algorithm to work with the DCW data. Upon examining the spectra of the residual interferograms at higher error, it was clear that the DEM errors were causing significant distortion in the residual interferogram phase estimates and hence the errors in the baseline estimates.

The processed data consist of an ERS-1 SAR image (orbit 21730) and a co-located ERS-2 SAR image (orbit 2057) collected on September 10^th, 1995 and September 11^th, 1995 respectively at a local time of approximately 10 AM PST. The SAR data were collected on track 199 at frame 2565 for both satellites. The center of the scene is approximately 35' 20" N, 59' 20" W. The SAR images are well correlated with a normal baseline of approximately ~42 meters [10]. The processed interferogram was of dimension 900 range samples by 2000 azimuth samples after filtering, or about 1/20 of a full SAR scene. See Figure 1.

The baseline was estimated using a linear model for both the baseline orientation and baseline magnitude for overlapping DCW data, a DTED-1 data set, and a TRIM data set. In addition, ERS precision orbit data also gives another basis of comparison. The estimated baselines are shown in Figure 2. The DTED-1 and TRIM based baseline estimates agree to approximately 1 meter with the precision orbit baseline while the DCW baseline is between 1 and 2 meters different from the precision baseline. The disagreement between the DTED and TRIM data was puzzling considering that they are both fairly good quality DEMs (DTED 90% significance level is 30 meters [11]). However, the difference of the two DEMs had linear and bi-linear trends which account for the disparity and an unexpectedly large standard deviation of approximately 44 meters.

We have reviewed an algorithm for estimating baseline geometry without phase unwrapping using existing DEM data. The performance limit of the algorithm comes from the influence of DEM errors on estimation of the error due to incorrect geometry (i.e. incorrect B and Q ). Noise is not a significant problem if the data is suitable for interferometry because of the large number of samples for estimating each term in the residual phase (see (4)). Processing larger data sets should partially alleviate the problems with trends found in the DTED DEM data used in the experiment. Note that even for large error, coarse DEMs, the baseline estimate provides a flattened interferogram suitable for further processing such as filtering, phase unwrapping and topographic estimation.