使用mossbSim軟件模擬穆斯堡爾實驗

王思廣，羅棱尹，賈春燕

(北京大學 物理學院 核物理與核技術國家重點實驗室，北京 100871)

穆斯堡爾譜儀分辨能力高，可精確探測原子核能級的變化，在物質結構研究領域有許多應用實例[1-7].北京大學開設穆斯堡爾實驗[8]課程多年.穆斯堡爾實驗涉及到的物理知識與常規的近代物理實驗所研究的射線與物質的相互作用有本質區別，其實驗設計非常巧妙，教學過程有些問題值得思考，例如：為什么選用薄的探測器用于穆斯堡爾實驗？將振動幅度增加1倍后穆斯堡爾譜會如何變化？如果上、下閾卡得不合適會如何影響最后的穆斯堡爾譜？放射源與樣品之間的距離對穆斯堡爾譜有何影響？由于實驗條件限制，實驗室未配備不同厚度的探測器進行實驗驗證，實驗課堂教學時間也不允許學生調節不同振動速度的幅度多次實驗進行檢驗.理論講解對于初步接觸γ射線與探測器的學生非常抽象.為了讓學生深入了解實驗，利用Geant4軟件包[9-10]開發了穆斯堡爾實驗的模擬軟件mossbSim.

1 穆斯堡爾實驗原理

1957年，穆斯堡爾研究191Ir的γ射線的共振散射現象時發現了穆斯堡爾效應：在固體中的核發射或吸收γ射線時存在部分原子核無反沖的現象.要理解原子核無反沖對于穆斯堡爾實驗的重要性，需要比對能級寬度與反沖能量的量級.

如果1個核的激發態壽命為τ，則其在退激發出的γ射線的能量不單一.能量分布I(E)與發出的γ射線的能量E之間的關系通常用洛倫茲函數描述：

(1)

(2)

對于靜止的原子核，通過能級差為E0的躍遷發出動量為pγ、能量為Eγ的γ射線后，原子核的反沖能量為

(3)

式中c是光在真空中的速度.

我校穆斯堡爾實驗用57Co放射源.根據以上公式及57Fe的第一激發態與基態之間的14.4 keV的能級差，可以計算出發生躍遷時的反沖能量為ER≈2×10-3eV.由于該能級的壽命τ≈0.1 μs，可以計算出對應的自然寬度Γ=4.9×10-9eV.因為所發出的γ射線的能量的改變ER遠大于共振的自然寬度Γ，故如果出現原子核反沖，就不可能觀察到14.4 keV的γ射線打在處于基態的57Fe原子組成的樣品上被共振吸收躍遷到第一激發態的現象.

發生共振吸收后，由于第一激發態的能級壽命τ僅約0.1 μs，故激發后的原子核將很快退激，再次放出14.4 keV的γ射線.與被吸收的γ射線相比，新放出的γ射線的出射方向各向同性，故有很大概率要偏離被吸收的射線的方向.如果放射源、樣品、探測器在1條軸線上，且樣品在放射源及探測器之間，這種吸收再發射的過程將導致探測器測量到的計數減少.

為了顯示出發生共振時探測器探測到的計數率與不發生共振時的計數率的不同，在穆斯堡爾實驗中，放射源相對樣品的運動速度將被逐漸改變，將探測器探測到的14.4 keV的射線峰的計數按照不同的補償運動速度段進行順序存儲.反復測量不同的補償速度對應的探測器探測到的14.4 keV能峰的計數(在能峰兩側設置上、下閾，記錄該峰內的計數).最后得到穆斯堡爾譜.如果補償速度掃描范圍包含各共振位置而且超出的范圍不算過大，就會看到清晰的補償速度滿足共振條件與不滿足共振條件時的計數率的差異：發生共振吸收時，在穆斯堡爾譜上對應位置出現吸收谷.

上述通過相對速度的改變而增減γ射線的能量的原理系多普勒效應.能量改變ED與放射源發出的γ射線的能量Eγ及相對于樣品的速度v成正比：

(4)

其中，θ是放射源運動方向與發出的γ射線的出射方向之間的夾角.

對于57Fe的14.4 keV的γ射線，如果放射源的運動速度為1 mm/s，通過多普勒效應對其改變的能量ED約為4.80×10-8eV，約為該能級10倍的自然寬度，α-Fe的穆斯堡爾譜各共振點所需要的補償速度的掃描范圍大致在±6 mm/s內.

基于以上理論，用Geant4軟件包開發了模擬穆斯堡爾實驗的軟件mossbSim.

2 模擬的核心方法

