1、1断层CT64层螺旋CT锥束CT首师大的锥束CT每次可扫描1536层、重建1900*1900*1500的三维CT图像，最高分辨率可达142微米（3.5pl/mm）一、计算机断层成像技术简介1.什么是计算机断层成像技术(Computer Tomography,CT)？广义上讲，计算机断层成像技术是：广义上讲，计算机断层成像技术是：(1)(1)检测一组被测物体的某一物理量检测一组被测物体的某一物理量A A；(2)(2)由该组物理量由该组物理量A,A,通过某种算法求解出通过某种算法求解出被测物体内部另一物理量被测物体内部另一物理量B B的分布的分布;(3)(3)由物理量由物理量B B的分布的分布,进
2、一步推知物体的进一步推知物体的内部结构。内部结构。o X射线射线(射线射线)断层成像技术断层成像技术(即狭义即狭义CT)o 激光激光CT成像成像(Optical Computer Tomography)o 核磁共振断层成像技术核磁共振断层成像技术(Magnetic Resonance Imaging,MRI)o 超声成像技术超声成像技术(Ultrasonic Tomography,UT)o 电阻抗成像技术电阻抗成像技术(Electrical Impedance Tomography)o 红外成像技术红外成像技术(Infrared Imaging)o 2.广义上的计算机断层成像技术 o透射式CT
3、（Transmission CT）衰减衬度CT成像（Attenuation contrast CT,即通常CT，重点讲）相位衬度CT成像（Phase contrast CT,简介）o发射式CT（Emission CT,简介）o散射式CT（Scattering CT，简介）X射线(射线)断层成像技术(即狭义CT)3.CT发展简史(狭义CT)o 1895年伦琴（年伦琴（Roentgen）发现发现X射线射线,1901年获年获Nobel prizes.o 1917年年 德国数学家德国数学家Radon提出提出“天线数学天线数学”，应用于应用于“无线电天文学的图像重建无线电天文学的图像重建”，为成像，为成
4、像理论奠基。理论奠基。o 1967年英国年英国Hounsfield研制成第一台用于临床的研制成第一台用于临床的医学医学CT机。机。o 1974年美国研制成功全人体年美国研制成功全人体CT。o 1979年年G.N.Hounsfield和和A.M.Cormack获获Nobel医医学奖。学奖。o 1980美国等研制成高精度工业美国等研制成高精度工业XCT机，机，90s美超美超大工业大工业CT机用于检测运载火箭发动机机用于检测运载火箭发动机(直径直径1.02.0m).o o 国内国内1980初期开始初期开始CT理论与应用研究，研究单理论与应用研究，研究单位主要有：清华大学、重庆大学、首都师范大位主要有
5、：清华大学、重庆大学、首都师范大学、中科院高能所、核九院电子所、上海交通学、中科院高能所、核九院电子所、上海交通大学，北京大学、北京航空航天大学大学，北京大学、北京航空航天大学4.CT成像特点 o 非接触性成像，适应检测领域宽非接触性成像，适应检测领域宽o 无影像重叠，密度和空间分辨率高无影像重叠，密度和空间分辨率高o 数字图像，易于处理、存储、传输，易于数字图像，易于处理、存储、传输，易于2维、维、3维应用维应用5.CT应用领域 o 医疗生物医疗生物(医疗诊断、骨密度医疗诊断、骨密度)o 工业工业(军工：导弹、炮弹、飞机发动机叶片、军工：导弹、炮弹、飞机发动机叶片、火箭推进器、气缸火箭推进器
6、、气缸.；民用：机械部件、陶；民用：机械部件、陶瓷瓷)o 安全检查安全检查(非金属武器、毒品、炸药、非金属武器、毒品、炸药、)o 材料材料(分子结构、岩心、管道、分子结构、岩心、管道、)o 地学地学(桥梁、隧道、水坝、物探、地震、桥梁、隧道、水坝、物探、地震、)o 天体天体(天体成像天体成像)6.检测对象尺度和分辨率 o 微焦点微焦点CT:CT:对象对象1-30mm,1-30mm,分辨率分辨率0.5-50.5-5微米微米o 医学医学CT:CT:对象对象(人人)300)300-800-800mm,mm,分辨率分辨率0.5mm-3mm0.5mm-3mmo 工业工业CT:CT:对象对象(构件构件)1
7、00-500mm,)100-500mm,分辨率分辨率0.1-0.5mm0.1-0.5mmo 大型工业大型工业CT:CT:构件构件300-1000mm,300-1000mm,分辨率分辨率0.5-1.0mm0.5-1.0mmo 超大型工业超大型工业CT:CT:构件构件1000-2000mm,1000-2000mm,分辨率分辨率1.0mm1.0mm左右左右o 工程工程CT:CT:桥梁、隧道等桥梁、隧道等,分辨率分辨率1 1米米-几米几米o 地学地学CT:CT:地壳、地质构造地壳、地质构造,.,.o 天体天体CT:.CT:.7.的扫描模式圆轨道锥束扫描方式圆轨道锥束扫描方式螺旋轨道锥束扫描方式螺旋轨道
8、锥束扫描方式Proj.image at angle 1Proj.image at angle 2扇束扇束扫描扫描 左图：肺部的左图：肺部的X X光影像（不能分辨深度）光影像（不能分辨深度）右图：肺部的右图：肺部的CTCT断层影像（能分辨深度）断层影像（能分辨深度）传统传统X X光影像与光影像与CTCT断层影像对比断层影像对比虚拟解剖（科研教学）虚拟内窥镜（诊断）虚拟手术（手术设计）由断层重构的三维影像由断层重构的三维影像工业工业CTCT的难点：的难点：(1)(1)涉及的检测对象复杂涉及的检测对象复杂(尺寸大小、材料种类、物尺寸大小、材料种类、物体形状、检测要求的分辨率差异大体形状、检测要求的分
9、辨率差异大)(2)(2)涉及的射线源种类多涉及的射线源种类多/差别大差别大(3)(3)涉及的探测器种类多涉及的探测器种类多/差别大差别大由探测器采集射线衰减量=投影数据(生数据)三维重构=表面、体显示(局部表面、体显示)2维/3维应用等（三维显示技术/图像处理/模式识别/.）通过特定的算法重建切片数据(或体数据)（我们研究的重点:图像重建技术）二、CT成像算法简介 基本原理基本原理 2.1.透射式CT的数学基础(TransmissionCT,TCT)o一个基本的事实一个基本的事实：X X射线穿过物质时射线穿过物质时，射线的强度会，射线的强度会发生衰减。发生衰减。o衰减系数衰减系数取决于取决于物
10、质密度和原子物质密度和原子序数，序数，X射线的射线的能量。能量。假定射线源是单色的(即所有光子具单一能量)，则1对均匀介质，有以下定律 )exp(0 xII (Beer Law)其中：0I是源的 X 射线的强度，I是检测到的 X 的强度；称为介质的线性衰减系数，即经过单位 长度的介质射线强度的衰减值；x是射线经过介质的长度。2对非均匀介质，有 ),(exp)(0dlyxILIL其中：),(yx是介质在点),(yx处的线性衰减系数；)(LI是检测到的 X 的强度；L是射线经过的直线，dl是直线的弧长。dlXfLIILRfLfL)()(ln)()(0Radon 变换变换)sin,(cos),(,)
11、(),(),(21xxXdlXfpRfpfpX其中或CTCT成像问题成像问题已知已知 ,其中其中 ，求求 的近似值。通常称的近似值。通常称 为为投影数据。投影数据。),(jkpfK,k21J,j21)(Xf),(jkpf模型1 投影(Phi=0)501001502002505101520253035投影(Phi=3/4)投影(Phi=)501001502002505101520253035501001502002505101520253035模型1投影数据图pdpdXppfXfp002),)(21)(Radon 反演公式反演公式(1917)实际成像中涉及到的问题实际成像中涉及到的问题o 投影数
12、据有噪声或有坏数据，需要去噪或修投影数据有噪声或有坏数据，需要去噪或修正。正。o 投影数据不完全。投影数据不完全。o 不同的数据噪声和数据缺少，导致影像伪影不同的数据噪声和数据缺少，导致影像伪影。o 射线硬化问题、散射问题。射线硬化问题、散射问题。o 重建速度问题。重建速度问题。o 空间分辨率、密度分辨率。空间分辨率、密度分辨率。o 与系统相关的问题与系统相关的问题(源的焦点尺寸、准直器尺源的焦点尺寸、准直器尺寸、探测器效率、转台偏心、寸、探测器效率、转台偏心、.).)o导算子放大噪导算子放大噪声声oHilbertHilbert变换变换，奇异积分，奇异积分o反投影算子导反投影算子导致伪影致伪影
13、dEdlEXExpESIpIXEpX),()(),(max00例如硬化问题：通常的X射线是多色的（即光子的能量服从一个分布）。在多色情况下，衰减服从1.非线性问题，不适定2.无精确解，需求逼近解3.迭代收敛问题 4.重建速度问题.变换类算法 o 滤波反投影算法：滤波反投影算法：Filter Back Projection(FBP)Algorithmo 直接傅立叶方法：直接傅立叶方法：Direct Fourier Transformation Algorithm迭代类算法迭代类算法o 代数重建算法代数重建算法 Algebraic Reconstruction Technique(ART)Algo
14、rithmo 同时代数重建算法同时代数重建算法Simultaneous Algebraic Reconstruction Technique(SART)o OSEM(Order Subset Expectation Maximization Algorithm)其它类型算法其它类型算法各有优缺点各有优缺点2.2 通常的重建算法通常的重建算法(Reconstruction Algorithms)l 扩展水平成像视野的扩展水平成像视野的RTRT扫描模式及直接重建方法扫描模式及直接重建方法l 扩展垂直成像视野的螺旋扫描模式及直接重建方法扩展垂直成像视野的螺旋扫描模式及直接重建方法l 针对感兴趣区域的
15、扫描模式及直接重建算法针对感兴趣区域的扫描模式及直接重建算法l 两种工业两种工业CTCT超分辨成像方法超分辨成像方法l 几种几种CTCT重建算法的重建算法的GPUGPU加速实现方法，使图像重建速度加速实现方法，使图像重建速度达到国际先进水平达到国际先进水平l CTCT扫描参数的自动获取方法扫描参数的自动获取方法l X X射线能量谱的间接测量方法射线能量谱的间接测量方法l 针对射线束流不稳定、探测器之间的非一致性、对能量针对射线束流不稳定、探测器之间的非一致性、对能量响应的非线性性，提出了数据校正方法响应的非线性性，提出了数据校正方法l 三、三、我们的一些工作我们的一些工作 3.1拓展水平方向视
16、野的研究拓展水平方向视野的研究RT 三次扫描模式三次扫描模式已有的算法需要将数据通过插值重排成平行束或扇形束，已有的算法需要将数据通过插值重排成平行束或扇形束，因此计算量大、分辨率低因此计算量大、分辨率低我们算法的特点我们算法的特点:a)直接重建直接重建,不需要重排不需要重排=速度快，分辨率高速度快，分辨率高b)具有并行计算的特点具有并行计算的特点动机：动机：a)需要用短的探测器对大的物体需要用短的探测器对大的物体扫描成像扫描成像b)长探测器价格高长探测器价格高c)加速器为窄视场源加速器为窄视场源 DBP图像CT 图像部分DBP图像面整列探测器我们方法重建结果我们方法重建结果重排方法重建结果重
17、排方法重建结果垂直成像视野拓展垂直成像视野拓展螺旋扫描模式及其断层直接重建的快螺旋扫描模式及其断层直接重建的快速方法速方法a3.拓展垂直方向视野的研究拓展垂直方向视野的研究Before 2D interpolationPhantomAfter 2D interpolation3.3 感兴趣区域扫描模式及直接重建算法感兴趣区域扫描模式及直接重建算法QxyxIOuOMSyBA转台平移的超分辨扫描模式转台平移的超分辨扫描模式探测器相叠的示意图探测器相叠的示意图平移旋转超分辨扫描模式平移旋转超分辨扫描模式 等价扫描模式等价扫描模式探测器相叠示意图探测器相叠示意图o 虚拟探测器的超分辨算法虚拟探测器的超
18、分辨算法特点：转台转动精度要求不高，近似特点：转台转动精度要求不高，近似算法算法特点：要求转台转动精度高，精特点：要求转台转动精度高，精确算法确算法3.5 CT超分辨成像方法超分辨成像方法o 虚拟焦点的超分辨算法-3-2-11230.20.40.60.81实际焦点和虚拟焦点的强度分布实际焦点和虚拟焦点的强度分布示意图示意图 线对卡线对卡The image reconstructed by FBP The image reconstructed by the virtual-focus-based method The image reconstructed by ART3.6扫描参数自动提取方
19、法（1）实际的扫描机械难以对准；（2）所需参数也无法直接测量。经典图像重建算法要求的对准关系不满足经典对准关系的一般情况12,R D h n n和 12,R D h n n 几何伪影的形态分析几何伪影的形态分析(a)(b)(c)NM95.5 96.0 96.5 97.0 97.51.00.50.51.01.5MN96.096.296.496.60.30.20.10.10.20.30.4MN96.096.196.296.396.496.50.30.20.10.10.20.3300.0300.1300.20.20.40.60.8 如美国如美国Varian公司型号公司型号PaxScan2520的面阵
20、探测器的每帧扫的面阵探测器的每帧扫描图像为描图像为1920*1536，A/D转换为转换为12位，数据输出为位，数据输出为16位。按位。按720个角度采样，数据量约为个角度采样，数据量约为4GB。相应的最大重建图像为。相应的最大重建图像为1920*1920*1536，数据量约为，数据量约为10GB。即使重建。即使重建1024*1024*800的三维图像，数据量仍有的三维图像，数据量仍有1.56GB3.8CT算法的加速问题算法的加速问题CTCT常用的硬件加速方法：常用的硬件加速方法：a)a)多多CPUCPU的并行机的并行机 价格昂贵价格昂贵b)b)SSESSE（Streaming SIMD Ext
21、ensionStreaming SIMD Extension）加速幅度小加速幅度小c)c)图形处理器（图形处理器（GPUGPU）发展速度快发展速度快 高带宽高带宽 高并行性高并行性 价格便宜价格便宜图形处理器（图形处理器（GPUGPU）加速的缺点：灵活性差）加速的缺点：灵活性差它是为图像显示它是为图像显示设计的，需要将算法设计成图形学的绘制形式设计的，需要将算法设计成图形学的绘制形式目前重建速度仍是制约三维锥束目前重建速度仍是制约三维锥束CT应用的主要因素之一应用的主要因素之一传统的方法：主要集中在传统的方法：主要集中在FDKFDK算法中算法中反投影部分的加速反投影部分的加速，其，其中中FFT
22、FFT和滤波等较复杂的计算则放在和滤波等较复杂的计算则放在CPUCPU上实现。纽约州立大上实现。纽约州立大学的学的MullerMuller研究小组的工作具有有代表性研究小组的工作具有有代表性.2007 2007年论文年论文“Real-time 3D computed Real-time 3D computed tomographictomographic reconstruction reconstruction using commodity graphics hardwareusing commodity graphics hardware，Phys.Med.Biol.Phys.Med.B
23、iol.52,34053419(2007)”52,34053419(2007)”中反投影时间是中反投影时间是8.9秒秒的反投影速的反投影速度度GPU加速锥束加速锥束FDK算法研究算法研究我们算法特点：在我们算法特点：在GPU上实现了上实现了FDK重建算法重建算法(包括滤波、包括滤波、反投影和累加全过程反投影和累加全过程)，并支持重建，并支持重建1024立方的立方的CT图像。用图像。用我们的方法由我们的方法由360个个512平方的投影数据重建平方的投影数据重建512立方反投影立方反投影需要需要商用锥束商用锥束CT机的重建速度机的重建速度滤波的GPU实现：实部数据纹理虚部数据纹理矩形光栅化蝶式算法
24、系数及索引纹理利用fragment Program shader多渲染目标技术实部结果纹理虚部结果纹理乒乓技术交换原始纹理和结果纹理FFT 的GPU实 现数据数组数据纹理FFT实部虚部0实部虚部频域滤波实部虚部滤波后的频域数据纹理逆FFT滤波后的数据纹理实部虚部GPUGPU反投影.滤波数据尺寸内存-显存传输耗时GPU滤波计算耗时显存-内存传输耗时GPU 滤波共耗时CPU滤波耗时256*256*3600.422秒1.172秒0.390秒1.984秒42.094秒512*512*3601.094秒4.187秒1.532秒6.813秒196.437秒1024*1024*3604.312秒15.719
25、秒6.564秒26.595秒673.564秒反投影部分的改进反投影部分的改进：渲染结果平行投影透视投影视点投影数据屏幕重建切片渲染结果平行投影R通道透视投影G通道透视投影B通道透视投影A通道透视投影重建切片投影纹理投影纹理改进的投影纹理改进的投影纹理2.4.Figure 2:Geometry of back projection.The slice under reconstruction has each filtered x-ray image projected onto it by projective texture mapping.GPU GPU加速实现的反投影算法加速实现的反投影
26、算法四个通道旋转后累加 角度的投影角度的投影角度的投影2角度的投影TexProj改进投影纹理数据纹理ImaProj重建结果纹理单通道GPU重建结果改进实现方法：角度的投影角度的投影角度的投影2角度的投影TexProj投影纹理数据纹理四通道ImaProj重建结果纹理四通道GPU重建结果数组重建结果现有的实现方法：扫描模式等价模式单层螺旋重建算法的GPU加速：“平行扇束平行扇束”投影模式投影模式：视点 投影平面zl视点 投影平面“平行扇束平行扇束”投影模式投影模式透视投影模式透视投影模式200000020000101nwfnfnfnfnP投影矩阵：投影矩阵：投影矩阵：投影矩阵：2000200020
27、00010nwnfnfnfnfnhP螺旋螺旋CT的的GPU加速，加速，重建体重建体大小大小投影数据投影数据尺寸尺寸内存内存-显存显存投影数投影数据求导据求导反投影反投影耗时耗时 显存显存-内内存存 GPU重建共重建共耗时耗时CPU重建时重建时间间256*256*256(256*256)*3600.406秒秒0.109秒秒1.641秒秒0.328秒秒2.484秒秒63.095秒秒512*512*256(512*256)*3600.600秒秒0.531秒秒5.094秒秒1.984秒秒8.209秒秒131.443秒秒1024*1024*256(1024*256)*3601.078秒秒0.953秒秒18.875秒秒9.359秒秒30.265秒秒411.950秒秒GPUGPU加速的螺旋加速的螺旋CTCT重建重建DrisDris模型模型58 o 锥锥束束CTCT应用软件应用软件 MagicEye-3dCT MagicEye-3dCT o 工业工业CTCT应用软件应用软件 MagicEyeMagicEye-ICT-ICTo X X线数字成像软件线数字成像软件 MagicEyeMagicEye-XDR-XDRo 三维可视化软件三维可视化软件 MagicEye-3dViewMagicEye-3dView
