SPECT醫(yī)學(xué)圖像重建及三維可視化研究

尹方平

廣東機(jī)電職業(yè)技術(shù)學(xué)院,廣東廣州510515

?

SPECT醫(yī)學(xué)圖像重建及三維可視化研究

尹方平

廣東機(jī)電職業(yè)技術(shù)學(xué)院,廣東廣州510515

摘要:針對醫(yī)學(xué)圖像的建模與仿真分析需要解決的幾何模型精度以及物理模型的正確性問題。本文以人體牙齒CT圖像為例實現(xiàn)建模過程，在牙齒分割方面提出了改進(jìn)。通過讀取DICOM圖像，對斷層圖像進(jìn)行裁剪和分割等預(yù)處理，使用移動立方體（MC）算法進(jìn)行表面重構(gòu)，同時使用二次誤差度量（QEM）方法消減和平滑三角面片，使用德洛奈（Delaunay）方法由面網(wǎng)格推出體網(wǎng)格模型。最后根據(jù)目標(biāo)組織選擇相應(yīng)的物理模型采用有限元（FEM）方法做分析，使用OpenGL顯示等效應(yīng)力的受力云圖。分析結(jié)果顯示，本文所進(jìn)行的SPECT醫(yī)學(xué)圖像重建及三維可視化研究結(jié)果與理論相符，通過SPECT圖像的重建，更加準(zhǔn)確的了解到人體組織器官的功能和代謝情況，提高了醫(yī)護(hù)人員對疾病診斷的準(zhǔn)確率。

關(guān)鍵詞:SPECT;圖像分割;三維建模;可視化

Study on SPECT Medical Image Reconstruction and 3D Visualization

YIN Fang-ping

Guangdong Mechanical and Electrical College, Guangzhou 510515, China

Abstract:Aiming at solving the problem of medical image accuracy of geometry model and physical model in the analysis of modeling and simulation. This paper took teeth CT image as an example to realize human teeth modeling process, and put forward to improve the tooth segmentation. By reading the DICOM image, the tomographic images were cropped and segmenting, and used the marching cubes (MC) algorithm for surface reconstruction. At the same time, we used two quadric error metric (QEM) method to reduce and smooth triangle piece, and used the method of Delaunay to deduce the surface mesh to the solid mesh. Finally, we analyzed the target tissue with the finite element method (FEM) according to the choice of the corresponding physical model, and used the OpenGL to display the stress cloud atlas of equivalent stress. The results showed that the reconstruction of medical image SPECT and 3D visualization were consistent with the theory. Doctors could understand the human tissue and organ more accurately through the reconstruction of SPECT images, and the correct rate of diagnosis on disease had been improved.Keywords: SPECT; image segmentation; 3D modeling; visualization

為了提高醫(yī)療診斷和治療規(guī)劃的準(zhǔn)確性與科學(xué)性，二維斷層圖像序列需要轉(zhuǎn)變成為具有直觀立體效果的圖像，以醫(yī)學(xué)圖像（SPECT、PET、CT、MRI等）三維可視化和計算機(jī)虛擬現(xiàn)實為基礎(chǔ)的生物建模與仿真技術(shù)應(yīng)用于現(xiàn)代醫(yī)學(xué)的各個領(lǐng)域，如虛擬手術(shù)，牙齒正畸，人工股骨置換、外科整形等，呈現(xiàn)出較好的發(fā)展前景[1]。

生物醫(yī)學(xué)建模仿真的研究有著很好的理論意義和應(yīng)用價值，因此已成為各國相關(guān)專家學(xué)者研究的一個熱點[2]。目前,基于醫(yī)學(xué)圖像可視化的生物力學(xué)仿真建模受到越來越廣泛的關(guān)注，成為生物力學(xué)研究的熱點領(lǐng)域之一[3]。近年來，國內(nèi)許多高校和科研院所相繼開展了醫(yī)學(xué)仿真相關(guān)問題的研究，其中在醫(yī)學(xué)圖像的幾何建模方面做了大量的研究并取得了突出的成果[4]。然而，目前虛擬人技術(shù)還僅僅處于起步階段，要想實現(xiàn)能夠完全模擬人體的各種生理過程還需要更多努力。

