本發(fā)明屬于生物力學(xué)仿真技術(shù)領(lǐng)域，尤其涉及外科手術(shù)軟組織變形的仿真方法及裝置。

背景技術(shù)：

近年來，盡管各種醫(yī)療手段都取得了巨大的進(jìn)步，但是對(duì)于情況較為嚴(yán)重的病患，手術(shù)依然是治愈率最高的一種。虛擬現(xiàn)實(shí)技術(shù)的手術(shù)仿真可以有效地縮短實(shí)習(xí)醫(yī)生進(jìn)行手術(shù)訓(xùn)練的周期，其主要應(yīng)用包括術(shù)前醫(yī)學(xué)影像處理、解剖結(jié)構(gòu)可視化、生物力學(xué)模擬等。

在虛擬現(xiàn)實(shí)技術(shù)的手術(shù)仿真領(lǐng)域，交互式生物力學(xué)仿真是一項(xiàng)關(guān)鍵技術(shù)，近年來發(fā)展起來的無網(wǎng)格方法闡述了不使用網(wǎng)格在傳統(tǒng)意義上的優(yōu)勢(shì)，這對(duì)于幾何形狀復(fù)雜的研究對(duì)象來說處理方式更為靈活。此外，由于拉格朗日框架對(duì)質(zhì)點(diǎn)分布的廣泛兼容性，無網(wǎng)格法可以較容易地控制有限應(yīng)力及大尺度形變，自五十年代有限元法出現(xiàn)之后，無網(wǎng)格方法已成為最具擴(kuò)展性的仿真方法。然而，無網(wǎng)格方法仍然存在著高精度幾何模型帶來的高計(jì)算復(fù)雜度的問題同時(shí)，在外科手術(shù)過程中軟組織變形的情況下，將不得不使用網(wǎng)格定義以控制可能出現(xiàn)的大尺度變形。

技術(shù)實(shí)現(xiàn)要素：

有鑒于此，本發(fā)明實(shí)施例提供了交互式生物力學(xué)的仿真方法及裝置，以解決無網(wǎng)格法無法處理外科手術(shù)仿真中軟組織大尺度變形的問題。

第一方面，提供了一種外科手術(shù)軟組織變形的仿真方法，包括：

在外科手術(shù)仿真過程中，對(duì)軟組織中的每個(gè)粒子i，計(jì)算其形變梯度其中，為粒子i的基本矩陣ai，所述j為所述粒子i的鄰域粒子，所述wij為權(quán)值，所述uij為初始時(shí)刻粒子i和鄰域粒子j的間距，所述xij為當(dāng)前時(shí)刻粒子i和鄰域粒子j的間距；

對(duì)所述形變梯度fi進(jìn)行奇異值分解得到奇異值矩陣并計(jì)算對(duì)角化后所述粒子i的第一應(yīng)力所述i為單位矩陣，所述λ和μ均為拉梅系數(shù)；

計(jì)算所述粒子i的每個(gè)鄰域粒子j對(duì)所述粒子i的作用力根據(jù)計(jì)算結(jié)果對(duì)所述軟組織進(jìn)行形變處理，所述vi為所述粒子i的體積。

第二方面，提供了一種外科手術(shù)軟組織變形的仿真裝置，包括：

第一計(jì)算單元，用于在外科手術(shù)仿真過程中，對(duì)軟組織中的每個(gè)粒子i，計(jì)算其形變梯度其中，為粒子i的基本矩陣ai，所述j為所述粒子i的鄰域粒子，所述wij為權(quán)值，所述uij為初始時(shí)刻粒子i和鄰域粒子j的間距，所述xij為當(dāng)前時(shí)刻粒子i和鄰域粒子j的間距；

第二計(jì)算單元，用于對(duì)所述形變梯度fi進(jìn)行奇異值分解得到奇異值矩陣并計(jì)算對(duì)角化后所述粒子i的第一應(yīng)力所述i為單位矩陣，所述λ和μ均為拉梅系數(shù)；

