基于改進(jìn)MC 和MT 算法的CT 圖像三維重建＊

吳建偉 ，王子牛

(1. 貴州大學(xué)計算機科學(xué)與信息學(xué)院，貴州貴陽550025;2. 貴州大學(xué)網(wǎng)絡(luò)與信息化管理中心中心，貴州貴陽550003)

醫(yī)學(xué)圖像三維重建就是利用二維醫(yī)學(xué)圖像序列(如CT 切片，即計算機斷層掃描切片)重建出三維圖像模型，圖像三維重建的常用算法主要包括體繪制和面繪制。體繪制能表現(xiàn)更多結(jié)構(gòu)上的細(xì)節(jié)，但要達(dá)到滿意的精度則需要很大的運算量，運行和存儲性能難以得到滿足。于是，面繪制成為醫(yī)學(xué)圖像三維重建的主流算法。MC(Marching Cubes)算法是面繪制中的經(jīng)典算法，由W.E.Lorenson 和H.E.Cline 在1987 年提出［1］。MC 算法較好的解決了在標(biāo)量數(shù)據(jù)場中抽取等值面的問題，適合在不規(guī)則的醫(yī)學(xué)體數(shù)據(jù)場中進(jìn)行等值面重建，并以大量的三角面片逼近真實的醫(yī)學(xué)三維圖像。但是，MC 算法也有一些缺點，主要在于抽取等值面時存在拓?fù)浣Y(jié)構(gòu)二義性［2-5］，處理的數(shù)據(jù)量大時對處理速度和存儲空間要求很高。隨著計算機硬件和并行計算的發(fā)展，后一問題得到了一定程度的緩和，而本文主要關(guān)注前一問題。作為MC 算法的發(fā)展，MT(Marching Tetrahedrons)算法［6］沒有等值面抽取的二義性問題，但計算量是MC 算法的好幾倍。為了在不明顯增加計算量的前提下，較好的解決MC 算法的二義性問題，本文結(jié)合兩種算法的優(yōu)點，提出一種改進(jìn)算法，在CT 圖像三維重建的實驗中進(jìn)行了驗證，并得到了較好的結(jié)果。

1 MC 算法

1.1 MC 算法的基本原理

假設(shè)有一個CT 切片數(shù)據(jù)集，可以看作是一個三維標(biāo)量數(shù)據(jù)場，其中每一個切片其實是一個二維數(shù)組，數(shù)組中的某一位置的值即該標(biāo)量數(shù)據(jù)場中該位置的標(biāo)量值。在相鄰兩層切片中分別各取4 個數(shù)據(jù)點，構(gòu)成一個小立方體，稱為立方體素，而整個標(biāo)量數(shù)據(jù)場可由這些立方體素的集合表示。MC算法就是對這些立方體素一個一個的進(jìn)行處理，抽取等值面。首先設(shè)定一個等值面值，對每一個立方體素的8 個頂點，若其頂點的標(biāo)量值:

(1)小于等值面值，則該頂點在等值面內(nèi)部;

(2)等于等值面值，則該頂點在等值面上;

(3)大于等值面值，則該頂點在等值面外部。

可將在等值面內(nèi)部的頂點標(biāo)記為1，其余頂點標(biāo)記為0，然后將8 個標(biāo)記值按一定順序排列，可由一個字節(jié)表示。由于每個頂點的標(biāo)記值有兩種可能，故每個立方體素的頂點狀態(tài)共有28=256 種可能。而每一個頂點狀態(tài)確定了立方體素各邊是否與等值面相交，從而確定立方體素三角面片的拓?fù)錉顟B(tài)。通過旋轉(zhuǎn)對稱性和互補對稱性可將一個立方體素的拓?fù)錉顟B(tài)歸納為以下15 種情形，如圖1 所示:

圖1 立方體素15 種基本拓?fù)錉顟B(tài)

1.2 立方體素的邊與等值面的交點

在確定立方體素的拓?fù)浣Y(jié)構(gòu)后，應(yīng)該計算出立方體素與等值面的交點?？梢约僭O(shè)數(shù)據(jù)場中各點的值沿立方體素邊界呈線性變化，這樣就可以通過線性插值的方法求出立方體素的邊與等值面的交點，即三角面片的頂點。設(shè)立方體素某邊與等值面相交于一點P，該邊的兩個端點坐標(biāo)為(x1，y1，z1)，(x2，y2，z2)，標(biāo)量值為V1，V2，等值面值為V，那么P點的坐標(biāo)(x，y，z)為:

