CN108804861A - 一种具有真实力反馈的脊柱微创手术培训系统及方法 - Google Patents

Abstract

本发明涉及一种具有真实力反馈的脊柱微创手术培训系统及方法，包括实验平台、手术器械和PC机，所述实验平台包括设置于底面的操作平台，在操作平台的面板中部设置脊柱物理模型夹持结构，用于夹持脊柱物理模型；在脊柱物理模型的正上方垂直于操作平台放置双目摄像机；所述双目摄像机通过双目摄像机支架设置于操作平台上，且通过电缆连接PC机；在所述手术器械的末端设置参考架。由于本发明的力反馈是真实的，可以给医生带来逼真的触觉感受，大大缩短培训医生的时间；使培训者多角度、多方位的观察椎弓根螺钉的植入情况，给培训者真实的视觉反馈，真正达到培训的目的；并且满足脊柱微创手术培训系统高精度的要求。

Description

一种具有真实力反馈的脊柱微创手术培训系统及方法

技术领域

本发明涉及医疗培训和计算机视觉领域，具体地说是一种具有真实力反馈的脊柱微创手术培训系统及方法。

背景技术

当前国内外绝大多数虚拟仿真培训系统主要是针对腹腔镜和内窥镜手术中软组织切割和变形的模拟和培训，对于骨组织的手术模拟和培训，只有少量国外学者进行了相关方面的研究。例如，D.Morris等开发了一套适用于颞骨手术培训的虚拟仿真系统。该系统集成了触觉、视觉和听觉等感官反馈，能够对一些简单的钻孔操作进行模拟和培训。Petersik等研究了一种基于多点碰撞检测的力触觉绘制算法，基于该算法他们所开发的岩骨手术仿真系统在钻孔模拟时能够提供非常逼真的振动感。J.Cordes等研发一套治疗脊柱疾病的虚拟仿真培训系统，通过人机交互界面和专家经验等，该系统实现了对脊柱手术操作的仿真。

以上培训系统存在一个共性的问题，无法向培训者提供真实准确的力反馈，其主要原因是骨组织手术操作的力学模型比较复杂，不仅与骨密度、骨厚度等骨组织属性有关，同时也会受到手术操作时切割速度、深度以及手术器械的微动频率等因素的影响，所以很难建立骨组织力学仿真的物理模型。并且脊柱手术对精度的要求更为严格，螺钉植入的宽度只有10mm～15mm，稍有偏差就会损伤脊柱神经，造成不可逆的后果。所以在手术培训系统中，获得精确、逼真的力反馈至关重要，这也一直都是限制虚拟手术培训系统质量提升的瓶颈所在。

现阶段动态跟踪算法无法同时满足准确性、鲁棒性、实时性要求。例如，传统的结合kalman滤波的meanshift算法，虽较好的解决了目标遮挡和目标快速移动的问题，但是在目标运动方向突然变化时，鲁棒性较差。结合传统不变矩和kalman的跟踪算法鲁棒性较强，但是实时性比较差且不适合目标遮挡的情况。

发明内容

针对现有技术的不足，本发明提供一种通过真实的手术器械和脊柱物理模型，结合仿射不变矩特征匹配和Kalman滤波的meanshift运动目标跟踪算法，获得真实的力反馈的脊柱微创手术培训系统及方法，给培训者提供更真实的手术体验。

本发明为实现上述目的所采用的技术方案是：

一种具有真实力反馈的脊柱微创手术培训系统，包括实验平台、手术器械6和PC机，所述实验平台包括设置于底面的操作平台1，在操作平台1的面板中部设置脊柱物理模型夹持结构2，用于夹持脊柱物理模型3；在脊柱物理模型3的正上方垂直于操作平台1放置双目摄像机4；所述双目摄像机4通过双目摄像机支架5设置于操作平台1上，且通过电缆连接PC机；在所述手术器械6的末端设置参考架7。

所述参考架7包括四个标志球和两个横梁，其中每个横梁两端设置两个颜色相同的标志球，且两个横梁上的标志球颜色不同，每个标志球到参考架中心的距离均相同。

