一种复杂射野的放疗剂量快速计算方法、设备和存储介质与流程

本发明属于放疗技术领域，涉及一种复杂射野的放射剂量快速计算方法、设备和存储介质。

背景技术：

在放疗技术领域，在用加速器肿瘤治疗时，通常要在射束路径上加楔形板过滤器以改变平野(即不加楔形板)的剂量分布。楔形板是用高密度材料做成的，它能使射野的输出剂量形成所需分布，以得到合适的处方剂量。由于在放射治疗中不可能在病人体内植入探测器来测定剂量，只能通过复杂的理论计算来确定病人体内的剂量分布。因此，在精确放射治疗技术的发展中，对剂量算法的研究在实现技术目标中具有特殊的重要地位。精确和快速的剂量算法和调强剂量算法的研究，已经成为适形和调强放射治疗技术的主攻内容和成功的关键。

《楔形野剂量计算中的误差分析和修正》(中国医学物理学杂志，第18卷第2期，2001年4月)提到的传统方法和改进方法都是采用公式计算的方法，其中没有改进的公式计算结果误差可以高达11％，而改进后的公式计算结果误差可以降低到1％，然而该方法需要有实验测量的多种参数，工作量较大，属于半经验的方法，对于垂直照射的情况、照射角度变化的情况的计算分析较为复杂。

同样，《一种放射治疗剂量的快速定量方法》(cn105288873a)，该方法采用的带楔形过滤器的射野的剂量计算是通过测量楔形因子得到；该方法的实验工作量大，分析复杂。

《关于x和γ射线束楔形板剂量的修正算法》(四川大学学报(自然科学版)，第40卷第5期，2003年10月)，提到的基于射线衰减修正的方法计算剂量，然而，并没有直接进行楔形野的剂量计算。而是通过衰减路径长度和衰减比例求解楔形因子，然后修正平野剂量，从而得到楔形野下的剂量分布。但是加速器的射束通常不是单能的，往往是一个能谱分布；而该文章并没有提及如何实现能谱的修正。

《在适形放疗中精确确定辐射野输出剂量的方法》(cn100431642c)和《动态楔形板控制点的生成方法及装置》(cn105105780a)也都只是提到有楔形过滤器的情况下采用质量衰减修正，采用网格划分的方法，依然没有提到存在能谱分布的情况如何修正。

另外，上述几种方法也都没有提到如何根据网格和线程数量进行并行计算的任务分配方法。

对于本发明涉及的复杂野，目前尚未见到快速的计算方法，现有技术的方法通常通过定义二维权重平面应用于三维剂量结果的修正，当进行三维剂量计算时，该二维权重需要反复微调，特别对于边缘处的权重需要采用梯度渐变权重模式，该方法费事费力，过度依赖于经验。

技术实现要素：

本发明的目的在于为克服上述现有技术的缺陷而提供一种复杂射野的放射治疗剂量快速计算方法、设备和存储介质。

如图1所示，本发明中的复杂射野包括为不规则射野、非对称射野、斜入射野或楔形野中的一种或一种以上的混合射野，其中非规则射野进一步包括带有挡块的射野和设置有多叶准直器的射野，本发明中的非对称射野为射野中心轴偏离射束中心轴的射野。

为实现上述目的，本发明采用以下技术方案：

一种复杂射野的放射治疗剂量快速计算方法，适于在计算设备中执行，包括如下步骤：

(1)定义限束装置的参数：

(2)射野网格化操作：

将射野均匀网格化，网格所在平面为射野平面；对权重相同的相邻网格进行网格合并；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给每个计算单元(线程)；

(3)结合单个网格对应的虚拟jaws坐标，对经过虚拟jaws的单个粒子进行剂量计算，得到单个线程内单个网格内单个粒子的剂量结果；

根据公式(1)，将单个网格内每个粒子产生的三维剂量按照体素进行叠加得到单个网格内所有粒子的三维剂量结果；

其中，d1(i,j,k)为单个网格产生的三维剂量，

m为进入该网格的粒子总数，m为该网格内的第m个粒子，

为该网格内单个粒子的剂量；

然后根据公式(2)，将单个线程下每个网格的三维剂量按照体素(用于统计剂量小体积单元或小质量单元)叠加，得到单个线程下的所有网格的总三维剂量；

其中，d2(i,j,k)为单个线程产生的三维剂量，

