杭仲斌,刘川凤,滕忠斌,耿 璇,胡凯漩,贾茗涵,李 俏,高 飞,王红玉,魏可新,宋明哲,刘蕴韬

(1.中国原子能科学研究院,计量校准与技术重点实验室,北京 102413;2.中国计量科学研究院,北京 100029)

125I粒子源由于其辐射特性满足粒子植入治疗手段的需求,已经广泛用于多类型实体肿瘤的放射治疗[1],治疗过程中粒子源周围组织剂量分布是否准确会直接影响治疗效果[2-3]。因此,粒子源剂量分布的准确性是临床放射治疗中急需解决的关键问题。美国医学物理师协会(AAPM)先后发布了TG 43报告[4]和TG 43-U1报告[5],明确了粒子源剂量学参数的数学模型和几种常见粒子源的相关参数,但是由于粒子源几何构造和生产工艺的差异,粒子源在病灶组织处的剂量分布也会不同,因此有必要针对国内生产的粒子源进行模拟研究。本研究选择在我国应用最广泛的6711型粒子源作为模型,采用MCNP软件对粒子源进行更加全面和细致的模拟研究,以便为新型粒子源的剂量参数研究和更加精准的临床放射治疗提供参考。

1 国产6711型粒子源

使用原子高科股份有限公司生产的6711型粒子源,其放射性活度为0.3～1.0 mCi,长轴剖面结构示于图1,整体结构为两端半椭球的圆柱体,侧壁由0.05 mm厚的钛包壳包围,钛管壁用来吸收衰变过程中发射的电子[6],整体长为4.5 mm,直径为0.8 mm,源芯为φ0.5 mm×3.0 mm吸附有1 μm厚的Br5125I2的圆柱形银棒,Ti包壳和源芯之间充满干燥空气。125I半衰期为(59.4±0.01) d,发射光子能量为27.202、27.472、30.98、31.71、35.492 keV,其中对应的发射率分别为0.259 5、0.498 1、0.155 6、0.034 7、0.052 1[7]。

图1 6711型125I型粒子源长轴剖面结构图Fig.1 Long axis cross-sectional structure of 6711 125I particle source

2 剂量学参数和剂量场分布

根据TG 43 U1报告,计算粒子源一维和二维剂量分布公式中,粒子源剂量率常数Λ为绝对量,几何函数GL(r,θ)、径向剂量函数gL(r)、一维各项异性函数φan(r)、二维各项异性函数F(r,θ)为相对量,且以上参数均可以通过理论计算和蒙特卡罗方法模拟得到。理论计算采用的极坐标系示于图2。参考点位置(r0,θ0)一般为距源中心1 cm,极角90°的位置,即r0=1 cm,θ0=90°。

图2 用于粒子源剂量学参数计算的极坐标Fig.2 Polar coordinates for particle source dosimetry parameter calculation

(1)

剂量率常数Λ使用参考点位置处的水吸收剂量率和空气比释动能强度的比值进行表征,如公式(2)所示。

(2)

基于线源的几何函数可通过离散点数据使用插值方法进行计算,高度简化的近似值为治疗规划提供了足够的准确性。基于线源的几何函数如公式(3)与(4)所示,放射性分布在圆柱形或环形空间时,圆柱体的长度L作为有效长度,即本次计算取L=3 mm。

(3)

GL(r,θ)=(r2-L2/4)-1(θ=0°)

(4)

式中:GL(r,θ)为几何点的几何函数;β为几何点距离圆柱形银棒两端的夹角,°;r为几何点的极径,mm;θ为几何点的极角,°。

径向剂量函数gL(r)考虑介质对射线的吸收和散射以及粒子源的外包壳对光子的吸收等影响,其计算见公式(5)。

(5)

一维各向异性函数φan(r)表示在距离r不变的情况下,4π范围内吸收剂量的积分平均值与θ=90°位置处的吸收剂量的比值,由公式(6)表示。

(6)

二维各向异性函数F(r,θ)表示当r相同角度不同时,粒子源的包壳厚度和周围介质(组织或水)对光子吸收、散射对剂量的影响,其物理含义表示以粒子源中垂线处剂量为参考,粒子源偏离中垂线的剂量相对于参考点处(r=1 cm,θ=90°)剂量的比值,其计算见公式(7)。

(7)

3 蒙特卡罗模拟和结果分析

3.1 蒙特卡罗模拟

利用MCNP软件对国产6711型粒子源结构进行精确建模。根据公式(1)可知,当源项相同时,即可用单位时间内的累积剂量来替代吸收剂量率和空气比释动能率。所以只需要分别模拟计算距离源1 cm处的水吸收剂量和1 m处的空气比释动能。某点的吸收剂量用质能吸收系数和能注量表征,比释动能则是用质能转移系数和能注量表征,查阅文献[8]可知,当介质的质量数Z和初始能量均很低时,由于韧致辐射影响很小,所以质能吸收系数和质能转移系数近似相等,在水和空气的介质中,当光子能量小于0.4 MeV时,质能吸收系数和质能转移系数推荐值相等,详细数据列于表1。

表1 水和空气中质能转移系数与质能吸收系数对照Table 1 Comparison of mass and energy transfer coefficient with mass and energy absorption coefficient in water and air

