CN105302972A - 一种基于元球模型的软组织形变方法 - Google Patents

Abstract

本发明提出了一种基于元球模型的软组织形变方法，包括四个步骤：元球模型拓扑结构的构建阶段，使用Bradshow？Gareth球树生成算法，根据原始网格模型生成维诺图，进而生成元球模型，通过设置阈值的方法构建其拓扑结构；元球模型形变计算阶段，使用Position？Based？Dynamic(PBD)算法结合Laplacian坐标约束来模拟软组织体模型的形变过程；软组织模型的蒙皮阶段，对体模型中的元球分别建立其距离场函数，将体模型和表皮网格模型建立映射关系，实现蒙皮过程；真实感绘制及实时触觉渲染阶段，根据软组织物理形变进行真实感绘制，并基于Geomagic？Touch力反馈设备进行实时触觉渲染。本发明可真实模拟虚拟手术中软组织的形变过程，具有物理真实感强，实时性好的特点。

Description

一种基于元球模型的软组织形变方法

技术领域

本发明涉及一种基于元球模型的软组织形变方法。

背景技术

随着虚拟现实(VirtualReality，VR)技术在医学领域应用的不断拓展和深入，虚拟手术(VirtualSurgery)及VR手术模拟器越来越受到人们的关注。虚拟手术技术具有手术环境及器械响应可控，可重复演练等多项优点，成为未来外科培训的发展趋势。医生借助VR手术模拟器进行手术训练及规划，可降低手术训练以及治疗的成本和风险；减少医生培训教学中对动物和尸体的依赖。

在虚拟手术中，对人体组织器官精确的物理建模是其关键技术。一些经典的力学模型，例如有限元和质子弹簧模型等，已被应用到各类手术中，但是这些算法有些缺乏实时运算能力，有的缺乏精确的表达能力，因而，研发更为有效的人体软组织物理模型就成为虚拟手术中的一大技术难题。

在这样的背景下，考虑到人体一些主要器官，例如：胃、肝脏、胆囊等通常具有连续、光滑、圆润的表面，而元球模型(Blobbymodel、Metaballs)在呈现和表达这样的形体时具有天然的优势。因此，本发明提出了一种基于元球模型的软组织形变方法。该算法包括软组织体模型形变阶段和蒙皮阶段，并在其中形变模拟中保证了其物理真实性和迭代实时性。

发明内容

本发明解决的技术问题是：克服了传统网格模型的真实性问题和传统物理方法的迭代周期性问题，提供了一种基于元球模型的软组织形变方法。该方法满足了软组织渲染的真实性和形变迭代的实时性。

本发明采用的技术方案为：一种基于元球模型的软组织形变方法，包括以下四个步骤：

步骤(1)、元球模型拓扑结构的构建：通过BradshowGareth球树生成算法生成软组织的元球模型，通过设置阈值的方法构建软组织体模型的拓扑结构。

步骤(2)、元球模型的形变计算：使用PositionBasedDynamic(PBD，位置动力学)算法与Laplacian坐标约束相结合的方法对软组织元球模型进行形变过程模拟。

步骤(3)、软组织体模型的蒙皮算法：对软组织体模型中的元球分别建立距离场函数，将体模型中的元球和表皮模型建立映射关系，在每次形变迭代中，表皮模型根据体模型的位置和映射关系更新其位置，完成蒙皮模拟。

步骤(4)、真实感绘制及实时触觉渲染：根据步骤(1)、(2)、(3)的算法，实现真实感绘制，并基于GeomagicTouch力反馈设备实现触觉渲染，并保证在此过程的实时性。

本发明的原理在于：

(1)为了实现软组织的连续、光滑、圆润的特性，本发明使用元球模型代替传统四面体模型作为物理模型来实现形变模拟。

(2)为了实现手术模拟中形变的物理真实性和实时性，又因为单纯PositionBasedDynamic(PBD)算法的局限性，本发明使用PBD算法和Laplacian坐标约束相结合的方法来进行体模型的形变模拟。

(3)为了实现软组织形变中的蒙皮模拟效果，本发明使用了生成元球模型距离场函数的方法，建立体模型与表皮网格模型的映射关系，并应用此算法来进行蒙皮的形变模拟。

(4)为了增加系统的物理真实性和触觉真实性，本发明使用了GeomagicTouch力反馈设备来进行手术器械的真实感绘制和触觉渲染。

本发明与现有技术相比的优点在于：

1、本发明采用的元球模型较传统的四面体模型对于软组织的连续、光滑、圆润的特性能够更好的表达出来。

