張益海，張 催，潘小東，商宏杰，駱巖紅，馬思漢，李公平

(蘭州大學 核科學與技術(shù)學院，蘭州 730000)

錐束工業(yè)CT射束硬化校正方法

張益海，張 催，潘小東，商宏杰，駱巖紅，馬思漢，李公平

(蘭州大學 核科學與技術(shù)學院，蘭州 730000)

錐束X射線CT系統(tǒng)中，X射線束具有寬能譜的特性，會導致X射線硬化效應的產(chǎn)生，從而使得重建圖像中出現(xiàn)杯狀和條狀偽影等，故有必要研究硬化校正方法以減少硬化效應帶來的影響。提出了基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正方法，基于試驗測量X射線初始能譜，利用蒙特卡羅模擬得到透射厚度與多能投影的關(guān)系，再利用多項式擬合得到多能投影與等效單能投影的轉(zhuǎn)換關(guān)系，從而實現(xiàn)硬化校正。分別采用該方法和傳統(tǒng)的濾波片法對一鋁制標準件進行射束硬化校正試驗。結(jié)果表明，所提出的方法比濾波片法能更好地消除錐束CT的硬化偽影。

錐束CT；硬化校正；多項式擬合法；濾波片法

工業(yè)CT(Industrial Computed Tomography，ICT) 能夠獲取被檢測物體內(nèi)部斷層的圖像信息，具有其他無損檢測方法無法替代的特點，在工業(yè)無損檢測領(lǐng)域具有廣闊的應用前景[1-6]。工業(yè)CT通常使用X射線管作為X射線源，其出射的X射線具有寬能譜的特性。由于低能光子比高能光子在透射物體過程中的衰減更多，透射物體后的X射線能譜發(fā)生變化，X射線譜峰向高能方向偏移，高能光子所占比例增大，射線變得更難衰減，探測器所探測到的高能光子所占比例增大，從而導致投影值偏小。在CT圖像重建常用的算法中，如濾波反投影算法(FBP)、代數(shù)迭代算法(ART)等，都是基于X射線能量單一的假設[7]。而實際上錐束工業(yè)CT的X射線是多色的，因此只能得到多能投影Pm(Multi-Energy Projection)數(shù)據(jù)，若直接由Pm數(shù)據(jù)進行重建圖像，射束硬化效應會導致本來密度均勻的被測物體在CT圖像上表現(xiàn)為亮度不一，灰度值不相等，通常直觀顯示為杯狀或條狀偽影等[8]。硬化偽影會降低CT圖像的分辨率，影響檢測的精度，因此有必要對其進行校正。多項式擬合硬化校正方法是CT系統(tǒng)中常用的校正方法，該方法是利用多項式擬合得到Pm和透射厚度L的關(guān)系曲線，再取曲線原點處的切線作為等效單能投影Ps(Single-Energy Projection)與L的關(guān)系曲線，最后通過多項式擬合得到Pm與Ps的轉(zhuǎn)換關(guān)系式。該方法的主要缺點是需要根據(jù)要求專門制備特定的吸收片，但吸收片的厚度很難做到足夠小，這就使得在原點處的切線變得不夠精確，并且對于由多種材料混合組成的樣品，很難找到與其材料相同的吸收片。因此，筆者采用蒙特卡羅方法模擬射線穿過吸收片的試驗，得到相應結(jié)果，相較于傳統(tǒng)的試驗測量，其優(yōu)點主要在于可以方便地定義吸收片的厚度和材料，并且減少了大量的測量時間和花費。但由于在利用蒙特卡羅方法模擬前沒有得到X射線管的初始能譜而使得結(jié)果誤差較大，因此，筆者基于前期對X射線管初始能譜測量的工作，提出了基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正方法；并且將該方法應用于一鋁制標準件開展試驗，并和傳統(tǒng)的濾波片法[9]做對比，試驗結(jié)果表明，相較于濾波片方法，提出的方法能給出更好的校正效果，且對射線利用率更高。

1 射線硬化校正原理分析

對于單能X射線，并且樣品材料均勻的情況，X射線的衰減量可由Lambert-Beer[10]定律來描述,見式(1)。

式中:I0為X射線的初始強度;I為X射線穿過樣品后的透射強度;L為射線透射過樣品的厚度，簡稱透射厚度;μ為樣品的線性吸收系數(shù)。

由式(1)變換得到ln(I0/I)=μL，通常令Ps=ln(I0/I)。則Beer定律可寫成：

