The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences. Vol. XXXVII. Part B7. Beijing 2008 
103 
process of function optimization, which is a really tough task in 
data processing of the conventional DInSAR. 
3.2 Parameter Adjustment by LS Solution 
With equation (3) we can compute the increments of elevation 
errors and linear deformation velocities along all the arcs in the 
network. By trials with simulated data, we have found that the 
arcs have an accurate solution for As and Av if y is larger 
than 0.45. The network can then be treated in a similar way as a 
leveling or GPS network; a weighted LS adjustment is applied 
to eliminate geometric inconsistency due to uncertainty in phase 
data, and thus obtaining the most probable values of the linear 
deformation velocities and elevation errors at PSs. 
Taking the adjustment of LOS linear deformation velocities as 
an example, we present here some basic ideas for computation. 
A prototype observation equation for an arc is expressed as 
Vp-v, =Av pl +r pl , p*l,Vp,l = \,2,K ,K (4) 
where v p = linear deformation velocity at PS p 
v, = linear deformation velocity at PS / 
K= total number of all the PSs 
Suppose there are Q arcs available in the network, we will have 
Q observation equations. The MC value of each arc can be used 
as a weight. With a weighted LS solution, the motion velocities 
at all the K PS can be eventually obtained. Such procedure can 
also be applied in a similar way onto the elevation- 
inconsistency network to estimate the elevation errors at all the 
K PSs. The Kriging interpolator can be applied to generate the 
deformation-velocity map and the elevation-error map. 
As a remark, we underline that a reference point without motion 
or elevation error should be selected according to a priori 
information to obtain a unique solution with LS adjustment, and 
thus making all the estimates be related to the benchmark. In 
addition, we solved the large sparse matrix system using the 
software package termed UMFPACK (Davis, 2002). Moreover, 
it should be noted that the FCN used here is much stronger in 
terms of reliability than the TIN. Our simulation study shows 
that the LS solution derived with the FCN is much more 
accurate than that derived with the TIN even though a small 
portion of measurements ( Ae , Av ) are intentionally set as 
outliers. This is because the redundancy number in the FCN is 
significantly larger than that in the TIN. 
3.3 Separation of Nonlinear Deformation and Atmospheric 
Effect 
The further analysis focuses on isolating nonlinear motion from 
atmospheric delay. For each interferometric pair, the residual 
phase increment (gradient) at each arc can be first derived. The 
integration of gradients (i.e., phase unwrapping) of all the arcs 
in the network is then performed by a weighted least squares 
method (Ghiglia & Pritt, 1998), and thus obtaining the residual 
phases in absolute sense at all the PS pixels for each pair. As 
mentioned before, the residual phase is contributed by nonlinear 
deformation, atmospheric delay and decorrelation noise. 
It is possible to separate nonlinear motion from atmospheric 
delay because the two terms have different spectral structure in 
space and time domain (Ferretti et al., 2000; Mora et al., 2003). 
In terms of atmospheric perturbation, a high correlation exhibits 
in space, but a significantly low correlation presents in time. In 
terms of nonlinear deformation, a strong correlation exists in 
space and a high correlation occurs in time. It is however not 
easy to discriminate the spectral bands between the nonlinear 
deformation and the atmospheric effect without availability of a 
priori information. Although an exact separation of the two 
terms is a challenged task, we try to achieve such purpose by 
introducing a strategy referred to as empirical mode 
decomposition (EMD) which was pioneered by Huang et al. 
(1998). This method is in principle different from the cascade 
filter applied by the previous studies (Ferretti et al., 2000, 2001; 
Mora et al., 2003). 
To decouple signatures at each PS, we first estimate the time 
series of unwrapped residual phases corresponding to all the 
SAR acquisition times. This can be done using a singular value 
decomposition (SVD) method (for details, see Berardino et al., 
2002). The separation between atmospheric delay and nonlinear 
motion is then conducted by EMD using the time series of 
unwrapped residual phases. 
The dataset being dealt with generally has several features: (a) 
the tens of samples are irregular; (b) the data is non-stationary 
(varying undulation); and (c) the data represents a nonlinear 
process. Theoretical study by Huang et al. (1998) indicated that 
the EMD approach is more advantageous for dealing with a 
nonlinear and nonstationary process than other signal analysis 
tools like Fourier transform. The EMD method separates the 
signal into a collection of intrinsic mode functions (IMF) that 
are useful for revealing some physical properties. For details of 
EMD operation, see Huang et al., 1998. 
By trials with real data, we can generally extract four IMFs 
from the time series of unwrapped residual phases. The IMFs 
with high frequencies correspond to atmospheric component, 
while the IMFs with low frequencies reflect nonlinear motion. 
It is not so direct to thoroughly discriminate them in the case of 
no a priori information available. However, reasonable results 
have been produced in our case when we consider the sum of 
the first and second level of IMFs as the atmospheric phases 
and treat the sum of the third and fourth level of IMFs as the 
nonlinear motions. The obtained two time series of phases for 
nonlinear motion and atmospheric delay at each PS are 
expressed as 
r ,def {t) = W n,def {t 0 ),r ,def {t x ),A (5) 
y,°‘ m {t) = \r ,m {t»\r ,m {t,),A ,r ,m {t N )} (6) 
where y/ nldef (t i ) = phase for nonlinear motion at (imaging time) 
y/ a,m (t i ) = phase for atmospheric delay at t,■ (imaging time) 
Finally, the total LOS deformation at any PS can be computed 
by summing the linear and nonlinear components, which can be 
written as a function of SAR acquisition time t t by 
def{t i ) = {t i -t Q )-v + ^-^ Mef ( ti ) (7)