摘 要：在三維表面建模技術(shù)中，Marching Cubes算法是應(yīng)用最為廣泛的方法之一。該算法簡(jiǎn)單高效，但與此同時(shí)，研究人員也發(fā)現(xiàn)它存在一些不足。在構(gòu)造等值面時(shí)，Marching Cubes算法要把所有體素全部檢測(cè)一遍，即使有些體素沒(méi)有和等值面相交，這影響了算法效率；此外在這個(gè)過(guò)程中，Marching Cubes算法還會(huì)忽略掉一些本來(lái)在等值面上的點(diǎn)，降低了表面重建的精度。針對(duì)這些問(wèn)題，本文對(duì)算法進(jìn)行了改進(jìn)。在構(gòu)造等值面時(shí)，不檢測(cè)空的體素以提高算法的速度，并且把一些被忽略的等值點(diǎn)添加進(jìn)來(lái)以提高算法的精度。

關(guān)鍵詞：三維表面建模；Marching Cubes算法；算法效率；逼近精度

中圖分類(lèi)號(hào)：TP311 文獻(xiàn)標(biāo)識(shí)碼：A

Abstract：Marching Cubes algorithm is one of the most widely used methods in 3D surface modeling technologies.Although the algorithm is simple and efficient，researchers have found some shortcomings in it.When constructing the isosurface，the Marching Cubes algorithm detects all voxels，even if some voxels do not intersect with the isosurface，which affects the efficiency of the algorithm.In addition，the Marching Cubes algorithm ignores some points that are originally on the isosurface in this process，which reduces the accuracy of surface reconstruction.In order to solve these problems，this paper improves the algorithm.When constructing isosurface，it does not detect empty voxels to increase the efficiency of the algorithm，and adds some neglected equivalents to improve the accuracy of algorithm.

Keywords：3D surface modeling；Marching Cubes algorithm；algorithmic efficiency；approximation accuracy

1 引言（Introduction）

隨著科學(xué)計(jì)算可視化和計(jì)算機(jī)硬件技術(shù)水平的發(fā)展，三維表面建模技術(shù)越來(lái)越多地應(yīng)用于生活、工程和科研領(lǐng)域[1]。在三維表面建模技術(shù)中，等值面抽取方法的重建效果比較好[2]，而Marching Cubes是其中比較經(jīng)典的算法[3]。在醫(yī)學(xué)領(lǐng)域，醫(yī)生需要在核磁共振圖像（MRI）的幫助下，掌握患者身體內(nèi)部的病變情況，對(duì)病人進(jìn)行治療[4，5]；在地質(zhì)領(lǐng)域，地質(zhì)專(zhuān)家需要通過(guò)地質(zhì)勘測(cè)圖像了解地層結(jié)構(gòu)，指導(dǎo)礦藏開(kāi)采或者預(yù)測(cè)地質(zhì)災(zāi)害[6，7]。這些圖像都可以由Marching Cubes算法重建獲得。然而Marching Cubes算法也存在一些不足之處，并且不斷有學(xué)者和研究人員進(jìn)行研究和改進(jìn)[8，9]。

本文以Marching Cubes算法為基礎(chǔ)，對(duì)算法進(jìn)行了一定的改進(jìn)，提高了算法的效率和逼近精度。利用三維礦體剖面輪廓線數(shù)據(jù)，對(duì)三維礦體表面進(jìn)行重建，還原了礦體原貌。在通過(guò)實(shí)驗(yàn)圖像驗(yàn)證了算法可行性之后，對(duì)算法進(jìn)行了時(shí)間復(fù)雜度分析，并且總結(jié)了算法的優(yōu)缺點(diǎn)。

2 基本原理（The basic principle）

本文的基礎(chǔ)數(shù)據(jù)是礦山剖面輪廓線數(shù)據(jù)，在進(jìn)行三維表面建模之前，先要對(duì)此數(shù)據(jù)進(jìn)行處理以生成體數(shù)據(jù)。生成體數(shù)據(jù)之后，用改進(jìn)的Marching Cubes算法構(gòu)造等值面，完成建模工作。主要步驟有兩個(gè)：（1）構(gòu)造體數(shù)據(jù)，這個(gè)步驟中排除了對(duì)空體素的檢測(cè)；（2）生成等值面，這個(gè)過(guò)程中增加了等值點(diǎn)。

2.1 構(gòu)造體數(shù)據(jù)

