CN111340127B - 基于材料聚类的能谱ct迭代材料分解方法和装置 - Google Patents

Abstract

本发明公开了一种基于材料聚类的能谱CT迭代材料分解方法和装置，其中方法包括：获取成像物体在不同能量射线下的能谱CT投影数据；利用所述能谱CT投影数据重建能谱CT图像；对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类；根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果。本发明能够有效地消除现有双能/能谱CT由于投影数据噪声或材料衰减系数接近带来的分解误差，大大提高能谱CT材料分解和物质识别的准确性，降低误报率和漏报率，具有重大的市场应用价值。

Description

基于材料聚类的能谱CT迭代材料分解方法和装置

技术领域

本发明涉及辐射成像技术领域，具体地，涉及一种基于材料聚类的能谱CT迭代材料分解方法和装置。

背景技术

X射线能谱CT通过采集两种及以上不同X射线能谱下的投影数据，然后通过专门的能谱CT材料分解算法，可以同时重建得到被扫描物体原子序数Z和电子密度ρ_e的分布图像，或几种基材料的分解系数，从而实现物质或材料识别。X射线能谱CT的材料分解算法可以分为三大类：投影域前处理方法、图像域后处理方法和迭代方法。无论是哪种方法，其核心都是认为任意物质的线衰减系数都可以分解为两个或更多个已知的、只与能量有关的基函数的线性组合。这个设定在很多情况下是合理的，也能够获得良好的材料分解结果。

以较为常用的投影域前处理方法为例，该方法认为任意物质的线衰减系数函数可以分解为两个已知的、只以能量为变量的基函数的线性组合，对每条射线路径，先用该路径上的能谱投影值求出线性组合系数在射线路径上的积分，然后利用CT图像重建算法得到任意像素/体素的组合系数，进而确定像素/体素的材料分解系数或者物质信息(原子序数和电子密度)。

但是，上述材料分解算法的效果比较依赖于X射线投影数据的噪声水平，因为材料分解算法过程存在一定的噪声放大，投影数据中过大的噪声会严重影响材料分解的准确性。另外，由于物质对X射线吸收的线衰减系数与X射线的能量有密切关系，并且随着X射线能量的增加，该系数迅速下降，当X射线能量在1～10MeV时，各类材料的线衰减系数差异变得很小，使得进行材料分解时对投影数据的噪声变得异常敏感。当投影数据噪声较为严重或各种材料线衰减系数较为接近时(例如MeV双能/能谱CT成像中)，会导致材料分解或物质识别的误差较大，比如低剂量CT成像，或者航空器、集装箱MeV高能CT成像。

发明内容

鉴于以上问题，本发明的目的是提供一种基于材料聚类的能谱CT迭代材料分解方法和装置，以解决现有技术中投影数据的噪声影响材料分解的准确性的问题。

为了实现上述目的，本发明的一个方面是提供一种基于材料聚类的能谱CT迭代材料分解方法，包括以下步骤：

获取成像物体在不同能量射线下的能谱CT投影数据；利用所述能谱CT投影数据重建能谱CT图像；对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类；根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果。

优选地，对重建得到的能谱CT图像进行聚类处理的步骤包括：计算不同能量的能谱CT图像对应的像素或体素点的能量值；根据所述像素或体素点的能量值得到多维数据集图像；确定多维数据集图像中的聚类中心；根据聚类中心对多维数据集图像中的所有点进行聚类。

优选地，确定多维数据集图像中的聚类中心的步骤包括：根据密度半径的预设值计算每个点的密度半径；根据计算得到的密度半径计算每个点的密度增量半径；根据每个点的密度半径和密度增量半径获取一个或多个密度中心，并将密度中心作为聚类中心。

优选地，根据聚类中心对多维数据集图像中的所有点进行聚类的步骤包括：以聚类中心为中心，将在密度半径的预设值范围内的所有点归类到所述聚类中心；将未归类到所述聚类中心的点归类到密度增量半径对应的点。

优选地，根据聚类结果进行能谱CT迭代材料分解的步骤包括：根据聚类结果获取聚类信息；将聚类信息和材料衰减系数的分解模型带入X射线的非线性投影公式中进行迭代求解，获取材料的分解系数。

优选地，所述材料的分解系数采用非线性最优化方法求解。

优选地，根据聚类结果进行能谱CT迭代材料分解的步骤还包括：根据材料的分解系数得到每个像素或体素点的电子密度和等效原子系数。

