CN102129576B - 全天空极光图像占空比参数的提取方法 - Google Patents

Abstract

本发明公开了一种全天空极光图像占空比参数的提取方法，主要解决现有技术只能对具有明显形状和边界的光斑进行分割导致占空比参数误差较大的缺陷。其特定是根据极光在天空所呈现的不同形态采用相应分割方法，即首先对全天空极光图像进行预处理；然后确定光斑区域最优分割阈值，对光斑区域进行分割；接着根据光线区域与背景天空区域纹理结构的差异，对光线区域进行分割；最后根据光斑和光线区域的分割结果，计算占空比参数。本发明通过灰度特征确定光斑分割阈值，对无明显形状或边界的光斑也能精确分割，并根据纹理特征实现光线区域的成功分割，有效降低了占空比参数的误差，可用于复杂形态的极光区域从天空背景区域的分割。

Description

全天空极光图像占空比参数的提取方法

技术领域

本发明属于图像处理技术领域，涉及特征提取和图像分割，可用于复杂形态的极光区域从天空背景区域的分割。

背景技术

空间科学的研究重点在于了解近地球周围环境的物理过程的语义，这些过程由太阳风等离子体、地球磁场和地球大气层相互之间的复杂作用所引起。极光就是由太阳风引起的变化而带来的能量所产生的主要现象之一，因此极光所引起的空间和时间上的变化在很大程度上反映了太阳风与地球磁场的变化，对极光的观察成为了研究地球科学的重要基础。早期主要通过一些物理属性，例如太阳风速度、电子密度等对极光进行研究，实际上，除了以上物理特征，极光在亮度、结构等方面的变化也是极光研究的一个重要方面。

随着数字全天空摄像仪的出现，全天空极光图像为极光现象的研究提供了新的研究途径。其中极光区域在全天空极光图像所占面积的比例，即占空比成为极光研究的重要参数。占空比参数在时间轴上的变化蕴含了极光变化的丰富信息，为了获得极光的占空比参数，首先要将极光区域从整个全天空图像中分割出来。为了准确的分割出极光区域，要对极光的形态进行深入的分析，对不同形态的极光区域采用不同的分割方法。

极光的形态与极光所属类型密切相关，中国极地中心的科学家以日侧冕状极光为研究重点，将极光分为四类：多弧状极光、热点冕状极光、辐射冕状极光和帷幔冕状极光。以上四类极光主要呈现出以下三种形态：

·光斑状：该类极光只包含光带或光斑，不含或几乎不含光线，通常是弧状和热点冕状极光；

·光线状：该类极光只包含光线，几乎没有光斑出现，通常以帷幔冕和辐射冕状极光为主；

·混合形态：该类极光不仅包含光斑还包含光线，多弧状极光、热点冕状极光、帷幔冕状和辐射冕状极光都有这种情况出现，在极光图像中比较普遍。

近年来图像分割领域非常活跃，有很多方法被提出并在相应的领域获得了良好的分割效果，但在极光图像的分割应用上却很少，这是因为极光在天空背景上是透明的、边缘模糊、形状不规则，而且变化速度快且变化形式多样，因此对极光的描述非常困难。

在1999年，Syrjasuo等人在文献“M.T.Syrjasuo，T.I.Pulkkinen.Determining theSkeletons of the Auroras.International Conference on Image Analysis and Processing，Venice，Italy pp.1063-1066，1999.”中为了勾画极光的骨架，在预处理阶段使用Fuzzy logic的方法对具有明显形状的极光区域进行了分割，他们通过对极光边缘的提取从而将极光区域与背景天空区域分割开，这项技术要求极光区域要有明显的边缘，对于具有明显光斑和光带的极光图像，该方法的效果较好，但对于边缘模糊的光斑尤其是没有明显边缘的极光光线部分，分割效果并不好。

在2004年，Syrjasuo等人在文献“M.T.Syrjasuo，E.F.Donovan and L.L.Cogger.Content-based retrieval of auroral images-thousands of irregular shapes.InternationalConference on Visualization，Imaging，and Image Processing，Marbella，Spain，pp.224-228，2004.”以及2005年在文献“M.T.Syrjasuo，and E.F.Donovan.Using relevance feedbackin retrieving auroral images.International Conference on Computational Intelligence，CalgaryAlberta，Canada，pp.420-425，2005.”中使用Isolable-contour map对一种亮度很高的特殊弧状极光即南北走向弧状极光进行分割，并提取分割结果的边缘特征该类极光的检索。Isolable-contour map方法要求被检测区域要有明确的形状，因此该方法适合有明显形状的且亮度较高的极光区域，对于亮度较低、前景与背景没有明显区域差别的光斑以及很难描述形状的光线区域，该方法无法获得理想的分割结果。

综上，已有的极光图像分割方法主要存在以下两个问题：

1)Fuzzy logic和Isolable-contour map适用于图像中明显光斑区域的分割，对于边缘模糊或没有明显形状的光斑区域，由于无法提取其边缘或形状特征，无法实现精确分割，导致光斑部分占空比参数不准确，从而使整个极光区域占空比参数的产生误差。

2)现有的极光分割技术只考虑了亮度较高的光斑或光带的分割，而多数极光是光斑和光线的混合体，甚至有些极光的光线区域在整个极光区域中占很大比例，现有技术无法实现光线区域的分割，从而导致极光光线区域分割的缺失，造成分割结果的不准确，进而导致占空比参数的不准确。

