光栅成像系统的信息提取方法及提取装置

摘要

本发明公开了一种光栅成像系统的信息提取方法，包括：获取来自光栅成像系统中的物体位移曲线与背景位移曲线；根据物体位移曲线与背景位移曲线，计算它们各自的积分或总和，并根据物体位移曲线积分或总和与背景位移曲线积分或总和的比值，确定散射角的零阶矩，从而确定吸收信息。本发明方法不需要对投影获得的位移曲线数据进行解卷积迭代获取散射角分布，然后计算其多阶矩进行信息提取，而是直接从投影数据计算图像信息，避免了解卷积迭代过程中花费的冗长计算时间及解卷积方法和停止标准的选择，精简了信息提取过程，同时还获得了更多的图像细节。

权利要求

1.一种光栅成像系统的信息提取方法，其特征在于包括：获取来自光栅成像系统中的物体位移曲线与背景位移曲线；

根据物体位移曲线与背景位移曲线，计算它们各自的积分或者总和，并根据物体位移曲线积分或者总和与背景位移曲线积分或者总和的比值确定散射角的零阶矩，从而确定吸收信息；

根据物体位移曲线与背景位移曲线，计算它们各自的积分或者总和，并根据物体位移曲线积分或者总和与背景位移曲线积分或者总和的比值确定散射角的零阶矩，包括：采用如下公式确定散射角α的零阶矩M0：

其中，T为物体位移曲线s(α)和背景位移曲线f(α)的周期，取变量所述的方法还包括：

根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分或总和，从而进行级数求和，确定一阶矩，从而确定相衬信息；

根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分或总和，从而进行级数求和，确定一阶矩，从而确定相衬信息，包括：采用如下公式确定散射角α的一阶矩M1：

其中：

其中，T为物体位移曲线s(α)和背景位移曲线f(α)的周期，取变量 d为成像系统中相位光栅和分析光栅的间距，n为任意自然数，p2为成像系统中分析光栅的周期；

所述的方法还包括：

根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分，从而计算级数求和，确定二阶矩，从而确定暗场信息；

根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分，从而计算级数求和，确定二阶矩，从而确定暗场信息，包括：采用如下公式确定散射角α的二阶矩M2：

其中，T为物体位移曲线s(α)和背景位移曲线f(α)的周期，取变量 d为成像系统中相位光栅和分析光栅的间距，p2为光栅成像系统中分析光栅的周期；

所述一阶矩M1的确定公式中，2阶及以上的部分被忽略不计，调整确定M1的公式为：其中，上式在运算过程中，

其中，k为任意大于等于2的自然数，φ为任意实数，d为成像系统中相位光栅和分析光栅的间距，p2为成像系统中分析光栅的周期；

所述二阶矩M2的确定公式中，2阶及以上的部分被忽略不计，调整确定M2的公式为：其中，上式在运算过程中，

其中，k为任意大于等于2的自然数，φ为任意实数，d为成像系统中相位光栅和分析光栅的间距，p2为成像系统中分析光栅的周期。

2.根据权利要求1所述的方法，其特征在于，还包括：根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分，从而计算级数求和，确定三阶以上阶矩，从而确定与之对应的物理信息。

3.根据权利要求2所述的方法，其特征在于，根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分，从而计算级数求和，确定三阶以上阶矩，从而确定与之对应的物理信息，包括：采用如下公式确定散射角α的n阶矩Mn：

其中，n为3以上的自然数，T为物体位移曲线s(α)和背景位移曲线f(α)的周期，取变量xn，x∈[-π，π]可以用傅里叶级数展开，从而获得上式的级数解析表达式。

4.一种光栅成像系统的信息提取装置，其特征在于，包括：处理器；

存储器，其上存储有可在所述处理器上运行的计算机程序；其中所述计算机程序被所述处理器执行时实现权利要求1-3所述的方法步骤。

5.一种计算机可读存储介质，其特征在于，所述计算机可读存储介质上存储有数据处理程序，所述数据处理程序被处理器执行时实现权利要求1-3所述的方法步骤。

说明书

光栅成像系统的信息提取方法及提取装置

技术领域

[0001] 本发明涉及辐射成像领域，进一步涉及一种光栅成像系统信息提取方法，以及光栅成像系统的信息提取装置。

背景技术

[0002] X射线光栅相衬成像技术，能够实现微米或亚微米级上的局部结构分辨，是对传统X射线成像技术的很好补充。该技术能够同时提取吸收、相衬和暗场像三种信息，适用于低原子序数、低密度物质间的区分，特别是针对包括乳腺在内的生物软组织结构。为了获得更好的图像质量，光栅成像技术通常通过相位步进的方式来获取数据。该方式指的是将其中一块光栅沿着垂直于光栅栅格线的方向在一个或多个周期内进行等距位移，从而获得放置物体前后的位移曲线。

