曠 敏，牟廉明，鄭云升，葉 薇，邱利瓊

（1.重慶理工大學(xué) 理學(xué)院，重慶 400054；2.內(nèi)江師范學(xué)院 數(shù)學(xué)與信息科學(xué)院，四川 內(nèi)江 641100；3.西華大學(xué) 理學(xué)院，成都 610039；4.顧縣中學(xué)，四川 廣安 638000）

1 背景

CT技術(shù)可在不破壞物體的前提下，通過(guò)斷層掃描清晰直觀地觀測(cè)到物體的結(jié)構(gòu)，其掃描部分由X射線、介質(zhì)樣品和一定數(shù)量的探測(cè)器組成，是一種利用X光對(duì)物體進(jìn)行掃描來(lái)檢測(cè)目標(biāo)斷層的投影技術(shù)［1－2］。Radon提出投影重建原理后，Hounsfield于1972年在英國(guó)放射學(xué)會(huì)上首次公布第1臺(tái)CT系統(tǒng)掃描裝置，從此CT技術(shù)得到迅猛發(fā)展［3－4］。而CT系統(tǒng)在安裝時(shí)往往存在誤差，導(dǎo)致成像質(zhì)量差，因此提高CT系統(tǒng)幾何參數(shù)的準(zhǔn)確性和精確度的研究顯得尤為重要。

目前，對(duì)CT系統(tǒng)幾何參數(shù)標(biāo)定方法的研究比較活躍。1999年，Gullberg等［5］利用質(zhì)點(diǎn)在探測(cè)器上運(yùn)動(dòng)形成的投影標(biāo)定出橢圓形參數(shù)，從而得到CT系統(tǒng)的幾何參數(shù)。2000年，Noo等［6］利用最小二乘法擬合出標(biāo)定模體的橢圓形投影軌跡，并與系統(tǒng)的幾何參數(shù)組合成關(guān)系式，從而確定CT系統(tǒng)的幾何參數(shù)。2006年，Liu等［7］比較了心圖法、對(duì)角插值法、曲線擬合和幾何法4種方法來(lái)測(cè)定CT系統(tǒng)旋轉(zhuǎn)中心（COR）的精確性。同年，孫怡等［8］選擇基于4個(gè)點(diǎn)狀目標(biāo)的金屬板作為定標(biāo)模型，利用標(biāo)準(zhǔn)模型的某一固定角度下的投影數(shù)據(jù)標(biāo)定幾何參數(shù)。2009年，李保磊等［9］提出利用正弦圖數(shù)據(jù)的冗余性采用總和均分的方式來(lái)確定COR的位置。2013年，孟凡勇等［10］利用中心投影具有π分割數(shù)據(jù)一致性的特點(diǎn)，提出利用原始數(shù)據(jù)精確確定COR的方法。2018年，龍建輝等［11］基于圖像輪廓曲線的提取，通過(guò)建立數(shù)學(xué)模型來(lái)標(biāo)定CT系統(tǒng)參數(shù)。2019年，夏彬珂等［12］和李江珊等［13］采用最小二乘法擬合優(yōu)化函數(shù)，并利用遺傳算法求幾何參數(shù)的值。以上方法側(cè)重提高COR的精確性，十分依靠算法，并在特定條件下才可以求出精確解，條件改變后，方法就不一定適用。據(jù)此，為提高X射線平行入射的二維CT系統(tǒng)幾何參數(shù)的精確性，并適用于不同模板，從多目標(biāo)規(guī)劃角度，提出一種多目標(biāo)優(yōu)化模型。當(dāng)條件增加或減少時(shí)，模型精確度高、適用性強(qiáng)且穩(wěn)定。

2 不同模板的幾何參數(shù)標(biāo)定優(yōu)化問(wèn)題及建模

對(duì)于典型的CT系統(tǒng)，平行入射的X射線垂直于探測(cè)器平面，將每個(gè)探測(cè)器單元看成一個(gè)接收點(diǎn)，且等距排列。X射線的發(fā)射器和探測(cè)器相對(duì)位置固定不變，整個(gè)發(fā)射－接收系統(tǒng)繞某固定的旋轉(zhuǎn)中心逆時(shí)針旋轉(zhuǎn)180次。對(duì)每個(gè)X射線方向，在具有512個(gè)等距單元的探測(cè)器上測(cè)量經(jīng)位置固定不動(dòng)的二維待檢測(cè)介質(zhì)吸收衰減后的射線能量，并經(jīng)過(guò)增益等處理后得到180組接收信息。

2.1 第1類(lèi)模板：正方形托盤(pán)上放置均勻介質(zhì)圓

