CN103177458A - 一种基于频域分析的可见光遥感图像感兴趣区域检测方法 - Google Patents

Abstract

本发明公开一种基于频域分析的可见光遥感图像感兴趣区域检测方法，属于遥感图像处理及图像识别技术领域。实施过程包括：1)对图像进行预处理；2)将预处理结果进行四元数傅里叶变换，获得图像的频域信息；3)保留相位谱并用巴特沃斯高通滤波器获取幅度谱高频信息；4)对相位谱与幅度谱高频信息做四元数傅里叶反变换，得到特征图；5)利用高斯金字塔对特征图滤波与降维，获得显著图；6)对显著图进行阈值分割与二值化，获得感兴趣区域模板；7)将模板升维并与原始图像做掩膜操作，得到最终的感兴趣区域。本发明实现了感兴趣区域的快速、准确定位，具有区域描述精度高、计算复杂度低等优点，可用于环境监测、土地利用与农业调查等领域。

Description

一种基于频域分析的可见光遥感图像感兴趣区域检测方法

技术领域

本发明属于遥感图像处理及图像识别技术领域，具体涉及一种基于频域分析的可见光遥感图像感兴趣区域检测方法。

背景技术

遥感图像的感兴趣区域检测是遥感与信息交叉学科的重点、难点问题，无论在军事领域还是民用领域，都具有重要的理论研究意义与实际应用价值。传统的遥感图像感兴趣区域检测方法往往都是知识驱动的，经过图像分割与特征提取等过程，先验知识库与分割阈值的好坏直接决定了后续图像处理效果，同时，先验知识库的建立本身又是一个比较复杂的问题，一般需要综合考虑各种专家知识库、多种区域目标特征模型、背景区域特点等。视觉注意模型为遥感图像感兴趣区检测提供了一个全新的视角，不同于传统的检测方法，视觉注意模型完全由数据驱动，不涉及知识库等外部因素的影响。将视觉注意模型引入遥感图像感兴趣区域的检测具有重大的意义。近些年来，已有一些研究与分析将视觉注意模型引入遥感图像感兴趣区域检测中，并且成为一项提高大规模数据分析实时性和准确性的重要技术方法。这些方法模拟了人类视觉系统(human visual system，HVS)的视觉注意，它可以提供观察者可能的感兴趣区域信息，继而在该区域中寻找具体目标，这样可以帮助制定合理的计算资源分配方案，从而大幅地提升已有图像处理系统的运行效率。这种遥感图像感兴趣区域检测方法在保证准确检测结果的前提下减少了时间消耗与资源占用，具有重要的应用价值和实际意义。

在基于低层视觉特征的视觉注意模型方面，1998年，Itti等人在文章“A Model ofSaliency-Based Visual Attention for Rapid Scene Analysis”中通过模仿人类视觉自底向上的注意机制，运用金字塔图像不同级交叉相减的方式，计算出图像颜色、方向和强度三个通道的图，通过叠加三个通道图，形成显著图。在显著图中，一些重要目标能够被显著性区域突出。Itti模型同人类视觉系统十分接近，利用了各种视觉特性产生显著图。张鹏等人在文章“基于视觉注意的遥感图像分析方法”中根据人类视觉系统的生理结构与认知特点，基于Itti的视觉注意模型提出了一个遥感图像分析模型，并在此基础上开发出了一种新的遥感图像分析方法。该方法在一定程度上具备了类似于人类视觉系统的选择性加工能力。将其应用于多种类型的真实遥感图像，获得了较为满意的实验结果。

在基于数学方法的视觉注意模型方面，Harel等人在文章“Graph-Based Visual Saliency”中提出了基于图论的算法(Graph-based visual saliency，GBVS)。这种算法在特征提取的过程中类似Itti算法模拟视觉原理，在显著图生成的过程引入马尔可夫链(Markow chains)，用纯数学计算得到显著值。王鑫等人在文章“基于图像显著性区域的遥感图像机场检测”中针对已有方法对图像逐像素进行分析的不足，基于GBVS模型将人眼的注意力选择计算模型引入到遥感图像的机场目标检测中，提出一种基于图像显著性区域的遥感图像中机场目标检测与识别的方法，以提高自动目标检测的效率。

