一种用于放射治疗的基于通量分布的蒙特卡罗机头源模型建立方法及系统与流程

本发明涉及放射治疗的计量计算，特别是一种用于放射治疗的基于通量分布的蒙特卡罗机头源模型建立方法及系统。

背景技术：

1、在放射治疗剂量计算领域，蒙特卡罗方法因其能够精确模拟粒子在介质中的输运过程，被公认为剂量计算的“金标准”。然而，蒙卡方法在临床实际应用中面临一个关键挑战：如何高效、准确地建立能够反映治疗机头和各通量平面状态的粒子源模型。

2、常见做法之一是预先使用蒙卡程序模拟治疗机头的靶、均整器、初级准直器、次级准直器及多叶准直器等部件，在指定平面记录大量粒子的位置、方向和能量等信息，形成相空间文件。在后续的患者剂量计算中，蒙卡程序反复读取该相空间文件，将其中的粒子记录作为入射粒子源。

3、然而，相空间文件的数据量极大。单个相空间文件通常包含数百万至数亿条粒子记录，文件体积可达数gb甚至数十gb。在多射野及容积旋转调强放疗等包含大量通量平面的治疗计划中，需要反复读取多个相空间文件，由此产生显著的存储压力和输入/输出瓶颈。同时，相空间文件与具体的机型、标称能量和机头配置深度绑定——不同加速器、不同能量档或不同机头配置均需单独制作对应的相空间文件，导致数据的制作、传输和版本管理成本居高不下。此外，相空间文件中所含粒子数量有限，直接使用时还存在潜在的统计方差问题。

4、另一种做法是在每次剂量计算时，显式建立完整的治疗机头几何模型，从靶开始逐粒子模拟粒子在机头内部的完整输运过程，直至粒子到达患者或模体。该做法虽然过程最为完整，不需要预先制作相空间文件，但存在严重的效率问题：机头内大量粒子在到达患者之前即被机头结构吸收或偏离有效射野，真正到达患者的有效粒子占比极低。对于包含数十至上百个通量平面的imrt或vmat计划，若对每个通量平面都重复模拟完整的机头输运，计算资源消耗巨大，难以满足临床计划验证对时间的要求。

5、卷积/叠加类剂量算法已广泛应用于放射治疗计划系统。该类算法通过解析方式计算剂量，速度较快，能够按通量平面输出主通量，其中已包含射野权重、机架和准直器状态、钨门及多叶准直器开野形状等计划信息。然而，cccs等解析算法的剂量计算精度在某些非均匀介质区域及小野条件下，与蒙卡方法相比仍存在一定偏差。更重要的是，cccs的剂量计算过程需要针对不同入射角分别执行通量投射和卷积计算，对于vmat计划通常需要依次处理数十至上百个通量平面并将结果累加，这一过程仍然涉及大量的方向相关计算和中间结果合并操作。

技术实现思路

1、为克服现有技术的上述缺点，本发明的目的是提供一种用于放射治疗的基于通量分布的蒙特卡罗机头源模型建立方法及系统，以解决上述背景技术提出的问题。

2、本发明解决其技术问题所采用的技术方案是：一种用于放射治疗的基于通量分布的蒙特卡罗机头源模型建立方法，包括以下步骤：

3、获取多个通量平面，每个通量平面关联源轴距、空间变换、能谱及粒子类型，计算各通量平面所有像素通量之和并乘以像素面积，得到总通量积分并保存；

4、将全部通量平面展平为一维通量序列，计算归一化累积分布；

5、对每个蒙卡粒子历史生成随机数，对所述归一化累积分布进行反演抽样，一次性确定被选中的通量平面及其像素，使通量平面及像素的抽样概率正比于其通量值；

6、在被选中像素的横向和纵向边界内连续抽样粒子位置；

7、根据被选中通量平面的源轴距和空间变换，将等中心平面上的像素坐标映射为患者坐标系中的粒子方向单位向量，按被选中通量平面的累积能谱抽样粒子能量并设置粒子类型；

8、在等中心平面完成粒子抽样后，沿粒子方向的反方向进行几何回退，将粒子平移至计算区域入射边界内侧后启动蒙卡输运，蒙卡输运过程中将各粒子历史的能量沉积累计到同一剂量网格；

