高分辨InSAR中的城市高层建筑特征提取

郭 睿，臧 博，彭树铭，邢孟道

(1. 西北工业大学 自动化学院， 陕西 西安 710129；2. 西安电子科技大学 电子工程学院， 陕西 西安 710071；3. 958899部队， 北京 100085)

现代的高分辨率星载和机载合成孔径雷达(Synthetic Aperture Radar, SAR)系统，为城市测绘提供了新的研究方向[1-7]。与中分辨率(10～30 m)和高分辨率(3～10 m)合成孔径雷达数据相比，空间分辨率高达1 m的合成孔径雷达数据能够在成像时更好地反映高层建筑的特性，如，相邻窗户的间距和楼层之间的距离(例如3～3.5 m)[3,7]。基于超高分辨合成孔径雷达数据，更多技术和方法被用于对城市地区建筑目标的信息进行反演。例如，先进干涉合成孔径雷达技术之一的层析合成孔径雷达(Tomography Synthetic Aperture Radar, TomoSAR)[2]能够通过对20～100幅图进行处理得到目标的三维散射点云分布。从减少层析合成孔径雷达反演所需图像数量的角度出发，预先从合成孔径雷达图像中提取出关于高层建筑的典型特征，如建筑掩模、方向与处于同一高度和具有相似方向的等高同向像素等，成为笔者研究的出发点。

作为城市遥感的热点，许多学者围绕基于多模式合成孔径雷达数据的建筑信息反演和重构这一主题开展研究。随着分辨率的不断提高，基于不同模式的合成孔径雷达幅度信息被用于建筑体检测、参数估计及个体建筑的高度反演[1-7]。不同结构的高分辨干涉合成孔径雷达数据，例如多时相和多基线干涉数据，被用于建筑投影模型的构建和参数的反演[2-5]。但由于可能的运动、大气扰动和时间去相关等，通过多基线干涉图像对建筑体进行高程反演，尤其是在密集城区内仍面临挑战[3,5]。与多航过技术相比，单航过干涉技术获取的高分辨干涉合成孔径雷达数据不受形变、大气扰动和时间去相干的影响，如TanDEM-X成功掀起了基于干涉合成孔径雷达技术的城市遥感应用尤其是建筑信息提取和反演研究的热潮[3-5,7]，尤其是意大利和德国等欧洲学者在相关方面进行了较多的研究，通常都是将非高层建筑作为研究对象。

文中以高层建筑的特征提取为主要研究目的，基于TanDEM-X获取的侧视干涉合成孔径雷达数据，对实验区域中的高层建筑在幅度、干涉相位图像中的投影特性进行分析，在此基础上提出高层建筑的等高同向像素特征定义，并提出一种自上而下逐步细化的高层建筑提取和特征估计方法，实现城市高层建筑特征信息的提取。

1 数据说明

1.1 测试区域

文中使用的单航过干涉合成孔径雷达数据入射角为36°，有效基线长度为169.4 m，空间分辨率为1.1 m×0.6 m (方位向×距离向)。选择的测试区域一和区域二都来自于美国拉斯维加斯的高层区域，如图1所示。图1(b)是两个区域的合成孔径雷达幅度图，相应的谷歌光学图为图1(a)。从光学图中可以看出，区域一中的高层建筑外形较为规则，而区域二中的高层建筑外形呈现有曲率的弧形。

图1 高层区域一(第一行)和高层区域二(第二行)

1.2 InSAR图像中的高层建筑投影

在侧视超高分辨合成孔径雷达图像中，高层建筑的墙体结构相较于顶部的散射特性更加明显和突出[6]，高强度的叠掩是合成孔径雷达幅度图中高层建筑的显著标志，因此文中主要研究高层建筑墙体部分的合成孔径雷达投影形成的叠掩。采用如图2所示的建筑散射投影模型，其中θ为入射角，建筑体宽为w、高为h，对于高层建筑h≫w·tanθ。R=L·sinθ是高层建筑的合成孔径雷达投影叠掩的长度。

图2 高层建筑的SAR几何投影模型对应图