发明内容

本发明的目的在于克服上述现有方法的不足，提出一种基于极光不同形态的全天空极光图像占空比参数的提取方法，以提高极光区域分割的准确率，减少占空比参数的误差。

实现本发明的技术思路是：首先采用改进的Otsu方法对光斑区域进行分割，得到光斑区域分割结果，然后对图像的光线区域进行分割，接着将上面两部分的分割结果进行组合形成极光整体区域的分割结果，最后计算极光区域在整个天空区域的比例得到极光占空比参数，具体实现步骤如下：

(1)对数据库中所有全天空极光图像进行预处理：将每一幅原始大小为512×512的全天空极光图像用半径为220的圆形区域进行掩模，去掉周围干扰灯光、山脉噪声，掩模成大小为440×440的掩模图像，所有掩模图像构成掩模图像集S_M；

(2)输入一幅掩模图像I_M，使用改进的Otsu方法对极光光斑区域进行分割，得到光斑区域分割结果R_P；

(3)从掩模图像集S_M中选择100幅不含极光或几乎不含极光的掩模图像构建纹理阈值训练集S_K，使用自适应局部二值模式ALBP特征向量提取算法，提取纹理阈值训练集S_K中每一幅图像的ALBP特征向量，构成ALBP特征向量集A＝{α_i}，其中α_i表示纹理阈值训练集S_K中第i幅掩模图像的ALBP特征向量；

(4)对纹理阈值训练集S_K的ALBP特征向量集A进行训练，得到参考向量R和纹理分割阈值Z；

(5)将掩模图像I_M分割成大小为16×16的图像块，构成图像块集X＝{x_i，j}，提取每一块的ALBP特征向量，构成ALBP特征向量集B＝{β_i，j}，其中x_i，j表示掩模图像I_M中第i行第j列的图像块，β_i，j表示图像块x_i，j的ALBP特征向量；

(6)计算ALBP特征向量β_i，j与参考向量R的相似度S_i，j：

$S_{i,j}=X^2({\mathrm{\β}}_{i,j},R)=\underset{l=1}{\overset{C}{\mathrm{\Σ}}}\frac{{({\mathrm{\β}}_{i,j}\left(l\right)-R\left(l\right))}^2}{{\mathrm{\β}}_{i,j}\left(l\right)+R\left(l\right)}$

其中X²(β_i，j，R)表示使用Chi-square统计量计算ALBP特征向量β_i，j与参考向量R之间的相似度，C表示ALBP特征向量的维数共有C维，l表示第l维，β_i，j(l)表示ALBP特征向量的第l维向量，R(l)表示参考向量R的第l维向量；

(7)将相似度S_i，j与纹理分割阈值Z进行比较，当相似度S_i，j大于或等于纹理分割阈值Z，则图像块x_i，j为极光光线区域，否则图像块x_i，j为背景天空区域，得到掩模图像I_M的光线区域分割结果R_T；

(8)对光斑区域分割结果R_P和光线区域分割结果R_T求并集，得到掩模图像I_M的完整分割结果R_A，求光斑区域分割结果R_P和光线区域分割结果R_T的交集，得到光斑区域分割结果R_P和光线区域分割结果R_T的重合部分R_X，将该重合部分R_X去掉，得到修正后的光线区域分割结果R′T＝R_T-R_X；

(9)根据极光区域的完整分割结果R_A、修正后的光线区域分割结果R′_T和光斑区域的分割结果R_P，计算掩模图像I_M的极光占空比参数：O_A，O_T，O_P；

