基于超像素分割与多信息融合的叠掩和阴影区域检测法

刘 鹏，李真芳，楼良盛，杨伟明，王 震

1. 西安电子科技大学雷达信号处理国家重点实验室，陕西 西安 710071； 2. 地理空间信息国家重点实验室，陕西 西安 710054； 3. 西安测绘研究所，陕西 西安 710054

干涉合成孔径雷达(interferometric synthetic aperture radar，InSAR)技术得益于其获取目标高程信息的能力，被广泛应用于地形测绘、海洋监测、灾害监测与预警等领域[1-4]。由于合成孔径雷达属于主动辐射测距成像模式的一种技术，使得叠掩和阴影成为SAR图像中普遍存在的现象。在InSAR信号处理中，叠掩和阴影区域的相位奇异性导致相位无法正确滤波和相位解缠，错误的相位信息导致相邻区域的相位展开错误，并造成高程反演出现错误。因此需要准确判断并标记SAR图像的阴影区和叠掩区，避免该区域影响整个相位解缠过程。

近年来，许多学者对InSAR中阴影和叠掩的检测方法进行了相关研究，并取得了一定成果，如基于幅度的阈值分割方法[5]、基于相干系数的阈值分割方法[6]、基于边缘特征和模糊理论的方法[7]、基于局部频率估计的方法[8]、基于特征值分解的方法[9]和基于恒虚警(CFAR)检测的方法[10]等。其中，局部频率估计和特征值分解更为有效，但其都有各自的局限性。一方面，上述方法都是逐像素进行判别，SAR图像实际上反映的是场景对电磁波的双站散射能力，目标的几何和物理特性通过影响场景电磁双站散射而被SAR系统捕获并成像，SAR图像上每个分辨单元是由完全随机分布的许多散射点矢量叠加组成，在SAR图像上形成相干斑噪声，进而使得简单阈值分割方法失效，逐像素的判别结果较为离散，即使做形态学处理，也很难反映SAR数据的纹理信息；另一方面，局部频率估计需要对每个像素周围指定窗口进行二维FFT操作，对SAR数据逐像素进行判别计算量较大，特征值分解方法需要获取同一区域的多幅SAR图像，通过构造图像阵列的方式判别叠掩和阴影，无法应用于大部分工程需求，恒虚警检测(constant false alarm，CFAR)基于目标与背景对比强烈的准则，需要估计目标的分布函数，对于较大的叠掩和阴影区域检测性能不佳。

针对逐像素检测方法在准确度和速度上的不足，考虑到在SAR图像中属于同一地物同等坡度的相邻像素具有相似信息，叠掩或阴影区域具有散射信息近似相等且距离上相近的特点，如果能够将具有相邻且相似信息的像素组合成一个基元块，则可以大大降低后续图像处理任务的复杂性。基于此，本文结合超像素分割方法提出一种基于多信息融合的超像素检测算法，实现了快速且准确的叠掩和阴影区域检测。近年来，常用的超像素分割方法是基于图论和梯度下降聚类的方法。基于图论的方法具体包括基于随机游走的熵率的方法[11]、伪布尔优化的方法[12]和基于图的贪心聚类算法[13]。基于梯度下降聚类的超像素分割包括简单线性迭代聚类(SLIC)方法[14]、迭代边缘细化方法[15]、多尺度分割方法[16]和分水岭算法[17]。这些方法由于其处理速度快、边界拟合度较好，被广泛应用于光学图像分割中。其中，SLIC方法由于在遵守图像边界的能力、速度、内存效率以及分割性能上的优越性而被广泛使用[18]。然而，由于SAR图像固有的相干斑噪声效应[19]，以及相对于光学图像而言SAR图像具有更多像素数量，通常达到106～108数量级，使得传统SLIC在超像素分割方面表现不佳。近年来多尺度分割和SLIC方法在SAR图像分割领域得到了广泛的应用。尽管过去已提出了许多的多尺度分割方法，但嵌入在eCognition软件中的多尺度分割(multiscale segmentation，MS)方法仍然是目前最有效的方法[20]。针对SAR图像存在的相干斑噪声问题以及SAR数据固有的数据特征，文献[21]基于概率密度函数的SLIC(PDF-based SLIC)方法，实现了L波段的SAR图像超像素分割。文献[22—23]通过改进SLIC超像素算法中的距离度量方法实现了PolSAR图像中目标的分类。然而SLIC算法应用于SAR图像分割中依旧有些问题没有很好地解决。首先，大部分研究都是集中在改进建立有效的距离度量以区分地物类别，并未考虑到初始化步骤对于SAR图像分割的重要性；其次，迭代过程中复杂的距离度量使得算法在大尺寸SAR图像中效率低，不利于实际工程应用。因此，需要考虑SAR图像中叠掩和阴影检测的实际情况，设计一种快速且准确的SLIC处理流程，以及将SLIC嵌入到后续处理中实现叠掩和阴影的快速检测。

