魏志晴,鄭文康,白艷萍,譚秀輝,程 蓉,胡紅萍,王 鵬

(中北大學(xué) 數(shù)學(xué)學(xué)院,山西 太原 030051)

0 引 言

傳統(tǒng)計(jì)算機(jī)斷層成像算法可分為解析法和迭代法。解析法由于投影數(shù)據(jù)等角度間隔缺失,重建圖像中會(huì)引入條紋狀偽影[1],且圖像細(xì)節(jié)嚴(yán)重模糊,重建圖像的噪聲較大,無(wú)法應(yīng)用于實(shí)際診斷。迭代法可以將先驗(yàn)信息融合到重建圖像過(guò)程當(dāng)中,具有較好的重建效果。斷層成像中應(yīng)用最廣泛的迭代算法是代數(shù)重建技術(shù)(Algebraic Reconstruction Technique,ART)、聯(lián)合代數(shù)重建技術(shù)(Simultaneous Algebraic Reconstruction Technique,SART)和期望最大化(Expectation Maximization,EM)。稀疏角度背景下,一般選取迭代重建法進(jìn)行圖像重建。近年來(lái)隨著壓縮感知理論的提出,基于稀疏先驗(yàn)的正則化方法得到廣泛關(guān)注,目前已經(jīng)成為主流[2]。其核心思想是利用重建圖像的稀疏先驗(yàn)信息來(lái)提高CT圖像重建質(zhì)量,在解決稀疏角度問(wèn)題中獲得了較好的效果。主要包括基于局部平滑的方法、基于字典的正則化方法和基于非局部相似的方法。

基于局部平滑的方法,主要是利用全變分(Total Variation,TV)模型對(duì)圖像進(jìn)行局部平滑。Candes等人[3]最先提出了TV正則化模型,此模型在圖像處理等領(lǐng)域可以有效抑制噪聲并保留圖像的邊緣、紋理等。Sidky等人[4]在Candes等人的研究成果上,設(shè)計(jì)了基于約束TV最小化的CT圖像重建模型,在稀疏角度采樣背景下獲得了較好的重建效果。但是TV正則化方法對(duì)于具有更多圖像紋理細(xì)節(jié)和復(fù)雜邊緣結(jié)構(gòu)的圖像,會(huì)使得重建的圖像過(guò)度平滑,導(dǎo)致部分細(xì)節(jié)丟失。

基于字典的正則化方法,如離散余弦變換等,其單一基函數(shù)無(wú)法滿(mǎn)足圖像信號(hào)的自適應(yīng)稀疏表示,因此在此基礎(chǔ)上出現(xiàn)了許多改進(jìn)方法。Chen等人[2]提出了一種基于TV與字典學(xué)習(xí)方法的磁共振成像算法,在保持精細(xì)結(jié)構(gòu)方面具有比較好的效果,進(jìn)一步提高了重建圖像質(zhì)量。Fan[5]提出了基于雙字典的稀疏角度MRI圖像重建方法,不僅提升了圖像質(zhì)量,而且收斂速度更快。但字典正則項(xiàng)的選取對(duì)重建圖像的質(zhì)量起到?jīng)Q定性作用,若選取不當(dāng),易造成某些特征被過(guò)度懲罰,不能保留微小的特征,使圖像細(xì)節(jié)模糊。

基于非局部相似方法主要是挖掘圖像在紋理等結(jié)構(gòu)信息具有重復(fù)性的特性。Buades等人[6]提出了著名的非局部均值(Nonlocal Mean,NLM)濾波算法。Huang[7]將NLM算法應(yīng)用在CT圖像ART算法的重建上,得到了較好的效果。崔樹(shù)輝[8]將BM3D (Block Matching and 3D collaborative filtering)加入到ART算法中,在稀疏角度情況下依舊取得了較高分辨率。Zhang[9]提出了群稀疏表示正則化(Group-Sparse Regularization,GSR)算法,具有良好的捕捉圖像奇異性的優(yōu)良特性,且單幅圖像中的稀疏性和非局部相似性是同時(shí)存在的,可以很好地達(dá)到對(duì)稀疏角度下SART重建圖像去除偽影的目的。

