CN110148160A - 一种正交x射线影像快速2d-3d医学图像配准方法 - Google Patents

Abstract

本发明公开了一种正交X射线影像快速2D‑3D医学图像配准方法，不断通过在两对正交角度下的DRR‑X射线影像间的2D‑2D配准来引导3D空间刚体变换。在3D空间刚体变换中对射线与DRR像素进行坐标反变换避免3D并且将计算量大的步骤DRR生成等步骤转移到GPU上进行并行计算，克服2D‑3D配准时间长的问题，达到快速进行2D‑3D医学图像配准的目的。

Description

一种正交X射线影像快速2D-3D医学图像配准方法

技术领域

本发明涉及一种正交X射线影像快速2D-3D医学图像配准方法， 属于医学图像数据处理技术领域，可应用于图像引导肿瘤放疗系统中 的患者姿态修正。

背景技术

图像引导系统主要用于辅助引导肿瘤放射治疗，通过放疗过程中 实时X射线影像与术前重构CT图像进行配准获取患者摆位误差,并修 正该误差使得放射线能够准确轰击到患者肿瘤靶灶，减小对患者正常 组织的损伤。因此在临床实际应用中为了缩短患者治疗过程中等待时 间并保证放疗效果，需要进行快速准确的2D-3D配准。

现有的2D-3D配准算法以基于灰度的迭代优化式配准算法为主， 其主要思想是由CT数据模拟X射线投影生成DRR(Digitally Reconstructured Radiograph，数字重建放射影像)，通过衡量DRR 与X射线影像的相似程度来来判断患者当前位置姿态与术前重构CT姿态是否存在差异，并且进一步通过优化算法寻找使得DRR与X射线 影像最相似的CT空间位置。在CPU的串行程序处理下需要不断遍历 CT体素、DRR以及X光像素等极其耗时的过程。由于DRR投影以及相 似性度量的计算本质上就是复杂度分别为O(n3)与O(n2)，这意味着一般的2D-3D配准算法是计算量巨大的算法，难以满足临床治疗中快 速配准的要求，需要对此不足进行改进。

发明内容

针对现有技术的不足，本发明提出了一种基于并行计算的正交X 射线影像快速2D-3D医学图像配准算法，其提供了一种采取GPU进行 并行计算，并且充分利用一对正交面板上2D坐标系统反映的空间位 置信息不断引导3D刚体变换进而实现快速配准的方法，输入一对射 线夹角为90°的两块探测面板获取的X射线影像与CT数据可输出术 前重构CT与当前X射线影像相对应的6个空间刚体变换参数。

本发明采取的技术方案为：一种正交X射线影像快速2D-3D医学 图像配准方法，充分利用在两个垂直方向进行拍摄的一对正交X射线 影像包含的空间位置信息，并通过在GPU上的并行计算实现快速精确 的2D-3D配准,主要包括以下步骤：

1)获取术前进行CT重构获得CT数据。并获取两幅在射线夹角 为90°下拍摄的X射线影像，保存两次X射线影像及其拍摄系统的 参数。系统参数主要包括射线源到等中心点距离，射线源到探测面板 距离，拍摄角度，X射线影像像素点的物理尺寸。

2)根据步骤1)中获取的拍摄参数，建立虚拟坐标系，设置射 线源点、等中心点、探测面板在虚拟坐标系中的位置，为虚拟探测面 板上每个像素点赋值一个虚拟坐标系坐标，并且将CT数据根据其体 素物理距离为每一个体素中心赋值一个虚拟坐标系的坐标，模拟X光 拍摄机械系统。

3)在步骤2)中所建立的虚拟坐标系中根据Siddon-Jacobs FastRay-Tracing算法进行DRR生成，但将每个DRR像素点的像素值 计算在GPU上并行化独立运行。

4)将步骤3)中生成DRR图像后与X射线影像在GPU上并行化 计算二值化局部加权均差相似度，以3D空间初始位置参数根据几何 结构正交分解得到的2D刚体参数为初始位置，对生成DRR图像进行 并行化的2D-2D配准。

5)将步骤4)中两对2D-2D图像间的配准结果，再次根据正交 分解关系反演到3D空间中，作为3D变换的初始位置，跟据这一初始 位置搜索优化得到新的刚体变换参数，对虚拟坐标系中的探测面板、 光源点进行反变换，并生成新的DRR图像，再次计算出归一化互相关 相似度。

6)重复步骤3)到步骤4)直到两对2D-2D图像间的归一化互相 关相似度不再增加，输出当前4)中最后的空间变换参数以及当前参 数对应的DRR图像。

其中局部加权均差相似度包含以下计算过程：

1)根据图像均值作为阈值将整幅图像二值化,对于大小为N*M的 2D图像Img，其第i行j列像素Img(i,j)二值化结果Img′(i,j)：

