CN105205816A - 多特征加权融合的高分辨率sar影像建筑区提取方法 - Google Patents

Abstract

本发明涉及一种多特征加权融合的高分辨率SAR影像建筑区提取方法，其主要特点在于应用距离公式确定特征权值的策略较为新颖，包括以下步骤：①对强度影像进行滤波等预处理；②分别基于灰度共生矩阵纹理分析方法和基于变差函数纹理分析方法进行特征提取；③依据巴士距离进行特征权值确定；④将步骤②得到的特征进行加权融合；⑤利用K均值聚类算法对融合后的特征图像进行非监督分类；⑥后处理分类结果去除小区域，填补空洞，并提取外部轮廓，得到建筑区。本发明充分综合了多种有利于建筑区提取的纹理特征，构造的融合特征图包含了比单特征图更加丰富的纹理信息，从而一定程度上提高了建筑区与非建筑区的分类精度，进而提高了建筑区提取精度，最终使得提取到的建筑区边界轮廓更加拟合真实建筑区轮廓。

Description

多特征加权融合的高分辨率SAR影像建筑区提取方法

技术领域

本发明属于遥感影像的合成孔径雷达图像解译领域，特别涉及机载合成孔径雷达遥感影像的特征提取和特征融合方法。

背景技术

遥感影像上的建筑区作为典型的人造结构，在影像解译中一直是研究的热点和重点。随着遥感技术的迅速发展，合成孔径雷达以其全天时全天候的优势在遥感技术中应用越来越多。在目前SAR影像解译领域，高分辨率SAR影像建筑区的提取研究已然成为一项备受关注的课题。高分辨率SAR影像建筑区提取是指通过一定的方法，提取有利于区分建筑区的纹理特征，通过纹理特征将建筑区和非建筑区区分开来。考虑到实际中，纹理分析方法的多样性，计算得到的纹理特征也有很多种，通常只用单特征或者几个经验特征进行分类计算。

纹理分类问题包含特征提取和特征分类两部分，其中特征提取的优劣是决定分类效果的关键因素。目前主流的纹理特征描述方法主要包括基于统计，基于结构，基于模型和基于滤波四种。由于SAR图像上建筑区往往表现出区别于其他地物的典型纹理特征，纹理分析用于SAR图像的信息提取中，纹理分类方法成为SAR图像上提取建筑区的主要手段。目前研究的方法有：(1)利用灰度共生矩阵计算高分辨率SAR图像纹理特征，以选取的特征矢量作为非监督分类的特征输入提取居民地。(2)基于变差函数纹理特征的高分辨率SAR图像建筑区提取方法，基于变程计算变差函数值作为分类特征，结合FCM分类器实现了建筑区和非建筑区两类分类问题。(3)综合灰度和纹理特征的高分辨率SAR图像建筑区提取方法，根据巴士距离对纹理特征进行选择，最后与原始灰度图像进行波段组合作为分类特征。(4)利用Markov随机场纹理分类方法对城区进行精细的分类。在上述研究中，多采用基于统计模型的单一纹理描述方法。Markov法在高分辨率条件下，对极不均匀的区域难以准确建模。基于结构信息的变差函数法在实验中表现了其优势性，但没有考虑高分辨率条件下建筑区的方向信息，忽视了变差函数的各向异性。综合多特征的方法，在特征选取上，采用了直接抛弃贡献较小特征的方法，难免会损失部分信息。

发明内容

本发明的目的在于提高现有建筑区提取方法的精度，提出一种基于巴氏距离权值计算方法，以此为基础，提出了多特征加权融合的高分辨率SAR影像建筑区提取方法。

本发明的技术方案相对现有方法，其主要改进在于在纹理特征选取策略上，引入一种距离测度来计算每种纹理特征值的权值，最后将多种纹理特征值融合为一种新特征。本文发明采用了应用距离公式确定特征权值的策略，充分考虑了多种统计纹理特征值和结构纹理特征值的贡献大小。

本发明提供了一种多特征加权融合的高分辨率SAR影像建筑区提取方法，该方法包括以下步骤：①对强度影像进行预处理；②分别基于灰度共生矩阵纹理分析方法和基于变差函数纹理分析方法进行纹理特征值提取；③依据巴氏距离进行特征权值的确定；④将步骤②得到的纹理特征值利用所述特征权值进行加权融合；⑤利用K均值聚类算法对融合后的特征图进行非监督分类；⑥后处理分类结果，去除小区域，填补空洞，并提取外部轮廓，得到建筑区。