n为进入该线程的网格总数，n为该线程的第n个粒子，

为该线程内单个网格的剂量；

最后，根据公式(3)，将所有单个线程的三维剂量按照体素叠加得到复杂射野的剂量计算结果，

其中，d3(i,j,k)为该复杂射野的放射治疗剂量，

p为进行复杂射野剂量计算的线程总数，p为第p个线程，

为该复杂射野内单个线程的剂量。

所述的复杂射野是由包括多叶准直器(mlc)的限束装置构造的，优选地，该复杂射野中的限束装置还包括钨门(jaws)、楔形板或挡块中的一种或一种以上。

步骤(1)中，所述的定义限束装置的参数是从dicom(digitalimagingandcommunicationsinmedicine，即医学数字成像和通信，是医学图像和相关信息的国际标准iso-12052，它定义了质量能满足临床需要的可用于数据交换的医学图像格式)读取限束装置的参数；所述的限束装置的参数包括mlc的控制点坐标、开口以及不同时刻的mlc的开口权重(快照)；优选地，所述的限束装置的参数还包括：钨门(jaws)的几何材料、坐标位置，楔形板或挡块的几何尺寸、材料等。

步骤(2)中，所述的网格化操作是将射野划分为若干网格，其中网格化操作的参数是射野大小以及每个小束流的大小，后者是自己根据需求定义的；优选地，划分的网格为大小均匀的网格；优选地，网格尺寸与探测器的尺寸一致。

步骤(2)中，网格的权重是将不同时刻的同一网格的坐标叠加求得的。

步骤(2)中，所述的计算单元为gpu计算单元或其他可编程计算芯片的计算单元，其中可编程计算芯片优选为现场可编程门阵列(field-programmablegatearray，fpga)。

步骤(3)中，采用蒙特卡罗计算模型对粒子进行剂量计算。

本发明还提供一种复杂射野的放射治疗剂量快速计算方法，适于在计算设备中执行，包括如下步骤：

(1)定义限束装置的参数；

(2)射野网格化：

对射野进行网格化操作，将射野划分为若干网格；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程)分配给每个计算单元(线程)；

根据放射源的位置和网格的划分尺寸和网格中心坐标，采用求交算法或投影定理计算出放射源粒子经过限束装置的路径长度；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

(3)根据限束装置的材料参数，限束装置的质量衰减系数(可从技术手册中查出)，结合放射源粒子经过限束装置的路径长度以及射野网格划分情况，计算出到达每个射束网格的粒子经过限束装置的衰减比；

(4)通过步骤(3)得到的衰减比作用于无限束装置时到达每个网格的粒子能谱(这些粒子是指从放射源，穿过限束装置到达网格的粒子)；结合虚拟jaws坐标，对射束能谱衰减后的射束粒子进行剂量计算，计算的剂量为每个线程中单个网格的剂量；将单个网格剂量进行衰减修正后求和得到单个线程的剂量结果；最后将所有线程的剂量求和得到复杂射野放射治疗剂量结果。

步骤(1)中，所述的限束装置选自楔形板或挡块中的一种或一种以上。所述的楔形板的几何尺寸能全部覆盖射野范围；所述的挡块的几何尺寸小于射野大小。优选地，所述挡块的材料为高原子序数的材料，例如pb(铅)、w(钨)等。

所述的限束装置的参数包括楔形板或挡块的材料、密度、放置位置或几何尺寸等参数中的一个或一个以上。

所述的步骤(2)中，网格大小可以根据实际测量所用探测器的尺寸来自定义；其中每个网格代表一个射束照射范围，网格所在平面为射野平面。

所述的步骤(3)中，针对到达每个网格的射线经过楔形板的路程不同，对每个网格入射粒子的能谱进行相应的衰减。

步骤(4)中，射束能谱衰减后基于蒙特卡罗计算模型进行剂量计算。

步骤(4)中，所述的衰减修正为三维剂量的二维衰减修正。

本发明还提供了一种复杂射野的放射治疗剂量快速计算方法，适于在计算设备中执行，包括如下步骤：

(1)定义非对称射野或斜入射野的参数；

(2)射野网格化操作：

对射野进行均匀网格化操作，将非对称射野或斜入射野划分为若干网格；其中每个网格代表一个射束照射范围，网格所在平面为射野平面；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程)分配给每个计算单元(线程)；

