1、数字图像处理学数字图像处理学 第7章 图像重建(Image Reconstruction)医学信息与图像教研室 李明彩?图像重建核心:由一系列投影重建一幅图像 理论基础:二维或者三维物体能够通过其无限多个投影确定。?二维图像重建 一个物体的多个轴向投影图重建目标图像?三维物体重建 由物体的图像重建三维物体模型?问题:能否从投影中恢复原图?答复是肯定的。引言引言 投影显示的仅仅是立体的“横断面”设想一个平头改锥(螺丝刀)的头部 投影只截取到物体所具有的信息中属于较低维度的一部分信息。例如“测量身高”,是将三维降低为一维的投影。?我们所见到的是否是真实的立体?位于眼球深处的“视网膜”是接受来自外界
2、光线的一个屏幕。左右眼球相隔一定的距离,同一物体映照在左眼和右眼视网膜上的二维图像并不相同。大脑会根据这种微小的差异补上远近信息。本章内容本章内容 计算机断层(CT)原理 2 傅里叶切片定理 4 投影和雷登变换 3 概述 3 1 使用平行射线束滤波反投影重建 5 1.图像重建原理的直观定性描述:假设右图是人体三维区域的一个横断面。图中背景表示均匀的软组织,它所环绕的物体是一个肿瘤,肿瘤也是均匀的,但有较高的吸收特性。7.1 概述概述 假设用一束细平行假设用一束细平行X 射线从左到右扫描图像平射线从左到右扫描图像平 面,并假设物体比背景吸 收的射线束能量多,利用收的射线束能量多,利用 另一端的另
3、一端的X射线吸收检测射线吸收检测 带将产生红框所示信号。该该图像幅度(亮度)与吸图像幅度(亮度)与吸 收成正比。信号中任意一点都是所穿过的相应空间点的该射线束中单一射线束吸收值的和(射线和)单一投影无法确定射 线路径处理的是单个物 体还是多个物体,基于 此开始重建,沿射线来 的方向把投影穿过该区 域反“涂抹”回去。沿垂直于射线方向复制 横穿图像的相同的一维 信号。这就是反投影 法。复制矩形图像中所有列中的一维信号创建 源源-检测器对旋转对旋转90度。重复刚才的步刚才的步骤,得到第二幅反投影图像。反投影相交处的亮度是各个反投影的亮度的两倍 将两次的反投影结果相加,我们感兴趣的物体已经包含在右图所
4、示的正方形中了,其幅度是单个投影幅度的两倍多。随着投影数量的增加,不相交反投影的强度相对于多个反投影相交区域的强度降低。0度和180度投影互为镜像 注意到“光晕”模糊效果 例1:包含两个不同吸收特性物体的简单平面的反投影 2.在三维重建中的数据形式有:透射模型;发射模型;反射模型 (1)透射断层重建成像透射断层重建成像(Transmission Computed Tomography,TCT):TCT是射线源的射线穿过物体或者人体组织,然后由接收器接收经过物体或人体组织吸收后的剩余射线,由于不同的物质或组织对射线的吸收不同,剩余射线反映了物体内部的状况和人体内部的不同组织,通过这些数据再重建图
5、像,从而发现物体内部欠缺或人体内部病变。如光,如光,x射线。(2)放射断层重建成像(Emission Computed Tomograhpy,ECT)ECT是在物体中注入放射性物质,然后从物体外部检测通过物体后放射出来的能量,由于物体内部不同的物质或人体不同的组织对放射能量有不同的吸收或衰减,由此可以获取不同的数据,用这些数据可以重建出所需要的图像,从而达到检测物体内部的分布情况或人体内部病变的目的。(3)反射断层重建成像反射断层重建成像(Reflection Computed Tomography,RCT)常用于雷达系统。雷达图像通常是由物体反射的回波产生的。例如光线、电子束、雷达,激光或超
6、声波等都可以用来进行这种测定。(4)核磁共振重建图像)核磁共振重建图像(Megnetic Resonance Imaging,MRI)MRI是由于具有奇数个质子或中子的原子核包含有一定的磁动量或旋量的质子,一般情况下质子在磁场中任意排列,当有适当强度和频率的共振场作用于物体时,质子吸收能量并转向与磁场相交的方向。如果此时把共振磁场去掉,质子吸收的能量会被释放并被检测器收集,根据检测的信号就可以确定质子的密度。通过控制共振磁场强度,可检测一条直线上的信号,从而通过该信号数据可以重建出物体图像。7.2 计算机断层原理计算机断层原理 1.传统的X线系统 胸部X线系统把物体放在对 X射线敏感的平板对
7、面并用圆锥形X射线束照射该个体。X射线平板产生 一幅图像,图像上每一点的亮度与 X射线通过该物 体后照射到该点上的 X射线能量成正比。通过多次反投影,并将多次反投影结果相加,可以产生胸腔结构的三维再现。2.计算机断层(Computerized Tomography,CT)原理 断层成像是一种将物体的每一片层完全隔离出 来进行观察的无损检测技术。它用于人体头部、腹部等内部器官的无损伤诊断,其基本方法就是 根据人体截面投影,经过计算机处理来重建截 面图象。通过身体产生的切片得到物体的三维 再现。三维描述就是堆积这些切片得到。主要在 放射医疗、工业检测设备中,显示人体(对象)各个部位断层图象。(1)
8、CT发展历史:CT 理论基础源于1917年,奥地利数学家约翰 雷登(Radon)所发表的论文。证明了二维或三维物体能够通过其无限多个投影来确定,但限于当时技术条件没能实现。45年后,Tufts大学的物理学者埃兰M考玛克重新发现了这些概念,并应于与 CT,在1963和1964年发表了他的最初发现,并说明了如何从不同角度方向得到的X射线图像重建人体横截面图像,给出了数学公式,并构建了一个用于展示其概念的实际 CT原型。英国EMI 公司中央研究所工程师高德弗里 N 豪斯菲尔德(Godfrey N.Hounsfield),经四年努力,在1972年研制成诊断头颅用的第一台计算机 X线断层摄影装置,这一新
9、设备在 1974年5月蒙特利尔(Montreal)召开第一次国际CT会议上,正式命名为“电子计算机断层摄影技术”简称 CT。1979年EMI 公司又研制出全身 CT。1979年为了表章这些学者的贡献,考玛克和豪斯菲尔德共同获得这项技术的诺贝尔医学奖。(2)几代产品简介(按扫描方式的改进):第一代(G1):单束扫描方式,”铅笔型”X射线管与单探测器沿所示线性方向增量式平移,0度-180度内所有期望的角度重复,生成完整投影,反投影获得一张切片图像,堆积切片生成人体截面三维体。G1时间长,一周需3-5分钟左右 第二代(G2):与第一代原理相同,窄角扇束射线束扫描,张角1020度,2030个探测器相配
10、合。扫描时间18秒左右。允许使用多个检测器,因此检测器的平移较少。第三代(G3):广角扇束扫描方式,张角为30度左右,探测器增加到2501000个,射线源和探测器同步旋转扫描,每个角度都产生完整的投影,消除了部分平移需要,扫描时间可缩短2.5秒。第四代第四代(G4):圆环检测器5000个左右,布满整个360。固定不动,X线源旋转扫描。约 2秒,仅射线源转动。G3,G4扫描速度快,主要缺点是造价高和较大的X射线散射。需要比G1,G2扫描器更高的剂量才能达到可以比拟的信噪比特性。第五代第五代(G5):电子束计算机断层:电子束计算机断层(EBCT)扫描器,排除了所有机械运动,使用电磁控制电子束,通过
11、触发环绕病人的钨极板,电子束产生X射线,然后然后X射线被整形为通过病人的扇形射线束,射线被整形为通过病人的扇形射线束,并激发G4扫描器那样的检测器环。第六代第六代(G6)螺旋螺旋CT:G3,G4扫描器使用一种所谓滑动环来配置,不需要在源/检测器和处理单元之间的电气和信号连接。然后,源/检测器对连续旋转对连续旋转360度,同时病人在垂直于扫描的度,同时病人在垂直于扫描的方向恒速移动,结果是连续的螺旋数值,这些数据经过处理后就可以得到各幅切片图像。第七代(G7)多切片CT扫描器:即将问世,采用“厚”的扇形射线束与平行检测器族相配合同时收集CT数据,即三维横截“厚片层”,而不是每个X射线脉冲产生单一
12、的横截面图像。优点是使用的X射线管更经济,从而降低了成本并降低了剂量。3.关于计算机断层成像 对于医学上的应用来说被计算的特性是组织的衰减系数,人体虽大部分软组织是水,但仍有足够的差异以产生不同的衰减系数,这样就可以给出一幅解剖横截面图像。计算机成像示意图如右侧:衰减系数的单位H(豪斯费尔德Hounsfield),豪斯费尔德将线性衰减系数分为 2000个单位,称为 CT值,一个豪斯费尔德等于水的衰减系数的 0.1%,标 度上选择H(水)=0 对于空气 H=-1000 衰减最小(不衰减)骨骼 H=+1000 衰减最大 1000)()()(?水水组织?H X射线经过物体时会发生衰减,不同的物质衰减
13、是不一样的。得到物体的图像最直接的方法是沿Y轴经衰减直接在胶片上成像。这与 X光透视是一样的,这样会造成图像的混叠。CT是把物体在Y轴方向划分成小的薄片,薄片的厚度是一个重要的参数,一般为1、2、3、4、5、8、10mm。每个薄片再划分为小的单元,即体素。设某一物体体素对 X射线的衰减系数为 ,体素厚度为d,和 为穿透物体前后的 X射线的 辐射强度。射线遵循如下的衰减定律:在断层扫描时,生成大量的数据,根据该数据再计算出每个体素的衰减系数,然后把这些衰减系数按一定的函数关系显示在屏幕上,这样,就产生了断层图像。de?0?01lnd?0?de101?假如一条直线上有 n个体素,第一个体素的衰减为
14、:第二个体素衰减为:d)(dddeeee212120012?对于第n个体素有:d)(nne?210nnlnd?0211?即:一一般般情情况况探探测测器器只只能能测测到到 ,而而不不能能测测到 ,因此,不能直接记录各个体素的衰减系数。但是,我们可以用数学方法求解衰减系数。学方法求解衰减系数。011lnniind?n?121,n?nnlnd?0211?假如某断层有2X22X2个体素,相应的衰减系数为 ,22211211,?分别从X和Z方向投影,测得的衰减系数为 A,B,C,D,即:DCBA?2212211122211211?从而可以解出 的值来。我们用一定的函数关系在屏幕上显示出来就可以得到相应的
15、断层图像。如果图像的分辨率为 512X512,则图像有262144个独立阵元,需要解 262144元的方程组,计算出值,重建出图像。22211211,?7.3 投影和雷登变换 下面详细推导根据 X射线计算机断层重建图像所需要的数学问题:?笛卡尔坐标系中一条直线可以用斜截式来表示,y=ax+b,也可以用法线来表示-极坐标 cossinxy?平行线束的投影可以由一组直线建模,投影信号中任意一点由沿着直线沿着直线 的射线和给出。工作在连续变量情况下,线求和变成了积分 cossinkkjxy?(,)(cossin)(,)kkjjkf x yxydxdgy?由冲击函数性质,积分只有沿着直线进行才会不为零
16、?给出沿xy平面中任意一条线的 f(x,y)的投影(线积分)的公式就是雷登变换,符号 都可以用来表示f的雷登变化,雷登变换是投影重建的基石。离散情况下:(,),(,)R f x yR fg?1100(,)(,)(cossin)MNxygf x yxy?,(,)x yf x y?是离散变量,如果固定 令 变化,投影简化为沿着这两个参数的特定值定义的线对的像素求和。通过增加覆盖图像的所要求的所有的值产生一个投影,改变 并重复这个过程产生另一个投影,类推.例2:使用雷登变换得到圆形区域的投影 圆的解析式 222,(,)0A xyrf x yAr?,其他常数,物体半径(,)(,)(cossin)(,)
17、()=(,)00gf x yxydxdyf x yxdxdyfy dyy?圆是对称的,只需求0度投影 2222222222(,)(,),(,0),(,)0,(,)(,)2,(,)()0,rrrrgfy dyLr gr gfy dyAdyA rrgg?沿直线积分其他不管 角是多少,投影都一样,说明 图像关于原点对称。?当雷登变换当雷登变换 以以 为直线坐标显示为一幅图像为直线坐标显示为一幅图像 时,结果称之为正弦图(sinogram),概念上类似所产生 的傅里叶谱。正弦图包括重建 f(x,y)所需要的数据。下 图右图是左图的正弦图,底部一行是矩形在水平方向 上的投影,中间一行是在垂直方向上的投影
18、,水平方 向较窄,物体是对称的,且平行于向较窄,物体是对称的,且平行于 x,y轴,图像平滑说轴,图像平滑说 明物体具 有均匀灰 度。度。(,)g?,?如下图是一幅如下图是一幅Shepp-Logan 幻影图像,是广泛使用的幻影图像,是广泛使用的 合成图像,设计用于模拟脑部主要区域的吸收,包括 几个小的肿瘤。解释这个正弦图是很困难的,其视觉 分析仅限于实际应用,但有时对于开发算法是有帮助 的。CT的关键目的:得到物体的三维表示 1)反投影每一个投影 2)反投影求和产生一个切片图像 3)堆积所有切片图像产生三维物体的再现 公式推导:单点开始:是全部投影 的一个固定 旋转值 的投影,反涂抹或者反投影后
19、,线 上每一点的值是 ,对投影信号中所有 重 复这个过程,角度不变,可得到 (,)kg?j(,)kg?k?(,)kL?j(,)kg?j?j(,)(,)(cossin,)kkkkkfx ygg xy?任意角度均成立,即 对所有反投影图像积分,得到最终图像 离散情况下(,)(,)(cossin,)kkkkkfx ygg xy?(,)(cossin,)fx yg xy?0(,)(,)f x yfx y?0(,)(,)f x yfx y?这种方法形成的反投影图像有时候称为层图(laminogram),是对产生投影的图像的一个近似。从正弦图得到反投影图像的例子,严重的模糊说明直接利用上述推导的公式在实际
20、中是不可接受的,必须在重建中做有意义的改善。傅里叶切片(投影切片)定理:一个投影的傅里叶 变换是得到投影区域的二维傅里叶变换的一个 切片。如图:任何一投影的一维傅里叶变换可以沿 着一个角度提取一条线的 F(u,v)来得到,该角度就 是产生投影所用的角度,得到所有一维投影的切片 组成二维傅里叶变化,反变换求原图 f(x,y)。7.4 傅里叶切片定理 关于 投影的一维傅里叶变换为 2(,)(,)jGged?(,)(,)(cossin)gf x yxydxdy?222(cossin)(,)(,)(cossin)(,)(co)ssi,n)(jjjxyxyGf x yxyedxdydf x yeddxd
21、yf x y edxdy?关于投影的一维傅里叶变换为 2()2()cos;sincos;sicossinn(,)(,)(,)(,)(cos,sin)jxyjux vyuvuvGf x y edxdyf x y edxdyF u vF?表示表示 的二维傅里叶变换。的二维傅里叶变换。(,)F u v(,)f x y傅里叶切片(投影切片)定理:傅里叶切片(投影切片)定理:一个投影的傅里叶 变换是得到投影区域的二维傅里叶变换的一个 切片。如图:任何一投影的一维傅里叶变换可以沿 着一个角度提取一条线的着一个角度提取一条线的 F(u,v)来得到,该角度就 是产生投影所用的角度,得到所有一维投影的切片 组成
22、二维傅里叶变化,反变换求原图 f(x,y)。7.5 使用平行射线束滤波反投影的重建 直接反投影会生成不可接受的模糊结果,如何 消除模糊反投影之前对投影进行简单的滤波。F(u,v)的二维反傅里叶变换为 极坐标表达 2()(,)(,)jux vyf x yF u v edudv?22(cossin)00(,)(cos,sin)cos,sin,jxyf x yFedduvdudvdd?令则 根据切片定理 22(cossin)0022(cossin)0022(cossin)2(cossin)0002(cossin)(,)(cos,sin)(,)(,)(,)(,)jxyjxyjxyjxyjxyf x y
23、ed ded dGed dGed dedFGGd?cossi2n(,)jxydGed?(,180)(,)GG?是一个一维傅里叶反变换的形式,附加 了一个 项,这是一个一维滤波函数,是一个斜坡滤波 器,由于它在两个积分方向上都是扩展到无穷,傅里叶反变 换无定义,所以对斜坡加窗,使得它在定义之外为零,这就 是一个加窗限制带宽的斜坡滤波器。cossn2i(,)(),xyjex yGdfd?2(,)jGed?a图:使用一个盒状滤波器限制带宽后的滤波器的频率域剖面线 b图:空间域表示 c图:汉明窗函数 d图:图a和图c乘积形成的加窗后的斜坡滤波器 E图:乘积的空间域表示(振铃效应)滤波反投影重建方法:滤波反投影重建方法:1)计算每一个投影的一维傅里叶变换 2)用滤波函数 乘以每一个傅里叶变换,就是乘以一个合适的窗(如海明窗)3)得到每一个滤波后的变换的一维反傅里叶变换 4)对步骤3得到的所有一维反变换积分(求和)得到?例3:滤波反投影重建 a图:使用斜坡滤波器的矩形滤波反投影 b图:用汉明窗斜坡滤波器的矩形反投影 c图:放大细节图 d图:放大细节图 E图:无滤波的模糊图 a图:使用斜坡滤波器的头部幻影滤波反投影 b图:用汉明窗斜坡滤波器的头部幻影反投影 c图:无滤波的反投影模糊图