能谱CT成像系统及数据采集和重建能谱CT图像的方法与流程

本公开涉及CT成像技术，具体而言涉及能谱CT成像系统及其数据采集方法，以及重建能谱CT图像的方法。
背景技术：
：X射线CT成像技术在医疗、安检、工业无损检测等领域有着非常广泛的应用。能谱CT是近年来受到普遍关注和迅速发展的研究方向。目前基于光子计数探测器的能谱CT系统是能谱CT的主要实现方式之一。国际上有多个研究单位和公司致力于研究和生产适用于X光CT的光子计数探测器。光子计数探测器对接收到的光子进行处理，按照能量阈值把接收到的光子归到相应的计数器，从而使得可以选择能量窗采集多个能量下的X射线穿过物体的信号数据。通过对应的能谱CT图像重建算法可以得到物体的三维物理信息，其中包括线性衰减系数、电子密度、等效原子序数、材料分布等等。因此，能谱CT相在功能上能够比传统CT提供被成像体更多层面的信息。另一方面，由于计数方式可以通过阈值方式消除电子学噪声，进一步提高数据的可靠性和信噪比。因此，能谱CT的图像可以具有更高的图像质量，或者在同等图像质量情况下降低剂量。光子计数探测器的性能是影响能谱CT成像效果的关键因素之一。目前应用于X光CT的像素化光子计数面阵探测器的像素分辨率主要在50-100微米之间，光子计数率在106-108/mm2之间。一方面，为了获得高质量的重建图像，需要采集足够多的光子以降低量子噪声的影响；另一方面，每个像素计数器的计数速度受电子学的硬件信号处理速度限制。由此，虽然减小像素大小有助于获得更高的每平方毫米计数率，但是减小像素使得像素之间的电荷共享问题严重化。而且，像素数目的增加对于探测器的读出速率要求也显著提高。因此，需要对能谱CT进行改进或开发信的能谱CT系统。在所述
背景技术：
部分公开的上述信息仅用于加强对本公开的背景的理解， 因此它可以包括不构成对本领域普通技术人员已知的现有技术的信息。技术实现要素：本申请公开一种能谱CT成像系统及其数据采集方法，以及重建能谱CT图像的方法，能够克服上述问题中的一个或多个。本公开的其他特性和优点将通过下面的详细描述变得显然，或部分地通过本公开的实践而习得。根据本公开的一个方面，提供一种能谱CT成像系统的数据采集方法，所述能谱CT成像系统包括具有多个像素的光子探测器,所述数据采集方法包括在多个投影角度从所述光子探测器对穿过被成像物的光子进行数据采集,对于每一投影角度,该数据采集方法包括：按照多个随机分布的像素地址编码从所述光子探测器获取数据。根据本公开的一些实施例，所述像素地址编码利用设定的像素地址分布的概率密度函数而生成。根据本公开的一些实施例，对于第一个投影角度，所述像素地址编码利用设定的像素地址分布的概率密度函数而生成；对于除第一个投影角度外的其余投影角度，所述像素地址编码或者利用为前一投影角度设定的像素地址分布的概率密度函数而生成、或者利用重新设定的像素地址分布的概率密度函数而生成、或者与用于前一投影角度的所述像素地址编码相同。根据本公开的一些实施例，所述按照多个随机分布的像素地址编码从所述光子探测器采集数据包括：从自所述光子探测器获得的全部数据中按照多个随机分布的像素地址编码选取数据。根据本公开的一些实施例，所述按照多个随机分布的像素地址编码从所述光子探测器采集数据包括：从所述光子探测器的全部像素中按照多个随机分布的像素地址编码实时选择像素来获取数据。根据本公开的另一发明，提供一种能谱CT成像系统的数据采集方法，所述能谱CT成像系统包括具有多个像素的光子探测器,所述数据采集方法包括在多个投影角度从所述光子探测器对穿过被成像物的光子进行数据采集,对于每一投影角度,该数据采集方法包括：从按照多个随机分布的像素地址编码设置的所述多个像素采集数据。根据本公开的一些实施例，所述光子探测器为面阵探 测器或环形探测器。根据本公开的另一方面，提供一种重建能谱CT图像的方法，包括：设置数据采集的多个能量窗；利用如前所述的任一数据采集方法从所述多个能量窗进行数据采集；利用所述像素地址编码生成对应的系统矩阵；将采集到的数据在投影域进行能谱信息分解，得到分解系数投影；通过迭代方法，利用所述分解系数投影进行分解系数的空间分布重建，得到分解系数；根据重建的分解系数进行合成，得到单能量的衰减系数和/或估算被成像物的电子密度和/或等效原子序数分布图。根据本公开的一些实施例，所述系统矩阵的维度是M×N，其中N是被重建的图像域像素数目，M是投影角度数与所述像素地址编码的数量之积。根据本公开的一些实施例，根据如下公式求分解系数投影：-ln∫wk(E)exp(-Στ∈[1,Γ]φτ(E)Aτ)dE=pk,Aτ={A1,τ,A2,τ,...,AL,τ}]]>其中τ为整数，对应第τ个分解项，wk(E)为归一化的第k个能量窗下的X光能谱分布，pk为对应于光子探测器上采集的数据的M维向量，φτ(E)为第τ个分解基函数，Aτ为第τ个分解系数投影，k﹦1,2,...K，K为能量窗数目，E是光子能量。根据本公开的一些实施例，根据如下公式通过迭代进行分解系数的空间分布重建：aτ=argminaτ||diag(Φ(i,j,z))▿aτ||1,s.t.||Haτ-Aτ||2≤ϵ]]>其中ε为根据数据噪声确定的阈值，Φ(i,j,z)为关于成像视野的像素位置(i,j,z)的先验函数，aτ为分解系数，H是所述系统矩阵。根据本公开的一些实施例，根据重建结果通过下式计算单能量的衰减系数：μ(E)=Στ∈[1,Γ]aτφτ(E).]]>根据本公开的一些实施例，在两个能量窗的情况下，K＝2，τ∈[1,2],φτ(E)分别取为光电效应系数、康普顿散射系数。根据本公开的一些实施例，等效原子序数分布图、电子密度分布图分别通过下式估算：Zeff=(diag(1a2)a1)1λ-1,λ≈4]]>ρe＝2a2其中a1、a2对应于τ分别为1和2时的aτ值。根据本公开的另一方面，提供一种重建能谱CT图像的方法，包括：设置数据采集的多个能量窗；利用如前所述的任一数据采集方法从所述多个能量窗进行数据采集；利用所述像素地址编码生成对应的系统矩阵；通过迭代方法，对每个能量窗下采集到的数据进行衰减系数的空间分布重建，得到衰减系数；利用衰减系数的重建结果，对被重建的图像域像素逐个进行能谱信息分解，得到分解系数；估算被成像物的电子密度和/或等效原子序数分布图。根据本公开的一些实施例，所述系统矩阵的维度是M×N，其中N是被重建的图像域像素数目，M是投影角度数与所述像素地址编码的数量之积。根据本公开的一些实施例，根据如下公式通过迭代进行衰减系数重建：gk=-log(pkpkair),]]>及μk=argminμk||diag(Φ(i,j,z))▿μk||1,s.t.||Hμk-gk||2≤ϵ]]>其中k表示第k个能量窗，gk为对采集到的数据进行归一化和负对数处理得到的投影值，μk为衰减系数，pk为对应于光子探测器上采集的数据的M维向量，ε为根据数据噪声确定的阈值，Φ(i,j,z)为关于成像视野的像素位置(i,j,z)的先验函数，H是所述系统矩阵。根据本公开的一些实施例，根据重建结果μk，通过求解如下线性方程组对被重建的图像域像素逐个进行能谱信息分解而得到分解系数：μ1(i,j,z)=Στ∈[1,Γ]aτ(i,j,z)φτ(E1)...μk(i,j,z)=Στ∈[1,Γ]aτ(i,j,z)φτ(Ek)...μk(i,j,z)=Στ∈[1,Γ]aτ(i,j,z)φτ(EK)]]>其中τ为整数，对应第τ个分解项，k＝1,...,K，K为能量窗数目，φτ(E)为第τ个分解基函数，aτ(i,j,z)为第τ个分解系数，Ek表示第k个能量窗的等效光子能量。根据本公开的一些实施例，在两个能量窗的情况下，K＝2，τ∈[1,2],φτ(E)分别取为光电效应系数、康普顿散射系数。根据本公开的一些实施例，等效原子序数分布图、电子密度分布图分别通过下式估算：Zeff=(diag(1a2)a1)1λ-1,λ≈4]]>ρe＝2a2其中a1、a2对应于τ分别为1和2时的aτ值。根据本公开的另一方面，提供一种能谱CT成像系统，包括：射线发生装置，包括射线源；光子探测器，包括多个像素；及数据采集系统，从所述光子探测器对穿过被成像物的光子进行数据采集，其特征在于，所述数据采集系统按照多个随机分布的像素地址编码从所述光子探测器获取数据。根据本公开的一些实施例，所述像素地址编码利用设定的像素地址分布的概率密度函数而生成。根据本公开的一些实施例，所述数据采集系统从所述光子探测器输出的全部数据中按照多个随机分布的像素地址编码选取数据。根据本公开的一些实施例，所述数据采集系统包括电子学系统，所述电子学系统配置为从所述光子探测器的全部像素中按照多个随机分布的像素地址编码实时选择像素来读取数据。根据本公开的一些实施例，能谱CT成像系统还包括数据处理系统，用于根据前述的任一项重建能谱CT图像的方法利用所述数据采集系统获取的数据重建能谱CT图像。根据本公开的一些实施例，所述光子探测器为面阵探测器或环形探测器。根据本公开的另一方面，提供一种能谱CT成像系统，包括：射线发生装置，包括射线源；光子探测器，包括多个像素；及数据采集系统，从所述光子探测器对穿过被成像物的光子进行数据采集。按照多个随机分布的像素地址编码设置所述光子探测器的所述多个像素。根据本公开的一些实施例，所述像素地址编码利用设定的像素地址分布的 概率密度函数而生成。根据本公开的一些实施例，能谱CT成像系统还包括数据处理系统，用于根据前述任一重建能谱CT图像的方法利用所述数据采集系统获取的数据重建能谱CT图像。根据本公开的一些实施例，所述光子探测器为面阵探测器或环形探测器。根据本公开的一些实施例的能谱CT成像系统及其数据采集方法以及重建能谱CT图像的方法，能够有效控制电荷共享效应，并提高CT系统处理速度，以及降低系统成本。附图说明通过参照附图详细描述其示例实施例，本公开的上述和其它特征及优点将变得更加明显。图1示意性示出根据本公开一些示例实施例的能谱CT成像系统的系统结构示意图；图2示意性示出根据本公开一些示例实施例的数据采集流程图；图3a和图3b示意性示出根据本公开一些示例实施例的随机分布或随机选择的探测器像素，其中图3a为随机均匀分布模式，图3b为按照特定概率密度例如中间概率大、边缘概率小的随机分布模式；图4a和图4b示意性示出根据本公开一些示例实施例的投影数据比较，其中图4a示出普通能谱CT采集到的一个能量窗的投影数据，图4b示出根据随机分布的像素地址编码采集前者十分之一数据得到的一个能量窗下的投影数据；图5a和图5b示出根据本公开一些示例实施例进行能谱CT重建的结果示例，其中图5a是电子密度分布图，图5b是等效原子序数分布图；及图6示出根据本公开一实施例的能谱CT结构示意图，其包括具有预定义的随机分布像素的探测器。具体实施例现在将参考附图更全面地描述示例实施例。然而，示例实施例能够以多种形式实施，且不应被理解为限于在此阐述的实施例；相反，提供这些实施例使 得本公开将全面和完整，并将示例实施例的构思全面地传达给本领域的技术人员。在图中相同的附图标记表示相同或类似的部分，因而将省略对它们的重复描述。此外，所描述的特征、结构或特性可以以任何合适的方式结合在一个或更多实施例中。在下面的描述中，提供许多具体细节从而给出对本公开的实施例的充分理解。然而，本领域技术人员将意识到，可以实践本公开的技术方案而没有所述特定细节中的一个或更多，或者可以采用其它的方法、组元、材料、装置、步骤等。在其它情况下，不详细示出或描述公知结构、方法、装置、实现、材料或者操作以避免模糊本公开的各方面。附图中所示的方框图仅仅是功能实体，不一定必须与物理上独立的实体相对应。即，可以采用软件形式来实现这些功能实体，或在一个或多个硬件模块或集成电路中实现这些功能实体，或在不同网络和/或处理器装置和/或微控制器装置中实现这些功能实体。附图中所示的流程图仅是示例性说明，不是必须包括所有的步骤。例如，有的步骤还可以分解，而有的步骤可以合并或部分合并，因此实际执行的顺序有可能根据实际情况改变。图1示意性示出根据本公开一些示例实施例的能谱CT成像系统100的系统结构示意图。如图1所示，能谱CT成像系统100可包括射线发生装置105、机械运动系统110、光子探测器120、以及数据采集系统130。根据本公开的系统可以通过圆轨道扫描实现，也可以通过螺旋轨迹扫描实现，可用于三维能谱CT成像。射线发生装置105可包括射线源，用于发射例如X射线。机械运动系统110用于使被成像物与所述射线源发生相对运动。机械运动系统110可例如包括机械运动装置和对应的控制系统(未示出)。可以是被成像物体115运动而射线源和/或探测器120保持静止(图1所示方式)，也可以是射线源和/或探测器运动，而物体保持静止。一般在医疗领域中避免转动病人，可通过转动源和/或探测器实现。在工业无损检测中，转动和平移物体的方式比较常见。对于CT成像，起作用的是相对运动，所以两种方式等效。光子探测器120可包括多个像素1202(如图3a所示)，用于接收并处理穿过 物体的光子。数据采集系统130从所述光子探测器120对穿过被成像物115的光子进行数据采集。根据本公开的发明构思，数据采集系统130按照多个随机分布的像素地址编码从光子探测器120获取数据，如下面将要详细描述的。根据本公开的一些实施例，光子探测器120包括电子学系统125。数据采集系统130从电子学系统125输出的全部数据中按照多个随机分布的像素地址编码选取数据。根据本公开的另一些实施例，数据采集系统130包括电子学系统125。电子学系统125配置为从光子探测器120的全部像素中按照多个随机分布的像素地址编码实时选择像素来读取数据。根据本公开的另一些实施例，如图3a、图3b和图6所示，可按照多个随机分布的像素地址编码设置光子探测器120的多个像素1205。另外，如图3a、图3b和图6所示，光子探测器120可为例如面阵探测器或环形探测器。能谱CT成像系统100还可包括主控制器140及数据处理装置135。主控制器140负责能谱CT系统运行过程的主控制，包括机械转动、电气控制、安全连锁控制等。数据处理装置135对由数据采集系统130获得的数据进行处理，获得能谱CT任意能量下的衰减系数图像。另外，也可由此计算等效原子序数和电子密度分布图。这些图可通过断层或者三维可视化方式在显示器上显示。主控制器140及数据处理装置135可以是单个PC，也可以是工作站或计算机集群。图2示意性示出根据本公开一些示例实施例的数据采集流程图。下面参照图2、图3a和图3b描述根据本公开的从光子探测器120获取数据的方法。根据本公开的发明构思的数据采集方法，对于每一投影角度，该方法包括：按照多个随机分布的像素地址编码从所述光子探测器获取数据。根据本公开的一些实施例，光子探测器120包括电子学系统125，该方法包括从自所述光子探测器获得的全部数据中按照多个随机分布的像素地址编码选取数据。根据本公开的另一些实施例，数据采集系统130包括电子学系统125。电子学系统125配置为从光子探测器120的全部像素中按照多个随机分布的像素地址编码实时选择像素来读取数据。该方法包括从所述光子探测器的全部像素中 按照多个随机分布的像素地址编码实时选择像素来获取数据。参照图2，在各个投影角度，可以按照一定的概率产生Ns个随机分布的像素地址编码，按照这些地址编码从探测器读取数值。这Ns个探测器像素的地址可以通过参数设置来控制其在整个探测器上的分布模式，即设置像素地址的概率密度函数。随机分布的被读取像素的具体示例如图3a和3b所示。虚线网格表示实际的探测器单元矩阵，深灰色小方格表示在某一角度数据采集过程中选择的数据输出像素。图3a为按照均匀分布随机采样。图3b为按照特定概率密度，例如中间概率大、边缘概率小的方式进行的随机采样。记探测器的总像素单元数目为Ndet，数据采集过程中每个视角读取的数据数目为Ns≦Ndet。根据本公开的技术方案，Ns可远小于Ndet，例如小于Ndet的50％，或更进一步地，小于Ndet的20％。如图2所示，根据本公开一些实施例的方法包括如下处理。数据采集开始205之后，判断是否重新设置地址分布概率210。如果为是，则设置地址分布的概率密度函数215，并进行生成地址编码225，然后按照地址编码读取探测器单元数据230。如果为否，则判断是否更新地址编码220。如果为是，则转到225，否则转到230。例如，对于除第一个投影角度外的其余投影角度，像素地址编码可以利用为前一投影角度设定的像素地址分布的概率密度函数而生成(转到225)；或者，与用于前一投影角度的所述像素地址编码相同(转到230)。在按照地址编码读取探测器单元数据230之后，判断是否完成所有角度扫描。如果为否，则对下一个角度判断是否重新设置地址分布概率210。如果为是，则完成数据采集240。下面说明根据本公开的利用通过前述方法采集的数据进行能谱CT图像重建的方法。设采集数据的多个能量窗(或谱)wk(E)，k﹦1,2,...K。K为能量窗数目，K＝2时为双能成像系统。各个wk(E)覆盖的能量范围可以有重叠部分也可以完全分开。多能CT的数据分别为用pk表示：-ln∫wk(E)exp(-Hμ(E))dE＝pk(1)其中，E是光子能量。wk(E)为归一化的第k个能量窗下的X光能谱分布，可以是多种方式产生，例如设置光子计数探测器的能窗阈值，或者在光源处使用不同的滤光片等。μ(E)是物体的衰减系数，H是M×N维投影矩阵，pk为M维向量，即从光子计数探测器的随机分布的像素上采集M条射线的数据。在本公开中，这样的光子计数探测器可称为随机地址编码的光子计数探测器Coded-PCD。能谱CT的重建从框架上可以有三种方式(参见Y.Xing,etal.,"AReconstructionMethodforDualHigh-EnergyCTWithMeVX-Rays,"NuclearScience,IEEETransactionson,vol.58,pp.537-546,2011)：(1)预处理方式：首先把采集到的数据在投影域进行能谱信息分解，然后通过迭代进行空间信息重建，再根据重建的空间信息进行合成；(2)后处理方式：首先对于每一个能谱下的数据进行空间信息重建，然后在图像域进行能谱信息分解；(3)综合处理方式：建立综合的数据模型，通过迭代重建，同时获得能谱相关和空间位置相关信息。根据本公开的能谱CT重建可以在预处理和后处理方法框架下完成。其中的能谱信息分解方式以及合成可以使用领域内现有的技术完成(参见G.Zhang,etal.,"Apracticalreconstructionmethodfordualenergycomputedtomography,"JournalofX-rayScienceandTechnology,vol.16,pp.67-882008以及Y.Xing,etal.,"AGeneralAdaptiveDecompositionmethodforMulti-EnergySpectralCT,"inNuclearScienceSymposiumandMedicalImagingConference(NSS/MIC),2011IEEE,Seul,Korea,2013,pp.M12-15)。空间信息重建方法可结合光子探测器的数据分布情况针对性地如下进行处理。记空间信息分布的投影数据为y，待重建的空间信息为x，根据光子探测器的数据采集方式，有：y＝Hx。(2)此处H的维度是M×N，M＝投影角度数×Ns，N是被重建的图像域像素数目。(2)式异于常规CT的空间信息重建之处在于系统矩阵H对应的是随机分布的探测器单元位置，而且可以是稀疏采样，即M<<N。即，每个角度下数据采集的探测器单元数目可以远小于总的探测器单元数目Ndet。当M<<N时，(2)是一个不定方程组求解问题。在本公开中，通过稀疏条件限制求解：x=argmin(|diag(Φ(i,j,z))▿x|),s.t.||y-Hx||2≤ϵ---(3)]]>ε为由数据噪声情况确定的常数，(i,j,z)为被重建像素的三维位置坐标，Φ(i,j,z) 是与成像视野的像素位置(i,j,z)有关的先验函数，可以预先设定或者通过探测器数据采样分布的概率密度函数计算得到(参见E.Y.Sidky,etal.,"Accurateimagereconstructionfromfew-viewsandlimited-angledataindivergent-beamCT,"JournalofX-RayScienceandTechnology,vol.14,pp.119-139,2006)。|·|为1阶范数：|x|＝Σ|xi|。||▿x||1=TV(x)=Σs,t(xs,t-xs,t-1)2+(xs,t-xs-1,t)2]]>根据拉格朗日乘子法，(3)可以通过下式求解：x=argmin(|diag(Φ(i,j,z))▿x|+λ||y-Hx||2)λ>0---(4)]]>此处λ是调节稀疏条件作用和数据保真性的因子，可以根据经验设置为一个常数。当Ns与Ndet接近时可以提高λ的取值。公式(3)和(4)的数值求解可以通过使用本领域公知的数值优化方法如ART-TV(参见Y.Sidky,etal.,"Accurateimagereconstructionfromfew-viewsandlimited-angledataindivergent-beamCT,"JournalofX-RayScienceandTechnology,vol.14,pp.119-139,2006)或者ADMM(参见J.YangandY.Zhang,"AlternatingDirectionAlgorithmsforl1-problemsinCompressiveSensing,"SIAMJ.Sci.Comput.,vol.33,pp.250-278,2011)的方法迭代完成。下面说明根据本公开一些实施例的利用通过前述方法采集的数据进行能谱CT图像重建的方法。根据一些实施例，采用图1所示的能谱CT系统架构，使用普通X光机，采集2个能量窗的数据，一圈采集Nθ个角度，每个角度使用同样参数的高斯分布的概率密度函数，但生成不同的样本，即不同的地址编码，进行探测器读数，且Ns﹦Ndet/10。图4a和图4b示意性示出根据本公开一些示例实施例的投影数据比较，其中图4a示出普通能谱CT采集到的一个能量窗的投影数据，图4b示出根据随机分布的像素地址编码采集前者十分之一数据得到的一个能量窗下的投影数据。图5a和图5b示出根据本公开一些示例实施例进行能谱CT重建的结果示例，其中图5a是电子密度分布图，图5b是等效原子序数分布图。下面分别以预处理能谱CT重建的方法框架、及后处理能谱CT重建的方法框架为例进行说明。预处理方式进行能谱CT重建1)按照概率密度函数生成地址编码，按照地址编码进行数据采集。记录地址编码，生成对应的系统矩阵H。但本公开不限于此，例如，也可以根据记录的地址编码，在重建过程中实时计算系统矩阵的元素值。系统矩阵H的生成为本领域所公知，此处不再赘述。2)计算Φ(i,j,z)。假设定义各个投影角度θ下探测器数据读取的采样点分布的概率密度函数为pdf(u,v,θ)，u,v表示探测器位置，θ表示投影角度，可以设(但是不限于)Φ(i,j,z)＝[H0T×pdf(u,v,θ)]-p,1/3＜p＜3，其中H0为采集全部探测器数据情况下的系统矩阵。此步骤通过一次反投影和结果求幂次运算完成。特殊情况下，例如无任何先验知识且均匀分布随机采样的情况下，可以设Φ(i,j,z)＝1。3)从2个能量窗的所有探测器数据进行解析得到分解系数投影。记得到的分解系数投影为Aτ。用公式表示为根据下式求Aτ：-ln∫wk(E)exp(-Στ∈[1,Γ]φτ(E)Aτ)dE=pk,]]>Aτ＝{A1,τ,A2,τ,…,AL,τ}。τ为整数，对应第τ个分解项，在2个能量窗的情况下τ∈[1,2]。φτ(E)为第τ个分解基函数，Aτ为第τ个分解系数投影。该步骤可以使用本领域内现有的方法完成，例如双效应分解或者基材料分解。此处可选择φ1,2(E)为光电效应系数、康普顿散射系数，但本公开不限于此。4)使用3)得到的结果进行分解系数重建，得到分解系数图像，即对所有Aτ,使 用如下方式重建获得aτ：aτ=argminaτ||diag(Φ(i,j,z))▿aτ||1,s.t.||Haτ-Aτ||2≤ϵ]]>此处ε为一个小的阈值。此处给出一种具体实现过程示例，但不排斥其它迭代方法实现这一步骤。因为对每个aτ均使用下面的过程得到，所以在下面的叙述中省略了下标τ。a.设定迭代初值为a0；b.进行保真项更新，即计算afitn+1=afitn+Am-Hmafitn||Hm||2(Hm)2,m=1,2,...,M]]>c.非负性约束更新(这一步根据能谱分解方法而定，例如基材料分解，也可以省略)aposn+1=max(0,afitn+1)]]>d.先验约束更新：i)acons0=aposn+1]]>ii)进行Q次全变分最小化迭代(此处可以使用梯度下降法，但不限于这个方法)，Q可为自行选择的整数，例如在5至100之间。aconsq+1=aconsq+α·▿(||diag(Φ(i,j,z))▿aconsq||1),q=1,...,Q]]>e.令再进行b-e步，直到满足收敛条件停止迭代。5)根据上述重建结果aτ，τ＝1,...,Γ，合成计算单能量的μ(E)，μ(E)=Στ∈[1,Γ]aτφτ(E),]]>也可由μ(E)、aτ，τ＝1,...,Γ估计被成像物的电子密度、等效原子序数分布图。在使用光电效应系数、康普顿散射系数作为φ(E)的情况下，Zeff=(diag(1a2)a1)1λ-1,λ≈4]]>ρe＝2a2。后处理方式进行能谱CT重建1)按照概率密度函数生成地址编码，按照地址编码进行数据采集。记录地址编码，生成对应的系统矩阵H。但本公开不限于此，例如，也可以根据记录的地址编码，在重建过程中实时计算系统矩阵的元素值。系统矩阵H的生成为本领域所公知，此处不再赘述。2)计算Φ(i,j,z)。例如定义各个投影角度θ下探测器数据读取的采样点分布的概率密度函数为pdf(u,v,θ)，可以设(但是不限于)Φ(i,j,z)＝[H0T×pdf(u,v,θ)]-p,1/3＜p＜3，其中H0为采集全部探测器数据情况下的系统矩阵。此步骤通过一次反投影和结果求幂次运算完成。特殊情况下，例如无任何先验知识且均匀分布随机采样的情况下，可以设Φ(i,j,z)＝1。3)从2个能量窗的所有探测器数据进行两个线衰减系数的空间分布重建。记对采集到的数据进行归一化和负对数处理得到的投影值为：gk=-log(pkpkair).]]>对所有gk,使用如下方式重建获得μk：μk=argminμk||diag(Φ(i,j,z))▿μk||1,s.t.||Hμk-gk||2≤ϵ]]>其中ε为一个小的阈值。此处给出一种具体实现过程示例，但不排斥其它迭代方法实现这一步骤。因为对每个μk均使用下面的过程得到，所以在下面的叙述中省略了下标k。a.设定迭代初值为μ0；b.进行保真项更新，即计算μfitn+1=μfitn+gm-Hmμfitn||Hm||2(Hm)T,m=1,2,...,M]]>c.非负性约束更新μposn+1=max(0,μfitn+1)]]>d.先验约束更新：i)μcons0=μposn+1]]>ii)进行Q次全变分最小化迭代(此处可以使用梯度下降法，但不限于这个方法)，Q可为自行选择的整数，例如在5至100之间。μconsq+1=μconsq+α·▿(||diag(Φ(i,j,z))▿μconsq||1),q=1,...,Q]]>e.令再进行b-e步，直到满足收敛条件停止迭代。4)根据上述重建结果μk，k＝1,...,K,对被重建的图像域像素逐个像素点进行能谱信息分解得到分解系数aτ(i,j,z)，即求解线性方程组：μ1(i,j,z)=Στ∈[1,Γ]aτ(i,j,z)φτ(E1)...μk(i,j,z)=Στ∈[1,Γ]aτ(i,j,z)φτ(Ek)...μk(i,j,z)=Στ∈[1,Γ]aτ(i,j,z)φτ(EK)]]>τ为整数，对应第τ个分解项。Ek表示第k个能量窗的等效光子能量。对所有像素点完成求解后，得到aτ。φτ(E)为第τ个分解基函数。在使用(但不限于)光电效应系数、康普顿散射系数作为φτ(E)的情况下，可以得到等效原子序数和电子密度分布图：Zeff=(diag(1a2)a1)1λ-1,λ≈4]]>ρe＝2a2。图6示出根据本公开一实施例的能谱CT结构示意图，其包括具有预定义的随机分布像素的探测器。根据本公开的发明构思，可以采用具有预定义的随机分布像素的探测器。例如，参照图3a和3b，虚线网格处可不设置探测器单元(或像素)，仅在深灰色小方格处设置探测器单元(或像素)。即，可按照多个随机分布的像素地址编码设置光子探测器的多个像素。像素地址编码可利用设定的像素地址分布的 概率密度函数而生成，如前所述。参照图6，该示例中没有采用一块完整的面阵光子计数探测器进行CT扫描实时编码采样的方式，而是采用另一种方式实现对投影数据的编码采集。如图6所示的CT系统结构，光子计数探测器以预先确定好的采样模式，例如利用设定的像素地址分布的概率密度函数而生成的采样模式，分布在一个环形探测器605的支架上。此外还有一与该环形探测器支架同轴的环形X射线源610。其中，环形X射线源610可以使用多个普通X射线源沿圆环排列的方式，也可以采用第五代电子束CT的环形靶组成。这种系统结构没有旋转结构，可以使用更高的速率采集投影数据，克服物体运动的模糊，提高CT重建图像的时间分辨率。对于这种结构，其数据处理方式基本与上述方式相同，只有在计算Φ(i,j,z)时，由于探测器环上光子计数像素已经是预先定义好的，所以可以提前计算好存储下来，在重建时直接使用。通过以上的详细描述，本领域的技术人员易于理解，根据本公开实施例的系统和方法具有以下优点中的一个或多个。本公开通过一种新的探测器设计方式和对应的能谱CT重建方法，提出了一种新的能谱CT成像系统。通过低采样率随机编码方式读取光子计数探测器的像素数据，降低探测器的电荷共享效应和探测器的数据读取速率需求。可以通过优化的能谱CT重建方法，获得常规系统固有的分辨率。可以同时实现高分辨率、小数据量和大成像视野。可以通过调整探测器数据采样分布的概率函数，灵活实现扫描视野内差异分辨率成像。通过以上的实施例的描述，本领域的技术人员易于理解，本公开实施例可以通过硬件实现，也可以通过软件结合必要的硬件的方式来实现。因此，本公开实施例的技术方案可以以软件产品的形式体现出来，该软件产品可以存储在一个非易失性存储介质(可以是CD-ROM，U盘，移动硬盘等)中，包括若干指令用以使得一台计算设备(可以是个人计算机、服务器、移动终端、或者网络设备等)执行根据本公开实施例的方法。本领域技术人员可以理解，附图只是示例实施例的示意图，附图中的模块或流程并不一定是实施本公开所必须的，因此不能用于限制本公开的保护范围。本领域技术人员可以理解上述各模块可以按照实施例的描述分布于装置中，也 可以进行相应变化位于不同于本实施例的一个或多个装置中。上述实施例的模块可以合并为一个模块，也可以进一步拆分成多个子模块。以上具体地示出和描述了本公开的示例性实施例。应该理解，本公开不限于所公开的实施例，相反，本公开意图涵盖包含在所附权利要求的精神和范围内的各种修改和等效布置。当前第1页1 2 3