在所述操作平台1上设置手术器械托架，用于放置手术器械6。

一种具有真实力反馈的脊柱微创手术培训方法，包括以下步骤：

步骤1：采用基于平面的棋盘格标定算法对双目摄像机进行标定，获得摄像机坐标系与图像像素坐标系之间的转换关系；

步骤2：通过双目摄像机对脊柱模型进行三维重建，获得实际空间脊柱模型的局部点云；

步骤3：对脊柱模型进行CT扫描，根据扫描得到的数据进行三维重建，获得脊柱的全局点云，将重建的脊柱导入虚拟界面；

步骤4：通过ICP空间配准算法将实际空间和虚拟空间进行配准，以获得实际空间和虚拟空间的转换关系；

步骤5：通过双目摄像机对手术器械进行实时动态跟踪，在虚拟界面实时显示手术器械相对于脊柱的位置。

所述实时动态跟踪包括以下过程：

步骤1：手动选定摄像机采集的视频中运动目标的初始化位置，确定kalman滤波器初始状态向量并初始化其他必要参数，计算目标模型的颜色直方图；

步骤2：利用初始化后的kalman滤波器预测当前时刻、当前帧的候选模型位置；

步骤3：以kalman滤波器预测的目标候选模型位置做为Mean Shift算法的迭代起始点，开始Mean Shift算法迭代过程，直至达到收敛条件，停止迭代，计算候选模型的颜色直方图；

步骤4：利用相似性度量函数进行相似性测定，通过Mean Shift算法获得最优候选模型位置y1，计算此时的相似度为P1，如果P1>T，将y1作为kalman滤波器的观测向量，对kalman滤波器进行更新；否则启动仿射变换不变矩算法，并与圆的标准仿射变换不变矩比较，得到新的候选模型位置y2；

步骤5：计算新候选模型的颜色直方图，并利用相似性度量函数计算P2；比较P1、P2的大小，将相似度大的候选模型位置，作为kalman滤波器的观测向量，对kalman滤波器进行更新。

所述其他必要参数包括：kalman滤波器中误差协方差矩阵初始值、kalman滤波器中状态转移矩阵、kalman滤波器中观测矩阵、kalman滤波器中系统噪声协方差矩阵和kalman滤波器中观测噪声协方差矩阵。

所述仿射变换不变矩算法为：

其中μ_pq为轮廓曲线的尺度规范化矩，其中p与q的和为轮廓曲线的尺度规范化矩的阶次，I₁、I₂、I₃为由前三阶轮廓曲线的尺度规范化矩构成的三个仿射不变矩。

所述轮廓曲线的尺度规范化矩为：

平面曲线L的(p+q)阶矩定义为：

m_p,q＝∫x^py^qds,(p,q＝0,1,2,…)

其中，ds为曲线L的弧微分，m_0,0就为轮廓周长；

平面曲线L的(p+q)阶中心矩定义为：

υ_p,q＝∫(x-x₀)^p(y-y₀)^qds,(p,q＝0,1,2,…)

其中，(x₀,y₀)是质心坐标，

轮廓曲线的尺度规范化矩为：

其中，m_0,0为零阶矩，υ_p,q为(p+q)阶中心矩。

所述与圆的标准仿射变换不变矩比较的过程如下：

步骤1：提取当前帧的图像进行高斯去噪，采用canny算法进行边缘检测；

步骤2：计算所有封闭边缘的仿射不变矩，并与圆的标准仿射不变矩比较，找到目标关键部位标志球的边缘；

步骤3：通过椭圆拟合算法找到标志球边缘的圆心，根据颜色信息，将同一种颜色的标志球圆心坐标的平均值作为新的候选模型位置y2。

所述收敛条件为上一帧目标的位置与Mean Shift算法搜索到的局部最优位置的欧式距离小于0.5。

本发明具有以下有益效果及优点：

