沐 芳，黃 歡 ，尚振宏,3，強(qiáng)振平

(1. 昆明理工大學(xué)信息工程與自動化學(xué)院，云南 昆明 650500；2. 西南林業(yè)大學(xué)大數(shù)據(jù)與智能工程學(xué)院，云南 昆明 650224；3. 昆明理工大學(xué)云南省人工智能重點(diǎn)實(shí)驗(yàn)室，云南 昆明 650500)

為了研究高度動態(tài)的太陽日冕結(jié)構(gòu)，人們一直在提高時空分辨率方面不斷探索。由于太陽日冕結(jié)構(gòu)的密度、位置和形狀幾秒鐘就會發(fā)生變化，較大的動態(tài)范圍造成經(jīng)過望遠(yuǎn)鏡成像的日冕結(jié)構(gòu)邊緣模糊且噪聲較大，圖像可視化困難。目前有幾種較為有效的增強(qiáng)圖像算法用于太陽日冕圖像增強(qiáng)?；谥狈綀D均衡的方法可以利用不同分布直方圖對不同圖像規(guī)定化處理[1]，也可以高斯濾波后，再獨(dú)立均衡各個分區(qū)，如保持亮度的動態(tài)直方圖均衡化[2](Brightness Preserving Dynamic Histogram Equalization, BPDHE)，這些方法能夠有效提高圖像整體的對比度。但由于直方圖的鄰域通常是正方形，導(dǎo)致嚴(yán)重的偽影，從而影響日冕環(huán)狀結(jié)構(gòu)；其次會放大圖像的低對比度噪聲，圖像細(xì)節(jié)模糊；亮度不同的區(qū)域之間也會失去有效結(jié)構(gòu)。基于曝光融合框架的對比度增強(qiáng)算法[3](Exposure Fusion Framework, EFF)可以解決直方圖均衡出現(xiàn)的對比度過高或過低的問題。基于多種經(jīng)典圖像增強(qiáng)算法融合的方法[4](Multi-deviation Fusion method, MF)對多種圖像增強(qiáng)算法取長補(bǔ)短，圖像增強(qiáng)時能夠兼顧細(xì)節(jié)增強(qiáng)、對比度提高和主觀感知效果等。這些方法的預(yù)處理大多采用對數(shù)變換，能夠降低計(jì)算復(fù)雜度。但由于對數(shù)變換抑制了亮區(qū)梯度幅度的變化，可能導(dǎo)致某些區(qū)域的精細(xì)結(jié)構(gòu)細(xì)節(jié)丟失。為解決這個問題，基于Retinex的方法在壓縮動態(tài)范圍的同時，增強(qiáng)邊緣并保持顏色穩(wěn)定，其中基于Retinex對噪聲魯棒的方法(Robust Retinex Model, rR)[5]可提高對噪聲的魯棒性，使用一種基于拉格朗日乘子的交替方向最小化算法代替對數(shù)變換，有效解決模型優(yōu)化過程中的問題，能有效保持細(xì)節(jié)結(jié)構(gòu)和邊緣穩(wěn)定。上述方法都僅為太陽日冕圖像提供了一種可視化方法，對于其中因?yàn)槿彰釀討B(tài)變化、望遠(yuǎn)鏡抖動、CCD欠采樣以及雜散光等因素造成的精細(xì)結(jié)構(gòu)模糊問題，無法有效增強(qiáng)。隨著數(shù)字圖像處理技術(shù)的發(fā)展，盲退卷積算法為恢復(fù)精細(xì)結(jié)構(gòu)的準(zhǔn)確性和真實(shí)性提供了科學(xué)的理論依據(jù)。文[6]根據(jù)望遠(yuǎn)鏡設(shè)計(jì)理論估計(jì)滿足大氣成像儀不同波段的點(diǎn)擴(kuò)散函數(shù)用于退卷積，但造成圖像模糊的原因除傳感器設(shè)備外，還有太陽活動位移、形狀變化等，所以理論估計(jì)的點(diǎn)擴(kuò)散函數(shù)并不能很好地退卷積得出清晰圖像。在使用空間望遠(yuǎn)鏡大氣成像儀觀測太陽日冕結(jié)構(gòu)時，對由于日冕位移、形狀改變、CCD欠采樣以及雜散光等因素造成成像模糊，可以建模為