建筑的方向角φ被定义为建筑墙体与方位向之间的夹角[7]，如图2所示。如果墙体外形不规则，则建筑方向角一般被定义为建筑物主要部分的方向。图中高层建筑叠掩处可认为处于相同的高度，同时位于同一方向上即具有相似的方向角，定义这些像素为等高同向像素，将其在合成孔径雷达幅度图中形成的直线或者曲线定义为等高同向线。文中研究的建筑特征就包含这些具有相同高度的像素及其形成的直线或曲线结构，从而为模拟建筑的合成孔径雷达投影模型[3]和从层析合成孔径雷达形成的点云[2]提取处于同一层的散射点提供重要参考。

干涉合成孔径雷达相干系数图中，建筑墙体部分尤其是墙角线对应的二次散射，具有很强的相干性[4,5]，而非人造设施区域的相干性则较低，如图1(d)所示。在图1(c)的未解缠干涉相位图中，高层建筑的叠掩区域呈现周期性的条纹特性，这是由于缠绕相位对应于建筑墙体的高度。叠掩长度越长，干涉相位的周期越多。

2 方法介绍

笔者采用自上而下的方法，通过图像、感兴趣区域、建筑体和像素四个层面的逐步细化处理，联合高分辨干涉合成孔径雷达数据的幅度、相干系数及未解缠干涉相位图来提取高层建筑的特征。

2.1 图像级别处理

在图像层面的处理中，包含高层建筑物或部分高层建筑的感兴趣区域(Region Of Interest, ROI)可以通过对干涉合成孔径雷达数据的两幅幅度图的图像处理进行提取，从而获得关于建筑位置的先验信息，具体步骤如下。

(1) 粗分割。作为一种有效的分割算法，均值漂移算法已成功应用于合成孔径雷达图像解译[7-8]。该算法不需要关于分类数量等先验信息，因此将两幅幅度图像作为输入，采用均值漂移算法进行全自动的粗略分割。

图3 建筑区域(左图)和非建筑区域(右图)的相干系数统计分布

(2) 分块检测。粗分割后得到不同分类的区域块为Sc={s1,s2,s3,…,si,…,sN},i∈[1,N]，其中N是分块总数，si是第i个分类块。这些分块区域有建筑区域，还有散射较强的非建筑区域。相比于非建筑区域，建筑区域有较高的相干系数。图3为建筑区域和非建筑区域的相干系数的统计分布。可以看出，建筑区域的相干系数明显高于非建筑区域，因此根据相干性的阈值来选择建筑区域。选择平均相干系数高于阈值Tco1(根据图5设置Tco1= 0.8) 的区段作为候选感兴趣区域，即

Scan={s1,s2,…,si,…,sk|E(Co(si))>Tco1},k

(1)

其中，E()表示样本均值，Co()表示相干系数。

(3) 感兴趣区域标记。由于高层建筑的合成孔径雷达投影叠掩区域相对较大，在候选感兴趣区域中去除区域小的Scan。通过连通性分析首先选择孤立的区域块，如果一个区域块内任一像素的邻域像素属于另一个区域块，则该区域块不是孤立的。考虑到高层建筑的严重叠掩特性，高层建筑区域块应包含较多的像素，设置最小像素数为NT，得到最终的感兴趣区域区域块为

SROI={s1,s2,…,si,…,sr|Nseg(sr)≥NT} ，

(2)

其中，SROI⊂Scan，Nseg(sr)表示sr中的总像素数。文中将NT设为Scan中候选感兴趣区域的平均大小。

2.2 ROI级别处理

经过图像级别处理，所提取的感兴趣区域与各高层建筑并不是一一对应的，单个建筑体可能被分割成多个感兴趣区域，而一个感兴趣区域中可能存在多个建筑体，因此对提取的感兴趣区域进行细化，使得感兴趣区域对应各个建筑体。

(1) 单个建筑体被分割成多个感兴趣区域的情况存在于空间上连接的区域块之间，因此依据区域块的可连接性将分割成多个感兴趣区域的单个建筑体进行区域合并，然后对存在多个建筑体的感兴趣区域进行分离。

