CN113125476B - 一种面阵工业ct散射校正方法 - Google Patents

Abstract

本发明涉及一种面阵工业CT散射校正方法，步骤1、使用X射线机对被检测样品进行周向DR扫描，获得I幅第一周向DR图像，并对其进行射束硬化校正，得到I幅第一图像；步骤2、在X射线机出束窗口前放置一块衰减板；并采用与步骤1中相同的扫描工艺对被检测样品进行周向DR扫描，获得I幅第二周向DR图像，并对其进行射束硬化校正，得到I幅第二图像；步骤3、将相同扫描角度的第一图像和第二图像进行散射校正，得到散射校正后的DR数据；步骤4、采用散射校正后的DR数据进行锥束CT重建，得到散射校正后的CT图像。该方法操作简便、衰减板加工要求低、算法步骤少、易于实现，大幅降低了后续工业CT重建后图像的散射伪像。

Description

一种面阵工业CT散射校正方法

技术领域

本发明涉及CT图像校正领域，特别涉及一种面阵工业CT散射校正方法。

背景技术

工业CT检测技术是一种在X射线检测技术上发展起来的实用无损检测手段，具有成像直观，定量、定位、定性准确等优势，广泛应用于工业无损探伤、医疗卫生等领域。工业CT中普遍存在散射现象，当X射线穿过物体时，康普顿效应使光子偏离了原来的运动方向，散射光子不仅降低探测器的信噪比，还给图像带来了很多伪影，为了压低散射光子，工业CT机通常在探测器前加准直器。随着X射线断层影像朝着快速成像、高能束的方向发展，业内普遍采用了发散的锥束X射线源以及面积较大的探测板，照射面积的增加不可避免会产生更多的散射光子，与此同时准直器的加工难度日益增加，设计人员想方设法减小准直器深度或摒弃准直器。因此，面阵工业CT的散射影响远大于线阵工业CT，限制了面阵工业CT系统发展。研究表明，散射强度不仅与X射线照射范围有关，而且还与成像空间几何构成、X射线光谱、物体几何结构以及物体材质有关。面阵工业CT中散射引入的伪影会导致图像的对比度降低且影响CT数的精确度，而图像对比度降低以及引入的杯状伪影、条状伪影随之带来的后果就是在CT扫描结果中目标定性、定量、定位变得更难。因此，研究面阵工业CT散射校正方法非常必要的。

国内外关于X射线图像散射校正的研究大致可以分为硬件校正和软件校正两大类。所谓的硬件校正方法是在X射线成像系统的各个仪器组件上添加一些校正工具，以达到散射校正的目的。比如对X射线进行预硬化，减少X射线能谱的低能部分来减小散射的产生。减少或者在探测器阵列前面加放后准直器或是金属栅格，阻挡通过物体产生的散射X射线进入探测器内。以及经典的空气隙(air gap)方法，即增加被照物体与探测器的距离,也能够减小散射成分。软件校正方法是指对已经得到的X射线投影图像，采用数字图像处理的方法，在计算机中通过对图像本身的分析和对被照物体的性质的估计，得出一个散射分布图，来进行散射校正。其中包括卷积法、反卷积法、蒙特卡罗模拟法、模型估计法等等。在面阵工业CT中，由于采用了面阵探测器，一些在二维CT中行之有效的散射校正方法，如后准直器等方法无法使用。因此一直都是利用软件校正的方法进行散射校正，但是软件校正方法算法复杂、计算量大、效率较低。近年来，国内外相关研究人员采用散射校正板法取得了较好的效果，通过二维铅块点阵列(二维孔阵列铅板)置于被检工件与探测器之前，通过扫描2次，获得物体+散射校正板图像a、单独物体图像b，利用图像a二维插值运算得到散射强度分布图，再用图像b减去散射强度分布图得到散射校正后图像，该方法对校正板加工要求较高。

发明内容

本发明所要解决的技术问题是针对现有技术的现状，提供一种简单易实现且能大幅降低了后续工业CT重建后图像的散射伪像的面阵工业CT散射校正方法。