B(x)=I(x)*k(x)+n(x)，

(1)

其中，*表示卷積運(yùn)算；B(x)表示獲得的模糊圖像，即大氣成像儀圖像；I(x)表示潛在的清晰圖像；k(x)表示未知的點(diǎn)擴(kuò)散函數(shù)；n(x)表示在每個像素x=(x,y)處相同且獨(dú)立分布的噪聲項(xiàng)。

根據(jù)求解方法的不同，盲退卷積算法大致可以分為兩類：(1)通過引入額外約束條件同時對目標(biāo)圖像和點(diǎn)擴(kuò)散函數(shù)估計(jì)進(jìn)行交替迭代，并基于最大后驗(yàn)概率估計(jì)(Maximum A Posterior, MAP)完成盲退卷積的算法；(2)先估計(jì)模糊核，再退卷積得到清晰圖像。第1類方法對噪聲有一定的抑制作用，但涉及目標(biāo)圖像和點(diǎn)擴(kuò)散函數(shù)估計(jì)的多次迭代計(jì)算，非常耗時，且嚴(yán)重依賴于清晰圖像和點(diǎn)擴(kuò)散函數(shù)的先驗(yàn)特性。對于大氣成像儀圖像而言，很難獲得清晰圖像的統(tǒng)計(jì)先驗(yàn)信息。第2類方法以基于顯著性邊緣[7]和功率譜的方法[8-9]為主，基于邊緣的方法是從原圖顯著邊緣估計(jì)模糊核，這類方法很大程度上取決于邊緣選擇，更適用于運(yùn)動模糊估計(jì)。大氣成像儀圖像結(jié)構(gòu)模糊，無法較好地選擇邊緣?；诠β首V的方法直接從輸入圖像功率譜中恢復(fù)點(diǎn)擴(kuò)散函數(shù)的功率譜，再使用相位提取算法求解點(diǎn)擴(kuò)散函數(shù)的相位。

因此，本文研究了一種利用模糊圖像功率譜估計(jì)點(diǎn)擴(kuò)散函數(shù)的單幀圖像盲退卷積算法，用于太陽大氣成像儀圖像增強(qiáng)。首先對數(shù)據(jù)預(yù)處理(包括去噪和灰度值變換處理)，從圖像中估計(jì)點(diǎn)擴(kuò)散函數(shù)的功率譜；再用相位提取算法恢復(fù)點(diǎn)擴(kuò)散函數(shù)，為了使算法對噪聲更魯棒，額外引入了權(quán)重參數(shù)α，限制使用點(diǎn)擴(kuò)散函數(shù)功率譜對相位估計(jì)的硬性約束，并對算法隨機(jī)初始化求解多次，避免陷入局部極小值。計(jì)算點(diǎn)擴(kuò)散函數(shù)所有可能解之間的L2范數(shù)距離，通過多次隨機(jī)初始化計(jì)算點(diǎn)擴(kuò)散函數(shù)，如果可能解之間的L2距離越小，可以認(rèn)為是最接近真解的點(diǎn)擴(kuò)散函數(shù)，最后對(1)式利用點(diǎn)擴(kuò)散函數(shù)退卷積得出目標(biāo)圖像。為進(jìn)一步減少噪聲影響，將待增強(qiáng)圖像分隔為4幅相鄰局部子圖，并按上述方法估計(jì)每幅子圖的點(diǎn)擴(kuò)散函數(shù)，對4個點(diǎn)擴(kuò)散函數(shù)進(jìn)行歸一化和求取均值得到整幅圖像的點(diǎn)擴(kuò)散函數(shù)。實(shí)驗(yàn)分析使用輪廓切片曲線、功率譜分析、點(diǎn)擴(kuò)散函數(shù)分析以及和其他較為先進(jìn)的增強(qiáng)算法結(jié)果對比分析，對重建的目標(biāo)圖像細(xì)節(jié)結(jié)構(gòu)進(jìn)行客觀評價(jià)，將高分辨率日冕成像儀(The High-Resolution Coronal Imager, Hi-C)圖像功率譜假設(shè)模型驗(yàn)證和精細(xì)結(jié)構(gòu)增強(qiáng)效果對比。結(jié)果顯示本文方法對大氣成像儀19.3 nm，21.1 nm和17.1 nm波段圖像的增強(qiáng)和復(fù)原有較好的效果，且對噪聲魯棒。此外，由于該方法不涉及目標(biāo)像的反復(fù)重建，在圖像恢復(fù)速度上也有一定的優(yōu)勢。

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