[0003] 光栅成像技术中的三种图像信息，一般是通过对位移曲线的后处理来获取的。主要方法包括余弦模型分析(cosine-model analysis，CMA)和小角散射(small angle X-ray scattering，SAXS)这两大类信息提取算法，其中前者中的三种图像信息解析表达为传统的基本定义，而后者则通过散射角的卷积假设及其多阶矩来进行信息提取。

发明内容

[0004] 有鉴于此，本发明的目的在于提供一种光栅成像系统的信息提取算法方法及装置，以至少部分解决上述的技术问题。

[0005] 根据本发明的一方面，提供一种光栅成像系统的信息提取方法，包括：

[0006] 获取来自光栅成像系统中的物体位移曲线与背景位移曲线；

[0007] 根据物体位移曲线与背景位移曲线，计算它们各自的积分或总和，并根据物体位移曲线积分或总和与背景位移曲线积分或总和的比值确定散射角的零阶矩，从而确定吸收信息。

[0008] 根据本发明又一方面，提供一种光栅成像系统的信息提取装置，包括：

[0009] 处理器；

[0010] 存储器，其上存储有可在所述处理器上运行的计算机程序；其中所述计算机程序被所述处理器执行时实现以上所述的方法步骤。

[0011] 根据本发明再一方面，提供一种计算机可读存储介质，所述计算机可读存储介质上存储有数据处理程序，所述数据处理程序被处理器执行时实现以上所述的方法步骤。

[0012] 本发明提出的光栅成像系统的信息提取方法及装置，不需要对投影获得的位移曲线数据进行解卷积迭代获取散射角分布，然后计算其多阶矩进行信息提取，而是直接从投影数据计算三种图像信息。

[0013] 本发明提出的光栅成像系统的信息提取方法及装置，避免了解卷积迭代过程中花费的冗长计算时间及解卷积方法和停止标准的选择，避免了迭代过程中部分结构细节的丢失，精简了信息提取过程，提高了图像信息的丰富度和对比度。

[0014] 本发明提出的光栅成像系统的信息提取方法及装置明显提高了计算效率，同时也获得了更多的图像信息。

附图说明

[0015] 图1a和图1b分别为本发明实施例中光栅成像系统及位移曲线示意图；

[0016] 图2为现有技术和本发明实施例的信息提取流程图对比；

[0017] 图3为现有技术和本发明实施例的信息提取方法的三种图像信息的比较；

[0018] 图4为现有技术和本发明实施例的信息提取分方法的三种图像信息重建图像侧视图的比较。

具体实施方式

[0019] 为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明作进一步的详细说明。

[0020] 本发明基于下述构想：现有技术中主要关注SAXS方法，SAXS方法通常需要通过解卷积迭代来获得卷积公式中的散射角分布g(α)，然后计算g(α)的零阶矩、一阶矩和二阶矩来分别获得吸收、相衬和暗场三种图像。由于采用了迭代的方式来解卷积计算散射角分布，通常SAXS方法会花费远比CMA方法要长的多的计算时间，并且解卷积方法和停止标准也需要预先设定，且会对迭代结果产生重大影响。本发明人通过推导SAXS方法中三种图像信息的解析表达式，从而可以直接从扫描获得的数据直接计算三种图像信息，避免了计算时间过于冗长的问题，同时还获得了更好的图像细节和对比度。

[0021] 本发明实施例中，提供一种光栅成像系统的信息提取方法，包括：

[0022] 获取来自光栅成像系统中的物体位移曲线与背景位移曲线；

[0023] 根据物体位移曲线与背景位移曲线，计算它们(它们是指物体位移曲线与背景位移曲线)各自积分或总和(离散的情况为总和，连续的情况为积分)，并根据物体位移曲线积分或总和与背景位移曲线积分或总和的比值确定散射角的零阶矩，从而确定吸收信息。一种具体的实施例为：根据物体位移曲线与背景位移曲线，确定散射角d的零阶矩M0，其中确定方式采用如下公式：

[0024]

[0025] 其中 为背景位移函数， 为物体位移函数，T为周期性函数 和 的周期。

