蒙特卡罗方法模拟137Cs的γ能谱

蒙特卡罗方法模拟137Cs的γ能谱

问题：蒙特卡罗方法模拟137Cs源的661keVγ射线在NaI(Tl)闪烁体中的输运过程。 一、建立坐标系

以闪烁体中轴线为z轴建立柱坐标系，闪烁体位于0二、编程思想

rur一个粒子的状态用(E,r,)表示。其中E为粒子rur能量，r(r,,z)为粒子的位置坐标，(,)为

粒子的运动方向。其中为粒子运动方向与 z 轴的夹角的余弦，为粒子运动方向在 x y 平面上投影的方位角。

对于一个粒子经历以下过程：

1． 源抽样：

由于是点源，能量和位置的分布均为δ函数，抽样得到 E661KeV,r0,z2。 运动方向：各向同性分布密度函数为

v1f(Ω)f1()f2()4

11f1(),f2()22抽样得：211,22。

2．到达闪烁体：

如果22，无法到达闪烁体，丢弃，返回源抽样重新产生粒子。粒子进入闪烁体

的瞬间，状态为：r212,,z0其他量不变。

- 1 -

3．输运过程：

rurlnL抽样得到到下次碰撞的距离，根据当前粒子状态中的r,算出下次碰撞

t(Em1)r的坐标r(r,,z)，如果不在闪烁体区域（0－>

碰撞有两种可能：光电效应和康普顿散射。

根据粒子当前的能量，（由NaI(Tl)闪烁体宏观界面数据）线性插值确定它的光电效应截面和康普顿散射界面。抽样得到本次反应的类型。

如果光电效应，E=0，输运过程结束。如果康普顿效应，抽样获取碰撞后的能量和运动方向（康普顿散射的能量分布密度函数知道，具体抽样方法参考讲义。）

如果E<1KeV输运过程结束，反之，重复本过程直到输运过程结束。

4．记录与统计：

记录末态能量Ef，计算沉积能量EdE0Ef，考虑到测量系统分辨率，多道记录能量为沉积能量的高斯展宽。

FWHM(Ed)0.010.05E0.4E2 记录能量E'Edx。其中，FWNM/22ln2FWNM/2.355。x由标准正态分布抽样得到。

- 2 -

三、流程图

开始 源抽样产生粒子 N 到达探测器？ Y 计算下一次碰撞的位置Y N闪烁体内？ Y 碰撞次数=０? Y 碰撞次数+1 光电效应 反应形式？ 康普顿散射 NE＝0 抽样决定碰撞后能量和运动方向 Y E<1KeV N 计算沉积能量，记录能量，多道寻址＋1 计算粒子数够？ Y 结束 N

程序见附件。 编程环境：Matlab6.5

- 3 -

四、模拟结果

计算了5×106个粒子。 1) 探测效率与峰总比

N732650――进入闪烁体的总粒子数

Nm270128――探测到的总粒子数

Np617――全能峰的总计数（全能峰半高宽以内道的计数之和）

探测效率：峰总比：

估计全能峰能量的相对误差： 对落在半高宽内的粒子，如下计算

NmN36.870%

NpNm61733.176%

2701281ENEi1Ni1.2465

21N21N2EEiEi1.0686865

Ni1Ni1E1.0338(chn)5.169KeV

误差估计取置信水平1－＝0.95,2

EN2.05.1690.03453KeV

617由于此处取全能峰半高宽以“内”的道未经过插值，123道14998个计数算做全能峰半高宽以“内”，而127道13159计数不算全能峰半高宽以“内”。所以E平均偏小，误差计算也不准。

2）能谱图

- 4 -

图2 全谱图

3）线性插值求得半高宽

3x 1042.521.510.50122124126128 图3 全能峰局部展宽

- 5 -

chn=126.93-122.87=4.06 E=4.065keV=20.3keV

在662KeV的能量分辨率：=

E20.33.07% E662- 6 -