当前位置：文档库 › 基于压缩传感和代数重建法的CT图像重建

基于压缩传感和代数重建法的CT图像重建

第35卷第3期2009年5月

光学技术OPTICAL TECHN IQU E Vol.35No.3

May 2009

文章编号:100221582(2009)0320422204基于压缩传感和代数重建法的CT 图像重建

练秋生,郝鹏鹏

(燕山大学信息科学与工程学院,河北秦皇岛　066004)

摘　要:代数重建法(ART )是一种重要的CT 图像重建方法,适合于不完全投影数据的图像重建,其缺点是重建速度慢。为提高图像重建的质量和速度,利用压缩传感理论提出了一种基于ART 的高质量图像重建算法。该算法将CT 图像的梯度稀疏性结合到ART 图像重建中,在每次迭代中的投影操作结束后用梯度下降法调整全变差,减小图像梯度的l 1范数。实验结果验证了该算法的有效性。

关

键

词:图像重建;代数重建法;压缩传感;全变差

中图分类号:TP391 文献标识码:A

Image reconstruction for CT based on compressed sensing and ART

LI AN Qiu 2sheng ,H AO Peng 2p eng

(Institute of Information Science and Technology ,Y anshan University ,Hebei Qinhuangdao 066004,China )Abstract :Algebraic reconstruction technique (ART )is an important method for the image reconstruction of computed to 2mography (CT ),which can recover the image from incomplete projection data.However ,ART suffers slow reconstruction speed.An algorithm based on ART is proposed to improve the quality of the recovered image and the speed of reconstruction ,making use of the theory of compressed sensing.The algorithm ,which combines the gradient sparsity of CT images and ART ,performs the projection on the corresponding hyperplane in each iteration ,and then adjusts the total variation with the gradient descend method to reduce the l 1norm of the image gradient.The experimental results show the effectiveness of the algorithm.

K ey w ords :image reconstruction ;algebraic reconstruction technique ;compressed sensing ;total variation

0　引　言

计算机断层扫描(Computed tomography ,简称CT )是利用射线穿过物体时,不同成分和厚度的物

体对射线的衰减不同来成像。它最早是由英国EM I 公司的hounsfield 发明的,并于1967至1970年间研制出第一台临床用CT 机[1]。目前CT 技术已广泛应用于医学、机械、航空航天、核工业等领域。

代数重建法(Algebraic reconstruction technique ,简称AR T )是CT 图像重建的一种重要方法。在CT 图像重建中,以滤波反投影算法(Filtered back projection ,简称FBP )为代表的变换法,具有重建速度快的优点,目前广泛应用于各领域,尤其是医学领域,但其对投影数据的完备性要求较高,而在实际应用中往往由于客观原因无法或很难获得完备的投影数据。

迭代重建的概念首先由G ordon 、Bender 、Her 2man 等人于1970年引入图像重建领域,第一台医用

CT 机所用的原理就是AR T 算法。AR T 算法适合

于不完全投影数据的图像重建[2],其最大的缺点是

重建速度慢,所需的内存空间大,难以用硬件实现。目前,对于提高该算法的重建速度,已经做了很多的研究,并取得了明显效果,主要包括投影数据的访问方式[3—5]、内存优化[6]和联合代数重建算法[7]等。随着计算机技术的发展,AR T 算法的缺点已经降为次要矛盾,它的优势更为突出,提高该算法的重建质量已经越来越引起人们的重视。

近三年来人们在图像重建领域有了重大突破。最近由美国斯坦福大学的Donoho 、Candes 等从信号分解和逼近理论基础上提出来的压缩传感(com 2pressed sensing ,简称CS )理论表明[8,9]:若图像在某个变换域具有稀疏表示,则可由少量的观测值重建出原始图像,观测值是通过对原始图像进行投影获得的。由于CT 图像具有稀疏表示(小波稀疏表示和梯度稀疏表示),利用压缩传感原理进行图像重建所需要的观测值远远少于奈圭斯特采样定理所要求

24Ξ收稿日期:2008206230;收到修改稿日期:2008212215 E 2m ail :lianqius @https://www.wendangku.net/doc/1f13657591.html,

基金项目:国家自然科学基金资助项目(60772079)

作者简介:练秋生(19692),男,燕山大学副教授,博士,从事多尺度几何分析、压缩传感、生物识别等研究。

的观测值。对于K 2稀疏图像(即在图像变换域中只有K 个有效值,其余为0或接近0),用3K ～5K 个观测值就可精确重建出原始图像