本发明解决上述技术问题所采用的技术方案为：一种面阵工业CT散射校正方法，其特征在于：包括以下步骤：

步骤1、使用X射线机对被检测样品进行周向DR扫描，获得I幅被检测样品的第一周向DR图像，并对I幅被检测样品的第一周向DR图像进行射束硬化校正，得到I幅第一图像；其中，将第i幅第一图像中像素点为(x,y)的像素值记为A_i(x,y)，i＝1、2…I；

步骤2、在X射线机出束窗口前放置一块厚度为S的衰减板，该衰减板采用的材料对应在该X射线能谱下的衰减系数为μ_S，μ_S和S之间满足以下关系：

并采用与步骤1中相同的扫描工艺对被检测样品进行周向DR扫描，获得I幅被检测样品的第二周向DR图像，并对I幅被检测样品的第二周向DR图像进行射束硬化校正，得到I幅第二图像，其中，将第i幅第二图像中像素点为(x,y)的像素值记为B_i(x,y)，i＝1、2…I；B_i(x,y)和A_i(x,y)的扫描角度相同；

步骤3、将相同扫描角度的第一图像A_i(x,y)和第二图像B_i(x,y)进行散射校正，得到散射校正后的DR数据C_i(x,y)；其中，

步骤4、采用散射校正后的DR数据C_i(x,y)进行CT重建，得到散射校正后的CT图像。

作为第一种方案，所述步骤1和步骤2中使用的射束硬化校正方法为阶梯试块或蒙特卡罗统计仿真方法。

作为第二种方案，所述步骤1和步骤2中使用的射束硬化校正方法包括以下步骤：

步骤a、加工N个不同穿透厚度的滤波片；其中，每个滤波片的穿透厚度均小于等于被检测样品的最大穿透厚度；N为正整数；

步骤b、将每个滤波片分次放置于步骤1中的X射线机出束窗口前，采用同一扫描工艺对每个滤波片进行DR扫描，获得每个滤波片的DR图像；每个滤波片的DR图像的尺寸均为c×d；

步骤c、根据步骤2中滤波片的DR图像的灰度值和滤波片对应的穿透厚度进行指数拟合，得到拟合函数q(t)，其中t为穿透厚度，q(t)为穿透厚度为t时拟合得到的灰度值；

步骤d、建立射束硬化校正函数T(e)，其中e为采集到的DR图像中各像素点的灰度值，T(e)为利用拟合函数q(t)＝e时反求获得的实际穿透厚度；

步骤e、利用步骤d中的射束硬化校正函数T(e)对步骤1中的I幅被检测样品的第一周向DR图像和步骤2中I幅被检测样品的第二周向DR图像进行射束硬化校正，得到I幅第一图像和I幅第二图像。

具体的，所述步骤c中的具体步骤为：

步骤c-1、将步骤b获取的每个滤波片的DR图像按照滤波片的穿透厚度从小到大的顺序进行排序，得到

其中，

为穿透厚度为z₁的滤波片所对应的DR图像，

为穿透厚度为z₂的滤波片所对应的DR图像，

为穿透厚度为z_n的滤波片所对应的DR图像，

为穿透厚度为z_N的滤波片所对应的DR图像，z₁＜z₂…＜z_n…<z_N；

步骤c-2、在上述

中选择同一个像素点(r,s)，1≤r≤c，1≤s≤d；并获取

中同一像素点(r,s)的灰度值；

步骤c-3、选取步骤c-1中排序后的前i₁张DR图像中同一像素点的灰度值，即：

中同一像素点(r,s)的灰度值，1≤i₁≤N，并结合穿透厚度对选取出的

中像素点(r,s)的灰度值进行指数函数拟合，得到第1次拟合后得到的指数函数g₁(t)；g₁(t)计算公式如下：

其中，t为穿透厚度，a₁为第1次拟合得到的幅值；b₁为第1次拟合得到的衰减系数，g₁(t)为穿透厚度为t时拟合得到的灰度值；

步骤c-4、计算出第1误差绝对值k₁；k₁的计算公式为：