根据放射源的位置、非对称射野或斜入射野的参数、网格的划分尺寸和网格中心坐标，采用求交算法或投影定理计算出每个网格对应的虚拟钨门jaws的坐标参数；

(3)采用蒙特卡罗计算模型对粒子进行剂量计算，得到每个线程中单个网格的剂量结果；

然后将单个网格的剂量进行叠加；将单个网格剂量求和得到单个线程的剂量；最后将所有单个线程的剂量求和得到复杂射野放射治疗剂量结果。

所述的非对称射野的参数包括非对称射野的几何中心、射野大小；斜入射野的参数包括斜入射野的角度。

本发明还提供了一种放射剂量的快速计算设备，包括：

一个或多个处理器；

存储器；以及

一个或多个程序，其实所述一个或多个程序存储在所述存储器中并被配置为由一个或多个处理器执行，所述的一个或多个程序包括用于上述复杂射野的放射治疗剂量快速计算方法的指令。

本发明还提供了一种存储一个或多个程序的计算机可读存储介质，所述的一个或多个程序包括指令，所述指令适于由存储器加载并执行上述复杂射野的放射治疗剂量快速计算方法。

本发明具有以下有益效果：

本发明基于并行计算方法和射束分配流程，实现了一种复杂射野的快速高效计算方法。通过将原来的单能或多能射束到达体模的或人体的射野平面划分多个网格，将网格等数量或不等数量的指派给每个并行的线程，实现了不损失精度的情况下，快速高效地进行复杂射野的模拟计算。

将网格化射野的方法以及衰减修正的方法应用于挡块和楔形板的计算，可以大大减少因为存在楔形板挡块导致的粒子输运深穿透问题。同样，将网格化计算剂量的方法应用于不对称射野和斜入射射野，以及mlc的动态或静态调强的剂量计算，经过实际案例计算验证表明，方法是有效且高效的。

本发明提供的复杂射野的快速高效计算方法使用了网格化操作将剂量计算化整为零，可以更方便的应用于多线程并行计算，例如本发明的技术方法可以应用于gpu、fpga或者其他并行计算平台。

附图说明

图1为本发明复杂射野分类示意图。

图2为本发明一个优选的实施例中带mlc的不规则射野示意图。

图3为本发明一个优选的实施例中带mlc的不规则射野的放射治疗剂量快速计算方法流程图。

图4为本发明一个优选的实施例中带mlc的不规则射野剂量计算网格合并图。

图5为本发明另一个优选的实施例中由mlc和钨门构造的不规则射野的示意图。

图6为本发明又一个优选的实施例中楔形野的放射治疗剂量快速计算方法流程图。

图7为本发明又一个优选的实施例中楔形板的结构与射束照射示意图。

图8为本发明又一个优选的实施例中射野网格划分示意图。

图9为本发明又一个优选的实施例中楔形野水箱中的离轴剂量分布对比图：

其中曲线(1)为实测的水箱深度为1.6cm的离轴剂量分布(点画线)；

曲线(1’)为利用本发明实施例1的计算方法模拟的水箱深度为1.6cm的离轴剂量分布(实线)；

曲线(2)为实测的水箱深度为5.0cm的离轴剂量分布(点画线)；

曲线(2’)为本发明实施例1的计算方法模拟的水箱深度为5.0cm的离轴剂量分布(实线)；

曲线(3)为实测的水箱深度为10.0cm的离轴剂量分布(点画线)；

曲线(3’)为利用本发明实施例1的计算方法模拟的水箱深度为10.0cm的离轴剂量分布(实线)；

曲线(4)为实测的水箱深度为20.0cm的离轴剂量分布(点画线)；

曲线(4’)为利用本发明实施例1的计算方法模拟的水箱深度为20.0cm的离轴剂量分布(实线)。

图10为本发明又一个优选的实施例中挡块构建的l形射野示意图，其中图(a)为l形射野的左视图；图(b)为l形射野的俯视图。

图11为本发明又一个优选的实施例中非对称射野示意图。

图12为本发明又一个优选的实施例中非对称野剂量计算结果示意图；

其中图(a)为离轴剂量分布图(左视图)；

图(b)为水箱三维剂量分布俯视图；

图(c)为水箱三维剂量分布左视图；

图(d)为水箱三维剂量分布主视图。

图13为本发明又一个优选的实施例中斜入射中机头位移示意图。