本文基于改进后的SLIC算法，结合相干系数、阈值幅度分割方法和局部频率估计方法，提出一种基于多信息融合的超像素检测算法，识别SAR图像中的叠掩和阴影区域，来保证InSAR后续处理的可靠性和有效性。由于缺少真实场景的数字高程模型，导致无法验证本文算法在实际应用中的准确性，因此本文提出一种基于电磁散射模型的快速回波仿真法，将先进积分方程近似模型(advanced integral equation model，AIEM)[24-25]和频域回波仿真[26]相结合，通过利用卫星工具包(satellite tool kit，STK)软件设计仿真所需的轨道参数，实现星载全极化双站SAR的原始回波数据仿真，回波能够反映多站、多角度、全极化、复杂地物的回波特性以及反映SAR图像的叠掩、阴影、相干斑噪声等物理特征，可用于检验本文所提出的基于多信息融合的超像素检测算法有效性。结合天绘二号卫星[27-28]的2019年9月获取的河北赤城精度检测场区域和2020年6月获取的哥伦比亚区域数据，证明本文算法的有效性。

1 SAR图像叠掩和阴影检测理论研究

1.1 叠掩和阴影的几何特性

工程上利用两个在俯仰向观测角度略有差别的接收天线，对同一目标进行各自成像处理，通过干涉的方法获取目标高程信息，称其为干涉合成孔径雷达(interferometric synthetic aperture radar，InSAR)技术。如图1所示，S1和S2代表装载雷达天线的两部卫星，当雷达观测具有陡峭地形的山区或城市时，具有相同斜距的目标点P1和P2在录取回波时会叠加到同一个像素内，将这一现象称为叠掩，将电磁波无法照射的区域成为阴影区域。由图1中的几何关系可以得到局部入射角θ的具体表达式为

θ=β-α

(1)

式中，β为雷达下视角；α为目标点的坡度角。当θ=0°时，斜坡的坡面与雷达视线垂直，在该目标点两侧的点将具有相同的斜距，使得不同层次高程的目标被叠加到同一距离单元内形成叠掩；当θ<0°时，此时α>β，底部的目标点到雷达的斜距大于顶部的目标点到雷达的斜距，使得SAR图像上出现顶底倒置的现象，进而形成叠掩现象。阴影区域一部分是由于电磁波被其他目标遮挡而形成，另一部分是由于θ>π/2，即α<β-π/2时电磁波镜面反射导致。

图1 几何关系Fig.1 Geometric relationship

1.2 基于信号特征的传统逐像素检测法

1.2.1 基于幅度阈值分割方法

1.2.2 基于局部频率估计方法

文献[8]指出InSAR系统形成的图像之间存在频率偏移Δf，并给出了具体表达式

(2)

式中，f0为载波的频率(单位Hz)；R为雷达到目标点的斜距(单位m)；k为常数，对于双星收发分置雷达系统取2；B为基线长度(单位m)；γ为基线倾角(单位rad)。由式(2)可以看出，当β-α<0或β-α>π/2时Δf取值为负，此时对应叠掩和阴影区域。

本文采用文献[8]中频谱分析方法估计局部频率，在指定窗口内做二维FFT求干涉相位图频谱，频谱峰值出现在二维局部频率对应的位置，通过判断频率符号确定可能的叠掩和阴影区域。

1.2.3 基于CFAR检测方法

均值类CFAR算法是目前主流的CFAR检测算法，其实现思路是以被检测单元为中心的周围单元作为参考窗，将窗口内的所有数据平均作为背景功率估值，将估计的背景功率与系数(门限因子)相乘得到决策阈值T，根据不同的背景功率和门限因子。文献[10]采用SO-CFAR算法进行检测，将参考窗分为对称的两个窗口，以较小的背景杂波数据来估计背景功率，其背景功率为

(3)

(4)