體數(shù)據(jù)是體素級(jí)表面建模的基礎(chǔ)。本算法將每一層剖面都定義為一個(gè)Nx×Ny的二維網(wǎng)格，所有剖面垂直疊加就形成一個(gè)規(guī)模為Nx×Ny×Nz的三維空間區(qū)域，算法過(guò)程就在這個(gè)區(qū)域中進(jìn)行。

在一系列二維平面上，定義場(chǎng)函數(shù)為f（x，y），對(duì)于平面上的某個(gè)坐標(biāo)為（x，y）的網(wǎng)格點(diǎn)P的場(chǎng)函數(shù)：（1）如果該點(diǎn)在剖面輪廓線內(nèi)部，其值為正；（2）如果該點(diǎn)位于剖面輪廓線的外部，其值為負(fù)；（3）如果該點(diǎn)落在剖面輪廓線上，其值為0。為了使得重建結(jié)果盡量平滑而且逼近物體真實(shí)表面，算法要選擇合適的場(chǎng)函數(shù)，否則會(huì)出現(xiàn)明顯的突變，造成繪制結(jié)果失真。本文采取歐氏距離函數(shù)作為場(chǎng)函數(shù)：

dis（x，y）表示在這個(gè)二維網(wǎng)格中，該網(wǎng)格點(diǎn)到該剖面上的輪廓線上點(diǎn)的最近距離。

在同樣大小的三維空間區(qū)域中，體數(shù)據(jù)的密度會(huì)影響Marching Cubes算法重建結(jié)果的逼近精度[10]，因此在構(gòu)造體數(shù)據(jù)的過(guò)程中，合理劃分網(wǎng)格單元的大小對(duì)算法結(jié)果有直接的影響。在對(duì)算法處理時(shí)間和存儲(chǔ)空間允許的情況下，盡量將網(wǎng)格劃分的精細(xì)一些，可以得到更為精細(xì)的重建結(jié)果。

在Marching Cubes算法中，等值面與體素棱邊的交點(diǎn)用線性插值求取。通常在計(jì)算場(chǎng)函數(shù)值時(shí)，要對(duì)所有體素端點(diǎn)都進(jìn)行計(jì)算，但是實(shí)際上與等值面相交的體素比起沒(méi)有與等值面相交的體素，數(shù)量要少很多[11]。本文是基于剖面數(shù)據(jù)進(jìn)行表面重建，因此在剖面上考慮此問(wèn)題，對(duì)原算法進(jìn)行了改進(jìn)。假設(shè)剖面是一個(gè)Nx×Ny的二維網(wǎng)格，在某一個(gè)剖面上，與等值線相交的網(wǎng)格數(shù)量占整個(gè)網(wǎng)格的比例很小。

如圖1所示，在一個(gè)剖面上，與等值線相交的網(wǎng)格數(shù)量所占比例很小，如果計(jì)算所有網(wǎng)格頂點(diǎn)的場(chǎng)函數(shù)值，無(wú)效計(jì)算很多，將會(huì)降低算法的效率。為了改善這一狀況，提高算法的效率，需要盡量減少無(wú)效的計(jì)算。為此要設(shè)置一個(gè)標(biāo)識(shí)符，一個(gè)表示各頂點(diǎn)狀態(tài)值的數(shù)組，以及一個(gè)存放場(chǎng)函數(shù)值的數(shù)組，首先判斷頂點(diǎn)與等值面的關(guān)系，若在等值線內(nèi)（包括在等值線上），則其狀態(tài)值為1，反之狀態(tài)值為0。當(dāng)確定了一個(gè)網(wǎng)格四個(gè)頂點(diǎn)的狀態(tài)值后，如果狀態(tài)值全為0或者全為1，則說(shuō)明該網(wǎng)格與等值面無(wú)交點(diǎn)，因此不需要計(jì)算其場(chǎng)函數(shù)值。步驟如下：

Step1：輸入剖面數(shù)據(jù)；

Step2：將剖面數(shù)據(jù)網(wǎng)格化；

Step3：計(jì)算每個(gè)網(wǎng)格四個(gè)頂點(diǎn)的狀態(tài)值，通過(guò)頂點(diǎn)狀態(tài)值確定是否與等值面相交，計(jì)算與等值面相交的網(wǎng)格的四個(gè)頂點(diǎn)的場(chǎng)函數(shù)值；

Step4：逐個(gè)處理每個(gè)剖面，構(gòu)造出體數(shù)據(jù)。

2.2 生成等值面

