伍建輝,劉郁紀(jì),李公平,湯振興,潘小東,林俏露

(蘭州大學(xué)核科學(xué)與技術(shù)學(xué)院,蘭州 730000)

CT檢測(cè)技術(shù)是射線掃描之后,利用不同算法求解出被測(cè)物體斷層的衰減系數(shù)分布,然后把衰減系數(shù)分布轉(zhuǎn)化為相應(yīng)的密度分布以重現(xiàn)物體斷層,從而達(dá)到檢測(cè)物體斷層內(nèi)材質(zhì)、幾何結(jié)構(gòu)和缺陷等目的。CT重建算法大致可分為變換法和代數(shù)重建法,變換法重建速度快,可以用硬件實(shí)現(xiàn),但要求數(shù)據(jù)的采集必須完全且均勻,而實(shí)際應(yīng)用中,由于一些客觀原因的存在,投影數(shù)據(jù)的采集往往不完全或不均勻[1]。出于這種情況的考慮,ART算法應(yīng)用而生,它很好地彌補(bǔ)了這一缺點(diǎn)。盡管ART算法重建速度較慢,但隨著計(jì)算機(jī)技術(shù)的發(fā)展,速度逐漸成為次要問(wèn)題,圖像質(zhì)量上升為主要問(wèn)題。近年來(lái),眾多學(xué)者對(duì)不同情況下ART算法重建圖像質(zhì)量進(jìn)行了研究[2-5]。在此基礎(chǔ)上,筆者基于 Visual C++6.0及Matlab 7.0,仿真CT圖像重建過(guò)程,對(duì)初始值的選取、先驗(yàn)知識(shí)如何影響重建圖像質(zhì)量、不同噪聲下不同投影數(shù)對(duì)重建圖像質(zhì)量的影響等做了詳細(xì)的探討,得出了一些有益的結(jié)論,以期為實(shí)際CT檢測(cè)提供一些有價(jià)值的參考。

1 ART算法

1.1 射線方程的建立

被測(cè)物體斷層離散化為J=n×n個(gè)像素,對(duì)每一像素進(jìn)行編號(hào),xj表示第j個(gè)像素的像素值。定義rij=l/σ,其中l(wèi)表示第i條射線通過(guò)第j個(gè)像素的長(zhǎng)度;σ表示每一像素的長(zhǎng)度或?qū)挾?。從而像素j對(duì)射線i的貢獻(xiàn)為：

射線經(jīng)斷層后,其投影值為：

其中,rij稱(chēng)為加權(quán)因子。上式可寫(xiě)為：

用矩陣形式表示為：

式中R,X,P分別為：

如式(3)所示,M條射線可以建立M個(gè)方程。由數(shù)學(xué)知識(shí)可知,要求出J個(gè)像素值,只有建立J個(gè)線性無(wú)關(guān)的方程才有可能精確地求出每一個(gè)像素值。在實(shí)際中,由于像素個(gè)數(shù)龐大,例如：128×128的重建區(qū)域,就需要建立128×128個(gè)方程,此外,由于噪聲的存在及其它一些特殊情況,投影數(shù)的個(gè)數(shù)不一定能達(dá)到方程的個(gè)數(shù)。所以,直接求解方程在實(shí)際應(yīng)用中是不可行的。ART算法通過(guò)迭代的方法很好地解決了這一問(wèn)題。

1.2 ART算法分析

ART算法的基本思路為：首先,任選一初始值X0,以這一初始值為基礎(chǔ),進(jìn)行投影運(yùn)算(式(3)中每一方程代表一超平面。初始值看作J維空間內(nèi)一點(diǎn),然后向一超平面投影,在該超平面上得到相應(yīng)的投影點(diǎn))得X1,使其滿足式(3)中方程1;然后以X1為基礎(chǔ),進(jìn)行類(lèi)似于上一步的運(yùn)算得X2,使其滿足式(3)中方程2;如此繼續(xù),直至所有射線方程計(jì)算完,最后得結(jié)果XM。然后以XM為初始值,進(jìn)行下一輪運(yùn)算。依此類(lèi)推,直至得到滿足要求的值。

ART算法中進(jìn)行相應(yīng)運(yùn)算的公式為：

式中k為迭代序號(hào);λ為松弛參數(shù)(其中λ應(yīng)滿足條件0＜λ＜2)。一般地,只有當(dāng)一輪計(jì)算完畢后,才稱(chēng)為一次完整的迭代[2]。

