基于快速魯棒性主成分分析的日冕噴流自動檢測方法*

耿成杰，李潤鑫，劉 輝，尚振宏

(1. 昆明理工大學信息工程與自動化學院，云南 昆明 650500；2. 昆明理工大學云南省計算機技術(shù)應用重點實驗室，云南 昆明 650500；3. 中國科學院云南天文臺，云南 昆明 650216；4. 昆明理工大學云南省人工智能重點實驗室，云南 昆明 650500)

日冕噴流是太陽外層大氣(日冕)中普遍存在的太陽瞬變現(xiàn)象[1]，是沿著開放磁力線運動的等離子體流，在不同位置和不同波段觀測到的噴流具有不同的形態(tài)學和動態(tài)學特征[2]。近年來的研究表明，日冕噴流在分析太陽磁場拓撲結(jié)構(gòu)、加熱局地日冕以及加速太陽風等方面可能起著十分重要的作用，同時噴流經(jīng)常伴隨著耀斑、暗條和日冕物質(zhì)拋射(Coronal Mass Ejection, CME)等太陽活動現(xiàn)象[3]。

文[4]通過隨機選擇20個X射線噴流和極紫外波段噴流分析得出，這20個噴流都起源于迷你暗條爆發(fā)，并且暗條的噴發(fā)始于噴流基部亮點附近。在噴流噴發(fā)過程中存在迷你暗條的上升運動，在暗條到達噴流頂端之前，測量得到的暗條平均長度等于X射線噴流的平均寬度，這與噴流噴發(fā)是由微暗條噴發(fā)驅(qū)動的觀點一致。文[5]首次對日冕爆裂噴流和由它引起的同時發(fā)生氣泡狀和噴流狀日冕物質(zhì)拋射的活動現(xiàn)象進行分析，提出了一個解釋日冕爆裂噴流的模型，對日冕爆裂噴流的冷熱分量和對應動態(tài)相關(guān)的日冕物質(zhì)拋射現(xiàn)象進行解釋：噴流冷分量是由噴流底部的暗條噴發(fā)引起，并進一步導致氣泡狀日冕物質(zhì)拋射；噴流狀日冕物質(zhì)拋射與噴流熱分量有關(guān)，熱分量是限制暗條的閉合場與周圍開放場磁力線重新聯(lián)接而產(chǎn)生的向外移動的熱等離子體。文[6]使用太陽動力學天文臺衛(wèi)星的大氣成像儀、日震學與磁場成像儀(Helioseismic and Magnetic Imager, HMI)的圖像數(shù)據(jù)和1 m新真空太陽望遠鏡(New Vacuum Solar Telescope, NVST)的Hα線心圖像數(shù)據(jù)對噴流的演化和噴流源區(qū)兩個相鄰暗條細絲間的變化進行分析發(fā)現(xiàn)，該雙向環(huán)狀噴流的爆發(fā)是由于噴流源區(qū)相鄰暗條細絲間的磁重聯(lián)導致，和目前普遍認為的雙向噴流的觸發(fā)機制有很大不同。

隨著太陽觀測設(shè)備的發(fā)展，高時空分辨率的日冕圖像數(shù)據(jù)越來越多，通過圖像處理技術(shù)從海量數(shù)據(jù)中自動檢測日冕噴流，有助于研究人員對噴流的形成過程、觸發(fā)機制以及噴流對日冕的加熱和太陽風加速背后的物理規(guī)律的研究，同時也能夠為分析噴流與暗條、日冕物質(zhì)拋射等伴隨活動現(xiàn)象之間的內(nèi)在聯(lián)系提供便利。

