基于频谱空间的高分辨率遥感图像城市道路绿地轮廓特征检测

赵好好

0 引言

城市绿地是城市生态系统的重要组成部分，对城市空气质量和居民健康具有重要影响[1].道路绿地属于城市绿地4大类别中的附属绿地，随着城市建设和规划管理对城市绿地结构信息的进一步需求，获取详细的道路绿地信息具有重要意义.当前，高分辨率遥感图像使得地物空间几何结构更加明显，地物纹理和形状信息更加精细，充分利用这些信息，使得基于空间结构特征的城市道路绿地信息提取成为可能[2-4].

轮廓特征检测是当前遥感图像分割、目标特征提取的核心技术和研究热点[5-6].传统的轮廓特征检测主要基于空域的梯度算子，如Roberts算子、Laplace算子、LOG算子等.由于高空间分辨率遥感图像具有灰度变化丰富、噪声复杂、边缘密度大等特点，采用传统的基于微分的边缘检测算子往往难以达到理想的效果[7].通过傅里叶变换可以将图像从空间域变换到频率域，频谱图中的高频信息可以反映图像中灰度变化剧烈的轮廓特征.频谱作为物质的能量特征之一，与波谱、能谱、重力、磁力等特征一样，都是用来区别物体属性的重要依据[8].频谱能量值在描述边缘形态特征上具有独特的优势，频谱特征能够简明、准确地表达地物空间结构特征[6-7]，基于频谱分析和频域滤波的图像特征提取，为分析结构复杂的城市地物提供了新的方法.

已有研究主要探讨了规则地物的频谱分析和滤波方法.如吴桂平等[6]利用分块傅里叶变换和频谱能量分段叠置的手段，对规则结构的农田边缘特征进行了快速检测；周立国等[8]分析了线状地物的频谱特点，设计滤波器进行了厂房结构的特征提取.规则地物遥感图像的方向性、周期性显著，表现在频谱图上则谐波特征明显，易于检测，而城市道路绿地属于自然纹理，轮廓特征复杂且不规则，频谱图中谐波特征不明显，传统的特征提取方法并不适用，有必要探讨一种更为通用的频谱分析和特征识别的方法.

本文以QuickBird图像为研究数据，以城市道路绿地为研究对象，进行频谱能量分析，确定道路绿地轮廓的截止频率，建立基于图像特征的较为通用的Gabor滤波器参数设置方法，同时利用道路绿地与道路的几何关系，提取道路绿地轮廓特征，最后评价了轮廓提取的精度和本方法的适用性，对具有复杂结构的城市地物特征识别和提取具有理论和应用价值.

1 研究数据

本文所用的高分辨率遥感图像数据为南京市某地区夏季7月的QuickBird图像，包含4个多光谱波段和1个全色波段，基于Pansharp方法融合后的空间分辨率为0.6 m.基于绿色植被在近红外波段具有最强的光谱反射特性，与道路色调特征差异明显，宜于进行道路绿地特征分析，本文实验选择近红外波段.如图1所示，研究区大小为2 276像素×2 144像素，道路绿地主要有5条，如图中数字编号，道路绿地①基本呈水平方向，道路绿地②沿图中竖直方向，道路绿地③和④为十字交叉型，道路绿地⑤近似为竖直方向.这5条道路绿地均没有完全遮盖道路，且树木生长情况良好.

道路可用其中心线抽象表达，道路绿地的轮廓可表示为闭合曲线，因此道路与道路绿地的几何关系可表示为线与面的关系.道路绿地位于道路红线以内，若道路中心线能够确定，通过中心线的缓冲区可以对道路绿地进行定位，从而获得道路绿地轮廓.

2 研究方法

2.1 频谱能量分析

遥感图像可以表示为空域的二维函数，对图像进行二维傅里叶变换，相当于将二维函数分解成不同频率的正弦函数或余弦函数在二维平面上的叠加[9].这些频率揭示了空域图像上的周期性成分，通过分析不同频率的能量分布特征，可区别空域图像上不同特征的地物[10].大小为M×N的图像的二维离散傅里叶变换如下式表示：

(1)