在100 mm×100 mm正方形托盤(pán)上放置1個(gè)均勻固體介質(zhì)（半徑為25 mm圓）組成的模板，模板的幾何信息如圖1所示。求CT系統(tǒng)的幾何參數(shù)，即X射線入射方向（IDX）、CT系統(tǒng)旋轉(zhuǎn)中心（COR）和探測(cè)器單元間距（DUS）。

圖1 第1類(lèi)模板示意圖

1）確定任意1條X射線方程

將模板的正方形托盤(pán)中心與圓中心重合，建立平面直角坐標(biāo)系。為進(jìn)一步確定COR與512條X射線之間的幾何關(guān)系，任取1個(gè)點(diǎn)作為旋轉(zhuǎn)中心點(diǎn)o′（x0，y0），任取X射線的1個(gè)初始方向?yàn)棣?。由于探測(cè)器單元是等距排列的，當(dāng)平行光入射時(shí)，旋轉(zhuǎn)中心點(diǎn)到每條射線的距離都可以表示，如圖2。假設(shè)X射線從不同方向入射的DUS不同，用Dj來(lái)表示，則任意1條X射線到旋轉(zhuǎn)中心o′的距離…，180。

圖2 X射線任意1個(gè)入射角度示意圖

根據(jù)圖3，作經(jīng)過(guò)點(diǎn)o′（x0，y0）的射線l0，設(shè)射線l0與y軸交點(diǎn)的縱坐標(biāo)為a0，則射線l0的斜截式為y0＝tanθjx0＋a0，即a0＝y(tǒng)0－tanθjx0。任意1條射線在線段l0的下方或上方，可以得到圖3所示的幾何關(guān)系，則構(gòu)成一個(gè)直角三角形△A0A′iAi，可以求得ai。設(shè)點(diǎn)Ai到A0的距離為a′i，則a′i＝因此任意1條X射線與y軸的相交點(diǎn)Ai的縱坐標(biāo)為ai＝a0－a′i。綜上，可以得到任意1條X射線的直線方程是y＝tanθjx＋ai。即：

圖3 任意1條射線與l0兩種位置關(guān)系示意圖

2）建立超定方程組

根據(jù)定義，圓的標(biāo)準(zhǔn)方程為x2＋y2＝252，單位為mm。聯(lián)立射線與圓的方程組：

3）建立第1類(lèi)模板的CT系統(tǒng)參數(shù)標(biāo)定優(yōu)化模型

需要求解該模板下CT系統(tǒng)的幾何參數(shù)，即IDX、COR和DUS。于是確定決策變量為IDX范圍為θj、COR的坐標(biāo)位置為o′（x0，y0）和DUS為Dj。

512條X射線穿過(guò)模板后，會(huì)在探測(cè)器接收器上顯示接收信號(hào)，因此DUS需要滿(mǎn)足一定的范圍，從而使得X射線從不同方向入射時(shí)DUS最小，因此確定目標(biāo)函數(shù)為180組方向下DUS的平均方差和最小，即。為減少旋轉(zhuǎn)中心的偏差，當(dāng)旋轉(zhuǎn)中心越接近坐標(biāo)原點(diǎn)時(shí)，512根X射線在180組方向下更容易檢測(cè)。于是選擇目標(biāo)函數(shù)為旋轉(zhuǎn)中心點(diǎn)o′到坐標(biāo)原點(diǎn)o的距離最小，即

有3類(lèi)條件對(duì)其約束。

約束1：X射線與圓不相交。當(dāng)X射線與物體相交時(shí)，X射線強(qiáng)度會(huì)衰減，因此探測(cè)器上接收的信號(hào)弱。滿(mǎn)足Δ≤0條件，即：

約束2：當(dāng)不同入射方向的X射線發(fā)生強(qiáng)度衰減時(shí)，確定探測(cè)器單元間距范圍。已知射入的X射線有512組，X射線入射角度有180組，因此不同的入射角度會(huì)影響DUS。為精確確定DUS范圍，根據(jù)幾何關(guān)系，用圖4中3種角度入射計(jì)算得DUS范圍為：

圖4 X射線3種入射角度示意圖

約束3：自身范圍限制。在平面上，平行束IDX范圍為一周，即0≤θj≤2π。由于圓在正方形托盤(pán)中，則COR一定在模板上，即

綜上，得到第1類(lèi)模板的CT系統(tǒng)多目標(biāo)參數(shù)優(yōu)化模型，i＝1，2，…，512，j＝1，2，…，180。

2.2 第2類(lèi)模板：正方形托盤(pán)中放置非均勻介質(zhì)橢圓和圓

