结合局部频率估计的小波域InSAR干涉相位滤波方法
【专利摘要】本发明提供了一种结合局部频率估计的小波域InSAR干涉相位滤波方法。该小波域InSAR干涉相位滤波方法利用局部频率估计实现了对复干涉相位的小波系数中有用信息子带和噪声子带的区分，利用通用阈值收缩和邻域阈值收缩两种方法分别具有去噪效果好和细节保持能力强的特点，对有用信息所在子带的小波系数进行邻域阈值收缩，而对其它子带的小波系数则进行通用阈值收缩，从而尽可能的滤除噪声，同时保持干涉条纹的细节信息不被破坏，实现高精度的干涉相位滤波，为高精度的干涉测量提供了条件。
【专利说明】结合局部频率估计的小波域InSAR干涉相位滤波方法
【技术领域】
[0001] 本发明涉及电子信息技术雷达【技术领域】，尤其涉及一种结合局部频率估计的小波 域InSAR干涉相位滤波方法。
【背景技术】
[0002] 干涉合成孔径雷达（Interferometric Synthetic Aperture Radar，InSAR)是利 用合成孔径雷达（SAR)两个通道的干涉相位信息提取地表的高程信息或变化信息，将SAR 的测量拓展到三维空间，具有全天时、全天候、高精度的特点，因此在地形测绘、冰川研究、 海洋测绘以及地面沉降监测等多个领域都有广泛的应用。
[0003] 干涉测量的精度和可靠性在很大程度上取决于干涉相位图的质量。然而，在实际 系统中，受热噪声去相干、时间去相干、基线去相干、配准误差等多种去相干因素的影响，干 涉相位图不可避免的存在相位噪声。低质量的干涉相位将会影响后续的干涉相位解缠及高 程反演的准确性。因此，在相位解缠前必须对干涉相位进行滤波，从而获取高质量的干涉相 位图。
[0004] 目前干涉相位的滤波方法可以大致分为空间域滤波和变换域滤波两类。圆周期均 值滤波或中值滤波是一种最基本的空间域滤波方法（见参考文献1)，它实现简单，但滤波 窗口大小不好确定，在条纹密集时容易破坏相位细节，降低分辨率。变换域滤波方法在实际 中应用更为广泛，如Goldstein滤波（见参考文献2)，然而，该方法受分块大小和滤波参数 的影响较大，在信噪比很低时，滤波效果较差；小波变换由于其良好的时频分析特性和多分 辨率特性，也可以应用于干涉相位滤波中（见参考文献3)，该方法通过增大小波系数中的 信号成分来实现滤波，能够很好地保持干涉条纹的细节信息，并在一定程度上提高了图像 的信噪比，但由于对噪声信息没有进行有效的抑制，使得去噪效果较差。因此，为满足InSAR 应用对干涉相位精度的要求，有必要进一步研究能有效去除噪声并保持相位细节的变换域 滤波方法。
[0005] 参考文献：
[0006] [1]R. Lanari. Generation of digital elevation models by using SIR-C/X-SAR multifrequency two-pass interferometry ：The Etna case study. IEEE Transactions on Geoscience and Remote Sensing,1996,34(5) :1097-1114.
[0007] [2]R. M. Goldstein, C. L. Werner. Radar Interferogram filtering for Geophysical Application. Geophysical Research Letters. 1998,25(21) :4035-4038.
[0008] [3] Lopez-Martinez C, Fabregas X. Modeling and reduction of SAR interferometric phase noise in the wavelet domain[J]. IEEE Trans, on Geoscience and Remote Sensing,2002,40(12) :2553-2566
【发明内容】
[0009] (一）要解决的技术问题
[0010] 为解决上述问题，本发明提供了一种结合局部频率估计的小波域InSAR干涉相位 滤波的方法，以克服现有的小波域干涉相位滤波方法不能兼顾去噪和细节保持的缺点，从 而提高干涉相位的精度。
[0011](二）技术方案
[0012] 根据本发明的一个方面，提供了一种结合局部频率估计的小波域InSAR干涉相位 滤波方法。该小波域InSAR干涉相位滤波方法包括：步骤A :将InSAR的干涉相位Φ变换 到复数域&_Φ，分别取复数域干涉相位&_φ的实部和虚部；步骤B :对干涉相位Φ进行局部 频率估计，得到干涉相位所在的频率范围；步骤C :对复数域干涉相位的实部和虚部，分别 进行尺度为s的小波分解，得到不同子带的小波系数的频率范围，其中，m， n为小波系数 的位置，i为分解的尺度，其范围为1?s ;步骤D :对复数域干涉相位的实部和虚部，分别根 据干涉相位所在的频率范围和不同子带小波系数的频率范围，分别确定有用信息和噪声所 在的子带；步骤E :对复数域干涉相位的实部和虚部，分别对噪声所在子带的小波系数分别 进行通用阈值收缩处理；步骤F :对复数域干涉相位的实部和虚部，分别对有用信息所在子 带的小波系数分别进行邻域阈值收缩处理；步骤G :对复数域干涉相位的实部和虚部，分别 将通用阈值收缩处理处理后的噪声所在子带的小波系数和邻域阈值收缩处理后的有用信 息所在子带的小波系数共同进行小波重构，得到滤波后的复数域干涉相位的实部和虚部； 以及步骤H :由滤波后的复数域干涉相位的实部和虚部得到InSAR滤波后的干涉相位。
[0013] (三）有益效果
[0014] 本发明利用局部频率估计实现了对复干涉相位的小波系数中有用信息子带和噪 声子带的区分，利用通用阈值收缩和邻域阈值收缩两种方法分别具有去噪效果好和细节保 持能力强的特点，对有用信息所在子带的小波系数进行邻域阈值收缩，而对其它子带的小 波系数则进行通用阈值收缩，从而尽可能的滤除噪声，同时保持干涉条纹的细节信息不被 破坏，实现高精度的干涉相位滤波，为高精度的干涉测量提供了条件。
【专利附图】
【附图说明】
[0015] 图1为根据本发明实施例结合局部频率估计的小波域InSAR干涉相位滤波方法的 流程图；
[0016] 图2为意大利Etna火山实测的星载InSAR干涉相位图；
[0017] 图3A?图3D为分别利用圆周期均值滤波,Goldstein滤波，Lopez-Martinez C提 出的小波滤波以及本实施例方法滤波后的干涉相位图
【具体实施方式】
[0018] 为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照 附图，对本发明进一步详细说明。需要说明的是，在附图或说明书描述中，相似或相同的部 分都使用相同的图号。附图中未绘示或描述的实现方式，为所属【技术领域】中普通技术人员 所知的形式。另外，虽然本文可提供包含特定值的参数的示范，但应了解，参数无需确切等 于相应的值，而是可在可接受的误差容限或设计约束内近似于相应的值。
[0019] 本发明将将干涉相位变换到小波域后，利用局部频率估计判断出小波系数中包含 有用信息的子带，对有用信息所在子带的小波系数利用邻域阈值进行收缩处理，而对其他 子带的小波系数利用通用阈值进行收缩处理，从而尽可能的滤除噪声，同时不损失干涉条 纹的细节信息，从而为InSAR干涉相位滤波提供了一种有效的方法。
[0020] 根据本发明的一个方面，提供了一种结合局部频率估计的小波域InSAR干涉相位 滤波方法。图1为根据本发明实施例结合局部频率估计的小波域InSAR干涉相位滤波方法 的流程图。如图1所示，本实施例结合局部频率估计的小波域InSAR干涉相位滤波方法包 括：
[0021] 步骤A :将InSAR预处理后得到的干涉相位Φ变换到复数域e# = cos<i)+jsin<i) 中，分别取复数域干涉相位6']<*>的实部1^{6' ]<*>}=〇08<1)，虚部11]1{6']<* >}=8；[11小；
[0022] 此处的预处理包括：双通道InSAR数据的成像处理，配准处理、复数据的共轭相乘 等等。这些均为本领域内的公知常识，此处不再详细说明。
[0023] 步骤B :对复数域干涉相位&_Φ进行局部频率估计，得到干涉相位所在的频率范 围；
[0024] 该步骤B具体包括：
[0025] 子步骤Bl :对干涉相位Φ的每一像素，对其进行局部频率估计，得到其方位向干 涉相位频率/3和距离向干涉相位频率，具体如下：
[0026] 子分步骤Bla :以该像素为中心取（2M+1) X (2N+1)的估计窗口，该估计窗口的相 位模型可表示为：
[0027]
【权利要求】
1. 一种结合局部频率估计的小波域InSAR干涉相位滤波方法，其特征在于，包括： 步骤A:将InSAR的干涉相位Φ变换到复数域&_φ，分别取复数域干涉相位^1*·的实部 和虚部； 步骤B:对干涉相位Φ进行局部频率估计，得到干涉相位所在的频率范围； 步骤C:对复数域干涉相位的实部和虚部，分别进行尺度为s的小波分解，得到不同子 带的小波系数4,,的频率范围，其中，m，n为小波系数的位置，i为分解的尺度，其范围为1? s； 步骤D:对复数域干涉相位的实部和虚部，分别根据干涉相位所在的频率范围和不同 子带小波系数的频率范围，分别确定有用信息和噪声所在的子带； 步骤E:对复数域干涉相位的实部和虚部，分别对噪声所在子带的小波系数分别进行 通用阈值收缩处理； 步骤F:对复数域干涉相位的实部和虚部，分别对有用信息所在子带的小波系数分别 进行邻域阈值收缩处理； 步骤G:对复数域干涉相位的实部和虚部，分别将通用阈值收缩处理处理后的噪声所 在子带的小波系数和邻域阈值收缩处理后的有用信息所在子带的小波系数共同进行小波 重构，得到滤波后的复数域干涉相位的实部和虚部；以及 步骤H:由滤波后的复数域干涉相位的实部和虚部得到InSAR滤波后的干涉相位。
2. 根据权利要求1所述的小波域InSAR干涉相位滤波方法，其特征在于，所述步骤B包 括： 子步骤Bl:对干涉相位Φ的每一像素，对其进行局部频率估计，得到其方位向干涉相 位频率尤和距离向干涉相位频率X; 子步骤B2:对所有像素估计得到的方位向干涉相位频率分别取最大值和最小值为 .元和，距尚向干涉相位频率/?分别取最大值和最小值为和，得到整幅干涉相 位图所在的频率范围为i.L,m..Lml·ι./;_.I*
3. 根据权利要求2所述的小波域InSAR干涉相位滤波方法，其特征在于，所述子步骤 Bl中，对干涉相位Φ的每一像素得到其方位向干涉相位频率￡和距离向干涉相位频率X 包括： 子分步骤Bla:以该像素为中心取（2M+1)X(2N+1)的估计窗口，该估计窗口的相位模 型为：
其中，炉(w,?)为估计窗口中心点的相位，k、1为窗口内像素相对窗口中心的位移，fa，f； 分别为该窗口沿方位向和距离向的干涉相位频率； 子分步骤Blb:通过最大化估计窗口的相位模型的代价函数来估计该像素沿方位向的 干涉相位频率Z和距离向的干涉相位频率X:
4. 根据权利要求1所述的小波域InSAR干涉相位滤波方法，其特征在于，所述步骤C 中，不同子带小波系数的频率范围为：
其中，< 为子带对应干涉相位图的低频部分，rff, 和分别为对应干涉相位图 在各个尺度的垂直方向、水平方向和对角线方向的高频部分。
5. 根据权利要求1所述的小波域InSAR干涉相位滤波方法，其特征在于，所述步骤C 中：用于小波分解的函数为Symlets函数、Daubechies函数或Coiflets函数，所述尺度s 介于2至8之间。
6. 根据权利要求1所述的小波域InSAR干涉相位滤波方法，其特征在于，所述步骤D根 据干涉相位所在的频率范围和不同子带小波系数的频率范围，分别确定有用信息和噪声所 在的子带包括： 分别判断不同子带小波系数的频率范围af? , 与 iUjxlUj是否有交集，如果有交集，则判断该子带为有用信息所在的子带；如果 交集为空，判断该子带是噪声所在子带。
7. 根据权利要求1所述的小波域InSAR干涉相位滤波方法，其特征在于，所述步骤E 中，按照下式对噪声所在子带的小波系数进行通用阈值收缩处理：
其中，为阈值收缩后的小波系数，sgn(·)为符号函数，下标+表示保持正值不变， 将负值置零；T为通用阈值。
8. 根据权利要求1所述的小波域InSAR干涉相位滤波方法，其特征在于，所述步骤F中 对有用信息所在子带的小波系数分别进行邻域阈值收缩处理包括： 以当前要处理的小波系数#为中心，选取大小合适的窗口W，令Σ 则阈值收缩后的小波系数P为：c,,=4,,，其中，下标+表示保持正值不变， 将负值置零，T为通用阈值。
9. 根据权利要求7或8所述的小波域InSAR干涉相位滤波方法，其特征在于，所述通用 阈值T按照下式计算： T =^2tog￥ 其中，N为干涉相位图包含的像素数，σ为噪声标准差，按下式进行估计： SJedian(yj) σ =-! L-〇 0.6745
10.根据权利要求1至8中任一项所述的小波域InSAR干涉相位滤波方法，其特征在 于，所述步骤H中，按照下式得到InSAR滤波后的干涉相位"
其中，&和?分别为滤波后的复数域干涉相位的实部和虚部。
【文档编号】G01S7/02GK104459633SQ201410717610
【公开日】2015年3月25日 申请日期:2014年12月1日 优先权日:2014年12月1日
【发明者】李芳芳, 丁赤飚, 雷斌, 林雪, 胡东辉, 仇晓兰 申请人:中国科学院电子学研究所