鐘海波， 余偉巍， 席 平， 何 飛

（1. 北京航空航天大學(xué)機(jī)械工程及自動(dòng)化學(xué)院，北京 100191；2. 昆明醫(yī)學(xué)院第一附屬醫(yī)院骨科，云南 昆明 650032）

隨著計(jì)算機(jī)技術(shù)和醫(yī)學(xué)三維圖像可視化技術(shù)的發(fā)展，計(jì)算機(jī)輔助虛擬外科手術(shù)在幫助醫(yī)師制訂周密的手術(shù)方案、術(shù)中導(dǎo)航及臨床教學(xué)中均具有重要的指導(dǎo)意義[1]。而在基于醫(yī)學(xué)圖像數(shù)據(jù)（CT, MRI等）所重建的三維模型中進(jìn)行虛擬切割是其中的重點(diǎn)及難點(diǎn)。

文中應(yīng)用醫(yī)學(xué)三維圖像可視化技術(shù)，結(jié)合開發(fā)工具VTK和MFC，建立了基于骨組織的計(jì)算機(jī)輔助模擬手術(shù)系統(tǒng)平臺(tái)，實(shí)現(xiàn)了基于常用繪制算法Marching Cube（移動(dòng)立方體法）和Ray Casting（光線投射法）的三維圖像模型重建。在此基礎(chǔ)上，給出虛擬切割的交互操作技術(shù)以及基于面切割和體切割的三維骨組織模型的切割仿真。

1 VTK簡(jiǎn)介

VTK是William J Schroeder, Kenneth M Martin, Willliam E Lorensen三人于1993年在OpenGL 的基礎(chǔ)上利用面向?qū)ο蠹夹g(shù)，設(shè)計(jì)和實(shí)現(xiàn)的可視化開發(fā)庫，由美國Kitware公司負(fù)責(zé)維護(hù)。

主要功能是進(jìn)行計(jì)算機(jī)圖形、圖像處理和科學(xué)計(jì)算可視化，尤其提供了強(qiáng)大的三維重建功能[2]。VTK的顯著特點(diǎn)就是管道(Pipeline)化，即一個(gè)VTK程序?qū)嶋H上就是一個(gè)完整的渲染管道，其具體流程如圖1所示。

圖1 VTK流程圖

2 切割交互

虛擬切割的交互操作是利用鼠標(biāo)模擬手術(shù)刀，在三維重構(gòu)模型上畫出切割線，選定切割區(qū)域，剖切位于其中的實(shí)體。

（1）坐標(biāo)系之間的轉(zhuǎn)換

在計(jì)算機(jī)圖形學(xué)里，通過二維圖像來表現(xiàn)三維立體圖形，如圖2所示三維立體投影到屏幕上得到平面圖像。VTK中提供了虛擬相機(jī)（位于視點(diǎn)），用作成像工具。

上述投影變換涉及到了多個(gè)坐標(biāo)系，即世界坐標(biāo)系，觀察坐標(biāo)系以及屏幕坐標(biāo)系，分別對(duì)應(yīng)于圖2中坐標(biāo)系OXYZ，O′X′Y′Z′，OeXeYeZe。

世界坐標(biāo)系是用戶所使用的基準(zhǔn)，為右手系，同時(shí)也用以確定相機(jī)及光源的位置和方向。

觀察坐標(biāo)系是虛擬相機(jī)用以成像的坐標(biāo)系，為左手系。坐標(biāo)原點(diǎn)O′位于視點(diǎn)上，也即相機(jī)的位置，視線方向?yàn)閆′軸。相機(jī)的特性可表示為一個(gè)4×4的變換矩陣，用以將世界坐標(biāo)轉(zhuǎn)換為觀察坐標(biāo)。

屏幕坐標(biāo)系是圖像平面上真正的像素位置，也是左手系。坐標(biāo)原點(diǎn)Oe位于視心（視線與屏幕畫面的交點(diǎn)），也即相機(jī)的焦點(diǎn)，Ze軸與Z′軸重合。而窗口在屏幕上的大小等因素決定了觀察坐標(biāo)如何映射為像素位置。

圖2 三維投影原理圖

（2）交互實(shí)現(xiàn)