本发明所述方法充分综合了多种有利于建筑区提取的纹理特征值，构造的融合特征图包含了比单特征图更加丰富的纹理信息，从而一定程度上提高了建筑区与非建筑区的分类精度，进而提高了建筑区提取精度，最终使得提取到的建筑区边界轮廓更加拟合真实建筑区轮廓。

进一步优选地，所述步骤①中，使用的强度影像是空间分辨率为1m的机载P波段影像块，并且通过执行滤波进行所述预处理。

进一步优选地，步骤②所述的，基于灰度共生矩阵纹理分析方法中，纹理灰度共生矩阵的计算表达式为P(i,j|d,θ)＝{P_k(i,j|d,θ)}_N×N，其中，i、j分别表示灰度共生矩阵的行列号，步长d、方向θ分别表示用来区分纹理的粗糙性与方向性，N表示图像的灰度级数，k为像元顺序号。本发明所述纹理特征值是由所述灰度共生矩阵P(i，j|d，θ)计算的以下5种参量：

均值 $MEA=\underset{i,j}{\Σ}i\×P(i,j|d,\mathrm{\θ}),$ 熵 $ENT=-\underset{i,j}{\Σ}P(i,j|d,\mathrm{\θ})\×log\left\{P(i,j|d,\mathrm{\θ})\right\},$

对比度 $CON=\underset{i,j}{\Σ}{(i-j)}^{2}P(i,j|d,\mathrm{\θ}),$ 方差 $Var=\underset{i}{\Σ}\underset{j}{\Σ}{(i-\mathrm{\μ})}^{2}P(i,j|d,\mathrm{\θ}),$

角二阶矩 $Energy=\underset{i,j}{\Σ}{P}^2(i,j|d,\mathrm{\θ}).$

进一步优选地，步骤②所述的基于变差函数的纹理分析方法中，用实验变差函数来估计无法直接获得的变差函数为：式中，N(h)表示观测数据中间距为h的点对数目，f(x_k)表示x_k处的像素值，将h取特定值时的变差函数值作为所述纹理特征值。

进一步优选地，所述步骤③中，确定特征权值确定的方法包括：

计算巴氏距离，计算公式为：

$BD=\frac{1}{4}\frac{{({\mathrm{\μ}}_1-{\mathrm{\μ}}_2)}^2}{{{\mathrm{\σ}}_1}^2+{{\mathrm{\σ}}_2}^2}+\frac{1}{2}log\left(\frac{{{\mathrm{\σ}}_1}^2+{{\mathrm{\σ}}_2}^2}{2{\mathrm{\σ}}_1{\mathrm{\σ}}_2}\right)$

式中，μ₁、μ₂、σ₁、σ₂分别表示在同一幅纹理特征图上建筑区和非建筑区的均值和方差，BD表示巴氏距离；

计算某一类型纹理特征值的特征权重，计算公式为：

$p_i=\frac{\stackrel{\‾}{{\mathrm{BD}}_i}}{\Σ{\mathrm{BD}}_i}$

式中，表示某一特征图上建筑区和非建筑区的巴氏距离值，p_i表示某一类型纹理特征值所占的权值。

进一优选步地，所述步骤④中，将步骤③中计算的特征权值分配给对应类型的纹理特征值，经过加权运算得到融合后的特征图。

进一优选步地，所述步骤⑤中，使用K-means聚类算法对融合后的特征图进行非监督分类。

进一优选步地，所述步骤⑥中，后处理分类结果去除小区域，填补空洞，并提取外部轮廓，得到建筑区。

进一优选步地，所述步骤②的基于灰度共生矩阵纹理分析方法中，计算纹理特征值是基于取纹理窗口来遍历灰度共生矩阵的，所述纹理窗口大小参数按如下方法选择：在3-45之间，采用以4为等差依次取值，作为窗口大小参数，绘制建筑区与非建筑区的纹理特征值随窗口大小参数变化而变化的曲线图，通过分析曲线图的变化趋势得出最佳纹理窗口大小参数；所述步骤②的基于变差函数纹理分析方法中，纹理间距按如下方法确定：将h作为自变量从1开始不断递增，计算不同h下建筑区与非建筑区样本均值，绘制变差函数曲线图，通过分析曲线图求得最佳纹理间距。

本发明所述方法充分综合了多种有利于建筑区提取的纹理特征值，构造的融合特征图包含了比单特征图更加丰富的纹理信息，从而一定程度上提高了建筑区与非建筑区的分类精度，进而提高了建筑区提取精度，最终使得提取到的建筑区边界轮廓更加拟合真实建筑区轮廓。

