标题

第 43 卷第 5 期应用科技 Vol.43.5 2016 年 10 月 Applied Science and Technology Oct. 2016 DOI:10.11991 / yykj.201512027 网络出版地址 :http: / / www.cnki.net / kcms / detail / 23.1191.U.20160914.1143.002.html 一种枝切法和质量图相结合的 InSAR 相位解缠算法 王霖郁, 李辉 摘 哈尔滨工程大学信息与通信工程学院, 黑龙江哈尔滨 150001 要 : 针对传统的基于残差点的路径积分相位展开方法对于相干性差 信噪比小的区域, 容易形成一个个独立的孤岛 而无法解缠的缺陷, 提出一种枝切法和质量图相结合的 InSAR 相位解缠的新算法 该算法以基于残差点的路径积分相 位展开方法为基础, 首先对大面积相位进行解缠, 然后通过包裹相位的质量图解缠剩余的孤岛相位以及残差点处的相 位, 其中质量图法采用堆排序算法提高运算速度 枝切法不能得到完全解和质量图法无法阻止相位解缠误差的扩散的 缺点得到了有效的克服, 通过该算法在实测和仿真数据的处理上, 也说明该算法即利用到了路径积分法的解缠速度快的 特点, 又利用到了质量图算法的相位解缠的有效性的特点, 能有效处理大规模缠绕相位的解缠需求, 相比于其他算法能 得到更高的精度 关键词 : 干涉 SAR; 枝切法 ; 质量图法 ; 大规模相位展开 ; 堆排序 中图分类号 :TP393 文献标志码 :A 文章编号 :1009 671X(2016)05 049 05 An algorithm based on the branch cut and quality map for InSAR phase unwrapping WANG Linyu, LI Hui College of Information and Communication Engineering, Harbin Engineering University, Harbin 150001, China Abstract:To overcome the shortcomings that it is easy to form a separate island without unwrapping caused by the traditional path integral phase unwrapping method, a new algorithm of InSAR phase unwrapping based on the tradi tional path integral phase unwrapping and quality map method was proposed. This new algorithm firstly applied branch cut method to unwrap large scale wrapped phase, then unwrap the remainder of the phase of separate islands by quality map, using heap sort algorithm to increase its operating efficiency. Simulation and experimental results show that the proposed algorithm can overcome the shortcomings of the branch cut method, which can not only get the complete phase unwrapping, but also surmount the difficulty in stopping the diffusion of phase unwrap ping error caused by the quality map method. In actual test and simulation, the method is proven to have both char acteristics of fastness in unwrapping of path integral and effectiveness in phase unwrapping of quality map method, and meet the requirement of unwrapping large scale wrapping phase, having higher precision by comparison with other algorithms. Keywords:interferometric SAR; branch cut method; quality map method; large scale phase unwrapping; heap sort 干涉合成孔径雷达技术 (InSAR) 是当前空间对 地观测领域的为获得高程信息的热门技术 它通过 两部天线对同一地区的多重的单景成像, 然后通过 对多重的单景成像进行相应的干涉处理, 然后通过 对缠绕的相位进行解缠的处理之后, 三维数字高程 图 (digital elevation models, DEM) 就可以通过解缠 收稿日期 :2015-12-31. 网络出版日期 :2016-09-14. 基金项目 : 国家自然科学基金项目 (61101141). 作者简介 : 王霖郁 (1977-), 女, 副教授 ; 李辉 (1990-), 男, 硕士研究生. 通信作者 : 李辉,E mail: lihui_nefu@ 163.com. 后的相位结合干涉合成孔径雷达的模型从而更进一 步得到地表形变等物理信息, 测量原理来源于杨氏 双缝干涉实验, 其中, 如何从缠绕的相位信息中获得 与高程直接对应的解缠相位是 InSAR 处理的关键 一步 [1-7] [1] 相位解缠算法主要有路径跟踪算法和最小 范数法 [4] 路径跟踪算法中的经典算法是枝切法, 它通过在残差点之间建立路径来分割极性不为零的 残差点, 从而避免解缠时的局部误差传播到别的待 解缠区域 质量图法也属于路径跟踪算法, 但是不

