CN104008269A - 一种基于人工标记物的手术导航系统空间自动注册方法 - Google Patents

Abstract

一种基于人工标记物的手术导航系统空间自动注册方法，该方法的步骤为：(1)在图像引导空间，对术前采集的医学影像图像进行表面重建，获得人头模型；采用不同的尺度空间因子，基于顶点预测器的原理，对模型进行光顺，构建模型的多尺度空间；(2)在多尺度空间中，计算相邻尺度的高斯差分，遍历当前尺度和相邻尺度的一环邻域寻找极值点；(3)对极值点进行标记物形状相关聚类，定位出标记物中心点；(4)在病人空间，使用双目相机实时获取病人头部图片，基于角点检测提取出标记物中心点，基于双目立体视觉原理恢复标记物中心点三维坐标；(5)基于随机迭代一致性（Ransac）方法匹配两空间标记物，计算欧式变换矩阵，完成手术导航系统的空间注册。

Description

一种基于人工标记物的手术导航系统空间自动注册方法

技术领域

本发明涉及计算机辅助手术领域，更具体地讲，涉及一种图像引导手术导航系统空间自动注册方法。 

背景技术

把病人在空间中的实际位置和三维模型中的虚拟位置进行精确配准的过程,在手术导航中称为注册。空间注册的精度直接影响着整个导航系统的精度，是手术前准备阶段中最重要的一步。基于点匹配的刚体配准法是目前手术导航系统临床应用中使用最多的空间注册法，其基本思想是分别在病人空间和三维图像空间中选取对应点集，通过数学原理计算出点集之间的变换关系，再扩展到病人空间中任意点在三维图像空间上的对应位置。 

目前基于点匹配的空间注册方法大都需要医生手动、顺序地在图像空间和图像引导空间选择人工标记物，这种方法获得的注册精度依赖于医生的经验，并且手动选取精度往往不稳定，通常需要多次选取然后求取平均值，增加了医生的工作量并延长了手术准备时间。 

本发明提出了一种基于人工标记物的手术导航系统空间自动注册方法，该方法设计了十字结构的圆形人工标记物，并利用标记物的形状特征，构建模型的多尺度空间，检测位于标记物圆周上的特征点，从而定位出标记物中心，在病人空间，利用角点检测和双目立体视觉原理，自动识别相同的标记物中心点，进而完成注册。这种方法不需要人工干预，稳定高效，大大节约了术前准备时间，提高了手术效率。 

发明内容

本发明的目的是：解决现有手术导航系统需要医生手动、顺序选择标记物 的费时费力的问题，提供了一种基于人工标记物的手术导航系统空间自动注册方法，该方法可以有效地完成病人术前医学影像空间与实际术中空间的注册，满足手术导航系统的速度和精度要求。 

为了实现上述目的，本发明的原理是：在图像引导空间，利用人工标记物和人头表面接触处不同的凹凸性，检测出曲率骤变的标记物圆周，根据几何特点计算标记物中心点。在病人空间，本发明设计了圆形十字标记物，以标记物中心点为圆点，以黑白两色填充，每90°变换颜色，这样十字标记物中心点为一稳定的角点。利于人工标记物的特殊设计，检测十字标记的中心角点。同时采用双目相机基于立体视觉的技术恢复这些标记物中心点的三维坐标，然后寻找两个空间的标记物对应关系，计算两空间的欧式变换矩阵，从而完成注册。 

本发明所采取的技术方案是：首先，在图像引导空间，对术前采集的医学影像图像进行表面重建，获得人头模型；采用不同的尺度空间因子，基于顶点预测器的原理，对模型进行光顺，构建模型的多尺度空间；在多尺度空间中，计算相邻尺度的高斯差分，遍历当前尺度和相邻尺度的一环邻域寻找极值点；然后对极值点进行标记物形状相关聚类，定位出标记物中心点。在病人空间，使用双目相机实时获取病人头部图片，基于角点检测提取出标记物中心点，基于双目立体视觉原理恢复标记物中心点三维坐标。最后，基于随机迭代一致性(Ransac)方法匹配两空间标记物，计算欧式变换矩阵，完成空间注册。 

本发明与现有技术相比的优点在于：充分利用了标记物的几何特征和灰度特征，并结合多尺度空间原理和立体视觉原理，自动定位出人工标记物，有效地提高了空间注册的精度和效率，减少了医生的工作量。 

附图说明

图1为本发明基于人工标记物的手术导航系统空间自动注册方法的流程图； 

图2为图像引导空间人工标记物定位的流程图。 

图3为定位出标记物示意图 

具体实施方法 

下面结合附图对本发明作进一步说明。 

本发明提出的基于人工标记物的手术导航系统空间自动注册方法的流程如图1所示，其步骤如下： 

(1)如图2所示，在图像引导空间，采用不同的尺度空间因子，基于顶点预测器的原理，对模型进行光顺，构建模型的多尺度空间；在多尺度空间中，计算相邻尺度的高斯差分，遍历当前尺度和相邻尺度的一环邻域寻找极值点；对极值点进行标记物形状相关聚类，定位出标记物中心点(如图3所示)： 

