摘要：郑南山,1,2, 张艳锁1,2, 李世金1,2, 杨化超1,2, 卞和方1,2, 张秋昭1,2, 张书毕1,2, 田雨1,2

本文内容来源于《测绘学报》2024年第10期（审图号GS京(2024)2165号）

基于相位质量融合估计与信息滤波的相位解缠方法

1.中国矿业大学自然资源部国土环境与灾害监测重点实验室，江苏 徐州 221116

2.

摘要

L波段差分干涉SAR卫星（陆探一号，LT-1）是我国首个L频段全极化民用SAR卫星星座，该组卫星将长期服务于自然资源监测监管，DEM生产是其主要任务之一。众所周知，相位解缠是影响DEM反演精度的关键因素之一，然而，高噪声及大梯度变化区域解缠结果精度低一直是相位解缠的瓶颈难题。对此，本文提出了一种相位质量融合估计策略引导的迭代信息滤波与解缠方法。首先，对干涉图像素的相干性、相位梯度、噪声分量和干涉条纹边缘进行估计与拟合，进而指导解缠路径布设，提高相位解缠结果的准确性；然后，采用中心差分信息滤波器实现解缠相位的准实时状态估计，进而提高解缠方法的稳健性；最后，通过TerraSAR-X/TanDEM-X双星及LT-1A/B双星SAR数据进行试验验证。结果表明，本文方法精度较统计费用流算法提高了约65.3%；较无迹信息滤波解缠算法提高了约13.6%。此外，试验结果不仅验证了本文方法的有效性，而且证明LT-1双星系统可达到高精度地形重建的预期效果，能为我国后续多云多雨区的DEM产品生产提供高精度的数据源。

关键词

基金项目

第一作者：高延东（1988—），男，博士，副教授，研究方向为单/多基线相位解缠、InSAR高精度DEM反演。E-mail：

通讯作者:郑南山 E-mail：znshcumt@cumt.edu.cn

本文引用格式

高延东, 郑南山, 张艳锁, 李世金, 杨化超, 卞和方, 张秋昭, 张书毕, 田雨.

GAO Yandong, ZHENG Nanshan, ZHANG Yansuo, LI Shijin, YANG Huachao, BIAN Hefang, ZHANG Qiuzhao, ZHANG Shubi, TIAN Yu.

阅读全文

InSAR技术已成为数字高程模型(DEM)产品获取的重要手段[1-2]。2022年我国发射了首个L波段差分干涉SAR卫星(陆探一号，LT-1)，该组卫星填补了我国差分干涉SAR卫星的空白，在轨运行后将服务于国土、地震、测绘、环境、减灾和林业等多个行业，其中，DEM产品获取是该组卫星的主要任务之一，LT-1为高精度DEM获取提供了发展机遇与挑战[3]。因此，充分挖掘LT-1地形重建能力，提高LT-1数据处理关键方法的可靠性，对高精度DEM获取具有极其重要的意义[4]。

相位解缠是影响InSAR技术DEM产品获取精度的关键因素之一[5]。相位解缠实质是一个病态问题，为得到解缠结果的唯一解，通常需要假设任意相邻像素相位梯度积分小于π(即满足相位连续性假设)[6]，因此，常规单基线相位解缠方法在大梯度变化区域难以得到理想的结果。此外，噪声也是影响相位解缠结果精度的主要因素之一，会增加相位解缠难度，甚至会导致相位解缠失败，进而影响最终DEM结果的质量[7]。现阶段常用相位解缠方法主要分为两类：一类是基于路径跟踪的局部最优相位解缠方法，主要以枝切法[8]和质量图引导算法[9-10]为代表；另一类是基于最小范数的全局最优相位解缠方法，主要以图割法[11]、最小费用流算法(minimum-cost flow algorithm，MCF)[12]和统计费用流相位解缠算法(statistical-cost network-flow algorithm for phase unwrapping, SNAPHU)[13]为代表。这些方法已经在实际工程中得到了广泛应用，然而在高噪声及大梯度变化区域仍会出现明显的解缠误差。鉴于此，国内外学者相继提出了一系列高效高精度相位解缠方法。

