基于相位一致性和方向?yàn)V波冕環(huán)震蕩識(shí)別與跟蹤*1

馮　松，王　洋，趙　燕，袁　泉，楊云飛，鄧　輝，王　鋒

(1. 昆明理工大學(xué)信息工程與自動(dòng)化學(xué)院計(jì)算機(jī)系，云南 昆明　650500; 2. 昆明理工大學(xué)云南省計(jì)算機(jī)技術(shù)應(yīng)用重點(diǎn)實(shí)驗(yàn)室，云南 昆明　650500)

?

基于相位一致性和方向?yàn)V波冕環(huán)震蕩識(shí)別與跟蹤*1

馮松1,2，王洋1,2，趙燕1,2，袁泉1,2，楊云飛1,2，鄧輝1,2，王鋒1,2

(1. 昆明理工大學(xué)信息工程與自動(dòng)化學(xué)院計(jì)算機(jī)系，云南 昆明650500; 2. 昆明理工大學(xué)云南省計(jì)算機(jī)技術(shù)應(yīng)用重點(diǎn)實(shí)驗(yàn)室，云南 昆明650500)

摘要：冕環(huán)結(jié)構(gòu)的準(zhǔn)確識(shí)別和跟蹤有利于更好地研究和分析冕環(huán)特征以及同耀斑爆發(fā)和日冕物質(zhì)拋射之間的關(guān)系。提出基于相位一致性和方向?yàn)V波的冕環(huán)震蕩識(shí)別和跟蹤方法。識(shí)別過(guò)程如下：首先，使用相位一致性技術(shù)對(duì)太陽(yáng)圖像中的冕環(huán)特征進(jìn)行識(shí)別；其次，對(duì)識(shí)別后的圖像進(jìn)行二值化處理；接著使用方向?yàn)V波去除非冕環(huán)的結(jié)構(gòu)特征；然后，對(duì)其進(jìn)行形態(tài)學(xué)處理；最后，再對(duì)識(shí)別結(jié)果進(jìn)行二次曲線擬合。運(yùn)用所提算法對(duì)2011年9月6日由太陽(yáng)動(dòng)力學(xué)觀測(cè)站(Solar Dynamics Observatory, SDO)在Fe IX 17.1 nm波段觀測(cè)到的冕環(huán)震蕩進(jìn)行了識(shí)別和跟蹤，結(jié)果驗(yàn)證了所提算法的可行性和準(zhǔn)確性。證明所提算法能被用于冕環(huán)特性以及震蕩過(guò)程的研究和分析。

關(guān)鍵詞：相位一致性；方向?yàn)V波；冕環(huán)震蕩；震蕩周期

日冕位于太陽(yáng)大氣的最外層，在上面可以觀測(cè)到一些亮的環(huán)形結(jié)構(gòu)，它們被稱為冕環(huán)。冕環(huán)由日冕磁場(chǎng)束縛的熱等離子體發(fā)亮而形成環(huán)狀結(jié)構(gòu)，是太陽(yáng)磁場(chǎng)結(jié)構(gòu)的一種表現(xiàn)形式。它們同日冕物質(zhì)拋射、太陽(yáng)耀斑、太陽(yáng)風(fēng)暴等有密切的關(guān)系[1]。研究冕環(huán)及其演化過(guò)程對(duì)于理解太陽(yáng)大氣層的磁場(chǎng)活動(dòng)和預(yù)報(bào)太陽(yáng)風(fēng)暴具有重要意義。

太陽(yáng)耀斑爆發(fā)和日冕物質(zhì)拋射釋放的大量能量引起附近局部磁場(chǎng)結(jié)構(gòu)的急劇變化，導(dǎo)致附近的冕環(huán)結(jié)構(gòu)產(chǎn)生震蕩[2-4]。這些震蕩不僅僅發(fā)生在一條冕環(huán)結(jié)構(gòu)上，也同時(shí)出現(xiàn)在多條冕環(huán)結(jié)構(gòu)上。最近的研究表明，這種震蕩是駐波的一種表現(xiàn)形式，但是每個(gè)環(huán)上由于震蕩的頻率、相位以及衰減系數(shù)的不同，往往呈現(xiàn)不同的震蕩特征[5-9]。研究這些特征對(duì)于更好地理解耀斑爆發(fā)和日冕物質(zhì)拋射對(duì)日冕磁場(chǎng)的影響具有很大的幫助。