50 应用科技第 43 卷 同于相位解缠的经典枝切算法, 它是直接利用能反 应缠绕相位质量的质量图, 然后从质量图的高质量 部分开始解缠缠绕相位, 依次向低质量缠绕相位区 域进行解缠 质量图法的明显缺陷是没有像经典枝 切法一样考虑相邻残差点引入的误差, 从而在低质 量区域, 解缠效果会非常不理想 此外, 其他解缠算 法在相位解缠过程中也都有一定的局限性 最小范 [7] 数法以全局误差作为算法的衡量指标, 容易造成 局部区域误差的偏移, 从而影响其他高质量区域的 解缠效果 本文提出一种枝切法和质量图法相结合的一种 合成算法 该算法不同于以往的解缠算法, 不仅利 用了枝切法中的残差点的概念, 以残差点的分布为 依据, 优化质量图中阀值的确定方法 [8], 还以此指 导质量图法的相位顺利解缠, 同时为了进一步的提 高算法的运行效率, 在质量图法中引入了堆排序的 概念, 使质量图算法的效率得到改进和提高 由于枝切法和质量图法的解缠都可以归属为路 径积分的一个类别, 所以这 2 种方法都很好地满足 相位解缠的一致性, 即展开相位重新缠绕后与原始 缠绕相位有很好的一致 所以提出的新的基于枝切 法和质量图法的融合算法就保证了相位解缠的一致 性问题 通常, 枝切法得不到相位解缠的完整解, 而 由于质量图法解缠过程没有考虑到 残差点 的概 念, 所以很容易在解缠过程中引入由残差点造成的 相位误差 所以针对这 2 种算法各自的优劣, 提出 一种枝切法和质量图法相结合的新型算法, 仿真和 实测数据处理结果表明该方法很好地克服了枝切法 和质量图法在相位解缠过程中的缺点, 同时保留了 [9-12] 两者的优点, 相比于其他方法具有更高的精度 1 堆排序算法 先简单的介绍下计算机的对象结构, 堆数据结 构 它是一棵完全二叉树 [13] 由于每一次的存取 分别对应堆的根节点和叶子节点, 所以这棵完全二 叉树有时候并不总是满的, 对于最后一层的叶子节 点, 往往跟数据的插入操作有关, 因此它的大小总是 保持着动态的变化 现在讨论堆 { r 1,r 2,,r n } 应该满足的性质, 并且做一些传统堆排序算法的简单 介绍和相关性质的整理 : r 2i r { i r 2i 或 r 2i +1 { (1) r 2i +1 式中 i = 1,2,,n / 2 该定义是 2 种常见的堆结构, 分别是最小堆和最大堆 本文使用的堆结构要满足 邻接列 的性质, 即取出其中最大的质量对应的相位的位置, 然后将该点邻接的相位的质量系数放入里面, 进行大小比较, 然后放到对应的位置上 通过需求描述的对比, 可以得知我们需要的堆结构应该是最大堆, 它的根节点的元素的值就是所有节点的最大的那个, 取出的操作也非常简单 2 一种基于堆排序的质量法快速实现 本节通过仿真意大利某火山的干涉数据分析验证基于堆排序的相位解缠算法的整体效率, 如图 1 所示 图 1 堆排序质量法算法示意图图 1(a) 为做过去平地处理和相位滤波后的缠绕相位, 该干涉数据的大小为 426 pixel 432 pixel, 图 1(b) 为通过所选用的相位导数变化图生成的相干系数, 然后通过堆结构的排序算法, 得到的解缠结果如图 1(c) 通过比较不同尺寸数据的解缠时间进一步的衡量加入堆排序的质量图法的效率问题, 如表 1 表 1 不同尺寸数据的解缠效率比较图像尺寸 / pixel 本文算法 / s 传统质量图法 / s 加速比 128 128 6.17 14.38 2.33 256 256 16.25 126.94 4.51 426 432 145.76 1138.38 7.81 表 1 中, 把传统的基于质量图的路径跟踪算法与新提出的基于堆排序的质量图算法做了对比和验证, 通过对不同大小的干涉图的解缠绕, 得出本文算法的运行时间, 可以很明显地看到所提出算法的优势 本次实验的环境为基于联想 Y460 ( 英特尔 Core i3 M 330 @ 2.13 GHz,6 GB 内存 ) 和 MATLAB R2014a 来进行的 3 合成相位展开算法 本文提出的方法主要是基于路径积分的思想 该方法首先通过使用传统的枝切法, 解缠出大范围的缠绕相位, 然后针对枝切法无法解缠的孤立区域以及残差点处的缠绕相位, 通过质量图法中的伪相干系数图为参考, 对剩余未解缠相位进行进一步的解缠 具体内容和步骤介绍如下

第 5 期 王霖郁, 等 : 一种枝切法和质量图相结合的 Insar 相位解缠算法 51 3.1 利用枝切法展开较大范围相位 由复变函数的相关理论可知闭合回路的积分满 足如下等式 : Δφ( r) dr = 2πN r (2) c 式中 N r 为闭合积分回路里面的残差点极性的代数 和 当式 (2) 里面的闭合路径的残差点的极性为 零, 或者等效极性为零时, 解缠的相位正确性与解缠 路径无关 但是这部分的解缠过程只针对的是残差点比较 少的部分 但是对于残差点过于密集的区域, 枝切 线无法完全放置下去, 并且容易形成孤立区域, 使相 位解缠的过程无法进行, 所以这部分的缠绕相位是 枝切法无能为力的, 这部分枝切法无法解缠的相位 可通过质量图法进一步解缠 3.2 利用质量图法展开剩余区域 利用枝切法解缠大范围缠绕区域后, 通过质量 图法解缠剩余的未解缠区域 图 2 为本文算法示意 图, 其中图 2( a) 中黑色表示枝切法解缠出来的相 位, 灰色为邻接相位 假如箭头指向为质量图法的 基准相位, 则解缠过程通过该基准相位向外蔓延 通过判断待解缠相位及其四周的解缠状况进行相应 处理 如果是已解缠相位, 则直接跳过该点并将其 标记为已解缠相位, 如果是未解缠相位, 则将其四周 各点放入邻接矩阵, 并扫描质量最高的相位是否为 已解缠相位, 是则将该相位视为参考相位对待解缠 点进行解缠, 否则重复将该点四周相位放入邻接矩 阵并重复操作 从图 2(a) 到图 2(b) 的过度是假设 图 2(b) 的箭头所指的相位对应的质量数在图 2(a) 的灰色相位中为最高, 因此该点先得到解缠, 并且将 该点及其周围未接触相位标记为灰色 图 2 基本操作步骤如下 : 本文算法示意图 1) 利用 3.1 中描述的枝切线法, 布置枝切线并 通过像元扩散法对缠绕相位进行积分来获得解缠相 位, 此时的解缠未完成, 因为对于大规模缠绕相位, 很容易形成闭合的枝切线而无法解缠 2) 利用伪相干系数图的定义得到相位质量图, 然后标记其中的已解缠相位, 并且把这些相位对应的 质量图中的质量数修改为 0, 因为对于伪相干系数图 来说, 数值越小, 则表示其相应的相位质量越高 3) 以步骤 1) 中基准相位为开始, 对其周围的 4 个相邻点分别进行遍历, 对其中质量最高的相位进 行解缠, 若其已在步骤 2) 中标记为解缠, 则跳过该 点, 否则标记其为已解缠相位 然后将该点周围的 4 个相邻点中未被处理的点标记到 邻接矩阵 中 4) 从 邻接矩阵 中找出质量值最小, 也即质量 最高的点, 对其重复步骤 3) 的处理, 这一过程是通 过堆排序实现的, 不在详述 5) 重复步骤 4), 直到 邻接矩阵 为空, 也即所 有像素点都得到处理 4 实验结果与分析 本节实验拟通过从仿真数据和实测数据两个维 度反应出本文算法在相位解缠过程中, 对解缠过程 指标的提高 这其中提高的指标反应在成像上的对 比和精确度上的对比, 从而验证本文算法的性能 4.1 仿真实验 这里的质量图选的是相位倒数变化图, 然后数据 是采用大小为像素点 128 128 的仿真相位, 它由 MATLAB 软件中 peaks 实现, 效果如图 3(a) 在仿真 过程中, 为了模拟现实数据的噪声的干扰, 手动加入 相干系数为 0.9 的均匀噪声 该相干系数对应的噪 声方差为 0.422 5 rad, 它的分布区间假设为 [-b,b), 则可得 b = 1.125 8 rad [14] 对加入噪声的待处理相位 进行缠绕处理后如图 3(b) 通过枝切法得到的残差 点如图 3(c) 所示 较高质量区域内缠绕相位的展开 如图 3(d) 引入质量图进行相位展开后如图 3(e), 从而得到最终解缠结果如图 3(f), 其与图 3(a) 的未 缠绕相位几乎一致 为了对本文算法进行定性地分析和比较, 我们 首先做出的是各类算法解缠后对应的解缠误差的二 维分布比较 其中解缠误差的定义是该点展开相位 与未加均匀噪声前的待处理的真实相位上对应点相 位值之差的绝对值 与常用的解缠算法对比, 效果 如图 4 可以看到本文方法误差没有很集中, 整体 的分布范围可以通过色调偏暗得出 并且可以看出 误差的幅度也不会太大 ; 最小二乘法的误差分布因 为算法本身的问题, 其误差范围最大, 并且在相位变 化较大缠绕相位处的解缠误差明显较大 同本文算 法主要是误差的范围较大, 这也是最小二乘法的局 限 ; 网络流法的误差也主要分布于相位变化较大的