1.由于本发明的力反馈是真实的，可以给医生带来逼真的触觉感受，大大缩短培训医生的时间；

2.本发明可以使培训者多角度、多方位的观察椎弓根螺钉的植入情况，给培训者真实的视觉反馈，真正达到培训的目的；

3.本发明满足脊柱微创手术培训系统高精度的要求。

4.本发明结合仿射不变矩特征匹配和Kalman滤波的meanshift运动目标跟踪算法，极大提高目标跟踪的准确性、实时性和鲁棒性。

附图说明

图1是本发明的系统结构连接图；

图2是本发明的手术器械示意图；

图3是本发明的方法流程图；

图4是提取实际空间脊柱模型的局部点云流程图；

图5是结合仿射不变矩特征匹配和kalman滤波的meanshift运动目标跟踪流程图；

图6是参考架上标志球三维测距的流程图。

具体实施方式

下面结合附图及实例对本发明做进一步的详细说明。

如图1所示为本发明的系统结构连接图。

脊柱微创手术培训系统，包括实验平台、手术器械6和PC机，实验平台包括设置于底面的操作平台1，用于盛放系统其他设备，在操作平台1的面板中间位置设置脊柱物理模型夹持结构2，脊柱物理模型夹持结构2用于夹持脊柱物理模型3，该脊柱物理模型夹持结构2可以调节夹持力度的松紧。

在脊柱物理模型3的正上方放置双目摄像机4，双目摄像机4的摄像头垂直于操作平台1，且双目摄像机4通过双目摄像机支架5设置于操作平台1上，双目摄像机4通过电缆连接PC机，与PC机进行通信。手术器械6包括尖端和末端，在手术器械6的末端设置参考架7。

在操作平台1上设置手术器械托架，用于放置手术器械6。

脊柱物理模型3是模拟人类真实的脊柱通过3D打印而成。

双目摄像机4的分辨率至少为1280×960。

如图2所示为本发明的手术器械示意图。

所述参考架7包括四个标志球和两个横梁，其中每个横梁两端设置两个颜色相同的标志球，且两个横梁上的标志球颜色不同，每个标志球到参考架中心的距离均相同。

如图3所示为本发明的方法流程图。包括以下几个步骤：

步骤1：采用基于平面的棋盘格标定算法对双目摄像机进行标定，获得摄像机坐标系与图像像素坐标系之间的转换关系；

步骤2：通过双目摄像机对脊柱模型进行三维重建，获得实际空间脊柱模型的局部点云；

步骤3：对脊柱模型进行CT扫描，根据扫描得到的数据进行三维重建，获得脊柱的全局点云，将重建的脊柱导入虚拟界面；

步骤4：通过ICP空间配准算法将实际空间和虚拟空间进行配准，以获得实际空间和虚拟空间的转换关系；

步骤5：通过双目摄像机对手术器械进行实时动态跟踪，在虚拟界面实时显示手术器械相对于脊柱的位置。

如上所述通过棋盘格平面模板对双目摄像机进行标定，获得实际空间坐标系及图像像素坐标系的转化关系，具体实施步骤如下：

步骤1：平面棋盘格标定板的制定。本系统采用的是定制棋盘格平面标定板，交点的个数为7×10、单位格大小为25mm×25mm。

步骤:2：图像采集。利用两台摄像机同时对标定板的不同角度进行采集，一般采集的图像越多，标定的结果越好，建议10-20组，本系统采集了15组。为了得到更好的标定效果，要求标定板画面占据整个图像的1/2以上，标定板的倾斜角度相对于操作平台小于45度。为了更好地建立畸变模型，使标定板尽可能的靠近左右摄像机公共视场的边缘。

步骤3：MATLAB标定过程。首先，读入左右摄像机图片，设定标定板的单位格的实际尺寸25mm。然后，设定畸变模型，可以设定只考虑径向畸变或径向畸变和切向畸变的模型。最后，利用MATLAB中的工具箱实现双目摄像机的标定。

本系统设定的畸变模型如下：

其中为校正后的新位置，(u，v)为原始位置。