利用鼠標(biāo)在屏幕上構(gòu)建一個(gè)封閉多邊形，并記錄邊界多邊形頂點(diǎn)的坐標(biāo)，由所獲取的平面封閉多邊形及視線方向OO′構(gòu)成一個(gè)無限長(zhǎng)柱體作為所定義的切割區(qū)域。為了方便交互操作，文中采用了橡皮筋畫線的方式，即隨著鼠標(biāo)的移動(dòng)，動(dòng)態(tài)顯示所畫直線的臨時(shí)圖元。

MFC與VTK均有自己獨(dú)立的顯示窗口和消息響應(yīng)機(jī)制，雖然在系統(tǒng)中已將MFC與VTK進(jìn)行了深層次的集成，構(gòu)建了統(tǒng)一界面和交互響應(yīng)，但是兩者的窗口坐標(biāo)系并不統(tǒng)一，所以在利用MFC鼠標(biāo)響應(yīng)獲取頂點(diǎn)坐標(biāo)構(gòu)建橡皮筋多邊形的同時(shí)，需要利用VTK的鼠標(biāo)響應(yīng)獲取多邊形頂點(diǎn)坐標(biāo)用于坐標(biāo)轉(zhuǎn)換。程序中關(guān)鍵函數(shù)代碼如下：

void OnLButtonDown(UINT nFlags, CPoint point)

{ //MFC的左鍵響應(yīng)函數(shù)，point為MFC窗口中屏幕坐標(biāo)值

……………………………

iren->GetEventPosition(pick[0],pick[1])

//獲取VTK中屏幕坐標(biāo)值

picker->Pick((double)pick[0], (double)pick[1],0.0, ren);

//到世界坐標(biāo)系的坐標(biāo)變換

picker->GetPickPosition( pickPos )

//獲取世界坐標(biāo)系中坐標(biāo)值

……………………………

}

3 VTK中模型切割實(shí)現(xiàn)

虛擬切割主要是模擬手術(shù)刀的效果或是剔除不感興趣部分，操作大多是直接作用于模型數(shù)據(jù)，通過破壞數(shù)據(jù)的結(jié)構(gòu)而達(dá)到切割效果。在VTK管道流各種數(shù)據(jù)結(jié)構(gòu)形式中，文中用到了其中兩種：

· Structured points（結(jié)構(gòu)化點(diǎn)集），其數(shù)據(jù)的拓?fù)鋵傩院蛶缀螌傩远际且?guī)則的。醫(yī)學(xué)圖像中DICOM格式即是這一種，對(duì)應(yīng)于VTK中類vtkImageData。

· Polygonal Data（多邊形體數(shù)據(jù)），由簡(jiǎn)單的圖形（三角形等）組成的復(fù)雜三維圖形通常是用這種數(shù)據(jù)結(jié)構(gòu)來表示，對(duì)應(yīng)于VTK中類vtkPolyData。

醫(yī)學(xué)圖像三維模型分別基于面繪制和體繪制算法重建，文中針對(duì)上述不同模型采用了不同的切割實(shí)現(xiàn)。

3.1 基于Marching Cube模型面切割實(shí)現(xiàn)

Marching Cube算法的基本思想是逐個(gè)處理數(shù)據(jù)場(chǎng)中的立方體（體素），分類出與等值面相交的立方體，并插值計(jì)算等值面與立方體邊的交點(diǎn)來構(gòu)造中間幾何圖元（一般為三角面片）。在VTK管道流中處理的數(shù)據(jù)結(jié)構(gòu)表示為多邊形數(shù)據(jù)(PolyData)。

（1）面切割原理

由切割區(qū)域邊界與構(gòu)成三維模型的三角形面片進(jìn)行求交運(yùn)算[4]。通過判斷三角面片頂點(diǎn)與切割區(qū)域邊界之間的相對(duì)位置來檢測(cè)三角面片與切割區(qū)域的關(guān)系，處在切割區(qū)域之內(nèi)的三角面片被去除，處在切割區(qū)域之外的三角面片保留在模型中，相交的三角面片將進(jìn)行切割運(yùn)算，并在交界處生成新的三角面片，如圖3所示。

結(jié)果改變了輸入多邊形數(shù)據(jù)的內(nèi)部幾何和拓?fù)浣Y(jié)構(gòu)，去除被切割的三角形單元，在交界處生成新的三角形單元，這樣就構(gòu)成了新的多邊形數(shù)據(jù)，實(shí)現(xiàn)了面切割的功能。

圖3 區(qū)域邊界與三角面片位置關(guān)系

（2）面切割實(shí)現(xiàn)