太陽動力學(xué)天文臺/大氣成像儀[10]主要目標(biāo)是了解太陽能量存儲與釋放機(jī)制，確定太陽如何以及為何發(fā)生變化。通過研究太陽周圍等離子體的相互作用，增強(qiáng)我們對太陽如何影響地球大氣和社會生產(chǎn)的理解。大氣成像儀于2010年2月11日正式啟動，由4臺望遠(yuǎn)鏡組成，可同時提供10個波段的全日面像，包括7個極紫外波段(9.4 nm，13.1 nm，17.1 nm，19.3 nm，21.1 nm，30.4 nm和33.5 nm)、2個紫外波段(160 nm和170 nm)和1個可見光波段(450 nm)；觀測范圍包括太陽過渡區(qū)和最高可達(dá)0.5R⊙的日冕，空間分辨率為1.5″，時間分辨率為12 s，視場達(dá)到41 × 41′。

高分辨率日冕成像儀[11-12]提供了日冕更精細(xì)結(jié)構(gòu)視角，獲得高空間分辨率(≈0.3～0.4″)，高時間分辨率(≈5 s)的太陽活動區(qū)域結(jié)構(gòu)圖像，是目前為止可獲得最高分辨率的太陽日冕儀，為了解小尺度日冕特征和發(fā)生在過渡區(qū)的活動現(xiàn)象，以及它們之間可能存在的相互關(guān)系提供依據(jù)。

由表1數(shù)據(jù)對比可知，高分辨率日冕成像儀雖然擁有更高的時空分辨率，但視場小，拍攝時長較短，無法對日冕進(jìn)行長期穩(wěn)定的觀測。相比而言，大氣成像儀能夠提供大視場、長期全日面觀測圖像，然而圖像質(zhì)量低于高分辨率日冕成像儀，因此對大氣成像儀圖像增強(qiáng)處理顯得尤為重要。

表1 大氣成像儀與高分辨率日冕成像儀數(shù)據(jù)特點(diǎn)對比Table 1 Comparison of SDO/AIA and Hi-C data characteristics

本文采用的大氣成像儀數(shù)據(jù)是2012年7月11日18時53分32秒可以觀察到熱耀斑的等離子體和日冕19.3 nm波段的局部像，和2014年2月2日23時40分11秒觀測活動區(qū)日冕21.1 nm波段的局部像，這兩組數(shù)據(jù)由同一望遠(yuǎn)鏡拍攝。此外，實(shí)驗(yàn)中還包括2014年10月4日10時13分11秒不同望遠(yuǎn)鏡拍攝的寧靜區(qū)日冕17.1 nm波段的局部像。其中19.3 nm圖像實(shí)驗(yàn)中，也采用了活動區(qū)域11520的2012年7月11日18時53分27秒高分辨率日冕成像儀圖像作為清晰圖像功率譜假設(shè)驗(yàn)證和實(shí)驗(yàn)對比，除此之外，實(shí)驗(yàn)對比圖還包括文[6]估計(jì)的點(diǎn)擴(kuò)散函數(shù)退卷積圖。數(shù)據(jù)預(yù)處理包括：根據(jù)大氣成像儀數(shù)據(jù)高動態(tài)的特點(diǎn)，對圖像對數(shù)變換，這一步能夠有效降低圖像的動態(tài)范圍；然后使用中值濾波對圖像去噪，該算法在去除部分噪聲的同時不會平滑過渡，能有效保持圖像的結(jié)構(gòu)。

2 盲退卷積算法理論

本文方法分為3部分：自然圖像功率譜假設(shè)模型、高分辨率日冕成像儀功率譜分析和大氣成像儀圖像點(diǎn)擴(kuò)散函數(shù)求解。