附图说明

下面结合附图对本发明进行进一步详细的描述，以使得本发明的上述优点更加明确。

图1是本发明所述的多特征加权融合的高分辨率SAR影像建筑区提取方法的流程图；

图2是本发明计算灰度共生矩阵纹理窗口参数分析曲线图；

图3是本发明计算变差函数纹理间距参数分析曲线图；

图4是本发明所述的特征权值计算和特征融合示意图；

图5是采用本发明方法和其他方法的分类结果对比图和后处理图；

图6是建筑区轮廓提取图。

具体实施方式

下面结合附图对本发明作进一步详细的描述。

图1表示本发明所述的多特征加权融合的高分辨率SAR影像建筑区提取方法的流程图。本发明提供的多特征加权融合的高分辨率SAR影像建筑区提取方法，针对建筑区的纹理特点，，分别采用经典的灰度共生矩阵方法提取统计纹理特征和采用变差函数方法提取结构纹理特征值，综合统计和结构纹理特征值，并考虑方向信息，然后利用巴氏距离特征权值计算方法计算特征权值，将所选纹理特征值进行加权融合，利用K均值聚类算法对融合后的特征图进行非监督分类，对分类图像进行后处理并提取外部轮廓。具体来说，本方法包括以下步骤：

步骤①：首先，对使用的强度影像执行滤波等预处理，使用的强度影像是执行滤波等预处理的空间分辨率为1m的机载P波段影像块。

步骤②：采用以下两种纹理分析方法，分别计算纹理特征值，具体方法如下：

a)基于灰度共生矩阵的纹理分析方法：计算纹理灰度共生矩阵的计算表达式为P(i,j|d,θ)＝{P_k(i,j|d,θ)}_N×N，其中，i、j分别表示灰度共生矩阵的行列号，步长d、方向θ分别表示用来区分纹理的粗糙性与方向性，N表示图像的灰度级数，k为像元顺序号。取图像(N×N)中任意一点(x，y)及偏离它的另一点(x+a，y+b)，设该点对的灰度值为(g1，g2)。令点(x，y)在整个画面上移动，则会得到各种(g1，g2)值，设灰度值的级数为k，则(g1，g2)的组合共有k的平方种。对于整个画面，统计出每一种(g1，g2)值出现的次数，然后排列成一个方阵，再用(g1，g2)出现的总次数将它们归一化为出现的概率P(g1，g2)。通过上面的表达式得到灰度共生矩阵是纹理特征提取的基础，作为图像纹理分析的纹理特征是由灰度共生矩阵计算获得的参量，下表中可计算获得的本发明用到的5种纹理特征值。

基于灰度共生矩阵的纹理分析方法中，计算纹理特征值是基于取纹理窗口遍历灰度共生矩阵，那么对纹理窗口的选择，采用以4为等差的从3-45之间的窗口大小，绘制建筑区与非建筑区的两种特征值随窗口变化而变化的曲线图(如图2)，通过分析曲线图的变化趋势可以得出最佳纹理窗口大小参数。参见图2所示的计算灰度共生矩阵纹理窗口分析曲线图，当窗口大小在29以上时，样本均值达到一个相对稳定值；建筑区和非建筑区的对比度特征值不会随窗口大小的改变而明显改变，而是保持一个稳定的值。考虑到大窗口会使边缘信息模糊，所以我们选择29*29的窗口大小，基于该窗口大小计算纹理特征值。

b)基于变差函数的纹理分析方法：变差函数用于度量区域化变量的空间相关性，在任意一点x∈Ω处，变差函数值表现为一个随机变量，在任意两点处，所对应的随机变量通常是不独立的。变差函数可定义为区域化变量在空间中相距为h(两点间距，包括距离和方向信息)的两个位置处的差值的方差，即：

2γ(x,h)＝Var[f(x+h)-f(x)]，其中f表示相应位置的像素值，γ(x，h)表示半变差函数值，Var表示进行方差运算。在实际应用中，，一般用实验变差函数来估计无法直接获得的变差函数：中，N(h)表示观测数据中间距为h的点对数目，在计算实验变差函数时一般会固定点对的方向，当h包含的标量h变化时，就可以计算一系列的变差函数值，从而可以建立变差函数曲线，帮助我们分析数据。

基于变差函数纹理分析方法，对于参数h的确定，我们可以将h作为自变量从1开始不断递增，计算不同h下建筑区与非建筑区样本均值，绘制变差函数曲线图(如图3)，通过分析曲线图求得最佳纹理间距。图3是本发明计算变差函数纹理间距参数分析曲线图，变差函数值会随着h的变化也展现出周期性。我们选择峰谷时，即h＝13对应的变差函数特征值作为纹理特征值。