$O_A=\underset{(i,j)\∈R_M}{\mathrm{\Σ}}{y}_A(i,j)/N^{\′},y_A(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\∈R_A\\ 0,& p(i,j)\∉R_A\end{array}\right.$

其中：O_A是完整极光区域的占空比参数，表示掩模图像I_M中以各种形态出现的极光在整个天空背景中所占比例，(i，j)是掩模图像I_M中第i行第j列的坐标，(i，j)∈R_M表示只选择在圆形区域R_M内的坐标(i，j)，p(i，j)是坐标(i，j)所对应的像素，y_A(i，j)表示像素p(i，j)是否属于极光区域的完整分割结果R_A，如果像素p(i，j)属于完整分割结果R_A，则y_A(i，j)的值为1，否则y_A(i，j)的值为0，N′是圆形区域R_M中像素总数；

$O_T=\underset{(i,j)\∈R_M}{\mathrm{\Σ}}{y}_T(i,j)/N^{\′},y_T(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\∈R_T^{\′}\\ 0,& p(i,j)\∉R_T^{\′}\end{array}\right.$

其中：O_T表示极光光线区域的占空比参数，是掩模图像I_M中以光线形态所呈现的极光在整个天空背景中所占比例，y_T(i，j)用于判断像素p(i，j)是否属于光线区域分割结果R′_T，如果像素p(i，j)属于光线区域分割结果R′_T，则y_T(i，j)的值为1，否则y_T(i，j)的值为0；

$O_P=\underset{(i,j)\∈R_M}{\mathrm{\Σ}}{y}_P(i,j)/N^{\′},y_P(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\∈R_P\\ 0,& p(i,j)\∉R_P\end{array}\right.$

其中：O_P是极光光斑区域的占空比参数，表示掩模图像I_M中以光斑或光带形态出现的极光在整个天空背景中所占比例，y_P(i，j)表示像素p(i，j)是否属于光斑区域分割结果R_P，如果像素p(i，j)属于光斑区域分割结果R_P，则y_P(i，j)的值为1，否则y_P(i，j)的值为0。

本发明与已有的极光图像分割算法相比，具有以下优点：

(1)本发明通过前后背景的灰度特征确定光斑区域分割阈值，无需光斑区域的形状特征或边界特征，对于无明显形状或者边界模糊的光斑，能够实现精确的分割，从而有效降低了光斑区域占空比参数的误差；

(2)本发明不仅考虑到极光光斑区域的分割，而且充分考虑到在多数极光图像中频繁出现的光线区域，针对光线区域特殊的形态特征，提出通过光线区域与背景天空区域纹理结构差异来实现极光光线区域的分割，弥补了现有方法在极光光线区域分割方法的缺失，提高了完整极光区域占空比参数的准确性；

(3)仿真结果表明，与现有极光分割技术相比，本发明根据极光在天空所呈现的不同形态采用相应的分割方法，能够更加细致准确的将极光光斑区域和光线区域从背景天空区域中分离出来，有效降低了占空比参数的误差。

本发明的技术过程和效果可结合以下附图详细说明。

附图说明

图1是本发明全天空极光图像占空比参数的提取总流程图；

图2是本发明中对光斑区域分割的子流程图；

图3是本发明中对光线区域分割的子流程图；

图4是本发明与现有方法对含光斑和光线的弧状极光仿真分割的结果比较图；

图5是本发明与现有方法对含光斑和光线的热点冕状极光仿真分割的结果比较图；

图6是本发明与现有方法对含光斑和光线的帷幔冕状极光仿真分割的结果比较图；

图7是本发明与现有方法对只含光线的帷幔冕状极光仿真分割的结果比较图；

图8是本发明与现有方法对含光斑和光线的辐射冕状极光仿真分割的结果比较图；

图9是本发明与现有方法对只含光线的辐射冕状极光仿真分割的结果比较图。

具体实施方式

参照图1，本发明全天空极光图像占空比参数的提取方法主要包括四个阶段：图像预处理、光斑区域分割、光线区域分割和占空比参数计算，具体步骤如下：

步骤1：对数据库中所有全天空极光图像进行预处理：将每一幅原始大小为512×512的全天空极光图像用半径为220的圆形区域进行掩模，去掉周围干扰灯光、山脉噪声，掩模成大小为440×440的掩模图像，所有掩模图像构成掩模图像集S_M。

步骤2：参照图2，计算掩模图像I_M中圆形区域R_M内灰度不同像素的发生概率，确定光斑区域最优分割阈值G^*，通过与最优光斑区域最优分割阈值G^*的比较，得到光斑区域的分割结果。

2.1)计算圆形区域R_M内灰度不同像素的发生概率

2.1.1)只保留掩模图像I_M中圆形区域R_M，计算该圆形区域内像素总数N′：

N′＝N-(c²-π(c/2)²)

其中：N是掩模图像I_M中像素总数，c是掩模图像I_M的边长，即圆形区域R_M的直径，c²-π(c/2)²是圆形区域R_M外四个角中像素总数，这四个角中像素灰度均为0；

2.1.2)计算圆形区域R_M中灰度为0的像素个数n′₀：

n′₀＝n₀-(c²-π(c/2)²)

其中：n₀是掩模图像I_M中灰度为0的像素个数；

2.1.3)根据圆形区域R_M内像素总数N′和灰度为0的像素个数n′₀，计算灰度为0的像素的发生概率g₀：

$g_0=\frac{n_0^{\′}}{N^{\′}}$

2.1.4)根据圆形区域R_M内像素总数N′和灰度为k的像素的个数n_k，计算灰度为k的像素的发生概率g_k：

$g_k=\frac{n_k}{N^{\′}},k=\mathrm{1,2}......M$

其中：M是圆形区域R_M内灰度的最大值；

2.2)根据不同灰度的像素在圆形区域R_M内发生的概率，使用Otsu方法计算前景类和背景类像素的最优分割阈值G^*：

2.2.1)假设G是圆形区域R_M中前景类像素和背景类像素的分割阈值，设G初值为0，假设max是前景类与背景类之间的最大类间方差，将max初始化为0；

2.2.2)在圆形区域R_M中，将所有像素分为两类，灰度大于或等于分割阈值G的像素，被划分为前景类像素，灰度小于分割阈值G的像素，则被划分为背景类像素，其属于前景类的像素和属于背景类的像素的发生概率分别为：

$p_B=\underset{k=0}{\overset{G}{\mathrm{\Σ}}}{g}_k$

$p_F=\underset{k=G+1}{\overset{M}{\mathrm{\Σ}}}{g}_k$

其中p_B是背景类像素的发生概率，p_F是前景类像素的发生概率，k表示像素的灰度，M是圆形区域R_M内灰度的最大值；

2.2.3)分别计算前景类像素的灰度均值：

，背景类像素的灰度均值：

和圆形区域R_M中像素的灰度均值：m＝p_Bm_B+p_Fm_F；

2.2.4)计算前景类和背景类之间的类间方差σ²：

σ²＝p_F(m_F-m)²+p_B(m_B-m)²

将类间方差σ²与最大类间方差max进行对比，如果σ²大于max，则将σ²的值赋值给最大类间方差max，并且将分割阈值G赋值给最优分割阈值G^*，否则，最大类间方差max和最优分割阈值G^*均保持不变；