第三計(jì)算單元，用于計(jì)計(jì)算所述粒子i的每個(gè)鄰域粒子j對(duì)所述粒子i的作用力根據(jù)計(jì)算結(jié)果對(duì)所述軟組織進(jìn)行形變處理，所述vi為所述粒子i的體積。

本發(fā)明實(shí)施例基于點(diǎn)的方法來重現(xiàn)外科手術(shù)交互過程中人體器官的生物力學(xué)反應(yīng)，在軟組織的形變中，其每個(gè)粒子的位置變化僅僅取決于材料性質(zhì)，能夠達(dá)到很好的模型精確度。

附圖說明

為了更清楚地說明本發(fā)明實(shí)施例中的技術(shù)方案，下面將對(duì)實(shí)施例或現(xiàn)有技術(shù)描述中所需要使用的附圖作簡單地介紹，顯而易見地，下面描述中的附圖僅僅是本發(fā)明的一些實(shí)施例，對(duì)于本領(lǐng)域普通技術(shù)人員來講，在不付出創(chuàng)造性勞動(dòng)性的前提下，還可以根據(jù)這些附圖獲得其他的附圖。

圖1是本發(fā)明實(shí)施例提供的外科手術(shù)軟組織變形的仿真方法的實(shí)現(xiàn)流程圖；

圖2是本發(fā)明另一實(shí)施例提供的外科手術(shù)軟組織變形的仿真方法的實(shí)現(xiàn)流程圖；

圖3是本發(fā)明實(shí)施例提供的實(shí)驗(yàn)中采用的帶有紋理的表面渲染模型圖；

圖4是本發(fā)明實(shí)施例提供的實(shí)驗(yàn)中采用的連續(xù)介質(zhì)粒子化表達(dá)示意圖；

圖5是有限元法和本發(fā)明實(shí)施例方案的變形效果圖；

圖6是基于點(diǎn)的方法與本發(fā)明實(shí)施例方案在體積上的比較示意圖；

圖7是基于點(diǎn)的方法與本發(fā)明實(shí)施例方案在運(yùn)行時(shí)間上的比較示意圖；

圖8是基于點(diǎn)的方法與本發(fā)明實(shí)施例方案在視覺效果上的對(duì)比示意圖；

圖9和圖10分別是肝和胃在大變形下的仿真效果示意圖；

圖11是本發(fā)明實(shí)施例提供的外科手術(shù)軟組織變形的仿真裝置的結(jié)構(gòu)框圖。

具體實(shí)施方式

以下描述中，為了說明而不是為了限定，提出了諸如特定系統(tǒng)結(jié)構(gòu)、技術(shù)之類的具體細(xì)節(jié)，以便透切理解本發(fā)明實(shí)施例。然而，本領(lǐng)域的技術(shù)人員應(yīng)當(dāng)清楚，在沒有這些具體細(xì)節(jié)的其它實(shí)施例中也可以實(shí)現(xiàn)本發(fā)明。在其它情況中，省略對(duì)眾所周知的系統(tǒng)、裝置、電路以及方法的詳細(xì)說明，以免不必要的細(xì)節(jié)妨礙本發(fā)明的描述。

圖1示出了本發(fā)明實(shí)施例提供的外科手術(shù)軟組織變形的仿真方法的實(shí)現(xiàn)流程，在本方案中，采用了基于點(diǎn)的方法來重現(xiàn)外科手術(shù)交互過程中人體器官的生物力學(xué)反應(yīng)，其實(shí)現(xiàn)流程詳述如下：

在s101，在外科手術(shù)仿真過程中，對(duì)軟組織中的每個(gè)粒子i，計(jì)算其形變梯度其中，為粒子i的基本矩陣ai，所述j為所述粒子i的鄰域粒子，所述wij為權(quán)值，所述uij為初始時(shí)刻粒子i和鄰域粒子j的間距，所述xij為當(dāng)前時(shí)刻粒子i和鄰域粒子j的間距。