步骤②中，基于灰度共生矩阵纹理分析选取的特征都是在建筑区提取中常用的特征；基于变差函数的纹理特征对于排列规则有序的建筑区表现优良。

步骤③：依据巴氏距离进行特征权值确定，详细方法包括：

图4是本发明所述的特征权值计算和特征融合示意图。

本发明提出的纹理特征值加权融合方法的具体流程中，在纹理特征值进行融合前，需要确定每种纹理特征值对应的权值。在现有技术的方法中，权值的确定通常根据经验，或者专家投票的方式人工分配，这加入了过多的人工干预，必然会影响其结果的可信性。权值是指某一部分对整体贡献大小的比例，本文从其本质出发，提出了结合巴氏距离的特征权值确定方法，其中巴氏距离公式为：

$BD=\frac{1}{4}\frac{{({\mathrm{\μ}}_1-{\mathrm{\μ}}_2)}^2}{{{\mathrm{\σ}}_1}^2+{{\mathrm{\σ}}_2}^2}+\frac{1}{2}\log \left(\frac{{{\mathrm{\σ}}_1}^2+{{\mathrm{\σ}}_2}^2}{2{\mathrm{\σ}}_1{\mathrm{\σ}}_2}\right)$

式中，μ₁、μ₂、σ₁、σ₂分别表示在同一幅纹理特征图上两个不同类别(建筑区和非建筑区)的均值和方差。BD值的大小决定了两个类别的区分能力，BD值越大区分能力越强，越小则区分能力越弱。

$p_i=\frac{\stackrel{\‾}{{\mathrm{BD}}_i}}{\Σ{\mathrm{BD}}_i}$

式中，表示某一特征图上建筑区与非建筑区的巴氏距离值，那么p_i表示某一特征所占的权值，BD值越大，p值越大。

步骤④：将上一步中得到的权值分配给对应的纹理特征值，经过加权运算得到融合后的特征图。

步骤⑤：选用K-means聚类算法对融合后的特征图进行非监督分类。

步骤⑥：最后，后处理分类结果去除小区域，填补空洞，并提取外部轮廓，得到建筑区。

本发明针对建筑区的纹理特点，分别采用经典的灰度共生矩阵方法提取统计纹理特征值和采用变差函数方法提取结构纹理特征值，综合统计和结构多种纹理特征值，并考虑方向信息，然后利用提出的巴氏距离特征权值计算方法，将所选纹理特征值进行加权融合，利用K均值聚类算法对融合后的特征图进行非监督分类，对分类图像进行后处理并提取外部轮廓。以国产机载P波段全极化SAR影像为数据源进行了实验，并对结果进行了定量分析，表明该方法能够有效高精度地提取高分辨率机载SAR影像中的建筑区。

图5是采用本发明方法和其他方法的分类结果对比图和后处理图，可以看出本发明较其他方法有更好的提取效果，图6是建筑区轮廓提取图，下表是权值计算数据表。

可见，本发明提出了以机载高分辨率SAR图像为数据源，综合基于统计的灰度共生矩阵和基于结构的变差函数纹理描述方法，并在特征提取过程中结合建筑物走向的先验信息充分考虑特征方向信息，提取多种分类纹理特征值，并对提取的多种纹理特征值进行加权融合作为非监督分类的特征输入进行分类，最后对分类结果进行后处理，并对结果进行了精度评定，证明该方法能够有效的高精度地对高分辨率SAR图像建筑区进行提取。

本发明所述方法充分综合了多种有利于建筑区提取的纹理特征值，构造的融合特征图包含了比单特征图更加丰富的纹理信息，从而一定程度上提高了建筑区与非建筑区的分类精度，进而提高了建筑区提取精度，最终使得提取到的建筑区边界轮廓更加拟合真实建筑区轮廓。

上述具体实例仅仅是示例性的，在本发明的上述教导下，本领域技术人员可以在上述实例的基础上进行各种改进和变形，而这些改进或者变形落在本发明的保护范围内。本领域技术人员应该明白，上面的具体描述只是为了解释本发明的目的，并非用于限制本发明。本发明的保护范围由权利要求及其等同物限定。

Claims (9)