文[7]將索貝爾(Sobel)邊緣檢測和霍夫(Hough)變換相結(jié)合對太陽動力學天文臺的大氣成像儀觀測數(shù)據(jù)中的日冕噴流進行檢測，并將檢測的日冕噴流相應的活動持續(xù)時間、長寬度和視線速度等參數(shù)存入日球物理學事件知識庫(Heliophysics Event Knowledgebase, HEK)?；舴蜃儞Q將檢測對象從目標空間轉(zhuǎn)換到參數(shù)空間，避免了在目標空間檢測時的分類運算，但是該方法比較依賴原始圖像空間中檢測對象的曲線表達式。而日冕噴流在形態(tài)上復雜多變，邊緣曲線的表達式難以確定，造成檢測效果不佳。文[8-9]使用幀間差分法對日出(Hinode)衛(wèi)星X射線望遠鏡(X-ray Telescope, XRT)觀測數(shù)據(jù)進行檢測，并對檢測的噴流長度、速度等參數(shù)進行計算和統(tǒng)計分析。幀間差分法原理簡單，易于實現(xiàn)，但是對噴流運動速度緩慢和強度較弱時的檢測結(jié)果不理想。文[2]首先對1 m新真空太陽望遠鏡觀測的多個Hα日浪事件進行研究，得到該類事件的觀測特性和主觀判斷結(jié)果，然后對觀測數(shù)據(jù)中比較暗的結(jié)構(gòu)(黑子、暗條和針狀體等)進行標記，再利用閾值分割方法分離觀測數(shù)據(jù)中的暗結(jié)構(gòu)，排除已經(jīng)標記的暗結(jié)構(gòu)，將剩余的暗結(jié)構(gòu)和主觀判斷結(jié)果進行比較，最終得到Hα日浪。該方法需要對黑子、暗條等結(jié)構(gòu)進行標記，比較費時。

事實上，日冕噴流的檢測可以看作在日冕序列圖像中檢測突變，或者是從復雜隨機變化的背景中檢測有規(guī)律變化的噴流活動。因此，如何實現(xiàn)復雜變化背景與稀疏前景的分離是檢測日冕噴流的關(guān)鍵?？焖亵敯粜灾鞒煞址治龇椒ㄔ谶\動目標檢測中有廣泛的應用。本文用該方法對太陽動力學天文臺的大氣成像儀觀測數(shù)據(jù)中不同時間段、不同波段、不同觀測位置的日冕噴流進行檢測。檢測結(jié)果表明，該方法與傳統(tǒng)的幀間差分法相比能夠檢出更加完整的噴流，檢出的噴流區(qū)域孔洞較少，比較連續(xù)。

1 太陽動力學天文臺/大氣成像儀實驗數(shù)據(jù)

太陽動力學天文臺是美國國家航空航天局(National Aeronautics and Space Administration, NASA)于2010年2月發(fā)射的一顆空間天氣衛(wèi)星，搭載的大氣成像儀具有較高的時空分辨率，能夠為進一步研究太陽爆發(fā)活動提供高質(zhì)量的觀測圖像數(shù)據(jù)。本文使用兩組不同時間段、不同波段和不同觀測位置的日冕噴流實驗數(shù)據(jù)集。數(shù)據(jù)集Ⅰ：波段為17.1 nm和19.3 nm，視場范圍449 × 179像元，像元分辨率為0.6″/pixel，視場中心坐標(205.2″，-253.8″)。17.1 nm波段的觀測時間為2014年10月4日10:10:00～10:41:00(UT)，共156幀圖像；19.3 nm波段的觀測時間為2014年10月4日10:10:08～10:41:08(UT)，共156幀圖像，時間分辨率都為12 s。數(shù)據(jù)集Ⅱ：波段為17.1 nm和19.3 nm，視場范圍186 × 313像元，像元分辨率為0.6″/pixel，視場中心坐標(-219″，-366.9″)。17.1 nm波段的觀測時間為2014年3月3日19:10:00～20:05:00(UT)，共276幀圖像；19.3 nm波段的觀測時間為2014年3月3日19:09:56～20:04:56(UT)，共276幀圖像，時間分辨率都為12 s。