基于體素級(jí)的三維表面建模，即基于等值面生成的表面建模，需要從大量的體數(shù)據(jù)中把近似表示物體表面的等值面提取出來(lái)。在Marching Cubes算法中，在一個(gè)體素中等值面的提取過(guò)程如下：首先計(jì)算體素每個(gè)頂點(diǎn)的狀態(tài)值，確定該體素是否是和等值面相交的體素；確定了體素和等值面的關(guān)系之后，對(duì)于和等值面相交的體素，計(jì)算其八個(gè)頂點(diǎn)的場(chǎng)函數(shù)值，然后根據(jù)已經(jīng)給定的等值面閾值求取交點(diǎn)（即等值點(diǎn)），求出體素內(nèi)所有的交點(diǎn)后，把這些交點(diǎn)按照一定的方式連接起來(lái)，再進(jìn)行三角剖分，就生成了該體素內(nèi)的等值面。Marching Cubes算法的前提是數(shù)據(jù)場(chǎng)沿著體素棱邊呈線性變化，求取等值點(diǎn)時(shí)根據(jù)其所在棱邊兩個(gè)端點(diǎn)的坐標(biāo)進(jìn)行插值即可求得。

體素的每個(gè)頂點(diǎn)都有兩種狀態(tài)，要么在等值面內(nèi)，要么在等值面外，而8個(gè)頂點(diǎn)，則共有28種狀態(tài)，也就是說(shuō)根據(jù)頂點(diǎn)狀態(tài)的不同，可以把體素類(lèi)型分成256種，記錄256種狀態(tài)比較復(fù)雜，而且這些狀態(tài)中很多是可以相互轉(zhuǎn)換的，根據(jù)旋轉(zhuǎn)對(duì)稱(chēng)性和互補(bǔ)對(duì)稱(chēng)性進(jìn)行簡(jiǎn)化，只需要記錄15種等值面連接構(gòu)型[12]。在求出體素與等值面的交點(diǎn)后，再根據(jù)這15種等值面構(gòu)型連接三角形組成最終的等值多邊形。

將數(shù)據(jù)場(chǎng)內(nèi)所有和等值面相交的體素處理完后，這些體素內(nèi)部的等值面組合起來(lái)就構(gòu)成了所要重現(xiàn)物體的表面網(wǎng)格模型，而可視化最終要將抽象數(shù)據(jù)轉(zhuǎn)變?yōu)槿菀桌斫獾膱D像信號(hào)，因此為了更好地逼真物體真實(shí)表面，需要對(duì)網(wǎng)格模型進(jìn)行光照處理。處理方法為在體素內(nèi)構(gòu)造等值面的同時(shí)，也要計(jì)算每個(gè)剖分三角形三個(gè)頂點(diǎn)的法向量，然后求取平均值作為該三角面片的法向量，利用此法向量就可以計(jì)算該面片上的光照強(qiáng)度，最終得到逼近真實(shí)物體的繪制結(jié)果[13]。

Marching Cubes算法在提取等值面時(shí)，通過(guò)插值計(jì)算出體素棱邊與等值面的交點(diǎn)，然后再通過(guò)一定的方式把插值點(diǎn)連接起來(lái)形成多邊形近似表示物體表面。在實(shí)際情況中，等值面和體素表面的交線是雙曲線，Marching Cubes算法以連接插值點(diǎn)的方法，用直線段近似表示曲線段生成物體表面。在這個(gè)過(guò)程中，原算法并不考慮是否有原輪廓線上的點(diǎn)落入體素之中，因此插值后形成的等值輪廓線和原輪廓線之間存在一定的誤差，而原輪廓線上的點(diǎn)是物體表面上的點(diǎn)，因此使用等值輪廓線提取出的等值面與物體原表面也存在一定的誤差。當(dāng)體數(shù)據(jù)場(chǎng)密度較大的時(shí)候，利用直線段來(lái)近似表示曲線段是可行的，因?yàn)楫?dāng)線段長(zhǎng)度非常小時(shí)，直線段和曲線段的逼近程度會(huì)很高，小到一定程度的時(shí)候肉眼是無(wú)法分辨的。但是如果體數(shù)據(jù)場(chǎng)密度較小，插值求取交點(diǎn)后連接形成的直線段長(zhǎng)度較大，此時(shí)其和實(shí)際的曲線段交線相比較，逼近程度就會(huì)差很多，這個(gè)逼近程度會(huì)隨著體數(shù)據(jù)場(chǎng)密度的減小而降低。