式中，εso为门限因子。SO-CFAR算法适用于多个相邻的目标检测，目标周围虚警率低，但是杂波边缘区域虚警概率提升。

2 SAR图像叠掩和阴影检测算法实现

2.1 SLIC超像素分割算法的改进

SAR图像反映的是场景反射电磁波的能力，因此SAR图像是强度图像，传统的SLIC算法适用于RGB图像的分割，鉴于SAR图像像素数量通常达到106～108数量级，需要在保证分割精度的同时尽可能减少传统SLIC算法运算量以达到工程上快速分割的目的，本节对SLIC算法进行改进，并给出了改进后的SLIC具体实现流程。预处理时需要将SAR图像转换成取值范围为[0,1]的灰度图像，然后放大100倍作为图像的亮度信息，得到亮度图L。

2.1.1 初始化聚类中心

图2 初始化种子点Fig.2 Initializing the seed point

为避免将超像素放置在图像边缘或噪声点上影响聚类性能，本文将聚类中心移动到与八连通域中的最低梯度位置相对应的种子位置。梯度度量方法为

(5)

2.1.2 像素点的聚类

SLIC算法通过检测每个像素点周围2S×2S范围内聚类中心之间的相似度，将相似度最高的聚类中心的标签赋予该像素。SLIC算法定义了最大的空间距离Ns=S和最大的颜色差异Nc来对空间距离和颜色距离进行差异整合，由于Nc的值不容易直接定义，所以SLIC算法将Nc取值为常量m∈[1,40]，最终给出了像素间的相似性公式为

d(i,Ck)=

(6)

d(i,Ck)描述了像素i和聚类中心Ck之间的相似度。可以将式(6)表示为向量的形式，特征向量f的加权L2范数形式为

(7)

式中，对角矩阵W=diag(1,λ,λ)，λ=m2/S2；f=[l,x,y]T。则d(i,Ck)可写成

(8)

(9)

进一步通过近似的方法将d(i,Ck)近似到ds(i,Ck)减少了计算量，将聚类问题转化为argminds(i,Ck)。

2.1.3 迭代更新

初始化标签图像l={l(i)=-1|i∈[1,N]}，距离图像D={d(i)=∞|i∈[1,N]}，循环k∈[1,K]计算Ck邻域2S×2S范围内ds(i,Ck)，若ds(i,Ck)

2.1.4 后处理增强连通性

经过上述迭代优化可能出现多连通情况、超像素尺寸过小以及出现孤立超像素等，这些情况可以通过增强连通性解决。本文具体处理步骤如下：

(1) 为确保不相邻区域的标记是不同的，将具有不同区域但相同标签的重新分配标签。

(2) 对标签图像进行形态学开操作处理，去除孤立的小点和毛刺。

(3) 上述操作完成后，某些区域可能会被关联，分成多个区域或被其他区域吸收，为确保不相邻区域的标记是不同的，需要对区域重新编号，以便他们从1开始顺序增加。

上述步骤可能导致超像素个数增多或减少，最终输出标签图像l={l(i)|i∈[1,N]}，其中l(i)表示SAR图像上每个像素对应的超像素编号，对应的超像素特征结构可表示为

(10)

2.2 基于多信息融合的超像素检测算法

基于改进后的SLIC超像素分割算法，结合InSAR相干系数、传统阈值分割和局部频率估计方法，本文提出了一种叠掩阴影判别算法。具体流程如图3所示。

图3 算法流程Fig.3 Algorithm flowchart

对于给定输入为R×A=N的SAR图像。首先，采用均值滤波的方式去除SAR图像中的斑点噪声；然后，对SAR图像中的轮廓特征进行精细化分割，获得标签图像l={l(i)|i∈[1,N]}和超像素特征结构SPCk,k∈[1,K′]参数。初始化叠掩、阴影和疑似叠掩阴影的标签图像Il、Is和Ils。计算SAR图像平均强度Lavr，遍历k∈[1,K′]，对每个超像素做局部频率估计得到叠掩和部分阴影区的区域集合Ils，通过强度和相干系数阈值分割法在区域Ils中识别叠掩区Il，则可得到部分阴影区Is=Ils-Il。估计出部分阴影区的平均强度avr(L(Is))，将2×avr(L(Is))和相干系数小于0.6作为最终判定阴影区的阈值进行阴影区域的检测。阴影区域通常没有回波能量，因此相干性较差，根据处理过程中的经验，一般认为相干系数小于0.6的区域包含阴影区域。在相干系数差的区域中通过设置阈值提取阴影区域。阈值分割得到的阴影区域中同样包含回波能量弱的区域，这是由于SAR图像每个像素值反映的是电磁波与目标点作用后被接收机接收的回波能量值，鉴于散射系数受到入射角影响较大，对于入射角大于70°的面目标而言，回波能量相对较弱，使得在数据分布特性方面与阴影区域相似，造成这部分区域相干系数较差。对于高程反演而言，相干系数较差的区域也会影响相位生成精度，因此回波能量弱的区域也应该被有效检测。

