CN101194845A - 基于双x射线机的同步实时动态跟踪立体定位方法及装置 - Google Patents

Abstract

本发明是基于双X射线机的同步实时动态跟踪立体定位方法及装置，方法步骤：采用摄象机对X射线机进行标定；同时采集两路X射线的图象；采用边缘检测的方法对图象进行实时处理，自动从图像中提取穿刺针和肺部小病灶信息，根据标定参数计算出穿刺针和病灶的三维位置坐标，从而完成定位。装置是在水平向C臂的二端分设水平向C臂X线球管和水平向影像增强器；在垂直向C臂的二端分设垂直向C臂X线球管和垂直向影像增强器。优点：在现有的小病灶穿刺过程中，对于穿刺针和病灶的位置关系的判断，由于只采用一个X射线机，只有一路信息，加之肺小病灶会随肺的呼吸而移动，在判断穿刺针是否穿中病灶时，需要靠经验，采用本发明，就比较容易判断。

Description

基于双X射线机的同步实时动态跟踪立体定位方法及装置

技术领域

本发明涉及的是一种基于双X射线机的同步实时动态跟踪立体定位方法及装置，属于医疗仪器技术领域。

背景技术

目前对于人体内部结构的医学三维定位，所采用的方式是，通过相关的成像设备，在不同时刻从不同角度获取图像再进行静态重建三维定位。该方式，对于静态目标，并且不需要实时定位的情形，是可以满足要求的，例如CT设备就是采用该方式进行定位的。但是对于运动目标，特别是肺部随呼吸上下移动的小病灶，以及要求进行实时定位的场合，上述方式就无法胜任。例如，在进行肺部穿刺取检的过程中，需要实时确定穿刺针和病灶之间的关系，即对这二者进行三维实时动态跟踪立体定位，以便随时调整进针的方向和深度。目前经皮肺穿刺所采用的定位方法是，即在穿刺前，采用金属标记物事先进行定位，以确定穿刺针进针的位置。在穿刺过程中，如使用CT设备进行定位，首先，CT设备无法做到同步实时动态跟踪定位穿刺活检，二则分辨率不高，通常对于1cm左右的肺部小病灶很难分辨，三是成本比较高。最常用的方法是依靠一台X射线机来观察穿刺针与病灶之间的位置关系，对于肺部小病灶，判断是否穿刺成功是就比较困难的。因为，肺的呼吸移动，肺部小病灶会随肺的呼吸而移动，当穿刺针位于病灶的正上方或正下方时，在X射线的图象上，所反映出来的效果都是两者和病灶重叠在一起，完全没有穿入肺部病灶。这时很大程度上需要靠医生的经验来进行判断。因此，这是目前经皮肺穿小病灶取检中的一个难点。

发明内容

本发明的目的在于提供一种能够为肺部小病灶的定位穿刺提供帮助，帮助医生能够快速准确地完成肺部病灶的动态跟踪三维立体定位穿刺取检工作的基于双X射线机的同步实时动态跟踪立体定位系统。

本发明的技术解决方案：基于双X射线机的同步实时动态跟踪立体定位装置，其结构是水平向C臂1的一端设水平向C臂X线球管4，另一端设水平向影像增强器6；垂直向C臂2的一端设垂直向C臂X线球管3，另一端设垂直向影像增强器5。

基于双X射线机的同步实时动态跟踪立体定位方法，其步骤分为，

(1)采用摄象机对X射线机进行标定；

(2)同时采集两路X射线的图象；

(3)采用边缘检测的方法对图象进行实时处理，自动从图像中提取穿刺针和肺部小病灶信息，根据标定参数计算出穿刺针和病灶的三维位置坐标，从而完成定位。

所述的步骤(1)采用摄象机对X射线机进行标定的具体方法：

1)确定图像中心，在X射线机的接收端和发射端之间，置一网格状的模板，使其和接收端平行，成像一次，然后，沿着垂直于接收端影像增强屏的方向，移动物体，改变模板与成像平面的距离，再次成像，计算成像中心；