圖1給出了單能投影Ps和多能投影Pm的投影值與透射厚度L的關(guān)系示意。由圖1可看出：對于單能X射線，投影值與L呈線性關(guān)系；對于多能X射線，投影值與L呈非線性關(guān)系，投影值隨著L的增加，增大的趨勢逐漸趨于平緩。這是由于出射X射線光強不但與L有關(guān)，還與射線的能量相關(guān)，如式(3)所示。

式中:E為射線的能量；S(E)為射線的能譜分布；μ(E)為樣品對應能譜的線性衰減系數(shù)。

則Beer定律可寫成多能投影Beer定律形式，如式(4)所示。

式中：Pm為過原點的曲線(見圖1)。

由于X射線管發(fā)出的X射線具有連續(xù)能譜分布，當X射線透射樣品時，射線束中低能光子比高能光子吸收得更快而導致射線硬化效應，使得X射線吸收不均勻，產(chǎn)生高頻信號；如將此線衰減系數(shù)μ(E)[11]代入濾波反投影重建算法經(jīng)過計算后，會產(chǎn)生條狀或杯狀偽影，這就降低了CT圖像的質(zhì)量，影響檢測的精度。因此，必須進行相應的校正：一方面可以通過找到Pm與Ps的轉(zhuǎn)換關(guān)系以實現(xiàn)硬化校正；另一方面可以降低X射線的多色性，使其近似為“單色光”以減小硬化效應。

圖1 多能投影、單能投影的投影值與透射厚度的關(guān)系

2 試驗過程

2.1 基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正方法

2.1.1 X射線管初始能譜測量

在用蒙特卡羅方法進行模擬試驗前，需要得到X射線管出射X射線能譜，而這可通過試驗測量或者計算機模擬得到。為了得到更準確的結(jié)果，筆者采用試驗測量的方法，采用瑞士COMET公司生產(chǎn)的MXR225/22 X射線管進行。主要試驗參數(shù)為：管電壓100 kV，管電流0.2 mA，小焦點d為1.0 mm，靶材為鎢，靶角20°，輻射范圍為40°，固有濾波為0.8 mm鈹窗，出射口有2 mm鋁濾波片。利用美國ORTEC公司GLP-10180/07P4型高純鍺探測器測量能譜。試驗結(jié)果如圖2所示，X射線由連續(xù)譜和特征譜構(gòu)成。特征峰主要包括鎢和銀元素產(chǎn)生的特征峰，這是由于靶材為鎢，而銀具有極佳的導熱性，可作為鎢靶與銅襯底的釬焊材料[12]。

圖2 試驗測量X射線管小焦點初始能譜(100 keV)

2.1.2 蒙特卡羅方法的軟件模擬計算

蒙特卡羅方法是一種隨機過程方法，其并不嚴格求解一些物理方程，而是通過模擬單粒子在介質(zhì)中的隨機運動過程，然后將大量粒子的平均結(jié)果作為系統(tǒng)結(jié)果的近似，給出粒子系統(tǒng)在介質(zhì)中運輸?shù)淖詈蠼Y(jié)果。筆者采用MCNP-4C軟件工具包建立X射線與不同厚度吸收片作用的模擬試驗模型(見圖3)。X射線源為點源，所出射光子能量服從如圖2所示的能譜分布。吸收片材料為Al，其厚度從0.1 mm增加到20 mm，每次增加的厚度為0.1 mm，共200組模擬試驗。探測器置于樣品后方，三者的幾何位置根據(jù)試驗測量所得，其中源-片距離為51 mm，源-探距離為815 mm。由于空氣對射線的吸收可忽略不計，為簡單起見，模擬環(huán)境設置為真空。

圖3 模擬試驗模型

所編寫的MCNP軟件輸入文件代碼，包括曲面塊設計，柵元塊設計和數(shù)據(jù)塊設計，模擬了射線穿過吸收片后射線強度發(fā)生變化的過程。在給定射線初始強度I0的情況下，通過探測器探測射線穿過不同吸收片后的透射強度I，得到200組試驗數(shù)據(jù)，利用投影值計算公式P= ln(I0/I)計算得到Pm，以L為橫坐標，Pm為縱坐標作圖，即可得到兩者的關(guān)系曲線，如圖4所示。

圖4 多能投影值與透射厚度的關(guān)系

2.1.3 校正公式計算

吸收片厚度很薄時，硬化效應很弱，投影值與L近似成線性關(guān)系，因此在圖4曲線原點處作其切線，該切線即代表Ps與L的關(guān)系，結(jié)果如圖5所示。由圖5可看到：Ps與Pm在L較小時，相差很小，但隨著L的增加，兩者的差距逐漸增大。這表明隨著樣品厚度的增大，硬化效應帶來的影響也將逐漸增大。

圖5 多能投影值、單能投影值與透射厚度的關(guān)系

