CN104463226A - 基于序列环扇形分块的弧状极光序列检测方法 - Google Patents

Abstract

本发明公开了一种基于序列环扇形分块局部方向模式SRLDP的弧状极光序列检测方法，主要解决现有技术忽略极光图像时间域信息检测不准确的问题。其实现步骤为：(1)利用序列局部方向模式编码方法提取极光序列局部方向模式特征，得到极光序列局部方向模式图；(2)将序列局部方向模式图分为环扇形块，统计每个环扇形块的像素点画出其直方图；(3)将一个序列的环扇形块的直方图转化为一个局部向量并串接，对所有序列重复操作，将串接后的局部向量作为极光序列局部方向模式特征矩阵的行向量；(4)利用自动调节谱聚类方法对极光序列局部方向模式特征矩阵进行聚类，得到分类结果。本发明提高了对弧状极光序列分类的准确率，可用于事件检测。

Description

基于序列环扇形分块的弧状极光序列检测方法

技术领域

本发明属于图像处理技术领域，涉及弧状极光序列检测方法，可用于极光序列特征提取、极光序列分类。

背景技术

极光是来自磁层的高能粒子沉降到高层大气并与中性成分碰撞激发的一种大气发光现象，是太阳风-磁层-电离层耦合造成的地球高纬活动的重要表现形式之一。所以人们通过极光形态及其演化的系统观测可以获得磁层和日地空间电磁活动的大量信息，有助于深入研究太阳活动对地球的影响方式与程度，对了解空间天气过程的变化规律具有重要意义。弧状极光作为日侧极光的一个重要形式，它的分布和发生揭示了地磁活动和磁层中的脉动现象。单幅图像的特征表示方法只能用来表征静态的极光图像，现在我们需要寻找一种能有效表示极光序列的方法。极光序列前后两帧是有关系的，应该利用这种时间域的关系来提取极光序列特征。

弧状极光作为日侧极光的一个重要形式，对它的形态特征和发生机理相对而言已经比较清楚。Wang Q等人在文章“Spatial texture based automatic classification of daysideaurora in all-sky images.Journal of Atmospheric and Solar-Terrestrial Physics,2010,72(5):498-508.”中通过研究静态图像得出了四种极光类型的统计分布规律，并初步验证了Hu等人在文章“Synoptic distribution of dayside aurora:Multiple-wavelength all-skyobservation at Yellow River Station inSvalbard.Journal of Atmospheric andSolar-Terrestrial Physics,2009,71(8):794-804.”中根据活动区域将日侧极光分为四类的合理性，还从形态学方向给出了这四种分类的解释。Han等人在文章“Aurora imagesclassification via features salient coding.Journal of xidian university,2013,40(6),180-186.”中提出的基于特征显著性编码的极光图像分类方法将极光图像分为四类。Yang Xi等人在文章“Wavelet hierarchical model for aurora images classification.Journal of xidian university,2013,40(2):18-24.”中提出的一种分层小波模型下的极光图像分类算法，利用支持向量机分类器把极光图像分成弧状和三种冕状类型。Han等人在文章“Aurora imageclassification based on LDA combining with saliency information.Ruan Jian XueBao/Journal of Software,2013,24(11):2758-2766.“中采用基于融合显著信息的潜在主题极光图像分类方法进行极光图像分类。然而，现有对于极光的研究多是对极光图像进行分类，研究单幅图像的特点，而没有把在时间上连续且类型相同的图像作为一个序列进行研究，从而影响了极光的准确分类。

发明内容

本发明的目的在于针对上述已有技术的不足，提出一种基于序列环扇形分块的弧状极光序列检测方法，提高对极光序列分类的准确率。

本发明的技术方案是：利用极光图像前后帧之间的时间域的关系来提取极光序列特征，对弧状极光序列进行检测，其实现步骤包括如下：

(1)提取极光图像序列局部方向模式图像I：

1a)输入极光序列图像后，提取单幅极光图像的局部方向模式特征e；