在s102中，對(duì)所述形變梯度fi進(jìn)行奇異值分解得到對(duì)角矩陣并計(jì)算所述粒子i的應(yīng)力其中，所述i為單位矩陣，所述λ和μ均為拉梅系數(shù)，由材料屬性決定。

在s103中，計(jì)算所述粒子i的每個(gè)鄰域粒子j對(duì)所述粒子i的作用力根據(jù)計(jì)算結(jié)果對(duì)所述軟組織進(jìn)行形變處理，其中，所述vi為所述粒子i的體積。

通過圖1所示方案，就可以基于點(diǎn)的方法來重現(xiàn)外科手術(shù)交互過程中人體器官的生物力學(xué)反應(yīng)，在軟組織的形變中，其每個(gè)粒子的位置變化僅僅取決于材料性質(zhì)，能夠達(dá)到很好的模型精確度。

此外，作為本發(fā)明的一個(gè)實(shí)施例，為了使現(xiàn)有模型高度可控和無條件穩(wěn)定，可以采用基于位置的方法去快速處理交互過程，如圖2所示，所述方法還包括：

在s201中，計(jì)算所述粒子i的位置修正值其中，所述m為質(zhì)量矩陣，所述ci(x)為所述粒子i周圍約束方程。

在s202中，根據(jù)所述位置修正值對(duì)所述粒子i的位置進(jìn)行修正。

在本發(fā)明實(shí)施例中，對(duì)手術(shù)器械表面網(wǎng)格進(jìn)行簡化采樣，利用這些點(diǎn)代替手術(shù)器械模型來實(shí)現(xiàn)快速碰撞檢測(cè)，在碰撞過程中，基于現(xiàn)有技術(shù)中解決保證物體不相互滲透問題的研究工作，在所述粒子i和所述粒子j的碰撞過程中，通過c(xi，xj)＝|xij|-r≥0來進(jìn)行碰撞限制，其中，所述r為在所述粒子i和所述粒子j的碰撞方向上，為所述粒子i和所述粒子j預(yù)設(shè)的碰撞限制距離，來防止物體間互相滲透現(xiàn)象的發(fā)生。此外，對(duì)于加速物體間的碰撞檢測(cè)，本方案中采用空間哈希網(wǎng)格來加速粒子鄰域信息的獲取速度。

作為本發(fā)明的一個(gè)實(shí)施例，為了防止相互滲透現(xiàn)象的發(fā)生，在粒子i和粒子j的碰撞過程中，還需要計(jì)算因摩擦力而產(chǎn)生的正切方向的位移：

其中，和分別為碰撞的粒子i和粒子j當(dāng)前可能到達(dá)的位置，xi和xj分別為碰撞開始時(shí)刻粒子i和粒子j所在的位置，是粒子i和粒子j接觸點(diǎn)的法線，所述

由此，對(duì)每個(gè)粒子i來說，因摩擦而產(chǎn)生的位移如下：

其中，uk和us分別是動(dòng)摩擦力系數(shù)和靜摩擦力系數(shù)，所述wi為粒子i權(quán)值，所述wj為粒子j權(quán)值，所述d為任意常量，當(dāng)粒子相對(duì)速度小于牽引力閾值，上式可以直接忽略正切方向位移。

對(duì)于與粒子i碰撞的粒子j來說，其因摩擦而產(chǎn)生的位移如下：

由于大多數(shù)的生物軟組織是不可壓縮的，因此，體積守恒在可變形物體的動(dòng)態(tài)仿真中有重要意義。對(duì)于四面體網(wǎng)格模型，在本發(fā)明實(shí)施例中，在每一個(gè)四面體網(wǎng)格上，對(duì)網(wǎng)格頂點(diǎn)(p1，p2，p3，p4)定義一個(gè)體積限制：

其中，pi，j＝pi-pj，所述v0為四面體的初始體積。

接下來，通過一系列實(shí)驗(yàn)來驗(yàn)證本發(fā)明實(shí)施例提供的生物力學(xué)仿真方法的準(zhǔn)確性、魯棒性以及不可壓縮性。首先，對(duì)實(shí)驗(yàn)環(huán)境的說明如下：