其中，f_t(r,s)为穿透厚度为t的滤波片所对应的DR图像中像素点(r,s)的灰度值；

步骤c-5、根据下述公式计算得到第1差值

Δf_t ⁽¹⁾(r,s)＝f_t(r,s)-g₁(t)/R；

其中，R为预设的常数；t＝z₁、z₂、…z_N；

步骤c-6、将t＝z₁、z2、…

代入到Δf_t ^(p-1)(r,s)中得到i_p个灰度值，其中T≥i_p＞i₁，并按照步骤c-3～步骤c-5中相同的方式，分别得到第p次拟合后得到的指数函数g_p(t)、第p误差绝对值k_p和第p差值Δf_t ^(p)(r,s)；

Δf_t ^(p)(r,s)＝Δf_t ^(p-1)(r,s)-g_p(t)/R；

其中，a_p为第p次拟合得到的幅值；b_p为第p次拟合得到的衰减系数；

p的初始值为2；

步骤c-7、并使p的值依次加1，直至p＝P，且每一次拟合时使用的灰度值个数满足以下条件：i_P>i_P-1…>i_p>i_p-1…>i₁，且i_P＝N；P为正整数；

按照步骤c-6中的方法，分别得到拟合后得到的P个指数函数g₁(t)、g₂(t)、…g_P(t)；P个误差绝对值k₁、k₂、…k_P和P个差值Δf_t ⁽¹⁾(r,s)、Δf_t ⁽²⁾(r,s)、…Δf_t ^(P)(r,s)；

步骤c-8、计算出步骤c-7中P个误差绝对值k₁、k₂、、、k_P中的最小值，则提取该误差绝对值所对应的次数M；M∈[1,P]；

当M＝1时，则

则M>1，则

优选的，所述步骤c-5和步骤c-6中R的取值范围为1＜R≤2。

本方案中，所述步骤2中的扫描工艺包括管电压、管电流、放大比、焦点尺寸和积分时间。

与现有技术相比，本发明的优点在于：利用不带衰减板扫描和带衰减板扫描获得周向360度DR图像，对所有DR图像进行射束硬化校正，通过差运算获得散射区域，最大程度降低了噪声干扰，另外可通过散射公式确定出衰减板的最优厚度，从而方便后续散射校正，最终通过最简单的计算公式对不带衰减板时扫描得到的DR图像进行散射校正。因此该方法操作简便、衰减板加工要求低、算法步骤少、易于实现，大幅降低了后续工业CT重建后图像的散射伪像。

附图说明

图1为本发明实施例一中被检测样品的示意图；

图2(a)为图1中被检测样品在无衰减板时对应射束硬化校正前的DR图像；图2(b)为图1中被检测样品在无衰减板时对应射束硬化校正后的DR图像；图2(c)为图1中被检测样品在有衰减板时对应射束硬化校正前的DR图像；图2(d)为图1中被检测样品在有衰减板时对应射束硬化校正后的DR图像；

图3为对图2(b)进行散射校正后的DR图像；

图4(a)为对图2(a)进行CT重建后得到的CT图像；图4(b)为对图2(b)进行CT重建后得到的CT图像；图4(c)为对图3进行CT重建后得到的CT图像；

图5为使用本发明实施例二中的射束硬化校正方法对图2(a)进行校正后得到的DR图像。

具体实施方式

以下结合附图实施例对本发明作进一步详细描述。

实施例一：

一种面阵工业CT散射校正方法，包括以下步骤：

步骤1、使用X射线机对被检测样品进行周向DR扫描，获得I幅被检测样品的第一周向DR图像，并对I幅被检测样品的第一周向DR图像进行射束硬化校正，得到I幅第一图像；其中，将第i幅第一图像中像素点为(x,y)的像素记为A_i(x,y)，i＝1、2…I；

本实施例中，步骤1和步骤2中使用的射束硬化校正方法为阶梯试块或蒙特卡罗统计仿真方法，这些射束硬化校正方法为现有方法；当然，也可以采用其他现有的射束硬化校正方法；