式中，(x,y)为空域坐标，(u,v)为频域坐标.图像的频谱值表示为|F(u,v)|.

如图2所示，遥感图像中有一条45°方向的道路绿地，从频谱图可以看到，图像的中高频丰富，且在135°方向存在显著谱线.频谱图中的亮线，反映了图像中与其垂直方向的地物轮廓特征.为了进一步定量分析频谱图，利用辐射扫描法计算频谱图的径向分布特征和角向分布特征.假设极坐标中的频谱函数为S(r,θ)，r和θ是坐标变量.径向能量S(r)和角向能量S(θ)可分别表示为

(2)

(3)

径向分布图主要揭示周期性频率成分的能量分布，频谱分析主要是通过最大径向能量确定图像的主频.如图2c所示，图像的径向分布曲线随着频率的增加迅速减小，在低频位置约12周期/图像处出现拐点，这部分属于低频区域，主要揭示图像色调特征；在频率为13～36周期/图像的中高频段，频谱能量下降速度明显减缓，此频段反映图像结构、轮廓特征，是轮廓提取的主要研究频段；在频率36周期/图像之后，频谱能量的下降速度最慢，此部分频段主要揭示图像中的细节信息及噪声.径向分布曲线没有出现显著的谐波特征，这是因为道路绿地属于自然纹理，周期性纹理特征不如建筑物等地物特征显著，表现在径向分布特征中无明显的主频谐波.从图2d的角向分布图中可以看到，图像的最大角向能量位于135°，且峰值特征显著，这揭示了图像中与之垂直方向存在显著的地物轮廓特征，即道路绿地的方向特征.

2.2 截止频率确定

图像的截止频率在频域滤波中起着至关重要的作用.截止频率通常会作为滤波器设计的中心频率，因此截止频率的准确性直接决定信息提取的质量.一些学者对图像纹理特征的截止频率进行了研究，如Bovik等[11]将频谱峰值作为滤波器的中心频率进行纹理分割，Tan[12]通过约束平方法寻找截止频率，Jernigan等[13]利用信息熵进行频率分析.这些研究对图像纹理所对应的截止频率进行了探讨，而对图像轮廓特征的截止频率探测方法则少有报道.图像轮廓特征对应频谱图的中高频部分.由于图像轮廓位置的灰度变化剧烈，造成中高频位置的频谱曲线出现显著波动，因此可对频谱曲线中的极值频率进行重点探讨.本文在前期研究中进行了大量实验，发现以频谱曲线的波谷位置为中心频率进行滤波，可获得轮廓特征，结合文献[14]，这里将主要探寻频谱曲线的波谷位置.本文的研究对象城市道路绿地通常具有显著的方向性特征，可通过分析频谱图最大角向能量方向的频谱能量分布，确定道路绿地轮廓的截止频率.图2d揭示道路绿地的方向为45°，此方向的频谱能量分布曲线如图2e所示.在频率36周期/图像处存在波谷成分，但谐波特征并不显著.

在反射光谱曲线分析中，一阶导数因可以在一定程度上削减大气和背景噪声影响，所以被广泛应用于土壤、植物等诸多地物的特征分析中[15-16].而在频谱曲线分析中，受噪声的影响，频谱曲线特征会出现规律不显著的情况，尤其是自然纹理，难以从频谱曲线中发现谐波特征.本文通过计算最大角向能量分布曲线的一阶导数，消除噪声对频谱分布规律的影响，从而寻找对应道路绿地轮廓特征的截止频率.如图2f所示，将最大角向能量一阶导数曲线的中高频段所反映的峰值能量频率作为对应图像轮廓特征的截止频率，这里截止频率为36周期/图像.

2.3 Gabor滤波器设计

Gabor滤波器与哺乳动物视网膜神经细胞的接收场模型相吻合，它能够同时在时域和频域获得最佳的局部化，因此在图像纹理分割和边缘提取中得到广泛应用[17-18].Gabor滤波器参数设置方法有搜索优化法和特征提取法[19-21].其中搜索优化法是采用一定的优化准则进行最佳参数的搜索，这种方法是数据无关的，但在数据特征不同时，采用相同的优化准则会导致结果不理想.特征提取法是采用实验的方法，根据研究数据特征，设计与之匹配的滤波器参数，这种方法受研究数据特征影响较大，往往通用性差.本文根据上述频谱特征分析，建立一种更为通用的参数设置方法.