如圖3所示，網(wǎng)格中的閉合實(shí)曲線是原輪廓線，虛線是經(jīng)過(guò)插值計(jì)算后連接插值點(diǎn)形成的等值輪廓線，實(shí)心黑點(diǎn)是原輪廓線上的數(shù)據(jù)點(diǎn)。從圖中可以看到原輪廓線上的點(diǎn)不全在等值線上，這樣就造成連接插值點(diǎn)形成的等值線和原輪廓線吻合程度較低。

為了降低逼近誤差，本文在構(gòu)造等值面時(shí)，除了計(jì)算體素棱邊和等值面的交點(diǎn)之外，還考慮了原輪廓線上的點(diǎn)，對(duì)于有原輪廓線上的點(diǎn)落入的體素，將插值點(diǎn)和原輪廓線上的點(diǎn)以一定的方式連接起來(lái)，然后再通過(guò)三角剖分生成等值面。通過(guò)這一方法，可以使得剖面數(shù)據(jù)中原輪廓線上的點(diǎn)包含在等值線上，而且在連接的時(shí)候加入原輪廓線上的點(diǎn)，那么生成的等值線段就不再是是線段，而是由兩條直線段組成的折線段。如圖4為一個(gè)有原輪廓線點(diǎn)落入的網(wǎng)格中等值線的連接方式，曲線段P1P2是等值面和該體素表面的實(shí)際交線，直線段P1P2是傳統(tǒng)Marching Cubes算法經(jīng)過(guò)插值計(jì)算出交點(diǎn)坐標(biāo)后連接形成的等值線段，Q是落入該表面上的原輪廓線上的點(diǎn)，連接的時(shí)候如果將點(diǎn)Q考慮在內(nèi)，那么最終以折線段P1QP2作為等值線，很明顯可以看到折線段P1QP2和點(diǎn)O圍城的多邊形面積更接近曲線段P1P2和點(diǎn)O圍城的扇形面積。

連接插值點(diǎn)形成等值線時(shí)，如果加入原輪廓線上的點(diǎn)，可能會(huì)改變?cè)瓉?lái)的連接方式。對(duì)于體素的某個(gè)面來(lái)說(shuō)，如果原輪廓線上的點(diǎn)落在該面的某條邊上，那么計(jì)算插值點(diǎn)的時(shí)候，插值點(diǎn)剛好就是這個(gè)原輪廓線上的點(diǎn)，此時(shí)即使將該原輪廓線上的點(diǎn)考慮在內(nèi)，也不會(huì)改變?cè)撁嫔喜逯迭c(diǎn)的連接方式。但是如果原輪廓線上的點(diǎn)完全落入該面內(nèi)，那么就會(huì)改變插值點(diǎn)的連接方式。如圖5所示，圖a和圖b分別是調(diào)整前后的等值點(diǎn)連接方式，調(diào)整后左邊和右邊兩個(gè)網(wǎng)格內(nèi)等值點(diǎn)的連接方式都發(fā)生了改變，中間的網(wǎng)格由于沒(méi)有原輪廓線上的點(diǎn)落入，因此連接方式?jīng)]有改變。

根據(jù)算法前提假設(shè)，場(chǎng)函數(shù)沿著體素棱邊呈線性變化，因此體素的一個(gè)表面四條邊與等值面的交點(diǎn)數(shù)量有0、2和4這三種情況。當(dāng)體素某個(gè)表面四條邊與等值面的交點(diǎn)數(shù)量為0的時(shí)候，必定沒(méi)有原輪廓線上的點(diǎn)落入；當(dāng)交點(diǎn)數(shù)量為2的時(shí)候，如果有原輪廓線上的點(diǎn)落入，那么改變連接方式比較簡(jiǎn)單，只需要依次連接三個(gè)點(diǎn)即可。如果有四個(gè)交點(diǎn)，該面是二義性面，在連接的同時(shí)也要考慮解決面二義性問(wèn)題，不同的連接方式可能會(huì)生成不同的等值面構(gòu)型，甚至可能會(huì)出現(xiàn)拓?fù)溴e(cuò)誤的情況[14]，比較經(jīng)典的解決方法有四面體剖分法[15]和雙曲線漸近線法[16]，不同的方法各有特點(diǎn)，本文采取了一種基于插值點(diǎn)連線交點(diǎn)的方法解決面二義性問(wèn)題[17]。

