付堉家 張津豪
摘要 CT (Computed Tomography)是一種可以在不破壞樣品的情況下獲取樣品的內(nèi)部結(jié)構(gòu)信息的一種技術(shù),目前在醫(yī)學(xué)、考古學(xué)等領(lǐng)域有廣泛應(yīng)用。我們根據(jù)已知的吸收率反推幾何模型時(shí),利用MATLAB和radon變換繪制出圖像,然后建立直角坐標(biāo)系,將這兩個(gè)圖像相結(jié)合進(jìn)行幾何分析,建立模型I——幾何分析模型。在射線垂直橢圓長軸入射、垂直短軸入射和射線穿過小球這三種情況下,分別計(jì)算探測器單元的間距,并對這些值取平均值得出探測器單元的間距為。利用托盤的中心和旋轉(zhuǎn)中心在X軸和Y軸方向的偏差求出旋轉(zhuǎn)中心的位置為。隨后尋找出圖像最寬和最窄處的角度為-90度和0度,進(jìn)行逐步類推得出射線的轉(zhuǎn)動(dòng)范圍。在已知未知介質(zhì)的吸收信息的情況下,我們采用“濾波反投影法”建立模型Ⅱ
濾波反投影模型,并對濾波后的投影值進(jìn)行內(nèi)插,進(jìn)行介質(zhì)幾何形狀的確定。通過畫出一個(gè)100*100的方框(正方形托盤)進(jìn)行介質(zhì)位置的確定。最后利用反投影得出十個(gè)位置處的吸收率。根據(jù)結(jié)果,由于物體內(nèi)部含有大量孔洞,會產(chǎn)生大量吸收率過低的薄層。薄層的存在會使得投影圖上少部分像素點(diǎn)缺失,從而導(dǎo)致反投影重建得到的原圖產(chǎn)生丟失微小結(jié)構(gòu)等后果。對于模型的改進(jìn)問題,本文著重考慮影響CT系統(tǒng)參數(shù)標(biāo)定的精度和穩(wěn)定性的兩個(gè)因素:橢圓長短軸測量次數(shù)少影響精度;噪聲和射束硬化效應(yīng)影響穩(wěn)定性。因此,本文設(shè)計(jì)了一個(gè)由四個(gè)等大等間距橢圓和一個(gè)小圓組成的全新的標(biāo)定模板以提高系統(tǒng)標(biāo)定的精度和穩(wěn)定性。
[關(guān)鍵詞]CT系統(tǒng) 濾波反投影 靈敏度分析radon變換 圖像處理
1 CT系統(tǒng)工作原理
X射線垂直于探測器平面平行入射,每個(gè)探測器單元看成一個(gè)接收點(diǎn),且等距排列。X射線的發(fā)射器和探測器相對位置固定不變,整個(gè)發(fā)射.接收系統(tǒng)繞某固定的旋轉(zhuǎn)中心旋轉(zhuǎn)若干次。探測器上的等距單元測量經(jīng)位置固定不動(dòng)的二維待檢測介質(zhì)吸收衰減后的射線能量,并經(jīng)過增益等處理后得到若干組接收信息并對物體進(jìn)行分析。
2 模型的假設(shè)
(1)射線在空氣中傳播過程中不會衰減;
(2)在傳播過程中不考慮干涉和衍射等光學(xué)問題;
(3)探測器旋轉(zhuǎn)中心在旋轉(zhuǎn)過程中不會發(fā)生偏移;
(4)光源不具有時(shí)間非均勻性和空間非均勻性。
3 模型的建立與求解
3.1 問題一的分析與求解
3.1.1 對問題的分析
我們利用MATLAB將附件一中的相關(guān)數(shù)據(jù)繪圖,對數(shù)據(jù)利用MATLAB進(jìn)行random變換處理,實(shí)現(xiàn)圖像的重建和復(fù)原。
3.1.2 模型1:幾何分析模型
模型的準(zhǔn)備:
(1) CT成像原理的進(jìn)一步分析。
有光學(xué)性質(zhì)可知,當(dāng)X射線穿過均勻材料的物質(zhì)時(shí),其強(qiáng)度的衰減率與強(qiáng)度本身成正比,即有:
射線穿過某物體時(shí),(1)式中的μ并非常數(shù),而是一個(gè)關(guān)于坐標(biāo)平面(x,y)的函數(shù)μ(x,y),當(dāng)射線沿平面內(nèi)某直線L穿行時(shí),當(dāng)射線沿平面內(nèi)某直線L穿行時(shí),(2)式變?yōu)椋?/p>
(2)利用Radon變換進(jìn)行投影重建。
一個(gè)己知平面內(nèi)沿不同直線對f(x,y)進(jìn)行線積分,得到的像F(d,a)就是函數(shù)f的Radon變換。
一個(gè)己知平面內(nèi)沿不同直線對f(x,y)進(jìn)行線積分,得到的像F(d,a)就是函數(shù)f的Radon變換。
即:
平行束重建采用的是平移加旋轉(zhuǎn)的掃描方式,如圖1所示,射線源在某一角度下水平移動(dòng),將物體全部照射后旋轉(zhuǎn)一角度,如此重復(fù),在這個(gè)過程中探測器相應(yīng)地運(yùn)動(dòng)以接收X射線。
平面中某一點(diǎn)密度可以看作是這一平面里所有經(jīng)過該點(diǎn)射線的投影值之和。整幅重建圖像可以看成是所有方向下的投影疊加而成。
因此有:
該式可當(dāng)做反投影重建算法的計(jì)算。其中Xk表示像素k的值,Pk,表示經(jīng)過像素k的第i條射線投影,np表示圖像內(nèi)的射線條數(shù)。
3.2 問題一的求解
以標(biāo)定模板中的橢圓中心作原點(diǎn),水平方向?yàn)閤軸,豎直方向?yàn)閥軸建系。
3.2.1 探測器單元之間的距離求解
(1)顯然,光源和探測器繞著旋轉(zhuǎn)中心旋轉(zhuǎn),當(dāng)射線正好垂直標(biāo)定模板的長軸入射時(shí),此時(shí)被介質(zhì)吸收的光源寬度最大,即圖像的高度最大。
利用MATLAB對180個(gè)角度中不為O的個(gè)數(shù)N(吸收系數(shù)不為O,即穿過介質(zhì))進(jìn)行計(jì)數(shù),發(fā)現(xiàn)個(gè)數(shù)N最多為298或288(由系統(tǒng)誤差引起),而橢圓實(shí)際長軸為L=80,所以探測器相鄰單元之間的間距為:
(2)同理,當(dāng)射線垂直于橢圓的短軸入射時(shí),此時(shí)被介質(zhì)吸收的光線的寬度最小,反映在圖1中為180個(gè)角度中不為O的個(gè)數(shù)N(吸收系數(shù)不為O,即穿過介質(zhì))最小。
(3)小圓球在投影圖中的體現(xiàn)為一條等距窄帶,求出其對應(yīng)的接收點(diǎn)個(gè)數(shù)N。
上述三種情況可得探測器單元間距。
3.2.2 CT系統(tǒng)旋轉(zhuǎn)中心在托盤中的位置
我們知道CT旋轉(zhuǎn)系統(tǒng)有其固定的旋轉(zhuǎn)中心,以下稱之為COR(Center of rotation)。理想狀態(tài)是COR與托盤幾何中心重合,在此情況下掃描可達(dá)最佳狀態(tài)。而實(shí)際安裝情況下往往因?yàn)榉N種外部因素,不能做到完全重合,因而存在一定的偏差。此偏差是影響成像質(zhì)量的重要參數(shù)。
(1)首先我們考慮oo的情況。因?yàn)樘綔y器有512個(gè)探測單元,因此該系統(tǒng)的旋轉(zhuǎn)中心位于第256和257個(gè)探測器單元中間的位置。
(2) 900的情況。同理旋轉(zhuǎn)中心縱坐標(biāo)Y1與托盤中心縱坐標(biāo)Yo的差值即為圖中所示YC。
點(diǎn)(xc,YC)即為所求CT系統(tǒng)旋轉(zhuǎn)中心與托盤中心的偏差。
由于探測器單元共有512個(gè),其中間單元為256和257?,F(xiàn)分別根據(jù)這兩個(gè)進(jìn)行計(jì)算:
3.2.3 X射線的180個(gè)方向
己知投影最窄出為0°方向,最寬處為-90°方向。其余的列數(shù)方案參照這兩個(gè)度數(shù)一次類推。
值得注意的是,這里-90°和0°之間的間隔并非剛好90,而是92個(gè),這也是因?yàn)橄到y(tǒng)存在一定的誤差所致,詳細(xì)的分析將在第4問的解答中體現(xiàn)。
3.3 問題二的分析與求解
3.3.1 模型的準(zhǔn)備
(1)反投影重建。平面中某一點(diǎn)的密度可看作這一平面內(nèi)所有經(jīng)過該點(diǎn)的射線投影之和的均值。
即“反投影重建”計(jì)算式為:
Xk表示像素k的值,pk,i表示經(jīng)過像素k的第j條射線投影,np表示圖像內(nèi)的射線條數(shù)。
(2)偽跡出現(xiàn)的原因。反投影的本質(zhì)即取自有限物體空間的投影均勻地反投影到射線到達(dá)無限空間各點(diǎn)上,包括原先像素為零的點(diǎn)。
(3)濾波器作用。減弱或消除傅里葉變換的高頻分量,從而減小噪聲。經(jīng)過查閱文獻(xiàn),為較完整保留目標(biāo)圖像的邊緣,我們選用中值濾波進(jìn)行去噪,可以用來減弱隨機(jī)干擾和脈沖干擾。
3.3.2 模型的求解
因此在本題中,為了去除“偽跡”,我們采用“濾波反投影法(FBP)”進(jìn)行本題的求解。
根據(jù)上問,我們己標(biāo)定該CT系統(tǒng)參數(shù),可以為具體測量做準(zhǔn)備。在本問中,我們根據(jù)標(biāo)定好的CT系統(tǒng)來求得物體的位置信息。
對每次測得的投影數(shù)據(jù)先進(jìn)行一維傅立葉變換。據(jù)中心切片定理,將此變換結(jié)果看成二維頻域中同樣角度下過原點(diǎn)的直線上的值:
步驟5重復(fù)180次,完成全部反投影。
由于實(shí)際中物體是連續(xù)的,但是采樣的投影數(shù)據(jù)是離散的,不一定能夠正好落在采樣點(diǎn)位置上,為了得到任意一點(diǎn)的值,需要對濾波后的投影值進(jìn)行內(nèi)插。
再根據(jù)matlab中自帶濾波反投影函數(shù)iradon對投影數(shù)據(jù)進(jìn)行反投影重建,再以旋轉(zhuǎn)中心為正方形幾何中心,畫出一個(gè)100*100的方框。
3.4 問題四的求解與分析
3.4.1 問題的分析
如果計(jì)算每一次探測所得到的吸收率之和,在理想情況下吸收率總和應(yīng)該為一個(gè)定值,但在問題(1)標(biāo)定摸板中,吸收率總和是隨探測次數(shù)而不斷波動(dòng)的,其中波動(dòng)最大的位置在第130次-170次之間,除了噪聲影響之外,這主要是射柬硬化效應(yīng)造成的。CT中的x射線是具有頻譜寬度的x線源,即x射線光子能量不同。在x射線穿過物體后,低能量部分易被吸收,高能部分較易穿過,因此在傳播過程中,射線的平均能量會變高,即射線會“變硬”,稱其為射束硬化效應(yīng)。本題中,在第130~170次測量中,射線所需要透過的物體厚度最大,則此范圍射線束較硬,導(dǎo)致所得吸收率誤差較大,精度和穩(wěn)定度降低。
3.4.2 問題的求解
為解決上述問題,我們設(shè)計(jì)了如圖1所示的標(biāo)定模板。
該模板采用四個(gè)大小相等的橢圓,如圖所示等距擺放在正方形托盤上,留足間距,保證X射線可以完整的照射到每個(gè)橢圓的長短軸。在左下角放置一個(gè)小圓形,用于確投影方向。相比于題目原來給定的一個(gè)橢圓和一個(gè)小圓的模板,該模板可以測量多次橢圓長短軸取均值,從而避免了單次測量的誤差和損失,是標(biāo)定更精確、更穩(wěn)定。
4模型的評價(jià)和改進(jìn)
4.1模型的優(yōu)點(diǎn)
(1)在解決問題時(shí),利用計(jì)算機(jī)軟件進(jìn)行計(jì)算較為方便,求解速度很快,工作效率高。
(2)在模型的推導(dǎo)中沒有十分艱深的數(shù)學(xué)概念。理解方便,應(yīng)用簡單。
(3)模型可廣泛應(yīng)用于醫(yī)學(xué),考古等現(xiàn)實(shí)問題中,與實(shí)際問題緊密聯(lián)系,通用性,推廣性較強(qiáng)。
4.2 模型的缺點(diǎn)
在反Radon變換過程中,對濾波器選取,規(guī)避信號損失及其造成的尾影過程優(yōu)化討論不足。
4.3 模型的推廣
系統(tǒng)優(yōu)化和參數(shù)的標(biāo)定在當(dāng)前具有廣泛的應(yīng)用,因而本文所用的模型建立方法不僅適用于CT系統(tǒng),還適用于其他在醫(yī)療、考古行業(yè)應(yīng)用的大型設(shè)備。我國是一個(gè)大國,對于各種高檔大型設(shè)備均具有很大的需求,所以從這個(gè)角度出發(fā),本文提出的模型具有廣闊的應(yīng)用的前景。