下述實(shí)驗(yàn)均在windows7操作平臺(tái)上完成，由c++實(shí)現(xiàn)并且通過opengl渲染。實(shí)驗(yàn)運(yùn)行在intelxeoncpu，2.40ghz，8gbmemory，nvidiagtx650ti的硬件平臺(tái)上。采用的模型為帶有紋理的表面渲染模型以及粒子化表達(dá)的連續(xù)介質(zhì)力學(xué)計(jì)算模型，其中，圖3為實(shí)驗(yàn)中采用的帶有紋理的表面渲染模型，圖4為實(shí)驗(yàn)中采用的連續(xù)介質(zhì)粒子化表達(dá)。

1、對(duì)模型仿真的準(zhǔn)確性驗(yàn)證：

通過對(duì)比sofa(simulationopenframeworkarchitecture)中有限元仿真方法同等條件下肝臟變形標(biāo)記點(diǎn)位移結(jié)果，來驗(yàn)證本方案所提模型在動(dòng)力學(xué)計(jì)算方面的準(zhǔn)確性。圖5左圖為有限元法的變形效果，右圖為本方案的變形效果，可以看出，本方案在模型仿真的精確度上逼迫有限元方法。

2、對(duì)模型的不可壓縮性的驗(yàn)證：

體積守恒可以用來保證手術(shù)對(duì)象的不可壓縮性。圖6是仿真過程中，基于點(diǎn)的方法與本方案在體積方面的比較，圖7是這兩種方法在運(yùn)行時(shí)間上的對(duì)比，通過對(duì)比可以看出，本方案提出的模型在沒有增加運(yùn)行時(shí)間的情況下，具有較好的體積守恒特點(diǎn)，從而更加逼近真實(shí)值。圖8為兩種方法視覺效果對(duì)比渲染圖。

3、對(duì)模型魯棒性的驗(yàn)證：

可變形固體的仿真方法，如質(zhì)點(diǎn)彈簧法，fem法等在大變形或大時(shí)間步長的條件下經(jīng)常會(huì)出現(xiàn)數(shù)值計(jì)算不穩(wěn)定的情況。在本方案的模型中，全部點(diǎn)的位置是由基于點(diǎn)的方法來計(jì)算的，這僅僅取決于材料性質(zhì)。此外，在本方案中引入基于位置的方法，來實(shí)現(xiàn)交互以及體積守恒，上述方案的采用使得本方案模型在大時(shí)間下以及大形變下也能夠達(dá)到穩(wěn)定性要求。圖9和圖10是肝和胃分別在大變形下仿真的情況。雖然超大變形很少在實(shí)際情況中出現(xiàn)，但這是一種檢驗(yàn)穩(wěn)定性的重要手段。

應(yīng)理解，上述實(shí)施例中各步驟的序號(hào)的大小并不意味著執(zhí)行順序的先后，各過程的執(zhí)行順序應(yīng)以其功能和內(nèi)在邏輯確定，而不應(yīng)對(duì)本發(fā)明實(shí)施例的實(shí)施過程構(gòu)成任何限定。

對(duì)應(yīng)于上文實(shí)施例所述的外科手術(shù)軟組織變形的仿真方法，圖11示出了本發(fā)明實(shí)施例提供的外科手術(shù)軟組織變形的仿真裝置的結(jié)構(gòu)框圖。為了便于說明，僅示出了與本實(shí)施例相關(guān)的部分。

參照?qǐng)D11，該裝置包括：

第一計(jì)算單元1101，在外科手術(shù)仿真過程中，對(duì)軟組織中的每個(gè)粒子i，計(jì)算其形變梯度其中，為粒子i的基本矩陣ai，所述j為所述粒子i的鄰域粒子，所述wij為權(quán)值，所述uij為初始時(shí)刻粒子i和鄰域粒子j的間距，所述xij為當(dāng)前時(shí)刻粒子i和鄰域粒子j的間距；