然而準(zhǔn)確地識(shí)別和提取冕環(huán)特征結(jié)構(gòu)面臨非常大的困難。這是由于：(1)冕環(huán)是由日冕磁場(chǎng)束縛的熱等離子體形成，這些等離子體分布并不均勻，雖然冕環(huán)結(jié)構(gòu)的強(qiáng)度較大，但強(qiáng)度變化(梯度)較小，導(dǎo)致冕環(huán)結(jié)構(gòu)的邊緣模糊，造成提取困難；(2)各個(gè)冕環(huán)之間的距離較近而且相互交錯(cuò)，導(dǎo)致產(chǎn)生的二維圖像環(huán)形結(jié)構(gòu)彼此交叉重疊；(3)由于磁場(chǎng)活動(dòng)的復(fù)雜性導(dǎo)致了冕環(huán)結(jié)構(gòu)呈現(xiàn)多種形態(tài)，如環(huán)形結(jié)構(gòu)出現(xiàn)斷裂等問(wèn)題；(4)在冕環(huán)附近還存在非環(huán)形的太陽(yáng)特征，如海綿狀亮斑，這些都增加了冕環(huán)提取的復(fù)雜性和準(zhǔn)確性。

目前冕環(huán)的識(shí)別方法主要分為兩類[10-13]：一類是運(yùn)用圖像空間域的強(qiáng)度和梯度閾值對(duì)特征結(jié)構(gòu)進(jìn)行識(shí)別和提取；另一類是利用頻域或小波域?qū)Σ煌念l率特征進(jìn)行濾波，識(shí)別和提取冕環(huán)特征。廣義上講，圖像的結(jié)構(gòu)特征，如特征邊緣、尺度和紋理等大量底層信息存儲(chǔ)在圖像的相位譜中。文[14-15]證明了相位一致性在表征特征結(jié)構(gòu)時(shí)對(duì)強(qiáng)度、對(duì)比度、照度具有不變性，而且通過(guò)計(jì)算圖像的局部能量實(shí)現(xiàn)了對(duì)圖像特征的穩(wěn)定提取。目前，相位一致性技術(shù)已經(jīng)應(yīng)用于黑子本影亮點(diǎn)和本影耀斑以及光球米粒等低對(duì)比度太陽(yáng)圖像特征的提取和識(shí)別[16-19]。

本文運(yùn)用相位一致性技術(shù)并結(jié)合方向?yàn)V波的方法對(duì)一個(gè)由太陽(yáng)動(dòng)力學(xué)觀測(cè)站搭載的太陽(yáng)大氣成像儀(Atmospheric Imaging Assembly, AIA)[20]在2011年9月6日由譜線Fe IX 17.1 nm波段觀測(cè)到的冕環(huán)震蕩進(jìn)行識(shí)別和分析，并同已有的文獻(xiàn)結(jié)果進(jìn)行比較來(lái)檢驗(yàn)所提方法的可行性和準(zhǔn)確性。

論文首先介紹所用的實(shí)驗(yàn)數(shù)據(jù)，其次簡(jiǎn)述冕環(huán)震蕩的識(shí)別過(guò)程，然后對(duì)識(shí)別的冕環(huán)震蕩進(jìn)行統(tǒng)計(jì)分析，并同已有文獻(xiàn)結(jié)果作比較，最后對(duì)本文進(jìn)行總結(jié)。

1實(shí)驗(yàn)數(shù)據(jù)