J=2計(jì)算過(guò)程如圖1所示。圖中超平面H1,H2交于F(方程解)。X0為任選初始值,X0向超平面H1投影得X1。X1向超平面H2投影得X2。如此反復(fù),由圖可知隨著投影次數(shù)增加,其投影值逐漸趨近于像素值F。此外,X′0在超平面H1上的投影值X′1較X4更接近于像素值F。所以選取合適的初始值可以在迭代次數(shù)少的情況下達(dá)到較高的精度。

圖1 二維計(jì)算過(guò)程

實(shí)際CT圖像重建中,由于噪聲的存在,超平面H1,H2將會(huì)有所移動(dòng),需加入松弛參數(shù)λ對(duì)投影值進(jìn)行適當(dāng)?shù)恼{(diào)整。當(dāng)0＜λ＜1時(shí),其投影值X1位于H1超平面右上方,并隨 λ增大,投影值X1逐漸趨于超平面H1;當(dāng)λ=1時(shí),其投影值X1恰好落在超平面H1上;當(dāng)1＜λ＜2時(shí),投影值X1位于超平面H1左下方,并隨λ增大,投影值X1逐漸遠(yuǎn)離超平面H1。依此類(lèi)推,Xk在Xk-1基礎(chǔ)上進(jìn)行相應(yīng)的投影得到[2]?？梢?jiàn),松弛參數(shù)的選取對(duì)圖像重建起著很重要的作用。

綜上所述,ART算法圖像重建過(guò)程中,松弛參數(shù)、初始值的選取等多方面因素對(duì)圖像重建具有重要影響。

2 多因素對(duì)圖像重建的影響

仿真試驗(yàn)中,選取128×128經(jīng)典Shepp-Logan[6]頭部模型為重建對(duì)象,采用每一投影角度下射線束為128的平行束掃描方式,權(quán)因子計(jì)算采用Siddon改進(jìn)型[7]。為了精確評(píng)價(jià)重建圖像質(zhì)量,采用歸一化均方距離判據(jù)d、歸一化平均絕對(duì)距離判據(jù)r來(lái)定量地評(píng)價(jià)圖像質(zhì)量。其具體定義為：

式中ti,j,ri,j分別為測(cè)試模型和重建后圖像中第i行、j列的像素密度;為測(cè)試模型密度的平均值;圖像的像素為N×N個(gè)。其中d較敏感地反映某幾點(diǎn)產(chǎn)生較大誤差的情況,而r則較敏感地反映許多點(diǎn)均有一些小誤差的情況。判據(jù)d,r越小,表示重建圖像質(zhì)量越好[1]。

2.1 先驗(yàn)知識(shí)對(duì)重建圖像質(zhì)量及松弛參數(shù)選取的影響