1b)根据极光序列图像中每幅图像的局部方向模式特征e，利用序列局部方向模式编码方法获取极光序列局部方向模式特征E；

1c)对极光序列局部方向模式特征E进行阀值设定和进制转换处理，获得该序列局部方向模式图像I；

(2)将序列局部方向模式图像I转换成极坐标图像I′，再将极坐标图像I′划分成X个环扇形块状区域，并对每个环扇形块状区域做直方图统计后，画出统计后的直方图SH_i,j，(i,j)对应不同方位环扇形块的位置；

(3)通过坐标轴转换将每个环扇形块状区域的直方图SH_i,j转换为局部行向量s′_1×Z，串接局部向量s′_1×Z获得全局向量s_1×v，将每一个全局向量s_1×v作为扇形分块局部方向模式特征矩阵S_hv的每一行，从而得到扇形分块局部方向模式特征矩阵S_hv，其中，h为输入极光序列总数，v＝Z×X；

(4)利用自动调节谱聚类方法对图像序列的环扇形分块局部方向模式特征矩阵S_hv进行聚类，将h个极光图像序列分为弧状极光序列h₁和非弧状极光序列h₂两类。

本发明具有如下优点：

1.本发明由于采用了序列局部方向模式编码方法对极光序列图像进行处理，充分利用了极光图像前后帧之间时间域的关系，从而提高了分类准确率；

2.本发明由于采用了环扇形分块的方式对极光序列局部方向模式图像进行处理，充分考虑了极光图像像素点的位置信息，从而更好的利用了极光图像的局部信息，使分类准确率更高。

附图说明

图1是本发明的实现流程图。

具体实施方式

参照图1，本发明的具体实现步骤如下：

步骤一：提取极光图像序列局部方向模式图像I。

1a)输入极光图像序列，并提取极光序列图像中每一幅图像的局部方向模式特征e：

1a1)对每一幅图像的每个像素点定义八个方向：0°、45°、90°、135°、180°、225°、270°、315°；

1a2)利用Kirsch算子与图像中的每个像素点进行卷积，得到该像素点的八个方向的边缘响应值：m_l＝I(x,y)*K_l(x,y)，l＝0,1,…,7，其中，m_l为该像素点的第l个边缘响应值，I(x,y)为像素点，K_l(x,y)为Kirsch算子的第l个边缘检测算子；

1a3)对每一幅图像的每个像素点进行步骤1a1)和1a2)的操作，得到每一幅图像的局部方向模式特征e；

1b)根据极光序列图像中每幅图像的局部方向模式特征e，在含有N幅极光图像的序列中，比较所有图像同一方位像素点的对应方向的边缘响应值m_l，取最大边缘响应值M_l＝max_N m_l，对序列图像中的每个像素点进行此操作，得到极光序列局部方向模式特征E；

1c)对极光序列局部方向模式特征E进行阀值设定和进制转换处理，获得该序列局部方向模式图像I：

1c1)设定8个最大边缘响应值M₀,M₁,…,M₇分别对应于8个方向0°、45°、90°、135°、180°、225°、270°、315°；

1c2)将8个方向的边缘响应值M_l按照从大到小排序，将M_l中前q个较大的值设置为M′_l＝1，其它8-q个值设置为M_l'＝0，用M_l'作为M_l的替代值，且M_l'的值仅为1或0；

1c3)将M'₀到M'₇从高位到低位依次排序，得到一串由1和0组成的数字串，即极光序列局部方向模式E的二进制码；

1c4)将极光序列局部方向模式E的二进制码转换为十进制码：其中，τ为极光序列局部方向模式特征E对应目标像素点的十进制值；

1c5)按照步骤1c1)～1c4)的方法计算序列图像的每个像素点的序列局部方向模式的十进制值，得到极光序列局部方向模式图像I。

步骤二：对序列局部方向模式图像I进行坐标转换，化分环扇形块状区域，画出直方图SH_i,j。