根據57Co衰變到57Fe的衰變綱圖[8]中提供的射線能量及分支比，利用能量分布函數I(E)抽樣出57Fe退激放出的γ射線的能量Eγ，抽樣定出γ射線的飛行方向，給出γ射線的四動量.然后根據衰變時刻的補償速度v計算出補償能量ED，得到調制后的γ射線新的四動量.

跟蹤模擬能量補償后的射線與物質的相互作用的任務由Geant4常規功能完成.但是樣品吸收γ射線再發射的過程需要用到Geant4的G4VFastSimulationModel特殊處理模塊.在mossbSim軟件中的做法是：在樣品的中間加了1層與樣品材質相同的α-Fe膜，膜的厚度比樣品的厚度至少小2個量級.當Geant4跟蹤模擬γ射線到這層膜時，利用自行設計的程序根據γ射線的截面σ(E)函數進行抽樣，決定發生共振吸收還是沿著原來的射線方向穿過這層膜.如果抽樣結果是穿過這層膜，將γ射線沿其動量方向移到膜的外邊界與樣品材質的交界處(在樣品內部)的點，并從交界點處沿著原來的方向且保持原來的動量繼續模擬與物質的相互作用.如果抽樣的結果是發生共振吸收，則通過隨機抽樣的方法確定共振吸收后退激放出新的γ射線的方向，該方向與原來的入射方向無關(模擬原子核退激各向同性隨機放出γ射線的物理現象)，新的γ射線的能量與被吸收的γ射線的能量、位置相同.將新的γ射線從其產生位置，沿其動量方向移到膜的外邊界與樣品材料的交界處(在樣品內部)，然后繼續模擬其與物質的相互作用，直到打在探測器(這里用的是NaI晶體)上將部分或全部能量沉積下來或者跟蹤到感興趣區域外.被吸收的γ射線在吸收發生時直接被殺死.

借助模擬時逐漸增大的Geant4事件序號(Event ID)，實現速度掃描：用當前的事件序號與總道數(這里用512道)計算得到余數i，定出該事件模擬得到的信號放在第i道.然后根據預設的速度變化規律，計算出第i道對應的掃描速度v，根據v計算出補償能量.最后抽樣給出這個事件產生的射線數量及對應的能量和角度，計算出考慮補償能量后的各射線的實際能量，模擬計算出該事件在探測器有效靈敏體積內沉積的能量總和及產生的總的光學光子數n，進而判斷n是否介于預設的上、下閾之間，如果是則在第i道加1.

設置各共振位對應的補償能量基于α-Fe的6條譜線v1，v2，…，v6的已知條件：

1)v6-v1=10.656 mm/s；

3)各譜線之間的能級間距已知[8].

3 模擬結果

首先模擬57Co放射源的能譜.如圖1所示，X為道數，N為每道計數.紅色豎線為所設置的14.4 keV的γ射線產生的信號峰左、右兩側的閾值.No.1能峰對應能量為6.4 keV的Fe的X射線，No.2能峰對應14.4 keV的γ射線，No.3能峰為123 keV及137 keV的γ射線在NaI中I元素產生的X射線逃逸后的剩余能量產生的信號峰，No.4為123 keV及137 keV的γ射線全能峰共同組成的信號.

圖1 模擬NaI探測器探測到的57Co放射源的 衰變能譜圖

模擬所用的NaI探測器晶體厚度為0.1 mm，模擬過程中考慮了57Fe的第二激發態到基態的躍遷發出的137 keV的γ射線，第二激發態到第一激發態的躍遷發出的123 keV的γ射線，第一激發態到基態的躍遷發出的14.4 keV的γ射線及第一激發態通過內轉換過程使得內層電子電離引發的Fe的K系標識X射線(模擬用的X射線能量是6.4 keV).

依據圖1所示的能譜，可以對14.4 keV的峰設置上、下閾，觀察落在上、下閾之間的信號所對應的補償速度的分布，即穆斯堡爾譜，如圖2所示.探測器厚度為0.1 mm，補償振動速度范圍為±8 mm/s.與實驗獲取的穆斯堡爾譜表示方法一致.其中X軸與補償速度v之間有轉換關系，該轉換關系是模擬程序預設的(mossbSim中預設了第1道對應于+v，中間第256道對應于-v，第512道對應于+v，其中+v到-v再到+v的過程系勻減速及勻加速變化).因在第1道至第512道之間對應的補償速度經歷了先減小再增大的過程，故如果設置的速度變化范圍±v內包括了這6個共振位置，則每個模擬周期會經歷2次共振點，故所模擬出的α-Fe的穆斯堡爾譜有12個共振位置，并且關于速度變化的轉折位置(這里是第256道)左右對稱.

