[生物学]生物医学工程博士论文答辩.ppt

浙江大学博士学位浙江大学博士学位论文答辩论文答辩答辩人：徐茂林徐茂林导师：汪元美：汪元美教授浙江大学生物医学工程与仪器科学学院2004.3.11 2004.3.11 下午下午2 2：0000论文题目论文题目透射及散射断层成像中若干反演透射及散射断层成像中若干反演算法的研究算法的研究ANewApproachtoReconstructionAlgorithmsforTransmissionandScatteringTomography说明说明：本课题承蒙国家自然科学基金（本课题承蒙国家自然科学基金（No.69972047）资助。资助。文中所用的实测数据由中国计量科学研究文中所用的实测数据由中国计量科学研究院提供。院提供。下图是中国计量科学研究院研制的工业下图是中国计量科学研究院研制的工业CTCT机。机。工业工业CTCT机机 本文立足于二维由投影重建图像理论，本文立足于二维由投影重建图像理论，着重讨论了三大块内容：透射断层成像中着重讨论了三大块内容：透射断层成像中变换法和级数展开法这两大基本算法的改变换法和级数展开法这两大基本算法的改进和创新，包括滤波反投影算法进和创新，包括滤波反投影算法(FBP)的的关键关键新滤波函数的设计；图像局部重新滤波函数的设计；图像局部重建的新方法以及建的新方法以及FBP算法的误差估计等；算法的误差估计等；对级数展开法中迭代算法做出改进，提出对级数展开法中迭代算法做出改进，提出新的对称块迭代算法。在散射断层成像中，新的对称块迭代算法。在散射断层成像中，讨论了类似于讨论了类似于CTCT中的滤波反投影算法，对中的滤波反投影算法，对S.NortonS.Norton提出的背散射成像的解析重建提出的背散射成像的解析重建论文主要内容公式，基于散射的物理过程重新作了推导，公式，基于散射的物理过程重新作了推导，使探测器响应函数有了明确的物理意义；针使探测器响应函数有了明确的物理意义；针对散射问题由于射线对散射问题由于射线“衰减衰减”所导致的问题所导致的问题的的“非线性性非线性性”及及“不适定性不适定性”，给出基于，给出基于散射能量谱的行扫描散射能量谱的行扫描ComptonCompton背散射成像模型背散射成像模型一个大型非线性方程组；用基于多目标一个大型非线性方程组；用基于多目标优化的神经网络算法求得重建问题的解。优化的神经网络算法求得重建问题的解。此外，用滤波反投影方法给出了此外，用滤波反投影方法给出了Laplace方程一般边值问题解的解析表达式，提供了方程一般边值问题解的解析表达式，提供了势函数的一种快速求解算法。势函数的一种快速求解算法。论文主要内容框架图：CT的数学基础：记被测物体对记被测物体对X射线的线性衰减系数为：射线的线性衰减系数为：根据根据Beer定理有：定理有：CT图像重建就是根据以上各个方向的线积分值图像重建就是根据以上各个方向的线积分值（即射线投影）来确定物体内部的组织结构；即（即射线投影）来确定物体内部的组织结构；即反演分布反演分布 。（A）反演算法变换法变换法 用给定的一系列解析式（算子或函数变换）反演公式(A)来重建图像。常用Randon求逆公式的“离散化”形式。级数展开法级数展开法 直接对公式(A)“离散化”，得到一大型稀疏线性方程组：其中：为投影数据；为投影数据；为待求图像向量。级数展开法就是由给定的投影数据级数展开法就是由给定的投影数据Y Y，估计估计出图像向量出图像向量X X；要求所估X X与一组个基图象b1,b2,bJ的组合能足够逼近所希望重建的图象。射线投影可由射线投影可由Radon变换变换:来描述来描述;其中其中（各参数如下图所示）（各参数如下图所示）Radon运用平均值的思想，建立运用平均值的思想，建立Abel积分方程获得积分方程获得求逆公式：求逆公式：Radon变换示意图：Radon逆变换公式中的广义积分是发散的，逆变换公式中的广义积分是发散的，不能直接用来进行图象重建，需要正则化处理不能直接用来进行图象重建，需要正则化处理。由由Fourier中心切片定理，中心切片定理，Radon逆变换公式亦可写为：对滤波器对滤波器作正则化处理作正则化处理，可获得可获得CT图像图像重建的变换法：重建的变换法：1.1.滤波反投影算法（滤波反投影算法（filterbackprojection简写成简写成FBP）：其中：中：称为滤波函数，称为滤波函数，W(R)称称为窗函数。为窗函数。2.2.卷积反投影卷积反投影算法算法（convolutionbackprojection简写成简写成CBP）滤波反投影算法：滤波反投影算法：采用先修正，后反投影的方法。对每个投采用先修正，后反投影的方法。对每个投影函数作滤波处理，得到修正的投影函数，然后影函数作滤波处理，得到修正的投影函数，然后在对修正的投影函数作反投影运算。关键问题是在对修正的投影函数作反投影运算。关键问题是投影函数的修正即：乘以权重因子投影函数的修正即：乘以权重因子|R|W(R)。卷积反投影算法：卷积反投影算法：用卷积方法修正投影函数，然后再作直接用卷积方法修正投影函数，然后再作直接反投影重建图像的方法。反投影重建图像的方法。优点：优点：能快速实现，精度高。一般要求：完全投影数据，均匀采样。一一.滤波反投影算法的改进滤波反投影算法的改进 及其它应用及其它应用新滤波函数及其时频特性分析新滤波函数及其时频特性分析1新滤波函数的表达式及时频特性曲线图新滤波函数的表达式及时频特性曲线图新滤波函数的频域曲线图：新滤波函数的频域曲线图：与几个常用滤波函数的频域曲线比较与几个常用滤波函数的频域曲线比较-A/4,A/4-A/2,A/2黑色曲线：矩形窗黑色曲线：矩形窗红色曲线：新滤波函数红色曲线：新滤波函数蓝色曲线：蓝色曲线：S_L滤波函数滤波函数绿色曲线：绿色曲线：Hamming滤波函数滤波函数时域特性曲线图：时域特性曲线图：新滤波函数新滤波函数矩形滤波函数矩形滤波函数S_L滤波函数滤波函数Hamming滤波函数滤波函数时域特性曲线比较：时域特性曲线比较：2.时域特性曲线图的数值分析：时域特性曲线图的数值分析：以下将精细地讨论新滤波函数的特性，给出其主瓣、以下将精细地讨论新滤波函数的特性，给出其主瓣、旁瓣的幅值和宽度旁瓣的幅值和宽度,并与常用滤波函数作比较。并与常用滤波函数作比较。求解上式，可依次得出主、旁瓣的宽度和幅求解上式，可依次得出主、旁瓣的宽度和幅值，列表如下：值，列表如下：各滤波函数主瓣、旁瓣的幅值和宽度各滤波函数主瓣、旁瓣的幅值和宽度 幅值和幅值和幅宽幅宽滤波函数滤波函数m0m2m4主瓣幅值主瓣幅值主瓣宽主瓣宽第一旁瓣第一旁瓣幅值幅值第一旁瓣第一旁瓣宽宽新滤波函数新滤波函数0.149080.01204340.001362750.149081.99762-0.01517010.300684RL滤波函数滤波函数0.250.031250.005208330.251.55939-0.08233620.322961Hamming滤波函数滤波函数0.08839230.006480620.0007958930.08839232.0577-0.02715680.75845S-L滤波函数滤波函数0.2026420.02343910.003752960.2026421.65596-0.05852090.273321注：注：mk为滤波函数的为滤波函数的k阶矩，定义为：阶矩，定义为：讨论：滤波函数对重建误差的影响主要有以下三滤波函数对重建误差的影响主要有以下三个因素：主瓣，临近旁瓣及远处旁瓣。主个因素：主瓣，临近旁瓣及远处旁瓣。主瓣越高越窄，空间分辨率越好；旁瓣越小，瓣越高越窄，空间分辨率越好；旁瓣越小，数值精度越高，通常密度分辨率就越好。数值精度越高，通常密度分辨率就越好。从时域特性图中及上表中可以看出，新滤从时域特性图中及上表中可以看出，新滤波函数的主瓣高而窄，尤其是远处旁瓣几波函数的主瓣高而窄，尤其是远处旁瓣几乎为零，能够在保证空间分辨率的同时，乎为零，能够在保证空间分辨率的同时，相对地提高密度分辨率。相对地提高密度分辨率。3.数字实验结果：数字实验结果：模拟数据模拟数据评价函数评价函数滤波函数滤波函数归一化均方根归一化均方根 距离测量值距离测量值归一化平均绝对归一化平均绝对距离测量值距离测量值最坏情况最坏情况 距离测量值距离测量值新滤波函数新滤波函数0.0814260.0137110.120148R_L滤波函数滤波函数0.0790580.0155800.129834S_L滤波函数滤波函数0.0787810.0138410.131357Hamming滤滤 波波 函函数数0.0895630.0137060.150031模拟数据重建图像评价测量值比较模拟数据重建图像评价测量值比较实测数据：实测数据：图图3.1.9：采用采用R_L滤波函数滤波函数图图3.1.10：采用采用S_L滤波函数滤波函数图图3.1.11：采用采用Hamming滤波函数滤波函数图图3.1.12：采用新滤波函数采用新滤波函数局部重建局部重建1.新滤波函数用于局部重建的特性新滤波函数用于局部重建的特性由于滤波函数是一偶函数，作分部积分后，有：由于滤波函数是一偶函数，作分部积分后，有：根据根据Fourier分析中分析中Riemann定理，有：定理，有：其中：其中：因此，当因此，当S较大时，成立如下的渐近近似式：较大时，成立如下的渐近近似式：对新滤波函数以及其对应的窗函数有以下结论：对新滤波函数以及其对应的窗函数有以下结论：并且：并且：（*）（*）于是得时域中新滤波函数的渐进展开为：于是得时域中新滤波函数的渐进展开为：根据（根据（*）式，对于任何滤波函数均有：）式，对于任何滤波函数均有：再由（再由（*）式）式不难看出，其他常用滤波函数的衰减不难看出，其他常用滤波函数的衰减都不会超过都不会超过2阶。综上，在阶。综上，在s s较大时，新滤波函数是较大时，新滤波函数是2阶阶衰减的，而其他滤波函数由于衰减的，而其他滤波函数由于H(A)不为零则是不为零则是1阶衰减阶衰减的，因此，新滤波函数远处旁瓣的衰减最快。的，因此，新滤波函数远处旁瓣的衰减最快。时域曲线局部区域放大比较图：时域曲线局部区域放大比较图：新滤波函数的时域特性曲线图局部放大新滤波函数的时域特性曲线图局部放大Hamming滤波函数的时域特性曲线图局部放大滤波函数的时域特性曲线图局部放大2.局部重建原理局部重建原理选择适当的滤波函数，使得滤波函数的旁瓣迅速衰减，选择适当的滤波函数，使得滤波函数的旁瓣迅速衰减，譬如在某一有限区间外很小，以至于几乎为零，则可譬如在某一有限区间外很小，以至于几乎为零，则可以通过如下的公式作重建：以通过如下的公式作重建：是示性函数。与滤波反投影公式作比较：是示性函数。与滤波反投影公式作比较：得出结论：采用满足如上条件的滤波函数，重建区域得出结论：采用满足如上条件的滤波函数，重建区域内任意一点的图像函数值，只需用到该点某一领域内内任意一点的图像函数值，只需用到该点某一领域内的的Radon投影数据，而不必用到全部投影数据。投影数据，而不必用到全部投影数据。局部重建采用数据范围示意图局部重建采用数据范围示意图：本文给出的新滤波函数即具有上述旁瓣迅速衰减本文给出的新滤波函数即具有上述旁瓣迅速衰减的特性，可以根据上图所示的投影数据采集方式作的特性，可以根据上图所示的投影数据采集方式作局部重建。局部重建。3.局部重建实验结果局部重建实验结果模拟数据模拟数据Shepp_Logan仿真图像仿真图像保持原图图像分辨率保持原图图像分辨率放大放大256*256图象分辨率图象分辨率注：注：数据量为数据量为301*200，局部重建区域半径为全局部重建区域半径为全局区域半径的局区域半径的0.5倍。倍。实测数据实测数据全局重建图像全局重建图像保持原图图像分辨率保持原图图像分辨率放大放大256*256图象分辨率图象分辨率注：注：这是用陶磁叶片的实测数据进行的图象重建实验这是用陶磁叶片的实测数据进行的图象重建实验。数据量为数据量为241*576，重建区域半径为全局区域半径的，重建区域半径为全局区域半径的0.36倍，倍，时间约为时间约为2 23 3s s。误差估计误差估计1.1.对相当广泛的目标函数，算法有如下误差估计对相当广泛的目标函数，算法有如下误差估计：分析上式可以得出，并计算得：分析上式可以得出，并计算得：因此，采用新滤波函数重建时，算法本身产生的误因此，采用新滤波函数重建时，算法本身产生的误差是高阶小量：差是高阶小量：2.2.局部重建的二次误差估计局部重建的二次误差估计 引进径向对称函数：引进径向对称函数：故采用上述方法作局部重建产生的二次误差为：故采用上述方法作局部重建产生的二次误差为：结论：CBPCBP算算法法是是CTCT中中标标准准的的方方法法，关关于于它它的的理理论论和和算算法法已已有有很很多多研研究究。同同时时，人人们们认认为为CBPCBP没没有有局局部部性性，似似乎乎不不可可能能直直接接用用来来作作局局部部重重建建。本本文文证证明明，无无论论从从理理论论或或数数据据实实验验来来看看，直直接接用用CBPCBP实实现现图图像像局局部部重重建建是是可可能能的的。所所谓谓直直接接用用CBPCBP是是指指，就就用用CBPCBP的的计计算算公公式式、算算法法，以以至至于于计计算算程程序序。但但是是，要要恰恰当当选选取取其其中中的的滤滤波波函函数数。本本文文选选取取一一种种新新的的滤滤波波函函数数，它它的的旁旁瓣瓣衰衰减减很很快快，经经证证明明在在可可能能选选取取的的滤滤波波函函数数中中，它它的的衰衰减减是最快的。是最快的。用新滤波函数，我们忽略远离计算点用新滤波函数，我们忽略远离计算点的数据，而实现图像局部重建，并给出了的数据，而实现图像局部重建，并给出了误差估计。本文建议的图像局部重建方法误差估计。本文建议的图像局部重建方法简单、快速，与全部数据用简单、快速，与全部数据用CBPCBP重建的图像重建的图像比较仍有较高的空间和密度分辨率。当比较仍有较高的空间和密度分辨率。当ROIROI取的足够大时，重建的正是整个取的足够大时，重建的正是整个上的图上的图像。像。滤波反投影算法的其它应用滤波反投影算法的其它应用考虑考虑Laplace方程的边值问题：方程的边值问题：其中其中D是是R2中有界单连通区域，中有界单连通区域，边界边界适当光滑：适当光滑：由于势函数由于势函数u生成的梯度场是一势场，其积分与路生成的梯度场是一势场，其积分与路径无关，因此，由边界条件则可以利用径无关，因此，由边界条件则可以利用Radon求逆公式求逆公式得得到其共轭函数到其共轭函数：利利用用Cauchy-Riemann条条件件，选选择择适适当当的的积积分分路路线线，得到势函数：得到势函数：M0为为D中任一固定点，中任一固定点，M为为D中的流动点。中的流动点。作数字模拟实验，取如下模拟函数：作数字模拟实验，取如下模拟函数：对应的共轭函数的精确解为对应的共轭函数的精确解为：在给定这些模拟函数满足的在给定这些模拟函数满足的Laplace方程边值条件后，方程边值条件后，利用滤波反投影方法求得其共轭函数，记解为利用滤波反投影方法求得其共轭函数，记解为V：共轭函数精确解共轭函数精确解v滤波反投影算法解滤波反投影算法解vv，v误差灰度图误差灰度图利用利用Cauchy-Riemann条件求得相应的势函数条件求得相应的势函数u，用三维立体图分别描述出势函数及其共轭函数解：用三维立体图分别描述出势函数及其共轭函数解：共轭函数解共轭函数解v空间坐标图空间坐标图势函数解势函数解u空间坐标图空间坐标图注：注：用滤波反投影法求得共轭函数的解较为精确用滤波反投影法求得共轭函数的解较为精确,误误差差r为为0.0018290.001829，计算所用时间也较短，求得空间分辨率，计算所用时间也较短，求得空间分辨率为为256*256的解的的解的BMP图大约需图大约需1秒。秒。二二.级数展开法的若干研究级数展开法的若干研究经典经典ART和和SIRT方法的迭代格式方法的迭代格式数字化图像的数字化图像的Radon变换可表示为：变换可表示为：则图像重建问题描述成向量的形式：则图像重建问题描述成向量的形式：为投影数据；为投影数据；为待求图像向量。为待求图像向量。目的：目的：寻找最优的图象向量寻找最优的图象向量优点：优点：适用于更普通的数据采集方法适用于更普通的数据采集方法,如：沿曲线或沿区域如：沿曲线或沿区域（带状区域）的积分（带状区域）的积分;容易利用先验条件；适于处理不完容易利用先验条件；适于处理不完全全数据问题。缺点：缺点：数据存储量大、重建速度慢、精度不高。数据存储量大、重建速度慢、精度不高。（*）（1 1）加法）加法ART算法迭代格式：算法迭代格式：（2）SIRT法采用逐点校正，其迭代格式为：法采用逐点校正，其迭代格式为：对称块迭代算法对称块迭代算法1.原理简介：利用几何对称性简化运算并改变迭代原理简介：利用几何对称性简化运算并改变迭代方式，这种对称性描述为：方式，这种对称性描述为：即在平面即在平面XOY上任给一条投影射线，总可以找到另上任给一条投影射线，总可以找到另外七条投影线而形成一几何对称图形；这外七条投影线而形成一几何对称图形；这一对称射线组一对称射线组 所穿过的八个对称象素所形成的图形成八卦状所穿过的八个对称象素所形成的图形成八卦状(见上图见上图)，结构十分简单。这一对称性可以用下面结构十分简单。这一对称性可以用下面的群来描述。的群来描述。记记 分别是的关于分别是的关于 轴和轴和 轴的翻转变轴的翻转变换，换，分别是关于直线分别是关于直线 和和 的的 翻转变换；翻转变换；是恒等变换。是恒等变换。构造如下集合构造如下集合:将将 中的乘法定义为变换的合成运算，中的乘法定义为变换的合成运算，那么由那么由翻转变换的性质可知，翻转变换的性质可知，是一个是一个AbelAbel群，群，而且是对称射线组而且是对称射线组 的变换群。的变换群。变变换换群群 不不仅仅是是保保距距的的，而而且且是是保保像像素素格格不不变变 的的(不不 计计 像像 素素 值值)，或或 者者 说说 像像 素素 格格 在在 变变 换换 群群 下下仍仍然然是是像像素素格格，只只是是位位置置不不同同。一一般般的的变变换换群群 例例如如旋旋转转变变换换群群是是保保距距的的，但但并并不不保保持持像像素素格格不不 变。变。若若将将在在变变换换群群 下下的的不不变变性性看看作作 平平面面上上射射线线之之间间的的一一种种关关系系，那那么么它它是是一一等等价价关关系系，它它将将平平面面上上所所有有的的射射线线分分成成了了等等价价类类，每每一一类类就就是是一一对对称称射射线线组组。于于是是，在在 平平面面上上，每每一一射射线线都都属属于于某某一一对对称称射射线线组组，即即所所有有对对称称射射线线组组覆覆盖盖了了 全全平平面面的的射射线线。而而且且两两个个对对称称射射线线组组(除除对对称称轴轴外外)不不含含同同一一射射线线，否否则则这这两两个个对对称称射射线线组组就就是是同同一一对对称称射射线线组组。这这样样，就就赋赋与与了了平平面面上上的的所所有有射射线线以以简简明明的的结结构构，便便与与计计算及分析。算及分析。这七条射线构造如下：这七条射线构造如下：若射线若射线L1穿过的像素点为穿过的像素点为j1，坐标位置为：坐标位置为：N N是图象向量是图象向量X的行的行(列列)象素个数，则其余对称象素象素个数，则其余对称象素的位置分别为：的位置分别为：利用这一对称性，同时可以得到八个利用这一对称性，同时可以得到八个Radon数据及数据及伪投影数据，则迭代式中所需的计算元素均已知，而计伪投影数据，则迭代式中所需的计算元素均已知，而计算量和时间均节省了约七倍，大大加快了代数重建的速算量和时间均节省了约七倍，大大加快了代数重建的速度。在这一基础上，提出对称块迭代算法：度。在这一基础上，提出对称块迭代算法：2.对称块迭代算法的迭代格式：对称块迭代算法的迭代格式：（3 3）对称）对称ARTART算法算法(SARTSART）：）：方程的循环方式变为：方程的循环方式变为：即迭代是按对称射线组的顺序进行的。即迭代是按对称射线组的顺序进行的。其中：其中：，是第是第k个对称射线组的指标集，个对称射线组的指标集，是对称射线组的个数，对称射线组的个数，M是总射线数。是总射线数。（4 4）对称块）对称块ARTART算法算法(SBIRTSBIRT）：）：这里 ，而 是关于对称块指标集 的控制序列，即控制对称块的迭代顺序；表示第 次迭代的权系数。3.数字实验数字实验Shepp-Logan仿真头像模型仿真头像模型 对称块迭代算法重建结果对称块迭代算法重建结果经典经典ART算法算法重建结果重建结果SIRT算法重建结果算法重建结果SART算算法法重建结果重建结果重建误差及时间比较：重建误差及时间比较：重建算法重建算法迭代次数迭代次数dre时间时间(秒秒)SART6o.1148970.07319719.44892925对称块迭代对称块迭代(SBIRT)60.0757810.0643908.36277841经典经典ART60.1168900.07835623.102844128SIRT60.2871150.24337541.630286253注：注：重建图像大小为重建图像大小为256*256，投影数据量为，投影数据量为301*200，即在，即在200个角度下取个角度下取301条投影。条投影。4.实测数据实测数据(a)SBIRTSBIRT(b)SARTSART(c)ARTART这这是是用用陶陶磁磁叶叶片片的的实实测测数数据据进进行行的的图图象象重重建建实实验验。所所给给的的平平行行投投影影数数据据量量是是576281281。图图 (a)、(b)、(c)是是分分别别由由算算法法SBIRTSBIRT、SARTSART及及ARTART重重建建的的256256256256T T图图象象；迭迭代代次次数数均均为为6 6次次，迭迭代代时时间间分分别别为为9595、5757及及301301秒秒。实实验验环环境境与与模模拟拟数数据据重重建建实实验验相相同同。结结果果同同样样表表明明，对对称称迭迭代代算算法法重重建建速速度度快快，而而且且重重建图象具有较高的分辨率。建图象具有较高的分辨率。结论：基于平面上射线和象素的几何对称性，提出了对称ART算法及对称块迭代算法（SBIRT）。模拟和实测数据实验表明，基于这种对称性结构的迭代算法（SART及SBIRT），使重建图象的各项性能指标都优于经典迭代算法，特别是速度快。值得注意的是这一对称性只是平面上射线和象素的几何性质，与重建图象函数的分布无关。因此，基于这一对称性的重建算法，对于所有重建问题都适用。传统ART算法迭代过程对称SART算法迭代过程对称块SBIRT算法迭代过程三三.背散射成像的多准则神经网络背散射成像的多准则神经网络优化重建优化重建Compton散射成像原理散射成像原理Compton散射成像原理示意图 上图是Compton散射成像示意图。设有一单能 光辐射源，其辐射强度为 ，发射的光子能量为E0，由该高能源辐射的光子经准直器后成为一束细光束照射到一样品上。光子与样品物质核外电子发生康普顿散射相互作用，散射后的光子经准直器，由探测器进行计数。探测器的光子计数率为：这里 是入射光线方向上的小线元，dS是入射光束的横截面积，d为探测器对散射点所张开的立体角，是Klein-Nishina微分散射截面，为测量点的电子密度，其中 的大小由放射源准直器的孔径决定，及 的大小由探测器准直器决定。，称为衰减因子：Ls和Ld分别为入射光线及散射光线的路径。此时，探测器的响应（即光子计数率）不仅与散射点的电子密度有关，而且还与入射光线和散射光线经过各点的电子密度有关。这样就使Compton散射成像问题成为一个“非线性问题”。基于能量谱的行扫描背散射成像模型基于能量谱的行扫描背散射成像模型 行扫描示意图 设光源为单色高能光,并选择广角理想点探测器，光源从左到右扫描。若不考虑多重散射效应及背景辐射，则入射光子经第i行物质电子散射后的能谱可写成如下形式：其中：是第i扫描行第k个能级的光子计数率，是探测器第k个能级的能量。表示第(i，j)个象素的电子密度。N是被测截面图象在一维方向的象素个数(设分成NN的网格)；是探测器在第i扫描行所能分辩的能级数。G是一个常数，它是由系统的几何形状和探测器的效率决定的。是Kronecker符号，表示第(i，j)个象素对第k个能级有贡献时为1，否则为0。是第(i，j)个象素的有效散射体积。是探测器对第(i，j)个象素所张的立体角；若 以表示探测器表面在散射方向的投影值，表示探测器到第(i，j)个象素的距离，则 /。现将第i扫描行的电子密度及能谱分别记为向量 及 ，则可写为如下矩阵形式：若记 ，即 是由()构成的准对角矩阵；则成像的数学模型可表为：记JNN是重建图象的象素总数，是能谱的能级总数，可见 是J维向量，是 维向量；而 是 K J维矩阵，其元素记为 。这是一个大型非线性方程组，重建问题即是要从所测能谱数据 由上式反演出电子密度来来。多准则神经网络优化重建多准则神经网络优化重建 如果考虑到多重散射效应及背景辐射等多种随机误差，成像模型应修正为：这里 表示随机误差向量。由于指数衰减因子的存在，上面方程式投影矩阵的条件数一般很大，因此是不适定的，任何测量误差或模型误差在重建过程中将会被放大。为了获得稳定的解，必须用一些适当的先验条件对解空间加以限制；也就是说，我们只能求出该方程在某些准则下的优化解来。上上述述图图像像重重建建问问题题的的正正则则化化求求解解可可表表为为在在如如下下定定义的可行域中义的可行域中 求使得下式求使得下式：成立的图像向量成立的图像向量 。于是图像重建问题转化为一个多目标优化问题。于是图像重建问题转化为一个多目标优化问题。该多目标优化问题可用线性加权法转化为如下单目标优该多目标优化问题可用线性加权法转化为如下单目标优化问题，即化问题，即：由由Kuhn-Tucker条件可知：此单目标优化问题条件可知：此单目标优化问题的最优解的最优解是上述是上述多目标优化问题的非劣解。多目标优化问题的非劣解。对图像重建问题来说，可选用以下三个目标函数：对图像重建问题来说，可选用以下三个目标函数：3.3.能量函数，它可使重建图象具有极小峰值。能量函数，它可使重建图象具有极小峰值。2.2.测量数据与再投影数据的加权平方误差函数：测量数据与再投影数据的加权平方误差函数：1.1.图像的负熵函数：图像的负熵函数：这里这里 是正则化常数。使是正则化常数。使 最小等价于使最小等价于使(-)最大，最大，它可保持重建图像的全局平滑性。它可保持重建图像的全局平滑性。是正则化常数是正则化常数是正则化常数，是正则化常数，是数据中噪声的方差是数据中噪声的方差 下面用人工神经网络（下面用人工神经网络（ArtificialNeuralNetworks简简记为记为ANN）方法求解这一优化问题。方法求解这一优化问题。人工神经网络是对人类大脑系统一阶特性的一种描述。人工神经网络是对人类大脑系统一阶特性的一种描述。它由大量的处理单元即人工神经元，通过适当的方式互连它由大量的处理单元即人工神经元，通过适当的方式互连构成，是一个大规模的非线性自适应系统。简单地讲，它构成，是一个大规模的非线性自适应系统。简单地讲，它是一个数学模型，可以用电子线路来实现，也可以用计算是一个数学模型，可以用电子线路来实现，也可以用计算机程序来模拟，是人工智能研究的一种方法。从数学的角机程序来模拟，是人工智能研究的一种方法。从数学的角度来看，这对应于一个微分动力学系统，即对应于一个常度来看，这对应于一个微分动力学系统，即对应于一个常微分方程组。如果神经网络系统的能量函数还满足微分方程组。如果神经网络系统的能量函数还满足Lyapunov函数条件，则该微分动力系统的平衡稳定点即函数条件，则该微分动力系统的平衡稳定点即为能量函数的极小值点。为能量函数的极小值点。现现在在将将优优化化问问题题映映射射到到一一神神经经网网络络，相相应应的的归归一一化化图图象象变变元元 映射为神经元输出变元映射为神经元输出变元 。这一神经网络的能量函数由下式给出：这一神经网络的能量函数由下式给出：这里这里J是神经元的总数；是神经元的总数；是增益函数。是增益函数。第第一一项项是是神神经经元元之之间间的的交交互互能能量量；第第二二项项与与约约束束惩惩罚罚有有关关，神神经经元元受受抑抑制制的的非非线线性性输输入入输输出出关关系系可可用用n元元函函数数 表表示示，而而 。这这里里，时时，；否否则则 ，是是一一个个较较大大的的惩惩罚罚参参数数，而而 ，。第第三三项是调节项。项是调节项。网络演化的动力学特性由非线性常微分方程组网络演化的动力学特性由非线性常微分方程组来描述，可以论证该网络系统是渐近稳定的，即来描述，可以论证该网络系统是渐近稳定的，即 是上是上述动力学方程组的述动力学方程组的Liapunov函数，亦即：函数，亦即：且当且仅当且当且仅当 时时于于是是该该网网络络系系统统是是渐渐近近稳稳定定的的，即即网网络络最最终终收收敛敛到到它它的的稳稳定定态态，它它对对应应于于 的的全全局局最最小小值值。这这里里 与与前前述述单单目目标标优优化化问问题题的的罚罚函函数数具具有有相相同同的的形形式式(忽忽略略调调节节项项)，因因此此当当 参参数数 较大时，网络的稳态解即是优化问题的非劣解。较大时，网络的稳态解即是优化问题的非劣解。实验结果：实验结果：注：注：假设该头像模型所代表的截面大小为假设该头像模型所代表的截面大小为50mm50mm，重建中被划分成重建中被划分成6464的网格。的网格。(b)b)是多准则神经网络算法两次迭代的重建图像；而是多准则神经网络算法两次迭代的重建图像；而(c)c)是最是最小二乘算法迭代三次的重建图像。小二乘算法迭代三次的重建图像。(a)(b)(c)Shepp-Logan头像模型头像模型(a)及多准则神经及多准则神经网络重建图像网络重建图像(b)与最小二乘重建图像与最小二乘重建图像(c)算法误差分析算法误差分析算法算法Algorithm多准则神经网络多准则神经网络Multicriterianeuralnetwork0.3152261.316479最小二乘法最小二乘法Least-squares0.6478531.101259归一化均方根误差归一化均方根误差 均方根信噪比均方根信噪比 行扫描背散射成像方法，具有扫描速度快，行扫描背散射成像方法，具有扫描速度快，几何位置灵活的特点。在噪声下的重建结果表明，几何位置灵活的特点。在噪声下的重建结果表明，所构造的神经网络对多准则背散射成像是可行的，所构造的神经网络对多准则背散射成像是可行的，适于处理不适定重建问题；性能指标明显优于最适于处理不适定重建问题；性能指标明显优于最小二乘法重建小二乘法重建。结论：结论：由投影重建图像是由投影重建图像是CT成像技术的理论基础，成像技术的理论基础，本文立足于二维由投影重建图像理论，本文立足于二维由投影重建图像理论，对透射及对透射及散射断层成像问题，从成像系统模型的角度，抽散射断层成像问题，从成像系统模型的角度，抽象出其数学的共同点，作了统一的研究和分类。象出其数学的共同点，作了统一的研究和分类。主要是：主要是：1提提出出了了一一种种新新的的滤滤波波函函数数，其其主主要要部部分分是是泛泛函函分分析析中中的的磨磨光光算算子子，它它有有许许多多优优异异的的性性质质，特特别别是是利利用用旁旁瓣瓣的的迅迅速速衰衰减减性性，用用滤滤波波反反投投影影算算法法简简便便快快速速地地进进行行了了图图像像局局部部重重建建，重重建建图图像像空空间间分分辨率高，且能较好地克服辨率高，且能较好地克服Gibbs效应。效应。四四.总结与展望总结与展望2 2利利用用投投影影矩矩阵阵的的几几何何对对称称性性，提提出出了了新新的的对对称称块块迭迭代代算算法法，对对经经典典ARTART及及SIRTSIRT算算法法作作了了较较大大改改进进，使使得得重重建建图图像像的的精精度度高高速速度度快快。由由于于这这一一几几何何对对称称性性与与重重建建图图像像函函数数的的分分布布无无关关，因因此此适适用用于于所所有有的的重建问题。重建问题。3 3针针对对ComptonCompton背背散散射射成成像像问问题题的的“非非线线形形性性”及及“不不适适定定性性”，提提出出了了多多准准则则神神经经网网络络优优化化算算法法，在在权权平平均均意意义义下下，克克服服了了单单准准则则优优化化算算法法偏偏重重单单个个指标的不足，适于处理这类指标的不足，适于处理这类“不适定性不适定性”问题。问题。展望展望对称块迭代对称块迭代局部重建局部重建由由Radon变换似变换似可进一步揭示一般区域上的可进一步揭示一般区域上的Green函数的形式和性质函数的形式和性质Compton背散射成像背散射成像三维三维CT技术：技术：螺旋螺旋CT、核磁共振（核磁共振（MRT）、）、三维正电子成像三维正电子成像（3DPET）、三维单电子反射成像（、三维单电子反射成像（SPECT）等。等。本文所提出的算法思想推广至三维成像，是否有本文所提出的算法思想推广至三维成像，是否有相应的理论和结果，同时是否会产生新的理论或相应的理论和结果，同时是否会产生新的理论或实验结果，这些都有待于展开进一步的研究工作。实验结果，这些都有待于展开进一步的研究工作。Thanks