2.3 基于回波模拟技术的算法有效性检验

为验证本文算法的正确性，需要明确叠掩和阴影的实际位置，由1.1节的分析可知，对于给定的观测场景数字高程模型(DEM)，通过坡度和下视角的关系可判断叠掩区域，通过射线追踪的方法可判断阴影区域。因此通过准确地模拟回波录取过程以获取真实场景的SAR数据，可以检验本文算法的有效性。SAR回波过程可被描述为发射信号与系统冲激响应的卷积[29]，将其扩展到双基，具体可表达为

SRx(τ,f,t)=STx(τ,f)⊗τhc(τ,t)+n0(τ,f,t)

(11)

式中，τ为快时间参数(单位s)；t为慢时间参数(单位s)；f为频率参数(单位Hz)；接收信号SRx(τ,f,t)被描述为式(11)中STx(τ,f)与通道脉冲响应函数hc(τ,t)的卷积以及0均值高斯白噪声n0(τ,f,t)的形式。在快时间的相干窗口内假设发射机位置PosTx(t)和接收机位置PosRx(t)相对时间是不变的。此时目标点的通道脉冲响应可被写为

(12)

式中，δ(τ-R(t,r)/c)为双基距离延时的δ函数；σ(t,r)为在采样时刻t对应的WGS-84坐标系位置r=(x，y，z)处目标点的反射率函数，反射率函数是由观测几何角、目标点材质等特性决定的；双基距离函数R(t,r)为

R(t,r)=|PosTx(t)-r|+|PosRx(t)-r|

(13)

在每个方位向慢时间采样时刻，4 dB波束宽度内对应的照射区域内的脉冲响应hc(τ,t)可被描述为

(14)

式中，W(t,r)是由雷达发射功率、收发天线增益以及雷达工作过程中由热损失等引起的系统损耗等因素共同决定的参数；St为慢时刻t时4 dB波束照射到的场景目标集合，去除其中的阴影区域。根据双站SAR雷达方程[29]可得到W(t,r)的具体表达为

(15)

式中，Pt为发射信号的功率；Gt(t,r)为慢时间t接收位置r处目标的发射天线方向图增益；Gr(t,r)为接收天线方向图增益。

SAR成像过程实质上是重构目标点反射率函数的过程，因此SAR数据仿真的准确性取决于反射率函数计算的准确性。文献[24]将AIEM模型与实测数据进行对比，证明了AIEM模型在计算随机粗糙面双站散射系数的可靠性。AIEM模型考虑基尔霍夫项、交叉项和补偿项，可以简化写成最终的简单形式，散射系数计算公式为

(16)

(17)

式中，K=sqrt((ksx-kx)+(ksy-ky)2)；ρ为相关长度(单位m)。

在AIEM模型的推导坐标系下φ=π，因此只要给出δ、ρ、θ、θs、φs这5个参数以及介电常数即可准确计算目标点的散射系数。通常对于目标场景采用三角面元建模，每个场景面元都有一个确定的观测几何角θ、θs、φs，若用A表示三角面元的面积，则每个面元接收电磁波的有效面积为Acosθ。由于SAR图像存在相干斑噪声，因此模拟的SAR图像也必须具有斑点噪声特性。文献[19]指出斑点噪声被建模为实部和虚部服从0均值高斯分布的乘性噪声。因此最终的反射率函数表达为

(18)

式中，N(0,1)表示服从均值为0方差为1的高斯分布特性的数。

根据上述建模分析，最终的双站SAR回波方程可被描述为

δ(τ-R(t,r)/c)dr+n0(τ,f,t)

(19)

3 试验结果及分析

3.1 改进SLIC算法性能分析

