莫礼东,徐伯庆
(上海理工大学 光电信息与计算机工程学院,上海 200093)
一种改进的TV图像重建算法
莫礼东,徐伯庆
(上海理工大学 光电信息与计算机工程学院,上海 200093)
针对当前稀疏角度下有限角图像重建过程中,边界部分出现伪影,降低了图像重建质量的缺陷。文中提出了一种新的ART+TV算法,该方法是在原始TV算法的基础上进行改进。原始TV梯度下降算法求解目标函数最小值时,使用固定函数作为目标函数,文中对其进行更改,采用带参数的目标函数,并对TV重建后的结果进行自适应步长修正,加速图像收敛。与传统的ART+TV算法相比,文中算法在不改变重建速度的基础上,且在少量迭代次数下,能重建出质量更高的图像,抑制图像伪影。
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成像已在医学和工业中有广泛应用[1]。但由于X射线的辐射剂量与X射线曝光次数成正比,而X射线的辐射对人体的身体健康有较大伤害。因此,减少X射线的辐射剂量成为了医学图像重建的研究方向。通过控制CT系统。对X射线投影角度进行稀疏采样和有限角度下的CT图像重建算法,是一种有效减少病人辐射剂量的技术手段,并被广泛应用。本文主要研究有限角度下的图像重建问题。
传统的图像重建算法主要包括滤波反投影(FBP)算法和代数迭代(ART)算法。其中FBP算法是从连续的成像模型推导出来的,并且间接认为在扫描角度范围内进行高密度的投影采样[2]。因此,该方法难以实现有限角度重建算法。ART算法是代数迭代重建,相比FBP算法,此算法更加实用于稀疏角度图像重建。文献[3~4]提出一种基于最小化和投影到凸集的迭代重建算法,且应用到稀疏角度CT重建,并取得了较好的重建结果。传统的ART+TV算法在扫描角度比较少的情况下,重建结果出现变形和伪影现象。为客服这一缺点,不少算法在原算法的基础上进行改进。文献[5~6]提出了一种基于先验图像约束的压缩感知重建算法。其主要思想是通过完整数据重建出来的图像作为先验条件,然后获取多个角度下的投影数据,将TV正则化思想用到待重建图像的基础上,还用到待重建图像与先验图像的差值图像。能够用少量的二维投影数据进行高质量的图像重建。但此算法对先验图像的信息具有一定的条件,给算法带来了局限性。本文在传统ART+TV算法的基础上,以ART算法为基础,加上全变差约束条件,对梯度图像的计算进行修改,采用可变的梯度计算函数,在保证重建图像速度的基础上,进行质量的优化。通过大量实验结果验证本文算法在原始算法的基础上确实有所改进。
1.1 ART算法原理
ART迭代算法最初是于1937年由Kaczmarz提出的[7-8],其旨在求解类似式(1)的方程组
(1)
在图像重建中,也是一个大稀疏矩阵的求解。因此,方程(1)的求解与图像重建求解方法类似。
假设一条强度为I0的单射线穿过一个内部分布未知的物体后,检测器上检测到的射线的强度为I,射线在物体中的衰减过程满足比尔定律[9]即
(2)
其中,f(x,y)表示物体在x-y平面内对单射线L的线性衰减系数。由式(2)可看出,f的值越大,衰减程度就越大。令穿过物体后射线投影数据为p,则有
(3)
因此p对应f(x,y)沿射线L方向的线积分。ART算法首先将待重建图像分解成n×n的小格,每一个小格对应一个像素值,如图1所示。
图1 ART重建示意图
在图1中,fi所表示的是第i个像素值,i=1,2,...,N,N=n2表示重建图像包含的总的像素数,因此f=(f1,f2,...,fN)表示一个N维空间的向量。在此我们定义扫描射线的宽度与像素网格的宽度一致。若在扫描时有P个方向,每个方向有R条射线,则总的射线条数 。图像像素在第j条射线上的积分记做Pj,j=1,2,...,M,所以图像向量f和投影向量p的关系可表示为
(4)
式(4)中,wij为加权因子,也即系统矩阵,反映了第j个像素对第i条射线投影值的贡献。每一条扫描射线的投影过程即为上述的一个方程,解得线性方程组的解也就得到了整个待重建截面的像素密度分布。
而由系数wij构成的矩阵是一个大型的稀疏矩阵,用常规的方法计算不仅计算过程复杂,且计算速度缓慢。在实际中,均用迭代方法进行计算,也即ART算法,其迭代公式为
(5)
式(5)中,k为迭代次数,wij为投影系数。由于噪声的影响,式(5)求解的结果就会存在较大的误差,于是引入松弛因子,来减少重建过程中的误差,其改进公式为
(6)
1.2 TV原理
在医疗图像重建中,由于X射线对人体的伤害以及人体器官组织的相似性。因此,稀疏角度图像重建被广泛运用。基于压缩感知的原理,可由少量的投影数据重建出原始图像,其被称为ART+TV算法。其重建图模型可表示为
min‖f‖TV,s.t.Wf=p
(7)
式(7)中,W是投影系数矩阵;p是投影值矩阵;f是灰度化的像素值。式(7)与式(1)是对应的,只是在式(1)的基础上加入了约束条件,其也是ART+TV与ART的主要区别。由式(7)可知,求解min‖f‖TV是关键。而并不能直接求出min‖f‖TV,一般均是通过用1范数L1来求解。
式(8)中▽fx,y=(Dxf,Dyf),Dxf和Dyf是离散微分算子,对于图像上的点(x,y),其梯度计算跟f(x,y)和f(x,y)附近的像素值有关
Dxfx,y=fx,y-fx-1,y;Dyfx,y=fx,y-fx,y-1
(8)
L1范数就是图像中每一个像素值的绝对值之和,其求导运算不方便。因此,大多采用L2范数来近似L1范数即
(9)
因此,对某个像素值f(x,y)求导,得
式中,ε=10-8其主要是防止‖f‖TV求导的时候出现分母为零的情况。通过上式的求解步骤,就容易求解方程(7)。
传统的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迭代次数;NTV是TV的迭代次数;pi是投影值;Wi是投影系统矩阵,本实验中NTV取20,Nview取10,K取20,a取0.2。
本文是在原始的ART+TV的算法上进行改进。原始的ART+TV算法,当用L1范数计算其最小值时,对于边界的方向,无法检测出来,由式(8)可知,L1范数求解最小值是求和的形式求解的,其中包含了方向不确定的边界的贡献,这样对重建图像会产生伪影。因此,本文提出了一种新的TV算法,对式(3)进行改进,采用式(10)进行替换
(10)
因此,此时对‖f‖TV求导公式变为
并对传统梯度下降的固定步长d进行修正,加速图像的收敛,其改动原理为
d2=‖fNview-fn+1‖
(11)
若d2>0.85×d,则d=0.85×d2。通过修改后的实验步骤为:(1)初始化图像及相关参数;(2)根据已经投影数据进行ART迭代,并进行非负修正,得到图像1;(3)利用本文提出的带参数求解的方法进行TV梯度下降算法求解目标函数最小值,得到重建图像2;(4)并对传统的TV算法进行动态步长d2修正,如式(11)所示,加速图像收敛,得图像3;(5)判断是否满足迭代终止条件,则退出循环;否则返回步骤(2),将图像3作为步骤(2)中ART迭代的初始图像,循环步骤(2)~(5)。
本实验采用分辨率为128×128的经典Shepp-Logan模型进行实验,分别在0~90;0~120;0~150范围,间隔1度进行仿真重建,对本文算法,传统ART+TV算法,ART算法进行对比,以说明改进后的算法的优势,其实验结果如图2所示。
图2 3组实验结果图
如图2所示,第1组图像实验是直接进行ART代数迭代算法的实验结果图;第2组图像实验是ART+TV即传统的TV算法实验;第3组是本文算法的实验结果图。
在实验过程中,为了客观地评价出重建图像的质量,本文采用峰值信噪比(PSNR)归一化均方距离判距d来定量评价重建图像与原始图像的误差。表1显示了3组实验的参数指标对比,其中参数PSNR值越大说明图像噪声越小,图像质量越好;归一化均方距离判距d的值越小代表重建误差越小。
表1 90个投影角度下重建指标对比
表2 120个投影角度下重建指标对比
表3 150个投影角度下重建指标对比
从表1~表3的实验结果可看出,本文算法在传统算法的基础上有所改进。
图3是在0°~90°范围内不同的迭代次数下的PSNR值对比图。
图3 不同算法的PSNR值的对比图
如图3可示,本文算法的PSNR值比其他两种算法的值要大。所以,本文算法对有限角度下图像重建算法的重建图像质量有所提高。
在稀疏角度图像重建中,本文针对有限角度图像
重建进行了研究。在有限角度重建过程中,受重建角度的限制,重建角度少的情况下边界出现伪影,为了避免这种情况,基于ART+TV算法提出了一种算法,其对梯度图像的求解过程进行修正,并且在梯度下降算法中加入动态修正因子,加速收敛。通过实验验证,本文算法在不改变重建的速度的基础上,比ART+TV算法在图像重建的质量上有所改进。
[1] 余维.不完备投影数据的CT重建算法研究[D].重庆:重庆大学,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].核电子学与探测技术,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的理论基础[M].北京:科学出版社,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