為了在平滑圖像的同時(shí)保留不同層次的結(jié)構(gòu),Zhang等人[10]提出一種用于保持邊緣的濾波器,即滾動(dòng)引導(dǎo)濾波(Rolling Guidance Filter,RGF)。裴佩佩等人[11]將RGF與卷積稀疏表示相結(jié)合,應(yīng)用在圖像處理中,克服了邊緣模糊,保留了對(duì)比度和邊緣紋理信息。

受上述研究的啟發(fā),該文提出了一種稀疏角度下基于滾動(dòng)引導(dǎo)濾波和群稀疏表示正則化的SART重建算法,稱(chēng)為SART-GSR-RGF。該算法首先對(duì)稀疏角度采樣下的圖像進(jìn)行FBP重建,作為SART算法的初始圖像;然后采用SART算法進(jìn)行重建,再經(jīng)過(guò)群稀疏正則化處理;最后采用滾動(dòng)引導(dǎo)濾波提升圖像的對(duì)比度,將所得圖像再次作為SART重建算法的輸入,繼續(xù)迭代重建,直到得到分辨率較高的圖像。這樣既利用了GSR的保護(hù)邊緣和平滑圖像的特點(diǎn),又利用了RGF保留不同層次結(jié)構(gòu)同時(shí)傳遞圖像特征的特點(diǎn)。經(jīng)過(guò)實(shí)驗(yàn)驗(yàn)證,所提算法對(duì)于圖像分辨率的提升具有很大作用。

1 理論依據(jù)

1.1 FBP重建和SART重建

FBP(Filtered Back-Projection)重建算法的依據(jù)是中心切片定理[12]。即空域進(jìn)行投影得到的一維傅里葉變換,是空域進(jìn)行二維傅里葉變換的一個(gè)切片。通過(guò)取遍所有投影值,再進(jìn)行傅里葉反變換得到空間域分布的圖像,就可以得到二維函數(shù)。FBP公式表示如下[13]：

(1)

上述公式表明,重建圖像f(x,y)在某一位置的值,是通過(guò)該點(diǎn)的所有濾波投影采樣的疊加。流程可以簡(jiǎn)化為圖1。

圖1 FBP重建示意圖

SART算法[14]是CT圖像重建中的一類(lèi)迭代重建算法,在ART和SIRT算法的基礎(chǔ)上,進(jìn)行一次更新只用到了一個(gè)投影采樣角度下的投影數(shù)據(jù),也即該投影角度下所對(duì)應(yīng)的系統(tǒng)方程。在修正某個(gè)像素值之前,需要計(jì)算一下經(jīng)過(guò)像素的所有射線(xiàn)的誤差,以對(duì)像素進(jìn)行修正,并進(jìn)行加權(quán)和歸一化處理,再把上述的結(jié)果更新到該像素上去。重復(fù)此過(guò)程,直至滿(mǎn)足收斂條件。算法公式為：

(2)

圖2 SART算法示意圖

其中灰色區(qū)域部分即為SART、SIRT和ART算法步驟中變化的部分。就圖像重建處理時(shí)間而言,解析法比迭代算法具有更好的性能。就降噪和有限數(shù)據(jù)而言,迭代法優(yōu)于解析法,因此該文先用FBP對(duì)平行透射數(shù)據(jù)進(jìn)行成像,之后作為SART迭代重建圖像的初值。

1.2 基于群的SR建模

圖像處理的稀疏表示(Sparse Representations,SR)就是將圖像信號(hào)變換到另一個(gè)變換域中,絕大部分變換系數(shù)的絕對(duì)值很小,所得到的變換向量是稀疏的或者近似稀疏的,可以將其作為圖像的簡(jiǎn)潔表達(dá)。SR模型可以表示為：

(3)