體素級(jí)表面建模在提取等值面的過(guò)程中，是以直線段近似代替曲線，而此處將原輪廓線上的點(diǎn)考慮進(jìn)來(lái)，實(shí)際上是以兩段折線段近似代替曲線，逼近程度要比之前高一些。在圖6（a）中可以看到，有兩個(gè)交點(diǎn)的情況有兩種，在每種情況中，考慮原輪廓線上的點(diǎn)之后，連接方式是唯一的，將原輪廓線上的點(diǎn)分別和兩個(gè)交點(diǎn)連接起來(lái)即可。但是在二義性面上，因?yàn)閮?nèi)部等值直線段有兩條，因此可能的連接方式有兩種，如圖6（b）所示，兩種不同的連接方式會(huì)導(dǎo)致體素內(nèi)的三角面片拓?fù)潢P(guān)系不同。在這種情況下，結(jié)合輪廓線和交點(diǎn)的關(guān)系，計(jì)算該輪廓線上點(diǎn)到兩個(gè)直線段的距離，取距離較近的那個(gè)直線段為需要調(diào)整的直線段，因?yàn)榫嚯x近的直線段其上面點(diǎn)的場(chǎng)函數(shù)值更接近等值面閾值，這樣更符合實(shí)際情況。

在調(diào)整了體素表面上的等值點(diǎn)連接方式之后，體素內(nèi)部三角剖分形成等值面的方式也隨著會(huì)發(fā)生改變，不能再完全按照之前介紹的15種基本構(gòu)型來(lái)構(gòu)造等值面了。但是如果全部重新調(diào)整等值面構(gòu)型會(huì)增加較大的計(jì)算量，而調(diào)整插值點(diǎn)連接方式的面上，調(diào)整方式都是相對(duì)簡(jiǎn)單的，因此可以在原來(lái)15種構(gòu)型的基礎(chǔ)上進(jìn)行改進(jìn)。

如圖7所示，以Marching Cubes算法中的一種等值面構(gòu)型為例。圖（a）是原來(lái)的等值面構(gòu)型，該體素中一個(gè)頂點(diǎn)在等值面內(nèi)（外），其余七個(gè)頂點(diǎn)都在等值面外（內(nèi)）。按照經(jīng)典Marching Cubes算法原理，在該體素中有三條棱邊和等值面相交，插值計(jì)算出交點(diǎn)坐標(biāo)后，連接這三個(gè)交點(diǎn)構(gòu)成一個(gè)三角形面片，最終以該三角形面片作為該體素內(nèi)的逼近等值面。但是當(dāng)加入原輪廓線上的點(diǎn)之后，上表面的連接方式發(fā)生了改變，此時(shí)就不能再按照?qǐng)D（a）的方式構(gòu)造等值面了。此時(shí)要在原來(lái)構(gòu)型的基礎(chǔ)上，如圖7（b）所示，點(diǎn)P1、P2和P3分別是插值計(jì)算出的交點(diǎn)，Q是原輪廓線上的點(diǎn)，做如下調(diào)整：首先確定該體素中插值點(diǎn)的連接方式，也就是說(shuō)要按照傳統(tǒng)Marching Cubes算法的方法確定等值面構(gòu)型，然后通過(guò)點(diǎn)P1和P2找到原構(gòu)型中和這兩個(gè)點(diǎn)處于同一個(gè)三角形面片中的另外一個(gè)頂點(diǎn)P3，之后分別連接QP1、QP2和QP3，最后去掉P1和P2之間的連線，這樣就完成了調(diào)整。由圖可知加入原輪廓線上的點(diǎn)Q后，通過(guò)調(diào)整等值面構(gòu)型，使用兩個(gè)三角形面片

△QP1P3和△QP2P3代替了原來(lái)的三角面片△P1P2P3，這樣就使得重建后的逼近精度有了一定的提高。

在原Marching Cubes算法15種基本構(gòu)型中，有些體素構(gòu)型中可能包含多個(gè)三角形等值面片，0號(hào)構(gòu)型中沒(méi)有三角形面片，1號(hào)構(gòu)型中有1個(gè)三角形面片，2、3、4、8號(hào)構(gòu)型中都有2個(gè)三角形面片，5、6、7號(hào)構(gòu)型中包含三個(gè)三角形面片，9—14號(hào)構(gòu)型中分別包含四個(gè)三角形面片。當(dāng)體素中三角形面片數(shù)量大于等于2的時(shí)候，就可能出現(xiàn)有公共邊的三角面片。但是通過(guò)觀察這些構(gòu)型，可以發(fā)現(xiàn)在所有存在公共邊面片的體素中，公共邊只在體素內(nèi)部，而在體素表面上不存在包含于兩個(gè)或兩個(gè)以上三角形面片中的等值直線段。本文是基于剖面輪廓線數(shù)據(jù)進(jìn)行三維表面重建，在構(gòu)造體數(shù)據(jù)的時(shí)候，要將相鄰兩層剖面進(jìn)行體素化，因此剖面上的原輪廓線數(shù)據(jù)只分布在體素的上下表面上，由于體素表面上的等值直線段不存在被幾個(gè)三角形面片共用的情況，因此在連接插值點(diǎn)和原輪廓線上的點(diǎn)時(shí)，不需要考慮圖7（b）中P1P2被幾個(gè)三角形面片共用的情況。