由于计算水吸收剂量时需要在粒子源径向距离1 cm处放置探测器,所以选用F5卡(点探测器)进行计算,邻域半径设为0,结合水的质能转移系数μtr/ρ,可得到单位时间内的吸收剂量;计算空气比释动能时,探测器结构使用PTB的GROVEX电离室[9],探测栅元形状为圆台型,极角为16°,厚度为1 cm,设置截止能量为5 keV,使用*F4卡用于确定探测器的能量通量,结合空气的质能吸收系数μen/ρ可得到在真空中距离粒子源1 m处的空气比释动能,乘以距离的平方即为单位时间的空气比释动能强度。

3.2 剂量率常数

根据剂量率常数的定义,分别模拟得出粒子源径向距离分别为1 cm水吸收剂量和1 m处的空气比释动能强度,代入公式(2),计算得出剂量率常数为0.957 cGy·h-1·U-1,TG 43-U1推荐值为0.965 cGy·h-1·U-1,相对误差为-0.83%。

3.3 径向剂量函数

将粒子源放入水中,在粒子源的径向不同距离上放置探测器,为进行更加细致的研究,将距离点分布更加密集,计算得出吸收剂量后代入公式(5),得出径向剂量函数,并与TG 43-U1推荐值和文献[10]推荐值进行比较,经整理,径向剂量函数的结果列于表2。

表2 径向剂量函数Table 2 Radial dose function

本次计算结果和TG 43-U1推荐值比较,相对误差范围为-0.49%～3.99%,与文献[10]比较,相对误差范围为-0.46%～4.00%,一致性较好。其中在r=(0.05～0.25) cm范围内,径向剂量函数出现上升的趋势,这是由于在这个范围内银棒的直径不可忽略,线源模型已经不能很好地适用于此类情况。目前很多商业治疗计划系统通过对距离r进行五阶多项式拟合,从而得到径向剂量函数公式。在对径向剂量函数关于距离r进行五阶多项式函数拟合时,多数文献均选用0.1 cm作为起始点[11-12]。本次拟合选用范围(0.25～10) cm,可得gL(r)的五阶多项式函数,其中a0=1.124 52,a1=8.806 6×10-2,a2=5.597 5×10-2,a3=1.147×10-2,a4=1.33×10-3,a5=4.379 48×10-5,相关系数R2=0.999 97,残差平方和为1.065 92×10-4,拟合精度好于0.1 cm～10 cm范围拟合的数据。

3.4 各向异性函数

将粒子源放置水中,在距离粒子源相同距离,在不同角度上放置探测器,距离选择0.25、0.5、0.75、1、1.5、2、2.5、3、4、5、6、7 cm,角度范围为(0～90)°,间隔5°,带入公式(6)和公式(7)计算各点的吸收剂量,分别整理得到一维和二维各向异性函数,列于表3。

表3 各向异性函数Table 3 Various anisotropic functions

将本研究中不同距离处的二维各向异性函数绘制成图,结果示于图3。通过分析可知,r=0.25 cm处的二维各向异性函数在θ=45°附近出现了明显的“鼓包”,这一现象与文献[10]、文献[12]模拟结果相同,其原因是粒子源内置银棒为直角型,在θ=45°附近,对125I产生的射线的衰减有所差异。通过比较发现,本次模拟二维各向异性函数的结果整体偏大,这是由于粒子源的构造引起的,文献[10]所用粒子源的银棒末端为平面倒角结构,Ti包壳厚度为0.11 cm,本研究所用粒子源则为直角结构,Ti包壳厚度为0.5 cm,出射射线所经历的衰减程度有所差异。这些原因导致银棒末端为直角型结构的粒子源二维各向异性大于银棒末端为斜切方向的二维各向异性函数的现象,这与文献[13]结论一致。

图3 国产6711型粒子源二维各向异性函数Fig.3 2D Anisotropy function of domestic 6711 particle source

将本研究与文献[10]的一维各向异性函数数据对比,示于图4,整体的变化趋势相同,不过整体数据大于文献[10]的推荐值,其原因也是两种粒子源的源芯底部结构差异导致。

图4 一维各向异性函数对比Fig.4 Comparison of 1D anisotropic functions

4 结论

本研究利用MCNP模拟软件对国产6711型粒子源的剂量学参数进行蒙特卡罗模拟研究,通过对粒子源进行精准建模,得到粒子源剂量率常数Λ、径向剂量函数gL(r)、一维各向异性函数φan(r)和二维各向异性函数F(r,θ),并将数据与TG 43-U1报告和相关文献推荐值进行比较,结果发现,剂量学常数与TG 43U1报告推荐值误差-0.83%,径向剂量函数结果分别与TG 43-U1报告推荐值比较,相对误差范围为-0.49%～3.99%,二维各向异性函数和一维各向异性函数结果整体比国外6711型粒子源偏大,这是由于国产6711型粒子源在源芯结构设计上与国外有细微区别导致。但由于TG 43-U1推荐值是蒙特卡罗模拟和实验测量的综合评价值,且蒙特卡罗软件自身存在的系统误差,结合本次模拟的数据和相关文献的相对误差值,本次模拟数据较为可信,并且本次模拟拓展了径向剂量函数和各向异性函数的数据范围。本研究可为国产6711型125I粒子源的临床治疗中所使用的剂量分布数值的准确性提供一定的参考。