1　醫(yī)學(xué)圖像的可視化

體繪制方法與面繪制相比，不僅可以顯示體數(shù)據(jù)的表面信息[5]，還可以體現(xiàn)其內(nèi)部信息，從而能實現(xiàn)三維醫(yī)學(xué)影像數(shù)據(jù)的真實感顯示，有利于醫(yī)生的全面理解與分析，因此在輔助醫(yī)生診斷和治療方面有著非常重要的作用。由于體繪制方法具有明顯的可以實施并行計算的特點，為了能夠可以實時顯示，本文采用基于GPU的統(tǒng)一計算設(shè)備架構(gòu)（Compute Unified Device Arehiteeture, CUDA）技術(shù)對光線投射法執(zhí)行并行計算。

由DICOM文件讀入的醫(yī)學(xué)圖像數(shù)據(jù)通常比較大，有上百M(fèi)，醫(yī)生感興趣的組織器官的信息只是圖像中的一部分，但是其他的背景信息還是要讀入，因此可以預(yù)先做粗分割，把感興趣的部位以立方體形式大致劃分出，這樣的圖像數(shù)據(jù)就大大減少，利于以后的圖像的處理[6]。

本文以牙齒CT圖像為例，分析牙齒的結(jié)構(gòu)特點（牙冠部分相接觸；牙根部分與牙槽骨離的很近，灰度值也相近；在牙根部分有分叉，拓?fù)浣Y(jié)構(gòu)產(chǎn)生變化），結(jié)合現(xiàn)已有的分割方法的特性，發(fā)現(xiàn)基于形變模型的水平集分割方法是最適合牙齒分割的。如圖1所示為牙齒的水平集分割程序框圖

圖1　水平集分割的程序框圖Fig.1 The block diagram of level set segmentation

其中，中間初始層使用水平集方法分割的參數(shù)定義如表1，除v的取值不同外，初始輪廓設(shè)定在外部與設(shè)定在內(nèi)部時的其他參數(shù)值均一樣。

表1　初始層輪廓參數(shù)設(shè)定Table 1 The parameters set for initial layer profile

分割中間初始層的結(jié)果如圖2所示，圖2(a)初始輪廓設(shè)定在外部，圖2(b)初始輪廓設(shè)定在內(nèi)部。左幅是手動設(shè)定的中間層（初始層）初始輪廓，右幅是中間層（初始層）的分割結(jié)果。

圖2　使用Li的方法進(jìn)行分割初始層Fig.2 The initial layer based on Li segmentation method

從圖2中可以看到，把初始輪廓放在內(nèi)部的效果相對好些。初始輪廓在外部，以圖2(a)所示，在尖部的灰度值與周圍的背景的灰度值相差比較大，因此輪廓向內(nèi)部收縮也就找到了背景與尖部區(qū)域之間的邊界；初始輪廓在內(nèi)部，以圖2(b)所示，在尖部的灰度值與內(nèi)部的目標(biāo)的灰度值相差也比較大，因此輪廓向外擴(kuò)張也就找到了目標(biāo)與尖部之間的邊界。

2　醫(yī)學(xué)圖像的三維建模

2.1面繪制-移動立方體算法

醫(yī)學(xué)圖像建立幾何模型是使用有限元進(jìn)行仿真分析的基礎(chǔ)。模型的表面實際上是連續(xù)的等值面[7]。假設(shè)醫(yī)學(xué)體數(shù)據(jù)中的某個體素的一個數(shù)據(jù)點的像素值小于給定閾值，則將該數(shù)據(jù)點標(biāo)記為1，認(rèn)為該數(shù)據(jù)點位于等值面的外部；相反，如果像素值大于(或等于)等值面的閾值，則將該數(shù)據(jù)點標(biāo)記為0，認(rèn)為該數(shù)據(jù)點位于等值面內(nèi)部。如果在某體素中它的一條邊的一個端點在等值面內(nèi)，而另一個端點在等值面之外，那么，該邊必定會和等值面相交。因此，通過對體素中數(shù)據(jù)點標(biāo)記的判斷結(jié)果就可以知道該體素是否和等值面相交，還可以知道體素中哪條邊上有交點。實際上只需要處理與等值面相交的體素就可以完成對醫(yī)學(xué)三維模型的繪制了，這樣能盡量的減少數(shù)據(jù)計算量。