优选地，所述方法应用于双能CT或多能CT成像的材料分解。

为了实现上述目的，本发明的另一个方面是提供一种基于材料聚类的能谱CT迭代材料分解装置，包括：数据采集模块，用于获取成像物体在不同能量射线下的能谱CT投影数据；图像重建模块，用于利用所述能谱CT投影数据重建能谱CT图像；聚类模块，用于对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类；分解模块，用于根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果。

优选地，所述聚类模块包括：能量值计算单元，用于计算不同能量的能谱CT图像对应的像素或体素点的能量值；图像获取单元，用于根据所述像素或体素点的能量值得到多维数据集图像；中心确定单元，用于确定多维数据集图像中的聚类中心；聚类单元，用于根据聚类中心对多维数据集图像中的所有点进行聚类。

与现有技术相比，本发明具有以下优点和有益效果：

本发明所述基于材料聚类的能谱CT迭代材料分解方法，通过材料聚类的方法将成像物体的像素或体素点聚类呈多个互不相交的类，再进行材料的迭代分解，得到最终的材料分解结果或物质识别结果。本发明能够有效地消除现有双能/能谱CT由于投影数据噪声或材料衰减系数接近带来的分解误差，大大提高能谱CT材料分解和物质识别的准确性，降低误报率和漏报率，具有重大的市场应用价值。

本发明对于临床医疗、安检、工业无损检测等领域的应用都有十分重要的应用价值，特别是针对材料选择跨度很大的高能双能/能谱CT具有更加重要的价值。

附图说明

图1是本发明所述基于材料聚类的能谱CT迭代材料分解方法的流程示意图；

图2是本发明中多维数据集图像示意图；

图3是本发明所述迭代材料分解方法在仿真模拟中的结果与现有材料分解方法的结果对比图；

图4是本发明所述基于材料聚类的能谱CT迭代材料分解装置的构成框图。

具体实施方式

下面将参考附图来描述本发明所述的实施例。本领域的普通技术人员可以认识到，在不偏离本发明的精神和范围的情况下，可以用各种不同的方式或其组合对所描述的实施例进行修正。因此，附图和描述在本质上是说明性的，而不是用于限制权利要求的保护范围。此外，在本说明书中，附图未按比例画出，并且相同的附图标记表示相同的部分。

图1是本发明所述基于材料聚类的能谱CT迭代材料分解方法的流程示意图，如图1所示，本发明所述基于材料聚类的能谱CT迭代材料分解方法，包括以下步骤：

步骤S1，获取成像物体在不同能量射线下的能谱CT投影数据；其中，采用能谱CT系统利用不同能量的X射线对成像物体进行扫描，获取能谱CT投影数据，所述能谱CT投影数据包括至少两个能谱的CT投影数据，既包括双能CT投影数据，也包括多能CT投影数据，例如，基于多能区光子计数探测器的能谱CT投影数据。

步骤S2，利用所述能谱CT投影数据重建能谱CT图像；

步骤S3，对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类，每个类代表一种未知的材料，具有一个衰减系数；

步骤S4，根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果。

本发明通过将重建的能谱CT图像聚类形成多个互不相交的类，在进行能谱CT迭代材料分解，得到的分解结果或物质识别结果更加准确，降低误报率和漏报率，具有重大的市场应用价值。

需要说明的是，本发明中，利用能谱CT投影数据重建能谱CT图像可以使用现有的图像重建方法获取，在此不再详细赘述。

在一个实施例中，对重建得到的能谱CT图像进行聚类处理的步骤包括：计算不同能量的能谱CT图像对应的像素或体素点的能量值，若为二维图像，则计算像素点的能量值，若为三维图像，则计算体素点的能量值；根据所述像素或体素点的能量值得到多维数据集图像；确定多维数据集图像中的聚类中心；根据聚类中心对多维数据集图像中的所有点进行聚类。其中，多维数据集图像的坐标为被扫描物体在不同能量下的线衰减系数或者线衰减系数的差值，例如，对于双能CT对应的数据集图像为二维数据集图像，坐标分别为μ_L,μ_L-μ_H，其中，μ_L,μ_H分别表示被扫描物体低能和高能下的线衰减系数；对于三能CT对应的数据集图像为三维数据集图像，坐标分别为μ₁,μ₁-μ₂，μ₂-μ₃，其中μ₁,μ₂,μ₃分别三种不同能量下被扫描物体的线衰减系数，若是对多能CT的能谱CT图像进行聚类时，则需要得到更多维度的数据集图像用于聚类。