ART算法圖像重建中,一般會(huì)采用先驗(yàn)知識(shí)來(lái)改善圖像質(zhì)量。如果已經(jīng)確定所有像素值都分布于某一范圍內(nèi),則在重建過(guò)程中,當(dāng)某像素值超出這一范圍時(shí),則強(qiáng)行把此像素值置為一特定值。例如,筆者所用重建模型的像素值都分布在0～1之間,當(dāng)某像素值f＜0時(shí),迭代過(guò)程中置此像素值f=0;當(dāng)某像素值f＞1時(shí),迭代過(guò)程中置此像素值f=1。圖2(a)和(b)分別表示在投影數(shù)(即投影角數(shù))為180,噪聲為0,初始值F(0)={0},加入先驗(yàn)知識(shí)與未加入先驗(yàn)知識(shí)情況下重建誤差分析圖,其中迭代次數(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 先驗(yàn)知識(shí)對(duì)重建誤差的影響

圖2中當(dāng)加入先驗(yàn)知識(shí)時(shí),其最佳松弛參數(shù)在0.7附近,同時(shí)松弛參數(shù)較大(如松弛參數(shù)＞0.4)的重建誤差隨著迭代次數(shù)的增加,其下降的速度很快。未加入先驗(yàn)知識(shí)時(shí),其最佳松弛參數(shù)在0.2附近,其較大松弛參數(shù)重建誤差下降緩慢。表1為松弛參數(shù)等于0.2,加入先驗(yàn)知識(shí)與未加入先驗(yàn)知識(shí)情況下,r和d值的分布情況。

表1 先驗(yàn)知識(shí)對(duì)歸一化絕對(duì)距離判據(jù)r及歸一化均方距離判據(jù)d的影響

由表1可看出,先驗(yàn)知識(shí)在松弛參數(shù)為0.2的情況下對(duì)判據(jù)d,r的影響很大。表1中,未加入先驗(yàn)知識(shí)迭代10次的r值較加入先驗(yàn)知識(shí)迭代3次的r值還大;同時(shí),未加入先驗(yàn)知識(shí)迭代10次的d值大于加入先驗(yàn)知識(shí)迭代5次的d值。

綜上所述,先驗(yàn)知識(shí)對(duì)大松弛參數(shù)重建誤差具有很好的抑制作用。當(dāng)加入先驗(yàn)知識(shí)時(shí),松弛參數(shù)可以適當(dāng)選大一點(diǎn),且經(jīng)過(guò)3～5次迭代便可以取得較好的重建結(jié)果。當(dāng)未加入先驗(yàn)知識(shí)時(shí),其松弛參數(shù)應(yīng)適當(dāng)選小一點(diǎn)。另外,相同松弛參數(shù)下,加入先驗(yàn)知識(shí),經(jīng)過(guò)3～5次迭代,就可以達(dá)到未加入先驗(yàn)知識(shí)迭代10次時(shí)的重建效果。所以在實(shí)際檢測(cè)中,應(yīng)盡量加入先驗(yàn)知識(shí)來(lái)提高重建圖像質(zhì)量。

2.2 初始值對(duì)圖像重建的影響

從圖1可看出,選取適當(dāng)?shù)某跏贾悼梢栽诘螖?shù)少的情況下重建出精度較高的圖像。圖3表示投影數(shù)為180,噪聲為0,松弛參數(shù)為0.4,加入先驗(yàn)知識(shí),初始值不同情況下重建誤差圖。其中所加入的誤差是以均值為0,標(biāo)準(zhǔn)差為某一值來(lái)定義(如誤差為0.01,表示均值為0,標(biāo)準(zhǔn)差為0.01);選取的隨機(jī)數(shù)在0～1之間服從均勻分布。

圖3 不同初始值重建誤差

圖3(a)中初始值誤差增大,其重建誤差也逐漸增大。同時(shí),初始值為0時(shí),迭代5～6次可以達(dá)到初始值誤差為0.07的重建誤差。而當(dāng)初始值為一隨機(jī)數(shù)時(shí),其重建誤差最大,大于初始值誤差為0.17的重建誤差。圖3(b)中,初始值取不同隨機(jī)數(shù)值時(shí),其重建圖像誤差大致相同。同時(shí),初始值取不同固定值(即使超出了重建區(qū)域像素值范圍)時(shí),其重建誤差也只在迭代次數(shù)少(如2～3次)時(shí)有細(xì)小差別,隨迭代次數(shù)增加,其差別逐漸消失。

綜上所述,選取誤差較小的初始值可以在迭代次數(shù)少的情況下重建出誤差較小的圖像,如當(dāng)初始值誤差＜0.07時(shí),迭代2～3次便可得出質(zhì)量較高的圖像。在完全不知真實(shí)值情況,但可以確定像素值分布范圍時(shí),初始值取一固定值為較佳選擇。

2.3 不同投影數(shù)在不同噪聲下對(duì)松弛參數(shù)選取及重建圖像的影響

式(3)中,如果線性無(wú)關(guān)方程個(gè)數(shù)越多,則M維空間中能確定精確解的范圍將會(huì)越小,由此可推知投影數(shù)對(duì)圖像重建具有重要的影響。此外,在實(shí)際CT檢測(cè)中,不可避免地存在噪聲,這對(duì)圖像重建的影響也是不可忽視的。在仿真投影數(shù)據(jù)中,分別加入均值為0,標(biāo)準(zhǔn)差為0.1,0.2,0.4的高斯白噪聲;采用投影數(shù)為180,128,90,36;為了考慮實(shí)際情況,選取松弛參數(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í),一些松弛參數(shù)重建誤差將發(fā)散,如圖6(a)和(b)中,松弛參數(shù)＞0.2的重建誤差不收斂;投影數(shù)較小時(shí),其松弛參數(shù)選取對(duì)噪聲的敏感度也較小,如投影數(shù)為36時(shí),在噪聲從0變?yōu)?.4的情況下,其較佳松弛參數(shù)的選取幾乎保持不變。同時(shí),從各圖中可看出,當(dāng)投影數(shù)減小及噪聲增加時(shí),都會(huì)增加重建圖像的誤差。

