l 獲取CT斷層圖像數(shù)據(jù)

　　葉輪斷層圖像的獲取是在江蘇大學附屬江濱醫(yī)院的多層螺旋CT機上完成的。該CT機為美國GE公司產(chǎn)品;與普通CT機相比，多層螺旋CT機具有高清晰度、高分辨率、高速掃描、性能優(yōu)越等特點，而且旋轉(zhuǎn)一圈就能采集出8層圖像，實現(xiàn)了真正的高速透視成像。

　　根據(jù)葉輪的結構特點，為了重建時能充分利用葉輪的結構特征，同時也為了減少三維圖像重建的難度與工作量，采集數(shù)據(jù)時盡可能使葉輪中心線與CT機水平軸心線重合(圖1)。采集精度采用512×512 x1.25，掃描后共獲得56幅斷層圖像，如圖2所示為第20、25幅圖像。從圖2中可以看出，經(jīng)多層螺旋CT透視掃描后，葉輪的斷層圖像包含了水泵葉輪的所有內(nèi)部結構，如葉片形狀、輪轂、軸孔、鍵槽等的結構信息，為后續(xù)的數(shù)字圖像處理與三維圖像重建提供了可靠的數(shù)據(jù)與信息基礎。

　　圖1 葉輪在CT中的放置圖

　　圖2 葉輪原始CT圖片

　　由于在CT斷層圖像獲取過程中，CT機硬件系統(tǒng)、外界干擾以及傳輸過程中的噪聲等都會對圖像的質(zhì)量產(chǎn)生影響，使圖像出現(xiàn)不同程度的質(zhì)量下降，必須進行濾波去噪;另一方面，為了消除因外界干擾以及斷層圖像重建過程中各種客觀因素的影響可能出現(xiàn)的偽跡現(xiàn)象，也很有必要通過灰度變換等措施使斷層圖像得到增強。

　　2.1 濾波去噪

　　一般在空間域內(nèi)采用鄰域平均即局部平均濾波來減少噪聲，而在頻率域內(nèi)，因噪聲頻譜多在高頻段，故一般采用中值濾波和高斯濾波等低通濾波的辦法。

　　中值濾波是一種非線性的信號處理方法，可以克服線性濾波器如均值濾波等帶來的圖像細節(jié)模糊，對濾除脈沖干擾及圖像掃描噪聲很有效。根據(jù)課題特點，采用中值濾波。

　　使用3×3模板分別對斷層數(shù)據(jù)第20、25斷層進行中值濾波，濾波后的圖像如圖3所示。可以看出，噪聲己基本消除，同時邊界也很好地得到保護。

　　圖3 中值濾波后CT圖片

　　2.2 灰度變換

　　灰度變換的目的是增強CT斷層圖像的對比度，擴大圖像灰度范圍，從而達到增強圖像的目的。為此采用局部分段線性變化法來實現(xiàn)灰度變化，該算法的數(shù)學模型為：

　　圖4是處理后的第20、25斷層圖像。需要注意的是，在處理時還可以通過小心調(diào)整折線拐點的位置并控制分段直線的斜率來動態(tài)地對任一灰度區(qū)間進行擴展或壓縮，直至達到滿意的效果。

　　圖4 灰度變換后得到增強的CT圖片

　　3 二維圖像處理

　　二維圖像處理主要是對前面處理過的圖像進行邊緣檢測并跟蹤、提取其輪廓，為二維幾何處理提供數(shù)據(jù)準備。邊緣檢測的方法很多，而且新算法還在不斷地出現(xiàn)。根據(jù)需要采用基于梯度的檢測算子來進行邊緣檢測，用八鄰域輪廓跟蹤來實現(xiàn)輪廓提取。

　　3.1 邊緣檢測

　　基于梯度的檢測算子有很多，本文采用Robert算子對灰度變換后的第20、25斷層圖像進行邊緣檢測。Robert邊緣檢測算子是根據(jù)任意一對互相垂直方向上的差分可用來計算梯度的原理，采用對角線方向相鄰像素之差。該算子的數(shù)學模型為：

　　圖5是處理后的圖像，可以看出：葉片的邊緣已比較完整、清晰地被檢測出來，說明使用Robert邊緣檢測算子的定位精度高，效果較好。

　　圖5 邊緣檢測后的CT圖片

　　3.2 輪廓跟蹤與提取

　　邊緣檢測后得到的邊緣輪廓是離散、獨立、互不相關的各邊界輪廓點的無序組合，需要利用輪廓跟蹤算法，沿著圖像的輪廓進行搜索，將搜索到的輪廓上點的坐標記錄在點列中存儲，供以后使用，其結果是，每一條輪廓線都表示成一個點列。

　　本文采用八鄰域輪廓跟蹤方法，在采用輪廓跟蹤算法時，邊緣灰度圖中背景像素顏色為白色，邊緣上的像素為黑色。

