PRECISE 2-D AND 3-D GROUND TARGET LOCALIZATION WITH TERRASAR-X

Previous studies have demonstrated the unprecedented absolute pixel localization accuracy of the German SAR (Synthetic Aperture Radar) satellites TerraSAR-X and TanDEM-X. Now, using global navigation satellite system data to locally correct for atmospheric delay, range accuracies of 1 cm are demonstrated to be attainable. Using globally available tropospheric and ionospheric delay information such accuracies can even be approached world-wide. When such accurate estimates are available for two or more SAR images from different orbits, techniques such as stereo SAR can be used to generate DEMs (Digital Elevation Models) using only the amplitude signal, completely bypassing interferometry. The first experiments of this type described here show that an absolute 3-D localization accuracy better than 10 centimeters is achievable. * Corresponding author.  INTRODUCTION Spaceborne SAR (Synthetic Aperture Radar) is known for its ability to provide weather and time of day independent observation of the earth’s surface and measurement of relative shifts based on the carrier phase (SAR interferometry). In contrast, our objective is the absolute pixel localization of a SAR image. Previous studies of the German SAR satellites TerraSAR-X (TSX-1) and TanDEMX (TDX-1) revealed an unprecedented localization accuracy at the centimeter level (Balss et al., 2011, Eineder et al., 2011, Cong et al., 2012, Schubert et al., 2011, 2012). However, we expect that by a thorough correction of SAR range and azimuth observations for all path delays and geodynamic effects this can be improved to the subcentimeter level (Balss et al., 2012). As these effects are likewise relevant for the global navigation satellite system (GNSS), GNSS is used as the primary data source for the respective corrections.  MEASUREMENT PRINCIPLE OF SAR Radar systems indirectly measure geometric distances by means of the two-way travel time of radar pulses from the radar transmitter to the ground and back to the radar receiver. In a focused SAR image, the instant of closest approach of the sensor and target as well as the signal travel time at this instant define the two radar time coordinates of azimuth and range. Usually, the conversion from range time to geometric distance refers to the vacuum velocity of light. However, electrons in the ionosphere, dry air and water vapor in the troposphere introduce additional signal delays which have to be taken into account. In addition, geodynamic effects including tides, loading and continental drift shift the true position of a ground target.  CORRECTION OF SIGNAL PATH DELAYS 3.1 Tropospheric Delay The amount of the tropospheric delay in SAR range measurements lies between 2.5 and 4 meters. In order to achieve centimeter or even millimeter level accuracy, precise tropospheric delay estimation from GNSS is required. In (Eineder et al., 2011), a height-dependent exponential model with hourly GNSS zenith tropospheric measurements was applied. In its most recent implementation, our GNSS based method uses IGS (International GNSS Service) zenith path delays (Byram et al., 2011) and the Vienna Mapping Function I (Kouba, 2007). Unfortunately, the permanent GNSS stations required for such accuracies are sparsely located on a global scale. Therefore, a general approach using the ECMWF (European Centre for Medium-Range Weather Forecasts) ERA Interim 3-D numerical weather model data (Dee et al., 2011) to estimate the tropospheric delay along the SAR signal path up to the mesopause (approx. 85 kilometers) was introduced in (Cong et al., 2012). 3.2 Ionospheric Delay The ionosphere ranges approximately from a height of 50 kilometers up to 1200 kilometers and it arises from the high energy fraction of solar insolation which releases ions and electrons (Kelly, 2009). The ionospheric effect on the signal International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 – 24 May 2013, Hannover, Germany


 INTRODUCTION
Spaceborne SAR (Synthetic Aperture Radar) is known for its ability to provide weather and time of day independent observation of the earth's surface and measurement of relative shifts based on the carrier phase (SAR interferometry).In contrast, our objective is the absolute pixel localization of a SAR image.Previous studies of the German SAR satellites TerraSAR-X (TSX-1) and TanDEM-X (TDX-1) revealed an unprecedented localization accuracy at the centimeter level (Balss et al., 2011, Eineder et al., 2011, Cong et al., 2012, Schubert et al., 2011, 2012).However, we expect that by a thorough correction of SAR range and azimuth observations for all path delays and geodynamic effects this can be improved to the subcentimeter level (Balss et al., 2012).As these effects are likewise relevant for the global navigation satellite system (GNSS), GNSS is used as the primary data source for the respective corrections.

 MEASUREMENT PRINCIPLE OF SAR
