《医学图像的三维重建》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
电子科技大学毕业设计(论文)任务书拟题单位计算机科学与工程学院审题人陈雷霆题目全称医学图像的三维重建主要任务:本文主要研究由医学图像(CT、MR1图像)二维断层序列来构建组织或器官的三维几何模型的技术。论文主要工作包括如下医学图像三维重建的关键技术:图像的输入与预处理,组织或器官的分割与提取,由二维轮廓线重构三维形体,基于规则体数据的三维表面重建。并主要分析MarchingCube算法抽取等值面的基本原理及实现方法,分析了MC算法中连接方式产生二义性的原因及消除二义性的方法。起止时间:2005年2月28日至2005年6月10日学生姓名张剑专业计算机科学与技术(软件技术)班次21082010班指导单位计算机科学与工程学院指导教师陈雷霆设计地点电子科技大学计算机学院成果形式论文2005年6月12日目录
1摘要:IABSTRACT:I第一章引言11.1课题研究背景及意义11.2国内外研究现状11.3本文的主要工作2第二章医学图像三维重建技术综述32.1弓1;;32.2各类三维可视化方法的典型算法42.1.1表面绘制方法42.1.2体绘制方法52.1.3混合绘制方法62.3各类算法的特点与应用7第三章医学图像的预处理与分割93.1概述93.2医学图像的获取及其数据类型93.2.1三维空间数据场的数据类型93.2.2医学图像的获取与输入113.3医学闇像的预处理123.3.1医学图像滤波123.3.2二维断层图像之间的插值183.4医学图像的分割193.1.1体数据的二值化203.1.2种子填充法进行组织提取223.1.3分割结果及分析233.5小结24第四章基于面绘制的医学图像三维重建264.1引言264.2切片级重建方法264.1.1轮廓的提取274.1.1.1梯度算子274.1.1.2Roberts算子294.1.1.3Sobel算子294.1.2平面轮廓重建314.3体素级重建方法324.1.1体素模型与等值面定义334.1.1.1体素模型334.1.1.2等值面定义354.1.2等值面的构造方法综述364.1.3移动立方体(MarchingCubes)算法抽取等值面374.1.3.1MC算法的基本原理374.1.3.2MC算法生成等值面的近似表达与连接的二义性434.1.3.3MC算法二义性的消除45
24.4MarchingCubes算法的实现过程474.1.1开发トロ及方法474.1.2体数据处理理程474.1.3算法具体实现484.1.4算法流程图504.1.5MarchingCubes算法主要功能模块504.5小结51第五章总结与展望525.1全文工作总结525.2课题研究发展方向53参考文献54英文原稿55中文翻译58致谢60
3医学图像的三维重建作者:张剑班级:21082010班指导老师:陈雷霆指导单位:计算机科学与工程学院摘要:医学图像三维可视化是目前的ー个研究热点问题,是一个多学科交叉的研究领域,是计算机图形学和图象处理在生物医学工程中的重要应用。它涉及数字图象处理、计算机图形学以及医学领域的相关知识。医学图像三可视化在诊断医学、手术规划及模拟仿真、整形及假肢外科、放射治疗规划、解剖教学等方面都有重要应用。因此,对医学图像三维可视化的研究,具有重要的学术意义和应用价值。医学图像三维可视化的主要研究内容包括医学图像的预处理,如插值、滤波等:组织或器官的分割与提取;面绘制与体绘制等。本文对医学CT,MRI图像三维可视化的关键技术进行了研究,研究了三维医学体数据的分割和体绘制算法。对组织或器官的分割与提取是保证重建模型准确表达其相应组织器官的前提。在三维医学体数据分割的研究中,本文将一些二维的图像处理方法,如边缘检测,区域生长等,推广到三维情况。关键词;三维医学图像,分割,三维可视化,模糊区域生长ABSTRACT:3Dvisualizationfrommedicalimagesisamulti-disciplinarysubject.Itisanimportantapplicationofcomputergraphicsandimageprocessinginbiomedicineengineering.Itrelatestothesubjectsofdigitalimageprocessing,computergraphicsandsomerelatedknowledgeofmedical.3Dvisualizationofmedicalimagesarewidelyusedindiagnostic,surgeryplanningandsimulating,plasticandartificiallimbsurgery,
4radiotherapyplanning,andteachinginanatomy.Studyon3Dreconstructionfrommedicalimageshasimportantsignificanceonscienceandworthinessinpracticalapplication.Themainresearchcontentsof3Dvisualizationfrommedicalimagesincludeimagepre-processing,suchasinterpolatingandfiltering,segmentingandextractingtissuesororgansofbody,surfacerenderingandvolumerendering.Inthisdissertation,keytechniquesfor3Dvisualizationfrommedicalimages(suchasCT,MRIimages)arestudied.Segmentingandextractingtissuesororgansfrommedicalimagesarepremisesof3Dvisualizationaccurately.Thesegmentationof3Ddatasetsisoftendifficult.Byextendingsomealgorithmstothreedimensions,weproposeanewhybridmethodfor3Dsegmentationofthewholebrainbasedonfuzzyregiongrowing,edgedetectionandmathematicalmorphology.Thethresholdcontrollingtheprocessofregiongrowingisdeterminedbyafuzzytechnique.Keywords:3Dmedicalimage,segmentation,3Dvisualization,fuzzyregiongrowing
5第一章引言1.1课题研究背景及意义现代各种医学成像技术如X射线断层成像(CT),磁共振成像(MRI),超声成像(US),单光子发射断层成像(SPECT)和正电子发射断层成像(PET),为人体不同部位提供了多种不同时空分辨率的二维医学图像,形成了现代医学影像诊断学的基础,极大的推动了医学诊断和临床治疗的发展。随着现代计算机辅助技术、虚拟现实等现代技术的发展,人们对直观的三维医学图像的需求也日益增强,医学图像三维可视化近来成为医学成像发展的热点。但是,这些医疗仪器只能提供人体内部的二维图像,二维断层图像只是表达某ー界面的解剖信息,医生们只能凭经验由多幅二维图像去估计病灶的大小及形状,“构思”病灶与其周围组织的三维儿何关系,这就给治疗带来了困难。在放射治疗应用中,仅由二维断层图像上某些解剖部位进行简单的坐标叠加,不能给出准确的三维影像,造成病变定位的失真和畸变。为提高医疗诊断和治疗规划的准确性与科学性,通过二维断层图像来构建人体器官、软组织和病变体等的三维模型并进行三维显示,可以辅助医生对病变体和周围组织进行分析,极大地提高医疗诊断的准确性和科学性,从而提高医疗诊断水平。同时还在矫形手术、放射治疗、手术规划与模拟、解剖教育和医学研究中发挥着重要作用。医学图像三维重建与可视化技术就是在这ー背景下提出的,医学图像三维重建在医疗中的应用主要有在医疗诊断中的应用、在手术规划及放射治疗规划中的应用、在整形与假肢外科中的应用和在虚拟手术及解剖教育中的应用等。1.2国内外研究现状国外现在己经有一些医学图像可视化系统,如美国的ANALYZE系统、加拿大的ViewWand系统、荷兰的COVIRA系统等。这些系统是根据不同目的而开发的,所以各个系统的侧重点不同,使用方法也不尽相同。但一般都基于高档工作站,能提供多模图像数据的可视化分析,但不能在普通微机上运行。国内在医学图像三维重建及可视化研究方面,浙江大学、清华大学、东南大学、中科院自动化所等均作了大量研究,开发了医学实验系统,但目前国内尚无成熟的商用系统。
61.1本文的主要工作本文主要研究由医学图像(CT、MRI图像)二维断层序列来构建组织或器官的三维儿何模型的技术。论文主要工作包括如下医学图像三维重建的关键技术:图像的输入与预处理,组织或器官的分割与提取,由二维轮廓线重构三维形体,基于规则体数据的三维表面重建。
7第二章医学图像三维重建技术综述2.1引言进行三维重建的主要任务就是实现三维可视化显示、操作及分析为诊断和治疗目的而提供的医学图像数据。三维可视化显示关心的是在显示设备上如何绘制出具有真实感的人体组织结构;在操作方面,完成交互式显示组织结构的改变,如进行手术导航等;而分析是对人体组织结构进行形态或功能上的定量处理。从80年代开始,就已经有许多可视化的体绘制方法被提出并成功的应用于医学领域了。目前图像三维重建的方法主要有两大类:ー类是通过几何单元拼接拟合物体表面来描述物体三维结构的,称为基于表面的三维面绘制方法(SurfaceFitting),又称为间接绘制方法;另ー类是直接将体素投影到显示平面的方法,称为基于体数据的体绘制方法(DirectvolumeRendering),又称为直接绘制方法。其中,表面绘制方法是基于二维图像边缘或轮廓线提取,并借助传统图形学技术及硬件实现的,而体绘制方法则是直接应用视觉原理,通过对体数据重新采样来合成产生三维图像。体绘制更能反应真实的人体结构,但是由于体绘制算法运算量太大,即使利用高性能的计算机,仍然无法满足实际应用中交互操作的需要,因此面绘制仍是目前的主流算法。基于表面的三维面绘制方法是首先在三维空间数据场中构造出中间儿何图元,然后再由传统的计算机图形学技术实现画面绘制。这种方法构造出的可视化图形虽不能反映整个原始医学图像数据场的全貌及细节,但是可以对感兴趣的等值面产生清晰的图像,而且可以利用现有的图形硬件实现绘制功能,速度较快,因而得到了广泛的应用。本文主要用移动立方体法(MarchingCubes)构造等值面的方法来重建组织器官的三维图像,也属于基于表面的三维面绘制方法。为了论文的完整性,以下对基于体数据的直接体绘制方法也作些简单介绍。另外,这种传统的划分,使得某些既以绘制表面为目标,又采用了体绘制原理、或者既以反映数据整体信息为目标又以儿何造型作为显示单元的算法无从归属,从而造成了长期以来算法讨论中的混乱。因此,可以将这一部分算法归为第三大类,即混合绘制方法(如图2-1)〇
8医学体数据三维可视化方法.图2-1医学体数据三维可视化方法分类示意图2.1各类三维可视化方法的典型算法2.2.!表面绘制方法由于表面可以简洁地反映复杂物体的三维结构,因此在医学图像中边界面轮廓是用于描述器官的最重要特征。表面绘制是ー种普遍应用的三维显示技术,其首先是从体数据中抽取一系列相关表面,并用多边形拟合近似后,再通过传统的图形学算法显示出来。而表面的提取通常是通过门限设定,必要时结合手工描制完成的。表面绘制方法的处理过程主要包括下面三部分:(1)体数据中待显示物体表面的分割;(2)通过儿何单元内插形成物体表面;(3)通过照明、浓淡处理、纹理映射等图形学算法来显示有真实感的图像,并突出特定信息。表面绘制有多种算法,但从重建过程处理的基本元素的级别上来分,可以把这些方法分成两大类:体素级重建方法和切片级重建方法。其具体形式有两种:边界轮廓线表示和表面曲面表示。表面绘制方法中,连接轮廓线法,是最早被用来进行表面绘制的方法。它是首先将每层图像的轮廓提取出来,然后用以轮廓线点为顶点的三角形将每层的轮廓线连接起来,从而拼接出物体表面。这种方法占用存储少,速度快,便于进行三维实时旋转操作,而且可以纠正由于分类不当导致的错误结果。但该算法中两相邻层轮廓线对应点的确定和连接常常是难以解决的问题。
9除了以轮廓线表示物体外,还可以由轮廓重建物体的表面来表示。最早的方法是基于多边形技术,主要用平面轮廓的三角形算法,根据在不同切片图像上抽取出的ー组轮廓线,用三角片拟合过这组轮廓线的曲面。Bussonnat提出了另外一种基于表面轮廓的Delaunay三角形方法,解决了系列表面轮廓的三维连通性问题。用三角形或多边形的小平面(或曲面)在相邻的边界轮廓线间填充形成物体的表面,所得出的只是分片光滑的表面,Lin采用从轮廓出发的B样条插值重建算法,得到了整体光滑的表面。Lorenoen提出了一种称为"MarchingCubes”法,这是ー种基于体素的表面重建方法,该方法先确定一个表面阈值,计算每一体素内的梯度值,并与表面阈值进行比较判断,找出那些含有表面的立方体,利用插值的方法求出这些表面。这种方法虽然也是用三角形拼接来形成表面,但与连接轮廓线法不同的是,这些三角形位于单个立方体元内。这种方法避免了相邻切片间等值线连接的困难,可以直接生成三维的等值曲面。该方法对于ー组体数据,可通过门限设定,将网格交点划分为在表面内和在表面外两种情况,这样,每个立方体元被表面切割的情况就可由8个顶点的值确定。然后用三角形将体元各边上的交点连接起来,从而构造出表面。另外,移动立方体法可以获得较高分辨率的图像,图像生成速度快,且易交互。但是,在相邻体元的共享面存在4个交点时,会产生歧义面,从而生成错误表面。为解决这ー问题,提出了分解立方体法(DividingCubesAlgorithm),即将立方体分解,直至像素大小,以直接绘制表面点。基于表面的方法主要优点是可以采用比较成熟的计算机图形学方法进行显示(如裁剪,隐藏面消除和浓淡计算等),计算量小,运行速度快,借助于专用硬件支持,可以实现实时交互显示。2.2.2体绘制方法近10年来,体绘制方法以其在体数据处理及特征信息表现方面的优势,已得到研究者越来越多的重视,被越来越广泛地应用于医学领域。这类方法依据视觉成像原理,首先构造出理想化的物理模型,即将每个体素都看成是能够接受或者发出光线的粒子,然后依据光照模型及体素的介质属性分配一定的光强和不透明度,并沿着视线观察方向积分,最后在像平面上就形成了半透明的投影图像。由于体绘制方法基于体绘制方程,在具体实现过程中表现出多种多样的思想,从而演绎出多种具体算法。现分别简略介绍如下:(1)按算法中数据处理顺序分类
10按照算法进程中数据输入及处理顺序,体绘制方法可以分为像空间序(Image-order)法、物体空间序(Object-order)法和混合序(Hybrid-order)法〇(2)按体数据域性质分类体绘制方法按数据域性质可以分为空间域方法和变换域方法。其中,空间域方法即直接对原始的空间数据进行三维显示。而变换域方法则是将体数据经过变换后再进行显示。目前讨论较多的变换域算法是傅立叶断层投影体绘制法和基于子波的体绘制法。此外还有基于离散余弦变换(DCT)的体绘制法,等等。2.2.2混合绘制方法混合绘制方法分为两种:ー种是表面的透明体素绘制法,它是以体绘制的原理来实现对ー个或多个表面的绘制;另ー种是体数据几何单元投影法,即将由体素集合构成的单元投影转化为儿何多边形显示。其中,表面的透明体素绘制法是将所关心的表面提取出来,并赋予其所在的体素相应的光强和不透明度,再运用体绘制方法来实现三维显示。而体数据几何单元投影法以往常被用于不规则网格体数据的三维显示。可是对于规则网格体数据,它需首先将数据分解成同性物质的长方体,再按深度划分,并将长方体的面扫描转换到像空间:最后在每个长方体的前后两个面之间,做体绘制积分,以计算出每个像素点的颜色和不透明度,再合成图像。2.3各类算法的特点与应用由于表面绘制方法所处理的数据通常仅是整个体数据的ー小部分,并且利用了计算机图形学多边形绘制技术,还借助图形硬件加速的支持,所以表面绘制法具有速度快的优点,而且可以决速灵活地进行旋转和变换光照效果。它适用于绘制表面特征分明的组织和器官(例如由CT数据生成骨骼三维图像),由于其形象清晰,可以ー・定程度上替代实物模型。但是因为其对表面分割精确程度要求高,所以表面绘制方法对于其它ー些应用,例如对形状特征不明显、有亮度变化特性的软组织,以及血管、细支气管等精细组织或器官的三维显示,常常效果不佳。而且,表面绘制方法不能保留数据的完整性,其物体仅显示为一个空壳,表面里面没有东西。目前一种补偿的办法是用二维平面切割表面,然后以原始体数据填充概面。此外,其单纯的表面模型不能提供触摸反馈,因为这种反馈是由物体内部结构产生的。这些将制约其广泛使用。与表面绘制不同,体绘制方法通常不要求对被显示物做精确的分割,而是
11对体数据场中每个体素分别进行处理,进而合成具有三维效果的图像。因此,对于形状特征模糊不清的组织和器官进行三维显示时,适合采用体绘制方法。更重要的是,在体绘制方法中,透明度的引入大大增强了数据整体显示效果。通过对不同的组织分配相位的透明度,可以同时将各组织器官的质地属性、形状特征及相互之间的层次关系表现出来,从而丰富了图像的信息。但是,在原始的体绘制过程中,一般要遍历体数据场中的每ー个体素,因而使得计算量大,图像生成速度慢,并且不能灵活地改变外部光照及视角,这样每一次变化都意味着整个绘制过程需重新开始。因此,体绘制更多地受到硬件技术发展的限制,因而在要求实时显示的应用场合难以胜任。这种缺点突出地体现在像空间序法,尽管能够有效地进行高质量的重采样,并且可以通过利用数据相关性及提前结束光线的方法来加速绘制,但总的来说,沿每条视线投影都涉及到儿乎整个体数据,因此不但计算采样点需花费很多时间,使得生成图像速度慢,而且其中庞大的内存开销也是ー个突出问题。为了获取较快的显示速度,通常采用的体绘制方法是物体空间序法,因为物体空间序法体素寻址简单,可以依存储顺序实现数据空间的无缝输入,因而能最优地利用高速缓存。因此不但绘制速度较快,而且可以实现图像渐进显示。此外还可以通过对体数据进行预处理,采用例如八叉树等特殊的数据结构来存储数据,或是采用最大强度体素投影法来简化投影过程,从而加速绘制。混合序法(即剪切形变法)完全继承了物体空间序法简化体素寻址计算的优越之处,同时又以二维图像变形替代了体数据从物体空间到像空间坐标变换的部分运算,因此具有更快的绘制速度。而且,该法对物体空间序法中的诸多加速技术也同样可以采用。傅立叶断层投影法不但可以快速地得到任意角度具有X光片效果的图像,而且某种程度上比直接用X光成像更具灵活性。此外,傅立叶断层投影法还可以被用于某些特殊场合,例如,在二维核磁共振成像机制中,图像是由接收线圈采集到的自由感应衰减(FID)信号经二维傅立叶变换得到的。因此可以考虑直接将FID信号作为原始输入的体数据,用傅立叶断层投影法进行MR!图像的三维显示。医学图像的数字化存储是可视化及其应用的基础,为了节省空间,必须对海量数据进行压缩。基于子波的体绘制方法充分利用了子波变换空域和频域局部化的优良特性,因而可以实现对体数据的压缩,以减小存储量和计算量,并且可以对体数据在不同分辨率上进行绘制,同时可对特定的区域选择特定的分辨率。因此,这种方法对于像远程医疗这样以网络传输为基础的课题更具有应用前景。混合绘制方法可以认为是基于对表面绘制和体绘制方法的折中考虑,而产
12生的一类算法。但表面的透明体素绘制法相对于表面绘制方法则降低了对数据分割的要求,其表面体素的分割可以是模糊的,并且可以显示多个不同属性的表面,因而增加了信息层次;而相对于体绘制方法则减少了计算量,并突出了表面信息。而且,这种方港便于进行手术开窗显示。由于体数据儿何单元投影法像表面绘制方法一样,以多边形作为过渡显示单元,因此可以借助图形硬件加速实现体绘制。因而对于不太复杂的组织结构,这种方法可以用相对较快的速度来获得体绘制的效果。第三章医学图像的预处理与分割3.1概述医学图像体数据主要源自各种影像设备,像CT、磁共振、超声等。按照体数据在三维空间上离散数据间的连接关系来分类,医学图像体数据属于规则网格或均匀网格结构化数据。经过图像预处理以后即可用于图像的分割和三维重建。预处理的过程则包括了断层图像的滤波、插值和规则体数据的封装等内容。图像分割是医学图像三维重建和可视化研究中一个非常重要的部分,体数据中蕴含着各种各样的物体,只有把他们区分开,才能有目的有选择的对它们进行分析和显示。体数据分割的结果是若干个区域,分别对应于不同物体,通过对其特征量的划分可以确定各个区域代表的是肌肉、骨骼还是其它组织。体数据的分割主要是借用计算机视觉中的分割方法,像门限法、边缘抽取和区域生长等。这些分割方法按输出结果的不同可以分为两类:二值分割(BinarySegmentation)和模糊分割(FuzzySegmentation)〇二值分割结果对每个体素做出ー个肯定的二值判断,即它属于或者不属于某个物体;而模糊分割则与此不同,它只给出每个体素属于某个物体的可能性或者某种物体占有率。对CT、MRI等图像数据均可转化为规则体数据的形式,亦称三维医学图像,因此,对其进行分割是在三维空间进行。对医学图像的分割,到目前为止还不存在ー个一致公认的最好方法,也没有一个判断分割是否成功的客观标准。3.2医学图像的获取及其数据类型3.2.1三维空间数据场的数据类型可视化技术的应用领域非常广泛,实现三维空间数据场可视化的算法与数据类型有极大的关系。因此,先简要介绍一下数据类型。
13数据类型可以按照两种方式进行分类:ー是按数据本身的物理特性进行的分类,二是按数据分布及连接关系进行的分类。1.数据本身的分类根据数据的物理特征,数据场可以分为标量场和向量场。标量场是指数据只有大小而无方向的数据场。比较常见的有密度场,温度场等。标量场的可视化主要是揭示各分类物质的空间分布。向量场的数据不仅有数值的大小,还有方向的变化。典型的有流体力学的数据场。向量场的可视化除了揭示各分类物质的空间分布外,还要反映其变化的趋势。2.数据分布及连接关系的分类根据数据场中样点分布的空间几何特征,数据场一般分为规则场与非规则场。所谓规则场,就是样点分布于直线型正交网格交点上,且这些网格点沿各轴向的间距是相等的。在三维空间,每一-个网格即为ー个小立方体,称为体元(voxel)或单元(cell)。非规则场可以进一步分成结构化的和非结构化的。前者的含义是数据场在逻辑上依然是较规整的,而每个单元的形状不是正交网格型的。后者又称为散乱的数据场。下面按照从简单到ー一般的顺序给出每ー类的定义。1)笛卡尔(Cartesian)均匀网格结构化数据:这类数据均匀的分布在三维网格点上,即在x,y,z三个方向上,网格点之间的距离均相等。高密度的人体CT扫描所得数据是其中一个例子。在这类数据中,无需给出各数据点的空间位置,只要给出三维网格某ー角点的空间位置和某一数据点的序号,即可根据网格间距所对应的距离求出该点的空间位置。2)规则的(Regular)网格结构化数据:这类数据分布在由长方体组成的三维网格点上,即在x,y,z三个方向上,网格点之间的距离互不相等,但在同一方向上则是相等的。与均匀网格一样,在这类数据中,也无需给出每个数据点的空间位置,可根据起始点的坐标、某一数据点的序号以及x,y,z三个方向的增量求出该点的空间位置。3)矩形网格结构化数据(Rectilinear):这类数据同样也分布在由长方体组成的三维网格点上。但是,长方体的大小可以各不相同,并无规则可循。4)不规则网格结构化数据(structured):这种类型也称为曲线型的(curvilinear),,尽管在逻辑上仍然被组织成三维数组,但是在空间位置的分布上却无规则可循,因此被称为不规则网格结构化数据。在这类网格数据中,对于每个数据点,除了给出它在三
14维数组中的下标以外,同样,也必须给出他所在的空间坐标。典型的有三维有限元网格,流体动力学网格。5)非结构化数据(Unstructured):与以上的类型不同,网格的儿何邻接关系需要另外提供,元素可以是四面体、五面体、六面体等。网格点之间用曲边连接,并且有时连线上有两个离散点。6)混合型的(Hybrid):数据场可由以上任何几种类型的元素共同组成。按照以上的分类,医学图像体数据是一个基于规则网格的标量数据场。3.1.1医学图像的获取与输入医学图像主要源自各种医学影像设备,像CT、磁共振(MRI)和超声等。这些影像设备的基本原理都是相似的,它们都采用傅立叶变换、卷积法或Radon反演运算从多个方向上的射线(或波)投影反算空间的某种物理属性,经图像重建(ImageReconstruction,它与物体三维重建是两个完全不同的概念)重构出ー一幅幅断层图像。图像重建是影像设备本身的计算机硬件系统完成的。CT图像的CT值反映组织对X射线吸收值(衰减系数U)I其单位为HounsfiehUnit(Hu),通常定义物质的CT值如下:CT值=常数X("ール.)/人.,式中,ル是水的X线衰减系数(对有效能量而言),对大多数CT机来说,常数取值为1000,CT值单位为享(HU),故CT值(HU)=1000义("ール,)/q,-一般认为空气对X线的吸收为0,故空气的CT值为ー1000:而骨组织的X线衰减是水的2-4倍(现取2倍),故相应CT值为+1(X)0,水的CT值为〇。由CT值知,人体的组织密度可分成2000个密度等级,在分析CT图像时,可以提出每个像素CT值的确切数据,以便比较与正常组织或周围组织的差异。但为了更快的做出诊断,人们还是习惯于直观的看CT图像,这就需要用荧光屏来显示。为此,CT机通常引入256级灰度级,即把ー100〇〜+1000的CT值用256级灰度来表示其黑白层次。•般我们所处理的图像是经过量化后的灰度图像,标准的CT灰度图像为12位灰度图像。MR工成像是根据电子、质子、中子等都具有自旋和磁距的特性,当机体
15置入磁场中之后,机体的质子都会被磁化,倾向于与磁场的方向一致或相反方向排列。磁化的强度是一个可测量的矢量。磁共振成像主要依赖的四个因素:质子密度、T1(T1弛豫时间反映质子置于磁场中产生磁化所需时间)、T2(T2弛豫时间表示在完全均匀的外磁场中横向磁化所维持的时间)和流空效应。应用不同的磁共振射频脉冲程序以得到反映这些因素不同侧重点的图像,即T1加强像、T2加强像和质子密度像。超声成像基本上是采用脉冲回波方式成像。即,用ー个短暂的电脉冲激励换能器晶片,使之振动产生超声波并射入体内,进入人体的超声波在遇到组织界面时,就会产生较强的回波信号。于是,根据接收到的回波信号就可以直接获取扫查平面上的人体结构图像。CT图像一般为512X512或1024X1024,扫描层距为2-5mm,新一代的螺旋CT可获得更小间距的断层序列。MRI图像一般为256X256,断层间距为2-10mm〇近年来医学成像设备的物理分辨率不断提高,MR!成像分辨率可达到0.2ルm。影像设备自身带的计算机一般是各生产厂家专用的,数据存储格式不尽相同,常见格式有DICOM、TIFF,BMPo图像的传输方式可用磁带、光盘和网络。医学图像的网络传输国际上已有通用标准,即DICOM协议。3.1医学图像的预处理3.1.1医学图像滤波一般情况下,各类图像系统中图像的传送和转换(如成像、复制、扫描、传输以及显示等)总要造成图像的某些品质下降。在医学图像的获取过程中,影像设备中各电子器件的随机扰动,不可避免地会带来噪声。对其进行滤波(Filtering)或平滑(Smoothing)就是要抑制噪声,增强图像特征,提高信噪比。通常改善的方法有两类:ー类是不考虑图像降质的原因,只将图像中感兴趣的特征由选择的突出,而衰减其次要信息;另ー类方法是针对图像降质的原因,设法去补偿降质因素,从而是改善后的图像尽可能的逼近原始图像。第一类方法通常就称为图像增强即滤波技术,能提高图像的可读性,但改善后的图像不一定逼近原始图像。第二类方法称为图像复原技术,能提高图像质量的逼真度。对于医学图像来说,设法补偿降质原因比较困难,所以我们主要利用第一类方法即图像的增强技术对其进行处理,突出感兴趣的特征,衰减其次要信息。当输入图像信号混入噪声,想要用滤波方法把噪声全部滤除而不损失原信
16号的强度几乎是不可能的。因此,对滤波处理所追求的目标有两条:最大限度地保持信号不受损失,不能损坏图像的轮廓及边缘等重要信息;同时尽可能多地滤除噪声,使图像清晰,视觉效果好。图像的滤波技术通常有两类方法:空间域法和频率域法。空间域法主要是在空间域中对图像像素灰度值直接进行运算处理。而频率域法就是在图像的某种变换域中(通常是频率域中)对图像的变换值进行某种运算处理,然后变换回空间域。以下我们主要介绍空间域中的滤波技术。空间域中的滤波方法常用的有以下三种:1.邻域平均法对于图像中的每ー个像素,取一个以它为中心的区域,用该区域内各像素灰度的加权平均值取代该像素的灰度值,这就是邻域平均法。具体的做法是取ー个方形区域,称为平滑窗口(Window)或掩膜(Mask),它是权值的二维阵列。滤波过程是用窗口在图像上滑动,窗口中心对着的像素根据下面式(3-1)来更新其灰度值。当每个像素都被扫描ー次之后,对ー・幅图像的平滑就完成了。这是空间域中对图像进行平滑的ー种最简单的方法。设f(i,j)是一幅待平滑的图像,平滑窗ロ为w的大小为(2N+l)X(2N+l),则平滑后的图像可以表示为NNyX%・ハ,+”,ノ+レ)gG,j)=……ッーーN(3-1)SZ%”1其中,吗,,是权值。我们可以对窗口内的权值进行归ー化,即使权值之和等于1。如式(3-2)〇
17NN2Z=1(3.2)u=-Nv=—N一般来说,窗口越大,平滑能力就越强。但是,噪声的消除程度和图像原有信号的衰减程度皆与窗口的大小成正比,因此并不是窗口越大越好。实际中常用的有3X3和5X5两种掩膜。3X3窗ロ又可分为4ー邻域和8ー邻域,如图3-2所示。4ー邻域8ー邻域图3-2滤波窗ロ1/162/16V102/164/162/161/16,J2/16ビ6J图3-3模板示意图
18在本课题开发的程序中设计了四种模板以供选择,•个平均模板,•个高斯模板,如图3-3所示,还有如上图3-2所示的两个模板。另外还设计了一个可自定义的模板,用它可实现任何〈5X5的模板操作。平均模板和高斯模板如图3-3所示,图3-4为模板设置对话框,由图所示,选择“自定义模板”,这样就可以自己定义模板的系数、模板的高度和宽度、模板中心元素的位置、以及模板元素的取值。图3-4模板设置对话框2.中值滤波法由于邻域平均法在消除噪声的同时会将图像中的一些细节模糊掉。如果既要消除噪声又要保持图像的细节,可以使用中值滤波。中值滤波是一种有效的非线性滤波,常用于消除随机脉冲噪声。在一定条件下,中值滤波可以克服线性滤波器如最小均方滤波,均值滤波等带来的图像细节模糊,而且对滤除脉冲干扰及图像扫描噪声最为有效。它的基本思想是:在图像上滑动ー个含有奇数个像素的窗口,对该窗口所覆盖像素的灰度按大小进行排序,处在灰度序列中间的那个灰度值称为中值,用它来代替窗口中心所对应像素的灰度。对于奇数个元素,中值是指按大小排序后,中间的数值;对于偶数个元素,中值是指排序后中间两个元素灰度值的平均值。中值滤波的关键在于选择合适的窗口大小及窗口形状,实际上一般多采用以像素(I,J)为中心的正方形窗口,若图像大小为MXN像素,窗口大小为KXL像素,则窗口处理时,总计算量为O(MXNXKXL)量级,也就是说,窗口处理的计算量不仅与被处理图像大小成正比,也与所用窗口大小成正比。比如可以设计了1X3、3X1和3X3三个模板窗口可供选择,同时还有一个自定义的模板,可以通过指定中值滤波器的高度和宽度,以及中心元素的X,Y坐标来选定滤波窗口。其对话框如图3-5所示。
19QiBtSi史・弊おを匿宀…L呼室目像中值^8¢¢%;.亘き匣凄遠透c1XJ3中值滤波:.ir-3ス3中值建波;璟定IL自定义中伍鸿波・波器高度,百ー[.酒该窃亮度.0—ーヨ便ー1」中心元素标,[ij-中心元素y坐标,n—图3-5中值滤波对话框2.一-种改进的增强边缘的中值滤波算法按照以上所述的传统的中值滤波算法,中值滤波后的图像中随机噪声将被有效的消除,这是因为领域中亮度值发生随机突变的像素经过排序后,要么排在队列的队首,要么排在队列的队尾,而中心像素的新值是取自队伍中位的那个像素值。但在对人体关节图像进行平滑处理中发现,传统的标准窗中值滤波尽管对图像的边缘有较少的损失,但由于扫描进入计算机后的关节图像本身有边缘模糊特性,因而达不到真正改善图像质量的目的,对于所需的有用信息,经过滤波处理后没有突出的表现出来,因而需要对传统的中值滤波算法进行改进,从而获得更为清晰、准确、无误的CT图像。以下介绍ー种改进的边缘增强的中值滤波算法。如何使中值滤波既可消除孤立噪声,又能保证边缘信息不被削弱,在实践中可以利用计算机编程技巧对中值滤波算法进行改进。现假设采用3X3的滤波窗口,其它滤波窗口的改进方法类似可推。在计算机中,ー个字节为8位二进制数,处理中所采用的灰度图像为256级,其中每一点像素的灰度值正好用一个字节储存。对3X3矩阵灰度值进行排序时,采用无符号字符型变量(UnsignedChar)来存储像素的灰度值。对于3X3的灰度矩阵H为:'%h1为ヽH=月h4h5其中Z1,为无符号数,/i,e[0128],ze[08]6る4ノ在传统的中值滤波方法中,就是直接对%为进行排序,然后取中间的值作为该像素的灰度值,这样边缘信息就会被削弱了。为了增强边缘信息,我们针对矩阵H构造ー个中间矩阵M:
20,加。m,ヽM=叫m4m5其中町.也为无符号数,M矩阵可由以下语句得来、团6叫恤ノ设定一个边缘附近的灰度阀值为F,贝リ:叫=F,ie[02];w,=h,ie[38]〇设灰度值小于阀值F的像素为背景像素,灰度值大于阀值F的像素为前景像素。按传统的中值滤波方法,对M矩阵进行排序,选取中位值。对于改进的中值滤波法,首先,在轮廓的边缘处,在巧外中,应既含有背景像素又含有前景像素。这样排序结果中值一定为飞吗中的一个值,即为F,这时,可以设定返回像素灰度值为0,即把边缘处得像素涂为黑色,这样就可以突出边缘信息了;其次,在其它地方,在吗恤中除孤立噪声像素外,要么同为前景色要么同为背景色像素,所以へ心对排序结果影响不大,效果同传统中值滤波,能去除孤立噪声。(a>延波前原始圏像华抚中体混波后的图像 21由图3-6的对比图可以看出,用传统的中值滤波得方法可以有效的去除孤立噪声,但是边缘信息被忽略了。而采用改进后的中值滤波算法,结果图中不仅有效的去除了噪声,而且使图像的边缘信息也得到了加强。3.3.2二维断层图像之间的插值当断层图像间的距离比断层图像内像素间的距离大得多时,就需要用图像插值方法在原来的断层图像之间再插值生成一些中间断层图像。因为许多体数据三维重建和三维显示方法要求体数据是各向同性的。在医学图像三维重建为了便于重建组织的剂量计算,要求每个体素为正方体。因而断层图像就要插值成各向同性,即经插值后的断层图像序列中断层间距等于断层图像内像素间距。然而,图像插值是ー个具有很大任意性的问题。为了使图像插值成为ー个确定的、可解的问题,通常引入下面三个约束条件:⑴插值图像要与原始断层图像相似;(2)插值图像与两个原始断层图像的相似度应该分别和它与这两个断层图像的距离成反比关系;(3)插值图像序列应该呈现出从一幅原始断层图像到另一幅原始断层图像的渐变过程。最简单的插值方法是对上下两个相邻的断层图像进行加权平均,产生一组插值图像。当断层间距与断层图像内像素间距相差不是很大时,采用这种插值方法是可行的。但是,当断层间距与像素间距相差很大时,这种方式生成的插值图像模糊不清。其原因是,两个相邻断层图像中处于同一位置的像素不一定对应同一物质,他们的加权平均没有什么意义。对断层间距较大的插值,ー种较好的方法是基于匹配的图像插值,这种方法能够更好的构造中间断层图像。关于断层图像的匹配是模式识别和计算机视觉领域所研究的主要问题之一。图像匹配的任务就是要在两幅图像之间建立对应关系,也就是寻找图像间的变换。传统的匹配方法一般包括特征选取、特征匹配和整幅像之间匹配映射关系的确定三大步骤,匹配的可靠性和准确性很大程度取决于是否能够找到鲁棒的方法去完成特征提取和特征匹配。基于匹配的断层图像插值就是在图像匹配的基础上对断层图像进行插值,这里只提出这种方法,详见文献。本文只简单介绍一下线性加权平均的图像插值方法。线性加权平均的图像插值方法可描述如下: 22设&+(•))分别是第k层和第k+1层切片图像。它们之间的插值图像可表示为:=(i,ノ)+え・S«m(3-3)其中,A=dl/(dl+d2),dl、d2分别是插值图像到第k、k+1层图像的距离。显然,当;I=0时,=Sk(i,j);当;I=1时,5式i,j)=Sk+l(i,j)〇给出ー组H邑€(0,1)「=1,2,...〃}就相应的得到ハ个插值图像。为了得到等间距的插值图像,参数序列L=应该取4=ユー3.4医学图像的分割正确的分割是重建模型能够正确表达重建的组织或器官的关键。由于分割是ー个不适定问题,至今仍没有通用的分割方法。只能根据具体的应用和图像的特点,设计各自的分割算法。本文医学图像的分割是针对三维体数据进行的,分割需要操作者根据图像特点及被分割重建的组织、器官的形态特征,来交互地输入适当的数值作为分割阈值,在采取形态学操作方法进行区域修整,然后由系统自动提取出被分割组织。分割的过程是事先在内存中拷贝ー份体数据,然后对拷贝体数据通过阀值法进行二值化,数学形态学区域修整,种子填充,最后提取并记录填充出的区域。待整个分割过程结束后,清除拷贝体数据,释放内存。3.4.1体数据的二值化在分割之前对体数据的灰度值进行分类处理,将三维图像处在要求提取组织的灰度范围内的像素(体素)置为1,其余置为。。这可以用阀值法进行处理,设定两个阀值hl和h2,然后设输入体数据的灰度值为f(x,y,z),输出体数据的灰度值为f'(x,y,z),则可以用公式表示如下:f1/rヽ,ク、_11若刀に,J丿ー[〇其余(3.4)由于在数字化的图像数据中,无用的背景数据和对象物的数据常常混在ー起,除此以外,在图像中还含有各种噪声。因此选择合适的阈值是ー个极为重要的问题,可以用灰度直方图来确定阈值。灰度直方图是数字图像处理中一个最简单、最有用的工具,它描述了一幅图像的灰度级内容。 23灰度直方图是灰度值的函数,描述的是图像中具有该灰度值的像素的个数,其横坐标表示像素的灰度级别,纵坐标是该灰度出现的频率(像素的个数)。直方图能给出该图像的概貌性描述,例如图像的灰度范围、每个灰度的频数和灰度的分布、整幅图像的亮度和平均明暗对比度等,由此可以得出进ー•步处理的重要依据。分析三维图像中的儿张断层二维图像,确定被提取组织或区域的大致灰度范围,然后根据三维图像的整体灰度直方图及其一些二维断层图像的灰度直方图,确定被提取组织或区域的灰度阈值的下限和上限。对于CT图像,像素灰度值的大小直接反映出组织密度的大小。对于MRI,由于成像机理与CT不同,图像不具有物理(密度)相似性,通常,皮肤和脂肪具有高灰度值而骨骼的灰度值儿乎与背景•样。根据图像的灰度直方图的峰值数基本能反映出图像包含组织的大致数量.谷值为不同组织间的分界限。图3-8为CT图像的一般情况,由图3-8的峰值个数可将图像中的物质(组织)分成四类:空气、脂肪、软组织和骨骼。谷值可视为四种物质(组织)间的分界限。 24软组织图3-8各种物质灰度值分布曲线对二维灰度图像进行滤波处理之后,编程绘制其二维断层图像的灰度直方图。图3-9(a)是64X64X93,256级灰度的CT图中的第10层和第40层,图3-9(b)是相应断层图像的灰度直方图。我们要提取的是骨骼组织。首先我们用光标在待提取区域扫描,预估出待提取区域的灰度范围,然后,分析整体灰度直方图和断层灰度直方图,骨骼组织灰度阈值在150至200之间,小于150为背景和脑组织,大于200为骨骼组织,在这之间的为ー•些其它软组织。按照这两个阀值可以分割出骨骼组织以及表面轮廓来。如图3-10所示。第府饪断层出像第4n层断,片灰衣直方图 25第10层红京图缘竺】。县も荒灰支史方图(a)断层图像(b)断层灰度直方图图3-9断层图像及灰度直方图(a)骨骼分割图(b)面部轮廓分割图图3-10分割二维结果图3.4.1种子填充法进行组织提取三维图形经过二值化操作后,会形成多个子区域,而我们分割的目的是分割出某ー待建组织或器官所在区域。为了在这些子区域中提取出待重建的区域,我们采用种子填充的方法来实现。我们处理的是三维医学图像,即填充在三维空间中进行,填充生成的区域是三维连通区域。可用二维图像的种子填充算法的基础上增加Z方向功能来实现。种子填充算法 261.开辟ー个堆栈:2.在要提取的区域中输入ー个点(可直接中键盘输入,也可以鼠标点取)作为种子点,将种子点压入堆栈中;3.从堆栈中取出ー个点,将该点置为0(涂为黑色),并记录下该点:4.按弹出点的左一上一前ー右ー下ー后的顺序查看它的六个相邻点,若为1(表示没有被填充),则将该点压进堆栈;5.检查堆栈是否为空,为空,清除堆栈,结束,否则,转3。3.4.3分割结果及分析图3-14为64X64X93的CT图像中的三幅断层图像。通过分析断层图像的灰度直方图,确定骨骼组织和皮肤组织的灰度阈值。骨骼组织灰度阈值在150至200之间,小于150为背景和脑组织,大于200为骨骼组织,在这之间的为ー些其它软组织。经过二值化处理,然后为切断与骨骼或皮肤组织具有弱连通的其它组织的联系,选用开启操作。经开启操作,种子填充提取出所要重建的组织。图3-14(a)、(b)、(c)分别为第10、40、60层的断层图像。图3-15为对应的骨骼组织的分割图像,图3-16为对应的皮肤组织的分割图像。图3-14断层图像 27图3-15分割后的骨骼组织 28值化后出现不同的组织器官连结在ー起,难以达到正确的分割。本文提出的分割过程是先分析待重建组织所处的灰度范围,给出阈值的上下限,对三维图像(体数据)二值化,然后根据被提取组织器官的形态特征进行区域修整,最后对重建区域用种子填充的方法提取出来,完成待建组织的分割。对于医学图像,要达到分割结果满足治疗的安全可靠,自动分割出的结果必须经医师或治疗计划师确认,对于重要组织或器官、敏感组织以及病变体必须由医师或治疗计划师手工轮廓勾画提取。要达到正确、满意的分割,分割一重建是-•个重复的过程,对分割出来的组织重建出几何模型后,对三维模型进行不同角度观察来判断其分割的正确性。如果对分割及重建结果不满意,可删除重建模型,调整分割阈值及选择适当的形态学操作方法及结构元素进行区域修整,再进行模型重建。如此重复,直到得到满意的分割与重建结果。 29第四章基于面绘制的医学图像三维重建4.1引言从医学图像中重建三维物体的方法很多,在第二章已经作了一个比较全面的介绍。基于面绘制的医学图像三维重建从重建过程处理的基本元素的级别上来分,可以把这些方法分成两大类:体素级重建方法和切片级重建方法。体素级重建方法首先要确定物体表面在每个体素内的小面片,然后将这些小面片连接起来就构成了物体的表面。切片级重建方法则是从・组平面轮廓重构通过这些轮廓的曲面。当原始图像的分辨率很高时,体素级重建方法比切片级重建方法更可靠、更有效;而当原始图像的分辨率很低时,体素级重建方法的重建精度也很低,这时切片级重建方法能够比较好的构造出光滑的表面。总体来说,体素级重建方法相对比较可靠,因为它只需要保证各个小面片之间的拓扑一致性,不需要考虑总体的拓扑关系,但重建的结果却产生大量的小面片,占用大量的存储空间,即使对于几何结构非常简单的物体也是如此。因此在保证一定精度的前提下减少小面片的数量就成为一个需要研究的问题。面切片级重建方法可以实现大幅度的数据压缩,但轮廓对应存在着多义性,特别是在出现分叉情况的时候使得轮廓对应问题的不确定性更加严重。这两种方法各有优缺点,本章将对这两种方法都进行讨论。主要介绍ー种体素级的重建方法及其具体实现,并给出部分实验结果。4.2切片级重建方法切片级重建方法就是由一组二维轮廓线重构三维形体,在医学数据的可视化中,这ー方法也有着广泛的应用前景。例如,如果医学断层图像两层切片之间的距离远远大于图像的分辨率,也就是Z方向的采样十分稀疏,那么,前面介绍的体绘制方法将难于应用,这是,需要在每层图像上就感兴趣的区域勾画出轮廓线,并由一个序列的二维轮廓线重构出由面模型表示的三维形体。又如,在人体断层解剖图像的可视化中,尽管相邻两片断层图像之间的距离很小(例如1mm左右),但是目前的物质分类技术尚不能将共存的多种不同物质加以分类,往往需要人工干预,由医学专家手工勾画出感兴趣的区域,如血管、神经等。因面也需要由一个序列的二维轮廓线重构三维形体。切片级重建的输入是ー组平行的平面轮廓。要使用切片级重建方法恢复体数据中所蕴含的物体结构,首先必须应用分割技术提取每个断层上的物体轮廓;然后由轮廓线重建三维形体,其中又包括四个基本问题:轮廓对应问题、轮廓拼接问题、分叉问题以及最后的曲面拟合。以下分别对这两部分进行讨论。 304.1.1轮廓的提取物体的边界在图像中是由灰度不连续性所反映的,人们希望直接从图中提取边界,但在成像条件不好,边界复杂时是很难做到的。局部边沿在图像的总像点中只占不大比例,那里局部以一定方式迅速改变。根据边缘的构型,可以将它分为阶跃型和屋脊型。基于边的分割方法是要找出灰度不连续的点作为边缘。阶跃型边缘点是灰度ー阶导数的极大值点,或者二阶导数的过零点;而屋脊型边缘点是灰度ー阶导数的过零点,或者二阶导数的极小值点。边缘检测算子是小空间范围的数学算子,常用来检测图中是否有局部边沿。最早使用的梯度算子是Roberts算子,它是ー种差分算子,由于窗口小,平滑噪声作用也小,改进的另外二种算子是prewitt算子和Sobel算子,它们对有噪声的阶跃形边沿检测能力比Roberts好,我们采用Sobel算子,根据可调的阈值进行二值化,最后采用八步跟踪的方法找出相应的边界。下面介绍儿种经典的边缘检测算子,用它们来进行医学图像的轮廓提取。4.1.1.1梯度算子图像的局部边缘定义为两个强度明显不同区域之间的过渡(由低灰度区域到高灰度区域的过渡或由高灰度区域到低灰度区域的过渡)点,图像的梯度函数即图像灰度变化的速度,将在这些过渡边界上存在极值(最大值),因此,我们可以首先得到灰度变化的速度(梯度算子),如果梯度值大于某个给定的阈值,则认为对应点为边缘点。若再将被确定的边缘点连接起来,就可以形成包围着区域的封闭曲线。对于二维图像中的每ー个像素(x,y),计算它的梯度。梯度定义为:图像函数f(x,y)在点(x,y)的梯度为: 31大小:(4-1)方向:Mx,y)=arctan(Gt/Gv)(4-2)对于数字图像的离散域,可以用图像的ー阶差分直接替代图像函数的偏导数:故-Iア”,f(x,y-V飲+な〃飲,初f(x+l,y)佟[ノ+り张J+1)任+,ノ+り二维离散图像函数,在x方向上的ー阶差分定义为:f(x+l,y)-f(x,y)在y方向上的ー阶差分定义为:f(x,y+l)-f(x,y)对于一阶微分算子,可以先求出一阶差分,然后代入g(x,y)的公式中,得到梯度的大小:g(…符十雷(4-3)梯度g(x,y)的大小反映了图像灰度局部变化的强弱,因此可以作为检测边缘点的依据。剃度算子可以有效的检测出阶跃型边缘,而对屋脊型边缘容易产生双边缘效果。4.1.1.1Roberts算子Roberts算子是Roberts提出的另ーー种类似于梯度算子的算子。是ー种利用局部差分算子寻找边缘的算子。它由下式给出:g(x,y)ay)—/(x+1,y+1)F+[/(x,y+l)-/(x+l,y)](4-4)其中f(x,y)是具有整数像素坐标的输入图像,平方根运算使该处理类似于在人类视觉系统中发生的过程。 32它是在ー个2X2的邻域上计算对角导数,即,g(x,y)的值是用方向差分的均方值R(x,y)来近似,R(x,y)又称为Roberts交叉算子。Roberts算子可以用模板进行计算,模板如下:。0)(0P(。T丿1-I。ノ两个模板分别实施相加。图4-2给出了部分CT断层图像的原始图像及其经过Roberts算子提取边缘后的边缘图像。从对断层图像提取的Roberts边缘图像可以看出,Roberts算子基本上检测出了边缘,但Roberts算子存在的问题是进行差分计算的过程中对噪声敏感,即有噪声影响的像素点可能被检测为边缘点。4.1.1.1Sobel算子Sobel算子是在一个3X3的邻域中,计算x和y方向的偏导数如下:4/■(切kD+Mx+LザザD+M%Ty)ザ(ア1时明A/Uy)ヵUTy+D+次v+Dザ・a+Ly+D]H/(xTy-D+ザ(えパザU+Ly-D]f(艾七I,y-Of(x-I,y)f(x,y)やりf(x,y^-1)7,y+2)实际上,Sobel算子是将差分运算与局部平均相结合,它应用f(x,y)邻域的灰度加权平均来求差值。g(x,y)«S=イ公J(x,y)-+△J(x,y), 33或:g(x,y)»S=|Ax/(x,y)\+|Av/(x,y)|(4-5)和其它的梯度算子ー样,Sobel算子也可以用卷积模板来实现即图像中的每个点都用这两个核做卷积,ー个核对通常的垂直边缘响应最大,而另ー个对水平边缘响应最大。两个卷积的最大值作为该点的输出。这…算子把重点放在接近于模板中心的像素点。Sobe!算子比梯度算子的抗噪声能力要强ー些,这是因为Sobe!算子是ー个8ー邻域算子:而梯度算子是ー•个4ー邻域算子。由图4-2的边缘图像可以看出,Sobel算子取得了较好的效果。Sobel算子具有计算简单,运算速度快的优点。但Sobel算子提取的边缘细节比较多,抑制噪声性能不强。 34GOあ星图堂原留ゴ9)R。加kr提取边簿图4-2各种算子提取边缘对比图4.1.1平面轮廓重建从ー组平行的平面轮廓重建三维形体必须解决下面四个问题:1.轮廓对应问题(correspondenceproblem)由于约束不足,轮廓对应存在很大的随意性。特别是当切片间距很大时,轮廓之间的错位也会很大,对应问题就变得愈发难以解决。到目前为止,还没有完全可靠的解决轮廓对应问题的全自动方法,因此人工干预在某些场合下是必不可少的。目前常用的手段有轮廓覆盖检查和基于全局的广义柱体方法。 351.轮廓拼接问题(Tillingproblem)因为通过空间中的两个轮廓的曲面可以有无穷多种形式,因此必须额外引入ー些约束才能使这个问题有解。常用的约束条件有体积最大、表面积最小、跨度最小或者是轮廓点匹配方向最大程度地与质心匹配方向一致。轮廓拼接就变成基于这些优化目标的优化问题。2.分叉问题(Branchingproblem)当ー个物体在一对相邻断层上的轮廓个数不相等时,就发生了分叉情况。在轮廓对应出现ー对多关系的情况下,确定形体表面的多边形拼接。分叉发生时的轮廓对应关系一般无法由分叉发生的局部信息来确定,分叉处理必须依赖于形体全局的拓扑和儿何结构。3.曲面拟合(surface-fittingproblem)根据上述步骤确定的点之间的连接关系和表面拓扑结构,确定最佳的拟合曲面。一般选用参数曲面拟和拼接得到的多边形网格,产生光滑的表面。在拟合过程中,网格点作为控制点,网格点之间的连接关系决定曲面的拓扑结构。在切片级重建中,轮廓对应和拼接以及对分叉情况的处理都是关键性的问题,它们的弱约束本性是造成它们具有非常大的随意性的根本原因。4.3体素级重建方法在第三章中我们对医学图像数据进行了分割,对于对用第三章的分割方法分割出的区域,分割出的待建组织区域均以体素形式记录。我们进行表面重建的目的就是对分割提取出的区域构建出其对应组织或器官的三维儿何模型。其方法是在提取的区域中,根据给定阈值,构建出表达重建组织或器官的等值面。等值面的构造是从体数据恢复物体三维几何描述的常用方法之一。如果我们把体数据看成是某个空间区域内关于某种物理属性的采样集合,非采样点上的值以其邻近采样点上的采样值的插值来估计,则该空间区域内所有具有某ー个相同值的点的集合将定义ー个或多个曲面,称之为等值面。因为不同的物质具有不同的物理属性,因此可以选定适当的值来定义等值面,该等值面表示不同物质的交界。也就是说,ー个适当值定义的等值面可以代表某种物质的表面。 36最早的体素级重建方法叫做立方块法,它是用边界体素的六个面拟合等值面,即把边界体素中相互重合的面去掉,只把不重合的面连接起来近似表示等值面。这种方法的特点是算法简单易行,便于并行处理,因为对每个体素的处理都是独立的:主要问题是出现严重的走样,显示图像给人ー种“块状”感觉,尤其在物体边界处锯齿形走样特别醒目,而且显示粗糙,不能很好地显示物体的细节。而Lorensen等人于1975年提出的移动立方体法(Marchingcubes),是最有影响的等值面构造方法,一直沿用至今。最初的Marchingcubes算法不能保证三角片所构成的等值面的拓扑一致性,会造成等值面上出现孔隙。MJ.Dorst首先提出了MC算法中的二义性,后来许多人在Lorensen方法的基础上做了许多改进。移动四面体法(MarchingTetrahedra)是在MC算法的基础上发展起来的,该方法首先将立方形体元剖分成四面体,然后在其中构造等值面。由于四面体是最简单的多面体,其它类型的多面体都可以剖分为四面体。在本章中,简要介绍体素模型及等值面的定义,讨论MC算法和MT算法的原理及实现方法,分析MC算法在连接方式上二义性的产生和消除,MT算法的剖分一致性等。对分割出的组织区域,给定等值面值,实现了MC算法和MT算法构造其三维表面模型。4.3.1体素模型与等值面定义4.3.1.1体素模型在三维空间的某一个区域内进行采样,若采样点在x,y,z三个方向上的分布是均匀的,采样间距分别为ハx,Ay,Az,则体数据可以用三维数字矩阵来表示,即s(i,j,k)o每八个相邻的采样点所定义的立方体区域就构成了一个体素,而这八个采样点称为该体素的角点,它们的坐标分别为(i,j,k),(i+1,j,k),(i,j+1,k),(i+1,j+1,k),(i,j,k+1),(i+1,j,k+1),(i,j+1,k+1)和(i+1,j+1,k+l),如图4-5所示。对于体素内的任一点,其值只能从体素的八个角点的采样值来估算。常用的体素模型有方向无关的三线性插值模型和方向有关的三线性插值模型两种。对于体素内任意一点P6=(x,y,z),其物理坐标可转换成图像坐标厶,女6),其中’6=丁,ノ6=ぎ,k6=3。AxAyAz当把方向无关的三线性插值作为体素模型时,其值可以表示为s(P6)=[s(P4)・(i+l-i6)+s(P5)・(i6ーり](4-6)而S(P4)=[S(Po)・(/+lーノ6)+5(。3)・(ノ6ーノ)](4-7) 37S(P5)=[s(P1)・(j+1ーノ6)+S(P2)•(ノ6-»](4-8)s(p“)=[s(f(〃)J(〃),幻・(た+1ーた6)+s(f(〃),/(〃),た+1)•/6-A)](4-9)(n=0,1,2,3)其中:プ(か)_(M〃=0』;)_(バ〃=0」;)1丿M+l(〃=2,3;)J1〃丿ー、ノ+1(〃=2,3;)公式(4-6)经过整理得:s(x,y,z)=。〇+a[x+a2y4-a3z+a4xy+a5xz+a6yz+^xyz(4-10)其中q(i=0,1,…,7)为常数,它们由体素的八个角点唯一决定。可以看出,三线性插值是ー维线性插值的组合。显然,其结果与这三个ー维线性插值的顺序无关。图4-5体素模型当把方向有关的三线性插值作为体素模型时,对于体素面上的点,和方向无关的三 38线性插值ー样,其值由该点所在的面上的四个角点值的双线性插值确定。而对于体素内的任意一点,经过它沿着某个给定方向所做的直线将与体素的表面交于两点,该点的值是关于这两个交点上的值的线性插值,如图4-5(b)所示。因此,这两种插值模型在体素的面上是相同的,而在体素内是不同的。三维医学图像是标量数据场,与方向无关,采用的体素模型为方向无关的三线性插值模型。4.3.1.1等值面定义等值面是空间中所有具有某个相同值的点的集合。它可以表示成{(x,y,z)Is(x,y,z)=c},c是常数。(4-11)并不是每个体素内都有等值面,当体素的八个角点都大于c或者都小于c时,其内不存在等值面。只有那些既有大于c的角点又有小于c的角点的体素オ含有等值面,我们称这样的体素为边界体素。等值面在ー个边界体素内的部分称为该体素内的等值面片。等值面是一个三次曲面,它与边界体素面的交线是ー条双曲线且这条双曲线仅由该面上的四个角点决定。这些等值面片之间具有拓扑一致性,即它们可以构成连续的无孔的无悬浮面的曲面(除非在体数据的边界处)。因为对于任何两个共面的边界体素,如果等值面与它们的公共面有交线,则该交线就是这两个边界体素中等值面片与公共面的交线,也就是说这两个等值面片完全吻合。所以,可以认为等值面是由许多个等值面片组成的连续曲面。由于等值面是三次代数曲面,构造等值面的计算复杂,也不便于显示,而多边形的显示则非常方便,所以,等值面的三角片拟合是常用的手段。我们本章论述的MC算法和MT算法都是在边界体素中生成三角面片,以三角面片拟合等值面。4.3.2等值面的构造方法综述等值面的构造方法主要有三种:立方体法、移动立方体法和分解立方体法。立方体法是Herman等人在1979年提出的,而移动立方体和分解立方体这两种方法都是由しorensen和cline两人合作,分别在1957年和1988年发表的,直沿用至今。本文主要采用移动立方体法进行等值面的构造,在这里先对其它两种方法简要介绍一下。1.方体法立方体法的执行过程由两步骤组成。第一步是边界体素的确定,对于某个 39设定的数值c,对每个体素遍历一遍,比较每个体素八个角点值与c的大小,并以此为依据判定他是不是边界体素。对于每个边界体素,它的六个面是六个多边形(而且是长方形)。第二步是对所有边界体素的六个由多边形表示的面进行显示。如果在具有硬件深度缓存(Z-buffer)功能的计算机上运行立方体法,可以把这组多边形不分次序的提交给硬件,由硬件完成隐藏面消除的任务。如果以软件方式执行立方体法,在算法中必须考虑多边形的遮挡问题。一个有效的方法是把遍历体素集合与显示两步骤合二为ー,遍历体素集合相对视点采用从后至前的次序。发现ー个边界体素,就立刻显示它的六个面。后显示到屏幕上去的多边形将覆盖先显示的多边形,这样就达到消除隐藏面的目的,这就是著名的画家算法的思想。立方体法的特点是算法简单易行,便于并行处理,因为对每个体素的处理都是相对独立的。其主要问题是出现严重的走样,显示图像给人ー种“块状”的感觉,尤其在物体边界处锯齿形走样特别醒目。而且显示粗糙,不能很好的显示物体的细节。立方体法的另ー个缺点是面的重叠冗余问题。两个相邻边界体素的公共面重复出现,实际上他们都不会出现在显示画面上,因为无论从哪个角度进行显示,他们都会被这两个体素的其它面遮挡。改进的立方体法删除了边界体素之间的公共面,减少了显示过程需要处理的多边形的数量。1.分解立方体法分解立方体法(Dividingcubes)是cline等人于1958年提出的另一-个等值面构造方法。与移动立方体法不同的是,分解立方体法是以点而不是三角片作为基本元素来构造等值面。只有当表示等值面的点的密度足够大时,其显示图像ォ会呈现“实”的物体。一般要求点的密度不低于显示图像的分辨率。Cline等人在提出分解立方体法时的一•个基本出发点是,点的显示ー一般要比多边形的显示快。当然,不同几何基元的显示速度决定于机器的图形处理硬件。有些系统己经带有多边形显示处理硬件,对多边形的显示速度不亚于对点的显示速度。此外,点的表示形式比多边形的表示形式更适合物体的剖切、合并等几何操作。然而,以ー组离散点集表示物体表面也存在ー些缺陷,ー个突出的问题是表面放大时点之间的不连续性会造成显示出现非实体感。因此,预先构造・定放大倍数的足够密的点集,或者动态地构造适合当前显示比例的点集是必要的。4.3.3移动立方体(MarchingCubes)算法抽取等值面在三维空间规则数据场中构造等值面的众多方法中,最有代表性的就是Marching 40cubes方法,它是由W.ELorenSon和H.Ecline在1987年提出来的。由于这一方法原理简单,易于实现,目前已经得到了较为广泛的应用。4.3.3.1MC算法的基本原理在MC算法中,假定原始数据是离散的三维空间规则数据场,用于医疗诊断的断层扫描(CT)及核磁共振成像(MRI)等产生的图像均属于这ー类型。MC算法的基本思想是逐个处理数据场中的立方体(体素),分类出与等值面相交的立方体,采用插值计算出等值面与立方体边的交点。根据立方体每ー顶点与等值面的相对位置,将等值面与立方体边的交点按一定方式连接生成等值面,作为等值面在该立方体内的一个逼近表示。因而,MC算法中每ー单元内等值面抽取的两个主要计算是:⑴体素中由三角片逼近的等值面计算; 41(2)三角片各顶点法向量计算。1.体素中等值面剖分方式的确定MC算法的基本假设是沿着立方体的边数据场呈连续线性变化,也就是说,如果一条边的两个顶点分别大于、小于等值面的值,则在该边上必有也仅有一点是这条边与等值面的交点。确定立方体体素中等值面的分布是该算法的基础。首先对立方体的八个顶点进行分类,以判定其顶点是位于等值面之外,还是位于等值面之内。再根据ハ个顶点的状态,确定等值面的剖分模式。顶点分类规则为:1:如立方体顶点的数据值?等值面的值,则定义该顶点位于等值面之外,记为“〇”;2Z如立方体顶点的数据值(等值面的值,则定义该顶点位于等值面之内,记为"1"。由于每一体素共有8个顶点,每个顶点共有2个状态,因此共有256种组合状态。根据互补对称性,即体素的顶点标记置反(0变为1,1变为0),不影响该体素内三角面片的拓扑结构,这样256种构型可以简化成128种。再根据旋转对称性,可将这128种构型进ー步简化成15种。图4-6给出了这15种基本构型的三角剖分,其中黑点表示标记为1的角点。对于8个角点的标记都为1或者都为〇的体素,它属于“。”号构型,没有等值面穿过该体素。当只有一个角点标记为1时,即“1”号构型,我们用ー个三角片代表体素内的等值面片,它将该角点与其它7个角点分成两部分。对于其余儿种构型,将产生多个三角面片。在Lorensen的处理中,首先建立一个“构型一三角剖分”查找表,它包含256个索引项,每个索引项包含索引,旋转以及指向15种三角剖分中的ー种的指针。索引旋转三角剖分模式(构型) 42其中,索弓l(Index)是体素8个角点标记的有序二进制编码,即匕V,V2匕匕匕匕匕Lorensen的方法是:对于每个体素,首先根据它的索引在“构型一三角剖分”查找表中确定其三角剖分形式,然后再根据相应索引项中的旋转参数确定最终的三角剖分。图4-6立方体中等值面的I5种构型1.求等值面与体素边界的交点在确定立方体内三角剖分模式后,就要计算三角片顶点位置。当三维离散数据场的密度较高时,即当体素很小时,可以假定函数值沿体素边界呈线性变化,这就是MC 43算法的基本假设。因此,根据这一基本假设,可以直接用线性插值计算等值面与体素边的交点。对于某棱边,如果它的两个端点匕、匕标记不同,那么等值面一定与此棱边相交。1)体素棱边与X轴平行时,设该边的两端点为K(i,j,k),V2(i+1,j,k),则交点为V(x,j,k)o其中:(4-13)2)体素棱边与丫轴平行时,设该边的两端点为・(i,j,k),V2(i,j+1,k),则交点为V(i,y,k)o其中:(4-14)3)体素棱边与Z轴平行时,设该边的两端点为匕(i,j,k),V2(i,j,k+1),则交点为丫(i,j,z)o其中:(4-15)ね[c—/(匕)]t/(v,)-/(v2)]求出了等值面与体素棱边的交点以后,根据索引我确定的三角剖分,即可将这些交点连接成三角片ー,得到该体素内的等值面片。1.等值面的法向计算为了利用图形硬件显示等值面图像,必须给出形成等值面的各三角面片的法向,选择适当的光照模型进行光照计算,生成真实感图形。对于等值面上的每一点,其沿面的切线方向的梯度分量应该是零,因此,该点的梯度矢量的方向也就代表了等值面在该点的法向。而且,等值面往往是两种具有不同密度的物质的分界面,因而其梯度矢量不为零值,即g(x,y,z)=Vf(x,y,z)(4-16) 44直接计算三角面片的法向是费时的,而且,为了消除各三角面片之间明暗度的不连续变化,只要给出三角面片各顶点处的法向并采用哥罗德(Gouraud)模型绘制各三角面片就行了。我们采用中心差分计算出体素各角点处的梯度,然后再一次通过体素棱边两个端点处梯度的线性插值求出三角面片各顶点的梯度,也就是各顶点处的法向,从而实现等值面的绘制。体素角点的中心差分公式如下:/(七+|,匕,Z*)一/(七メ刀く)8x~2Ar__/(七,力+1,z*)-f(xt,yH,zk)あ2国/(七,匕,Z«+])—/(X"匕,)Sz=(4-17)其中,Ax,Ay,Az分别是体素的边长。1.MC算法抽取等值面的算法流程MC算法求等值面的步骤:(1)将三维离散规则数据场分层读入;(2)扫描两层数据,逐个构造体素,每个体素中的8个角点取自相邻的两层;(3)将体素每个角点的函数值与给定的等值面值c做比较,根据比较结果,构造该体素的索引表;(4)根据索引表得出将与等值面有交点的边界体素;(5)通过线性插值方法计算出体素棱边与等值面的交点;(6)利用中心差分方法,求出体素各角点处的法向量,再通过线性插值方法,求出三角面片各顶点处的法向;(7)根据各三角面片各顶点的坐标及法向量绘制等值面图像。MC算法描述述:(算法4-1) 45for(k=l;k<;k++){读入k-l,k,k+l,k+2四层数据点值{for(j=l;j 46433.2MC算法生成等值面的近似表达与连接的二义性ー、MC算法构造的三角面片是待求等值面的近似表示在MC算法中,等值面与体素边界的交点是基于函数值在体素边界上作线性变化这ー假设而求出来的。当数据场密度高、体素很小时,这ー假设接近于实际情况。例如,由CT或MRI得到的医学图像就属于这种情况。但是,在稀疏的数据场中,体素较大,如果仍然认为函数值在体素边界上具有线性变化,将会产生较大误差。这时,需要根据不同的应用背景对函数值沿体素边界的变化作其它适当的假设,才能较准确地求出等值面。另外,即使函数值沿体素边界作线性变化这ー假设符合实际,通过线性插值求出的交点位置是准确的。但是,将体素中同一个面上两条相邻边上的交点简单地用直线连接起来也是ー种近似,如图4-7所示。⑻(b)图4-7MC算法生成等值面的近似表达由421节分析可知,以方向无关的三线性插值模型生成的等值面为三次曲面。由式(4-10)、(4-11)可以方便地求出等值面与体素边界面的交线方程。不失一般性,设某边界面所在平面的方程为:z=z。,代入式(4-11)得:b0+btx+b2y+byXy=c(4-18)其中毎=%+。5こ0仇=4+ム2〇 47b2=a2+a5z0b3=a4+a7z0显然,方程(4-18)表示的是ー对双曲线,如果用一条直线来表示这条双曲线,则会引起误差(如图4-7所示)。对于体素很小时,如CT、MR【等医学图像,当Z方向间距较大时,可以用断层间插值方法来减小Z间距,这种近似引起的误差是可以忽略的。二、等值面连接方式上的二义性在MC算法中,在体素的ー个面上,如果值为1的角点和值为〇的角点,分别位于对角线的两端,那么就会有两种可能的连接方式,因而存在着二义性,如图4-8所示。这样的面称为二义性面,包含1个以上的二义性面的体素即为具有二义性的体素。在图3.2中的15种构型中,第3,6,7,10,12,13等6种构型是具有二义性的。MC算法存在着连接方式上的二义性问题在Lorensen提出此算法后不久,就由MJ.Durst在ー篇文章中提了出来。(“)ぽ抜庁ズ二义隼的:小ぎボ 48433.2MC算法二义性的消除从上面分析可知,由于体素间存在二义性的面上连接不当可能出现拓扑不一致,会产生等值面的拓扑不一致,从而使等值面出现孔隙。如构型3和构型6是两个相邻体素,它们的公共面是二义的。如果它们的三角剖分分别采用3a和6a(或者3b和6b),那么在它们的公共面上就会导致矛盾的拓扑流型。所构造的等值面将会出现孔隙。如果它们分别采用3a和6b或者3b和6a三角剖分方式,则会消除孔隙。(a)(b)图4-9拓扑不・致造成空隙只有解决二义性面上的不确定性,才能保证所构造的等值面具有拓扑ー•致性。目前,人们在Lorensen算法基础上做了不同的改进,以消除等值面构造中的二义性。渐近线判别法是G.M.Nielson等人提出的,该方法最为常用。本文的重建系统中亦采用该方法来消除等值面连接方式的二义性。下面详细论述渐近线判别和消除二义性的方法。图4-10双曲线与体素边界面的相互位置关系正如式(4-18)所表示的,在一般情况下,等值面与体素边界面所在平面的交线是双曲线。该双曲线的两支及其渐近线与体素的一个边界的相互位置关系可用图4-10来表示。在该图所列的4种状态中,当双曲线的两支均与某边界面相交时,就产生了连接方式的二义性。此时,双曲线的两支将边界面划分为3个区域,可见,双曲线中两条 49渐近线的交点必然与边界面中位于对角线上的ー对交点落在同一个区域内。方程4-19所表示的双曲线的两条渐近线的交点坐标为a2+a5Z0J«4+«7Z0y=_i_"6二0_(4-19)ゝa4+。7z0当出现二义性时,需要计算f(x,y,z。)的值。如果f(x,y,q)",则渐近线的交点应与其函数值大于c的ー对角点落在同一区域内。如果f(x,y,z0) 504.4MarchingCubes算法的实现过程4.4.1开发平台及方法ANSIC语言和OpenGL三维图形软件包为主要编程语言。因为OpenGL软件包已经成为三维领域中的工业标准。•些基本操作如难度很大的消隐由软件包本身帮助完成。你只需调用函数完成某项操作。可供选择的开发途径主要有二条:1.处理二维断层图像生成某・ー标准文件格式,或转换成某ー标准文件格式然后利用现有的可视化系统显示三维图像。这么做工作量小,功能丰富强大,但是必须购买软件。2.处理二维断层图像生成某一标准文件格式,并提供格式转换工具。利用图形语言OpenGL自我编写显示程序,特点是功能不是很强大,需要不断完善功能。这条正是我们选择的。4.4.2体数据处理过程体数据处理的过程分为如下4个部分:1.数据获取使用CT机的脑部CT断层图像,共10层,间距为10mm。对于MarchingCubes算法来说,10mm的间距太宽.必须插值。-般线性插值即可.更深ー步的插值对数据的影响不。MarchingCubes算法对不同的数据精度将产生不同重建效果。精度越高,效果越好。2.图象处理常规处理有平滑、去噪、剪切等。CT的二维图像中有机床的影子,在处理头部图像时如果不剪裁掉.将影响最终的图像。3.表面构造选用MarchingCubes算法生成三维数据集合。生成的数据主要有faceline.dat 51点坐标集合;index.dat点的索引关系集合:normal.dat点的法向量集合1.显ホ选用VC++4.0为开发平台,OpenGL三维图形软件包。4.4.1算法具体实现体元大小的一般为实际物体的1/30-1/50,所以如果是256个单位长度的话,那么CUBE的大小以4X4X2为宜。slice.c具体负责读人断层图像,并提供dividing-slice按2的次基减少数据量。由于各CT生产厂家的CT和MRI文件格式不同,所以读入断层图像的例程提供了对不同厂家数据文件的支持。读入的图像可以是有文件头的CT数据或无文件头的CT数据。poly.c输出生成的三个数据集合;faceling.dat,index.dat,normal.dat〇在输出中index.dat的格式是ANSH码形式,每一行有三个点的索引号,这样存储数据比较大,可以存储为更紧凑的格式。现说明如下:如原来的格式是659,697,1086,10611,可存储为6957,106811,只要读写时相应修改一下读写例程。如果上下两行没有相关联项,则用某ー约定的特殊数值开头,表明无关联的数据。三个输出文件的格式如下:第一行:数据值表明有多少行数据。第二行:实际的数据值。一直到文件结尾。main.c是MarchingCubes算法的具体实现。计算法向量时提供了6点和26点两种方法,视需要的精度而定。command.c可以在屏幕上输出帮助内容.处理从命令行输入的参 52数。cube.c计算CUBE的索引值。设计时指定了高低两个阈值,这样就可以把不感兴趣的高端的值去掉。在实际选取低端阈值时取得值是!75,如果高端阈值取得过低.可能导致重建后的图像有漏洞。需要不断调整阈值以取得最佳效果。cube.h预定义MarchingCubes的索引表,通过查索引表,计算速度将加快。在解决算法的缺陷时由基表的I5种情况通过循环在每次执行程序时生成256种情况的具体的等势面通过CUBE时的边和三角片的情况。hmc.h预定义数据结构及全局变量。在command.c中的需要指定的参数有:两个高低阈值、数据的路径、层间距、输出文件名、减少数据量参数、缩放数据的输出范围参数、数据的扩展名、需要处理的层数、起始层。在cube.h中预定义的索引表有:edge-rotation[24][12](ー个CUBE有六个面,每个面有上下左右四个旋转方向,所以有24个入口。12表示12条边的情况);init-facet-data[15][12](15种情况中有交叉点通过的12条边的情况);init-facet-index[15][8](15种情况中8个顶点的情况);facet-data[256][12]及edge-data[256][l2]在主程序中生成。4.4.1算法流程图算法的数据流图如图4-12所示。生成(多面体的)面genfacetdat() 53姑命令行处理commandline()<0t初始化边数据initedgedate()〇读人数据read_any_slice()Do_marching_cubes()写入数据化writeanypoly()图4-12算法的数据流图4.4.1MarchingCubes算法主要功能模块应用程序中用到的主要功能模块如下:FUNCTIONNAMES 54voidpaint_muscles//displaysthemuscles.voidpaint_polyline//paintthepolyline.paint_polygons//displaythepolygons.calculate_polygon_vertex_normal//calculatethevertexnormalcalcnormal//calculateanormal4.4小结本章对基于面绘制的三维重建方法进行了讨论,论述了切片级重建和体素级重建的基本原理和实现方法。对于切片级重建给出了提取轮廓线的各种实现方法,并对结果进行了比较。对从ー组平行的平面轮廓重建三维形体的原理及实现方法进行了论述,着重讨论了轮廓拼接问题。对于体素级重建简要给出了规则数据场的可视化中体素模型的定义及等值面的定义。论述了MarchingCubes算法抽取等值面的基本原理及实现方法。分析了MC算法中连接方式产生二义性的原因,给出了用渐近线判别法消除连接方式二义性的具体实现方法。运用MC算法对医学图像分割出的组织或器官重建出其表面儿何模型,并对重建结果进行了分析。本章论述了MT算法的基本原理及实现方法。最后给出了采用MC算法重建的结果图。 55第五章总结与展望5.1全文工作总结本文主要研究了基于CT或MRI断层序列图像构建组织或器官的三维几何模型,对医学图像三维重建的・些关键技术进行了研究。在第二章中对三维重建的各种方法进行了分类和评价。然后对基于表面绘制的三维重建方法进行了详细论述,并实现了一个基于MC算法的三维重建,现对全文的工作总结如下。I.通过查阅文献,对医学图像三维可视化的各种方法进行了分类和讨论。目前图像三维重建的方法主要有两大类:ー类是通过几何单元拼接拟合物体表面来描述物体三维结构的,称为基于表面的三维面绘制方法(SurfaceFitting),又称为间接绘制方法;另ー类是直接将体素投影到显示平面的方法,称为基于体数据的体绘制方法(DirectVolumeRendering),又称为直接绘制方法。然后对这两大类方法又可细分处很多重建算法,对这些重建方法的优缺点和适用条件都进行了讨论。2.由于从医疗设备获取的图像不可避免的存在噪声和失真,对输入图像需要进行滤波处理。本文在第三章中设计了邻域平均法、中值滤波法、改进的中值滤波法及保持边缘滤波法的几种模板,可根据输入图像的特征及分割与重建的需要选择使用。设计了两种断层插值方法,将断层图像序列插值成各向同性,然后封装成规则体数据。定义了体数据的内存记录方式及外存文件格式,压缩存储空间。3.对组织或器官的正确分割是重建模型能准确表达组织或器官的关键。本文提出并实现了三维医学图像分割的方法,分割的技术路线是:先分析断层图像,交互给定分割闭值,对图像二值化,然后选择适当操作进行区域修整,最后用种子填充的方法填充出所要区域。此分割可获得满意的分割结果。4.实现了MarchingCubes算法构造表面模型。由分割提取的区域给定等值面值构建的表面模型能准确表达组织或器官。最后给出了采用MC算法重建的结果。 565.2课题研究发展方向1.由MC算法抽取等值面构造的三维模型包含数量巨大的三角面片,可以继续通过模型网格简化来减少模型的三角面片数量,提高实时绘制能力。2.由于辅助诊断及治疗规划的需要,不仅要能对所重建的组织进行旋转、缩放等操作,施以外表的重建组织隐藏显示或半透明效果处理,使内部的组织能清楚地显示出来,以确定内部组织所处的三维空间位置;还要对重建组织施以剖切、开窗等处理,以确定给定切面对选定组织实施剖切的切面形状及切面大小,或通过对外表组织的开窗观察内部组织。因此,将来的系统还要把对重建组织的剖切开窗加进去。3.本文主要讨论的都是基于表面绘制的重建方法。在将来的继续研究中。可以把表面绘制与直接体绘制相结合,以表面特征表现组织或器官的层次结构,同时保留内部细节信息。利用混合绘制来绘制出表面特征较清晰、层次较分明的组织器官。4.虚拟手术工作环境,未来的医学图像可视化将与虚拟现实技术相结合,创造ー个虚拟的环境,让操作者像置身于真实世界中一样。在医学实践中,医生可通过三维输入设备直接对病人的模型实施各种手术方案,整个模拟过程和实际的手术过程非常接近,这将会帮助医生制定最有效、最安全的手术方案。 57参考文献1.石教英,蔡文立。科学计算可视化算法与系统[M]。北京:科学出版社,19962.荆仁杰,叶秀清,徐胜荣,陈存椿。计算机图像处理[M]。浙江:浙江大学出版社,19903.李晓梅,黄朝晖。科学计算可视化导论。国防科技大学出版社,19964.钟尔杰。实用数值计算方法。电子科技大学出版社,20035.DonaldHearn,M.PaulineBaker〇计算机图形学。电子工业出版社20036.郭钧锋,蔡元龙,“基于对象内插的形态学及描术及实现算法”,见:张作生主编,《生物医学工程前沿》,中国科学技术大学出版社,1993年9月,85-88.7.童劲松.科学计算可视化现状及发展趋势.计算机科学,1995,22(4)8.薛强,蔡文立,石教英.MarchingBoxes:一个多精度等值面抽取算法.计算机辅助设计与图形学学报.1998,10(1):8-14 58英文原稿Wepresentanewalgorithm,calledmarchingcubes,thatcreatestrianglemodelsofconstantdensitysurfacesfrom3Dmedicaldata.Usingadivide-and-conquerapproachtogenerateinter-sliceconnectivity,wecreateacasetablethatdefinestriangletopology.Thealgorithmprocessesthe3Dmedicaldatainscan-lineorderandcalculatestriangleverticesusing1inearinterpolation.Wefindthegradientoftheoriginaldata,normalizeit,anduseitasabasisforshadingtheimagesproducedfromthegeneratedsurfacemodelsistheresultofmaintainingtheinter-sliceconnectivity,surfacedata,andgradientinformationpresentintheoriginal3Ddatamodels.Resultsfromcomputedtomography(CT),magneticresonance(MR),andsingle-photonemissioncomputedtomography(SPECT)illustratethequalityandfunctionalityofmarchingcubes.Wealsodiscusesimprovementsthatdecreaseprocessingtimeandaddsolidmodelingcapabilities.Three-dimensionalsurfacesoftheanatomyofferavaluablemedicaltool.Imagesofthesesurfaces,constructedfrommultiple2Dslicesofcomputedtomography(CT),magneticresonance(MR),andsingle-photonemissioncomputedtomography(SPECT),helpphysicianstounderstandthecomplexanatomypresentintheslices.Interpretationof2Dmedicalimagesrequiresspecialtraining,andalthoughradiologistshavetheseskills,theymustoftencommunicatetheirinterpretationstothereferringphysicians,whosometimeshavedifficultyvisualizingthe3Danatomy.Researchershavereportedtheapplicationof3Dmedicalimagesinavarietyofareas.Thevisualizationofcomplexacetabularfractures,craniofacialabnormalities,andintracranialstructureillustrate3D'spotentialforthestudyofcomplexbonestructures.Applicationsinradiationtherapyandsurgicalplanningshowinteractive3Dtechniquescombinedwith3Dsurfaceimages.Cardiacapplicationsincludearteryvisualizationandnongraphicmodelingapplicationstocalculatesurfaceareaandvolume.Existing3Dalgorithmslackdetailandsometimesintroduceartifacts.Wepresentanew,high-resolution3Dsurfaceconstructionalgorithmthatproducesmodelswithunprecedenteddetail.Thisnewalgorithm,calledmarchingcubes,createsapolygonalrepresentationofconstantdensitysurfacesfroma3Darrayofdata.Theresultingmodelcanbedisplayedwithconventionalgraphics-renderingalgorithmsimplementedinsoftwareorhardware. 59Afterdescribingtheinformationflowfor3Dmedicalapplications,wedescriberelatedworkanddiscussthedrawbacksofthatwork.Thenwedescribethealgorithmaswellasefficiencyandfunctionalenhancements,followedbycasestudiesusingthreedifferentmedicalimagingtechniquestoillustratethenewalgorithm*scapabilities.Medicalapplicationsof3Dconsistoffoursteps.Althoughonecancombinethelastthreestepsintoonealgorithm,welogicallydecomposetheprocessasfollows:1.Dataacquisition.Thisfirststep,performedbythemedicalimaginghardware,samplessomepropertyinapatientandproducesmultiple2Dslicesofinformation.Thedatasampledependsonthedataacquisitiontechnique.2.Imageprocessing.Somealgorithmsuseimageprocessingtechniquestofindstructureswithinthe3Ddataortofiltertheoriginaldata.MRdata,inparticular,needsimageprocessingtoselectappropriatestructure.3.Surfaceconstruction.Surfaceconstruction,thetopicofthispaper,involvesthecreationofasurfacemodelfromthe3Ddata.Themodelusuallyconsistsof3Dvolumeelements(voxels)orpolygons.Usersselectthedesiredsurfacebyspecifyingadensityvalue,Thisstepcanalsoincludethecreationofcutorcappedsurfaces4.Display.Havingcreatedthesurface,thefinalstepdisplaysthatsurfaceusingdisplaytechniquesthatincluderaycasting,depthshading,andcolorshading. 60中文翻译我们提出了一种叫做移动立方体的新的算法,它能由3D连续密集的表面生成三角形模型。通过采用分割法使内部切片连通,我们创作了一种定义三角形拓扑结构的体格。这种算法通过顺序扫描和线性插值计算三角形顶点来处理医学三维数据。我们得到原始数据的梯度,把它初始化,并以之作为描绘模型的基础。Thedetailin从已生成的表面模型中产生的图像细节是保持内部切片连通性,表面数据,以及由初始三维数据体现的梯度信息等的结果。来自X断层,磁共振,单光子发射的计算结果阐明了移动立方体算法的特点和功能。同时我们也讨论了其改进来减少处理时间,并增强立体建模的能力。解剖学中的三维表面方法提供了一种有价值的医学工具。这些表面图像,通过诸如CT,MR,SPECT等多种二维切片的重建,帮助内科医生了解复杂的切片解剖情形。解释这二维医学图像需要专门的训练,尽管放射线学者有这些技能,不过他们必须常常和相关的内科医生交流他们的见解,因为内科医生不时难以用肉眼观察三维的解剖。研究者报道说医学三维图像运用在广泛的领域。复杂关节的碎裂,颅面的畸形,颅骨内结构等的可视化阐述了三维方法在复杂骨骼结构研究中的潜カ。而在放射理论和外科移植方面的运用显示了交互式的三维技术结合了三维表面图像。心脏方面的运用则包括动脉可视化和用以计算表面与体积的非图形建模运用。现有的3D算法缺乏细节并且有时候会引入虚假。我们提出了一种新的,高效解决的3D表面重建算法来生成带有以前所未有的细节。这种称为移动立方体的新算法,生成多边形,来表示出来自三维数组的连续密集的表面。所得的模型可以通过由软件或硬件实现的传统图形学描述算法来体现。在描述了三维医学运用的数据流之后,我们开始描述与之相关的工作和讨论工作中的缺点。然后我们描述算法以及其效率,功能的增强,之后用研究的事例,采用三种不同的医学图像技术去阐述新算法的能力。三维图像的医学应用由四个步骤组成。尽管我们可以把后三部结合成ー块,但我们仍把处理过程作如下分解:1.数据的获取这是第一步,由医学图像硬件完成。数据样本依赖于数据获取技术。2.图像处理ー些算法采用图像处理技术来得到3D数据结构或过滤原始数据。特别是像MR 61数据,需要图像处理来选择适当的结构。1.表面重建表面重建,是这篇论文的主题,涉及到用3D数据来生成表面模型。这类模型通常由3D体元素或多边形组成。.用户依靠指定一个亮度值来选择需要的表面,同时这ー步也包括表面的分割和遮掩。2.显示生成表面,最后一步是显示。利用显示技术来体现表面,包括光线投射,深度以及颜色的描绘。 62在论文即将完成之际,我首先要感谢我的老师陈雷霆副教授。感谢陈老师在学术上给予我的悉心指导,提供给我ー个很好的学习环境以及实践锻炼的机会。在此,谨向培养、帮助和关怀我的陈老师致以由衷的敬意和诚挚的谢意。感谢实验室研究生何晓乾对本课题研究工作的帮助和支持,为本论文的研究提供了大量的资料和设备。还要感谢所有教育过,帮助过和指点过我的教研室老师、同学,感谢他们的培养和教诲,感谢同学们在朝タ相处的日子里给予的友谊与支持。张剑2005年6月
此文档下载收益归作者所有