SKA中CLEAN算法的優(yōu)化研究

錢慎一，劉慧慧

（鄭州輕工業(yè)大學(xué)計(jì)算機(jī)與通信工程學(xué)院，鄭州 450000）

0 引言

平方公里陣列（Square Kilometer Array，SKA）是一項(xiàng)大型國(guó)際科研合作項(xiàng)目，比目前最大的射電望遠(yuǎn)鏡JVLA的靈敏度提高約50倍[1]、巡天速度提高10000倍，由2500個(gè)反射面天線[2]、250組中頻和低頻孔徑陣列組成[3]，分布范圍大于3000公里[4]，形成旋臂陣列望遠(yuǎn)鏡，其接收面積達(dá)一平方公里，能夠探測(cè)距離地球數(shù)百萬(wàn)甚至數(shù)十億光年的物體的無(wú)線電波。隨著SKA項(xiàng)目的實(shí)施，海量的天文數(shù)據(jù)迅速增加，如何有效、快速地處理這些數(shù)據(jù)成為了一個(gè)研究重點(diǎn)。ARL算法參考庫(kù)是SKA的候選算法庫(kù)，數(shù)據(jù)處理模塊是基于Python語(yǔ)言實(shí)現(xiàn)的，主要功能是創(chuàng)建校準(zhǔn)模型、對(duì)數(shù)據(jù)進(jìn)行預(yù)處理、成像等功能。

CLEAN算法是ARL中的一個(gè)重要算法，在射電天文成像中起著至關(guān)重要的作用。目前的射電望遠(yuǎn)鏡天線數(shù)量有限，導(dǎo)致了空間頻率覆蓋不完整[5]，從而對(duì)最終的圖像的構(gòu)建造成影響。由于CLEAN算法可以消除該影響，因此很多學(xué)者對(duì)CLEAN算法做了相關(guān)研究。1974年，Hogbom[6]首次提出了CLEAN算法，一種非線性的迭代方法，用以消除旁瓣干擾。該算法雖然能消除旁瓣的干擾，但需要消耗大量的時(shí)間，運(yùn)算效率低。2004年，Bhatnagar[7]提出了一種用于無(wú)線電干涉圖像的尺度敏感的反卷積算法（Asp-clean）算法，將圖像建模為自適應(yīng)尺度像素的集合，可以更準(zhǔn)確地重建非對(duì)稱結(jié)構(gòu)的天空?qǐng)D像，但是增加了算法復(fù)雜性和計(jì)算成本，其計(jì)算時(shí)間是CLEAN算法的三倍多。2008年，Tim J.Cornwell[8]提出了一種Multiscale CLEAN算法，可以更好地處理擴(kuò)展源，提高圖像質(zhì)量，但是該算法比Hogbom CLEAN算法需要更長(zhǎng)的運(yùn)行時(shí)間。2011年，U.Rau[9]在Multiscale CLEAN算法的基礎(chǔ)上，結(jié)合了multi-frequency算法，可以在更高的靈敏度和采樣頻率的情況下，重建天空?qǐng)D像。2011年，Sarod Yatawatta[10]將扁長(zhǎng)橢球波函數(shù)應(yīng)用到CLEAN算法中，減少去除擴(kuò)展源的誤差。2015年，A.La Camera[11]提出了將圖像域劃分成等平面區(qū)域，并對(duì)每個(gè)區(qū)域應(yīng)用具有邊界效應(yīng)校正的CLEAN方法。將該方法應(yīng)用于以空間可變PSF為特征的恒星系統(tǒng)的模擬圖像，通過(guò)對(duì)去模糊的圖像進(jìn)行光圈測(cè)光，獲得了良好的圖像質(zhì)量。2017年，Jun Cheng[12]提出了一種小波CLEAN算法，以解決傳統(tǒng)CLEAN算法出現(xiàn)的問(wèn)題，并對(duì)小波濾波器的參數(shù)進(jìn)行了優(yōu)化，提高了圖像質(zhì)量，但沒(méi)有縮短運(yùn)行時(shí)間，未提高算法運(yùn)行效率。大部分學(xué)者都是針對(duì)圖像質(zhì)量問(wèn)題進(jìn)行的研究，但是在運(yùn)算效率方面進(jìn)行的研究很少。如今隨著海量天文數(shù)據(jù)的產(chǎn)生，CLEAN算法的運(yùn)算效率的提高變得刻不容緩。

1 CLEAN算法分析

