基于移動頂點的肝臟模型漸進切割算法

張峰峰， 黃軻， 于凌濤， 詹蔚

(1.蘇州大學 機電工程學院，江蘇 蘇州 215006； 2.蘇州大學 蘇州納米科技協同創新中心，江蘇 蘇州 215123； 3.哈爾濱工程大學 機電工程學院，黑龍江 哈爾濱 150001； 4.蘇州大學附屬第一醫院，江蘇 蘇州 215000)

原發性肝癌是常見的消化系統惡性腫瘤之一，其中約85%為肝細胞癌。目前，肝癌已成為全球第二大癌癥相關死亡原因，每年新增肝癌患者30萬～100萬例[1-2]。肝癌的治療復雜，部分肝臟切除手術是避免患者死亡和肝癌轉移最好的方法，但手術過程中醫生無法直觀地觀察到肝臟內部的組織結構。應用增強現實(augment reality，AR)技術能夠將虛擬肝臟模型疊加到真實圖像或視頻中，為醫生提供更多手術信息[3]。通過圖像配準可以將肝臟虛擬模型與真實肝臟的位置重疊，疊加的虛擬模型能夠清晰地顯示肝臟的內部結構，為醫生提供直觀的指導。為了保證增強現實指導的準確性，疊加的虛擬模型必須與真實場景的肝臟保持一致，虛擬模型需要實時進行相應的切割和變形。

肝臟模型的切割既要滿足一定的準確性,又要滿足一定的實時性。實時性、穩定性和真實性是衡量切割效果的關鍵因素[4]，直接決定了切割的效果。實時性與切割的刷新率有關，為了達到流暢的切割效果，通常要求刷新率不低于25 Hz。這就需要算法滿足一定的時間復雜度和空間復雜度要求[5-7]。

用于增強現實的虛擬模型需要由患者CT圖像或MRI圖像重建得到，目前三維重建的方法主要有2大類：體繪制、面繪制。體繪制效果更好，能更好地反映患者體內真實情況，但體繪制模型計算量巨大，很難滿足實際使用的需求。面繪制模型計算量小，便于操作，是目前主流的三維重建方法[8]。

針對面繪制的模型，切割方法又主要有面模型切割法和體模型切割法。面模型切割法由頂點連接成大量三角形構成模型的表面網格。切割過程中主要通過對三角形單元進行操作，改變模型網格的拓撲結構。體模型切割法需要將網格模型進行多面體劃分，使模型內部以多面體單元填充，切割過程中主要對多面體單元進行操作，常用的多面體有四面體、六面體和超六面體[9]。目前大多數學者集中于四面體切割法的研究，根據不同的剖分策略可將四面體切割的方法分為去除法和體元剖分法[10]。去除法是由Bronielsen等[11]提出的，去除法的基本思想是通過碰撞檢測找到與手術器械發生碰撞的四面體單元，將其從模型中刪除。但其真實感差，切口邊緣會產生嚴重的鋸齒。體元剖分法真實感強，切口細膩，但其計算量巨大，隨著切割的進行，單元體數量急劇增加，很難滿足實時性要求[12]。另外Mor[13]提出了一種累進切割法，其基本思想是在手術刀完全掠過某四面體單元之前提前將其分裂。存在計算量大的問題。Nienhuys[14]和Serby等[15]提出了一種移動頂點的方法，在切割過程中搜索切割路徑附近的頂點，將其移動到切割線上。但這種移動模型頂點容易產生退化單元。

目前針對面模型切割的研究較少，各種方法分別存在計算量大和網格質量差的問題，且大多數方法只對表面切割進行了研究，沒有對切割面進行研究。本文針對上述問題對面模型切割的頂點移動法進行了研究和優化，并按照Delaunay三角剖分準則對局部三角面片進行重構，最后對切割面的構造進行了研究。

1 切割算法

如圖1所示，本文算法的核心思路為：

1)進行碰撞檢測，獲得切割點具體位置；

2)搜索最近頂點，將搜索到的最近點移動到切割點的位置；

3)優化局部面片，檢測最近點附近的退化三角形并進行優化；

4)分裂網格，將切割線兩側的頂點分別向不同的方向移動一定的距離以形成切口；

5)構造切割面，連接相應頂點形成切割面。

圖1 移動頂點算法思路Fig.1 Block diagram of moving vertex algorithm

以上所有計算在一幀中完成，稱為一個計算幀，即每一個計算幀從檢測到切割操作開始，切割面構造完結束。

1.1 碰撞檢測

