专利名称：在存在移动时用于在3d空间中的x射线标志定位的方法
技术领域：
本发明的实施例涉及图像处理。具体而言，实施例提供用于对位于X射线图像中的标志关联点进行检测、定位及分类的方法和系统。
背景技术：
现代成像系统在单个过程或扫描期间获取许多图像。使用已知为重建的过程来组合图像，以制作表示穿过患者的有限宽度“切片”的一个或多个图像。因为在为医疗或牙科过程做准备时或者为了适当地校准成像设备，经常要对所得图像加以分析，所以重要的是图像尽可能地清楚。然而，扫描期间患者移动会导致图像中的误差或模糊。
在成像过程之前常将标志放置于患者上。这些标志在所得图像上显现为标志点，并用以识别图像中的特定特征。标志点也可以用以针对患者移动来修正图像。在典型的运动修正情形中，必须于在扫描期间捕捉的数百个(例如，600个)或更多图像中的每个中识别相对较小数量的标志(例如，7个)。因而,使用标志来检测及修正患者运动仅在能够在图像中自动识别数千个标志点(例如，使用以上的示例为4200个标志点)的时候才是可行的。标志点识别由于X射线的非线性特性以及X射线图像的相对较大的动态范围而进一步被复杂化。又一复杂之处在于，标志点识别和所得运动修正需要快速执行，这是因为其必须在图像的重建可以开始之前完成。现有的用于识别图像上的标志相似点的方法依赖于这些点的一般为圆形的特性，并采用诸如图像梯度处理或者霍夫(Hough)变换的一些变型的点识别技术，这些技术每图像要花费30秒之久才能完成。即使最简单的方法也要花费每图像1/10秒或更多，且经常要花费每图像5秒之久。这些运动修正方法的准确度与提取标志点位置的准确度直接相关。例如，在大多数修正方法中，标志必须定位至像素的1/4至1/8的精度。这个程度的准确度必须在存在宽的背景动态范围值、对图像背景的非线性依赖性、以及噪声的情况下实现。许多现有的修正方法不能正确处理图像背景的非线性特性，导致图像边缘处的位置误差。这些现有的方法当图像中的标志彼此相对较为靠近时也不能很好地工作。一旦标志点在图像中被识别出，这些标志点也就需要被映射回其所源自于的物理标志。物理标志通过它们的三维(“3D”)位置坐标而被唯一识别。因而，为了制作与特定标志相关联的每个投影图像的点的列表，一个步骤是识别每个标志的物理位置。这可以通过识别与由标志所产生的点的位置如何随着台架(gantry)转动而在投影图像之间的变化最为一致的每个标志的3D位置来进行。然而，因为患者移动或者计划外的系统或设备移动，所以每个标志的关联点位置会与其预期的位置明显改变。因此，识别标志的三维坐标的准确度在存在噪声或患者运动或者当多个标志位置彼此靠近时会受到损害。现有的用于识别标志的三维位置的技术可使用依赖于位置的电磁场，其由远程传感器检测并转换成位置。除该发射器和传感器系统的复杂度以外，系统所使用的技术也依赖于多个传感器，以确定标志的三维位置。其他定位技术包括光学技术，这类技术依赖于对多个被检测部位之间的固定关系的了解，从而增加了技术的复杂度。即使在与图像点相关联的物理标志被识别，并且位置被分配给每个之后，将每个图像上的每个点分配给其原标志时仍然存在问题。在此完成之前，所有图像中的所有点仅作为单个列表存在。点分配(标志跟踪)试图将每个图像上的每个标志与实际生成它的标志相关联。在完成点分配之后，每个标志将处理与该标志相关联的标志点的分尚列表。每个标志的分离列表可具有每个投影图像上的仅一个或零个标志点。重要的是，标志与其关联的标志点之间的正确对应关系即使在存在运动和测量误差的情况下也要被维持。
发明内容
实施例提供用于快速识别图像中的标志关联点的系统和方法，其然后可以用于运动修正系统中。一个方法包括处理每个投影图像，以减小除圆形物体(即，标志点的形状)之 外所有物体的强度，识别每个图像中的局部极大值，生长局部极大值以生成候选区域，并对每个候选区域应用圆度和投影连续性标准。该方法以与其他方法相比将总图像处理时间减小至少一个数量级的方式，使用一系列应用的简单均值滤波器、差分以及除法/减法操作，来减小除圆形物体以外所有物体的强度。一个方法包括获得图像、应用2D高通滤波器，以生成背景抑制后的图像、以及从所述背景图像中提取标志关联点。实施例还提供用于对图像中的标志点位置进行准确的子像素定位的系统和方法。一个方法包括在两个维度中每个上对图像中所识别的每个标志点的位置独立地进行评价。该评价通过对标志点的一侧或两侧采样，从而确定每个标志点的背景轮廓来进行。采用基线轮廓和标志点轮廓之比来修正背景影响，同时解决X射线衰减的非线性性质。残留的背景偏移通过减去标志点轮廓的拖影的平均值来去除。最后，对所得标志点轮廓进行评价，以确定其一阶矩(moment)，该一阶矩表示标志点的定位。替代地，可以通过使用位置拟合参数，将理论标志点轮廓拟合到该标志点轮廓来处理轮廓。该方法以明显比投影图像像素大小的分辨率高的分辨率，确定在计算机断层摄影扫描的投影图像中的标志点的位置。此外，该方法在存在明显不均质的背景和相邻结构时执行标志识别或定位。一个方法包括获得包含标志点的图像、得出背景修正后的标志点轮廓、得出背景和基线修正后的标志点轮廓、以及定义标志点的子像素位置以作为背景修正后的标志点轮廓的尖峰的中心。实施例也提供用于快速地识别并分配近似的3D位置给与图像中的标志点相关联的标志的系统和方法。一个方法包括获取观察到的投影图像之间的标志点的位置变化，并将位置变化拟合到通过使用V、Z和R柱坐标的不同的可能值所定义的理论投影轨迹。该步骤包括为柱坐标维度中的一个或多个来定义一组可能值，并从该组中选择这些值中的一个。该方法然后在每个投影图像中确定落在由从所述组中选择的值所定义的轨迹的距离范围内的标志点。然后通过找出每个标志点距轨迹范围的边界的均方距离来给落在该范围内的每个标志点分配度量。对于落在该范围内的所有标志点对度量求和，以得出与特定实验坐标值集合相关联的单个度量。对可能坐标值的集合中的每个值重复这些步骤，并评价度量的所得值(作为坐标值的函数)以确定函数尖峰，这些函数尖峰均假定为表示一个或多个标志的群集(cluster)。将与尖峰相关联的一个或多个柱坐标的这些坐标值分配给对应的得出的物理标志。对其余的柱坐标重复这些步骤，直到给标志分配了所有三个标志为止。因为度量拟合的性质，所以在存在噪音和患者运动时可以找出平均标志定位，并且提取并分开位置彼此靠近的多个标志。一个方法包括获得图像序列，每个图像均表示扫描器的转动角度，并具有第一值和第二值。该方法还包括贯穿该图像序列来分析标志点的第一值的性能，以确定标志的径向距离以及标志的角度位置。该方法还包括贯穿该图像序列分析标志点的第二值的性能，以确定标志的轴向位置。实施例也提供用于通过称为标志跟踪的处理而将每个投影图像中的每个标志点分配给最可能与该标志点相关联的单个物理标志的系统和方法。一个方法包括获得带有多个标志点的图像序列，并计算每个图像中的物理标志的预期U和V位置及范围。该方法还包括选择物理标志，并生成候选点列表，所述候选点列表由落在为该物理标志从理论上确定的每个投影图像中的预期标志点定位的范围内的标志点所构成。该方法进一步包括，通过将该列表中每个标志点的实际图像定位减去物理标志的理论轨迹，而将与理论的误差分配给每个标志点，从而来处理该列表。该方法然后逐个投影地遍历(traverse)该列表，并且基于U和V方向上的候选列表中的连续点的误差之间的一阶时间导数(derivative)以及其他标准，将该列表缩窄到每投影图像零或一个“优良点”。通过对具体说明书以及附图加以考虑，本发明的其他方面将会变得明显。
图I示意性示出锥形束计算机断层摄影系统。图2示意性示出图I的计算机。图3示出包括大量标志点的图像。图4是示出标志点识别方法的流程图。图5是示出图4的标志点识别方法的背景减小过程的流程图。图6是作为图5的背景减小过程的一部分的、用以在一个维度上对图像进行快速均值滤波的过程的图解说明。图7示出在应用了图5的背景减小过程之后的图3的图像。图8是示出图4的标志点识别方法的标志点提取过程的流程图。图9和图10示出在执行了图8的标志点提取过程的步骤之后的图7的图像。图11是示出标志点子像素定位方法的流程图。图12是示出图11的标志点子像素定位方法的标志点轮廓得出过程的流程图。图13是在图12的标志点轮廓得出过程中进行分析的围绕标志点的示例性像素列的图解说明。图14是示出在图12的标志点轮廓得出过程中所使用的标志点轮廓和基线轮廓的示例性值的曲线图。图15是示出图11的标志点子像素定位方法的标志点位置确定过程的流程图。图16是示出使用图15的积分(integral)标志中心确定过程的示例性标志点轮廓曲线的子像素标志点中心位置的图示。图17是示出使用图15的曲线拟合标志点中心确定过程,对标志点轮廓曲线的示例性最佳拟合高斯曲线的子像素标志点中心位置的图示。图18是示出物理标志三维定位方法的流程图。图19A-B是示出对标志点在U维度上的性能(behavior)进行分析的过程的流程图。图19C-D是示出对标志点在V维度上的性能进行分析的过程的流程图。图20是示出在U方向上存在运动时使Cp的试验值拟合到通过投影图像序列的大量标志点的多个U位置的过程的示例的曲线图。图21是示出R取固定值作为函数f的试验R和f值对的累积拟合度量的示例性一维绘图的曲线图。图22是示出使Z的试验值拟合到贯穿投影图像序列的三个标志的V轨迹的过程的示例的曲线图。图23是示出作为函数Z的累积拟合度量的示例性一维绘图的曲线图。 图24A-E是示出标志点到物理标志映射方法的流程图。图25是示出关于贯穿投影图像序列而识别的所有标志点的示例性物理标志和示例性U位置的预期U范围的曲线图。图26是示出具有在图25中所示U范围内的U位置的图25中所示标志点的曲线图。图27是示例性候选标志点列表的图解说明。图28是示例性候选标志点列表的变平过程的图解说明。图29是示例性候选标志点列表的同步化过程的图解说明。图30是当“下一个”投影图像包含两个或更多候选标志点时处理“下一个”投影图像的图解说明。图31是当“下一个”投影图像包含单个候选标志点时处理“下一个”投影图像的图解说明。图32是当“下一个”投影图像不包含候选标志点时处理“下一个”投影图像的图解说明。图33是当优良点列表包含优良点和可疑点时处理“下一个”投影图像的图解说明。
具体实施例方式在详细解释本发明的任何实施例之前，要理解的是，本发明在其应用时不限于在下面的描述中所阐述的或者在附图中示出的构造的细节以及部件的布置。本发明能够具有其他实施例，并且能够以各种方式实践或实施。还应予以注意的是，可利用多个基于硬件和软件的装置、以及多个不同结构的部件来实现本发明。此外，如随后的篇章中所描述的，图中所示特定配置旨在例示本发明的实施例。替代配置也是可能的。X射线是一种形式的射线，其类似于光线，但可以穿透人体的射线。然而，当X射线穿过人体时，X射线的强度减弱。强度的减弱量与X射线穿过的组织的密度有关。X射线的这一特性被用以产生人体以及其他物体内部的结构的图像。X射线系统传统上包括产生X射线的发散束的X射线管。人或物体置于X射线管与一张膜之间，这张膜对在每个定位处射在其上的X射线的强度敏感。开关X射线管，对该膜进行显影，并且所得图像呈现X射线路径中的结构的“阴影(shadow)图像”，且其中，结构越密集则在图像中显得越白。现代X射线过程使用一个或多个检测器，而非使用膜。检测器以电子方式检测并量化到达检测器的X射线的量(即，密度)。使用检测器，建立X射线计算机断层摄影(“CT”)系统，其使X射线源绕着患者转动，并且利用在患者的与X射线源相反的一侧的单个宽条状(single-wide strip)检测器元件，以电子方式检测所得X射线。针对所有不同的X射线位置收集来自检测器的数据，然后在已知为重建的过程中对其进行组合。组合后的图像表示穿过患者的单个有限宽度“切片”，其中每个点处的图像密度表示特定物理定位处的组织的X射线密度。X射线源、检测器、以及允许X射线源转动的机械结构的组合已知为“CT台想”
o通过在图像采集或扫描之间移动患者、X射线源、或者两者的同时重复上述过程，可以获取多个切片。例如，移动支承患者的平台且同时也转动X射线源，代替数据切片而产生“螺线(helix)”。另外，将检测器条或环的大小或宽度从单行检测器增加到多行(例如， 达到256行)，允许对更多数据进行更加快速的采集。此外，用较大的二维检测器替换检测器条，获取每个X射线源位置处的整个检测器面板图像，而非仅单条数据。对在数量上达到600或者更多的这些的图像的收集已知为投影图像。每个投影图像表示随着X射线源和检测器绕着患者同步地转动从不同的透视或角度对患者的X射线快照。因为锥形X射线束需要覆盖二维检测器，所以这种CT成像已知为锥形束(“CB”)CT成像。图I示意性示出了根据本发明的一个实施例的CB CT成像系统10。CB CT成像系统10包括扫描器12和计算机14。扫描器12包括CT台架13，其包括X射线源16、检测器18、和转动载台20。X射线源16和检测器18在转动载台20上在彼此对面对准，该转动载台20使X射线源16和检测器18绕着患者22移动。患者22由座位24所支承。图I中所示的成像系统10为牙科成像系统。因此，患者22坐在座位24中并将他或她的下颚放置到托架(rest) 26中。当台架13转动以完成对患者的头部的扫描时，托架26保持患者的头部相对静止。扫描器12将图像数据输出至计算机14。图像数据表示在扫描期间由检测器18所检测到的X射线的强度水平。计算机14连接至控制台30，其包括显示器32和一个或多个输入和/或输出装置，例如键盘34。使用者使用控制台30以与计算机14交互。例如，用户可以使用控制台30以从计算机14请求图像或其他数据。计算机14将所请求的信息提供至控制台30，并且控制台30将该信息显示在显示器32上、将该信息打印至打印机(未示出)，和/或将该信息保存至计算机可读存储器模块(未示出)。图2示意性示出了根据本发明的一个实施例的图I的计算机14。计算机14包括输入/输出接口 40、电子处理单元(“EPU”)42、以及一个或多个存储器模块，诸如随机访问存储器(“RAM”)模块44和只读存储器(“ROM”)模块46。输入/输出接口 40从扫描器12接收图像数据，并将该图像数据提供给EPU 42。在一些实施例中，输入/输出接口 40将图像数据储存至存储器模块，诸如RAM模块44，并且EPU 42从存储器模块获得图像数据用于处理。输入/输出接口 40也可将数据传输至扫描器12，诸如校准和操作数据或指令。EPU 42接收图像数据，并通过执行一个或多个应用或模块来处理信息。在一些实施例中，应用或模块储存在存储器中，诸如ROM模块46。在其他实施例中，应用或模块可储存在硬盘存储单元上并被传送至RAM用于执行。如图2所示，ROM模块46储存标志点识别模块50、标志点子像素定位模块52、物理标志三维(“3D”)定位模块54、标志点映射模块56、患者运动识别模块58、和图像重建模块59。EPU42检索并执行标志点识别模块50，以识别所接收的图像数据中的标志点(参见图4-10)。EPU 42检索并执行标志点子像素定位模块52，以提供对标志点位置的准确的子像素定位(参见图11-17)。EPU 42检索并执行物理标志3D定位模块54，以确定通过物理标志在投影图像中的标志点识别的所述物理标志的三维位置(参见图18-23)。EPU 42检索并执行标志点映射模块56，以将在投影图像中识别的标志点映射到特定的物理标志(参见图24-33)。EPU 42检索并执行患者运动识别模块58，以得出患者运动信息。EPU 42检索并执行图像重建模块59，以将投影图像与运动信息结合并生成一组运动修正轴图像。完成台架旋转需要大约8到40秒。在此期间，患者可能移动，导致所得的图像的模糊。典型的图像分辨率为0.25毫米的量级。因此，这个量级的患者移动通常导致图像模糊，并且大范围的患者移动可以使得所得图像对于其所期望的临床目的而言变得不可接 受。计算机14 (参见图12)可以通过跟踪所得图像中刚性物体的移动来针对患者移动对图像进行修正。例如，在其中没有患者移动的理想条件下，被成像的刚性物体随着台架绕着患者转动而以明确定义的方式改变在投影图像的二维上的位置。患者移动导致图像中物体的期望位置之间的偏离。通过测量明确定义的物体与其期望位置的偏离，可以测量并修正患者运动量。具体而言，如果图像中有至少三个物体，则可以组合这些物体与其期望位置的所测量的偏离来确定六参数患者位置误差值(例如，患者运动矢量)，其可以应用于图像以修正患者移动。为了确保在图像中存在所期望的数量的明确定义的刚性物体，可以将一个或多个基准标志置于患者上。这些标志一般包括铅或钢BB，其较为致密并可防止或限制X射线透过。这些标志也可以由其他材料制成，并且构造成其他的形式或形状，其在扫描期间生成的投影图像中以高比例可见。每个标志或多个标志可位于多层粘合剂之间，并且该粘合剂可以涂覆于患者，以确保标志在手术期间不移动。将标志置于患者上，使得由CB CT成像系统所生成的每个视野或图像均包括至少三个标志点。例如，在患者上可放置七到九个标志，以确保至少三个标志点在由CB CT成像系统所生成的每个图像内，以减少位置测量噪声，并允许结果的统计组合。在一些实施例中，标志均匀地分隔，以相对于患者具有最大的空间分布，并避免对于正常的图像判读(interpretation)的干扰。识别图像中的标志在将标志置于患者上之后，对患者进行扫描，诸如，对患者的头部进行扫描，并将表示投影图像序列的所得图像数据传输至计算机14 (参见图2)。图3示出了包括八个标志点102的示例性投影图像100。作为针对患者移动而修正所生成的投影图像时的初始步骤，计算机14处理所生成的投影图像，以识别每个图像中的标志点。图4示出了根据本发明的一个实施例的标志点识别方法60。当EPU 42执行标志点识别模块50时，由计算机14的EPU 42来执行标志点识别方法60。如图4所示，方法60的第一步通常包括减小每个投影图像的背景的强度(步骤62)。此背景抑制步骤使图像的背景部分最小化。图像的背景部分包括除了标志点之外的图像中的所有事物。因为基准标志具有明确定义的形状(即，球形或圆形)，所以步骤62使图像中未具有此明确定义的形状的所有事物最小化。在一些实施例中，通过对每个图像应用线性高通滤波器来实现背景抑制。为了获得可接受的速度，执行三个优化。首先，通过将原始图像与低通滤波后的图像组合来实现高通操作。其次，代替在两个图像方向上同时执行滤波器，而在两个或更多个维度上(即，以0°、45°、90°和135° )对图像依次应用一维高通操作。最后，利用通过下述的累积像素方法所实施的简单均值滤波器来实现低通滤波器操作。一维高通滤波器操作将趋于加强在滤波器的方向上快速变化的区域。在不同方向上执行的一系列这类操作将趋于相对于所有其他形状而加强圆形物体。由于标志所产生的标志点为圆形，所以这一系列滤波将加强标志点，且减弱任何周围的结构。 图5示出了背景减小步骤62中所涉及的过程或步骤。如图5所示，通过首先对原始图像A进行简单的一维的、累积均值低通图像滤波(步骤69)，以生成滤波后的图像B来执行多维高通滤波操作。然后，如下所述，将滤波后的图像B与原始图像A组合以生成图像C(步骤70)。使用来自于先前的滤波步骤的输出图像作为对下一个滤波步骤的输入，对每个滤波方向重复这两个步骤(步骤71-73)。图5所示过程允许快速的、多维高通滤波，以减小图像的背景。图6示出了一个方向上的一维低通滤波中所涉及的这些步骤的示例。代替使用更典型的基于内核的滤波方法，通过使用累积均值滤波器(图5A中所绘)能更加快速地实现步骤69和71的滤波。此快速滤波通过生成新图像(图像BI)来实现(步骤74)，其中，新图像中的各个像素均被设定为在滤波器的方向上从该特定像素向下直至第一像素的、来自原始图像A的所有像素值的累积和。例如，如图6所示，当计算机14 一开始从扫描器12获得原始图像(图像A)时，该方法考虑包括七个像素的图像A的一行(其中，处理方向为水平的)，然后通过将该行中的各个像素设定为该像素左侧的像素值们加上该像素的值之和来生成累积图像(图像BI)。一经在步骤74生成累积图像，图像BI就移位与低通滤波因子相等的量(例如，2)，以生成移位后的图像(图像Cl)(步骤75)。通过对移位后的边缘像素简单地调零来处理图像边缘。例如，如图6所示，累积图像(图像BI)中的每个像素向左移位两个位置，以生成移位后的图像(图像Cl)。因为前两个像素(即，边缘像素)的左侧两个位置没有像素，所以这些像素设定为O。然后将累积图像(图像BI)减去移位后的图像(图像Cl)以生成差值图像(图像Dl)(步骤76)。如图6所示，为了减去图像，将未移位的累积图像(图像BI)中的各个像素值减去移位后的图像(图像Cl)中的对应像素值。最后的步骤是通过将差值图像中的每个像素除以滤波器因子(移位值)以归一化差值图像(图像D1)，从而产生低通滤波后的图像B (步骤 77)。步骤74-77的计算时间不依赖于低通滤波器的内核大小。这些步骤中的计算粗略涉及2XN个加法和减法，而不是(内核大小)XN个操作，其一般需要普通的基于内核的滤波器。对于512X512个像素的图像和典型的平滑因子15，这表示计算时间相对于传统的I-D内核方法实现了 7. 5:1的减少。与2-D内核相比，(相对于两个ID累积均值滤波器操作)对于滤波器因子15，速度提高为56:1。在生成一维低通滤波后图像之后，通过将低通图像(图像B)与用于生成其的图像(图像A)组合而生成高通滤波后图像(图像C)(步骤73)。这可通过两个图像之间的简单差值操作来进行。然而，除法操作用于第一方向滤波步骤70中的处理。使用除法以补偿X射线衰减的指数性质和归一化图像(图像B)，以生成图像C (步骤70)。对差值图像执行除法是有效的，这是因为X射线衰减是以材料密度的负指数的形式出现的。因此，由于除法等同于对数相减然后求指数，所以对差值图像执行除法避免了计算开销高的取对数和求指数操作。在顺序的维度中的每个上重复应用一维的、累积均值低通滤波器，并组合输入和输出图像的过程(步骤71-73)。然而，在处理第一维度之后，在步骤72可以使用减法代替除法，这是因为对X射线衰减的指数性质的修正通过在处理第一维度时执行除法(步骤70)来完成。一般来讲，0、45、90和135度处的四个通道(pass)足以将球状标志点从宽广范围的背景中适当地隔离。图7示出了在应用图5的背景减小方法之后的图3的图像100。如图7所示，背景的大部分(即，图像中除标志点102外的部分)已被减小或暗化，允许标志点 102更容易识别。回到图4，一旦背景通过上述步骤减小，则从高通滤波后的图像中提取标志点(步骤80)。图8示出了从高通图像中提取标志点的过程或步骤。如图8所示，该过程开始于在每个高通图像中定义局部极大值(步骤82)。此步骤包括识别每个图像中的每个其所有8个相邻像素都具有比其小的值的像素。接着，使用所识别的局部极大值作为种子点以生长区域(步骤84)。此过程将特定区域内所包含的局部极大值组合为单个候选物体或区域。在一些实施例中，所产生的区域的大小或范围受根据经验指定的下阈值的限制。图9示出了在区域生成局部极大值(步骤82和84)之后的图7的图像100，其生成候选区域104。在步骤86，对每个候选区域应用形状标准，以识别哪个区域最有可能代表标志点。在一些实施例中，此过程通过取在(I)候选区域(a)的像素面积与(2)通过候选区域的平均半径的平方乘以n ((a的平均半径)2X Ji)限定的面积之间的比率来执行。该比值用于确定候选区域是否具有合适的形状以作为标志点(例如，候选区域是否具有圆形形状)。最后，在步骤88，通过比较相邻投影图像中的候选区域的定位并应用接近度(nearness)标准来进一步处理候选区域，以丢弃不可能的候选者。例如，如果候选区域在顺序的投影图像之间其位置移动太多，则该候选区域不可能代表标志点。类似地，如果两个候选区域太靠近在一起，则这两个区域不可能代表标志点。图10示出了在对候选区域104应用形状标准并丢弃不可能的候选区域104 (步骤86和88)之后的图9的图像100。如图10所示，在这些步骤之后，仅留下与八个标志点102 (参见图7和8)相对应的候选区域104。在一些实施例中，对图像应用二维高通滤波器包括首先对该图像应用二维低通滤波器以生成低通图像。之后，将原始图像与低通图像组合。在一些实施例中，标志点识别方法60可以在时间上以每图像少于50毫秒的量级来执行。另外，当使用球形标志时，标志点识别方法60表现出对标志的高选择率，以及对非标志相关联结构的高丢弃率。例如，在一些实施例中，标志点识别方法60可以具有100%的灵敏度，并且在实际情形中具有比40:1更好的特异性(低于2. 5%的错误检测)。确定标志的子像素中心点即使在使用标志点识别方法60在图像中识别标志点之后，在这些标志点可以用于执行运动修正之前也需要对其进一步地或者更准确地限定。具体地，必需以比由标志点识别方法60所提供的像素分辨率(子像素)好得多的分辨率来知晓标志点的有效二维位置。可能需要小到像素分辨率的1/8的分辨率。再者，这必须在图像内以及图像之间都存在极大变化的背景中完成。图11示出了根据本发明一个实施例的标志点子像素定位方法110。当计算机14的EPU 42执行标志点子像素定位模块52时，由EPU 42执行标志点子像素定位方法110。在一些实施例中，标志点子像素定位方法110如下所述使用由标志点识别方法60所识别的标志点作为起始点,来确定准确的子像素标志定位。如图11所示,标志点子像素定位方法110通常包括标志点轮廓(profile)得出过程(步骤112)和标志点中心确定过程(步骤114)，对投影图像中的两个维度中的每个重复 执行这两个步骤(步骤116)。如以下更详细描述的(参见图12-14)，标志点轮廓得出过程112应用一系列步骤，以去除背景中的大部分，并在每个期望的维度上生成多像素标志点轮廓。在标志点中心确定过程114中，将所生成的标志点轮廓与标志的预期形状相比较，并得出一维子像素标志点中心定位(参见图15-16)。图12示出了根据本发明一个实施例的标志点轮廓得出过程112。考虑得出标志点的垂直位置的操作，其中，垂直方向上的像素序列称为列，且水平方向上的像素序列称为行。用于得出标志点的水平位置的过程是相似的，但所有对行的提及都由列替换，而对列的提及都由行替换。像素的行和像素的列均可称为像素线(Pixel line)。该过程的第一步涉及通过在围绕标志的背景列的范围内对像素值求平均值来得出垂直背景轮廓。图13示出了围绕标志点122并以其为中心的像素的示例性列121。步骤120对围绕标志点122的背景列121a的像素值求平均值，并排除以标志点122为中心的标志列121b的像素值。对背景列121a求平均值包括对跨列121a的每一行的像素值求和并将和除以列数。对背景列121a求平均值生成单列背景值，其为标志点122的背景轮廓。接着，通过对标志列121b求平均值来得出标志点轮廓，其中，所述标志列121b是以标志点122为中心的列(步骤124)。再者，对标志列121b求平均值包括对跨列121b的每行的像素值求和并将和除以列数。对标志列121b求平均值生成单列标志值，其为标志点122的标志点轮廓。在步骤126，以逐行为基础，将背景轮廓除以标志点轮廓。例如，将组成背景轮廓的平均值的列中的每一行中的值除以组成标志点轮廓的平均值的列中的对应的行的值。按此顺序对行值执行除法确保出现图像反转(例如，暗区变亮而亮区变暗)。另外，步骤126使用除法(例如，代替减法)以补偿X射线衰减的负指数性质(这是因为，如所指出的，除法在数学上等同于对两个轮廓的对数之差求指数)。作为替代实施例，像素值的对数可以在找出标志点轮廓和背景轮廓之前计算。然后，通过以逐个元素为基础将两个轮廓相减并对每个元素求指数来实现对标志点轮廓的背景补偿。此方法较之于以上方法略微精确一些，这是因为其能更准确地处理在求平均值的行上的图像变化，但其比简单的求平均值和除法的方法更耗时间。此时，已获得标志点的合理的基线修正后的轮廓。然而，该轮廓的拖影(tail)与零相比仍会存在一些残留的偏移，如果拖影是非对称的则尤其如此。该轮廓通过下列程序而达到零。平均基线偏移值通过对离开标志点的尖峰(例如，具有最高值的标志点内的像素)的拖影中的列上的像素值求平均值来确定，例如在其中像素值不再随着其从标志点中心移动离开而减小的列上对像素值求平均值来确定(步骤128)。关注于从此平均像素中排除与和此标志点相邻的任何结构，例如图像中的其他标志或其他结构相关联的拖影。图14用曲线图示出了示例性标志点轮廓129和示例性基线轮廓130的值，以及对标志点轮廓133进行完全基线修正之后的结果。回到图11，在生成基线修正后的标志点轮廓之后，定位方法110的第二步包括标志点位置确定过程114。图15更详细地示出了此过程。如图15所示，存在至少两个替代项或选项以用于执行此过程。第一选项包括确定标志点轮廓的积分(步骤140)，以及将标志点位置定义为表示该积分的一半的沿着积分的位置(步骤142)。此标志点位置有效地为一阶矩的积分，或者是在其中积分在该点左侧的部分等于积分在该点右侧的部分的位置。图16示出了示例性标志点轮廓曲线146的示例性子像素标志点位置定位144。回到图15，用于确定标志点位置的第二选项使用最小二乘法或其他拟合将适当的曲线函数拟合到标志点轮廓(步骤146)，并且将标志点位置定义为表示曲线函数的中心的参数(步骤148)。图17示出了拟合到示例性标志点轮廓曲线151的最佳拟合高斯曲线150 的示例性子像素标志点位置定位149。步骤146中所使用的曲线函数可以是标志点的高斯曲线或预期曲线，该曲线将与成像物理学有关的问题，例如部分体积效应(partial volume effect)和指数衰减特性纳入考虑。通常，步骤146所使用的曲线函数可以是其中曲线定义参数之一表示曲线中心的任何曲线。曲线函数也可以处理非常量的基线轮廓。在一些实施例中，如果选择处理非常量基线轮廓的曲线函数，则标志点轮廓得出步骤112中的步骤128和132可以去除。此外，可以基于可测量的图像特性，例如已知的标志维度、像素大小、以及其他维度上的标志的近似子像素定位，来针对每个标志、图像、或者图像组而动态地修改曲线函数。确定标志的3D梓坐标即使在图像中识别标志点和标志点中心(例如,使用标志点识别方法60和标志点子像素定位方法110)之后，标志点的积累也需要分类成使得每个标志点均与负责生成该点的物理标志适当地相关联。为了做到这一点，第一步是唯一地识别每个物理标志点。物理标志点均通过其三维(3D)坐标而唯一识别。因此，将标志点与物理标志相关联的第一步是识别物理标志并且给每个分配其3D坐标。然而，此过程由于患者移动以及在投影图像中确定标志点位置的误差而会复杂化。另外，由于存在多个标志，所以在一些投影图像中，一些标志点会重叠在其他标志点之上。图18示出了用于识别物理标志中的每个并在3D空间中给每个物理标志分配位置或定位的物理标志3D定位方法160。当计算机14的EPU 42执行物理标志3D定位模块54时，由EPU 42来执行物理标志3D定位方法160。物理标志3D定位方法160通常贯穿(through)随着CT台架转动而生成的图像序列来跟踪由标志产生的点的位置，以将标志映射到3D空间中的位置。现有的跟踪技术在多个被检部位之间使用多个传感器和多个固定关系来跟踪患者在3D空间中的位置。然而，方法160中所应用的跟踪技术或方法使用隔离的单点检测部位。因此，代替比较同时收集的多个部位的特性，该跟踪方法收集单个点随时间的特性，并将观察到的随时间的位置变化匹配到3D位置，其将产生相同的随时间的位置变化。此跟踪方法除CT成像台架外无需额外的设备，并且容许相对较大量的图像噪音和患者移动。如图18所示，方法160开始于得到每个标志点在每个投影图像上在标为“U”和“V”的两个图像维度上的位置(例如，通过标志子像素定位方法110而确定的每个标志点的中心)(步骤162)。U是与CT台架的转动圆形相切的维度，且V是与CT台架的转动轴线相平行的维度。标志点的这些维度位置可以通过标志点识别方法60、标志点子像素定位方法110、或其组合来确定。如果标志在3D空间中具有固定的位置，则其在每个投影图像上的U和V位置(或
值)都可由下列公式表示
权利要求
1.一种确定在投影图像序列中所识别的标志的径向距离(R)、角位置Ccp)和轴向位置(Z)的方法，所述方法由包括扫描器、带有电子处理单元的计算机、和存储有可由所述电子处理单元执行的标志三维定位模块的存储器模块的成像系统执行，所述方法包括 在所述计算机处，基于由所述扫描器所产生的图像数据获得图像序列，所述图像序列中的每个图像表示所述扫描器转动的角度，且包含标志点位置，所述标志点位置具有在与所述扫描器的转动的圆形相切的第一维度上的第一值以及在与所述扫描器的转动轴线相平行的第二维度上的第二值； 利用所述电子处理单元，贯穿所述图像序列，分析所述标志点位置的所述第一值的性能，以确定所述标志的所述径向距离和所述标志的所述角位置；以及 利用所述电子处理单元，贯穿所述图像序列，分析所述标志点位置的所述第二值的性能，以确定所述标志的所述轴向位置。
2.根据权利要求I所述的方法，其中，分析所述标志点位置的所述第一值的性能包括 Ca)选择径向距离和角位置值对的子采样第一集合， (b)从所述子采样第一集合中选择径向距离和角位置值对， (C)基于所选择的径向距离和角位置值对，确定用于所述图像序列中每个图像的值的第一范围， (d)如果所述标志点位置的所述第一值在用于各个图像的值的所述第一范围内，则为所述图像序列中相应图像中的所述标志点位置确定拟合度量， Ce)对所述度量拟合求和，以生成累积拟合度量，并将所述累积拟合度量分配给所选择的径向距离和角位置值对， (f)对所述子采样第一集合中的每个径向距离和角位置值对，重复步骤(b)到(e)， (g)对分配给所述子采样第一集合中的每个径向距离和角位置值对的所述累积拟合度量进行评价，作为所述子采样第一集合中的径向距离和角位置值对的函数， (h)检测所评价的函数的尖峰，所述尖峰对应于尖峰径向距离值和尖峰角位置值，以及 (i)将所述尖峰径向距离值分配给所述标志的所述径向距离，并将所述尖峰角位置值分配给所述标志的所述角位置。
3.根据权利要求2所述的方法，其中，步骤(c)包括基于所选择的径向距离和角位置以及由图像表示的转动角度，确定用于每个图像的理论第一值，并以所述理论第一值为所述值的第一范围的中心。
4.根据权利要求2所述的方法，其中，步骤(d)包括确定所述标志点位置的所述第一值与所述第一范围中的最小范围值和所述第一范围中的最大范围值中至少一者之差，并对所述差求平方。
5.根据权利要求I所述的方法，其中，分析所述标志点位置的所述第一值的性能包括 (a)选择假定的径向距离值， (b)选择角位置值的子采样第一集合， (c)从所述子采样第一集合中选择角位置值， (d)基于所假定的径向距离值和所选择的角位置值，确定所述图像序列中每个图像的值的第一范围， (e)如果所述标志点位置的所述第一值在用于各个图像的所述值的所述第一范围内，则为所述图像序列中相应图像中的所述标志点位置确定拟合度量， (f)对所述度量拟合求和以生成累积拟合度量，并将所述累积拟合度量分配给所选择的角位置值， (g)对所述子采样第一集合中的每个角位置值，重复步骤(C)到(f)， (h)对分配给所述子采样第一集合中的每个角位置值的所述累积拟合度量进行评价，作为所述子采样第一集合中的角位置值的函数， (i)检测所评价的函数的尖峰，所述尖峰对应于尖峰角位置值，以及 (j)将所假定的径向距离值分配给所述标志的所述径向距离，并将所述尖峰角位置值分配给所述标志的所述角位置。
6.根据权利要求5所述的方法，其中，步骤(c)包括基于所假定的径向距离值、所选择的角位置以及由图像表示的转动角度，确定用于每个图像的理论第一值，并以所述理论第一值为所述值的第一范围的中心。
7.根据权利要求5所述的方法，其中，步骤(d)包括确定所述标志点位置的所述第一值与所述第一范围中的最小范围值和所述第一范围中的最大范围值中至少一者之差，并对所述差求平方。
8.根据权利要求I所述的方法，其中，贯穿所述图像序列，分析所述标志点位置的所述第二值的性能，以确定所述标志的所述轴向位置包括 Ca)选择轴向位置值的子采样第二集合， (b)从所述子采样第二集合中选择轴向位置值， (C)基于所述标志的径向距离、所述标志的角位置以及所选择的轴向位置值，确定用于所述图像序列中每个图像的值的第二范围， (d)如果所述标志点位置的所述第二值在用于各个图像的值的所述第二范围内，则为在所述图像序列中的相应图像中的所述标志点位置确定拟合度量， (e)对所述度量拟合求和以生成累积拟合度量，并将所述累积拟合度量分配给所选择的轴向位置值， (e)对所述子采样第二集合中的每个轴向位置值，重复步骤(b)到(e)， (f)对分配给所述子采样第二集合中的每个轴向位置值的累积拟合度量进行评价，作为所述子采样第二集合中的轴向位置值的函数， (g)检测所评价的函数的尖峰，所述尖峰对应于尖峰轴向位置值，以及 (h)将所述尖峰轴向位置值分配给所述标志的所述轴向位置。
9.根据权利要求8所述的方法，其中，步骤(c)包括基于所述标志的径向距离、所述标志的角位置、所选择的轴向位置值以及由所述图像表示的转动角度，确定用于每个图像的理论第二值，并以第二第一值为所述值的第二范围的中心。
10.根据权利要求8所述的方法，其中，步骤(d)包括确定所述标志点位置的所述第二值与所述第二范围中的最小范围值和所述第二范围中的最大范围值中至少一者之差，并对所述差求平方。
11.一种成像系统，所述成像系统包括扫描器、带有电子处理单元的计算机、和存储有可由所述电子处理单元执行的用于确定在投影图像序列中所识别的标志的径向距离(R)、角位置(cp)和轴向位置(Z)的标志三维定位模块的存储器模块，其中，标志点子像素定位模块当执行时被配置为 在所述计算机处，基于由所述扫描器所产生的图像数据获得图像序列，所述图像序列中的每个图像表示所述扫描器转动的角度，且包含标志点位置，所述标志点位置具有在与所述扫描器的转动的圆形相切的第一维度上的第一值，以及在与所述扫描器的转动轴线相平行的第二维度上的第二值； 利用所述电子处理单元，贯穿所述图像序列，分析所述标志点位置的所述第一值的性能，以确定所述标志的所述径向距离和所述标志的所述角位置；以及 利用所述电子处理单元，贯穿所述图像序列，分析所述标志点位置的所述第二值的性能，以确定所述标志的所述轴向位置。
12.根据权利要求11所述的成像系统，其中，为了分析所述标志点位置的所述第一值的性能，所述标志点子像素定位模块当执行时被配置为用于 (a)选择径向距离和角位置值对的子采样第一集合， (b)从所述子采样第一集合中选择径向距离和角位置值对， (C)基于所选择的径向距离和角位置值对，确定用于所述图像序列中每个图像的值的第一范围， (d)如果所述标志点位置的所述第一值在用于各个图像的所述值的第一范围内，则为所述图像序列中相应图像中的所述标志点位置确定拟合度量， Ce)对所述度量拟合求和以生成累积拟合度量，并将所述累积拟合度量分配给所选择的径向距离和角位置值对， (f)对所述子采样第一集合中的每个径向距离和角位置值对，重复步骤(b)到(e)， (g)对分配给所述子采样第一集合中的每个径向距离和角位置值对的所述累积拟合度量进行评价，作为所述子采样第一集合中的径向距离和角位置值对的函数， (h)检测所评价的函数的尖峰，所述尖峰对应于尖峰径向距离值和尖峰角位置值，以及 (i)将所述尖峰径向距离值分配给所述标志的所述径向距离，并将所述尖峰角位置值分配给所述标志的所述角位置。
13.根据权利要求12所述的成像系统，其中，步骤(c)包括基于所选择的径向距离和角位置以及由所述图像表示的转动角度，确定用于每个图像的理论第一值，并以所述理论第一值为所述值的第一范围的中心。
14.根据权利要求12所述的成像系统，其中，步骤(d)包括确定所述标志点位置的所述第一值与所述第一范围中的最小范围值和所述第一范围中的最大范围值中至少一者之差，并对所述差求平方。
15.根据权利要求11所述的成像系统，其中，为了分析所述标志点位置的所述第一值的性能，所述标志点子像素定位模块当执行时被配置为用于 (a)选择假定的径向距离值， (b)选择角位置值的子采样第一集合， (c)从所述子采样第一集合中选择角位置值， (d)基于所假定的径向距离值和所选择的角位置值，确定用于所述图像序列中每个图像的值的第一范围， (e)如果所述标志点位置的所述第一值在用于各个图像的所述值的第一范围内，则为所述图像序列中相应图像中的所述标志点位置确定拟合度量， (f)对所述度量拟合求和以生成累积拟合度量，并将所述累积拟合度量分配给所选择的角位置值， (g)对所述子采样第一集合中的每个角位置值，重复步骤(C)到(f)， (h)对分配给所述子采样第一集合中的每个角位置值的所述累积拟合度量进行评价，作为所述子采样第一集合中的角位置值的函数， (i)检测所评价的函数的尖峰，所述尖峰对应于尖峰角位置值，以及 (j)将所假定的径向距离值分配给所述标志的所述径向距离，并将所述尖峰角位置值分配给所述标志的所述角位置。
16.根据权利要求15所述的成像系统，其中，步骤(c)包括基于所假定的径向距离值、所选择的角位置以及由所述图像表示的转动角度，确定用于每个图像的理论第一值，并以所述理论第一值为所述值的第一范围的中心。
17.根据权利要求15所述的成像系统，其中，步骤(d)包括确定所述标志点位置的所述第一值与所述第一范围中的最小范围值和所述第一范围中的最大范围值中至少一者之差，并对所述差求平方。
18.根据权利要求11所述的成像系统，其中，为了贯穿所述图像序列分析所述标志点位置的所述第二值的性能以确定所述标志的所述轴向位置，所述标志点子像素定位模块当执行时被进一步配置为用于 Ca)选择轴向位置值的子采样第二集合， (b)从所述子采样第二集合中选择轴向位置值， (C)基于所述标志的径向距离、所述标志的角位置以及所选择的轴向位置值，确定用于所述图像序列中每个图像的值的第二范围， (d)如果所述标志点位置的所述第二值在用于各个图像的值的所述第二范围内，则为所述图像序列中的相应图像中的所述标志点位置确定拟合度量， Ce)对所述度量拟合求和以生成累积拟合度量，并将所述累积拟合度量分配给所选择的轴向位置值， (e)对所述子采样第二集合中的每个轴向位置值，重复步骤(b)到(e)， (f)对分配给所述子采样第二集合中的每个轴向位置值的累积拟合度量进行评价，作为所述子采样第二集合中的轴向位置值的函数， (g)检测所评价的函数的尖峰，所述尖峰对应于尖峰轴向位置值，以及 (h)将所述尖峰轴向位置值分配给所述标志的所述轴向位置。
19.根据权利要求18所述的成像系统，其中，步骤(C)包括基于所述标志的径向距离、所述标志的角位置、所选择的轴向位置值以及由所述图像表示的转动角度，确定用于每个图像的理论第二值，并以第二第一值为所述值的第二范围的中心。
20.根据权利要求18所述的方法，其中，步骤(d)包括确定所述标志点位置的所述第二值与所述第二范围中的最小范围值和所述第二范围中的最大范围值中至少一者之差，并对所述差求平方。
全文摘要
一种确定在投影图像序列中所识别的标志的径向距离(R)、角位置和轴向位置(Z)的方法和系统。本发明的实施例允许可由电子处理单元执行的标志三维定位模块基于由扫描器所生成的图像数据而获得图像序列。该投影图像序列中的每个图像均表示扫描器转动的角度，并且包含标志位置。贯穿所述图像序列，对标志点位置的第一值和第二值的性能进行分析，以确定标志的径向距离、标志的角位置、和标志的轴向位置。因此，本发明的实施例允许在投影图像中对放置于患者上的外部标志快速地进行检测和定位。
文档编号A61B6/00GK102753097SQ201080062135
公开日2012年10月24日 申请日期2010年8月26日 优先权日2009年11月25日
发明者大卫·谢伯克 申请人:登塔尔图像科技公司