1.一种多特征加权融合的高分辨率SAR影像建筑区提取方法，该方法包括以下步骤：①对强度影像进行预处理；②分别基于灰度共生矩阵纹理分析方法和基于变差函数纹理分析方法进行纹理特征值提取；③依据巴氏距离进行特征权值的确定；④将步骤②得到的纹理特征值利用所述特征权值进行加权融合；⑤利用K均值聚类算法对融合后的特征图进行非监督分类；⑥后处理分类结果，去除小区域，填补空洞，并提取外部轮廓，得到建筑区。

2.根据权利要求1所述的方法，其特征在于，所述步骤①中，使用的强度影像是空间分辨率为1m的机载P波段影像块，并且通过执行滤波进行所述预处理。

3.根据权利要求1所述的方法，其特征在于，步骤②所述的，基于灰度共生矩阵纹理分析方法中，纹理灰度共生矩阵的计算表达式为P(i,j|d,θ)＝{P_k(i,j|d,θ)}_N×N，其中，i、j分别表示灰度共生矩阵的行列号，步长d、方向θ分别表示用来区分纹理的粗糙性与方向性，N表示图像的灰度级数，k为像元顺序号；

所述纹理特征值是由所述灰度共生矩阵P(i，j|d，θ)计算的以下5种参量：

均值 $MEA=\underset{i,j}{\Σ}i\×P(i,j|d,\mathrm{\θ}),$ 熵 $ENT=-\underset{i,j}{\Σ}P(i,j|d,\mathrm{\θ})\×log\left\{P(i,j|d,\mathrm{\θ})\right\},$

对比度 $CON=\underset{i,j}{\Σ}{(i-j)}^{2}P(i,j|d,\mathrm{\θ}),$ 方差 $Var=\underset{i}{\Σ}\underset{j}{\Σ}{(i-\mathrm{\μ})}^{2}P(i,j|d,\mathrm{\θ}),$

角二阶矩 $Energy=\underset{i,j}{\Σ}{P}^2(i,j|d,\mathrm{\θ}).$

4.根据权利要求1所述的方法，其特征在于，步骤②所述的基于变差函数的纹理分析方法中，用实验变差函数来估计无法直接获得的变差函数为：式中，N(h)表示观测数据中间距为h的点对数目，f(x_k)表示x_k处的像素值，将h取特定值时的变差函数值作为所述纹理特征值。

5.根据权利要求4所述的方法，其特征在于，

所述步骤③中，确定特征权值确定的方法包括：

计算巴氏距离，计算公式为：

$BD=\frac{1}{4}\frac{{({\mathrm{\μ}}_1-{\mathrm{\μ}}_2)}^2}{{{\mathrm{\σ}}_1}^2+{{\mathrm{\σ}}_2}^2}+\frac{1}{2}log\left(\frac{{{\mathrm{\σ}}_1}^2+{{\mathrm{\σ}}_2}^2}{2{\mathrm{\σ}}_1{\mathrm{\σ}}_2}\right)$

式中，μ₁、μ₂、σ₁、σ₂分别表示在同一幅纹理特征图上建筑区和非建筑区的均值和方差，BD表示巴氏距离；

计算某一类型纹理特征值的特征权重，计算公式为：

$p_i=\stackrel{\‾}{\frac{{\mathrm{BD}}_i}{{\mathrm{\ΣBD}}_i}}$

式中，表示某一特征图上建筑区和非建筑区的巴氏距离值，p_i表示某一类型纹理特征值所占的权值。

6.根据权利要求5所述的方法，其特征在于，所述步骤④中，将步骤③中计算的特征权值分配给对应类型的纹理特征值，经过加权运算得到融合后的特征图。

7.根据权利要求6所述的方法，其特征在于，所述步骤⑤中，使用K-means聚类算法对融合后的特征图进行非监督分类。

8.根据权利要求7所述的方法，其特征在于，所述步骤⑥中，后处理分类结果去除小区域，填补空洞，并提取外部轮廓，得到建筑区。

9.根据权利要求8所述的方法，其特征在于，所述步骤②的基于灰度共生矩阵纹理分析方法中，计算纹理特征值是基于取纹理窗口来遍历灰度共生矩阵的，所述纹理窗口大小参数按如下方法选择：在3-45之间，采用以4为等差依次取值，作为窗口大小参数，绘制建筑区与非建筑区的纹理特征值随窗口大小参数变化而变化的曲线图，通过分析曲线图的变化趋势得出最佳纹理窗口大小参数；所述步骤②的基于变差函数纹理分析方法中，纹理间距按如下方法确定：将h作为自变量从1开始不断递增，计算不同h下建筑区与非建筑区样本均值，绘制变差函数曲线图，通过分析曲线图求得最佳纹理间距。