9、以所述总通量积分除以实际模拟的蒙卡粒子历史数作为每个蒙卡历史所代表的单位控制量粒子数，对各体素累计的蒙卡沉积能量按该比值进行缩放换算，输出绝对吸收剂量。

10、作为本发明的进一步改进：所述多个通量平面包括n个主源通量平面和n个散射源通量平面共2n个平面，其中前n个通量平面按填入主源通量，后n个通量平面按填入散射源通量；

11、其中，k=0,...,n-1为通量平面索引，为第k个通量平面在像素(i,j)处的通量值，为第k个通量平面的主源粒子数通量，为第k个通量平面对应的散射源通量。

12、作为本发明的进一步改进：所述主源通量平面接收卷积/叠加剂量算法已生成的各通量平面主通量，每个主源通量平面配置一个对应的散射源通量平面，主源源轴距记为，散射源源轴距记为；

13、所述散射源通量平面复制对应主源通量平面的机架角、准直器角、治疗床角和等中心，配置不同的散射源源轴距和散射能谱；

14、总通量积分的计算方式为：

15、，其中，δa为单个像素面积；

16、全部2n个通量平面共同进入所述归一化累积分布的构建和反演抽样。

17、作为本发明的进一步改进：根据治疗机头配置均整块采用ff模型：

18、对于ff机头，主源粒子数通量按进行离轴修正，主源平均能量按进行离轴软化，其中，控制离轴粒子数变化率，控制离轴能量软化率，设通量网格像素中心坐标为，为像素离轴半径；

19、散射源按构造均整块散射高斯分布，其中，为第k个通量平面的均整块散射分布，为散射强度系数，为第k个主源切片最大通量，、为均整块散射高斯分布宽度；

20、并按计算准直器开野形状的过滤因子，其中，为第k个通量平面的孔径过滤因子，为双线性插值，为限幅函数；

21、散射源通量按生成，其中，为第k个通量平面的散射源通量，为以、为宽度的高斯卷积操作，采用散射源源轴距<主源源轴距并经准直器和多叶准直器开野形状投影后过滤，再经二维高斯卷积模拟次级准直器和多叶准直器散射。

22、作为本发明的进一步改进：根据治疗机头未配置均整块采用fff模型：

23、对于fff机头，主源粒子数通量通过进行离轴修正，max(0,·)确保修正后通量非负；

24、主源平均能量按进行离轴软化且b取较小值或零；散射源通量按生成，用于等效描述有限源尺寸及次级准直器和多叶准直器散射。

25、作为本发明的进一步改进：所述总通量积分，其中δa为单个像素面积，为第k个通量平面在像素(i,j)处的通量值，该φ在归一化累积分布构建前单独计算并保存，φ不参与概率归一化；

26、所述归一化累积分布仅保留各通量网格的通量占比用于决定抽样概率，整体等比例缩放全部通量值不改变抽样概率而仅等比例改变φ和最终绝对剂量；

27、所述绝对吸收剂量按下式输出：

28、；

29、其中，为体素v的吸收剂量单位为gy/控制量，为实际模拟的蒙卡粒子历史数，为所有粒子历史在体素v的加权沉积能量总和，单位为mev，为第h个粒子历史的统计权重，为第h个粒子历史在体素v的沉积能量，为体素密度单位为g/cm³，为体素体积单位为cm³，1.60218×10-¹³为mev到j的换算系数，ρ_v·v_v×10-³为体素质量单位为kg。

30、作为本发明的进一步改进：所述实际模拟的蒙卡粒子历史数通过以下方式控制：

31、设基础批次数为、粒子数倍率为、有效野面积为、参考面积为，每个批次使用个线程块、每块个线程各发射一个初级粒子，则实际批次数为：

32、；

33、其中，floor(·)为向下取整函数；有效野面积根据高于预设通量阈值的像素面积确定，多通量平面任务按通量平面数取平均有效野面积；

34、作为本发明的进一步改进：所述方法还支持目标不确定度驱动的自适应粒子数控制：

35、在高剂量区域统计批次间剂量方差得到当前相对不确定度ε，当ε高于目标时，按估计下一轮所需粒子数，其中，n为当前已模拟粒子数，分批追加个粒子对应的粒子历史，确保追加量为非负值。