一般情况下，二维Gabor函数h(x,y)可表示为

h(x,y)=g(x′,y′)exp[2πj(Ux,Vy)],

(4)

式中，(Ux,Vy)表示特定的空间中心频率，g(x,y)为高斯函数，(x′,y′)是(x,y) 旋转θ角.

二维Gabor函数的傅里叶变换为

(5)

式中，(U′,V′)和(u′,v′)分别表示角度旋转θ后的中心频率和空间频率.

Gabor函数实际上是复数形式，将其分解可表示为

h(x,y)=Re(x,y)+jIm(x,y),

(6)

式中：

Re(x,y)=g(x,y)cos[ω(xcosθ+ysinθ)],

(7)

Im(x,y)=g(x,y)sin[ω(xcosθ+ysinθ)].

(8)

Gabor函数的实部是个偶函数，用实部对图像滤波仅起到图像平滑的作用.为提高检测精度并减少运算量，仅使用Gabor函数的虚部进行道路绿地轮廓的初步提取.一个θ方向的二维奇Gabor函数可表示为

(9)

其中ω为角频率，ω=2πf，f为频率，σ为空间常数.

奇Gabor函数中的角频率ω由频谱分析的截止频率来确定，方向θ则取决于道路绿地的方向.为选择最佳的空间常数σ，本文计算奇Gabor函数与式(10)所示的理想阶跃边缘函数的卷积，固定角频率ω和方向θ，在卷积函数随着空间常数σ的变化曲线中，取最大值处的σ作为滤波器设计参数.

(10)

然后，本文提取道路色调特征，进一步提取其中心脊线作为道路的中心线.色调信息集中在频谱图的低频中心[22]，可通过低通滤波提取.在近红外波段上，道路的方向特征显著，但色调最暗.因此，本文对原图像进行反变换，突出道路色调信息，然后基于低通滤波提取道路色调信息.用式(5)所示的Gabor滤波器进行低通滤波，滤波器中心设置为(0,0)，方向选择最大角向能量方向.

3 研究结果及评价

3.1 研究结果

研究区遥感图像频谱的径向分布曲线如图3a所示，频谱能量随着频率的增加迅速减小，在低频约60周期/图像处出现拐点，而在中高频都没有出现峰值谐波.从图3b的角向分布曲线中可以看到，最大角向能量位于0°(180°)和90°，且峰值特征显著，这揭示了图像中竖直和水平两条道路绿地(道路绿地②和道路绿地①)的方向特征；第3个显著的波峰位于45°，对应图像中135°的道路绿地(道路绿地③)方向特征；同时在105°、135°和173°也存在较为显著的谐波成分.其中135°的谐波成分反映了图像中45°的道路绿地(道路绿地④)方向特征；道路绿地⑤并非完全呈竖直直线，173°谐波则恰能反映此方向特征；研究区图像中与105°垂直的方向并不存在道路绿地，而是存在大量的建筑物，因此对105°方向的频谱能量不进行详细分析.

0°、45°、90°、135°、173°各方向的频谱能量分布曲线及一阶导数曲线如图4所示.将一阶导数曲线在中高频段的峰值频率作为截止频率，设计提取道路绿地轮廓的奇Gabor滤波器.结合一阶导数曲线的低频区域较高峰值，判断最大角向能量分布曲线的低频区域的拐点频率，设计提取道路信息的滤波器参数.具体滤波器参数如表1所示.

表1 滤波器参数Table 1 Parameters of filters

对研究区图像进行低通滤波，获得各方向的道路色调信息，然后对道路信息进行细化，提取其脊线，得到道路中心线.为了使结果显示更加直观，这里对提取结果进行反相处理，如图5a所示.整体来说，道路中心线保持了较好的连续性，没有出现较多断点、中断等现象.道路②在与道路③的交叉口处出现了部分不连续，这是因为此交叉口是属于环形交通岛，此处的道路不在竖直线上；同时道路⑤的最下端有部分没有被检测出来，这是因为此处道路方向有改变的趋势，所以在方向滤波时受到了影响.