基于低层视觉特征的视觉注意模型较好地模拟了人眼视觉的关注方式，但是没有充分考虑图像的频域特征，同时计算速度慢、效率低，难以达到实时应用的要求。基于数学方法的视觉注意模型计算速度快，但是仅考虑了图像的自身特征，容易出现漏检、误检的现象。近些年Ell和Sangwine等人将四元数傅里叶变换应用到彩色图像中。其优势在于把彩色图像视为一个矢量整体，不仅获得了彩色图像频域信息的完整描述，而且保留了彩色图像各通道间的相互关系，可以用于可见光遥感图像的感兴趣区域检测中。

在感兴趣区域描述方面，传统方法中，用一个固定半径圆来描述感兴趣区域。它在识别随机区域时会带来大量冗余信息，不利于图像的压缩。使用单一阈值的速度非常快，但是感兴趣区域会有很多小碎片，区域中会有一些黑斑，区域的描述不够准确。最大类间方差法(Ostu方法)是一种自动的非参数非监督的阈值选择法，该方法可以说是自适应计算单阈值的简单高效方法。该方法具有计算简单、自适应强等优点。

发明内容

本发明的目的在于提供了一种基于频域分析的可见光遥感图像感兴趣区域检测方法，该方法用于对可见光遥感图像的感兴趣区域进行精确检测。现有的感兴趣区域检测方法主要依靠图像在空域中表现出来的各种特征完成检测，没有充分考虑图像的频域特性，针对可见光遥感图像的感兴趣区域进行检测时，对细节丰富且显著的地物信息无法完成精确检测。所以本发明方法主要关注两个方面：

1)提升可见光遥感图像感兴趣区域检测精度，获得更为准确的感兴趣区域信息；

2)降低可见光遥感图像感兴趣区域检测的计算复杂度。

本发明所使用的技术方案包括可见光遥感图像的显著图生成，感兴趣区域模板生成以及感兴趣区域生成三个主要过程，具体包括以下步骤：

步骤一：利用HSI(hue色调，saturation色饱和度，intensity亮度)彩色空间变换对可见光遥感图像进行预处理，即分别提取图像的R(red，红色)、G(green，绿色)、B(blue，蓝色)三个分量，利用HSI彩色空间变换将所获得的R、G、B三个分量转换到HSI彩色空间，获取可见光遥感图像的H、S、I三个分量。之所以进行彩色空间变换是由于HSI彩色空间与人类的颜色感知系统是一致的并且其优于RGB颜色空间，所以本发明将可见光遥感图像从RGB空间变换到了HSI空间。

步骤二：将步骤一得到的H、S、I三个分量作为输入变量，输入到四元数傅里叶变换中，计算后获得可见光遥感图像频域的幅度谱与相位谱信息；

步骤三：将频谱中的相位谱保留，采用巴特沃斯高通滤波器对可见光遥感图像的幅度谱进行高通滤波，得到幅度谱的高频信息；

步骤四：将相位谱与幅度谱的高频信息作为输入变量，输入到四元数傅里叶反变换中，计算后得到可见光遥感图像空间域的特征图。具体过程为：

1)将相位谱与幅度谱的高频信息作为输入变量，输入到四元数傅里叶反变换中；

2)计算四元数傅里叶反变换，在空间域构成输出图像；

3)将输出图像转化为灰度图像，得到可见光遥感图像空间域的特征图。

步骤五：采用高斯金字塔对步骤四得到的特征图进行高斯低通滤波，并进行空间降维，获得显著图；

步骤六：采用最大类间方差法对显著图进行阈值分割并对分割结果完成二值化操作，获得可见光遥感图像的感兴趣区域模板。具体过程为：

1)采用最大类间方差法得到显著图的最优分割阈值；

2)根据得到的阈值对图像进行二值化分割，感兴趣区域为1，非感兴趣区域为0；

2)得到二值化感兴趣区域模板。

步骤七：将感兴趣区域模板升维到与原始可见光遥感图像大小一致，然后基于该模板与原始图像进行掩膜操作，得到可见光遥感图像的感兴趣区域，即完成对原始可见光遥感图像的感兴趣区域检测。

本发明提出的方法具有如下的优点：