在VTK中，提供了類vtkClipPolyData來裁剪多邊形數(shù)據(jù)，用戶可以構(gòu)建該類對(duì)象并設(shè)置相關(guān)屬性來實(shí)現(xiàn)對(duì)三維重建模型的任意面切割。另外，通過調(diào)用函數(shù)GenerateClippedOutputOn來獲取并保存被切割部分的多邊形數(shù)據(jù)，以此實(shí)現(xiàn)后續(xù)一系列交互功能。如圖4所示為一面切割過程示意圖。

圖4 面切割過程圖

3.2 基于RayCasting模型體切割實(shí)現(xiàn)

RayCasting算法的基本思想是對(duì)于圖像平面上的每個(gè)像素向數(shù)據(jù)場(chǎng)投射光線，在光線上采樣并沿線段積分計(jì)算光亮度和透明度，按采樣順序進(jìn)行圖像合成，得到結(jié)果圖像?？梢陨筛哔|(zhì)量的顯示圖像，由于體繪制計(jì)算數(shù)據(jù)量大，所以相比面繪制，速度較慢。在VTK管道流中處理的數(shù)據(jù)結(jié)構(gòu)表示為結(jié)構(gòu)化點(diǎn)集(Structured points)。

（1）體切割原理

體繪制模型的切割方法一般有兩種：分別基于投射光線或三維體數(shù)據(jù)進(jìn)行處理來達(dá)到切割效果。

· 基于投射光線的切割[5]通過投射光線與切割區(qū)域邊界面的求交，確定新的光線采樣區(qū)間，僅在區(qū)間上進(jìn)行重采樣計(jì)算并合成最終圖像。這種方法的優(yōu)點(diǎn)是保存了數(shù)據(jù)的完整性，并且由于減少了RayCasting算法采樣點(diǎn)的計(jì)算，極大地提高了繪制速度。但是，這種方法不能對(duì)被切割掉的部分進(jìn)行操作。

· 基于三維體數(shù)據(jù)的切割[6]通過直接對(duì)體數(shù)據(jù)進(jìn)行操作，檢測(cè)體素與切割區(qū)域的位置關(guān)系，將處于切割區(qū)域內(nèi)的體素值置為所設(shè)定閾值，使體數(shù)據(jù)丟掉該部分灰度信息，從而使其在最終圖像中不可見來達(dá)到切割效果。由于三維體數(shù)據(jù)量比較大，相比面切割，切割速度會(huì)比較慢。

由于在虛擬切割過程中，往往需要恢復(fù)已被前面切除的部分并進(jìn)行重新定位切割，以及在切割完成之后，可能需要保留被切割部分并對(duì)之進(jìn)行交互操作。所以文中采用第二種方法，可以直接對(duì)體數(shù)據(jù)進(jìn)行操作來完成上述功能。

（2）體切割實(shí)現(xiàn)

由于VTK中體繪制管道流數(shù)據(jù)結(jié)構(gòu)Structured points的特殊性，不能改變其數(shù)據(jù)的拓?fù)浣Y(jié)構(gòu)，并不能像多邊形數(shù)據(jù)那樣去除掉數(shù)據(jù)單元來實(shí)現(xiàn)切割，所以只能改變每個(gè)體素所包含的灰度值信息。

體數(shù)據(jù)量非常巨大，達(dá)到千萬級(jí)（文中CT數(shù)據(jù)512×512×233），如果遍歷體數(shù)據(jù)中的每個(gè)體素并進(jìn)行位置關(guān)系判斷，計(jì)算量將非常大，可以在預(yù)處理階段通過設(shè)定閾值threshold構(gòu)建一vtkIdList鏈表，來記錄有效數(shù)據(jù)點(diǎn)的Id，這樣極大地減少了所需處理的體素?cái)?shù)，提高了切割速度。

空間點(diǎn)與切割區(qū)域之間的位置關(guān)系，一般可以通過點(diǎn)與切割區(qū)域邊界（多邊形柱體）之間的距離dis來判斷。若dis<0，則處于切割區(qū)域內(nèi)，將會(huì)被“切割”。但是，空間點(diǎn)與面之間三維距離的計(jì)算包含大量的乘除法運(yùn)算，將會(huì)影響計(jì)算速度。文中通過投影將三維距離計(jì)算轉(zhuǎn)為點(diǎn)與直線二維距離的計(jì)算，并在判斷過程中加入了包圍盒的檢測(cè)，大大優(yōu)化了計(jì)算速度。