選擇震蕩事件發(fā)生在編號(hào)為11283的活動(dòng)區(qū)為實(shí)驗(yàn)數(shù)據(jù)，時(shí)間從2011年9月6日22時(shí)17分到2011年9月6日22時(shí)34分。觀測(cè)資料由太陽(yáng)大氣成像儀在譜線Fe IX 17.1 nm波段記錄。冕環(huán)震蕩由一個(gè)附近爆發(fā)的X2.1級(jí)耀斑引起。該組數(shù)據(jù)共包含65幅圖像，間隔時(shí)間為12 s，圖像的像元分辨率是0.6″。圖1(a)展示了2011年9月6日22時(shí)17分獲得的圖像，其中白色箭頭所示的位置為耀斑爆發(fā)的位置，它引起了圖1(a)中白色方框標(biāo)記區(qū)(即圖1(b))冕環(huán)產(chǎn)生了震蕩。圖1(b)展示了這個(gè)區(qū)域的細(xì)節(jié)。它的大小為210 × 144″。在圖1(b)中，由數(shù)字1和2指示的位置為耀斑爆發(fā)能量釋放所產(chǎn)生的圖像結(jié)構(gòu)特征，它導(dǎo)致了附近的冕環(huán)產(chǎn)生劇烈震蕩，尤其是箭頭LoopA和LoopB指示的兩個(gè)冕環(huán)結(jié)構(gòu)震蕩尤為明顯。為了同文[9]做比較，本文主要針對(duì)冕環(huán)A和冕環(huán)B的震蕩運(yùn)用所提算法進(jìn)行識(shí)別和跟蹤，并對(duì)其結(jié)果進(jìn)行統(tǒng)計(jì)分析來(lái)檢驗(yàn)算法的可行性和準(zhǔn)確性。

圖1冕環(huán)

Fig.1Loops

2識(shí)別和跟蹤方法

從圖2(b)可以看出，那些海綿狀亮斑和由于耀斑爆發(fā)產(chǎn)生的近似垂直的一些特征(如圖中上方中間位置被白色方框標(biāo)記)干擾冕環(huán)A和B的提取。所要識(shí)別和提取的冕環(huán)A和冕環(huán)B方向近似水平。因此算法選用方向?yàn)V波思想識(shí)別近似垂直的特征并去除。算法選用sobel邊緣檢測(cè)算子g(g=[1, 0, -1; 2, 0, -2; 1, 0, -1])作為方向?yàn)V波算子檢測(cè)干擾線段的垂直邊緣并去除它們，處理結(jié)果如圖2(c)。

從圖2(c)可以發(fā)現(xiàn)，那些近似垂直的干擾特征通過(guò)使用方向?yàn)V波的方式已經(jīng)被很準(zhǔn)確地去除。然而所需提取的冕環(huán)結(jié)構(gòu)仍然存在不連續(xù)以及冕環(huán)上存在小的孔洞且邊緣不規(guī)則、噪聲點(diǎn)大和冕環(huán)粗細(xì)不均勻等問(wèn)題。算法采用形態(tài)學(xué)的去除小面積、填充、膨脹、閉運(yùn)算、細(xì)化等一系列操作對(duì)上述問(wèn)題進(jìn)行處理，處理結(jié)果如圖2(d)。

從圖2(d)箭頭處可以看出，識(shí)別的冕環(huán)結(jié)構(gòu)不光滑，存在許多分支結(jié)構(gòu)，并且圖中除了所要識(shí)別的冕環(huán)A和冕環(huán)B外，還有一些短的非目標(biāo)線。因此算法根據(jù)節(jié)點(diǎn)處各分枝的長(zhǎng)短，保留最長(zhǎng)的分枝結(jié)構(gòu)，去除短的分枝雜點(diǎn)，并根據(jù)環(huán)的長(zhǎng)度，提取了冕環(huán)A和冕環(huán)B，結(jié)果如圖2(e)。

最后對(duì)所識(shí)別的65幅圖像中的冕環(huán)A和冕環(huán)B進(jìn)行了二次方程的曲線擬合，獲取一個(gè)光滑的弧線，其二次方程定義如下：

y=ax2+bx+c,

(1)

