一種改進的TV圖像重建算法

莫禮東，徐伯慶

(上海理工大學(xué) 光電信息與計算機工程學(xué)院，上海 200093)

?

一種改進的TV圖像重建算法

莫禮東，徐伯慶

(上海理工大學(xué) 光電信息與計算機工程學(xué)院，上海 200093)

針對當(dāng)前稀疏角度下有限角圖像重建過程中，邊界部分出現(xiàn)偽影，降低了圖像重建質(zhì)量的缺陷。文中提出了一種新的ART+TV算法，該方法是在原始TV算法的基礎(chǔ)上進行改進。原始TV梯度下降算法求解目標(biāo)函數(shù)最小值時，使用固定函數(shù)作為目標(biāo)函數(shù)，文中對其進行更改，采用帶參數(shù)的目標(biāo)函數(shù)，并對TV重建后的結(jié)果進行自適應(yīng)步長修正，加速圖像收斂。與傳統(tǒng)的ART+TV算法相比，文中算法在不改變重建速度的基礎(chǔ)上，且在少量迭代次數(shù)下，能重建出質(zhì)量更高的圖像，抑制圖像偽影。

ART+TV算法；迭代重建；稀疏角度CT

MO Lidong, XU Boqing

(School of Optical-Electrical and Computer Engineering, University of Shanghai for Science and Technology,Shanghai 200093, China)

X射線CT成像已在醫(yī)學(xué)和工業(yè)中有廣泛應(yīng)用[1]。但由于X射線的輻射劑量與X射線曝光次數(shù)成正比，而X射線的輻射對人體的身體健康有較大傷害。因此，減少X射線的輻射劑量成為了醫(yī)學(xué)圖像重建的研究方向。通過控制CT系統(tǒng)。對X射線投影角度進行稀疏采樣和有限角度下的CT圖像重建算法，是一種有效減少病人輻射劑量的技術(shù)手段，并被廣泛應(yīng)用。本文主要研究有限角度下的圖像重建問題。

傳統(tǒng)的圖像重建算法主要包括濾波反投影(FBP)算法和代數(shù)迭代(ART)算法。其中FBP算法是從連續(xù)的成像模型推導(dǎo)出來的，并且間接認為在掃描角度范圍內(nèi)進行高密度的投影采樣[2]。因此，該方法難以實現(xiàn)有限角度重建算法。ART算法是代數(shù)迭代重建，相比FBP算法，此算法更加實用于稀疏角度圖像重建。文獻[3～4]提出一種基于最小化和投影到凸集的迭代重建算法，且應(yīng)用到稀疏角度CT重建，并取得了較好的重建結(jié)果。傳統(tǒng)的ART+TV算法在掃描角度比較少的情況下，重建結(jié)果出現(xiàn)變形和偽影現(xiàn)象。為客服這一缺點，不少算法在原算法的基礎(chǔ)上進行改進。文獻[5～6]提出了一種基于先驗圖像約束的壓縮感知重建算法。其主要思想是通過完整數(shù)據(jù)重建出來的圖像作為先驗條件，然后獲取多個角度下的投影數(shù)據(jù)，將TV正則化思想用到待重建圖像的基礎(chǔ)上，還用到待重建圖像與先驗圖像的差值圖像。能夠用少量的二維投影數(shù)據(jù)進行高質(zhì)量的圖像重建。但此算法對先驗圖像的信息具有一定的條件，給算法帶來了局限性。本文在傳統(tǒng)ART+TV算法的基礎(chǔ)上，以ART算法為基礎(chǔ)，加上全變差約束條件，對梯度圖像的計算進行修改，采用可變的梯度計算函數(shù)，在保證重建圖像速度的基礎(chǔ)上，進行質(zhì)量的優(yōu)化。通過大量實驗結(jié)果驗證本文算法在原始算法的基礎(chǔ)上確實有所改進。

1 ART+TV算法原理

1.1 ART算法原理

ART迭代算法最初是于1937年由Kaczmarz提出的[7-8]，其旨在求解類似式(1)的方程組

(1)

在圖像重建中，也是一個大稀疏矩陣的求解。因此，方程(1)的求解與圖像重建求解方法類似。

假設(shè)一條強度為I0的單射線穿過一個內(nèi)部分布未知的物體后，檢測器上檢測到的射線的強度為I，射線在物體中的衰減過程滿足比爾定律[9]即

(2)

其中，f(x,y)表示物體在x-y平面內(nèi)對單射線L的線性衰減系數(shù)。由式(2)可看出,f的值越大，衰減程度就越大。令穿過物體后射線投影數(shù)據(jù)為p，則有

(3)

因此p對應(yīng)f(x,y)沿射線L方向的線積分。ART算法首先將待重建圖像分解成n×n的小格，每一個小格對應(yīng)一個像素值，如圖1所示。

圖1 ART重建示意圖