2.1 自然圖像功率譜假設(shè)模型

文[8]的研究表明，清晰圖像的功率譜(或自相關(guān))是各向同性的，模糊過程中圖像的中高頻分量降低，模糊圖像不同方向微分濾波后再計(jì)算自相關(guān)，當(dāng)模糊方向與微分方向一致時，損失最大，此時能量絕對值之和最小的方向?yàn)槟：较?。根?jù)以上結(jié)論，我們對清晰圖像和模糊圖像分別進(jìn)行不同角度(θ: 1～180°)的傅里葉中心投影變換，得到不同投影角度對應(yīng)的一維切片后，再對這些切片一維微分濾波, 最后計(jì)算其自相關(guān)得到

(Rd*Ρθ(I))(x)≈cθ·δ(x)，

(2)

(Rd*Ρθ(B))(x)≈cθ·RΡθ(k)(x)，

(3)

其中，d為一維微分濾波器；cθ為一個隨角度變化的乘性因子；R為對圖像的自相關(guān)計(jì)算；Ρ為投影算子；δ為沖激函數(shù)，即

(4)

2.2 高分辨率日冕成像儀圖像功率譜分析

高分辨率日冕成像儀圖像的功率譜特性分析如圖1。將圖1(a)的一個局部子圖的Rd*Ρθ(I)繪制得到圖1(b)，圖1(b)45°，90°，135°和180°方向功率譜截面得到圖1(c)，將圖1(c)中4條曲線對應(yīng)數(shù)值相除得到圖1(d)。通過觀察圖1(c)，清晰圖像圖1(a)經(jīng)過(2)式計(jì)算得到的功率譜并不是一個精確的δ函數(shù)。由圖1(d)可以看出，不同角度的功率譜相除得到不同常數(shù)，cθ隨著投影角度的變化而變化，即太陽清晰圖像功率譜隨著θ不同有不同乘性因子的變化量。進(jìn)而證明(2)式和(3)式用于大氣成像儀圖像點(diǎn)擴(kuò)散函數(shù)估計(jì)是可靠的。圖中多個尖峰是由于除數(shù)存在0值造成的。

圖1 高分辨率日冕成像儀功率譜模型說明Fig.1 Explanation of Hi-C power spectrum model

為了不增加圖像的總能量，濾波器d均值為0，在進(jìn)行d*Ρθ(B)計(jì)算時，會導(dǎo)致Ρθ(B)的均值缺失(為0)，但圖像本身功率譜并非為0，所以在計(jì)算過程中，引入?yún)?shù)mθ代替0值。即RΡθ(k)(x)的均值使用mθ/cθ代替。令fθ(x)=Rd*Ρθ(B)(x)，由(3)式可以得到

RΡθ(k)(x)=[fθ(x)+mθ]/cθ.

(5)

2.3 點(diǎn)擴(kuò)散函數(shù)估計(jì)

(5)式是一個欠定方程，為了成功恢復(fù)RΡθ(k)(x)，需要引入額外約束條件估計(jì)cθ和mθ。

2.3.1 估計(jì)cθ

我們采用關(guān)于點(diǎn)擴(kuò)散函數(shù)的3個假設(shè)[13]對(5)式進(jìn)行等價(jià)計(jì)算：

(1)相機(jī)曝光，光子在傳感器上累積圖像，由于光子能量非負(fù)，故點(diǎn)擴(kuò)散函數(shù)是非負(fù)的，即k≥0，所以關(guān)于k的投影Ρθ(k)和自相關(guān)RΡθ(k)都是非負(fù)的，根據(jù)(5)式可以得到fθ(x)+mθ≥0。

綜上，(5)式可以等價(jià)于

(6)

2.3.2 估計(jì)支持范圍sθ

通過在點(diǎn)擴(kuò)散函數(shù)估計(jì)和sθ估計(jì)之間進(jìn)行最大期望(Expectation-maximization, EM)迭代，從而估計(jì)sθ和點(diǎn)擴(kuò)散函數(shù)。最大期望算法的每次迭代是對點(diǎn)擴(kuò)散函數(shù)的再一次修正，能夠更準(zhǔn)確估計(jì)點(diǎn)擴(kuò)散函數(shù)，使其在盲復(fù)原工作中有效抑制振鈴和偽影。具體步驟為