一方面，通过将叠掩和阴影区域检测之后的结果作为掩膜图屏蔽叠掩和阴影区域进行干涉相位解缠，可以大大缩短解缠时间并提高解缠精度。对检测出的叠掩和阴影区域进行标注，以便后续通过升降轨融合方法或已有的DEM对叠掩和阴影区域进行修正，生成高精度的DEM产品。另一方面，TH-2致力于构建全球DEM库，因此必须尽可能提高叠掩和阴影的检测效率，而超像素分割过程作为叠掩和阴影检测的重要步骤，需要提高其运算效率和分割性能。为了证明本文提出的改进SLIC算法的有效性，对TH-2获取的SAR数据集进行了性能和运算效率分析，并将本文算法与目前已被开发能够用于SRA图像分割的方法(传统SLIC[18]、MR方法[20]及PDF-based SLIC方法[21])进行对比。其中，MR方法是集成在商业软件eCognition v8.9中的方法，在该版本中MR方法并未使用多核并行计算技术进行加速。为对比4种算法的运算效率，其余方法均采用C++编程语言实现且不进行并行加速。

为了评估不同分割方法的性能，本文采用3个指标来评估分割结果，分别为边界契合度(boundary fit，BF)、平均变异系数(mean coefficient of variation，MCV)和运算效率。

BF参数用来衡量超像素分割边界与真实边界的一致性，定义为

(20)

式中，ES∩T用来描述通过算法分割的边界和地物真实边界的重合像素数；ET表示地物真实边界的像素数。因此BF的值域为[0,1]之间，BF的值越接近于1表示超像素块与图像边缘的一致性更好。本文中BF参数选取局部分割区域进行统计求解。

通常不同分割算法分割出的超像素块之间差异较大，因此每个超像素块的平均值略有差异，采用方差难评估多种算法之间的性能，CV参数可以消除由于不同超像素块的平均水平不同对离散程度比较的影响。MCV参数定义为

(21)

如图4所示，试验数据选自TH-2号2020年6月获取的哥伦比亚区域实测数据，大小为18 214×23 566像素，截取其中大小为3000×3000像素的SAR图像作为输入。在MR方法中，尺度参数设置为50，形状参数设置为0.4，颜色参数设置为0.5。在SLIC方法中，初始超像素个数K取值均为10 000，颜色和空间重要性度量参数m取值15，其余所需参数均参考对应文献。由图4结果可以看出，SAR图像中通过分割得到的超像素区域中的颜色是均匀的，即属于同一类的像素在回波能量上具有相似的表现。从细节上可以看出不同算法在性能上有一定的差异，主要体现在边界的定位以及分割效率上。对图4中矩形的区域进行评估，由图4中黄色框和表1中BF和MCV参数结果可以看出MR方法和本文算法都能较好地分割出不同区域的边界，而PDF-based SLIC方法和传统SLIC方法分割性能相对较差。本次试验中采用MR方法进行分割用时18.28 s，本文算法用时1.99 s，结果表明本文算法在更短的时间内分割出的超像素更加均匀和紧凑，这为后续叠掩和阴影区域检测提供了重要的保障。

表1 性能比较

图4 TH-2数据超像素分割结果Fig.4 Superpixel segmentation results for the TH-2

在TH-2号2020年6月获取的哥伦比亚区域数据集中截取[500,1000,2000,…,10 000]2数据量大小的图像作为输入，在每次试验中超像素个数设置为图像行数的10倍，颜色和空间重要性度量参数m固定为15。具体结果如图5所示。

图5 运算效率Fig.5 Operating efficiency

随着输入数据量的增大，MR方法需要根据相异性准则对空间上非相邻的区域进行合并使得运算量增大，PDF-based SLIC需要估计像素周围的概率密度函数进行聚类使得计算量增加，由于本文算法在计算上的有效近似以及后处理的优化使得计算效率明显优于其他方法。

3.2 仿真验证

本文利用STK软件生成观测河北赤城精度检测场区域的轨道数据，通过构造500 m基线得到发射机和接收机对应的天线相位中心文件，其余具体仿真参数见表2。

表2 仿真参数

观测区域的DEM见图6(a)，本文采用西安电子科技大学InSAR课题组研发的RDSpace软件对仿真回波数据进行成像，本文给出发射星的成像结果见图6(b)，对比图6(d)成像后地理编码的结果和图6(a)可以证明本文回波仿真方法的有效性。

图6 回波仿真SAR成像结果Fig.6 SAR results of echo simulation