(2) 采用已成功用于合成孔径雷达图像解译的形态学方法[8]，首先通过闭运算对感兴趣区域进行扩展，使得周边属于同一建筑的像素也被包含在感兴趣区域块内；之后执行开运算将错连起来的建筑体分开。

(3) 对各高层建筑对应的掩模进行标记。

2.3 建筑级别处理

通过对干涉相位图处理，估计出高层建筑的方向角，在此基础上提取建筑体对应的平行四边形模型。由于高层建筑周边的非人造区域(即低相干区域)的相位噪声会影响建筑方向角估计，因此采用如图1(b)～(d)所示的处理，通过干涉相位与相干系数相乘，将低相干区域的干涉相位消减到接近零。

(1) 建筑方向角估计。基于提取的建筑掩模，首先确定局域干涉条纹频率估计采用的窗口。采用极大似然估计方法[9]在建筑叠掩所在窗口内进行干涉条纹的二维频率估计。通过最大化二维离散傅里叶变化，找到其峰值分布的频率作为建筑叠掩干涉条纹的频率，从而得到这一局部频率的方向角作为建筑的方向角φ=arctan(fy/fx)。当估计得到的方向角小于0时，则给该角度加上π作为建筑方向角。

(2)建筑块提取。合成孔径雷达图像中建筑物更直观的外形特征是平行线[5-7]，因此将高层建筑叠掩对应于不同的平行四边形模型。比较容易确定的是平行四边形沿距离向的两条边，然后结合提取的建筑掩模和建筑方向角确定平行四边形的其余两边。通过确定方位向上建筑掩模边界的最远点和最近点，构建以建筑方向角为斜角的连接两点的线段，从而完成建筑块的提取。

2.4 像素级别处理

在像素级别处理中，将对2.2节中定义的建筑特征进行估计。如图4所示，根据相似的干涉相位、相干系数和幅度值等共同特点对等高同向像素进行提取。干涉合成孔径雷达的有效基线造成了高层建筑叠掩的非解缠干涉相位呈现出干涉条纹，通过干涉条纹的方向，可以提取等高同向像素的粗略信息；超高分辨合成孔径雷达幅度图像中不同楼层对应的强角反射规则结构被用来估计等高同向像素的精确信息；最后通过干涉幅度和相位值对提取的等高同向像素进行改进。

(1) 提取初始等高同向像素。考虑到噪声对方法和实验的影响，使用非局域滤波方法实现干涉相位的噪声抑制[10-11]。然后，提取具有相似干涉相位值的像素作为候选等高同向像素。在文中，非解缠相位值位于[-π, π)，选择相位值大于0.96π的像素。建筑体沿墙体方向存在一条以上的等高同向线(取决于建筑物高度和分辨率)，因此再次采用均值漂移算法对提取的像素进行自动分类，将像素数最多类别中的所有像素作为初始等高像素。

(2) 精细化等高同向像素。由于沿距离向通常存在若干个已选择的初始等高同向像素，因此通过如式(4)的计算，使得在同一个方位向上只存在一个固定的等高同向像素，即精细等高同向像素的距离坐标r为

图4 图像处理

(3)

其中，αi表示固定方位向坐标上第i个初始等高同向像素的干涉相位，Nround()表示四舍五入取整。然后再次检验精细等高同向像素处的干涉相位值，排除绝对相位值小于0.96π的像素。图4(a)和(c)分别对应提取的初始等高同向像素和按式(3)计算后得到的等高同向像素示例。图4(b)为按式(3)得到的每个像素的加权系数。图4(d)是最终的精细等高像素。

(3) 估计等高线。根据图4(d)所示的精细等高像素，用包括傅立叶，高斯，谐波和多项式等4种曲线拟合方法来拟合等高线。通过比较发现，多项式模型足以获得准确的结果。

(4) 最终的等高像素选择。通过设置幅度和相干系数的阈值为

Tam=E(Am)+2×ψ(Am),Tco2=E(Co)+2×ψ(Co) ，

(4)