Radar systems indirectly measure geometric distances by means of the two-way travel time of radar pulses from the radar transmitter to the ground and back to the radar receiver.In a focused SAR image, the instant of closest approach of the sensor and target as well as the signal travel time at this instant define the two radar time coordinates of azimuth and range.Usually, the conversion from range time to geometric distance refers to the vacuum velocity of light.However, electrons in the ionosphere, dry air and water vapor in the troposphere introduce additional signal delays which have to be taken into account.In addition, geodynamic effects including tides, loading and continental drift shift the true position of a ground target.

Tropospheric Delay
The amount of the tropospheric delay in SAR range measurements lies between 2.5 and 4 meters.In order to achieve centimeter or even millimeter level accuracy, precise tropospheric delay estimation from GNSS is required.In (Eineder et al., 2011), a height-dependent exponential model with hourly GNSS zenith tropospheric measurements was applied.In its most recent implementation, our GNSS based method uses IGS (International GNSS Service) zenith path delays (Byram et al., 2011) and the Vienna Mapping Function I (Kouba, 2007).
Unfortunately, the permanent GNSS stations required for such accuracies are sparsely located on a global scale.Therefore, a general approach using the ECMWF (European Centre for Medium-Range Weather Forecasts) ERA Interim 3-D numerical weather model data (Dee et al., 2011) to estimate the tropospheric delay along the SAR signal path up to the mesopause (approx.85 kilometers) was introduced in (Cong et al., 2012).

Ionospheric Delay
The ionosphere ranges approximately from a height of 50 kilometers up to 1200 kilometers and it arises from the high energy fraction of solar insolation which releases ions and electrons (Kelly, 2009).The ionospheric effect on the signal International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany path delay depends on the radio wavelength.In order to compensate for the ionospheric delay, we first determine the vTEC (vertical total electron content) based on GNSS (Hofmann-Wellenhof et al., 2008).The vTEC is then converted to the slant range signal path delay based on the radar frequency and the measurement geometry of the individual scene (Gisinger, 2012).
Our method makes use of the well-established SLM (single layer model) that assumes the ionosphere is a vanishingly thin spherical layer at an altitude of 450 kilometers (Schaer, 1999).This allows a consistent remapping of the slant ionospheric delays that are obtained from the geometry-free linear combination of dual-frequency GNSS observations for every GNSS satellite within an epoch.Aside from the vTEC, the geometry-free equation contains the DCB (differential code bias) of the observed satellite as well as the DCB of the GNSS receiver.The satellite and receiver biases can be decorrelated from the vTEC signal and separated from each other by introducing them as daily constants and imposing an additional zero-mean constraint on the satellite biases.For our particular experiment, we rely on the local Wettzell GNSS receivers that are part of the global IGS tracking station network.Their DCBs are known as they are estimated along with the daily computation of the CODE (Center for Orbit Determination Europe) GIMs (Global Ionosphere Maps) (Schaer, 1999).The daily GIMs include the satellite biases as well as those of the stations.With both sets of biases known, we can perform a straight forward computation of the vTEC observed by any GNSS satellite every 30 seconds.
In order to remove outliers from the dual-frequency measurements and obtain a much more homogenous set of code observations by combining phase and code, the code smoothing procedure of the Bernese GNSS software is applied (Dach et al., 2007).Depending on the observed GNSS satellite constellation, 8 to 10 vTEC values are obtained per epoch that correspond to the SLM line of sight pierce point at an altitude of 450 kilometers.After interpolating every GNSS vTEC series at the actual time at which the corner reflector (the pixel of interest) was imaged in the SAR acquisition, a 2-D plane modeling the vTEC gradient during the SAR acquisition is estimated by least squares.The vTEC evaluated at the SAR satellite pierce point combined with the corresponding SLM mapping function yields the ionospheric signal delay in the radar range measurement.The above procedure is carried out individually for all three IGS GNSS receivers at Wettzell (WTZA, WTZS, WTZZ) and the final ionospheric delay for the SAR signal is taken as the average.The standard deviation of the 3 receiver values can be interpreted as an overall consistency indicator for the final value.For the 30 datatakes acquired to date, almost all standard deviations of the 3 receiver values lie around 2 millimeters with a maximum of 5 millimeters.Hence, we conclude that the method delivers stable results.
As the orbits of TSX-1 and TDX-1 are still within the upper ionosphere, the upper portion of the ionosphere contributes to the GNSS based measurement value but not to the path delay of the radar signal.For this reason, we apply a 75% weighting factor which was derived from an analysis of the IRI-2007 (International Reference Ionosphere) model (Balss et al., 2012, Bilitza, 2012).
As for the troposphere, we investigated two different methods for determining the electron content.The first was based on local measurements provided by a GNSS reference station as outlined above.The second was based on CODE GIMs (Schaer, 1999) which enable a straight forward computation of the vTEC.CODE GIMs are provided for the entire globe every 2 hours on a spatial grid of 5.0º longitude x 2.5º latitude for a resolution on the equator of 560 km East-West x 280 km North-South which is large compared to the dimensions of a SAR image.In the vicinity of the reference station the former local method yields better results when one considers the final localization accuracy although the corrections themselves are very similarthe average (standard deviation) of the GNSS and CODE values to date are 6.3 (2.2) cm and 6.5 (2.1) cm respectively.The latter global method has of course the advantage of world-wide applicability.

 CORRECTION OF GEODYNAMIC EFFECTS