图14为本发明又一个优选的实施例中带楔形板的斜入射复杂野剂量计算结果示意图；

其中图(a)为离轴剂量分布图(左视图)；

图(b)为水箱三维剂量分布俯视图；

图(c)为水箱三维剂量分布左视图；

图(d)为水箱三维剂量分布主视图。

具体实施方式

以下结合实施例和附图进一步说明本发明。

实施例1mlc构造的不规则野(不带jaws)

本实施例示出了带有mlc的不规则射野(如图2)的放射治疗剂量快速计算方法，该不规则射野为等中心射野，包括如下步骤(如图3所示)：

(1)定义限束装置mlc的参数110：

从dicom读取mlc的参数，包括mlc的控制点坐标、开口以及不同时刻的mlc的开口权重(快照)；两个相对的mlc叶片之间的开口权重是通过tps(放射治疗计划系统)优化时形成的，本实施例中不考虑片间漏射的情况；

(2)射野网格化操作：

将包含mlc开口形状的矩形区域(如图2所示)均匀网格化121，网格所在平面为射野平面；

对权重相同的相邻网格进行网格合并122(如图4所示)，其中网格的权重是将不同时刻的同一网格的坐标叠加求得的；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标123；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给每个gpu计算单元(或cpu线程)124；

(3)结合单个网格对应的虚拟jaws坐标，对经过虚拟jaws的单个粒子进行剂量计算，得到单个线程内单个网格内单个粒子的剂量结果131；

根据公式(1)，将单个网格内每个粒子产生的三维剂量按照体素进行叠加得到单个网格内所有粒子的三维剂量结果132，

其中，d1(i,j,k)为单个网格产生的三维剂量，

m为进入该网格的粒子总数，m为该网格内的第m个粒子，

为该网格内单个粒子的剂量；

然后根据公式(2)，将单个线程下每个网格的三维剂量按照体素(用于统计剂量小体积单元或小质量单元)叠加，得到单个线程下的所有网格的总三维剂量133；

其中，d2(i,j,k)为单个线程产生的剂量，

n为进入该线程的网格总数，n为该线程的第n个粒子，

为该线程内单个网格的剂量；

最后，根据公式(3)，将所有单个线程的剂量按照体素叠加，得到复杂射野的剂量计算结果134，

其中，d3(i,j,k)为该复杂射野的放射治疗剂量，

p为进行复杂射野剂量计算的线程总数，p为第p个线程，

为该复杂射野内单个线程的剂量。

实施例2mlc和jaws构造的不规则野

本实施例示出了由mlc和钨门(jaws)共同构造的不规则射野(如图5)的放射治疗剂量快速计算方法，包括如下步骤：

(1)定义限束装置mlc和钨门的参数：

从dicom读取mlc的参数，包括mlc的控制点坐标、开口以及不同时刻的mlc的开口权重(快照)，两个相对的mlc叶片之间的开口权重是通过tps(放射治疗计划系统)优化时形成的；钨门的参数包括钨门的坐标；本实施例中不考虑片间漏射的情况。

(2)射野网格化操作：

将jaws与mlc共同形成的开口形状均匀网格化，网格所在平面为射野平面；对权重相同的相邻网格进行网格合并；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给每个gpu计算单元(或cpu线程)；

(3)结合单个网格对应的虚拟jaws坐标，对经过虚拟jaws的单个粒子进行剂量计算，得到单个线程内单个网格内单个粒子的剂量结果；

根据公式(1)，将单个网格内每个粒子产生的三维剂量按照体素进行叠加得到单个网格内所有粒子的三维剂量结果；

其中，d1(i,j,k)为单个网格产生的三维剂量，

m为进入该网格的粒子总数，m为该网格内的第m个粒子，

为该网格内单个粒子的剂量；

然后根据公式(2)，将单个线程下每个网格的三维剂量按照体素(用于统计剂量小体积单元或小质量单元)叠加，得到单个线程下的所有网格的总三维剂量；

其中，d2(i,j,k)为单个线程产生的剂量，

n为进入该线程的网格总数，n为该线程的第n个粒子，

为该线程内单个网格的剂量；

最后，根据公式(3)，将所有单个线程的剂量按照体素叠加得到复杂射野的剂量计算结果，

其中，d3(i,j,k)为该复杂射野的放射治疗剂量，