經(jīng)過(guò)上面的處理，在構(gòu)造等值面的時(shí)候既考慮了體素棱邊與等值面的交點(diǎn)，也考慮了原輪廓線上點(diǎn)的分布，因此在一個(gè)剖面上，最終形成的等值輪廓線與原輪廓線之間的誤差相比之前有了一定的減小，一定程度上提高了算法的逼近精度，圖8是對(duì)圖3進(jìn)行處理前后的對(duì)比，可以明顯看到二者的區(qū)別。

綜合上述步驟，得出構(gòu)造等值面的步驟如下：

Step1：讀入一個(gè)體素?cái)?shù)據(jù)；

Step2：根據(jù)體素八個(gè)頂點(diǎn)的場(chǎng)函數(shù)值插值計(jì)算棱邊上的交點(diǎn)；

Step3：根據(jù)頂點(diǎn)狀態(tài)索引值查找對(duì)應(yīng)的邊狀態(tài)索引值，再根據(jù)邊狀態(tài)索引值查找該體素內(nèi)的等值面構(gòu)型，確定等值面連接方式；

Step4：確定該體素表面是否包含原輪廓線上的點(diǎn)，如果包含原輪廓線上的點(diǎn)，則執(zhí)行Step5，否則執(zhí)行Step8；

Step5：查找包含原輪廓線點(diǎn)Q的面，計(jì)算該面上插值求出的等值點(diǎn)個(gè)數(shù)n，如果n=4，則執(zhí)行Step6，如果n=2則執(zhí)行Step7；

Step6：按照Step3中已經(jīng)確定的連接方式，計(jì)算點(diǎn)Q到該面上兩條等值線段的距離，取距離較近的等值線段進(jìn)行調(diào)整；

Step7：查找等值線段P1P2所在三角形面片的另外一個(gè)頂點(diǎn)P3，分別連接QP1、QP2和QP3，并且去掉等值線段P1P2；

Step8：算法結(jié)束。

2.3 基于剖面數(shù)據(jù)進(jìn)行三維表面重建的過(guò)程

根據(jù)之前對(duì)算法的介紹，基于剖面數(shù)據(jù)的改進(jìn)MC表面建模流程如下：

Step1：讀入剖面數(shù)據(jù)；

Step2：將相鄰兩層剖面數(shù)據(jù)網(wǎng)格化，得到體素；

Step3：選中體素，求取體素頂點(diǎn)狀態(tài)值，判斷體素是否與等值面相交；

Step4：對(duì)于和等值面相交的體素，求取每個(gè)頂點(diǎn)的場(chǎng)函數(shù)值；

Step5：檢查體素是否還有二義性面，對(duì)二義性面進(jìn)行處理；

Step6：根據(jù)頂點(diǎn)狀態(tài)值找到相應(yīng)的邊狀態(tài)值，根據(jù)邊狀態(tài)值查找對(duì)應(yīng)的等值面連接方式，對(duì)于有原輪廓線上的點(diǎn)落入的體素，對(duì)其中的等值面連接方式進(jìn)行調(diào)整；

Step7：計(jì)算剖分三角形的法向量，計(jì)算光照，繪制該三角形對(duì)應(yīng)的等值面片；

Step8：選擇新的體素轉(zhuǎn)入Step3繼續(xù)處理，直到所有體素都處理完畢；

Step9：結(jié)束。

3 實(shí)驗(yàn)結(jié)果（Experimental result）

本文在對(duì)剖面數(shù)據(jù)進(jìn)行體素化構(gòu)造體數(shù)據(jù)的過(guò)程中，將剖面網(wǎng)格劃分為規(guī)模大小為Nx×Ny的網(wǎng)格，假設(shè)每個(gè)剖面上平均有m條礦體輪廓線剖面，每條輪廓線上有n個(gè)點(diǎn)。下面對(duì)算法效率和結(jié)果進(jìn)行分析說(shuō)明。