36、作为本发明的进一步改进：所述主源通量平面和散射源通量平面的模型参数通过标准方野的水箱测量数据联合标定获得：建立一组覆盖小野、中等野和大野的标准方野，对每个标准方野生成蒙卡源并在与水箱测量一致的几何条件下计算剂量，提取中心轴百分深度剂量pdd、代表深度处的横向及纵向剂量分布profile、以及各射野相对于参考野的总散射因子scp；

37、标定过程中以目标函数χ²=[(d_calc-d_meas)/σ_meas]²最小化为准则，其中χ²为卡方统计量，遍历所有标准方野，遍历每个射野的所有测量点，为计算剂量，为水箱测量剂量，为测量不确定度；

38、将各标准方野的计算结果与水箱测量数据逐项比较并迭代调节模型参数，直至全部标准方野在预设验收标准内同时吻合；

39、标定完成后将参数按机型、标称能量及ff/fff模式保存为固定机头模型，总通量积分作为绝对通量标度一并保存。

40、另一方面，本发明还提供一种用于放射治疗的基于通量分布的蒙特卡罗机头源模型建立系统，包括：

41、数据输入模块，配置为接收多个通量平面及各平面关联的源轴距、空间变换、能谱和粒子类型；

42、散射源生成模块，配置为根据治疗机头是否存在均整块分别采用ff模型或fff模型生成对应的散射源通量平面；

43、累积分布构建模块，配置为按计算总通量积分φ并单独保存，以及将全部通量平面展平后按构建归一化累积分布；

44、粒子抽样生成模块，配置为基于所述归一化累积分布反演抽样确定通量平面和像素，在像素内连续抽样位置，按构造方向并按能谱抽样能量；

45、几何回退模块，配置为将粒子沿反方向几何回退至计算区域入射边界内侧；

46、蒙卡输运与剂量累计模块，配置为对有效粒子执行蒙卡输运并将全部通量平面能量沉积直接累加到同一剂量数组；

47、粒子数控制模块，配置为按有效野面积控制或按目标不确定度控制以确定实际模拟粒子数；

48、绝对剂量换算模块，配置为以作为每历史粒子数，按输出gy/控制量的绝对吸收剂量；

49、参数管理模块，配置为按机型、能量和ff/fff模式存储已标定的机头模型参数的复用。

50、与现有技术相比，本发明的有益效果是：

51、1、本发明通过将全部通量平面展平为一维序列并构建统一的归一化累积分布，使每个蒙卡粒子仅需一次随机数反演即可同时确定被选中的通量平面及其像素，抽样概率自动正比于通量分布，无需预先为各通量平面人工分配粒子数，也无需像相空间源那样反复读取数gb的逐粒子记录文件。同时，将多个不同入射角、不同通量平面的通量平面置于同一个源对象中，在同一次蒙卡输运任务中完成全部通量平面的粒子输运，患者体素、材料、随机数及剂量数组只需建立一次，全部能量沉积直接累加到同一剂量网格，避免逐通量平面反复初始化计算环境和逐角度累加中间结果的开销，更重要的是，多个通量平面只要一次蒙卡任务。此外，在等中心平面完成粒子抽样后，沿粒子方向反方向进行几何回退，将粒子直接平移至计算区域入射边界内侧再启动蒙卡输运，跳过从源平面到患者边界之间不产生有效剂量的真空或空气传播段。

52、2、本发明将主源与散射源分别建模，为每个主源通量平面配置对应的散射源通量平面，散射源通量平面复制对应主源的空间几何但配置不同的源轴距和能谱，主源与散射源作为的通量平面共同参与统一抽样，从而在蒙卡输运中完整保留主源与散射粒子的差异。针对ff机头，本发明依据均整块中心厚、边缘薄的结构特征，采用离轴参数a描述主源粒子数中心少、边缘多的线性分布趋势，采用离轴参数b描述主源平均能量随离轴半径增加而软化的规律，同时以中心高斯分布描述均整块产生的散射并经准直器开野过滤和卷积展宽后形成散射源通量；针对fff机头，本发明采用不同的离轴修正方式，以高斯卷积将修正后的主源投影到散射平面来等效描述有限源尺寸及准直器散射。上述差异化建模具有明确的机头依据，避免统一经验修正造成的系统偏差。进一步地，本发明将总通量积分在概率归一化前单独保存，以作为每个蒙卡历史代表的单位控制量真实粒子数，对蒙卡沉积能量进行换算后直接输出绝对吸收剂量。