步骤2、在X射线机出束窗口前放置一块厚度为S的衰减板，该衰减板采用的材料对应在该X射线能谱下的衰减系数为μ_S，μ_S和S之间满足以下关系：

并采用与步骤1中相同的扫描工艺对被检测样品进行周向DR扫描，获得I幅被检测样品的第二周向DR图像，并对I幅被检测样品的第二周向DR图像进行射束硬化校正，得到I幅第二图像，其中，将第i幅第二图像中像素点为(x,y)的像素记为B_i(x,y)，i＝1、2…I；B_i(x,y)和A_i(x,y)的扫描角度相同；其中，扫描工艺包括管电压、管电流、放大比、焦点尺寸和积分时间；

步骤3、将相同扫描角度的第一图像A_i(x,y)和第二图像B_i(x,y)进行散射校正，得到散射校正后的DR数据C_i(x,y)；其中，

步骤4、采用散射校正后的DR数据C_i(x,y)进行CT重建，得到散射校正后的CT图像。本实施例中，通过锥束CT重建即可得到CT图像；

上述步骤2中衰减板的衰减系数和厚度之间的关系式以及步骤3中的散射校正后的DR数据C_i(x,y)的计算公式是本发明中提出的改进点，其推导过程为：

X射线机与探测器之间没有被检测样品的情况下，探测器各个位置接收的射线能量均为E₀；设探测器上任意一点的位置(x1,y1)，由于探测器的位置呈正方形，因此x1、y1的范围均为[0,n]，取被检测样品上的其中一点(a,b)为例，探测器上点接收到的穿透被检测样品的能量包括以下两部分：

(1)透射能量：

根据朗伯-比尔定律：射线衰减基本规律，得到入射光子垂直透过被检测样品位置为(a,b)处后到达探测器的透射能量E_t的计算公式为：

E_t＝E₀·e^-μT(a,b)； 公式1

其中，μ为被检测样品材料在该X射线能谱下的衰减系数；T(a,b)为被检测样品上位置为(a,b)处的厚度；

(2)散射能量：

空间上所有位置(x1,y1)的入射光子透过对应位置的被检测样品厚度T(x1,y1)，产生散射，散射角度正好经过被检测样品上的位置(a,b)，则所有散射能量叠加可得到散射能量E_S：

其中，假设一次散射光子在探测器(x1,y1)处接收到的能量为E_S(x1,y1)，E_S与E₀之间遵循Klein-Nishina(K-N)公式：

其中，M为电子能量；θ_(x,y)为散射极角，即入射光子的方向与散射光子的方向之间的夹角，通过该角度保证散射角度正好经过被检测样品上的位置(a,b)；

对公式3进行化简，得到：

最终得到：

将公式5代入到公式2中得到：

由于X射线穿透被检测样品的能量，能量通过探测器接收形成的是图像。能量就是物理量，图像就是这个能量的量化，从而：

其中，A_i(a,b)为第i幅第一图像中像素点为(a,b)的能量，等同于第i幅第一图像中像素点为(a,b)的像素值；

根据衰减板的衰减原理得到：

其中，B_i(a,b)为第i幅第二图像中像素点为(a,b)的能量，等同于第i幅第二图像中像素点为(a,b)的像素值；

为了能去除B_i(a,b)中的透射衰减能量，即：将公式8和公式7代入到下述公式中，得到：

其中，根据公式7得到散射衰减能量计算公式为：

为了减少散射衰减能量的计算量且保证散射校正易实现，从而使：

在X射线机的当前工艺满足透射要求时，公式9中

远大于1，因此，上述公式11中需要满足：

对公式12求解，得到：

此时，则将

作为散射的计算公式；

因此得到散射校正后的计算公式为：