Kalman滤波相位解缠被认为是一种噪声稳健性较好的解缠算法，由文献[14]首次提出，同时给出了该方法的梯度估计模型[15]，随后将扩展Kalman滤波模型引入解缠领域中，在高噪声区域得到了较为理想的结果[16]。然而，扩展Kalman滤波模型在处理线性化问题时会损失高阶信息，故在低质量相位区域该方法并不能得到满意的结果；对此，文献[17]将无迹Kalman滤波引入解缠技术中，将其与路径跟踪策略和极大似然梯度估计相结合，进一步改善了Kalman滤波解缠算法的效率与精度；但该类方法模型的准确性及解缠路径的可靠性仍需要进一步研究[18]。大梯度变化区域解缠精度低一直是研究的热点问题，多基线相位解缠被认为是解决大梯度变化区域解缠的有效方法之一。文献[19]提出一种极具代表性的两阶段规划多基线解缠模型，该模型可以将单基线解缠方法的优点融入多基线解缠中，极大提高了大梯度变化区域的DEM结果可靠性。然而多基线相位解缠方法仍然存在噪声稳健性差、数据源缺失及基线配比难等瓶颈，严重限制了该类方法的应用范围[20-21]。因此，如何通过单干涉对解决大梯度变化引起的解缠问题成为研究热点。文献[22]提出一种距离向分频干涉测量辅助解缠算法，利用干涉相位减去分频干涉技术模拟相位获得残余相位，进而降低了常规干涉图解缠难度，有效减少解缠误差。文献[23]提出利用地理编码干涉图及其相位噪声相干估计来创建解缠路径，并通过残差检测出解缠误差，在每个孤立分量上增加±2kπ来校正解缠相位。文献[24]提出基于间断自适应马尔可夫随机场模型的相位解缠算法，建立的能量函数能够在不连续相位处阻断像素连接，有利于大梯度变化场景的地形重建。虽然这些研究从外部信息辅助、多技术融合、函数模型优化等不同角度提高了相位解缠结果精度，但同时面对地表复杂环境、高噪声及大梯度变化仍然无法获得理想的解缠结果。近年来，PhaseNet[25]、PGNet[26]和GAUNet[10]等半网络半模型的解缠方法相继被提出，这些方法利用优化的深度学习网络可以直接从缠绕相位中预测出梯度模糊数，从而不再依赖相位连续性假设，比较适合大梯度变化区域的相位梯度计算，但网络模型的泛化性能差仍然限制了该类方法的应用范围。

综上所述，本文基于Kalman滤波相位解缠技术，提出一种基于相位质量融合估计与信息滤波的相位解缠方法。首先，对干涉图像素的相干性、相位梯度、噪声分量和干涉条纹边缘进行估计与拟合，提出一种更为准确的解缠路径策略，进而指导解缠路径布设，提高相位解缠结果的精度；然后，采用中心差分信息滤波器实现解缠相位的状态估计，提高解缠模型的稳健性；最终，获取高精度DEM产品。经过TerraSAR-X/TanDEM-X双星及LT-1 A/B双星SAR数据试验，结果证明本文方法可以获得较已有方法更高精度的结果，可以为LT-1卫星高精度DEM生产提供可靠的技术方案。

质量图是影响解缠结果精度的主要因素之一[27]，相位梯度的二范数信息通常能够相对准确地标定像素质量。单一质量图反映的仅是一类干涉图信息，而影响相位质量的因素是多方面的，因此有必要融合多种相位质量估计方式。本文充分挖掘相干性、相位梯度、噪声分布及干涉条纹边缘等一系列相位质量的评价方式，优选相干系数、相位导数标准差、罗伯茨梯度幅值标准差、缠绕相位的统计标准差和噪声分布进行融合。其中相位导数标准差与罗伯茨梯度幅值标准差均反映的是相位梯度二范数信息，主导了质量融合估计结果；缠绕相位统计标准差反映的是相位二范数信息，能够表达干涉条纹边缘特征；噪声分布反映的是相位一范数信息，辅助识别高噪声像素。各种质量评价方式之间具有一定差异性，冗余性低，融合结果能够更加全面综合地表征像素间的质量相对关系。相干系数的计算公式为

(6)