2a)取极光序列局部方向模式图像I的中心点为直角坐标系的原点O，以平行于图像I的两个边记为横轴x和纵轴y；

2b)将直角坐标系中的原点O作为极坐标中的极点O'，对直角坐标系中的序列局部方向模式图像I进行极坐标系转换：x＝ρ cos θ,y＝ρ sin θ，获得极坐标图像I′，其中(x,y)为直角坐标，(θ,ρ)为对应的极坐标，图像I中的极光图像有效圆形区域转换为图像I′后仍是以圆形出现，且其极点O'为图像I′内圆形区域的中心点；

2c)对极坐标图像I′进行扇形分区处理：以极坐标图像I′内有极光的圆形区域为基准得到最大的圆O″和其半径R，从圆O″的中心点O'引出一条射线作为极轴，按逆时针方向，以45°为单位将圆O″分为8个等分扇形区域i，i＝1,2,…,8，其中，点O'为圆O″的中心点；

2d)在进行扇形分区后图像I′的基础上进行环形划分处理：对圆O″内的一条半径R进行4等分，获得3个等分点w_p，p＝1,2,3，以w_p到O'的长度为半径作以O'为圆心的同心圆，获得到4个圆环区域j，j＝1,2,3,4，得到X＝32个环扇形块状区域，其中，(i,j)对应不同的环扇形块状区域；

2e)对第(i,j)个环扇形块状区域做直方图统计：

2e1)确定直方图横轴所代表的像素点的十进制值τ，τ的所有可能值为含有k个1，8-k个0的二进制码进行全排列后求得的十进制值，τ的总数为：个；

2e2)对直方图纵轴所代表的像素点的十进制值τ出现的次数SH_i,j(τ)做统计，其公式如下：

${\mathrm{SH}}_{i,j}\left(\mathrm{\τ}\right)=\underset{i=1}{\overset{8}{\mathrm{\Σ}}}\underset{j=1}{\overset{4}{\mathrm{\Σ}}}f(\mathrm{SLDP}(\mathrm{\θ},\mathrm{\ρ}),\mathrm{\τ})$    0°≤θ≤360°，0≤ρ≤R，

其中，SLDP(θ，ρ)为(θ，ρ)对应像素点的十进制值；函数f(SLDP(θ，ρ),τ)所代表的映射关系为：当f中的两个参数值相等时，f(SLDP(θ，ρ),τ)＝1，当其中的两个参数值不相等时，f(SLDP(θ，ρ),τ)＝0；

2f)将横轴用Z个τ的不同的值按从小到大的顺序在横轴上标记，纵轴为每个不同的τ出现的次数SH_i,j(τ)，得到环扇形块状区域的直方图SH_i,j。

步骤三：获取环扇形分块局部方向模式特征矩阵S_hv。

3a)将一个图像序列中的每个环扇形块状区域的直方图SH_i,j转换为局部向量s′_1×Z，局部向量s′_1×Z中含有Z个值，这些值分别对应于直方图SH_i,j上的标记点τ的值从小到大出现的次数SH_i,j(τ)；

3b)将X个局部向量s′_1×Z进行串接获得一个全局向量s_1×v，其中，v＝Z×X；

3c)按照步骤3a)～3b)的方法对所有图像序列的环扇形块状区域的直方图SH_i,j进行转换，获得h个全局向量s_1×v，其中，h为输入极光序列总数；

3d)将每一个全局向量s_1×v作为环扇形分块局部方向模式特征矩阵S_hv的每一行，从而得到环扇形分块局部方向模式特征矩阵S_hv。

步骤四：利用自动调节谱聚类方法对环扇形分块局部方向模式特征矩阵S_hv进行聚类处理。

4a)将环扇形分块局部方向模式特征矩阵S_hv的行向量由第一行到第h行分别记为s₁,s₂,…s_h，用s_ξ和s_η来表示s₁,s₂,…s_h中两个任意的两个行向量，其中，ξ,η∈{1,…,h}；