圖2(a)為支持域sθ迭代路徑，圖2(b)為根據(jù)最終sθ迭代結(jié)果恢復(fù)的點(diǎn)擴(kuò)散函數(shù)功率譜。圖2(a)中藍(lán)色曲線表示初始化的sθ，根據(jù)最大期望迭代更新sθ(紅色曲線)，達(dá)到終止條件后sθ(綠色曲線)的數(shù)值用于作為圖2(b)圖橫軸的支持域數(shù)值。

圖2 (a)sθ迭代路徑；(b)一維點(diǎn)擴(kuò)散函數(shù)功率譜Fig.2 sθ Iterative path (a) and one dimensional PSF power spectrum (b)

算法1(求解k)的步驟可歸納為

輸入：模糊圖像B；

(1)fθ(x)=Rd*Ρθ(B)。

(3)外循環(huán)開始，令次數(shù)t=1。

(4)內(nèi)循環(huán)開始：

①使用給定的sθ，根據(jù)(6)式計(jì)算模糊核功率譜|k^|2(^ 表示對信號做傅里葉變換)；

②通過算法2計(jì)算模糊核kt；

③最大期望迭代更新sθ。

(5)當(dāng)滿足sθ>0.1max(RΡθ(k))時內(nèi)循環(huán)結(jié)束。

(6)t=Touter(T表示剪切圖像張數(shù))時外循環(huán)結(jié)束。

圖3是本文算法流程，首先對原圖(a)去噪和對數(shù)變化得到圖3(b)，計(jì)算圖3(b)的一維自相關(guān)圖3(c)，根據(jù)圖3(c)迭代支持域sθ圖3(d)，最后使用該sθ校準(zhǔn)點(diǎn)擴(kuò)散函數(shù)一維功率譜圖3(e)，復(fù)原為二維功率譜圖3(f)，相位提取算法計(jì)算相位后恢復(fù)點(diǎn)擴(kuò)散函數(shù)圖3(g)，最后非盲退卷積得到目標(biāo)圖像圖3(h)。

圖3 本文算法的執(zhí)行過程：(a)輸入大氣成像儀觀測的局部太陽模糊圖；(b)對數(shù)增強(qiáng)后的圖；(c)一維自相關(guān)投影；(d)sθ迭代路徑；(e)sθ優(yōu)化后的一維自相關(guān)投影；(f)恢復(fù)后的二維模糊核功率譜；(g)相位恢復(fù)后的模糊核；(h)本文算法增強(qiáng)后的太陽圖

2.4 點(diǎn)擴(kuò)散函數(shù)的相位恢復(fù)

求出|k^|2之后，我們使用相位提取算法，通過添加約束條件求解模糊核k的相位分量φ(ω)，最終得到模糊核k。相位提取算法是一種欠定算法，該算法的解存在平凡模糊的問題，通常會收斂到平凡解或局部極小值。典型的相位恢復(fù)方法以誤差下降算法(Error-Reduction, ER)為例，在迭代過程中降低誤差，但將不滿足約束的信號設(shè)為0值的方式常常使迭代停滯，混合輸入輸出(Hybrid Input-Output, HIO)算法[14]改進(jìn)了此問題，對解設(shè)置一個小擾動，可以使迭代跳出停滯現(xiàn)象快速收斂。文[15]提出松弛平均交替反射(Relaxed Averaged Alternating Reflection, RAAR)算法，該方法能夠使收斂更快更平穩(wěn)，但將輸入的功率譜作為一種硬性約束，即需要恢復(fù)的信號完全滿足輸入的功率譜，當(dāng)輸入信號含有噪聲時，會使收斂震蕩，魯棒性欠佳。大氣成像儀圖像噪聲復(fù)雜，并不能在數(shù)據(jù)預(yù)處理階段完全消除，故本文算法放寬了模糊核功率譜的硬性約束要求(通過調(diào)整權(quán)重參數(shù)α)，可提高算法對噪聲的魯棒性。初始階段采用隨機(jī)初始化相位分量，多次重復(fù)此隨機(jī)過程，可以避免收斂到局部極小值；其次，引入混合輸入輸出算法對解的擾動方法(即算法2中β項(xiàng))，可以使解快速收斂。最后計(jì)算所有可能解的L2距離，求出最優(yōu)解。