[9]

。

本文在AR T 算法的基础上,利用CT 图像具有梯度稀疏表示的先验知识,提出了一种结合压缩传感的代数重建算法,并通过实验对该算法的有效性进行了验证。

1　ART 算法

代数重建法(AR T )是一个迭代的过程。AR T 算法认为穿过物体的射线是有一定宽度的,首先把一幅n ×n 的图像f (x ,y )看作是由一个n ×n 的正方形网格组成的,在每个网格内f (x ,y )被视为常量。根据成像的物理过程和相应的数学模型,建立重建图像和投影数据之间关系的代数方程组

w 11f 1+w 12f 2+…+w 1N f N =p 1w 21f 1+w 22f 2+…+w 2N f N =p 2

…

w M 1f 1+w M 2f 2+…+w M N f N =p

(1)

即

p i =

∑

j =1

w ij f j ;i =1,2,…,M (2)

式中p i 对应于第i 条射线的投影值;M 是投影总数,它等于投影角度个数与某角度投影射线数的乘积;f i 代表第j 个网格(像素)内的常数值;N =n ×

n 是像素总数;w ij 代表权因子,反映了第j 个像素

对第i 条射线投影值的贡献,它的计算直接关系到重建速度和重建精度。在实际应用中,为了提高重建速度,往往对其进行简化,认为像素值集中在像素的中心,射线宽度等于像素宽度,当射线穿过像素中心时,权因子定义为1,即

w ij =

1,i 号射线经过j 号像素中心0,其他

(3)

式(1)可用矩阵表示为

P =W F

(4)

式中P 为M 维观测值向量;F 为N 维图像向量;W 为M ×N 维投影矩阵。图像重建的任务就是如何由

W 和P 求出F ,从而把图像重建问题转化为求解一

个大型线性方程组的问题。该方程组有着N 个未知数,M 个方程,在一般情况下M

[10]

。该算法首先给出一个初始解F 0,然后

进行如下迭代

f (k +1)j =f (k )j +λp i -

∑N

n =1

in f (k )n

∑N

n =1

w 2

w ij

(5)

式中λ为松弛因子(0<λ<2),λ根据投影数的多少、是否有噪声以及噪声的大小而有所不同,在迭代过程中要寻求合适的λ来加快收敛速度。

方程组(1)中的每一个方程代表一个超平面,每次按式(5)迭代一次,表示将当前重建出来的图像投影到相应的超平面上,对所有M 个超平面迭代完成后为一次完整的迭代,如果没有收敛可以再次开始迭代,直到满足一定的收敛条件。

2　基于压缩传感原理的图像重建

压缩传感是利用图像具有稀疏性表示的先验知识,由投影获得的少量观测值来进行图像重建。令Ψ代表某一稀疏性变换,Φ代表投影矩阵,Ψ和Φ是不相干的。对于图像F ,y =ΦF 是经过投影获得的观测值,则通过求解下面的非线性方程就可重建出原始图像

min u

‖Ψu ‖0,　

s.t.　Φu =y (6)其中零范数‖?‖n 代表非零元素的个数。式(6)为

非凸优化问题,是典型的NP 2hard 问题,计算复杂度非常高。因此Candes 和Donoho 提出用l 1范数来代替l 0范数[8,11]

min u

‖Ψu ‖1,　

s.t.　Φu =y (7)

其中‖u ‖1=

∑

i =1

|u i |,从而把问题(6)转化为

一个凸优化问题来求解,此时原始图像仍可精确重

建。在图像具有稀疏表示的前提下,式(6)和式(7)的解是等价的[8]。

在对图像进行稀疏变换时,常用的是DCT (Discrete cosine transform ),小波变换和有限差分等。用有限差分做稀疏变换时,常用全变差(Total variation ,简称TV )来衡量。由于CT 图像具有梯度稀疏性,本文采用图像梯度的l 1范数即全变差作为图像的稀疏性先验模型,此时图像重建方程成为

min u

‖u ‖TV ,　s.t.　

Φu =y (8)式中

‖u ‖TV =

∑

s.t

|u (s +1,t )-u (s ,t )|2+|u (s ,t +1)-u (s ,t )|2

(9)

3　基于ART 和压缩传感的重建算法

传统的AR T 重建算法没有利用图像的任何先

24第3期练秋生,等:　基于压缩传感和代数重建法的CT 图像重建

图2　shepp 2logan 图像不同投影角下重建的PSNR 变化情况

