（计算数学专业论文）图象重建块迭代算法的研究.pdf

e巫窑适厶兰亟翌位丝塞空塞趟耍 中文摘要 摘要：在图象重建算法中，最主要的两种重建算法是解析重建算法和基于迭代的重建算法， l a n d w e b e r 迭代算法是图象重建算法中基于迭代算法的重要图象重建算法。本文将针对 l a n d w e b e r 分块迭代算法中松弛参数的选取进行研究，在重建过程中采取对投影矩阵按投 影角度分块的方法，选取特定的松弛参数。按照本文采集数据模式进行重建图像。通过数 值实验得出结论：对于按角度分块的块迭代算法，松弛参数选取为 乘以块矩阵与其共轭 转置矩阵乘积的最大特征值分之，当采集完全投影数据，且凡接近l 6 - 1 7 时效果最好。 另外，在本文中，按角度分块的做法和松弛参数的选取方法对于有限角度图象重建问题也 是可行的，我们也做了相应的实验证实了这种可能性，并进一步证明：当图像大小和每个 投影角度下的射线条数一定时，投影角度增加，丑应适当变小，可以使得平均误差最小； 投影角度减小，a 应适当变大，可以使得平均误差最小。 关键词：代数重建算法：松弛参数；块矩阵：投影矩阵；有限角 分类号：t p 3 9 1 4 1 e 夏至亟厶：j 兰亟上堂焦论塞垦s ! 基i a b s t r a c t a b s t r a c t ：t h et w om o s ti m p o r t a n ti m a g er e c o n s t r u c t i o nm e t h o d sa r ea n a l y t i c a l g o r i t h ma n di t e r a t i v ea l g o r i t h m l a n d w e b e r si t e r a t i v ea l g o r i t h mi sa ni m p o r t a n t m e t h o db a s e do nt h ei t e r a t i v e a l g o r i t h m i nt h i s a r t i c l ew es t u d yt h er e l a x a t i o n p a r a m e t e r ss e l e c t i o ni nt h el a n d w e b c r sb l o c k i t e r a t i v ea l g o r i t h m s w es e l e c tt h e s p e c i a lr e l a x a t i o np a r a n m t e ra n dp u tt h ep r o j e c t i o nm a t r i xi n t om a n yb l o c k sa c c o r d i n g t ot h ep r o j e c t i o n a n g l e si nt h er e c o n s t r u c t i o np r o c e d u r e ，a n dr e c o n s t r u c tt h ei m a g e a c c o r d i n gt ot h ef o r mw eg e tt h ep r o j e c t i o nd a t as t r i c t l yw es e l e c tt h ep r o d u c to f k a n d o n er a t i ot h em a x i m a le i g e n v a l u e so ft h em u l t i p l e so ft h eb l o c k m a t r i x e sa n dt h e i r c o n j u g a t et r a n s p o s em a t r i x e sa st h es p e c i a lr e l a x a t i o np a r a m e t e r s b yp e r f o r m i n g n u m e r i c a le x p e r i m e n t ，w ec o n c l u d et h a ti f 入a p p r o a c ht o1 1 6 - 1 7 w ew i l lg e tt h eb e s t r e s u l tw h e nw ec a no b t a i nt h ec o m p l e t ed a t u m ；o t h e r w i s e ，t h a tw es e l e c tt h es p e c i a l r e l a x a t i o np a r a m e t e ra n dp u tt h ep r o j e c t i o nm a t r i xi n t om a n yb l o c k sa c c o r d i n gt ot h e p r o j e c t i o na n g l e si sa v a i l a b l et ot h el i m i t e da n g l ep m b l v m i nt h i sa r t i c l e ，w eh a v e c e r t i f i e dt h ep o s s i b i l i t yo ft h ea b o v eb yp e r f o r m i n gt h ec o n c r e t en u m e r i c a le x p e r i m e n t f u r t h e r m o r e ，w eh a v ef o u n d e dt h a to nc o n d i t i o nt h a tt h es i z eo ft h ei m a g ea n dt h e n u m b e ro ft h el i n e sa te a c ho fp r o j e c t i o na n g l ea r ec e r t a i n ，w e c a n m a k e t h ea v e r a g ee r r o rs m a l l e rb yr e d u c i n gt h e 旯p r o p e r l yw h e nt h en u m b e ro ft h e p r o j e c t i o na n g l e si n c r e a s eg r a d u a l l yo r i n c r e a s i n gt h eaw h e nt h en u m b e ro f t h ep r o j e c t i o na n g l e sd e c r e a s e k e y w o r d s ：a r t ；r e l a x a t i o np a r a m e t e r ；b l o c k m a t r i x ；p r o j e c t i o nm a t r i x ； l i m i t e da n g l e c l a s s n o ：1 1 p 3 9 1 4 1 学位论文版权使用授权书 本学位论文作者完全了解北京交通大学有关保留、使用学位论文的规定。特 授权北京交通大学可以将学位论文的全部或部分内容编入有关数据库进行检索， 并采用影印、缩印或扫描等复制手段保存、汇编以供查阅和借阅。同意学校向国 家有关部门或机构送交论文的复印件和磁盘。 ( 保密的学位论文在解密后适用本授权说明) 学位论文作者签名： 导师签名： 签字日期：年 月 日 签字日期：年月日 e 立窑擅厶生亟堂位熊塞 丝剑鲑岜唆 独创性声明 本人卢明所薯交的学仲论文是本人在导师指导卜进行的研究i ：作和取得的研究成果，除 了文中特别加以标注和i 致谢之处外，论文中不包含其他人已经发表或撰q 过的研究成果也 不包含为获得j 匕京交通人学或其他教育机构的学位或证i ；而使圳过的材料。与我一同l ：作的 同。占对本研究所做的任何贡献均己在论文中作了明确的说明并表示了谢意。 学位论文作者签名：签字日期：年月日 致谢 本论文的工作是在我的导师渠刚荣教授的悉心指导下完成的，渠刚荣教授严 谨的治学态度和科学的工作方法给了我极大的帮助和影响。在此衷心感谢两年半 来渠刚荣老师对我的关心和指导。 渠刚荣教授悉心指导我们完成了实验的科研工作，在学习上和生活上都给予 了我很大的关心和帮助，在此向渠刚荣老师表示衷心的谢意。 渠刚荣教授对于我的科研工作和论文都提出了许多的宝贵意见，在此表示衷 心的感谢。 在工作及撰写论文期自j ，张玉川、毕清华、马有为等同学对我论文中的一些 细节问题的研究工作给予了热情帮助，在此向他们表达我的感激之情。另外也感 谢我的母亲和我的姐姐，他们的理解和支持使我能够在学校专心完成我的学业。 e立窑埋厶 生亟= ；三位迨塞 虚 序 从理论上讲，c t 是一个从投影重建图像的反问题有很广的适用性。广泛用 于医学、工业无损检测、地球物理勘探和农林业等方面。c t 技术的医学应用与非 医学应用各有其自身的特点。但是，总的来说图象重建算法有三类；解析重建算 法，如反投影重建算法和f o u r i e r 算法等：系列重建算法，如代数重建算法等：另 外，还有其它重建算法。在得到确切完全的投影数据时，解析重建算法( 如卷积 反投影重建算法) 重建质量高、并且重建速度快，是一种理想的图像重建算法。 但是，正是因为此种算法对投影数据的采集要求较高，系列重建算法才又一次引 起了人们的重视。 在系列重建算法中，迭代重建算法历来是重要的图像重建算法，但是由于迭 代重建中数据运算量大，运算速度慢等特点，限制了这种算法在实际中的应用。 在另外一方面来说，迭代重建算法也有其自身的优点，比如：如果在投影数据不 是均匀分布的，或者无法得到足够多的投影数据等情况下，反投影重建算法就不 如迭代重建算法具有优势。所以我们有必要对迭代重建算法做更加深入的研究， 充分利用它的优点，努力克服它的缺点，充分发挥迭代重建算法在实际应用中的 优势，将有着重要的现实意义。很重要的是，在迭代重建过程中技巧的应用将会 给重建质量带来很大的改进。本文针对松弛参数的选取问题进行研究，并得到了 一定的结果。 本论文结构分布：第2 章图像重建中的解析算法，主要介绍解析算法的起源 以及解析重建算法的主要内容；第3 章图像重建中的迭代重建算法，主要介绍迭 代算法的由来和几种主要的迭代算法内容：第4 章松弛参数的选取问题，主要介 绍松弛参数选取中的一般方法和问题；第5 章松弛参数选取的研究，主要介绍本 文的工作，以及试验结果：最后是结论总结。 在论文写作的过程中，渠刚荣教授悉心指导我们完成了实验的科研工作，在 学习上和生活上都给予了我很大的关心和帮助，在此向渠刚荣老师表示衷心的谢 意。渠刚荣教授对于我的科研工作和论文都提出了许多的宝贵意见，在此表示衷 心的感谢。 另外在本人做论文期问，感谢张玉川，毕清华、马有为等同学的热情帮助。 感谢我的母亲和我的姐姐，他们的理解和支持使我能够在学校专心完成我的学业。 e瘟銮堑厶堂亟堂位监室i l矗 l 引言 目前，图象重建算法中最主要的两种是解析重建算法和基于迭代的重建算法， 解析重建算法有滤波反投影重建算法和f o u r i e r 重建算法等。解析重建算法在一 定条件下，重建速度快、重建质量高，是当下比较流行的算法。 解析重建算法( 如滤波反投影算法) 对数据的采集有一定的要求，如果投影数 据不是均匀分布的，或者无法得到足够多的数据等情况下，滤波反投影重建算法 就不如代数重建算法具有优势。文献 2 研究了迭代松弛参数的选取在一般 l a n d w e b e r 迭代算法，包括同时迭代算法( 简称s a r t 算法) 和块迭代算法中的收 敛性问题。文献i s 用矩阵的奇异值定理研究了s a r t 迭代算法对更广泛松弛参数 选取的收敛性。怎样选取松弛参数加快收敛速度，并改进图像的质量是人们所关 注的。但是，现在我们没有好的方法去选择松弛参数，可靠的方法是通过实验验证， 这个问题在理论方面需要进一步的研究【4 】。 本文在所选松弛参数满足文献 2 和 3 中块迭代算法和s a r t 算法收敛性的条 件下，通过数值实验，研究了投影矩阵按角度分块的块迭代算法中松弛参数的不 同选取对迭代算法的收敛速度和成像质量的影响。实验表明当按本文采集模式采 集完全投影数据时，对于标准s h e p p l o g a n 图或其他图像，选取松弛参数为x ( 0 1 ，x r ) 乘以块矩阵与其共轭转置矩阵乘积的最大特征值分之一的重建效 果比直接选取为最大特征值分之一要好，且a 接近1 6 - 1 7 时效果最好，重建图像 与原始图像之间的平均误差最小。按角度分块的方法和松弛参数的选取方法，对 有限角度图象重建问题有启发作用，我们也做了相关数字实验，并且得到可以接受 的结果。另外，本文通过试验发现：当图像大小和每个投影角度下的射线条数一 定时，投影角度增加，五应适当变小，可以使得平均误差最小；投影角度减小，五 应适当变大，可以使得平均误差最小。 e 塞窑道厶堂亟堂位迨塞幽堡重建主的鲤蚯簋壁 2 图像重建中的解析算法 2 1拉东变换及其逆变换 层析技术的图像重建计算首先是拉东( r a d o n ) 在1 9 1 7 年完成的，后人称之为拉 东变换。在二维空刚中，若用极坐标( r ，矿) 代替直角坐标( 石，y ) ，如图卜l 所示： v 翁一， 这时有工= r e o s 口o ) ，= r s i n q p = a r c l g ( y x ) 在x o y 平面内，某被测物理量 厂( r ，妒) 沿直线s 的积分为 p ( 1 ，口) = c ，( 扩了，0 + n ，c f g ( 训凼 ( 2 1 ) 记为观厂( ，秽) ，它是厂( ，妒) 在8 方向上的投影，称为函数厂( r 伊) 的拉东变换，其 拉东逆变换l 孵“pl ( ，妒) 为 h ) = 订1re 篙筹争口 旺z ， 即 巾纠= 专r c 篙筹争口 眩， 理论上证明：只要投影采样间隔无限小，并能得到所有扫描积分测量的投影值，就 能精确地得到场分布，( ，妒) ”1 。 实现拉东逆变换的困难： 1 ，在c t 中，我们仅能得到有限的测量值。从有限数的线积分不足以唯一地决定 一个图像。也不可能准确，即有数据不完全的困难。 2 ，c t 中得到的测量值，只能用以估算线积分值。这些估算中的不准确，是由于x 2 e 立窑垣厶翌亟堂位迨塞笪傻重建g ! 曲鲣蚯簋选 一射线束是有宽度、散射、束的硬化、光子统计及检测器的不准确等造成的，拉 东逆变换对这些不准确是敏感的。 3 对拉东逆变换不宜于做数值计算，我们需要一个有效的、不烦琐和速度上允许 的算法【6 1 。 2 2 解析重建算法 为解决逆拉东变换不宜作数值计算的困难，我们对逆拉东变换进行变换，产 生了卷积反投影算法。平行束的卷积重建算法近似公式可用两步来实现： l 计算关于第一个变量卷积：记岛( ，) = p ( ，口) ，一( ) = _ 4 石e ”“巴( “) c o s ( 2 # u l s ) 出，则关于第一个变量的卷积为： 岛一) = e p ( ，护减( p - p ) d p ( 2 4 ) 2 计算图像近似值： 厂( r ，v ) = - ( 1 2 , 0 r d 口e p ( ，口减( ，c o s ( o - v ) - p ) ( 2 5 ) 卷积法具有计算简单，其运算速度快、所需存储量小以及重建质量高的特点， 所需费用较小。但是，只有在投影数据确切、数据量丰富的情况下，重建质量才 可以和其它重建算法取得的图像相媲美。 些盛窑迪厶芏亟芏位硷塞幽塑重建巾盥丝瑾堡建珏篮 3 图象重建中的迭代重建算法 计算机层析成像算法大致可以分为：精确算法、近似算法和迭代算法。迭代 算法中，代数重建算法是提出最早，最为人们所熟知的算法。它有许多优点，但 由于计算量大、重建时间长，在很长一段时间内限制了其在医学和工业领域的应 用。提高迭代型算法的计算速度一直是人们关注的问题。近年来人们提出了不少 提高迭代算法计算速度的方法，加上计算机计算速度的提高，迭代算法重新受到 人们的亲睐“】。 3 1迭代重建算法 迭代重建算法是图象重建算法中求解线性方程组的重要方法，大致有两种形 式：逐线进行修正的代数重建方法( a r t ) ；逐点进行修正的联合迭代算法( s i r t ) 。 另外，每一种迭代算法又主要包括两种修正方法：加法和乘法【6 j 。 3 1 1系列重建算法的离散重建问题 系列重建算法开始就将需要重建的图像离散化。如果假设密度在像素内取平 均值，那么，估计函数值被转化为寻找一个有限数集：假设我们固定j 个基图像 岛，岛，t ，它的线性组合使我们能够逼近所希望的任何图像，。此处，用系列 逼近r ，故称为此列重建法。 在此处，我们取，= 月2 ( j 为象素个数) 的数字化图像，并且将像素从1 到， 进行编号，并定义： 惭) = 船勰扛凳价镰内 ( 3 - 1 ) 作为基图像。图像厂的n 数字化是图像厂，并定义： ，( 葺_ ，) = e 妒，( z ，y ) ( 3 2 ) 其中，x ，是图像，在篼j 个像素内的平均值，并且记x = 而乇，。，勺) 为图像向量。 基图像的选取还有其它方法，一旦基图像 乜，也，6 j 选定之后，就能用基图 像的线性组合来表达任何图像r ，因此，将一般的“估计图像，”的问题限制为 “寻找一个图像向量x ”，使得( 3 2 ) 式所定义的，( 使用给定的基图像) 尽可能 e 丛窑迪厶堂亟堂位途塞幽叁里建主曲丝岱重建簋选 接近原图像厂。 由拉东变换胄和冠( i = l 2 ，在此，假设射线束是由一条条平行的具有一 定宽度的射线组成的，为假定的射线条数) 的定义，可得 ， r i f - 兰弓，= j ，r i b y ( x ，y ) ， ( 3 3 ) ，。l 其中，f = l ，2 ，j ，j = l ，2 ，j ，6 f ( j ，y ) 的定义如 ( 3 1 ) ，r , b j ( j ，y ) 刚巧就是光源一检测器对的第f 个位 置上的直线和第个像素交线的长度，如图1 ，令 = 鹃( x , y ) ，假设y 是投影数据，y j 是y 的第f ( f = 1 ，2 ，i ) 个分量，则由( 3 3 ) 式： y i 兰x j a 0 ( 3 4 ) i - 1 其中，f = 1 ，2 ，1 ，_ ，= l ，2 ，j 。 我们称矩阵r = i ) 为投影矩阵，由于每条射线仅 图1 口角度下的投影射线 与很少的像素相交，可知对于任意固定的第i 个位置，口。将会有很多为零，易知r 为稀疏矩阵。令e 为，维向量，的第i 个分量岛是( 3 4 ) 式中左端与右端之差， 称e 为误差向量。于是，( 3 4 ) 式可以改写为： y = r x + e ( 3 5 ) 因此，代数重建方法可以表述为方程( 3 5 ) 式，给定数据y ，在误差向量e 允 许的范围内估计出图像向量x 。如果离散重建问题解的估计式是向量r ，则被重 建图像的估计值是厂： ， f = e ( 3 6 ) 3 1 2 将图像重建问题划归为解线方程组的方法 直接离散方法 如3 1 1 所述，我们假设图像向量x ，根据每个角度和每一条射线的投影值建 立方程( 3 4 ) 。这样，把由所有的射线得到的方程组合起来就组成了线性方程组， 从而，图象重建问题转化为解线性方程组问题。 从“投影“角度引入 所谓一个投影值既是一条射线的线积分，现在我们假设第i 条射线的线积分己 知( 1s i s i ) 即 e 峦窑堑厶堂亟：兰位途塞幽塑重建生趁丝岱重建簋遮 量膨= 咒 ( 3 _ 7 ) ， 此时，我们离散等式的左端可得：咒兰石，口f ( 1 f ，= 1 2 ，j ) ，即得同直 接离散方法相同的等式( 3 4 ) ，以下同恐1 1 ) 分析。于是，图象重建问题可以 转化为解线性方程组的问题。 3 1 3 按角度分块 在计算投影矩阵时，按投影角度分块可做如下解释( 如图1 所示) ：假设投影 角度固定，则在此投影角度下，任何一条射线与所有象素的交线长可以作为一 行元素。那么，所有条射线和所有象素的交线长就组成了一个矩阵，就是所指的 块矩阵。当角度口变化时，块矩阵也会随之变化。这样，在不同的角度下得到了 不同的块矩阵，作为每次迭代中的投影矩阵。 3 1 4 迭代重建算法的讨论 迭代算法具有解析重建算法不具备的优点。实际应用过程中，有时无法测到 大量的投影；有时投影不是均匀分布在1 8 0 或3 6 0 范围内。例如，为了避免物体 ( 如脏器) 运动，或为了减少剂量，投影数据采集不足；又如某一物体在某一方 向上特别长，在此方向上无法或较难测得数据。再如用成像法勘探地球资源时， 用孔穴法采集的数据，数量不足又分配不均匀。此外，在有衍射和折射的场合， 射线不是沿直线行进，滤波反投影方法则不如迭代重建方法优越川。 迭代法的最大缺点是速度较馒，不能通过硬件实现。数据量较大，例如a r t 重建方法中，对于一幅n * n 的图像，取m 个投影，每个方向的投影有n 条射线， 如果直接用a r t 方法重建，则系数矩阵的元素个数约为0 ( n 4 ) 。对于一副2 5 6 * 2 5 6 的图像或者更大的图像，这个数据是很大的。如何节省内存，提高运算速度是十 分关键的问题1 引。但是，a r t 方法有它自身的优点所需数据较少，比较灵活的 特点将更加突出，这将会打破滤波反投影算法一统天下的局面。 3 2 代数重建算法( a r t ) 3 2 1a r t 算法 自1 9 7 0 年以来，众多学者对a r t 算法进行了大量的研究和实践，分析t a r t 算 6 e 盛銮道厶堂亟：i 兰位论塞崮塞堡建虫丝丝i l 重建簋垄 法并予以改进。影响a r t 算法的因素有： ( 1 ) 松弛因子的选择， ( 2 ) 最优迭代次数的估计， ( 3 ) 初值的选取， ( 4 ) 先验知识的多少。 所有图象重建的a r t 方法都是迭代过程，即他们产生一系列向量 x ，工”，一”，使得序列收敛于工。这就意味着，对l ，如果k 选得足够大， 叠鼬保证任意趋近于工：。a r t 的第k 步迭代能用函数来描述，它的变量是两个 ，维向量和一个实数，它的值是一个，维向量。于是， 一“= 口x ( j “，版，) ( 3 8 ) 其中石和是两个，维向量，是一个实数值。不同的a r t 方法所选用的 各不相同，但是，图象重建问题就是要解线性方程组( 3 5 ) 。 3 2 2a r t 算法一般步骤 解线性方程组，a r t 算法的一般步骤： l ，对x 的各个分量设初值， 工f = 砰，= l ，2 ，乃： 2 ，计算第i 个方程的估计值， 即一= 吩掣； 3 ，计算 ，= y i 一一： ”1 4 ，计算第f 个方程的修正值： 厂、 c o = a ，x 1 名( ) 2 ， 5 ，对各个工；值进行修正， 彤枷= + 掣 6 ，将对第一个方程修正后的各z ，= 工：0 ) + g ，代入第二个方程求废，重复2 5 步 骤，即由第二条射线对各个z ，值进行修正，然后再带入第三个方程求硝，。 如此下去，直到第，个方程； 7 ，经过第l 到，个方程的修正，得到的各z ，值完成了第一轮迭代。一般经过足轮 迭代，各工；的值和前一轮中各x j 的值之间的误差均小于指定的正数，则可以 认为迭代达到了收敛的要求，可以接受最后的迭代值。否则，继续迭代。 对于a r t 算法，在迭代过程中，为了得到第( k + 1 ) 次迭代值，我们应用口。到 第k 次迭代值0 ”，投影矩阵a 的行及测量向量y 的第个分量上。这种算法只 7 e 巫窑堑厶主亟：羔位盈塞笪叁重建生啦堡毡重建簋选 需要存储上次迭代得到的一的值，一次迭代完毕，可以用原来存储一。，的存储位 来存储最新得到的工球“的值，所以a r t 算法被称为存储器效率较高的算法。另外， 当获取的投影数据不完全或是存在噪声时，a r t 方法较卷积反投影方法要好。 3 3 联合迭代算法( s i r t ) 讨论 3 3 1 联合迭代算法 联合迭代重建法又称为同时迭代重建法，它和a r t 最大的区别是：a r t 每次 仅考虑一条射线，而s i r t 是利用一个像素内通过的所有射线的修正值来确定对这 一像素的平均修正值。另外，s 玎是所有射线通过方格网后，才算完成一次迭代， 这样可以压制一些干扰因素，而且计算结果与射线的使用次序无关。而a r t 算法 则与射线的使用次序有关，并且上一次计算的误差可能影响到下一次的计算结果。 联合迭代算法基本公式： o 砖“= 砖+ 0 n a 碟 ( 3 9 ) i = l 厂jj、 其中c 掣= i ( 乃一_ ) ( 勺) 2i x r , x ，i = l ，2 ，i ，j = l ，2 ，j 。虬是通过第足 k卢ij - l， n r 个像素的射线条数，叮是迭代次数，( 1 虬) x 艺c 掣是对第k 个像素的平均修正值。 i = l 联合迭代重建算法虽然可以抑制一些噪声影响和其它一些干扰因素，但是由 于此算法迭代速度相比之下较慢。只有在测量数据特别不准确时，s i r t 才显示出 重建的优越性。除此之外，无论在重建质量或是在重建速度上均不如解析重建算 法和a r t 算法优越。所以，此算法的广泛应用受到了限制。 3 3 2 联合迭代算法原理 s i i 仆一同时迭代重建法是与a r t 并行的又一种迭代重建算法，而在实际重 建过程中不可避免地会引入误差，也即需修改式： g = a f + e( 3 1 0 ) 其中e 为误差向量，可以认为e 是包含在测量和计算过程中所产生的所有误差。s i r t 算法的提出旨在使重建图像对测量误差不敏感。我们知道a r t 算法每次迭代只用 到一条射线的线积分值( 或射线投影值) 。如果这一射线投影值包含误差，则所得 的解跟着引入误差。要减少误差，办法之一是使上式的解满足最j 、, - - - 乘准则，即 e 塞塞堑厶堂亟上堂位迨塞幽筮重建主的丝丛重建簋鎏 ( f ) = ( g a f ) 7 ( g - a f ) 最，j 、， 其条件是： a 7 g = a 7 a f( 3 1 1 ) 解出方程中的f ，即是所求满足最小二乘准则的解。( 3 11 ) 式的物理意义是，a 7 g 为g 的反投影，g 为测量矢量。因a f 表示伪射线和，所以a 7 a f 表示伪射线和的反 投影。如果求解满足最小二乘准则，那么，这两种情况下的反投影重建应该相等。 或者说估计图像所得射线和的反投影值等于实测投影的反投影值，计算出图像矢 量f 即可】。 3 4块迭代算法与其它迭代算法 3 4 1 块迭代算法 在分块大小上，s a r t 块迭代算法介于前两者之间的一种算法。他将投影矩阵 分块，既继承了a r t 算法的优势：重建质量高、重建速度较s 琢t 快的特点，又 继承了s i r t 算法的优势：抑制了噪声和其它一些干扰因素的影响。其中，l a n d w e b e r 分块迭代算法是最重要的一种分块迭代算法。 l a n d w e b e r 迭代方法是图象重建中的一般迭代算法。迭代方程为【2 】： ( “”= x t + 无v 一1 a g ( s a x ( ) ， ( 3 1 2 ) 其中，k = 0 1 ，2 ，s ，膏为迭代次数，k 为松弛因子，v ，w 分别为，阶和 m 阶方阵，且其元素非负的对角矩阵，a 为 彳行，列矩阵，a 。为4 共轭转置矩阵， j p 为图象向量，口是投影矩阵【2 1 。 假设6 = l ，2 ，m - u 包，将m 分成r 份。任意两个岛的交集可能非空， i $ d g t 并且任意一个6 ，都不为空， 岛= ，- - ti ，) ，任意n o ，【n = n ( m o d t ) + l ， 并且 丘= 形= 4 饥l 矿 9 ，( ，p 肼i ，) e 塞窑垣厶堂亟堂位盈塞幽盈重建虫的丝重建簋选 口= 一， k ”f 。) 分别是4 ，和口对应于6 的分块。a 。是a 的第f 行，是的第i 个对角线上的 块矩阵和b 是口的第i 个元素。根据以上分块，可以把迭代公式变换为： x “= x ”+ 五矿。禾】嘶i 】x l k + j ) ( 岛 1 4 t l x ) ( 3 1 3 ) 其中，k = 0 ，l ，2 ，3 ，s 。从公式可以看出，当t = i 时，公式( 3 1 3 ) 转化为公 式( 3 1 2 ) 。 3 4 2 a r t 、块迭代算法、s i r t 三者之间的比较 a r t 、块迭代算法、s i r t 三者之间有很多相似的地方：如都属于迭代算法； 都具有卷积反投影所不具备的优点，如：对于投影数据缺失重建效果良好和可以 运用迭代技巧的优点等。另外，他们还有很多不同之处： 1 通过实验发现，从重建速度上来讲，三者之中，a r t 最快，块迭代算法次 之，s i r t 最慢； 2 从存储量上来看，a r t 所需存储空间最少，块迭代算法次之，s i r t 所需存 储空间最多： 3 从重建质量上来看，a r t 至少与卷积反投影方法相当，块迭代算法在选取 恰当松弛参数时较a r t 好，而s i r t 算法，虽然总是收敛的，但是只有在测量数 据特别不准确时才显出它在重建质量上胜过其它算法。 在技巧的运用上来看，a r t 、块迭代算法、s i r t 三者均可以运用卷积反投影 方法不可使用的策略，如：可以根据先验知识选择初值，来加快迭代速度；可以 通过选择恰当的松弛参数来改善重建的质量等：在应用上，a r t 、块迭代算法、 s i r t 三者适用于投影数据缺失、噪声较大、射线弯曲等情况，而卷积反投影方法 最易用于投影数据完全、重建图像基本固定的情况下。 1 0 e 峦銮适厶堂亟堂位迨塞丝丝叁筮的垄躯四壁 4 松弛参数的选取问题 在使用迭代算法重建过程中，最重要的技巧之一有松弛参数的选取。如何选 取最恰当的松弛参数以得到最优的重建结果，是非常关键的。但是由于对松弛参 数的选择没有理论根据，只有一些实践的经验值，因此松弛参数的选择主要是通 过试验的方法来证明其合理性。我们有必要对松弛参数的选取从理论上作进一步 的研究。 4 1 松弛参数选取的一般方法 我们可以用形如 乙”西sq i 的不等式组来代替爿程( 前面要解的线性方程组，p 个方程) ，改写为内积形式 “，x ) - q i( 1 f p ) 其中珥是给定的j 向量，q i 是给定的实数。 我们假定f l i 中至少有一个像素为非零向量，如果珥的所有像素都为零，那么， 其对应方程对重建结果不影响。 引入某个向量集m 和 川= x l ( 吩，x ) q j ( 1 i p ) ， = n m i = 1 寻找 厂的一个元素。给出一个a r t 型的方法，它由下面函数序列定义，称为 不等式松弛法：任意 x 似+ 1 ) = j + a 。( ( ：二j ；2 三：；嫉耄，气，) 。，其它 h ， 其中a 是一个实数，即通常被称为松弛参数，= k ( m o d p ) + 1 ，即 f 0 = 1 ，= 2 ，i p l = p ，i p = l ，。 定理：如果松弛参数对所有k 满足弱条件 0 占a k g ， 2( 4 2 ) 那么，不等式松弛参数法产生一个序列x 仰，x m ，x m ，仅当n 是非空时，这个 序列收敛于n 中的一个向量。 e 基至适厶堂亟上堂位迨塞丝塾墨墼的选抠四星 注意l 由公式可知a r t 方法就是松弛参数选取为l 时的一种等式松弛 算法，有不等式收敛性定理，可知a r t 算法是收敛的。 注意2a r t 算法中，我们可以利用松弛参数的选取来改善重建质量。一般情况 下，选取小松弛参数会得到好的重建结果，另外，初值x 的选取也关系到最后的 重建质量。 松弛参数的选取是一个复杂的过程，如文献 6 】所述，现在我们没有任何好的 方法去选择松弛参数，仅仅可靠的方法是通过实验验证，这个问题在理论方面需 要进一步的研究。 一般，选取松弛参数丑的方法有： l = 2 ( + ) ，其中和分别是投影矩阵的最大和最小特征值。 2 五= 旯( v ) ，其中是投影矩阵的最大特征值，五是o 旯 1 的实数。 3 五= 2 + 几。一( 一p m ) c o s 石( 2 p ( 七) 一0 8 ，当我们知道迭代的 总次数时适用。 4 2松弛参数选取的一些观点 通过对松弛参数的恰当选取可以大大改善重建图像的质量，但是我们至今还 没有确切的方式去选择合适的松弛参数，大多数只是在经验上、实践上证明如何 选取松弛参数会更理想。所以说，如何从理论上选取恰当的松弛参数成了松弛参 数选取的基本问题。 迭代算法实际上是一个解超大规模线性方程组的带松弛参数的迭代过程，松 弛参数的选择对算法的执行速度有很大影响，松弛参数选择的一个基本标准为式 ( 4 2 ) 。在这个区间内a r t 算法才能按照某种标准保持其收敛性。一些学者认为松 弛参数的选择可以根据不同的目的，通过试验的方法来选择最佳的松弛参数。其 基本做法是：根据不同的重建任务，用实验的方法来模拟实际环境，选取不同的松 弛参数重建图像，从重建结果中选择被认为最好的松弛参数，然后在实际的医疗 科研中采用这个松弛参数进行图像重建。尽管这种方法在有些时候能取到较好的 图像效果，但是实验室是不可能模拟出千变万化的实际环境的。因此这种方法有 很大的局限性。 另一种松弛参数的选择方法是选择较小的松弛参数，一般的选择标准是使松 弛参数处于区间( 0 0 2 5 0 2 5 ) 之间，特别是在投影数据噪声比较大情况下。按照 这种方式选择松弛参数，重建后的图像很平滑且噪声较小。但是为什么会这样， 没有一个较好的理论证明。一个合理的解释是：小的松弛参数起到低通滤波作用， 2 巫窒适厶堂亟堂位途塞控熊叁熬的垂熟四壁 从而有效地避免了高频成分的出现，由于高频成分往往含有许多噪声，因而选择 较小的松弛参数就有效地抑制了噪声。但是这种方法的缺点是图像中本应有的高 频成分被去掉了，从而影响了原始图像的真实性。前两种松弛参数的选择方法都 是常数的形式，这种常量的形式易于实现，但却有一个最大的不足：它不能反映投 影数据的变化”。 在本论文中，我们首先找到分块矩阵与其共轭转置矩阵的乘积的最大的特征 值，然后选取 = a ( 1 ) ( 其中o o ，并且0 p 2 五2k = o ，l ，2 ，3 ，s ，那么，定理l 的结论 对于任意投影矩阵和及其相应投影数据的分块都是成立的1 2 】口 松弛参数 选取为l “( 以是1 4 。】的最大特征值) ，那么由于 1 1 4 。4 ，j i ：= 4 l = 以，则当0 冬。i i ：= p 时，9 0 _ 2 4 2 2 4 6 1 9 朱扬明，庄大戈不完全投影数据图像重建的最小交叉墒算法上海交通人学学报v 0 1 2 6 ， n o 3 ，1 9 9 2 ：p 2 9 3 5 2 0 m i n a r b og m e n t ：hm a x i m u me n t r o p yh l g o r i t h mf o rr e c o n s t r u c t i n gas o u r c ef r o m p r o j e c t i o nd a t a c o m p u t e rg r a p h i c sa n di m a g ep r o c c e s s i n g 1 9 7 9 ：( 1 0 ) ：4 8 6 8 2 1 m i b r a h i ms e z a n t o mg r a p h i ci m a g er e c o n s t r u c t i o nf r o mi nc o m p l e t ev i e wd a t ab y c o n v e xp r o j e c t i o n sa n dd i r e c tf o u r i e ri n v e r s i o ni e e et r a n s o nm e d i c a l