在圖1中,fi所表示的是第i個像素值，i=1,2,...,N,N=n2表示重建圖像包含的總的像素數(shù)，因此f=(f1,f2,...,fN)表示一個N維空間的向量。在此我們定義掃描射線的寬度與像素網(wǎng)格的寬度一致。若在掃描時有P個方向，每個方向有R條射線，則總的射線條數(shù) 。圖像像素在第j條射線上的積分記做Pj,j=1,2,...,M,所以圖像向量f和投影向量p的關(guān)系可表示為

(4)

式(4)中，wij為加權(quán)因子，也即系統(tǒng)矩陣，反映了第j個像素對第i條射線投影值的貢獻。每一條掃描射線的投影過程即為上述的一個方程，解得線性方程組的解也就得到了整個待重建截面的像素密度分布。

而由系數(shù)wij構(gòu)成的矩陣是一個大型的稀疏矩陣，用常規(guī)的方法計算不僅計算過程復(fù)雜，且計算速度緩慢。在實際中，均用迭代方法進行計算，也即ART算法，其迭代公式為

(5)

式(5)中，k為迭代次數(shù)，wij為投影系數(shù)。由于噪聲的影響，式(5)求解的結(jié)果就會存在較大的誤差，于是引入松弛因子，來減少重建過程中的誤差，其改進公式為

(6)

1.2 TV原理

在醫(yī)療圖像重建中，由于X射線對人體的傷害以及人體器官組織的相似性。因此，稀疏角度圖像重建被廣泛運用?；趬嚎s感知的原理，可由少量的投影數(shù)據(jù)重建出原始圖像，其被稱為ART+TV算法。其重建圖模型可表示為

min‖f‖TV,s.t.Wf=p

(7)

式(7)中,W是投影系數(shù)矩陣;p是投影值矩陣;f是灰度化的像素值。式(7)與式(1)是對應(yīng)的，只是在式(1)的基礎(chǔ)上加入了約束條件，其也是ART+TV與ART的主要區(qū)別。由式(7)可知，求解min‖f‖TV是關(guān)鍵。而并不能直接求出min‖f‖TV，一般均是通過用1范數(shù)L1來求解。

式(8)中▽fx,y=(Dxf,Dyf)，Dxf和Dyf是離散微分算子，對于圖像上的點(x,y)，其梯度計算跟f(x,y)和f(x,y)附近的像素值有關(guān)

Dxfx,y=fx,y-fx-1,y;Dyfx,y=fx,y-fx,y-1

(8)

L1范數(shù)就是圖像中每一個像素值的絕對值之和，其求導(dǎo)運算不方便。因此，大多采用L2范數(shù)來近似L1范數(shù)即

(9)

因此，對某個像素值f(x,y)求導(dǎo)，得

式中，ε=10-8其主要是防止‖f‖TV求導(dǎo)的時候出現(xiàn)分母為零的情況。通過上式的求解步驟，就容易求解方程(7)。

傳統(tǒng)的ART+TV算法步驟如下[10]

Givenpi,Nview,Wi(i=1,2,3,…Nview)andNTVInitialization:f0=0.

fork=1,2,3,…do

ART Updating:

f0=fk

fori=1,2,3,…Nviewdo

end

TV Minimization:

forn=1,2,3,…,NTV-1 do

d=‖fn-f0‖2,

fn+1=fn+a·d·v

end

其中，Nview是ART迭代次數(shù);NTV是TV的迭代次數(shù);pi是投影值;Wi是投影系統(tǒng)矩陣，本實驗中NTV取20，Nview取10，K取20，a取0.2。

2 本文算法

本文是在原始的ART+TV的算法上進行改進。原始的ART+TV算法，當(dāng)用L1范數(shù)計算其最小值時，對于邊界的方向，無法檢測出來，由式(8)可知，L1范數(shù)求解最小值是求和的形式求解的，其中包含了方向不確定的邊界的貢獻，這樣對重建圖像會產(chǎn)生偽影。因此，本文提出了一種新的TV算法，對式(3)進行改進，采用式(10)進行替換

(10)

因此，此時對‖f‖TV求導(dǎo)公式變?yōu)?/p>

并對傳統(tǒng)梯度下降的固定步長d進行修正，加速圖像的收斂，其改動原理為

d2=‖fNview-fn+1‖

(11)

若d2>0.85×d,則d=0.85×d2。通過修改后的實驗步驟為：(1)初始化圖像及相關(guān)參數(shù)；(2)根據(jù)已經(jīng)投影數(shù)據(jù)進行ART迭代，并進行非負修正，得到圖像1；(3)利用本文提出的帶參數(shù)求解的方法進行TV梯度下降算法求解目標(biāo)函數(shù)最小值，得到重建圖像2；(4)并對傳統(tǒng)的TV算法進行動態(tài)步長d2修正，如式(11)所示，加速圖像收斂，得圖像3；(5)判斷是否滿足迭代終止條件，則退出循環(huán)；否則返回步驟(2)，將圖像3作為步驟(2)中ART迭代的初始圖像，循環(huán)步驟(2)～(5)。

3 實驗結(jié)果

