使用截短的投影和在先采集的3D CT图像的锥形束CT设备

摘要

根据用x－射线装置采集的截短锥形束投影数据和代表对象较大区域的在先CT图像计算对象的一个区域的3D图像。外插截短投影数据以获得与探测器之外的投影方向相关的伪投影数据，并且基于由伪投影数据完成的截短投影数据重建中间CT图像。然后将在先CT图像与中间CT图像配准。根据截短投影数据和已配准的在先CT图像计算与探测器之外的投影方向相关的向前投影数据。最后基于由向前投影数据完成的截短投影数据重建3D图像。

说明书

本发明涉及三维图像处理，尤其在3D x-射线医学成像领域中的三维图像处理。

在本发明的上下文中，测量设备典型是装配有图像增强器或平板探测器的C形臂系统。C形臂系统是普通类型的x-射线成像装置。其特征在于臂形状像“C”，并且在其端部载有诸如x-射线管的x-射线源，以及x-射线探测器。x-射线源具有小的、像点的焦点。源的焦点相应于锥形的顶点。

普通类型的探测器是x-射线图像增强器和CCD照相机的组合。图像增强器具有圆形灵敏区域，其典型直径约为20至40cm。入射的x-射线光子转化成由CCD照相机读出的中间光学图像。探测器还可以是具有矩形灵敏区域的平板探测器，通常尺寸从20×20到40×40cm²。平板探测器的灵敏区域划分成直接读出的小探测器元件的2D阵列。使用这两种类型的探测器时，探测器的最终输出是数字数据的2D阵列，其代表当进行测量时照射到灵敏区域的x-射线光子的空间强度分布的采样版本。

同时，源的焦点和探测器的灵敏区域定义了x-射线的锥形束。锥形束可以被隔板或准直仪变窄。如果x-射线管工作时，在源和探测器之间存在对象，就获得所谓的对象的锥形束投影。这样的锥形束投影提供了3D对象的2D投影图像，并且携带了关于对象内x-射线衰减系数的空间分布的有限信息。

如果对象没有完全地被锥形束照射，那么该投影被称为截短的。截短量取决于探测器的尺寸、对象的尺寸、投影方向以及其他几何因素。如果对象是人体，测量非截短锥形束投影需要无法实施之大的探测器。因而，在医学应用中，测得的锥形束总是截短的。

在由C形臂系统的机械设计强加的限制内，可以将源移动到同心球面上的任意点。这使得可能从定位在同心上的对象的任意方向上采集锥形束投影。源还可以沿着局限于所提及的同心球的一些轨迹以连续的方式移动。这使得可能以连续改变的方向采集定位在同心上的对象的一系列锥形束投影。如果移动仅包括围绕单一轴的旋转，那么产生的轨迹是同心圆弧。通过组合至少两种旋转运动，可以产生非平面源轨迹。

通常以所谓的荧光透视模式，在介入手术期间广泛使用C形臂系统。在该模式操作下，C形臂可能用改变的投影方向，产生待成像对象的时序锥形束投影。产生的投影图像协助医生，例如，移动导管或放置支架。在介入手术中，仅投影相对小的3D感兴趣区域(ROI)常常就足够了。这可以用比大探测器更便宜的小探测器完成。而且，小探测器允许锥形束变窄，因而减少释放到患者的辐射剂量。

最近，已经建议使用C形臂用于真实3D成像或体积成像。通过首先采集对象的一系列锥形束投影实现体积成像，其中源沿着围绕对象的一些轨迹移动，并且然后使用执行重建算法的计算机根据这些投影重建对象。

重建图像包括体积元素或体素的3D阵列。其代表了对象3D区域中x-射线衰减系数空间分布的离散近似。精确已知体素相对于附于C形臂系统的参照系的坐标。

给定采集的几何配置，如由源轨迹和锥形束的轨迹所定义的，投影区域(ROP)定义为沿着源轨迹的所有锥形束中包含的体积。例如，如果源轨迹是完整的圆并且探测器的灵敏区域是圆盘，那么ROP是同心球，其直径(通常为15至25cm)部分取决于探测器的尺寸。如果探测器的灵敏区域不是圆盘或者源轨迹不是完整的同心圆，ROP的形状更加复杂，但是ROP仍然定位在同心上，并且其尺寸具有与前述例子相同的数量级。

