田應貴，程鵬

（1.四川大學計算機學院，成都 610065；2.四川大學空天科學與工程學院，成都 610065）

0 引言

近年來，計算機技術(shù)飛速發(fā)展，隨著交叉學科的興起，三維重建技術(shù)被廣泛應用到各個行業(yè)中，特別是在醫(yī)療行業(yè)。三維重建算法可分為體繪制和面繪制兩個大類，體繪制愈加真實地反映物體構(gòu)造，然而運算量大，導致實時性較差，而面繪制則比體繪制運算速度快，實時性相對較強[1]，所以在實踐應用中，目前圖像三維重建的主流方法仍舊是面繪制。

Marching Cubes（MC）算法[2]是面繪制算法中的經(jīng)典算法，它是W.Lorensen等人于1987年提出來的一種體素級重建方法[3]；MC算法因原理好理解并且算法易實現(xiàn)而得到了廣泛的應用，但是算法在拓撲一致性、繪制精度以及重建速度上仍有許多需要改進和完善的地方。

本文以保持好的視覺效果同時提高重建實時性為出發(fā)點，展示了一種改進的MC算法，最終在Visual Studio 2015（v140）環(huán)境下，使用OpenGL可編程管線基于GPU繪制出重建后的三維圖像，有效地實現(xiàn)了等值面快速提取，提高了醫(yī)學圖像可視化效果以及三維重建的效率。

1 移動立方體算法

移動立方體算法（下稱MC算法）是在三維數(shù)據(jù)場中生成等值面的經(jīng)典算法，它使用體素立方體抽取等值面，一般處理三維正交數(shù)據(jù)場[4]，表示：

定義1：

設(shè) f(x,y,z)為三維連續(xù)圖像，頻域上占有限帶寬，覆蓋空間區(qū)域V，空間取樣函數(shù)S(x,y,z)是理想的δ脈沖陣列[5]：

其中：

得到圖像：

表示為三維立體數(shù)據(jù)。

1.1 等值面剖分模型

MC算法以體元作為基本處理單元處理正交場中的數(shù)據(jù)[6]，體元是由上下兩張切片中8個相鄰像素點構(gòu)成的立方體，如圖1。

圖1 體元模型

在使用MC算法處理體元數(shù)據(jù)過程中，一般從左到右，從上到下，從前到后依次處理每個體元，從而通過判定體元中各個頂點是在等值面外還是在等值面內(nèi)對體元的8個頂點進行分類，根據(jù)頂點位置分布即可確定當前立方體等值面所對應的三角剖分狀態(tài)。因為體元的8個頂點各有兩個狀態(tài)，由此我們可以得到256種索引狀態(tài)，再由立方體對稱性、旋轉(zhuǎn)性、共軛性可將這256種剖分狀態(tài)簡化為15種，如圖2。

圖2 15種簡化剖分模型

圖中實心標記的頂點即為處于等值面內(nèi)的點，有了這15種三角剖分模型便可建立一個包含256個索引項的三角剖分構(gòu)型查找表，使用它來記錄所有等值面的連接方式。此時根據(jù)編號進行線性插值運算[7-9]便可得到等值面的頂點坐標，公式如下：

其中等值點坐標為P，邊上兩個端點的坐標分別為P1、P2，灰度值分別為V1、V2，給定的閾值為 u0。

1.2 等值面法向量的計算

在三維重建的可視化過程中，得到體元基本數(shù)據(jù)后，還要保證正確的鏈接各個等值點以形成等值面，由于每個點在等值面內(nèi)沿面切線方向上的梯度分量為0，所以此時該點的梯度矢量方向就是要求的法向量[10]。但是等值面通常是兩種不同密度物質(zhì)的分解面，所以導致處于等值面上的點的梯度矢量不為0，傳統(tǒng)MC算法使用線性插值計算法向量[11]：

其中，等值點法向量坐標為N，邊上兩個端點的法向量坐標分別為 N1、N2，灰度值為 V1、V2，u0代表閾值。

2 中點法簡化計算

三維重建中使用線性插值來計算等值點和法向量顯得太過復雜，對于體元之間共有的棱邊還會進行重復計算，一條棱邊重復計算了4次，這增加了重建的時間代價和復雜性，因此本文使用中點法[12]簡化計算，使用棱邊中點代替線性插值點，三角面頂點坐標和法向量坐標表示如下：

其中頂點坐標為 P，頂點法向量為 N，P1、P2，V1、V2分別為棱邊的兩個角點坐標和法向量。

3 OpenGL可編程管線

使用MC算法進行三維重建過程中會產(chǎn)生大量的三角片面，若是使用傳統(tǒng)的固定OpenGL渲染管線，這將會增加CPU的負擔，占用大量的時間空間片，大大降低了三維重建的效率，因此本文采用最新的OpenGL可編程管線進行三維重建，通過頂點著色器傳遞頂點數(shù)據(jù)，片段著色器處理法向量和顏色值進行可視化三維重建。

3.1 頂點處理

在使用MC算法進行三維重建過程中，會產(chǎn)生大量的三角片面，其中包括大量的頂點數(shù)據(jù)需要處理，如果每次得到一個三角片或者一個頂點數(shù)據(jù)就從CPU發(fā)送頂點數(shù)據(jù)到顯卡，CPU會進行大量的重復工作，同時消耗大量的時間，因此本文提出使用頂點緩沖對象（Vertex Buffer Objects，VBO）[13]來管理內(nèi)存分配，它會將大量的頂點數(shù)據(jù)存儲在GPU中，我們可以一次性將大量的頂點數(shù)據(jù)發(fā)送到顯卡，此時頂點著色器幾乎能做到立即訪問接受到的頂點數(shù)據(jù)，這個過程是非?？斓?。頂點緩沖對象的定義代碼如下：

