基于光子計數(shù)探測器的CT系統(tǒng)在蒙特卡羅模擬平臺的實現(xiàn)

孫永剛 孔慧華 張海嬌

摘 要：為了解決利用Matlab，C++等軟件對投影圖像、重建算法的模擬無法真實反映CT系統(tǒng)的成像過程的問題，在研究模擬理論的基礎上，首先利用Geant4和Gate數(shù)據(jù)包搭建蒙特卡羅模擬平臺；然后通過光子計數(shù)探測器模擬多能譜CT成像，同時分析在不同標準下的閾值劃分對重建圖像的影響，利用光子計數(shù)探測器同時掃描出多個不同能段下的投影圖像，并對比出各個能段下的投影效果；其次，根據(jù)不同能段的劃分，選取每個能段的360度下的投影數(shù)據(jù)，利用FBP算法重建投影切片圖像；最后，對所得到的切片圖像進行材料區(qū)分，得到硫和鋁的衰減系數(shù)曲線。實驗結(jié)果表明，模擬方法對光子計數(shù)探測器的成像過程做出了模擬，得到了較好的模擬圖像，并給出了相關閾值的劃分依據(jù)，可為其他CT成像研究提供參考。

關鍵詞：圖像處理；多能譜CT；能譜劃分；Geant4；Gate

中圖分類號：TP391.4 文獻標志碼：A

文章編號：1008-1534（2018）05-0335-06

可見光的差別在于波長和頻率，在各個能量下波長和頻率是不同的，因此可以利用這一現(xiàn)象（即光譜探測原理），設計具有X射線能量分辨能力的多能譜CT投影系統(tǒng)[1]。陳平等[2]基于能譜濾波分離獲取多個能譜通道下的投影。楊曉飛[3]設置2個射線源得到2個能譜投影圖像，并對光子計數(shù)型探測器設定閾值獲取多能譜投影。能譜分離一是通過加裝有濾波片的探測器，識別不同能量下的光子，以此實現(xiàn)多能譜的獲取；二是設置多個探測源，可以利用多個入射能量獲得不同能量的投影。當前多數(shù)情況下主要利用雙源雙探測器實現(xiàn)對2個不同能譜的分辨。但上述辦法效率不高，投影過程相對繁瑣；三是先設定好需要的能量閾值，然后探測器可以區(qū)分X射線的光子脈沖，對選定的不同能譜區(qū)域內(nèi)的光子數(shù)進行累計計數(shù)，從而獲取多個能譜通道下的投影。但光子計數(shù)探測器的模擬成本高昂，步驟復雜，因此搭建計算機模擬平臺成為理想的選擇。上述3種能譜獲取方式如圖1所示。

筆者通過搭建蒙特卡羅平臺來模擬光子計數(shù)探測器的成像過程，通過Geant4數(shù)據(jù)包的數(shù)據(jù)算法程序計算探測器接收到的光子個數(shù)，然后通過Gate平臺的模擬程序，設計出光子計數(shù)探測器，在進行實驗后，得到了每個能譜閾值之間的投影測量數(shù)據(jù)；最后利用投影測量數(shù)據(jù)用FBP算法重建出切片圖像[4]，以此分析不同能譜通道的劃分對圖像的影響。

1 光子計數(shù)探測器成像模型

光子計數(shù)探測器能對X射線中的不同能量的光子分別記錄，可獲得多個不同能譜通道下的投影，根據(jù)蘭博比爾（Lambert-Beers）定理，一個閾值S下的光子強度為

2 蒙特卡羅模擬平臺

Geant4可對X射線中不同粒子的運動過程進行仿真模擬，具有開放且可編輯的源代碼，可根據(jù)具體的實驗需求構(gòu)建出CT仿真系統(tǒng)。Gate是用于對CT成像的模擬軟件，可以設置對應的探測器的參數(shù)，比如探測器像素、材料及厚度。同時設置好模體的材料及形狀，最后取得對模體掃描后的投影圖像，進而完美地實現(xiàn)對CT系統(tǒng)的模擬仿真[5]。

2.1 蒙特卡羅模擬平臺的搭建