第二計(jì)算單元1102，對(duì)所述形變梯度fi進(jìn)行奇異值分解得到奇異值矩陣并計(jì)算對(duì)角化后所述粒子i的第一應(yīng)力所述i為單位矩陣，所述λ和μ均為拉梅系數(shù)；

第三計(jì)算單元1103，計(jì)算所述粒子i的每個(gè)鄰域粒子j對(duì)所述粒子i的作用力根據(jù)計(jì)算結(jié)果對(duì)所述軟組織進(jìn)行形變處理，所述vi為所述粒子i的體積。

可選地，所述裝置還包括：

第四計(jì)算單元，計(jì)算所述粒子i的位置修正值其中，

修正單元，根據(jù)所述位置修正值對(duì)所述粒子i的位置進(jìn)行修正。

可選地，所述裝置還包括：

碰撞限制單元，在所述粒子i和所述粒子j的碰撞過程中，通過|xij|-r≥0進(jìn)行碰撞限制，其中，所述r為在所述粒子i和所述粒子j的碰撞方向上，為所述粒子i和所述粒子j預(yù)設(shè)的碰撞限制距離。

可選地，所述裝置還包括：

第五計(jì)算單元，在所述粒子i和所述粒子j的碰撞過程中，計(jì)算所述粒子i因摩擦而產(chǎn)生的位移：

并計(jì)算所述粒子j因摩擦而產(chǎn)生的位移：

其中，uk和us分別是動(dòng)摩擦力系數(shù)和靜摩擦力系數(shù)，xi*和xj*分別為所述粒子i和所述粒子j當(dāng)前可能到達(dá)的位置，xi和xj分別為碰撞開始時(shí)刻所述粒子i和所述粒子j所在的位置，n＝xij*/xij*是所述粒子i和所述粒子j接觸點(diǎn)的法線，所述

可選地，所述裝置還包括：

體積限制單元，對(duì)于四面體網(wǎng)格模型，在每一個(gè)四面體網(wǎng)格上，對(duì)網(wǎng)格頂點(diǎn)(p1，p2，p3，p4)定義體積限制其中，pi，j＝pi-pj，所述v0為四面體的初始體積。

所屬領(lǐng)域的技術(shù)人員可以清楚地了解到，為了描述的方便和簡潔，僅以上述各功能單元、模塊的劃分進(jìn)行舉例說明，實(shí)際應(yīng)用中，可以根據(jù)需要而將上述功能分配由不同的功能單元、模塊完成，即將所述裝置的內(nèi)部結(jié)構(gòu)劃分成不同的功能單元或模塊，以完成以上描述的全部或者部分功能。實(shí)施例中的各功能單元、模塊可以集成在一個(gè)處理單元中，也可以是各個(gè)單元單獨(dú)物理存在，也可以兩個(gè)或兩個(gè)以上單元集成在一個(gè)單元中，上述集成的單元既可以采用硬件的形式實(shí)現(xiàn)，也可以采用軟件功能單元的形式實(shí)現(xiàn)。另外，各功能單元、模塊的具體名稱也只是為了便于相互區(qū)分，并不用于限制本申請(qǐng)的保護(hù)范圍。上述系統(tǒng)中單元、模塊的具體工作過程，可以參考前述方法實(shí)施例中的對(duì)應(yīng)過程，在此不再贅述。

本領(lǐng)域普通技術(shù)人員可以意識(shí)到，結(jié)合本文中所公開的實(shí)施例描述的各示例的單元及算法步驟，能夠以電子硬件、或者計(jì)算機(jī)軟件和電子硬件的結(jié)合來實(shí)現(xiàn)。這些功能究竟以硬件還是軟件方式來執(zhí)行，取決于技術(shù)方案的特定應(yīng)用和設(shè)計(jì)約束條件。專業(yè)技術(shù)人員可以對(duì)每個(gè)特定的應(yīng)用來使用不同方法來實(shí)現(xiàn)所描述的功能，但是這種實(shí)現(xiàn)不應(yīng)認(rèn)為超出本發(fā)明的范圍。