這兩組數(shù)據(jù)中噴流變化尺度較大，像素強度增加比較明顯，同時還存在一些亮度增強的非日冕噴流變化，導致視場范圍內(nèi)的圖像強度發(fā)生復雜變化。圖1(a)顯示了數(shù)據(jù)集Ⅰ在10:25:00(UT)17.1 nm波段對數(shù)對比度增強后的全日面圖像；圖1(b)顯示了(a)中矩形區(qū)域視場范圍內(nèi)17.1 nm波段在該時刻對應的對數(shù)對比度增強后的日冕圖像；圖1(c)顯示了19.3 nm波段在10:24:56(UT)對應的對數(shù)對比度增強后的日冕圖像。(b)和(c)中的日冕噴流在像素強度上分布不均勻，且(b)中的圖像內(nèi)容更復雜，日冕噴流與圖像背景像素強度差異大，易于區(qū)別。19.3 nm波段在10:27:20(UT)的日冕噴流強度較弱，與圖像背景強度差異小，此時17.1 nm波段中的日冕噴流像素強度與圖像背景依然有較大的區(qū)別。圖2顯示了數(shù)據(jù)集Ⅱ在2014年3月3日19:10:00～20:05:00(UT)時間段17.1 nm/19.3 nm波段的對數(shù)對比度增強后的日冕圖像，圖2(b)中的圖像內(nèi)容比圖2(c)中的復雜，但(c)中的日冕噴流在直觀上易于辨別。

圖1 (a)數(shù)據(jù)集Ⅰ在10:25:00(UT)時刻17.1 nm波段對數(shù)對比度增強后的全日面日冕圖像；(b)為(a)中紅色矩形框標記的日冕圖像；(c)為(a)中紅色矩形框在10:24:56(UT)時刻19.3 nm波段標記的日冕圖像

圖2 (a)數(shù)據(jù)集Ⅱ在19:31:49(UT)時刻17.1 nm波段對數(shù)對比度增強后的全日面日冕圖像；(b)為(a)中紅色矩形框標記的日冕圖像；(c)為(a)中紅色矩形框在19:31:56(UT)時刻19.3 nm波段標記的日冕圖像

2 快速魯棒性主成分分析方法

魯棒性主成分分析[10]方法的數(shù)學模型為

(1)

(1)式是一個非凸的優(yōu)化問題，目前還沒有針對該問題有效的解決方案，可以通過松弛技術(shù)，將(1)式中矩陣S的l0范數(shù)使用l1范數(shù)代替，使用核范數(shù)‖·‖*代替rank(·)近似地逼近矩陣L的秩，由此可以將(1)式松弛為一個凸優(yōu)化問題進行求解，數(shù)學表達式為

(2)

其中，D為觀測矩陣；‖L‖*為矩陣L的核范數(shù)；‖S‖1為矩陣S的l1范數(shù)；λ(λ> 0)為權(quán)重參數(shù)，用于調(diào)節(jié)低秩矩陣和稀疏矩陣的占比。假設(shè)觀測矩陣中含有觀測到的m個數(shù)據(jù)，每一個數(shù)據(jù)列向量化為矩陣D∈Rm×n的一列，矩陣D可以分解為低秩矩陣L和稀疏矩陣S，即D=L+S。其中L和S未知，L滿足低秩特性，S滿足稀疏特性，數(shù)值可以為任意大小，這也是主成分分析魯棒性的體現(xiàn)。

由于魯棒性主成分分析是通過代價高昂的優(yōu)化來計算的，在數(shù)據(jù)量較大時，實時性較低，因此魯棒性主成分分析快速算法的開發(fā)顯得尤為重要。文[11]對(2)式的懲罰約束進行改進，得到在視頻背景建模/前景檢測下比非精確增廣拉格朗日乘子法(Inexcat Augmented Lagrangian Multiplier, IALM)快一個數(shù)量級的快速魯棒性主成分分析方法，數(shù)學模型為