Solid earth tides and continental drift are the most prominent geodynamic effects.While the former causes a shift of up to a few decimeters over the course of a day, the latter slowly increases over years.However, in order to achieve millimeter localization accuracy, smaller geophysical effects must be compensated for.Atmospheric pressure loading and ocean tidal loading weigh on the tectonic plate and their variations cause a shift of up to several millimeters each.An even weaker effect is caused by pole tides.
The correction of these effects follows the conventions issued by IERS (International Earth Rotation Service) 2010 edition (Petit and Luzum, 2010).An exception is atmospheric pressure loading that is not yet included in the IERS conventions since the accuracy of modeling such non-tidal displacements for correcting reference points is still in debate (Petit and Luzum, 2010).However, in the field of VLBI (Very Long Baseline Interferometry), accounting for atmospheric pressure loading gave significant improvements and thus we decided to apply the same data from the NASA atmospheric pressure loading service (Petrov and Boy, 2004).Aside from this, the IERS conventions include all effects outlined in the previous paragraph as well as the submillimeter signals induced by atmospheric tidal loading and ocean pole tide loading.While these two signals may not be significant for the Wettzell experiment, they could have millimeter impact on other observations sites.In particular, ocean pole tide loading will become a more important effect with the planned expansion of this study to O'Higgins station in Antarctica as described on section 7.
All of the effects described above were computed for the Wettzell location in local North, East, Up coordinates at the time of the datatakes and were then transformed into radar range and azimuth time coordinates (Eineder et al., 2011).

 MEASUREMENT OF THE 2-D LOCALIZATION ACCURACY
In order to verify the pixel localization accuracy of a SAR system, the range and azimuth times of the corner reflector in the focused SAR images corrected for the estimated propagation delay were compared with the values obtained from precise on ground measurements of their positions.The conversion of the spatial geodetic coordinates into radar time International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XL-1/W1, ISPRS Hannover Workshop 2013, 21 -24 May 2013, Hannover, Germany coordinates is based on the zero-Doppler equation (Cumming and Wong, 2005) and interpolation of the satellite's position.
Our recent measurement series is based on a trihedral corner reflector with edge length 1.5 meters which was installed in the summer of 2011 at the Geodetic Observatory Wettzell, Germany.In this way, we benefit from the very close distance (approx.240 meters) to the local GNSS reference station.
The ground position of the corner reflector which is defined by its phase center (point number 401 in Figure 1) is known to within 1 centimeter relative to the reference station as obtained from a terrestrial geodetic survey.A strong temporal correlation between the final range location errors was observed as can be seen in Figure 4.A visual inspection of the temporal progression reveals that a good portion of the standard deviation results from a slow variation, perhaps due to seasonal variation of the measured range offset.In contrast, almost all immediately neighboring measurement values typically differ by millimeters.Thus, on a short-term scale, the measurements suggest that the localization accuracy for TSX-1 could be further improved.This effect is the subject of an ongoing investigation.

 MEASUREMENT OF THE 3-D LOCALIZATION ACCURACY