在本發(fā)明所提供的實(shí)施例中，應(yīng)該理解到，所揭露的裝置和方法，可以通過其它的方式實(shí)現(xiàn)。例如，以上所描述的系統(tǒng)實(shí)施例僅僅是示意性的，例如，所述模塊或單元的劃分，僅僅為一種邏輯功能劃分，實(shí)際實(shí)現(xiàn)時(shí)可以有另外的劃分方式，例如多個(gè)單元或組件可以結(jié)合或者可以集成到另一個(gè)系統(tǒng)，或一些特征可以忽略，或不執(zhí)行。另一點(diǎn)，所顯示或討論的相互之間的耦合或直接耦合或通訊連接可以是通過一些接口，裝置或單元的間接耦合或通訊連接，可以是電性，機(jī)械或其它的形式。

所述作為分離部件說明的單元可以是或者也可以不是物理上分開的，作為單元顯示的部件可以是或者也可以不是物理單元，即可以位于一個(gè)地方，或者也可以分布到多個(gè)網(wǎng)絡(luò)單元上?？梢愿鶕?jù)實(shí)際的需要選擇其中的部分或者全部單元來實(shí)現(xiàn)本實(shí)施例方案的目的。

另外，在本發(fā)明各個(gè)實(shí)施例中的各功能單元可以集成在一個(gè)處理單元中，也可以是各個(gè)單元單獨(dú)物理存在，也可以兩個(gè)或兩個(gè)以上單元集成在一個(gè)單元中。上述集成的單元既可以采用硬件的形式實(shí)現(xiàn)，也可以采用軟件功能單元的形式實(shí)現(xiàn)。

所述集成的單元如果以軟件功能單元的形式實(shí)現(xiàn)并作為獨(dú)立的產(chǎn)品銷售或使用時(shí)，可以存儲(chǔ)在一個(gè)計(jì)算機(jī)可讀取存儲(chǔ)介質(zhì)中?；谶@樣的理解，本發(fā)明實(shí)施例的技術(shù)方案本質(zhì)上或者說對(duì)現(xiàn)有技術(shù)做出貢獻(xiàn)的部分或者該技術(shù)方案的全部或部分可以以軟件產(chǎn)品的形式體現(xiàn)出來，該計(jì)算機(jī)軟件產(chǎn)品存儲(chǔ)在一個(gè)存儲(chǔ)介質(zhì)中，包括若干指令用以使得一臺(tái)計(jì)算機(jī)設(shè)備(可以是個(gè)人計(jì)算機(jī)，服務(wù)器，或者網(wǎng)絡(luò)設(shè)備等)或處理器(processor)執(zhí)行本發(fā)明實(shí)施例各個(gè)實(shí)施例所述方法的全部或部分步驟。而前述的存儲(chǔ)介質(zhì)包括：u盤、移動(dòng)硬盤、只讀存儲(chǔ)器(rom，read-onlymemory)、隨機(jī)存取存儲(chǔ)器(ram，randomaccessmemory)、磁碟或者光盤等各種可以存儲(chǔ)程序代碼的介質(zhì)。

以上所述實(shí)施例僅用以說明本發(fā)明的技術(shù)方案，而非對(duì)其限制；盡管參照前述實(shí)施例對(duì)本發(fā)明進(jìn)行了詳細(xì)的說明，本領(lǐng)域的普通技術(shù)人員應(yīng)當(dāng)理解：其依然可以對(duì)前述各實(shí)施例所記載的技術(shù)方案進(jìn)行修改，或者對(duì)其中部分技術(shù)特征進(jìn)行等同替換；而這些修改或者替換，并不使相應(yīng)技術(shù)方案的本質(zhì)脫離本發(fā)明實(shí)施例各實(shí)施例技術(shù)方案的精神和范圍。

以上所述僅為本發(fā)明的較佳實(shí)施例而已，并不用以限制本發(fā)明，凡在本發(fā)明的精神和原則之內(nèi)所作的任何修改、等同替換和改進(jìn)等，均應(yīng)包含在本發(fā)明的保護(hù)范圍之內(nèi)。