算法2(求解kt)迭代包括以下步驟：

輸入：ρ(ω)=|k^(ω)|，點(diǎn)擴(kuò)散函數(shù)尺寸s；

輸出：kt。

(1)外循環(huán)n=1：

①隨機(jī)初始化相位φ(ω)([-π,π] 之間)；

③內(nèi)循環(huán)m=1：

(a)施加傅里葉域約束(即是對幅值ρ放寬至α倍)：

2=-1{[αρ+(1-α)|^|]ei·φ^}；

(b)施加空間域約束(即點(diǎn)擴(kuò)散函數(shù)非負(fù)性約束和支持域約束)：

(7)

其中，Ω={x: 2g2(x)-g(x)<0}∪{x:x?[0,s]×[0,s]}，為使m次迭代中，信號能夠平穩(wěn)快速地收斂，本文采用混合輸入輸出算法中對β的設(shè)置，即β=β0+(1-β0)[1-e(-m/7)3]。本實(shí)驗(yàn)中β0=0.75，當(dāng)出現(xiàn)x∈Ω時，為賦予一個較小的值β(x)+(1-2β)2(x)，以便能夠逃出迭代過程中的停滯現(xiàn)象；

④m=Mouter時，內(nèi)循環(huán)結(jié)束；

(2)n=Nouter時，外循環(huán)結(jié)束，kn兩兩計(jì)算L2距離，得到最優(yōu)解kt。

算法2包括兩個循環(huán)結(jié)構(gòu)，設(shè)置外循環(huán)為了避免解陷入局部極小值，初始化值越接近真值，越容易收斂到全局最小值。多次相位隨機(jī)初始化求解k(x)，最終通過求不同解之間的范數(shù)距離，范數(shù)距離最小的解可認(rèn)為是最優(yōu)解k(x)。設(shè)置內(nèi)循環(huán)能夠使算法平穩(wěn)快速地收斂，在迭代求解過程中交替引入頻域約束(幅值約束)和空間域約束(點(diǎn)擴(kuò)散函數(shù)的非負(fù)性約束和支持域約束)，當(dāng)出現(xiàn)x∈Ω時，為賦予一個較小的值β(x)+(1-2β)2(x)，以便能夠避免迭代過程中的停滯現(xiàn)象。文[15]指出，通過逐步增大β值能夠引導(dǎo)解穩(wěn)定收斂。故本文通過內(nèi)循環(huán)將β值設(shè)置為一個不斷接近1的動態(tài)值。

如(1)式，已知模糊核k(x)和模糊圖像B(x)，求解目標(biāo)清晰圖像的方法叫作非盲退卷積算法。實(shí)驗(yàn)使用文[16]中交替最小化的方法，即半二次分裂方法來求解目標(biāo)圖像I(x)。

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

為了評估該方法的有效性和準(zhǔn)確性，本文實(shí)驗(yàn)中對大氣成像儀19.3 nm，21.1 nm和17.1 nm 3個波段不同區(qū)域進(jìn)行增強(qiáng)，圖像均通過對數(shù)顯示。圖4中3個波段的圖像在增強(qiáng)后，細(xì)節(jié)結(jié)構(gòu)更銳，尤其冕環(huán)之間的相對模糊結(jié)構(gòu)增強(qiáng)后能較好地分離，邊界更為清晰。由第3列本文算法估計(jì)的點(diǎn)擴(kuò)散函數(shù)尺寸大小可知，19.3 nm和21.1 nm波段點(diǎn)擴(kuò)散函數(shù)較17.1 nm更大，這說明前兩個波段圖像模糊程度更高，由于19.3 nm和21.1 nm波段的圖像主要是拍攝太陽活動區(qū)得到的，活動區(qū)更大的動態(tài)范圍和冕環(huán)相對望遠(yuǎn)鏡的運(yùn)動造成圖像模糊程度更高。在第4列功率譜圖中，3個波段圖像增強(qiáng)后的功率譜曲線相對大氣成像儀圖像曲線在中頻都有一定程度的提高。中頻曲線的提高在增強(qiáng)圖像中表現(xiàn)為已有結(jié)構(gòu)邊界更銳，以及由于模糊沒觀察到的細(xì)節(jié)結(jié)構(gòu)也得以恢復(fù)。