碰撞檢測用于決定何時發生切割并得到切割位置，是后續一切操作的基礎。如圖2所示，常用的碰撞檢測方法為包圍盒法，其中包圍盒又分為包圍球、軸對齊包圍盒(AABB包圍盒)、方向包圍盒(OBB包圍盒)和離散方向凸包圍盒(k-Dops包圍盒)[16-17]。包圍盒算法簡單，計算量小，但無論哪種包圍盒都無法完全貼合模型。為了能夠精確地檢測切割位置，本文針對面模型切割提出了一種新的碰撞檢測方法。

圖2 包圍盒算法對比示意Fig.2 Comparison of bounding box algorithms

在模型切割仿真中，切割工具的簡化形式有3種：1)將切割工具視為與網格表面的一系列碰撞點；2)使用簡單的平面或直線表示切割工具;3)渲染出完整的切割工具，但在切割計算時使用其中軸線進行計算[18]。本文使用第3種方法(如圖3)，在手術刀刀刃上固定2個點Se和St，使用2點之間的連線進行切割計算。

圖3 本文所采用手術刀具示意Fig.3 Schematic diagram of surgical knife used in this article

在切割過程中從點Se向點St投射一條光線，如果光線被肝臟模型遮擋無法到達St即發生切割。檢測到切割后查找遮擋光線的三角面片，記為切割三角形ΔT1T2T3。該三角面片所在平面可以由三點式平面方程得到：

(1)

過點Se和點St的直線由兩點式直線方程等到：

(2)

由式(1)、(2)即可得到切割點的精確位置Ci(xi,yi,zi)。

1.2 搜索最近點

在本文方法中準確地獲取到最近點是算法的關鍵，其結果決定了切割能否成功以及切割的最終效果。本文使用與傳統模型切割不同的計算頻率，本文方法中手術刀每移動一定距離，進行一次碰撞檢測和切割計算，該方法不僅具有更好的延遲反饋，還能通過設置不同的步長來調整計算頻率，在精度和效率之間取得平衡。

通過計算ΔT1T2T3的頂點與切割點Ci距離的方式搜索最近點，但最近點并不總是取計算得出距離最近的點，在切割過程中應該避免圖4所示的特殊情況。

圖4 特殊的最近點Fig.4 Special closest point

在圖4中，當相臨2次計算得到的最近點Mi-1和Mi為一個四邊形的2個相對頂點時，Mi-1和Mi沒有直接連接，無法形成切割線，三角面片不能正常分裂。故每次計算得到臨時最近點Mi后，需要對其和上一個最近點Mi-1進行判斷，判斷條件為：點Mi-1屬于當前切割三角形ΔT1T2T3或點Mi屬于上一個切割三角形ΔT1T2T3。滿足該條件說明點Mi和點Mi-1屬于同一個三角形，2個點由三角形的一條邊相連，能夠正常形成切割線。若不滿足判斷條件，取距離次之的點做為臨時最近點再進行判斷，如果如圖5所示點Mi-1和Mi仍然不能直接連接,則取第3個點做為臨時最近點。如果仍然不滿足條件，表示步長設置過大，2次切割檢測跨過了若干三角形，導致點Mi和點Mi-1所屬的2個三角形沒有公共邊或公共點，此時應適當減小檢測步長。

圖5 特殊的次最近點Fig.5 Special second closest point

1.3 優化面片

狹長或小面積的退化三角形會影響后續的切割計算，還會導致模型顯示質量下降。切割使用的原始模型遵循Delaunay三角剖分原則，不包含退化三角形。本文的移動頂點法每次只移動一個頂點，故不會出現小面積三角形，但并不能避免出現狹長三角形。因為本文的方法每次只移動一個最近點Mi，狹長三角形也只會出現在包含點Mi的三角形中。為了消除狹長三角形，只需要對包含點Mi的三角形進行檢測優化即可。

Delaunay三角剖分是一種由給定點集V形成優質表面網格的算法。Delaunay 三角剖分包含2個重要準則，一個是空圓特性：在Delaunay三角網中任一三角形的外接圓中不存在V中的其他頂點；另一個是最大化最小角特性：在給定點集V可能形成的三角剖分中，Delaunay三角剖分所形成的三角形的最小角最大[19]。

圖6 Delaunay三角剖分Fig.6 Delaunay triangulation

在本文方法中，模型頂點點集V并不是固定不變的，隨著切割的進行，總會有部分頂點坐標發生變化，不能直接使用Delaunay三角剖分。但在4個頂點構成的2個相臨三角形中，如果2個三角形的外接圓范圍內均不包含第4個頂點，則這4個頂點滿足空圓特性，否則總是可以連接四邊形的另一個對角線使其滿足空圓特性。且滿足空圓特性連接方式的最小角大于不滿足空圓特性連接方式的最小角。即滿足空圓特性的連接方式也滿足最大化最小角特性。

圖7 交換對角線Fig.7 Exchange diagonal