2.2.5)给分割阈值G加1，即G′＝G+1，将G′与圆形区域R_M中灰度最大值M相比较，如果G′小于或等于灰度最大值M，则将G′作为新的分割阈值带入步骤2.2.2)；如果G′大于灰度最大值M，则循环完毕，得到最优分割阈值G^*；

2.3)根据最优分割阈值G^*，完成极光光斑区域分割：将掩模图像I_M中每个像素与光斑分割阈值G^*进行比较，如果大于或等于G^*，则将该像素划分为前景类，如果小于G^*则该像素被划分为背景类，得到掩模图像I_M的光斑分割结果R_P。

步骤3：参照图3，对特征训练集S_T进行训练，构建主要概率模式集P，根据主要概率模式集P，提取ALBP特征向量，确定参考向量R和纹理分割阈值Z，根据纹理分割阈值Z，得到光线区域的分割结果R_T。

3.1)构建主要概率模式集P：

3.1.1)从掩模图像集S_M中选择100幅含有大量光线区域的掩模图像构成特征训练集S_T，对特征训练集S_T中每一幅掩模图像从左上角开始，按从左到右、从上到下的顺序，取大小为5×5的正方形邻域作为处理单元，对于每个处理单元，以中心像素的灰度作为阈值，邻域中像素的灰度与中心像素灰度进行对比，如果邻域像素的灰度大于中心像素的灰度，则该邻域像素被置为1，否则置0，邻域的所有像素转换为一串0～1交替的二进制序列，为二进制序列中每个元素分配一个权重因子2ⁱ，则该二进制序列转换为一个确切的十进制数，即为中心像素的LBP特征向量：

$\mathrm{LBP}(Q,D)=\underset{i=0}{\overset{Q-1}{\mathrm{\&Sigma;}}}{b}_i{2}^i,b_i=\left\{\begin{array}{cc}1,& p_i\&GreaterEqual;p_c\\ 0,& p_i<p_c\end{array}\right.$

其中：Q表示邻域中元素个数，取值为16，D表示方形邻域的边长，取值为5，i表示邻域中第i个像素，p_i表示邻域中第i个像素的灰度，p_c表示中心像素灰度值，b_i表示邻域中第i个像素的灰度p_i是否大于或等于中心像素灰度p_c，当b_i为1时表示邻域中第i个像素的灰度p_i大于或等于中心像素灰度p_c，b_i为0表示邻域中第i个像素的灰度p_i小于中心像素灰度p_c；

根据LBP特征向量，提取旋转不变LBP特征向量：

LBP^ri(Q，D)＝f_m(f_c(LBP(Q，D)，j))，j＝0，1，...Q-1

其中：f_c(LBP(Q，D)，j)表示二进制序列LBP(Q，D)顺时针旋转j次，j取值为从0到Q-1的整数，f_n(f_c(LBP(Q，D)，j))表示从所有取值中选择最小值；

3.1.2)将特征训练集S_T中所有掩模图像的旋转不变LBP特征向量取均值，得到平均旋转不变LBP特征向量V；

3.1.3)按照模式发生概率从大到小的顺序，对平均旋转不变LBP特征向量V进行排序得到排序后的平均旋转不变LBP特征向量V′，选择该向量V′中的前C个模式，这C个模式的概率值之和Sum满足大于或等于90％，用这C个模式组成主要概率模式集P；

3.2)对于输入的掩模图像I，提取I的旋转不变LBP特征向量r_M，根据主要概率模式集P，保留r_M中与其对应的C个模式及其概率值，即为掩模图像I的ALBP特征向量；

3.3)确定参考向量R和光线区域分割阈值Z：

3.3.1)掩模图像集S_M中选择100幅不含极光或几乎不含极光的掩模图像构建纹理阈值训练集S_K，使用ALBP特征向量提取算法，提取纹理阈值训练集S_K中每一幅图像的ALBP特征向量，构成ALBP特征向量集A＝{α_i}，其中α_i表示纹理阈值训练集S_K中第i幅掩模图像的ALBP特征向量；

对该训练集S_K中所有掩模图像的ALBP特征求平均值，得到参考向量R：

$R=\underset{i=1}{\overset{N_R}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i/N_R;$

其中：N_R表示纹理阈值训练集S_K中掩模图像的数量，N_R为100；

3.3.2)将训练集中每一幅掩模图像分割为16×16大小的块，对每个块提取ALBP特征向量，计算每个块ALBP特征向量与参考向量R之间相似度：

$L_{i,j}=X^2(W_{i,j},R)=\underset{l=1}{\overset{C}{\mathrm{\&Sigma;}}}\frac{{(W_{i,j}\left(l\right)-R\left(l\right))}^2}{W_{i,j}\left(l\right)+R\left(l\right)}$

其中：W_i，j表示第i幅图像第j个块的ALBP特征向量，l是ALBP特征向量中的第l维，W_i，j(l)表示ALBP特征向量W_i，j的第l维向量，R(l)表示参考向量R的第l维向量，C是ALBP特征向量的维数；

3.3.3)从所有块的ALBP特征向量与参考向量R之间的相似度L_i，j中找到最大值，即为纹理分割阈值：

$Z=\underset{i,j}{\max }\left(L_{i,j}\right).$

3.4)根据光线区域分割阈值Z，得到光线区域分割结果R_T：