4b)对环扇形分块局部方向模式特征矩阵S_hv的行向量s_ξ和s_η进行计算值a_ξη，当ξ≠η时，当ξ＝η时，a_ξη＝0，把a_ξη作为相似度矩阵A中第ξ行、第η列的元素，得到一个相似度矩阵A∈R^h×h，其中，d(s_ξ,s_η)为欧氏距离函数，σ_ξ和σ_η为环扇形分块局部方向模式特征矩阵S_hv中的局部尺度；

4c)给定相似度矩阵B∈R^h×h，b_ξη为对角矩阵B中第ξ行、第η列的元素，利用相似度矩阵A中的各元素a_ξη计算对角矩阵B中元素b_ξη，当ξ＝η时，当ξ≠η时，b_ξη＝0；

4d)根据相似度矩阵A和对角矩阵B，得到归一化矩阵：L＝B^-1/2 AB^-1/2，其中，L∈R^h×h；

4e)计算归一化矩阵L的前c个较大的特征值所对应的特征向量矩阵U∈R^h×c，u_ωγ为特征向量矩阵U中第ω行、第γ列的元素，其中，ω∈{1,…,h}，γ∈{1,…,c}，c为自主设定的聚类分组数；

4f)给定二次对角矩阵y_ωγ为二次对角矩阵Y中第ω行、第γ列的元素，利用相似度矩阵U中的各元素u_ωγ计算二次对角矩阵Y中元素y_ωγ， $y_{\mathrm{\ω\γ}}=\frac{u_{\mathrm{\ω\γ}}}{{\left({\mathrm{\Σ}}_{\mathrm{\γ}}{u_{\mathrm{\ω\γ}}}^2\right)}^{\frac{1}{2}}};$

4g)利用梯度下降法对二次对角矩阵Y做分类处理，将极光序列分为了弧状极光序列和非弧状极光序列两类。

对环扇形分块局部方向模式特征矩阵S_hv进行聚类处理不局限于自动调节谱聚类方法，还可以用网络聚类、分层聚类等方法进行处理。

本发明的效果通过下面的实验仿真进一步说明：

1.仿真条件与方法：

硬件平台为：Intel Core i3、2.93GHz、3.45GB RAM；

软件平台为：Windows7操作系统下的MATLAB R2011b；

实验的数据：是由我国北极黄河站的全天空极光数据中已经建立好的弧状极光数据库和非弧状极光数据库，这些数据来自2003年到2007年，其中弧状极光序列有222个共11150幅图像，非弧状极光序列有338个共12561幅图像，序列中的极光图像均为512*512的灰度图像。

2.仿真内容与结果

用本发明和现有基于局部方向模式的弧状极光序列检测方法对上述实验数据同时进行分类处理，其分出的弧状极光序列的查全率和查准率结果如表1。

表1 弧状极光序列查全率和查准率

从表1可见，本发明在弧状极光序列分类上有很好的效果，分类准确率上有显著提高。

Claims (7)

1.一种基于序列环扇形分块的弧状极光序列检测方法，包括如下步骤：

(1)提取极光图像序列局部方向模式图像I：

1a)输入极光序列图像后，提取单幅极光图像的局部方向模式特征e；

1b)根据极光序列图像中每幅图像的局部方向模式特征e，利用序列局部方向模式编码方法获取极光序列局部方向模式特征E；

1c)对极光序列局部方向模式特征E进行阀值设定和进制转换处理，获得该序列局部方向模式图像I；

(2)将序列局部方向模式图像I转换成极坐标图像I′，再将极坐标图像I′划分成X个环扇形块状区域，并对每个环扇形块状区域做直方图统计后，画出统计后的直方图SH_i,j，(i,j)对应不同方位环扇形块的位置；

