陳鋒 鄭娜 許海波

1)(中國(guó)工程物理研究院研究生部,北京 100088)

2)(北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所,北京 100094)

提出了一種質(zhì)子能量在中高能時(shí)利用能量損失進(jìn)行密度重建的方法,并利用Bethe-Bolch公式給出了利用能量損失進(jìn)行密度重建的方程及條件.針對(duì)1.6 GeV的質(zhì)子能量,通過(guò)定量計(jì)算常見(jiàn)材料的阻止本領(lǐng),得出質(zhì)子能量在1.45—1.6 GeV范圍內(nèi)時(shí),材料的阻止本領(lǐng)的變化率小于1%,可近似為常數(shù).最后,通過(guò)理論計(jì)算和Geant 4模擬,得出質(zhì)子能量在1.6 GeV時(shí),可以對(duì)面密度為113 g/cm2的縮比法國(guó)實(shí)驗(yàn)客體進(jìn)行密度重建.

1 引 言

1996年,美國(guó)洛斯·阿拉莫斯國(guó)家實(shí)驗(yàn)室(LANL)的科學(xué)家Gavron等提出用質(zhì)子代替X光進(jìn)行輻射照相,為流體動(dòng)力學(xué)實(shí)驗(yàn)提供了一種先進(jìn)的診斷方法.自此,美國(guó)[1,2]和俄羅斯[3]均在高能質(zhì)子照相上投入了大量的人力物力,針對(duì)流體動(dòng)力學(xué)的一些關(guān)鍵科學(xué)問(wèn)題,開(kāi)展相關(guān)理論和實(shí)驗(yàn)研究,結(jié)果論證了高能質(zhì)子照相技術(shù)在致密物質(zhì)內(nèi)部幾何結(jié)構(gòu)和物理特性的診斷方面具有革命性的進(jìn)步.與傳統(tǒng)的X射線照相相比,質(zhì)子照相在穿透能力、精細(xì)結(jié)構(gòu)分辨、物質(zhì)混合診斷、密度和空間分辨、多時(shí)刻和多角度成像等方面有明顯的優(yōu)勢(shì).美國(guó)LANL利用洛斯·阿拉莫斯中子散射中心(LANSCE)800 MeV的質(zhì)子照相裝置對(duì)小尺寸內(nèi)爆做了大量研究,論證了質(zhì)子照相技術(shù)對(duì)動(dòng)態(tài)樣品的連續(xù)成像能力和材料分辨能力;2011年LANL利用布魯克海文實(shí)驗(yàn)室的交變同步質(zhì)子加速器(AGS)對(duì)法國(guó)實(shí)驗(yàn)客體(FTO)進(jìn)行了質(zhì)子照相,其能量為24 GeV,分析得出質(zhì)子照相對(duì)于客體具有1%的密度分辨能力和200μm的空間分辨能力.俄羅斯理論物理所和高能所也分別利用800 MeV的TWAC加速器和70 GeV的質(zhì)子U70加速器裝置,開(kāi)展了多種爆轟實(shí)驗(yàn)和分解實(shí)驗(yàn)[4].德國(guó)重離子研究中心(GSI)利用反質(zhì)子和離子研究裝置(FAIR)進(jìn)行了質(zhì)子照相研究[5].質(zhì)子照相在國(guó)內(nèi)主要處于理論分析和模擬階段[6?10].實(shí)驗(yàn)方面,中國(guó)原子能研究院設(shè)計(jì)了100 MeV的質(zhì)子照相裝置[11].蘭州近代物理研究所在重離子加速器(HIRFL)設(shè)計(jì)了600 MeV/u的碳離子照相裝置[12,13].中國(guó)工程物理研究院流體物理研究所開(kāi)展了低能11 MeV的設(shè)計(jì)研究[14],2017年,在神光II升級(jí)裝置上開(kāi)展了首輪激光加速質(zhì)子對(duì)間接驅(qū)動(dòng)快點(diǎn)火靶內(nèi)爆過(guò)程的照相實(shí)驗(yàn)研究,經(jīng)過(guò)優(yōu)化激光與靶參數(shù),獲得了能量高于18 MeV的質(zhì)子束[15].高能質(zhì)子照相中的成像部件主要是Zumbro磁透鏡[16],它的特點(diǎn)是可以使得客體平面不同位置出發(fā)而散射角度相同的質(zhì)子在中心平面上匯聚,因而可以在中心平面上設(shè)置一個(gè)角度準(zhǔn)直器,將大角度散射的粒子屏蔽掉.LANL科學(xué)家利用核衰減和庫(kù)侖散射兩種作用的信息來(lái)重建圖像.