3.4.1)将掩模图像I_M分割成大小为16×16的图像块，构成图像块集X＝{x_i，j}，提取每一块的ALBP特征向量，构成ALBP特征向量集B＝{β_i，j}，其中x_i，j表示掩模图像I_M中第i行第j列的图像块，β_i，j表示图像块x_i，j的ALBP特征向量；

3.4.2)计算ALBP特征向量β_i，j与参考向量R的相似度S_i，j：

$S_{i,j}=X^2({\mathrm{\&beta;}}_{i,j},R)=\underset{l=1}{\overset{C}{\mathrm{\&Sigma;}}}\frac{{({\mathrm{\&beta;}}_{i,j}\left(l\right)-R\left(l\right))}^2}{{\mathrm{\&beta;}}_{i,j}\left(l\right)+R\left(l\right)}$

其中：X²(β_i，j，R)表示使用Chi-square统计量计算ALBP特征向量β_i，j与参考向量R之间的相似度，C表示ALBP特征向量的维数共有C维，l表示第l维，β_i，j(l)表示ALBP特征向量的第l维向量，R(l)表示参考向量R的第l维向量；

3.4.3)将相似度S_i，j与纹理分割阈值Z进行比较，当相似度S_i，j大于或等于纹理分割阈值Z，则图像块x_i,j为极光光线区域，否则图像块x_i，j为背景天空区域，得到掩模图像I_M的光线区域分割结果R_T。

步骤4.根据光斑区域、光线区域和完整极光区域的分割结果，计算占空比参数。

4.1)对光斑区域分割结果R_P和光线区域分割结果R_T求并集，得到掩模图像I_M的完整分割结果R_A，求光斑区域分割结果R_P和光线区域分割结果R_T的交集，得到光斑区域分割结果R_P和光线区域分割结果R_T的重合部分R_X，将该重合部分R_X去掉，得到修正后的光线区域分割结果：R’_T＝R_T-R_X；

4.2)计算占空比参数

根据极光区域的完整分割结果R_A、修正后的光线区域分割结果R′_T和光斑区域的分割结果R_P，计算掩模图像I_M的极光占空比参数：O_A，O_T，O_P；