在100 mm×100 mm正方形托盤(pán)上放置由2個(gè)非均勻固體介質(zhì)（長(zhǎng)半軸為15 mm、短半軸為40 mm的橢圓和半徑為4 mm的圓）組成的模板，模板的幾何信息如圖5所示。其他條件相同，求CT系統(tǒng)的幾何參數(shù)。

圖5 第2類(lèi)模板示意圖

1）確定任意1條X射線方程

在本模板中固體介質(zhì)為非均勻的，且從1個(gè)上升為2個(gè)，即增加了橢圓。為確定模板下的多目標(biāo)優(yōu)化模型，利用相同的建模思路。X射線不論射入橢圓還是圓都會(huì)發(fā)生衰減，根據(jù)圖7的幾何關(guān)系確定任意1條X射線方程過(guò)程相同，則可以確定任意1條X射線方程為y＝tanθjx＋bi，其中2，…，512，j＝1，2，…，180。

2）建立超定方程組

圖6 X射線任意1個(gè)入射角度示意圖

圖7 射線和l0兩種位置關(guān)系示意圖

同樣，利用弦長(zhǎng)公式組成超定方程組，則橢圓和圓的弦長(zhǎng)分別為：

3）建立第2類(lèi)模板的CT系統(tǒng)參數(shù)標(biāo)定優(yōu)化模型

其中：k＝2，3，i＝1，2，…，512，j＝1，2，…，180。

3 適用于不同模板下平行束CT系統(tǒng)參數(shù)標(biāo)定的多目標(biāo)優(yōu)化模型

為提高CT系統(tǒng)幾何參數(shù)的準(zhǔn)確性和精確度，并保證在樣品數(shù)量改變的情況下精確度相對(duì)較高，根據(jù)式（4）和式（8）兩個(gè)模型，在構(gòu)成模板樣品數(shù)量不同、形狀不同、介質(zhì)不同時(shí)，提出平行束CT系統(tǒng)參數(shù)標(biāo)定的多目標(biāo)優(yōu)化模型。則利用射線與模板上樣品之間的幾何關(guān)系，以COR的坐標(biāo)（x0，y0）、探測(cè)器單元間距（DUS）和X射線的入射方向（IDX）為決策變量，以X射線從180個(gè)方向入射時(shí)DUS方差和的平均值最小和旋轉(zhuǎn)中心點(diǎn)o′到坐標(biāo)原點(diǎn)的距離最小為目標(biāo)函數(shù)，以X射線與介質(zhì)樣品不相交、3種方向入射時(shí)（即水平、豎直、東北45°方向）確定IDX范圍和自身限制為約束條件，建立適用于不同模板的平行束CT系統(tǒng)參數(shù)標(biāo)定的多目標(biāo)優(yōu)化模型。

其中：用Dj表示DUS，n表示整個(gè)發(fā)射－接收系統(tǒng)繞旋轉(zhuǎn)中心的次數(shù)，j＝1，2，…，n。探測(cè)器的個(gè)數(shù)為m，i＝1，2，…，m。利用判別式Δrij表示X射線是否射入樣品中，X射線被截?cái)嗟南议L(zhǎng)用lrij表示，r為X射線射入不同樣品的數(shù)量，r＝1，2，3，…。探測(cè)器間距與正方形托盤(pán)和樣品大小有關(guān)，因此最大間距和最小間距分別用c和d表示。COR在正方形托盤(pán)內(nèi)，因此橫縱坐標(biāo)確定。用θj表示IDX，考慮發(fā)射－接收系統(tǒng)按逆時(shí)針旋轉(zhuǎn)是在1周內(nèi)，即θj∈［0，2π］。

4 定量分析

參考文獻(xiàn)［21］的部分?jǐn)?shù)據(jù)，擬合并多次迭代求解第2類(lèi)模板的幾何參數(shù)，得到結(jié)果：CT系統(tǒng)的旋轉(zhuǎn)中心（COR）為（－9.213 4，6.035 1），探測(cè)器單元間距（DUS）為0.276 2 mm，X射線入射初始方向（IDX）為120.210 0°，每次逆時(shí)針旋轉(zhuǎn)約0.997 1°，并得到180組方向。

為驗(yàn)證模型的精度，以第2類(lèi)模板為例，將幾何參數(shù)結(jié)果和9個(gè)樣本［12－20］的初始IDX、COR和DUS與樣本1的數(shù)據(jù)［22］進(jìn)行誤差分析。分別得到10個(gè)樣本初始IDX、COR和DUS的誤差數(shù)據(jù)，具體見(jiàn)表1～3。

表1 IDX誤差分析

表2 COR誤差分析

表3 DUS誤差分析