(3)通过坐标轴转换将每个环扇形块状区域的直方图SH_i,j转换为局部行向量s′_1×Z，串接局部向量s′_1×Z获得全局向量s_1×v，将每一个全局向量s_1×v作为扇形分块局部方向模式特征矩阵S_hv的每一行，从而得到扇形分块局部方向模式特征矩阵S_hv，其中，h为输入极光序列总数，v＝Z×X；

(4)利用自动调节谱聚类方法对图像序列的环扇形分块局部方向模式特征矩阵S_hv进行聚类，将h个极光图像序列分为弧状极光序列h₁和非弧状极光序列h₂两类。

2.根据权利要求1所述的方法，其中步骤1a)提取单幅图像的局部方向模式特征e，按如下步骤进行：

1a1)对单幅图像每个像素点定义八个方向，以图像中3×3的九个像素点的块为准，中间位置为目标像素点，其它为该目标像素点对应的八个不同方向；

1a2)利用Kirsch算子分别与该像素点的八个方向对应卷积，得到该像素点的八个方向的边缘响应值：m_l＝I(x,y)*K_l(x,y)，l＝0,1,…,7，其中，m_l为图像对应像素点的8个方向的卷积值，K_l(x,y)为八个方向的Kirsch边缘检测算子；

1a3)对单幅图像每个像素点进行1a1)和1a2)的操作，得到单幅图像的局部方向模式特征e。

3.根据权利要求1所述的方法，其中步骤1b)所述的利用序列局部方向模式编码方法获取极光序列局部方向模式特征E，是在含有N幅极光图像的序列中，比较所有图像同一方位像素点的对应方向的卷积值m_l，取最大值M_l＝max_N m_l，对序列图像中的每个像素点进行此操作，得到极光序列局部方向模式特征E。

4.根据权利要求1所述的方法，其中步骤1c)的对极光序列局部方向模式特征E进行阀值设定和进制转换处理，获得该序列局部方向模式图像I，按如下步骤进行：

1c1)设定目标像素点8个卷积值M_l在3×3的块中的分布位置为：卷积值M₀为3×3的块中右下角所处位置的值，按顺时针方向依次将8个方向记为M₀,M₁,…,M₇；

1c2)将目标像素点的8个方向的卷积值M_l按照从大到小排序，将M_l中前q个较大的值设置为M′_l＝1，其它8-q个值设置为M′_l＝0，用M′_l作为M_l的替代值，即M′_i在3×3的块中的位置对应于M_l，且M′_l的值仅为1或0；

1c3)将M'₀到M'₇从高位到低位依次排序，得到一串由1和0组成的数字串，即极光序列局部方向模式E的二进制码；

1c4)将极光序列局部方向模式E的二进制码转换为十进制码：其中，τ为极光序列局部方向模式特征E对应目标像素点的十进制值；

1c5)按照步骤1c1)、1c2)、1c3)和1c4)的方法计算序列图像的每个像素点的序列局部方向模式的十进制值，得到极光序列局部方向模式图像I。

5.根据权利要求1所述的方法，其中所述步骤(2)中将极坐标图像I′划分成X个环扇形块状区域，按如下步骤进行：

2a)对极坐标图像I′进行扇形分区处理：以极坐标图像I′内有极光的圆形区域为基准得到最大的圆O″和其半径R，从圆O″的中心点O'引出一条射线作为极轴，按逆时针方向，以45°为单位将圆O″分为8个等分扇形区域i，i＝1,2,…,8，其中，点O'为圆O″的中心点；

2b)在进行扇形分区后图像I′的基础上进行环形划分处理：对圆O″内的一条半径R进行4等分，获得3个等分点w_p，p＝1,2,3，以w_p到O'的长度为半径作以O'为圆心的同心圆，获得4个圆环区域j，j＝1,2,3,4，得到X＝32个环扇形块状区域，其中，(i,j)对应不同的环扇形块状区域。

6.根据权利要求1所述的方法，其中所述步骤(2)中对每个环扇形块状区域做直方图统计后，画出统计后的直方图SH_i,j，按如下步骤进行：

2c)对第(i,j)个环扇形块状区域做直方图统计：