进一步地，确定多维数据集图像中的聚类中心的步骤包括：根据密度半径的预设值计算每个点的密度半径；根据计算得到的密度半径计算每个点的密度增量半径；根据每个点的密度半径和密度增量半径获取一个或多个密度中心，并将密度中心作为聚类中心，其中，密度中心的数量根据实际应用中可能的材料的种类数量进行事先预估得到，比如，医学成像一般选择3～5种，安检应用则更多一些。例如，如图2所示，以双能CT得到的多维数据集图像为例，该图像中，横坐标为被扫描物体低能下的线衰减系数，纵坐标为被扫描物体在低能和高能下的线衰减系数之差，图中的各聚焦的散点中心的黑色圆点则表示聚类中心，多维数据集图像中包括多个聚类中心。

进一步地，在一个可选的实施例中，根据聚类中心对多维数据集图像中的所有点进行聚类的步骤包括：以聚类中心为中心，将在密度半径的预设值范围内的所有点归类到所述聚类中心；将未归类到所述聚类中心的点归类到密度增量半径对应的点，从而实现所有点的聚类。

在一个实施例中，根据聚类结果进行能谱CT迭代材料分解的步骤包括：根据聚类结果获取聚类信息，聚类信息指的是聚类得到的互不相交的类；将聚类信息和材料衰减系数的分解模型带入X射线的非线性投影公式中进行迭代求解，获取材料的分解系数。优选地，所述材料的分解系数采用非线性最优化方法求解，例如Levenberg–Marquardt算法。

进一步地，根据聚类结果进行能谱CT迭代材料分解的步骤还包括：根据材料的分解系数得到每个像素或体素点的电子密度和等效原子系数。

需要说明的是，上述基于材料聚类的能谱CT迭代材料分解方法不仅可以应用于X射线双能CT成像的材料分解，当可以采集到更多个能谱的投影数据，例如三能CT数据或者基于多能区光子计数探测器的能谱CT数据时，还可以应用于多能CT成像的材料分解。

下面以双能CT成像为例进行详细说明。

采用双能CT系统利用低能和高能的X射线对成像物体进行扫描，得到双能CT投影数据，利用双能CT投影数据重建能谱CT图像，并对重建得到的CT图像进行聚类处理。具体地，分别计算低能CT图像和高能CT图像对应各像素/体素点的(μ_L,μ_L-μ_H)，其中，μ_L,μ_H分别表示被扫描物体低能和高能下的线衰减系数；然后画出如图2所示的(μ_L,μ_L-μ_H)数据集图像。

多维数据集图像中包括多个散点，根据密度半径的预设值计算每个散点的密度半径，具体地：预设密度半径的数值为ρ₀，例如0.001，以及N个P维向量点

这些向量点组成了集合Φ；图2中的各个散点，根据下式(1)计算每个点的密度半径：

其中，

表示散点

的密度半径，θ表示集合Φ内的任意散点，

是散点

和散点θ在P维实空间R^P的欧几里得距离。

通过下式(2)计算每个点的密度增量半径：

其中，

表示散点

的密度增量半径，

表示散点

的密度半径，ρ(θ)表示散点θ的密度半径。

根据密度半径和密度增量半径获取密度中心，并将密度中心作为聚类中心。具体地，在(ρ,τ)散点图上，将

和

都较大的点

作为密度中心，可以使用从大到小排序的算法，分别对

和

进行排序，获取这两个数值均排在前列的一个或多个点

作为密度中心，密度中心的个数可以根据实际应用中可能的材料的种类数量进行事先预估，以每个密度中心作为一个聚类中心，如图1中的各聚焦的散点中心的黑色圆点即为密度中心。根据聚类中心对所有点进行归类时，先将密度半径预设值ρ₀范围内的所有点归类到该聚类中心；对于其他未归类的点

将其归类到其密度增量半径对应的点θ，即θ满足

通过这样的归类方式实现所有点的聚类。待扫描成像的物体可以通过双能CT各能量下的CT重建图像聚类成N个互不相交的类Ψ_t,t＝1,2..,N，每个类代表一种未知的材料，即N个衰减系数μ(t,E),t＝1,2..,N。

