CN104750992A - 模拟粒子输运和计算放疗中人体剂量的方法、装置及系统 - Google Patents

Abstract

本发明涉及模拟粒子输运和计算放疗中人体剂量的方法、装置及系统。所述模拟粒子输运的方法包括：记录所产生输入粒子的输运径迹；基于每批次运行粒子的径迹计算每个栅元的不确定度，若栅元的不确定度不超过第一阈值，则该栅元为达标栅元；获取感兴趣区域的达标率，所述感兴趣区域至少包括一个栅元，所述感兴趣区域的达标率为该区域达标栅元占该区域所有栅元的比例；若所述感兴趣区域的达标率超过第二阈值，则停止继续输入粒子，并输出历史输入粒子的输运径迹。本发明能够提高对放射治疗过程所涉及粒子输运过程的模拟效率。

Description

模拟粒子输运和计算放疗中人体剂量的方法、装置及系统

技术领域

本发明涉及放射治疗技术领域，特别涉及一种模拟粒子输运的方法及装置、一种计算放射治疗中人体剂量的方法及一种放射治疗系统。

背景技术

在放射治疗技术领域中，人体组织剂量分布的计算方法大致可以分为两大类：半经验解析方法及蒙特卡罗方法。

半经验解析方法包括基于离轴比(OAR，Off Axis Ratio)经验公式的方法，以及基于笔束核(Pencil Beam Kernel)和点核(Point Kernel)的卷积/叠加(Convolution/Superposition)方法等。但是半经验解析方法的精确度有限。

而蒙特卡罗方法，因为其处理复杂问题(复杂几何、复杂的放射源布置等)的能力而成为一种不可替代的方法。蒙特卡罗方法可以精确地对放射治疗过程所涉及到的物理过程进行建模，而使用较少的近似。蒙特卡罗方法的最大不足之处在于计算强度较大，耗时较长。

发明内容

本发明技术方案所解决的技术问题为，如何提高对放射治疗过程所涉及粒子输运过程的模拟效率。

为了解决上述技术问题，本发明技术方案提供了一种模拟粒子输运的方法，适于模拟粒子在栅元中的能量分布，包括：

估算所需入射总粒子数，产生入射粒子并分批输入；

记录所输入粒子的输运径迹；

基于每批次运行粒子的径迹计算每个栅元的不确定度，若栅元的不确定度不超过第一阈值，则该栅元为达标栅元；

获取感兴趣区域中栅元的达标率，所述感兴趣区域至少包括一个栅元，所述感兴趣区域的达标率为该区域中达标栅元占该区域所有栅元的比例；

若所述感兴趣区域中栅元的达标率超过第二阈值，则停止继续输入粒子，并输出历史输入粒子的输运径迹，否则继续输入粒子，直到总粒子数运行完毕。

可选的，分批输入的每一批入射粒子的能量接近，或类型相同，或能量接近并且类型相同。

可选的，所述分批输入包括：

根据所述入射粒子的类型不同交错地分批输入。

可选的，所述输运径迹包括：入射粒子的能量信息、速度信息及其他径迹信息，所述速度信息包括所述入射粒子的入射方向信息，所述记录所入射粒子的输运径迹包括：

若入射粒子与已记录粒子的能量信息与入射方向信息接近，则将该已记录粒子的输运径迹赋值为该入射粒子的输运径迹。

可选的，所述其他径迹信息包括：所述入射粒子的类型信息、入射位置信息、权重信息及所经栅元信息；

所述所经栅元信息包括：所经栅元的能量分布以及不确定度。

可选的，所述方法还包括：当每一批粒子输入时：

若其中一个粒子从第一重要性栅元进入第二重要性栅元，则所进入的粒子以第一概率进行分裂，分裂后的粒子降低权重以使该批粒子的总权重不变，所述第一重要性栅元的重要性低于第二重要性栅元的重要性；

若其中一个粒子从第三重要性栅元进入第四重要性栅元，则所进入的粒子中以第二概率被杀死，未被杀死的粒子增加权重以使该批粒子的总权重不变，所述第三重要性栅元的重要性高于第四重要性栅元的重要性。

可选的，所述第一概率为第一重要性栅元与第二重要性栅元的重要性比值，所述第二概率为第四重要性栅元与第三重要性栅元的重要性比值。

可选的，所述重要性栅元的重要性为预先手动设定，或根据栅元信息自动设定；所述栅元信息包括栅元的不确定度，或者栅元的物理属性。

可选的，所述方法还包括：

基于历史输入粒子的输运径迹对入射粒子的剂量分布的不确定度做动态降噪处理。

可选的，所述动态降噪处理通过如下方式实现：

得到所述栅元中粒子的三维的剂量分布和不确定度；