(1)本发明使用的四元数傅里叶变换可以采取快速傅里叶变换形式，解决了可见光遥感图像的计算效率问题。它可以满足实时计算的要求；

(2)本发明保留了四元数傅里叶变换后频域的相位谱，很好地凸显了图像中物体的位置信息，提升了感兴趣区域的定位准确度；

(3)本发明对四元数傅里叶变换后频域的幅度谱进行了巴特沃斯高通滤波，很好地保留了细节信息的能量，使阈值分割的结果更为精确；

(4)本发明利用最大类间方差法进行阈值分割。该方法能够在低计算复杂度下快速获得显著图的分割阈值，从而快速、高效地获得感兴趣区域边界信息的精确描述。

本发明的研究工作得到了国家自然科学基金项目(61071103)与中央高校基本科研业务费专项资金(2012LYB50)的资助。

附图说明

图1为本发明的流程图；

图2为本发明所使用的可见光遥感图像示例图片；

图3为本发明所使用示例图片的H、S、I三个分量的结果。(a)为H分量，(b)为S分量，(c)为I分量；

图4为本发明的显著图和感兴趣区域图。(a)为示例图片的显著图，(b)为示例图片的感兴趣区域图；

图5为示例图片采用本发明方法和其他方法生成的显著图的比较。(a)为Itti方法生成的显著图，(b)为GBVS方法生成的显著图，(c)为本发明方法生成的显著图；

图6为示例图片采用本发明方法和其他方法生成的感兴趣区域图的比较。(a)为Itti方法生成的感兴趣区域图，(b)为GBVS方法生成的感兴趣区域图，(c)为本发明方法生成的感兴趣区域图；

图7为示例图片的地面实况(Ground-Truth)显著度图；

图8为ROC曲线图，评价了Itti方法、GBVS方法与本发明方法在可见光遥感图像感兴趣区域检测中的性能。

具体实施方式

下面结合附图对本发明做进一步详细说明。本发明的总体框架如图1所示，现介绍每一步实现细节。

步骤一：利用HSI(hue色调，saturation色饱和度，intensity亮度)彩色空间变换对可见光遥感图像(如图2所示)进行预处理，即分别提取图像的R(red，红色)、G(green，绿色)、B(blue，蓝色)三个分量，利用HSI彩色空间变换将所获得的R、G、B三个分量转换到HSI彩色空间，获取可见光遥感图像的H、S、I三个分量：

