CN104063057A - 一种支持增强力触觉反馈的涡卷弹簧仿真系统 - Google Patents

Abstract

本发明公开了一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征是采用一种支持增强力触觉反馈的涡卷弹簧模型，在交互过程中，虚拟柔性体产生实时变形力触觉信息，该力触觉信息与有色状态噪声信号相加，作为反馈输出信号。本发明系统采用的涡卷弹簧每层的拉伸或压缩变形量计算方法相同，计算简单，加快了变形计算速度，能够保证变形仿真具有较高精度，通过调节涡卷弹簧第一层半径，任意两层弹簧之间的间距和弹簧丝的直径等，就可模拟不同类型的柔性体；加入的有色状态噪声信号能提高人的力触觉感知能力；为开发基于力触交互的虚拟手术平台提供了解决思路。

Description

一种支持增强力触觉反馈的涡卷弹簧仿真系统

技术领域

本发明涉及一种支持增强力触觉反馈的涡卷弹簧仿真系统，属于虚拟现实人机交互领域。

背景技术

现代医疗技术的日新月异，对外科手术的精细程度要求越来越苛刻，但一名优秀的外科医生需经过大量训练、实验，花费8-10年才能成长起来。和传统以动物和尸体作为实验对象不同，虚拟手术在模拟虚拟柔性体的各种物理、生理特性的同时，还可得到实验数据，为分析总结提供依据，极大地改善了手术训练效果，同时还可根据医生的需要提供不同的可以反复利用的实验对象，大大降低了外科医生手术训练成本。

虚拟手术仿真除在视觉上给使用者营造一种沉浸感、真实感之外，还要提供一种力触觉反馈机制，使人们在看见虚拟对象的同时，还可触摸到它的几何形状和物理特征，视觉和力触觉的有机结合使虚拟手术仿真更加真实可信。用于力触觉交互的形变模型是虚拟手术仿真系统的关键技术，国内外学者在这方面已做了很多工作。弹簧-质点模型虽然不需要连续参数化，且在静态和动态都可进行分析，但当变形较大时不能准确描述柔性体的形变，系统稳定低。有限元模型虽然计算精度较高，算法适应性也比较广泛，但模型复杂度高、计算时耗多，实时性较差。边界元模型通过对模型边界进行离散来降低问题的维数，减少自由度数目使得计算简化，但由于考虑内部结点位移使得计算复杂度高。以上说明，这些常用的柔性体变形物理建模方法均存在计算较为繁杂和仿真精度不高等问题，继而影响了计算的实时性和有效性。因此，在保证形变建模方法准确性的同时，简化算法，提高计算实时性是当前虚拟力触觉交互系统亟待解决的首要问题。

人机工程学的大量研究表明：人的力触觉感知存在着信号与噪声之间的协同行为，即适量噪声的加入会提高人的力触觉感知能力，这主要原因是因为人手皮肤对一个持续作用力具有很强的适应性，而对噪声信号比较敏感，所以当加入合适的噪声信号，会产生比较真实的感觉。根据这一原理，为了提高虚拟现实人机力触觉交互过程的稳定性和感知真实性，提出了一种支持增强力触觉反馈的涡卷弹簧仿真系统。

现代随机信号的大量研究表明：实际的观测和数据处理中，观测误差和动力学模型往往是具有一定时间相关或空间相关性的有色噪声。因此我们加入适量的有色的状态噪声，对于探索建模效果的真实性以及精确性的提高，形变模型的好坏是手术成功的关键，同时柔性体变形所需实时性和准确性直接决定了虚拟手术能否成功进行。

发明内容

本发明所要解决的技术问题在于克服现有技术的不足，提供一种支持增强力触觉反馈的涡卷弹簧仿真系统。

为解决上述技术问题，本发明所采用的技术方案是：

一种支持增强力触觉反馈的涡卷弹簧仿真系统，包括：计算机主机和力触觉交互设备，在计算机主机上连接有显示器，所述计算机主机至少包括硬盘和1394卡，其特征在于，所述硬盘至少包括：用于检测虚拟代理位置的位置检测模块、用于检测虚拟代理与虚拟柔性体是否发生碰撞的碰撞检测模块、力触觉计算变形模块、根据力触觉计算变形模块输出的最新力触觉变形信息进行图形显示的图形刷新模块及力触觉信息反馈模块，所述硬盘通过1394卡及1394连接线与力触觉交互设备连接，用于将力触觉计算变形模块产生的柔性体表面的变形信息传输至力触觉交互设备，

所述力触觉计算变形模块243用于在检测到虚拟代理碰撞到虚拟柔性体表面上任何一点时，在给定虚拟外力F作用下，虚拟代理与虚拟柔性体交互的局部区域内部，填充支持力触觉反馈的涡卷弹簧模型；在交互过程中，虚拟柔性体产生实时变形力触觉信息，该力触觉信息与有色状态噪声信号η(t)相加，即采用支持增强力触觉反馈的涡卷弹簧模型，产生的力触觉信息作为反馈输出信号；所述涡卷弹簧模型的建模步骤如下：