2)dx和dy的确定，在求取dx和dy的比值以及乘积的时，利用已知的四个点和四条边所能提供的组合，算出所有的dx*dy以及dx/dy，再分别加以平均，作为最终的dx*dy以及dx/dy的值；

3)外参数和畸变系数的标定，采用共面点的方法进行X射线机的标定，利用图像点和实际点的对应关系，可以计算出旋转矩阵中的r1，r2，r4，r5以及tx和ty，然后再利用正交性，计算出矩阵R中其余的所有元素，采用最速下降法，同时对焦距、tz以及畸变系数进行优化得到了所需的参数。

所述的步骤(2)同时采集两路X射线图象的具体方法：采用一个可以两路同时输入的视频采集卡，把两个X线机的视频输出信号直接接到该卡上，通过该采集卡，把模拟信号转为数字信号。

所述的步骤(3)采用边缘检测的方法对图象进行实时处理的具体方法是根据病灶的位置，自动选择能够包含穿刺针和病灶的矩形区域，在该区域内，采用边缘检测算子进行边缘检测，然后采用Hough变换进行直线检测，检测出穿刺针，在此基础上，根据灰度信息检测出针头的位置，对于病灶，根据上一帧的位置，在其周围进行搜索，利用灰度信息和形态信息。

本发明的优点：本发明采用两台X线机从两个不同角度，同时采集图象，然后自动分析图象中的目标和病灶，然后根据事先标定好的成像参数，计算出穿刺针和和病灶的三维坐标，从而完成同步实时动态跟踪的三维立体定位工作。通过事先对两台X射线机进行标定，确定X线成像的相关参数，在穿刺过程中，两台X线机的图象同时输入计算机，通过计算机的实时自动处理，自动从图象中提取穿刺针和病灶，再根据事先标定好的参数值，计算两者的三维坐标之间的关系，从而完成三维实时动态跟踪定位。在现有的小病灶穿刺过程中，对于穿刺针和病灶的位置关系的判断，由于只采用一个X射线机，只有一路信息，加之肺小病灶会随肺的呼吸而移动，因此，在判断穿刺针是否穿中病灶的时候，很多时候需要靠经验，采用本发明，就比较容易判断。

附图说明

图1是本发明的立体定位系统，

图2是本发明的水平向X线投射与检查床中轴关系示意图，

图3是本发明的线板和缩光器等准直器开口设计。

图4是原始模板的形状，

图5是模板所成的像。

图中的1是水平向C臂，2是垂直向C臂，3是垂直向C臂X线球管，5是垂直向影像增强器，4是水平向C臂X线球管，6是水平向影像增强器，7是穿刺检查床，8是X线图像处理器及计算机辅助系统。

具体实施方式

对照附图1，其结构是水平向C臂1的一端设水平向C臂X线球管4，另一端设水平向影像增强器6；垂直向C臂2的一端设垂直向C臂X线球管3，另一端设垂直向影像增强器5。

所述的水平向影像增强器6和垂直向影像增强器5属现有技术。

实施例1，

水平向C臂的跨度为130cm，拱深为140cm。由于该C臂的X线球管位于医生操作的一侧，选择α夹角为30。～45。，这样可观察到穿刺针长度并给医生穿刺操作留出较为适宜的空间。

垂直向C臂，主要用来完成肺部病灶垂直方向X线成像(其C臂位于水平向C臂内，并且由床的头侧伸进，C臂的跨度为110cm，拱深为150cm)，观察穿刺针的进针深度及左右偏移程度，与水平向C臂共同完成引导穿刺针进入病灶内任务。

根据几何关系，可计算进针角度和穿刺针可观测到的长度。

对照图2，图中的矩形为检查床，CD为检查床宽，S为X线发射端，AB为接收屏直径，F为接收屏中心点，SF为水平向C臂X线球管到影像增强器的距离，O为肺部病灶点，也是C臂旋转中心点并位于检查床的中轴线上，EO为穿刺进针方向，OP为水平向X线视野内穿刺针进入病灶的可视长度，ΔASB为经准直后的水平向X射线成像区，α为水平向C臂X射线中心线与检查床中轴线的旋转夹角。