圖4 大氣成像儀不同波段增強(qiáng)圖Fig.4 Enhancement images of AIA in different bands

為了更好地對比細(xì)節(jié)結(jié)構(gòu)的增強(qiáng)情況，我們在圖5中引入結(jié)構(gòu)的輪廓切片圖。文[6]根據(jù)大氣成像儀的成像特點(diǎn)估計(jì)不同波段點(diǎn)擴(kuò)散函數(shù)用于增強(qiáng)圖像，太陽一個波段的所有圖像使用同一個點(diǎn)擴(kuò)散函數(shù)退卷積得到，這是不準(zhǔn)確的，圖像模糊原因也只考慮望遠(yuǎn)鏡設(shè)備的因素，所以恢復(fù)結(jié)果不佳(如圖5第2列)。圖5中，我們對不同波段不同的區(qū)域分別增強(qiáng)，并將本文估計(jì)的點(diǎn)擴(kuò)散函數(shù)與文[6]點(diǎn)擴(kuò)散函數(shù)作對比。如19.3 nm波段的兩幅局部圖，本文估計(jì)的點(diǎn)擴(kuò)散函數(shù)尺寸相當(dāng)，但形狀不同，對于結(jié)構(gòu)A和結(jié)構(gòu)B，冕環(huán)能較好地分離，如輪廓切片圖所示。同樣，在21.1 nm和17.1 nm波段的增強(qiáng)圖中，C，D和E結(jié)構(gòu)的模糊方向與點(diǎn)擴(kuò)散函數(shù)形狀一致，不同的區(qū)域，太陽活動不一致，冕環(huán)位置和形狀的變化，相互干擾是不可控的。這也是必須分區(qū)域單獨(dú)估計(jì)點(diǎn)擴(kuò)散函數(shù)退卷積圖像的原因。

圖5 不同波段不同區(qū)域細(xì)節(jié)增強(qiáng)圖Fig.5 Detail enhancement images of different bands and regions

如圖6，通過引入更高分辨率的日冕成像儀圖像進(jìn)行細(xì)節(jié)結(jié)構(gòu)對比，評估本文增強(qiáng)方法的準(zhǔn)確性。為了便于與高分辨率日冕成像儀圖對比，圖6(a)，(b)經(jīng)過放大和對數(shù)處理。觀察三者的能量圖可知，大氣成像儀圖像細(xì)節(jié)區(qū)域顏色相近，模糊明顯，幾乎看不到細(xì)節(jié)結(jié)構(gòu)，增強(qiáng)后的圖像與高分辨率日冕成像儀圖像相近。通過三者的功率譜圖圖6(d)可以看出，大氣成像儀幾乎沒有高頻細(xì)節(jié)，中頻部分也很弱。經(jīng)過本文方法增強(qiáng)后的圖像，中高頻都有一定的提高，略弱于高分辨率日冕成像儀圖像，說明恢復(fù)細(xì)節(jié)結(jié)構(gòu)的數(shù)量少于高分辨率日冕成像儀圖像，中頻處有一點(diǎn)略高于高分辨率日冕成像儀圖像。在細(xì)節(jié)圖6(e)中可觀察到，本文方法不僅能較好地分離冕環(huán)的相對模糊結(jié)構(gòu)，且對于已存在的冕環(huán)，能更好地保持邊緣。細(xì)節(jié)結(jié)構(gòu)邊緣較高分辨率日冕成像儀尖銳，所以中頻曲線略高于高分辨率日冕成像儀。圖中的A，B和C都是由于模糊無法分辨的結(jié)構(gòu)，增強(qiáng)后，與高分辨率日冕成像儀對比可知，結(jié)構(gòu)保持穩(wěn)定，且恢復(fù)情況較為準(zhǔn)確。如對應(yīng)的輪廓切片圖6(f)A結(jié)構(gòu)，大氣成像儀中表現(xiàn)為一個冕環(huán)的結(jié)構(gòu)，成功地分離為兩條冕環(huán)，與高分辨率日冕成像儀保持一致。