为了证明本发明中散射校正方法的有效性，本实施例中以直径为40mm的不锈钢圆片作为应用对象，如图1所示，根据被检测样品的密度、尺寸、结构，确定检测工艺(管电压400kV、管电流1800uA、放大比1.215、焦点尺寸0.4mm、积分时间1200ms)，按照步骤1中的处理，得到一幅被检测样品的第一周向DR图像，如图2(a)所示，对图2(a)进行射束硬化校正，得到图2(b)中所示的第一图像；采用紫铜作为衰减板的材料，在当前X射线能谱下，经过射束硬化校正，紫铜衰减板的衰减系数μ_S＝0.392746，对应的最佳衰减板厚度为

则按照步骤2中的方式得到一幅被检测样品的第二周向DR图像，如图2(c)所示，对图2(c)进行射束硬化校正，得到图2(d)中所示的第二图像；

并将相同扫描角度的放置衰减板前后的射束硬化校正后的DR数据A_i(x,y)和B_i(x,y)，进行处理得到散射校正后的DR数据C_i(x,y)，如图3所示。

另外对图2(a)、图2(b)和图3分别进行CT重建，得到4(a)、4(b)和4(c)，分别对应为对原始DR图像、原始DR图像进行射束硬化校正后得到的DR图像以及原始DR图像进行射束硬化校正和散射校正后的DR图像重建后获得CT图像，因此很清楚的看到了本方法大幅降低了后续工业CT重建后图像的散射伪像。

实施例二：

与实施例一不同的是，实施例二中使用的射束硬化校正方法包括以下步骤：

步骤a、加工N个不同穿透厚度的滤波片；其中，每个滤波片的穿透厚度均小于等于被检测样品的最大穿透厚度；N为正整数；

步骤b、将每个滤波片分次放置于X射线机出束窗口前，采用同一扫描工艺对每个滤波片进行DR扫描，获得每个滤波片的DR图像；每个滤波片的DR图像的尺寸均为c×d；其中，本实施例中扫描工艺至少包括管电压、管电流、放大比、焦点尺寸和积分时间这些参数；这些参数均为X射线机进行DR扫描时采用的工艺参数；

步骤c、根据步骤2中滤波片的DR图像的灰度值和滤波片对应的穿透厚度进行指数拟合，得到拟合函数q(t)，其中t为穿透厚度，q(t)为穿透厚度为t时拟合得到的灰度值；

其中，具体步骤为：

步骤c-1、将步骤b获取的每个滤波片的DR图像按照滤波片的穿透厚度从小到大的顺序进行排序，得到

其中，

为穿透厚度为z₁的滤波片所对应的DR图像，

为穿透厚度为z₂的滤波片所对应的DR图像，

为穿透厚度为z_n的滤波片所对应的DR图像，

为穿透厚度为z_N的滤波片所对应的DR图像，z₁＜z₂…＜z_n…<z_N；

步骤c-2、在上述

中选择同一个像素点(r,s)，1≤x≤c，1≤y≤d；并获取

中同一像素点(r,s)的灰度值；

步骤c-3、选取步骤c-1中排序后的前i₁张DR图像中同一像素点的灰度值，即：

中同一像素点(r,s)的灰度值，1≤i₁≤N，并结合穿透厚度对选取出的

中像素点(r,s)的灰度值进行指数函数拟合，得到第1次拟合后得到的指数函数g₁(t)；g₁(t)计算公式如下：

其中，t为穿透厚度，a₁为第1次拟合得到的幅值；b₁为第1次拟合得到的衰减系数，g₁(t)为穿透厚度为t时拟合得到的灰度值；

步骤c-4、计算出第1误差绝对值k₁；k₁的计算公式为：

其中，f_t(r,s)为穿透厚度为t的滤波片所对应的DR图像中像素点(r,s)的灰度值；

步骤c-5、根据下述公式计算得到第1差值Δf_t ⁽¹⁾(r,s)；

Δf_t ⁽¹⁾(r,s)＝f_t(r,s)-g₁(t)/R；

其中，R为预设的常数；t＝z₁、z₂、…z_N；

R的取值范围为1＜R≤2；本实施例中，R＝1.5；上述R的取值是根据实际实验确定出的，为了能对准确率进行微调；

步骤c-6、将t＝z₁、z₂、…