圖2 將補償振動速度范圍設置為±8 mm/s后 模擬α-Fe的穆斯堡爾譜

考慮到補償能量與補償速度的一一對應關系，如果將振動速度的范圍增加或減小，則共振位置將因所對應的補償速度位于坐標軸上的相對位置變化而變化.為了驗證這點，將振動速度范圍由原來的±8 mm/s提高到±20 mm/s后，所得到的穆斯堡爾譜如圖3所示.

圖3 將補償振動速度范圍設置為±20 mm/s后模擬的α-Fe的穆斯堡爾譜

圖2和圖3的明顯差異在于圖3較圖2中共振點之間的距離收縮(以道數表示)，但實際上如果以補償速度描述共振點的距離，圖2與圖3上相鄰的共振點所對應的補償速度應該完全相同，當然共振點間的補償速度之差也應該完全相同.

以上模擬所用的探測器的厚度為0.1 mm.選用如此薄的探測器的原因主要是為了壓低放射源發射的123 keV及137 keV的γ射線對實驗結果產生的干擾：實驗選用的放射源為57Co，其發出的γ射線中除用來觀察穆斯堡爾現象的14.4 keV的γ射線，還有123 keV和137 keV的γ射線，并且兩射線的分支比分別為91%及9%，因從第一激發態退激到基態發出14.4 keV的γ射線有內轉換過程與其進行競爭，內轉換系數α=15(即發出內轉換電子的概率與發出γ射線的概率之比為15∶1)，故14.4 keV的γ射線的實際分支比非常小，約為5.7%，具體參閱參考文獻[8]提供的衰變綱圖.123 keV的γ射線與14.4 keV的γ射線有級聯衰變效應，即從第二激發態躍遷到第一激發態發出123 keV的γ射線后，因第一激發態的能級壽命很短，故很快繼續躍遷到基態發射出14.4 keV的γ射線或通過內轉換過程發射出內轉換電子.如果這2條前后時間相差很小的γ射線幾乎同時發出并都被探測器探測到，而探測器的響應時間較慢，無法區分2條γ射線，則探測器給出的總的能量將為2條γ射線在探測器上沉積的能量相加，這種真符合效應將減少14.4 keV信號數.對于14.4 keV的γ射線，很薄的1層NaI晶體就能夠將其能量完全沉積下.而對于能量比其高的123 keV的γ射線，因穿透能力較強，故被探測到的效率低，通過真符合效應對14.4 keV的信號影響小.

同樣的原因，對于薄的探測器，137 keV的γ射線通過偶然符合對14.4 keV的信號的影響也小.但是，如果選用厚的探測器，因為對123 keV以及137 keV的γ射線的效率的提高而導致對14.4 keV的γ射線的符合效應的影響增大.

另外，如果選用厚的探測器，123 keV以及137 keV的γ射線在14.4 keV能峰附近產生的康普頓散射平臺將比較薄的探測器上產生的康普頓散射平臺高，從而感興趣的14.4 keV信號的信噪比低.

為了驗證以上說法，將NaI的厚度分別設為0.1，0.5，10.0mm進行模擬，得到57Co的能譜如圖4所示，對應的穆斯堡爾譜如圖5所示.從圖4可以看出:隨著探測器厚度的增加，123 keV以及137 keV組成的復合峰包的γ射線信號的幅度較14.4 keV的γ射線峰的幅度迅速增加.從圖5也可以看出：隨著探測器厚度的增加，模擬同樣的事件數(每種厚度條件下模擬的事件總數均為1×1010)得到的最后譜圖中的總計數也增加.

圖4 模擬NaI探測器晶體的厚度分別為0.1，0.5，10.0 mm時57Co放射源的衰變能譜圖

圖5 模擬NaI探測器晶體的厚度分別為0.1，0.5，10.0 mm時所測得的穆斯堡爾譜

為了便于比較圖5中共振點對應的谷的深度與本底基線的比例(對于穆斯堡爾譜，本文定義該比值為信噪比)，將各穆斯堡爾譜進行歸一化，如圖6所示.

圖6 模擬NaI探測器晶體的厚度分別為0.1，0.5，10.0 mm時的歸一化穆斯堡爾譜

從圖6可以看出：0.1 mm厚度的探測器所探測到穆斯堡爾譜上的共振信號的信噪比好于10.0 mm厚的探測器得到的信噪比.另外，因電子學處理每個信號都需要一定的時間，故實驗中123 keV及137 keV的γ射線在厚的探測器中因探測效率的增加將增加譜儀的死時間率.故選用較薄的探測器.