為了利用圖形硬件顯示等值面圖像，必須給出形成等值面的各三角面片的法向（或三角形頂點的法向），并選擇適當(dāng)?shù)墓庹漳Ｐ秃筒馁|(zhì)進(jìn)行光照計算，才可生成真實感圖形顯示。對于等值面上的每一點，其沿面的切線方向的梯度分量應(yīng)該是零，因此，該點的梯度矢量的方向也就代表了等值面在該點的法向。設(shè)醫(yī)學(xué)體數(shù)據(jù)中某空間坐標(biāo)為(I,j,k)的數(shù)據(jù)點像素值用I(I,j,k)表示，則該數(shù)據(jù)點處的梯度值為g=(gx, gy, gz)。使用中心差分法來計算體素中該數(shù)據(jù)點處的梯度：

同樣利用線性插值，可以求得構(gòu)成等值面的三角面片在某體素邊上的交點的法向量，即三角面片各個頂點的法向量：N=N1+(value- V1)′(N2- N1)/(V2- V1)(3)其中，N1、N2是三角形頂點所在邊的兩端點的法向量，V1、V2分別為兩端點的像素值，value為該等值面的閾值[8]。最后利用求得的三角面片的三個頂點的法向量就可以實現(xiàn)等值面的真實繪制。分割得到的牙齒組織執(zhí)行MC算法后，顯示為圖3。

圖3　MC結(jié)果顯示Fig.3 MC results show

由結(jié)果可以看到，在考慮二義性的情況下得到的結(jié)果與沒有考慮二義性的結(jié)果相比，單從顯示上看不出有什么區(qū)別，但二者生成的三角面片數(shù)，頂點數(shù)不同。此處沒有對MC的結(jié)果進(jìn)行平滑操作，因此從顯示效果看有比較明顯的棱角。

2.2基于Delaunay準(zhǔn)則的建模算法

在進(jìn)行建模的過程中，考慮到Delaunay算法的執(zhí)行效率，選用三角面片經(jīng)過消減98%之后的數(shù)據(jù)。為了檢測點是否在多面體內(nèi)部，針對檢測點Q，隨機(jī)生成不同方向的射線，如果找到一條不與多面體相交的射線，則說明該點在多面體外。以這個思路為基礎(chǔ)，重心射線法可以很好的解決內(nèi)部冗余網(wǎng)格的問題。以三角面片的重心為射線的端點，只要存在一條射線除了該三角面片本身外，不與三角網(wǎng)格中的其他面片相交，則可以知道該三角面片是表面的，否則該三角面片是在內(nèi)部的。圖4左幅圖像是重構(gòu)表面前的結(jié)果，右幅是重構(gòu)表面后的。

圖4　重構(gòu)表面前后的對比圖Fig.4 Comparison chart for reconstruction surface

MC得到的結(jié)果與簡化后的結(jié)果均存儲為V-F表，即頂點-面表，因此以V-F表中的點集作為輸入，使用三維Delanay增量算法進(jìn)行四面體剖分。在此步驟中，使用了標(biāo)準(zhǔn)模板庫（Standard Template Library ,STL）容器類操作，在Delauany增量算法中，使用到了不少的插入和刪除操作，而容器恰好比較適合大量的插入和刪除步驟。

經(jīng)過Delaunay增量算法的操作，得到的結(jié)果圖5(a)，很明顯圖5(a)與圖4相差太大，原來的表面沒有恢復(fù)。解決方法就是刪除重心在重構(gòu)表面外部的四面體。鑒于已經(jīng)得到了沒有內(nèi)部冗余網(wǎng)格的表面，即重構(gòu)后的表面可以看做是一般的多面體。本文判斷點在一般多面體的方法：射線的端點位于檢測點Q，射線的方向為d (1,0,0)，計算射線與多面體的交點。假設(shè)該射線僅與多面體相交于面的內(nèi)點，交點數(shù)量的奇偶性說面了點在多面體內(nèi)部還是外部。如果是奇數(shù)個點，則點位于內(nèi)部；如果是偶數(shù)個點，則點位于外部。刪除重心在重構(gòu)表面外部的四面體后，顯示如圖5(b)右幅。