(1.1)设原始的人头模型为M，p为M中的一个顶点。对于尺度空间因子σ，M经过以σ为参数的高斯函数光顺后的结果为M_σ，p在M_σ中对应的顶点是p_σ，定义顶点p处的信号强度l(p,σ)为p在p_σ法向方向n_σ(p)上的位移： 

l(p,σ)＝＜n_σ(p),p-p_σ＞ 

(1.2)在以σ为因子的尺度空间中，估算p点对应顶点p_σ的法向n_σ(p)。 (是数学符号，表示对于属于M的任意一点p)，设N(p)是p的邻域，q是邻域内的一点，即n_q是q点的法向量。定义两点空间距离权值为||q-p||，同时考虑到网格M可能是非均匀采样的曲面，引入q的Voronoi面积作为顶点的面积权值Aq，利用如下公式估算出n_σ(p)： 

$n_{\mathrm{\σ}}\left(p\right)=\frac{{\mathrm{\Σ}}_{q\∈N\left(p\right)}{n}_q{A}_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}{{\mathrm{\Σ}}_{q\∈N\left(p\right)}{A}_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}$

其中G_σ(x)是高斯函数： 

$G_{\mathrm{\σ}}\left(x\right)=\frac{1}{2\mathrm{\π}{\mathrm{\σ}}^2}{e}^{-x^2/2{\mathrm{\σ}}^2}$

(1.3)计算p点的信号强度l(p,σ)。已知n_σ(p)，利用p在预测器∏q下的预测信号强度S(∏q(p))的加权平均来计算l(p,σ)，如下式： 

$l(p,\mathrm{\σ})=\frac{{\mathrm{\Σ}}_{q\∈N\left(p\right)}S\left({\mathrm{\Π}}_q\left(p\right)\right)A_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}{{\mathrm{\Σ}}_{q\∈N\left(p\right)}{A}_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}$

其中p点在预测器∏q下的预测信号强度S(∏q(p))为： 

S(∏_q(p))＝n_σ(p)·(q-p) 

(1.4)选取不同的尺度空间因子:σ₁，σ₂，….，σ_n，其中σ_i+1＝kσ_i，i＝1,2,…,n-1，分别做上述处理，就可以得到三维模型M在不同尺度空间中的表示L(M,σ_i),i＝1,2,…,n； 

(1.5)计算出高斯差分函数的值D(M,σ_i)： 

D(M,σ_i)＝L(M,σ_i+1)-L(M,σ_i) 

(1.6)在p点所在尺度的1环邻域，以及相邻尺度的对应顶点和该顶点的1环领域内，寻找具有最大和最小DOG值的顶点作特征点。 

(1.7)从检测到的特征点中使用Kmeans++的方法选取k个(k大于实际粘贴的标记点数目)特征点做初始聚类中心，然后扫描所有特征点并计算它到每一个类中心的距离，并将该点加入到距离最近的类中。 

(1.8)迭代循环所有的类，如果发现某两个类之间的距离(类间距离定义为两个类中心之间的距离)小于一个阈值d(例如20mm，约为标记点直径)，则将这两个类合并，直至所有类间距离都大于这个阈值为止。遍历所得的聚类，如果类中包含的特征点数少于阈值t，则删除该类，否则保留该类。 

(1.9)计算类内特征点到类中心的距离，如果距离在某一范围内(如7～13mm，约为标记点直径±3mm)，则保留该特征点，否则删除。处理完所有类后，重复(1.7)、(1.8)。 

(1.10)由于对应于标记物的类中的特征点都分布在标记物圆周上，任意取三个特征点，根据三角形三边垂直平分线的交点即三角形外接圆圆心的原理，计算外接圆圆心，圆心到三点中任一点的距离即外接圆半径。然后计算类内特征点到圆心的距离，若与半径差值小于一定阈值(标记物尺寸是已知的)，则认为该点合格，统计合格的特征点的比例，如果大于2/3，则认为该圆心合格。 

(2)在病人空间，使用双目相机实时获取病人头部图片，基于角点检测提取出标记物中心点，基于双目立体视觉原理恢复标记物中心点三维坐标。 

(3)基于随机迭代一致性(Ransac)方法匹配病人空间和图像引导空间标记物，计算欧式变换矩阵，完成空间注册： 

(3.1)假设在图像引导空间和病人空间得到的标记物中心点集分别是P、Q，从P、Q中分别任意选三个点，组成3个点对。 

(3.2)根据这3个点对计算出刚体变换矩阵T。 

(3.3)用矩阵T对P进行变换，遍历变换后P中的点，在Q中找其最近点，如果小于预定义的阈值，则认为(Pi,Qj)是正确匹配的点对，否则(Pi,Qj)是错误匹配的点对。正确匹配的点对记为局内点(Inlier)，错误的点对记为局外点(Outlier)。 