对所述三维的剂量分布进行滤波处理，以使所述剂量分布在三个维度上连续可导；

重新计算得到所述剂量分布对应的不确定度。

可选的，所述方法还包括：导入几何模型，所述几何模型包括：模拟对象的栅元、物理材料、栅元权重或/和几何虚拟截面，其中，所述几何虚拟截面适于定义所述栅元对应物理材料，以使所述栅元指向的模拟对象具有均匀化材料，所输入粒子的输运径迹与所述几何虚拟截面有关。

为了解决上述技术问题，本发明技术方案还提供了一种计算放射治疗中人体剂量的方法，包括：根据如上所述的模拟粒子输运的方法获得的粒子在栅元中的能量分布以计算放射治疗中的人体剂量。

为了解决上述技术问题，本发明技术方案还提供了一种模拟粒子输运的装置，适于模拟粒子在栅元中的能量分布，包括：源处理模块、输运处理模块、噪声处理模块及输出模块；

所述源处理模块，适于估算所需入射总粒子数、产生入射粒子并分批输入；

所述输运处理模块，适于记录所输入粒子的输运径迹；

所述噪声处理模块，适于实现如下步骤：

基于每批次运行粒子的径迹计算每个栅元的不确定度，若栅元的不确定度不超过第一阈值，则该栅元为达标栅元；

获取感兴趣区域的达标率，所述感兴趣区域至少包括一个栅元，所述感兴趣区域的达标率为该区域达标栅元占该区域所有栅元的比例；

若所述感兴趣区域中栅元的达标率超过第二阈值，则停止继续输入粒子，并通过所述输出模块输出历史输入粒子的输运径迹，否则使继续输入粒子，直到总粒子数运行完毕。

为了解决上述技术问题，本发明技术方案还提供了一种放射治疗系统，包括：

如上所述的模拟粒子输运的装置，适于模拟粒子在栅元中的能量分布；

计算人体剂量的装置，适于根据所述模拟粒子输运的装置获得的粒子在栅元中的能量分布计算放射治疗中的人体剂量。

上本发明技术方案的有益效果至少包括：

本发明技术方案估算所需入射的总粒子数，通过用户对不确定度的要求，估算出具体模体内需要输运的总粒子数，作为总的计算目标，在能够在达到用户目标的情况下，减少不必要的粒子输运，并在一些无法截断的情况下，减少输运粒子数，提高模拟效率。

本发明技术方案能够提高放射治疗过程中模拟粒子输运的速度，通过对运行粒子不确定度的处理和评估，可快速达到计算目标，并在达到目标后对模拟过程进行截止，从而在减少了计算用时的同时达到用户对不确定度的要求，大大提高模拟效率。

本发明技术方案在感兴趣区域，可以根据运行粒子的不确定度动态地调整抽样的粒子数，能够在保证全局不确定度均衡的前提下对粒子进行分批输运，可减少粒子抽样，减少不必要的粒子输运，节约了大量计算用时。

在本发明技术方案的可选方案中，还采用了粒子全同化并行处理的技术手段，对入射粒子进行分类分别处理；本发明技术方案对入射粒子按能量与类型进行分类，并将能量接近和类型相同的粒子作为同一批次，可有利于并行计算单元在接近的时间内完成计算，从而加快并行计算的模拟速度。

在本发明技术方案的可选方案中，还在每批粒子的模拟后对粒子的剂量分布和不确定度进行滤波处理，实现动态降噪，提高不确定度的均衡性，使所有感兴趣区域的计算点的不确定度降到可接受范围。

在本发明技术方案的可选方案中，还增加了入射粒子的均匀性，对入射粒子的方向进行均匀化处理，降低源入射粒子的不确定度，从而有利于粒子动态降噪。

在本发明技术方案的可选方案中，还在粒子的输运过程中，对其模拟过程中的虚碰撞反应(虚拟反应)与真实物理反应(真实反应)进行了区别抽样，对虚拟反应下的粒子仅抽样输运程度而不重新抽样粒子方向及能量，仅在真实反应下才重新对粒子方向及能量进行抽样，从而减少总体粒子的抽样数目。本发明技术方案采用了如下方式减少总体粒子的抽样数目：对虚拟截面减少粒子输运过程中发生能量抽样和方向抽样数目；对相同的入射粒子的径迹进行直接复制，避免重复抽样；对粒子进行能量截断以减少低能粒子的输运；对粒子进行权截断以减少重要性低的粒子的输运。

附图说明

图1为本发明技术方案提供的一种模拟粒子输运的方法流程示意图；

图2为本发明技术方案提供的另一种模拟粒子输运的方法流程示意图；