对于3500×3588像素大小的SAR图像，分别设置初始超像素个数K=[3,10,3]×10 000，颜色和空间重要性度量参数m=[15,15,5]。为方便呈现分割后的结果，本文选取图7(a)红色框出来的区域进行分析，图7(a)为红色框中的SAR图像，3种分割参数均实现了对SAR图像的分割，并很好地保证了SAR图像中目标的轮廓信息。对比图7(b)、(c)可以看出当m取值相同时，超像素个数越多分割的结果越细，分割完成后的超像素个数接近初始的超像素个数，初始超像素越多适应目标边界的能力越强。对比图7(b)、(d)可以看出当初始超像素取值相同时，m取值越小超像素形状和尺寸越不规则，但是分割出的超像素更加均匀和紧凑，边界粘附性更好。从图7(d)结果可以看出，生成的超像素个数比输入的超像素多，这是因为参数m取值越小颜色权值越大，此时算法更倾向于像素强度值上的邻近性，使得在初次分割时大量相同标签的区域并未接壤，在后续处理中为确保不相邻区域的标记是不同的，对具有不同区域但相同标签的重新分配标签，进而使得超像素个数增多。在叠掩和阴影检测时，更关注不同类别区域地物能量的变化，因此通常m取值较小，为保证结果的可靠性，一般m取值范围为1～40。

超像素形状和尺寸的规则性并不能提高叠掩和阴影检测正确性，对于叠掩和阴影的检测通常更关心目标能量上的变化。通常希望分割的结果能够更紧密地黏附到目标边界，因此设置初始超像素个数K=30 000，颜色和空间相对重要性度量指标m=5进行叠掩和阴影的检测。并将本文算法与文献[5—6，10]以及根据几何关系判断的真实叠掩和阴影区域进行对比验证，详细结果见图8。图8(a)、(e)是根据DEM和轨道信息几何关系准确判断的叠掩和阴影位置，由于SAR图像反映的是目标点的反射率函数被脉冲响应函数调制后的点散布函数形式，在该试验中SAR图像一个像素上的点是5 m×4 m范围内的DEM点的叠加，因此首先需要将DEM反投影到SAR图像上再进行插值，才能确保高程点与SAR像素对应。可以看出CFAR检测方法出现很多漏检区域和无效区域，叠掩的检测正确率不足50%。本文算法几乎检测出全部的叠掩区域，未检测的部分大多是离散的孤立点区域，这些区域通常只占1～2个像素在形态学处理时被去除，对后续InSAR处理影响可忽略。对于阴影区域而言，阴影区域多为噪声，回波能量较弱，因此阴影区域几乎都可以根据强度值检测出来，由图8(h)中可以看出，检测出的阴影区域比真实阴影区域较多，这是因为这些区域普遍集中在背坡面，具有较大的入射角，几乎没有回波信号，导致在图像上呈现与阴影区域相似的特征。

图7 SLIC超像素分割结果Fig.7 The results of SLIC superpixel segmentation

图8 不同方法检测结果Fig.8 The results of different methods

为验证检测结果的有效性，定义叠掩和阴影检测正确率为检测出的真实叠掩和阴影数占准确的叠掩阴影总数比例，定义叠掩检错率为检测错误的总像素占检测出的总像素的比例，为验证本文算法对DEM产品误差精度的提升效果，通过计算未掩膜处理时DEM生成的高程误差均值以及经过叠掩和阴影区域掩膜处理后的高程误差均值，并计算相位解缠速度提升指标和DEM产品误差精度提升指标。具体指标结果见表3。

表3 叠掩阴影检测正确率

由表3可以看出，CFAR检测方法对于复杂地形的SAR图像表现效果较差，相干系数联合幅度阈值检测的方法检测正确率相对于本文算法较低，具体是因为相干系数是由周围信息估计得到的，幅度阈值是按照像素点去检测的，所以实际的相干系数和像素点的幅度会有不一致性导致正确率降低，错检率增加，这一现象在边界处更为明显。通过将本文算法检测出的叠掩和阴影区域作为掩膜图，在相位解缠过程中进行掩膜使得相位解缠时间从16.64 s降低到13.13 s，相位解缠速度提升17.31%，生成的DEM产品高程误差均值从2.32 m降低到1.91 m，DEM产品误差精度提升17.58%，优于其他方法。图9给出了仿真数据InSAR处理后生成的DEM产品高程误差值，误差分析试验区域为图9中红色线段位置，图9(b)为未掩膜的DEM生成误差，图9(c)给出了掩膜后的DEM生成误差，DEM生成误差是通过InSAR处理得到的DEM产品和原始回波仿真输入的观测场景DEM做差得到。图9中对比红色方框中的结果可以明显看出，由于解缠错误导致误差较大的区域，通过被有效地掩膜使得高程反演误差降低，证明本文算法的有效性。