One particular application that can benefit from an improved TSX-1 geolocation capability is stereo SAR.In general this method can be used for constructing a DEM (Digital Elevation Model) from just two acquisitions taken from different orbits (Leberl, 1990) but due to the higher sensitivity of the phase, the interferometric approach is favored for DEM generation as in the TanDEM mission (Krieger et al., 2007).
However, dual orbit stereo simulations carried out for the Wettzell corner reflector using the achieved pixel localization accuracy suggest that precise, absolute 3-D scatterer reconstruction is possible.These simulations are based on a rigorous solution of the range and azimuth equation system at the zero-Doppler location and results in the most general case of parameter estimation, i.e. an adjustment of conditions with parameters.Since the system is highly non-linear, linearization and iteration is required, but the investigations performed have shown a stable and reliable convergence.
If the geolocation accuracy of the second orbit has the same quality level as the first, two datatakes from different orbits allow 3-D localization with centimeter level accuracy without the need for external information other than the corrections.Two datatakes over Wettzell can be acquired using TSX-1 and TDX-1 within just 5 days and by increasing the number of acquisitions for each orbit, the accuracy is also improved since the system becomes more and more over-determined.
Compared to established coherent SAR scatterer methods such as PSI (Persistent Scatterer Interferometry) which make use of the phase component and rely on precise coregistration, our refined stereo method exploits solely the magnitude and is therefore far less dependent on viewing angle.This allows larger baselines of several hundred kilometers which in turn provides a more stable geometrical configuration for absolute 3-D geolocation.
First practical results with 2 datatakes from the second Wettzell orbit with 47 degree incidence angle are very promising.By combining these measurements with 2 datatakes from the original 34 degree orbit, one obtains a stereo baseline of 168 kilometers.The actual passes used are from 3.2.2013,8.3.2013(34 degree, both TSX-1), 8.2.2013 (47 degree, TSX-1) and 2.3.2013(47 degree, TDX-1).Third degree polynomials obtained by least squares are used to model the TSX-1/TDX-1 phase center as a 3-D trajectory in ITRF2008 coordinates during a pass.These computations are based on the L1B product annotation for TSX-1/TDX-1 (Fritz et al., 2007).The range and azimuth radar times of the corner reflector are extracted from focused SAR images in the same way as for the 2-D case as described in (Balss et al. 2012).Finally, all the effects described in Sections 3 and 4 are transformed to radar range and azimuth times for each of the 4 datatakes and can thus be corrected in every observation.
The remaining geometrical inconsistencies are handled by the least squares estimate of the corner reflector position.In theory, the perpendicular circles described by the pairs of azimuth and range equations intersect at a single location as shown in Figure 5, but in practice they do not because of error in the along track (azimuth time) and across track (range time) components.By utilizing the a priori accuracies known from the 2-D experiment (±1 cm range, ±4 cm azimuth) in a weighted least squares procedure, the estimate is forced to rely more on the more precise range observations.The procedure is initialized with an approximate value for the corner reflector (simply a coordinate somewhere within the radar scene) and in successive iterations, the remaining inconsistencies are distributed between the conditions (range equation, zero-Doppler) and the parameter estimation (corner reflector coordinates) using the L2-norm.Convergence is reached after 4 to 5 iterations.The solution is the absolute 3-D position in the ITRF2008 reference frame at the epoch of the first datatake (3.2.2013) since continental drift is referenced to the first datatake of the stereo stack.The remaining residuals with respect to the "true" corner coordinates obtained from a terrestrial survey were -0.9, -5.8 and -3.5 cm in the X, Y and Z directions respectively.From the residuals and variance-covariance matrix of the least squares estimate, the standard deviations of the estimated corner reflector coordinates were 4.9, 13.3 and 10.2 cm.These values are reasonable since the system of equations is hardly over-determined.In fact it is close to its minimal configuration of only one data take per orbit.
These initial results illustrate the high potential of such accurate, absolute SAR range and azimuth observations, although one has to emphasize that these results were obtained with a 1.5 meter corner reflector.When applied to the persistent scatterers present in urban areas we expect a somewhat reduced accuracy.However, considering that for some locations stacks of 40 or more images are available, a highly over-determined least squares estimate can be set up.Investigations using such configurations are planned for the TSX-1 datasets of Munich.

 CONCLUSIONS AND OUTLOOK