图3为本发明技术方案提供的又一种模拟粒子输运的方法流程示意图；

图4为本发明技术方案涉及的粒子在输运过程中发生截面反应的示意图。

具体实施方式

在放射治疗领域中，关注并需要模拟各类治疗机产生的射线束，比如加速器产生的高能电子束、质子束、重离子和光子束、钴-60(Co)光子束以及X射线治疗机产生的X线束等。模拟治疗机产生的射线束主要输出的是粒子的模拟信息，可记录所有粒子到达或穿过用户所定义的几何空间的相关信息，包括粒子所带电荷、能量、位置、方向和粒子所经过材料的历程标记等。

基于蒙特卡罗模拟粒子输运方法需要输入射野内粒子能谱分布、角分布和空间分布等有关治疗射线束特性的详细信息，按照所获取的信息可利用蒙特卡罗算法模拟大量粒子的输运过程。

一种利用蒙特卡罗算法模拟粒子输运的方法，实现对加速器治疗头进行模拟，其主要步骤如下：

首先，建立模拟加速器治疗头的用户程序和所需的粒子与介质相互作用的截面数据；

其次，完成用户的输入，主要包括加速器治疗头各部件模块的几何定义，入射粒子束定义以及程序运行控制参数的选择；

再次，进行模拟计算以及分析计算的结果，并把模拟计算结果作为体模吸收剂量计算的源项输入。

使用蒙特卡罗算法进行粒子模拟是一种采用随机抽样的精确方法，但是，现有技术普遍认为，蒙特卡罗的抽样特性使该算法具有固有耗时性，需要大量的随机数抽样才能提高模拟精度。

本发明技术方案使用蒙特卡罗模拟粒子输运方法，能够结合粒子数估算以及粒子不确定度的算法，减少粒子抽样，从而大大提高了基于蒙特卡罗算法模拟粒子输运过程的速度及效率。

为了使本发明的目的、特征和效果能够更加明显易懂，下面结合附图对本发明的具体实施方式做详细说明。

在下面的描述中阐述了很多具体细节以便于充分理解本发明，但是本发明还可以采用其他不同于在此描述的方式来实施，因此本发明不受下面公开的具体实施例的限制。

实施例一

使用蒙特卡罗算法模拟粒子输运过程，是一种采用大量随机抽样来达到精确计算的方法，尤其，为了满足用户对粒子输运不确定度的要求，往往需要模拟大量粒子，具有耗时的特性。

本实施例提供了一种模拟粒子输运的方法，适于模拟粒子在栅元中的能量分布，可基于蒙特卡罗模拟算法，其通过计算几何栅元的不确定度、区分栅元的重要性及均衡感兴趣区域内栅元的不确定度，能够减少粒子的抽样，节省了大量计算用时。

参考图1，本实施例提供的模拟粒子输运的方法包括：

步骤S100，估算入射粒子数、产生入射粒子并分批输入。

所述入射粒子是基于源产生的，所述源是蒙特卡罗工具对一已知放射源的模拟结果，其表现形式为相空间源，基于相空间文件。

估算入射粒子数方法通过预先模拟不同粒子数在均匀模体(例如水模体)，或者仿真人体，或者参考人体，获得粒子数与不确定度关系，并通过插值或者拟合获得粒子数与不确定度的一一映射关系，从而得到估算结果。估算结果可以是一次模拟粒子输运过程的入射粒子总数。

所输入的粒子是基于不同类型的源产生的，比如，光子类型的源、电子类型的源或质子类型的源，在本实施例中，是对入射粒子输入手段进行了限定的，具体是：

将类型相同，或能量接近，或类型相同并且能量接近的入射粒子分为一类，并基于每一类的入射粒子进行该类粒子的分批处理；

基于入射粒子的类别及批次，将入射粒子按照其所产生的源的分布情况交错地分批输入。

例如，以光子类型与电子类型的混合源为例，将该混合源分为M批次，对于每一批次，都有光子与电子，其比例不变，交错运行，直到完成所有批次计算；每一批次的相同类型中，还可以按照能量，分为N子批次。

本实施例对入射粒子按能量与类型进行分类，将能量接近及类型相同的粒子作为同一子批次，可减少并行计算过程中多个线程之间互相等待的时间，从而加快并行计算的模拟速度。

本实施例还对不同粒子类型按照源的分布情况交错运行，能够降低模拟计算的不确定度，在提高并行计算速度的前提下，有利于降低引入误差。

当然，在其他实施例中，也可以在上述输入手段上做适应性变更，如仅采用将能量接近及类型相同的粒子作为同一批次，由此直接分批输入；或者，直接不同类型的源产生粒子交错输入，或者部分粒子直接输入，部分粒子按照本实施例的上述手段进行优化输入，其都是可被实施的。