其中,x表示觀測(cè)到的信號(hào)向量,D是字典,α表示信號(hào)的系數(shù),λ是正則化參數(shù),‖·‖0表示l0范數(shù)。SR的目標(biāo)是為訓(xùn)練好的D尋找一個(gè)稀疏向量α來(lái)表示x。為了更好地用α表示x,必須選擇一個(gè)有效的字典D。已經(jīng)提出了一些近似算法來(lái)交替優(yōu)化D和α,如MOD,k-SVD和在線(xiàn)學(xué)習(xí)。

1.3 滾動(dòng)引導(dǎo)濾波

滾動(dòng)引導(dǎo)濾波器具有尺度感知和邊緣恢復(fù)的優(yōu)良特性[15],能夠保留圖像邊緣,有效避免模糊現(xiàn)象。RGF首先對(duì)物體的尺度進(jìn)行了篩選,小尺度圖像信息主要指的是細(xì)節(jié)、紋理、小物體和一些噪聲。大尺度圖像信息主要指的是邊界、顏色過(guò)渡和變化不明顯的區(qū)域。小尺度的邊緣被平滑掉,僅僅保留了大尺度物體的邊緣。大多數(shù)邊緣保持濾波器沒(méi)有考慮到在濾波器設(shè)計(jì)中融入尺度,不能很好地實(shí)現(xiàn)小結(jié)構(gòu)與大結(jié)構(gòu)之間的分離。因此,RGF對(duì)于雙邊濾波和非線(xiàn)性擴(kuò)散濾波具有很大優(yōu)勢(shì),可以用來(lái)解決偽影噪聲問(wèn)題。

2 SART-GSR-RGF算法

受群稀疏正則化圖像重建[16]研究的啟發(fā),該文提出一種基于群稀疏性正則化與滾動(dòng)引導(dǎo)濾波算法相結(jié)合的SART-GSR-RGF重建方法。

圖3 相似群構(gòu)造示意圖

GSR算法對(duì)每個(gè)相似群fGk尋找合適的稀疏字典DG,在確定稀疏字典后,

(4)

其中,αG稱(chēng)為結(jié)構(gòu)GSR,在GSR算法中把αG作為正則項(xiàng),則求解GSR的無(wú)約束優(yōu)化模型為：

(5)

此正則化模型利用文獻(xiàn)[9]來(lái)求解。

最后,利用具有尺度感知和邊緣恢復(fù)的優(yōu)良特性的RGF進(jìn)行處理,先要進(jìn)行小結(jié)構(gòu)去除。

當(dāng)用小尺度的高斯濾波去除小結(jié)構(gòu)的邊緣時(shí),會(huì)使得圖像邊緣被嚴(yán)重模糊。相反,當(dāng)經(jīng)過(guò)相同方式用同樣尺度的高斯濾波去除大結(jié)構(gòu)的邊緣時(shí),只是達(dá)到了模糊的效果,仍能復(fù)原出為視覺(jué)所能接受的圖像。利用加權(quán)平均的形式來(lái)表示高斯濾波器算子,可得到：

(6)

其次進(jìn)行邊緣恢復(fù)。在雙邊濾波器基礎(chǔ)上將Ft作為引導(dǎo)圖像,迭代過(guò)程如式7所示：

(7)

算法1：SART-GSR-RGF

重復(fù)

SART重建：

forp=1,2,…,P

ifp=1

up=SART(u0,ω)

else

up=SART(up-1,ω)

end if

end for

ifup<0

up=0,u=up

GSR正則化：

fort=0,1,…,T

updateft+1,et+1=ft+1-bt

for eachfGk

end for

end for

u0=RGF(DG*αG)

滿(mǎn)足條件時(shí)停止

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