图3　实际CT 图像不同投影角下重建的PSNR 变化情况验知识。本文将CT 图像的梯度稀疏性先验知识加入到图像重建中,结合压缩传感的思想提出了一种新的CT 图像重建算法,算法的具体实现步骤如下:

(1)获得投影数据(观测值)P ;

(2)对投影数据进行排序,使得在投影过程中彼此相邻的超平面尽可能的正交,从而加快收敛速度;

(3)初始化:对重建图像

赋初值:F 0←0,设定初始迭代值:j ←1和迭代次数的最大值maxiter (本文中,maxiter

=10);

(4)AR T 迭代:对F j -1按式(5)在每一个超平

面上执行投影操作,执行完一次完整的迭代后,得到

重建图像F j ;

(5)全变差调整:采用梯度下降法减小F j 的全

变差‖F j ‖TV ,即F ′j =F j -a[5(‖F j ‖TV )/5F j ]。

式中a 是步长因子,它的取值与投影角度数有关;

(6)F j ←F ′j ,j ←j +1,如果j ≤maxiter ,返回

第(4)步继续进行迭代,否则停止迭代。

在获得投影数据和进行AR T 迭代时,都要用

到权因子矩阵W ,但W 是M ×N 的矩阵,M 和N 通常都很大,若事先计算并存储到硬盘中,使用时再调用,需要很长的时间和很大的存储空间。本文并不对W 进行存储而是直接获得投影值。在每次AR T 迭代时,采用文献[12]中提到的快速网格遍历算法来实时计算第i 个方程的权因子。该算法利用了一条射束穿过重建区域时只有极少数的网格包含在射束内,且数目不超过2n 的性质,有效减少了算法所需的内存空间

,加快了重建速度。

4　实验结果与分析

(a )shepp 2logan (b )实际CT 图像

图1　标准CT 图像

首先选用shepp 2logan 标准CT 模型进行实验(如图1(a )所示),图像尺寸大小为256×256像素。在此实验中对本文算法(AR T +TV )和AR T 算法分别在180、100、50个投影角度时的重建情况进行了比较。每个角度下的射线束为256条。不同投影角度数下松弛因子和步长因子的选取如表1所示。

表1　shepp 2logan 模型不同投影角度数下的参数设置

投影角度数松弛因子λ步长因子a

1800.80.02/j

100 1.050.03/j

1.2

0.06/j

图2是shepp 2logan 标准CT 模型在180、100和

50个投影角度数下迭代10次重建的PSNR 比较情况。在180、100和50个投影角度数下迭代10次时,本文算法重建图像的质量分别比AR T 算法高12.58dB 、11.69dB 和6.97dB 。

从图2可以看出,本文算法迭代10次时重建图像的PSNR 远远高于AR T 算法。AR T 算法在迭代3～4次后,PSNR 的增长趋于平缓,收敛于线性方

程组的解,而本文算法在迭代3～4次后的PSNR 仍

然保持较快的增长速度。在相同投影角度数下达到相同PSNR 时,本文算法所需要的迭代次数远远少于AR T 算法。

为进一步比较两种算法的性能,在实验中对实际的CT 图像(如图1(b )所示)进行了重建,该图像尺寸大小为256×256像素。投影角度数仍为180、100和50,每个角度下的射束数也是256条,不同投影角度数下松弛因子和步长因子的设置如表2所

24光　学　技　术第35卷

示。

表2　实际CT 图像在不同投影角度下的参数设置

投影角度数

松弛因子λ

步长因子a

1800.5 3.5/j 1000.7 5.5/j 50

1.0

7.5/j

图3是实际CT 图像在不同投影角度数下迭代10次的PSNR 变化情况。在180、100和50个投影

角度下迭代10次时,本文算法的PSNR 值分别比

AR T 高1.86dB 、2.19dB 和1.

87dB 。图4是实际CT 图像在50个投影角下本文算法和AR T 算法的重建图像。从图4可知本文算法主观视觉效果明显优于AR T 算法。AR T 算法的重建图像存在明显的伪吉布斯效应。

5　结　论

(a )ART :PSNR =27.88dB (b )ART +TV :PSNR =29.75dB

图4　实际CT 图像在50个投影角下迭代10次的重建图像

本文将CT 图像的梯度稀疏性先验知识加入到AR T 算法中,结合压缩传感的思想提出了一种新的CT 图像重建算法。对shepp 2logan 标准CT 模型和实际CT 图像进行重建的实验结果表明,本文算法重建图像的峰值信噪比和主观视觉效果均优于AR T 算法。