对照图3，采用限线板和缩光器等准直器开口设计来控制X线的视野和散射，使双C臂双球管X射线同步清晰成像。准直器最大开口(I)按公式I＝L×d/D计算。

X线机的标定

首先从计算机视觉的角度出发，借鉴摄像机标定技术，对两台X射线机进行标定和图像矫正。按照X线机成像的几何模型，求解出在该成像过程中的内参数、外参数等有关参数，在此基础上利用径向畸变模型，对X射线机的成像畸变进行矫正处理。

成像的几何过程：

1)从世界坐标系的坐标(X_w，Y_w，Z_w)到X射线机坐标系的坐标(X_c，Y_c，Z_c)的变换

世界坐标系的坐标(X_w，Y_w，Z_w)到X射线机坐标系的坐标(X_c，Y_c，Z_c)的变换

$\left[\begin{array}{c}{X}_c\\ Y_c\\ Z_c\end{array}\right]=R\left[\begin{array}{c}{X}_w\\ Y_w\\ Z_w\end{array}\right]+T$

2)从X射线机坐标系的坐标(X_c，Y_c，Z_c)到理想图像坐标系的坐标(X_u，Y_u)的变换

$X_u=\frac{f{X}_c}{Z_c}$

$Y_u=\frac{f{Y}_c}{Z_c}$

3)从理想图像坐标系的坐标(X_u，Y_u)到实际的图像坐标(X_d，Y_d)的变换，即径向畸变

X_d＝X_u(1+k₁r²)

Y_d＝Y_u(1+k₂r²)

r²＝X_u ²+Y_u ²

4)从实际成像的图像坐标(X_d，Y_d)到计算机帧存中的像素坐标(X_f，Y_f)变换

X_f＝X_d/d_x+cx

为了标定外参数，成像中心(cx，cy)，X方向和Y方向单位像素之间的间隔dx和dy等内参数必须是已知的。但是X射线机的这些内参数，是未知的。因此，必须首先标定这些内参数。再根据图像点和实际点的对应关系，利用现有的摄像机标定算法，进行外参数的标定，计算出旋转矩阵R和平移向量T及畸变系数k。为了实现图像点和物体点的对应，需设计特定形状的模板，对模板成像，然后提取模板上的特征点，完成图像点和实际点的对应。本方法采用的是网格状的模板。

确定图像中心

在X射线机的接收端和发射端之间，置一网格状的模板，使其和接收端平行，成像一次。然后，可以沿着垂直于接收端影像增强屏的方向，移动物体，改变模板与成像平面的距离，再次成像，这样就可以计算成像中心了。假设模板上有两个点，Pi和Pj，Pi1(xi1，yi1)，Pj1(xj1，yj1)，Pi2(xi2，yi2)，Pj2(xj2，yj2)分别表示这两个点在第一次和第二次所成像的位置，那么成像中心可以表示成下面的形式：