$H=\left\{\begin{array}{cc}\mathrm{\θ}& (\mathrm{ifB}\≤G)\\ 360-\mathrm{\θ}& (\mathrm{ifB}>G)\end{array}\right.$

其中 $\mathrm{\θ}=\arccos \left\{\frac{\frac{1}{2}[(R-G)+(R-B)]}{{[{(R-G)}^2+(R-B)(G-B)]}^{1/2}}\right\}$

色饱和度分量S为：

$S=1-\frac{3}{(R+G+B)}\left[\min (R,G,B)\right]$

最后，亮度分量I为：

$I=\frac{1}{3}(R+G+B)$

假定RGB值已被归一化至[0，1]范围内，角度θ根据HSI空间的红轴来度量。可见光遥感图像的H、S、I三个分量的结果如图3所示，(a)为H分量，(b)为S分量，(c)为I分量。

步骤二：将步骤一得到的H、S、I三个分量作为输入变量，输入到四元数傅里叶变换中，计算后获得可见光遥感图像频域的幅度谱与相位谱信息。利用纯四元数，一幅可见光遥感图像可以用四元数形式表示为：

f(n，m)＝H(n，m)μ₁+S(n，m)μ₂+I(n，m)μ₃

其中(n，m)是每个像素点的位置。H(n，m)，S(n，m)和I(n，m)分别为每个像素的色调、色饱和度和亮度。μ_i，i＝1，2，3满足μ_i ²＝-1，μ₃＝μ₁μ₂，μ₁⊥μ₂，μ₂⊥μ₃，μ₁⊥μ₃。只要满足T.Ell和S.Sangwin在文章“Hypercomplex Fourier Transforms of Color Images”提出的四元数条件，选择任意的μ_i，i＝1，2，3，不会影响显著图的结果。

f(n，m)可以表示为以下形式(sympletic form)：

$\left\{\begin{array}{c}f(n,m)=f_1(n,m)+f_2(n,m){\mathrm{\μ}}_2\\ f_1(n,m)=H(n,m){\mathrm{\μ}}_1\\ f_2(n,m)=S(n,m)+I(n,m){\mathrm{\μ}}_1\end{array}\right.$

四元数傅里叶变换概括如下：

F[u，v]＝F₁[u，v]+F₂[u，v]μ₂

$F_i(u,v)=\frac{1}{\sqrt{\mathrm{MN}}}\underset{m=0}{\overset{M-1}{\mathrm{\Σ}}}\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}{e}^{-{\mathrm{\μ}}_{1}2\mathrm{\π}((\mathrm{mv}/M)+(\mathrm{nu}/N))}{f}_i(n,m)$

其中(n，m)和(u，v)表示每个像素的空间域和频率域的位置，N和M是图像的高度和宽度。f_i，i∈{1，2}是通过上述f(n，m)的表示形式计算得到的。

步骤三：将频谱中的相位谱保留，采用巴特沃斯高通滤波器对可见光遥感图像的幅度谱进行高通滤波，得到幅度谱的高频信息。利用极坐标，图像的频域表示F[u，v]可以写为：

$F[u,v]=\left|F\right[u,v\left]\right|e^{\mathrm{\μ\Φ}[u,v]}$

其中，e^μΦ＝cosΦ+μsinΦ(欧拉公式)，μ为一个单位纯四元数(μ²＝-1)。|F[u，v]|称为傅里叶变换的幅度谱或频率谱，Φ[u，v]称为相角或者相位谱。对幅度谱进行巴特沃斯高通滤波，保留相位谱信息。n阶且截至频率距原点为D₀的巴特沃斯型高通滤波器的传递函数由下式给出：

$H(u,v)=\frac{1}{1+{[D_0/D(u,v)]}^{2n}}$

其中D(u，v)＝[(u-M/2)²+(v-N/2)²]^1/2，是点(u，v)到傅里叶变换中心(原点)的距离。

步骤四：将相位谱与幅度谱的高频信息作为输入变量，输入到四元数傅里叶反变换中，计算后得到可见光遥感图像空间域的特征图。四元数傅里叶反变换为：

$f_i(n,m)=\frac{1}{\sqrt{\mathrm{MN}}}\underset{v=0}{\overset{M-1}{\mathrm{\Σ}}}\underset{u=0}{\overset{N-1}{\mathrm{\Σ}}}{e}^{{\mathrm{\μ}}_{1}2\mathrm{\π}((\mathrm{mv}/M)+(\mathrm{nu}/N))}{F}_i(u,v)$

然后在空间域构成输出图像：

f′(n，m)＝x(n，m)μ₁+y(n，m)μ₂+z(n，m)μ₃

将该图像转化为灰度图像，即为特征图。

步骤五：采用高斯金字塔对步骤四得到的特征图进行高斯低通滤波，并进行空间降维，获得显著图，如图4(a)所示。高斯金字塔分解会产生一系列分辨率逐步降低的图像。该低通滤波是通过与高斯滤波器内核卷积和一个下采样运算完成的。设原图为g₀，一共有C行R列像素。那么它的一阶金字塔滤波图像g₁的分辨率和采样密度都会比原图低。通过迭代该算法可以得到g₁的滤波图像g₂。以此类推，得到一系列图像。高斯金字塔算法可描述为：

g_k＝REDUCE(g_k-1)

式中，层数满足0＜k＜N，N表示高斯金字塔的总层数，点(i，j)满足0≤i＜C_k，0≤j＜R_k，C_k和R_k分别表示高斯金字塔结构中第k层图像的水平宽度和垂直高度。第k层图像在5×5窗口下可以认为是：

$g_k(i,j)=\underset{m=-2}{\overset{2}{\mathrm{\Σ}}}\underset{n=-2}{\overset{2}{\mathrm{\Σ}}}w(m,n)g_{k-1}(2i+m,2j+n)$

其中w(m，n)称为生成核，将其分离变量写作

需要满足归一性、对称性和均匀分布性。在这三个约束下，可以用一个参数a来表示：

$\left\{\begin{array}{c}\hat{w}\left(0\right)=a\\ \hat{w}(-1)=\hat{w}\left(1\right)=\frac{1}{4}\\ \widehat{w(-2)}=\hat{w}\left(2\right)=\frac{1}{4}-\frac{a}{2}\end{array}\right.$

当参数a＝0.4时，生成核的等效权值函数同高斯函数特别相似。所以这个多尺度分解算法也称为高斯金字塔分解。

步骤六：采用最大类间方差法对显著图进行阈值分割并对分割结果完成二值化操作，获得可见光遥感图像的感兴趣区域模板。设图像像素总量为N，灰度范围为[1，2，…，L]。灰度为i的像素总数为n_i。则其出现的概率为：

p_i＝n_i/N(i＝1，2，...，L)

$\underset{i=1}{\overset{L}{\mathrm{\Σ}}}{p}_i=1$

通过一个灰度阈值t，这些像素被分为A、B两类(前景和背景)。A的灰度范围为[1，…，k]，B的灰度范围为[k+1，…，L]。它们的概率分别为：

$\left\{\begin{array}{c}{\mathrm{\ω}}_A=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{p}_i=\mathrm{\ω}\left(k\right)\\ {\mathrm{\ω}}_B=\underset{i=k+1}{\overset{L}{\mathrm{\Σ}}}{p}_i=1-\mathrm{\ω}\left(k\right)\end{array}\right.$

这两类的灰度均值为：

$\left\{\begin{array}{c}{\mathrm{\λ}}_A=\underset{i=i}{\overset{k}{\mathrm{\Σ}}}i{p}_i/{\mathrm{\ω}}_A=\frac{\mathrm{\λ}\left(k\right)}{\mathrm{\ω}\left(k\right)}\\ {\mathrm{\λ}}_B=\underset{i=k+1}{\overset{L}{\mathrm{\Σ}}}i{p}_i/{\mathrm{\ω}}_B=\frac{{\mathrm{\λ}}_T-\mathrm{\λ}\left(k\right)}{1-\mathrm{\ω}\left(k\right)}\end{array}\right.$

其中

λ_T是整幅图像的平均灰度值，它满足λ_T＝ω_Aλ_A+ω_Bλ_B

这两个类的方差为：

${\mathrm{\σ}}^2\left(k\right)==\frac{{[{\mathrm{\λ}}_T\mathrm{\ω}\left(k\right)-\mathrm{\λ}\left(k\right)]}^2}{\mathrm{\ω}\left(k\right)[1-\mathrm{\ω}\left(k\right)]}$

使用最大类间方差的准则，让k从1到L范围依次取值，使σ²最大的k值即为最大类间方差法的最优阈值，则最优阈值k^*为：

${\mathrm{\&sigma;}}^2\left(k^{*}\right)=\underset{1\&le;k<L}{\max }{\mathrm{\&sigma;}}^2\left(k\right)$

利用最大类间方差法(Ostu法)找到适合的分割阈值后，根据得到的阈值对图像进行二值化分割，感兴趣区域为1，非感兴趣区域为0。

步骤七：将感兴趣区域模板升维到与原始可见光遥感图像大小一致，然后基于该模板与原始图像进行掩膜操作，得到可见光遥感图像的感兴趣区域，如图4(b)所示，即完成对原始可见光遥感图像的感兴趣区域检测。

本发明的效果可通过以下实验结果与分析进一步说明：

1.实验数据

本发明从SPOT5卫星源图中选取了一幅北京某郊区的可见光遥感图像，并分别从中截取生成了大小为2048×2048的一幅图作为本文实验源图，如图2所示。

2.对比实验及实验评价指标

本发明方法对比了传统的Itti方法和GBVS方法。从主观上分别对比了不同方法生成的显著图和感兴趣区域图，如图5和图6所示。图5中，(a)为Itti方法生成的显著图，(b)为GBVS方法生成的显著图，(c)为本发明方法生成的显著图。图6中，(a)为Itti方法生成的感兴趣区域图，(b)为GBVS方法生成的感兴趣区域图，(c)为本发明方法生成的感兴趣区域图。

本发明还使用了ROC(Receiver Operating Characteristic)曲线(又称受试者工作特征曲线)从客观上评价了这几种感兴趣区域检测方法。ROC常被用来衡量二值分类器在不同判断标准下的效果，是很多算法中用来评价视觉注意模型的常用指标。其常用形式有ROC曲线和ROC曲线下的面积。ROC曲线表征了分类器在不同的判断标准下产生的真阳性率(TruePositive Rate，简称：TPR)和假阳性率(False Positive Rate，简称：FPR)。其基本思想是使用不同的阈值，从预测的显著图中选出感兴趣对象，并使用地面实况(Ground-Truth)显著度图(即标定好的真实数据)进行验证。假设已确定一个阈值，预测显著度大于阈值的对象记为感兴趣对象，预测显著度小于这个值的对象为非感兴趣对象，进而得到以下四个值：

真阳性＝预测为感兴趣对象&实际为注视点

真阴性＝预测为非感兴趣对象&实际为非注视点

假阳性＝预测为感兴趣对象&实际为非注视点

假阴性＝预测为非感兴趣对象&实际为注视点

由此可以计算真阳性率(TPR)和假阳性率(FPR)，计算公式为：

当阈值很大时，三种方法都无法识别出实际的感兴趣对象，因此TPR和FPR在阈值很大时基本一致为0；当阈值非常小时，三种方法将显著图中的对象均预测为感兴趣对象，因此不存在预测为非感兴趣对象的情况，所以TPR和FPR在阈值很小时也基本一致为1；在0-1期间，在相同FPR值下，TPR越大越好；在相同TPR值下，FPR越小越好。图7为地面实况(Ground-Truth)显著度图。图8为ROC曲线图。从图中可以看出本发明的方法性能接近且略优于Itti方法，明显优于GBVS方法。

Claims (3)

1.一种基于频域分析的可见光遥感图像感兴趣区域检测方法，该方法用于检测可见光遥感图像得到精确的感兴趣区域，其特征在于，包括以下步骤：

步骤一：利用HSI(hue色调，saturation色饱和度，intensity亮度)彩色空间变换对可见光遥感图像进行预处理，即分别提取图像的R(red，红色)、G(green，绿色)、B(blue，蓝色)三个分量，利用HSI彩色空间变换将所获得的R、G、B三个分量转换到HSI彩色空间，获取可见光遥感图像的H、S、I三个分量；

步骤二：将步骤一得到的H、S、I三个分量作为输入变量，输入到四元数傅里叶变换中，计算后获得可见光遥感图像频域的幅度谱与相位谱信息；

步骤三：将频谱中的相位谱保留，采用巴特沃斯高通滤波器对可见光遥感图像的幅度谱进行高通滤波，得到幅度谱的高频信息；

步骤四：将相位谱与幅度谱的高频信息作为输入变量，输入到四元数傅里叶反变换中，计算后得到可见光遥感图像空间域的特征图；

步骤五：采用高斯金字塔对步骤四得到的特征图进行高斯低通滤波，并进行空间降维，获得显著图；

步骤六：采用最大类间方差法对显著图进行阈值分割并对分割结果完成二值化操作，获得可见光遥感图像的感兴趣区域模板；

步骤七：将感兴趣区域模板升维到与原始可见光遥感图像大小一致，然后基于该模板与原始图像进行掩膜操作，得到最终的可见光遥感图像的感兴趣区域，即完成对原始可见光遥感图像的感兴趣区域检测。

2.根据权利要求1所述的一种基于频域分析的可见光遥感图像感兴趣区域检测方法，其特征在于，所述步骤四的具体过程为：

1)将相位谱与幅度谱的高频信息作为输入变量，输入到四元数傅里叶反变换中；

2)计算四元数傅里叶反变换，在空间域构成输出图像；

3)将输出图像转化为灰度图像，得到可见光遥感图像空间域的特征图。

3.根据权利要求1所述的一种基于频域分析的可见光遥感图像感兴趣区域检测方法，其特征在于，所述步骤六的具体过程为：

1)采用最大类间方差法得到显著图的最优分割阈值；

2)根据得到的阈值对图像进行二值化分割，感兴趣区域为1，非感兴趣区域为0；

2)得到二值化感兴趣区域模板。

Images