代入到Δf_t ^(p-1)(r,s)中得到i_p个灰度值，其中T≥i_p＞i₁，并按照步骤c-3～步骤c-5中相同的方式，分别得到第p次拟合后得到的指数函数g_p(t)、第p误差绝对值k_p和第p差值Δf_t ^(p)(r,s)；

其中，a_p为第p次拟合得到的幅值；b_p为第p次拟合得到的衰减系数；

p的初始值为2；

步骤c-7、并使p的值依次加1，直至p＝P，且每一次拟合时使用的灰度值个数满足以下条件：i_P>i_P-1…>i_p>i_p-1…>i₁，且i_P＝N；P为正整数；优选的，每一次迭代拟合时使用的灰度值个数按照相同的个数增加，i_p-i_p-1＝i_P-i_P-1；

按照步骤c-6中的方法，分别得到拟合后得到的P个指数函数g₁(t)、g₂(t)、…g_P(t)；P个误差绝对值k₁、k₂、…k_P和P个差值Δf_t ⁽¹⁾(r,s)、Δf_t ⁽²⁾(r,s)、…Δf_t ^(P)r,s)；

步骤c-8、计算出步骤c-7中P个误差绝对值k₁、k₂、、、k_P中的最小值，则提取该误差绝对值所对应的次数M；M∈[1,P]；并得到最优多项指数函数，并将最优多项指数函数作为拟合函数q(t)；

当M＝1时，则

当M>1时，则

步骤d、建立射束硬化校正函数T(e)，其中e为采集到的DR图像中各像素点的灰度值，T(e)为利用拟合函数q(t)＝e时反求获得的实际穿透厚度；即：实际穿透厚度t＝q^-1(e)；

获取步骤b中DR图像的位数s，利用拟合函数q(t)计算得到DR图像中像素点的灰度值为[0,2^s]时一一对应的穿透厚度t；

步骤e、利用步骤d中的射束硬化校正函数T(e)对步骤1中的I幅被检测样品的第一周向DR图像和步骤2中I幅被检测样品的第二周向DR图像进行射束硬化校正，得到I幅第一图像和I幅第二图像。其中，利用射束硬化校正函数T(e)对被检测样品的第一周向DR图像和第二周向DR图像进行射束硬化校正为现有的方法，在此不进行赘述；

本发明中以下述实验来对本方法进行说明，以直径为40mm的不锈钢圆片作为应用对象，如图1所示；采用紫铜作为滤波片的材料(不仅限于紫铜)，设置6种紫铜滤波片的穿透厚度(1mm、2mm、4mm、6mm、8mm和12mm)，将6种滤波片分次放置于X射线机出束窗口前，采用固定扫描工艺(管电压400kV、管电流1800uA、放大比1.215、焦点尺寸0.4mm、积分时间1200ms)进行DR扫描，采集不同穿透厚度的滤波片的DR图像，如图2所示；

对DR图像中每一点都进行射束硬化拟合校正，任意选取上述DR图像中的同一点(x,y)，在不同滤波片穿透厚度下的图像灰度值见表1。

表1滤波片穿透厚度与DR图像灰度对应表

序号	滤波片穿透厚度(mm)	DR图像灰度值
			1	1	48585.8
2	2	27492
			3	4	14052.6
4	6	8956.57
			5	8	6164.52
6	12	3225.98

第1次拟合时选取穿透厚度为1mm、2mm、4mm的滤波片的DR图像灰度值；第2次拟合时选取穿透厚度为1mm、2mm、4mm和6mm的滤波片的DR图像灰度值；第3次拟合时选取穿透厚度为1mm、2mm、4mm、6mm和8mm的滤波片的DR图像灰度值；第4次拟合时选取穿透厚度为1mm、2mm、4mm、6mm、8mm和12mm的滤波片的DR图像灰度值；设定R＝1.5，并按照上述步骤3中的具体方法得到每次拟合时对应的幅值、衰减和误差绝对值，具体结果见表2所示。

表2迭代过程及结果数据

根据表2中误差绝对值中的最小值，即第4次迭代时误差绝对值为最小值，建立对应的拟合函数为：