圖2～6都基于圖1所示的上、下閾區間挑選的信號給出的穆斯堡爾譜.如果該上、下閾區間選取不合適，例如誤將6.4 keV的X射線峰及其下本底包含在閾值內，因這些射線不會發生共振吸收，故它們抬高了整個基線的高度而對共振信號谷的深度的增大沒有實質貢獻，最終將導致信噪比拉低.

圖7給出對于同樣厚度的探測器設置不同上、下閾后的差異.圖7中標記為“Wrong range”譜數據系將圖1所示的上、下閾區間范圍由原來的35～65道擴展為5～65道，擴展后的上、下閾包含6.4 keV的峰及峰下本底計數.標記為“Right range”譜數據與圖2中的穆斯堡爾譜數據相同.2種譜數據各自歸一化后畫在一起.顯然標記為“Wrong range”的穆斯堡爾譜信號的信噪比要比標記為“Right range”的穆斯堡爾譜信號的信噪比差.

圖7 設置不同上、下閾所得的穆斯堡爾譜圖

利用mossbSim模擬軟件，也可以研究放射源與樣品之間距離對穆斯堡爾譜的影響.圖8所示系放射源距離樣品5 mm及100 mm所對應的穆斯堡爾譜.顯然距離較近時的穆斯堡爾譜分辨變差.可以引導學生從距離短時大角度入射到樣品中的射線份額比例較大進行思考原因.

圖8中12個谷對應著共振點的位置.左側1～6號吸收谷在模擬程序中被設置的初始相對強度分別為1.0，0.8，0.3，0.3，0.8，1.0,其中1號與6號所需要的補償速度較2號和5號大.D=5 mm時，因當入射角較小與較大的射線不能同時滿足共振條件，故1號和6號對應的分辨率較后者差，吸收谷變寬，導致深度反而不如2號和5號谷.D=100 mm時入射角差異不大，故分辨率較好.

(a)D=5 mm

除以上列舉的改變相關參量進行穆斯堡爾譜形狀的研究外，學生也可以從源代碼上做多樣化的探究，比如：為研究振動系統不完美的狀況下對穆斯堡爾譜的影響，可將振動速度由現有的勻加速和勻減速改變換成曲線變化的速度增減；加入振動的空程差(加減速度時有偏移量，即道址與振動速度不是一一對應，而是人為隨機加入錯位)，觀察有、無空程差的穆斯堡爾譜之間的差異；改變探測器的分辨率，即將沉積的能量與轉換成的光子數對應的系數改變，研究探測器的分辨率與共振信號的寬度、信噪比之間的關系.

mossbSim軟件模擬結果存儲為ROOT文件形式[11]，其優點是可以很方便地進行數據分析并畫圖比較.也為學生提供了可以將ROOT數據格式文件轉換為普通文本格式文件的程序，供不愿學習ROOT軟件的學生利用分析實驗數據的方法分析模擬數據，以便將實驗數據的分析結果與理論模擬數據的分析結果進行比較.

作為比較，圖9為實驗穆斯堡爾譜圖.模擬的穆斯堡爾譜(圖8)與實驗譜相似.細微差異主要來自于探測器各組成的尺寸參量及掃描速度的差異.模擬中各吸收谷的相對強度是臨時設置的1組數據，與實際也有差異.如果將實驗裝置的各參量輸入準確，可以通過比較模擬與實際測得的穆斯堡爾譜給出α-Fe共振谷的相對強度.模擬的輸出、輸入比值可作為數據處理過程所需的探測效率修正因子.故本模擬軟件可以作為細致研究穆斯堡爾譜的起點.

圖9 α-Fe樣品的實驗穆斯堡爾譜圖

4 結束語

利用Geant4開發了模擬穆斯堡爾實驗的軟件mossbSim.該軟件的核心是利用Geant4的G4VFastSimulationModel功能模塊處理樣品共振吸收再發射的過程，實現了穆斯堡爾實驗的蒙特卡羅模擬.利用mossbSim軟件，學生可以調節探測器的厚度、速度振動范圍，對14.4 keV的峰設置不同的上、下閾，模擬出α-Fe的穆斯堡爾譜.通過模擬，學生將會深入理解穆斯堡爾效應的原理；通過對模擬代碼的研究，將有助于了解和掌握Geant4這一通用的開源模擬軟件包的基本功能，同時該軟件也為學生提供接觸ROOT分析軟件的機會.將該軟件應用于教學課堂的優點是學生可以學到很多模擬及數據分析相關的知識，缺點是對任課教師要求比較高(需要Linux相關的基本知識).為了減輕教師負擔，我們也準備了教學錄像，學生可以通過多媒體演示了解程序的運行方式.全部源代碼見http://wlsy.nenu.edu.cn/doku.php?id=wenzhang:2101:wsg.