对研究区图像进行奇Gabor滤波，获得轮廓图像.结合图像中的道路特征以及实地调查，5条道路宽度为24～30 m，以图5a所示的道路中心线为中心，以5条道路宽度分别为40、50、50、40、40像素为缓冲区，提取道路绿地轮廓，并进行数学形态学开运算，以去除内部噪声和边缘毛刺，得到道路绿地轮廓图像，其反相图像如图5b所示.可见，沿着道路绿地方向上的轮廓曲线基本完整，整体受噪声影响较小，基本不受建筑物遮挡等影响.

3.2 结果评价

将所提取的道路绿地轮廓与原图像叠加，如图6所示.所提取的图像轮廓与原图像中道路绿地轮廓边缘符合情况良好，没有出现大片漏检或过检等现象.但在道路②的下半部分，有少量道路绿地轮廓没有检测完整，这是因为此处的道路中心线受交叉路口的影响，出现了不连续性.在道路③、④的交叉处，有小部分道路绿地被重复检测，道路④、⑤的连接处，也有小部分被重复检测.

为了进一步评价研究结果，本文将讨论滤波器参数设置对道路绿地提取的影响，并利用Canny定量评价准则，与Canny算子的提取结果进行对比分析.

Mehrotra等[23]对奇Gabor滤波器的参数设置研究中，指出不同的ωσ取值，可以得到不同尺度的边缘检测效果.同时指出，对一维信号来说，当ωσ=1.0时，可以取得较好的边缘检测结果.傅一平等[24]在此基础上进行了二维图像边缘检测研究，指出当ωσ=1.8，且θ=π/2时可获得较好的边缘检测结果.但在图像特征检测时，如何选择适应识别目标的角频率ω或者空间频率f，目前还未见讨论.这里从研究区道路绿地②和③分别裁取两幅500像素×500像素的图像进行验证比较，分别令ω=1.00,σ=1.00及ω=1.80,σ=1.00，结果如图7所示，与本文轮廓提取结果差异较大.由此可知，图像的特征提取中，准确的滤波器参数设置会显著提高滤波效果.

进一步利用Canny准则进行精度评价，使用信噪比准则∑与定位精度准则∨的乘积∑∨以及单边缘响应准则C来评价边缘检测质量[25]，所得结果如表2所示.本文的信噪比和定位精度与一阶微分算子基本相当，与Mehrotra等[23]所提出的Gabor滤波器参数设置方法相比，信噪比与定位精度有所提高，而本文的单边缘响应准则明显优于一阶微分算子及文献[23]提出的参数设置准则.

表2 边缘检测精度比较Table 2 Precision comparision of edge detection

Canny边缘检测结果如图8所示，在检测出的道路绿地轮廓中，靠近道路一侧的轮廓比较清晰，但靠近建筑物一侧的道路绿地轮廓极不完整，多数无法形成闭合的曲线.Canny算子为全局检测算子，检测结果中没有显著突出道路绿地的响应特征，也无法突出道路的线状特征，不能利用道路中心线信息进行道路绿地的定位.

4 结论

本文结合高分辨率遥感图像的频谱空间特征，提出了一种基于Gabor滤波的城市道路绿地轮廓提取方法.通过频谱分析，寻找道路绿地轮廓对应的方向和截止频率，建立奇Gabor滤波器提取绿地轮廓信息，并利用低通滤波提取道路中心线，计算道路中心线的缓冲区获得道路绿地轮廓信息.本文利用最大角向能量一阶导数曲线的峰值确定截止频率，建立基于地物特征的遥感图像信息提取方法，优化了Gabor滤波器的参数，对城市地物的特征提取具有一定通用性.

本文研究区的道路绿地主要呈直线或近似直线型的分布，方向是Gabor滤波的重要参数之一，本文结果对方向性较强的城市道路绿地信息提取具有一定代表意义.对于城市道路方向多变的地区，在滤波器设计时需要适当调整方向参数.