1.1 Multi-Scale CLEAN算法

2008年，Cornwell提出 Multi-Scale CLEAN 算法[13]，是對(duì)Hogbom CLEAN算法的一種改進(jìn)，能夠較好地處理展源。Multi-Scale CLEAN算法步驟如下所示：

Step1：初始化模型等于0，殘差圖等于“臟”圖。模型：IM=0；殘差圖：IR=ID。

Step2：對(duì)于每個(gè)尺度αq，計(jì)算尺度卷積殘余和尺度偏差 S（α）=1-0.6α/α；對(duì)于每

max對(duì)尺度αp，αq，計(jì)算交叉項(xiàng) B*m（αp）*m（αq）。

Step3：對(duì)于每個(gè)尺度，找到殘差圖的峰值及其位置。

Step4：在乘以比例相關(guān)的偏差項(xiàng)之后，選擇殘差最大的尺度。

Step5：將此組件添加到當(dāng)前模型中，按循環(huán)增益縮放。

Step7：通過(guò)潔束對(duì)當(dāng)前模型進(jìn)行卷積：BG*IM

Step8：添加殘差圖，得到恢復(fù)的圖像：BG*IM+IR

在整個(gè)算法流程中，尺度的處理至關(guān)重要。對(duì)于尺度范圍的選擇，沒(méi)有明確的定義，但這通常不是太重要。像等差數(shù)列這樣的尺度范圍會(huì)浪費(fèi)計(jì)算時(shí)間，太大的范圍又會(huì)導(dǎo)致收斂性差。所以我們通常選取的尺度范圍一般是等比數(shù)列的形式，如 0、2、4、8、16、32，或者是 0、3、10、30 這樣的尺度范圍。

1.2 Clark CLEAN算法

CLEAN算法是在天文圖像復(fù)原中最常用的去卷積算法之一，但是由于其計(jì)算量較大，因此Clark提出了一種改進(jìn)的CLEAN算法，可以降低數(shù)據(jù)處理過(guò)程的計(jì)算工作量[14]。CLEAN算法中的大部分計(jì)算在于移動(dòng)和縮放“臟”束，而Clark CLEAN算法用二維傅里葉變換（Fast Fourier Transformation，F(xiàn)FT）可以更有效地完成任務(wù)，并且只使用一部分“臟”束就可以找到點(diǎn)源的近似位置和強(qiáng)度。

Clark CLEAN算法提出了兩種簡(jiǎn)化方式，第一種是只考慮圖像的最大值，第二種是只考慮“臟”束的中心部分。由于“臟”束是對(duì)稱的，所以只需要存儲(chǔ)一個(gè)邊長(zhǎng)為2：1的矩形區(qū)域。為了定位組件，將圖像之外的“臟”束計(jì)為零。這些簡(jiǎn)化大大減少了CLEAN算法的計(jì)算工作量。

詳細(xì)來(lái)說(shuō)，Clark CLEAN算法有兩個(gè)周期，稱為次周期和主周期。次周期的過(guò)程如下：

Step1：如果圖像峰值點(diǎn)的強(qiáng)度大于“臟”束外部的最大旁瓣，則從“臟”圖中選擇該點(diǎn)，并記錄其位置坐標(biāo)。

Step2：把挑選出來(lái)的部分“臟”束與“臟”圖中找到的點(diǎn)進(jìn)行潔化操作。

Step3：當(dāng)殘差圖像的點(diǎn)小于外部的最大旁瓣，則停止循環(huán)；否則，跳轉(zhuǎn)到Step1繼續(xù)執(zhí)行循環(huán)。

算法的主周期是首先找到“臟”束的主瓣和外部的最大旁瓣；其次將次循環(huán)找到的點(diǎn)進(jìn)行快速傅立葉變換，并與經(jīng)過(guò)快速傅里葉逆變換（IFFT）的“臟”束相乘，再經(jīng)過(guò)IFFT得到亮度分布；接著從“臟”圖中減去該亮度分布；最后判斷是否達(dá)到閾值，如果達(dá)到閾值，則跳出主循環(huán)，否則跳轉(zhuǎn)到主循環(huán)的第一步繼續(xù)執(zhí)行循環(huán)。在一定程度上，光束在次周期引入的誤差可以在隨后的次周期中得到修正。

影響算法計(jì)算時(shí)間的因素有圖像大小、增益值、圖像的復(fù)雜性和旁瓣等。由于這些因素的存在，計(jì)算時(shí)間的變化很大，改變其中一個(gè)因素，都會(huì)對(duì)時(shí)間造成一定的影響。