p为进行复杂射野剂量计算的线程总数，p为第p个线程，

为该复杂射野内单个线程的剂量。

实施例3楔形野

一种复杂射野的放射治疗剂量快速计算方法，如图6-9所示(本实施例中的限束装置包括楔形板和jaws)，包括如下步骤：

(1)定义楔形板的参数210：

定义楔形板的材料、密度、放置位置和几何尺寸等参数，楔形板的几何尺寸要保证全部覆盖射野范围；楔形野的中心轴与射束束流的中心轴对齐重合(如图7所示)；

(2)射野网格化：

对射野进行网格化操作，根据自定义将射野划分为若干网格221；如图8所示为大小为10cm*10cm的射野，将该射野均分为0.5cm*0.5cm的网格，其中图8中的“+”表示每个网格的中心坐标；其中每个网格代表一个射束照射范围，网格所在平面为射野平面，通常指的是等中心平面。

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给每个gpu计算单元(线程)222；

根据放射源的位置、网格的划分尺寸和网格中心坐标，采用求交算法计算出放射源粒子经过楔形板的路径长度223；

其中，网格大小可以根据实际测量用的探测器的尺寸来定义，每个网格代表一个射束照射范围，网格所在平面为射野平面，通常指的是等中心平面；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标224；

综上，通过楔形板几何模块和射野网格化操作，实现了非均匀射野的构造。

(3)计算出到达每个射束网格的粒子经过楔形板的衰减比230：

根据楔形板的材料参数，楔形板的质量衰减系数(物质的质量衰减系数可从下列书籍中查得：《防护辐射导论》，原子能出版社，附表1。)，结合放射源粒子经过楔形板的路径长度以及射野网格划分情况，根据公式(4)计算出到达每个射束网格的粒子经过楔形野(限束装置)的衰减比i/i0；

i＝i0e^(-μρd)式(4)

其中，i为衰减后的没有改变能量和方向的粒子数目；

i0为不含限束装置时到达每个网格的某能量和入射方向的粒子数目；

μ为限束装置的质量衰减系数，与限束装置的材料成分和粒子能量有关，cm²/g；

ρ为限束装置的密度,g/cm³；

d为限束装置的厚度,cm；

i/i0每个射束网格的粒子经过限束装置的衰减比；

(4)通过步骤(3)得到的衰减比，结合无楔形板时到达每个网格的射束粒子能谱(这些粒子是指从放射源，穿过限束装置到达网格的粒子)计算得到经过楔形板衰减后的射束粒子能谱241；

结合虚拟jaws坐标，对射束能谱衰减后的射束粒子进行剂量计算，该计算采用蒙特卡罗计算模型，计算的剂量为每个线程中单个网格的剂量242；

将单个网格剂量进行衰减修正后求和得到单个线程的剂量243，其中衰减修正为三维剂量的二维衰减修正；

最后将所有单个线程的剂量求和得到复杂射野放射治疗剂量结果244。图9为楔形野的水箱中的离轴剂量分布对比图，其中实线为通过本发明模拟得到的离轴剂量分布示意图；点画线为水箱中实测的离轴剂量分布图。通过对比，本发明提供的模拟方法与实测值高度吻合。

实施例4挡块构造l形野

一种复杂射野的放射治疗剂量快速计算方法，本实施例中的限束装置为包括挡块和jaws，通过挡块构建l形射野，该l形射野是在一个16cm×16cm开放野的一角用锥形合金挡块去掉一个12cm*12cm部分形成的，计算平面垂直l形射野的其中一段并且通过射束中心轴(参见中华人民共和国医药行业标准yy0775-201o测试例7，如附图10所示)。包括如下步骤：

(1)定义挡块的参数：

定义挡块的材料、密度、放置位置和几何尺寸等参数，挡块的几何尺寸要小于射野范围；

(2)射野网格化：

对射野进行网格化操作，根据自定义将射野划分为若干网格；将该l型射野均分为0.5*0.5平方厘米的网格，其中每个网格代表一个射束照射范围，网格所在平面为射野平面，通常指的是等中心平面。

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给现场可编程门阵列(field-programmablegatearray，fpga)；

根据放射源的位置、网格的划分尺寸和网格中心坐标，采用求交算法计算出放射源粒子经过挡块的路径长度；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

综上，通过挡块几何模块和射野网格化操作，实现了非均匀射野的构造。