步骤1，建立空间直角坐标系，确定柔性体表面上任意点处安放的涡卷弹簧模型，其过程如下：

步骤1-1，建立空间直角坐标系，

在给定虚拟外力F作用下，在柔性体表面任意点O₁处安放一涡卷弹簧，在距离任意点O₁为r处，取涡卷弹簧第一层的一个外侧点P₀，即r为涡卷弹簧第一层的半径，以任意点O₁为原点，射线O₁P₀所在方向为X轴正方向，建立XYZ空间坐标系；

步骤1-2，按层设置涡卷弹簧，涡卷弹簧任意两层的间距均为d，涡卷弹簧每层弹簧丝的直径均为h；

所述涡卷弹簧的总层数包含两种，分别是奇数和偶数，令n＝1,2,3……S,S为自然数；则所述涡卷弹簧具体描述如下：

A.当该涡卷弹簧总层数为奇数时，设层数为2n+1，在Z轴上等间距取2n+1个点，分别为O₁,O₂,…,O_i,…,O_2n,O_2n+1，间距为d；以O_i为圆心，设置涡卷弹簧第i层，且该涡卷弹簧第i层的半径为R_i，其表达式如下：

$R_i=\left\{\begin{array}{cc}r+(i-1)\mathrm{\ρ}& i\≤n+1\\ r+(2n-i+1)\mathrm{\ρ}& i>n+1\end{array}\right.$

也即前n+1层涡卷弹簧的半径R_i组成以r为首项，以ρ为公差的等差数列，ρ为实数；

B.当该涡卷弹簧总层数为偶数时，设层数为2n，则在Z轴上等间距取2n个点，分别为O₁,O₂,…,O_i,…,O_2n-1,O_2n，间距为d；以O_i为圆心，设置涡卷弹簧第i层，且该涡卷弹簧第i层的半径为R_i，其表达式如下：

$R_i=\left\{\begin{array}{cc}r+(i-1)\mathrm{\ρ}& i\≤n\\ r+(2n-i)\mathrm{\ρ}& i>n\end{array}\right.$

也即前n层涡卷弹簧的半径R_i组成以r为首项，以ρ为公差的等差数列；

涡卷弹簧每层的截面惯性矩均为涡卷弹簧各层的挠度值相同，挠度值为X′，弹性模量E取决于柔性体的材质，且柔性体的材质均相同；

步骤2，确定涡卷弹簧任一层所消耗的外力，

假定给定虚拟外力F的作用线与涡卷弹簧中心线一致，且在给定虚拟外力F作用下，如果涡卷弹簧共有前C层产生变形，则第C层称为变形截止层；被拉伸或压缩时，变形截止层产生的变形量不大于挠度值X′，涡卷弹簧前C-1层中，任意一层的变形量均等于挠度值X′；

步骤2-1，确定涡卷弹簧前C-1层消耗的外力，

涡卷弹簧第一层消耗的外力F₁为：

$F_1=\left\{\begin{array}{cc}{X}^{\′}\·K_1& C>1\\ F& C=1\end{array}\right.$

其中，K₁表示涡卷弹簧第一层的刚度系数；

除第一层和变形截止层第C层外，涡卷弹簧第i层消耗的外力F_i为：

$F_i=\frac{{\mathrm{Gh}}^4}{8{R_i}^3}(d-h)$

其中，G表示该涡卷弹簧的切变模量，2≤i≤C-1；

步骤2-2，确定涡卷弹簧第一层和变形截止层第C层的刚度系数，分别表示为K₁和K_C；其表达式分别如下：

$K_1=\frac{\mathrm{EI}}{l_1}$

$K_C=\frac{\mathrm{EI}}{l_C}$

其中，l₁、l_C分别为该涡卷弹簧第一层和变形截止层的工作长度，其计算式分别如下：

l₁＝2πR₁

$l_C=2\mathrm{\π}\underset{i=1}{\overset{C}{\mathrm{\Σ}}}{R}_i$

步骤3，涡卷弹簧变形截止层产生拉伸或压缩的变形量为X_C，其计算公式式为：

$X_C=\frac{F-\underset{i=1}{\overset{C-1}{\mathrm{\Σ}}}{F}_i}{K_C}$

步骤4，将涡卷弹簧变形量相加，所得变形量之和X等效为虚拟柔性体表面的变形：

步骤5，所述有色状态噪声信号η(t)，其表达式为：

η(t)＝D(q^-1)C(q^-1)ξ(t)

式中，

D(q^-1)、C(q^-1)分别表示单位滞后算子q^-1的多项式，其表达式分别为：

D(q^-1)＝1+D₁q^-1+D₂q^-2+D₃q^-3

C(q^-1)＝1+C₁q^-1+C₂q^-2

其中，D₁、D₂、D₃、C₁、C₂皆是常数；

ξ(t)是零均值白噪声，t为时间。

步骤2中，所述常数D₁、D₂、D₃、C₁、C₂的取值分别为：

D₁＝-1.5,D₂＝0.7,D₃＝0.1,C₁＝0.5,C₂＝0.2。

所述涡卷弹簧产生拉伸或压缩变形总计需要的时延时间不大于1ms，也即满足刷新频率不低于1000Hz的要求。

所有涡卷弹簧层数中，每层产生拉伸或压缩变形需要的时延时间构成等比数列，即满足：

t_i＝p^i-1t₁，

其中，t_i表示涡卷弹簧第i层产生拉伸或压缩变形需要的时延时间，p为等比数列的公比，t₁为涡卷弹簧第一层产生拉伸或压缩变形需要的时延时间，1≤i≤C。

在所述涡卷弹簧上任一层产生拉伸或压缩变形达到给定的挠度值后，其下一层才开始产生拉伸或压缩变形。

力触觉交互设备采用PHANTOM OMNI力触觉交互设备，所述PHANTOMOMNI力触觉交互设备通过APS开关电源供电，通过APS开关电源之后，将供电电源220V的交流电转化为18V，2.22A所需的电压，显示卡、1394卡、硬盘均通过总线实现双向传输，显示器通过VGA连接线与显示卡连接。

本发明的有益效果是：本发明公开了一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征是采用一种支持增强力触觉反馈的涡卷弹簧模型，在交互过程中，虚拟柔性体产生实时变形力触觉信息，该力触觉信息与有色状态噪声信号相加，作为反馈输出信号。该系统采用的涡卷弹簧模型每层涡卷弹簧的拉伸或压缩变形量计算方法相同，计算简单，加快了变形计算速度，能够保证变形仿真具有较高精度，通过调节涡卷弹簧第一层半径，任意两层弹簧之间的间距和弹簧丝的直径等，就可模拟不同类型的柔性体；加入的有色状态噪声信号能提高人的力触觉感知能力。该系统为开发基于力触交互的虚拟手术平台提供了解决思路。

附图说明

图1是系统电路原理图；

图2是涡卷弹簧模型示意图；

图3是柔性体变形仿真流程图；

图4是涡卷弹簧模型构建方法流程图；

图5是涡卷弹簧模型变形、变形量与时延时间关系示意图。

具体实施方式

下面结合附图所示流程对本发明提出的一种支持增强力触觉反馈的涡卷弹簧仿真系统进行详细说明：

如图1所示的系统电路原理图。一种虚拟柔性体的力触觉反馈仿真系统，包括：计算机主机2和力触觉交互设备3，在计算机主机2上连接有显示器1，所述计算机主机2至少包括硬盘24和1394卡23，其特征在于，所述硬盘24至少包括：用于检测虚拟代理位置的位置检测模块241、用于检测虚拟代理与虚拟柔性体是否发生碰撞的碰撞检测模块242、力触觉计算变形模块243、根据力触觉计算变形模块243发出的最新力触觉变形的信息进行图形显示的图形刷新模块244及力触觉信息反馈模块245，所述硬盘24通过1394卡23及1394连接线与力触觉交互设备3连接，用于将力触觉计算变形模块243产生的柔性体表面的变形信息传输至力触觉交互设备3，

所述力触觉计算变形模块243用于在检测到虚拟代理碰撞到虚拟柔性体表面上任何一点时，在给定虚拟外力F作用下，虚拟代理与虚拟柔性体交互的局部区域内部，填充支持力触觉反馈的涡卷弹簧模型，如图2所示的涡卷弹簧模型示意图；在交互过程中，虚拟柔性体产生实时变形力触觉信息，该力触觉信息与有色状态噪声信号η(t)相加，即采用支持增强力触觉反馈的涡卷弹簧模型，产生的力触觉信息作为反馈输出信号；所述涡卷弹簧模型的建模步骤如下：

步骤1，建立空间直角坐标系，确定柔性体表面上任意点处安放的涡卷弹簧模型，其过程如下：

步骤1-1，建立空间直角坐标系，

在给定虚拟外力F作用下，在柔性体表面任意点O₁处安放一涡卷弹簧，在距离任意点O₁为r处，取涡卷弹簧第一层的一个外侧点P₀，即r为涡卷弹簧第一层的半径，以任意点O₁为原点，射线O₁P₀所在方向为X轴正方向，建立XYZ空间坐标系；

步骤1-2，按层设置涡卷弹簧，涡卷弹簧任意两层的间距均为d，涡卷弹簧每层弹簧丝的直径均为h；

所述涡卷弹簧的总层数包含两种，分别是奇数和偶数，令n＝1,2,3……S,S为自然数；则所述涡卷弹簧具体描述如下：

A.当该涡卷弹簧总层数为奇数时，设层数为2n+1，在Z轴上等间距取2n+1个点，分别为O₁,O₂,…,O_i,…,O_2n,O_2n+1，间距为d；以O_i为圆心，设置涡卷弹簧第i层，且该涡卷弹簧第i层的半径为R_i，其表达式如下：

$R_i=\left\{\begin{array}{cc}r+(i-1)\mathrm{\ρ}& i\≤n+1\\ r+(2n-i+1)\mathrm{\ρ}& i>n+1\end{array}\right.$

也即前n+1层涡卷弹簧的半径R_i组成以r为首项，以ρ为公差的等差数列，ρ为实数；

B.当该涡卷弹簧总层数为偶数时，设层数为2n，则在Z轴上等间距取2n个点，分别为O₁,O₂,…,O_i,…,O_2n-1,O_2n，间距为d；以O_i为圆心，设置涡卷弹簧第i层，且该涡卷弹簧第i层的半径为R_i，其表达式如下：

$R_i=\left\{\begin{array}{cc}r+(i-1)\mathrm{\ρ}& i\≤n\\ r+(2n-i)\mathrm{\ρ}& i>n\end{array}\right.$

也即前n层涡卷弹簧的半径R_i组成以r为首项，以ρ为公差的等差数列；

涡卷弹簧每层的截面惯性矩均为涡卷弹簧各层的挠度值相同，挠度值为X′，弹性模量E取决于柔性体的材质，且柔性体的材质均相同；

步骤2，确定涡卷弹簧任一层所消耗的外力，

假定给定虚拟外力F的作用线与涡卷弹簧中心线一致，且在给定虚拟外力F作用下，如果涡卷弹簧共有前C层产生变形，则第C层称为变形截止层；被拉伸或压缩时，变形截止层产生的变形量不大于挠度值X′，涡卷弹簧前C-1层中，任意一层的变形量均等于挠度值X′；

步骤2-1，确定涡卷弹簧前C-1层消耗的外力，

涡卷弹簧第一层消耗的外力F₁为：

$F_1=\left\{\begin{array}{cc}{X}^{\′}\·K_1& C>1\\ F& C=1\end{array}\right.$

其中，K₁表示涡卷弹簧第一层的刚度系数；

除第一层和变形截止层第C层外，涡卷弹簧第i层消耗的外力F_i为：

$F_i=\frac{{\mathrm{Gh}}^4}{8{R_i}^3}(d-h)$

其中，G表示该涡卷弹簧的切变模量，2≤i≤C-1；

步骤2-2，确定涡卷弹簧第一层和变形截止层第C层的刚度系数，分别表示为K₁和K_C；其表达式分别如下：

$K_1=\frac{\mathrm{EI}}{l_1}$

$K_C=\frac{\mathrm{EI}}{l_C}$

其中，l₁、l_C分别为该涡卷弹簧第一层和变形截止层的工作长度，其计算式分别如下：

l₁＝2πR₁

$l_C=2\mathrm{\π}\underset{i=1}{\overset{C}{\mathrm{\Σ}}}{R}_i$

步骤3，涡卷弹簧变形截止层产生拉伸或压缩的变形量为X_C，其计算公式为：

$X_C=\frac{F-\underset{i=1}{\overset{C-1}{\mathrm{\Σ}}}{F}_i}{K_C}$

步骤4，将涡卷弹簧变形量相加，所得变形量之和X等效为虚拟柔性体表面的变形：

步骤5，所述有色状态噪声信号η(t)，其表达式为：

η(t)＝D(q^-1)C(q^-1)ξ(t)

式中，

D(q^-1)、C(q^-1)分别表示单位滞后算子q^-1的多项式，其表达式分别为：

D(q^-1)＝1+D₁q^-1+D₂q^-2+D₃q^-3

C(q^-1)＝1+C₁q^-1+C₂q^-2

其中，D₁、D₂、D₃、C₁、C₂皆是常数；

ξ(t)是零均值白噪声，t为时间。

步骤2中，所述常数D₁、D₂、D₃、C₁、C₂的取值分别为：

D₁＝-1.5,D₂＝0.7,D₃＝0.1,C₁＝0.5,C₂＝0.2。

所述涡卷弹簧产生拉伸或压缩变形总计需要的时延时间不大于1ms，也即满足刷新频率不低于1000Hz的要求。

所有涡卷弹簧层数中，每层产生拉伸或压缩变形需要的时延时间构成等比数列，即满足：

t_i＝p^i-1t₁，

其中，t_i表示涡卷弹簧第i层产生拉伸或压缩变形需要的时延时间，p为等比数列的公比，t₁为涡卷弹簧第一层产生拉伸或压缩变形需要的时延时间，1≤i≤C。

在所述涡卷弹簧上任一层产生拉伸或压缩变形达到给定的挠度值后，其下一层才开始产生拉伸或压缩变形。

力触觉交互设备采用PHANTOM OMNI力触觉交互设备，所述PHANTOMOMNI力触觉交互设备通过APS开关电源供电，通过APS开关电源之后，将供电电源220V的交流电转化为18V，2.22A所需的电压，显示卡、1394卡、硬盘均通过总线实现双向传输，显示器通过VGA连接线与显示卡连接。

下面以虚拟手和虚拟肝脏模型为例，列举本发明技术方案的具体实施方式。

构建虚拟手和虚拟肝脏模型，实现虚拟场景的初始化。

本实例中所有虚拟手和虚拟肝脏模型都直接采用从3DS MAX2013软件中导出的OBJ格式，以1795个质点，3587个三角网格构成的虚拟手和3845个质点，7685个三角网格构成的虚拟肝脏模型为例来进行变形仿真，实验过程中模型获取和修改非常方便；操作系统为Windows2000，以3DS MAX2013、OpenGL图形库为基础，在Microsoft Visual C++2012软件开发平台上进行仿真。

当检测到虚拟手碰撞到虚拟肝脏表面上任何一点时，在给定虚拟外力F＝7×10^-4N作用下，虚拟手与虚拟肝脏交互的局部区域内部，填充支持力触觉反馈的涡卷弹簧模型；在交互过程中，虚拟肝脏产生实时变形力触觉信息，该力触觉信息与有色状态噪声信号相加，作为反馈输出信号，如图3所示的柔性体变形仿真流程图；所述涡卷弹簧模型的建模步骤为：

在给定虚拟外力F＝7×10^-4N作用下，在柔性体表面任意点O₁处安放一涡卷弹簧，在距离任意点O₁为r＝1×10^-3m处，取涡卷弹簧第一层的一个外侧点P₀，即r为第一层涡卷弹簧的半径，以任意点O₁为原点，射线O₁P₀所在方向为X轴正方向，建立XYZ空间坐标系；

按层设置涡卷弹簧，分别如图2所示的涡卷弹簧模型示意图和如图4所示的涡卷弹簧模型构建方法流程图；涡卷弹簧任意两层的间距均为d＝1×10^-3m，涡卷弹簧每层弹簧丝的直径均为h＝5×10^-4m；所述涡卷弹簧的总层数包含两种，分别是奇数和偶数，也即前n+1层或者前n层涡卷弹簧的半径R_i组成以r＝1×10^-3m为首项，以ρ＝0.2×10^-3为公差的等差数列；

以R₁＝r＝1×10^-3m为半径，形成该涡卷弹簧第一层；

以R₂＝r+ρ＝1×10^-3+0.2×10^-3＝1.2×10^-3m为半径，形成该涡卷弹簧第二层；

以R₃＝r+2ρ＝1×10^-3+2×0.2×10^-3＝1.4×10^-3m为半径，形成该涡卷弹簧第三层；

依此类推，形成该涡卷弹簧的任意层；

计算的中间过程、最后的数据均按四舍五入法保留小数点后3位。

涡卷弹簧每层的截面惯性矩均为涡卷弹簧各层的挠度值相同，假定挠度值为X′＝1.295×10^-4m，弹性模量E＝3.09×10¹¹Pa，取决于柔性体的材质，且柔性体的材质均相同；

假定给定虚拟外力F＝7×10^-4N的作用线与涡卷弹簧中心线一致，且在给定虚拟外力F＝7×10^-4N作用下，如果涡卷弹簧共有前C层产生变形，则第C层称为变形截止层；被拉伸或压缩时，变形截止层产生的变形量不大于挠度值X′＝1.295×10^-4m，在涡卷弹簧前C-1层中，任意一层被拉伸或压缩时，产生的变形量均等于挠度值X′＝1.295×10^-4m；

所有涡卷弹簧层数中，每层产生拉伸或压缩变形需要的时延时间构成以t₁＝10^-5s为首项、以p＝1.2为公比的等比数列；

如图5所示的涡卷弹簧模型变形、变形量与时延时间关系示意图，假定力触觉再现刷新频率为1200Hz，则力触觉再现刷新频率的倒数

涡卷弹簧的切变模量G＝8×10⁴Pa；

根据涡卷弹簧特性设定：

涡卷弹簧第一层的工作长度为：l₁＝2πR₁，

确定涡卷弹簧第一层的刚度系数K₁为：

$K_1=\frac{\mathrm{EI}}{l_1}=\frac{3.09\×{10}^{11}\×4.909\×{10}^{-14}}{2\mathrm{\π}{\×10}^{-3}}=2.414N/m$

假定在给定虚拟外力F作用下，能使涡卷弹簧第一层产生拉伸或压缩变形达到给定的挠度值，则：涡卷弹簧第一层消耗的外力F₁为：

F₁＝X′·K₁＝1.295×10^-4×2.414＝3.126×10^-4N<7×10^-4N

因此，涡卷弹簧第一层消耗的外力F₁<F，且涡卷弹簧第一层产生拉伸或者压缩变形总计需要的时延时间T₁＝t₁＝10^-5s<T，满足刷新频率的要求；

只有当涡卷弹簧第一层产生拉伸或压缩变形达到给定的挠度值后，涡卷弹簧第二层才开始产生拉伸或压缩变形。

若涡卷弹簧第二层，在给定虚拟外力F＝7×10^-4N作用下，被拉伸或压缩时，产生的变形量等于挠度值X′＝1.295×10^-4m时，则涡卷弹簧第二层消耗的外力F₂为：

$F_2=\frac{G{h}^4}{8R_2^3}(d-h)=\frac{8\&times;{10}^4{(5\&times;{10}^{-4})}^4}{8\&times;{(1.2\&times;{10}^{-3})}^3}\&times;(1\&times;{10}^{-3}-5\&times;{10}^{-4})=1.808\&times;{10}^{-4}N$

涡卷弹簧前二层共消耗的外力之和为：

F₁+F₂＝3.126×10^-4+1.808×10^-4＝4.934×10^-4N<F＝7×10^-4N，

涡卷弹簧前二层产生拉伸或压缩变形总计需要的时延时间为：

T₂＝t₁+t₂＝(1+q)t₁＝(1+1.2)×10^-5＝2.2×10^-5s<T，这里T为力触觉再现刷新频率的倒数，满足刷新频率的要求；

只有当涡卷弹簧第二层产生拉伸或压缩变形达到给定的挠度值后，涡卷弹簧第三层才开始产生拉伸或压缩变形。

若涡卷弹簧第三层，在给定虚拟外力F＝7×10^-4N作用下，被拉伸或压缩时，产生的变形量等于挠度值X′＝1.295×10^-4m时，则涡卷弹簧第三层消耗的外力F₃为：

$F_3=\frac{G{h}^4}{8R_3^3}(d-h)=\frac{8\&times;{10}^4{(5\&times;{10}^{-4})}^4}{8\&times;{(1.4\&times;{10}^{-3})}^3}\&times;({10}^{-3}-5\&times;{10}^{-4})=1.139\&times;{10}^{-4}N,$

涡卷弹簧前三层共消耗的外力之和为：

F₁+F₂+F₃＝3.126×10^-4+1.808×10^-4+1.139×10^-4＝6.073×10^-4N<F＝7×10^-4N，

涡卷弹簧前三层产生拉伸或压缩变形总计需要的时延时间为：

T₃＝t₁+t₂+t₃＝(1+q+q²)t₁＝(1+1.2+1.2²)×10^-5＝3.64×10^-5s<T，这里T为力触觉再现刷新频率的倒数，满足刷新频率的要求；

只有当涡卷弹簧第三层产生拉伸或压缩变形达到给定的挠度值后，涡卷弹簧第四层才开始产生拉伸或压缩变形。

若涡卷弹簧第四层，在给定虚拟外力F＝7×10^-4N作用下，被拉伸或压缩时，产生的变形量等于挠度值X′＝1.295×10^-4m时，则涡卷弹簧第四层消耗的外力F₄为：

$F_4=\frac{G{h}^4}{8R_4^3}(d-h)=\frac{8\&times;{10}^4{(5\&times;{10}^{-4})}^4}{8\&times;{(1.6\&times;{10}^{-3})}^3}\&times;({10}^{-3}-5\&times;{10}^{-4})=7.629\&times;{10}^{-5}N,$

涡卷弹簧前四层共消耗的外力之和为：

F₁+F₂+F₃+F₄＝3.126×10^-4+1.808×10^-4+1.139×10^-4+7.629×10^-5＝6.836×10^-4N<F＝7×10^-4N，

涡卷弹簧前四层产生拉伸或压缩变形总计需要的时延时间为：

T₄＝t₁+t₂+t₃+t₄＝(1+q+q²+q³)t₁＝5.368×10^-5s<T，

这里T为力触觉再现刷新频率的倒数，满足刷新频率的要求；

只有当涡卷弹簧第四层产生拉伸或压缩变形达到给定的挠度值后，涡卷弹簧第五层才开始产生拉伸或压缩变形。

若涡卷弹簧第五层，在给定虚拟外力F＝7×10^-4N作用下，被拉伸或压缩时，产生的变形量等于挠度值X′＝1.295×10^-4m时，则涡卷弹簧第五层消耗的外力F₅为：

$F_5=\frac{G{h}^4}{8R_5^3}(d-h)=\frac{8\&times;{10}^4{(5\&times;{10}^{-4})}^4}{8\&times;{(1.8\&times;{10}^{-3})}^3}\&times;({10}^{-3}-5\&times;{10}^{-4})=5.358\&times;{10}^{-5}N,$

涡卷弹簧前五层共消耗的外力之和为：

F₁+F₂+F₃+F₄+F₅

＝3.126×10^-4+1.808×10^-4+1.139×10^-4+7.629×10^-5+5.358×10^-5

＝7.372×10^-4N>F＝7×10^-4N

因此，涡卷弹簧前五层共消耗的外力之和不小于给定虚拟外力F，故涡卷弹簧第五层为变形截止层，不需要再判断是否满足刷新频率的要求。

涡卷弹簧第五层即变形截止层实际消耗的外力F₅为：

F₅＝F-(F₁+F₂+F₃+F₄)

＝7×10^-4-6.836×10^-4，

＝0.164×10^-4N

涡卷弹簧第五层的工作长度为： $l_5=l_C=2\mathrm{\&pi;}\underset{i=1}{\overset{5}{\mathrm{\&Sigma;}}}{R}_5=2\mathrm{\&pi;}(R_1+R_2+R_3+R_4+R_5),$

涡卷弹簧第五层的刚度系数为：

$\begin{array}{c}{K}_5=\frac{\mathrm{EI}}{l_5}=\frac{3.09\&times;{10}^{11}\&times;4.909\&times;{10}^{-14}}{2\mathrm{\&pi;}({10}^{-3}+1.2\&times;{10}^{-3}+1.4\&times;{10}^{-3}+1.6\&times;{10}^{-3}+1.8\&times;{10}^{-3})}\\ =0.345N/m\end{array}$

涡卷弹簧第五层即变形截止层产生拉伸或压缩的变形量为X₅，其为：

$X_5=\frac{F-(F_1+F_2+F_3+F_4)}{K_5}=\frac{0.164\&times;{10}^{-4}}{0.345}=0.475\&times;{10}^{-4}m$

将涡卷弹簧变形量相加，所得变形量之和X等效为虚拟肝脏表面的变形，

$X=\underset{i=1}{\overset{4}{\mathrm{\&Sigma;}}}{X}^{\&prime;}+X_5=4X^{\&prime;}+X_5=4\&times;1.295\&times;{10}^{-4}+0.475\&times;{10}^{-4}=5.655\&times;{10}^{-4}m.$

所述有色状态噪声信号η(t)，其表达式为：

η(t)＝D(q^-1)/C(q^-1)ξ(t)

式中，

D(q^-1)、C(q^-1)分别表示单位滞后算子q^-1的多项式，其表达式分别为：

D(q^-1)＝1-1.5q^-1+0.7q^-2+0.1q^-3

C(q^-1)＝1+0.5q^-1+0.2q^-2

ξ(t)是零均值白噪声，t为时间；

注意：在计算在给定虚拟外力F作用下，虚拟柔性体产生实时变形仿真的过程中，若r、d、h、ρ等这些参数选取的过大，则支持增强力触觉反馈的涡卷弹簧模型中，涡卷弹簧变形截止层数值就少，计算量小，实时性好，但变形仿真效果不佳；若r、d、h、ρ等这些参数选取的过小，则支持增强力触觉反馈的涡卷弹簧模型中，涡卷弹簧变形截止层数值就大，计算量大，实时性不佳，但变形仿真效果较好；另外在设置t₁和t_i之间的关系时，要应考虑到程序运行时计算机本身的硬件配置，故在调试整个程序的过程中，要折中选择这些参数，不断反复调试，从而使变形效果更加逼真。

为验证本发明的实施效果，操作者通过PHANTOM OMNI手控器端部的手柄来触摸、感知和控制虚拟手对虚拟肝脏进行旋转的变形仿真，并将交互过程中产生的力触觉信息实时反馈给操作者。实验结果表明：该系统模拟效果逼真，无在力触觉感知方面，能使操作者产生身临其境的感觉，增强了交互的感知准确性和真实性。

力触觉交互设备3采用PHANTOM OMNI力触觉交互设备3，所述PHANTOM OMNI力触觉交互设备3通过APS开关电源4供电，通过APS开关电源之后，将供电电源220V的交流电转化为18V，2.22A所需的电压，显示卡21、1394卡23、硬盘24均通过总线22实现双向传输，显示器1通过VGA连接线与显示卡21连接。

Claims (5)

1.一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征在于，包括：计算机主机(2)和力触觉交互设备(3)，在计算机主机(2)上连接有显示器(1)，所述计算机主机(2)至少包括硬盘(24)和1394卡(23)，其特征在于，所述硬盘(24)至少包括：用于检测虚拟代理位置的位置检测模块(241)、用于检测虚拟代理与虚拟柔性体碰撞的碰撞检测模块(242)、力触觉计算变形模块(243)、根据力触觉计算变形模块(243)输出的最新力触觉变形信息进行图形显示的图形刷新模块(244)及力触觉信息反馈模块(245)，所述硬盘(24)通过1394卡(23)及1394连接线与力触觉交互设备(3)连接，用于将力触觉计算变形模块(243)产生的柔性体表面的变形信息传输至力触觉交互设备(3)，

所述力触觉计算变形模块(243)用于在检测到虚拟代理碰撞到虚拟柔性体表面上任何一点时，在给定虚拟外力F作用作用下，虚拟代理与虚拟柔性体交互的局部区域内部，填充支持力触觉反馈的涡卷弹簧模型；在交互过程中，虚拟柔性体产生实时变形力触觉信息，该力触觉信息与有色状态噪声信号η(t)相加，作为反馈输出信号；所述涡卷弹簧模型的建模步骤如下：

步骤1，建立空间直角坐标系，确定柔性体表面上任意点处安放的涡卷弹簧模型，其过程如下：

步骤1-1，建立空间直角坐标系，

在给定虚拟外力F作用下，在柔性体表面任意点O₁处安放一涡卷弹簧，在距离任意点O₁为r处，取涡卷弹簧第一层的一个外侧点P₀，即r为涡卷弹簧第一层的半径，以任意点O₁为原点，射线O₁P₀所在方向为X轴正方向，建立XYZ空间坐标系；

步骤1-2，按层设置涡卷弹簧，涡卷弹簧任意两层的间距均为d，涡卷弹簧每层弹簧丝的直径均为h；

所述涡卷弹簧的总层数包含两种，分别是奇数和偶数，令n＝1,2,3……S,S为自然数；则所述涡卷弹簧具体描述如下：

A.当该涡卷弹簧总层数为奇数时，设层数为2n+1，在Z轴上等间距取2n+1个点，分别为O₁,O₂,…,O_i，…,O_2n,O_2n+1，间距为d；以O_i为圆心，设置涡卷弹簧第i层，且该涡卷弹簧第i层的半径为R_i，其表达式如下：

$R_i=\left\{\begin{array}{cc}r+(i-1)\mathrm{\&rho;}& i\&le;n+1\\ r+(2n-i+1)\mathrm{\&rho;}& i>n+1\end{array}\right.$

B.当该涡卷弹簧总层数为偶数时，设层数为2n，则在Z轴上等间距取2n个点，分别为O₁,O₂,…,O_i,…,O_2n-1,O_2n，间距为d；以O_i为圆心，设置涡卷弹簧第i层，且该涡卷弹簧第i层的半径为R_i，其表达式如下：

$R_i=\left\{\begin{array}{cc}r+(i-1)\mathrm{\&rho;}& i\&le;n\\ r+(2n-i)\mathrm{\&rho;}& i>n\end{array}\right.$

其中，ρ是等差数列的公差，为实数；

涡卷弹簧每层的截面惯性矩均为涡卷弹簧各层的挠度值相同，挠度值为X′，弹性模量E取决于柔性体的材质，且柔性体的材质均相同；

步骤2，确定涡卷弹簧任一层所消耗的外力，

假定给定虚拟外力F的作用线与涡卷弹簧中心线一致，且在给定虚拟外力F作用下，如果涡卷弹簧共有前C层产生变形，则第C层称为变形截止层；被拉伸或压缩时，变形截止层产生的变形量不大于挠度值X′，涡卷弹簧前C-1层中，任意一层的变形量均等于挠度值X′；

步骤2-1，确定涡卷弹簧前C-1层消耗的外力，

涡卷弹簧第一层消耗的外力F₁为：

$F_1=\left\{\begin{array}{cc}{X}^{\&prime;}\&CenterDot;K_1& C>1\\ F& C=1\end{array}\right.$

其中，K₁表示涡卷弹簧第一层的刚度系数；

除第一层和变形截止层第C层外，涡卷弹簧第i层消耗的外力F_i为：

$F_i=\frac{{\mathrm{Gh}}^4}{8{R_i}^3}(d-h)$

其中，G表示该涡卷弹簧的切变模量，2≤i≤C-1；

步骤2-2，确定涡卷弹簧第一层和变形截止层第C层的刚度系数，分别表示为K₁和K_C；其表达式分别如下：

$K_1=\frac{\mathrm{EI}}{l_1}$

$K_C=\frac{\mathrm{EI}}{l_C}$

其中，l₁、l_C分别为该涡卷弹簧第一层和变形截止层的工作长度，其计算式分别如下：

l₁＝2πR₁

$l_C=2\mathrm{\&pi;}\underset{i=1}{\overset{C}{\mathrm{\&Sigma;}}}{R}_i$

步骤3，涡卷弹簧变形截止层产生拉伸或压缩的变形量为X_C，其计算公式为：

$X_C=\frac{F-\underset{i=1}{\overset{C-1}{\mathrm{\&Sigma;}}}{F}_i}{K_C}$

步骤4，将涡卷弹簧变形量相加，所得变形量之和X等效为虚拟柔性体表面的变形：

步骤5，所述有色状态噪声信号η(t)，其表达式为：

η(t)＝D(q^-1)/C(q^-1)ξ(t)

式中，

D(q^-1)、C(q^-1)分别表示单位滞后算子q^-1的多项式，其表达式分别为：

D(q^-1)＝1+D₁q^-1+D₂q^-2+D₃q^-3

C(q^-1)＝1+C₁q^-1+C₂q^-2

其中，D₁、D₂、D₃、C₁、C₂皆是常数；

ξ(t)是零均值白噪声，t为时间。

2.根据权利要求1所述的一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征在于，步骤2中，所述常数D₁、D₂、D₃、C₁、C₂的取值分别为：

D₁＝-1.5,D₂＝0.7,D₃＝0.1,C₁＝0.5,C₂＝0.2。

3.根据权利要求1所述的一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征在于，所述涡卷弹簧产生拉伸或压缩变形总计需要的时延时间不大于1ms。

4.根据权利要求1所述的一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征在于，所有涡卷弹簧层数中，每层产生拉伸或压缩变形需要的时延时间构成等比数列，即满足：

t_i＝p^i-1t₁，

其中，t_i表示涡卷弹簧第i层产生拉伸或压缩变形需要的时延时间，p为等比数列的公比，t₁为涡卷弹簧第一层产生拉伸或压缩变形需要的时延时间，1≤i≤C。

5.根据权利要求1或2或3或4所述的一种支持增强力触觉反馈的涡卷弹簧仿真系统，其特征在于，在所述涡卷弹簧上任一层产生拉伸或压缩变形达到给定的挠度值后，其下一层才开始产生拉伸或压缩变形。