其中，x為圖像中冕環(huán)的橫坐標(biāo)；y為圖像中冕環(huán)的縱坐標(biāo)；a、b、c為參數(shù)。圖2(e)中冕環(huán)A和冕環(huán)B二次擬合參數(shù)a、b、c的值如表1。圖2(e)環(huán)的擬合結(jié)果被疊加到原圖上并顯示在圖2(f)中。由于冕環(huán)震蕩導(dǎo)致每一幅圖像中冕環(huán)A和冕環(huán)B的形態(tài)均不相同，因此本文分別對(duì)其進(jìn)行了逐一擬合，其它的擬合參數(shù)并未在此列出，這些擬合線為本文要進(jìn)行統(tǒng)計(jì)分析的特征。

圖2識(shí)別過(guò)程

Fig.2The process of Identification

3實(shí)驗(yàn)結(jié)果與分析

為了驗(yàn)證所識(shí)別冕環(huán)的結(jié)果，分別選取冕環(huán)A和冕環(huán)B上在X軸方向上10個(gè)位置，然后根據(jù)每幅圖的二次擬合方程(1)得到Y(jié)值，并對(duì)這10個(gè)Y值分別取它們?cè)诿?幅圖像上的均值作為這一幅圖像冕環(huán)震蕩的位置。圖3(a)和(b)中星號(hào)標(biāo)記的位置即為冕環(huán)A和冕環(huán)B的Y均值位置隨時(shí)間T的變化情況。在圖3中，橫軸表示時(shí)間(單位為分鐘)，縱軸表示冕環(huán)移動(dòng)距離(單位為千米)。從圖中可以看出這些星號(hào)呈現(xiàn)一個(gè)類似衰減的正弦函數(shù)，因此采用如下函數(shù)擬合這些位移，具體的衰減正弦擬合函數(shù)定義如下：

表1　對(duì)圖1(b)中冕環(huán)A和B的二次擬合參數(shù)

(2)

其中，A為振幅；τ為衰減時(shí)間；t為時(shí)間(單位為分鐘)；T為周期；φ是初相位(單位為度)。

圖3　震蕩

對(duì)于這兩個(gè)曲線的震蕩周期，是采用傅里葉變換的頻率分析方法獲得。這兩個(gè)振蕩信號(hào)的頻譜分析圖如圖4。依據(jù)傅里葉頻譜分析，這兩個(gè)環(huán)的震蕩頻率約為8.4 MHz，周期為118 s。運(yùn)用這個(gè)周期值T(即118 s)和擬合方程(2)擬合圖3中的星號(hào)位置。擬合的曲線如圖3(a)和(b)，它們的衰減系數(shù)約為122 s，但兩個(gè)震蕩序列的相位差近似41°，冕環(huán)B先于冕環(huán)A，而且B環(huán)的振幅也高于A環(huán)。這表明由耀斑所引起的冕環(huán)震蕩不是從最近位置的冕環(huán)震蕩開(kāi)始，而是導(dǎo)致較遠(yuǎn)位置的冕環(huán)首先震蕩。不過(guò)這種位置的遠(yuǎn)近也可能是由成像的視角導(dǎo)致的。這個(gè)分析結(jié)果同文[9]的結(jié)論一致，進(jìn)一步驗(yàn)證了本文所提的方法可以運(yùn)用于冕環(huán)特征的研究和分析。

圖4頻譜分析圖

Fig.4Spectrum analysis diagram

4討論與總結(jié)

大多數(shù)的冕環(huán)震蕩是由附近的耀斑爆發(fā)和日冕物質(zhì)拋射引起的，對(duì)這些冕環(huán)的震蕩進(jìn)行分析可以更好地理解冕環(huán)的磁場(chǎng)特征與太陽(yáng)耀斑爆發(fā)和日冕物質(zhì)拋射的能量釋放之間的關(guān)系，準(zhǔn)確地識(shí)別和提取冕環(huán)特征至關(guān)重要。