参考文献:

[1]庄天戈.CT 原理与算法[M ].上海:上海交通大学出版社,

1992.

[2]Muller K.Fast and accurate three 2dimensional reconstruction from

cone 2beam projection data using algebraic methods[D ].The Ohio State University ,1998:32—43.

[3]王宏钧,路宏年,傅健.代数重建技术中投影序列选择次序的研

究[J ].光学技术,2006,32(3):389—391.

[4]Guan H ,G ordon R A.Projection access order for speedy conver 2

gence of ART :a multilevel scheme for computed tomography[J ].Physics in Medicine and Biology ,1994,39(11):2005—2022.[5]Herman G T ,Meyer L B.Algebraic reconstruction can be made

computationally efficient [J ].IEEE Trans Med Image ,1993,12(3):600—609.

[6]秦中元,牟轩沁,王平,等.一种内存优化的代数重建算法及其

快速实现[J ].电子学报,2003,31(9):1327—1329.

[7]王旭,陈志强,熊华,等.联合代数重建算法中基于像素的投影

计算方法[J ].核电子学与探测技术,2005,25(6):784—788.[8]Donoho https://www.wendangku.net/doc/1f13657591.html,pressed sensing[J ].IEEE Transactions on Infor 2

mation Theory ,2006,52(4):1289—1306.

[9]Candes E J ,Romberg J.Practical signal recovery from random pro 2

jections[C ].SPIE Computational Imaging III ,San Jose ,2005.5674.

[10]G ordon R ,Bender R ,Herman G T.Algebraic reconstruction

techniques (ART )for three 2dimensional electron microscopy and X 2ray photography [J ].Theoretical Biology ,1970,29(3):471—481.

[11]Candes E J ,Wakin M.An introduction to compressive sampling

[J ].IEEE Signal Processing Magazine ,2008,25(2):1—30.[12]张顺利,张定华,王凯,等.一种基于ART 算法的快速图像重

建技术[J ].核电子学与探测技术,2007,27(3):479—483.

(上接第421页)

由设计结果可以看出,在整个视场范围内,点列

图弥散圆RMS 半径值小于2

μm ;传函已经达到衍射极限,接近理想状态。三条谱线(014μm 、0165

μm 、019

μm )在017孔径附近接近相交,实现了复消色差。各个透镜薄厚、间隔和光学工艺性合理(玻璃材料除了N 2KZFS11外,都采用国产常规玻璃),能够满足实际加工的要求。

3　结　论

本文结合摄远型物镜的光学和像差特性,通过

建立的能够分配三分离透镜组光焦度的复消色差方程组,并在此基础上,较好的利用了光学设计软件,在适当的光焦度搭配范围之内,采用摄远型的初始结构设计了一种宽光谱长焦距折射式物镜。设计过程中能够更加有目标和快速的达到理想的成像质量,通过实例验证了设计理论的可行性。而且在设计过程中具有较强的方法性,对于宽光谱长焦距折

射式物镜的设计具有一定的指导意义。参考文献:

[1]王庆有.CCD 应用技术[M ].天津:天津大学出版社,2000.[2]汪明强,李林,黄一帆.三反射镜空间遥感器的光学设计[J ].

光学技术,2007,33(2):170—172.[3]李晓彤,岑兆丰.几何光学?像差?光学设计[M ].杭州:浙江大

学出版社,2003.

[4]安连生,李林,李全臣.应用光学(第三版)[M ].北京:北京理

工大学出版社,2002.

[5]K inSlake R.Lens design fundamentals[M ].Beijing :Science Pub 2

lishing Company ,1978.

[6]Sigler R D.Designing apochromatic telescope objectives with liquid

lenses[J ].SPIE ,1991,1535:89—112.

[7]Morian H F.New glasses for optics and optoelectronics[J ].SPIE ,

1990,1400:146—147.

[8]Mercado R I.Correction of secondary and higher 2order spectrum

using special materials[J ].SPIE ,1992,1535:184—198.[9]萧泽新.工程光学设计[M ].北京:电子工业出版社,2003.[10]李润顺,王鹏,张爱红,等.波差法设计长焦距复消色差平行光管物镜[J ].哈尔滨工业大学学报,1996,28(2):44—49.[11]G eary J M.Introduction to lens design :with practical ZEMAX

examples[M ].Virginia ,Willmann 2Bell ,Inc ,2002.

24第3期练秋生,等:　基于压缩传感和代数重建法的CT 图像重建