继续参考图1，本实施例模拟粒子输运的方法还包括：

步骤S101，记录所输入粒子的输运径迹。

在采用蒙特卡罗方法对输入粒子进行模拟时，其是根据不同粒子的反应类型特征来抽样粒子输运情况的，本实施例将一个粒子被抽样的输运情况称为所述输运径迹，输运径迹是描述该粒子被抽样的物理反应类型的信息集合，由于粒子可能会被抽样的物理反应类型的基本物理参数及模型已预先存入蒙特卡罗工具，因此，在具体的执行过程中，需要根据实际粒子的类型、粒子能量、粒子速度、粒子所处位置的材料性质来调用对应被抽样的模型及基本物理参数，从而得到粒子的输运径迹。

在上述论述中，所述基本物理参数可以是指所述物理反应的微分散射截面、平均自由程等类似参数，所述模型可以是指描述光子的光电效应、康普顿散射、对反应等类似模型。

基于上述，在本实施例中，所述输运径迹包括粒子被抽样物理反应类型的一切信息，如：被抽样的物理反应类型的基本物理参数及模型、入射粒子的能量信息、速度信息及其他径迹信息。其中，所述速度信息包括所述入射粒子的入射方向信息，所述其他径迹信息包括：所述入射粒子的类型信息、入射位置信息、权重信息、所经栅元信息及相对于所经栅元的不确定度信息。

基于粒子被抽样的情况来看，本实施例步骤S101的记录手段包括粒子的抽样及粒子输运径迹的存储，由于上述粒子被抽样时，主要涉及物理反应类型的基本物理参数及模型的抽样，在上述信息被确定的时候，本实施例还对其记录手段进行了简化：

若入射粒子与已记录粒子的能量信息与入射方向信息接近，则复制已记录粒子的输运径迹，赋值为该入射粒子的输运径迹。

当然，在其他实施例中，也可以仅凭输入粒子与已记录粒子的能量信息与入射方向信息接近，直接复制物理反应类型的基本物理参数及模型的抽样结果及已记录粒子的输运径迹，赋值为该入射粒子的输运径迹。

本实施例涉及的简化记录手段能够基于相同的入射粒子(相同入射粒子指能量信息与输入方向信息接近的粒子)，对重复径迹的粒子进行复制而避免重复抽样，能够减少总体粒子的抽样数目，进一步提高蒙特卡罗工具的模拟算法效率。

对于非直接复制输运径迹的入射粒子，本实施例的对该粒子的抽样及输运径迹的存储，都是基于随机数确定粒子的上述与抽样物理反应有关的信息，并用于基于蒙特卡罗的输运模拟，从而生成该粒子的输运径迹。

基于现有技术还可知，蒙特卡罗模拟粒子输运的过程需要模拟大量粒子，以满足模拟算法对不确定度的要求，故而算法本身非常耗时。继续参考图1，本实施例还提供了如下步骤，能够基于模拟过程的全局不确定度减少抽样粒子数，从而确保算法全局不确定度的均衡：

步骤S102，基于每批次运行粒子的径迹计算每个栅元的不确定度，若栅元的不确定度不超过栅元阈值，则该栅元为达标栅元。

步骤S103，获取感兴趣区域中栅元的达标率，所述感兴趣区域至少包括一个栅元，所述感兴趣区域的达标率为该区域达标栅元占该区域所有栅元的比例。

所述栅元阈值为栅元达标率的评估，使栅元的达标率符合于预定要求。

在每一批入射粒子完成抽样及完成存储粒子输运径迹后可执行步骤S102及S103，其中，栅元不确定度是对与粒子发生截面反应的几何栅元的不确定度的评估。由于考虑到蒙特卡罗模拟时，其对于不确定度的要求可能有所不同，故而，设置感兴趣区域的达标率的评估，使上述不确定度的均衡可基于各感兴趣区域的达标率的均衡。

当各感兴趣区域的达标率都超过预先设定的阈值(图1中为了区别，将该阈值称为感兴趣区域阈值)时，则评估各感兴趣区域的达标率符合要求，可以对蒙特卡罗模拟过程进行截断，即停止继续输入粒子，并输出本地模拟结果(包括已完成记录的粒子输运径迹)。上述感兴趣区域的不确定度的评估，能够对全局不确定度进行均衡，并在满足全局不确定度的要求下，减少了不必要的粒子抽样输运，可节省大量耗时，大大提升蒙特卡罗模拟效率。

另外，本实施例还涉及对一个栅元的不确定度进行计算，评估一个栅元的不确定度基于不同的不确定度标准可能是具有多种定义的，本实施例公开如下两种评估方式，以供参考：