Through an in-depth analysis of TSX-1 carried out as part of the Wettzell corner reflector experiment, the pixel geolocation accuracy has been improved to approximately 1 centimeter in range.This was accomplished by introducing precise local GNSS-based corrections for propagation delays and state of the art modeling of geophysical effects at the target.In addition, alternative model-based path delay corrections were implemented to ensure global applicability for when local GNSS data is not available.A potential application that can exploit this increased measurement accuracy is stereo SAR with prominent scatterers.Simulations and first practical results for the Wettzell corner reflector indicate an absolute 3-D geolocalization accuracy better than 10 centimeters.
The next stage of project will investigation the angular dependency of the corrections by installing a second corner reflector with different geometry at Wettzell as well as making acquisitions at an additional test site at O'Higgins Station (Figure 6), Antarctica, to evaluate location independence.Furthermore, TSX-1 has been scheduled to continue acquiring datatakes of the Wettzell corner reflector from the second 47 degree orbit in order to further investigate the 3-D performance of stereo SAR.In addition, the first practical experiments using this method in an urban environment are planned for existing datasets of Munich.

Figure 1 .
Figure 1.Corner reflector at the geodetic observatory at Wettzell, Germany.At present we have 30 TSX-1 and 4 TDX-1 acquisitions, both with a 34 degree incidence angle, of the Wettzell test site.In the first stage of our project, we intentionally focused on only one sensor (TSX-1) in order to avoid possible unknown systematic sensor differences.The study of such sensor dependencies is the subject of ongoing investigations.After correcting the measured radar time coordinates for signal propagation delays and geodynamic effects, we obtain a pixel localization accuracy (1σ) for TSX-1 of 38 millimeters in azimuth and either 11 or 14 millimeters in range, depending upon whether the atmospheric delay corrections are based on local GNSS measurements (Figure 2, top) or the globally applicable model-based concept using ECMWF and CODE data (Figure 2, bottom).Within the constraints of the sparse statistics of TDX-1 datatakes which were recorded up to now, we did not find any coarse deviation between the localization results of both sensors.It should be noted that the attained accuracy in range is at the level of the orbit accuracy which was determined by satellite laser ranging to have a mean error of -2 mm and a standard deviation of 15 mm.

Figure 2 .
Figure 2. Difference between SAR and GNSS coordinate estimates in radar geometry for the Wettzell corner reflector TSX-1 dataset (blue).1σ error bars in azimuth and range for the entire dataset (red).Local IGS corrections (top) and global ECMWF and CODE corrections based (bottom).A breakdown of the contributions of the geodynamic effects and signal propagation delays to the localization accuracy is shown in Figure3.The largest contributions are in range from the troposphere and solid earth tides, both with 1σ values of 6 cm.Next comes the ionosphere in range and the solid earth tides in azimuth, both being approximately 2 cm.The remaining contributions are between 0.2 and 0.5 cm.

Figure 4 .
Figure 4. Temporal progression of the pixel localization error in range for the Wettzell corner reflector TSX-1 dataset using local IGS corrections.

International
Figure 5. Intersection of the circles (red) described by the azimuth and range equations for two orbits (blue) in the ITRF2008 reference frame at a point in the SAR image (center of green ellipse).