圖5　刪除重心在重構(gòu)表面外部的四面體前后的對比Fig.5 Comparison of the external surface of the reconstruction without the focus

3　質(zhì)量評價與仿真分析

3.1模型的質(zhì)量評價

高質(zhì)量的四面體網(wǎng)格可以提高有限元數(shù)值計算的精度和效率，如果有很多劣質(zhì)單元，就會嚴(yán)重影響計算的性能，甚至可能導(dǎo)致計算失去意義。所以對四面體網(wǎng)格做質(zhì)量評價是十分有必要的。長期以來，四面體單元質(zhì)量衡量準(zhǔn)則并沒有一個公認(rèn)的標(biāo)準(zhǔn)，在本文中使用一種計算簡單，也經(jīng)常使用的一種標(biāo)準(zhǔn)：

其中r和R分別是四面體內(nèi)切球和外接球的半徑(0,1)。如果一個四面體為正四面體則它的質(zhì)量參數(shù)為1，如果四面體四點共面則質(zhì)量參數(shù)為0。如果一個四面體單元的質(zhì)量參數(shù)小于0.01則可以認(rèn)定此四面體單元是薄元，即認(rèn)定為質(zhì)量差的四面體模型[7]。四面體剖分后的質(zhì)量系數(shù)分布如表2。

表2　四面體的質(zhì)量系數(shù)Table 2 Quality factor of the tetrahedron

由表中可以看到，消減98%后的三角面片和頂點數(shù)本身較少，Delaunay四面體化后大部分的四面體質(zhì)量系數(shù)集中在0.5到1之間，這說明生成的四面體網(wǎng)格質(zhì)量是比較高的?；謴?fù)率體現(xiàn)的是重構(gòu)后的表面三角面片有多少的比例在四面體網(wǎng)格結(jié)構(gòu)中存在。本文對消減90%的數(shù)據(jù)同樣使用上述步驟，結(jié)果四面體的整體質(zhì)量不高，大部分都在0.1到0.5之間，這就需要在使用空間分解法在內(nèi)部添加點，提高整體質(zhì)量。經(jīng)過加內(nèi)部點后，質(zhì)量在0.5到1之間的達(dá)到65%以上，質(zhì)量得到提升，但是恢復(fù)率略有下降，但還在97%左右。

3.2模型的仿真分析

對于人體軟組織模型的研究難度較大，通常采用軟組織的彈性模量、阻尼系數(shù)、密度等物理量表征軟組織的粘彈性、各向異性、非均勻性等[8,9]。在本文中，針對牙齒的材料和性質(zhì)，選擇了線彈性的物理模型，將模型的相關(guān)參量設(shè)置為：彈性模量為1.86 e10 pa，泊松比為0.31。這里設(shè)定施加載荷于牙冠面，方向為沿坐標(biāo)軸Z軸的負(fù)方向，大小為50 N，這個數(shù)值是牙齒正常咀嚼時的大小。設(shè)定的邊界條件是位移邊界條件，在牙齒的根部的節(jié)點的位移向量均為0。通過對模型進(jìn)行結(jié)構(gòu)離散、單元分析、整體分析、引入邊界條件和方程求解等步驟實現(xiàn)對模型的基本有限元分析過程[10]。

在云圖的繪制中，如何選擇從物理量到顏色模式的映射方式是非常重要的。顏色映射的好壞直接影響圖像生成的質(zhì)量，也對物理信息連續(xù)變化的體現(xiàn)產(chǎn)生影響。本文使用的云圖的顯示方法，首先先構(gòu)建一個顏色表，將等效應(yīng)力按照大小分區(qū)間，每個區(qū)間使用不同的顏色值表示。

