伍建輝,劉郁紀,李公平,湯振興,潘小東,林俏露
(蘭州大學(xué)核科學(xué)與技術(shù)學(xué)院,蘭州 730000)
CT檢測技術(shù)是射線掃描之后,利用不同算法求解出被測物體斷層的衰減系數(shù)分布,然后把衰減系數(shù)分布轉(zhuǎn)化為相應(yīng)的密度分布以重現(xiàn)物體斷層,從而達到檢測物體斷層內(nèi)材質(zhì)、幾何結(jié)構(gòu)和缺陷等目的。CT重建算法大致可分為變換法和代數(shù)重建法,變換法重建速度快,可以用硬件實現(xiàn),但要求數(shù)據(jù)的采集必須完全且均勻,而實際應(yīng)用中,由于一些客觀原因的存在,投影數(shù)據(jù)的采集往往不完全或不均勻[1]。出于這種情況的考慮,ART算法應(yīng)用而生,它很好地彌補了這一缺點。盡管ART算法重建速度較慢,但隨著計算機技術(shù)的發(fā)展,速度逐漸成為次要問題,圖像質(zhì)量上升為主要問題。近年來,眾多學(xué)者對不同情況下ART算法重建圖像質(zhì)量進行了研究[2-5]。在此基礎(chǔ)上,筆者基于 Visual C++6.0及Matlab 7.0,仿真CT圖像重建過程,對初始值的選取、先驗知識如何影響重建圖像質(zhì)量、不同噪聲下不同投影數(shù)對重建圖像質(zhì)量的影響等做了詳細的探討,得出了一些有益的結(jié)論,以期為實際CT檢測提供一些有價值的參考。
被測物體斷層離散化為J=n×n個像素,對每一像素進行編號,xj表示第j個像素的像素值。定義rij=l/σ,其中l(wèi)表示第i條射線通過第j個像素的長度;σ表示每一像素的長度或?qū)挾取亩袼豭對射線i的貢獻為:
射線經(jīng)斷層后,其投影值為:
其中,rij稱為加權(quán)因子。上式可寫為:
用矩陣形式表示為:
式中R,X,P分別為:
如式(3)所示,M條射線可以建立M個方程。由數(shù)學(xué)知識可知,要求出J個像素值,只有建立J個線性無關(guān)的方程才有可能精確地求出每一個像素值。在實際中,由于像素個數(shù)龐大,例如:128×128的重建區(qū)域,就需要建立128×128個方程,此外,由于噪聲的存在及其它一些特殊情況,投影數(shù)的個數(shù)不一定能達到方程的個數(shù)。所以,直接求解方程在實際應(yīng)用中是不可行的。ART算法通過迭代的方法很好地解決了這一問題。
ART算法的基本思路為:首先,任選一初始值X0,以這一初始值為基礎(chǔ),進行投影運算(式(3)中每一方程代表一超平面。初始值看作J維空間內(nèi)一點,然后向一超平面投影,在該超平面上得到相應(yīng)的投影點)得X1,使其滿足式(3)中方程1;然后以X1為基礎(chǔ),進行類似于上一步的運算得X2,使其滿足式(3)中方程2;如此繼續(xù),直至所有射線方程計算完,最后得結(jié)果XM。然后以XM為初始值,進行下一輪運算。依此類推,直至得到滿足要求的值。
ART算法中進行相應(yīng)運算的公式為:
式中k為迭代序號;λ為松弛參數(shù)(其中λ應(yīng)滿足條件0<λ<2)。一般地,只有當一輪計算完畢后,才稱為一次完整的迭代[2]。
J=2計算過程如圖1所示。圖中超平面H1,H2交于F(方程解)。X0為任選初始值,X0向超平面H1投影得X1。X1向超平面H2投影得X2。如此反復(fù),由圖可知隨著投影次數(shù)增加,其投影值逐漸趨近于像素值F。此外,X′0在超平面H1上的投影值X′1較X4更接近于像素值F。所以選取合適的初始值可以在迭代次數(shù)少的情況下達到較高的精度。
圖1 二維計算過程
實際CT圖像重建中,由于噪聲的存在,超平面H1,H2將會有所移動,需加入松弛參數(shù)λ對投影值進行適當?shù)恼{(diào)整。當0<λ<1時,其投影值X1位于H1超平面右上方,并隨 λ增大,投影值X1逐漸趨于超平面H1;當λ=1時,其投影值X1恰好落在超平面H1上;當1<λ<2時,投影值X1位于超平面H1左下方,并隨λ增大,投影值X1逐漸遠離超平面H1。依此類推,Xk在Xk-1基礎(chǔ)上進行相應(yīng)的投影得到[2]??梢?松弛參數(shù)的選取對圖像重建起著很重要的作用。
綜上所述,ART算法圖像重建過程中,松弛參數(shù)、初始值的選取等多方面因素對圖像重建具有重要影響。
仿真試驗中,選取128×128經(jīng)典Shepp-Logan[6]頭部模型為重建對象,采用每一投影角度下射線束為128的平行束掃描方式,權(quán)因子計算采用Siddon改進型[7]。為了精確評價重建圖像質(zhì)量,采用歸一化均方距離判據(jù)d、歸一化平均絕對距離判據(jù)r來定量地評價圖像質(zhì)量。其具體定義為:
式中ti,j,ri,j分別為測試模型和重建后圖像中第i行、j列的像素密度;為測試模型密度的平均值;圖像的像素為N×N個。其中d較敏感地反映某幾點產(chǎn)生較大誤差的情況,而r則較敏感地反映許多點均有一些小誤差的情況。判據(jù)d,r越小,表示重建圖像質(zhì)量越好[1]。
ART算法圖像重建中,一般會采用先驗知識來改善圖像質(zhì)量。如果已經(jīng)確定所有像素值都分布于某一范圍內(nèi),則在重建過程中,當某像素值超出這一范圍時,則強行把此像素值置為一特定值。例如,筆者所用重建模型的像素值都分布在0~1之間,當某像素值f<0時,迭代過程中置此像素值f=0;當某像素值f>1時,迭代過程中置此像素值f=1。圖2(a)和(b)分別表示在投影數(shù)(即投影角數(shù))為180,噪聲為0,初始值F(0)={0},加入先驗知識與未加入先驗知識情況下重建誤差分析圖,其中迭代次數(shù)為1~10,松弛參數(shù)分別選取λ=0.03,0.06,0.1,0.2,0.4,0.7,1.0,1.2,1.5和1.9。
圖2 先驗知識對重建誤差的影響
圖2中當加入先驗知識時,其最佳松弛參數(shù)在0.7附近,同時松弛參數(shù)較大(如松弛參數(shù)>0.4)的重建誤差隨著迭代次數(shù)的增加,其下降的速度很快。未加入先驗知識時,其最佳松弛參數(shù)在0.2附近,其較大松弛參數(shù)重建誤差下降緩慢。表1為松弛參數(shù)等于0.2,加入先驗知識與未加入先驗知識情況下,r和d值的分布情況。
表1 先驗知識對歸一化絕對距離判據(jù)r及歸一化均方距離判據(jù)d的影響
由表1可看出,先驗知識在松弛參數(shù)為0.2的情況下對判據(jù)d,r的影響很大。表1中,未加入先驗知識迭代10次的r值較加入先驗知識迭代3次的r值還大;同時,未加入先驗知識迭代10次的d值大于加入先驗知識迭代5次的d值。
綜上所述,先驗知識對大松弛參數(shù)重建誤差具有很好的抑制作用。當加入先驗知識時,松弛參數(shù)可以適當選大一點,且經(jīng)過3~5次迭代便可以取得較好的重建結(jié)果。當未加入先驗知識時,其松弛參數(shù)應(yīng)適當選小一點。另外,相同松弛參數(shù)下,加入先驗知識,經(jīng)過3~5次迭代,就可以達到未加入先驗知識迭代10次時的重建效果。所以在實際檢測中,應(yīng)盡量加入先驗知識來提高重建圖像質(zhì)量。
從圖1可看出,選取適當?shù)某跏贾悼梢栽诘螖?shù)少的情況下重建出精度較高的圖像。圖3表示投影數(shù)為180,噪聲為0,松弛參數(shù)為0.4,加入先驗知識,初始值不同情況下重建誤差圖。其中所加入的誤差是以均值為0,標準差為某一值來定義(如誤差為0.01,表示均值為0,標準差為0.01);選取的隨機數(shù)在0~1之間服從均勻分布。
圖3 不同初始值重建誤差
圖3(a)中初始值誤差增大,其重建誤差也逐漸增大。同時,初始值為0時,迭代5~6次可以達到初始值誤差為0.07的重建誤差。而當初始值為一隨機數(shù)時,其重建誤差最大,大于初始值誤差為0.17的重建誤差。圖3(b)中,初始值取不同隨機數(shù)值時,其重建圖像誤差大致相同。同時,初始值取不同固定值(即使超出了重建區(qū)域像素值范圍)時,其重建誤差也只在迭代次數(shù)少(如2~3次)時有細小差別,隨迭代次數(shù)增加,其差別逐漸消失。
綜上所述,選取誤差較小的初始值可以在迭代次數(shù)少的情況下重建出誤差較小的圖像,如當初始值誤差<0.07時,迭代2~3次便可得出質(zhì)量較高的圖像。在完全不知真實值情況,但可以確定像素值分布范圍時,初始值取一固定值為較佳選擇。
式(3)中,如果線性無關(guān)方程個數(shù)越多,則M維空間中能確定精確解的范圍將會越小,由此可推知投影數(shù)對圖像重建具有重要的影響。此外,在實際CT檢測中,不可避免地存在噪聲,這對圖像重建的影響也是不可忽視的。在仿真投影數(shù)據(jù)中,分別加入均值為0,標準差為0.1,0.2,0.4的高斯白噪聲;采用投影數(shù)為180,128,90,36;為了考慮實際情況,選取松弛參數(shù)為0.01,0.03,0.06,0.1,0.2,0.4,0.8和1.0。重建誤差如圖4,5和6所示。
比較圖4,5和6可看出,相同噪聲下,投影數(shù)增大,其最佳松弛參數(shù)逐漸減小。相同投影數(shù)下,噪聲增大,其最佳松弛參數(shù)逐漸減小。此外,投影數(shù)、噪聲都較大時,一些松弛參數(shù)重建誤差將發(fā)散,如圖6(a)和(b)中,松弛參數(shù)>0.2的重建誤差不收斂;投影數(shù)較小時,其松弛參數(shù)選取對噪聲的敏感度也較小,如投影數(shù)為36時,在噪聲從0變?yōu)?.4的情況下,其較佳松弛參數(shù)的選取幾乎保持不變。同時,從各圖中可看出,當投影數(shù)減小及噪聲增加時,都會增加重建圖像的誤差。
圖5 噪聲為0.2的重建誤差
圖6 噪聲為0.4的重建誤差
上述仿真過程為加先驗知識時的情況。但未加先驗知識時,其噪聲、投影數(shù)對松弛參數(shù)選取及圖像重建的影響具有相同的趨勢,只是其最佳松弛參數(shù)較加先驗知識時要小得多,如投影數(shù)為180時,其最佳松弛參數(shù)隨噪聲的變化,其值在0.2~0.06之間變化,限于篇幅,此不贅述。
實際CT檢測在一些特殊情況下,可能不能掃描某些角度,此時,可以用ART算法來重建圖像。圖7描述了松弛參數(shù)為0.4、迭代次數(shù)為8、在均勻投影數(shù)為180的基礎(chǔ)上缺失4,8,12,…,60°的重建誤差。圖8為松弛參數(shù)等于0.4、迭代 8次、在均勻投影數(shù)為80,100,120,140,160,185,210,230,250,270和290基礎(chǔ)上缺失20°與40°的誤差分析圖。圖9為均勻投影數(shù)為80,185和290基礎(chǔ)上缺20°的重建圖及第80行重建像素值分析圖(注:某均勻投影數(shù)基礎(chǔ)上缺角度是指某投影數(shù)在180°范圍內(nèi)均勻分布后缺失角度數(shù))。
由圖7可看出,其重建誤差與所缺度數(shù)幾乎成線性關(guān)系。在圖 8中,不同投影數(shù)下缺 20°或 40°,都可以通過均勻增加投影數(shù)減小其重建誤差。從圖9可以看出,通過均勻增加投影數(shù),可以很好地抑制重建過程中產(chǎn)生的誤差。綜上所述,實際掃描過程中,當缺失一定的度數(shù)時,可以通過均勻增加投影數(shù)來減小其重建誤差。
圖 9 不同投影數(shù)基礎(chǔ)上缺20°的重建分析圖
初始值、先驗知識、投影數(shù)和噪聲等對圖像重建具有很大影響。在實際CT檢測中,如能預(yù)先獲得一些被測物像素值信息(如要檢測某已知產(chǎn)品一些細小裂縫或雜質(zhì),則可以把接近于被測物的像素值作為初始值)時,則迭代2~3次便可以達到較高的分辨率,這樣就可以縮短檢測時間。此外,要盡可能利用被測物先驗知識,這樣可以在很大程度上改善重建圖像質(zhì)量,同時在加入先驗知識后,應(yīng)適當?shù)靥岣咚沙趨?shù)值。一般地,當加入先驗知識時,在噪聲較大、投影數(shù)較多的情況下,應(yīng)選擇較小的松弛參數(shù)(如0.06~0.1之間);當噪聲較小投影數(shù)較少時,松弛參數(shù)可以適當?shù)厝〈笠稽c(如0.4~0.7之間)。ART算法作為一種有效的迭代算法,可以在缺失角度情況下重建出較好的圖像,當缺失一定角度時,可以通過均勻增加投影數(shù)的方式來改善重建圖像質(zhì)量。另外,通過進一步研究發(fā)現(xiàn),在一定投影數(shù)下,隨著缺失度數(shù)的增加,其最佳松弛參數(shù)逐漸增大。
[1]莊天戈.CT原理與算法[M].上海:上海交通大學(xué)出版社,1992.
[2]張順利,張定華,熬波,等.ART算法高質(zhì)量重建研究[J].核電子學(xué)與探測技術(shù),2007,27(4):719-720.
[3]張順利,張定華,王成,等.投影數(shù)對ART算法重建質(zhì)量的影響[J].無損檢測,2008,30(12):889-891.
[4]冷帥,張麗,陳志強,等.迭代算法用于投影數(shù)據(jù)不完備時重建斷層圖像的研究[J].核電子學(xué)與探測技術(shù),2003,23(3):216-219.
[5]高欣,夏順仁,汪元美,等.一種不完全投影圖像重建的快速迭代算法[J].浙江大學(xué)學(xué)報(工學(xué)版),2004,38(9):1108-1111.
[6]Shepp L,Logan B P.The Fourier reconstruction of a head section[J].IEEE Trans Nucl SCI,1974(NS-21):21-43.
[7]張順利,張定華,李山,等.ART算法快速圖像重建研究[J].計算機工程與應(yīng)用,2006,42(24):2-3.