基于GPU的MC算法

李丹1 彭海欣2

（1.山東科技大學(xué)計算機學(xué)院，山東 青島 266590；2.山東青島經(jīng)濟技術(shù)開發(fā)區(qū)第一中學(xué)，山東 青島 266555）

基于GPU的MC算法

李丹1 彭海欣2

（1.山東科技大學(xué)計算機學(xué)院，山東 青島 266590；2.山東青島經(jīng)濟技術(shù)開發(fā)區(qū)第一中學(xué)，山東 青島 266555）

MarchingCubes算法是醫(yī)學(xué)圖像三維重建獲取等值面常用的方法，該方法原理簡單，但在遍歷立方體的時候會遍歷許多空立方體，浪費時間。本文提出一種基于頂點狀態(tài)獲取相鄰非空體元的方法，根據(jù)頂點的狀態(tài)判斷某一面是否含有等值邊，含有等值邊的那面相鄰體元有很大可能含有等值邊，減少空立方體的遍歷。此外，本文還將算法在基于GPU的基礎(chǔ)上進(jìn)行并行化處理。

MarchingCubes；醫(yī)學(xué)圖像；相鄰體元；GPU

近年來，三維重建技術(shù)迅速發(fā)展，其應(yīng)用領(lǐng)域十分廣泛，其中在醫(yī)學(xué)方面有著重要的學(xué)術(shù)意義及實踐意義。移動立方體MarchingCubes算法［1］（簡稱MC算法）原理簡單，實現(xiàn)容易，逼近精度較高［2］。但是其仍存在不足之處，遍歷立方體時會存在空立方體，浪費大量精力，避免對空立方體的遍歷，可以提高三維重建速度。由高峰等［3］提出的選取種子體元后根據(jù)種子體元衍生出整個器官的等值面，王旭初［4］提出建立接查找子表約束體元搜索路徑并根據(jù)鄰接查找子表特點設(shè)計堆棧結(jié)構(gòu)實現(xiàn)搜索算法。王溪［5］提出了表面判斷查找表尋找等值面，不需要判斷體數(shù)據(jù)中的每個立方體，但需耗費存儲空間存儲記錄256中每個狀態(tài)下每個三角形的等值邊與面的相對情況，以及計算每條邊的情況。本文在此基礎(chǔ)上直接通過頂點狀態(tài)實現(xiàn)判斷，而不需要建立查找表及計算等值邊。GPU從出現(xiàn)后發(fā)展飛速，并且被利用到并行處理中，而GPU強大的計算能力及高效處理效率，可以更快提高M(jìn)C方法的處理速度。因此，基于GPU的MC方法可以實現(xiàn)更快的重建速度。

1 MarchingCubes算法基本原理

MC算法也就是移動立方體算法，就是利用一個個立方體來獲取等值面信息，而等值面則被化為三角面片，MC算法獲得的就是三角面片信息。立方體也稱體元，是由相鄰8個頂點構(gòu)成，掃描每個體元的8個頂點信息，根據(jù)給定的閾值，若包含閾值范圍，則體元中必定存在等值面，等值面與體元相交的點成為等值點，相交的邊稱為等值邊。將等值點連接則生成所需的等值面，具體過程如下。

1.1 確定含有等值面的立方體

將三維離散數(shù)據(jù)場看作立體柵格系統(tǒng)，圖像信息讀入后，相鄰兩幅圖像構(gòu)成一層，各自相對應(yīng)相近的4個點構(gòu)成一個立方體。由于數(shù)據(jù)場沿著立方體是線性變化的，所以，若一條邊上的2個頂點的灰度值分別大于、小于所規(guī)定的閾值，則該邊上有且僅有一個與等值面相交的點，即等值點。立方體8個頂點的灰度值分別與規(guī)定的閾值比較，以此來判斷等值面是否與該立方體相交。根據(jù)閾值可以將頂點標(biāo)記為標(biāo)記和不標(biāo)記2種狀態(tài)：若頂點的閾值大于規(guī)定的閾值，則頂點在等值面外部，不標(biāo)記該頂點；若頂點的閾值小于規(guī)定的閾值，則頂點在等值面內(nèi)部，標(biāo)記該頂點。每個頂點有2種狀態(tài)，則1個立方體的8個頂點有256種狀態(tài)。