利用相同的扫描工艺(管电压400kV、管电流1800uA、放大比1.215、焦点尺寸0.4mm、积分时间1200ms)采集图1被检测样品的周向DR图像，得到图2(a)，利用射束硬化校正函数T(e)对图2(a)进行射束硬化校正，获得校正以后的DR图像，如图5，用于CT图像重建。

将本发明中方法与单一指数函数拟合方法、分段三次样条插值拟合方法进行比较；分别采用穿透厚度为1mm、2mm、3mm、4mm、6mm、8mm、10mm、12mm的滤波片采集获得DR图像，任意选取一点图像灰度值，见表3所示。选择1mm、2mm、4mm、6mm、8mm、12mm作为已知量，采用不同拟合方法计算3mm、10mm两点的灰度值，计算拟合结果与实测结果的误差。

表3滤波片穿透厚度与DR图像灰度对应表

序号	滤波片穿透厚度(mm)	DR图像灰度值
			1	1	48585.8
2	2	27492
			3	3	18764.5
4	4	14052.6
			5	6	8956.57
6	8	6164.52
			7	10	4395.67
8	12	3225.98

根据1mm、2mm、4mm、6mm、8mm、12mm对应的图像灰度值，利用单一指数函数拟合方法建立最优拟合结果p(t)为:

p(t)＝67537.2e^-0.3806t；

本发明中方法获得的最优多项指数函数q(t)为：

计算3mm、10mm两点的灰度值，结果见表4所示。

表4单一指数拟合、分段三次样条拟合和本发明中拟合方法的结果对比

可见，本专利方法拟合结果与实测结果误差最小。本方法中通过多次进行指数拟合的方法，使拟合结果更加接近于实测结果，而现有单一指数拟合方法通过一次数据的拟合而出现误差较大的问题，因此本方法的拟合结果更加精确。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (4)

1.一种面阵工业CT散射校正方法，其特征在于：包括以下步骤：

步骤1、使用X射线机对被检测样品进行周向DR扫描，获得I幅被检测样品的第一周向DR图像，并对I幅被检测样品的第一周向DR图像进行射束硬化校正，得到I幅第一图像；其中，将第i幅第一图像中像素点为(x,y)的像素值记为A_i(x,y)，i＝1、2…I；

步骤2、在X射线机出束窗口前放置一块厚度为S的衰减板，该衰减板采用的材料对应在X射线能谱下的衰减系数为μ_S，μ_S和S之间满足以下关系：

并采用与步骤1中相同的扫描工艺对被检测样品进行周向DR扫描，获得I幅被检测样品的第二周向DR图像，并对I幅被检测样品的第二周向DR图像进行射束硬化校正，得到I幅第二图像，其中，将第i幅第二图像中像素点为(x,y)的像素值记为B_i(x,y)，i＝1、2…I；B_i(x,y)和A_i(x,y)的扫描角度相同；

其中步骤1和步骤2中使用的射束硬化校正方法包括以下步骤：

步骤a、加工N个不同穿透厚度的滤波片；其中，每个滤波片的穿透厚度均小于等于被检测样品的最大穿透厚度；N为正整数；

步骤b、将每个滤波片分次放置于步骤1中的X射线机出束窗口前，采用同一扫描工艺对每个滤波片进行DR扫描，获得每个滤波片的DR图像；每个滤波片的DR图像的尺寸均为c×d；

步骤c、根据步骤2中滤波片的DR图像的灰度值和滤波片对应的穿透厚度进行指数拟合，得到拟合函数q(t)，其中t为穿透厚度，q(t)为穿透厚度为t时拟合得到的灰度值；

步骤c中的具体步骤为：

步骤c-1、将步骤b获取的每个滤波片的DR图像按照滤波片的穿透厚度从小到大的顺序进行排序，得到

其中，

为穿透厚度为z₁的滤波片所对应的DR图像，

为穿透厚度为z₂的滤波片所对应的DR图像，

为穿透厚度为z_n的滤波片所对应的DR图像，