$O_A=\underset{(i,j)\&Element;R_M}{\mathrm{\&Sigma;}}{y}_A(i,j)/N^{\&prime;},y_A(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\&Element;R_A\\ 0,& p(i,j)\&NotElement;R_A\end{array}\right.$

其中：O_A是完整极光区域的占空比参数，表示掩模图像I_M中以各种形态出现的极光在整个天空背景中所占比例，(i，j)是掩模图像I_M中第i行第j列的坐标，(i，j)∈R_M表示只选择在圆形区域R_M内的坐标(i，j)，p(i，j)是坐标(i，j)所对应的像素，y_A(i，j)表示像素p(i，j)是否属于极光区域的完整分割结果R_A，如果像素p(i，j)属于完整分割结果R_A，则y_A(i，j)的值为1，否则y_A(i，j)的值为0，N′是圆形区域R_M中像素总数；

$O_T=\underset{(i,j)\&Element;R_M}{\mathrm{\&Sigma;}}{y}_T(i,j)/N^{\&prime;},y_T(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\&Element;R_T^{\&prime;}\\ 0,& p(i,j)\&NotElement;R_T^{\&prime;}\end{array}\right.$

其中：O_T表示极光光线区域的占空比参数，是掩模图像I_M中以光线形态所呈现的极光在整个天空背景中所占比例，y_T(i，j)用于判断像素p(i，j)是否属于光线区域分割结果R′_T，如果像素p(i，j)属于光线区域分割结果R′_T，则y_T(i，j)的值为1，否则y_T(i，j)的值为0；

$O_P=\underset{(i,j)\&Element;R_M}{\mathrm{\&Sigma;}}{y}_P(i,j)/N^{\&prime;},y_P(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\&Element;R_P\\ 0,& p(i,j)\&NotElement;R_P\end{array}\right.$

其中：O_P是极光光斑区域的占空比参数，表示掩模图像I_M中以光斑或光带形态出现的极光在整个天空背景中所占比例，y_P(i，j)表示像素p(i，j)是否属于光斑区域分割结果R_P，如果像素p(i，j)属于光斑区域分割结果R_P，则y_P(i，j)的值为1，否则y_P(i，j)的值为0。

本发明的优点可通过以下仿真实验进一步说明：

本发明所使用的数据来自2003年12月到2004年1月由北极黄河站所采集的全天空极光图像。为了验证本发明的有效性，仿真实验所采用的数据均经过中国极地中心的专家进行人工标记和分割，便于本发明分割仿真的结果与专家分割结果进行比较。

本发明在大量全天空极光图像上进行了测试实验，其中包括弧状极光、帷幔冕状极光、热点冕状极光和辐射冕状极光，并且每一种类型都选择多幅含有不同形态极光的图像，以对各种类型极光的不同形态进行分割仿真，来验证本发明的有效性。

本发明与现有Fuzzy logic和Isolable-contour map方法在相同的数据集上将极光区域从天空背景区域中分割出来，对三种方法分割仿真的结果分别进行：以人眼感知为评分依据的主观视觉评价实验和以错误率γ为评价依据的客观数据评价。

实验1：主观视觉评价

由于人眼非常敏感，能够观察到两幅图之间最细微的差别，因此人眼对分割结果的评价非常重要。在本实验中，邀请12人对全天空极光图像分割仿真的结果进行评分，由于极光图像不同于自然图像，评分人员是否具有专业背景知识对分割仿真结果的评价有一定影响，因此实验由两组人员完成，专家组由6位具有极光背景知识的专家组成，无经验组由6位没有专业背景知识的普通人组成。对分割仿真的结果的评价采用评分的形式来实现：

其中：U表示对分割仿真的结果的评价标准，分数越高表示对分割仿真的结果越满意。

对多弧状极光、热点冕状极光、帷幔冕状极光和辐射冕状极光四种类型，每类100幅图像，分别使用本发明方法与现有Fuzzy logic和Isolable-contour map方法进行分割仿真。分割仿真的结果如图4所示，其中：图4是本发明与现有方法对含光斑和光线的弧状极光分割仿真的结果，图5是本发明与现有方法对含光斑和光线的热点冕状极光分割仿真的结果，图6是本发明与现有方法对含光斑和光线的帷幔冕状极光分割仿真的结果，图7是本发明与现有方法对只含光线的帷幔冕状极光分割仿真的结果，图8是本发明与现有方法对含光斑和光线的辐射冕状极光分割仿真的结果，图9是本发明与现有方法对只含光线的辐射冕状极光分割仿真的结果。

表1.专家组评价结果

表2.无经验组评价结果

两组人员对分割仿真的结果按评价标准U进行评分，专家组评价结果如表1，无经验组评价结果如表2。

从表1和表2可看出，对于多弧状极光仿真分割的效果，两组人员对本发明的评分略比Fuzzy logic和Isolable-contour map评分高，这是因为多弧状极光中所出现的极光形态主要以明亮的光斑和光带为主，三种方法均能很好的将其分割出来，但是有些多弧状极光中还会出现少量光线区域，Fuzzy logic和Isolable-contour map无法将这部分分离出来，而本发明可将这部分成功的分离出来，分割仿真的结果如图4所示，其中：图4(a)是待分割的含有光斑和光线区域的弧状极光，图4(b)是Fuzzy logic方法分割仿真的结果，图4(c)是Isolable-contour map方法分割仿真的结果，图4(d)是本发明分割仿真的结果。

对于热点冕状极光、辐射冕状极光和帷幔冕状极光，无论是专家组还是无经验组，对三种方法的分割仿真结果的评价中，本发明的评分都明显高于Fuzzy logic和Isolable-contour map两种方法，这是由于在三种冕状极光中光斑区域的边界较为模糊，而且出现了大量的光线区域，本发明使用前景像素与背景像素的灰度差异作为分割依据，对光斑区域的边缘和形状没有特殊要求，对于边缘模糊或没有明显形状的光斑，本发明也能够很好的进行分割，本发明还实现了现有方法没有考虑到的极光光线区域的分割，因此本发明仿真分割的效果较之于Fuzzy logic和Isolable-contour map更加精确。对三种冕状极光中所包含的典型极光形态的分割仿真结果如图5-9所示：

图5是含光斑和光线的热点冕状极光分割仿真的结果，其中：图5(a)是待分割的含有光斑和光线区域的热点冕状极光，图5(b)是Fuzzy logic方法分割仿真的结果，图5(c)是Isolable-contour map方法分割仿真的结果，图5(d)是本发明分割仿真的结果。

图6是含光斑和光线的帷幔冕状极光分割仿真的结果，其中：图6(a)是待分割的含有光斑和光线区域的帷幔冕状极光，图6(b)是Fuzzy logic方法分割仿真的结果，图6(c)是Isolable-contour map方法分割仿真的结果，图6(d)是本发明分割仿真的结果。

图7是只含光线的帷幔冕状极光分割仿真的结果，其中：图7(a)是待分割的只含光线的辐射冕状极光，图7(b)是Fuzzy logic方法分割仿真的结果，图7(c)是Isolable-contourmap方法分割仿真的结果，图7(d)是本发明分割仿真的结果。

图8是含光斑和光线的辐射冕状极光分割仿真的结果，其中：图8(a)是待分割的含有光斑和光线区域的辐射冕状极光，图8(b)是Fuzzy logic方法分割仿真的结果，图8(c)是Isolable-contour map方法分割仿真的结果，图8(d)是本发明分割仿真的结果。

图9是只含光线的辐射冕状极光分割仿真的结果，其中：图9(a)是待分割的只含光线的辐射冕状极光，图9(b)是Fuzzy logic方法分割仿真的结果，图9(c)是Isolable-contourmap方法分割仿真的结果，图9(d)是本发明分割仿真的结果。

实验二：客观数据评价

为了客观的衡量分割仿真的结果，将本发明对四类极光的分割仿真的结果与专家手工分割的结果进行对比。

由于极光本身所具有的透明且边缘不清晰的特征，即使是专家手工分割也不可能做到完全精确，为了更加客观的对分割效果进行评价做出以下约定：将专家分割的结果和分割仿真的结果都划分成2×2大小的块，在每个小块中，只要有一个白点就认为该区域是极光区域，如果小块中一个白点都没有，则认为该区域是天空区域。对分割仿真的结果进行客观评价采用下面公式表示：

$\mathrm{\&gamma;}=\frac{{\mathrm{\&gamma;}}_a+{\mathrm{\&gamma;}}_s}{n}$

其中：γ是分割仿真的结果的错误率，γ_a是专家划分为极光区域但仿真分割的结果划分为天空区域的错误块个数，γ_s是专家划分为天空区域但仿真分割的结果划分为极光区域的错误块个数，n是一幅图像中所有的块个数。最终的实验结果见表3。

表3.错误率(％)

从表3可看出，本发明仿真分割的结果的错误率明显低于Fuzzy logic和Isolable-contour map两种方法。其中：多弧状极光主要由明显的光带和少量光线组成，因此Fuzzy logic和Isolable-contour map能够检测到大部分极光区域，只有少量光线区域检测不完全，因此本发明仿真分割的结果的错误率略低于Fuzzy logic和Isolable-contour map；但对于热点、帷幔和辐射冕状极光，光线在整个极光区域中占了很大部分，因此Fuzzylogic和Isolable-contour map未检测到的极光区域增多，错误率也随之增加，本发明使用与极光光线区域形态相对应的方法，通过纹理结构的差异来进行分割，仿真分割结果的错误率远远低于Fuzzy logic和Isolable-contour map两种方法。

综上，本发明方法对极光在全天空图像上所呈现的不同形态采用相应的方法进行分割，提高了极光分割结果的准确性，降低了极光占空比参数的误差。

Claims (3)

1.一种全天空极光图像占空比参数的提取方法，包括如下步骤：

(1)对数据库中所有全天空极光图像进行预处理：将每一幅原始大小为512×512的全天空极光图像用半径为220的圆形区域进行掩模，去掉周围干扰灯光、山脉噪声，掩模成大小为440×440的掩模图像，所有掩模图像构成掩模图像集S_M；

(2)输入一幅掩模图像I_M，使用改进的Otsu方法对极光光斑区域进行分割，得到光斑区域分割结果R_P：

(2a)只保留掩模图像I_M中圆形区域R_M，计算该圆形区域内像素总数：

N'=N-(c²-π(c/2)²)

其中：N是掩模图像I_M中像素总数，c是掩模图像I_M的边长，即圆形区域R_M的直径，c²-π(c/2)²是圆形区域R_M外四个角中像素总数，这四个角中像素灰度均为0；

(2b)计算圆形区域R_M中灰度为0的像素个数n′₀：

n′₀=n₀-(c²-π(c/2)²)

其中：n₀是掩模图像I_M中灰度为0的像素个数；

(2c)根据圆形区域R_M内像素总数N'和灰度为0的像素个数n′₀，计算灰度为0的像素的发生概率：

$g_0=\frac{n_0^{\&prime;}}{N^{\&prime;}}$

(2d)根据圆形区域R_M内像素总数N'和灰度为k的像素的个数n_k，计算灰度为k的像素的发生概率：

$g_k=\frac{n_k}{N^{\&prime;}},k=\mathrm{1,2},......M$

其中：M是圆形区域R_M内灰度的最大值；

(2e)根据灰度为0的像素发生的概率g₀和灰度为k的像素的发生概率g_k，使用Otsu方法计算前景类和背景类像素的最优分割阈值G^*；

(2f)将掩模图像I_M中每个像素与光斑分割阈值G^*进行比较，如果大于或等于G^*，将该像素划分为极光光斑区域，否则将该像素划分为背景天空区域，得到掩模图像I_M的光斑分割结果R_P；

(3)从掩模图像集S_M中选择100幅不含极光或几乎不含极光的掩模图像构建纹理阈值训练集S_K，使用自适应局部二值模式ALBP特征向量提取算法，提取纹理阈值训练集S_K中每一幅图像的ALBP特征向量，构成ALBP特征向量集A＝{α_i}，其中α_i表示纹理阈值训练集S_K中第i幅掩模图像的ALBP特征向量；

(4)对纹理阈值训练集S_K的ALBP特征向量集A进行训练，得到参考向量R和纹理分割阈值Z；

(5)将掩模图像I_M分割成大小为16×16的图像块，构成图像块集X＝{x_i，j}，提取每一块的ALBP特征向量，构成ALBP特征向量集B＝{β_i，j}，其中x_i，j表示掩模图像I_M中第i行第j列的图像块，β_i，j表示图像块x_i，j的ALBP特征向量；

(6)计算ALBP特征向量β_i，j与参考向量R的相似度S_i，j：

$S_{i,j}=X^2({\mathrm{\&beta;}}_{i,j},R)=\underset{l=1}{\overset{C}{\mathrm{\&Sigma;}}}\frac{{({\mathrm{\&beta;}}_{i,j}\left(l\right)-R\left(l\right))}^2}{{\mathrm{\&beta;}}_{i,j}\left(l\right)+R\left(l\right)}$

其中X²(β_i，j，R)表示使用Chi-square统计量计算ALBP特征向量β_i，j与参考向量R之间的相似度，C表示ALBP特征向量的维数共有C维，l表示第l维，β_i，j(l)表示ALBP特征向量的第l维向量，R(l)表示参考向量R的第l维向量；

(7)将相似度S_i，j与纹理分割阈值Z进行比较，当相似度S_i，j大于或等于纹理分割阈值Z，则图像块x_i，j为极光光线区域，否则图像块x_i，j为背景天空区域，得到掩模图像I_M的光线区域分割结果R_T；

(8)对光斑区域分割结果R_P和光线区域分割结果R_T求并集，得到掩模图像I_M的完整分割结果R_A，求光斑区域分割结果R_P和光线区域分割结果R_T的交集，得到光斑区域分割结果R_P和光线区域分割结果R_T的重合部分R_X，将该重合部分R_X去掉，得到修正后的光线区域分割结果R′_T=R_T-R_X；

(9)根据极光区域的完整分割结果R_A、修正后的光线区域分割结果R′_T和光斑区域的分割结果R_P，计算掩模图像I_M的极光占空比参数：O_A，O_T，O_P；

$O_A=\underset{(i,j)\&Element;R_M}{\mathrm{\&Sigma;}}{y}_A(i,j)/N^{\&prime;},$ $y_A(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\&Element;R_A\\ 0,& p(i,j)\&NotElement;R_A\end{array}\right.$

其中：O_A是完整极光区域的占空比参数，表示掩模图像I_M中以各种形态出现的极光在整个天空背景中所占比例，(i，j)是掩模图像I_M中第i行第j列的坐标，(i，j)∈R_M表示只选择在圆形区域R_M内的坐标(i，j)，p(i，j)是坐标(i，j)所对应的像素，y_A(i，j)表示像素p(i，j)是否属于极光区域的完整分割结果R_A，如果像素p(i，j)属于完整分割结果R_A，则y_A(i，j)的值为1，否则y_A(i，j)的值为0，N′是圆形区域R_M中像素总数；

$O_T=\underset{(i,j)\&Element;R_M}{\mathrm{\&Sigma;}}{y}_T(i,j)/N^{\&prime;},$ $y_T(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\&Element;{R^{\&prime;}}_T\\ 0,& p(i,j)\&NotElement;{R^{\&prime;}}_T\end{array}\right.$

其中：O_T表示极光光线区域的占空比参数，是掩模图像I_M中以光线形态所呈现的极光在整个天空背景中所占比例，y_T(i，j)用于判断像素p(i，j)是否属于光线区域分割结果R′_T，如果像素p(i，j)属于光线区域分割结果R′_T，则y_T(i，j)的值为1，否则y_T(i，j)的值为0；

$O_P=\underset{(i,j)\&Element;R_M}{\mathrm{\&Sigma;}}{y}_P(i,j)/N^{\&prime;},$ $y_P(i,j)=\left\{\begin{array}{cc}1,& p(i,j)\&Element;R_P\\ 0,& p(i,j)\&NotElement;R_P\end{array}\right.$

其中：O_P是极光光斑区域的占空比参数，表示掩模图像I_M中以光斑或光带形态出现的极光在整个天空背景中所占比例，y_P(i，j)表示像素p(i，j)是否属于光斑区域分割结果R_P，如果像素p(i，j)属于光斑区域分割结果R_P，则y_P(i，j)的值为1，否则y_P(i，j)的值为0。

2.根据权利要求1所述的全天空极光图像占空比参数的提取方法，其中步骤(3)所述的使用自适应局部二值模式ALBP特征向量提取算法，提取纹理阈值训练集S_K中每一幅图像的ALBP特征向量，按如下步骤进行：

(3a)构建主要概率模式集P：

首先，从掩模图像集S_M中选择100幅含有大量光线区域的掩模图像构成特征训练集S_T，对特征训练集S_T中每一幅掩模图像提取旋转不变LBP特征；

其次，将特征训练集S_T中所有掩模图像的旋转不变LBP特征向量取均值，得到平均旋转不变LBP特征向量V；

最后，按照模式发生概率从大到小的顺序，对平均旋转不变LBP特征向量V进行排序得到排序后的平均旋转不变LBP特征向量V′，选择该向量V'中前C个模式，这C个模式的概率值之和满足大于或等于90％，用这C个模式组成主要概率模式集P；

(3b)对于输入的掩模图像I，提取I的旋转不变LBP特征向量r_M，根据主要概率模式集P，保留r_M中与其对应的C个模式及其概率值，即为掩模图像I的ALBP特征向量。

3.根据权利要求1所述的全天空极光图像占空比参数的提取方法，其中步骤(4)所述的对纹理阈值训练集S_K的ALBP特征向量集A进行训练，按如下步骤进行：

(4a)提取纹理阈值训练集S_K中每幅掩模图像的ALBP特征α_i，并对该训练集S_K中掩模图像的ALBP特征求平均值，得到参考向量R：

$R=\underset{i=1}{\overset{N_R}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i/N_R$

其中：N_R表示纹理阈值训练集S_K中掩模图像的数量，N_R为100；

(4b)将训练集中每一幅掩模图像分割为16×16大小的块，对每个块提取ALBP特征向量，计算每个块的ALBP特征向量与参考向量R之间相似度：

$L_{i,j}=X^2(W_{i,j},R)=\underset{l=1}{\overset{C}{\mathrm{\&Sigma;}}}\frac{{(W_{i,j}\left(l\right)-R\left(l\right))}^2}{W_{i,j}\left(l\right)+R\left(l\right)}$

其中：W_i，j表示第i幅图像第j个块的ALBP特征向量，l是ALBP特征向量中的第l维，W_i，j(l)表示ALBP特征向量W_i，j的第l维向量，R(l)表示参考向量R的第l维向量，C是ALBP特征向量的维数；

(4c)从所有块的ALBP特征向量与参考向量R之间的相似度L_i，j中找到最大值，即为纹理分割阈值：

$Z=\underset{i,j}{\max }\left(L_{i,j}\right).$

Images