能量損失是阻止本領(lǐng)在路徑上的積分,能量損失的密度重建只在醫(yī)學(xué)上有廣泛的應(yīng)用[17,18].人體大多數(shù)是水,并且體量元素大多是較輕的元素,醫(yī)學(xué)上用水近似代替其他物質(zhì)組分(比如脂肪、軟骨、肺等).不同物質(zhì)等價(jià)于不同密度的水,所以可以用水的平均電離勢(shì)代替其他物質(zhì),這意味著阻止本領(lǐng)正比于實(shí)際密度,因而可以利用能量損失進(jìn)行密度重建.當(dāng)客體所含物質(zhì)原子序數(shù)變化較大時(shí),除了Z/A值的變化(當(dāng)Z大于29時(shí),隨著原子序數(shù)增大而減小),Z的不同也使得平均電離勢(shì)有很大的不同,不能看作常數(shù),這種情況下,不能直接利用能量損失對(duì)多材料客體進(jìn)行密度重建.對(duì)于中高能的質(zhì)子,材料的阻止本領(lǐng)隨能量的變化比較緩慢,在一定近似條件下可以認(rèn)為阻止本領(lǐng)只是材料密度的函數(shù).可以先通過(guò)能量損失得到阻止本領(lǐng),再?gòu)淖柚贡绢I(lǐng)得到密度,從而達(dá)到診斷密度的目的.本文將對(duì)阻止本領(lǐng)隨能量的變化進(jìn)行研究,在適當(dāng)?shù)哪軈^(qū)范圍內(nèi),通過(guò)阻止本領(lǐng)實(shí)現(xiàn)對(duì)客體的密度重建.

2 基于能量損失的密度重建原理

2.1 理論分析

當(dāng)質(zhì)子入射到靶物質(zhì)時(shí),與靶原子發(fā)生相互作用而損失能量,可以通過(guò)修正后的Bethe-Bloch公式描述[19]:

在非常高的質(zhì)子能量下,密度修正因子可寫為

從(1)式可知,阻止本領(lǐng)與入射質(zhì)子能量有關(guān).為了方便討論,以S表示阻止本領(lǐng),引入?yún)⒘縂=βγ=P/M,則入(1)式可得

以U代替(4)式中除密度以外的多項(xiàng)式,則

因?yàn)槊芏葍H由材料決定,因此與入射質(zhì)子相關(guān)的量是U,對(duì)G求導(dǎo)有

將(2)式代入(6)式,整理并計(jì)算當(dāng)G → ∞時(shí)的極限

因此,對(duì)于G非常大的能量區(qū)間內(nèi),U基本為常數(shù).

圖1 材料的值與G 的關(guān)系 (a)鋁、銅、鉛;(b)鎂、鐵、鎢Fig.1.The relationship betweenand G for di ff erent materials:(a)Al,Cu,Pb;(b)Mg,Fe,W.

表1 質(zhì)子能量分別為1.45 GeV和1.6 GeV時(shí),材料的?U/?G值Table 1.The ?U/?G value of the material,the energy of the proton is 1.45 GeV and 1.6 GeV.

從表1可以看出,?當(dāng)質(zhì)子能?量在1.45—1.6 GeV范圍內(nèi)時(shí),各材料的??U/?G?值在零附近,并隨著原子序數(shù)的增大有遞減的?趨勢(shì);?從表2可以看出,材料的U值相對(duì)變化率??U/U?小于1%,并且隨著原子序數(shù)增大有逐漸減小的趨勢(shì).因而可以認(rèn)為,在此能量區(qū)間內(nèi),各材料的U值隨能量的變化可以忽略,因此,在此能區(qū)內(nèi)?U/?G≈0,即各材料的U值近似為常數(shù),U值由入射質(zhì)子能量決定.

將密度函數(shù)沿每一條射線的積分離散化,由(1)式和(4)式可得多材料客體的能量損失公式為

其中?Ei是第i條射線對(duì)應(yīng)的能量損失值,Umk是體素k中材料m的U值,lik是幾何矩陣的第mk個(gè)元素,它與沿射線i通過(guò)體素k的光程相關(guān).可以通過(guò)(8)式進(jìn)行密度重建.

表2 質(zhì)子能量分別為1.45 GeV和1.6 GeV時(shí),材料的U值Table 2.The U value of the material,the energy of the proton is 1.45 GeV and 1.6 GeV.

2.2 能量歧離對(duì)密度重建的影響

單能的質(zhì)子束通過(guò)靶材后出射質(zhì)子將不是單能的,而是近似為高斯分布,其中能量損失的漲落稱為能量歧離,能量岐離直接導(dǎo)致重建密度出現(xiàn)漲落.能量損失值相對(duì)于初始能量不超過(guò)20%時(shí),能量歧離值由Bohr理論[20]計(jì)算:

其中E(Ein,x)是能量為Ein的入射質(zhì)子通過(guò)靶厚x時(shí)的平均能量值;Z,A,ρ,K和β的定義同(1)式.

本文討論的質(zhì)子能量區(qū)間為1.45—1.6 GeV,質(zhì)子初始能量為1.6 GeV時(shí),最大能量損失相對(duì)于初始值將不超過(guò)10%,因此可以用Bohr理論計(jì)算能量歧離.通過(guò)(8)式和(9)式可以得出能量歧離造成的密度誤差為