为穿透厚度为z_N的滤波片所对应的DR图像，z₁＜z₂…＜z_n…＜z_N；

步骤c-2、在上述

中选择同一个像素点(r,s)，1≤r≤c，1≤s≤d；并获取

中同一像素点(r,s)的灰度值；

步骤c-3、选取步骤c-1中排序后的前i₁张DR图像中同一像素点的灰度值，即：

中同一像素点(r,s)的灰度值，1≤i₁≤N，并结合穿透厚度对选取出的

中像素点(r,s)的灰度值进行指数函数拟合，得到第1次拟合后得到的指数函数g₁(t)；g₁(t)计算公式如下：

其中，t为穿透厚度，a₁为第1次拟合得到的幅值；b₁为第1次拟合得到的衰减系数，g₁(t)为穿透厚度为t时拟合得到的灰度值；

步骤c-4、计算出第1误差绝对值k₁；k₁的计算公式为：

其中，f_t(r,s)为穿透厚度为t的滤波片所对应的DR图像中像素点(r,s)的灰度值；

步骤c-5、根据下述公式计算得到第1差值Δf_t ⁽¹⁾(r,s)；

Δf_t ⁽¹⁾(r,s)＝f_t(r,s)-g₁(t)/R；

其中，R为预设的常数；t＝z₁、z₂、…z_N；

步骤c-6、将

代入到Δf_t ^(p-1)(r,s)中得到i_p个灰度值，其中T≥i_p＞i₁，并按照步骤c-3～步骤c-5中相同的方式，分别得到第p次拟合后得到的指数函数g_p(t)、第p误差绝对值k_p和第p差值Δf_t ^(p)(r,s)；

Δf_t ^(p)(r,s)＝Δf_t ^(p-1)(r,s)-g_p(t)/R；

其中，a_p为第p次拟合得到的幅值；b_p为第p次拟合得到的衰减系数；

p的初始值为2；

步骤c-7、并使p的值依次加1，直至p＝P，且每一次拟合时使用的灰度值个数满足以下条件：i_P＞i_P-1…＞i_p＞i_p-1…＞i₁，且i_P＝N；P为正整数；

按照步骤c-6中的方法，分别得到拟合后得到的P个指数函数g₁(t)、g₂(t)、…g_P(t)；P个误差绝对值k₁、k₂、…k_P和P个差值Δf_t ⁽¹⁾(r,s)、Δf_t ⁽²⁾(r,s)、…Δf_t ^(P)(r,s)；

步骤c-8、计算出步骤c-7中P个误差绝对值k₁、k₂、、、k_P中的最小值，则提取该误差绝对值所对应的次数M；M∈[1,P]；

当M＝1时，则

则M＞1，则

步骤d、建立射束硬化校正函数T(e)，其中e为采集到的DR图像中各像素点的灰度值，T(e)为利用拟合函数q(t)＝e时反求获得的实际穿透厚度；

步骤e、利用步骤d中的射束硬化校正函数T(e)对步骤1中的I幅被检测样品的第一周向DR图像和步骤2中I幅被检测样品的第二周向DR图像进行射束硬化校正，得到I幅第一图像和I幅第二图像；

步骤3、将相同扫描角度的第一图像A_i(x,y)和第二图像B_i(x,y)进行散射校正，得到散射校正后的DR数据C_i(x,y)；其中，

步骤4、采用散射校正后的DR数据C_i(x,y)进行CT重建，得到散射校正后的CT图像。

2.根据权利要求1所述的面阵工业CT散射校正方法，其特征在于：所述步骤1和步骤2中使用的射束硬化校正方法为阶梯试块或蒙特卡罗统计仿真方法。

3.根据权利要求1所述的面阵工业CT散射校正方法，其特征在于：述步骤c-5和步骤c-6中R的取值范围为1＜R≤2。

4.根据权利要求1～3任一项所述的面阵工业CT散射校正方法，其特征在于：所述步骤2中的扫描工艺包括管电压、管电流、放大比、焦点尺寸和积分时间。

Images