将幅度和相干系数均大于阈值的精细等高像素识别为最终的等高像素。其中，Am和Co分别是精细等高像素的幅度和相干系数，ψ()表示方差。

3 实验结果与分析

对应于上节中所提的从上至下4个级别的处理，对测试区域一和区域二中的高层建筑进行特征提取的分步结果进行讨论和分析。

图5 (a)是高层区域一和区域二中提取的感兴趣区域结果，不同颜色表示不同的感兴趣区域标记。可以看出，经过图像级别的处理，所提取的感兴趣区域内存在单个建筑体被分割成多个感兴趣区域，和一个感兴趣区域中存在多个建筑体的情况。因此，需要进行感兴趣区域级别的处理，将感兴趣区域与高层建筑一一对应。经过感兴趣区域层面的处理，在区域一和区域二中提取到的高层建筑掩模标记如图5(b)所示。需要说明的是，在区域二中3号楼被标记为两部分，这一点从图1中区域二的光学图像和幅度独享可以看出，3号楼的主墙体与方位向接近垂直，因此对应的建筑叠掩部分侧墙体分布占主要分布。

区域一和区域二中提取到的高层建筑的方向角估值如表1。需要说明的是，在区域一中方向角较小的高层建筑叠掩强度更强，例如3号、4号和5号楼，证明了具有较小方向角的建筑体在SAR幅度图中显示出更强的二次反射[3]。在表1中，区域二中3号楼两部分掩模的方向角分别估计为77.47°和147.53°；1号楼和2号楼的方向角很接近，这点在光学图像也可以看出。

表1 区域一和区域二中的高层建筑方向角估值

基于估计的方向角和提取的建筑掩模，构建的平行四边形建筑块如图5 (c)所示。在后续的等高线估计时，将区域二中3号楼的两个建筑块合并，同时提取整个建筑体的等高线。

图5 区域一和区域二的中间结果

图6 等高同向像素提取结果及等高线估计结果

如图6所示，将精细化等高同向像素(实心点)、最终选择的等高同向像素(空心点)、等高线(曲线)叠加于幅度图中的高层建筑叠掩。在区域一中，等高线的主体部分平行于平行四边形的边，这表明反射主要来自平坦的主墙体。然而，等高线的方向在平行四边形的左侧或右侧附近略有变化，表明这部分反射来自侧墙体。在区域二中等高线不是直线型，主要是由于复杂的曲线型主墙体表面结构，如图1 中的光学图像所示。

表2列出两个区域中精细化的等高同向像素和最终的等高同向像素的相干系数均值和归一化幅度均值。可以看出，最终的等高同向像素的相干系数均值和归一化幅度均值都大于精细化等高像素，即有更高的信噪比。

表2 精细化及最终的等高同向像素的相干系数均值及归一化幅度均值

表3 精细化及最终等高同向像素与所估计等高线的均方误差

为了评估等高同向像素的精度，表3中列出了精细化等高像素和最终等高像素相对于所估计的等高线的均方误差。最终的等高像素距离等高线的均方误差更小，远低于距离方向上的像素尺寸。同时可以看出，相比于区域二，区域一中估计出的等高线对于等高同向像素拟合得更精确，这是由于区域一中的高层建筑主墙体更平坦。

4 结束语

文中依据高层建筑在高分辨干涉合成孔径雷达图像中投影的结构与特性，提出建筑方向角和等高同向像素等高层建筑的特征参数定义；联合干涉合成孔径雷达幅度、相干系数及未解缠干涉相位图像的信息，通过图像、感兴趣区域、建筑体和像素四个层面的逐步细化处理，实现了高层建筑的特征提取。提出了基于未解缠干涉相位的建筑方向角估计方法；并在平行四边形假设下，对建筑块进行检测提取；对所提出的高层建筑特征进行提取估计。通过TanDEM-X系统获取的拉斯维加斯高层区域的实验验证了所提方法的有效性。论文方法较好地实现了对个体高层建筑的检测及建筑个体的等高同向像素等特征的提取，为进一步的高层建筑三维信息反演等其他信息特征的估计以及层析合成孔径雷达的应用提供了参考信息。