2、本发明提出的将PBD算法和Laplacian坐标约束相结合的方法来模拟形变过程,克服了单纯PBD算法的局限性，并且物理真实性不亚于传统物理方法，又能很好地满足系统实时性要求。

3、本发明提出的将元球建立距离场，进而建立体模型和表皮模型映射关系的蒙皮算法，以及使用GeomagicTouch力反馈设备实现触觉渲染，大大增加了形变模拟过程中的真实性。

附图说明

图1：基于元球模型的软组织形变方法的处理流程图；

图2(a)：用于模拟的肝脏原始网格模型；

图2(b)：使用Bradshow算法生成的肝脏元球模型图；

图2(c)：肝脏的表面网格模型和元球模型图；

图3(a)：肝脏半透明表皮模型和元球模型图；

图3(b)：在图3-(a)的基础上添加元球模型的拓扑结构图；

图3(c)：肝脏的半透明表皮模型和元球模型的拓扑结构图；

图4：PBD算法中的拉伸约束图示；

图5(a)：肝脏形变模拟的实验结果图，渲染半透明元球模性和其拓扑结构；

图5(b)：肝脏形变模拟的纹理表面的实验结果图；

图6(a)：胆囊形变模拟的实验结果图；

图6(b)：小肠形变模拟的实验结果图；

图6(c)：胃形变模拟的实验结果图；

图7(a)：基于有限元算法对肝脏形变模拟的实验结果图；

图7(b)：基于弹簧质子模型对肝脏形变模拟的实验结果图；

图7(c)：基于传统PBD算法对肝脏形变模拟的实验结果图；

图7(d)：基于发明对肝脏形变模拟的实验结果图。

具体实施方式

图1给出了基于元球模型的软组织形变方法的处理流程，下面结合其他附图及具体实施方式进一步说明本发明。

本发明提供一种基于元球模型的软组织形变方法，主要步骤如下：

1、元球模型拓扑结构的构建及其优化法

首先使用BradshowGareth球树生成算法，根据原始三角形网格模型(图2a)生成维诺图，根据维诺图，生成需要的元球模型，如图2b所示；图2c显示了网格模型和元球模型在一起的比较结果。

在已得到元球模型的基础上，对其构建拓扑结构，以下是其构建算法。

对于元球i，中心是c_i，与其拓扑相连的点函数定义为：