$\left\{\begin{array}{c}\mathrm{cx}=\frac{2}{n\·(n-1)}\underset{i=1}{\overset{n-1}{\mathrm{\Σ}}}\underset{j=i+1}{\overset{n}{\mathrm{\Σ}}}\frac{x_{i1}\·x_{j2}-x_{j1}\·x_{i2}}{x_{i1}-x_{i2}+x_{j2}-x_{j1}}\\ \mathrm{cy}=\frac{2}{n\·(n-1)}\underset{i=1}{\overset{n-1}{\mathrm{\Σ}}}\underset{j=i+1}{\overset{n}{\mathrm{\Σ}}}\frac{y_{i1}\·y_{i2}-y_{j1}\·y_{i2}}{y_{i1}-y_{i2}+y_{j2}-y_{j1}}\end{array}\right.$

单位像素之间的间隔dx和dy的确定

dx和dy表示的是单位像素之间的间隔，这对于计算最终的像素位置是很重要的。确定的方法为：用模板成像，把模板放置于X射线机的接收端，更准确地说是置于影像增强器的输入屏上，这样所成的像，除了受到畸变参数的影响外，不再受其他任何参数的影响。如果没有畸变的话，按照成像几何，物体位于该位置所成像的大小应该等于实际物体的大小。

由于实际所成的像存在有一定的畸变，而且是越远离中心越严重，因此在选择图像点时，采用中心附近的点。通常采用半径为图像边长五分之一到六分之一的区域，在这样的区域里，畸变通常是可以忽略的。我们采用构成一个正方形的四个顶点，该正方形的边长已知，面积也已知。为了求解dx和dy，至少需要两个方程。如果单纯利用边长作为已知条件来建立方程，得到的方程往往是病态的，它的解不可靠。可利用正方形边长相等的属性来建立一个方程，再利用其它属性得到另外一个方程。设任意相邻两边的三个顶点坐标分别为(x1，y1)，(x2，y2)，(x3，y3)，则所构造的方程为：

((x1-x2)*dx)2+((y1-y2)*dy)2＝((x1-x3)*dx)2+((y1-y3)*dy)2化简上式后，得到的是dx和dy的比例关系。

另外一个方程可以通过下列方式得到。设构成该区域的四个点为p1，p2，p3和p4，每个点的坐标是已知的。假设平面上任意三个点的坐标分别为(x1，y1)，(x2，y2)，(x3，y3)，三点围成的面积可以表示为：

$\left|\begin{array}{ccc}x1& y1& 1\\ x2& y2& 1\\ x3& y3& 1\end{array}\right|$

假设p1，p2，p3三点的像素坐标分别为(x1，y1)，(x2，y2)，(x3，y3)，那么三点围成的实际面积的表达式为：

$\left|\begin{array}{ccc}x1*\mathrm{dx}& y1*\mathrm{dy}& 1\\ x2*\mathrm{dx}& y2*\mathrm{dy}& 1\\ x3*\mathrm{dx}& y3*\mathrm{dy}& 1\end{array}\right|$

把上述行列式展开后得到的结果如下：

(x1*y2-x2*y1+x2*y3-y2*x3+x3*y1-y1*x3)*dx*dy。

由于三点形成的区域的实际面积是正方形面积的一半，从而得到了另外一个方程，这个方程所提供的是dx和dy的乘积。有了这两个方程，就可以解得dx和dy了。为了进一步减小误差，可以采用平均的方法。即在求取dx和dy的比值以及乘积的时候，利用已知的四个点和四条边所能提供的组合，算出所有的dx*dy以及dx/dy，再分别加以平均，作为最终的dx*dy以及dx/dy的值。

得到了内参数(cx，cy)和(dx，dy)后，就可以进行X射线机外参数以及畸变系数的标定。

外参数和畸变系数的标定

外参数和畸变系数的标定通常采用的方法是：先假设畸变系数为0，在此前提下，计算外参数，然后再利用有关信息，进行优化计算，得到畸变系数。

内参数标定之后，X射线机的标定问题就可以采用一般的摄像机标定方法来解决了。本文基本上采用Tsai的标定方法，来计算后面的有关参数，但是由于X射线机的特点，有个别步骤与之有所差异。

在Tsai的方法中，主要是提出了共面点和不共面点的两种标定方法。为了简化问题，我们采用共面点的方法进行X射线机的标定。利用图像点和实际点的对应关系，可以计算出旋转矩阵中的r1，r2，r4，r5以及tx和ty，然后再利用正交性，计算出矩阵R中其余的所有元素。不过，与Tsai的方法不同的是，因X射线机的畸变程度比较严重，估算出来的f和tz误差较大。所以不对焦距f进行估计，而是把它作为已知值。尽管该值可能会存在一定的误差，但是要远小于同时对f和tz进行估算所得到的f的误差。通过上述步骤，得到了除畸变系数以外的所有参数。摄像机如果不是采用广角镜头，且对于精度的要求不是很高的话，则完全可以忽略畸变系数；但是对于X射线机，则不能忽略，因为它的畸变是很明显的。求取畸变系数是一个非线性问题，无法通过求解一般的线性方程得到，通常采用优化策略加以解决。可采用最速下降法，同时对焦距、tz以及畸变系数进行优化。这样，就得到了所需的参数。

在上述系统配置的基础上，把两台X射线机的图象同时输入计算机，应用计算机图像处理和计算机视觉技术图像的矫正和实时动态跟踪算法设计。实时动态跟踪法的设计：首先通过边缘检测等技术，从图像中得到所需的目标点，包括肺部病灶和穿刺针的信息，再采用HOUGH变换得到穿刺针所在的位置，根据穿刺针的灰度信息以及图像序列信息，到穿刺针头在图像中的位置；其次，穿刺前利用链码跟踪技术，计算出肺部病灶区域的形态和灰度等特征，在穿刺过程中根据前面计算所得的特征信息及位置，根据前一幅图中病灶的位置，实时动态地跟踪病灶，并在前一幅图中病灶位置的周围区域进行搜索，从而可以获取两幅图象中病灶和针的信息，再利用标定的结果，得到病灶和针的三维信息，从而可以确定穿刺针、病灶之间的三维空间位置关系，达到计算机辅助动态跟踪导引穿刺的目的。

原始模板的形状如图4所示，模板成像后的图像如图5所示。

Claims (5)

1.基于双X射线机的同步实时动态跟踪立体定位方法，其特征是该方法的步骤分为，

(1)采用摄象机对X射线机进行标定；

(2)同时采集两路X射线的图象；

(3)采用边缘检测的方法对图象进行实时处理，自动从图像中提取穿刺针和肺部小病灶信息，根据标定参数计算出穿刺针和病灶的三维位置坐标，从而完成定位。

2.根据权利要求2所述的基于双X线机的实时立体定位方法，其特征是步骤(1)采用摄象机对X射线机进行标定的具体方法：

1)确定图像中心，在X射线机的接收端和发射端之间，置一网格状的模板，使其和接收端平行，成像一次，然后，沿着垂直于接收端影像增强屏的方向，移动物体，改变模板与成像平面的距离，再次成像，计算成像中心；