将上述双能/能谱CT各能量下重建的CT图像聚类结果，Ψ_t,t＝1,2..,N，以及材料衰减系数的分解模型μ(t,E)＝b₁(t)μ₁(E)+b₂(t)μ₂(E)，其中，μ(t,E)表示衰减系数函数，b₁(t)表示分解系数，μ₁(E)表示衰减系数，代入到X射线的非线性投影公式中得到基于图像域材料聚类的迭代求解问题，如下式(3)所示：

其中，

表示高低能投影数据，m表示不同的入射能谱的索引，例如m＝L时表示低能X射线，m＝H时表示高能X射线，S_m(E)表示不同的入射X射线等效能谱,

表示高低能不同能谱下的最高能量，ray_i表示第i个投影路径，最内层关于dl的积分表示沿X射线路径上对线衰减系数的线积分，而外层则表示对能量E的积分。

为了求解分解系数b_k(t,E),k＝1,2,t＝1,2...,N，上式可以写成加权最小二乘的形式，如下式(4)所示：

其中，

是不同投影路径的权重，

是在投影路径ray_i上X射线束的归一化能谱。

公式(4)可以通过常规的非线性最优化方法求解的，例如Levenberg–Marquardt算法，然后，可以用下面公式(5)求出每个像素/体素的电子密度ρ_e(t)和等效原子系数Z_eff(t)：

其中，ρ_k,Z_k,A_k,k＝1,2分别表示每个像素/体素点上对应两种基材料的密度、原子序数、原子量。

图3是本发明所述迭代材料分解方法在仿真模拟中的结果与现有材料分解方法的结果对比图，其中，利用蒙特卡洛方法模拟9MeV和6MeV双能CT的材料分解结果，模型中不同的圆内填充不同的单质材料，原子序数从13增加到50，大圆柱模型内为水。图3中的第一行为原子序数分解结果，第二行为电子密度分解结果；(a)和(d)是理论的标准值；(b)和(e)是现有双能CT投影域材料分解算法结果，由于投影数据中的噪声导致原子序数分解的误差较大；(c)和(f)是本发明的迭代分解结果。根据图3对比分析本发明和现有分解结果可知，本发明能够有效抑制噪声，且分解的数值也更准确。对比分析本发明和仿真结果可知，本发明能够有效地提高现有X射线双能/能谱CT的材料分解准确性，对于提升当前各类双能/能谱CT的成像质量和物质识别精度具有重要现实意义和应用价值。

本发明所述基于材料分类的能谱CT迭代材料分解方法通过双能/能谱CT图像重建，由不同能量的重建图像进行聚类分析，将成像物体聚类成N个互不相交的类Ψ_t,t＝1,2..,N，每个类代表一种未知的材料，即N个衰减系数函数；然后利用材料衰减系数的分解模型对衰减系数函数进行分解，并由此建立起基于图像域材料聚类的迭代分解目标函数，求解该目标函数即得到最终的双能/能谱CT材料分解或物质识别结果。本发明能够有效地消除现有双能/能谱CT在高噪声或MeV高能双能/能谱成像时的分解误差，大大提高双能/能谱CT材料分解和物质识别的准确性。

图4是本发明所述基于材料聚类的能谱CT迭代材料分解装置的构成框图，如图4所示，本发明所述基于材料聚类的能谱CT迭代材料分解装置，包括：数据采集模块1，用于获取成像物体在不同能量射线下的能谱CT投影数据，可以是双能CT投影数据，也可以是多能CT投影数据，通过不同能量的X射线对成像物体进行扫描得到；图像重建模块2，用于利用所述能谱CT投影数据重建能谱CT图像；聚类模块3，用于对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类，每个类代表一种未知的材料；分解模块4，用于根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果。

在一个实施例中，所述聚类模块包括：能量值计算单元，用于计算不同能量的能谱CT图像对应的像素或体素点的能量值；图像获取单元，用于根据所述像素或体素点的能量值得到多维数据集图像；中心确定单元，用于确定多维数据集图像中的聚类中心；聚类单元，用于根据聚类中心对多维数据集图像中的所有点进行聚类。

