本发明涉及放射治疗系统技术领域,尤其涉及一种卷积叠加剂量计算系统。
背景技术:
放射治疗是目前治疗恶性肿瘤的主要手段之一。剂量计算是放射治疗计划的核心,剂量计算的速度与精度,对放射治疗计划制定的效率和质量具有重要影响。研究表明,照射剂量的准确性提高l%,治愈率可提高2%。一般把icru(internationalcommissionradiationunits&measurements)24号报告中推荐的±5%作为照射剂量误差的允许范围。在调强放射治疗的逆向计划中,优化过程需要进行多次的剂量计算(大约10到1000次),因此对计算速度的要求也非常苛刻。一个具备临床实用性的剂量计算模型应在1分钟以内完成单野、低精度的剂量计算;在1小时内完成多野、高精度或优化剂量计算。
计算剂量分布的模型可分为3大类:经验模型、半解析模型和解析模型。为了满足临床放射治疗计划的质量要求,一般使用半解析模型进行逆向计划的剂量计算,如基于核(笔形束核、点核)模型的卷积叠加剂量计算方法。虽然解析模型的剂量计算精度最高,但是所需的计算量非常大,不能用于逆向计划的剂量计算,一般只用来计算最终治疗计划的剂量分布。相对于解析模型,基于半解析模型的剂量计算方法的计算量相对较小,但是在逆向计划中多次计算剂量分布,其总的计算量也是非常可观的。一些硬件加速的方法被用来加速基于核模型的卷积叠加剂量计算方法,如fpga和gpu等。
能够满足临床放射治疗精确要求的半解析模型是点核剂量计算方法。点核剂量计算方法具有较大的计算复杂度。临床中,射线源被认为是一个点源,射线以点源为中心发散射出照射肿瘤。在发散射线倾斜入射模体表面的情况下,剂量计算过程中每个碰撞点处的点核都要进行旋转,与过该碰撞点的射线平行。在直角坐标系下旋转每个点核,又增加了点核剂量计算方法的计算复杂度,计算时间增加了2-3倍。卷积叠加剂量计算方法是一种半解析算法,其计算精度和计算速度都有待进一步提高。
鉴于此,现有技术有待改进,本案由此产生。
技术实现要素:
本发明目的在于提供一种卷积叠加剂量计算系统,以提高剂量计算精度和剂量计算速度。
为实现本发明目的,采用如下技术方案。
一种卷积叠加剂量计算系统,包括:
信息输入模块,用以输入剂量计算所需数据信息,所需数据信息包括模体的三维密度信息、器官勾画信息、治疗头信息和射野信息;
点核能量分布模拟模块,用以根据输入的治疗头信息,利用蒙特卡罗算法模拟点核能量分布;
点核模型参数提取模块,用以提取点核在各立体角方向的能量分布,进行参数拟合,得到点核模型参数;
点核查找表生成模块,用以构建新的剂量沉积点处采样模型,利用新模型的轴对称特性和旋转不变特性,将剂量沉积点到碰撞点的线段与碰撞点处入射射线的夹角存储,生成查找表;
term(totalenergyreleasedperunitmass)值计算模块,用以在直角坐标系下,根据治疗头信息计算模体表面的二维注量分布,利用射线追踪算法计算入射射束与各体素的相交长度,从而计算各体素的term值;
剂量计算模块,用以计算剂量沉积点的绝对位置,并根据该位置信息读取查找表;以剂量沉积点射线入射方向作为极角方向采样的起始方向,计算各立体角轴线方向上的碰撞点位置;根据碰撞点的位置信息从查找表中读取该碰撞点处剂量沉积点到碰撞点的线段与碰撞点处入射射线的夹角值;进行剂量计算。
信息输出模块,在直角坐标系下输出三维剂量分布,以及统计各器官的剂量-体积曲线。
进一步的,所述点核查找表生成模块构建新的剂量沉积点处采样模型,以剂量沉积点处射线入射方向作为极角方向采样的起始方向。
进一步的,所述新模型在同一剂量沉积点处,极角相同、方位角不同的两立体角轴线上,相对位置相同的两碰撞点处的剂量沉积点到碰撞点的线段与碰撞点处入射射线的夹角值是相同的。
进一步的,所述新模型以射线源为中心,射线源到剂量沉积点的距离为半径进行旋转,在新的剂量沉积点处,剂量沉积点到碰撞点的线段与碰撞点处入射射线的夹角值与原剂量沉积处相对位置相同的碰撞点处的夹角值相同。
进一步的,所述查找表对于相同射线源到剂量沉积点距离下,存储一份不同极角下碰撞点处的夹角值。
进一步的,所述卷积叠加剂量计算系统采用fpga和/或gpu硬件加速。
本发明具有的有益效果:本发明一方面构建剂量沉积点处的采样模型,进一步提高剂量计算精度;另一方面利用新模型的结构特点,将点核的旋转信息存储生成查找表,以提高剂量计算速度;由于本发明算法的特殊结构,非常适合于fpga和gpu硬件加速。
附图说明
图1是本发明一种卷积叠加剂量计算系统的模块框图。
图2是本发明一种卷积叠加剂量计算系统的射线平行入射示意图。
图3是本发明一种卷积叠加剂量计算系统的射线倾斜入射示意图。
图4是本发明一种卷积叠加剂量计算系统的新模型的示意图。
图5是本发明一种卷积叠加剂量计算系统的新模型的轴对称性示意图。
图6是本发明一种卷积叠加剂量计算系统的新模型的新模型的旋转不变性示意图。
具体实施方式
以下结合附图及具体实施例对本发明做详细的说明。
如图1所示,本发明提供一种卷积叠加剂量计算系统,包括信息输入模块10、点核能量分布模拟模块20、点核模型参数提取模块30、点核查找表生成模块40、term值计算模块50、剂量计算模块60和信息输出模块70。
所述信息输入模块10,用以输入剂量计算所需数据信息,所需数据信息包括病人的三维密度信息、器官勾画信息、治疗头信息和射野信息。其中病人的三维密度信息可以是ct图像、mr图像或其它方式获取的病人密度信息。器官勾画信息可以通过物理师在三维密度信息上进行勾画而获得的,也可以通过自动勾画软件进行自动勾画而获取。治疗头信息包括完整的治疗头形状结构和参数信息,以及照射方向,等中心的位置等。本实施例中,使用螺旋ct设备获取ct数据,该ct数据可以用来表示病人的密度信息。ct数据输入到器官勾画软件中,由物理师为每个病人勾画器官形状,从而获得病人的各器官信息。
所述点核能量分布模拟模块20,用以根据输入的治疗头信息,利用蒙特卡罗算法模拟点核模型的能量分布。在球坐标系下模拟点核的能量分布,极角方向的采样间隔为3.75°,采样数为48;方位角采样间隔为360°,采样数为1;径向使用非等间隔采样,采样数为24,最大范围为60cm。
所述点核模型参数提取模块30,用以提取点核在各立体角方向的能量分布,进行参数拟合,得到点核模型参数。本实施中,点核模型的能量扩散函数表示为:
其中aθ、aθ、bθ和bθ是与极角方向有关的参数值。提取每个立体角方向沿径向的能量分布,通过拟合方法计算点核模型中(公式1)的参数值。
所述点核查找表生成模块40,用以生成查找表。根据交互理论,剂量计算可以表示为term值和核模型的卷积形式:
其中t(s)为s点处来自于原射线的term值。严格意义上讲,只有在单能射线平行入射无限大均匀模体的情况下才能满足核的空间不变性,公式(2)才能够成立。这里我们只针对射线倾斜入射对公式(2)计算精度的影响展开讨论。尽管对方位角φ和极角θ采用等间隔的采样方式,不同采样位置处立体角
图3为射线倾斜入射的示意图。以z轴方向为极角的采样起始点,相对于r点,s点所在极角位置为θ,方位角为φ;以射线入射方向为极角的采样起始点,相对于s点,r点所在极角位置为β≠π-θ,方位角为φ+π。由于θ与β不在互补,s点所在立体角ωθ,φ与r点所在立体角ωβ,φ+π≠ωπ-θ,φ+π也不相等。因而,在射线倾斜入射的情况下,公式(2)的剂量计算卷积公式不在成立,不仅仅需要对s点处的点核模型h(s)沿该点入射射线方向进行旋转,还需要考虑相同位置处旋转前后立体角采样间隔一致对剂量计算精度所造成的影响。将r点处极角采样起始点与r点处保持一致,类似于射线平行入射,则可以减少由于s点处点核模型旋转引起的采样不一致而导致的剂量计算误差。然而,各碰撞点s处的射线入射方向是不一致的,根据其中一个碰撞点s处的射线入射方向确定r点处极角采样的起始方向,必然导致其它碰撞点处产生剂量计算误差。为了减小由于极角采样起始方向而引入的剂量总误差,本发明以剂量沉积点r处的入射射线方向作为极角的采样起点,如图4所示。由于射线源距离模体表面的距离较大,距离r点较近的碰撞点处的入射射线近似平行于r点处的入射射线,因而近似满足交互定理,且这些碰撞点对r点处剂量的贡献远大于其它碰撞点;距离r点较远的碰撞点处的入射射线不平行于r点处的入射射线,因而不满足交互定理,但是这些碰撞点对r点处剂量的贡献较小。
利用新模型进行剂量计算,需要计算剂量沉积点处极角采样的起始点,不可避免地增加了剂量计算的算法复杂度。然而,新模型具有轴对称性和旋转不变性,利用这两个特性,可以进一步减小卷积/叠加剂量计算的算法复杂度。图5所示为新模型的轴对称性。如图5所示,相对于剂量沉积点r,碰撞点s1和s2具有相同的极角θ,且碰撞点s1和s2到剂量沉积点r的距离相同。在这两碰撞点处,剂量沉积点到碰撞点线段与过碰撞点的入射射线的夹角α是相等的。因此,在r点处相同的极角下,碰撞点处的夹角α计算一次即可。该夹角值被用来确定公式(1)中各参数的值。图6所示为新模型的旋转不变性。如图6所示,点r1和点r2是距离射线源长度相同的两个剂量沉积点,点s1和点s2是相对点r1和点r2具有相同的极角、方位角和距离处的两个碰撞点。在这两碰撞点处,夹角α也是相等的。因此,距离射线源长度相同的剂量沉积点,碰撞点处的夹角α计算一次即可。
在卷积叠加剂量计算方法中,计算碰撞点处的夹角值所需计算量非常大。在本发明中,利用新模型的轴对称性和旋转不变性,将碰撞点处的夹角值存入查找表中,大大减少了剂量计算中计算夹角值所需计算量。在本实施中,极角采样数为48,径向采样数为60,沿模体深度方向采样间隔为0.01cm,采样数为10000,则总查找表的大小为109.86mb。点核查找表生成模块对于同一治疗头,只需生成一次查找表。
所述term值计算模块50,用以在直角坐标系下,根据治疗头信息计算模体表面的二维注量分布,利用射线追踪算法计算入射射束与各体素的相交长度,从而计算各体素的term值。
所述剂量计算模块60,用以计算剂量沉积点的绝对位置,并根据该位置信息读取查找表;以剂量沉积点射线入射方向作为极角方向采样的起始方向,计算各立体角轴线方向上的碰撞点位置;根据碰撞点的位置信息从查找表中读取该碰撞点处剂量沉积点到碰撞点的线段与碰撞点处入射射线的夹角值;进行剂量计算。碰撞点s释放的能量在r点处沉积的剂量可以写为:
其中ωmn为相对于r点,s点所在立体角,ηrmn和ρrmn为r点处的相对密度值和密度值,t(s)、σ(s)和ds为s点处的term值、密度值和径向微元长度,点r到点s之间的距离为rl,该距离被分为l段,每一段的长度为δri,每一段的相对密度为ηimn。已知s点的位置,就可以计算t(s)和σ(s)值;已知r点的位置,就可以计算ηrmn和ρrmn的值;已知s点和r点的位置,就可以计算r点到s点的线段与s点处入射射线的夹角,从而确定参数am、am、bm和bm的值;以及知道点r到点s之间各采样间隔的长度及相对密度,就可以计算∑ηimnδri的值。本实施中,以剂量沉积点为中心进行剂量计算,极角采样间隔为3.75°,采样数为48;方位角采样间隔为45°,采样数为8;径向使用非等间隔采样,采样数为60,最大半径为60cm。因而,径向的采样间隔值δri固定的,每个采样间隔中的相对密度值、密度值和term值近似等于该采样间隔中心点处的相对密度、密度值和term值。综上所述,在已知s点和r点的位置,以及这两点间各采样间隔的中心位置,就可以利用公式(3)计算碰撞点s处释放的能量在r点处的剂量沉积。
所述信息输出模块70,用以在直角坐标系下输出三维剂量分布,以及统计各器官的剂量-体积曲线。
为了验证本发明中算法的有效性,在水模体中对新算法进行了测试。所使用水模体大小为128×128×128,x、y和z三个方向的采样间隔都为0.25cm;使用6mev射线源照射水模体,射线源到模体表层的距离为100cm;模体表层二维注量图为均匀注量图,其大小为40×40,x和y两个方向的采样间隔都为0.5cm。新方法与射线平行入射方法和射线倾斜入射方法进行了比较。从实验结果中可以看出,新方法避免了剖面剂量射线边缘处的剂量“隆起”现象,与直接旋转点核的剂量计算方法的结果更为接近。射线平行入射方法所需计算时间为108.8s;射线倾斜入射方法所需计算时间为305.4s;新方法所需计算时间为200.5s。相比射线倾斜入射方法,新方法所需计算时间缩短了34%。
通过以上描述可知,本发明提供的一种卷积叠加剂量计算系统,一方面构建新的剂量沉积点处采样模型,提高了剂量计算精度;另一方面利用新模型的轴对称性和旋转不变性,将点核的旋转信息存储生成查找表,提高了剂量计算速度。
以上显示和描述了本发明的基本原理和主要特征以及本发明的优点。本行业的技术人员应该了解,本发明不受上述实施例的限制,上述实施例和说明书中描述的只是说明本发明的原理,在不脱离本发明精神和范围的前提下,本发明还会有各种变化和改进,这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。