一种干涉相位滤波方法
【技术领域】
[0001] 本发明属于遥感和信号处理的交叉技术领域，特别涉及一种干涉相位的滤波方 法。
【背景技术】
[0002] 对两幅合成孔径雷达图像共辄相乘，并对相乘结果取相位便得到干涉相位。干涉 相位是合成孔径雷达干涉测量中重要的物理量，其质量的好坏将决定最终产品一一数字高 程模型或地形形变量的精度。然而，受去相关因素的影响，干涉相位中总是存在严重的空变 噪声。空变噪声不仅会引入残差点还会破坏干涉相位的分布，从而增加了后续相位解缠的 难度，最终导致产品精度的降低，因此必须通过干涉相位滤波过程予以滤除。
[0003] 现有的干涉相位滤波方法主要分为空域方法和变换域方法。空域方法主要根据干 涉相位中干涉条纹的方向，采用自适应窗口来滤波，以达到滤除噪声并保持细节的效果，这 类方法主要包括均值方法、Lee方法、旋方法等。但这类方法对于低相干区的效果不是很好， 并且容易破坏干涉相位变化剧烈区域的细节信息。
[0004] 变换域方法可以分为频域方法和小波域方法。变换域方法的主要依据是:干涉相 位和噪声在变换域可以得到很好的分离。频域方法主要包括Goldstein方法及其改进方法， 小波域方法主要包括离散小波变换方法等。由于干涉相位和噪声不可能做到完全分离，变 换域方法可能导致过滤波或欠滤波。
【发明内容】
[0005] 本发明提出一种干涉相位滤波方法，能够在有效滤除噪声的同时保持干涉相位的 细节信息。
[0006] 本发明技术方案是:一种干涉相位滤波方法，其特征在于，包括下述步骤：
[0007] 第一步，频率滤波：
[0008] 将干涉相位图转换到复数域得到复干涉相位图，将复干涉相位图分成有重叠区域 的正方形图像块，下述图像块均指正方形图像快；
[0009] 对每一个图像块，分别进行下述处理：
[0010] 用二维傅里叶正变换计算图像块的二维频谱；
[0011]对图像块的二维频谱进行取绝对值、平滑、归一化和幂指数操作得到权重，其中， 幂指数是利用复干涉相位图的噪声标准差和干涉相位图的归一化倒置相位导数方差图确 定；
[0012] 用得到的权重对图像块的二维频谱进行加权，从而实现滤波；
[0013] 对加权后的图像块二维频谱进行二维傅里叶逆变换，得到频域滤波后的图像块；
[0014] 对频域滤波后的图像块，用取均值的方法对重叠区域进行平滑，从而得到频域滤 波后的复干涉相位图；
[0015] 第二步，空域滤波：
[0016] 用复干涉相位图减去频域滤波后的复干涉相位图，得到残余复干涉相位图，用均 值滤波对残余复干涉相位图进行滤波，得到空域滤波后的残余复干涉相位图；用频域滤波 后的复干涉相位图和空域滤波后的残余复干涉相位图相加，得到滤波后的复干涉相位图， 对其取相位得到滤波后的干涉相位图。
[0017] 采用本发明可取得以下技术效果：
[0018] 本方法包含一步频域滤波和一步空域滤波，频域滤波的优势在于滤除噪声的同时 对密集干涉条纹的保持，而空域滤波的优势在于进一步滤除噪声和对部分细节信息的保 持，本方法集中上述两种滤波过程的优势;此外，本发明的另一创新点在于，第一步频域滤 波时采用复干涉相位图的噪声标准差和干涉相位图的归一化倒置相位导数方差图来确定 幂指数，提高了滤波精度;同时通过均值滤波保持了干涉相位的细节信息。
【附图说明】
[0019] 图1为本发明提供的干涉相位滤波方法流程示意图；
[0020] 图2为待滤波的干涉相位图；
[0021 ]图3为Goldstein方法得到的滤波后的干涉相位图；
[0022]图4为Lee方法得到的滤波后的干涉相位图；
[0023]图5为本发明得到的滤波后的干涉相位图；
[0024] 图6为图3,图4和图5的性能对比图。
【具体实施方式】
[0025]下面结合图1对本发明进行详细说明。
[0026] 步骤①:将干涉相位图转换成复干涉相位图并对其分块：
[0027]
[0028] 其中，供是干涉相位图，Φ是复干涉相位图。对复干涉相位图进行分块，图像块根据 实际情况确定，本实施方式中图像块大小为32 X 32像素的正方形，分块与分块之间有重叠， 每个方向重叠14个像素点，然后对于每个图像块完成步骤②。
[0029]步骤②:图像块滤波
[0030] 第（1)步，计算图像块的二维频谱：
[0031] B=FFT2(P)
[0032] 其中，P是图像块，B是图像块的二维频谱，FFT2( ·)二维傅里叶正变换。
[0033]第(2)步:对图像块的二维频谱进行取绝对值、平滑、归一化和幂指数操作得到权 重。
[0034]计算图像块二维频谱的绝对值A:
[0035] A = abs(B)
[0036] 其中，abs( ·)表示取绝对值。然后对A值进行平滑：
[0037] S = smooth(A)
[0038] 其中，S是平滑后的图像块二维频谱的绝对值，sm〇〇th( ·)是二维高斯滤波器且其 窗口大小根据实际情况确定，本发明实施方式中为7X7像素。接着对S进行归一化：
[0039] N=normalization(S)
[0040] 其中，normalization!； ·)是归一化算子。最后，对N取幂指数得到权重W:
[0041 ] ff=Na
[0042] 其中：
[0043] 0 = (1.10)^
[0044]
[0045] 其中，σ是复干涉相位图像块的噪声标准差，f是滤波窗口对应的干涉相位图中归 一化倒置相位导数方差图的均值，计算公式为：
[0046] σ = 1.4826median( | VP-median( VP) | )
[0047]
[0048]
[0049] 其中，PDV是滤波窗口对应的干涉相位图的相位导数方差图，mean( ·)是取均值算 子。
[0050] 第(3)步:对图像块的二维频谱加权：
[0051 ] Bw=B · ff
[0052] Bw是加权后的图像块二维频谱。
[0053] 步骤③，计算频域滤波后的复干涉相位图
[0054] 得到所有图像块加权后的二维频谱后，将它们通过二维傅里叶逆变换计算加权后 的复干涉图像块Pw:
[0055] Pw=IFFT2(Bw)
[0056] 其中，IFFT2( ·)是二维傅里叶逆变换。重叠区域的频域滤波结果来自不同的图像 块，因此采用对所有重叠区域的频域滤波结果进行取均值的方法，从而得到频率滤波后的 复干涉相位图如 requency 〇
[0057] 步骤④:计算空域滤波后的残余复干涉相位图 [0058]第（1)步，计算残余复干涉相位图itresiduai:
[0059] l^residual - Φ-^frequency
[0060] 第(2)步，计算空域滤波后的残余复干涉相位图Itspatial:
[0061 ] l^spatial - IllGanf ? 1 ?θΓ (l^residual)
[0062]其中，meanfilter( ·)是均值滤波器，其窗口大小根据实际情况确定，本发明实施 方式中为5X5像素。
[0063]步骤⑤:计算滤波后的干涉相位图
[0064]第（1 )步，计算滤波后的复干涉相位图filtered :
[0065] Φ? iltered - l^spatial+ltf requency
[0066] 第(2)步，计算滤波后的干涉相位图：
[0067]
j：
[0068]其中，angle( ·)是取相位算子。
[0069]关于图2-图6的详细说明：
[0070] 图2为待滤波的干涉相位图。该干涉相位截自意大利Etna火山的实测数据，大小为 500 X 500像素，其中横坐标是列索引，单位是像素;纵坐标是行索引，单位是像素。
[0071] 图3为Goldstein方法得到的滤波后的干涉相位图。可以看到，虽然干涉条纹得到 了很好保持，但是部分噪声仍然存在。
[0072] 图4为Lee方法得到的滤波后的干涉相位图。可以看到，噪声得到了很好的滤除，但 是干涉条纹得到了破坏，尤其在干涉条纹密集的区域。
[0073] 图5为本发明得到的滤波后的干涉相位图。可以看到，噪声不但得到了有效的滤 除，干涉条纹也得到了很好的保持。
[0074] 图6为图3,图4和图5的性能对比图。该性能对比图包含了一个性能指标:残差点。 残差点的计算采用下述方法：对于滤波后的干涉相位图中任意的4个相邻的像素点： 供(/,./：)，M/ + 1,./)，以/,./ + 1)，+ + ，我们分别计算两两像素点的缠绕梯度值： ' ζ . X / \ "V
[0075]
[0076] 其中，(：{>}是缠绕算子，其取值范围为[-31，31)，&1，&2，& 3，&4是缠绕梯度值。然 后对4个缠绕梯度值求和：
[0077] Δ total = Δ 1+ Δ 2+ Δ 3+ Δ 4
[0078] 如果Δ tQtal矣〇,称像素点（i，j)为残差点，残差点数量越少，滤波后的干涉相位图 质量越好。可以看出，相比Goldstein方法和Lee方法，本发明极大的减少了干涉相位图中的 残差点数量，更有利于后续的相位解缠。
【主权项】
1. 一种干涉相位滤波方法，其特征在于，包括下述步骤： 第一步，频率滤波： 将干涉相位图转换到复数域得到复干涉相位图，将复干涉相位图分成有重叠区域的正 方形图像块，下述图像块均指正方形图像快； 对每一个图像块，分别进行下述处理： 用二维傅里叶正变换计算图像块的二维频谱； 对图像块的二维频谱进行取绝对值、平滑、归一化和幂指数操作得到权重，其中，幂指 数是利用复干涉相位图的噪声标准差和干涉相位图的归一化倒置相位导数方差图确定； 用得到的权重对图像块的二维频谱进行加权； 对加权后的图像块二维频谱进行二维傅里叶逆变换，得到频域滤波后的图像块； 对频域滤波后的图像块，用取均值的方法对重叠区域进行平滑，得到频域滤波后的复 干涉相位图； 第二步，空域滤波： 用复干涉相位图减去频域滤波后的复干涉相位图，得到残余复干涉相位图，用均值滤 波对残余复干涉相位图进行滤波，得到空域滤波后的残余复干涉相位图；用频域滤波后的 复干涉相位图和空域滤波后的残余复干涉相位图相加，得到滤波后的复干涉相位图，对其 取相位得到滤波后的干涉相位图。
【专利摘要】本发明提出一种干涉相位滤波方法。技术方案包括下述步骤：第一步，频率滤波：将干涉相位图转换到复数域得到复干涉相位图，将复干涉相位图分成有重叠区域的正方形图像块，对每个图像块分别进行频域处理后得到频域滤波后的复干涉相位图；第二步，空域滤波：用复干涉相位图减去频域滤波后的复干涉相位图，得到残余复干涉相位图，用均值滤波对残余复干涉相位图进行滤波，得到空域滤波后的残余复干涉相位图；用频域滤波后的复干涉相位图和空域滤波后的残余复干涉相位图相加，得到滤波后的复干涉相位图，对其取相位得到滤波后的干涉相位图。本发明能够在有效滤除噪声的同时保持干涉相位的细节信息。
【IPC分类】G06T5/00
【公开号】CN105469368
【申请号】CN201510853164
【发明人】黄海风, 汪洋, 董臻, 张永胜, 孙造宇, 何志华, 张启雷, 杜湘瑜, 刘奇, 陈筠力, 陈国忠, 陈重华, 李威, 童庆为, 景海涛
【申请人】中国人民解放军国防科学技术大学, 上海卫星工程研究所
【公开日】2016年4月6日
【申请日】2015年11月30日