1.2 提取等值面

分析立方體的256種狀態(tài)可以發(fā)現(xiàn)，立方體頂點狀態(tài)具有反對稱性和旋轉(zhuǎn)性，即旋轉(zhuǎn)不會改變等值面的拓?fù)浣Y(jié)構(gòu)，反對稱性就是所有頂點狀態(tài)相反不會改變等值面的連接方式。根據(jù)這種情況可以將256種情況化簡為15種情況，如圖1所示。

圖1 立方體15種基本拓?fù)浣Y(jié)構(gòu)

可以建立邊索引表，將每條邊進(jìn)行編號，邊索引中代表這256種狀態(tài)下每條邊的狀態(tài)，根據(jù)8個頂點的狀態(tài)可以得到一個頂點索引值，根據(jù)該索引值可以查看邊索引表，得到等值邊情況。同時，還可以建立三角形索引記錄256種狀態(tài)下每個三角面片的頂點情況（用頂點編號記錄頂點，每3個代表一個三角面片），用得到的邊索引值獲取三角形索引的三角面片情況。根據(jù)8個頂點的狀態(tài)獲取頂點索引值，根據(jù)索引值查找邊索引，即現(xiàn)在狀態(tài)下的12條邊情況，獲取相對應(yīng)邊索引值。再根據(jù)2個索引值查找三角索引，獲取該立方體的等值面情況。

1.3 計算等值點信息

獲取立方體等值面位置后，就可以計算等值點坐標(biāo)及法向量等信息。等值點計算采取插值計算，常用的有線性插值法和中點選擇法。中點選擇法是取該等值邊的中點，若用p1、p2分別代表2個頂點的坐標(biāo)，n1、n2代表2個頂點的法向量，則公式可表示為：

兩頂點的灰度值分別為I1、I2，等值面閾值為T，則線性插值法公式為：

MC算法的算法步驟為：①逐步掃描兩層數(shù)據(jù)圖像，構(gòu)建體元；②將體元的所有頂點的灰度值與規(guī)定的閾值做比較，建立頂點索引；③若體元含有等值面，利用灰度差分法計算立方體各頂點的梯度；④根據(jù)邊頂點索引查找在邊索引表，獲得存在等值點所在的邊；⑤根據(jù)該相交的邊的2個頂點，用線性插值法獲取等值點的坐標(biāo)信息及法向量；⑥根據(jù)索引值查找三角面片索引表，確定該立方體內(nèi)三角面片的頂點連接方式；⑦由三角面片構(gòu)成等值面；⑧重復(fù)以上步驟直到所有的體元遍歷結(jié)束。

2 基于GPU的改進(jìn)的MC算法

2.1 改進(jìn)的MC算法算法原理

一些MC算法中為提高計算速度，通常采用中點法來計算等值點，其計算量比線性插值法少。此外，大多數(shù)是建立在減少空立方體的遍歷上。在獲取圖像信息時多數(shù)區(qū)域信息是無效的，因此會建立許多空立方體，若避免對這些立方體的遍歷，基本會減少1/2的時間。本文算法主要思想為確定一個種子體元，以此體元為為起點，遍歷相鄰體元。

若某一立方體的某一面存在著等值邊，又因為模型是相連通的，則該面相鄰的立方體的同一面必定存在等值邊，除了該立方體正好是邊界的情況外，含有等值邊的立方體一定存在著等值面。因此，通過判斷立方體面上是否有等值邊，可判斷相鄰立方體中一定存在等值面的立方體。逐步遍歷這些非空立方體，就可避免對非空立方體的處理。又因為模型是連通的，所以基本不會遺漏主要的特征的三角網(wǎng)格。

在觀察立方體的狀態(tài)時可以發(fā)現(xiàn)，若一條邊的2個頂點狀態(tài)相反時必定存在一個等值點，當(dāng)立方體的一個面存在2個等值點時，必定存在等值邊。由此可以發(fā)現(xiàn)，當(dāng)某個面的4個頂點存在2個標(biāo)記、2個不標(biāo)記或者3個頂點狀態(tài)相同、另一頂點狀態(tài)相反時，該面必定存在等值邊，如圖2所示。

圖2 面的等值邊情況