1.3 三角面片的頂點法向量

在OpenGL 中繪制三角面片需要提供三角面片頂點的法向量，從而可以使用光照效果。由于數(shù)據(jù)場中各數(shù)據(jù)點的標(biāo)量值的梯度垂直于等值面，所以數(shù)據(jù)場中某點P 的法向量可用該點處的梯度代替。設(shè)P 點處的法向量為n(x，y，z)，其計算公式為:

其中，V(i，j，k)代表點(i，j，k)處的標(biāo)量值，Δx，Δy，Δz 為三個方向上的立方體素的邊長。

三角面片的頂點是通過對邊上兩端點坐標(biāo)進(jìn)行插補得到的，如式(1)。同樣，在已經(jīng)求得各數(shù)據(jù)點的法向量后，可以通過對邊上兩端點的法向量進(jìn)行插補來得到三角片頂點的法向量。頂點法向量的計算公式仍為式(1)，但(x1，y1，z1)，(x2，y2，z2)分別表示兩端點的法向量，(x，y，z)為要求的頂點法向量。

1.4 MC 算法的二義性

在圖1 中，有些頂點狀態(tài)下，可以有兩種不同的交點連接方式，也就會得到不同的三角面片，即存在二義性問題。MC 算法的二義性問題包括立方體素面二義性和體二義性。

(1)立方體素面二義性

在立方體某一面的兩個對角頂點的標(biāo)量值小于等值面值，而另外兩個對角頂點的標(biāo)量值大于等值面值時，就會有兩種交點連接方式，如圖2 所示:

圖2 立方體素表面的不同交點連接方式

(2)立方體素體二義性

正如立方體素表面存在不同的交點連接方式，立方體素內(nèi)部也有類似情況，比如當(dāng)立方體素的兩個體對角頂點在等值面內(nèi)部，而其余頂點均在等值面外部時，有如圖3 所示兩種連接方式:

圖3 立方體素內(nèi)部的不同交點連接方式

2 MT 算法

MT 算法與MC 算法的不同在于它把一個立方體再細(xì)分為若干個四面體，常見的劃分方法將一個立方體劃分為5 個，6 個或24 個四面體。四面體是最簡單的多面體，其等值面拓?fù)浣Y(jié)構(gòu)可歸結(jié)為3種，如圖4 所示:

圖4 四面體素的3 種拓?fù)錉顟B(tài)

類似于MC 算法，只要確定了四面體的4 個頂點的標(biāo)記狀態(tài)，就可以確定等值面與四面體邊的交點的連接方式，也就得到了用來逼近等值面的三角面片。重要的是，對于一個四面體來說，其邊上的交點的連接方式是唯一的，即對于單個四面體來說其交點的連接不存在二義性。但是若相鄰兩個立方體的剖分為四面體的方式不一致，則會產(chǎn)生剖分二義性。如圖5 所示，這兩個立方體剖分后，在其相接的面上兩條剖分線并不重合，這可能使得兩個立方體中生成的三角面片不能以一致的方式進(jìn)行連接，導(dǎo)致空洞的產(chǎn)生。

圖5 相鄰兩立方體剖分不一致

3 基于MC 和MT 算法的改進(jìn)算法

為了避免MC 算法的面和體二義性，需先判斷立方體是否具有面和體二義性，這可以用一個二義性表來枚舉所有256 種情況。若立方體不具有面和體二義性，可以用MC 算法進(jìn)行處理;否則，對立方體執(zhí)行MT 算法，把立方體剖分為多個四面體進(jìn)行處理［7］。然而，對立方體進(jìn)行剖分會產(chǎn)生剖分二義性問題:一方面，相鄰的兩個剖分立方體之間可能拓?fù)洳灰恢?，如前面的圖5 所示;另一方面，剖分立方體也可能與相鄰的未剖分的立方體拓?fù)洳灰恢隆Ｈ鐖D6 所示，兩個相鄰的立方體中的三角面片沒有嚴(yán)密的連接起來，從而形成了中間白色區(qū)域的空洞。

圖6 剖分立方體與相鄰的未剖分立方體拓?fù)洳灰恢?/p>

下面將分別討論這兩種情形。

3.1 相鄰兩剖分立方體之間拓?fù)洳灰恢?/h3>

圖7 立方體6 -剖分的一種方式