2)根据图像二值化结果Img′，计算第i行j列处的局加权均值 M(i,j)：

其中参数m,n为局部矩形区域的长宽，alpha,beta计算如下：

alpha＝(2n+1)(2m+1)-1,

根据两幅图像加权均值M₁(i,j)、M₂(i,j)，计算相似度 Metric(M₁,M₂)：

3)对于图像I₁,I₂归一化互相关相似度Ncc(I₁,I₂)包含以下计算过 程：

Cvar(I₁,I₂)计算方式如下：

对于图像I，Avg(I)计算方式如下：

对于图像I，Var(I)计算方式如下：

本发明改进之处如下：

本发明与已有技术相比，采用正交角度下拍摄一对X射线影像作 为输入图像，通过分析正交拍摄系统中空间刚体变换参数与探测面板 上平面刚体变换之间的几何关系，不断采用2D-2D配准来引导2D-3D 配准过程中的参数优化，减少优化过程中的迭代次数。并且通过将时 间复杂度高的DRR生成、2D双线性插值、3D空间坐标变换以及两种 相似度计算进行基于GPU的并行化计算，缩短每一次迭代时间。通过 本发明，可在肿瘤放射临床治疗中，不损失精度的情况下更快速的进 行2D-3D配准，修正患者摆位误差。

附图说明

图1-a是本发明中设立坐标系间几何结构及其几何参数关系分析 图；

图2-a是本发明简要流程图；

图3-a是本发明实际案例中输入的一对头颅模体的人工合成正 交X射线影像；

图3-b是本发明实际案例中输入的头颅模体CT视图；

图4-a是本发明实际案例中相对于图3-a的输出结果图像；

图4-b是本发明中输入图像与输出配准结果图像的棋盘显示效 果。

具体实施方式

为了体现本发明的技术特点，具体实施一次2D-3D配准案例并结 合附图说明本发明具体流程。下面以一次对患者进行术中X射线影像 拍摄完毕开始至2D-3D配准结束完成姿态位置偏差的输出为实施案 例，X光图像以及CT数据采用公开的头颅CT数据及其模拟X光图像， 结合附图和技术方案进一步说明本发明的具体实施方式。具体步骤如 下：

1)采用一个如附图3-b所示头颅模体CT作为CT输入其物理大 小为400*400*284，体素物理尺寸为2*2*2mm³，分别在0°和90° 角度，无旋转，无平移的条件下人工合成模拟一对正交X射线影像图 像作为配准目标参数的对应的2D图像，如附图3-a，即配准目标位 置为Tx,Ty,Tz,Rx,Ry,Rz都是0。保存合成X射线影像采取X射线影 像像素物理尺寸，射线源到探测面板距离，射线源到等中心点距离， CT体素物理尺寸，预设配准中使用512*512的配准图像，设定配准 的初始位置为沿着X,Y,Z三个方向平移-5mm,-5mm,-5mm,旋转 5°,5°,5°。

2)根据步骤1)中获取的机械参数，建立虚拟坐标系如附图1-a。 设置射线源点、等中心点、探测面板在虚拟坐标系中的位置，为虚拟 探测面板上每个像素点赋值一个虚拟坐标系坐标，并且将CT数据根 据其体素物理距离为每一个体素中心赋值一个虚拟坐标系的坐标，模 拟X光拍摄机械系统。

3)在步骤2)中所建立的虚拟坐标系中,计算所有射线源到DRR 像素中心点物理坐标间的光路。将每一条光路上累积经过的CT体素 值进行加权累积，权值由光路经过的体素内部长度占虚拟光路在CT 内部总长度之比计算得到。将DRR像素点的像素值计算分别放在在 GPU的一个计算单元上并行化独立运行，限于本实验中硬件限制，所 有的DRR像素值分块进行并行计算后再汇总成DRR图像。在此步骤中 生成与0°、90°下X射线影像相对应的两幅DRR图像。

4)将步骤3)中生成两幅DRR图像分别与其对应的X射线影像 以3D空间初始位置参数根据附图1-a几何结构关系带入1)中得到 的几何参数数据，计算出2D刚体参数为初始位置。在GPU上计算对 生成DRR图像计算局部加权均差相似度，其中delta取值为0。

5)从步骤4)中初始2D位置及其对应局部加权均差相似度出发， 由powell搜索优化方法搜索下一步2D变换参数。根据搜索得到的新 的参数，在GPU上对DRR图像进行双线性插值的2D刚题变换，得到 新的DRR图像。

6)重复步骤4)到步骤5)这一过程直到新的DRR图像与X射线 影像间的相似度在相邻两次变换中的差值小于阈值0.001。对0°、 90°下的DRR-X射线影像对，分别独立进行上述2D-2D配准。

