基于克里金插值和Chan-Vese模型边界探测的InSAR相位解缠方法_王升.pdf

收稿日期2019-09-05 基金项目国家自然科学基金项目 （编号 41464001） ， 江西省自然科学基金项目 （编号 20151BAB203042） 。 作者简介王升 （1995） ， 男， 硕士研究生。通讯作者朱煜峰 （1981） ， 男， 副教授， 博士， 硕士研究生导师。 总第 521 期 2019 年第 11 期 金属矿山 METAL MINE 基于克里金插值和Chan-Vese模型边界探测的 InSAR相位解缠方法 王升朱煜峰 1 （东华理工大学测绘工程学院， 江西 南昌 330000） 摘要相位解缠是InSAR技术测量高程中的重要步骤， 相位解缠的精度直接决定InSAR测量的精度。目前 常用的相位解缠方法的解缠精度易受低相干区域的干涉噪声影响， 针对在低相干区域相位解缠不准确的问题， 提 出了基于克里金插值和Chan-Vese模型的相位解缠方法。该方法在低相干区域进行克里金插值， 重新计算低相干 区域的干涉相位值， 高相干区的干涉相位值不变， 生成新干涉图， 使用Chan-Vese法提取出新干涉图的干涉条纹边 界， 按照一定的方法计算出跃变数， 达到相位解缠的目的。选取西藏自治区班戈县某矿山为试验区， 采用C波段的 Sentinel-1A合成孔径雷达数据对的干涉相位图为试验数据， 对所提出的相位解缠方法与Goldstein枝切法等其他5 种方法的相位解缠结果进行了对比分析。研究表明 所提方法具有一定的适用条件， 在相干性较高的区域， 该方法 与其他方法的相位解缠结果差异很小； 在低相干区域跨越不超过2个干涉条纹块的情况下， 相比于其他方法， 由所 提方法相位解缠结果反演出的高程与SRTM1数字高程模型之间的差异在合理范围内， 说明该方法仍然适用， 而在 低相干区跨越2个以上干涉条纹块的情况下， 该方法不适用。 关键词InSAR相位解缠克里金插值Chan-Vese模型相干性 中图分类号TD17， P236文献标志码A文章编号1001-1250 （2019） -11-142-09 DOI10.19614/ki.jsks.201911024 InSAR Phase Unwrapping Based on Kriging Interpolation and Chan-Vese Model Boundary Detection Wang ShengZhu Yufeng2 （Faculty of Geomatics， East China University of Technology， Nanchang 330000， China） AbstractPhase unwrapping is an important step in InSAR elevation measurement， the accuracy of phase unwrapping directly determines the accuracy of InSAR measurement.The unwrapping accuracy of commonly used phase unwrapping meth- ods is susceptible to interferometric noise in low coherence region. In order to solve the problem that the result of phase un- wrapping is inaccurate in low coherence region， a phase unwrapping based on Kriging interpolation and Chan-Vese model is proposed.In this ， Kriging interpolation is pered in the low coherence region，and the interferometric phase values in the low coherence region are recalculated，the interferometric phase values in the high coherence region are unchanged， and a new interferogram is generated.The interferometric fringe boundary of the new interferogram is extracted by Chan-Vese ， and the jump coefficient is calculated according to certain s， so as to achieve the purpose of phase unwrapping.A mining area in Bange County of Tibet was selected as the test area， the interference phase images which are pro- duced by the pair of sentinel-1A synthetic aperture radar data C-band are used as experimental data， the phase unwrapping results of proposed unwrapping are compared with those of other unwrapping such as branch cut pro- posed by Goldstein.The study results show that proposed has some applicable conditions， there are little differences be- tween the phase unwrapping results of this and other s in the high coherence region.In the condition of low con- herence region spanning no more than two interferometric fringes， compared with other s，the difference between the uation inverted from the unwrapping result of proposed and SRTM1 digital elevation model is within a reasonable range， thus， the proposed is still suitable.However， in the condition of low coherence region spanning more than two in- Series No. 521 November2019 142 ChaoXing terferometric fringes， the proposed is not suitable. KeywordsInSAR， Phase unwrapping， Kriging interpolation， Chan-Vese model， Coherence 星载合成孔径干涉雷达 （InSAR） 技术是一种当 前极有潜力的遥感测绘技术， 具有高精度、 高空间分 辨率和不易受天气条件影响等优点 ［1］， 目前已经被广 泛应用于地表地形测绘、 地震变形、 地面沉降、 滑坡 监测等领域 ［2］。相位解缠是InSAR技术中的关键步 骤， 该步骤是在测量出的相位中增加合适的整数周 期， 用以恢复出真实的干涉相位 ［3］。 InSAR相位解缠结果的精度直接影响反演出的 DEM的可靠性， 部分学者针对相位解缠方法进行了 探索， 提出了一系列的相位解缠算法， 主要分成三大 类， 即 路径积分法、 最小范数法和网络流法 ［4-6］。最 经典的路径积分法是Goldstein等 ［7］于1988年提出的 枝切法， 具体思路是计算干涉相位图中的正残差点 和负残差点， 根据一定的规则连接正负残差点， 使所 有的残差点达到电荷平衡状态， 计算路径积分并在 相位解缠时避开这些连接线， 从而减少干涉噪声对 解缠过程的影响。枝切法也有一定的缺陷， 解缠结 果易形成没有解缠的 “孤岛” 区域， 为此， Xu等 ［8］提出 质量图指导法， 它是在枝切法基础上的一种改进， 不 同于枝切法的相位解缠路径， 质量图指导法按照干 涉相位质量高低的顺序进行相位解缠。最小范数法 最常用的范数是L2范数， 也称为最小二乘法， 是通过 解算泊松方程使得相位解缠结果在全局上最优， 由 于最小二乘法在相位解缠过程中需要通过迭代求解 方程， 所以解缠效率较低 ［9］。网络流法是基于图论的 一种相位解缠方法， 最著名的是最小费流量法， 它利 用有向网络最优化的思想， 在一些情况下， 能够取得 较好的相位解缠结果 ［10］。总体而言， 目前存在的相 位解缠算法都各有特点， 没有一种完全适用于各种 情况的相位解缠算法。 在整幅干涉图的干涉相位质量较高的情况下， 采用上述方法能够进行正确的相位解缠， 然而， 干涉 相位质量易受到研究区域及其周围区域复杂环境的 影响， 如高大的山脉对微波信号产生阻挡效应， 形成 雷达阴影区域， 阴影区内的相干性较低， 形成低相干 区域， 测量出的干涉相位常带有无规律的噪声毛刺 现象 ［11］。上述3类相位解缠方法难以克服干涉噪声 的影响， 从而易造成相位解缠结果错误， 进一步造成 InSAR高程测量不精确。因此， 本研究提出基于克里 金插值和Chan-Vese模型的相位解缠方法， 目的是在 低相干区域实现正确的相位解缠。该方法在低相干 区域进行克里金插值， 尽可能减少低相干区域的干 涉噪声， 使用Chan-Vese法提取出干涉条纹块边界， 由邻接矩阵计算出跃变数并进行相位解缠。本研究 提出的相位解缠方法能够在特定的低相干区域实现 正确的相位解缠， 即在低相干区域跨越不超过2个干 涉条纹的情况下， 该方法能够正确解缠并且可提高 InSAR高程测量精度。 1InSAR相位解缠算法 1. 1克里金插值法 克里金插值是一种根据已知点数值计算未知点 数值的一种插值方法， 由于其具有线性、 最优估计、 无偏的特点， 而被广泛应用于地形数据重构、 降水预 测、 地质矿产勘查等领域 ［12-13］。克里金插值的核心是 半变异函数理论， 半变异函数用来描述随着空间位 置的变化而产生的空间随机变量的变化， 且随着距 离的增大而减小 ［13］。半变异函数表达式为 γi,j 1 2 Dzi- zj σ2- Covzi,zj，（1） 式中，Dzi,zj是第i个与第j个测量点的测量值之差 的方差。式 （1） 满足以下假设 区域内的每个测量点 的方差都是σ2， 有了这个前提假设， 半变异函数γi,j可 与两测量点zi与zj之间的协方差Covzi,zj产生联系。 Covzi,zj通常按照步长分组近似解算出， 根据已知 值， 计算出离散的半变异函数值， 而后使用拟合模型 拟合出不同距离的半变异函数值， 常用的模型有球 状模型、 指数模型和高斯模型 ［14］。未知点插值数值 是其周围已知数值的加权和， 为满足插值数值无偏 性， 有如下关系式成立 E ∑ i 1 n λizi E z0，（2） 式中，λi为每个zi所对应的权重；z0为未知点的数值； zi为已知的n个点中第i个点的数值。当式 （2） 满足内 蕴平稳假设时， 可推导出 ∑ i 1 n λi 1，（3） 式 （3） 表示所有的权重之和为1。据此建立最优估计 量限定条件， 未知点估计量与其真实值之差的方差 最小， 可表示为 D ∑ i 1 n λizi- z0 min.（4） 联立式 （1） 、 式 （3） 与式 （4） 可得关于权重λ和半 变异函数γ取最小值的表达式 2019年第11期王升等 基于克里金插值和Chan-Vese模型边界探测的InSAR相位解缠方法 143 ChaoXing 2∑ i 1 n λiγi,0- γ0,0-∑ i 1 n ∑ j 1 n λiλjγi,j min.（5） 式 （3） 和式 （5） 构成了条件极值问题， 通过拉格 朗日乘数法解算出未知点周围的已知点权重， 通过 对未知点进行克里金插值得到的数值即为其周围数 值的加权和， 即未知点z0的克里金插值的数值为 ∑ i 1 n λizi。 1. 2Chan-Vese模型 Chan-Vese模型由Chan ［15］等于2001年提出。该 模型是在Mumford-Shah图像分割模型的基础上改进 而来， 并且引入水平集方法进行图像分割， 相比于 Mumford-Shah图像分割模型， Chan-Vese模型能够减 少图像噪声对图像分割的影响 ［16］。在使用 Chan- Vese模型提取物体轮廓时不依赖于局部梯度， 通过 初始轮廓曲线不断演化， 逐渐接近物体的真实轮廓， 具体方法是通过最小泛函原理计算出图像中的轮廓 曲线。Chan-Vese模型的能量函数公式为 ［17］ EC,c1,c2 vAC μLC λ1∫ Ω1 ||u x,y - c1 2 H[]ϕ x,y dxdy λ2∫ Ω2 ||u x,y - c2 2{ }1 - H[]ϕ x,y dxdy， （6） 式中，A C为演化的轮廓曲线C所包围的区域面积； Ω1和Ω2分别表示被轮廓曲线C分割的内部区域与外 部区域；L C为轮廓曲线的长度；ux,y为轮廓曲线 C所包围区域内的所有像素点的灰度值；c1、c2分别为 区域内和区域外像素灰度值的平均值；v、μ、λ1、λ2都 为能量项的系数；ϕx,y为水平集函数；H[ ]⋅为 Heaviside函数； 能量函数公式EC,c1,c2是关于曲线 的泛函， 使用变分法， 根据能量最小准则， 初始的轮 廓不断演化， 轮廓曲线根据如下梯度下降公式进行 演化 ［18］ ∂ϕ ∂t δ ϕ ■ ■ ■ ■ ■ ■ ■ ■ μdiv ∇ϕ ||∇ϕ - v - λ1 u - c12 λ2 u - c22, 6 式中，δ ⋅和div ⋅分别为Dirac函数和散度函数；t为 轮廓曲线的演化时间； 为书写方便， 代表图像灰度值 的ux,y用u代替，ϕx,y用ϕ代替。 在应用Chan-Vese模型提取干涉条纹边界时， 首 先需要手动确定出大致边界， 然后根据变分公式迭 代并计算干涉条纹边界参数， 当能量泛函达到最小 值时， 此时由边界参数确定的边界即为精确的干涉 条纹边界。 1. 3基于克里金插值和Chan-Vese模型边界探测 解缠方法 SAR影像对经过配准和干涉处理后得到的干涉 相位是缠绕相位， 缠绕相位无法直接用于InSAR高程 测量， 相位解缠的目标是恢复出真实的干涉相位值， 由于缠绕相位的范围在为-π～π， 其与真实干涉相位 值相差整数倍的周期， 即2nπ， 其中，n也被称为跃变 数， 所以相位解缠的核心是确定跃变数n的值。在 InSAR获取的地形干涉相位图中， 完全相干区域的干 涉条纹清晰， 首先手动设置大致的干涉条纹边界， 然 后使用Chan-Vese模型经过不断迭代， 逼近精确的干 涉条纹边界， 最后根据精确的干涉条纹边界计算出 邻接矩阵， 并根据邻接矩阵计算出各个干涉条纹块 的跳跃系数k， 从而完成干涉条纹相位解缠。 实际干涉相位图中的每个像素不可能完全相 干， 绝大多数的相干系数值都小于1， 有些区域的相 干系数普遍低于其周围的相干系数， 形成低相干区 域。低相干区域的干涉相位存在较多噪声， 该类区 域会对相位解缠过程产生一定的干扰， 并影响InSAR 高程测量精度 ［19］。为此， 本研究将克里金插值法引 入相位解缠方法中， 用于减少低相干区域的干涉相 位噪声。本研究相位解缠过程为 ①从左到右、 从上 到下依次扫描InSAR干涉相位图和相干系数图， 并设 置相干系数阈值， 高于该阈值的像素保留原来的干 涉相位值， 低于该阈值像素的干涉相位值被认为含 有较大的噪声， 根据克里金插值法计算出新值替换 原来的干涉相位值； ②在经过插值的干涉相位图中， 手动设置干涉条纹大致的边界轮廓， 采用Chan-Vese 模型提取干涉条纹精确的边界轮廓， 进一步提取出 干涉条纹块； ③根据各个干涉条纹块之间的位置关 系， 生成邻接矩阵， 通过邻接矩阵计算各个干涉条纹 块的跃变数n， 每个干涉条纹块的相位值加上2nπ， 完 成相位解缠， 如图1所示。 2试验结果与分析 2. 1试验数据 试验采用Sentinel-1A卫星的SAR单视复视数据 对， 成像时间分别为2017年10月28日和2017年12 月3日， 极化方式VV， 空间基线206.484 m， 时间基线 36 d， SAR 成像区域为西藏自治区班戈县境内某矿 区， 在对SAR数据对进行配准、 干涉、 去平地相位和 Goldstein滤波处理后， 得到含有地形相位的干涉图。 从干涉图中裁剪出3块区域作为验证本研究方法有 效性的数据， 并分析适用条件， 第一块区域相干性较 高， 第二块区域相干性较低且干涉条纹较稀疏， 第三 块区域相干性较低且干涉条纹较密集。 金属矿山2019年第11期总第521期 144 ChaoXing 2. 2相干性较高区域 选取一块相干性较高的区域作为试验区域， 图2 （a） 为相干系数图， 总体相干系数值较高， 与其对应的 同一区域干涉相位图 （图2 （b） ） 的干涉条纹边界较为 清晰， 且每个条纹内的连续性较好， 干涉随机噪声较 少， 信噪比较高。 使用本研究方法进行相位解缠试验， 相干系数 阈值为0.6， 相干系数低于0.6的干涉相位被认为是低 相干区域， 干涉相位值的误差较大， 使用克里金插值 法对该类低相干区域的相位值进行重新计算， 克里 金插值的距离-半变异函数采用高斯函数拟合。经 过克里金插值后生成新干涉相位图， 在新干涉相位 图中使用Chan-Vese模型分割出干涉条纹块， 共计分 割成12块， 每块干涉条纹标上序号， 如图3所示， 曲 线所包围的区域是通过Chan-Vese模型分割出的干 涉条纹块。根据这12块干涉条纹的位置相邻关系， 写成12行12列的邻接矩阵 （图4） ， 邻接矩阵按照如 下规则生成 若序号为i与序号为j的干涉条纹块具 有公共边界线， 且边界线两侧的干涉相位值之差为 正数， 则邻接矩阵的第i行第j列的邻接矩阵元素值 为1， 若为负数， 则元素值为-1， 其他情况时， 邻接矩 阵的元素值为0， 邻接矩阵的对角线与下三角的元素 全为0。从邻接矩阵的第1行第1列出发， 按图4中的 箭头所指方向作邻接矩阵元素的累加和， 若累加最 后的元素值不为0， 则累加和就是干涉条纹块序号， 即为该元素列号的干涉条纹块的跃变数。按照上述 方法计算出12块干涉条纹的跃变数n分别为0、 1、 1、 0、 0、 1、 0、 1、 1、 1、 1、 1， 每块干涉条纹块的相位值增加 2nπ得到解缠相位值， 解缠结果如图5 （a） 所示。分别 采用Goldstein枝切法、 最小费流量法、 Delaunay最小 费流量法、 质量图指导法、 区域生长法进行相位解缠 试验， 并与本研究方法相同区域的解缠结果进行对 比， 本研究方法与其他5种方法的解缠结果如图5所 示。 由图5分析可知 各方法的相位解缠结果十分相 似， 解缠相位数值范围差异很小， 其中Delaunay最小 2019年第11期王升等 基于克里金插值和Chan-Vese模型边界探测的InSAR相位解缠方法 145 ChaoXing 由表1可知 本研究方法与其他5种方法解缠结 果之间的差异均值都很小， 说明解缠结果具有高度一 致性， 其中Delaunay最小费流量法的差异均值最大， 达到0.008 rad。以最小费流量法为评价标准， 差异均 值仅有0.007 rad， 可见在相干性较高的情况下， 本研 究方法能够得到较高精度的解缠结果。 2. 3相干性较弱区域 2. 3. 1干涉条纹较稀疏的情况 为进一步验证本研究提出的相位解缠算法性 能， 选取相干性较弱的区域作为试验区， 该区域中的 干涉条纹较稀疏， 且存在高大山峰， 阻隔了山峰附近 山谷区域的回波信号， 形成雷达阴影。图6 （a） 所示 的相干系数图的中间位置呈现出形状近似字母 “Y” 的深灰色区域， 显示出该区域存在大面积低相干的 雷达阴影区域， 在干涉相位图 （图6 （b） ） 中与其对应 的低相干区内含有较多噪声， 干涉条纹边缘出现了 多处断裂现象。 使用本研究方法对缠绕的干涉相位图进行相位 解缠， 首先对低相干区域的干涉相位进行克里金插 值计算， 相干系数阈值设为0.6， 使用高斯模型拟合距 离-半变异函数， 经过克里金插值处理后重新生成的 干涉图如图7所示， 由图7可以看出 首先恢复出了 干涉条纹， 并且在一定程度上减少了干涉相位噪声； 然后使用Chan-Vese模型提取干涉条纹块； 最后根据 费流量法解缠结果的数值范围为-1.894～7.792 rad， 与其他方法解缠相位的数值范围只有0.001 rad左右 的差异。解缠结果之间的差异主要集中于右上角与 左下角的较低相干区域， 是因为在该类区域内的干 涉条纹边界较模糊。 为定量描述本研究方法与其他5种方法解缠结 果之间的差异程度， 采用差异均值来衡量， 差异计算 公式为 ［20］ k 1 mn∑ i 1 m ∑ j 1 n ||Ti,j- Ui,j，（8） 其中，m、n分别为图像的行数与列数；T为本研究方 法的解缠结果图矩阵；U为其他方法的解缠结果图矩 阵；k为本研究方法与其他方法解缠结果之间的差异 均值， 计算结果如表1所示。 金属矿山2019年第11期总第521期 146 ChaoXing 邻接矩阵计算出每块干涉条纹块的跃变数， 得到的 解缠结果如图8 （a） 所示。为比较相位解缠效果， 使 用本研究相位解缠算法与其他常用的解缠算法进行 对比， 对相同区域分别采用Goldstein枝切法、 最小费 流量法、 Delaunay最小费流量法、 质量图指导法、 区域 生长法5种方法进行相位解缠， 结果如图8 （b） ～8 （f） 所示。 对比图8中对同一低相干区使用不同相位解缠 方法得到的结果， 可以看到， 在 “Y” 形低相干区域噪 声强烈干扰下， 解缠结果各不相同。Goldstein枝切法 的解缠结果 （图8 （b） ） 中出现一些块状的孤岛现象， 特别是低相干区域有多处不合理的低相位值区， 表 明该类区域解缠时发生错误。对于最小费流量法 （图8 （c） ） 与Delaunay最小费流量法 （图8 （d） ） 的解缠 结果， 由于设定的初始相位区域不同， 两种相位解缠 结果在相同区域相差一个周期， 除去这种差别， 两者 解缠结果基本相同。最小费流量法与Delaunay最小 费流量法的解缠结果同样存在低相干区域解缠错误 的问题， 在图8 （c） 与图8 （d） 右下侧都出现明显的断 裂现象， 与实际情况不相符。质量图指导法的解缠 结果 （图8 （e） ） 右上方有明显的明暗区分界线， 分界 线左右侧的相位值不连续， 说明解缠结果不正确。 区域生长法的相位解缠结果 （图8 （f） ） 右下方变化较 剧烈， 低相干区域的解缠相位值发生错误。相比较 而言， 本研究方法的解缠结果 （图8 （a） ） 的相位变化 范围最小， 不同区域之间的变化较为平缓， 低相干区 域没有明显的断裂与突变现象。 通过相位转高程步骤将解缠后的相位值反演成 数字高程模型， 由本研究方法与Goldstein枝切法、 最 小费流量法、 Delaunay最小费流量法、 质量图指导法 和区域生长法5种方法的解缠结果反演的数字高程 模型的第418行剖面曲线如图9所示。 由图9分析可知 由不同方法解缠出的相位图反 演的数字高程模型形态各不相同。由Goldstein枝切 法解缠相位得到的数字高程模型剖面 （图 9 （b） ） 的 450～550列有不连续现象， 说明相位解缠结果出现了 孤岛效应。由最小费流量法和Delaunay最小费流量 法得到的数字高程模型剖面 （图9 （c） 和图9 （d） ） 450～ 550列的高程值均产生了无规则的剧烈变化。由质 量图指导法解缠相位得到的数字高程模型剖面 （图9 （e） ） 的450～550列的高程值忽然剧烈升高之后剧烈 下降， 表明由于噪声的影响， 反演出的高程值不稳定。 区域生长法相位解缠得到的数字高程模型的剖面图 （图9 （f） ） 形态明显有别于其他方法， 其最右侧峰顶处 的高程值偏高。通过本研究方法得到的解缠相位值 反演出的数字高程 （图9 （a） ） 坡度变化较缓， 雷达阴影 形成的低相干区域 （450～550列） 没有突变现象。 图9中的数字高程剖面图反映的是山谷地形， 山 谷两侧是高大的山峰， 其中右侧比左侧的峰顶高， 采 用空间分辨率为30 m的SRTM1数字高程模型验证不 同相位解缠方法的测量结果。SRTM1数字高程模型 显示右侧与左侧峰顶之间高差为111 m， 该数字高程 模型的高程精度为16 m。因此InSAR技术测量出 的两峰高差与SRTM1的高差之间的差异在理论上应 为-32～32 m， 如果超出这个范围， 可认为相位解缠发 生错误。本研究计算了不同相位解缠方法得到的数 字高程模型中的两峰高差， 并计算了与SRTM1两峰 高差之间差异值， 结果如表2所示。 分析表2可知 只有本研究方法和质量图引导法 2019年第11期王升等 基于克里金插值和Chan-Vese模型边界探测的InSAR相位解缠方法 147 ChaoXing 金属矿山2019年第11期总第521期 148 ChaoXing 的高差差异值处于-32～32 m范围内， 表明受到雷达 阴影区噪声的影响， 其他4种解缠方法在相位解缠过 程中产生了错误， 其中区域生长解缠方法得到的结 果与SRTM1的高差差异最大， 达到100.533 m， Gold- stein枝切法、 最小费流量法和Delaunay最小费流量法 的高差差异值也远远超出理论值范围。进一步分析 可知 虽然本研究方法和质量图指导法的高差差异 值均在合理范围内， 但是质量图指导法的相位解缠 结果存在多处突变， 相比之下， 由本研究方法得到的 解缠相位整体变化较平缓， 更符合实际情况， 因此， 在干涉条纹较稀疏的情况下， 本研究方法的相位解 缠效果比其他5种方法更优越。 2. 3. 2干涉条纹较密集的情况 为验证在干涉条纹较密集的情况下本研究提出 的相位解缠方法的有效性， 选取一幅只有地形相位 的干涉图 （图10 （b） ） 进行试验， 该干涉图的相干系数 图 （图10 （a） ） 显示， 干涉图上方的相干系数值较高， 下方的相干系数值很低。仔细观察相干系数图， 可 以发现干涉条纹较密集， 低相干区域横跨至少3个干 涉条纹块， 可明显看出低相干区域将上方的两条干 涉条纹边界割裂开。若需验证本研究方法在干涉条 纹较密集的情况下是否仍然有效， 就需要验证经过 克里金插值后， Chan-Vese模型是否能够正确提取出 干涉条纹块。 对相干系数值低于 0.6的干涉相位进行了克里 金插值计算， 结果如图11所示。从图11 （a） 中可以看 到， 由于干涉相位图大面积失相干， 导致断裂的干涉 条纹块在经过克里金插值后仍然没有有效连接起 来， 插值后的干涉条纹之间的边界也较模糊， Chan- Vese模型提取的干涉条纹块 （图11 （b） ） 也发生了错 误。由此看来， 在干涉条纹较密集的低相干区域本 研究方法失效， 对比上一组干涉条纹较稀疏的低相 干区域的相位解缠试验结果， 通过分析解缠失效原 因， 不难发现由于本组试验区域的低相干区跨越3个 以上的干涉条纹块， 导致克里金插值法无法完全恢 复出干涉条纹块之间的边界， 而上一组试验中虽然 低相干区域的面积较大， 但是干涉条纹较稀疏， 低相 干区域仅仅跨越2个干涉条纹块， 从而干涉条纹块之 间的边界也较易恢复， Chan-Vese模型能够正确提取 出干涉条纹块， 因此能够得到正确的相位解缠结果。 3结语 选用3组InSAR干涉图数据验证了本研究提出 的基于克里金插值和Chan-Vese模型的相位解缠方 法的有效性， 分析了其适用条件， 并且与Goldstein枝 切法、 最小费流量法、 Delaunay最小费流量法、 质量图 指导法、 区域生长法5种相位解缠方法进行了对比分 析。研究表明 在相干性较高的区域， 本研究方法与 其余5种方法解缠结果的差异均值较小， 说明解缠结 果差别很小， 并且本研究方法与最小费流量法的解 缠结果高度一致。在低相干区跨越2个干涉条纹的 情况下， 由本研究方法的解缠结果反演出的数字高 程模型与SRTM1 DEM之间的差异在合理区间， 解缠 结果仍然有效， 而在低相干区跨越2个以上干涉条纹 块的情况下， 本研究方法的相位解缠结果发生错误 并失效。综合试验结果， 认为本研究提出的相位解 2019年第11期王升等 基于克里金插值和Chan-Vese模型边界探测的InSAR相位解缠方法 149 ChaoXing ［1］ ［2］ ［3］ ［4］ ［5］ ［6］ ［7］ ［8］ ［9］ ［10］ ［11］ ［12］ ［13］ ［14］ ［15］ ［16］ ［17］ ［18］ ［19］ ［20］ 缠方法适用于高相干区或低相干区跨越不多于2个 干涉条纹的区域。然而， 该方法也存在不足， 如在使 用Chan-Vese模型提取干涉条纹块之前需要手动给 定大致轮廓， 自动化水平较低。尽管如此， 该方法在 解决特定InSAR失相干的情况下常规相位解缠方法 易发生解缠错误的问题方面具有一定的优势。 参 考 文 献 朱建军， 李志伟， 胡俊.InSAR变形监测方法与研究进展 ［J］ . 测绘学报， 2017， 46 （10） 1717-1733. Zhu Jianjun， Li Zhiwei， Hu Jun.Research progress and s of InSAR for deation monitoring［J］ . Acta Geodaetica et Carto- graphicaSinica， 2017， 46 （10） 1717-1733. 郝华东， 刘国林， 陈贤雷， 等.基于DEM的低相干区SAR干涉图 卡尔曼滤波相位解缠算法 ［J］ .国土资源遥感， 2013， 25 （1） 50- 55. Hao Huadong， Liu Guolin， Chen Xianlei， et al.Kalman filter phase unwrapping algorithm of SAR interferogram in low coherence re- gion based on DEM［J］ .Remote Sensing for Land Resources， 2013， 25 （1） 50-55. 李振叶， 陈星彤.InSAR数据处理中相位解缠算法综述 ［J］ .矿山 测量， 2014 （3） 78-80. Li Zhenye， Chen Xingtong. Overview of phase unwrapping algo- rithms in InSAR data processing ［J］ .Mine Surveying， 2014 （3） 78- 80. Zebker H A， Lu Y P.Phase unwrapping algorithms for radar interfer- ometryresidue-cut， least-squares， and synthesis algorithms ［J］ .Jour- nal of the Optical Society of America A， 1997， 5 （13） 586-598. 向茂生， 李树楷.用狄氏边界的泊松方程进行InSAR相位解缠 的研究 ［J］ .测绘学报， 1998， 27 （3） 204-210. Xiang Maosheng， Li Shukai.InSAR phase unwrapping using Pois- sons equation with Dirichlet boundary conditions ［J］ .Acta Geodaet- ica et Cartographica Sinica， 1998， 27 （3） 204-210. 许才军， 王华.InSAR相位解缠算法比较及误差分析 ［J］ .武汉 大学学报 信息科学版， 2004， 29 （1） 67-71. Xu Caijun， Wang Hua.Comparison of InSAR phase unwrapping al- gorithms and error analysis ［J］ .Geomatics and Ination Science of Wuhan University， 2004， 29 （1） 67-71. Goldstein R M， Zebker H A， Werner C L.Satellite radar interferome- trytwo-dimensional phase unwrapping ［J］ .Radio Science， 1988， 23 （4） 713-720. Xu W， Cumming L.A region growing algorithm for InSAR phase un- wrapping［J］ .IEEE Transaction on Geoscience and Remote Sens- ing， 1999， 37 （1） 124-134. 王世超.矿区大量级沉降InSAR相位解缠算法研究 ［D］ .青岛 山东科技大学， 2017. Wang Shichao.Research of InSAR Phase Unwrapping Algorithm at High Volume Subsidence Mining Area ［D］ .QingdaoShandong Uni- versity of Science and Technology， 2017. Chen C W， Zebker H A. Network approaches to two-dimensional phase unwrappingintractability and two new algorithmserratum ［J］ . Journal of the Optical Society of America A， 2000， 17 （3） 401-414. 韩松涛， 向茂生.一种干涉雷达阴影区的处理方法 ［J］ .电子测量 技术， 2008， 31 （3） 4-6. Han Songtao， Xiang Maosheng.Approach for shadow in int