(3)

由于性能的原因，文中主要考慮rank(L)=t時的情況，

(4)

對(4)式使用交替最小化方式進行求解，具體求解過程為

(5)

(6)

(5)式通過對D-Sk進行部分奇異值分解來計算矩陣L，t為奇異值的個數(shù)，當t很小時，算法有非常高的計算效率。(6)式通過對D-Lk+1進行逐像素的軟閾值操作[shrink(D-Lk+1,λ)]計算矩陣S。軟閾值函數(shù)定義為

shrink(x,ε)=sign(x)max{0,|x|-ε}.

(7)

快速魯棒性主成分分析方法是根據(jù)低秩和稀疏分解的思想進行運動物體的檢測，基本思路是假設(shè)由序列圖像組成的觀測矩陣D可以分解為低秩背景矩陣L和稀疏前景矩陣S。由于太陽日冕時刻發(fā)生復雜隨機的變化，所以日冕圖像中的內(nèi)容復雜多變。從數(shù)據(jù)集Ⅰ和Ⅱ中的序列圖像來看，日冕噴流在時間序列上變化尺度大，而序列圖像中日冕變化尺度比較小的部分在時間序列上差異不大，且在序列圖像中有比較高的比例，因此變化尺度比較小的這部分日冕可以近似認為是低秩、相關(guān)的。日冕噴流序列圖像數(shù)據(jù)的這一特點與快速魯棒性主成分分析的低秩、稀疏分解思想相似，因此，本文嘗試將快速魯棒性主成分分析方法用于日冕噴流檢測，并在實驗對比分析中驗證方法的可行性和優(yōu)勢。

3 日冕噴流檢測與結(jié)果分析

在使用快速魯棒性主成分分析方法檢測日冕噴流運動之前，我們先對Level 1.0的原始圖像數(shù)據(jù)進行預處理。預處理包括FITS(Flexible Image Transport System)文件標準化、太陽自轉(zhuǎn)校正、曝光時間校正和圖像對比度增強。對于FITS文件的標準化和太陽自轉(zhuǎn)校正，太陽軟件包中提供了相應的解決方法(aia_prep.pro和drot_map.pro)。FITS文件標準化主要實現(xiàn)文件日面中心對齊，像元分辨率統(tǒng)一為0.6″并消除側(cè)傾角。太陽自轉(zhuǎn)校正主要消除由太陽自轉(zhuǎn)在投影平面造成的位移，減小對日冕噴流演化分析時的影響。FITS頭文件包含的圖像數(shù)據(jù)曝光時間大小不一定相同，導致圖像序列的亮度差異大，把圖像數(shù)據(jù)除以曝光時間得到單位時間內(nèi)曝光量的圖像數(shù)據(jù)，然后再使用標準差歸一化方法，使圖像序列的強度統(tǒng)一。原始FITS文件的圖像數(shù)據(jù)對比度較低，這里使用對數(shù)函數(shù)對圖像進行對比度增強，方便進行圖像內(nèi)容觀察和噴流檢測。

對圖像進行預處理之后，本文使用快速魯棒性主成分分析方法對序列圖像進行噴流檢測。我們使用快速魯棒性主成分分析方法分別對數(shù)據(jù)集Ⅰ和數(shù)據(jù)集Ⅱ進行實驗，對檢測結(jié)果進行塊匹配和3D濾波(Block-Matching and 3D filtering, BM3D)[12]去噪后，使用文[13]提出的最大類間方差法分割圖像。并將二值化后連通域面積小于20的區(qū)域剔除，得到最終檢測結(jié)果如圖3、圖4。

