（信号与信息处理专业论文）基于灰色系统的ct不完备投影重建算法的研究.pdf

大连理工大学硕士学位论文 摘要 计算机断层成像( c t ) 是一门用来获取观测目标断层图像的技术，它广泛地应用于医 疗诊断和工业无损检测等领域。具有无损、无接触、高分辨率等优点，被国际无损检测 技术界称为最佳无损检测手段，随着c t 技术在多个领域的应用，c t 重建算法也因此成 为研究的热点。 本文主要对c t 的重建算法进行研究，其内容从两个方面展开，即不完备投影重建算 法和多机a i 盯重建算法。 现阶段的c t 重建算法一般采用完全投影数据重建的，而实际中由于诸多因素的影 响，难以获得完全的投影数据，给图像重建带来了困难。针对这种情况，本文主要讨论 了c t 重建算法中经常出现的三种不完备投影问题：内部问题、外部问题和有限角度问 题。投影数据不完备可能由设备引起，也可能由于被测物体的自身密度过高致使x 射线 难以穿透引起。那么探索一种通用的方法来解决投影数据不完备的情况就显得至关重 要。本文引入了一种新方法，利用灰色系统理论中的g m ( 1 ，1 ) 灰色模型，建立动态微分 方程对未知数据进行预测，从而构建完全的投影数据，最后利用滤波反投影算法进行重 建。同其他的局部重建方法相比，这种方法减少了辐射量，可以更有效、更便捷地解决 不完备投影问题，而且对外部缺失、内部缺失和角度缺失情况都适用。因此，利用g m ( i ，1 ) 灰色模型解决不完备投影问题，具有理论研究价值与实际应用意义。 本文的另一个工作就是局域网中多机迭代重建( a r t ) 算法，a r t 算法计算量巨大， 为了达到理想的重建时间，提高重建速度很有必要，但由于a r t 算法的局限性，仅通 过算法的优化很难再有较高的提升空间。本文提出一种提高a r t 重建速度的方法，利 用操作系统提供的数据共享传输功能，通过在主机上开辟共享缓存区，在其中设立数据 文件，其他从机分别调用共享区的数据进行重建，并将重建结果传输给主机，最后由主 机对来自不同从机的结果选取不同的融合方法进行图像融合，得到理想的c t 断层重建 结果。该方法原理简单，易于实现，重建时问成倍减少，具有很高的性价比。 关键词：c t 重建；不完备数据；灰色系统；a r t 并行算法；图像融合 大连理工大学硕士学位论文 t h es t u d yo f i n c o m p l e t ep r o j e c t i o nc t r e c o n s t r u c t i o n b a s e do ng r a ys y s t e m a b s t r a c t c o m p u t e dt o m o g r a p h y ( c di sat e c h n i q u ef o ri m a g i n gc r o s s s e c t i o n so fa no b j e c tu s i n g as e r i e so fx - r a ym e a s u r e m e n t s i th a sb e e nw i d e l ya p p l i e di nd i a g n o s t i cm e d i c i n ea n d i n d u s t r i a ln o n - d e s t r u c t i v et e s t i n g c o m p u t e dt o m o g r a p h yi so n eo ft h em v a n c e dm e a n sf o r n o n - d e s t r u c t i v et e s t i n g , w h i c hh a st h ev i r t u e so fn o n - d e s t r u c t i v e ，n o n - c o n m c ta n dh i g h r e s o l u t i o n w i t l lt h ew i d ea p p l i c a t i o no f c t r e c o n s t r u c t i o na l g o r i t h mh a sb e c o m et h ef o c u so f m a n ys t u d i e s 1 1 1 er e s e a r c h e si nt h i sp a p e rf a l li n t ot w op a r t s ：i n c o m p l e t ep r o j e c t i o nr e c o n s t r u c t i o na n d m u l t i - a r t ( a l g o r i t h mr e c o n s t r u c t i o nt e c h n i q u e ) a tp r e s e n t ，m o s to ft h ec ts y s t e m sn e e dc o m p l e t ep r o j e c t i o n so fd a t a , b u ti ti sd i f f i c u l t t oo b t a i nc o m p l e t ep r o j e c t i o n sf o rm a n yr e a s o n s i nv i e wo f t h i sp o i n t ，t h i sp a p e rs t u d i e st h r e e i n c o m p l e t ed a t ap r o b l e m s ( t h ee x t e r i o rp r o b l e m ，t h ei n t e r i o rp r o b l e ma n dt h el i m i t e da n g l e p r o b l e m ) o fc o m p u t e dt o m o g r a p h y 1 1 1 ei n c o m p l e t ed a t ap r o b l e m sm a yb ec a u s e db ym a n y f a c t o r s ，s u c ha s ，t h ee q u i p m e n ta n dt h ed e n s i t yr e g i o no fd e t e c t e do b j e c t s t h e r e f o r e , i ti s e s s e n t i a lt oe x p l o r ean e wu n i v e r s a lm e t h o da p p r o a c h i n gt ot h ei n c o m p l e t ed a t ap r o b l e m s t h i sp a p e ri n t r o d u c e san e wm e t h o do fc o m p l e t i o nw h i c hb u i l d sd y n a m i cd i f f e r e n t i a l e q u a t i o n st op r e d i c to t h e ru n k n o w n d a t ab a s e do nt h eg r a ys y s t e mt h e o r y , a n dw ec a ne x p l o i t f b pt os o l v et h ei n c o m p l e t ed a t ap r o b l e m s c o m p a r e dt oo t h e rl o c a lt o m o g r a p h ya l g o r i t h m s ， t h i sm e t h o dh a sr e d u c e dt h ea m o u n to fr a d i a t i o n m o r e o v e r , i t sm o r ee f f i c i e n t l yt os o l v e a b o v ep r o b l e i n s h e n c et h em e t h o do fc o m p l e t i o nb a s e do ng m ( 1 ，1 ) m o d e lc a nb ew i d e l y a p p l i e di np r a c t i c e n eo t h e rw o r ko ft h ep a p e ri sm u l t i a r tr e c o n s t r u c t i o no nt h el o c a la r e an e t w o r k ( l a n ) a st h ec o s to f a r t i sv e r yb i g , t oo b t a i na c c e p t a b l er e c o n s t r u c t i o nt i m e , o p t i m i z a t i o n i sv e r yn e c e s s a r y b u tt h e r ea r es o m el i m i t s an e wm e t h o db a s e do nt h el a ni sp r e s e n t e dt o i m p r o v et h es p e e di nt h i sp a p e r t h es h a r e db u f f e r ，i nw h i c ht h ed a t af i l ei su s e d ，i sc r e a t e di n t h em a s t e rc o m p u t e r , a n dt h em a s t e rc o m p u t e ra n dt h es l a v ec o m p u t e rc a l le x c h a n g e i n f o r m a t i o n a f t e rt h es l a v ec o m p u t e r ss e n dt h er e s u l t st ot h es h a r e db u f f e r , t h ef i n a l r e c o n s t r u c t i o ni sg o tb yt h em a s t e rc o m p u t e ru s i n gd i f f e r e n ti m a g ef u s i o na l g o r i t h m s i t s s i m p l ea n de a s yt oi m p l e m e n t ，a n dt h er e c o n s t r u c t i o nt i m ec a n b er e d u c e dg r e a t l y k e yw o r d s ：c tr e c o n s t r u c t i o n ；i n c o m p l e t ep r o j e c t i o n ；g r a ys y s t e m ；p a r a l l e la r t ；i m a g e f u s i o n 独创性说明 作者郑重声明：本硕士学位论文是我个人在导师指导下进行的研究工 作及取得研究成果。尽我所知，除了文中特别加以标注和致谢的地方外， 论文中不包含其他人已经发表或撰写的研究成果，也不包含为获得大连理 工大学或者其他单位的学位或证书所使用过的材料。与我一同工作的同志 对本研究所做的贡献均已在论文中做了明确的说明并表示了谢意。 作者签名： 训7 1 2 。衫 大连理工大学硕士研究生学位论文 大连理工大学学位论文版权使用授权书 本学位论文作者及指导教师完全了解“大连理工大学硕士、博士学位 论文版权使用规定”，同意大连理工大学保留并向国家有关部门或机构送 交学位论文的复印件和电子版，允许论文被查阅和借阅。本人授权大连理 工大学可以将本学位论文的全部或部分内容编入有关数据库进行检索，也 可采用影印、缩印或扫描等复制手段保存和汇编学位论文。 作者签名：塾煎 导师签名一了嘞 盟年上月互日 大连理工大学硕士学位论文 1 绪论 1 1c t 的产生和发展 计算机断层成像技术【l 】( c o m p u t e r i z e dt o m o g r a p h y ，美国放射学界统称为c o m p u t e d t o m o g r a p h y ，简称c t ) 是通过对物体进行不同角度的射线投影测量而获取物体横截面信 息的成像技术，涉及到放射物理学、数学、计算机图形图像学和机械学等多个学科领域。 c t 技术从理论上讲是一个从投影重建图像的理论，有其普遍性，在数学界已经引起了 广泛的重视。计算机断层扫描作为一种技术，既有坚实的数学理论为依据，又有现代计 算机技术相支持，使其在其它领域得到广泛的应用。 1 8 9 6 年，伦琴在n a t u r e 杂志上发表了关于x 射线的第一篇论文，并因此成为诺贝 尔物理学奖的第一位获得者。2 0 世纪初，x 射线透射成像开始应用于临床医学诊断。在 1 9 1 7 年，奥地利数学家j r a d o n 在其发表的论文中提出了c t 图像重建的基本数学理论， 他指出：任何物体均可用无限多个投影来表示；反之，如果知道无限多个投影，便可重 建该物体对象【2 】。在该文中，r a d o n 系统地论证了由积分值( 即投影) 确定被积函数( 即重 建图像) 的整套理论方法，为c t 技术的形成和发展提供了可靠的理论依据p j ，但是限于 当时的技术条件而未能实现。直到1 9 5 7 年，美国物理学家a m c o r m a c k 教授在哈佛大 学电子加速实验室里，开始研究图像重建的数学理论，自五十年代以来，他发表了一系 列论文，不仅证明了医学领域中从x 射线投影数据重建图像的可能性，而且提出了相应 的实现方法并完成了仿真实验研究，并于1 9 6 3 年首先提出用断层的多个方向的投影重 建断层图像的代数计算方法【4 1 。1 9 7 2 年英国e m i 公司中央研究所工程师g n h o u n s f i e l d ，从实验技术角度解决了吸收值的问题，成功地研制了诊断头颅用的第一台 电子计算机x 射线断层摄影装置，并申请了英国专利。从此之后，人们才真正实现了人 体断层成像。h o u n s f i e l d 和c o r m a c k 两人因此而共同获得1 9 7 9 年诺贝尔医学奖。x c t 断层摄影装置的研究成功可以说是在x 射线被发现后的七八十年中放射医学领域中最 重要的突破性进展，它是二十世纪科学技术的重大成就之一。至此，c t 技术在医学上 的成功应用已轰动全球，其在成像方面无可争辩的优越性引起世人瞩目。1 9 8 0 年h e r m a n 教授在其专著中系统深入地阐述了c t 的理论基础1 5 】。 c t 技术自诞生以来，技术发展r 新月异，各个公司几乎每年都推出两种以上型号 的产品，由于c t 产品与临床应用相互结合、相互促进，因此c t 产品不断进行着技术 更新。下面我们就简要回顾一下c t 技术的发展过程。 基于灰色系统的c t 不完备投影重建算法的研究 传统意义上的c t 分为五代 6 1 。1 9 7 1 年，在a t k i n s o nm o d e ) , 医院，安装了第一台头 部c t ；1 9 7 3 年出现第一台全身c t ( m a r k1e m i ) ，是传统意义上的第一代c t ；1 9 7 4 年 在o h i on u c l e a rd e l t a 医院安装了传统意义上的第二代c t ；同年，在a r t o n i x 医院安装 了传统意义上的第三代c t ；1 9 8 5 年出现了第一台螺旋c t ；1 9 9 8 年，出现了第一台真 正的多排c t 。 从扫描方式角度可以分为：第一代平行束平移加旋转扫描c t ；第二代扇形束平移 加旋转扫描c t ；第三代扇形束旋转扫描c t ：第四代锥形束扫描体积c t ；第五代螺旋 扫描c t 【“。前三代c t 均存在着扫描速度慢、重建效率低、射线剂量大、设备复杂等缺 点。而第四代体积c t 采用锥形束旋转扫描方式，对物体进行一周扫描得到投影数据后 即可重建该物体的所有断层图像，具有扫描速度快、重建效率高等优点。自滑环技术的 引入，第五代螺旋c t 诞生了，其包含两方面的基本内容：一是x 射线管及探测器连续 3 6 0 度旋转；二是病人同时随检查床匀速推进。螺旋c t 已成为2 1 世纪c t 的主流，是 c t 发展史上的里程碑。 在c t 技术中，重建算法一直处于非常重要的地位。不同的算法不仅会产生不同重 建质量的图像，而且它们对计算量和数据量的要求也有很大不同。在过去的三十年里， 尽管x 射线c t 成像技术得到了长足的发展，但直到现在c t 重建算法仍然是一个非常 重要的研究领域。c t 重建算法大致可分为以下两类。 变换法( 即解析法) ：该方法以r a d o n 变换理论为基础，可分为直接r a d o n 逆变换法 和滤波反投影法。前者以s m i t h c s 】和g r a n g e a t 分别提出的算法最为著名，后者则以 f d d k a r n p 、d a v i d 和k r e s s 9 z 人联合提出的f d k 算法最为著名。 迭代法：该方法通过对初始假定图像进行数学迭代运算重建图像，包括代数迭代算 法【1 0 1 、期望最大化算法( e m ) 、同时代数迭代算法等。 滤波反投影算法是c t 重建方法中最重要的一种近似重建算法，它具有重建精度高， 重建速度较快等优点。因而在c t 设备中得到了广泛应用。在实际应用中，由于各种原因 可能造成投影数据的不完备，在这种情况下，直接利用传统的解析重建方法( f b p ) 或者 统计重建方法( e m ) 都不能获得理想的重建结果，要克服投影数据缺失为重建算法带来 的局限，就必须对缺失的投影数据进行补偿，使其完备化。 1 2 不完备投影数据c t 重建及迭代算法的背景及意义和发展研究情况 在c t 的检测过程中，成像结构与待测物体的不理想均会产生投影数据不完备的问 题，而投影数据的不完备包含多种情况，其中研究比较多的有三种：内部问题、外部问 大连理工大学硕士学位论文 题和有限角度问题。例如，很多情况下我们为了减少x 射线剂量仅对感兴趣的区域进行 扫描，或者在检测火箭外壳时，x 射线不能够穿过整个火箭断层，仅需对火箭的外层进 行扫描，或者勘探中由于受到地理条件的限制，不能够获取所有角度下的投影数据。 内部问题和外部问题主要是划分出一部分区域，仅对该区域部分感兴趣，称为局部 重建，也叫做感兴趣区域( r o i ) 重建，这方面的研究比较多，kt s m i t h 1 于1 9 8 4 年首 先提出了l a m d at o m o g r a p h y 算法。这也是研究局部重建领域的第一种算法，该算法的 优点是具有严格的局部性，适用于重建任意维的空间信息，缺点是重建图像的质量仅在 密度函数不连续处较好，其他区域较差。a d df a r i d a n i ，e r i kl r i t m a n 和k e n n a n t s m i t h 1 2 - 1 4 1 等于1 9 9 2 年对l a m d a t o m o g r a p h y 算法进行了详细的推导和计算。a n d r e a s r i e d e r ，r a i n e rd i e t z ，t h o m a ss c h u s t e r ”1 于2 0 0 0 年提出一种近似的转置方法，其核心思 想就是首先对密度函数厂进行变形，以获得更好的局部性，然后进行局部重建。f a r r o k h r a s h i d f a 玎o l d 一嘲，j o ed e s t e f a n o ，t i mo l s o n t l 7 1 ，f r a s h i d f a r r o k h i ，k j r “u ， c a b e r e n s t e i n 【”】等分别对雷登数据利用小波变换进行局部化，以获得局部信息的重建结 果。 解决有限角度问题的方法是角度缺失重建，又叫有限角度c t 重建。如果扫描角度 小于1 8 0 度，则对其投影数据的重建就属于有限角度重建。 有限角度图像重建算法很多，大致可分为两类：基于变换的迭代一解析重建算法和 基于级数展开的迭代一代数，统计重建算法【”之o l 。在有限角度重建中，约束条件、先验知 识以及正则化的引入非常重要，比如图像非负有界，投影数据对称等。i r r ( i t e r a t i v e r e c o n s l m c t i o nr e p r o j e e t i o n ) 方法1 2 1 - 2 4 1 和g p ( g e r c h b e r g p a p o u l i s ) 方法1 2 5 - 2 8 1 是最常用的迭 代重建方法。1 9 8 2 年，m n a s s i 和m e d o f f 提出了用i r r 算法对心脏成像【2 1 - 2 2 1 ，1 9 8 5 年，k i m 提出可以直接在投影数据间迭代p s i r r ( p r o j e c f i o ns p a c em r ) 。a r t 、s i r t 和 p w l s 等代数迭代方法 2 9 1 与m e 、e m 、m l 以及b a y e s i a n 估计等优化准则例，以及基 于图像平滑、边缘保持等的正则化条件都可以应用到有限角度的图像重建中【3 1 1 。 现有的c t 重建算法主要可分为解析重建法和迭代重建法。其中迭代算法的缺点主 要是计算量大，收敛速度慢，这就对计算机的处理能力有了非常高的要求，近些年来计算 机硬件有了突飞猛进的发展，但还不能够完全满足人们对图像重建速度的要求。因此， 如何减少迭代算法的重建时间成为研究的重点。王小璞，张朋等在文献 3 2 1 q b 将图像逐 级分块重建，张顺利等人在文献1 3 3 q b 提出了一种快速、实时的投影系数计算方法，在 保证图像质量不受损失的酊提下取得了6 倍以上的重建加速比。此外，并行算法的引入 也大大减少了重建时问，很多学者将不同的并行方法引入c t 当中，在文献 3 4 1 q b 应用 v l s i ( v e r y - l a r g e s c a l ei n t e g r a t i o n ) 平台实现并行算法。文献 3 5 1 中l a t t a r d 提出一种基于 基于灰色系统的c r 不完备投影重建算法的研究 多处理器的计算机来进行重建，每个处理器对其中一些数据进行计算然后与其它的处理 器交换处理结果。在m i m d ( m u l t i p l e i n s t r u c t i o ns t r e a m ，m u l t i p l e d a t as t r e a m ) 计算机上， 并行算法得到很好的应用，c l a u r e n t 等人【3 6 1 提出了两种并行化的方法，在重建效率上 有明显提高。 1 3 本文的主要工作 本文的工作主要分为两部分：基于灰色系统的c t 不完备投影重建和局域网迭代融 合方法的研究。 ( 1 ) 本文主要研究了不完备投影数据重建的三个重要方向：内部问题、外部问题和 有限角度问题。首先采用灰色系统理论建立一阶灰微分方程，对缺失数据进行灰色预测， 然后采用滤波反投影重建算法进行重建。并给出了建模方法和仿真结果，理论分析和重 建的仿真实验说明了采取本文的方法能够有效解决上述问题。 ( 2 ) 本文的另一个工作就是研究局域网内多机迭代融合方法。文中提出一种提高 a r t 重建速度的方法，利用操作系统提供的数据共享传输功能，通过在主机上开辟共享 缓存区，在其中设立数据文件，其他从机分别调用共享区的数据进行重建，并将重建结 果传输给主机，最后由主机对来自不同从机的结果选取不同的融合方法进行图像融合， 得到理想的c t 断层重建结果。该方法原理简单，易于实现，重建时间成倍减少，具有 很高的性价比。 1 4 章节安排 本文共分五章：第一章为绪论，简单介绍了c t 检测技术的应用背景和发展，不完 备投影数据c t 重建及迭代融合算法的背景、意义和发展研究情况以及本文所做的工作。 第二章首先介绍了关于c t 技术理论及c t 检测技术的基础知识，为后文做理论上的铺 垫，然后介绍了c t 重建的两种基本算法：滤波反投影重建算法和迭代重建算法。第三 章首先介绍了灰色系统理论的基础知识，然后给出灰微分方程的建模过程，最后利用 g m ( 1 ，1 ) 模型对缺失投影数据进行预测，用来解决内部问题、外部问题和有限角度问题， 并给出仿真条件和计算机模拟结果。第四章主要介绍了a r t 算法通过局域网多机实现 并采取图像融合技术获得最终重建结果的方法，而后给出了实验结果及性能分析。最后， 在第五章对本文内容进行了总结，并提出了下一步工作的努力方向。 大连理工大学硕士学位论文 2c t 图像重建的基本原理 2 1c t 成像原理 c t 技术的物理原理是基于射线与物质的相互作用原理。当射线( t 和x 射线) 穿越 物体时，由于产生光电效应、康普顿效应及电子对效应等物理过程，射线( 即入射光子) 将被物质吸收，使得射线强度发生衰减。如图2 1 所示，一束x 射线沿直线穿过一个物 体并被物体各部分吸收，强度被衰减的x 射线透过物体，最终被探测器记录下来。如果 x 射线覆盖了整个样品断层，就可以得到该投影角度下，样品在扫描平面内的一维投影。 再将样品相对x 射线和探测器进行旋转，就可以得到各个方向上经衰减的x 射线信号， 这些投影数据反映了物体在这一截面不同位置对x 射线吸收系数的分布情况，如果我们 利用一定的数学方法进行运算，就可以得到物体这个截面上的吸收系数分布函数。随后 我们再将物体相对于x 射线和探测器沿着与扫描截面相垂直的方向移动一个位移，重新 进行上面的步骤，就可以得到物体另一个断层截面上的吸收系数分布。通过一系列这样 的扫描，就可以获得整个物体各个部分的吸收系数分布，从而获得物体吸收系数的三维 分布图。 x 射 线 物体 探测器 图2 1x 射线c t 成像原理 f i g 2 1x - r a yc ti m a g i n gp r i n c i p l e 具体的射线强度衰减3 _ 7 】情况如图2 2 所示，假设入射强度为s o 的单色射线穿过厚度 为1 、线性衰减系数为z 的均匀物体后，探测器检测到的射线强度，服从比尔定理( b e e r s l a w ) ，即 i = l o e 一肛( 2 1 ) 或肛= i n ( s o s ) ( 2 2 ) 基丁灰色系统的c t 不完备投影重建算法的研究 ，dr 习， 惨 陵= 五。匕翻 图2 2 比尔定理 f i g 2 2 b e e r sl a w 在实际的射线检测中，沿射线束路径方向上的物体密度是不均匀的，比尔定理可写成： ，：如一脚埘 ( 2 3 ) 经变换可得： i ( z ) d l = i n f i ( 2 4 ) 其中，( z ) 是线性衰减系数沿射线路径变化的函数，反映了被测物体内部的物质密度 分布情况，因此它在空间上就形成了物体内部的图像，( f ) 实际上反映了图像的灰度分 布。在c t 检测中，j 。是探测器测得的射线原始能量，是射线经被测物体衰减后探测 器测得的出射能量。由此，可定义 州) = l n 手= i 州) 讲 ( 2 5 ) 这里将p ( ，) 称为投影数据。因此，图像重建问题的实质就可以归结为从投影数据 p ( ，) 来计算函数( ，) 值的反演问题。 2 2c t 重建算法的分类 在c t 技术发展的数十年中，人们提出了各种各样的重建算法，这些算法大致分为 两类：解析法和迭代法。解析法以r a d o n 变换理论为基础。是对横截剖面的直接数学反 计算，包括滤波反投影算法( f b p ) 、卷积反投影算法( c b p ) 和r a d o n 逆变换算法。迭代 法是通过对初始假定图像进行数学迭代运算来重建图像，常用的有代数重建法( a r t ) 、 极大熵迭代算法和根据估计理论得到的期望最大化重建算法如m l e m 、o s e m 、 b a y e s i a n 重建算法等。各种方法有着自己的优缺点，其中滤波反投影算法简单、重建速 度较快，是目前c t 重建技术中最常用的算法，其应用非常广泛，但该算法要求采集投 影角度必须完全，而且该算法重建图像中易产生星状伪迹，这在一定程度上影响了其图 像重建质量。在实际应用中如果采集的投影数据不完备，用解析法重建往往得不到满意 的结果，而迭代算法可以弥补这一缺陷，通常迭代算法可以采用较少的投影重建出较高 6 大连理工大学硕士学位论文 质量的图像。当c t 系统在数据采集过程中受到严重的噪声干扰，用确定性重建方法也 将得不到满意的重建结果，而采用以估计理论为基础的e m 迭代方法进行重建，将使图 像重建质量有很大改善。不过迭代方法重建质量虽好，但运算量较大，计算速度较慢。 2 ，3 滤波反投影重建算法 2 3 1 反投影重建算法 反投影重建算法的基本思想： 断层平面中某一点的密度值可以看作是这一平面内所有经过该点的射线的投影值 之和( 或其平均值) 。整幅重建图像可以看作是所有方向下的投影累加而成。 因而反投影算法的基本做法是把每次测得的投影数据按“原路”反投影到投影射线经 过的各个像素上。也就是说，指定投影射线上所有各点的值等于测得的投影值。 为了更好的说明反投影重建算法的物理概念，特举例如下： 图2 3 断层的像素值和射线 f i g 2 3r a ya n dp i x e lv a l u eo f c r o s ss e c t i o n 图2 - 3 为一与纸平行的断层，具有四个像素，像素值( 代表密度) 分别为而，而，而，- 。 赋值如下：五= 5 ，屯= 0 ，而= 2 ，= 1 8 。射线标号示于图中根据射线 ( ( 1 ) ，( 2 ) ，( 3 ) ，( 4 ) ，( 5 ) ，( 6 ) ，( 7 ) ) 投影的定义有： p 1 2 舶+ x 2 = 5 p 2 2 x 3 + 斟= 2 0 p 3 2 触+ 石3 = 7 p 4 2 x 2 + x 4 = 1 8 p s 。x 3 = 2 p 6 = 朋+ 工4 = 2 3 p 7 = x 2 = 0 根据反投影重建算法的基本思想，重建得到图像中各像素值如下： 基于灰色系统的c r 不完备投影重建算法的研究 而= p l + p 3 + p 6 = 3 5 x 2 = p l + p 4 + p 7 = 2 3 b = p 2 + p 3 + p 5 = 2 9 工4 = p 2 + p 4 + p 6 = 6 1 重建后图像如图2 4 ( b ) 所示。可以看到原图中像素不为零的点反投影重建后比较突 出。但原图中为零的像素，经反投影后不再为零。即有伪迹。为了使重建后图像的像素 值更接近于原图的像素值，在求反投影时，把数据除以投影数目n p ，这样有； 吆= 去和 c z e ， 该式可作为反投影重建算法的计算式。其中表示像素七的值，p k ，表示经过像素七 的第f 条射线投影，n ，表示图像内的射线条数。按式( 2 6 ) 重建后图像如图2 4 ( c ) 所示。 圈圈圈 ( a ) 原图( b ) 反投影重建图( c ) 求平均后图像 图2 4 反投影重建示例 f i g 2 4e x a m p l eo f b a c k - p m j e c d o nr e c o n s t r u c t i o n 从反投影重建图像可以看到，原图中像素值为零的点经反投影重建后其像素值不等 于零，而是等于一个相对较小的数值。在整幅图像的重建中，这种现象表现为灰雾背景， 被称为星状伪迹。产生星状伪迹的原因在于：反投影的本质是把取自有限物体空间的投 影均匀的回抹( 反投影) 到射线所及的无限空间的各点上，包括原先像素值为零的点。 反投影重建算法的“取投影”，“反投影重建”这些环节可以看作一个以原图为输入， 重建后图像为输出的成像系统，该系统的点扩展函数( p s f ) 可以表示为： 坶2 崭南 亿z ， 其中( ，口) 为极坐标表示，o ，y ) 为直角坐标表示，由式( 2 7 ) 可见，相应于反投影重 建算法的系统，它的点扩展函数不是占函数，系统不是完美的。虽然在，= 0 处能反映原 图是点源的情况，但在r 0 处，像素值不为零，式( 2 7 ) 定量地描述了反投影重建算法 星状伪迹的本质。 一8 一 大连理工大学硕士学位论文 若原像为y ) ，则将原像取投影后再按反投影算法得到的重建图像为： ，( 五y ) = ( 五y ) ( 五，) ( 2 8 ) 要去除反投影算法的星状伪迹，可以在输出端加一滤波器，使加了滤波器后的反投 影重建等效成像系统的点扩展函数为8 ( x ，y ) 。设滤波器的p s f 为q ( x ，y ) ，相应的传递函 数为烈参叮) 。这里q ( 善，玎) = 最【g ( x ，y ) 】，要求： 1 了告+ g ( x ，y ) = 8 ( x , y ) ( 2 9 ) 露n x + y 做二维傅立叶变换，得： q ( f ，7 ) = 石孝2 + ，7 2 = 巧口 ( 2 1 0 ) 这是一个二维滤波器，实现起来较麻烦。若p 的变化范围扩至。o ，则根本不能实现。 但它提供了去除星状伪迹的一个方向。 2 3 2 滤波反投影重建算法 反投影重建算法的缺点是引入星状伪迹，即原来图像中密度为零的点，重建后不一 定为零，从而使图像失真。去除伪迹的方法通常有以下两种： 箍龛黔二伍豳 屯垂亟园一巨亘垂亘囝一。篙袁冒像) ( a ) 毪龛黔匝匣卜巨亟哥 豆受垂习一( 篙豢豁 ( b ) 图2 5 去除星状伪迹的方案 f i g 2 5 n es c h eo f e 蚍t m gs t a rm i f a c t s 由于二维滤波器实现麻烦，实际应用中通常采用方案( b ) ，这一方案被称为滤波反 投影重建算法。 j i 弧 基于灰色系统的c t 不完备投影重建算法的研究 图2 6 滤波反投影算法所用的坐标系 f i g 2 6 t h ec o o r d i n a t eu s i n gb yf i l t e rb a c k - p r o j e c t i o na l g o r i t h m 图2 6 中x y 为固定坐标系；x ，一) ，为旋转坐标系，它绕原点转动，使得投影线总 是沿着y ，方向；( r , o ) 为极坐标。x r y ，坐标的原点与x - y 坐标的原点重合，两者的夹 角为，不同的代表不同的投影视角，射线位置可由坐标( ，矿) 完全确定。空间任一 点的位置可用( 工，y ) ，( x r ，y r ) ，或用极坐标( ，p ) 表示。投影数据用p ( ，妒) 表示。 滤波公式为： g ( ，妒) = p ( ，) + ( ) ( 2 1 1 ) 公式( 2 1 1 ) 实现对投影数据的滤波，这里用g ( t ，) 表示滤波后的投影数据，式中 | i l ( _ ) 为一维滤波器，该滤波器的傅立叶变换为l d 。 反投影重建公式为： a ( r ，口) = f g 【，c o s ( o 一矿) ，矿 却 ( 2 1 2 ) 式中g r c o s ( o 一声) ，庐 的物理意义是投影p ( x ，庐) 经过传递函数为i p i 的滤波器后得 到的修正后投影g ( ，) 在满足x r = r e o s ( o 一妒) 时的值。 公式( 2 1 2 ) 称为滤波反投影方程，其物理意义是经过给定点( ，口) 的所有滤波后的投 影在= 0 石范围内的累加反投影重建，得出( r ，0 ) 点的像素值。 滤波( 卷积) 反投影算法包含以下三大步： ( 1 ) 把在固定视角以下测得的投影p ( ，矿) 经过滤波，得到滤波后的投影g ( ，) ； ( 2 ) 对每一个或，把g ( ，办) 反投影于满足= r c o s ( o - # 1 ) 的射线上的所有各点 ( r ，口) ； ( 3 ) 将步骤( 2 ) 中的反投影值对所有0 矿丌进行累加( 积分) ，得到重建后的图像。 2 4 迭代算法 2 4 1 迭代算法的基本思想 迭代重建算法的思想与变换法出发点不同，迭代重建算法一开始就把连续的图像， 离散化，把整个图像区域划分为j = n 一有限个像素，每个像素内部的灰度值为常数。 射线的定义与变换法稍有不同，这里把射线看作宽为f 的粗线( 如图2 7 所示) ，这条粗线 与所经过的像素的重叠区域就是该像素对这条射线投影的贡献3 0 l 【3 1 1 。图2 7 中，j 号像 1 0 大连理工大学硕士学位论文 素的像素值为_ ，f 号射线与，号像素重叠的面积为鲋曰c 的面积，它与像素面积6 2 之 比为0 = a a b c 6 2s l 称为加权因子。 - 2 等m x 图2 7 迭代重建算法所用的射线 f i g 2 7 t h e r a yu s e di ni t c r a t i v ea l g o r i t h m 于是，像素_ ，对射线i 的射线投影的贡献为 p f = 白o ( 2 1 3 ) i 号射线还经过其他像素，其总的射线投影为 p l = p 口= 白_ ( 2 1 4 ) 把p ；称为( 第i 号射线的) 射线和。 第二种表示射线和及定义屹的方法是认为像素值集中在像素的中心。此时 ，。；f l 得射线( 宽为订经过号像素中心7( 2 铲，具i d - 他i l 。、15)io _ q g - - 种定义嘞的方法相似的第三种定义0 的方法( 图2 8 ) 是令各射线的宽度为 零。它们的间距为f 。于是 ：j 1 ，堂碧线通过丹像素内的一点； ( 2 1 6 ) q lo ，其他。 不难看出，射线宽度f 及吒的定义将影响式( 2 1 4 ) 中射线和的值，从而会影响重建 后图像的精度。采用第三种方法来选取勺并求射线和，可得： p f = t a x ，i = 1 , 2 ，8 ( 2 1 7 ) 下oi譬上 基丁：灰色系统的c r 不完备投影重建算法的研究 或用矩阵表示为 图2 8 利用第三种定义的方法求射线和 f i g 2 8r a y - s u mb yt h et h i r dd e f i n i t i o n p = r x 其中p = b 见热r ，x = k 而b r r = r t ir t 2 r 1 9 r 2 1r 2 2 r 2 9 ir s 2 r s 9 ( 2 。1 8 ) ( 2 1 9 ) 把其推广至，个像素，条射线的一般情况。此时： x 为，维矢量，称图像矢量，p 为维矢量，称测量矢量，r 为i ，维矩阵，称投 影矩阵。现在的任务是根据测量得到的p 和已知的矩阵r 求图像矢量x 。 实际应用中，测量误差不可避免，噪声影响也不可忽视，应对方程2 1 8 作如下修正 p = r x + 。 ( 2 2 0 ) 这里，口为误差矢量。 因此，从投影重建图像问题可归结为：根据测量矢量p ，利用式( 2 2 0 ) 估计图像矢 量工。 2 4 2a r r 算法 a r t 重建算法的步骤是： 先假设一初始图像x ( “，根据x ( o 求近似图像x ( n ，再根据x ( 1 求近似图像工( ”，如 此继续，直至满足预定条件后结束。在根据x 求x ( t + l 时需加一校正值出( ”。出m 只 大连理 ：大学硕士学位论文 设矢量工的维数为j ，一；为给定的，维矢量；q j 为实数。欲求满足下列方程的x ： + ”乙气q ；k - - i p t x ( g ) 。，警 眨2 3 ) x ( “”2 b “t ) 广其他 。 l1 1 2 “ 。 式中，k 为迭代序号，= t ( m o d ，) + 1 ，即乇= l ， = 2 ，“= ，= l ，“= 2 ，。名”称松 只= 扛i n j 苫= q t ， ( 2 2 5 ) 一“1 ) 一l 。：。一，r 。御 r i k x ( k ) 2 两式相减得石0 ( _ j ) 一x o ( j 一1 ) = 口( x o ( 七一1 ) 一工1 ( | i 一2 ) ) 整理得 工o ( i ) = x o ( _ j 一1 ) + 口工1 ( 七一1 ) 一g x ( 1 ) ( _ j 一2 ) 将( 后一2 ) = 算1 ( 七一1 ) 一x o ( | j 一1 ) 代入上式得 x o ( 后) = ( 1 一a ) x o ( | j 一1 ) 方法3 ：将g m ( 1 ，1 ) 模型x o ( 露) + 1 ( 七) = b 可以转化为 川七) = i l + - 0 0 5 5 口( z x ( o ) 一1 ) ，| j = 3 ，聘 一o ( 2 ) = 一口o ( 1 ) 其中= 上l + 0 5 a ，口= 志 方法4 ：g m ( 1 ，1 ) 模型工o ( 七) + 船1 ( 后) = 6 可以转化为 x i 。( t ) = ：x i 1 i ) ( j k i - 二j 1 ) i - ：而0 5 b x 。( 一1 ) ，t = 3 ，h o ( 2 ) = 一盯算o ( 1 ) 其中= 丽b 五，口=