Long Y , Fessler J A , Balter J M . 3D Forward and Back-Projection for X-Ray CT Using Separable Footprints[J]. IEEE Transactions on Medical Imaging, 2010, 29(11):1839-1850.
摘要
在x射线计算机断层扫描(CT)中，迭代重建三维图像的方法比传统的滤波反投影(FBP)具有提高图像质量的潜力。然而，三维锥束前后投影仪的计算负担是迭代方法在x射线CT实际应用中面临的最大挑战之一。此外，投影精度对迭代法也很重要。本文提出了两种新的可分离足迹投影方法，将体素足迹函数近似为二维可分离函数。由于这些足迹函数的可分离性，在探测器单元上计算它们的积分大大简化，并且可以有效地实现。SF-TR投影机在跨轴方向上使用梯形函数，在轴向上使用矩形函数，而SF-TT投影机在两个方向上都使用梯形函数。仿真和实验表明，两种SF投影方法都比目前最先进的距离驱动(DD)投影方法精度更高。SF-TT投影仪比SF-TR投影仪对与大锥角相关的射线更精确。SF-TR投影机的计算速度与DD投影机相近，SF-TT投影机的计算速度约为DD投影机的2倍。
1、介绍
与滤过反投影(FBP)[4]等传统方法相比，用于三维层析成像重建的迭代统计方法[1]-[3]具有许多优点，如可能改善图像质量和减少剂量。它们基于测量统计和物理模型，可以很容易地结合先验信息、系统几何和探测器响应。
统计重构方法的主要缺点是迭代算法的计算时间较长，通常需要最小化某些代价函数。对于大多数迭代重建方法，每次迭代都需要一个前向投影和一个后向投影，其中前向投影大致是Radon变换的离散化求值，后向投影是前向投影的伴随。这些操作是迭代重构方法的主要计算瓶颈，特别是在三维图像重建。前向投影仪方法也用于制作数字渲染x光片(DRR)[5]，[6]。
传统的前向和后向投影仪计算每个层析射线与每个图像基函数之间的相交长度。许多加速这一过程的方法已经被提出，例如[7]-[13]。由于探测器单元的尺寸有限，平均每个探测器单元的交叉长度被认为是一个更精确的建模[14]-[19]。在数学上，它类似于计算每个基函数和一些检测器模糊(如二维矩形函数)的足迹的卷积。
任何投影仪方法都必须考虑到成像系统的几何形状。锥形束几何图形是轴向和螺旋锥形束x射线计算机断层扫描(CT)所需要的。在三维平行光束几何投影空间，有四个独立的指标。射线方向由射线的方位角和极坐标指定，并分别表示射线的方位角和极坐标，并表示a上的局部坐标二维面积检测器。相反，轴向锥束投影空间由三个独立的指标和两个距离参数表征，其中表示光源点与轴逆时针方向的夹角，表示探测器坐标，表示光源到旋转中心的距离，表示等中心到探测器的距离(见图1)。轴向锥束几何是螺旋螺距为零的螺旋锥束几何的一种特殊情况。
锥束几何中层析射线的发散导致像基函数的深度依赖性放大，即接近x射线源的体素比接近探测器的体素在探测器上投下更大的阴影。这种复杂性在平行光束几何中没有出现。因此，设计了许多现有的投影和反投影方法对于三维平行梁几何[16]-[18]，[20]，[21]不直接适用于锥梁几何。
图1 轴向锥束平面探测器几何结构
各种各样的三维圆锥梁几何投影方法已经被提出[5]，[14]，[15]，[22]-[25]。
所有的方法都在计算的复杂性和准确性之间做出了一定的妥协。其中，球对称基函数[15]，[22]在图像表示方面比简单立方体素或其他基函数有许多优点，例如，它们的外观与观看角度无关。然而，计算它们的足迹函数的积分需要大量的计算。Ziegler等人将这些积分存储在一个查找表中。如果使用优化的blob，并且需要较高的精度，由于加载一个大的表，并且blob与更多的层析射线相交的事实，向前和向后投影的计算仍然是昂贵的。
为了加快圆锥光束前后投影的计算，引入了校正技术[24]。Riddell等人[24]将原始数据重新采样到与两个重建体积主轴对齐的平面上，这样原始的锥束几何结构就可以被一个只涉及连续平面放大的更简单几何结构所取代。在迭代方法中，重采样测量可以简化每次迭代的前向和后向投影。然而，重采样涉及插值，可能会轻微降低空间分辨率。该方法的另一个缺点是，原始投影数据样本的统计独立性通常假设在校正后不再成立，因为插值引入了统计相关性。
距离驱动(DD)投影仪[14]是目前最先进的方法。它将图像体素和检测器单元格的水平和垂直边界映射到一个公共平面(如或平面)上，并用矩形逼近它们的形状。(这一步类似于整改。)它计算沿(或)方向和沿方向的重叠长度，然后将它们相乘得到重叠面积。DD投影仪在x射线源的方位角为的奇数倍左右时误差最大，因为在这些角上跨轴足迹近似为三角形而不是矩形。
本文介绍了两种用于三维前向投影和后向投影的新方法，即SF- tr[26]和SF- tt[27]投影机。他们将体素足迹函数近似为二维可分离函数。这种近似对于典型的轴向或螺旋锥束CT几何形状是合理的。这些足迹函数的可分离性极大地简化了它们在探测器单元上的积分计算，并允许SF投影的高效实现。SF-TR投影机在跨轴方向上使用梯形函数，在轴向上使用矩形函数，而SF-TT投影机在两个方向上都使用梯形函数。对于锥面角小的锥束几何，如多层探测器几何，在轴向使用矩形近似是准确的;对于锥面角大的CT系统，如平板探测器几何，使用梯形近似是准确的。
我们的研究表明，两种SF投影机方法都比距离驱动(DD)投影机更准确。特别地，SF方法减少了数的奇数倍左右的误差，SF-TT投影仪比DD投影仪更精确，SF-TR投影仪用于与大锥角相关的体素。SF-TR投影仪的计算速度与传统投影仪相当DD投影机和SF-TT投影机大约慢2倍。
为了平衡计算和准确性，人们可以结合SF-TR和SF-TT投影机，即使用SF-TR投影机来处理与小锥角相关的体素，例如x射线源平面附近的体素，在那里矩形近似是足够的，而使用SF-TT投影机来处理与大锥角相关的体素。
本文的组织结构如下。第二节回顾了锥梁的几何和投影，描述了锥梁三维系统模型。并给出体素基函数锥束投影的解析公式。第三部分介绍了SF投影机，并对SF投影机与DD投影机进行了对比。第四节给出仿真结果，包括SF-TR的精度和速度比较，SF-TT，和DD投影仪作为独立模块和迭代重建。最后，第V得出结论部分。
2、 锥束投影
A.锥束几何
为了表示简单，我们关注平面探测器的轴向锥束几何(见图1)。这些方法很容易推广到电弧探测器和螺旋几何。
源位于以平面旋转中心为圆心的半径圆上的点上。源位置可参数化如下:
其中Dso为光源到旋转中心的距离，Beta为光源点与轴逆时针方向的夹角。为了简单起见，我们提出了x射线理想点源的情况。为了部分解释非理想x射线源，可以修改下面(20)和(26)中的足迹函数。
让(s,t)表示二维探测器平面上的局部坐标，其中s轴垂直于z轴，t轴平行于z轴。二维探测器上的一个点可以表示为：
其中旋转中心到探测器的距离为Dod=Dsd-Dso。射线从P0到P1的方向向量可以表示为:
两个式子分别表示射线从P0到P1的方位角和极角。
积分是沿着线段的(下面是一条射线的定义)：
对于源与探测器之间的一点x(x,y,z)，其投影坐标S为：
投影坐标t为:
连接光源和X的光线方位角和极角为:
B.锥梁三维系统模型
在迭代图像重建的实践中，我们不是对一个连续的对象F(x)进行操作，而是将一个由公共基函数表示的离散化对象投影到笛卡尔网格上，如下所示:
(公式内的参数和一些前提)
大多数投影/反投影方法使用线性模型，忽略了由比尔定律[28]，[29]的非线性引起的“指数边缘梯度效应”。我们在这里采用相同类型的近似。假设检测器模糊坐标h(s,t)是独立于β，不随移动而改变的，并且仅沿s和t坐标作用。那么理想的无噪声投影满足：
p(s,t,β)由公式（7）给出，为锥束投影值，（sk,tl）表示由索引（k,l）指定的检测器单元格的中心，我们提出的方法适用于任意样本（sk,tl），但为了表示和实现的简单性，我们将重点放在均匀间隔样本的情况下：
deta s和deta t分别为s和t的样本间距，用户可选择的参数cs和ct表示检测器的偏移量，例如，对应于cs=1/4对应检测器的偏移量为1/4 [30]，[31]。
将对象的基展开模型(14)代入(15)并使用(7)得到线性模型：
其中系统矩阵A的元素是以下单基函数的圆锥束投影的样本，其中心在c[n]，
“模糊足迹”功能是：
计算体素的足迹也称为“溅射”[32]。
前向投影仪的目标是快速而准确地计算(17)。尽管系统矩阵A是稀疏的，但对于锥束CT中感兴趣的问题大小，即使是预计算和存储非零系统矩阵值也是不切实际的，因此实用的方法(包括我们提出的方法)基本上是动态计算这些值。
我们专注于探测器模糊的一个简单的可分离模型
rs和rt分别表示沿s和t的宽度。该模型考虑了探测器元件的有限尺寸。注意，rs和rt可以不同于样品间距sk-sk-1和tl-tl-1，以说明探测器间隙。
C.体素基函数的足迹
我们将重点讨论三次体素基函数，但也可以推导出其他基函数的足迹的解析公式。三次体素基函数由：
其中
表示指示函数。
将(22)代入(20)，则第n个基函数的锥束投影足迹解析式为：
对于典型的锥束几何形状，射线的极坐标角要比90小得多，因此不需要考虑。结合(18)，(19)和(23)得到锥束CT中立方体素的“理想”投影仪。
3、可分离足迹(SF)投影
计算体素基函数的动态真实足迹(23)和“模糊足迹”(19)是非常昂贵的，因此需要适当的近似“模糊足迹”(19)来简化二重积分计算。
为了探索替代方案，我们用Dso=541mm和Dsd=949mm模拟了平面探测器锥束几何。我们用（23）的解析方法在以下样本位置计算体素的锥束投影，
图2的左列显示了以β=30°为原点时体素deta1=deta2=deta3=1MM的精确足迹函数及其轮廓。
图2的中间列表示当β=0°以（100，150，15）为中心的那些体素，连接光源和这个体素中心的光线的方位角和极角分别为14.3°和2.1°。典型的64层锥束CT几何锥角约为2。图2的右列显示β=0°时以（93，93，93）mm为中心的体素。连接光源和这个体素中心的光线的方位角和极角分别为11.7°和11.5°。典型的40×40cm2平板探测器锥束CT几何的锥角约为12°。前两个真实足迹看起来像二维可分离函数。第三个足迹除了左上角和右下角的小区域外，几乎是可分离的。
受真实脚印形状的启发(见图2)，我们将其近似如下:
表示具有单位最大振幅的二维可分离函数
右边两个式子分别表示s和t中的近似函数。式（25）中l()为qsf的振幅。
对于小基函数和窄模糊h(s,t)，每个探测器单元内与每个基函数相交的射线角度非常相似，因此L比h(s,t)和q平滑得多。将(25)代入(19)得：
不等式利用的事实是L在每个探测器单元上近似一个常数。该值表示检测器单元（sk,tl）的常数，*表示二维卷积。
如果探测器模糊也被建模为可分离的，即，
那么模糊足迹函数(27)有以下可分离近似:
A.振幅近似方法
B.带有梯形/矩形功能的SF投影(SF- TR)
C.带有梯形/梯形功能的Sf投影(SF- TT)
以上A、B、C略，都是算法公式和推倒过程，需要的时候看原文公式即可。
D. SF投影仪的应用
我们使用系统矩阵模型(18)和可分离足迹方法(29)来进行前投影和后投影，这确保了SF前投影和后投影互为精确伴随算子。略…
E. SF与DD的比较
DD方法基本上是在一个公共平面上如xz何yz平面，使用两个方向的矩形逼近体素足迹。它还在检测平面上使用了可分离和移不变的检测器模糊(21)。然而，基于原始检测器模糊映射边界在公共平面上近似可分离检测器模糊不再是移不变的。这似乎防止使用内循环，以提高SF方法的效率。
4、结果
为了评估我们提出的SF-TR和ST-TT投影机，我们将它们与DD投影机进行了比较，DD投影机是目前最先进的方法。我们比较了它们作为单个模块和迭代重建方法的准确性和速度。
A.前后投影仪作为单个模块
B.迭代重建中的前向和后向投影仪
A、B略
5、总结
我们提出了两个新的3D向前和向后投影仪x线CT: SF-TR和SF-TT。仿真结果表明，SF-TR投影仪比DD投影仪精度更高，计算速度相近，SF-TT投影仪精度更高，但计算速度较SF-TR投影仪慢。DD投影仪是特别有利的相对于其他以前发表的投影仪在速度和准确性之间的平衡。SF-TR方法在跨轴方向上使用梯形函数，在轴向上使用矩形函数，而SF-TT方法在两个方向上都使用梯形函数。轴向的矩形近似适用于锥角较小的连续油管系统，例如多层几何结构。梯形近似对于具有大锥角的几何图形更实际，例如平板探测器几何图形。为了平衡精度和计算量，我们建议将SF-TR和SF-TT方法结合使用，即小锥角对应的体素使用SF-TR投影仪，大锥角对应的体素使用SF-TT投影仪。
本文的模型和模拟都考虑了理想点源。对于有限大小的x射线源，会有更多的模糊，SF和DD方法之间的差异可能会更小。将足迹函数近似为二维可分离函数是该方法的关键贡献。由于可分性大大简化了足迹函数积分的计算，在跨轴和轴向使用更精确的函数是可能的，而不会使计算明显复杂化。
SF方法的计算效率依赖于探测器平面的垂直轴与旋转轴平行的假设。如果探测器平面轻微旋转，则需要稍微插补以重新采样到与旋转轴平行的坐标上。
虽然我们在本文中关注的是体素基函数，但二维可分离足迹近似的思想也可以应用到其他轴向和跨轴方向具有可分离性的基函数上，只要选择适当的函数即可。
进一步的研究将解决实现基于图形处理单元(GPU)的SF投影仪[6]、[44]编程技术来提高速度的问题。