CN110533100A - 一种基于机器学习进行cme检测和跟踪的方法 - Google Patents

Abstract

本发明公开了一种基于机器学习进行CME检测和跟踪的方法，以白光日冕仪图像为基础，通过预处理、分类、检测和跟踪四个阶段实现CME(日冕物质抛射)的检测和跟踪。在预处理阶段，对图像进行下采样和差分图像计算。在分类阶段，使用有监督图像分类方法将给定的白光日冕仪图像分类为有、无CME类别，再从已训练好的CNN中提取的中层特征用于检测CME的区域。在检测阶段，使用无监督图像协同定位的方法挖掘出CME区域。在跟踪阶段，将图像转换为极坐标系，对不同位置的CME进行划分，并定义了一些规则清洗过滤CME，计算所有检测到的CME的基本参数。

Description

一种基于机器学习进行CME检测和跟踪的方法

技术领域

本发明属于计算机图像处理和计算机视觉领域，尤其涉及一种基于机器学习进行CME检测和跟踪的方法。

背景技术

对日冕物质抛射(coronal mass ejection，CME)进行观测的太空任务可以追溯到20世纪70年代。SOHO卫星上的日冕仪对CME观测做出了巨大贡献，如大角度和光谱日冕仪(Large Angle and Spectrometric Coronagraph Experiment，LASCO)可以追踪1.1至约30Rs的CME。自从日地关系卫星(Solar TErrestrial RElations Observatory，STEREO)发射升空后，可以从日地关联日冕和太阳风层探测器的仪器包中两个不同的望远镜COR1、COR2观测CME。随着日冕图像数据的大量积累，能够自动检测、跟踪不同特征并构建相应的事件目录(尤其是CME)的能力变得越来越重要。一方面，这可以为CME关键参数的统计研究提供了更便捷的数据访问。另一方面，通过自动检测具有CME标记的日冕仪图像，可以进行实时天气预测。

目前已经开发了不同的CME目录，这些目录中包含长期运行的日冕仪观测资料，它们被分类为手动或自动目录。一般主要使用的数据是由LASCO日冕仪观测并创建的CME目录，这些数据在CDAW数据中心进行维护，CDAW提供了LASCO和其他相关仪器的观测的事件视频以及每个CME的关键参数。虽然CDAW目录已被广泛采用，但CME检测和跟踪是通过人工完成的，这显然是主观和耗时的。因为不同操作员的经验不同，可能会得到不同的检测结果和物理参数。当太阳接近其最大活动时，CME的检测和跟踪更是会需要消耗大量的人力，以上手动CME目录的缺点促使自动目录的开发。

发明内容

发明目的：本发明所要解决的技术问题是针对现有技术的不足，提供一种新的CME检测和追踪的方法。

技术方案：本发明公开了一种基于机器学习进行CME检测和跟踪的方法，该方法用于从连续的白光日冕仪图像中检测出CME并进行跟踪。包括以下步骤：

步骤1，采集原始日冕图像，对图像进行预处理，得到日冕差分图像；

步骤2，图像分类：训练LeNet-5卷积神经网络模型进行日冕差分图像二分类，判定单张日冕差分图像是否有CME发生，得到有CME类别的差分图像和无CME类别的差分图像，有CME发生的图像记为有CME类别的差分图像；

步骤3，CME区域检测：将判断为有CME类别的差分图像进行LeNet-5卷积神经网络特征图提取，使用无监督的图像协同定位方法DDT进行CME区域的粗定位，再利用图像切割方法进行分段区域平滑，得到带有精细CME区域注释的日冕差分图像；

步骤4，进行CME跟踪并计算CME事件的基本参数，。

步骤1包括如下步骤：

步骤1-1，从等级为0.5的LASCO C2数据中提取原始日冕图像，将分辨率为1024×1024的原始日冕图像下采样到512×512分辨率并根据太阳中心的坐标对齐，再将图像通过噪声滤波器，所述噪声滤波器为带有尺寸为3×3的滑动窗口的标准化盒式滤波器；

步骤1-2，使用公式u_i＝n_i-n_i-1计算得到日冕差分图像u_i，其中n_i是当前图像，n_i-1是前一帧图像。

步骤2包括如下步骤：

步骤2-1，将日冕差分图像分辨率下采样为112×112；