52 应用科技第 43 卷 区域, 虽然误差范围较小, 但是误差的幅度非常明 显, 这也是跟本文相比的不足之处 图 3 本文方法的仿真结果 图 4 本文算法与其他算法的二维比较 最后, 为了进一步的检验本文算法相较于最小二乘法和网络流算法的稳健性, 我们通过改变加入 4.2 实测数据实验实测数据采用的是文献 [2] 提供的某山地干涉 均匀噪声的方差的幅值, 然后计算出各类算法展开后所得到的均方根误差, 见表 2 从表 2 中可以看到 : 本文方法的误差值要小于其他的 2 种算法, 而这充分说明了本文算法比最小二乘法和网络流算法有更好的稳健性 相位数据, 该数据为 512 pixel 512 pixel 的缠绕相位, 如图 5(a) 所示 ; 其对应的残差点分布如图 5(b); 部分未孤立的相位通过枝切法展开得到的枝切线放置如图 5(c); 部分展开相位如图 5(d); 剩余的没有展开的区域留待进一步通过质量图法进行解缠, 在图 表 2 本文算法与其他算法对应均方根误差比较 5(e) 中用黑色表示 ; 利用相位倒数变化法得到的质量 噪声方差 对应的相干系数 本文算法误差 最小二乘法误差 网络流法误差 图如图 5(f); 通过质量图法展开图 5(e) 中未展开的区域 ; 图 5(g) 为最后的展开相位 为了直观地比较 0.422 5 0.9 2.498 9 8.523 5 3.562 4 和分析本文算法对缠绕相位解缠的有效性和一致性, 把最终得到的展开相位重新缠绕如图 5(h), 该缠绕 0.810 0 0.8 4.084 4 9.945 6 7.432 3 相位与原始缠绕相位非常吻合, 充分说明本文方法能 1.210 0 0.7 5.972 7 11.790 8 10.659 9 够恢复出完整且较为平滑的相位曲面, 具有一定的现实可推广性