若存在上述狀態(tài)即只要一個面的頂點狀態(tài)不全相同，就可判定該面含有等值邊，同時也就可以判定該面相鄰立方體為非空立方體，則可以遍歷相鄰該方向上的立方體。可以對立方體的頂點和邊進(jìn)行編號，以及根據(jù)立方體6個面定義6個方向，如圖3所示。

圖3 立方體方向及立方體編號情況

若某一面存在等值邊，則該方向標(biāo)記為1，否則為0，并將該方向的體元壓入棧中，而下一個遍歷的體元就從棧中提取。為了防止同一體元壓入棧中或者壓入處理過的體元，應(yīng)對體元設(shè)置標(biāo)記1代表不需要再壓入體元，否則為0，在每次壓入棧前都需查看該體元狀態(tài)標(biāo)記。每個體元在處理過程中不僅需要獲取等值面信息，還需獲取相鄰未處理過的存在等值面的體元信息并將其壓入棧中，然后下一個體元不需要再去遍歷判斷是否是非空體元，直接從棧中讀取非空體元，在讀取體院同時刪除棧中該因素，同時該體元狀態(tài)更新為1。每遍歷一個體元，都可能壓入新的非空體元，需要不斷從棧中提取體元，直至棧中為空。

該算法的算法步驟為：①逐步掃面圖像數(shù)據(jù)構(gòu)造體元；②從中間體元開始，按原MC算法找到第一個非空體元；③該體元獲取完等值面時，獲取6個面的方向標(biāo)記索引；④若某個方向標(biāo)記為1，并且該方向的體元的狀態(tài)0，則將該方向體元坐標(biāo)壓入棧中；⑤下個處理體元從棧中讀取，讀取完后將該體元信息從棧中刪除,并將該體元的狀態(tài)更改為1；⑥以此讀取棧中體元，直到棧為空。

2.2 基于GPU的并行化處理

GPU（Graphics Processing Unit）是圖形處理器的簡稱，這個概念是由NVIDIA公司在發(fā)布GeForce256繪圖處理芯片時首先提出。GPU與CPU類似，為復(fù)雜的數(shù)學(xué)計算和幾何計算而制造，在浮點計算及并行計算等方面有著比CPU高幾十倍的性能。因此，對于GPU對于浮點計算的高性能，能對MC算法的插值計算等計算步驟進(jìn)行更快的計算，并且利用GPU的海量線程進(jìn)行并行計算，在基于GPU的MC方法上具有更高的效率。將原始的MC方法及本文中改進(jìn)的MC方法在基于GPU的環(huán)境下進(jìn)行并行。

2.2.1 原MC算法并行化。原始MC算法主要借用GPU的浮點計算能力，其算法流程如下：①并行的獲取醫(yī)學(xué)圖像數(shù)據(jù)；②根據(jù)獲取的醫(yī)學(xué)圖像構(gòu)建體元；③根據(jù)頂點狀態(tài)判斷是否含有等值點，若沒有返回上一步，繼續(xù)遍歷下一個體元，若有等值面則進(jìn)行下一步；④獲取該體元內(nèi)的等值面信息；⑤GPU中計算，通過線性插值方法，計算出體元邊界與等值而的交點，然后利用中心差分法，求出體元各角點處的法向，再通過線性插值方法，求出三角形各頂點處的法向；⑥遍歷體元直至獲取完所有體元等值面；⑦利用三角面片信息繪制模型。

2.2.2 基于頂點狀態(tài)的MC算法并行化。改進(jìn)的在讀入離散數(shù)據(jù)構(gòu)造立方體后，首先找出第一個非空立方體后，分別對6個面進(jìn)行是否有等值邊的狀態(tài)判斷，并進(jìn)行標(biāo)記，獲取完該體元的等值面后根據(jù)等值邊標(biāo)記將標(biāo)記為1的方向的相鄰體元壓入棧中。計算同樣是在GPU中進(jìn)行的。其算法流程如下：①并行獲取兩層數(shù)據(jù)構(gòu)建立方體；②并行的從立方體中找出第1個非空立方體，從中間開始查找；③對立方體按照正常MC方法進(jìn)行并行獲取等值面信息；④對6個面的狀態(tài)進(jìn)行判斷，獲取相鄰含等值面的體元信息并壓入棧中；⑤并行從棧中獲取待處理的體元，同時將該體元從棧中刪除以及將該體元狀態(tài)改為1；⑥多線程處理體元，重復(fù)從步驟③開始處理，在步驟④時，壓入棧前先判斷體元狀態(tài)是否為0，若為1，則不壓入棧中；⑦并行的重復(fù)處理棧中體元直至為空為止。

