专利名称：基于极化干涉合成孔径雷达数据优化干涉相干系数的方法
技术领域：
本发明涉及合成孔径雷达遥感信息处理技术领域，是一种基于单基线 极化干涉合成孔径雷达数据的优化干涉相干系数的方法，它能简单、快速 实现干涉相干系数的优化并提取对应于散射中心的干涉相位。
背景技术：
极化干涉合成孔径雷达技术将全极化数据引入干涉处理，不但可利用 极化信息对散射体形状、指向、介电系数等参数的敏感性[l]，而且还可基
于干涉相位获得高度、运动、形变信息，因此它比单独的干涉或极化SAR 在信息提取、分类等应用中具有更突出的优势，现已成为SAR领域的研 究热点之一。 参考文献 Papathanassiou LP.and Cloude S.R. (2001)， Single-Baseline Polarimetric SAR Interferometry， 7hg/75". o"
6"eoscie/ ce朋d j e膨te 5"e/ S7.恥 39(11): 2352-2363. Cloude S. R. and Papathanassiou K. P. (1998)， Polarimetric SAR Interferometry [J]， IEEE Trans, on Geoscience and Remote Sensing, 36(5) : 1551—1565. Maxim Neumann, et al. (2006)， PolInSAR Coherence Set Theory and Application, ^7)S4y Dresden: Germany. Papathanassiou K. P. (1999), Polarimetric SAR Interferometry, Doctor Dissertation, Institut filr Hochfrequenztechnik des DLR， Oberpfaffenhofen.[5] Colin C. ， et al. (2003)， Investigation on Different Coherence Optimization Methods, t力e _T￡ i77tezv7aJtio/7sJ /for^s/ o;
J"ter/er譜etrj"， Frascati: ESA. Sagues L, et al. (2000)， Indoor Experiments on Polariraetric SAR Interferometry, Zfi^ 7]r朋51. o/ ( eo cj'e/7ce朋c/ 7 柳Me Se/^i/^; 38(2): 671-684. Colin C. ， et al. (2006), An interferometric coherence opmtimization method in radar polarimetry for high-resolution imagery, 7!ra/751.6"eoscj'e77ce a/ J j e膨te 5"e/7si恥 44(1):
167-175. LuBai， . et al. (2006)， An Improved Coherence Optimization MethodinPolarimetricSARInterferometry, ihte^raaf 6b/7/ez"e/7ce a/ yfec/az; Shanghai: China. TaoXiong, et al. (2006)， A New Approach for DEM Generation Based on Polarimetric SAR Interferometry, 7/7tei77atio朋J C0/7/are77ce or 7 aokr， Shanghai: China.
作为极化干涉SAR数据处理中一项关键技术，优化相干性的方法利
用干涉相干系数对极化状态的依赖性[2]，通过改变极化状态优化千涉相干
系数、提高干涉相位的质量。从数学的角度讲，优化相干性是在所有复矢 量或者极化状态的空间中^,^)l5^scv,、哨、叫，寻找使干涉相干系数达
到最大的一组投影矢量或极化状态的过程[3]:
优化相干性在改善干涉性能的同时，还可提取出最接近于点状散射体 的多个确定性散射机理，从而将高度上分布的有效散射相位中心等效为点 处理，抑制体散射去相干。该方法通过选择与散射中心散射机理有关的最
优投影矢量A与A，补偿时间去相干、优化信噪比、改善噪声去相干[4]。因此，基于单基线极化干涉数据的优化相干性的方法在分类、参数反演应 用中具有巨大的应用潜力。
Papathanassiou于1997年引入Lagrange乘数法求解最优投影矢量，该 方法可采用两次特征分解或一次奇异值分解即可求解除最优投影矢量，然 而该算法应用前提为三维相干矩阵t,'、 t。可逆[2，4]。由于特征矢量绝对相
位不唯一，因而还需引入额外的限定条件。rgKdO-()才能避免主、副天线
选用最优投影矢量对干涉相位估计的影响。由于该方法未对投影矢量作严 格的限定，因而造成临近的位置的最优矢量间变化程度较大，导致干涉相
位不连续，说明其提取的干涉相位不一定准确[5]。
由于主要散射中心的散射机理通常是稳定不变的，因此应选用相同的 极化状态或投影矢量计算干涉相干系数，从而出现了约束投影矢量对内两 矢量相同的优化相干性的方法。与之前的极化全局空间的优化方法相比，
约束类方法解决的是极化子空间的优化问题。Sagues从不同极化基中搜索 出具有最高干涉相干性的极化状态作为优化结果，该方法计算量大且其估 计精度受搜索步长的影响[6]。
Colin从Lagrange函数求解出发，通过引入完全相干假设将最优投影
矢量的求解过程简化为一次特征矢量的求解过程。该算法要求矩阵(t,'+t^
可逆，且受到完全相干假设的限制，其适用于处理基线去相关、时间去相 关影响较小的数据[5]。另一种约束类相干性优化方法，利用近似关系 '1^122-(71| + ^/2在矩阵 +1;2)可逆条件下将优化问题转换为数值距离的计 算；该算法同样受到^,, + 122)可逆的限制，且由于该算法利用数值逼近求解， 因而其估计结果的准确性与初始估计、搜索步长密切相关[7]。此外，Lu Bai尝试利用特征分解结果作为初始估计，并利用干涉相干性能否改善作为结
果收敛条件，但是该方法仍要求矩阵(1"11 + 122)可逆[8]。 Tao Xiong比较HH、 HV、 VV三种极化状态的干涉结果中具有最大干涉相干性的干涉相位作为 最优结果[9]，虽然该算法运算效率高，但其干涉相干性的改善程度有限。 本发明方法通过分析一个极化干涉数据矩阵提取最优投影矢量对，且 不要求该极化干涉数据矩阵可逆，因而普遍适用于各类极化干涉数据。本 发明方法无需补偿最优投影矢量间的绝对相位差，简化了处理；并通过检 验优化结果的有效性，可确保优化后干涉相干系数以及干涉相位结果的可 靠性。
发明内容
本发明的目的在公开一种基于极化干涉合成孔径雷达数据优化干涉 相干系数的方法，该方法计算速度快、相干性改善程度大、适用于对不同 类型的极化干涉数据的优化相干系数。
为达到上述目的，本发明的技术解决方案是
基于极化干涉合成孔径雷达数据优化干涉相干系数的方法，其包括步 骤如下
A) 利用快速傅里叶变换计算各极化状态的平地相位；
B) 去除各极化状态的平地相位对极化干涉数据的影响；
C) 对极化干涉数据平滑滤波；
D) 利用特征分解提取最优投影矢量对；
E) 计算最优干涉相干系数与干涉相位，检验最优估计的有效性；F)将平地相位补偿至最优干涉相位中，按干涉相干系数大小顺序排 列最优千涉相干系数以及最优干涉相位。
所述的方法，其所述A)步先利用快速傅里叶变换计算各极化状态的 平地相位，其各极化状态的平地相位计算包含如下步骤
Al、计算指定范围内^^^2，各极化状态下的干涉相位&、 &、 fe、
<formula>formula see original document page 12</formula>
其中,,表示主天线接收到的复图像数据；^表示副天线接收到的复图像数 据；O",n)指示像素在图像中的行、列位置；m,、 ^表示用于计算平地相位 的数据所在行的起始位置、终止位置；xx、 、 m、 n表示复图像数据的极 化方式，符号.表示复数乘运算，'表示共轭运算，^表示复数的辐角主值。 A2、利用快速傅立叶变换计算指定范围内m,"^2，各极化状态的复 干涉相位序列的频谱lw、义 "、 Am &，
Lrr =尸
其中^(加,..)表示由&—,1) 、 &(m,2)、 ... (^(m，W)构成的序列；A3、计算指定范围内^^^2各极化状态复相位序列频谱能量的平均
值<formula>formula see original document page 13</formula>A4、确定各极化状态复相位序列频谱能量的平均值的最大值位置 A5、计算各极化状态平地相位序列仏w、 6。,"、 ,"
<formula>formula see original document page 13</formula>其中W表示快速傅里叶变换长度。
所述的方法，其所述B)步去除各极化状态的平地相位对极化干涉数
据的影响，是去除平地相位的影响可减小相邻位置干涉相位梯度，增强相
邻位置的极化干涉数据的一致性，具体以下式计算<formula>formula see original document page 14</formula>所述的方法，其所述C)步极化干涉数据平滑滤波，包含以下步骤 Cl、以待滤波位置为窗口中心，大小为^Xw;的滤波窗口内，生成包 含待滤波位置的，沿行、列以及两对角线方向的，中心位置不同的滤波子
C2、计算各个子窗口内极化干涉数据总能量的方差^—,n):<formula>formula see original document page 14</formula>
其中^—， )表示以— )中心的第L个子窗口内极化干涉数据总能量的方
差；符号ll表示复数的模；符号〈〉(一以—,n)中心的第L个子窗口内极化干
涉数据的算术平均；
C3、以极化千涉数据总能量方差最小为原则确定滤波子窗口；； C4、计算选定的滤波子窗口内极化干涉数据的算术平均T,一,n)、<formula>formula see original document page 14</formula>所述的方法，其所述D)步利用特征分解提取最优投影矢量对，是基于 对矩阵^(，)的分析，提取三组最优投影矢量对R,W、 R，W、 K,d5j，
具体步骤为
Dl、 计算矢巨阵￡1,2(/ ,")￡1「2—,")禾口 flf2( i,")￡l12(m," )的单位特征矢量
ft「2(m，")lll2 (m，w) = J]A,(i52 32 w
其中，4、 a2、 4是矩阵n,2(m,")ti「2("!,")和si「2(/n,")n,2("!,")的三个特征矢量；A、
^分别是矩阵il,2—,")0,")和0，")il,2(m,")对应于特征值义,的单位化特征矢
D2、按照A4^顺序排列单位特征矢量A、 52,，生成三组最优投影 所述的方法，其所述E)步计算最优干涉相干系数与干涉相位，检验
最优估计的有效性，包含以下步骤
El 、利用最优投影矢量对计算最优干涉相干系数"，)与干涉相位
V(《Tn(^nK,)(《T22(m，")^)
E2、检验干涉相干系数值的合理性，设置错误估计标记 当干涉相干系数r,(琴)大于1或(《1,|(罕)(5|,)(《1;2(，)&) = 0时，最优干涉
15相干系数以及干涉相位无效，并设置错误估计标记y,(M,"卜Q，
所述的方法，其所述F)步将平地相位补偿至最优干涉相位中，是按 干涉相干系数大小顺序排列最优干涉相干系数以及最优干涉相位，具体包 含如下步骤
Fl、将平地相位仏"")补偿至最优干涉相位f,—, )中，获得干涉相位
F2、按照+,,)"2(一)^>,力的顺序，排列生成三组最优干涉相干系数
以及干涉相位结果^(m,nU(M,4 、 {,2(/ ,")} 、 {r3(M,")A(m, )};其中第一最
优结果的干涉相干系数y,0^)最高，明显优于优化前的XX、 XY、 YX、 YY 极化的干涉相干系数。
本发明的方法计算速度快、相干性改善程度大、适用于对不同类型的 极化干涉数据的优化相干系数。
图1本发明中所使用的极化干涉数据示意图2本发明方法的处理流程图3a设定5x5窗口时统一滤波窗口示意图3b设定5x5窗口时Lee自适应滤波时采用的滤波子窗口示意图； 图3c设定5x5窗口时本发明中在所选用的滤波子窗口示意图； 图4a SIR-C系统L波段数据去除平地相位影响前HH极化干涉相位量 化图；图4b SIR-C系统L波段数据去除平地相位影响前HV极化干涉相位量化图4c SIR-C系统L波段数据去除平地相位影响前VH极化干涉相位量化图4d SIR-C系统L波段数据去除平地相位影响前VV极化干涉相位量化图4e图4a-图4d中干涉相位的量化方式示意图；图5a SIR-C系统L波段数据计算的HH极化平地相位量化图；图5b SIR-C系统L波段数据计算的HV极化平地相位量化图；图5c SIR-C系统L波段数据计算的VH极化平地相位量化图；图5d SIR-C系统L波段数据计算的VV极化平地相位量化图；图5e图5a-图5d中干涉相位的量化方式示意图；图6a SIR-C系统L波段数据去除平地相位影响后HH极化干涉相位量化图6b SIR-C系统L波段数据去除平地相位影响后HV极化干涉相位量化图6c SIR-C系统L波段数据去除平地相位影响后VH极化干涉相位量化图6d SIR-C系统L波段数据去除平地相位影响后VV极化干涉相位量化图6e图6a-图6d中干涉相位的量化方式示意图7a SIR-C系统L波段数据经本发明方法处理获得第一最优干涉相
17干系数的量化图7b SIR-C系统L波段数据经本发明方法处理获得第二最优干涉相干系数的量化图7c SIR-C系统L波段数据经本发明方法处理获得第三最优干涉相
干系数的量化图7d图7a-图7c中干涉相干系数的量化方式示意图8a SIR-C系统L波段数据经本发明方法处理获得第一最优干涉相
干系数的量化图8b SIR-C系统L波段数据经本发明方法处理获得第二最优干涉相干系数的量化图8c SIR-C系统L波段数据经本发明方法处理获得第三最优干涉相
干系数的量化图8d图8a-图8c中干涉相干系数的量化方式示意图9a SIR-C系统L波段数据用5 X 5滤波窗口计算的HH极化干涉相
千系数的量化图9b SIR-C系统L波段数据采用5X5滤波窗口计算的HV极化干涉相干系数的量化图9c SIR-C系统L波段数据采用5 X 5滤波窗口计算的VH极化干涉相干系数的量化图9d SIR-C系统L波段数据采用5X5滤波窗口计算的VV极化干涉相干系数的量化图9e图9a-图9b中干涉相干系数的量化方式示意图；图9d图9a-图9b最优干涉相位量化方式示意图；图10 SIR-C系统L波段数据的经本方法优化干涉相干系数与未经本方法优化的各极化干涉相干系数的对比图。
具体实施例方式
如图2所示，基于极化干涉合成孔径雷达数据优化干涉相干系数的方法包括步骤如下
A) 利用快速傅里叶变换计算各极化状态的平地相位；
B) 去除各极化状态的平地相位对极化干涉数据的影响；
C) 对极化干涉数据平滑滤波；
D) 利用特征分解提取最优投影矢量对；
E) 计算最优干涉相干系数与干涉相位，检验最优估计的有效性；
F) 将平地相位补偿至最优干涉相位中，按干涉相干系数大小顺序排列最优干涉相干系数以及最优干涉相位。
步骤A)利用快速傅里叶变换计算各极化状态的平地相位，是为去除平地相位造成的相邻位置上数据的波动。考虑到不同极化状态下干涉相位的差异，因此分别基于各极化下干涉相位估计平地相位；为提高计算效率，仅选用干涉相位变化规律的部分干涉相位数据计算即可。具体步骤如下
Al、计算指定范围内^^^2，各极化状态下的干涉相位fe、 &、 &、
At(附，")=々g(』('"，")■ AlT (W，"))
其中A表示主天线接收到的复图像数据；^表示副天线接收到的复图像数据；O",力指示像素在图像中的行、列位置；m,、 ^表示用于计算平地相位的数据所在行的起始位置、终止位置；^、 w、 ^、 〃表示复图像数据的极化方式，符号.表示复数乘运算，'表示共轭运算，^g表示复数的辐角主值。A2、利用快速傅立叶变换计算指定范围内"^^2，各极化状态的复干涉相位序列的频谱z^、 x, AT、 x . 、
其中!^卜,..)表示由&—,1)、 ,2)、 ... (^0",A0构成的序列。
A3、计算指定范围内""^2各极化状态复相位序列频谱能量的平均
值|4|、 |&|、 |&|、 |iYY|:
w2 — w, 1 'm2 —
1
肌2 — 1
S卜:.
20A4、确定各极化状态复相位序列频谱能量的平均值的最大值位置
A5、计算各极化状态平地相位序列仏w、
W灯(")=
2;r(w --i)
-i)
.0
2;r(n -i)
其中A表示快速傅里叶变换长度。
步骤B)去除各极化状态的平地相位对极化干涉数据的影响，是为减小 相邻位置干涉相位梯度，增强相邻位置的极化干涉数据的一致性。为了避 免平地相位去除对主天线极化数据的影响，采用去除平地相位对副天线的 全极化数据的处理方法。由于不同极化下平地相位存在差异，因而需对各
极化状态的数据分别处理。具体步骤如下
步骤C)极化干涉数据平滑滤波，是为抑制噪声影响。若采用如图3所示的统一矩形窗口滤波，能造成非同质、不均匀数据之间的平滑。为保 证用于平滑滤波的极化干涉数据是均匀的，因而需根据极化干涉数据的波 动性自适应选择具有方向性的如图3b所示的滤波子窗口。考虑到相邻位 置数据的一致性高的特点，使得滤波子窗口能包含尽量多的与滤波位置临 近的数据，本方法通过平移图3b所示的子窗口，提供了更丰富的滤波子 窗口。由于在极化千涉数据滤波结果中补偿平地相位比较繁琐，因而暂不 补偿平地相位。极化干涉数据平滑滤波包含步骤如下
Cl、以待滤波位置为窗口中心大小为1^X^滤波窗口内，生成包含待
滤波位置的，沿行、列以及两对角线方向的，中心位置不同的如图3b所
示的滤波子窗口；
C2、计算各个子窗口内极化干涉数据总能量的方差^ —n): 0"2,.—,") = 〈|、 「 +h |2 +h 「 +k'2;a|2 +>'2 |2 +|s;「 +k'2Kr「〉，、，
-〈(l'lXvl "H'ml "H"^Wl +|"^2^>￡"| +|^2"1 +|^2灯| +|lf2>x| )〉
其中一,.—，力表示以—,力中心的第L个子窗口内极化干涉数据总能量的方 差；符号M表示复数的模；符号〈〉(—以—,")中心的第L个子窗口内极化干
涉数据的算术平均。
C3、以极化干涉数据总能量方差最小为原则确定滤波子窗口 ；。 C4、计算选定的滤波子窗口内极化干涉数据的算术平均1,1—,")、
22T22 (m，w)
A议
(m,H)，iF
步骤D)利用特征分解提取最优投影矢量对。本方法基于包含极化、 干涉信息的矩阵&2，利用特征分解即可提取最优投影矢量对KA,J、 ^Aj、 一^J。本方法不要求矩阵1,,、 T。可逆，可用于各类极化干涉数据
的处。利用特征分解提取最优投影矢量对的具体步骤为
Dl、计算矩阵 fll2(w，")￡^(w,")禾口 Q「2(m，")ft,2(m，")的单"f立牛争征矢量
ft^(m，")￡212(m,") = Z^.^2 62 w
其中，V A2、 A2是矩阵il,2(M,")il「2(，")和i2「2(m,")fl,2(;n，")的三个特征矢量；<5,,、 ^分别是矩阵il,2(m,")J^(M,")和O,n)iM一)对应于特征值;i,的单位化特征矢
D2、按照W^^顺序排列单位特征矢量A、 A,，生成三组最优投影 矢量对{^,,<52,} 、 {^,(52」、(A，^,卜
步骤E)计算最优干涉相干系数与干涉相位，检验最优估计的有效性。
已有的优化干涉相干系数的方法在矩阵7|,、 122不可逆时，无法获得正确估
23计，却并未对优化结果进行有效性检验。本方法添加了对最优估计的有效
性检验，其具体步骤如下
El、利用最优投影矢量对计算最优干涉相干系数y,如,")与干涉相位
E2、检验干涉相干系数值的合理性，设置错误估计标记 当干涉相干系数/,如，n)大于1或(《1|1(，)51,,2%(^)^) = 0时，最优千涉
相干系数以及干涉相位无效，并设置错误估计标记,,如,")=0， r如，力-o。
步骤F)将平地相位补偿至最优干涉相位中，按干涉相干系数大小顺 序排列最优干涉相干系数以及最优干涉相位。为避免将平地相位补偿至极
化干涉数据滤波结果的复杂过程，将平地相位补偿至最优干涉相位中；该
平地相位补偿方法可保证最优干涉相位结果的连续性。具体包含如下步
骤
Fl、将平地相位 。,^W补偿至最优干涉相位f,(—中，获得干涉相位 F2、按照+,力"2—的顺序，排列生成三组最优干涉相干系数以及
干涉相k结果(^(m，")，A(m,")) 、 {y2(w，/i),A(m，")} 、 {;K3(m,rt)，A(m>")}。
下面结合附图详细说明本发明技术方案中所涉及的各个细节问题。应指出的是，所描述的实施例仅旨在便于对本发明的理解，而对其不起任何 限定作用。
首先介绍一下极化干涉数据及其表示方法
如图1所示，全极化合成孔径雷达数据是由采用某两种正交极化状态
X、 Y发送、接收雷达信号的天线系统，采集的到的XX、 XY、 YX、 YY
四种极化状态下的雷达信号分别经过成像处理获得的这四种极化状态下 的复数据矩阵&、 &、 &、；。而单基线极化干涉合成孔径雷达数据是 由不同时间采集的或者同时在不同位置采集的两组全极化合成径雷达数 据构成，且每组的各个极化状态下的雷达信号要分别经过成像处理后，两 组全极化合成孔径雷达数据还需要经过极化定标处理、图像配准和重采样 处理才能获得两组全极化复数据矩阵sw、 s^、 s,M、 s，与s^、 & 、 & 、 s, 。其中，标号l表示参考图像数据，本发明中叫做主天线全极化数据，
标号2表示非参考图像数据，本发明中叫做副天线全极化数据。
本发明实例中使用的数据是美国喷气推进实验室的SIR-C系统于 1994年10月7日、8日采集自中国天山山麓的L波段全极化数据构成的 单基线极化干涉数据，这两组全极化数据的编号为122.20和138.20。图像 尺寸为6000X900 (行方位位置，歹U:斜距位置，单位像素)。SIR-C 的全极化系统采用H、 V极化状态。
如图2所示，基于极化干涉合成孔径雷达数据优化干涉相干系数的方
法包含如下步骤
A)利用快速傅里叶变换计算各极化状态的平地相位；B) 去除各极化状态的平地相位对极化干涉数据的影响；
C) 对极化干涉数据平滑滤波；
D) 利用特征分解提取最优投影矢量对；
E) 计算最优干涉相干系数与干涉相位，检验最优估计的有效性；
F) 将平地相位补偿至最优干涉相位中，按干涉相干系数大小顺序排 列最优干涉相干系数以及最优干涉相位。
下面介绍各个步骤具体的实施过程
Al、利用行标在3601-4600范围内的数据，计算HH、 HV、 VH、 VV 极化下干涉相位
^ (m,w)=(m，") /2fT (w，w))
将各位置的干涉相位&(^)、 ^(m,力、(M^)分别按照如图4d所
示的方式量化为灰度值，生成如图4a、图4b、图4c、图4d所示未去去平 地相位前HH、 HV、 VH、 VV极化干涉相位量化。观察图4a-图4d发现， 干涉相位梯度较大、干涉条纹密集。
A2、将行标在3601-4600范围内的数据补零为长度1024的复数序列，
利用快速傅立叶变换计算频谱
26<formula>formula see original document page 27</formula>
其中l""表示行标号为w的复干涉相位序列快速傅立叶变换结果的第K 个位置的值。
A3、计算指定范围内^^"2各极化状态复相位序列频谱能量的平均
<formula>formula see original document page 27</formula>
A4、确定各极化状态复相位序列频谱能量的平均值的最大值位置
<formula>formula see original document page 27</formula>
通常情况下，不同方位位置的干涉条纹距离向频谱除了受到平地相位的影 响还受到地形的影响。为了能准确估计平地相位，本发明方法采用一定范 围内的干涉相位距离向平均频谱估计平地相位。由于各个范围位置的地形变化不同，因而在平均过程中将地形变化引起的频谱当作噪声抑制了，而 强化了该方位范围内共同的平地相位引起的频谱，因而平均频谱获得的平 地相位估计更准确。
A4、由于平均频谱主要受到平地相位的影响，因此可直接利用能量峰 值位置计算出平地相位
-1)
512
85""-1)
512
85;r("-1)
512
85+-1)
512(5)
"=1,2,...,卿
其中(V^W表示XX极化状态下第n个距离向像素位置的平地相位。其中 平地相位表示当观测区域为平坦地面时干涉相位，平地相位随着距离向位 置而变化。图5a、图5b、图5c、图5d是HH、 HV、 VH、 VV极化的平 地相位结果的量化图。
步骤B:将副天线全极化图像中各个相位位置的复数据乘以对应平地 相位
2册(m,")"2冊(附,")e '"")
为了对比说明去除平地相位前后各个极化状态下干涉相位的情况，利 用去除平地相位后的副天线全极化数据与主天线全极化干涉数据计算HH、 HV、 VH、 VV极化状态下的干涉相位
(*'朋("',")=Ag(J, v,細
图6a、图6b、图6c、图6d是去平地相位影响后的HH、 HV、 VH、 VV
极化的干涉相位量化图。与图4a-图4d相比较，平地相位去处后干涉相位 梯度明显减小，而且相邻像素之间干涉相位的差异也明显减小，从而保证 了后续滤波窗口内极化干涉数据的一致性。 步骤C的具体实施过程为
Cl、以待滤波位置为窗口中心大小为5x5滤波窗口内，生成包含待滤 波位置的，沿行、列以及两对角线方向的，中心位置不同的如图3c所示 的滤波子窗口；
C2、计算各个子窗口内极化干涉数据总能量的方差AO^):
a2丄(附,")-〈l"W卩+K^「+l^r+kwT+kwl2+k'2Ay「+ku7l2+k;l2〉，、，
-〈(一ml +klJlT| "H^JItI +^2议| 十l^2^yl +11 21^| )〉
其中<0^)表示以—《)中心的第L个子窗口内极化干涉数据总能量的方 差；符号ll表示复数的模；符号〈〉(一以On,n)中心的第L个子窗口内极化干
涉数据的算术平均。
C3、以极化干涉数据总能量方差最小为原则确定滤波子窗口^。 C4、计算选定的滤波子窗口内极化干涉数据的算术平均^0^)、该处理的目的是在适合方向的子窗口内进行平滑滤波，尽量保证滤波区域 内的极化干涉数据是同质均匀的。该自适应滤波处理在抑制噪声影响的同 时，尽可能的保留各个极化的幅度和相对相位关系，从而保证后续提取到 的主要散射中心信息的准确性。
步骤D的具体实施过程为
Dl、计算矩阵 )和 )的单位特征矢量
其中，A、 A2、 A2是矩阵n,2(m,")i1「2—,")和Q「2(/n,")&2(/n,")的三个特征矢量；d5,、 ^分别是矩阵iU^)0,力和^>，一12(^)对应于特征值2,的单位化特征矢
D2、按照PA^顺序排列单位特征矢量A、 A,，生成三组最优投影
(10)
步骤E的具体实施过程为:El、利用最优投影矢量对计算最优干涉相干系数y,—,")与干涉相位
f,.(m，")
y.(附,w) = ~ - ■ ■ ',A/ T',(m,n)A)(《T22—,"
(11)
E2、检验干涉相干系数值的合理性，设置错误估计标记
当干涉相干系数r,—,")大于1或(《i;,(m,"H)(《T22("wH,)"时，最优干涉
相干系数以及干涉相位无效，并设置错误估计标记K—，力-Q， O,n) = 0。 步骤F的具体实施过程为
Fl、将平地相位仏."")补偿至最优干涉相位^,—，》)中，获得干涉相位
F2、按照,,—)^2(—)的顺序，排列生成三组最优干涉相干系数 以及干涉相位结果&,0 ，") ,")} 、 、 &3(w,")A(/n,")}。
F3、重复执行步骤C-步骤F2，求得各个位置的最优干涉相干系数以 及干涉相位，并构成最终优化结果
、,(1,1) ... y,(l,W)—
..."风W)」 (13) 、(U) ... )]...
图7a、图7b、图7c分别是第一、第二、第三最优干涉相干系数r,、 r2、 r3 的量化图。图8a、图8b、图8c分别是第一、第二、第三最优干涉相位ct),、 *2、 。3的量化图。本发明的一个实例，是用美国喷气推进实验室的SIR-C系统重复飞行 采集的中国天山山麓的L波段单基线极化干涉数据的优化干涉相干系数 的结果。
图4a-图4d是步骤Al中利用未去除平地相位影响的行标号在3601 -4000的极化干涉数据计算的HH、 HV、 VH、 VV极化下干涉相位量化图。 图4e指明了生成图4a-图4d的量化方式。从图4a-图4d中可发现，干涉 相位在列方向形成明暗相间的条纹，且整体的干涉条纹十分密集，这说明 干涉相位在距离向上变化很快、相位梯度大，有必要去除平地相位的影响。 如图4a-图4b所示，干涉相位形成的条纹十分规律，适用于计算平地相位。
本实例采用方位向标号从3601-4600的干涉相位的距离向频谱估计平 地相位。图5a-图5d为计算出的HH、 HV、 VH、 VV极化的平地干涉相位 主值的量化图。图5e指明了生成图5a-图5d的量化方式。从图5a-图5b 可发现，平地相位在列方向形成平行的明暗相间的条纹。
本实例数据经去除平地相位处理后，行标号从3601-4600数据计算 HH、 HV、 VH、 VV极化下的干涉相位。图6a-图6d为计算出的HH、 HV、 VH、 VV极化的平地干涉相位主值的量化图。图6e指明了生成图6a-图 6d的量化方式。与图4a-图4d比较发现，平地相位去除后干涉相位梯度明 显降低，相邻像素之间的干涉相位差异明显减小。
本实例数据经本发明方法处理，最终获得了三组最优干涉相干系数以 及最优干涉相位结果。图7a、图7b、图7c分别是经过本方明的优化相干系数方法处理得到的第一、第二、第三最优干涉相干系数结果量化图；图
7d指示了图7a-图7c量化方式。图8a、图8b、图8c分别是利用本方明优 化干涉相千系数的方法计算得到的第一、第二、第三最优干涉相位的量化 图。图8d指示了图8a-图8c的量化方式。从图8a-图8c发现，三组最优 的干涉相位具有较好的连续性。
为对比应用本发明方法后干涉相干系数的改善程度，计算5X5矩形 窗口内HH、 HV、 VH、 VV极化的干涉相干系数^—,"、 ^vOn,"、 ^—, )、
V(m,")=
rw(m,")=
'朋—:)l卜冊(m,")'/2冊(m,")1
F(m,M)'
如r—,")./,肿(—〖l卜2抓(附，")./服
w(附，")l;卜卵(m,力'/鹿(m，")1
(仍，rt) ■ "^ikk (附，")*^ (附，")■ Z2rr (W，")
(14)
图9a、图9b、图9c、图9d分别是HH、 HV、 VV极化状态下的干涉相干 系数均匀量化图。图9e是图9a-图9d中干涉相干系数的量化方式。从图 9a-图9d发现，HH、 W极化下的干涉相干系数图比HV、 VH极化下干涉 相干系数图亮，这说明HH、 VV极化下干涉相干系数整体优于HV、 VH 极化下干涉相干系数，并且还证实了干涉相干系数对极化状态的依赖性以 及利用极化状态优化干涉相干系数的可行性。将图9a-图9d与图7a-图7c 对比后发现，第一最优干涉相干系数图明显亮于HH、 W极化的干涉相干系数图，这证实了本发明方法确实可以明显优化干涉相干系数。
为了更加清晰的描述本发明提出相干性优化方法提高干涉相干性的
能力，图10展示了 HH、 HV、 VV的干涉相干系数与三组最优干涉相干系 数的统计直方图。具体的统计方法为将
区间划分为均匀50小区间，
统计图像中干涉相干系数落入各个小区间的频数，从而计算各区间内干涉
相干系数的出现频率，从而绘制出统计直方图。从图10发现，第一最优 干涉相干系数总体分布更接近于1，其总体干涉相干性明显优于HH、 VV 的结果；验证了本发明方法的有效性。
通过对本实例结果的分析可证实，本发明提出的基于单基线极化干涉 数据优化干涉相干系数的方法确实可以明显提高干涉相干系数、改善干涉 相位质量。因而，本发明方法在极化干涉的分类、参数反演的数据处理中 有巨大的应用潜力。
权利要求
1、基于极化干涉合成孔径雷达数据优化干涉相干系数的方法，其特征在于包括步骤如下A)利用快速傅里叶变换计算各极化状态的平地相位；B)去除各极化状态的平地相位对极化干涉数据的影响；C)对极化干涉数据平滑滤波；D)利用特征分解提取最优投影矢量对；E)计算最优干涉相干系数与干涉相位，检验最优估计的有效性；F)将平地相位补偿至最优干涉相位中，按干涉相干系数大小顺序排列最优干涉相干系数以及最优干涉相位。
2、如权利要求1所述的方法，其特征在于，所述A)步先利用快速 傅里叶变换计算各极化状态的平地相位，其各极化状态的平地相位计算包 含如下步骤Al、计算指定范围内",,^"2，各极化状态下的干涉相位&、 &、 fe、<formula>formula see original document page 2</formula>其中A表示主天线接收到的复图像数据；^表示副天线接收到的复图像数 据；—一指示像素在图像中的行、列位置；^、 ^表示用于计算平地相位的数据所在行的起始位置、终止位置；w、 w、 ^、 w表示复图像数据的极 化方式，符号.表示复数乘运算，'表示共轭运算，^表示复数的辐角主值。A2、利用快速傅立叶变换计算指定范围内m,^^,，各极化状态的复干涉相位序列的频谱义^、 & 、^ 其中&—，..)表示由&(^1)、 &—,2)、…fc一,W)构成的序列；A3、计算指定范围内《^ ^2各极化状态复相位序列频谱能量的平均 值M、 l&l、 |^卜附2 -/n1W2 —1m2 —1"=m2,A4、确定各极化状态复相位序列频谱能量的平均值的最大值位置 A5、计算各极化状态平地相位序列^。^、6。m ("): 0麵 w2-i)02;r(w -…《f" 4。其中K表示快速傅里叶变换长度。
3、 如权利要求l所述的方法，其特征在于所述B)步去除各极化 状态的平地相位对极化干涉数据的影响，是去除平地相位的影响可减小相 邻位置干涉相位梯度，增强相邻位置的极化干涉数据的一致性，具体以下 式计算;(仍,")=j2JX (m, "y"(") y; (m,") = ^ (仿，一'。'(")
4、 如权利要求1所述的方法，其特征在于所述C)步极化干涉数 据平滑滤波，包含以下步骤Cl、以待滤波位置为窗口中心，大小为^X^的滤波窗口内，生成包含待滤波位置的，沿行、列以及两对角线方向的，中心位置不同的滤波子 窗口；C2、计算各个子窗口内极化干涉数据总能量的方差^—力°" - 〈(>si;ty| +|A;fyl +|isi>t1 +l52;cf| "Hs'2at1 +|s'2;\t| )、 / (/",")，丄-〈(—ml十kml 十^2Xvl +|J'2Jkr| +|;2灯| )〉其中A—,力表示以—,")中心的第L个子窗口内极化干涉数据总能量的方差；符号ll表示复数的模；符号〈〉一以(m,")中心的第L个子窗口内极化干涉数据的算术平均；C3、以极化干涉数据总能量方差最小为原则确定滤波子窗口；； C4、计算选定的滤波子窗口内极化干涉数据的算术平均1;,—，。、
5、如权利要求l所述的方法，其特征在于所述D)步利用特征分解 提取最优投影矢量对，是基于对矩阵^On,力的分析，提取三组最优投影矢具体步骤为Dl、计算矩阵 )的单位特征矢量其中，4、 A2、义2是矩阵Q,2(,",")ftf2(m,")和ft「2(m,")S1,2—,")的三个特征矢量；<5,,、&分别是矩阵￡1,2—,力^2—, )和，")iM,)对应于特征值义,的单位化特征矢D2、按照;^44顺序排列单位特征矢量* 、 a，生成三组最优投影
6、如权利要求1所述的方法，其特征在于所述E)步计算最优干 涉相干系数与干涉相位，检验最优估计的有效性，包含以下步骤El、利用最优投影矢量对计算最优干涉相干系数"吣)与干涉相位<formula>formula see original document page 6</formula>E2、检验干涉相干系数值的合理性，设置错误估计标记 当干涉相干系数n(w，")大于1或051%(一)<51,)(<52%( ,+2,) = 0时，最优干涉相干系数以及干涉相位无效，并设置错误估计标记^^;) = ()， O,n) = 0。
7、如权利要求1所述的方法，其特征在于所述F)步将平地相位补偿至最优干涉相位中，是按干涉相干系数大小顺序排列最优干涉相干系数以及最优干涉相位，具体包含如下步骤Fl、将平地相位^。,^(n)补偿至最优干涉相位0,n)中，获得干涉相位<formula>formula see original document page 6</formula>F2、按照^(，)^(")">, )的顺序，排列生成三组最优干涉相干系数 以及干涉相位结果！y,(/n,";uK"" 、 b("i,";u(m,"》、{^(/ ,"))3(^")};其中第一最优结果的干涉相干系数;^z，")最高，明显优于优化前的XX、 XY、 YX、 YY 极化的干涉相干系数。
全文摘要
本发明提供一种基于极化干涉合成孔径雷达数据优化干涉相干系数的方法，涉及合成孔径雷达技术，该方法先利用快速傅里叶变换计算各极化状态的平地相位；去除各极化状态的平地相位对极化干涉数据的影响；对极化干涉数据平滑滤波；利用特征分解提取最优投影矢量对；计算最优干涉相干系数与干涉相位，检验最优估计的有效性；将平地相位补偿至最优干涉相位中，按干涉相干系数大小顺序排列最优干涉相干系数以及最优干涉相位。本发明提供方法不受极化干涉数据矩阵可逆的条件限制，无需补偿最优投影矢量间的绝对相位差，并可检验优化结果的有效性。本发明的方法原理简单、计算方便，可用于极化干涉分类、参数反演等数据处理中。
文档编号G01S7/48GK101685155SQ200810223609
公开日2010年3月31日 申请日期2008年9月27日 优先权日2008年9月27日
发明者芳 曹, 文 洪, 璐 白 申请人:中国科学院电子学研究所