本文以兩組不同時間段、不同波段和不同觀測位置的日冕圖像序列作為檢測對象，用快速魯棒性主成分分析的檢測結(jié)果與幀間差分法的結(jié)果進行對比分析，說明快速魯棒性主成分分析方法的可行性和優(yōu)勢。從圖3、圖4中第2行(使用MATLAB的Image labeler應用程序中的像素標記工具，將有噴流的地方直接標記，標記的地方值為1，其余地方為0，將標記結(jié)果作為標簽使用)和第3行(快速魯棒性主成分分析檢測結(jié)果)圖像序列對比，快速魯棒性主成分分析能夠檢出噴流的絕大部分，且與幀間差分法(第4行，當前時間幀與前12 s相減)檢出的噴流區(qū)域相比，空洞較少，且檢出的區(qū)域比較連續(xù)，說明魯棒性主成分分析方法能夠用來檢測日冕噴流。

對于數(shù)據(jù)集Ⅰ，從圖3(a)中矩形框標記的日冕噴流的檢測結(jié)果來看，在10:27:12時刻，快速魯棒性主成分分析方法的檢測結(jié)果明顯優(yōu)于差分法的，且檢測結(jié)果中離散點較少，檢出的噴流區(qū)域更加完整；圖3(b)中矩形框標記的日冕噴流的檢測結(jié)果顯示，快速魯棒性主成分分析方法能夠檢出噴流的微小區(qū)域，差分法則沒能檢出。

圖3 數(shù)據(jù)集Ⅰ。(a)17.1 nm波段檢測對比結(jié)果；(b)19.3 nm波段檢測對比結(jié)果Fig.3 Detection and comparison results of data Ⅰ. (a) 17.1 nm images; (b) 19.3 nm images

對于數(shù)據(jù)集Ⅱ，從圖4中矩形框標記的日冕噴流的檢測結(jié)果來看，與圖3矩形框標記的日冕噴流的檢測結(jié)果基本一致。圖4(a)中矩形框標記的日冕噴流的檢測結(jié)果顯示，快速魯棒性主成分分析方法能夠檢出更多的噴流，檢測結(jié)果中不存在大量的非噴流區(qū)域。圖4(b)中矩形框標記的日冕噴流的檢測結(jié)果顯示，快速魯棒性主成分分析方法與差分法相比能夠檢出日冕噴流，且日冕噴流檢測結(jié)果從直觀上容易辨認，不存在誤檢的情況。

從圖3、圖4日冕噴流的檢測結(jié)果來看，快速魯棒性主成分分析方法能夠在復雜多變的日冕背景中檢出運動緩慢且強度較弱的日冕噴流。從原理上講，幀間差分法對圖像內(nèi)容單一、背景穩(wěn)定的圖像數(shù)據(jù)時具有較好的適用性，但是對包含復雜隨機變化的日冕圖像卻存在漏檢(圖4中矩形框4)的情況。

但是快速魯棒性主成分分析方法和幀間差分法的檢測結(jié)果有時相差不大(如圖3中(a)和(b)，圖4中(a)和(b))，且差分法在檢出噴流區(qū)域的同時，也檢出許多不是噴流的區(qū)域。為了對快速魯棒性主成分分析方法和幀間差分法的檢測結(jié)果進行客觀分析，同時也為快速魯棒性主成分分析方法中的t和λ的設(shè)置提供客觀依據(jù)，本文采用F1指標[14]對兩者的檢測效果進行綜合評估，公式為

圖4 數(shù)據(jù)集Ⅱ。(a)17.1 nm波段檢測對比結(jié)果；(b)19.3 nm波段檢測對比結(jié)果Fig.4 Detection and comparison results of data Ⅱ. (a) 17.1 nm mages; (b) 19.3 nm images

(8)

(9)