进一步地，所述中心确定单元包括：密度半径计算子单元，用于根据密度半径的预设值计算每个点的密度半径；密度增量半径计算子单元，用于根据计算得到的密度半径计算每个点的密度增量半径；密度中心确定子单元，用于根据每个点的密度半径和密度增量半径获取一个或多个密度中心，并将密度中心作为聚类中心，其中，密度中心的数量根据预先估计的材料种类数量确定。例如，如图2所示，以双能CT得到的多维数据集图像为例，该图像中，横坐标为被扫描物体低能下的线衰减系数，纵坐标为被扫描物体在低能和高能下的线衰减系数之差，图中的各聚焦的散点中心的黑色圆点则表示聚类中心，多维数据集图像中包括多个聚类中心。

进一步地，在一个可选的实施例中，聚类单元通过下述方式对多维数据集图像中的所有点进行聚类，包括：以聚类中心为中心，将在密度半径的预设值范围内的所有点归类到所述聚类中心；将未归类到所述聚类中心的点归类到密度增量半径对应的点，从而实现所有点的聚类。

需要说明的是，本发明所述基于材料聚类的能谱CT迭代材料分解装置的具体实施方式与上述能谱CT迭代材料分解方法的具体实施方式大致相同，在此不再一一赘述。

以上所述仅为本发明的优选实施例，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (9)

1.一种基于材料聚类的能谱CT迭代材料分解方法，其特征在于，包括以下步骤：

获取成像物体在不同能量射线下的能谱CT投影数据；

利用所述能谱CT投影数据重建能谱CT图像；

对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类；

根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果；

其中，根据聚类结果进行能谱CT迭代材料分解的步骤包括：

根据聚类结果获取聚类信息；

将聚类信息和材料衰减系数的分解模型带入X射线的非线性投影公式中进行迭代求解，获取材料的分解系数。

2.根据权利要求1所述的基于材料聚类的能谱CT迭代材料分解方法，其特征在于，对重建得到的能谱CT图像进行聚类处理的步骤包括：

计算不同能量的能谱CT图像对应的像素或体素点的能量值；

根据所述像素或体素点的能量值得到多维数据集图像；

确定多维数据集图像中的聚类中心；

根据聚类中心对多维数据集图像中的所有点进行聚类。

3.根据权利要求2所述的基于材料聚类的能谱CT迭代材料分解方法，其特征在于，确定多维数据集图像中的聚类中心的步骤包括：

根据密度半径的预设值计算每个点的密度半径；

根据计算得到的密度半径计算每个点的密度增量半径；

根据每个点的密度半径和密度增量半径获取一个或多个密度中心，并将密度中心作为聚类中心。

4.根据权利要求3所述的基于材料聚类的能谱CT迭代材料分解方法，其特征在于，根据聚类中心对多维数据集图像中的所有点进行聚类的步骤包括：

以聚类中心为中心，将在密度半径的预设值范围内的所有点归类到所述聚类中心；

将未归类到所述聚类中心的点归类到密度增量半径对应的点。

5.根据权利要求1所述的基于材料聚类的能谱CT迭代材料分解方法，其特征在于，所述材料的分解系数采用非线性最优化方法求解。

6.根据权利要求1所述的基于材料聚类的能谱CT迭代材料分解方法，其特征在于，根据聚类结果进行能谱CT迭代材料分解的步骤还包括：

根据材料的分解系数得到每个像素或体素点的电子密度和等效原子系数。

7.根据权利要求1所述的基于材料聚类的能谱CT迭代材料分解方法，其特征在于，所述方法应用于双能CT或多能CT成像的材料分解。

8.一种基于材料聚类的能谱CT迭代材料分解装置，其特征在于，包括：

数据采集模块，用于获取成像物体在不同能量射线下的能谱CT投影数据；

图像重建模块，用于利用所述能谱CT投影数据重建能谱CT图像；

聚类模块，用于对重建得到的能谱CT图像进行聚类处理，形成多个互不相交的类；

分解模块，用于根据聚类结果进行能谱CT迭代材料分解，得到材料分解或物质识别结果；

其中，所述分解模块根据聚类结果获取聚类信息，将聚类信息和材料衰减系数的分解模型带入X射线的非线性投影公式中进行迭代求解，获取材料的分解系数。

9.根据权利要求8所述的基于材料聚类的能谱CT迭代材料分解装置，其特征在于，所述聚类模块包括：

能量值计算单元，用于计算不同能量的能谱CT图像对应的像素或体素点的能量值；

图像获取单元，用于根据所述像素或体素点的能量值得到多维数据集图像；

中心确定单元，用于确定多维数据集图像中的聚类中心；

聚类单元，用于根据聚类中心对多维数据集图像中的所有点进行聚类。

Images