本文運(yùn)用相位一致性和方向?yàn)V波的方法分析了一個(gè)由太陽(yáng)大氣成像儀在17.1 nm波段觀測(cè)到的冕環(huán)震蕩。然而，由于冕環(huán)結(jié)構(gòu)的復(fù)雜性，本文僅選取了震蕩區(qū)域中兩個(gè)比較明顯的冕環(huán)結(jié)構(gòu)運(yùn)用相位一致性和方向?yàn)V波的圖像特征提取方法識(shí)別和追蹤它們，并根據(jù)識(shí)別和追蹤的結(jié)果對(duì)這兩條冕環(huán)的周期、相位差以及衰減系數(shù)進(jìn)行統(tǒng)計(jì)分析，并同已有文獻(xiàn)的結(jié)果進(jìn)行比較，證明算法可行，識(shí)別結(jié)果準(zhǔn)確。

盡管相位一致性可以很好地解決低信噪比條件下圖像特征提取的問(wèn)題，但由于冕環(huán)特征存在于復(fù)雜的背景之中以及冕環(huán)結(jié)構(gòu)邊緣模糊和信噪比低等因素，使得提取的冕環(huán)結(jié)構(gòu)不太完整，一些比較細(xì)、暗的冕環(huán)結(jié)構(gòu)未能識(shí)別成功。對(duì)于這些情況僅僅從圖像處理的角度識(shí)別和提取是不夠的，需要結(jié)合其他波段的數(shù)據(jù)及磁場(chǎng)反演等物理過(guò)程分析，提高冕環(huán)識(shí)別的準(zhǔn)確性，這也是下一階段研究工作的重點(diǎn)。

參考文獻(xiàn)：

[1]Vaiana G S, Rosner R. Recent advances in coronal physics[J]. Annual Review of Astronomy and Astrophysics, 1978, 16: 393-428.

[2]Aschwanden M J. The 3D geometry, motion, and hydrodynamic aspects of oscillating coronal loops[J]. Space Science Reviews, 2009, 149(1-4): 31-64.

[3]White R S, Verwichte E. Transverse coronal loop oscillations seen in unprecedented detail by AIA/SDO[J]. Astronomy & Astrophysics, 2012, 537: A49-A58.

[4]Wang T J, Ofman L, Davila J M, et al. Growing transverse oscillations of a multistranded loop observed by SDO/AIA[J]. The Astrophysical Journal Letters, 2012, 751(2): L27-L32.

[5]Aschwanden M J, Fletcher L, Schrijver C J, et al. Coronal loop oscillations observed with the transition region and coronal explorer[J]. The Astrophysical Journal, 1999, 520(2): 880-894.

[6]Luna M, Terradas J, Oliver R, et al. Transverse oscillations of two coronal loops[J]. The Astrophysical Journal, 2008, 676(1): 717-727.

[7]Luna M, Terradas J, Oliver R, et al. Transverse oscillations of a multi-stranded loop[J]. The Astrophysical Journal, 2010, 716(2): 1371-1380.

[8]Terradas J, Arregui I, Oliver R, et al. Resonant absorption in complicated plasma configurations: applications to multistranded coronal loop oscillations[J]. The Astrophysical Journal, 2008, 679(2): 1611-1620.

[9]Jain R, Maurya1 R A, Hindman B W. Fundamental-mode oscillations of two coronal loops within a solar magnetic arcade[J]. The Astrophysical Journal Letters, 2015, 804(1): L19-L27.

[10]Lee J K, Newman T S, Gary G A. Oriented connectivity-based method for segmenting solar loops[J]. Pattern Recognition, 2006, 39(2): 246-259.

[11]Carcedo L, Brown D S, Hood A W, et al. A quantitative method to optimise magnetic field line fitting of observed coronal loops[J]. Solar Physics, 2003, 218(1): 29-40.