步骤2-2，从线上CME人工检测数据库CDAW目录中获得日冕差分图像标签，用于LeNet-5卷积神经网络的训练；LeNet-5卷积神经网络根据日冕差分图像标签学习CME图像和非CME图像的特征，得到判别CME的模型[参考文献：Wei X S，Zhang C L，Wu J，eta1.Unsupervised object discovery and co-localization by deep descriptortransformation[J].Pattern Recognition，2019，88：113-126.]；

步骤2-3，LeNet卷积神经网络架构包括特征提取器模块和分类器模块，特征提取器模块包括两个卷积层、两个非线性激活层和两个下采样层，分类模块包括两个全连接层，在卷积层和下采样层之后，每个图像的特征映射被下采样到25×25分辨率，再通过两个全连接层获得高级语义知识，并输出最终的CME发生概率。

LeNet卷积神经网络第l卷积层(此处l取值为2)的第j个特征图通过如下公式获得：

f(x)＝max(0，x)，

N表示第l-1卷积层的特征图的数量，x是特征图，表示第l-1卷积层的第i个特征图，表示第l卷积层的第i个特征图对应的卷积核，表示第l卷积层的第j个特征图对应的偏置项，f表示非线性激活函数。

步骤2-3包括如下步骤：

步骤2-3-1，将最初的LeNet卷积神经网络架构的输出层改为2个输出单元，对输出单元进行归一化得到CME发生的概率P_CME，具体架构如下式所示：

其中x_CME和x_non-CME是输出层的输出单元，也称输出得分。x_CME得分越高，CME发生的概率越高，x_non-CME得分越高，CME未发生的概率越高。概率值P_CME大于0.5的图像被视为检测到有CME类别的差分图像；

步骤2-3-2，在LeNet卷积神经网络训练阶段，使用随机梯度下降训练模型，目标损失函数L定义如下：

其中M表示训练数据的数量；是从线上CME人工检测数据库CDAW目录中获得的图像标签值，等于0或1；y_i是LeNet卷积神经网络输出概率，大于0且小于1；LeNet卷积神经网络训练的目标是使得关于权重和偏差的函数L最小化，权重和偏见的更新规则如下：

其中i是迭代索引，η是学习率，k_i，b_i是第i次迭代更新后的权重和偏差，LeNet卷积神经网络每层中的权重初始化为均值为零的高斯分布，标准偏差为0.01，并且神经元偏差在每个卷积层和完全连接层中初始化为零；

步骤2-3-3，在LeNet卷积神经网络测试阶段，对一组连续的日冕差分图像进行分类，判定单张图像是否有CME发生，将一组检测到发生CME的连续帧差分图像视为CME演化的图像序列，用于CME协同定位和跟踪。

步骤3中，使用无监督图像协同定位方法DDT进行图像粗定位[参考文献：Wei X S，Zhang C L，Wu J，et al.Unsupervised object discovery and co-localization bydeep descriptor transformation[J].Pattern Recognition，2019，88：113-126.]，包括如下步骤：

步骤3-1，将一组检测到发生CME的连续帧差分图像进行LeNet-5卷积神经网络特征图提取，即从LeNet-5卷积神经网络的第二个卷积层提取特征图，每张图像的特征图为25×25×50维向量，所有图像的特征图构成N×25×25×50的特征向量；

步骤3-2，提取N×25×25×50的特征向量进行主成分分析降维，获得降维后的特征向量为N×25×25×1维，降维后，对图像每个位置进行特征值计算；

步骤3-3，图像所有位置的特征值形成指示矩阵F，将特征值大于设定阈值的像素位置设为CME区域(阈值设为30)，从而得到25×25分辨率的CME粗略定位图像；