尽管ROP内容的“精确”重建通常是不可能的(在医学应用中)，但是ROP中相当好的重建仍然可能。这是因为锥形束投影的缺失部分对ROP中精确结果仅仅有小的作用。而且，当相关积分线与ROP的距离增加时，缺失数据的影响迅速减少。对缺失部分的粗略猜测已经使得能在ROP中实现满意的重建，并且在重建过程中，可以安全地忽略沿着以大距离穿过ROP的线的线积分。尽管如此，必须明智地选择源轨迹。这些观察提示延伸截短的投影，使得它们表现为稍大于ROP的对象的非截短投影，并且提示将延伸的投影反馈到一种用于非截短投影的已知重建算法中。产生的图像与真实图像将区别仅在于，一些未知但微弱且平稳变化的幻像，其不混乱小的解剖细节。在ROP的边缘附近，幻象最强并且朝着该区域的内部迅速减少。投影截短得越少并且投影的延伸部分越精确，重建图像将越好。

延伸截短的投影的简单方法是对对象的形状和其中的x-射线衰减系数进行似乎合理的假设，并且根据这些假设对投影的测得部分进行外推。即使非常简单地假设对象是具有恒定x-射线衰减系数的球或者圆柱，其也运行意外得好。如果根本的假设是非常错误的，则结果可能不够满意。在医学应用中，当ROP处于人脑中时，出现这种情况。于是，用被骨头包围的软组织填充ROP，但是骨头对x-射线衰减贡献相当大。

上面概述的体积成像类型，是锥形束计算机断层摄影(CBCT)的一种形式。严格地说，并不强制使用C形臂进行数据采集。还可以使用其他类型的台架。源和探测器甚至可以连接到机器人臂并由机器人臂移动。然而，理想的是台架是“开放的”，从而医生可以自由地接近患者。C形臂台架符合这个要求。

基于C形臂系统的CBCT提供了高的、各向同性的空间分辨率，并且成为介入手术中标准荧光透视模式的受欢迎的助手。例如，在神经系放射学中，CBCT操作模式可以用于确认介入的成功或者如果在介入期间出现复杂情况时，辅助医生做出决定。患者可以保持在平台上并且不需要移动到另一个扫描器。释放给患者的剂量优选为尽可能低。由于这个原因，ROP应当尽可能小。此外，小ROP可以用比大探测器更便宜的小探测器实现。

患者的ROI中的x-射线衰减系数的3D分布还可以使用具有旋转台架的标准CT系统确定。用较老的CT系统时，患者平台并且患者与其一起，沿着CT系统的长轴逐步地平移。在每一步，根据当x-射线源围绕系统(患者)的长轴旋转时获取的大量扇形束投影，重建“当前”切片的2D图像。不出现截短的扇形束投影的问题(除非，可能有肥胖患者)。最后堆叠产生的2D图像以形成所需3D图像。较新的CT系统使用患者的连续平移从而源沿着围绕患者的螺线移动。无论如何，最终3D图像的体素以附于CT系统的坐标系表达。这样一种CT图像的分辨率是各项异性的并且沿着患者的长轴相对较差。常规CT台架是“封闭的”并且严格地限制医生接近患者。

在神经学应用中，出于诊断和计划的目的，通常在介入之前数小时或数天进行常规3D CT扫描。对于介入本身，代替地使用C形臂系统。如在上指出的，那么探测器应当相对小，并且当C形臂系统在CBCT模式运转时，小探测器可能导致有害的伪影。

在“Tomographic Reconstruction for Truncated Cone Beam DataUsing Prior CT Information”(R.E.Ellis和T.M.Peters eds.，MICCAI2003，LNCS2879，第134-141页，2003)中，K.Ramamurthi和J.Prince建议使用先前采集的3D CT图像来合成用C形臂系统采集的截短的锥形束投影的缺失部分。这实际上是沿着缺失的积分线数值积分3D CT图像中的x-射线衰减系数，如果C形臂系统具有大探测器，则其就可以采集到所述缺失的积分线。假设先前采集的CT图像在C形臂系统的坐标系中可利用。实际上，不存在这种情况。基于C形臂系统的CBCT图像并不是事前就知道的，从而不能用于与先前采集的CT图像配准。本文留下这个问题用于进一步研究以发展用于此目的的可靠的2D-3D配准算法，但是由于锥形束投影的复杂几何形状和这些投影不覆盖与先前CT图像相同信息的事实，这是困难的任务。

本发明的一个目的是提高CBCT图像的质量，所述CBCT图像依靠横向延伸相对小的探测器，诸如在C形臂系统中使用的那些探测器而获得。

因而，本发明提出一种如权利要求1所述的用于重建3D图像的设备。在权利要求2-7中引用该设备的优选特征。本发明的另一方面涉及权利要求8和从属权利要求9-13定义的计算机程序产品。

本发明利用中间CBCT图像，其通过使用诸如上面提到的简单延伸方法延伸可获得的截短的锥形束投影而获得。产生的图像不如最终图像好，但是足以用于配准目的。这样的配准可以在3D域中执行，优选基于该中间图像和先前CT图像之间的局部相关性。一旦已经完成配准，就使用先前的3D CT图像合成截短投影的缺失部分。