3 試驗分析

在visual studio 2013環(huán)境下，用OpenGL實現(xiàn)三維繪制，脊椎骨和頭顱閾值為60。由試驗結(jié)果可以看出，改進(jìn)的MC算法與原MC算法相比，模型樣子基本不變，只是頭顱頸椎部分少了一部分邊緣瑣碎部分，這部分不影響模型本身，而脊椎圖像基本沒有什么變化（見圖4）。

圖4 試驗結(jié)果顯示

表1 試驗結(jié)果對比

表2 并行化試驗結(jié)果對比

由表1結(jié)果可以看出，該算法在原算法的基礎(chǔ)上節(jié)省了大約50%的時間，雖然可能忽略了一些三角面片，但原模型基本特征沒有發(fā)生變化。

并行化試驗結(jié)果見表2，并行化后的算法比原算法快了將近一倍的速度，很明顯比串行速度快；而改進(jìn)后的MC算法時間并行化后比最初時間快了將近3倍，速度有了很大的提高。

4 結(jié)語

對于醫(yī)學(xué)圖像三維重建，本文在原來的MC算法基礎(chǔ)上設(shè)置立方體方向來描述相鄰立方體，用頂點狀態(tài)確定含有等值面的立方體方向，將其壓入棧中待處理，以此避免空立方體的遍歷。本文對其進(jìn)行實驗研究，結(jié)果表明該算法比原算法提高了將近一半的速度。同時，還將算法進(jìn)行并行化處理，可以發(fā)現(xiàn)并行化的繪制速度明顯提高。基于頂點狀態(tài)的MC算法是在圖像模型具有連通性的情況下，在一個模型中可能存在互相不連通的情況，需要建立多個種子節(jié)點。同時，算法需要三維數(shù)組記錄每個立方體的是否處理情況，由于三角面片數(shù)量多，所以三維數(shù)組占用空間內(nèi)存較大。

［1］Willian E Lorense，Harvey E Cline.Marching Cubes：A High Resolution 3d Surface Construction Algorithm［J］.Computer Graphics，1987（4）：163-169

［2］王旭初，王贊.基于最近鄰Marching Cubes的醫(yī)學(xué)圖像三維重建［J］.計算機工程與應(yīng)用，2012（18）：154-158.

［3］高峰，付忠良.基于改進(jìn)移動立方體的醫(yī)學(xué)圖像三維重建算法［J］.計算機應(yīng)用，2013（S1）：201-203，213.

［4］］王旭初，王贊，牛彥敏，等.融合構(gòu)型查找表與鄰接查找子表的改進(jìn)MC方法［J］.重慶大學(xué)學(xué)報，2012（12）：68-77，83.

［5］王溪，秦新強，黨發(fā)寧，等.醫(yī)學(xué)圖像三維重建的規(guī)則移動立方體法［J］.2009（4）：477-481.

MC Algorithm Based on GPU

Li Dan1Peng Haixin2
（1.School of Computer Science，Shandong University of Science and Technology，Qingdao Shandong 266590；2.Shandong
Qingdao Economic and Technological Development Zone the First Middle School，Qingdao Shandong 266555）

MarchingCubes algorithm is a commonly used method for 3D reconstruction of medical images.The meth?od is simple,but a lot of empty cubes will be traversed when traversing the cube,which wastes time.This paper pre?sented a method based on vertex state to obtain adjacent non-empty body elements and determine whether a side contains the equivalent edge according to the state of the vertex to.It is very likely that the adjacent surface element with the same edge can contain the equivalent edge,which reduces the ergodic of the empty cube.In addition,the al?gorithm was parallelized based on GPU.

MarchingCubes；medical image；adjacent voxel；GPU

TP391.41

：A

：1003-5168（2017）03-0057-04

2017-02-16

李丹（1992-），女，碩士，研究方向：CAD與圖形圖像處理。