式中，Qm，n表示干涉图(m，n)位置像素根据多个相位质量估计值拟合的三次函数曲线。在曲线中选择与相位导数标准差拥有同样横向坐标的点作为最终相位质量估计结果。

相位质量融合估计过程可简化为图1。由图1可知，本文方法可以兼顾不同质量图的特点，因此基于相位质量融合估计结果进行整体排序，可创建出更加可靠的解缠路径，避免或延迟错误的梯度估值过早出现在解缠路径上；相较于传统质量图引导的相位解缠算法，本文提出的融合估计质量图增加了对大梯度变化区域边缘的识别能力，解缠路径更加合理，泛化性能更强。

图1 三阶函数拟合的相位质量融合估计

Fig.1 Phase quality fusion estimation based on fitted third-order function

信息滤波是Kalman滤波的一种信息形式，利用信息矩阵和信息向量来描述高斯分布，状态更新过程计算复杂度小，因此无迹信息滤波相位解缠(unscented information filter phase unwrapping，UIFPU)方法较Kalman滤波解缠方法效率更高[28]。假设解缠相位变量ψ的均值为

，协方差为

，则基于无迹变换的采样点计算为

(7)

式中，λ=α2(κ+1)表示比例因子，用来降低总的预测误差，其中α和κ控制采样点散布范围。将采样点集通过状态方程进行传递，计算待解缠像素的解缠相位一步预测均值和预测误差协方差。均值与协方差的权重为

(8)

式中，β用于引入解缠相位变量概率分布的高阶矩信息。容易看出，无迹变换的核心在于选择一组确定的采样点作非线性变换，采样点值及其权重是由一组(α，β，κ)参数确定。然而在实际数据处理中，使用特殊参数值甚至可以得到更精确的解缠相位，其中α的设定对无迹Kalman滤波解缠影响较为显著。这在一定程度上也说明推荐的参数值能够保证解缠算法稳定性，但并不一定是最优数值。

中心差分变换思想核心是利用一阶、二阶差分算子代替非线性函数泰勒展开式的一阶项与二阶项。中心差分信息滤波器相比于无迹信息滤波器具有稍高的理论精度，其采样点的计算为

(10)

通过对比能够看出，中心差分变换过程更加简单高效，且只受一个中心差分步长参数调控。当确定采样点后，可根据信息滤波器进行解缠相位的预测与更新，增强解缠模型的精确性与稳健性，这一过程的伪代码如下。

Algorithm 迭代信息滤波与解缠

Input：解缠起点相位ψ(1)，相位梯度估计误差方差Q，量测相位误差方差R，观测方程Φ(σ)

Output: ψ(i)，Pψ(i)

初始化均方根误差：P1=Q(1)

For i∈[2，M×N] do //M与N分别表示干涉图的行数和列数

； //解缠相位的一步预测值

； //解缠相位预测误差方差

； //优化的解缠相位预测误差协方差，u为调节因子

； //解缠相位预测信息矩阵

； //解缠相位预测信息量

； //状态变量的量测预测值

； //状态变量与其量测预测值之间的互协方差

； //后验估计信息矩阵

； //解缠相位

； //解缠相位的估计误差方差

End for

综上所述，本文提出一种基于相位质量融合估计与信息滤波的相位解缠方法，首先使用相位质量融合估计策略对待解缠干涉图每一像素进行质量标定，然后采用文献[17]所述局部频率估计算法获取待解缠干涉图的相位梯度结果，随后对相位质量融合结果依据最大堆进行排序，确定解缠路径的布设。在路径迭代解缠过程中，采用中心差分信息滤波器对解缠相位进行实时估计与修正，最终获得高精度解缠相位。

为了验证本文方法的性能，分别使用陕西省渭南市TerraSAR-X/TanDEM-X及山西省晋城市LT-1数据进行试验，并将本文方法与MCF、SNAPHU及高精度UIFPU算法进行试验对比，为了验证不同方法的解缠精度，文中采用SRTM3对不同解缠方法获取的DEM进行精度评价。值得注意的是，虽然SRTM3的精度为90 m左右，但是SRTM3可以验证不同解缠方法在大梯度变化区域的解缠斑块现象，并可以通过统计解缠误差的均方根误差(root mean squared error，RMSE)，评价不同解缠方法结果的平均误差分布情况。