圖5 噪聲為0.2的重建誤差

圖6 噪聲為0.4的重建誤差

上述仿真過(guò)程為加先驗(yàn)知識(shí)時(shí)的情況。但未加先驗(yàn)知識(shí)時(shí),其噪聲、投影數(shù)對(duì)松弛參數(shù)選取及圖像重建的影響具有相同的趨勢(shì),只是其最佳松弛參數(shù)較加先驗(yàn)知識(shí)時(shí)要小得多,如投影數(shù)為180時(shí),其最佳松弛參數(shù)隨噪聲的變化,其值在0.2～0.06之間變化,限于篇幅,此不贅述。

2.4 缺角度對(duì)圖像重建的影響

實(shí)際CT檢測(cè)在一些特殊情況下,可能不能掃描某些角度,此時(shí),可以用ART算法來(lái)重建圖像。圖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°,都可以通過(guò)均勻增加投影數(shù)減小其重建誤差。從圖9可以看出,通過(guò)均勻增加投影數(shù),可以很好地抑制重建過(guò)程中產(chǎn)生的誤差。綜上所述,實(shí)際掃描過(guò)程中,當(dāng)缺失一定的度數(shù)時(shí),可以通過(guò)均勻增加投影數(shù)來(lái)減小其重建誤差。

圖 9 不同投影數(shù)基礎(chǔ)上缺20°的重建分析圖

3 結(jié)論

初始值、先驗(yàn)知識(shí)、投影數(shù)和噪聲等對(duì)圖像重建具有很大影響。在實(shí)際CT檢測(cè)中,如能預(yù)先獲得一些被測(cè)物像素值信息(如要檢測(cè)某已知產(chǎn)品一些細(xì)小裂縫或雜質(zhì),則可以把接近于被測(cè)物的像素值作為初始值)時(shí),則迭代2～3次便可以達(dá)到較高的分辨率,這樣就可以縮短檢測(cè)時(shí)間。此外,要盡可能利用被測(cè)物先驗(yàn)知識(shí),這樣可以在很大程度上改善重建圖像質(zhì)量,同時(shí)在加入先驗(yàn)知識(shí)后,應(yīng)適當(dāng)?shù)靥岣咚沙趨?shù)值。一般地,當(dāng)加入先驗(yàn)知識(shí)時(shí),在噪聲較大、投影數(shù)較多的情況下,應(yīng)選擇較小的松弛參數(shù)(如0.06～0.1之間);當(dāng)噪聲較小投影數(shù)較少時(shí),松弛參數(shù)可以適當(dāng)?shù)厝〈笠稽c(diǎn)(如0.4～0.7之間)。ART算法作為一種有效的迭代算法,可以在缺失角度情況下重建出較好的圖像,當(dāng)缺失一定角度時(shí),可以通過(guò)均勻增加投影數(shù)的方式來(lái)改善重建圖像質(zhì)量。另外,通過(guò)進(jìn)一步研究發(fā)現(xiàn),在一定投影數(shù)下,隨著缺失度數(shù)的增加,其最佳松弛參數(shù)逐漸增大。

[1]莊天戈.CT原理與算法[M].上海：上海交通大學(xué)出版社,1992.

[2]張順利,張定華,熬波,等.ART算法高質(zhì)量重建研究[J].核電子學(xué)與探測(cè)技術(shù),2007,27(4)：719-720.

[3]張順利,張定華,王成,等.投影數(shù)對(duì)ART算法重建質(zhì)量的影響[J].無(wú)損檢測(cè),2008,30(12)：889-891.

[4]冷帥,張麗,陳志強(qiáng),等.迭代算法用于投影數(shù)據(jù)不完備時(shí)重建斷層圖像的研究[J].核電子學(xué)與探測(cè)技術(shù),2003,23(3)：216-219.

[5]高欣,夏順仁,汪元美,等.一種不完全投影圖像重建的快速迭代算法[J].浙江大學(xué)學(xué)報(bào)(工學(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].計(jì)算機(jī)工程與應(yīng)用,2006,42(24)：2-3.