基于卷积核的锥束CT散射伪影校正方法与流程

本发明涉及一种锥束ct散射伪影校正方法，特别是一种基于卷积核的锥束ct散射伪影校正方法。

背景技术：

锥束ct从出现以来就凭借着其精度高、成像速度快、剂量小等优势被广泛的应用于临床医学和工业无损检测中。但是受其大视野和使用大面积探测器的影响，散射问题已经成为了影响cbct成像质量的主要原因，受到了越来越多的关注。为了解决cbct系统中的散射引起的图像质量下降的问题，学者们提出了许多方法。根据方法中使用手段的不同，主要可以分为三类方法：硬件校正法、软件校正法、软硬件混合校正法。其中硬件校正法是在ct系统中添加额外的硬件，通过减少x射线来达到减少散射的目的，例如准直器、抗散射光栅等。这些方法对于散射的削弱均有一定的效果，但是存在散射校正不完全、与cbct空间结构不相适应等问题。

软件法是在获得全部的投影数据后，利用数字图像处理算法估计散射光子的分布，并在投影数据中将其去除。常见的有卷积法、反卷积法和蒙特卡洛模拟法等。软件法不需要添加额外的硬件，实现简单，但是校正精度与卷积核的选取有很大关系。shaw等人1987年提出卷积法，通过将投影数据与核函数卷积、加权后得到散射估值，再从投影数据中将其减去实现散射校正。现在核函数主流的求解方法是通过蒙特卡洛模拟笔形束穿过不同厚度的体膜，每个笔形束都用点扩展函数来表征，散射分布就是每个点扩展函数的累积响应。谢世鹏、msun、esthermeyer、nanminabhatia等人从2010年到2016年均有相关的文章。

混合法是硬件法与软件法的结合，在向ct系统中添加一些特制的硬件如bsa、调制板等后经过一定的处理达到散射校正的目的。混合法效果一般较好，但是需要特制的硬件和扫描次数。

技术实现要素：

鉴于此，本发明提供一种基于卷积核的锥束ct散射伪影校正方法，提出了一种新的卷积核，其中最关键的部分是散射核函数的求解，针对同材质均匀物体，分析x射线与物质相互作用的每个过程，并用数学公式加以描述，从而求得整个平面的散射分布。

为了达到上述目的，本发明是通过以下技术方案实现的：

本发明提供一种基于卷积核的锥束ct散射伪影校正方法，包括：

对被测物进行ct扫描，得到投影数据；

估算x射线的初始的光子数n0；

利用初始光子数n0计算探测器上的光子数分布；

计算x射线散射分布；

从原始投影数据中扣除散射分布。

进一步地，估算x射线的初始的光子数n0，包括：

利用对被测物体进行ct扫描时的电压、电流和x射线能谱估算x射线的初始的光子数n0。

进一步地，利用初始光子数n0计算探测器上的光子数分布，包括：

计算被测物体的入射点到散射点的光子衰减；

计算经过散射点后光子散射角的概率分布；

计算从散射点到探测器的光子衰减。

进一步地，计算被测物体的入射点到散射点的光子衰减，包括：

根据比尔定律，被测物体的入射点到散射点的光子衰减为：

其中，n1即经过物体对x射线的衰减后光子数，μ为物体对x射线的衰减系数。

进一步地，计算经过散射点后光子散射角的概率分布，包括：

x射线在散射点发生散射，x射线方向发生偏转，偏转方向与原方向的夹角θ根据初始能量的大小呈现一定的概率分布，公式如下：

其中，e0为入射光子能量,e1为入射光子能量，α＝e0/m0c²，m0为静态电子质量，c为真空中光速，re为经典电子半径；

以探测器为中心点建立坐标系，散射点在探测器上对应的坐标为(x’,y’)，则有

根据上述公式可以得到散射光子偏转θ时的概率

进一步地，计算从散射点到探测器的光子衰减，包括：

计算散射前后光子数的变化：其中n2为散射前后的光子数；

根据比尔定律，计算被测物体的入射点到散射点的光子衰减为：

其中，n即经过被测物体对x射线的衰减后光子数，μ为物体对x射线的衰减系数。

进一步地，计算x射线散射分布，包括：

计算散射点对探测器上任一点造成的散射光子数

nf＝ρ×n3×p×δω，

其中，ρ为一调节常数，表示电子密度等对光子数的影响；n3表示x射线与物体作用过程中光子数；p表示散射光子偏转θ的概率；δω表示探测器对a点的张角。

本发明提供一种基于卷积核的锥束ct散射伪影校正方法，具有如下有益效果：

本发明基于锥束ct散射伪影矫正方法的卷积法提出了一种新的卷积核，其中最关键的部分是散射核函数的求解，针对同材质均匀物体，分析x射线与物质相互作用的每个过程，并用数学公式加以描述，从而求得整个平面的散射分布。

附图说明

图1为本发明实施例所提供的基于卷积核的锥束ct散射伪影校正方法的流程示意图；

图2为本发明实施例所提供的基于卷积核的锥束ct散射伪影校正方法另一种实施方式的流程示意图；

图3为x射线与被测物体相互作用示意图；

图4为以探测器中心点为极坐标下散射光子的角分布。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