一种评估方式为：

基于所述每批次运行粒子的径迹计算粒子与每个栅元相互作用的截面数据，所述截面数据包括粒子在该栅元中实际发生反应的概率；

比较所述栅元中粒子实际发生反应的概率与预期发生反应的概率，以得到该栅元的不确定度。

在上述栅元不确定度的评估方式中，是以粒子截面反应为考量标准的，从现有技术可知，粒子截面反应是指一个粒子与栅元截面发生反应的概率，所述反应的发生具体是一个粒子被所经过的栅元吸收能量，因吸收能力的不同，则发生的反应不同；根据不同的粒子反应，可以比较相对于一个栅元而言，所预期的粒子发生一种反应的概率与实际粒子发生该种反应的概率，从而得到该种反应发生概率的不确定度(比如上述预期概率与实际概率之差)，基于栅元可能会发生的多种反应发生概率不确定度的评估函数，可获得该栅元的不确定度。所述评估函数可以是和函数或均值函数等。

另一种评估方式为：

基于所述每批次运行粒子的径迹计算粒子在每个栅元中的粒子数密度；

根据所计算栅元及其粒子数密度序列拟合得到有关栅元中粒子数密度的分布曲线；

基于所述分布曲线确定各个栅元的不确定度。

上述栅元不确定度的评估方式中，是以经过每个栅元的粒子数密度为考虑标准的，所述栅元粒子数密度指的是经过该栅元的入射粒子个数。在一批输入粒子运行完毕后，可从该批粒子的输运径迹来获得各栅元的粒子数密度，并得到栅元与粒子数密度的对应关系(表现在二维坐标上则为一串离散序列)，基于该对应关系可拟合得到栅元粒子数密度的分布曲线，比较所得分布曲线上的栅元粒子数密度值与实际输运过程中栅元的粒子数密度，即可得到栅元的粒子数密度差值，由此确定各个栅元的不确定度。

当各感兴趣区域的达标率还没有符合所述预先设定的阈值时，则继续产生并输入下一批入射粒子，继续进行重复步骤S101～S103，直到各感兴趣区域的达标率都符合预先设定的阈值，或者依据估算入射粒子总数，所有粒子已经基于蒙特卡罗模拟运行完毕。在后者的情况下，也会最终输出本地模拟结果(包括全部完成记录的粒子输运径迹)。

实施例二

本实施例基于实施例一，提供了一种如图2所示的模拟粒子输运的方法，可根据模拟过程全局不确定度，对不确定度进行均衡，且根据不确定度的分布动态调整粒子抽样数目。该方法具体包括如下步骤：

步骤S100～S103，与实施例一一致。

当下一批入射粒子输入时，基于入射粒子的输运过程，对输入粒子做如下处理：

步骤S104，若入射粒子的输运过程显示，该批粒子中有一个入射粒子从重要性低的栅元进入重要性高的栅元，则对进入该重要性高的栅元的入射粒子以第一概率进行分裂，分裂后的粒子的输运轨迹可复制分裂前的粒子，但其粒子权重降低以使该批粒子的总权重不变。

步骤S105，若入射粒子的输运过程显示，该批粒子中有一个入射粒子从重要性高的栅元进入重要性低的栅元，则进入该重要性低的栅元的入射粒子以第二概率被杀死，未被杀死的粒子的输运轨迹不变但增加其粒子权重以使该批粒子的总权重不变。

为了简单起见，在图2中，将步骤S104中所述重要性低的栅元定义为第一重要性栅元，将重要性高的栅元定义为第二重要性栅元；而步骤S105中所述重要性高的栅元定义为第三重要性栅元，重要性低的栅元定义为第四重要性栅元，但可以理解的是，上述第一至第四重要性栅元并非一定是不同的重要性栅元，其仅是为了区分栅元重要性的高低而采用的一种相对概念，第一至第四重要性栅元所指重要性栅元的定义范围在实际执行中可能重叠。

所述重要性栅元和所述重要性可以为预先手动设定，或根据栅元信息自动设定；所述栅元信息包括栅元的不确定度，或者栅元的物理属性。

在具体实施中，所述重要性栅元和所述重要性可以在进行放射治疗前由操作人员通过放射治疗系统附带的应用软件预先手动设定，例如：可以将感兴趣区域中的若干个可能有肿瘤的区域设定为重要性栅元，并分别设定这些重要性栅元的重要性(具体的值)。

所述重要性栅元的重要性也可以根据栅元的不确定度自动设定，一般情况下，栅元的不确定度越高，则所述栅元的重要性也越高。

步骤S104及S105中，涉及的第一概率为第一重要性栅元与第二重要性栅元的重要性比值，所述第二概率为第四重要性栅元与第三重要性栅元的重要性比值。