unsigned int VBO;

glGenBuffers（1,&VBO）;

//綁定獨一無二的ID

glBindBuffer（GL_ARRAY_BUFFER,VBO）;

//綁定緩沖類型

glBufferData（GL_ARRAY_BUFFER,sizeof（vertices）,GL_STATIC_DRAW）;

//將頂點數(shù)據(jù)復制到緩沖內(nèi)存中

3.2 頂點數(shù)據(jù)解析

在OpenGL可編程管線中，只要是以頂點屬性的形式進行輸入的數(shù)據(jù)，我們都可以將其傳遞給頂點著色器，這就要求在配置頂點時要具有很強的靈活性，因此我們要手動配置輸入數(shù)據(jù)與頂點著色器頂點屬性的對應關(guān)系，以保證在渲染前告訴OpenGL該如何處理接收到的頂點數(shù)據(jù)，OpenGL中頂點數(shù)據(jù)解析如圖3所示：

圖3 頂點數(shù)據(jù)解析

解析后位置數(shù)據(jù)以32位（4字節(jié)）浮點值的形式存儲，每個位置有三個這樣的數(shù)據(jù)，其中第一個值處于緩沖的起始位置，得到這些信息后OpenGL使用glVertex?AttribPointer函數(shù)配置頂點數(shù)據(jù)的解析方式：

glVertexAttribPointer（0,3,GL_FLOAT,GL_FALSE,3*sizeof（float）,（void*）0）;

glEnableVertexAttribArray（0）;

3.3 著色器狀態(tài)配置

當我們繪制多個物體時，就需要重復的進行數(shù)據(jù)的綁定和鏈接以保證綁定正確的緩沖對象，從而重建過程中配置正確的頂點屬性就顯得極為麻煩，因此本文提出使用頂點數(shù)組對象（Vertex Array Objects，VAO）來把這些狀態(tài)存儲到一個對象中去，在每次需要配置頂點屬性指針時，只需要將那些調(diào)用執(zhí)行一次，之后在繪制物體時將其綁定到相應的VAO即可，頂點數(shù)組對象會存儲glEnableVertexAttribArray和glDisableVertex?AttribArray的調(diào)用，通過glEnableVertexAttribArray配置頂點屬性，通過glVertexAttribPointer調(diào)用相關(guān)聯(lián)的頂點緩沖對象。

圖4 頂點數(shù)組對象（VAO）綁定VBO

創(chuàng)建好VAO后，類似的可以使用glBindVertexAr?ray綁定VAO，再對屬性指針和對應的VBO進行綁定，完成之后解綁VAO以供使用，此后每次繪制物體時，將VAO綁定到特定的設(shè)定上即可，代碼如下：

//1.綁定頂點VAO

glBindVertexArray（VAO）;

//2.把頂點數(shù)組復制緩沖中

glBindBuffer（GL_ARRAY_BUFFER,VBO）;

glBufferData（GL_ARRAY_BUFFER,sizeof（vertices）,ver?tices,GL_STATIC_DRAW）;

//3.設(shè)定頂點屬性指針

glVertexAttribPointer（0,3,GL_FLOAT,GL_FALSE,3*sizeof（float）,（void*）0）;glEnableVertexAttribArray（0）;

//..::繪制代碼（渲染循環(huán)中）::..

//4.繪制代碼

glUseProgram（shaderProgram）;

glBindVertexArray（VAO）;

glDrawArrays（GL_TRIANGLES,0,3）;

4 實驗與分析

采用傳統(tǒng)MC算法與改進的MC算法基于OpenGL可編程管線對醫(yī)學CT掃描數(shù)據(jù)進行三維重建，實驗環(huán)境為AMD A8-4500M 1.90GHz，4GB內(nèi)存，64位Windows 10操作系統(tǒng)，數(shù)據(jù)為200×160×160的人腦CT數(shù)據(jù)，實驗結(jié)果如圖5所示。表1、表2為實驗結(jié)果對比。

圖5 大腦表面重建實驗結(jié)果

表1 OpenGL可編程管線下的實驗對比

表2 改進MC算法在固定管線和可編程管線下的實驗對比

實驗對比結(jié)果可見，本文改進的算法相較于傳統(tǒng)MC算法在生成的三角片數(shù)量上略有增加，重建圖像質(zhì)量上有一定提升，同時在使用了中點法改進算法之后，重建時間提升了21.8%；由表2實驗對比可以看出，在醫(yī)學圖像重建中，使用MC算法會產(chǎn)生大量三角片和頂點坐標數(shù)據(jù)，本文使用OpenGL可編程管線中的頂點緩沖對象和頂點數(shù)組對象處理三角片和頂點數(shù)據(jù)，重建時間提升了55.1%，綜上本文提出的改進方法在提升了重建圖像效果的同時，運行時間取得了明顯的改進。

5 結(jié)語

針對MC算法進行三維重建耗時過長、重建效果不佳等問題，本文使用OpenGL可編程管線結(jié)合改進的MC算法，使用中點法代替線性插值、頂點緩沖對象和頂點數(shù)組對象處理大量的三角片頂點數(shù)據(jù)，經(jīng)過實驗表面該方法在提升圖像重建效果的同時，在算法耗時上得到大大提升，增強了MC算法的實時性，為算法在臨床醫(yī)學領(lǐng)域的應用做出了一定的貢獻。