试验数据及结果如图2所示。图2(a)是该试验区的干涉相位图，图2(b)是根据SRTM3 DEM模拟的解缠相位，图2(c)是依据相位质量融合估计方法计算得到的质量图。由图2(a)可知，该区域条纹较为密集，该区域存在明显的大梯度变化，加大了相位解缠难度。图2(d)—(k)为不同方法的解缠结果，可以看出由于地形大梯度变化，MCF方法产生了明显的解缠误差，其解缠误差绝对值大部分在20 rad以上；SNAPHU方法的解缠结果与误差分布与MCF类似，也产生了严重的解缠误差，严重影响最终DEM产品的反演精度；UIFPU方法的解缠结果得到了显著改善，没有出现大范围误差传递现象，然而在局部区域仍然存在明显的解缠误差；相较而言，本文方法获得了最为理想的解缠结果。对比图2(j)与(k)可发现本文方法显著减少了两处局部解缠误差，这是由于相位质量融合估计策略能够引导出更优的解缠路径。为进一步验证不同方法结果精度，本文将不同解缠结果获取了DEM，并进行了精度分析，结果见表1，本文方法的高程RMSE为6.617 3 m，在4种方法中精度最高，其精度相较于MCF提升了约81.9%，较SNAPHU提升了约82.0%，较UIFPU提升了约26.9%，显著提高了大梯度地形变化区域InSAR相位解缠的可靠性。

图2 基于TerraSAR-X/TanDEM-X双星干涉数据的不同解缠方法结果及解缠误差分布

Fig.2 Unwrapped phase and error distribution of different PU approaches based on TerraSAR-X/TanDEM-X bistatic interferometric data

表1基于TerraSAR-X/TanDEM-X双星干涉数据的不同相位解缠方法性能对比

Tab.1

新窗口打开| 下载CSV

为了进一步验证本文方法的性能，选取山西省晋城市LT-1数据进行试验分析，图3(a)为试验区域及其干涉相位图。图3(b)是该区域的先验SRTM3 DEM，由其可获得图3(c)所示的模拟解缠相位。图3(d)是依据本文方法计算得到的质量图，可以看出该区域干涉像素相位质量普遍较高，证实了LT-1双星系统具有较强的干涉能力。使用枝切法、SNAPHU、UIFPU和本文方法对该干涉图进行解缠。各种方法的解缠相位、解缠误差分布和相高转换结果如图4所示，解缠相位转换为高程后统计的RMSE见表2。

图3 试验区域及其对应的模拟解缠相位和相位质量融合结果

Fig.3 The SRTM3 DEM of the experimental area and its corresponding simulated unwrapped phase and fused phase quality results

图4 基于LT-1A/B双星干涉数据的不同解缠方法结果及解缠误差分布和相高转换结果

Fig.4 Unwrapped phase, error distribution and elevation results of different PU approaches based on LT-1 bistatic interferometric data

表2基于LT-1 A/B双星干涉数据的不同相位解缠方法性能对比

Tab.2

新窗口打开| 下载CSV

对比解缠相位相高转换后的高程结果可以看出，4种方法整体上较完整地恢复出了地形特征，由于所用雷达影像分辨率较高，所以基于LT-1双星干涉数据反演的地形细节较SRTM3 DEM更加丰富。得益于LT-1双星准确的量测相位以及系统噪声少的优势，枝切法并未出现大量解缠空洞。通过解缠误差分布图可以看出，枝切法与SNAPHU方法均存在相对较多的局部解缠偏差，而UIFPU与本文方法的局部解缠误差相对较少，且本文方法比UIFPU在左上角区域少了一处明显的解缠误差斑块，差异的本质原因是本文提出相位质量融合估计方法能够延缓解缠梯度估计误差大的像素，避免造成解缠误差过早传递。由解缠相位转换为高程后的统计均方根误差可得，本文方法具有最优的地形反演表现，高程精度较SNAPHU提高了约65.3%，较UIFPU方法提高了约13.6%。通常MCF与SNAPHU的解缠性能相当，然而针对此幅地形干涉图，MCF以及图割法的解缠结果RMSE均在15 rad以上，基本可以认定为解缠失败，因此并未列举对比。为了突出本文方法在面对不同地形、不同噪声情况下模型的稳健性，基于模拟的SAR数据对本文方法模型的稳健性进行了测试，同时为了表明本文方法对于解缠路径布设更为合理，在传统UIFPU方法基础上采用不同相位质量估计策略进行引导处理，得到的已解缠像素绝对相位均方根误差随解缠路径走向变化如图5所示。