在本实施例中还需要强调的是，在入射粒子从重要性低的栅元进入重要性高的栅元，设置粒子以第一概率分裂，使本发明技术方案能够动态增加输运粒子的数目，降低重要性较高的栅元的不确定度，使各栅元的不确定度动态地变化至均衡，从而提高粒子输运模拟的精度，上述手段也可以使各感兴趣区域在实施例一的方案中尽快达到均衡，故而可减少粒子运输的批次，也进而提高模拟效率；而入射粒子从重要性高的栅元进入重要性低的栅元，设置粒子以第二概率分裂，使本发明技术方案能够在保证粒子模拟精度的前提下动态减少输运粒子的数目，进一步提高模拟效率。

实施例三

本实施例基于实施例一，提供了一种如图3所示的模拟粒子输运的方法，可根据模拟过程对粒子的分批输运，在模拟完成预定批次粒子的输运径迹后，基于历史批次粒子的输运径迹对全局不确定度进行动态降噪，从而降低模拟粒子输运过程的不确定度，有助于减少粒子抽样数目及提高模拟效率。该方法具体包括如下步骤：

步骤S100～S103，与实施例一一致。

在模拟粒子输运的过程中，若粒子输运的批次达到预定数目，则执行：

步骤S106，对入射粒子的剂量分布做动态降噪处理。

在本实施例中，所述动态降噪处理具体通过如下方式实现：

得到每个栅元中粒子剂量分布的三维曲线(也即所述栅元中粒子三维的剂量分布)及所述剂量分布对应的不确定度，所述每个栅元中粒子的剂量分布可以根据模拟获得的粒子在栅元中的能量分布获得，具体为单位质量粒子的能量分布；

对上述三维曲线进行滤波处理，以使所述三维曲线在三个维度上连续可导；

重新计算所述剂量分布对应的不确定度。

步骤S107，输出动态降噪处理后得到的每个栅元剂量分布的不确定度。

在具体实施中，根据滤波处理后的剂量分布的三维曲线，可计算得到每个栅元剂量分布的不确定度。

上述滤波处理可以消除上述剂量分布三维曲线上的毛刺等噪声，平滑三维曲线；再根据平滑后的三维曲线计算栅元中粒子的剂量分布的不确定度，可对模拟过程进行整体降噪，从而有助于下一次粒子输运的模拟。

实施例四

本实施例提供了一种模拟粒子输运的方法，其对模拟过程中模拟对象(所述模拟对象可为放射治疗领域中的各类治疗机，比如加速器治疗头等)的几何模型进行了特别限定，本实施例包括步骤：

导入几何模型；以及其他模拟步骤。其中，所述其他模拟步骤可参考实施例一至实施例三任一种技术方案中所记载的步骤流程。

所述几何模型包括：模拟对象的栅元、物理材料、栅元权重或/和几何虚拟截面，其中，所述几何虚拟截面适于定义所述栅元对应物理材料，以使所述栅元指向的模拟对象具有均匀化材料，所输入粒子的输运径迹与所述几何虚拟截面有关。

本实施例关注的是，上述模拟中的几何虚拟截面，所述几何虚拟截面为粒子进行截面反应的抽样概率，所述截面反应包括真实反应与虚拟反应，对于一粒子的输运，设所定义几何虚拟截面为∑_max，在一次抽样中粒子发生真实反应的概率为∑_r、虚拟反应的概率为∑_r′，r＝1～R；其中，r为自然数，代表该粒子被抽样的次数，R为大于或等于1的自然数，代表粒子被有限次抽样的具体数目；所述几何虚拟截面∑_mxa满足如下条件：

∑_max＝max(∑₁，∑₂，…，∑_R-1，∑_R)；

∑_max＝∑₁+∑′₁＝∑′₂+∑′₂＝…＝∑_R+∑′_R；

从上述条件来看，所述截面反应的抽样概率为粒子有限次被抽样经过所述栅元的最大概率值。在满足上述条件时，对整个模型进行粒子输运过程的模拟下，所有材料均采用∑_max的几何虚拟截面进行输运长度抽样，每种材料发生真实反应概率为∑_r/∑_max，发生虚拟反应的概率为∑′_r/∑_max，因此采用虚拟截面输运时能够与真实物理输运过程保持高度一致性。

结合图4，图4为粒子输运过程中发生截面反应的示意图，其中，箭头方向为粒子的入射方向，白色圆形代表一个几何栅元，黑色圆形代表粒子在该处发生虚拟反应(虚碰撞)，而阴影圆形代表粒子在该处发生真实反应(真实物理反应，比如光电效应、康普顿效应等)，粒子的路径为1-2-3-4-5-6-7-8，其中，1、5、8处代表真实反应，2、3、4、6、7代表虚拟反应，可以看到发生虚碰撞反应后，只需重新抽样输运程度即可，不需要重新抽样粒子方向和能量，只有在真实反应时才需要重新进行粒子方向和能量抽样。