本文方法對每個包含最近點Mi的三角形和其相臨三角形組成的點集進行空圓特性判斷，將不滿足空圓特性的點集按照上述方式重新連接。保證每次移動最近點后三角面片中沒有退化三角形。

1.4 構造切割面

在網格模型的切割研究中更多的研究集中于模型表面的切割。不同于先生成完整的切割線再計算切口的方法，漸進式切割的計算是隨著切割動作逐漸進行的，切割面的構造必須也是逐漸進行的。對于本文的方法即在每次頂點移動后對切割面進行計算。另外由于在漸進式切割中不能確定當前計算幀切割是否結束，所以在構造切割面時每一個計算幀都要按照切割結束來生成閉合的切割面。圖8以切割線左側部分為例說明切割面的構造方法。

圖8 切割面的構造過程Fig.8 Cutting surface construction process

為了使切割真實性更好，通過在每一個計算幀中記錄手術刀尖端Sti的坐標的方式得到切割深度。然后依次連接點Mi、Sti-1、Mi-1和點Mi-1、Sti-1、Sti-2構成ΔMiSti-1Mi-1和ΔMi-1Sti-1Sti-2即可形成完整的切割面。

由于每次碰撞檢測為一個計算幀，因此手術刀接觸模型的瞬間即為第一個計算幀，此時表面最近M1與手術刀尖端點St1距離很近，可認為點M1與點St1完全重合，ΔM2St1M1為狹長的退化三角形，故只保留點M1，刪除點St1和ΔM2St1M1。另外，在每一個計算幀中，切割線左右部分共用點Mi、點Sti和點Sti-1，即共用ΔMiStiSti-1，所以不需要構造ΔMiStiSti-1。但在下一個切割循環中點Mi變為Mi-1并在三角面片的分裂過程中被復制一份連接切割線右側三角形，左右部分不再共用ΔMiStiSti-1，所以該三角形需要在下一個計算幀中構造(如圖9)。

圖9 切割面上的特殊三角形Fig.9 Special triangle on the cutting surface

2 實驗結果與分析

使用上述針對面模型切割的移動頂點法對肝臟模型進行切割實驗。實驗硬件環境為：Intel Core i3-8100(主頻3.6 GHz)，8 G內存，NVIDIA Geforce GTX1050Ti顯卡，在Unity 2018.2和Visual Studio 2015中使用C#語言進行編程實現。使用的肝臟模型是經患者CT掃描分割后三維重建得到，由8 706個頂點和15 872個三角形單元組成。

圖10展示了針對面模型切割的移動頂點法的切割效果，圖10(a)為直接將移動頂點法應用于面模型切割的效果。在點M1移動后，出現退化三角形概率較大，雖然隨著切割的進行，隨后的頂點移動有可能使退化三角形正?；匦伦兂烧Ｈ切?，但不能保證每一個退化三角形都能正常化。圖10(b)為使用本文提出的方法優化的結果。在原移動頂點法出現退化三角形的位置經過交換對角線的操作后成功消除了退化三角形，每次頂點移動后都進行檢測優化可以確保每一個退化三角形都能正常化。圖10(c)為本文方法構造出的切割面。本文致力于肝臟面模型切割的拓撲結構變化研究，不涉及生物力學和變形，肝臟切口大小通過一個參數控制。為了更好的顯示切割面，圖10(c)將切口大小參數進行了適當的放大。當切口較小時切口邊緣頂點移動距離小不會出現退化三角形，如圖10(c)切口增大時切口邊緣將出現退化三角形。在實際操作中不能預測切口大小，故需要對切口邊緣的三角形進行優化，使用與最近點相關三角形優化相同的方法對切口邊緣三角形優化結果如圖10(d)。

圖10 模型切割效果Fig.10 Model cutting effect

圖11 多路徑、交叉切割示意Fig.11 Schematic diagram of multi-path, cross cutting

圖12 交叉切割實驗刷新率Fig.12 Cross cutting experiment refresh rate

另外本文對切割過程中出現切割路徑交叉和多次切割的情況進行了實驗，如圖11，并記錄了實驗過程的刷新率。實驗顯示本文的方法具有良好的穩定性。切割過程中，圖像的平均刷新率達到了119 Hz(如圖12)，為模擬組織變形的計算留出了足夠的空間。

3 結論

1)本文提出光線投射算法進行精確的切割檢根據Delaunay三角剖分準則對頂點移動后的局部面片進行了優化和重構以消除原移動頂點法可能出現的退化的畸形三角形。

2)本文還提出了一種針對面模型切割的切割面構造方法，能夠在漸進切割過程中逐步生成完整閉合的切割面，切面真實度大大提高。

3)本文提出的方法很好地解決表面切割模型的實時性問題并具有良好的穩定性，且實時刷新率能夠穩定在119±2 Hz，具有很好的實時性。