(3)计算出到达每个射束网格的粒子经过挡块的衰减比：

根据挡块的材料参数，挡块的质量衰减系数(《防护辐射导论》，原子能出版社)，结合放射源粒子经过挡块的路径长度以及射野网格划分情况，根据公式(4)计算出到达每个射束网格的粒子经过挡块的衰减比i/i0；

(4)针对到达每个网格的射线经过挡块的路程不同，对每个网格入射粒子的能谱进行相应的衰减，通过步骤(3)得到的衰减比，结合无限束装置时到达每个网格的射束粒子能谱，计算得到经过楔形板衰减后的射束粒子能谱；

结合虚拟jaws坐标，对射束能谱衰减后的粒子进行剂量计算，该计算采用蒙特卡罗计算模型，计算的剂量为每个线程中单个网格的剂量；

将单个网格剂量进行衰减修正后求和得到单个线程的剂量，其中衰减修正为三维剂量的二维衰减修正；

最后将所有单个线程的剂量求和得到复杂射野放射治疗剂量结果。

实施例5非对称野

一种非对称射野的放射治疗剂量快速计算方法(如图11-12所示)，包括如下步骤：

(1)定义非对称射野的参数：

本实施例中示出的复杂射野为非对称射野，其中非对称野是指射野中心轴偏离且平行于射束中心轴的情景，如图10所示当准直器的角度为0时，图中射野区域的长为x1-x2＝10cm，宽为y1-y2＝10cm：偏心坐标点offsetx0＝-1/2(x1-x2)＝-2.5cm，offsety0＝-1/2(y1-y2)＝-3.0cm。

定义非对称射野的几何中心(偏心坐标点的位置为射野的几何中心)，射野大小等参数；

(2)射野网格化：

对射野进行均匀网格化操作，根据自定义将射野划分为若干大小相同的网格；其中每个网格代表一个射束照射范围，网格所在平面为射野平面；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给每个gpu计算单元(线程)；

根据放射源位置、网格坐标等参数，利用求交算法(或投影定理)求出每个小网格对应的虚拟钨门jaws的坐标参数；

(3)根据步骤(2)中每个小网格对应的虚拟钨门jaws的坐标参数，采用蒙特卡罗计算模型，对射束粒子进行剂量计算，得到的剂量为每个线程中单个网格的剂量；

将单个网格剂量求和得到单个线程的剂量；

最后将所有单个线程的剂量求和得到复杂射野放射治疗剂量结果。图12为本实施例非对称野剂量计算结果示意图。

实施例6斜入射野

一种带楔形板的斜入射复杂野的放射治疗剂量快速计算方法，本实施例中的射束入射方式为斜入射，限束装置包括楔形板和jaws(如图13所示，其中机头从a位置偏移到b位置)，包括如下步骤：

(1)定义斜入射野和楔形板的参数：

本实施例中示出的复杂射野为斜入射射野机头旋转角度(楔形板同时旋转)。

定义楔形板的材料、密度、放置位置和几何尺寸等参数，楔形板的几何尺寸要保证全部覆盖射野范围；

(2)射野网格化：

对射野进行网格化操作，根据自定义将射野划分为若干大小相同的网格；其中每个网格代表一个射束照射范围，网格所在平面为射野平面；

启动多线程模块，将每个网格的计算任务(穿过网格的粒子能量、方向、数量、位置等及辐射输运过程等)分配给每个gpu计算单元(线程)；

根据放射源位置、网格坐标等参数，利用求交算法(或投影定理)求出每个小网格对应的虚拟钨门jaws的坐标参数和每个网格内对应的放射源射线粒子经过楔形板的路径长度；

(3)计算出到达每个射束网格的粒子经过楔形板的衰减比：

根据楔形板的材料参数，楔形板的质量衰减系数(物质的质量衰减系数可从下列书籍中查得：《防护辐射导论》，原子能出版社，附表1。)，结合放射源粒子经过楔形板的路径长度以及射野网格划分情况，根据公式(4)计算出到达每个射束网格的粒子经过楔形野(限束装置)的衰减比i/i0；；

(4)通过步骤(3)得到的衰减比和步骤(2)中每个小网格对应的虚拟钨门jaws的坐标参数，对经过楔形板的射束粒子进行能谱衰减；

采用蒙特卡罗计算模型计算这些射束粒子在模体中的三维剂量分布，得到的剂量为每个线程中单个网格的剂量；