圖6 與高分辨率日冕成像儀細(xì)節(jié)對比圖Fig.6 Detail comparison with Hi-C

圖7給出了經(jīng)典的基于直方圖的BPDHE算法，以及目前較為先進(jìn)的EFF，MF和rR算法的實(shí)驗(yàn)結(jié)果。從第1行局部放大細(xì)節(jié)中可以看出，BPDHE方法使圖像局部過度增強(qiáng)且細(xì)節(jié)模糊，本文算法的增強(qiáng)圖對細(xì)節(jié)增強(qiáng)效果更好，邊緣結(jié)構(gòu)也較為穩(wěn)定。第2行為各個算法對應(yīng)的直方圖，根據(jù)直方圖像素值的分布情況可以看出，MF算法的可視化效果最好。對比其他增強(qiáng)算法，本文算法在保證對比度和亮度增強(qiáng)的同時，能復(fù)原更精細(xì)的細(xì)節(jié)結(jié)構(gòu)。如表2，我們對不同算法分別計(jì)算了標(biāo)準(zhǔn)差和信息熵，這兩者可以有效地反映圖像的對比度和亮度信息，其值越大，表示圖像信息豐富。亮度誤差[4](Lightness Order Error, LOE)能夠衡量圖像的亮度失真情況，值越小代表亮度越好。結(jié)構(gòu)相似度(Structural Similarity, SSIM)可以用來衡量圖像恢復(fù)質(zhì)量，是圖像亮度、對比度和結(jié)構(gòu)對比三者評估的指標(biāo)，其值越接近1，圖像質(zhì)量越好。表中的方法均與高分辨率日冕成像儀圖像作對比計(jì)算結(jié)構(gòu)相似度。相比而言，雖然本文算法運(yùn)行時間較長，但對細(xì)節(jié)結(jié)構(gòu)的增強(qiáng)更為準(zhǔn)確。

4 結(jié) 論

本文提出了一種基于盲退卷積的太陽大氣成像儀圖像質(zhì)量增強(qiáng)方法，基于圖像功率譜特性的點(diǎn)擴(kuò)散函數(shù)估計(jì)模型，從大氣成像儀圖像中估計(jì)點(diǎn)擴(kuò)散函數(shù)，最后退卷積得出增強(qiáng)圖像。文中采用輪廓切片分析評估方法的準(zhǔn)確性，功率譜分析評估方法的有效性，不同區(qū)域點(diǎn)擴(kuò)散函數(shù)分析評估方法的合理性。實(shí)驗(yàn)結(jié)果表明，本文方法對大氣成像儀19.3 nm，17.1 nm和21.1 nm波段圖像增強(qiáng)效果較好，尤其紋理細(xì)節(jié)和邊緣結(jié)構(gòu)的增強(qiáng)效果最好。其次，本文算法也有一定的噪聲抑制效果，得益于點(diǎn)擴(kuò)散函數(shù)相位估計(jì)過程的松弛設(shè)計(jì)，恢復(fù)新細(xì)節(jié)結(jié)構(gòu)也較為準(zhǔn)確，因?yàn)檫@三個波段的圖像噪聲較小，日冕結(jié)構(gòu)特征也更豐富和明顯，能夠更為準(zhǔn)確地估計(jì)模糊核。對比其他圖像增強(qiáng)方法，本文方法的實(shí)時性較差，點(diǎn)擴(kuò)散函數(shù)估計(jì)和退卷積的最優(yōu)化迭代過程是主要的耗時原因。將來可以通過提高計(jì)算設(shè)備的性能和采用較為快速且有效的最優(yōu)化方法來提高算法的實(shí)時性。在大氣成像儀的其他波段圖像中，9.4 nm，13.1 nm和33.5 nm波段圖像噪聲很大，異常點(diǎn)較多，結(jié)構(gòu)相對模糊導(dǎo)致增強(qiáng)效果不佳，將來能夠通過研究更好的去噪算法和對噪聲更加魯棒的退卷積算法對其他波段的圖像進(jìn)行有效增強(qiáng)。