本實驗采用分辨率為128×128的經(jīng)典Shepp-Logan模型進行實驗，分別在0～90；0～120；0～150范圍，間隔1度進行仿真重建，對本文算法，傳統(tǒng)ART+TV算法，ART算法進行對比，以說明改進后的算法的優(yōu)勢，其實驗結(jié)果如圖2所示。

圖2 3組實驗結(jié)果圖

如圖2所示，第1組圖像實驗是直接進行ART代數(shù)迭代算法的實驗結(jié)果圖；第2組圖像實驗是ART+TV即傳統(tǒng)的TV算法實驗；第3組是本文算法的實驗結(jié)果圖。

在實驗過程中，為了客觀地評價出重建圖像的質(zhì)量，本文采用峰值信噪比(PSNR)歸一化均方距離判距d來定量評價重建圖像與原始圖像的誤差。表1顯示了3組實驗的參數(shù)指標(biāo)對比，其中參數(shù)PSNR值越大說明圖像噪聲越小，圖像質(zhì)量越好；歸一化均方距離判距d的值越小代表重建誤差越小。

表1 90個投影角度下重建指標(biāo)對比

表2 120個投影角度下重建指標(biāo)對比

表3 150個投影角度下重建指標(biāo)對比

從表1～表3的實驗結(jié)果可看出，本文算法在傳統(tǒng)算法的基礎(chǔ)上有所改進。

圖3是在0°～90°范圍內(nèi)不同的迭代次數(shù)下的PSNR值對比圖。

圖3 不同算法的PSNR值的對比圖

如圖3可示，本文算法的PSNR值比其他兩種算法的值要大。所以，本文算法對有限角度下圖像重建算法的重建圖像質(zhì)量有所提高。

4 結(jié)束語

在稀疏角度圖像重建中，本文針對有限角度圖像

重建進行了研究。在有限角度重建過程中，受重建角度的限制，重建角度少的情況下邊界出現(xiàn)偽影，為了避免這種情況，基于ART+TV算法提出了一種算法，其對梯度圖像的求解過程進行修正，并且在梯度下降算法中加入動態(tài)修正因子，加速收斂。通過實驗驗證，本文算法在不改變重建的速度的基礎(chǔ)上，比ART+TV算法在圖像重建的質(zhì)量上有所改進。

[1] 余維.不完備投影數(shù)據(jù)的CT重建算法研究[D].重慶:重慶大學(xué),2014.

[2] Pan X,Sidky E Y,Vannier M.Why do commercial CT scanners still employ traditional, filtered back-projection for image reconstruction[J].Inverse Problems,2009,25(12):1541-1548.

[3] Gopi V P,Fayiz T K,Palanisamy P. Regularization based CT image reconstruction using algebraic techniques[C].Washington,USA:International Conference on Electronics and Communication Systems, IEEE,2014.

[4] Sidky E Y,Kao C M,Pan X.Accurate image reconstruction from few-views and limited-angle data in divergent-beam CT[J].Journal of X-ray science and technology,2009,14(2):119-139.

[5] Chen G J,Leng S.Prior image constrained compressed sensing (PICCS):a method to accurately reconstruct dynamic CT images from highly undersampled projection data sets[J].Medical Physics,2008,35(2):660-663.

[6] 郭靜鈺,齊宏亮,袁媛,等.先驗圖像約束的有限角度CT圖像重建算法[J].核電子學(xué)與探測技術(shù),2014(12):1421-1424.

[7] Kaczmarz S.Angenaherte auflosung von systemen lincarer glerichungen[J].Bull Acad P0l Sciencever Letter A,1937(6-8A):355-357.

[8] Tanabe K.Projection method for solving a singular system of linear equations and its applications[J].Numer Math,1971,17(3):203-214.

[9] G T赫爾曼.由投影重建圖像-CT的理論基礎(chǔ)[M].北京:科學(xué)出版社,1985.

[10] Jin X,Li L,Chen Z,et al.Anisotropic total variation for limited-angle CT reconstruction[C].Tianjin:IEEE Nuclear Science Symposium Conference Record, Nuclear Science Symposium,2010.

An Improved TV Image Reconstruction Algorithm

For the current limited angle image reconstruction process sparse angle, the boundary part artifacts appear to reduce the image quality of the reconstructed defects. This paper proposes a new algorithm of ART + TV, which is the basis of the original TV algorithms, to improve it. In solving the minimum objective function by the original TV gradient descent algorithm, the fixed function is used as the objective function. While in this paper, the text to change it, the objective function with parameters is employed, and the results of TV reconstruction are adaptive step corrected to accelerate the graphics convergence. With compare the proposed algorithm offers better quality of reconstructed images with a few iterations than the conventional ART + TV algorithms at the same reconstruction speed and suppresses image artifacts

ART+TV algorithm; iterative reconstruction; sparse angle CT

2016- 01- 14

莫禮東(1991-)，男，碩士研究生。研究方向：信號與信息處理。徐伯慶(1958-)，男，博士，副教授。研究方向：通信及圖像處理。

10.16180/j.cnki.issn1007-7820.2016.10.014

TP391.41

A

1007-7820(2016)10-047-04