在该合成中，粗略的采样通常是足够的，其加速了手术。然后通过插值法上采样(up-sample)缺失数据，并且使用合成的和插值的数据完成用于第二(最终)重建的截短锥形束投影。产生的CBCT图像几乎没有截短伪影。

参考附图，在下列非限制性例子的描述中，本发明的其他特征和优点将变得清楚，其中：

图1是根据本发明的设备的结构框图；以及

图2是该设备的配准模块的结构框图。

图1示出了可用于根据截短的锥形束(CB)投影重建3D图像的图像重建设备，依靠具有相对小的探测器的C形臂x-射线装置获得所述截短的锥形束(CB)投影。因此，x-射线装置具有相对小的投影区域(ROP)，其例如在介入手术中覆盖感兴趣的区域(ROI)。这样小的探测器允许减少辐射剂量。其还比大探测器更便宜。

设备10具有两个输入，一个用于接收截短的CB投影，以及描述采集几何配置(c形臂外形，源轨迹，......)的数据，另一个用于接收先前计算的CT图像。在先CT图像可以是整个对象的图像，其在截短CB投影中部分可见。在这种情况下，在先CT图像例如得自在用大探测器执行的扫描中预先采集的锥形束投影。或者，在先CT图像可以仅仅覆盖对象的部分，只要该部分包括ROP同时显著大于ROP。例如，这种在先CT图像可以作为在此描述的重建设备10的输出获得。

在结构简图形式的图1中示出了设备10。在典型实施方式中，图1中所示的模块11-15是软件模块，即，包括写成执行下面详细描述的数据处理步骤的计算机程序代码。然后，设备10包括诸如工作站的计算平台，所述工作站与用于获得截短的CB投影的C形臂x-射线装置、存储在先CT图像的存储器、用于显示重建的CBCT图像的可视切片的显示器相联。还可以载入包括数据处理模块11-15的计算机程序以升级与C形臂x-射线装置相联的现有工作站。模块11-15还可以实施为硬件电路。

在重建程序的第一步骤中，外推模块11使用简单延伸方法外推探测器的视场之外的线积分。该延伸方法在此可以是粗略的，例如在US专利6,542,573中公开的简单的“放射状”延伸方法。通过将椭圆弧加到其两端来延伸每排的每个截短投影。延伸的横向范围由延伸因数控制，所述延伸因数是延伸投影的横向范围与截短投影的横向范围的比。在探测器之外，外推的线积分称作伪投影数据。

然后，由3D重建模块12将合适的重建算法应用到利用伪投影数据完成的截短CB投影，合适的重建算法例如是在1994年IEEE Trans.Med.Img.第13卷第186-195页，M.Defrise和R.Clark的“A cone-beam reconstruction algorithm using shift-variant filtering and cone-beam back-projection”中公开的算法。如果源轨迹是圆弧，对于模块12而言，甚至更简单的“切向”延伸方法就足够，例如在1984年J.opt.Soc.Amer.A第1卷第612-619页，L.A.Feldkamp、L.C.Davis和W.J.Kress的“Practical cone-beam algorithm”中公开的FDK算法。产生的中间CBCT图像在C形臂系统的坐标系中表达。

由模块13执行的接下来的步骤在于，使用合适的图像配准算法，将在先CT图像与在先前步骤中获得的ROP的中间CBCT图像配准。这样的算法能够处理不同的体素和图像尺寸，以及两幅图像之间适度的解剖学差异。该步骤产生了在C形臂系统的坐标系中表达的在先CT图像。

向前投影模块14使用由模块13输出的配准的在先CT图像，估算缺失的x-射线衰减系数的线积分，如果C形臂系统具有大探测器，则其已经测量了所述缺失的线积分。

有利地，模块14以粗略的空间分辨率执行X-射线衰减系数的实际积分，即，仅对以相同像素分辨率横向延伸的由真实C形臂探测器表示的实际探测器的一些像素。例如，可以每四个体素明确地计算积分，例如索引为偶数行和列的那些。如果投影数量大，仅对每第n个投影进行积分也是足够的，其中n是小的正整数，例如n＝2。粗略的采样节省了计算时间。然后模块14使用内插值法估算未明确估算的线积分。这比计算所有的缺失数据更快，并且仍然足够精确。

另外，模块14可以对向前投影数据应用线性变换，以便使这些数据在探测器边沿与截短CB数据匹配。例如，依靠在沿着边沿的交叠边缘中向前投影数据和截短CB数据之间的最小二乘法拟合，来计算增益和偏移量。这些参数用于补偿不同扫描中可能不同的束质量和散射贡献。

