针对SAR数据高噪声及大梯度变化区域的相位解缠方法与流程

针对sar数据高噪声及大梯度变化区域的相位解缠方法
技术领域
1.本发明涉及一种相位解缠方法，尤其适用于处理雷达卫星sar数据使用的一种针对sar数据高噪声及大梯度变化区域的相位解缠方法。


背景技术：

2.合成孔径雷达干涉(interferometric syntheticaperture radar,insar)测量技术近年来已广泛应用于地形三维重建和由地震活动的物理模型反演，尽管目前已经存在几十种相位解缠方法，但是在噪声量大的区域相位解缠方法的鲁棒性一般较差，在地形梯度大的区域相位解缠结果误差一般较大。这几十种相位解缠方法可以大致分为两类：一类是单基线相位解缠方法，另一类是多基线相位解缠方法。由于单基线相位解缠方法受相位连续性假设的限制，在地形梯度大的区域难以获得理想的结果，虽然多基线相位解缠方法可以摆脱相位连续性假设的限制，但是在解缠效率和基线配比等方面仍然存在一些困难。
3.由于修正矩阵束模型具有一定的频率估计精度和效率，所以本发明保留了修正矩阵束模型的框架，在其奇异值处理过程中引入了低通滤波技术，并且可以根据干涉图像素均方差值自适应选择局部频率估计窗口的大小，尽可能地实现方位向和距离向相位梯度的准确估计。这种低通滤波技术是由butterworth滤波器来实现的，其特点是通频带内的频率响应曲线最大限度平坦，没有纹波，而在阻频带则逐渐下降为零，这意味着butterworth低通滤波器可以减少高频噪声。此外，我们也引入了自适应平方根无迹kalman技术到相位解缠问题中，平方根无迹kalman算法只需要误差协方差的平方根，而不需要完整的协方差阵，其采用误差协方差阵的平方根形式替代协方差阵进行递推计算，保证了误差协方差阵的非负定性，从而提高了数值的稳定性，并且在计算中只需要存储和运算平方根因数，这样可以降低计算负担，以实现对无迹kalman相位解缠技术性能的提高。


技术实现要素：