圖6是使用FEM分析方法對牙齒體網(wǎng)格進(jìn)行計算后的顯示的等效應(yīng)力分布圖，本文同時把體網(wǎng)格導(dǎo)入ABAQUS軟件中，在同樣的位置加載載荷并設(shè)置邊界條件，等效應(yīng)力分布結(jié)果如圖7所示。

圖6　牙齒等效應(yīng)力分布Fig.6 Stress distribution of teeth equivalent

圖7　ABAQUS中牙齒的應(yīng)力分布 Fig.7 Stress distribution of ABAQUS in teeth

圖6、7可見，圖像的應(yīng)力分布很近似。圖6顯示在舌側(cè)齒尖加載荷時，牙頸與牙根部分等效應(yīng)力較大，這個結(jié)論與文獻(xiàn)一致并驗證了MC操作、面片簡化QEM算法、Delaunay剖分的正確性。

4　結(jié)論

SPECT成像技術(shù)廣泛應(yīng)用在核醫(yī)學(xué)臨床診斷中，本文目的是進(jìn)一步實現(xiàn)SPECT成像技術(shù)的可視化程度，進(jìn)而基于有限元方法對醫(yī)學(xué)圖像進(jìn)行建模仿真分析，把建模過程與分析過程整合到一個平臺中。通過應(yīng)用面繪制中的MC算法進(jìn)行三維表面重建并進(jìn)行顯示，基于Delaunay準(zhǔn)則的增量方法建立四面體網(wǎng)格。把有限元法應(yīng)用到牙齒數(shù)據(jù)上，最后求得的等效應(yīng)力使用OpenGL顯示其分布云圖。有力的提升了SPECT圖像重建技術(shù)的使用效果，為該技術(shù)的進(jìn)一步推廣奠定了理論基礎(chǔ)。

參考文獻(xiàn)

[1] Lacroix Damien, Pati ?o Juan Fernando Ramírez. Finite element analysis of donning procedure of a prosthetic transfemoral socket[J].Annals of biomedical engineering, 2011,39(12):2972-2983

[2] Erdemir Ahmet, Guess Trent M, Halloran Jason , et al. Considerations for reporting finite element analysis studies in biomechanics[J]. Journal of biomechanics , 2012,45(4):625-633

[3]周筠,樊曉平,廖志芳,等.醫(yī)學(xué)體數(shù)據(jù)中四面體化方法的研究進(jìn)展[J].計算機(jī)應(yīng)用研究,2011,28(10):3615-3622

[4]呂曉琪,吳建帥,張明,等.基于擬蒙特卡羅方法的三維醫(yī)學(xué)重建模型體積測量方法研究[J],計算機(jī)應(yīng)用研究,2014,31(2):612-614,618

[5]李艷波,印桂生,張菁,等.Delaunay四面體軟組織建模方法[J].計算機(jī)輔助設(shè)計與圖形學(xué)學(xué)報,2010,22(12):2119-2124

[6]朱代輝,林時苗,楊育彬.醫(yī)學(xué)三維影像體數(shù)據(jù)閾值分割方法[J].計算機(jī)科學(xué),2013,40(1):269-272

[7]周筠.面向生物醫(yī)學(xué)仿真的表面重建和四面體化研究[D].長沙:中南大學(xué),2012

[8]丁麗娟,程杞元.數(shù)值計算方法[M].北京:高等教育出版社,2011:15-86

[9]任靖.基于水平集主動輪廓模型的醫(yī)學(xué)圖像分割方法的研究[D].合肥:合肥工業(yè)大學(xué),2011

[10] Gao H, Chae O. Individual tooth segmentation from CT images using level set method with shape and intensity prior[J]. Pattern Recognition, 2010,43(7):2406-2417

作者簡介:尹方平(1980-),男,廣東惠州人,副教授,碩士,主要從事應(yīng)用數(shù)學(xué)，圖像處理，模式識別、智能交通研究. E-mail:styfp@tom.com

收稿日期:2013-05-10修回日期: 2013-05-22

中圖法分類號:R318.6

文獻(xiàn)標(biāo)識碼:A

文章編號:1000-2324(2015)04-0544-05