将单个网格剂量求和得到单个线程的剂量；

最后将所有单个线程的剂量求和得到斜入射野放射治疗剂量计算结果。图14为本实施例复杂野剂量计算结果示意图。

实施例7

本实施例提供了一种放射剂量的快速计算设备，包括：

一个或多个处理器；

存储器；以及

一个或多个程序，其实一个或多个程序存储在所述存储器中并被配置为由一个或多个处理器执行，其中一个或多个程序包括用于上述放射剂量的快速计算方法的指令，包括如下步骤：

(1)定义限束装置的参数：

(2)射野网格化操作：

将射野均匀网格化，网格所在平面为射野平面；

对权重相同的相邻网格进行网格合并；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

启动多线程模块，将每个网格的计算任务分配给每个计算单元或线程；

(3)结合单个网格对应的虚拟jaws坐标，对经过虚拟jaws的单个粒子进行剂量计算，得到单个线程内单个网格内单个粒子的剂量结果；

根据公式(1)，将单个网格内每个粒子产生的三维剂量按照体素进行叠加得到单个网格内所有粒子的三维剂量结果；

其中，d1(i,j,k)为单个网格产生的三维剂量，

m为进入该网格的粒子总数，m为该网格内的第m个粒子，

为该网格内单个粒子的剂量；

然后根据公式(2)，将单个线程下每个网格的三维剂量按照体素叠加，得到单个线程下的所有网格的总三维剂量；

其中，d2(i,j,k)为单个线程产生的三维剂量，

n为进入该线程的网格总数，n为该线程的第n个粒子，

为该线程内单个网格的剂量；

最后，根据公式(3)，将所有单个线程的三维剂量按照体素叠加得到复杂射野的剂量计算结果，

其中，d3(i,j,k)为该复杂射野的放射治疗剂量，

p为进行复杂射野剂量计算的线程总数，p为第p个线程，

为该复杂射野内单个线程的剂量；

或

(1)定义限束装置的参数；

(2)射野网格化：

对射野进行网格化操作，将射野划分为若干网格；

启动多线程模块，将每个网格的计算任务分配给每个计算单元或线程；

根据放射源的位置和网格的划分尺寸和网格中心坐标，采用求交算法或投影定理计算出放射源粒子经过限束装置的路径长度；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

(3)根据限束装置的材料参数，限束装置的质量衰减系数，结合放射源粒子经过限束装置的路径长度以及射野网格划分情况，计算出到达每个射束网格的粒子经过限束装置的衰减比；

(4)通过步骤(3)得到的衰减比作用于无限束装置时到达每个网格的粒子能谱；

结合虚拟jaws的坐标，对射束能谱衰减后的射束粒子进行剂量计算，计算的剂量为每个线程中单个网格的剂量；

将单个网格剂量进行衰减修正后求和得到单个线程的剂量结果；

最后将所有线程的剂量求和得到复杂射野放射治疗剂量结果；

或

(1)定义非对称射野或斜入射野的参数；

(2)射野网格化操作：

对射野进行均匀网格化操作，将非对称射野或斜入射野划分为若干网格；其中每个网格代表一个射束照射范围，网格所在平面为射野平面；

启动多线程模块，将每个网格的计算任务分配给每个计算单元或线程；

根据放射源的位置、非对称射野或斜入射野的参数、网格的划分尺寸和网格中心坐标，采用求交算法或投影定理计算出每个网格对应的虚拟钨门jaws的坐标参数；

(3)结合虚拟jaws的坐标参数采用蒙特卡罗计算模型对粒子进行剂量计算，得到每个线程中单个网格的剂量结果；

然后将单个网格的剂量进行叠加；将单个网格剂量求和得到单个线程的剂量；最后将所有单个线程的剂量求和得到复杂射野放射治疗剂量结果。

实施例8

一种存储一个或多个程序的计算机可读存储介质，其中，一个或多个程序包括指令，该指令适于由存储器加载并执行上述放射剂量的快速计算方法，包括如下步骤：

(1)定义限束装置的参数：

(2)射野网格化操作：

将射野均匀网格化，网格所在平面为射野平面；

对权重相同的相邻网格进行网格合并；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

启动多线程模块，将每个网格的计算任务分配给每个计算单元或线程；