首先，确定直方图横轴所代表的像素点的十进制值τ，τ的所有可能值为含有k个1，8-k个0的二进制码进行全排列后求得的十进制值，τ的总数为：个；

然后，对直方图纵轴所代表的像素点的十进制值τ出现的次数SH_i,j(τ)做统计，其公式如下：

其中，(θ，ρ)为各个像素点的极坐标，SLDP(θ，ρ)为(θ，ρ)对应像素点的十进制值；函数f(SLDP(θ，ρ),τ)所代表的映射关系为：当f中的两个参数值相等时，f(SLDP(θ，ρ),τ)＝1，当其中的两个参数值不相等时，f(SLDP(θ，ρ),τ)＝0；

2d)将横轴用Z个τ的不同的值按从小到大的顺序在横轴上标记，纵轴为每个不同的τ出现的次数SH_i,j(τ)，得到环扇形块状区域的直方图SH_i,j。

7.根据权利要求1所述的方法，其中步骤(4)所述对图像序列的环扇形分块局部方向模式特征矩阵S_hv进行聚类，按如下步骤进行：

4a)将图像序列的环扇形分块局部方向模式特征矩阵S_hv的行向量由第一行到第h行分别记为s₁,s₂,…s_h，用s_ξ和s_η来表示s₁,s₂,…s_h中两个任意的两个行向量，其中，ξ,η∈{1,…,h}；

4b)对特征矩阵S_hv的行向量s_ξ和s_η进行计算，得到一个相似度矩阵 $A=\left(\begin{array}{ccc}{a}_{11}& ...& a_{1h}\\ .& & .\\ .& a_{\mathrm{\ξ\η}}& .\\ .& & .\\ a_{h1}& ...& a_{\mathrm{hh}}\end{array}\right),$ A∈R^h×h，a_ξη为相似度矩阵A中第ξ行、第η列的元素，当ξ≠η时，当ξ＝η时，a_ξη＝0，其中，d(s_ξ,s_η)为欧氏距离函数，σ_ξ和σ_η为特征矩阵S_hv中的局部尺度；

4c)用相似度矩阵A的各元素a_ξη进行计算得到对角矩阵 $B=\left(\begin{array}{ccc}{b}_{11}& ...& b_{1h}\\ .& & .\\ .& b_{\mathrm{\ξ\η}}& .\\ .& & .\\ b_{h1}& ...& b_{\mathrm{hh}}\end{array}\right),$ B∈R^h×h，b_ξη为对角矩阵B中第ξ行、第η列的元素，当ξ＝η时，当ξ≠η时，b_ξη＝0；

4d)根据相似度矩阵A和对角矩阵B，得到归一化矩阵：L＝B^-1/2AB^-1/2，其中，L∈R^h×h；

4e)计算归一化矩阵L的前c个较大的特征值得到一个特征向量矩阵 $U=\left(\begin{array}{ccc}{u}_{11}& ...& u_{c1}\\ .& & .\\ .& u_{\mathrm{\ω\γ}}& .\\ .& & .\\ u_{1h}& ...& u_{\mathrm{ch}}\end{array}\right),$ U∈R^h×c，u_ωγ为特征向量矩阵U中第γ行、第ω列的元素，其中，c为自主设定的聚类分组数ω∈{1,…,c}，γ∈{1,…,h}；

4f)对特征向量矩阵U中各元素u_ωγ进行计算，得到二次对角矩阵： $Y=\left(\begin{array}{ccc}{y}_{11}& ...& y_{c1}\\ .& & .\\ .& y_{\mathrm{\ω\γ}}& .\\ .& & .\\ y_{1h}& ...& y_{\mathrm{ch}}\end{array}\right),$ Y∈R^h×c，y_ωγ为二次对角矩阵Y中第γ行、第ω列的元素，

4g)利用梯度下降法对二次对角矩阵Y做分类处理，将极光序列分为了弧状极光序列和非弧状极光序列两类。