平臺主要利用Geant4和Gate數(shù)據(jù)包進行模擬實驗，2種數(shù)據(jù)包需要在Linux系統(tǒng)上運行。為了便于仿真模擬，筆者在Windows系統(tǒng)環(huán)境下，通過Red Hat公司開發(fā)的Cygwin虛擬鏡像軟件，配置運行環(huán)境，搭建出模擬平臺的外部框架。安裝VC++2010進行環(huán)境變量的設置，并在安裝完成后運行實例程序進行檢測。

至此平臺的搭建工作基本完成，在進行模擬實驗的時候，通過Gate調(diào)用Geant4數(shù)據(jù)包，實現(xiàn)對CT成像過程的模擬。在Gate中，可以分別設計CT系統(tǒng)的點源、模體及探測器，同時可以根據(jù)具體實驗需求來設定所需的入射能量、物理過程及數(shù)據(jù)輸出的種類。

2.2 平臺對X射線散射過程的模擬

在一個CT系統(tǒng)里，首先要設置相應的X射線源，在蒙特卡羅模擬平臺里，可以通過點源文件設置入射光子強度、入射角度以及點源數(shù)量。為了提高模擬效率，可以設置一個錐束X射線源，射線源發(fā)出的X射線到達模體后，會產(chǎn)生一系列物理現(xiàn)象，例如光電效應、韌致輻射、康普頓散射等。

在模擬康普頓散射及韌致輻射現(xiàn)象時，Geant4數(shù)據(jù)包具有多種散射模型，例如所包含的G4-KleinNishinaCompton，G4-KleinNishinaMode模型。第1種模型的模擬機制不同，其模擬速度較快，第2種模型中考慮了粒子的各種物理現(xiàn)象，包括原子殼效應、光電效應，對比之前的方法更為客觀，但模擬速度也隨之變慢，因此筆者采用第2種散射模型。所模擬的散射過程符合“克來茵-仁科”截面公式：

平臺中的韌致輻射模型可以對低能時的入射截面進行計算，對應的X射線成像過程也可獲得準確的模擬。每次仿真模擬時，首先要確定每個光子的入射能量，入射角度及光子個數(shù)，然后對光子的所有物理過程進行模擬，記錄下光子的運動信息。當入射的光子產(chǎn)生物理反應時，根據(jù)康普頓、韌致輻射截面公式，可以計算出光子的隨機概率分布，同時根據(jù)一定規(guī)則對入射光子進行篩選，從而得到了到達探測器的光子個數(shù)N及光子能量E，最后由此生成投影的測量數(shù)據(jù)。

2.3 掃描模體的設置

筆者選取一個包含了4個不同半徑且不同材料的子圓柱體的有機圓柱體，4個子圓柱體的材料分別為金屬物質(zhì)鋁和鐵，非金屬物質(zhì)碳和硫，外部的母圓柱體材料選取人體骨骼有機化合物。其橫截面切片圖如圖2所示。模體中的具體參數(shù)如表1所示。

2.4 光子計數(shù)探測器的設計

當前光子計數(shù)探測器選取的晶體材料大多為硅（Si），砷化鎵（GaAs），碲化鎘（CdTe）和碲鋅鎘（CZT），不同探測材料的能量分辨率和探測范圍也不盡相同?？紤]材料的性能及成本，設置厚度1.5 mm，大小1.0 mm×1.0 mm的碲鋅鎘（CZT）材料為探測器的晶體單元。這里構(gòu)建的面陣探測器包含256×256個探測單元，設置圓柱體載物臺，系統(tǒng)旋轉(zhuǎn)中心軸為載物臺中心軸。射線源高度與載物臺中心軸處于同一平面，這樣就構(gòu)成一個錐束CT成像系統(tǒng)，見圖3。

3 模擬實驗

3.1 基于光子計數(shù)型的多能譜CT模擬

平臺通過根據(jù)CT系統(tǒng)的各部分來設置對應的模塊文件，以此完成模擬實驗，模塊文件主要有X射線源模塊Source.mac，掃描模體模塊Phantom.mac以及探測器主體模塊CTScanner.mac。