步骤4：提高标定精度。首先，删除投影误差大、模糊不清、角点提取不准确和倾斜角度过大的图像。其次，在图像少于10对、模板没有充分的覆盖摄像机视场的每一个角落和模板相对于摄像机没有足够的角度变化时，添加图片。通过以上处理后重新标定。

如图4所示为提取实际空间脊柱模型的局部点云的流程图。具体实施过程如下：

通过左右摄像机分别采集图片，对左右图像进行立体校正，采用canny算法进行边缘提取，通过灰度差绝对值和(SAD)算法进行立体匹配，计算初始视差图，然后采用左右一致性滤波器，置信滤波器和唯一性滤波器对初始化视差图进行优化，对三维脊柱模型进行重建，提取实际空间的脊柱物理模型局部点云。

如上所述实际空间坐标系(摄像机坐标系)与虚拟空间的坐标系进行空间配准。具体实施过程如下：

首先采用实际空间和虚拟空间的脊柱的特征点的三维坐标，进行粗配准，然后采用ICP算法对实际空间脊柱的局部三维点云与虚拟空间脊柱的全局三维点云进行精配准，得到实际空间与虚拟空间的空间转换。

如图5所示为结合仿射不变矩特征匹配和kalman滤波的meanshift运动目标跟踪流程图。具体实施步骤如下：

步骤1：手动选定摄像机采集的视频中运动目标的初始化位置，计算目标模型的颜色直方图，并确定kalman滤波器初始状态向量并初始化其他必要参数。

Kalman滤波器初始状态向量设定为：[x₀ y₀ 0 0]，其中x₀和y₀表示初始目标模型的中心坐标。

Kalman滤波器其他参数设定如下：

(1)Kalman滤波器中误差协方差矩阵初始值p₀＝0；

(2)Kalman滤波器中状态转移矩阵其中Δt为相邻两帧图像间的时间间隔，本系统Δt等于0.1s；

(3)Kalman滤波器中观测矩阵为

(4)Kalman滤波器中系统噪声协方差矩阵和观测噪声协方差矩阵分别设定为四阶和二阶单位阵。

步骤2：利用初始化后的kalman滤波器预测当前时刻、当前帧的候选模型位置；

步骤3：以kalman滤波器预测的目标候选模型位置为meanshift算法的迭代起始点，开始meanshift算法迭代过程，直至上一帧目标的位置与meanshift算法搜索到的局部最优位置的欧式距离小于0.5，停止迭代，计算候选模型的颜色直方图；

步骤4：利用相似性度量函数进行相似性测定，通过meanshift算法获得最优候选模型位置y1，计算此时的相似度为P1，如果P1>0.4，将y1作为kalman滤波器的观测向量，对kalman滤波器进行更新；否则启动仿射变换不变矩算法，并与圆的标准仿射变换不变矩比较，得到新的候选模型位置y2；

步骤5：计算新候选模型的颜色直方图，并利用相似性度量函数计算P2；比较P1、P2的大小，将相似度大的候选模型位置，作为kalman滤波器的观测向量，对kalman滤波器进行更新。

如上所述meanshift算法的具体实施过程如下：

步骤1：初始帧目标模型和当前帧候选模型的建立；

以一定间隔的颜色值为单位，将取值为像素颜色值的特征空间分为多个特征值。在包含目标的搜索窗口中，根据第u个特征值的概率为目标和候选模型建模。

目标模型建模定义如下：

候选模型建模定义如下：

式中x₀和y₀是目标和候选模型搜索窗口(n个像素)的中心像素坐标，x_i是第i个像素的坐标；k(||x||²)是核函数，h表示核函数的带宽，一般等于窗口宽度的一半；函数b和δ的作用是判断x_i的颜色值是否属于特征值u；C和C_h是一个标准化的常量，使得所有的特征值的概率和为1。

步骤2：目标模型与候选模型相似度计算；

相似性函数描述初始帧目标模型和当前帧候选模型的相似性度量，定义如下：