其中，P為準確率(Precision)；R為召回率(Re-call)；F1為綜合性能指標；TP(True Positives)，F(xiàn)P(False Positives)和FN(False Negatives)分別表示正確檢測為日冕噴流的前景像素點、錯誤檢測為日冕噴流的背景像素點和錯誤檢測為背景的前景像素點。綜合性能指標F1越高，表明方法的檢測效果越好，如表1。對于數(shù)據(jù)集Ⅰ和數(shù)據(jù)集Ⅱ，快速魯棒性主成分分析方法在綜合性能指標上均有提高，且日冕噴流檢測的準確率明顯提高。經(jīng)過多次實驗，根據(jù)F1指標高低設(shè)置快速魯棒性主成分分析模型中的t和λ值，數(shù)據(jù)集Ⅰ中17.1 nm波段對應的t為30，λ為0.01，數(shù)據(jù)集Ⅱ中17.1 nm波段對應的t為7，λ為0.01。19.3 nm波段采用和17.1 nm波段一樣的值。

表1 快速魯棒性主成分分析與幀間差分檢測結(jié)果綜合比較

此外，在對日冕噴流進行統(tǒng)計研究時，日冕噴流的速度、長度和寬度等參數(shù)經(jīng)常作為統(tǒng)計對象，這里使用快速魯棒性主成分分析方法的日冕噴流檢測結(jié)果二值圖像的最小面積外接矩形估算噴流的長寬比，繼而估計噴流的大致形態(tài)。最后根據(jù)檢測結(jié)果的二值圖像凸殼質(zhì)心的位移疊加，估算噴流的總體速度。

根據(jù)估算結(jié)果得到結(jié)論：日冕噴流的長寬比都大于1，形態(tài)上具有較明顯的朝向性。在數(shù)據(jù)集Ⅰ 17.1 nm波段測得的速度為174 km·s-1，長寬比為3.19，在10:25:00時達到最大長度約8.5 × 104km，活動時長為10:14:00～10:28:00；在19.3 nm波段測得的速度為194 km·s-1，長寬比為3.25，在19:33:44時達到最大長度約9.7 × 104km，活動時長為10:13:31～10:26:03。在數(shù)據(jù)集Ⅱ的17.1 nm波段測得的速度為174 km·s-1，長寬比為4.74，在19:33:12時達到最大長度約9.3 × 104km，活動時長為19:25:24～19:34:48；在19.3 nm波段測得的速度為157 km·s-1，長寬比為3.15，在19:32:20時達到最大長度約8.03 × 104km，活動時長為19:22:44～19:34:44。

4 結(jié) 論

本文將快速魯棒性主成分分析方法應用于日冕噴流檢測，通過低秩和稀疏分解的思想對日冕的動態(tài)背景進行建模，從而檢出作為前景變化的日冕噴流。實驗結(jié)果表明，快速魯棒性主成分分析方法在噴流變化尺度較大、像素強度增加明顯，同時還存在亮度增強的小尺度非噴流變化的日冕序列圖像中檢測日冕噴流時，檢測結(jié)果與幀間差分法相比，檢出的日冕噴流區(qū)域孔洞較少，區(qū)域連續(xù)，且檢測的準確度有所提高。差分法原理簡單，易于實現(xiàn)，但是在針對日冕噴流檢測這種除了日冕噴流的變化之外，還存在其他小的、亮度增加的非噴流變化的動態(tài)場景時效果不理想，而快速魯棒性主成分分析方法能夠檢出幀間差分法不能檢測到的強度較弱的日冕噴流。

目前，快速魯棒性主成分分析方法存在的問題是模型中參數(shù)t和λ的設(shè)置，這些參數(shù)還不能實現(xiàn)自動設(shè)置，需要通過多次試驗人工選擇，而如何通過客觀指標分析對模型參數(shù)進行自動設(shè)置和對可行性進行客觀評價是改進的重點。日冕噴流的觀測特征參數(shù)(如速度、長度和寬度等)還需要更加精準的測量。此外，快速魯棒性主成分分析方法在數(shù)據(jù)量比較大的情況下運算量大，耗費時間長，如何對全日面上的日冕噴流特別是弱的日冕噴流進行檢測是下一步的重點內(nèi)容。