[0026] 对于相位步进式光栅成像系统，如图1a和图1b所示，本发明实施例的信息提取方法是基于相位步进式的光栅成像系统而提出的，该系统及位移曲线示意图见附图1a和1b.该系统主要包括X射线球管、探测器、源光栅G0、相位光栅G1及分析光栅G2部分。被检查物体按照附图1a中的方向进行相位步进扫描，从而获得附图1b中的放置物体前后的位移曲线。需要注意的是，光栅成像系统通常通过相位步进方式来进行扫描成像，相位步进过程中移动三个光栅的任意一个均可完成扫描，最常见的是移动分析光栅G2，以下说明内容均以移动G2为例。

[0027] 其中，当采用SAXS进行信息提取时，基本假设是：物体位移曲线(放置物体)s(α)为背景位移曲线(不放置物体)f(α)和散射角分布g(α)的卷积：

[0028]

[0029] 其中，α为散射角，定义是为 Δl为相位步进的距离，d为相位光栅G1和分析光栅G2的间距。

[0030] 通常，吸收(A)、相衬(P)和暗场(P)分别对应着散射角分布g(α)的零阶矩、一阶矩和二阶矩：

[0031]

[0032]

[0033]

[0034] SAXS方法中，首先根据公式(1)对物体位移曲线和背景位移曲线的投影数据进行解卷积迭代来获取散射角分布g(α)，采用的迭代方法一般为Richard-Lucy；然后根据公式(2)-(4)来计算三种图像信息。该方法优势在于不需要对位移曲线进行模型假设即可进行信息提取，且获取的图像噪声水平较低，限制因素在于解卷积迭代过程计算时间过长，以及迭代可能导致部分结构细节的损失。

[0035] 为了克服SAXS方法中的限制因素，本申请发明人提出了一种信息提取方法，该方法可称之为ASAXS(analytic small angle X-ray scattering，解析式小角散射)方法，以下对该方法进行介绍。

[0036] 在光栅成像技术中，位移曲线s(α)和f(α)为周期性信号，周期为 p2为分析光栅G2的周期，d为相位光栅G1和分析光栅G2的间距。为了方便以下推导，采用新变量来替代散射角α，且有x∈[-π，π]。根据公式(1)的傅里叶变换性质，可以得到[0037]#(5)

[0038] 其中n为任意自然数。

[0039] 当n＝0时，从公式(5)可以得到g(α)的零阶矩：

[0040]

[0041] 一些实施例中，还根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分，从而计算特定级数求和，确定一阶矩，从而确定相衬信息。以及根据物体位移曲线与背景位移曲线，计算它们与多阶正弦、余弦函数乘积的积分，从而计算特定级数求和，确定二阶矩，从而确定暗场信息。具体的实施例可以为：

[0042] 当n≥1时，从公式(5)可以得到

[0043]

[0044]

[0045] 其中，

[0046]

[0047]

[0048]

[0049]

[0050] 在计算公式(3)和(4)中的一阶矩和二阶矩前，本申请发明人先将函数x和x2写成傅里叶级数形式(其中x∈[-π，π])：

[0051]

[0052]

[0053] 将公式(7)-(9)代入公式(3)和(4)中，可以得到

[0054]

[0055]

[0056] 考虑到位移曲线可以近似看成一阶余弦曲线，且有：

[0057]

[0058]

[0059] 其中k为任意大于等于2的自然数，φ为任意实数。由此，公式(10)和(11)可以近似为：

[0060]

[0061]

[0062] 以上两个公式中2阶及以上的部分被忽略不计。

[0063] 总而言之，通过公式(6)、(13)和(14)，本申请发明人可以获得SAXS方法中三种信息的直接解析表达式，从而可以从原始的扫描数据中直接进行信息提取，而不需要SAXS中的解卷积过程，该算法流程如附图2所示。

[0064] 由此可以进一步引申到更高阶次的矩(Mn，n≥3)，其表达式为

[0065]

[0066] 同样，公式(15)中，α用 x∈[-π，π]来替代，得到

[0067]

[0068] 考虑到函数xn，x∈[-π，π]可以用傅里叶级数展开，从而获得公式(16)的级数表达式，最后结合公式(12)和近似余弦函数的位移曲线，从而获得高阶矩的级数解析表达式。

[0069] 此外，计算散射角分布g(α)的n阶矩时，公式(16)将会表示为无穷级数求和的形式(类似公式(10)和(11))。考虑到位移曲线并非是完美的一阶余弦曲线，阶次n≥2的级数还可能存在图像信息，因此本方法也适用于对前N项级数求和，其中N＝1，2，...，+∞。