在Source.mac模塊文件中，設置好入射點源位置、入射強度及粒子種類；在Phantom.mac模塊中設置掃描模體的幾何形狀、模體材料和模體坐標；在CTScanner.mac中設置探測器材料及像素大小。在平臺上模擬出CT成像過程，得到對應的投影數(shù)據(jù)文件[7]。

設置好X射線源、掃描模體、光子計數(shù)探測器及輸出文件后，通過SpectrumGUI軟件進行能譜數(shù)值模擬，獲取能譜為0～140 keV、電壓在125 kVP時的離散能譜數(shù)據(jù)，然后加入到source.mac文件中，即可得到模擬所需要的X射線源的能譜分布數(shù)據(jù)文件。首先模擬截取能譜為10～20 keV，25～35 keV，60～90 keV，95～125 keV下的投影過程，每個通道下進行360度的掃描，并得到所有能譜通道的投影數(shù)據(jù)；然后選取角度在137°時的投影數(shù)據(jù)，并生成投影圖像；最后利用Matlab軟件，由FBP算法重建出投影圖上第69層的切片圖像。能譜圖像及對應結(jié)果如圖4和圖5所示。

分析重建好的切片圖像，通過對比可以看出不同能譜通道下的投影數(shù)據(jù)也有不同，根據(jù)徐品等[8]的研究理論，計數(shù)值在最低時的閾值對應較高的入射電壓，所以在高電壓下光子計數(shù)探測器的光子計數(shù)值在較高能量下會相對較少，探測器在高能量下接收到的光子數(shù)減少，使得部分物質(zhì)的重建亮度降低。為進一步驗證實驗結(jié)論，需要做出每個能譜通道下所對應的衰減系數(shù)曲線圖像，可以選取部分材料區(qū)域，以便對衰減系數(shù)與理論數(shù)據(jù)作出對比。其中深色是實驗數(shù)據(jù)，淺色是理論數(shù)據(jù)（該數(shù)據(jù)選取平臺模擬不添加生成噪聲的物理現(xiàn)象，但受其他因素諸如探測器的探測效率及模擬硬件的限制，噪聲無法完全避免，只能選取噪聲最低時的效果做對比），結(jié)果如圖6所示。

可以看出不同的掃描有不同的衰減系數(shù)曲線，在重建出的切片圖像中，圖像顏色越淺，該區(qū)域的衰減曲線越接近0。圖像的噪聲越小，曲線波動也就越小，而密度較高的元素變化越大。低能譜時的光子強度比高能譜小，輪廓分明。對比硫和鋁區(qū)域的衰減系數(shù)曲線圖可以得到，在25～35 keV之間，曲線更為平緩，與理論數(shù)值的擬合較好，因此圖像質(zhì)量更好。但分別在25～35 keV，60～90 keV和95～125 keV之間時，當光子數(shù)減少時，圖像對比度在60～90 keV，95～125 keV之間比在25～35 keV之間要小，衰減曲線較之前波動大，與理論數(shù)值擬合較差，物質(zhì)成分的區(qū)分度降低。

3.2 能譜劃分對投影的影響

在得到多個能譜下的投影后，可以看出，在不同能段下，探測器對光子的吸收效率不同。而在相同的能譜范圍里，無法保證探測器在每一段能譜下也能接受相同數(shù)量的光子數(shù)，這樣會對所生成的投影圖像有一定影響。

為此，設置不同的閾值，但要使探測器能接收數(shù)量大致相同的光子數(shù)；同樣再設置范圍相同的閾值做對比。閾值設定范圍在16 keV左右，投影、重建結(jié)果所對應的衰減系數(shù)如圖7所示[9]。

在探測器接收的光子數(shù)接近時，投影結(jié)果在高能和低能時比較接近，重建出來的切片圖像中，金屬與非金屬物質(zhì)區(qū)別較大，不同材料的投影對比度也有較大的差異，但衰減系數(shù)的曲線更加平穩(wěn)。在選取相同閾值范圍時，低能譜下的投影對比度有較大差異；在高能時，探測器能接收到的光子數(shù)減少，圖像物質(zhì)識別度降低。