$Topo\left(i\right)=\left\{\begin{array}{c}L\left(c_i\right),num\left(i\right)\≤N\\ L^{\′}\left(c_i\right),num\left(i\right)>N\end{array}\right.---\left(1\right)$

其中，N是规定阈值，num(i)是与元球重叠的元球数量。L函数指与元球i重叠的元球全部定义为拓扑连接，L′是指距离i最近的N个元球定义为其拓扑连接。由此算法构建的拓扑结构如图3所示：

其中，图3a是表皮模型和元球模型；图3b是在图3a的基础上添加拓扑结构，图3c是只有表皮模型的情况下的拓扑结构。

2、元球模型形变计算

元球模型的形变模拟使用位置动力学：PositionBasedDynamic(PBD)算法与Laplacian坐标约束相结合的方法。

首先是PBD算法。PBD(PositionBasedDynamic)是一种直接由位置计算形变的启发式算法。在进行形变迭代时，不需要计算单元相互之间的加速度，只需要根据基于位置的约束函数将每个顶点投射到合适的位置，而从当前位置到最终位置的该变量可以用约束函数的梯度表示。

由于PBD的BendingConstraints和VolumeConstraints的局限性，本文方法只能采用其中的StretchConstraints(拉伸约束)。

图4给出了StretchConstraints的示例。距离约束函数为：

C_stretch(p₁,p₂)＝|p₁-p₂|-d(3)

其中，d为顶点p₁和p₂之间最初的原始距离，有式子可以得到， ${\▿}_{p_2}C(p_1,p_2)=-\frac{p_1-p_2}{\left|p_1-p_2\right|},$ 最终可以得到：

${\mathrm{\Δp}}_1=-\frac{w_1}{w_1+w_2}\left(\right|p_1-p_2|-d)\frac{p_1-p_2}{\left|p_1-p_2\right|}---\left(4\right)$

${\mathrm{\Δp}}_2=\frac{w_2}{w_1+w_2}\left(\right|p_1-p_2|-d)\frac{p_1-p_2}{|p_1-p_2|}---\left(5\right)$

其中，w₁与w₂为两个顶点的权重。

在形变过程中，只有PBD算法对其进行约束是不够的，因为其中的StretchConstraint只能约束到2维，因此，还需要添加一个3维约束，这里我们就引进了LaplacianCoordinatesConstraint(Laplacian坐标约束)。它的算法如下。

对于任一元球，设为m，中心设为c_m，假设其拓扑结构上相邻元球有n个，而这些元球中心设为c_i，那么元球i的拓扑中心为c_center，即：

$c_{center}=\underset{i}{\overset{n}{\Σ}}\frac{c_i}{n}---\left(6\right)$

然后，对元球模型的形变模拟进行预处理，对于元球m，定义其Laplacian坐标为：

$L_m=c_m-c_{center}=c_m-\underset{i}{\overset{n}{\Σ}}\frac{c_i}{n}---\left(7\right)$

根据元球模型的Laplacian坐标约束，即L_m，它是一个固定矢量，在每次形变迭代中，元球中心位置更新为：

c′_m＝L_m+c′_center(8)

式中c′_m和c′_center为c_m和c_center更新之后的新的位置。

3、软组织模型的蒙皮算法

为了增加软组织模拟中的真实性，需要添加表皮模型，而表皮模型的形变是由蒙皮算法决定的。如下是本发明所使用的蒙皮算法。

首先，需要对元球模型进行预处理。即对于体模型中任一元球，建立其距离场函数：

$f(d,r)=c\·{\exp }^{-\frac{{(d-r)}^2}{r^2}}---\left(9\right)$

在上式中，r是其半径，d是到球心的距离，c是常量参数，它的值由实验决定。

然后，对于表皮模型中的任一点，设为v，假设满足条件的元球有n个，则这n个元球构成一个虚拟拓扑结构，假设此拓扑的中心为p_center，定义其式为：

$p_{center}=\underset{i}{\overset{n}{\Σ}}\frac{f(d_{v,c_i},r_i)}{{\mathrm{\Σ}}_i^{n}f(d_{v,c_i},r_i)}{c}_i---\left(10\right)$

在此式中，c_i为元球i的球心。

预处理结束后，开始形变模拟。在形变迭代过程中，首先建立表皮顶点v与其对应的虚拟拓扑结构中心p_center的关系式：

disp＝v-p_center(11)

然后，在每次迭代计算过程中，表皮点v的位置更新后的位置如式(12)所示：

v′＝p′_center+disp(12)

其中v′和p′_center为迭代后更新的v和p_center。

4、真实感绘制和触觉渲染

根据上述算法，完成对软组织的真实感绘制。视觉渲染使用OpenGLAPI，触觉渲染使用GeomagicTough力反馈设备(图5)及其提供的OpenHapticAPI。

实现本次算法的实验设备为NVIDIAGeForceGTX460,Intel(R)Core(TM)2QuadCPU(2.66GHz,4cores),and4GRAM运行在Windows764位系统上。

图(5)显示了本次算法中的肝脏部分的实验结果，图5a是肝脏手术模拟中元球模型形变的效果，而图5b表示肝脏表面网格模型的形变效果。

图(6)是使用本实验中除了肝脏外的其他软组织器官的实验结果。其中，图(a)(b)(c)分别代表着胆囊、小肠和胃。然后，第一列是原始状态中的表皮模型和元球体模型；第二列显示了表皮模型和拓扑结构；第三列和第四列分别表示对添加纹理后的软组织不同位置拉伸之后形变的实验结果。

表1不同软组织模型实验数据表

表1是本文方法的实验结果数据表。其中，四行代表着四种软组织器官；第一列是软组织器官的表皮模型中的顶点数，第二列是软组织器官体模型中的元球数，第三列是元球模型中元球之间构成的拓扑连接线数，第四列是生成元球模型的时间，最后一列是形变过程中每次迭代所需的时间。

图7给出了本文方法与传统物理方法的比较。其中(a)图是FEM方法，(b)图是弹簧质子模型的实验结果，(c)图是基于传统PBD算法的实验结果，(d)图是本次算法的实验结果。从中可以看出，本算法与FEM算法相较于弹簧质子和传统PBD算法有较好的实验效果，但是本算法的计算迭代周期比FEM算法短得多。

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

尽管上面对本发明说明性的具体实施方式进行了描述，以便于本技术领的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。

Claims (4)

1.一种基于元球模型的软组织形变方法，其特征在于包括以下四个步骤：

步骤(1)、元球模型拓扑结构的构建：通过BradshowGareth球树生成算法生成软组织的元球模型，通过设置阈值的方法构建软组织体模型的拓扑结构；

步骤(2)、元球模型的形变计算：使用PositionBasedDynamic(PBD),位置动力学算法,结合Laplacian坐标约束来对软组织元球模型进行形变过程模拟；

步骤(3)、软组织体模型的蒙皮算法：对软组织体模型中的元球分别建立距离场函数，将体模型中的元球和表皮模型建立映射关系，在每次形变迭代中，表皮模型根据体模型的位置和映射关系更新其位置，最终完成表面网格的变形；

步骤(4)、真实感绘制及实时触觉渲染：根据步骤(1)、(2)、(3)的算法，实现真实感绘制，并基于GeomagicTouch力反馈设备实现触觉渲染，并保证此过程的实时性。

2.根据权利要求1所述的基于元球模型的软组织形变方法，其特征在于：步骤(1)中所述的元球模型拓扑结构的构建中，在通过使用Bradshow算法生成软组织的元球模型的基础上，使用如下算法构建其拓扑结构：

对于元球i，中心是c_i，与其拓扑相连的点函数定义为：

$Topo\left(i\right)=\left\{\begin{array}{c}L\left(c_i\right),num\left(i\right)\≤N\\ L^{\′}\left(c_i\right),num\left(i\right)>N\end{array}\right.---\left(1\right)$

其中，N是规定阈值，num(i)是与元球i重叠的元球数量，L函数指与元球i重叠的元球全部定义为拓扑连接，L′是指距离i最近的N个元球定义为其拓扑连接。

3.根据权利要求1所述的基于元球模型的软组织形变方法，其特征在于：步骤(2)中所述的元球模型形变计算中，使用PositionBasedDynamic(PBD)算法与Laplacian坐标约束算法相结合的方法来对其进行形变模拟：

首先，使用PBD算法中的拉伸约束，拉伸约束函数如(2)式：

C_stretch(p₁,p₂)＝|p₁-p₂|-d(2)

其中，d为顶点p₁和p₂之间最初的原始距离， ${\▿}_{p_2}C(p_1,p_2)=-\frac{p_1-p_2}{|p_1-p_2|},$ 最终可以得到：

${\mathrm{\Δp}}_1=-\frac{w_1}{w_1+w_2}\left(\right|p_1-p_2|-d)\frac{p_1-p_2}{|p_1-p_2|}---\left(3\right)$

${\mathrm{\Δp}}_2=\frac{w_2}{w_1+w_2}\left(\right|p_1-p_2|-d)\frac{p_1-p_2}{|p_1-p_2|}---\left(4\right)$

其中，w₁与w₂为两个顶点的权重；

然后，Laplacian坐标约束算法如下：

对于任一元球，设为m，中心设为c_m，假设其拓扑结构上相邻元球有n个，而这些元球中心设为c_i，那么元球i的拓扑中心为c_center，即：

$c_{center}=\underset{i}{\overset{n}{\Σ}}\frac{c_i}{n}---\left(5\right)$

然后，对元球模型的形变模拟进行预处理，对于元球m，定义其Laplacian坐标：

$L_m=c_m-c_{center}=c_m-\underset{i}{\overset{n}{\Σ}}\frac{c_i}{n}---\left(6\right)$

由于L_m它是一个固定矢量，在每次形变迭代中，元球中心位置更新为：

c′_m＝L_m+c′_center(7)

式中c′_m和c′_center为c_m和c_center更新之后的新的位置。

4.根据权利要求1所述的基于元球模型的软组织形变方法，其特征在于：步骤(3)中所述的软组织体模型的蒙皮算法中，首先，需要对元球模型进行预处理：即对于体模型中任一元球，建立它的距离场函数，如(8)式所示：

$f(d,r)=c\·{\exp }^{-\frac{{(d-r)}^2}{r^2}}---\left(8\right)$

其中，r是其半径，d是表皮模型中的任一点到该球心的距离，c是常量参数；

然后，对于表皮模型中的任一点，设为v，假设满足条件的元球有n个，ri为第i个元球的半径，T是阈值,则这n个元球构成一个虚拟拓扑结构，假设c_i为元球i的球心，则定义拓扑的中心p_center为：

$p_{center}=\underset{i}{\overset{n}{\Σ}}\frac{f(d_{v,c_i},r_i)}{{\mathrm{\Σ}}_i^{n}f(d_{v,c_i},r_i)}{c}_i---\left(9\right)$

预处理结束后，开始形变模拟迭代，在形变迭代过程中，首先建立表皮顶点v与其对应的虚拟拓扑结构中心p_center的关系式：

disp＝v-p_center(10)

然后，在每次迭代计算过程中，表皮点v的位置更新后的位置如式(11)所示,其中v′和p′为迭代后更新的v和p：

v′＝p′+disp(11)。