(3.4)如果局内点的数目达到一定数量(与总点对的比例达到一定阈值)，对所有的局内点采用基于单位四元数的参数估计方法重新计算M，否则迭代(3.1)-(3.3)步。 

本发明中未详细阐述的部分属于本领域技术人员的公知技术。 

Claims (2)

1.一种基于人工标记物的手术导航系统空间自动注册方法，其特征在于包括如下步骤：

(1)在图像引导空间，采用不同的尺度空间因子，基于顶点预测器的原理，对通过术前采集的医学影像图像进行表面重建的人头模型进行光顺，构建模型的多尺度空间；在多尺度空间中，计算相邻尺度的高斯差分，遍历当前尺度和相邻尺度的一环邻域寻找极值点；对极值点进行标记物形状相关聚类，定位出标记物中心点；

(2)在病人空间，使用双目相机实时获取病人头部图片，基于角点检测提取出标记物中心点，基于双目立体视觉原理恢复标记物中心点三维坐标；

(3)基于随机迭代一致性（Ransac）方法匹配两空间标记物，计算欧式变换矩阵，完成空间注册。

2.根据权利要求1所述的基于人工标记物的手术导航系统空间自动注册方法，其特征在于：步骤(1)中基于多尺度空间的标记物中心点定位方法，包括步骤如下：

(2.1)设原始的三维模型为M，p为M中的一个顶点;对于尺度空间因子σ，M经过以σ为参数的高斯函数光顺后的结果为M_σ，p在M_σ中对应的顶点是p_σ，定义顶点p处的信号强度l(p,σ)为p在p_σ法向方向n_σ(p)上的位移：

l(p,σ)＝＜n_σ(p),p-p_σ＞

(2.2)在以σ为因子的尺度空间中，估算p点对应顶点p_σ的法向n_σ(p);设N(p)是p的领域，q是邻域内的一点，即n_q是q点的法向量;定义两点空间距离权值为||q-p||，同时考虑到网格M可能是非均匀采样的曲面，引入q的Voronoi面积作为顶点的面积权值Aq，利用如下公式估算出n_σ(p)：

$n_{\mathrm{\σ}}\left(p\right)=\frac{{\mathrm{\Σ}}_{q\∈N\left(p\right)}{n}_q{A}_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}{{\mathrm{\Σ}}_{q\∈N\left(p\right)}{A}_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}$

其中G_σ(x)是高斯函数：

$G_{\mathrm{\σ}}\left(x\right)=\frac{-}{{2\mathrm{\π\σ}}^2}{e}^{-x^2/{2\mathrm{\σ}}^2}$

(2.3)计算p点的信号强度l(p,σ)；已知n_σ(p)，，利用p在预测器∏q下的预测信号强度S(∏q(p))的加权平均来计算l(p,σ)，如下式：

$l(p,\mathrm{\σ})=\frac{{\mathrm{\Σ}}_{q\∈N\left(p\right)}S\left({\mathrm{\Π}}_q\left(p\right)\right)A_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}{{\mathrm{\Σ}}_{q\∈N\left(p\right)}{A}_q{G}_{\mathrm{\σ}}\left(\right||q-p|\left|\right)}$

其中p点在预测器∏q下的预测信号强度S(∏q(p))为：

S(Π_q(p))＝n_σ(p)·(q-p)

(2.4)选取不同的尺度空间因子:σ₁，σ₂，….，σ_n，其中σ_i+1=kσ_i，i=1,2,…,n-1，分别做上述处理，就得到三维模型M在不同尺度空间中的表示L(M,σ_i),i=1,2,…,n；n为选取的不同尺度因子的个数;

(2.5)计算出高斯差分函数的值D(M,σ_i)：

D(M,σ_i)＝L(M,σ_i+1)-L(M,σ_i)

(2.6)在p点所在尺度的1环领域，以及相邻尺度的对应顶点和该顶点的1环领域内，寻找具有最大和最小DOG值的顶点作特征点;

(2.7)从检测到的特征点中使用Kmeans++的方法选取k个特征点做初始聚类中心，k大于实际粘贴的标记点数目，然后扫描所有特征点并计算它到每一个类中心的距离，并将该点加入到距离最近的类中;

(2.8)迭代循环所有的类，如果发现某两个类之间的距离小于一个阈值d，则将这两个类合并，直至所有类间距离都大于这个阈值为止;所述类间距离定义为两个类中心之间的距离;遍历所得的聚类，如果类中包含的特征点数少于阈值t，则删除该类，否则保留该类;

(2.9)计算类内特征点到类中心的距离，如果距离在某一阈值范围内，则保留该特征点，否则删除;处理完所有类后，再重复一遍（2.7）、（2.8）步;

(2.10)由于对应于标记物的类中的特征点都分布在标记物圆周上，任意取三个特征点，根据三角形三边垂直平分线的交点即三角形外接圆圆心的原理，计算外接圆圆心，圆心到三点中任一点的距离即外接圆半径;然后计算类内特征点到圆心的距离，若与半径差值小于一定阈值，则认为该点合格，统计合格的特征点的比例，如果大于2/3，则认为该圆心合格。