其中为当前帧的候选模型，为初始目标模型，m为特征值u的个数。

步骤3：meanshift算法反复迭代，最后得到在当前帧目标的最优位置。

通过对相似性函数求最大值，可以推导出meanshift向量

Meanshift向量所指的方向就是相似度增大的方向，反复迭代最后得到在当前帧目标的最优位置。

传统的不变矩仅仅在目标旋转、平移和尺度变化条件下有效，但是在实际的三维空间中，在双目摄像机不同的视角下，目标成像发生仿射变换的情况下是不起作用的，于是本文在传统不变矩的基础上提出仿射不变矩。本系统采用在仿射变换条件下保持不变的三个仿射不变矩：

其中μ_pq为轮廓曲线的尺度规范化矩。

轮廓曲线的尺度规范化矩μ_pq的定义如下：

平面曲线L的(p+q)阶矩定义为：

m_p,q＝∫x^py^qds,(p,q＝0,1,2,…)

其中，ds为曲线L的弧微分，m_0,0为轮廓周长；

平面曲线L的(p+q)阶中心矩定义为：

υ_p,q＝∫(x-x₀)^p(y-y₀)^qds,(p,q＝0,1,2,…)

其中，(x₀,y₀)是质心坐标，

轮廓曲线的尺度规范化矩为：

其中，m_0,0为零阶矩，υ_p,q为(p+q)阶中心矩。

如上所述系统动态跟踪算法中候选模型位置y2的确定步骤如下：

步骤1：提取当前帧的图像进行高斯去噪，采用canny算法进行边缘检测；

步骤2：计算所有封闭边缘的仿射不变矩，仿射不变矩与圆的标准仿射不变矩比较，找到目标关键部位标志球的边缘。圆的标准仿射不变矩满足|I₁-0.006332|≤0.0003，|I₂|<0.0000001,|I₃|<0.0000001，所以只要满足这三个条件的边缘就为关键部位；

步骤3：通过椭圆拟合算法找到标志球边缘的圆心，根据颜色信息，将同一种颜色的标志球圆心坐标的平均值作为y2。

如图6所示为参考架上标志球三维测距的流程图。在动态跟踪算法的基础上，提取图像中目标模型中标志球的像素坐标，通过双目摄像机对标志球进行三维测量。首先，采用自适应二值化，高斯去噪及形态学算法对摄像机采集的图片进行图像处理。然后，采用Hough变换算法提取标志球的圆心坐标，如果提取的精度不高，可以修改Hough变换算法的参数，从新提取。最后，通过双目摄像机的重建原理对标志球进行三维测量，获得参考架上四个标志球在实际空间的位置，进而求出手术器械在实际空间的位置。

Claims (10)

1.一种具有真实力反馈的脊柱微创手术培训系统，包括实验平台、手术器械(6)和PC机，其特征在于，所述实验平台包括设置于底面的操作平台(1)，在操作平台(1)的面板中部设置脊柱物理模型夹持结构(2)，用于夹持脊柱物理模型(3)；在脊柱物理模型(3)的正上方垂直于操作平台(1)放置双目摄像机(4)；所述双目摄像机(4)通过双目摄像机支架(5)设置于操作平台(1)上，且通过电缆连接PC机；在所述手术器械(6)的末端设置参考架(7)。

2.根据权利要求1所述的具有真实力反馈的脊柱微创手术培训系统，其特征在于：所述参考架(7)包括四个标志球和两个横梁，其中每个横梁两端设置两个颜色相同的标志球，且两个横梁上的标志球颜色不同，每个标志球到参考架中心的距离均相同。

3.根据权利要求1所述的具有真实力反馈的脊柱微创手术培训系统，其特征在于：在所述操作平台(1)上设置手术器械托架，用于放置手术器械(6)。

4.一种根据权利要求1～3所述系统的脊柱微创手术培训方法，其特征在于：包括以下步骤：

步骤1：采用基于平面的棋盘格标定算法对双目摄像机进行标定，获得摄像机坐标系与图像像素坐标系之间的转换关系；