4.针对现有技术的不足之处，提供一种步骤简单，使用效果好的一种针对sar数据高噪声及大梯度变化区域的相位解缠方法
5.为实现上述技术目的，本发明的针对高噪声及大梯度变化区域的相位解缠方法，所包括以下步骤：
6.步骤s1，收集需要提取高程量或微小形变量地区的sar数据，对sar数据进行预处理得到干涉相位图；通过改进的修正矩阵束模型对干涉相位图所有像元的局部频率进行估计，改进的修正矩阵束是在修正矩阵束模型的奇异值处理步骤中替用一阶butterworth低通滤波器对分解后得到的每一个奇异值赋予不同的权值来实现信号重构，从而获取准确的局部频率，以此对干涉相位图中所有像元的方位向和距离向的相位梯度进行估计：局部频率估计窗口的大小根据干涉图像元相位均方差自适应选择，窗口包括小的窗口和大的窗口，像元的干涉相位均方差反映条纹的密集程度，条纹密集则均方差大，利用小的窗口进行局部频率估计以免破坏局部条纹平稳的假设；条纹稀疏则均方差小，利用大的窗口进行局
部频率估计，既满足局部条纹平稳的假设，又增加了样本点数量，提高了频率估计结果的准确性；
7.步骤s2，由于局部频率估计结果不存在大的跳变，当局部频率估计结果出现明显的不连续时，则认定该局部像元的频率估计结果存在误差，需要进行相应地修正，利用下式从距离向和方位向上分别对频率估计结果的连续性进行计算和检测；
[0008][0009][0010]
式中m为与判定连续性窗口大小相关的参数，例如m＝3，即表示判定连续性窗口大小为(2m 1)
×
(2m 1)＝7
×
7；u和v表示局部窗口内像元偏离k位置中心像元的程度，和分别表示干涉相位图中k位置中心像元的频率估计结果在距离向和方位向的可信值，f
x
和f
y
分别是距离向和方位向上的频率估计结果；
[0011]
步骤s3，根据求得的可信值和设置阈值，分别设置为和其中和为的最大值和的最大值，利用下式将干涉相位图中k位置中心像元的频率估计结果优化为以k位置像元为中心的窗口(2m 1)
×
(2m 1)范围内所有像元频率的平均值：
[0012][0013][0014]
其中和表示修正后干涉相位图中k位置中心像元在方位向和距离向的窗口范围内所有像元频率的平均值，从而解决了局部频率估计结果存在异常值的问题；
[0015]
步骤s4，建立无迹kalman滤波相位解缠观测方程的自适应平方根模型，并利用无迹kalman滤波相位解缠观测方程的自适应平方根模型获得相位解缠起始像元；首先计算自适应平方根无迹kalman相位解缠模型的sigma采样点及其权值系数；将sigma点通过无迹kalman滤波相位解缠状态方程计算传播后均值和误差协方差的平方根，利用误差协方差阵的平方根形式替代协方差阵进行递推计算，从而保证了误差协方差阵的非负定性，提高了数值的稳定性，并且在计算中只需要存储和运算平方根因数，从而降低计算负担；
[0016]
步骤s5，计算无迹kalman滤波相位解缠观测方程的自适应平方根模型的量测预测值、量测预测值协方差平方根和状态与量测间的互协方差；其中采用测量噪声等价协方差平方根阵原理，引入自适应因子，实现对测量噪声协方差平方根的异常值分离和自适应修正，进而提高平方根无迹kalman相位解缠模型的精度；
[0017]
步骤s6，计算无迹kalman滤波相位解缠观测方程的自适应平方根模型的滤波增益；利用求得的kalman滤波增益矩阵对状态估计及其协方差平方根进行更新，获取鲁棒性
更好的无迹kalman相位解缠模型；
[0018]
步骤s7，通过基于堆排序的有效质量指导策略对干涉相位图进行逐个相位展开，将干涉复数信号所有像元的相位质量，即质量图作为指导策略，以最大堆作为数据结构，通过不断更新最大堆的节点来实现对质量图的排序，保证了相位解缠过程一直从高质量像元向低质量像元计算；干涉相位图中每一个像元的质量是由相干系数表示，相干系数位于区间[0，1]，相干系数越大则代表像元的质量越高，相反若相干系数越小则代表像元的质量越低，所有像元的相干系数组合在一起即为质量图；重复步骤s4
‑
步骤s6，直到干涉相位图中所有像元都被解缠。
[0019]
进一步，使用改进的修正矩阵束模型对干涉图距离向和方位向梯度进行估计：
[0020]
步骤a1，计算干涉图中所有像元干涉相位的均方差ξ，公式为
[0021][0022]
其中l
×
l为计算像元干涉相位均方差的一个窗口大小，即判断条纹密集程度的一个窗口大小，一般可设置l＝5，a
s
表示局部窗口内第s个像元的干涉相位复数的振幅；为a的平均值；
[0023]
步骤a2，使用公式(1)求得的干涉相位均方差ξ，根据ξ值自适应选择待估像元的局部频率估计窗口大小，其中7
×
7或9
×
9为小的窗口，17
×
17或19
×
19为大的窗口，公式为：
[0024][0025]
其中，h
×
h表示干涉相位图中k位置中心像元进行局部频率估计计算的窗口大小；
[0026]
步骤a3，假设i
k
是以k像元为中心，窗口大小为h
×
h的干涉相位样本，利用公式：i
k
＝udv
t
对信号矩阵i
k
进行奇异值分解，其中u和v表示h阶正交(酉)矩阵；d表示对应的特征值矩阵，d＝diag(σ1,σ2,
…
σ
h
)，且σ1≥σ2≥
…
≥σ
h
≥0，σ
i
,i＝1,2,
…
,h表示信号矩阵i
k
奇异值分解后对角线上的特征值；
[0027]
步骤a4，由一阶butterworth低通滤波器计算传递函数，即对信号矩阵i
k
分解后对角线上的每个特征值赋予不同的权值以降低高频噪声对信号重构的影响，具体计算过程为：
[0028][0029]
式中w
t
表示特征值矩阵d对角线上第t个奇异值即将被分配的权值系数，σ
a
表示前t个奇异值中第a个奇异值，σ
t
表示特征值矩阵d对角线上第t个奇异值；
[0030]
步骤a5，利用公式(4)计算处理后的干涉相位样本矩阵即尽量剔除信号中高频噪声成分，尽量完整地保留有用的相位信息：
[0031]
[0032]
步骤a6，利用公式(5)根据处理后的干涉相位样本矩阵构建母hankel矩阵以及两个子hankel矩阵：
[0033][0034]
其中及分别为母hankel矩阵以及两个子hankel矩阵；
[0035]
步骤a7，利用公式(6)对母hankel矩阵进行奇异值分解：
[0036][0037]
其中和仅具有与有用信号相关的主特征值对应的信息；
[0038]
步骤a8，通过奇异值分解，并结合公式(5)和公式(6)得到降维处理后的相位数据信号，由于相位数据矩阵具有低秩结构，经过降维处理后以抑制噪声，推导获得如下公式：
[0039][0040]
其中及分别是通过母hankel矩阵以及两个子hankel矩阵和重构后的信号矩阵；
[0041]
步骤a9，利用公式(8)表示以k像元为中心窗口内样本干涉相位的距离向和方位向的局部频率：
[0042][0043]
其中和分别为和的伪逆矩阵，conj(
·
)为矩阵的共轭相乘因子；和分别是k位置像元在距离向和方位向上的频率估计结果；
[0044]
步骤a10，由于局部频率估计结果不存在大的跳变，当频率估计结果出现明显的不连续时，则认定该局部像元的频率估计结果存在误差，需要进行相应地修正，利用下式从距离向和方位向上分别对频率估计结果的连续性进行计算和检测；
[0045][0046]
式中m为与判定连续性窗口大小相关的参数，例如m＝3，即表示判定连续性窗口大小为(2m 1)
×
(2m 1)＝7
×
7；u和v表示局部窗口内像元偏离k位置中心像元的程度，和分别表示干涉相位图中k位置中心像元的频率估计结果在方位向和距离向的可信值，f
x
和f
y
分别是像元在距离向和方位向上的频率估计结果；
[0047]
步骤a11，根据求得的可信值和设置阈值，分别设置为和其中和为的最大值和的最大值，利用下式计算该像元的频率估计结果优化为以k位置像元为中心的窗口(2m 1)
×
(2m 1)范围内所有像元频率的平均值：
[0048][0049]
其中f
x
和f
y
分别是像元在距离向和方位向上的频率估计结果；和表示修正后像元k在距离向和方位向的频率估计结果，以解决局部频率估计结果存在异常值的问题；
[0050]
步骤a12，根据步骤a11获得局部频率可以获得像素k的距离向和方位向梯度，公式如下：
[0051][0052]
其中和分别为像素k的距离向和方位向梯度；
[0053]
重复步骤a2
‑
步骤a12，直至估计出干涉相位图中所有像元距离向和方位向的梯度。
[0054]
进一步，建立自适应平方根无迹kalman滤波相位解缠模型，获得相位解缠起始像元的步骤为：
[0055]
步骤b1，利用下式计算自适应平方根无迹kalman相位解缠模型解缠起始像元的sigma采样点及其权值系数，公式为：
[0056][0057]
其中ρ
i
为起始像元的sigma采样点，是第k
‑
1位置解缠相位的状态估计值，s
k
‑1是第k
‑
1位置像元的误差协方差的平方根，n表示状态变量的维数(n＝1)；{
·
}
i
表示括号内矩阵的第i列，η表示尺度函数；
[0058][0059]
其中和为权值系数，η＝τ2(n α)
‑
n表示尺度函数，τ，α，ζ用来调节sigma点，通常τ＝0.01,ζ＝2,α＝0，τ决定sigma点在均值附近的散布程度，α为冗余量；
[0060]
步骤b2，将sigma点通过无迹kalman滤波相位解缠状态方程计算传递后均值和误差协方差的平方根，公式为：
[0061][0062]
其中f[
·
]表示状态方程转移函数；符号qr{
·
}和cholupdate{
·
}表示qr分解运算和cholesky因子一阶更新运算，表示第k个像元sigma点预测值；表示干涉相位图第k个像元的状态估计；s
ψψ
表示干涉相位图第k个像元状态预测误差协方差平方根；q(k
‑
1)是干涉相位图第k个像元的相位梯度估计误差方差。
[0063]
进一步，计算自适应平方根无迹kalman滤波相位解缠模型的量测预测值、量测预测值协方差平方根和状态与量测间的互协方差；采用测量噪声等价协方差平方根阵原理对测量误差方差的平方根进行自适应修正的具体步骤为：
[0064]
步骤c1，计算待解缠像元的量测预测值，公式为：
[0065][0066]
其中，表示k
‑
1像元的sigma点经传递后得到的k像元的sigma点预测值；表示待解缠像元的量测预测值；h[
·
]为观测方程系数。
[0067]
步骤c2，利用公式(16)和公式(17)计算待解缠像元的量测预测值协方差的平方根和状态与量测间的互协方差：
[0068][0069][0070]
其中表示待解缠像元k的量测预测值协方差的平方根；表示待解缠像元k的状态与量测间的互协方差；
[0071]
步骤c3，利用下式计算待解缠像元的预测残差信息：
[0072][0073]
式中为k位置像元的预测残差，表示k位置像元的真实缠绕相位，为待解缠像元的量测预测值；
[0074]
步骤c4，利用下式计算待解缠像元测量噪声协方差的等价协方差阵的对角线元
素：
[0075][0076]
其中r
i
表示测量噪声协方差阵r的第i个对角线元素；表示测量噪声协方差的等价协方差阵第i个对角线元素；u0和u1为小数值的常数阈值，这里可设置u0＝0.45,u1＝3；v
i
为标准化预测残差。表示预测观测值协方差阵的第i个对角线元素；median(
·
)表示取中位数算子。
[0077]
步骤c5，将公式(19)求得的对角线元素组成的测量噪声协方差的等价协方差阵带入式(16)中替换测量噪声协方差阵r重新更新计算量测预测值的协方差平方根。
[0078]
进一步，计算自适应平方根无迹kalman滤波相位解缠模型的滤波增益，对状态估计及其协方差平方根进行更新，公式为
[0079][0080]
其中υ
k
表示干涉相位图第k个像元的增益矩阵；表示干涉相位图第k个像元的状态估计值；s(k)表示干涉相位图第k个像元状态估计的协方差平方根。
[0081]
有益效果：
[0082]
本发明通过在修正矩阵束模型中嵌入一个butterworth低通滤波器，并可以根据干涉图的相位均方差自适应调整局部频率估计分析窗的大小，然后对频率估计结果异常值进行修正以获取高精度的距离向和方位向梯度；通过自适应平方根无迹kalman相位解缠方法结合最大堆解缠策略得到解缠相位。相比于其它常规相位解缠方法，本方法对于高噪声区域和地形大梯度变化区域的相位解缠具有更好的解缠精度和模型鲁棒性，有利于提高最终insar数据处理结果的精度。
附图说明
[0083]
图1是本发明的针对高噪声及大梯度变化区域的相位解缠方法流程图流程图；
[0084]
图2为本发明研究区的地理位置区域图；
[0085]
图3为本发明研究区的干涉相位图；
[0086]
图4为本发明距离向的相位梯度估计结果示意图；
[0087]
图5为本发明方位向的相位梯度估计结果示意图；
[0088]
图6为本发明修正后的距离向相位梯度估计结果示意图；
[0089]
图7为本发明修正后的方位向相位梯度估计结果示意图；
[0090]
图8为本发明针对高噪声及大梯度变化区域的相位解缠方法的解缠结果示意图；
[0091]
图9为mcf方法的相位解缠结果示意图；
[0092]
图10为snaphu方法的相位解缠结果示意图；
[0093]
图11为无迹kalman相位解缠方法的相位解缠结果示意图；
[0094]
图12为由高精度dem生成的模拟解缠相位图。
具体实施方式
[0095]
下面结合附图，对本发明作进一步的阐述：
[0096]
如图1所示，本发明的一种针对sar数据高噪声及大梯度变化区域的相位解缠方法，输入信息为研究区的地理位置区域图的干涉相位图，其中实际地理位置区域图如图2所示，输入的干涉相位图如图3所示，包括一个干涉像对经过粗配准、精配准、去平/去地形、复共轭相乘等处理过程，所述方法包括以下步骤：
[0097]
步骤s1，使用改进的修正矩阵束模型对干涉图所有像元的距离向和方位向梯度进行估计。
[0098]
步骤s1.1，计算干涉图中所有像元干涉相位的均方差ξ，公式为：
[0099][0100]
其中l
×
l为计算像元干涉相位均方差的一个窗口大小，即判断条纹密集程度的一个窗口大小，一般可设置l＝5；a
s
表示局部窗口内第s个像元的干涉相位复数的振幅；为a的平均值；
[0101]
步骤s1.2，使用公式(1)求得的干涉相位均方差ξ，根据ξ值自适应选择待估像元的局部频率估计窗口大小，其中7
×
7或9
×
9为小的窗口，17
×
17或19
×
19为大的窗口，公式为：
[0102][0103]
其中，h
×
h表示干涉相位图中k位置中心像元进行局部频率估计计算的窗口大小；
[0104]
步骤s1.3，假设i
k
是以k像元为中心，窗口大小为h
×
h的干涉相位样本，利用公式：i
k
＝udv
t
对信号矩阵i
k
进行奇异值分解，其中u和v表示h阶正交(酉)矩阵；d表示对应的特征值矩阵，d＝diag(σ1,σ2,
…
σ
h
)，且σ1≥σ2≥
…
≥σ
h
≥0，σ
i
,i＝1,2,
…
,h表示信号矩阵i
k
奇异值分解后对角线上的特征值。
[0105]
步骤s1.4，由一阶butterworth低通滤波器计算传递函数，即对信号矩阵i
k
分解后对角线上的每个特征值赋予不同的权值以降低高频噪声对信号重构的影响，具体计算过程为：
[0106][0107]
式中w
t
表示特征值矩阵d对角线上第t个奇异值即将被分配的权值系数，σ
a
表示前t个奇异值中第a个奇异值，σ
t
表示特征值矩阵d对角线上第t个奇异值；
[0108]
步骤s1.5，利用公式(4)计算处理后的干涉相位样本矩阵即尽量剔除信号中高频噪声成分，尽量完整地保留有用的相位信息：
[0109][0110]
步骤s1.6，利用公式(5)根据处理后的干涉相位样本矩阵构建母hankel矩阵以及两个子hankel矩阵：
[0111][0112]
其中及分别为母hankel矩阵以及两个子hankel矩阵；
[0113]
步骤s1.7，利用公式(6)对母hankel矩阵进行奇异值分解：
[0114][0115]
其中和仅具有与有用信号相关的主特征值对应得信息；
[0116]
步骤s1.8，通过奇异值分解，结合公式(5)和公式(6)可以得到降维处理后的相位数据信号，由于相位数据矩阵具有低秩结构，经过降维处理后可以抑制噪声，推导获得如下公式：
[0117][0118]
其中及分别是通过母hankel矩阵以及两个子hankel矩阵和重构后的信号矩阵；
[0119]
步骤s1.9，利用公式(8)表示以k像元为中心窗口内样本干涉相位的距离向和方位向的局部频率：
[0120][0121]
其中和分别为和的伪逆矩阵，conj(
·
)为矩阵的共轭相乘因子；和分别是k位置像元在距离向和方位向上的频率估计结果；
[0122]
步骤s1.10，根据步骤s1.9获得局部频率可以获得像素k的距离向和方位向梯度，公式如下：
[0123][0124]
其中和分别为像素k的距离向和方位向梯度；
[0125]
重复步骤s1.2
‑
s1.10，直至估计出所有像元距离向和方位向的梯度。方位向和距离向的相位梯度估计结果如图4和图5所示。
[0126]
步骤s2，根据公式(9)，从距离向和方位向上分别对频率估计结果的连续性进行计算和检测。阈值分别设置为和其中和为的最大值和的最大值。当可信值超过阈值时，认为该像元的频率估计结果存在异常，需要进行修正。
[0127][0128]
式中m为与判定连续性窗口大小相关的参数，例如m＝3，即表示判定连续性窗口大小为(2m 1)
×
(2m 1)＝7
×
7；u和v表示局部窗口内像元偏离k位置中心像元的程度，和分别表示干涉相位图中k位置中心像元的频率估计结果在方位向和距离向的可信值，f
x
和f
y
分别是像元在距离向和方位向上的频率估计结果；
[0129]
步骤s3，从距离向和方位向上分别对频率估计结果的异常值进行修正，异常像元的频率估计结果优化为以k位置像元为中心(2m 1)
×
(2m 1)范围内所有像元频率的平均，公式为
[0130][0131]
其中f
x
和f
y
分别是像元在距离向和方位向上的频率估计结果；和表示修正后像元在方位向和距离向的频率估计结果。由公式(10)的频率结果重新根据公式(28)计算修正后的方位向和距离向相位梯度。本发明实例干涉图在方位向和距离向上的相位梯度修正后的结果如图6和图7所示。
[0132]
步骤s4，建立自适应平方根无迹kalman相位解缠模型；获得相位解缠起始像元，假设解缠起始像元在干涉图中的位置为k
‑
1，待解缠像元的位置为k。
[0133]
步骤s4.1，计算解缠起始像元的sigma采样点及其权值系数，公式为
[0134][0135]
其中ρ
i
为起始像元的sigma采样点，是第k
‑
1位置解缠相位的状态估计值，
s
k
‑1是第k
‑
1位置像元的误差协方差的平方根，n表示状态变量的维数(n＝1)；{
·
}
i
表示括号内矩阵的第i列，η表示尺度函数；
[0136][0137]
其中和为权值系数，，η＝τ2(n α)
‑
n表示尺度函数，τ，α，ζ用来调节sigma点，通常τ＝0.01，α＝0，ζ＝2，τ决定sigma点在均值附近的散布程度，α为冗余量。
[0138]
步骤s4.2，将解缠起始像元的sigma点通过自适应平方根无迹kalman滤波相位解缠模型的状态方程计算传播后的均值和误差协方差的平方根，公式为
[0139][0140]
其中符号qr{
·
}和cholupdate{
·
}表示qr分解运算和cholesky因子一阶更新运算，表示第k个像元sigma点预测值；表示干涉相位图第k个像元的状态估计；s
ψψ
表示干涉相位图第k个像元状态预测误差协方差平方根；q(k
‑
1)是干涉相位图第k个像元的相位梯度估计误差方差。
[0141]
步骤s5，计算待解缠像元的量测预测值、量测预测值协方差平方根和状态与量测间的互协方差；采用测量噪声等价协方差平方根阵原理对测量误差方差的平方根进行自适应修正。
[0142]
步骤s5.1，计算待解缠像元的量测预测值，公式为：
[0143][0144]
其中，表示k
‑
1像元的sigma点经传递后得到的k像元的sigma点预测值；表示待解缠像元的量测预测值；h[
·
]为观测方程系数。
[0145]
步骤s5.2，计算待解缠像元的量测预测值协方差的平方根和状态与量测间的互协方差，公式为
[0146]
[0147][0148]
其中表示待解缠像元k的量测预测值协方差的平方根；表示待解缠像元k的状态与量测间的互协方差；
[0149]
步骤s5.3，计算待解缠像元的预测残差信息，公式为
[0150][0151]
式中为k位置像元的预测残差，表示k位置像元的真实缠绕相位，为待解缠像元的量测预测值；
[0152]
步骤s5.4，计算待解缠像元的测量噪声协方差的等价协方差阵的对角线元素，公式为
[0153][0154]
其中r
i
表示测量噪声协方差阵r的第i个对角线元素；表示测量噪声协方差的等价协方差阵第i个对角线元素；u0和u1为小数值的常数阈值，这里可设置u0＝0.45,u1＝3；v
i
为标准化预测残差。表示预测观测值协方差阵的第i个对角线元素；median(
·
)表示取中位数算子。
[0155]
步骤s5.5，将公式(19)求得的对角线元素组成的测量噪声协方差的等价协方差阵带入式(16)中替换测量噪声协方差阵r重新更新计算量测预测值的协方差平方根。
[0156]
步骤s6，计算自适应平方根无迹kalman滤波相位解缠模型的滤波增益，对状态估计及其协方差平方根进行更新，公式为
[0157][0158]
其中υ
k
表示干涉相位图第k个像元的增益矩阵；表示干涉相位图第k个像元的状态估计值；s(k)表示干涉相位图第k个像元状态估计的协方差平方根。
[0159]
步骤s7，通过基于堆排序有效质量指导策略对干涉图进行逐个相位展开；重复步骤s4
‑
s6，直到所有像元都被解缠。最终得到干涉图的相位解缠结果。相位解缠的结果如图8所示。
[0160]
为本发明针对sar数据高噪声及大梯度变化区域的相位解缠方法的示例性实施例更加有说服力，我们对此同样的干涉图分别采用mcf、snaphu、无迹kalman相位解缠等方法进行解缠处理。mcf方法的解缠结果如图9所示；snaphu方法的解缠结果如图10所示；无迹kalman相位解缠方法的解缠结果如图11所示。从图中可以看出，本专利提出的方法可以获得较其它方法更好的解缠结果；为了定量描述相位解缠的质量，我们将不同方法的解缠结果与由高精度dem生成的模拟干涉图(图12)进行差分处理，计算各种相位解缠方法的解缠结果的平均误差，其中，mcf方法的解缠结果平均误差为0.8411rad；snaphu方法的解缠结果的平均误差为0.7822rad，无迹kalman相位解缠方法的解缠结果的平均误差为0.6981rad，本发明的一种针对sar数据高噪声及大梯度变化区域的相位解缠方法的解缠结果的平均误差为0.5948rad。结果表明，本发明的一种针对sar数据高噪声及大梯度变化区域的相位解缠方法具有更好的相位解缠精度。