本節(jié)將對(duì)提出算法的重建性能進(jìn)行驗(yàn)證,測(cè)試平臺(tái)如下：64位Intel Core i7-10700處理器,主機(jī)頻率2.90 GHz,8G運(yùn)行內(nèi)存的Windows10操作系統(tǒng),數(shù)學(xué)軟件為Matlab R2021a。實(shí)驗(yàn)所用模型為Shepp-Logan。在比較不同算法的重建性能之前,首先利用所提算法進(jìn)行消融實(shí)驗(yàn),即對(duì)SART算法、SART-GSR算法SART-RGF算法以及所提算法進(jìn)行對(duì)比分析。同時(shí)驗(yàn)證了在不同稀疏角度下各算法的性能,觀察重建效果。圖4和圖5展示了在256個(gè)換能器背景下32角度和64角度時(shí)不同算法的重建效果;圖6和圖7展示了在512換能器背景下32角度和64角度時(shí)不同算法的重建效果,相對(duì)應(yīng)地也展示了同一放大部分的重建結(jié)果。

圖4 256換能器32角度下圖像重建效果比較

圖5 256換能器64角度下圖像重建效果比較

圖6 512換能器32角度下圖像重建效果比較

圖7 512換能器64角度下圖像重建效果比較

通過(guò)圖4至圖7可以看出,所提算法具有較好的重建效果,但僅通過(guò)視覺(jué)效果進(jìn)行評(píng)估算法的重建效果時(shí),評(píng)估結(jié)果具有一定的主觀性。為了更客觀地評(píng)價(jià)所提算法的有效性,使用具體的評(píng)價(jià)指標(biāo)進(jìn)行評(píng)價(jià)。選擇PSNR[17],MSE以及FSIM進(jìn)行客觀評(píng)價(jià),其中圖8、圖10和圖12展示了256換能器設(shè)置下32/64稀疏角度的20次迭代過(guò)程中的PSNR值、MSE值和FSIM值的比較結(jié)果及趨勢(shì);圖9、圖11和圖13展示了512換能器設(shè)置下32/64稀疏角度的20次迭代過(guò)程中的PSNR值、MSE值和FSIM值的比較結(jié)果及趨勢(shì)。

圖8 256換能器圖像重建效果PSNR比較

圖9 512換能器圖像重建效果PSNR比較

圖10 256換能器圖像重建效果MSE比較

圖11 512換能器圖像重建效果MSE比較

圖12 256換能器圖像重建效果FSIM比較

圖13 512換能器圖像重建效果FSIM比較

根據(jù)圖8至圖13可以很清晰地看到,所提算法SART-GSR-RGF具有最高的PSNR值和FSIM值,以及最低的MSE值。其中PSNR值越大,意味著圖像重建后的保真效果越好,圖像越清晰;FSIM值越高,意味著當(dāng)前所測(cè)試的圖像越接近參考圖像;MSE值越小,意味著與原始重建圖像相差越低,具有較好的重建性能。并且從圖中可以看到,所提算法在迭代開(kāi)始前期就已經(jīng)達(dá)到了比較好的性能,說(shuō)明所提算法具有很好的收斂性能,即算法收斂速度快。

為比較各類(lèi)方法的優(yōu)劣性,驗(yàn)證所提算法的有效性,表1給出了SART算法、TV-POCS算法、SART-NLM算法以及SART-GSR-RGF算法的PSNR,RMSE以及FSIM的數(shù)據(jù)。結(jié)果表明,所提算法在重建結(jié)果上均優(yōu)于其他4個(gè)對(duì)比算法。

表1 不同算法重建性能比較

4 結(jié)束語(yǔ)

該文提出了一種稀疏角度下基于滾動(dòng)引導(dǎo)濾波和群稀疏表示正則化的SART重建算法,稱(chēng)為SART-GSR-RGF。由于群稀疏正則化在去掉偽影的同時(shí),可能將邊緣或細(xì)節(jié)過(guò)度平滑,使得對(duì)比度降低,無(wú)法獲得符合人類(lèi)視覺(jué)效果的高分辨率圖像,因此利用滾動(dòng)引導(dǎo)濾波進(jìn)行一定的對(duì)比度提升,再次作為SART的輸入進(jìn)行迭代。經(jīng)過(guò)算法的分析以及實(shí)驗(yàn)的驗(yàn)證表明,所提算法具有良好的重建效果。接下來(lái),將對(duì)此算法進(jìn)行進(jìn)一步研究,以期獲得三維頭骨模型的較高質(zhì)量的重建。