由表1～3可知：第2類(lèi)模板的初始IDX為0.563 7°＜0.617 7°；COR的誤差距離0.248 5 mm＜0.640 1 mm；DUS的誤差為0.000 6 mm＜0.000 8mm。3部分都低于平均值，且COR的誤差距離僅次于樣本3、4，而本文中DUS的誤差僅為0.000 6 mm，誤差非常小。從整體上看，利用提出的優(yōu)化模型求解出的CT系統(tǒng)的參數(shù)比現(xiàn)有大部分結(jié)果更精確，因此所提出的優(yōu)化模型可行度高，精確度也較高。

5 模型應(yīng)用

為考察模型的實(shí)際使用情況，檢驗(yàn)?zāi)Ｐ偷姆€(wěn)定性，選擇第2類(lèi)模板進(jìn)行應(yīng)用驗(yàn)證。根據(jù)計(jì)算出的幾何參數(shù)和文獻(xiàn)數(shù)據(jù)［21］重建圖像，并與精度較高的樣本3的圖像進(jìn)行對(duì)比，檢查成像的質(zhì)量，進(jìn)而說(shuō)明模型的穩(wěn)定性。

由于模板上的樣品對(duì)X射線有吸收和散射作用，導(dǎo)致X射線發(fā)生衰減。將模板上介質(zhì)樣品的各單位近似均勻處理，得到X射線衰減系數(shù)公式根據(jù)文獻(xiàn)［23］，基于Radon變換和反變換原理，利用Matlab內(nèi)置圖像重構(gòu)函數(shù)iradon，采用濾波反投影算法重建圖像。此時(shí)，圖像需要校正，采用簡(jiǎn)單的平移、旋轉(zhuǎn)，平移前后坐標(biāo)用（x，y）與（x′，y′）表示，在水平和豎直方向的偏移量分別用d x，d y表示，旋轉(zhuǎn)角度用α表示。

由圖8、9可知，校正后的圖像不僅消除了偽影，而且圖像的邊緣更光滑，輪廓更清晰。將本文的校正圖像與樣本3的校正圖像對(duì)比，雖然圖像重建的算法不同、模型不同，但兩者差別甚微。由此可見(jiàn)，所提出的多目標(biāo)優(yōu)化模型重建的圖像能達(dá)到其他方法重建的效果。換言之，說(shuō)明本文中模型在保證成像質(zhì)量的同時(shí)精度高，表明模型的穩(wěn)定性強(qiáng)。

圖8 第1個(gè)圖校正前、本文方法校正后和 樣本3校正后的圖像

圖9 第2個(gè)圖校正前、本文方法校正后和 樣本3校正后的圖像

6 結(jié)論

為提高平行束CT系統(tǒng)幾何參數(shù)的精確性，從多目標(biāo)規(guī)劃角度，提出一種平行束CT系統(tǒng)參數(shù)標(biāo)定的多目標(biāo)優(yōu)化模型。該模型適用于不同數(shù)量、不同形狀、不同介質(zhì)構(gòu)成的模板。利用超定方程組方程多、求解值少且便于迭代優(yōu)化的特點(diǎn)和幾何關(guān)系確定相應(yīng)模板的多目標(biāo)優(yōu)化模型。通過(guò)分析發(fā)現(xiàn)，在不同數(shù)量、不同形狀、不同介質(zhì)構(gòu)成模板的條件下，兩類(lèi)模板的模型結(jié)構(gòu)類(lèi)似，因此總結(jié)提出了一種平行束CT系統(tǒng)參數(shù)標(biāo)定的多目標(biāo)優(yōu)化模型。通過(guò)對(duì)第2類(lèi)模板的計(jì)算，并將9組樣本分別與樣本1數(shù)據(jù)進(jìn)行定量分析發(fā)現(xiàn)，本文結(jié)果優(yōu)于大部分樣本，精度較高。另外，將第2類(lèi)模板的幾何參數(shù)結(jié)合相關(guān)數(shù)據(jù)重構(gòu)圖像，發(fā)現(xiàn)相較于樣本3所重構(gòu)的圖像差別不大。因此，本文中所提出的模型不但適用于該類(lèi)問(wèn)題，還能保證精度較高、穩(wěn)定性較強(qiáng)。

由于未對(duì)圖像重建算法做深入研究，在多目標(biāo)規(guī)劃領(lǐng)域，第2類(lèi)模板的計(jì)算結(jié)果是全局最優(yōu)解，而實(shí)際問(wèn)題中可能不存在這樣的解。因此，需要進(jìn)一步從多目標(biāo)規(guī)劃問(wèn)題出發(fā)求有效解，確保最優(yōu)解在實(shí)際中一定存在。