[0070] 本发明实施例的方法针对光栅相衬成像技术采用相位步进方式采集数据的特点，在不进行解卷积迭代求取散射角分布的情况下，直接对原始数据进行计算，高效率地获得了吸收、相衬和暗场三种图像信息。本方法与常规SAXS方法相比，极大地提高了计算的效率，降低了SAXS算法对局部细节的损耗，提高了图像对比度和分辨能力。

[0071] 采用本发明实施例在相位步进式光栅成像系统中，ASAXS信息提取流程如图2所示。

[0072] 首先，通过公式(6)、(13)和(14)直接计算吸收A、相衬P和暗场D三种图像信息。

[0073] A＝-ln M0#(17)

[0074] P＝M1#(18)

[0075] D＝M2#(19)

[0076] 然后，通过滤波反投影等算法重建CT图像。

[0077] 部分实验结果如图3和图4所示，实验样本为石蜡固定的人体切除乳腺标本。光栅成像系统中，源光栅到相位光栅间距为1020mm，源光栅到分析光栅间距为1280mm，源光栅、相位光栅和分析光栅周期分别为16.8μm、4.2μm和2.4μm，X射线球管参数为35kVp和35mA，探测器像素为74.8μm×74.8μm，相位步进数为8。

[0078] 图3为SAXS和ASAXS两种方法提取的三种图像信息的比较，其中图上比例尺为5mm，图像尺寸为516×591像素。图4为SAXS和ASAXS两种方法提取的三种图像信息重建图像侧视图的比较，其中图上比例尺为5mm，建像素尺寸为70μm×70μm，图像尺寸为526×472像素。实验中共有2,880个投影图像，SAXS方法提取图像花费了1，469.50秒，而ASAXS方法则只花了164.28秒，效率提高了将近9倍。两种方法提取的吸收像及其重建图像相互符合得很好，对应SSIM(structural similarity，结构相似性)值均为1.00；相衬像及其重建图像比较接近，对应SSIM值分别为0.87和0.94；暗场像及其重建图像差异较大，对应SSIM值分别为0.64和0.72。如图3中实线框所示，在相衬和暗场像中，ASAXS方法提供了更多的局部结构细节，同时对应的重建图像信号强度也更高。这可能是由于SAXS方法中解卷积迭代过程的平滑效果所致。此外，还对两种方法提取的对应图像的CNR(contrast-to-noise ratio，对比噪声比)值进行了比较。图3和图4中的虚线框表示计算CNR选择的区域，可以得到，相衬像及其重建图像中ASAXS对应的图像CNR值稍微高一些(ASAXS：1.97、6.19；SAXS：1.90、5.78)；而暗场信息中ASAXS对应的图像CNR值则明显高于SAXS(ASAXS：1.42、0.98；SAXS：0.99、0.52)。最后，需要指出的是，两种方法提取的相衬重建图像中，均可以清楚分辨出肿瘤区域(图4中黑色箭头所示)。

[0079] 总而言之，该方法不仅能够大幅度降低信息提取计算时间，还能够提供更多的图像细节。需要指出的是，该方法也会导致图像噪声水平的提升，不过与图像对比度、分辨能力的提升，是影响不大的。

[0080] 根据本发明实施例的另一方面，提供一种光栅成像系统信息提取装置，其特征在于，包括处理器和存储器：其中存储器上存储有可在所述处理器上运行的计算机程序；其中所述计算机程序被所述处理器执行时实现以上实施例所述的方法步骤。本发明实施例中的光栅成像系统信息提取装置包括处理器(例如，微处理器(μP)、数字信号处理器(DSP)等)。处理器可以是用于执行本发明实施例描述的方法流程的不同动作的单一处理单元或者是多个处理单元。此外，装置可以包括具有非易失性或易失性存储器形式的至少一个计算机可读存储介质，以下将具体介绍。

[0081] 根据本发明实施例的又一方面，还提供一种计算机可读存储介质，所述计算机可读存储介质上存储有数据处理程序，所述数据处理程序被处理器执行时实现以上实施例所述的方法。本发明实施例的信息提取方法可以采取存储有指令的计算机可读存储介质上的计算机程序产品的形式，该计算机程序产品可供指令执行系统(例如，一个或多个处理器)使用或者结合指令执行系统使用。在本发明实施例的上下文中，计算机可读存储介质可以是能够包含、存储、传送、传播或传输指令的任意介质。例如，计算机可读存储介质可以包括但不限于电、磁、光、电磁、红外或半导体系统、装置、器件或传播介质。计算机可读存储介质的具体示例包括：磁存储装置，如磁带或硬盘(HDD)；光存储装置，如光盘(CD-ROM)；存储器，如随机存取存储器(RAM)或闪存；和/或有线/无线通信链路。

[0082] 以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。