利用MATLAB實現(xiàn)CT斷層圖像的三維重建

1、文章編號:1004-4140(200402-0024-06利用MATLAB實現(xiàn)CT斷層圖像的三維重建曾箏1,董芳華1,陳曉2,周宏2,周建中1(1 華中科技大學水電與數(shù)字化工程學院,武漢430074; 2 總后軍需裝備研究所,北京100088摘要:醫(yī)學圖像三維重建為人體結構提供了真實、直觀的反映,便于醫(yī)學人員對病灶的觀察及手術的進行。但圖像三維重建編程實現(xiàn)困難,不易被非計算機專業(yè)人士所掌握。本文運用MA TLAB6.5軟件中的圖像處理工具箱實現(xiàn)了CT斷層圖像的三維表面重建及體重建,原理簡單,編程實現(xiàn)方便。在對頭部CT圖片進行的三維表面重建及體重建實驗中,重建速度快,顯示效果良好,便于各類非計算

2、機專業(yè)人士推廣應用。關鍵詞:CT圖像;表面重建;體重建;MA TLAB中圖分類號:TP391 文獻標識碼:AThree Dimensions Reconstruction of CT Image by MATLAB ZENG Zheng1,DONG Fang-hua1,CHEN Xiao2,ZHOU Hong2,ZHOU Jian-zhong1(1. College of Hydroelectricity & Digital Engineering, HUST, Wuhan 430074, China;2. The Quartermaster Research Institute of

3、 the General Logistics Department of CPLA, Beijing 100088,ChinaAbstract: Three dimensions reconstruction of CT image gives a intuitionistic image of body. It is convenient to observation and operation. But it is difficult to those people not major in computer. In order to solve this problem, how to

4、reconstruct of CT image using image processing toolbox of MATLAB6.5 is discussed. MABLAB gives us a convenient method to reconstruct and acquires a good result. It can be mastered by anybody.Key words: CT image; surface reconstruction; volume reconstruction; MABLAB1引言隨著計算機科學技術的發(fā)展及其與其他各學科間的相互交叉滲透,特別是

5、計算機斷層掃描(CT、核磁共振掃描(MIT等技術逐漸應用于臨床醫(yī)學,使得獲取各種二維醫(yī)學圖像的能力得到了空前的提高1。如何將這些二維數(shù)據(jù)轉化為更容易被接受和處理的信息,正是醫(yī)學數(shù)據(jù)三維可視化所要解決的問題。三維重建技術的實現(xiàn)方法包括兩種2:一種是通過幾何單元拼接擬合物體表面來描述物體的三維結構,稱為表面重建;另一種是直接將體像素以一定的顏色和透明度投影到顯示平面的方法,稱為體重建。比較而言,表面重建運算量小,表面顯示清晰,但對邊緣檢測的要求比較高;而體重建直接基于體數(shù)據(jù)進行顯示,避免了重建過程中所造成的偽像痕跡,但運算量較大3。無論采用哪種重建方法,對于非計算機專業(yè)的人士,編程實現(xiàn)都是很困難的

6、,2期曾箏等:利用MA TLAB實現(xiàn)CT斷層圖像的三維重建 25 因此,本文介紹了基于MATLAB進行三維表面重建及體重建的方法,實現(xiàn)起來靈活、方便,便于掌握。三維重建步驟如圖1所示。 圖1 三維重建流程圖為了有利于從CT圖像中準確地提取出有用的信息,需要對原始圖像進行預處理,以突出有效的圖像信息,消除或減少噪聲的干擾。(1 CT圖像格式的轉換與讀寫。原始的CT圖像是采用DICOM3.0標準進行存儲的,不能被MATLAB所識別,因此必須進行圖像格式的轉換。在這里,使用Visual C+6.0自行開發(fā)轉換程序,在正確讀取DICOM圖像后,通過選擇合適的窗寬、窗位,將窗寬范圍內的值通過線性或非線性

7、變換轉換為小于256的值,將CT圖像轉換為256色BMP圖像4。使用MATLAB中的圖像讀入函數(shù)imread( ,可以讀取BMP圖像,同時可以使用圖像寫出函數(shù)imwrite( 及圖像顯示函數(shù)image( 、imshow( 對圖像進行寫出和顯示5。(2 圖像增強。圖像增強就是根據(jù)某種應用的需要,人為地突出輸入圖像中的某些信息,從而抑制或消除另一些信息的處理過程6。使輸入圖像具有更好的圖像質量,有利于分析及識別。圖像增強技術主要包括直方圖修改、圖像平滑、圖像邊緣銳化及偽彩色增強等。分別采用以下函數(shù)實現(xiàn):灰度直方圖均衡化。均勻量化的自然圖像的灰度直方圖通常在低灰度區(qū)間上頻率較大,使得圖像中較暗區(qū)域中

8、的細節(jié)看不清楚,采用直方圖修整可使原圖像灰度集中的區(qū)域拉開或使灰度分布均勻,從而增大反差,使圖像的細節(jié)清晰。直方圖均衡化在MATLAB 中使用histeq( 函數(shù)實現(xiàn)?；叶茸儞Q法。照片或電子方法得到的圖像,常表現(xiàn)出低對比度即整個圖像偏亮或偏暗,為此需要對圖像中的每一個像素的灰度級進行標度變換,擴大圖像灰度范圍,以達到改善圖像質量的目的。這一灰度調整過程可以用imadjust( 函數(shù)實現(xiàn)。平滑與銳化濾波。平滑技術用于平滑圖像中的噪聲,基本采用在空間域上的求平均值或中值,或在頻域上采取低通濾波。在MATLAB中,各種濾波方法都是在空間域中通過不同的卷積模板即濾波算子實現(xiàn),可用fspecial( 函

9、數(shù)創(chuàng)建預定義的濾波算子,然后用filter2( 或26 CT 理論與應用研究 13卷conv2( 函數(shù)在實現(xiàn)卷積運算的基礎上進行濾波。中值濾波是一種基于排序統(tǒng)計理論的抑制噪聲的非線性信號處理技術,其在除去圖像中的孤立點、線的噪聲的同時,很好地保護了圖像的邊緣信息,適用于一些線性濾波器無法勝任地場合。在MATLAB 中使用medfilt( 函數(shù)實現(xiàn)中值濾波7。適當運用上述方法對原始圖像進行處理,將使原始圖像變得較清晰,能夠較真實地反映圖像的結構特征,便于三維重建的處理及顯示。2CT 圖像三維表面重建計算機三維表面重建是指首先運用圖像技術從二維圖像中分割出興趣區(qū)的輪廓曲線,然后經(jīng)圖形處理,得到其三

10、維結構,從而再現(xiàn)原物體的空間結構。因此,對于三維表面重建而言,邊界輪廓的提取尤為重要8。為了便于面部邊界的提取,先對各CT 圖片進行顏色處理,去掉非有效區(qū),如頭發(fā)、支架等部分,并使其色素盡量減少,圖2a 為經(jīng)過格式轉換的頭部斷層圖像。在提取邊界時,首先采用逐行掃描圖片的辦法,通過比較相鄰點的像素值,找到圖片邊界上的一個點,作為切片邊界的起點。然后從邊界起點開始,逐點判斷與之相鄰的八個點,如果某點為圖片的邊界點則記錄下,并開始下一步判斷,直到獲得所有的邊界點。圖2b 通過上述方法得到的面部邊界輪廓曲線。通過這種方法得到的輪廓曲線是以大量的點坐標形式存儲,在保證擬合精度的情況下,希望精簡數(shù)據(jù)點,以

11、減少存儲空間和提高計算機處理速度。因此將圖2b 中的面部邊緣輪廓曲線以中心點為極點轉換到極坐標系中。然后,將圖像畫入以極角0,2ð為橫坐標,極徑為縱坐標的直角坐標系中,如圖3a 所示。采用12階傅立葉級數(shù)對該圖像進行擬合,得到擬合后的圖像如圖3b 所示。 圖2面部邊界提取(a原始CT 圖像(b面部邊界輪廓y0 50 100 150 200 250 300 350 400 75808590 95100 105 110 115(b 擬合后的邊界曲線 (a 擬合前的邊界曲線2期 曾 箏等:利用MA TLAB 實現(xiàn)CT 斷層圖像的三維重建 27在這里,面部模型可以采用以下數(shù)學表示:=+=12

12、1,0,sin cos (2(n n i n i i i M n b n a a 式中n i a ,n i b ,為各層邊界曲線傅立葉級數(shù)擬合的系數(shù):d n f a i n i =cos (1,d n f b i n i sin (1,=。這樣,某一層斷層面部邊界曲線僅使用25個數(shù)據(jù)就可以表示,所占用的存儲空間很小。 下面,用這些表征面部邊界輪廓曲線的數(shù)據(jù)進行三維表面重建。 (1 重建數(shù)據(jù)的采集。運用上述傅立葉級數(shù)的系數(shù),求出邊界上若干個點x ,y 向坐標值,并為其加上適當?shù)膠 坐標值。xo=0: pi/180:2*pi; % x 的值在0,2中選取yo=yo+a(i*cos(i-1*xo+b

13、(i*sin(i-1*xo; % 通過傅立葉系數(shù)求y 值,其中yo 初始值為a 0 consx=consx;yo.*cos(xo; %將x ,y 值從極坐標系轉換到直角坐標系 consy=consy;yo.*sin(xo;consz=consz;ones(1,length(xo*iLayer*(-4.0; %為每一切片層賦予z 坐標值,iLayer 為層數(shù)(2 邊界輪廓曲線表面繪制。surf(consx,consy,consz; %利用surf( 函數(shù)進行三維表面繪制。(3 設置圖像的顏色及陰影效果。colormap(gray; %利用colormap( 函數(shù)為圖像定義顏色集 shading

14、flat; %利用shading 定義顯示圖像的顏色陰影 (4 設置圖像光照效果。light('Position',-80,-262,-200,'style','infinite' %利用light( 函數(shù)為圖像設置光照效果 light('Position',-500,-0,-4500,'style','infinite' light('Position',5000,100,-300,'style','infinite'(5 設置圖像的顯示效果view

15、(-144,20; %利用view( 函數(shù)定義觀察者視角 lighting gouraud; %利用lighting 定義顯示圖像的光線陰影 axis equal; %利用axis 定義顯示圖像的軸圖4 頭部CT 圖像三維表面重建(a view(-144,20(b view(3(c view(-70,3028 CT理論與應用研究 13卷當函數(shù)view(AZ,EL取不同的值時,可以得到圖像不同的視角,其中AZ表示圖像水平方向旋轉的角度,EL表示圖像垂直方向的高度,MATLAB6.5自行對view(2和view(3的方位進行了定義。圖4a為view(-144,20得到的圖像;圖4b為view(3得

16、到的圖像;圖4c為view(-70,30得到的圖像。采用view( 函數(shù)可以從不同視角對三維重構圖像進行觀察。根據(jù)上述表面重建步驟,對頭部切片進行表面重建,結果如圖4所示。由于采用上述的傅立葉級數(shù)擬合含有耳朵的CT圖像層時,耳朵的輪廓信息無法很好地表達,所以對這些層進行處理時,首先去除耳朵。從圖4中可以看出,運用MATLAB程序在進行CT圖像邊界輪廓提取的基礎上得到三維表面重建圖像。重建速度快、效果好;但是面繪制的缺點是信息的丟失比較大,運算量與景物和物體形狀有關。3CT圖像三維體重建體繪制通過計算所有體素對光線的作用得到二維投影圖像,基于體繪制的三維體重建方法計算量不依賴于景物的復雜程度和物

17、體形狀的復雜程度,也不需要對切片的邊界輪廓進行提取,其計算過程不依賴于視點,處理三維采樣信號方便,便于顯示物體的內部結構。但是,三維體重建所需數(shù)據(jù)量大,運算速度較慢。下面,介紹利用MATLAB進行頭部CT圖像三維體重建的過程。(1 重建數(shù)據(jù)的采集。對現(xiàn)有的n幅頭部CT圖像數(shù)據(jù)進行三維數(shù)據(jù)集D的構造,得到的數(shù)據(jù)集D為一個x×y×n的矩陣。image1 = imread('01.bmp' %使用imread( 函數(shù)讀入現(xiàn)有的n幅圖像image2 = imread('02.bmp'imagen = imread('n.bmp'D=c

18、at(3,image1,image2,image3,imagen; %使用cat( 函數(shù)創(chuàng)建三維矩陣D(2 重建數(shù)據(jù)預處理采用上述方法構造的三維數(shù)據(jù)集D,數(shù)據(jù)量大,在體重建中速度慢,并且可能在計算中超出內存。因而,可以根據(jù)實際情況,對數(shù)據(jù)集D進行預處理,減少數(shù)據(jù)量。x y z D = reducevolume(D,a b c; %使用reducevolume( 函數(shù)減少數(shù)據(jù)量,其中a,b,c為x,y和z軸數(shù)據(jù)抽取的比例,根據(jù)數(shù)據(jù)情況自行定義。D = smooth3(D; %使用smooth( 函數(shù)對數(shù)據(jù)進行平滑處理(3 計算數(shù)據(jù)集在顯示平面累計投影。fv = isosurface(x,y,z,

19、D,isovalue; %使用isosurface( 函數(shù)計算數(shù)據(jù)集在顯示平面累計投影,isovalue根據(jù)實際情況自行定義(4 構造三維體重建碎片p = patch(fv, FaceColor', 'yellow', 'EdgeColor', 'none' %使用patch( 函數(shù)對碎片進行構造,并對圖像的顏色,光線進行定義,其中fv是第(3步中得到的。(5 設置圖像的顏色、陰影及顯示效果。2期 曾 箏等：利用 MATLAB 實現(xiàn) CT 斷層圖像的三維重建 colormap(gray; 29 %利用 colormap( 函數(shù)為圖像定義顏

20、色集 view(3; %利用 view( 函數(shù)定義觀察者視角 lighting gouraud; %利用 lighting 定義顯示圖像的光線陰影 axis equal; %利用 axis 定義顯示圖像的軸 daspect(x y z; %使用 daspect( 定義 x、y、z 軸的顯示比例 根據(jù)上述體重建步驟， 對頭部切片 進行體重建，結果如圖 5 所示。 從圖 5 中可以看出，運用 MATLAB 程序在構造體重建碎片的基礎上，實現(xiàn)體三維重建。 重建速度較表面重建慢，但是體重建盡量保護了更多的表面信息，耳朵部分也不需要任何特 殊處理，就能得到很好的重建及顯示效果。 圖5 頭部 CT 圖像三維體重建 結論 本文介紹了運用 MATLAB6.5 圖像處理工具箱進行 CT 斷層圖像的三維重建的原理及實 現(xiàn)方法。并結合實際應用，實現(xiàn)了頭部 CT 圖像的三維表面重建及體重建，編程實現(xiàn)簡單， 顯示效果理想，大大提高了實驗效率，