　　結果如圖6所示。可以看出，葉輪斷面的輪廓是封閉的，同一切片上的多輪廓也互不相交，實現(xiàn)輪廓跟蹤后的每一條輪廓線是單像素寬的互不相交的曲線。同時也得到了每個輪廓的輪廓點序列，為矢量化處理打下了良好的基礎。

　　圖6 提取輪廓后的CT圖片

　　4 二維幾何處理

　　經(jīng)上述處理后，得到的是葉輪斷層上每一個輪廓的離散點序列，由于存在大量的冗余信息，不能直觀地描述斷層上輪廓的幾何拓撲結構，在商用CAD或自編軟件中進行編輯處理時的工作量很大，必須進一步用直線、圓弧來擬合這些用離散點序列表達的輪廓，實現(xiàn)用幾何的方式來描述這些輪廓，即對輪廓進行矢量化處理，通過坐標變換將以像素為單位的矢量數(shù)據(jù)變換為笛卡爾坐標系等以國際長度單位為單位的矢量數(shù)據(jù)。

　　采用直線段和圓弧對提取的輪廓進行擬合實現(xiàn)了二維輪廓的矢量化。基本步驟是：對跟蹤得到的離散輪廓點序列，先提取輪廓上的特征點，利用特征點作為分界點將平面輪廓劃分為若干分段，并判斷各線段的曲線類型，再應用分段一合并算法實現(xiàn)輪廓的直線段擬合，最后應用垂直平分線求交法實現(xiàn)輪廓的直線圓弧擬合。

　　坐標變換的變換公式為：

　　式中(Xi，yi)為圖像像素位置，(a，b)為笛卡爾坐標的原點相對于搜索起點的坐標，(Xi，Yi)為圖像像素位置在設定的笛卡爾坐標系下的坐標，坐標的單位均為像素;S為系統(tǒng)標定參數(shù)，其含義為圖像中一個像素長度代表的實際長度為Smm。

　　5 三維模型重建

　　三維重建主要有3種途徑：基于輪廓的表面重建、基于體素的表面重建以及基于體素的直接體視法。前兩種重建方法都屬于表面重建，其中基于輪廓表面重建法是從一組平面輪廓重構通過這些輪廓的表面;基于體素的表面重建法首先要確定物體表面在每個小體素內(nèi)的小面片，然后將這些小面片連接起來構成物體的表面。這兩種算法可采用現(xiàn)有的圖形學顯示方法。直接體視法直接對體數(shù)據(jù)進行顯示，而不是構造物體表面的幾何描述，可顯示體數(shù)據(jù)中所包含的所有物體，但由于不包含物體的幾何拓撲結構關系，因而不能很好地表現(xiàn)空間層次。本文先用第一種重建算法實現(xiàn)了水泵葉輪的三維重建，然后再通過一系列算法將該模型轉(zhuǎn)換成實體模型。

　　5.1 表面三維重建

　　首先在前面已提取并矢量化的輪廓線上尋找特征點，試圖通過這些點的連線來近似地表達輪廓線

　　圖7 輪廓點所圍面積計算

　　由于特征點數(shù)目遠小于矢量化的輪廓線上的點數(shù)，故構造三角面片的工作量大為減小。本文采用了一種基于Delaunay三角剖分思想的表面重建的算法，算法考慮了組成最佳重建表面的三角片的形態(tài)特點，按照剖分三角片集中的最小內(nèi)角在整體上盡可能大的優(yōu)化準則進行剖分，算法簡單，重建表面質(zhì)量較好。

　　5.2 表面模型轉(zhuǎn)換為實體模型

　　由于實體模型的任意曲匾可用多邊形網(wǎng)格來逼近，而邊界模型(B-Rep模型)雖然是用實體的表面來表示實體的形狀，但在具體處理時是將實體的表面形狀表示為多面體，即用一組三維多邊形所圍成的形體來表示實體。所以可以在前面已經(jīng)構造出的由三角面片形成的表面模型的基礎上用邊界模型來構造和表達重建的三維實體模型。重建后的葉輪三維模型如圖8所示。可看出，該模型在基本形狀、結構葉輪三維模型尺寸與葉輪原型基本一樣，但是在光順程度上仍有欠缺，需要在今后對算法與重建過程繼續(xù)加以改進。

　　圖8 重建后的水泵葉輪三維模型

　　6 結束語

　　利用多層螺旋CT對水泵葉輪進行透視掃描，采集得到了含有葉輪設計及內(nèi)部構造信息的CT斷層圖像;在VC++6.0開發(fā)平臺上，對這些斷層圖像進行一系列數(shù)字圖像處理，實現(xiàn)了水泵葉輪的三維表面模型重建以及實體模型的轉(zhuǎn)換，為水泵葉輪的逆向工程提供了一種新的思路。