步骤2：通过双目摄像机对脊柱模型进行三维重建，获得实际空间脊柱模型的局部点云；

步骤3：对脊柱模型进行CT扫描，根据扫描得到的数据进行三维重建，获得脊柱的全局点云，将重建的脊柱导入虚拟界面；

步骤4：通过ICP空间配准算法将实际空间和虚拟空间进行配准，以获得实际空间和虚拟空间的转换关系；

步骤5：通过双目摄像机对手术器械进行实时动态跟踪，在虚拟界面实时显示手术器械相对于脊柱的位置。

5.根据权利要求4所述的脊柱微创手术培训方法，其特征在于：所述实时动态跟踪包括以下过程：

步骤1：手动选定摄像机采集的视频中运动目标的初始化位置，确定kalman滤波器初始状态向量并初始化其他必要参数，计算目标模型的颜色直方图；

步骤2：利用初始化后的kalman滤波器预测当前时刻、当前帧的候选模型位置；

步骤3：以kalman滤波器预测的目标候选模型位置做为Mean Shift算法的迭代起始点，开始Mean Shift算法迭代过程，直至达到收敛条件，停止迭代，计算候选模型的颜色直方图；

步骤4：利用相似性度量函数进行相似性测定，通过Mean Shift算法获得最优候选模型位置y1，计算此时的相似度为P1，如果P1>T，将y1作为kalman滤波器的观测向量，对kalman滤波器进行更新；否则启动仿射变换不变矩算法，并与圆的标准仿射变换不变矩比较，得到新的候选模型位置y2；

步骤5：计算新候选模型的颜色直方图，并利用相似性度量函数计算P2；比较P1、P2的大小，将相似度大的候选模型位置，作为kalman滤波器的观测向量，对kalman滤波器进行更新。

6.根据权利要求5所述的脊柱微创手术培训方法，其特征在于，所述其他必要参数包括：kalman滤波器中误差协方差矩阵初始值、kalman滤波器中状态转移矩阵、kalman滤波器中观测矩阵、kalman滤波器中系统噪声协方差矩阵和kalman滤波器中观测噪声协方差矩阵。

7.根据权利要求5所述的脊柱微创手术培训方法，其特征在于：所述仿射变换不变矩算法为：

其中μ_pq为轮廓曲线的尺度规范化矩，其中p与q的和为轮廓曲线的尺度规范化矩的阶次，I₁、I₂、I₃为由前三阶轮廓曲线的尺度规范化矩构成的三个仿射不变矩。

8.根据权利要求7所述的脊柱微创手术培训方法，其特征在于，所述轮廓曲线的尺度规范化矩为：

平面曲线L的(p+q)阶矩定义为：

m_p,q＝∫x^py^qds,(p,q＝0,1,2,…)

其中，ds为曲线L的弧微分，m_0,0就为轮廓周长；

平面曲线L的(p+q)阶中心矩定义为：

υ_p,q＝∫(x-x₀)^p(y-y₀)^qds,(p,q＝0,1,2,…)

其中，(x₀,y₀)是质心坐标，

轮廓曲线的尺度规范化矩为：

其中，m_0,0为零阶矩，υ_p,q为(p+q)阶中心矩。

9.根据权利要求5所述的脊柱微创手术培训方法，其特征在于，所述与圆的标准仿射变换不变矩比较的过程如下：

步骤1：提取当前帧的图像进行高斯去噪，采用canny算法进行边缘检测；

步骤2：计算所有封闭边缘的仿射不变矩，并与圆的标准仿射不变矩比较，找到目标关键部位标志球的边缘；

步骤3：通过椭圆拟合算法找到标志球边缘的圆心，根据颜色信息，将同一种颜色的标志球圆心坐标的平均值作为新的候选模型位置y2。

10.根据权利要求5所述的脊柱微创手术培训方法，其特征在于，所述收敛条件为上一帧目标的位置与Mean Shift算法搜索到的局部最优位置的欧式距离小于0.5。