用模块14计算的向前投影数据完成测得的截短锥形束投影，以在第二3D重建模块15的输入提供完成的CB投影。该模块使用诸如在上面提及的M.Defrise等人的论文中公开的合适的算法，根据完成的锥形束投影计算ROP的最终3D图像。如果源轨迹是圆弧，可以使用FDK算法代替。产生的CBCT图像几乎没有截短伪影。

模块13应用严格的3D图像配准技术以将在先CT图像与中间CBCT图像配准。使用诸如交叉相关或共有信息的相似性手段的基于强度的配准方法，不容易应用，这是因为由于截短，参考图像和中间图像不共有相同的解剖信息，尤其对骨结构。优选使用这样的配准方法，其基于对中间CBCT图像中的许多小邻域进行的局部相关系数计算。在中间CBCT图像的骨的附近选择所述邻域，其首先通过应用阈值来确定。这样一种选择策略具有这样的优点，即仅在提供相关灰度值信息的区域内估算图像相似性。

该3D图像配准技术由图2中示出的配准模块13的结构简图示出。阈值应用模块20借助于确定的阈值分割中间CBCT图像，以产生可能代表图像中骨的体素。从分割的图像中，模块20可以进一步移除位于每个2D切片中的圆之外的所有体素，所述圆的直径定义为ROP在x和y方向上的最小延伸。模块20的输出是中间CBCT图像中的相关骨体素的3D坐标的列表。

模块13还包括常规边缘滤波器21，例如，Sobel或Canny滤波器，中间CBCT图像应用于其上。滤波器21向ROP的每个体素赋予与3D空间或2D切片中的x-射线衰减系数的最大梯度对应的边缘强度值。分类模块22将体素分类成按照它们的边缘强度值排序(首先最大值)的列表。该列表中的体素由选择模块23顺序处理，对于列表中的每个体素i，选择模块23对在某一局部邻域中的骨体素的数目k(i)计数(如模块20所指示)。这样的邻域是例如以正在考虑的列表的体素为中心的7×7×7立方体。

模块23建立了模板，该模板由在它们相应的局部邻域中计数得至少K个骨体素的列表中的体素构成。整数K可以取得相对小，例如K＝3，以便于消除最不可能代表骨结构的很小的吸收区域。该模板可以包括在列表顶部取得的并且满足k(i)≥K标准的固定数目N个体素。可选地，首先筛选该列表以消除不满足k(i)≥K标准的体素，并且模块23在该模板中保持已筛选列表的前Q％(例如，Q％＝3％)，因而在模板中产生可变数目N个体素。在列表顶部选择有限数目的模板体素，确保了这N个体素代表最显著灰度值转变，并且因此代表存在于中间CBCT图像中的最可能的骨的轮廓。

在后续步骤中，模块24使用由图2中的索引j指示的模板体素，作为基于局部相关相似性执行配准的锚定点。模块24确定在先CT图像的空间变换，其最大化在模板体素处中间CT图像和变化的在先CT图像之间的局部相关性度量M(T)。在优化程序中测试的每个变换例如由6个参数的参数集T确定，所述6个参数对应于平移向量(3个坐标)和3个旋转角度。

局部相关性度量可以具有这样形式：

$>>M>>(>T>)>>=>1>N>>>\Sigma >j>>>>>[>>\Sigma >>i>\in >n>>(>j>)>>>>>(>>b>i>>->ver>>b>‾>>j>>)>>>(>>t>i>>->ver>>t>‾>>j>>)>>]>2>>>>\Sigma >>i>\in >n>>(>j>)>>>>>>(>>b>i>>->ver>>b>‾>>j>>)>2>>>\Sigma >>i>\in >n>>(>j>)>>>>>>(>>t>i>>->ver>>t>‾>>j>>)>2>>>>>s>$

其中，在模板的N个体素j之上计算总和，n(j)指示体素j的局部邻域(例如，在两幅图像的体素分辨率下，以体素j为中心的5×5×5立方体)，bi是中间CT图像中的n(j)的体素i的x-射线衰减值，ti是通过T变换的在先CT图像中的n(j)的体素i的x-射线衰减值，和分别是n(j)内的bi和ti的平均值。

在2001年T.Netsch等人在International Conference on ComputerVision(ICCV’01)(Vancouver，BC，第501-508页)的“Towardsreal-time multi-modality 3-D medical image registration”中，公开了可由模块24用来最大化M(T)的有效算法，将参考所述文献。最佳参数集T构成了由模块24确定并且应用于在先CT图像的配准数据，以提供由向前投影模块14处理的配准的图像。

虽然已经参考其的优选实施方式描述了本发明，将认识到，不脱离本发明的精神和范围，可以对那些实施方式产生各种各样的修改方案。