(3)结合单个网格对应的虚拟jaws坐标，对经过虚拟jaws的单个粒子进行剂量计算，得到单个线程内单个网格内单个粒子的剂量结果；

根据公式(1)，将单个网格内每个粒子产生的三维剂量按照体素进行叠加得到单个网格内所有粒子的三维剂量结果；

其中，d1(i,j,k)为单个网格产生的三维剂量，

m为进入该网格的粒子总数，m为该网格内的第m个粒子，

为该网格内单个粒子的剂量；

然后根据公式(2)，将单个线程下每个网格的三维剂量按照体素叠加，得到单个线程下的所有网格的总三维剂量；

其中，d2(i,j,k)为单个线程产生的三维剂量，

n为进入该线程的网格总数，n为该线程的第n个粒子，

为该线程内单个网格的剂量；

最后，根据公式(3)，将所有单个线程的三维剂量按照体素叠加得到复杂射野的剂量计算结果，

其中，d3(i,j,k)为该复杂射野的放射治疗剂量，

p为进行复杂射野剂量计算的线程总数，p为第p个线程，

为该复杂射野内单个线程的剂量；

或

(1)定义限束装置的参数；

(2)射野网格化：

对射野进行网格化操作，将射野划分为若干网格；

启动多线程模块，将每个网格的计算任务分配给每个计算单元或线程；

根据放射源的位置和网格的划分尺寸和网格中心坐标，采用求交算法或投影定理计算出放射源粒子经过限束装置的路径长度；

根据网格划分情况，计算出每个网格坐标投影到真实jaws所在平面的坐标，得到虚拟jaws坐标；

(3)根据限束装置的材料参数，限束装置的质量衰减系数，结合放射源粒子经过限束装置的路径长度以及射野网格划分情况，计算出到达每个射束网格的粒子经过限束装置的衰减比；

(4)通过步骤(3)得到的衰减比作用于无限束装置时到达每个网格的粒子能谱；

结合虚拟jaws的坐标，对射束能谱衰减后的射束粒子进行剂量计算，计算的剂量为每个线程中单个网格的剂量；

将单个网格剂量进行衰减修正后求和得到单个线程的剂量结果；

最后将所有线程的剂量求和得到复杂射野放射治疗剂量结果；

或

(1)定义非对称射野或斜入射野的参数；

(2)射野网格化操作：

对射野进行均匀网格化操作，将非对称射野或斜入射野划分为若干网格；其中每个网格代表一个射束照射范围，网格所在平面为射野平面；

启动多线程模块，将每个网格的计算任务分配给每个计算单元或线程；

根据放射源的位置、非对称射野或斜入射野的参数、网格的划分尺寸和网格中心坐标，采用求交算法或投影定理计算出每个网格对应的虚拟钨门jaws的坐标参数；

(3)结合虚拟jaws的坐标参数采用蒙特卡罗计算模型对粒子进行剂量计算，得到每个线程中单个网格的剂量结果；

然后将单个网格的剂量进行叠加；将单个网格剂量求和得到单个线程的剂量；最后将所有单个线程的剂量求和得到复杂射野放射治疗剂量结果。

上述实施例1-6中，基于并行计算方法和射束分配流程，实现了一种复杂射野的快速高效计算方法。通过将原来的单能或多能射束到达体模或人体的射野平面划分多个网格，将网格等数量或不等数量的指派给每个并行的线程，实现了不损失精度的情况下，快速高效地进行复杂射野的模拟计算。另外，本发明提供的复杂射野的快速高效计算方法还能够移植到gpu或者其他并行计算平台。

应当理解，这里描述的各种技术可结合硬件或软件，或者它们的组合一起实现。从而，本发明的方法和设备，或者本发明的方法和设备的某些方面或部分可采取嵌入有形媒介，例如软盘、cd-rom、硬盘驱动器或者其它任意机器可读的存储介质中的程序代码(即指令)的形式，其中当程序被载入诸如计算机之类的机器，并被该机器执行时，该机器变成实践本发明的设备。

以示例而非限制的方式，计算机可读介质包括计算机存储介质和通信介质。计算机存储介质存储诸如计算机可读指令、数据结构、程序模块或其它数据等信息。通信介质一般以诸如载波或其它传输机制等已调制数据信号来体现计算机可读指令、数据结构、程序模块或其它数据，并且包括任何信息传递介质。以上的任一种的组合也包括在计算机可读介质的范围之内。