MT 算法常用的立方體剖分方式有5 -剖分，6-剖分和24 -剖分。為了避免相鄰兩個剖分立方體之間的拓?fù)洳灰恢滦?，可以使用一種6 -剖分方式，如圖7 所示:在這種剖分方式下，立方體的左右兩面的剖分線方向相同，因此在橫向方向上的兩個相鄰剖分立方體之間是拓?fù)湟恢碌摹Ｍ?，在另外兩個方向上，兩個相鄰剖分立方體之間也是拓?fù)湟恢碌?。所以，只要相鄰的剖分立方體的剖分都使用如圖7 所示的方式，就可以避免它們之間的剖分二義性問題。

3.2 剖分立方體與相鄰的未剖分立方體拓?fù)洳灰恢?/h3>

對于相鄰兩個立方體中一個是有二義性需要剖分成四面體，另一個是沒有二義性不需要剖分成四面體的情形?？梢詫]有二義性的立方體在每個面上作一條對角線，該對角線的方向需與圖7 中的一致。這樣，剖分立方體和未剖分立方體在其相接的表面上的拓?fù)涫峭耆恢碌模簿筒淮嬖谄史侄x性了，如圖8 所示:

圖8 兩個不同立方體的拓?fù)湟恢滦?/p>

需要注意的是，當(dāng)一個立方體在每個面上各增加一條對角線之后，立方體的邊由12 條增加到18條。因此，MC 算法原有的相交邊索引表不再適用，需要重新生成。同樣，原有的三角面片索引表也需要根據(jù)新的立方體結(jié)構(gòu)重新生成。

4 實驗結(jié)果及分析

本文的實驗環(huán)境為:CPU 為AMD Athlon II X4 640 Processor，800 MHz × 4;內(nèi) 存2G;顯 卡 為NVIDIA GeForce GTS 450;操作系統(tǒng)為32 位Ubuntu14.04;編譯器為gcc 4.6。

使用的CT 數(shù)據(jù)為:大腦，200 ×160 ×160，8位;腹部及骨盆，512 ×512 ×174，8 位。傳統(tǒng)的MC算法，MT 算法與改進(jìn)后的算法的三維重建效果如圖9，圖10 所示。

實驗結(jié)果數(shù)據(jù)的比較見表1。

圖9 大腦表面三維重建圖像

圖10 腹部及骨盆部位三維重建圖像

表1 傳統(tǒng)算法與改進(jìn)算法的實驗結(jié)果數(shù)據(jù)

對比圖9 右邊的局部放大圖可以看出，MC 算法生成了不正確的三角面片，而MT 算法和本文的改進(jìn)算法均生成了正確的三角面片。另外，MT 算法和本文改進(jìn)算法的重建圖像顯示了更多的細(xì)節(jié)，在重建精度上比MC 算法高。從表1 可知，本文改進(jìn)算法的執(zhí)行時間介于MC 算法和MT 算法之間，并且更接近于MC 算法。

5 結(jié)語

針對MC 算法二義性和MT 算法計算量大的問題，提出了一種結(jié)合兩種算法優(yōu)點的改進(jìn)算法。該算法解決了MC 算法的二義性問題，三維重建效果更接近于MT 算法，而計算量更接近于MC 算法，從整體上優(yōu)于這兩種算法。本文實現(xiàn)了傳統(tǒng)的MC 算法，MT 算法和本文提出的改進(jìn)算法，實驗結(jié)果與預(yù)想的效果一致。

［1］Lerensen W E，Cline H E. Marching cubes:A high resolution 3D surface construction algorithm［J］. Computer Graphics，1987，21(3):163 -169.

［2］T S Newman，H Li.A survey of the marching cubes algorithm［J］.Computers Graphics，2006，30(5):854 -879.

［3］王曉，何建英.考慮MC 算法二義性后的三角化方法［J］.湖北工業(yè)大學(xué)學(xué)報，2007，22(6):65 -68.

［4］楊吉宏，薛凌燕，張民，等.基于Double Marching Cubes 的表面重建算法［J］.計算機工程與設(shè)計，2010，31(4):795 -797.

［5］A Gueziec，R Hummel.Exploiting Triangulated Surface Extraction using Tetrahedral Decomposition［J］.IEEE Transactions on Visualisation and Computer Graphics，1995(4):328 -342.

［6］Nielson G M，Hamann B. The asymptotic decider:Resolving the ambiguity in marching cubes［C］//Proceedings of Visualization ’91，San Diego:IEEE，1991:83 -91.

［7］Huiyan Jiang，Ye Zhang，Yudong Zhao.A medical image 3D reconstruction method based on improved MC and MT algorithms［J］.Industrial Electronics and Applications，2009(25 -27):3765 -3768.