因此，本实施例的几何虚拟截面不同于现有技术，其是通过如下方式定义的：

所述粒子以所述截面反应的抽样概率被抽样，以经过所述栅元；

所述截面反应的抽样概率为所述真实反应的抽样概率与虚拟反应的抽样概率之和；

粒子与栅元进行真实反应时的抽样包括输运程度、方向及能量的抽样，而粒子与栅元进行虚拟反应时的抽样仅包括输运程度的抽样。

实施例五

本实施例提供了一种计算放射治疗中人体剂量的方法，包括步骤：

模拟粒子输运以获得的粒子在栅元中的能量分布；

基于粒子在栅元中的能量分布计算放射治疗中的人体剂量。

其中，所述模拟粒子输运的方式可采用实施例一至实施例四中任一种方法。

实施例六

本实施例提供了一种模拟粒子输运的装置，对应于实施例一，适于模拟粒子在栅元中的能量分布，包括：源处理模块、输运处理模块、噪声处理模块及输出模块。

所述源处理模块，适于执行步骤S100；

所述输运处理模块，适于执行步骤S101；

所述噪声处理模块，适于执行步骤S102及S103，并通过所述输出模块输出本地模拟结果。

在其他实施例中，上述模拟粒子输运的装置还可以对应实施例二，不同于本实施例之处在于，所述噪声处理模块还适于执行步骤S104及S105。

在其他实施例中，上述模拟粒子输运的装置还可以对应实施例三，不同于本实施例之处在于，所述噪声处理模块还适于执行步骤S106及步骤S107。

实施例七

本实施例提供了一种模拟粒子输运的装置，对应于实施例四，适于模拟粒子在栅元中的能量分布，包括：源处理模块、输运处理模块、噪声处理模块、输出模块、输入模块及几何处理模块。

所述源处理模块，适于执行步骤S100；

所述输运处理模块，适于执行步骤S101；

所述噪声处理模块，适于执行步骤S102及S103，并通过所述输出模块输出本地模拟结果；

所述输入模块适于导入实施例四所述几何模型至所述几何处理模块。

实施例八

本实施例提供了一种放射治疗系统，包括：

模拟粒子输运的装置，其适于模拟粒子在栅元中的能量分布，其具体结构可参考实施例六及实施例七中的论述内容。

计算人体剂量的装置，其适于根据所述模拟粒子输运的装置获得的粒子在栅元中的能量分布计算放射治疗中的人体剂量。

本发明虽然已以较佳实施例公开如上，但其并不是用来限定本发明，任何本领域技术人员在不脱离本发明的精神和范围内，都可以利用上述揭示的方法和技术内容对本发明技术方案做出可能的变动和修改，因此，凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化及修饰，均属于本发明技术方案的保护范围。

Claims (16)

1.一种模拟粒子输运的方法，适于模拟粒子在栅元中的能量分布，其特征在于，包括：

估算所需入射总粒子数，产生入射粒子并分批输入；

记录所输入粒子的输运径迹；

基于每批次运行粒子的径迹计算每个栅元的不确定度，若栅元的不确定度不超过第一阈值，则该栅元为达标栅元；

获取感兴趣区域中栅元的达标率，所述感兴趣区域至少包括一个栅元，所述感兴趣区域的达标率为该区域中达标栅元占该区域所有栅元的比例；

若所述感兴趣区域中栅元的达标率超过第二阈值，则停止继续输入粒子，并输出历史输入粒子的输运径迹，否则继续输入粒子，直到总粒子数运行完毕。

2.如权利要求1所述的方法，其特征在于，分批输入的每一批入射粒子的能量接近，或类型相同，或能量接近并且类型相同。

3.如权利要求2所述的方法，其特征在于，所述分批输入包括：

根据所述入射粒子的类型不同交错地分批输入。

4.如权利要求1所述的方法，其特征在于，所述输运径迹包括：入射粒子的能量信息、速度信息及其他径迹信息，所述速度信息包括所述入射粒子的入射方向信息，所述记录所入射粒子的输运径迹包括：

若入射粒子与已记录粒子的能量信息与入射方向信息接近，则将该已记录粒子的输运径迹赋值为该入射粒子的输运径迹。

5.如权利要求4所述的方法，其特征在于，所述其他径迹信息包括：所述入射粒子的类型信息、入射位置信息、权重信息及所经栅元信息；