图9 DEM产品高程误差Fig.9 Elevation errors of DEM product

3.3 天绘二号数据试验

将本文所提算法应用在天绘二号卫星在2019年9月获取的河北赤城精度检测场SAR图像中，对于6046×6971像素大小的SAR图像，设置初始超像素个数K=30 000，颜色和空间重要性度量参数m=5。考虑到原始接收到的SAR数据中由于地物信息复杂存在强散射点数据，导致如果直接对SAR图像量化后显示，造成SAR图像显示较暗问题，若不进行处理会使得后续阴影检测失效，因此需要去除SAR数据中后5%强度值数据再进行后续试验。

如图10所示，图10(a)、(c)为本文算法检测出的叠掩区域以及将叠掩区域边界标记到SAR图像上结果，图10(b)、(d)为本文算法检测出的阴影区域以及将阴影区域边界标记到SAR图像上结果。为证明本文算法的有效性，选取图10中绿色框中的区域的SAR图像以及对应的相干系数图，分别采用本文算法、CFAR检测方法和相干系数联合幅度阈值检测方法进行试验。试验结果如图11所示，从中可以看出相干系数联合幅度阈值检测的方法检测出的叠掩区域较真实叠掩区域更大，无法适应叠掩的边界，检测出的阴影区域较小，且叠掩和阴影漏检严重。相干系数联合幅度阈值检测方法受到SAR图像固有的斑点噪声影响检测结果离散，检测性能较差。实际SAR图像相比于仿真的SAR图像而言斑点噪声更严重，使得相干系数联合幅度阈值检测方法失效。而本文算法检测的叠掩和阴影区域边界一致性较好，检测性能更优。

图10 河北赤城精度检测场SAR图像试验结果Fig.10 The results of the SAR Image of Hebei Chicheng calibration field

图11 不同方法检测结果Fig.11 The results of different methods

为进一步验证本文算法应用于实测数据的有效性，将检测出的叠掩和阴影区域在相位解缠过程中进行掩膜。用于计算高程误差的辅助DEM选自天绘中心测量得到的5 m精度的DEM数据，未掩膜时相位解缠用时23.77 s，DEM生成的高程误差均值为3.79 m，掩膜后相位解缠用时16.93 s，DEM生成的高程误差均值为2.74 m，相位解缠速度提升了28.78%，DEM产品精度提升了27.7%。选取图12中红色线段区域进行分析，其对应的经纬度范围从(115.623 683°E,41.210 212°N)到(115.666 799°E,41.179 376°N)，图12(a)为未掩膜时DEM生成后的高程误差值，图12(b)为掩膜后DEM生成后的高程误差值。为便于对比掩膜前后高程误差变化结果，采用红色框标记同一区域，可以看出经过掩膜处理后，图12中红色框区域内的高程误差明显降低，证明了本文提出的基于超像素分割与多信息融合的叠掩和阴影区域检测方法的有效性。

图12 TH-2实测数据DEM产品高程误差Fig.12 Elevation errors of DEM products based on TH-2 measured data

4 结 论

试验结果表明，本文提出的基于多信息融合的超像素检测算法能够快速准确地提取SAR图像中的叠掩和阴影区域，在叠掩区域识别上具有显著的优势，在阴影区域识别上，由于低回波能量区域和实际地物遮挡形成的阴影区域在数据特征上相似，在没有实际精确DEM辅助的情况下无法准确区分低回波能量区域和实际地物遮挡形成的阴影区域，但是在大部分实际需求中这些区域都会对后续处理造成影响，因此阴影区域的准确区分并不是必要的。本文所提算法检测的结果可作为相位解缠环节中的掩膜产品，提升相位解缠的速度和DEM产品的精度。同时提取的超像素信息还可用于SAR目标识别、升降轨融合等领域。本文提出的基于AIEM模型的频域回波仿真算法，实现了星载全极化双站SAR的原始回波数据仿真，模拟的回波能够反映SAR图像的叠掩、阴影、相干斑噪声等物理特征，回波仿真结果可用于成像算法验证和系统设计参数验证。