因此，只有合理設定能譜劃分閾值，使探測器能在每個能譜通道接受數(shù)量相同的光子數(shù)，才能獲得更好的投影效果[10]。

4 結(jié) 論

通過研究多能光子計數(shù)探測器的成像過程，首先驗證了在蒙特卡羅模擬平臺可以準確模擬CT成像；然后對相應的物理過程、能譜閾值的選擇及成像效果有了更深層次的研究，得到的投影圖像真實可靠，更加容易識別出不同材料的物質(zhì)；最后根據(jù)實驗得到的圖像，分析物質(zhì)的衰減變化，充分驗證了實驗結(jié)果的準確性。

利用計算機平臺模擬CT系統(tǒng)的成像過程，豐富了相關的研究方法，進一步簡化了實驗過程，降低了實驗所需要的人力、物力等成本，可為其他的CT成像研究提供模擬手段及研究基礎。

參考文獻/References：

[1] 莊天戈.CT原理和算法[M].上海：上海交通大學出版社，1992.

[2] 陳平，郭蓉，潘晉孝， 等. 基于Geant4的能譜濾波分離虛擬平臺研究[J]. 中國測試，2016，42（6）：65-69.

CHEN Ping，GUO Rong，PAN Jinxiao，et al. Study on virtual platform based on Geant4 energy spectrum filtering separation[J].China Measurement & Test， 2016，42（6）：65-69.

[3] 楊曉飛. 基于GATE平臺的雙能CT成像系統(tǒng)模擬設計[D].沈陽：東北大學， 2014.

YANG Xiaofei. The Simulation Design of Dual Energy CT Imaging System Based on GATE Platform[D]. Shenyang： Northeastern University， 2014.

[4] 郭蓉，潘晉孝，陳平，等. X射線譜的Geant4模擬及其濾波效應研究[J]. 核電子學與探測技術，2015，35（8）：767-769.

GUO Rong， PAN Jinxiao， CHEN Ping， et al. Geant4 simulation of X-ray spectrum and the study of filtering effect[J]. Nuclear Electronics & Detection Technology，2015，35（8）：767-769.

[5] 廖玉婷，王永波，許可欣，等. 基于光子計數(shù)探測器的X線CT蒙特卡羅仿真[J]. 中國醫(yī)學物理學雜志 2016，33（2）：122-127.

LIAO Yuting， WANG Yongbo， XU Kexin， et al. X-ray computed tomography Monte Carlo simulation based on photon countingdetector[J]. Chinese Journal of Medical Physics，2016，33（2）：122-127.

[6] Geant4 Collaboration.Geant4 User's Guide for Application Developers（Version：geant4 10.0）[R].[S.l.]：CERN，1996：159-163.

[7] BAUMER C， MARTENS G， MENSER B， et al. Testing an energy-dispersive counting-mode detector with hard X-rays from a synchrotron source[J]. IEEE Transactions on Nuclear Science， 2008，55（3）： 1785-1790.

[8] 徐品，陳奭，袁剛，等.基于光子計數(shù)探測器的能譜CT的研究[J].CT理論與應用研究， 2015，24（6）：822-823.

XU Pin， CHEN Shi， YUAN Gang， et al. Study of spectral CT based on a photon counting detector[J]. CT Theory and Applications， 2015，24（6）：822-823.

[9] 徐亞東， 介萬奇， 查剛強， 等. CdZnTe平面探測器對低能X/γ射線的光譜響應[J]. 光學學報， 2009， 29（11）： 3072-3077.

XU Yadong， JIE Wanqi， ZHA Gangqiang， et al. A study on the low energy X/γ-ray spectral resonse of CdZnTe planar detectors[J]. Acta Optic Sinica， 2009， 29（11）： 3072-3077.

[10] 陳平， 潘晉孝， 劉賓. 連續(xù)能譜X-CT投影仿真算法[J]. 無損檢測， 2009， 31（2）： 102-104.

CHEN Ping， PAN Jinxiao， LIU Bin. Simulation arithmetic of X-CT projection based on consecutive spectrum[J]. Nondestructive Testing， 2009， 31（2）： 102-104.