（1）算法的時(shí)間復(fù)雜度

算法主要時(shí)間用于構(gòu)造體數(shù)據(jù)和生成等值面兩部分，這兩個(gè)部分中有大量的計(jì)算和比較。體數(shù)據(jù)的構(gòu)造由判斷體素頂點(diǎn)與等值面關(guān)系、計(jì)算頂點(diǎn)場(chǎng)函數(shù)值兩部分組成。生成等值面時(shí)需要查詢邊狀態(tài)表和等值面結(jié)構(gòu)表。

在構(gòu)造體數(shù)據(jù)時(shí)，計(jì)算量主要集中網(wǎng)格頂點(diǎn)狀態(tài)值和場(chǎng)函數(shù)值的計(jì)算，計(jì)算頂點(diǎn)狀態(tài)值的次數(shù)為Nx×Ny，所以時(shí)間復(fù)雜度為O（NxNy）；計(jì)算場(chǎng)函數(shù)的次數(shù)為4×m×n，由于m<

O（NxNy）+O（N）=O（NxNy）

生成等值面的過(guò)程中，需要在每個(gè)體素內(nèi)提取等值面，計(jì)算次數(shù)Nx×Ny；而查詢邊狀態(tài)表和等值面結(jié)構(gòu)表的時(shí)間復(fù)雜度都是是線性的O（N），因此二者相加得到生成等值面等值面的總體時(shí)間復(fù)雜度：

O（NxNy）+O（N）+O（N）=O（NxNy）

算法的總體時(shí)間復(fù)雜度由構(gòu)造體數(shù)據(jù)和生成等值面兩部分的時(shí)間復(fù)雜度相加所得，近似為O（NxNy）。

（2）算法實(shí)驗(yàn)效果

本文以礦體剖面輪廓線數(shù)據(jù)為實(shí)驗(yàn)數(shù)據(jù)，通過(guò)實(shí)驗(yàn)對(duì)算法進(jìn)行了仿真和驗(yàn)證，編程環(huán)境為VS2010，結(jié)果如下。

使用本文算法對(duì)礦體剖面數(shù)據(jù)進(jìn)行處理，正確重建了礦體的三維表面模型，消除了二義性帶來(lái)的孔洞現(xiàn)象，效果如圖9所示為兩種礦體剖面數(shù)據(jù)重建后的線框模型和表面模型。

在重建過(guò)程中，隨著網(wǎng)格規(guī)模的增加，圖形的逼近程度也會(huì)增加，但是由此產(chǎn)生的計(jì)算量增大，因此所耗時(shí)間也會(huì)增加。圖10是對(duì)同一種礦體的剖面數(shù)據(jù)采取不同規(guī)模的網(wǎng)格重建后的結(jié)果。

如表1所示，同樣的剖面數(shù)據(jù)，采用不同的網(wǎng)格規(guī)模進(jìn)行重建，當(dāng)網(wǎng)格規(guī)模為20×20時(shí)所需時(shí)間為0.823s，當(dāng)網(wǎng)格規(guī)模為40×40時(shí)所需時(shí)間為1.452s，因此要得到效果較好的重建圖形并且盡量節(jié)省時(shí)間，需要合理選用網(wǎng)格規(guī)模。

4 結(jié)論（Conclusion）

本文在Marching Cubes算法的基礎(chǔ)上，改進(jìn)了原算法，不再對(duì)空體素頂點(diǎn)場(chǎng)函數(shù)值進(jìn)行計(jì)算，提升了算法效率，并且通過(guò)增加交點(diǎn)的方式，提高了算法的逼近精度，此外也存在一些不足，歸納如下：

減少了無(wú)效的場(chǎng)函數(shù)值計(jì)算，傳統(tǒng)Marching Cubes算法在尋找等值點(diǎn)時(shí)要計(jì)算所有的體素頂點(diǎn)，包括很多本來(lái)和等值面不相交的體素也進(jìn)行了計(jì)算，本算法通過(guò)設(shè)置標(biāo)志位只對(duì)和等值面相交的體素的頂點(diǎn)計(jì)算了場(chǎng)函數(shù)值，減少了無(wú)效計(jì)算，提升了算法效率。

在構(gòu)造等值面的過(guò)程中，傳統(tǒng)Marching Cubes算法只考慮插值計(jì)算求出的交點(diǎn)，將這些交點(diǎn)按照體素對(duì)應(yīng)的等值面構(gòu)型連接起來(lái)表示等值面。為了提高重建結(jié)果的逼近精度，本文在構(gòu)造等值面時(shí)，將原輪廓線上的點(diǎn)和插值點(diǎn)一起連接起來(lái)，使得原輪廓線上的點(diǎn)包含于生成的等值面之中，并且由于使用折線段代替直線段近似表示等值線，使得算法重建的逼近精度有了一定的提高。