第 5 期 王霖郁, 等 : 一种枝切法和质量图相结合的 Insar 相位解缠算法 53 5 结束语 为了充分利用经典枝切法和质量图法的一致性和有效性的两大优点, 本文提出了一种即结合了枝切法又兼顾到质量图法的新算法, 该方法对于下面的问题有针对性的诠释 1) 根据枝切法的算法特点, 有效地避免了质量图法误差传递的缺陷, 并利用质量图法的特点有效地解决了枝切法解缠不完全以及容易形成无法解缠的孤立区域的缺陷 2) 结合堆排序算法进一步提高算法整体的运行效率 3) 能够恢复出完整且较为平滑的三维曲面 4) 提高了整体的缠绕相位的解缠精度 参考文献 : [1]GOLDSTEIN R M, ZEBKER H A, WERNER C L. Satellite radanterferometry two dimensional phase unwrapping[ J]. Radio science, 1988, 23(4): 713-720. [ 2] GHIGLIA D C, PRITT M D. Two dimensional phase un wrapping: theory, algorithms, and software [ M ]. New York: Wiley, 1998. [3]WANG C, ZHANG H, LIU Z. Spaceborne synthetic ap erture radanterferometry[ M]. Beijing: Science Press, 2002. [4] PRITT M D, SHIPMAN J S. Least squares two dimensional phase unwrapping using FFT s[ J]. IEEE transactions on ge oscience and remote sensing, 1994, 32(3): 706-708. 图 5 本文算法的实测数据结果 [5] COSTANTINI M. A novel phase unwrapping method based on network programming[ J]. IEEE transactions on geosci ence and remote sensing, 1998, 36(3): 813-821. [6] 刘双亚, 李世强, 冯锦. 改进的基于频谱分割的 InSAR 绝对相位确定方法 [ J]. 国外电子测量技术, 2016(7): 27-33. [7] 张金芝, 黄海军, 毕海波, 等. SBAS 时序分析技术监测现代黄河三角洲地面沉降 [ J]. 武汉大学学报 : 信息科学版, 2016, 42(2): 242-248. [8] 余应淮, 谢仕义, 梅其祥. 基于核回归修正的上采样相位相关精确运动估计算法 [ J]. 计算机应用, 2016, 36 (8): 2316-2321, 2326. [9] 陈强, 杨莹辉, 刘国祥, 等. 基于边界探测的 InSAR 最小二乘整周相位解缠方法 [J]. 测绘学报, 2012, 41(3): 410-416. [10] 钟何平, 唐劲松, 张森. 一种基于质量引导和最小不连续合成的 InSAR 相位解缠算法 [ J]. 电子与信息学报, 2011, 33(2): 369-374. [11] 李芳芳, 占毅, 胡东辉, 等. 一种基于质量指导的 In SAR 相位解缠快速实现方法 [ J]. 雷达学报, 2012, 1 (2): 196-202. [12] 张小京, 刘冰. 数据共轭复用策略 InSAR 干涉相位估计 [J]. 测绘科学, 2016, 41(8): 28-32, 18. [13] 肖枫, 伍吉仓, 刘朝功, 等. 不同质量图在相位解缠算法中的比较分析 [ J]. 大地测量与地球动力学, 2010, 30(2): 80-85. [14]LEE J S, HOPPEL K W, MANGO S, et al. Intensity and phase statistics of multilook polarimetric and interferometric SAR imagery[ J]. IEEE transactions on geoscience and re mote sensing, 1994, 32(5): 1017-1028. 本文引用格式 : 王霖郁, 李辉. 一种枝切法和质量图相结合的 Insar 相位解缠算法 [J]. 应用科技, 2016, 43(5): 49-53. WANG Linyu, LI Hui. An algorithm based on the branch cut and quality map for InSAR phase unwrapping[ J]. Applied Science and Technology, 2016, 43(5): 49-53.