圖2是1.6 GeV質(zhì)子穿過(guò)面密度為115.5 g/cm2的鎢靶時(shí)的出射能量分布,得出能量平均值為1.468 GeV,能譜滿足高斯分布,其半高全寬(FWHM)為12 MeV,左右各偏離6 MeV.造成的密度誤差為0.86 g/cm3,相對(duì)誤差4.47%.

圖2 出射質(zhì)子的能量分布Fig.2.The energy distribution of the emitted proton.

2.3 物理建模

以縮比FTO為研究對(duì)象,開(kāi)展密度重建研究.FTO是美、英、法等國(guó)為了研究初級(jí)武器的動(dòng)態(tài)行為而設(shè)計(jì)的與之等效的靜態(tài)客體,最大面密度約200 g/cm2.FTO結(jié)構(gòu)如圖3所示,由一組同心球?qū)咏M成,中心區(qū)是空區(qū).空球半徑為1 cm,第二層和第三層分別是鎢和銅,球半徑分別是4.5和6.5 cm,入射質(zhì)子為1.6 GeV,進(jìn)行密度重建的能區(qū)為1.45—1.0 GeV,因此最大能損值為150 MeV.本文將FTO按比例縮小,使其對(duì)應(yīng)的最大能損值為150 MeV,此時(shí)球半徑從外到內(nèi)依次為4.0,2.7和0.6 cm,對(duì)應(yīng)的面密度為113 g/cm2.

圖3 FTO客體Fig.3.French test object.

2.4 蒙特卡羅模擬

掃描質(zhì)子成像系統(tǒng)可以避免MCS在接受平面處造成的圖像模糊,如圖4所示,質(zhì)子細(xì)束沿著x方向?qū)腕w進(jìn)行等間隔掃描,并且兩個(gè)探測(cè)器分別記錄入射質(zhì)子能量和出射質(zhì)子能量,將客體以等角度間隔進(jìn)行旋轉(zhuǎn),并且重復(fù)此過(guò)程,直到旋轉(zhuǎn)180?或360?,從而可以通過(guò)入射能量值和出射能量值進(jìn)行圖像重建.影響質(zhì)子輻射照相掃描系統(tǒng)的密度分辨率的因素有能量歧離、入射質(zhì)子能量的分散和能量探測(cè)器的噪聲,其中能量歧離是限制密度分辨率的主要因素.

圖4 掃描質(zhì)子成像示意圖Fig.4.The schematic diagram of a proton CT scanner.

本文利用Geant 4[21,22]軟件,基于入射能量為1.6 GeV的質(zhì)子實(shí)現(xiàn)掃描成像.入射質(zhì)子為1.6 GeV的單能質(zhì)子束,質(zhì)子數(shù)為106,掃描間隔為0.5 mm,旋轉(zhuǎn)角度間隔為0.9?(如果是球?qū)ΨQ客體則僅測(cè)1次投影值),探測(cè)器為薄層介質(zhì)板.

3 數(shù)值模擬與密度重建

通過(guò)數(shù)值模擬可以得到能量損失隨位置的分布,從而可以利用(8)式通過(guò)求解線性方程組,重建得到密度分布.圖5(a)和圖6(a)分別是理論計(jì)算和Geant 4軟件模擬的對(duì)稱軸上能量損失隨位置的分布,圖5(b)和圖6(b)分別是采用兩種方法利用(8)式重建的對(duì)稱軸上密度隨位置的分布.由圖5(a)和圖6(a)均可以看出能量損失值在材料界面處有明顯變化,是典型的有利于界面確定和密度重建的分布圖像.

4 結(jié) 論

通過(guò)理論分析得出利用能量損失進(jìn)行密度重建的方程.當(dāng)質(zhì)子能量區(qū)間在1.45—1.6 GeV時(shí),理論計(jì)算和Geant 4模擬結(jié)果表明,材料阻止本領(lǐng)的變化率小于1%,可以認(rèn)為各材料的阻止本領(lǐng)不隨能量變化,僅由入射質(zhì)子能量決定,因此可以直接對(duì)多材料客體進(jìn)行密度重建.醫(yī)學(xué)上主要用不同密度的水近似人體組織,即主要是對(duì)不同密度的水進(jìn)行密度重建,并且醫(yī)用質(zhì)子照相能量約在200 MeV的低能區(qū),本文的研究可以將利用能量損失進(jìn)行密度重建的方法從醫(yī)學(xué)低能區(qū)推廣到GeV量級(jí)的中高能區(qū).

圖5 理論計(jì)算結(jié)果 (a)對(duì)稱軸上能量損失值隨位置的分布;(b)對(duì)稱軸上的密度分布Fig.5.Theoretical calculation results:(a)The distribution of the energy loss value with the position on a symmetric axis;(b)density distribution on the symmetric axis.

圖6 Geant 4軟件模擬結(jié)果 (a)對(duì)稱軸上能量損失值隨位置的分布;(b)對(duì)稱軸上的密度分布Fig.6.Geant 4 simulation results:(a)The distribution of the energy loss value with the position on as symmetric axis;(b)density distribution on the symmetric axis.