核磁引导下的三维蒙特卡罗剂量独立验证方法及系统与流程

1.本发明涉及放射治疗剂量计算和剂量验证技术领域，尤其涉及核磁引导下的三维蒙特卡罗剂量独立验证方法及系统。


背景技术：

2.自适应放疗（art）通过图像引导来获得患者摆位精度、解剖生理变化，进行治疗计划的修改或重新设计，可以最大限度地减少放疗误差，减少靶区外扩，更精准地保护正常组织。自适应放疗常用cbct引导，与cbct相比，mri有成像速度快、软组织对比度高、生物信息更丰富、无附加剂量等优点，更加适合自适应放疗。
3.将mri与放疗加速器相结合，形成一种新的放疗技术—基于核磁引导的放疗方法。相比于常规加速器，核磁引导加速器由于其特有的磁场洛伦兹力和复杂的核磁成像设备，需要一种不同于常规加速器的剂量验证方法。美国医学物理学家协会（aapm）的tg219报告充分表述了第三方剂量验证软件的重要性，国外商用tps已有用于临床核磁引导放疗的剂量计算方法，但还没有一种精确的第三方剂量验证方法。


技术实现要素：

4.针对现有技术的不足，本发明提供了核磁引导下的三维蒙特卡罗剂量独立验证方法及系统，解决了现有技术中缺少一种精确的第三方剂量验证方法的技术问题。
5.为解决上述技术问题，本发明提供了如下技术方案：一种核磁引导下的三维蒙特卡罗剂量独立验证方法，该方法包括以下步骤：s1、从dicom文件中获取关键信息，dicom文件包括ct文件、rtstruct文件、rtplan文件和rtdose文件，关键信息包括图像hu矩阵、放疗计划的实施参数、剂量网格和感兴趣器官信息；s2、根据已知参数搭建虚拟核磁加速器模型；s3、将关键信息输入虚拟核磁加速器模型中进行蒙特卡罗剂量计算得到蒙卡剂量；s4、通过3d伽马测试对比三维蒙卡剂量与从rtdose文件中获取的待验证剂量之间的差异大小，并生成剂量验证报告和dvh图；若三维蒙卡剂量与待验证剂量差异小，则验证通过并进入步骤s6；若三维蒙卡剂量与待验证剂量差异大，则验证不通过，并进入步骤s5；s5、引入实测剂量作为金标准，将待验证剂量和实测剂量作进一步对比；若实测剂量与待验证剂量差异小，则验证通过并进入步骤s6；若实测剂量与待验证剂量差异大，则验证不通过并结束；s6、治疗计划验证通过，并执行照射。
6.进一步地，在步骤s1中，从dicom文件中获取关键信息包括以下过程：读取ct文件，从ct文件中获取患者的三维ct图像，以及三维ct图像的起始点
origin、间隔space、维度size；读取rtstruct文件，从rtstruct文件中获取患者外轮廓范围、各感兴趣区域位置和辅助结构的位置、密度信息，以及文件内包含的support结构信息；读取rtplan文件，从rtplan文件中获取放疗计划的实施参数，包括源轴距sad、射束个数、射束角度、照射时间，以及每个控制点control point的准直结构位置、准直结构开口大小信息；读取rtdose文件，从rtdose文件中获取剂量网格的信息与放疗计划系统tps计算的剂量作为待验证剂量，剂量网格的信息包括剂量网格的起始点origin、间隔space、维度size。
7.进一步地，在步骤s2中，根据已知参数搭建虚拟核磁加速器模型包括以下步骤：s21、搭建虚拟核磁加速器模型的基础模型，基础模型包括靶、初级准直器、电离室、反射镜与患者无关的结构，包含mlc模型和钨门jaw的准直结构，以及包含恒温冷却器、线圈、治疗床的mr-linac特殊结构；s22、模拟基础模型中第一相空间平面处的粒子生成，第一相空间平面包括靶、初级准直器、电离室、反射镜；s23、模拟粒子在钨门jaw准直结构中的反应；s24、模拟粒子在多叶准直器mlc模型中的反应；s25、模拟粒子在恒温冷却器中的粒子散射。
8.进一步地，在步骤s23中，模拟粒子在钨门jaw准直结构中的反应，具体包括以下步骤：s231、通过钨门jaw的上表面、下表面、倾斜面的各个端点来确定钨门jaw的纵向位置，根据rtplan文件确定钨门jaw的横向开口大小；s232、通过点到平面的距离方程，判断粒子在钨门jaw内穿行的长度；若粒子与钨门jaw不相交，穿行长度为0，粒子进入mlc模型中；若粒子与钨门jaw相交，若粒子是电子，将粒子直接杀死；若粒子是光子，将粒子沿当前运动方向移动至钨门jaw的某一表面，并进入下一步；s233、找到粒子穿出钨门jaw所需要的路程，求得粒子在钨门jaw内发生康普顿散射的射程，比较粒子射程与路程的大小；若射程大于路程，则粒子直接穿出钨门jaw；若射程小于路程，则发生康普顿散射反应，改变粒子的速度、能量后重复步骤s233，直到粒子沉积或穿出钨门jaw。
9.进一步地，在步骤s24中，模拟粒子在多叶准直器mlc模型中的反应，具体包括以下步骤：s241、通过已知参数确定mlc模型中各叶片位置，已知参数包括mlc模型的叶片长度、叶片半径、叶片间距、叶片旋转角参数；s242、根据粒子沿前进方向与叶片上表面、下表面、侧面、圆形端面的交点，计算粒子在叶片内穿行的长度；若粒子与叶片不相交，即穿行长度为0，粒子存活，进入恒温冷却器中；若粒子与叶片相交，若粒子是电子，将粒子直接杀死；若粒子是光子，将粒子沿当
前运动方向移动至mlc模型的某一表面，并进入下一步；s243、找到粒子穿出mlc模型所需要的路程，求得粒子在mlc模型内发生康普顿散射的射程，比较粒子射程与路程的大小；若射程大于路程，则粒子直接穿出mlc模型；若射程小于路程，则发生康普顿散射反应，改变粒子的速度、能量后重复步骤s243，直到粒子沉积或穿出mlc模型。
10.进一步地，在步骤s25中，模拟粒子在恒温冷却器中的粒子散射，具体包括以下步骤：s251、将恒温冷却器的外圆外、圆环内、内圆内分别定义为区域1、区域2和区域3；s252、求粒子在恒温冷却器中的运动距离length；s253、根据粒子的运动距离length判断粒子的所在位置；当粒子在区域1且将与区域2相交时，粒子沿原方向运动length长度到达圆环上；当粒子在区域1但不会与区域2相交时，直接杀死粒子；当粒子在区域2时，根据材料和粒子能量求出粒子平均自由程，判断自由程和length的大小，若自由程大于length，则粒子直接跑出圆环，进入区域3；若自由程小于运动长度，则粒子在自由程末端发生反应；当粒子在区域3时，判断粒子是否与ct图像的hu矩阵相交，如果相交则退出循环，进行模体内的剂量计算；如果不相交，则粒子移动length长度，重新进入区域2内；s254、更新粒子位置，重复上述过程，直到粒子计算被终止或退出循环。
11.进一步地，在步骤s3中，计算三维蒙卡剂量包括以下步骤：s31、根据关键信息改变准直结构开口大小，模拟粒子在加速器机头的反应，得到打出加速器机头的粒子信息，粒子信息包括粒子位置坐标、粒子速度方向、粒子能量、粒子类型、粒子权重；s32、通过坐标转换将加速器机头打出的粒子转换到dicom图像坐标系，加速器机头打出的粒子包含粒子信息；s33、根据材料矩阵和密度矩阵判断粒子所在体素的材料和密度；s34、进行若干粒子的输运；s35、若干粒子的输运完成后，计算各个体素的三维蒙卡剂量得到蒙卡剂量。
12.进一步地，在步骤s34中，进行粒子运输具体包括以下步骤：s341、根据粒子所在体素的材料、密度和粒子的能量求得粒子在当前体素的各反应截面，并采用接受拒绝抽样方法判断粒子的反应类型；s342、根据粒子的反应类型改变粒子的能量、位置和方向，并将次级粒子存储为新粒子进入步骤s341循环；当粒子能量小于截止能量，粒子直接沉积在当前体素内；当粒子离开剂量网格范围时，则粒子死亡；s343、重复上述步骤直到所有的粒子沉积或死亡结束。
13.该技术方案还提供了一种用于实现上述剂量独立验证方法的系统，该验证系统包括：前处理模块、机头输运模块、坐标系转换模块、剂量计算模块和后处理模块；所述前处理模块用于提取dicom文件中的关键信息，并对ct图像进行重采样；
所述机头输运模块用于完成虚拟核磁加速器的建模，并模拟粒子在与患者无关部分的粒子生成、在准直结构里的粒子准直、在恒温冷却器的粒子散射和衰减；所述坐标系转换模块用于在机头输运模块中将准直结构参数从iec坐标系转换为机头坐标系，并在机头输运模块完成后将粒子从机头坐标系转换为dicom坐标系，用于剂量计算模块计算蒙卡剂量；所述剂量计算模块用于模拟粒子在患者体内或ct矩阵内的反应过程，并根据反应过程计算蒙卡剂量，反应过程包括带电粒子在磁场下的偏转过程；所述后处理模块用于分析剂量差异，并给出剂量验证结果报告。
14.进一步地，所述前处理模块包括dicom文件预处理模块、截面库导入模块和重采样模块，所述dicom文件预处理模块用于，读取ct文件，得到患者的三维hu矩阵，读取rtstruct文件，得到患者外轮廓范围、各感兴趣区域位置和辅助结构的位置、密度信息，尤其是文件内包含的support结构信息，读取rtplan文件，得到放疗计划的实施参数，读取rtdose文件，得到剂量网格的信息并将剂量网格中所包含的三维蒙卡剂量存储以供与待验证剂量对比；所述截面库导入模块用于预先将常用的关键信息的截面库提前制表并存在本地存储截面库中；所述重采样模块用于将ct文件中的三维hu矩阵向rtdose剂量网格进行重采样。
15.借由上述技术方案，本发明提供了核磁引导下的三维蒙特卡罗剂量独立验证方法及系统，至少具备以下有益效果：1、本发明提供了一种完整的独立剂量验证方法，同时可作为剂量的计算使用，该方法充分考虑核磁加速器的特殊部件并对其精确建模，在充分考虑粒子在磁场下的偏转，取带电粒子到体素几何边界的距离、因发生反应而停下所走的距离、因速度变化限制而停下所走的距离之中的最小值作为带电粒子真正的步长，并且对不同材料采用不同的速度变化限制，在满足临床速度的条件下，大大提高了对肺部等低密度区域的剂量计算精度，总的来说，该方法能在较短时间实现磁场下的基于蒙卡的高精度剂量计算，进而实现独立于计划系统的二次剂量验证，具有较高的临床价值。
16.2、本发明本发明充分考虑核磁加速器的特殊部件，如恒温冷却器、线圈、治疗床以及磁场下带电粒子偏转过程，并采用自适应的步长限制方法提高剂量计算精度，通过图形处理器gpu加速在两分钟内完成核磁引导计划的三维剂量分布的蒙特卡罗剂量验证，由此对治疗计划进行验证，提高治疗计划的精准性。
附图说明
17.此处所说明的附图用来提供对本技术的进一步理解，构成本技术的一部分，本技术的示意性实施例及其说明用于解释本技术，并不构成对本技术的不当限定。在附图中：图1为本发明剂量独立验证方法的流程图；图2为本发明虚拟核磁加速器的结构示意图；图3为本发明钨门jaw结构的示意图；
图4为本发明 mlc模型结构的示意图；图5为本发明恒温冷却器中粒子输运的示意图；图6为本发明粒子的坐标转换的示意图；图7为本发明剂量独立验证系统的结构框图；图8为本发明计量网格的示意图；图9为本发明示例性的dvh对比图。
18.图中：10、前处理模块；20、机头输运模块；30、坐标系转换模块；40、剂量计算模块；50、后处理模块；101、dicom文件预处理模块；102、截面库导入模块；103、重采样模块。
具体实施方式
[0019] 为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图 和具体实施方式对本发明作进一步详细的说明。借此对本技术如何应用技术手段来解决技术问题并达成技术功效的实现过程能充分理解并据以实施。
[0020]
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分步骤是可以通过程序来指令相关的硬件来完成，因此，本技术可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本技术可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、cd-rom、光学存储器等)上实施的计算机程序产品的形式。
[0021]
请参照图1-图9，示出了本实施例的一种具体实施方式，本实施例能在较短时间实现磁场下的基于蒙卡的高精度剂量计算，进而实现独立于计划系统的二次剂量验证，具有较高的临床价值。
[0022]
请参照图1，一种核磁引导下的三维蒙特卡罗剂量独立验证方法，该方法包括以下步骤：s1、从dicom文件中获取关键信息，dicom文件包括ct文件、rtstruct文件、rtplan文件和rtdose文件，关键信息包括图像hu矩阵、放疗计划的实施参数、剂量网格和感兴趣器官信息；在步骤s1中，从dicom文件中获取关键信息包括以下过程：读取ct文件，从ct文件中获取患者的三维ct图像，以及三维ct图像的起始点origin、间隔space、维度size；读取rtstruct文件，从rtstruct文件中获取患者外轮廓范围（包含人体行的所有部位）、各感兴趣区域位置（人体重要器官（如肺、脑、脊髓等器官）和辅助结构（如治疗床、成像线圈、冷却装置和面罩等）的位置、密度信息，尤其是文件内包含的support结构信息（如床、床垫等）；读取rtplan文件，从rtplan文件中获取放疗计划的实施参数，包括源轴距sad、射束个数、射束角度、照射时间，以及每个控制点control point的准直结构位置、准直结构开口大小信息；读取rtdose文件，从rtdose文件中获取剂量网格的信息与放疗计划系统tps计算的剂量作为待验证剂量，剂量网格的信息包括剂量网格的起始点origin、间隔space、维度size，待验证剂量为已知的数据，可以从rtdose文件中直接获取。
[0023]
若图8所示，剂量网格的信息包括剂量网格的起始点（origin
x
，originy，originz）、间隔（space
x
，spacey，spacez）、维度（size
x
，sizey，sizez）。剂量网格的起始点可以视为由x方向size
x
个体素，y方向spacey个体素，z方向sizez个体素，按照一定的顺序排列的三维网格，x轴最小值、y轴最小值、z轴最小值所在的点是起始点，各体素的大小由间隔决定，网格各方向的体素数目由维度决定，当粒子从网格内打出网格时，粒子死亡。
[0024]
使用python软件包simpleitk将三维ct图像采样成与剂量网格的信息相同维度、间隔、起始点的图像，即将ct文件的hu矩阵向rtdose剂量网格进行重采样，并根据需要填充空气。重采样完成后，从三维ct图像中提取具体数值（hu矩阵），按照ct机校准后的电子密度曲线将hu矩阵转换为密度矩阵。根据使用的材料种类设定材料阈值，将hu矩阵转换为材料矩阵。
[0025]
例如使用四种材料，材料的阈值为a1、a2、a3、a4、a5，当ct值在和之间时，材料为种类i（i= 1，2，3，4）。即hu矩阵的每一个数值代表患者ct序列第k张ct、第j列、第i行体素的hu值。相同索引下的材料矩阵数值和密度矩阵数值，代表此体素的材料序号和密度大小。
[0026]
s2、根据已知参数搭建虚拟核磁加速器模型；已知参数为根据现有资料和技术文献查找，得到核磁加速器模型各部分结构的详细参数（材料、密度、角度等），根据各部分结构的详细参数搭建虚拟的核磁加速器模型。然后根据实际测量数据调节虚拟核磁加速器模型的参数，完成充分的临床验证，保证虚拟核磁加速器建模的准确性。
[0027]
在步骤s2中，根据已知参数搭建虚拟核磁加速器模型包括以下步骤：s21、搭建虚拟核磁加速器模型的基础模型，如图2所示，基础模型包括靶、初级准直器、电离室、反射镜与患者无关的结构，包含mlc模型和钨门jaw的准直结构，以及包含恒温冷却器、线圈、治疗床的mr-linac特殊结构；与患者无关的结构由靶、初级准直器、电离室、反射镜组成，初始打靶电子呈高斯分布，平均能量在7mev左右。
[0028]
在本实例中，钨门jaw是斜面末端，沿y轴向移动；mlc模型是80对圆形末端叶片模拟elekta的agility叶片，沿x轴向移动。
[0029]
治疗床位于患者模体下方，由多个复杂精细的部件组成，并在冠状轴方向关于等中心点左右对称。对于每例患者都要对治疗床位置进行微调，无法建立统一的模型。因此本实例中采用对hu矩阵赋值的方法，从rtstruct文件提取治疗床各结构的位置、密度等信息，将密度转成hu值，添加到hu矩阵对应位置处，作为ct图像的一部分进行计算。
[0030]
线圈分为上下两部分，上线圈位于病人体表的正上方，由矩形和矩形两端的圆弧组成。下线圈架在治疗床bridge结构的上表面，为矩形结构，如图2所示。在本实例中对线圈的建模也是根据密度转换成hu值，添加到hu矩阵对应位置处，作为ct图像的一部分进行计算，相应的剂量可以忽略，主要是改变离子状态，主要应用在粒子输运方法中。
[0031]
恒温冷却器结构呈圆环状，圆心放置在等中心点处，患者体模、支撑部件（治疗床等）、线圈均放置在内圆内部。外圆距等中心93.7cm，内圆距等中心46.5cm。
[0032]
在步骤s2中，虚拟核磁加速器模型如图2所示，包括靶、初级准直器、电离室、反射镜等与病人无关的结构，和mlc模型和钨门jaw等准直结构，以及恒温冷却器、线圈、治疗床等mr-linac的特殊结构。采用蒙特卡罗程序构建与病人无关的机头部件，并模拟粒子在机头中的输运过程。另外，将线圈、治疗床添加到病人ct图像上，模拟粒子打入人体之前的衰减。最后，通过一组衰减系数来衡量冷却器的不均匀性。
[0033]
s22、模拟基础模型中第一相空间平面处的粒子生成，需要模拟粒子在靶、初级准直器、电离室、反射镜中的反应，第一相空间平面包括靶、初级准直器、电离室、反射镜，本实施例采用egsnrc的beamnrc模块来模拟第一相空间平面以上的结构，其中，采用source19模拟入射电子束，采用slab、cons3r、chamber、mirror结构分别模拟靶、初级准直器、电离室、反射镜，生成iaea格式的相空间文件存储粒子。
[0034]
s23、模拟粒子在钨门jaw准直结构中的反应；钨门jaw的结构图如图3所示，在步骤s23中，模拟粒子在钨门jaw准直结构中的反应，具体包括以下步骤：s231、通过钨门jaw的上表面、下表面、倾斜面的各个端点来确定钨门jaw的纵向位置，根据rtplan文件确定钨门jaw的横向开口大小；入射粒子为source19模拟的入射电子束，在经过与病人无关结构反应后的粒子后，生成iaea格式的相空间文件存储的粒子，从rtplan文件中能够直接获取源轴距（sad）、射束个数、射束角度、照射时间，以及每个控制点（control point）的准直结构位置、准直结构开口大小信息。
[0035]
s232、通过点到平面的距离方程，判断粒子在钨门jaw内穿行的长度；若粒子与钨门jaw不相交，穿行长度为0，粒子进入mlc模型中；若粒子与钨门jaw相交，若粒子是电子，将粒子直接杀死；若粒子是光子，将粒子沿当前运动方向移动至钨门jaw的某一表面，并进入下一步；s233、找到粒子穿出钨门jaw所需要的路程，求得粒子在钨门jaw内发生康普顿散射的射程，比较粒子射程与路程的大小；若射程大于路程，则粒子直接穿出钨门jaw；若射程小于路程，则发生康普顿散射反应，改变粒子的速度、能量后重复步骤s233，直到粒子沉积或穿出钨门jaw。
[0036]
s24、模拟粒子在多叶准直器mlc模型中的反应；mlc模型的结构如图4所示，在步骤s24中，模拟粒子在多叶准直器mlc模型中的反应，具体包括以下步骤：s241、通过已知参数确定mlc模型中各叶片位置，已知参数包括mlc模型的叶片长度、叶片半径、叶片间距、叶片旋转角参数；s242、根据粒子沿前进方向与叶片上表面、下表面、侧面、圆形端面的交点，计算粒子在叶片内穿行的长度；若粒子与叶片不相交，即穿行长度为0，粒子存活，进入恒温冷却器中；若粒子与叶片相交，若粒子是电子，将粒子直接杀死；若粒子是光子，将粒子沿当前运动方向移动至mlc模型的某一表面，并进入下一步；具体的，用点到平面的距离公式计算粒子沿当前方向，到各个表面的距离，取最小
的距离，就是粒子在叶片内的穿行长度，其中点到平面的距离公式如下：根据叶片各端点，将叶片某一表面写成ax+by+cz+d=0的形式，a、b、c为常规系数，假设粒子的坐标为，则粒子到该平面的距离d为：在本实施例中，点到平面的距离均采用上述公式进行计算，其他部分所涉及到该处的内容为了避免重复，便不再详细赘述，相互参见即可。
[0037]
s243、找到粒子穿出mlc模型所需要的路程，求得粒子在mlc模型内发生康普顿散射的射程，比较粒子射程与路程的大小；若射程大于路程，则粒子直接穿出mlc模型；若射程小于路程，则发生康普顿散射反应，改变粒子的速度、能量后重复步骤s243，直到粒子沉积或穿出mlc模型。
[0038]
粒子在叶片中的康普顿散射反应和在钨门jaw内的相同，在准直结构中主要考虑康普顿散射，预先将粒子反应截面库保存到本地，使用时直接读入，能节省大量的截面计算时间。
[0039]
通过步骤s23和步骤s24的模拟，能够得到粒子在第二相空间平面上的粒子分布。
[0040]
s25、模拟粒子在恒温冷却器中的粒子散射。
[0041]
在步骤s25中，模拟粒子在恒温冷却器中的粒子散射，具体包括以下步骤：s251、将恒温冷却器的外圆外、圆环内、内圆内分别定义为区域1、区域2和区域3；s252、求粒子在恒温冷却器中的运动距离length；具体的，请参照图5，假设粒子沿当前方向运动，与外圆交于点a和点d，与内圆交于点b和点c，其中粒子到各点的距离分别为disa、disb、disc和disd，且disa＞disd，disb＞disc。距离为负代表粒子沿反方向运动才能与圆相交，取最小的正的dis作为粒子下一步的运动距离length。若所有dis都是负值或任一dis不存在，即粒子不与圆环相交，则终止粒子的计算。
[0042]
s253、根据粒子的运动距离length判断粒子的所在位置；当粒子在区域1且将与区域2相交时，粒子沿原方向运动length长度到达圆环上；当粒子在区域1但不会与区域2相交时，直接杀死粒子；当粒子在区域2时，根据材料和粒子能量求出粒子平均自由程，判断自由程和length的大小，若自由程大于length，则粒子直接跑出圆环，进入区域3；若自由程小于运动长度，则粒子在自由程末端发生反应；根据蒙特卡罗方法采样粒子反应，若粒子发生康普顿散射，则改变光子的方向和能量，忽略次级电子，若抽样到其他反应，则终止粒子的计算。
[0043]
当粒子在区域3时，判断粒子是否与ct图像的hu矩阵相交，如果相交则退出循环，进行模体内的剂量计算；如果不相交，则粒子移动length长度，重新进入区域2内。
[0044]
s254、更新粒子位置，重复上述过程，直到粒子计算被终止或退出循环。另外，实际的恒温冷却器在各角度是不均匀的，各射束的剂量会根据机架角度乘以一组衰减系数来表征恒温冷却器在各角度的不均匀性。
[0045]
在步骤s2中，对于搭建的基础模型还需要按照加速器验收的标准对模型进行验收
测试，不断调整模型各参数，如部件离源距离、倾斜角、材料、密度等，直到满足验收要求。一旦验收完成与病人无关的结构便不再改变。临床上用于剂量验证时，通过rtplan文件中读入与病人相关的参数，改变准直开口、照射时间等，算出剂量分布。
[0046]
s3、将关键信息输入虚拟核磁加速器模型中进行蒙特卡罗剂量计算得到蒙卡剂量；在步骤s3中，计算三维蒙卡剂量包括以下步骤：s31、根据关键信息改变准直结构开口大小，模拟粒子在加速器机头的反应，得到打出加速器机头的粒子信息，粒子信息包括粒子位置坐标（x，y，z）、粒子速度方向（u，v，w）、粒子能量（e）、粒子类型（光子或电子）、粒子权重（wt），加速器机头为采用蒙特卡罗程序构建与病人无关的机头部件。
[0047]
s32、通过坐标转换将加速器机头打出的粒子转换到dicom图像坐标系，加速器机头打出的粒子包含粒子信息，请参照图6，加速器机头打出的粒子的初始状态为图中的机头坐标系，而将该坐标系下的粒子转换到dicom图像坐标系，只需要改变x、y、z的坐标即可实现，转换到dicom图像坐标系时，机头坐标系下的x轴不变，y轴的朝向转变为y1，z轴的朝向转变为z1。
[0048]
s33、根据材料矩阵和密度矩阵判断粒子所在体素的材料和密度，从ct文件中获取患者的三维hu矩阵，根据电子密度曲线将患者的三维hu矩阵转为材料矩阵和密度矩阵；s34、进行若干粒子的输运；在步骤s34中，进行粒子运输具体包括以下步骤：s341、根据粒子所在体素的材料、密度和粒子的能量求得粒子在当前体素的各反应截面，并采用接受拒绝抽样方法判断粒子的反应类型。
[0049]
在本实施例中，接受拒绝抽样方法为常规的现有技术手段，本实施例对该方法可以直接使用，并不需要做改变，同样的，以粒子在当前体素的各反应截面来判断粒子的反应类型，判断的方法为接受拒绝抽样方法，本实施例对该方法不再详细赘述。
[0050]
s342、根据粒子的反应类型改变粒子的能量、位置和方向，并将次级粒子存储为新粒子进入步骤s341循环，粒子由加速器机头打出后，入射到人体内部的是初级粒子，在人体内部，由初级粒子产生的粒子是次级粒子；当粒子能量小于截止能量，粒子直接沉积在当前体素内；具体的，当粒子能量较小时，无法穿过当前体素进入下一体素，直接沉积在当地，设定某一能量阈值，当粒子能量小于阈值时，认为粒子沉积在当地，不需要考虑该粒子的具体反应，该阈值叫截止能量，当改变后粒子的能量小于截止能量时，沉积在当前体素内。
[0051]
当粒子离开剂量网格范围时，则粒子死亡，具体的可以参照图8，给出了剂量网格的示意图。
[0052]
s343、重复上述步骤直到所有的粒子沉积或死亡结束。
[0053]
因为考虑光子的光电效应、康普顿散射和电子对效应。电子采用连续慢化近似方法，考虑电子的弹性碰撞、moller散射、韧致辐射等，因此对于若干粒子的输运方法采用光子-电子耦合输运方法，能够完美规避上述的因素所带来的影响。
[0054]
与常规核磁加速器相比，本实施例所提供的虚拟核磁加速器的剂量计算还需要特别考虑磁场下带电粒子的偏转作用。本实施例创新性地采用一种自适应步长的磁场下带电
粒子的输运方法。其中带电粒子加速度f的公式如（1）所示。
[0055]
根据相对论和一阶近似等，可以得到磁场下的速度变化量的单位矢量，如公式（2）所示。
[0056]
在一阶近似中，需要保证在步长内粒子的速度方向和速度大小变化很小，如公式（3）和公式（4）所示。
[0057]
电子到体素几何边界之间的距离记为。电子因发生反应而停下，所走的距离记为；电子因公式（3）和（4）的限制而停下，所走的距离记为。取、、之中的最小值作为带电粒子真正的步长。
[0058]
在本实施例中针对采用自适应调节的方法，在本实施案例中采用：低密度的材料中，如空气，为0.02，在较低密度的材料中（如肺），为0.1，在高密度的材料中（如骨），忽略此限制，即为无穷大。该方法可以在较小的精度损失下大幅提高计算速度，满足在线自适应放疗的临床要求。
[0059]
上式中，其中m是粒子质量，c表示光速，表示电场强度，表示粒子当前的速度向量，表示磁场强度，表示动量对时间t的微分，q表示粒子所带电荷量，表示粒子的初始速度向量，s表示粒子移动步长，表示初始光子能量，表示速度方向变化量，表示粒子能量变化量，代表光子能量变化量，代表粒子速度单位矢量的变化量，在各个材料中，和的取值不同的，根据材料类型而改变。
[0060]
s35、若干粒子的输运完成后，计算各个体素的三维蒙卡剂量得到蒙卡剂量。
[0061]
蒙卡剂量d的计算公式为：上式中，e是在单位体素内沉积的剂量，代表当前体素的密度，v是当前体素的体积，n是模拟的粒子数。of是归一化因子，用来将模拟剂量转换为实测剂量，具体做法是，在水箱中模拟1mu的粒子剂量，找到纵向剂量最大值处的值，则of可以通过下述公式求得，
最后还可以根据加速器实际情况进行微调。
[0062]
在此步骤中，预先将常用的各种材料的截面库制表并存在本地，计算前从本地读取截面库到内存中，计算时只需做简单的查询工作即可，可以节省大量时间。
[0063]
截面库包含光子的光电效应、康普顿散射、电子对效应、瑞利散射，以及电子的韧致辐射、弹性散射、穆勒散射等。
[0064]
制表是指从国际原子能机构（iaea）官网下载各元素在不同能量段下对各反应的微观截面。
[0065]
然后针对每一种常用材料，根据其组成元素和元素占比求得该材料在不同能量段下对各反应的微观截面。将此截面存在本地，计算时只需根据粒子所在材料和粒子能量，就可查得粒子发生各反应的微观截面，进行粒子运动模拟。将上百种元素的查询、求和转换为几种材料的查询，节省了大量时间。在本实施例中，选用四种常用材料，空气、肺、水和骨，其电子能量范围为200kev-21mev，光子的能量范围为10kev-21mev。
[0066] s4、通过3d伽马测试对比三维蒙卡剂量与从rtdose文件中获取的待验证剂量之间的差异大小，并生成剂量验证报告和dvh图，供临床参考。
[0067]
在步骤s4中，三维蒙卡剂量与待验证剂量之间的差异大小在对比后，会同步生成剂量验证报告和dvh图，如图9所示，为本实施例示例性的dvh对比示意图，图中给出了几个示例性的剂量验证结果对比，其中，（a）为靶区、（b）为肋骨、（c）为胃、（d）为小肠、（e）为左肺，而图中的1代表验证剂量，2代表三维蒙卡剂量，因此根据验证报告输出dvh图对比能够更直观的得到对第三方剂量验证的结果。
[0068]
剂量验证报告包含病人基本信息（姓名、病例号），计划信息（计划名，治疗次数，处方剂量，加速器名称，治疗技术等），验证结果（是否通过），待验证剂量与三维蒙卡剂量的对比结果（整体伽马通过率结果、重要器官的伽马通过率结果、重要器官的dvh对比、病人横断面剂量对比），时间节点（接收病例时间、检测时间、生成报告时间）。
[0069]
若三维蒙卡剂量与待验证剂量差异小，则验证通过并进入步骤s6；若三维蒙卡剂量与待验证剂量差异大，则验证不通过并进入步骤s5；具体的，当3d伽马测试结果的全局剂量大于95%时，则验证通过，若3d伽马测试结果的全局剂量小于95%时，则验证不通过，3d伽马测试结果的全局剂量为3%/2mm，是指对于三维蒙卡剂量（三维蒙卡剂量最大值记为）中每一个体素的剂量d，如果在待验证剂量中，在距离该体素2mm范围内，能找到与三维蒙卡剂量中每一个体素的剂量d的剂量差异在以内的体素，就视为这个体素验证通过。验证通过的体素与总的体素个数的比值，就是伽马测试结果。
[0070]
s5、引入实测剂量作为金标准，将待验证剂量和实测剂量作进一步对比；若实测剂量与待验证剂量差异小，则验证通过并进入步骤s6；此处的实测剂量相当于步骤s4的三维蒙卡剂量，差异的大小对比与步骤s4相同，均为3%/2mm，此处不再纤细赘述，若多次实测剂量与待验证剂量差异小，但三维蒙卡剂量与待验证剂量差异大，需要重新对虚拟核磁加速器模型进行校正，对于虚拟核磁加速器模型
进行校正则是对其参数进行调整，也就是在步骤s2构建虚拟核磁加速器模型的过程中校正各项参数。
[0071]
若实测剂量与待验证剂量差异大，则独立验证计划不通过并结束或在调整治疗计划后返回步骤s1重新进行剂量验证。
[0072]
实测剂量为病人qa计划在arccheck上的测量结果，金标准具有权威可靠性，待验证剂量与金标准的实测剂量差异大，则此次的独立待验证计划可信性存疑，不进行实际照射。
[0073]
s6、治疗计划验证通过，并执行照射。
[0074]
本实施例提供了一种完整的独立剂量验证方法，同时可作为剂量的计算使用，该方法充分考虑核磁加速器的特殊部件并对其精确建模，在充分考虑粒子在磁场下的偏转，取带电粒子到体素几何边界的距离、因发生反应而停下所走的距离、因速度变化限制而停下所走的距离之中的最小值作为带电粒子真正的步长，并且对不同材料采用不同的速度变化限制，在满足临床速度的条件下，大大提高了对肺部等低密度区域的剂量计算精度，总的来说，该方法能在较短时间实现磁场下的基于蒙卡的高精度剂量计算，进而实现独立于计划系统的二次剂量验证，具有较高的临床价值。
[0075]
与上述实施例提供的剂量独立验证方法相对应，本实施例还提供剂量独立验证方法的系统，由于本实施例提供的剂量独立验证系统与上述实施例提供的剂量独立验证方法相对应，因此前述剂量独立验证方法的实施方式也适用于本实施例提供的剂量独立验证系统，在本实施例中不再详细描述。
[0076]
请参阅图7，其所示为本实施例提供的剂量独立验证系统的结构框图，该剂量独立验证系统包括：前处理模块10、机头输运模块20、坐标系转换模块30、剂量计算模块40和后处理模块50；前处理模块10用于提取dicom文件中的关键信息，并对ct图像进行重采样；机头输运模块20用于完成虚拟核磁加速器的建模，并模拟粒子在与患者无关部分的粒子生成、在准直结构里的粒子准直、在恒温冷却器的粒子散射和衰减；坐标系转换模块30用于在机头输运模块20中将准直结构参数从iec坐标系转换为机头坐标系，并在机头输运模块20完成后将粒子从机头坐标系转换为dicom坐标系，用于剂量计算模块40计算蒙卡剂量；剂量计算模块40用于模拟粒子在患者体内或ct矩阵内的反应过程，并根据反应过程计算蒙卡剂量，反应过程包括带电粒子在磁场下的偏转过程；后处理模块50用于分析剂量差异，并给出剂量验证结果报告。
[0077]
前处理模块10包括dicom文件预处理模块101、截面库导入模块102和重采样模块103，dicom文件预处理模块101用于读取ct文件，得到患者的三维hu矩阵，读取rtstruct文件，得到患者外轮廓范围（包含人体行的所有部位）、各感兴趣区域位置（人体重要器官（如肺、脑、脊髓等器官）和辅助结构（如治疗床、成像线圈、冷却装置和面罩等）的位置、密度信息，尤其是文件内包含的support结构信息（如床、床垫等），读取rtplan文件，得到放疗计划的实施参数，读取rtdose文件，得到剂量网格的信息并将剂量网格中所包含的三维蒙卡剂量存储以供与待验证剂量对比；
截面库导入模块102用于预先将常用的关键信息的截面库提前制表并存在本地存储截面库中，计算前从本地读取截面库，此方法可以节省大量时间，加速程序运行；重采样模块103用于将ct文件中的三维hu矩阵向rtdose剂量网格进行重采样。
[0078]
需要说明的是，上述实施例提供的系统，在实现其功能时，仅以上述各功能模块的划分进行举例说明，实际应用中，可以根据需要而将上述功能分配由不同的功能模块完成，即将设备的内部结构划分成不同的功能模块，以完成以上描述的全部或者部分功能。另外，上述实施例提供的系统与方法实施例属于同一构思，其具体实现过程详见方法实施例，这里不再赘述。
[0079]
以上实施方式对本发明进行了详细介绍，本文中应用了具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制。