在相同透射厚度下，對應得到200組(Ps,Pm)值，并以多能投影值與單能投影值的關(guān)系Pm為橫坐標，Ps為縱坐標，得到兩者的關(guān)系曲線，如圖6所示。從圖6可看出：Ps隨Pm的增大而增大，且兩者逐漸趨于線性關(guān)系。

圖6 多能投影值與單能投影值的關(guān)系

再對該試驗數(shù)據(jù)點進行多項式擬合。由于4階多項式擬合就能得到很好的擬合效果，由此得到Ps與Pm的轉(zhuǎn)換關(guān)系，如式(5)所示。

成像過程中，只需將試驗采集得到的Pm數(shù)據(jù)先代入式(5)轉(zhuǎn)換為Ps數(shù)據(jù)，再進行圖像重建，即可實現(xiàn)射束硬化校正的目的。

2.2 濾波片硬化校正方法

濾波片硬化校正方法，就是在CT系統(tǒng)的X射線源出束口處加裝濾波片，以對X射線進行預硬化，濾除低能X射線，使盡量多的高能射線通過。圖7為初始透射能譜與濾波后能譜的對比，其中兩個譜的總計數(shù)均為400 000±2 000。筆者選擇0.46 mm銅片作為濾波片，可看出射線的低能部分被很好地過濾掉了，X射線的平均能量得到了增大，從而降低了硬化效應。該方法最大的不足在于，其硬化校正質(zhì)量與濾波片厚度成正相關(guān)，但隨著濾波片的增厚，射線強度損失嚴重，因此必須增加輻照時間或者增大X射線管的管電流以保證探測器擁有足夠的曝光。

圖7 初始能譜與過濾能譜的對比圖

3 校正結(jié)果及對比

3.1 多項式擬合硬化校正結(jié)果 試驗使用的是COMET公司的X射線管，型號為MXR225/22。試驗主要參數(shù)為：管電壓100 kV，管電流0.2 mA，采用小焦點(d=1.0 mm)。使用美國VARIAN公司的非晶硅平板探測器，型號為PaxScan1313DX，該探測器的探測面積為13 cm×13 cm，像元矩陣總數(shù)為1 024×1 024，像元尺寸為127 μm，極限分辨率為3.94 lp·mm-1，動態(tài)范圍為16位。對一鋁制標準件進行錐束X射線CT試驗，得到360個角度下的多能投影數(shù)據(jù)，選擇第512層數(shù)據(jù)進行斷層圖像重建，重建圖像像素為400像素×400像素，采用FBP重建算法。

圖8 多項式擬合硬化校正前后重建切片效果對比

得到基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正前后切片重建效果的對比，如圖8所示。從圖8(a)、(b)可看到經(jīng)過校正以后，投影正弦圖變得更加清晰。從圖8(c)、(d)可看到：校正前，樣品周圍比較模糊，同時樣品中間相對較暗；校正后，樣品更為清晰，也更為均勻；從圖8(e)、(f)中能更加直觀地看出，未校正的重建圖像存在明顯的杯狀偽影且噪聲起伏較大，而經(jīng)過校正后的重建圖像中的杯狀偽影得到了有效消除，切片的灰度值幾乎均勻，同時噪聲起伏減小，重建圖像質(zhì)量得到明顯提高。

圖9 濾波片硬化校正前后重建切片效果對比

3.2 濾波片法硬化校正結(jié)果

采用濾波片硬化校正方法的結(jié)果如圖9所示。試驗條件同3.1節(jié)，此外加裝0.46 mm銅片，管電流0.5 mA。通過校正前后結(jié)果的對比可看出，該校正方法能在一定程度上減少硬化偽影，使圖像質(zhì)量得到改善，但圖像仍存在硬化偽影。

圖10 三維CT圖像硬化校正前后對比

3.3 兩種校正方法對比

由3.1節(jié)和3.2節(jié)的結(jié)果可看到，傳統(tǒng)的濾波片硬化校正方法只能在一定程度上抑制硬化偽影現(xiàn)象，不能完全解決射束硬化帶來的問題，并且由于濾波片對X射線的吸收強烈，探測器曝光減少，因此不得不增大管電流以保證足夠的曝光。筆者提出的基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正方法，能在不降低射線利用率的情況下，得到更好的硬化校正結(jié)果。

此外，對于由較高原子序數(shù)組成的樣品，例如鐵和銅樣品，濾波片硬化校正方法開始不太適用，而基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正方法仍有較好的效果，圖10給出了將該方法應用于錐束工業(yè)CT系統(tǒng)三維CT圖像的結(jié)果，樣品為鋼卷尺。由圖10可看到，經(jīng)過校正后的三維CT圖像在細節(jié)處更為清晰，偽影更少。試驗條件同3.1節(jié)，對于鐵樣品的硬化校正，Ps與Pm的轉(zhuǎn)換關(guān)系如式(6)所示。