图5 解缠结果均方根误差随解缠路径走向变化

Fig.5 The mean square error of the PU results vary with the unwrapping path

对于路径跟踪解缠方法而言，相位梯度的二范数通常较相位梯度一范数、相位一范数及相位二范数等相位质量评价方式能够引导出更为理想的解缠结果，原因在于其能够衡量相位梯度波动变化，而相位梯度又是质量图引导算法的最为关键因素。因此，主要统计了3种含有相位梯度二范数信息的质量估计策略引导的解缠路径布设。由图5可以看出，相位导数标准差引导的路径在解缠到第2.5万与第5万左右像素时出现了两处显著的局部误差，而罗伯茨梯度幅值标准差与相位质量融合估计方法均可避免此两处误差传递。虽然罗伯茨梯度幅值标准差引导的路径在解缠到第9万左右像素时出现了明显的局部误差，但最终解缠精度仍略优于相位导数标准差。由于相位质量融合估计策略能够延迟甚至避免一些解缠误差出现，因此解缠效果最为理想，实质是因为相位质量融合估计策略在相位导数标准差的基础上还考虑到了相干性、噪声以及条纹边缘等信息，更为准确地标记出低质量像素，提高了局部乃至整体像素解缠的可靠性。

为了进一步示范LT-1 A/B双星干涉数据及本文方法在地形重建方面的能力，选取数据尺寸为3000×3000像素的较大研究区域进行试验分析，该区域实际面积约为523 km2。通过上述分析，该部分本文仅对UIFPU与本文方法两种高精度算法进行处理，最终得到的高程结果及其误差分布如图6所示。采用SRTM3 DEM对数据处理结果进行精度分析。UIFPU方法得到的高程均方根误差为12.904 5 m，本文方法则为9.221 0 m。针对较大规模的LT-1 A/B双星干涉数据，信息滤波类解缠方法仍可以实现较精准的地形重建。然而，随着地形干涉数据规模增大，梯度变化、噪声等因素会降低相位解缠性能，造成解缠误差局部传递，如图6(c)中凹凸的误差斑块。与之相比，本文方法通过更加全面准确的像素质量估计以优化解缠路径，最终尽可能避免了解缠误差局部传递，提升了LT-1双星系统的地形重建效果。随着陆探一号01组卫星的正式投入使用，必将更有效支撑我国多云多雨多植被地区高精度L波段大规模测绘，本文方法可以为LT-1 SAR卫星数据处理提供有效的技术支撑。

图6 基于大范围LT-1A/B双星干涉数据应用不同解缠方法得到的高程结果及高程误差分布

Fig.6 Elevation results and error distribution based on large-scale LT-1A/B bistatic interferometric data using different PU approaches

本文提出了一种相位质量融合估计策略引导的迭代中心差分信息滤波与解缠方法，通过优化解缠路径布设方法实现减少或延迟解缠误差局部累积；在质量排序解缠过程中，采用高效的中心差分信息滤波实现解缠相位的实时估计与修正。基于TerraSAR-X/TanDEM-X双星干涉数据的试验结果表明，本文方法的解缠结果相高转换后的精度较MCF提升了约81.9%，较UIFPU提升了约26.9%，尤其在大梯度地形区域能够得到更接近真实场景的绝对相位。基于LT-1 A/B双星干涉数据的试验结果表明，本文方法的解缠结果相高转换后的精度较SNAPHU提升了约65.3%，较UIFPU提升了约13.6%，得到的DEM精度远优于SRTM3 DEM，在大梯度变化区域能够获取到更可靠的高程观测结果，同时也表明LT-1双星系统在地形测绘领域具有较好的应用成效。LT-1的高时空分辨率使它反演的地形细节更为丰富，因此，陆探双星系统能够为我国高精度地形数据的生产与更新提供支撑。

来源：测绘学报