步骤3-4，将25×25分辨率的CME粗略定位图像进行上采样到512×512分辨率，结合步骤1得到的差分日冕图像，使用图像切割方法Graph-cut进行CME区域精细分割[参考文献：Boykov Y，Veksler O，Zabih R.Fast approximate energy minimization via graphcuts[C]//Proceedings of the Seventh IEEE International Conference on ComputerVision.IEEE，1999，1：377-384.]。

步骤3-2包括：提取的特征图具有25×25个单元，并且每个单元包含一个d维特征向量，主成分分析使用正交变换将d维相关变量转换为一组线性不相关变量，这些不相关变量通过协方差矩阵的特征分解得到主成分，所述协方差矩阵Cov(x)通过以下公式计算：

其中参数K＝h×w×S，S表示具有h×w分辨率的输入特征图的数量，h表示行像素数，w表示列像素数，表示第s个图像(即具有h×w分辨率的输入特征图)在像素位置(i，j)处的CNN特征，表示S个图像的CNN特征平均值；在特征分解之后，得到协方差矩阵的50维特征向量ξ₍₁₎，...，ξ_(d)，其与降序排序的特征值λ₁≥...λ_d≥0对应，其中ξ_(d)，λ_d分别表示协方差矩阵的第d个特征向量和第d个特征值，选取最大特征值对应的特征向量作为主投影方向计算第n个图像在像素位置(i，j)的特征值

其中n取值为s～S；通过上述方式，将特征向量尺寸从50减小到1，并且变换后的特征值被视为每个像素位置处的共同对象的出现可能性。

步骤3-3中，每个像素位置的特征值组成指示矩阵F，特征值大于30的像素位置判定为有CME存在，从而得到CME区域的粗略位置，所述指示矩阵F为：

其中，f_(h，w)表示像素位置(h，w)处的特征值，指示矩阵F的尺寸为h×w。

步骤3-4包括如下步骤：

步骤3-4-1，对步骤3-3得到的每个像素位置的特征值进行归一化，得到每个像素位置发生CME的概率；

步骤3-4-2，计算步骤1得到的差分日冕图像每个像素位置与其周围临近像素的亮度差值；

步骤3-4-3，构造如下能量函数E(l)：

E(l)＝λ_sE_smooth(l)+λ_dE_ddata(l)，

其中l类表示图像是CME类别，λ_s和λ_d是非负常数，(λ_s，λ_d一般取值分别为0.7和1.8)，E_smooth(l)表示通过相邻邻域位置的亮度差异来测量l的类一致性，E_data(l)测量了l和预测数据之间的不一致，设置E_smooth(l)和E_data(l)如下：

E_data(l)＝∑_p-log(pr(l_p))，

其中p、q表示像素位置，l_p表示位置p被指定为l类，l_q表示位置q被指定为l类，pr(l_p)表示像素位置p被指定为l类的概率，N₈是特征图中的8个周围像素，I_p、I_q分别表示位置p处的亮度、q处的亮度；

使用图像切割方法求解能量函数E(l)对应图结构的最小割集，进而获得能量函数的最小解，从而获得带有精细CME区域注释的日冕差分图像[参考文献：BoykovY，VekslerO，Zabih R.Fast approximate energy minimization via graph cuts[C]//Proceedingsof the Seventh IEEE International Conference on Computer Vision.IEEE，1999，1：377-384.]。

步骤4包括如下步骤：

步骤4-1，将步骤3得到的带有精细CME区域注释的所有日冕差分图像转换为具有360×360分辨率的极坐标系，太阳北极与图像北极一致；

步骤4-2，计算给定时间范围内日冕差分图像中每个角度的CME区域的最大高度，对出现CME区域的角度按照位置连通性进行合并，合并后组成若干不同宽度大小的CME事件，一个CME事件代表固定角度区间的CME区域随时间演化的结果，检查每个CME事件是否满足如下两个条件：

条件1，存在两帧以上；

条件2，CME的最大高度能够超过日冕仪观测的视场边界；

如果任一条件不满足，则丢弃包含所述CME事件的日冕差分图像；否则保留；

步骤4-3，计算CME事件的基本参数，包括中心角位置，角宽度和速度，其中，中心角位置定义为末端位置角度和起始位置角度的中间位置，角宽度等于末端位置角度减去起始位置角度，通过高度时间线性拟合测量得到每个位置角度的CME传播速度，所有CME检测到的位置角度的速度分布的中值为CME的代表速度。

目前，机器学习技术在许多不同的研究领域得到了越来越广泛的应用。它将计算机科学与和日光层物理学跨学科结合在一起。针对不同的太阳特征和太空天气，机器学习技术已经有了相当多的应用。例如，Huang等(2018)将深度学习方法应用于耀斑预测。Camporeale等(2017)和Delouille等(2018)分别使用机器学习技术对太阳风和日冕洞进行分类。最近，Galvez等人(2019)甚至为太阳动力学观测站(SDO)任务编制了一个适合机器学习研究的数据集。

在计算机视觉领域，机器学习在图像分类，特征检测和跟踪方面表现出色。鉴于其巨大的成功以及对CME快速检测和跟踪的需求，本发明基于LASCO C2数据采用机器学习技术CAMEL进行自动CME检测和跟踪。

有益效果：为了评估本发明的机器学习技术的性能，将本发明的结果与其他一些现有的自动检测和跟踪技术，即CACTus、CORIMP和SEEDS的结果进行了比较。本发明的方法不仅能够检测和跟踪明亮的信号，而且还能够跟踪微弱的信号，输出更完整的CME结构。同时具有较好的性能，可以尽早捕捉到LASCO C2中CMEs的出现时间，得到更精确的CME形态和时间信息，可以用来推导更精确的CME运动。本发明的方法记录了CMEs的像素级位置。因此，本发明不仅可以检测和跟踪CME前缘，还可以在检测到的CME区域中检测和跟踪任何其他感兴趣的结构。本发明的方法在计算上既便宜又快速，经过适当的训练，在一台普通的PC机上，对一张LASCO C2图像中检测和跟踪CME只需要几秒钟。基于机器学习进行CME检测和跟踪是图像分类的副产品，对CME标记图像进行分类，将CME数据及时有效地分布到不同的空间天气预报中心，对预测CME的到达具有重要意义。

附图说明

下面结合附图和具体实施方式对本发明做更进一步的具体说明，本发明的上述或其他方面的优点将会变得更加清楚。

图1是LeNet架构示意图，包含2个卷积层，2个下采样层和2个全连接层。

图2是图像共定位流程，从训练好的LeNet模型中提取每个CME检测图像的特征，并使用PCA作为投影方向来转换这些特征以评估它们的相关性。

图3是CME区域改进示例，使用图形切割方法获得具有良好注释的CME精确区域的图像。

图4a是跟踪模块的输入：CME演变的图像序列，包括一组CME类的灰色图像和一组用红色标注的像素级CME区域的图像。选定的时间范围是2012年2月4日。

图4b上方是19时35分转换为极坐标系的CME图像，下方是用于计算CME角位置的清洗后的CME图像。

图4c是最大速度角位置处的CME高度-时间图，其中虚线表示所得到的被跟踪CME的时间范围。

图5是本发明流程图。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

如图5所示，对于一组给定的LASCO C2图像，本发明经过如下步骤，完成对CME的检测和追踪：

步骤(1)，图像预处理阶段。对1024×1024分辨率的输入图像LASCO C2进行下采样得到分辨率为512x512的图像并降噪，得到差分图像(后一帧减前一帧)。

步骤(2)，图像分类阶段。从CDAW目录中获得具有或不具有CME的图像标签，用于训练模型。测试阶段，判断是否有CME事件发生。使用LeNet提取图像的中级特征，用于步骤(3)的图像检测。

步骤(3)，图像检测阶段。使用无监督的图像共定位方法——DDT进行CME区域粗定位，再利用图像切割方法(Graph-Cut)方法进行分段区域平滑，得到有精细调整的CME区域的图像。

步骤(4)，图像跟踪阶段。定义规则来识别一系列CME，首先一个CME事件必须存在超过2帧。其次，CME的最大高度必须超出视场(FOV)。计算所有检测到的CME事件的中心角位置、角宽度、角速度。

具体来说，本发明针对用户输入的一组图像，首先对图像进行下采样/降噪等预处理获得高质量的图片，再经过分类、检测、追踪三个算法模块获得精确的CME区域和CME事件的具体参数。分类模块确认该图像是否有CME发生，检测模块检测CME在整个图像上的具体区域，追踪模块则是追踪一个CME事件的发生以及发生的具体参数。

实施例

下面根据实施例说明本发明的各个步骤。

步骤(1)，图像预处理阶段。对1024×1024分辨率的输入图像LASCO C2进行下采样和降噪，并计算得到差分图像。

步骤(11)，使用太阳物理软件包(SSW)中的lasco_readfits.pro读取下载0.5级LASCO C2文件，然后使用SSW中的reduce_level_1.pro处理1级数据，该处理包括暗电流，杂散光，失真，渐晕，测光，时间和位置校正的校准。经过处理后，太阳的北极方向和图像的北极方向一致。所有具有1024×1024分辨率的输入LASCO C2图像首先被下采样到512×512分辨率并根据太阳中心的坐标对齐。然后，所有下采样图像都通过噪声滤波器，以抑制一些尖锐的噪声特征。此处使用带有尺寸为3×3的滑动窗口的标准化盒式滤波器。标准化的盒式滤波是一种基本的线性图像滤波，它计算周围像素的平均值。

步骤(12)，使用如下公式计算得到日冕差分图像ui，

u_i＝n_i-n_i-1

其中n_i是当前图像，n_i-1是前一帧图像。根据LASCO图像的曝光时间比，采用亮度归一化方法。对于部分含有缺失块的LASCO图像，根据之前的图像计算缺失块掩膜：如果之前图像中像素的值为零，那么差分图像中相同像素的值也为零。最后的差分图像乘以缺失的块掩膜。

步骤(2)，获取图像标签，用于训练模型，如果具有CME则图像标签为1，否则为0。测试阶段，测试单张图片是否有CME事件发生。使用LeNet提取图像的中级特征，用于图像检测。具体步骤为：

步骤(21)，将日冕差分图像分辨率下采样为112x112。

步骤(22)，从CDAW数据目录中获取图像标签，用于LeNet-5网络的训练。

步骤(23)，LeNet网络架构分为特征提取器模块和分类器模块。特征提取器模块由卷积层，非线性激活层和下采样层组成，分类模块由两个全连接层构成。在卷积层和下采样层之后，每个图像的特征映射被下采样到25×25分辨率，最后两个全连接的层构成分类模块。卷积层可以看作是局部连接的网络，其中每个隐藏单元将仅连接到输入图像的小的连续区域，并在每个位置获得不同的特征激活值。卷积核在工作时，会有规律地扫过输入特征，在感受野内对输入特征做矩阵元素乘法求和并叠加偏差量。第1层的第j个特征图如下获得：

f(x)＝max(0，x)

N表示第l-1层的特征图的数量，k表示卷积核，b表示偏置项，f表示非线性激活函数，这里f使用的是Rectified Linear Units(ReLUs)激活函数，它可以使CNN训练速度提高几倍。下采样层有助于扩大各个位置的感受野并且聚合特征，此层分别处理每个特征图，它计算每个特征图中邻域的最大值。

步骤(231)，最初的LeNet架构设计用于手写数字识别，输出层输出10个单位，代表每个类别(0-9)的概率。本发明中修改了输出层，更改后输出层输出2单元，x_CME和x_non-CME是最终输出层的输出单元，也称输出得分，对二者进行归一化处理得到CME的发生概率。LeNet架构如图1所示。

x_CME得分越高，CME发生的概率越高，x_non-CME得分越高，CME未发生的概率越高。输出概率值大于0.5的的图像被视为检测到CME的图像。图1显示了本发明使用的LeNet架构。

步骤(232)，在训练阶段，使用随机梯度下降训练模型，作为分类问题，目标损失函数可以定义如下：

其中M表示训练数据的数量，是从人工检测数据库CDAW目录中获得的图像标签值，等于0或1，y_i是LeNet卷积神经网络输出概率，大于0且小于1。可以看出L是非负的，所以LeNet网络训练的目标是使得关于权重和偏差的函数L最小化。使用随机梯度下降训练我们的模型，批量大小为128个示例。权重和偏见的更新规则如下：

其中i是迭代索引，η是学习率，学习率初始化为0.0001，并在终止前减少三次。只有一批训练样例用于更新每次迭代中的权重和偏差。每层中的权重初始化为均值为零的高斯分布，标准偏差为0.01，并且神经元偏差在每个卷积层和完全连接层中初始化为零。

步骤(233)，在测试阶段，连续运行差分图像按时间顺序分类。可以将一组CME检测到的连续帧视为CME演化的图像序列，其用于CME共定位和跟踪。

步骤(3)，图像检测阶段。使用无监督的图像共定位方法——DDT挖掘用于分割CME区域的高级信息，再利用图像切割方法(Graph-Cut)方法进行分段区域平滑，得到有精细调整的CME区域的图像。

步骤(31)，从CNN模型的最后一个卷积层提取卷积特征图。使用深度描述符变换(DDT)挖掘用于分割CME区域的隐藏信息。DDT其利用主成分分析(PCA，Pearson 1901)来分析CNN特征。

步骤(32)，提取的特征图具有25×25个单元，并且每个单元包含一个d维特征向量。PCA使用正交变换将d维相关变量转换为一组线性不相关变量，这些变量通过协方差矩阵的特征分解得到主成分。输入数据的协方差矩阵Cov(x)通过以下公式计算：

其中K＝h×w×S。S表示具有h×w分辨率的输入特征图的数量，并且表示图像s在像素位置(i，j)处的CNN特征。对协方差矩阵进行特征分解得到50维特征向量ξ₍₁₎，...，ξ_(d)，其与降序排序的特征值λ₁≥...λ_d≥0对应。选取最大特征值对应的特征向量作为主投影方向。对于图像n的CNN特征向量的特定位置(i，j)进行特征计算，其主要特征计算如下：

通过上述方式，将其特征尺寸从50减小到1，并且变换后的特征值可以被视为每个像素处的共同对象的出现可能性。

步骤(33)，所有位置的特征值形成指示矩阵F，其尺寸为h×w：

共同定位的流程可以在图2中找到。从训练的CNN模型中获得的CME演化的图像序列包括一组CME图像，其通过DDT算法直接处理以用于CME区域共定位。图像共定位的最终输出是一组CME区域掩模图像，使用最近插值法，使其分辨率与输入特征图的分辨率相同。

步骤(34)，在计算机视觉中使用图形切割方法来进行分段区域平滑。显然，指示矩阵只能粗略计算像素位置属于CME或非CME类的概率，而且相邻像素之间也可能出现类一致性的问题。

步骤(341)，对步骤(33)得到的每个像素位置的特征值进行归一化，得到每个像素位置发生CME的概率P。

步骤(342)，计算步骤(1)得到的差分图像每个像素位置与其周围临近像素的亮度差值。

步骤(343)，为了解决相邻像素之间可能出现类一致性的问题，制定了能量最小化的框架。在这个框架中，最小化能量的图像像素标记为l：

E(l)＝λ_sE_smooth(l)+λ_dE_data(l)

其中λ_s，λ_d是非负常数，以平衡每项的影响。E_smooth(l)通过相邻邻域位置的亮度差异来测量l的类一致性。E_data(l)测量了l和预测数据之间的不一致，这些数据主要基于步骤(2)计算的概率进行优化。设置E_smooth(l)和E_data(l)如下：

其中p、q表示像素位置，l_p、l_q分别表示位置p被指定为l类，位置q被指定为l类，pr(l_p)表示像素位置p被指定为l类的概率，N₈是特征图中的8个周围像素，I_p、I_q表示位置p、q处的亮度。然后可以采用图切割优化来有效地解决能量最小化问题。图切割算法生成标记问题的相关图，求解该图的最小割集，然后获得能量函数的最小解。图3展示了优化之前和之后的比较结果的一个例子。

步骤(4)，图像追踪阶段。定义规则来识别一系列CME，计算所有检测到的CME事件的基本参数。

步骤(41)，将带有精细CME区域注释的图像转换为具有360×360分辨率的极坐标系。

步骤(42)，定义规则来识别一系列CME，首先一个CME事件必须存在超过2帧。其次，CME的最大高度必须超出视场(FOV)。

步骤(43)，计算CME事件的基本参数：中心角位置，角宽度，速度。中心角位置定义在末端PA(Position Angle位置角度)和起始PA的中间位置。每个CME的高度可以通过中心角位置的峰值高度来定义CME的角宽度等于末端PA减去起始PA。确定每个CME的开始时间和结束时间，分析CME的中心位置角度的高度随时间变化，所有时间段CME的高度都在上升。接下来，检查每个时间段中的CME是否满足定义的两个条件。不满足条件的分段(不存在至少两帧或该时间段内的CME未到达FOV)被丢弃。为了得出CME速度，计算所有检测到CME位置的速度分布中值，通过线性拟合高度——时间来计算每个位置的传播速度。图4a～图4c显示了位置角度与最大速度的CME高度——时间图(由于说明书附图为灰度图，图中红色部分显示为颜色较浅的部分)，图4a是跟踪模块的输入：在2012年2月4日18时23分31秒和2012年2月4日19时35分31秒时的CME演变的图像序列，包括一组CME类的灰色图像和一组像素级CME区域的图像。图4b上方是2012年2月4日19时35分31秒转换为极坐标系的CME图像，下方是清洗后用于计算CME角位置的CME图像。图4c是最大速度角位置处的CME高度-时间图，其中横轴表示时间，纵轴表示CME事件的高度与太阳高度的比值(R/Rsun，一般取值2.2-6.2)，图中虚线表示所得到的被跟踪CME事件的时间范围。

本发明提供了一种基于机器学习进行CME检测和跟踪的方法，具体实现该技术方案的方法和途径很多，以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。本实施例中未明确的各组成部分均可用现有技术加以实现。

Claims (10)

1.一种基于机器学习进行CME检测和跟踪的方法，其特征在于，包括以下步骤：

步骤1，采集原始日冕图像，对图像进行预处理，得到日冕差分图像；

步骤2，图像分类：训练LeNet-5卷积神经网络模型进行日冕差分图像二分类，判定单张日冕差分图像是否有CME发生，有CME发生的图像记为有CME类别的差分图像；

步骤3，CME区域检测：使用LeNet-5卷积神经网络对有CME类别的差分图像进行特征图提取，使用无监督的图像协同定位方法进行CME区域的粗定位，再利用图像切割方法进行分段区域平滑，得到带有精细CME区域注释的日冕差分图像；

步骤4，进行CME跟踪并计算CME事件的基本参数。

2.根据权利要求1所述的方法，其特征在于，步骤1包括如下步骤：

步骤1-1，从等级为0.5的LASCO C2数据中提取原始日冕图像，将分辨率为1024×1024的原始日冕图像下采样到512×512分辨率并根据太阳中心的坐标对齐，再将图像通过噪声滤波器，所述噪声滤波器为带有尺寸为3×3的滑动窗口的标准化盒式滤波器；

步骤1-2，使用公式u_i＝n_i-n_i-1计算得到日冕差分图像u_i，其中n_i是当前图像，n_i-1是前一帧图像。

3.根据权利要求2所述的方法，其特征在于，步骤2包括如下步骤：

步骤2-1，将日冕差分图像分辨率下采样为112×112；

步骤2-2，从线上CME人工检测数据库CDAW目录中获得日冕差分图像标签，用于LeNet-5卷积神经网络的训练；LeNet-5卷积神经网络根据日冕差分图像标签学习CME图像和非CME图像的特征，得到判别CME的模型；

步骤2-3，LeNet卷积神经网络架构包括特征提取器模块和分类器模块，特征提取器模块包括两个卷积层、两个非线性激活层和两个下采样层，分类模块包括两个全连接层，在卷积层和下采样层之后，每个图像的特征映射被下采样到25×25分辨率，再通过两个全连接层获得高级语义知识，并输出最终的CME发生概率。

4.根据权利要求3所述的方法，其特征在于，LeNet卷积神经网络第l卷积层的第j个特征图通过如下公式获得：

f(x)＝max(0，x)，

N表示第l-1卷积层的特征图的数量，x是特征图，表示第l-1卷积层的第i个特征图，表示第l卷积层的第i个特征图对应的卷积核，表示第l卷积层的第j个特征图对应的偏置项，f表示非线性激活函数。

5.根据权利要求4所述的方法，其特征在于，步骤2-3包括如下步骤：

步骤2-3-1，将最初的LeNet卷积神经网络架构的输出层改为2个输出单元，对输出单元进行归一化得到CME的发生概率P_CME，具体架构如下式所示：

其中x_CME和x_non-CME是输出层的2个输出单元，也称输出得分，x_CME得分越高，CME发生的概率越高；x_non-CME得分越高，CME未发生的概率越高，概率值P_CME大于0.5的图像被视为检测到有CME类别的差分图像；

步骤2-3-2，在LeNet卷积神经网络训练阶段，使用随机梯度下降训练模型，目标损失函数L定义如下：

其中M表示训练数据的数量；是从人工检测数据库CDAW目录中获得的图像标签值，等于0或1；y_i是LeNet卷积神经网络输出概率，大于0且小于1；LeNet卷积神经网络训练的目标是使得关于权重和偏差的函数L最小化，权重和偏差的更新规则如下：

其中i是迭代索引，η是学习率，k_i，b_i分别是第i次迭代更新后的权重和偏差，LeNet卷积神经网络每层中的权重初始化为均值为零的高斯分布，标准偏差为0.01，并且神经元偏差在每个卷积层和完全连接层中初始化为零；

步骤2-3-3，在LeNet卷积神经网络测试阶段，对一组连续的日冕差分图像进行分类，判定单张图像是否有CME发生，将一组检测到发生CME的连续帧差分图像视为CME演化的图像序列，用于CME协同定位和跟踪。

6.根据权利要求5方法，其特征在于，步骤3中，使用无监督图像协同定位方法进行图像粗定位，包括如下步骤：

步骤3-1，将一组检测到发生CME的连续帧差分图像进行LeNet-5卷积神经网络特征图提取，即从LeNet-5卷积神经网络的第二个卷积层提取特征图，每张图像的特征图为25×25×50维向量，所有图像的特征图构成N×25×25×50的特征向量；

步骤3-2，提取N×25×25×50的特征向量进行主成分分析降维，获得降维后的特征向量为N×25×25×1维，降维后，对图像每个位置进行特征值计算；

步骤3-3，图像所有位置的特征值形成指示矩阵F，将特征值大于设定阈值的像素位置设为CME区域，从而得到25×25分辨率的CME粗略定位图像；

步骤3-4，将25×25分辨率的CME粗略定位图像进行上采样到512×512分辨率，结合步骤1得到的差分日冕图像，使用图像切割方法Graph-cut进行CME区域精细分割。

7.根据权利要求6所述的方法，其特征在于，步骤3-2包括：提取的特征图具有25×25个单元，并且每个单元包含一个d维特征向量，主成分分析使用正交变换将d维相关变量转换为一组线性不相关变量，这些不相关变量通过协方差矩阵的特征分解得到主成分，所述协方差矩阵Cov(x)通过以下公式计算：

其中参数K＝h×w×S，S表示具有h×w分辨率的输入特征的数量，h表示行像素数，w表示列像素数，表示第s个图像在像素位置(i，j)处的CNN特征，表示S个图像的CNN特征平均值；在特征分解之后，得到协方差矩阵的特征向量ξ₍₁₎，...，ξ_(d)，其与降序排序的特征值λ₁≥...λ_d≥0对应，其中ξ_(d)，λ_d分别表示协方差矩阵的第d个特征向量和第d个特征值，选取最大特征值对应的特征向量作为主投影方向，计算第n个图像在像素位置(i，j)的特征值

通过上述方式，将特征向量尺寸从50减小到1，并且变换后的特征值被视为每个像素位置处的共同对象的出现可能性。

8.根据权利要求7所述的方法，其特征在于，步骤3-3中，每个像素位置的特征值组成指示矩阵F，特征值大于30的像素位置判定为有CME存在，从而得到CME区域的粗略位置，所述指示矩阵F为：

其中，f_(h，w)表示像素位置(h，w)处的特征值，指示矩阵F的尺寸为h×w。

9.根据权利要求8所述的方法，其特征在于，步骤3-4包括如下步骤：

步骤3-4-1，对步骤3-3得到的每个像素位置的特征值进行归一化，得到每个像素位置发生CME的概率；

步骤3-4-2，计算步骤1得到的差分日冕图像每个像素位置与其周围临近像素的亮度差值；

步骤3-4-3，构造如下能量函数E(l)：

E(l)＝λ_sE_smooth(l)+λ_dE_data(l)，

其中l类表示图像是CME类别，λ_s和λ_d是非负常数，，E_smooth(l)表示通过相邻邻域位置的亮度差异来测量l的类一致性，E_data(l)测量了l和预测数据之间的不一致性，设置E_smooth(l)和E_data(l)如下：

E_data(l)＝∑_p-log(pr(l_p))，

其中p、q表示像素位置，l_p表示位置p被指定为l类，l_q表示位置q被指定为l类，pr(l_p)表示像素位置p被指定为l类的概率，N₈是特征图中的8个周围像素，I_p、I_q分别表示位置p处的亮度、q处的亮度；

使用图像切割方法求解能量函数E(l)对应图结构的最小割集，进而获得能量函数的最小解，从而获得带有精细CME区域注释的日冕差分图像。

10.根据权利要求9所述的方法，其特征在于，步骤4包括如下步骤：

步骤4-1，将步骤3得到的带有精细CME区域注释的所有日冕差分图像转换为具有360×360分辨率的极坐标系，太阳北极与图像北极一致；

步骤4-2，计算给定时间范围内日冕差分图像中每个角度的CME区域的最大高度，对出现CME区域的角度按照位置连通性进行合并，合并后组成不同宽度大小的CME事件，检查每个CME事件是否满足如下两个条件：

条件1，存在两帧以上；

条件2，CME的最大高度能够超过日冕仪观测的视场边界；

如果任一条件不满足，则丢弃包含所述CME事件的日冕差分图像；否则保留；

步骤4-3，计算CME事件的基本参数，包括中心角位置，角宽度和速度，其中，中心角位置定义为末端位置角度和起始位置角度的中间位置，角宽度等于末端位置角度减去起始位置角度，通过高度时间线性拟合测量得到每个位置角度的CME传播速度，所有CME检测到的位置角度的速度分布的中值为CME的代表速度。