4 結(jié)論

提出了基于試驗能譜和蒙特卡羅模擬的多項式擬合硬化校正方法，并將其應用于錐束CT系統(tǒng)的硬化校正工作中。結(jié)果表明，與傳統(tǒng)的濾波片硬化校正方法相比，所提出的方法不但有更好的硬化校正效果，還能提高射線的利用率，可成為錐束工業(yè)CT射束硬化校正的一種有效方法。由于研究有限，此方法暫時僅應用于單材料樣品的硬化校正，還未應用于多材料混合樣品的硬化校正工作，后續(xù)將做進一步的工作。

[1] KASTNER J, HARRER B, REQUENA G, et al. A comparative study of high resolution cone beam X-ray tomography and synchrotron tomography applied to Fe- and Al-alloys[J]. Ndt& E International, 2010, 43(7):599-605.

[2] AGRAWAL A K, SARKAR P S, KASHYAP Y S, et al. Application of X-Ray CT for non-destructive characterization of graphite fuel-tube[J]. Journal of Nondestructive Evaluation, 2016, 35(3):1-8.

[3] LIU T. Cone-beam CT reconstruction for planarobject[J]. Ndt& E International, 2012, 45(1):9-15.

[4] IZUMI S,KAMATA S, SATOH K, et al. High energy X-ray computed tomography for industrial applications[J]. IEEE Transactions on Nuclear Science, 1993, 40(2):158-161.

[5] 張豐收,張定華,鄧效忠,等. 工業(yè)CT射束硬化研究[J]. 機械強度,2007,29(6):1022-1025.

[6] 吳志宏, 叢鵬, 劉錫明. 基于重投影的CT圖像硬化偽影校正[J]. 原子能科學技術(shù), 2015, 49(5):935-938.

[7] JIANG H. Computed tomography principles, design, artifacts, and recent advances[M]. Bellingham:SPIE Press, 2009.

[8] 曾鋼, 劉力, 闕介民,等. 工業(yè)CT圖像均勻性校正[J]. 高能物理與核物理, 2004, 28(10):1131-1134.

[9] 吳世法. 近代成像技術(shù)與圖像處理[M].北京:國防工業(yè)出版社,1997:143-144.

[10] HANKE R, BOBEL F. Determination of material flaw size by intensity evaluation of polychromatic X-ray transmission[J]. NDT & E International,1992,25(2):87-93.

[11] KROHN B R, SILVER M D. Optimizing contrast resolution in industrial computed tomographic scanners[J]. Materials Evaluation,1990,49(10):1296-1331.

[12] 張催,李公平,潘小東,等. 工業(yè)CT用X射線管射線強度分布研究[J]. 核技術(shù), 2015, 38(9): 8-12.

The Methods for Beam Hardening Correction of Cone-beam Industrial Computed Tomography

ZHANG Yihai, ZHANG Cui, PAN Xiaodong, SHANG Hongjie, LUO Yanhong, MA Sihan, LI Gongping

(School of Nuclear Science and Technology, Lanzhou University, Lanzhou 730000，China)

Since the X-ray beam has the characteristics of wide spectrum in the cone beam X-ray computed tomography imaging system, this leads to the X-ray hardening effect, so as to make the reconstruction image in the goblet artifact and strip artifact.It is necessary to study methods of the hardening correction in order to reduce the effects of hardening effect.A method of beam hardening correction is put forward, which is based on the experimental spectrum and Monte Carlo simulation of polynomial fitting.Based on the initial X-ray energy spectrum by experimental measurement, the relationship between transmission thickness and pluripotent projection can be obtained by exploiting Monte Carlo simulation. And then using polynomial fitting can obtain relationship between pluripotent projection and equivalent single projection, so as to realize hardening correction. This method and the traditional filter method for X-ray hardening correction on the aluminum workpiece were respectively used. The experimental results show that the proposed method is better than the filter method which can decrease cone-beam CT beam hardening artifacts effectively.It provides an effective correction method for cone beam industrial computed tomography beam hardening correction.

cone-beam CT; beam hardening correction; polynomial fitting method; filtering method

2016-10-12

蘭州大學中央高?；究蒲袠I(yè)務費專項資金資助項目(lzujbky-2016-32, lzujbky-2016-208)

張益海(1988-)，男，碩士，主要從事X射線成像研究

李公平， ligp@lzu.edu.cn

10.11973/wsjc201706002

TL99;TG115.28

A

1000-6656(2017)06-0008-05