所述所经栅元信息包括：所经栅元的能量分布以及不确定度。

6.如权利要求1所述的方法，其特征在于，还包括：当每一批粒子输入时：

若其中一个粒子从第一重要性栅元进入第二重要性栅元，则所进入的粒子以第一概率进行分裂，分裂后的粒子降低权重以使该批粒子的总权重不变，所述第一重要性栅元的重要性低于第二重要性栅元的重要性；

若其中一个粒子从第三重要性栅元进入第四重要性栅元，则所进入的粒子中以第二概率被杀死，未被杀死的粒子增加权重以使该批粒子的总权重不变，所述第三重要性栅元的重要性高于第四重要性栅元的重要性。

7.如权利要求6所述的方法，其特征在于，所述第一概率为第一重要性栅元与第二重要性栅元的重要性比值，所述第二概率为第四重要性栅元与第三重要性栅元的重要性比值。

8.如权利要求6或7所述的方法，其特征在于，所述重要性栅元的重要性为预先手动设定，或根据栅元信息自动设定；所述栅元信息包括栅元的不确定度，或者栅元的物理属性。

9.如权利要求1所述的方法，其特征在于，还包括：

基于历史输入粒子的输运径迹对入射粒子的剂量分布的不确定度做动态降噪处理。

10.如权利要求9所述的方法，其特征在于，所述动态降噪处理通过如下方式实现：

得到所述栅元中粒子的三维的剂量分布和不确定度；

对所述三维的剂量分布进行滤波处理，以使所述剂量分布在三个维度上连续可导；

重新计算得到所述剂量分布对应的不确定度。

11.如权利要求1所述的方法，其特征在于，还包括：导入几何模型，所述几何模型包括：模拟对象的栅元、物理材料、栅元权重或/和几何虚拟截面，其中，所述几何虚拟截面适于定义所述栅元对应物理材料，以使所述栅元指向的模拟对象具有均匀化材料，所输入粒子的输运径迹与所述几何虚拟截面有关。

12.一种计算放射治疗中人体剂量的方法，其特征在于，包括：根据如权利要求1～11任一项所述的模拟粒子输运的方法获得的粒子在栅元中的能量分布以计算放射治疗中的人体剂量。

13.一种模拟粒子输运的装置，适于模拟粒子在栅元中的能量分布，其特征在于，包括：源处理模块、输运处理模块、噪声处理模块及输出模块；

所述源处理模块，适于估算所需入射总粒子数、产生入射粒子并分批输入；

所述输运处理模块，适于记录所输入粒子的输运径迹；

所述噪声处理模块，适于实现如下步骤：

基于每批次运行粒子的径迹计算每个栅元的不确定度，若栅元的不确定度不超过第一阈值，则该栅元为达标栅元；

获取感兴趣区域的达标率，所述感兴趣区域至少包括一个栅元，所述感兴趣区域的达标率为该区域达标栅元占该区域所有栅元的比例；

若所述感兴趣区域中栅元的达标率超过第二阈值，则停止继续输入粒子，并通过所述输出模块输出历史输入粒子的输运径迹，否则使继续输入粒子，直到总粒子数运行完毕。

14.如权利要求13所述的装置，其特征在于，所述噪声处理模块还适于实现如下步骤：所述当下一批粒子输入时：

若其中一个粒子从第一重要性栅元进入第二重要性栅元，则所进入的粒子以第一概率进行分裂，分裂后的粒子降低权重以使该批粒子的总权重不变，所述第一重要性栅元的重要性低于第二重要性栅元的重要性；

若其中一个粒子从第三重要性栅元进入第四重要性栅元，则所进入的粒子中以第二概率被杀死，未被杀死的粒子增加权重以使该批粒子的总权重不变，所述第三重要性栅元的重要性高于第四重要性栅元的重要性；

所述重要性栅元的重要性为预先手动设定，或根据栅元信息自动设定；所述栅元信息包括栅元的不确定度，或者栅元的物理属性。

15.如权利要求13所述的装置，其特征在于，还包括：输入模块及几何处理模块；所述输入模块适于导入几何模型至所述几何处理模块；

所述几何模型包括模拟对象的栅元、物理材料、栅元权重及几何虚拟截面，其中，所述几何虚拟截面适于定义所述栅元对应物理材料，以使所述栅元指向的模拟对象具有均匀化材料，所输入粒子的输运径迹与所述几何虚拟截面有关。

16.一种放射治疗系统，其特征在于，包括：

如权利要求13～15任一项所述的模拟粒子输运的装置，适于模拟粒子在栅元中的能量分布；

计算人体剂量的装置，适于根据所述模拟粒子输运的装置获得的粒子在栅元中的能量分布计算放射治疗中的人体剂量。