[12]Steger C. An unbiased detector of curvilinear structures[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1998, 20(2): 113-125.

[13]James McAteer R T, Kestener P, Arneodo A, et al. Automated detection of coronal loops using a wavelet transform modulus maxima method [J]. Solar Physics, 2010, 262(2): 387-397.

[14]Kovesi P. Image features from phase congruency[J]. Journal of Computer Vision Research, 1999, 1(3): 1-27.

[15]Kovesi P. Phase congruency: a low-level image invariant[J]. Psychological Research, 2000, 64(2): 136-148.

[16]Feng Song, Xu Zhi, Wang Feng, et al. Automated detection of low-contrast solar features using the phase-congruency algorithm[J]. Solar Physics, 2014, 289(10): 3985-3994.

[17]Feng Song, Yu Lan, Yang Yunfei, et al. Identification of emission sources of umbral flashes using phase congruency[J].Research in Astronomy and Astrophysics, 2014, 14(8): 1001-1010.

[18]Feng Song, Zhao Yan, Yang Yunfei, et al. Identifying and tracking of peripheral and central umbral dots[J].Solar Physics, 2015, 290(4): 1119-1133.

[19]韓翠翠, 蔣霞, 楊云飛, 等. 基于NVST觀測(cè)的米粒識(shí)別和形態(tài)特征分析[J]. 科學(xué)通報(bào), 2016, 61(8): 881-890.

Han Cuicui, Jiang Xia, Yang Yunfei, et al. Identification and morphology analysis of granules using the New Vacuum Solar Telephone Observations[J]. Chinese Science Bulletin, 2016, 61(8): 881-890.

[20]Lemen J R, Title A M, Akin D J, et al. The Atmospheric Imaging Assembly (AIA) on the Solar Dynamics Observatory (SDO) [M]// The solar dynamics observatory. New York: Springer US, 2012: 17-40.

*基金項(xiàng)目：國(guó)家自然科學(xué)基金 (U1531132, 11463003, 11573012, 11303011) 資助.

收稿日期：2015-09-01；

修訂日期：2015-10-11

作者簡(jiǎn)介：馮松，男，副教授. 研究方向：天文技術(shù)，數(shù)字圖像處理. Email: ynkmfs@escience.cn

中圖分類號(hào)：TP391.41

文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào):1672-7673(2016)03-0345-06

Identification and Tracking of Coronal Loop Oscillation Using Phase Congruency and Directional Filter

Feng Song1,2, Wang Yang1,2, Zhao Yan1,2, Yuan Quan1,2, Yang Yunfei1,2, Deng Hui1,2, Wang Feng1,2

(1. Faculty of Information Engineering and Automation, Kunming University of Science and Technology, Kunming 650500, China; 2. Computer Technology Application Key Laboratory of Yunnan Province, Kunming University of Science and Technology, Kunming 650500, China, Email: ynkmfs@escience.cn)

Abstract:Solar flares occurring usually cause oscillations of nearby coronal loops. A precondition of studying these oscillations is to explore how to identify coronal loops. Meanwhile, accurate identification is also crucial to study and analyze the statistical features of the coronal loops and the relationship to coronal magnetic field. However, the traditional identification methods based on intensity and/or gradient threshold have some limitations, because coronal loops hold a relatively low contrast, blurred boundary, and the image noise is close to coronal loops. So we propose a new method to identify coronal loops using phase congruency and directional filter. The identification procedures are as follows: firstly, phase congruency technique is used to identify linear shape in a coronal image; secondly, a binary image is obtained from the image of phase congruency for detecting the loop features; thirdly, a directional filter based on Sobel Operator is applied to remove the redundant features that are nearly perpendicular to those target loops; fourthly, morphological filters are used to extract the shape of those loops; finally, quadratic fitting is used to smooth the identified loop features. For illustrating the performance of the method, the coronal loop oscillations taken with the Atmospheric Imaging Assembly on aboard the Solar Dynamic Observatory on September 11, 2011 are used to evaluate the identification and tracking process. The coronal loop oscillations were caused by a nearby X2.1 flare emission. We analyzed the identified results, and calculated their oscillation period, phase difference and decay coefficient. They are in good agreement with those of previous works. This implies that the proposed method can accurately identify the shape of a loop in an image; furthermore, the results can also be used to study the physical properties of coronal loop oscillations and the relationship to others solar phenomena.

Key words:Phase congruency; Directional filter; Coronal loop oscillation; Oscillation period

CN 53-1189/PISSN 1672-7673