5)将步骤4)中2对DRR-X射线影像像间的配准结果，根据附 图1-a中正交分解关系反演到3D空间中，作为3D变换的初始位置。 由这一初始位置生成新的一对正交DRR图像，分别与对应的X射线影 像在GPU上计算NCC相似度，并且取两对DRR-X射线影像NCC度量均 值作为初始迭代的相似度。以该初始相似度，结合由2D-2D配准结果 根据附图1-a所示几何关系公式反演的3D变换初始位置，对6个维 度分别进行一维搜索，计算出新的3D空间变换参数。

6)由该新的3D空间变换参数，对DRR像素物理坐标及其对应的 射线源坐标进行反变换。该变换结果等价于对CT体素进行正变换。 得到新的光源点与DRR像素物理坐标点，进入2)生成新的一对正交 DRR图像。

7)重复步骤2)到步骤6)的过程不断进行2D-2D配准，引导 3D空间变换，直到两次循环间的步骤5)计算得到NCC相似度差值 小于0.001，停止迭代。输出步骤5)中计算的NCC相似度值0.996358， 参数Tx,Ty,Tz,Rx,Ry,Rz分别为[-0.022,-0.081,-0.106,0.0,-0.215,-0.165]， 耗时6.62秒,同时输出一对正交DRR图像如附图4-a。为直观展现配 准结果，将输入、输出的两对2D-2D图像进行棋盘状拼接，如附图 4-b。

表1-a中显示了在输入相同数据的情况下经典医学图像处理库 ITK与本发明的运行结果对比，其中Tx,Ty,Tz单位为毫米，Rx,Ry,Rz 单位为度，采用角度制。相对ITK库，在非常短的时间内取得更高的 配准精度。

表1-a

本发明采用正交角度下拍摄一对X射线影像作为输入图像，通过 分析正交拍摄系统中空间刚体变换参数与探测面板上平面刚体变换 之间的几何关系，不断采用2D-2D配准来引导2D-3D配准过程中的参 数优化，减少优化过程中的迭代次数。并且通过将时间复杂度高的 DRR生成、2D双线性插值、3D空间坐标变换以及两种相似度计算进 行基于GPU的并行化计算，缩短每一次迭代时间。通过本发明，可在 肿瘤放射临床治疗中，进行快速高精度的2D-3D配准，修正患者摆位 误差。

以上所述只是本发明的优选实施方式，对于本技术领域的普通技 术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和 润饰，这些改进和润饰也被视为本发明的保护范围。

Claims (4)

1.一种正交X射线影像快速2D-3D医学图像配准方法，其特征在于，该方法包括下述步骤：

步骤一：重复利用2D-2D近似刚体配准结果由几何关系反演出3D空间参数；

步骤二：在2D-2D配准过程中采用一种局部加权均值误差来衡量2D图像间的相似程度。

2.根据权利要求1所述一种正交X射线影像快速2D-3D医学图像配准方法，其特征在于，步骤一中2D-2D近似刚体配准的具体配准过程如下：

1)先根据3D初始参数，利用正交分解将6个刚体参数中的5个分解到2个正交面板坐标系中；

2)以步骤1)中分解的参数作为初始参数，在两个正交面板坐标系中分别进行2D-2D配准；

3)将步骤2)中的2D-2D配准结果反演至3D刚体变换参数的5个，计算NCC相似度并以此为起点优化搜索下一组3D参数；

4)重复步骤1)到步骤3)的过程，直到3)中相邻两次迭代所计算的NCC相似度的增加量小于预设阈值。

3.根据权利要求2所述的一种正交X射线影像快速2D-3D医学图像配准方法，其特征在于，步骤2)中包含的2D-2D配准过程中，对于并非完全一致的两幅2D图像，鲁棒的进行相似度计算；具体计算步骤如下：

1)根据图像均值作为阈值将整幅图像二值化，对于大小为N*M的2D图像Img，其第i行j列像素Img(i,j)二值化结果Img′(i,j)：

2)根据图像二值化结果Img′，计算第i行j列处的局加权均值M(i,j)：

其中参数alpha，beta计算方式如下：

alpha＝(2n+1)(2m+1)-1,

3)根据两幅图像加权均值M₁(i,j)、M₂(i,j)，计算相似度Metric(M₁,M₂)：

4.根据权利要求2所述的一种正交X射线影像快速2D-3D医学图像配准方法，其特征在于，步骤1)中3D刚体参数搜索完成后不是对CT体素进行3D刚体变换与插值，而是对射线源坐标、DRR像素坐标在GPU上进行反变换，反变化具体步骤如下：

1)3D空间参数取相反数；

2)对两个虚拟射线源点及其对应DRR像素按照输入的X射线影像拍摄过程的机械参数在GPU上并行化赋值物理坐标；

3)对步骤2)中的物理坐标根据步骤1)中得到的变换参数在GPU上对虚拟射线源点以及对应DRR像素物理坐标并行化计算空间变换。