2 CLEAN算法優(yōu)化

2.1 算法描述與實(shí)現(xiàn)

Multi-Scale CLEAN算法是CLEAN算法的一個(gè)很好的改進(jìn)算法，可以更好地處理展源，使最終的復(fù)原圖像更加的清晰，可以獲得比CLEAN算法更好的效果。但是由于其算法實(shí)現(xiàn)較為復(fù)雜，比CLEAN算法多了尺度的計(jì)算部分，需要更多的處理步驟，因此算法的運(yùn)行時(shí)間較長(zhǎng)。為了減少算法的運(yùn)行時(shí)間，本文提出了一種Multi-Scale CLEAN算法與Clark CLEAN算法相結(jié)合的算法，在文中稱為MS-Clark CLEAN算法。

MS-Clark CLEAN算法的主要思想就是結(jié)合Multi-Scale CLEAN算法與Clark CLEAN算法的優(yōu)點(diǎn)，來(lái)實(shí)現(xiàn)一種運(yùn)算時(shí)間快、圖像質(zhì)量高的算法。將Clark CLEAN算法中選取部分“臟”束的思想用于Multi-Scale CLEAN算法。在保證圖像質(zhì)量的同時(shí)減少運(yùn)算時(shí)間。

MS-Clark CLEAN算法的部分偽代碼如表1所示，該部分偽代碼表示的是選取部分“臟”束。根據(jù)分析Clark CLEAN算法可知，可以將“臟”束的主瓣作為有效束代替“臟”束進(jìn)行迭代。

MS-Clark CLEAN算法的部分偽代碼：

MS-Clark CLEAN算法的流程圖如圖1所示。首先是輸入“臟”圖，“臟”束的主瓣；接著，計(jì)算每個(gè)尺度的偏差以及交叉項(xiàng)；其次，在“臟”圖中找到最大值點(diǎn)及其位置坐標(biāo)；然后，判斷殘差圖像的最大值是否小于給定的噪聲水平，若滿足條件，則跳出循環(huán)，繼續(xù)向下執(zhí)行，否則，繼續(xù)找殘差圖像中大于最大旁瓣的點(diǎn)；接著，對(duì)滿足條件的點(diǎn)組成的函數(shù)圖和潔束進(jìn)行卷積；最后加上殘差圖像，就得到了最終的圖像，也就是復(fù)原的清晰圖像。

圖1 MS-Clark CLEAN算法流程圖

2.2 參數(shù)選取

為了達(dá)到最好的復(fù)原效果，需要選取合適的參數(shù)，通過(guò)分析MS-Clark CLEAN算法可知，增益值、“臟”束的大小、尺度以及迭代次數(shù)都會(huì)對(duì)最終的實(shí)驗(yàn)效果造成一定的影響。本小節(jié)主要針對(duì)以上4個(gè)參數(shù)分別進(jìn)行了實(shí)驗(yàn)，最終選取了最合適的參數(shù)值。在實(shí)驗(yàn)中采用峰值信噪比[15]（Peak Signal to Noise Ratio，PSNR）來(lái)衡量算法的效果。

在選取合適參數(shù)的實(shí)驗(yàn)中，首先進(jìn)行的是增益值的選取，通過(guò)分析比較發(fā)現(xiàn)，在0-1范圍內(nèi)，增益值為0.2左右時(shí)，圖像的PSNR值最小，此時(shí)的圖像相對(duì)來(lái)說(shuō)最清晰，圖像質(zhì)量最好，因此本文選取增益值為0.2。如果增益值選擇的不合適，反而會(huì)得不到質(zhì)量最好的圖像。

第二個(gè)調(diào)整的參數(shù)是“臟”束的大小，如圖2（a）所示是“臟”束的大小與PSNR值的關(guān)系，如圖2（b）所示是“臟”束的大小對(duì)時(shí)間的影響。從圖中可以看出，時(shí)間隨著“臟”束的變大而增加，在“臟”束大小為300×300時(shí)，變化巨大，其后，趨于平衡。結(jié)合“臟”束的大小與PSNR值的關(guān)系，可以發(fā)現(xiàn)在“臟”束大小為300×300后，PSNR值變化不大，時(shí)間的變化也不大，都是趨于平緩的。在“臟”束大小為200×200到300×300之間，PSNR值的變化不是很大，可是時(shí)間的變化很大。因此“臟”束大小可以選擇在200×200到250×250之間的尺寸，本文選取的大小為200×200。