另外，在體切割過程中需要保存被切割的體數(shù)據(jù)，可以通過構(gòu)造一與初始數(shù)據(jù)具有相同拓?fù)浣Y(jié)構(gòu)的空vtkImageData數(shù)據(jù)ImageData_T，來臨時(shí)儲(chǔ)存體數(shù)據(jù)，并在切割完成后，獲取被切割體數(shù)據(jù)的最小包圍盒來優(yōu)化ImageData_T。

綜上所述，通過遍歷vtkIdList中體素，體繪制切割算法的具體描述如下：

1）獲取當(dāng)前體素的空間坐標(biāo)CurrentPoint；

2）將CurrentPoint按視線方向投影到觀察坐標(biāo)平面，獲取投影坐標(biāo)ProjectPoint；

3）將ProjectPoint與切割多邊形包圍盒Bounds進(jìn)行位置關(guān)系的初步判斷。若在區(qū)域內(nèi)，則轉(zhuǎn)至 4)，否則轉(zhuǎn)至 1)；

4）將ProjectPoint與切割多邊形繼續(xù)進(jìn)行位置關(guān)系判斷，若在區(qū)域內(nèi)，則轉(zhuǎn)至 5)，否則轉(zhuǎn)至1)；

5）將當(dāng)前體素灰度值按對(duì)應(yīng)位置存入ImageData_T中，并修改當(dāng)前體素灰度值為所設(shè)定閾值。然后轉(zhuǎn)至 1)；

6）獲取被切割體數(shù)據(jù)最小包圍盒，優(yōu)化ImageData_T。

圖5為一體切割過程示意圖。

圖5 體切割過程圖

4 討 論

本文利用兩種三維重建方法構(gòu)建了基于醫(yī)學(xué)CT數(shù)據(jù)的三維骨組織可視化模型，并分別成功地在三維重建模型完成了虛擬切骨。對(duì)于所述兩種切割方法各有優(yōu)缺點(diǎn)，由于可視化算法和處理數(shù)據(jù)的不同，面切割的速度要明顯優(yōu)于體切割，兩者的切割速度還受切割區(qū)域大小的影響；但面切割是對(duì)骨面的切割，切割后只是空洞的面，而體切割后能看到骨組織的內(nèi)部結(jié)構(gòu)。另外，在利用上述方法進(jìn)行切割過程中，可以把切割后和被切割的數(shù)據(jù)分別保存，這樣就可以實(shí)現(xiàn)重復(fù)切割，撤銷切割以及顯示和操作被切割的數(shù)據(jù)等功能。上述方法均已在輔助模擬手術(shù)系統(tǒng)中實(shí)現(xiàn)。

結(jié)合MFC和VTK實(shí)現(xiàn)了虛擬任意切割功能。采用源代碼公開的VTK類庫作為開發(fā)醫(yī)學(xué)系統(tǒng)的基礎(chǔ)，避免了過多介入計(jì)算機(jī)圖形圖像學(xué)底層的算法，也可以根據(jù)自己的需要，改進(jìn)和增加醫(yī)學(xué)圖像處理所需要的特殊算法，因而可以把更多時(shí)間和精力用于模擬手術(shù)系統(tǒng)功能的開發(fā)，可大大降低開發(fā)系統(tǒng)的難度，減少開發(fā)時(shí)間。

[1]Rhodes Michael L, Roberstson Douglas D. Applications in surgery and therapy [J]. Computer Graphics and Applications, 1996, 16(1): 28-29.

[2]Kitware, Inc. The VTK user’s guide [EB/OL]. USA,http://www.vtk.org.com.

[3]祁俐娜, 羅述謙. 基于VTK的醫(yī)學(xué)圖像三維重建[J].北京生物醫(yī)學(xué)工程, 2006, 25(1): 1-5.

[4]秦緒佳, 歐宗瑛, 侯建華. 醫(yī)學(xué)圖像三維重建模型的剖切與立體視窗剪裁[J]. 計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2002, 14(3): 275-279.

[5]位 軍, 齊 華, 王 毅, 等. 一種實(shí)時(shí)體切割繪制的光線分段預(yù)處理方法[J]. 計(jì)算機(jī)工程與應(yīng)用,2006, (34): 225-227.

[6]王滿寧, 韓正之, 宋志堅(jiān). 三維醫(yī)學(xué)圖像的快速重建與任意切割[J]. 復(fù)旦學(xué)報(bào)(醫(yī)學(xué)版), 2002, 29(2):142-144.