以下结合具体情况说明本发明的示例性实施例：

请参考图1，图1为本发明实施例所提供的基于卷积核的锥束ct散射伪影校正方法的流程示意图；本实施例提供一种基于卷积核的锥束ct散射伪影校正方法，包括以下步骤：

步骤s101、对被测物进行ct扫描，得到投影数据；

步骤s102、估算x射线的初始的光子数n0；

步骤s103、利用初始光子数n0计算探测器上的光子数分布；

步骤s104、计算x射线散射分布；

步骤s105、从原始投影数据中扣除散射分布。

请参考图2，图2为本发明实施例所提供的基于卷积核的锥束ct散射伪影校正方法另一种实施方式的流程示意图；本实施例提供一种基于卷积核的锥束ct散射伪影校正方法，包括：

步骤s201、对被测物进行ct扫描，得到投影数据。

步骤s202、利用对被测物体进行ct扫描时的电压、电流和x射线能谱估算x射线的初始的光子数n0。

步骤s203、计算被测物体的入射点到散射点的光子衰减；

根据比尔定律，被测物体的入射点到散射点的光子衰减为：

其中，n1即经过物体对x射线的衰减后光子数，μ为物体对x射线的衰减系数。

步骤s204、计算经过散射点后光子散射角的概率分布；

其中，x射线在散射点发生散射，x射线方向发生偏转，偏转方向与原方向的夹角θ根据初始能量的大小呈现一定的概率分布，公式如下：

其中，e0为入射光子能量,e1为入射光子能量，α＝e0/m0c²，m0为静态电子质量，c为真空中光速，re为经典电子半径；

以探测器为中心点建立坐标系，散射点在探测器上对应的坐标为(x’,y’)，则有

根据上述公式可以得到散射光子偏转θ时的概率

步骤s205、计算从散射点到探测器的光子衰减；

在本实施例中，首先，计算散射前后光子数的变化：其中n2为散射前后的光子数。

然后，根据比尔定律，计算被测物体的入射点到散射点的光子衰减为：

其中，n即经过被测物体对x射线的衰减后光子数，μ为物体对x射线的衰减系数。

步骤s206、计算x射线散射分布；

在本实施例中，计算散射点对探测器上任一点造成的散射光子数

nf＝ρ×n3×p×δω，

其中，ρ为一调节常数，表示电子密度等对光子数的影响；n3表示x射线与物体作用过程中光子数；p表示散射光子偏转θ的概率；δω表示探测器对a点的张角。

步骤s207、从原始投影数据中扣除散射分布。

该方法的具体操作过程如下：

首先对被测物体进行ct扫描，利用得到的数据粗略的估计出物体的轮廓信息，然后利用扫描物体时的电压和电流以及相关的能谱估算出初始的光子数n0，并用求得物体轮廓信息和光子数依据所给公式计算出探测器上对应的光子分布。结合投影数据，建立线性方程组求解初始光子数，利用初始光子数解得散射分布，从初始的投影数据中口出散射分布。

x射线与物体相互作用的过程可以分为三个部分，具体过程如图3，图中a为散射点，b为所选取的截面。首先是从入射点到发生散射点a之前的部分l1段，此时物体对x射线起衰减作用，该部分遵从比尔定律，假设x射线初始的光子数为n0，那么有：

n1即经过物体对x射线的衰减后光子数，μ为物体对x射线的衰减系数。在同材质的均物体中，μ保持恒定不变。

其次x射线在a点发生散射，此时x射线方向发生偏转，能量也随之变化。x射线偏转方向与原方向的夹角根据初始能量的大小呈现一定的概率分布。klein-nishina公式描述的是实际散射物理过程中单个光子发生散射其康普顿散射微分截面，即这一概率分布。其具体公式如下：

其中e0为入射光子能量，α＝e0/m0c²，m0为静态电子质量，c为真空中光速，re为经典电子半径，te＝e/m0c²＝2.8×10^-13cm。图4给出了极坐标表示的单个电子康普顿散射微分截面与入射光子能量与散射角度的关系，从图中可以看出，不同能量下散射发生时，散射光子偏转不同方向的概率，并且随着入射光子能量的增加，散射光子偏转方向趋向于前向分布。

k-n公式描述的就是在入射能量为e0的条件下，散射角度发生变化的概率。

在只考虑一次散射的情况下，在a点发生散射后，角度发生偏转，以探测器中心点建立坐标系，a点在探测器上对应的坐标为(x’,y’)则偏转角：

根据k-n公式，可以计算出散射光子偏转θ时的概率p：

同时发生散射前后时光子数的变化可表示为：

发生散射后，物体对散射光子的衰减同样服从比尔定律：

n＝n2e^-μl；

散射点a对探测器上任一点造成的散射光子数可以表示为：

nf＝ρ×n3×p×δω；

其中ρ为一调节常数，表示电子密度等对光子数的影响。n3表示x射线与物体作用过程中光子数。p表示散射光子偏转θ的概率。δω表示探测器对a点的张角，可以认为是一个常数。

最后，还需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。

以上通过具体个例对本发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想；同时，对于本领域的一般技术人员，依据本发明的思想，在具体实施方式及应用范围上均会有改变之处，综上所述，本说明书内容不应理解为对本发明的限制。