當(dāng)數(shù)據(jù)場(chǎng)數(shù)據(jù)密度特別小時(shí)，沒(méi)有解決逼近精度低的問(wèn)題，傳統(tǒng)Marching Cubes算法在較高密度的數(shù)據(jù)場(chǎng)中可以得到很好的表面繪制結(jié)果，但是隨著數(shù)據(jù)場(chǎng)數(shù)據(jù)密度的降低，算法重建結(jié)果的逼近程度也隨著降低，本文在數(shù)據(jù)場(chǎng)數(shù)據(jù)密度特別小時(shí)沒(méi)有做相應(yīng)的處理，將作為后續(xù)問(wèn)題繼續(xù)研究。

參考文獻(xiàn)（References）

[1] Minho Chang，et al.Interactive marching cubes algorithm for intraoral scanners[J].The International Journal of Advanced Manufacturing Technology，2017，89（5-8）：2053-2062.

[2] 黃偉.三維表面建模方法的研究與實(shí)現(xiàn)[D].長(zhǎng)沙：中南大學(xué)，2010.

[3] William E.Lorensen，Harvey E.Cline.Marching cubes：A high resolution 3D surface construction algorithm[J].ACM SIGGRAPH Computer Graphics，1987，21（4）：163-169.

[4] 于洋.肺部CT血管分割及三維重建[D].哈爾濱工業(yè)大學(xué)，2016.

[5] 姚翠萍，周湘連，王晶，等.納米金標(biāo)記細(xì)胞的熒光壽命成像及其三維重建[J]. 西安交通大學(xué)學(xué)報(bào)，2016，50（04）：153-158.

[6] 劉致寧，宋承云，李志勇，等.基于多地震屬性融合分割的地質(zhì)異常體三維模型構(gòu)建[J].Applied Geophysics，2016（3）：519-528.

[7] 鄒艷紅，何建春.移動(dòng)立方體算法的地質(zhì)體三維空間形態(tài)模擬[J].測(cè)繪學(xué)報(bào)，2012，41（06）：910-917.

[8] Seungki Kim.A Novel Interpolation Scheme for Dual Marching Cubes on Octree Volume Fraction Data[J].Computers & Graphics，2017，66（8）：169-178.

[9] 畢碩本，陸源，曾曉文，等.Marching Cubes改進(jìn)算法及其氣象三維模擬[J]系統(tǒng)仿真學(xué)報(bào)，2017，29（07）：1405-1410；1418.

[10] 孫偉，張彩明，楊興強(qiáng).Marching Cubes算法研究現(xiàn)狀[J].計(jì)算機(jī)輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào)，2007，19（7）：947-952.

[11] 錢(qián)峰，馬秀麗，楊勝齊，等.移動(dòng)立方體算法的研究和改進(jìn)[J].計(jì)算機(jī)工程與應(yīng)用，2010，46（34）：177-180.

[12] Adriano Lopes and Ken Brodlie.Improving the Robustness and Accuracy of the Marching Cubes Algorithm for Isosurface[J].IEEE Transactions on Visualization and Computer Graphics，2003，9（1）：l6-29.

[13] 顧耀林，呂理偉.移動(dòng)立方體算法中的三角剖分[J].計(jì)算機(jī)工程與設(shè)計(jì)，2006，27（1）：120-123.

[14] 梁秀霞，張彩明.拓?fù)浣Y(jié)構(gòu)正確的三線性插值曲面的三角片逼近[J].計(jì)算機(jī)研究與發(fā)展，2006，43（3）：528-535.

[15] Hummel R. Exploiting Triangulated Surface Extraction Using Tetrahedral Decomposition[M].IEEE Educational Activities Department，1995，1（4）：328-342.

[16] Nielson G，Hamann B.The asymptotic decider：resolving the ambiguity in marching cubes[C].Proceedings of Visualization' 91，Los Alamitos CA，1991：83-91；413.

[17] 王錚，李瑞明.移動(dòng)立方體算法面二義性問(wèn)題研究[J].軟件工程，2017，20（09）：6-8；5.

作者簡(jiǎn)介：

王 錚（1985-），男，碩士，助教.研究領(lǐng)域：計(jì)算機(jī)圖形學(xué).