圖2

接下來(lái)對(duì)尺度這個(gè)參數(shù)進(jìn)行實(shí)驗(yàn)，根據(jù)具體數(shù)據(jù)發(fā)現(xiàn)，0,3,10,30這組尺度PSNR值最小，圖像較為清晰，圖像質(zhì)量較高。因此，本文選取的尺度范圍為0,3,10,30。

最后對(duì)迭代次數(shù)進(jìn)行分析實(shí)驗(yàn)，在迭代次數(shù)為1000左右時(shí)，可以得到較低的PSNR值，此時(shí)的圖像較為清晰，圖像質(zhì)量較好。迭代次數(shù)太少，找出的“臟”圖中的點(diǎn)源就少，影響最終的圖像質(zhì)量。迭代次數(shù)大于1000次時(shí)，PSNR值并沒(méi)有較大的改變，但迭代次數(shù)增多會(huì)使計(jì)算時(shí)間增長(zhǎng)，計(jì)算代價(jià)變大。因此，本文選擇的迭代次數(shù)為1000。

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

實(shí)驗(yàn)通過(guò)對(duì)MS-Clark CLEAN算法進(jìn)行測(cè)試，同時(shí)與優(yōu)化之前的Multi-Scale算法進(jìn)行對(duì)比。主要采用以主觀觀察評(píng)價(jià)為主，客觀量化評(píng)價(jià)標(biāo)準(zhǔn)為輔的圖像質(zhì)量評(píng)價(jià)方法?？陀^量化評(píng)價(jià)標(biāo)準(zhǔn)采用PSNR來(lái)評(píng)價(jià)。

為了測(cè)試算法的性能，進(jìn)行了如下的圖像實(shí)驗(yàn)。選取的圖像大小為 512×512，圖 3（a）為“臟”圖，圖 3（b）為CLEAN算法的復(fù)原結(jié)果，圖3（c）為Multi-Scale算法的復(fù)原結(jié)果，圖3（d）為MS-Clark CLEAN算法的復(fù)原結(jié)果。通過(guò)主觀觀察可以看出，CLEAN算法對(duì)于圖像的邊緣區(qū)域處理的效果不好，復(fù)原的點(diǎn)源較少。Multi-Scale算法比CLEAN算法對(duì)于邊緣區(qū)域的處理效果更好，復(fù)原的點(diǎn)源較多。MS-Clark CLEAN算法的處理效果比CLEAN算法的好，但是對(duì)于邊緣區(qū)域的處理，沒(méi)有Multi-Scale算法的效果好。

圖3

圖4 兩種算法在不同圖像大小下的PSNR值

在圖像質(zhì)量方面，Multi-Scale算法和MS-Clark CLEAN算法的結(jié)果相似，但是在時(shí)間方面，MS-Clark CLEAN算法的運(yùn)行時(shí)間要低于Multi-Scale算法的運(yùn)行時(shí)間。兩種算法的運(yùn)行時(shí)間如表1所示。從表中可以看出，Multi-Scale算法的運(yùn)行時(shí)間隨著圖像大小的增大而急劇增多，在圖像大小為8192×8192時(shí)，運(yùn)行時(shí)間達(dá)到了一個(gè)小時(shí)十五分鐘左右，每運(yùn)行一次耗時(shí)巨大，代價(jià)太大，不利于進(jìn)行大圖像的測(cè)試。而MSClark CLEAN算法在圖像大小為8192×8192時(shí)，運(yùn)行時(shí)間為33分鐘左右，大大降低了算法的運(yùn)行時(shí)間。

表1 兩種算法的運(yùn)行時(shí)間

4 結(jié)語(yǔ)

本文提出了一種MS-Clark CLEAN算法，該算法在Multi-Scale CLEAN算法的基礎(chǔ)上，結(jié)合了Clark CLEAN算法的優(yōu)點(diǎn)，將“臟”束的主瓣作為有效束代替“臟”束進(jìn)行迭代，并通過(guò)調(diào)整參數(shù)使圖像復(fù)原得到較好的結(jié)果。由于選取了部分“臟”束，使得在保證圖像質(zhì)量的情況下，提高了運(yùn)算速度。實(shí)驗(yàn)結(jié)果證實(shí)了改進(jìn)后的算法的可行性和有效性。