2)dx和dy的确定，在求取dx和dy的比值以及乘积的时，利用已知的四个点和四条边所能提供的组合，算出所有的dx*dy以及dx/dy，再分别加以平均，作为最终的dx*dy以及dx/dy的值；

3)外参数和畸变系数的标定，采用共面点的方法进行X射线机的标定，利用图像点和实际点的对应关系，可以计算出旋转矩阵中的r1，r2，r4，r5以及tx和ty，然后再利用正交性，计算出矩阵R中其余的所有元素，采用最速下降法，同时对焦距、tz以及畸变系数进行优化得到了所需的参数。

3.根据权利要求2所述的基于双X线机的实时立体定位方法，其特征是步骤(2)同时采集两路X射线图象的具体方法：采用一个可以两路同时输入的视频采集卡，把两个X线机的视频输出信号直接接到该卡上，通过该采集卡，把模拟信号转为数字信号。

4.根据权利要求2所述的基于双X线机的实时立体定位方法，其特征是步骤(3)采用边缘检测的方法对图象进行实时处理的具体方法是根据病灶的位置，自动选择能够包含穿刺针和病灶的矩形区域，在该区域内，采用边缘检测算子进行边缘检测，然后采用Hough变换进行直线检测，检测出穿刺针，在此基础上，根据灰度信息检测出针头的位置，对于病灶，根据上一帧的位置，在其周围进行搜索，利用灰度信息和形态信息。

5.基于双X射线机的同步实时动态跟踪立体定位装置，其特征是在水平向C臂的一端设水平向C臂X线球管，另一端设水平向影像增强器；在垂直向C臂的一端设垂直向C臂X线球管，另一端设垂直向影像增强器。

Images