CN111643104A - Ct散射校正方法及系统 - Google Patents

Ct散射校正方法及系统 Download PDF

Info

Publication number
CN111643104A
CN111643104A CN202010133896.6A CN202010133896A CN111643104A CN 111643104 A CN111643104 A CN 111643104A CN 202010133896 A CN202010133896 A CN 202010133896A CN 111643104 A CN111643104 A CN 111643104A
Authority
CN
China
Prior art keywords
scattering
value
energy spectrum
correction
detector
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010133896.6A
Other languages
English (en)
Other versions
CN111643104B (zh
Inventor
高河伟
邢宇翔
陈志强
张丽
李亮
张涛
周浩
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Tsinghua University
Original Assignee
Tsinghua University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Tsinghua University filed Critical Tsinghua University
Priority to CN202010133896.6A priority Critical patent/CN111643104B/zh
Publication of CN111643104A publication Critical patent/CN111643104A/zh
Application granted granted Critical
Publication of CN111643104B publication Critical patent/CN111643104B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
    • A61B6/5258—Devices using data or image processing specially adapted for radiation diagnosis involving detection or reduction of artifacts or noise
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/02—Arrangements for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
    • A61B6/03—Computed tomography [CT]
    • A61B6/032—Transmission computed tomography [CT]
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/44—Constructional features of apparatus for radiation diagnosis
    • A61B6/4411—Constructional features of apparatus for radiation diagnosis the apparatus being modular

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Medical Informatics (AREA)
  • Radiology & Medical Imaging (AREA)
  • Biomedical Technology (AREA)
  • Biophysics (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Veterinary Medicine (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Optics & Photonics (AREA)
  • Pathology (AREA)
  • Public Health (AREA)
  • Physics & Mathematics (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • General Health & Medical Sciences (AREA)
  • Pulmonology (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Analysing Materials By The Use Of Radiation (AREA)
  • Apparatus For Radiation Diagnosis (AREA)

Abstract

一种CT散射校正方法及系统，该CT散射校正方法包括：基于预先建立的投影数据值与散射值映射关系查找表和测量的投影值得到散射估计值，映射关系中的所述散射值为经过能谱校正后的散射值；或者基于能谱物理对散射值的分布函数进行校正后迭代求解得到散射估计值。上述方法通过在散射估计过程中加入能谱校正的因素，基于能谱硬化对于散射分布影响的不同假设，按照迭代求解的方式或者基于散射遮挡硬件测量的方式可以得到能谱校正后的散射估计值，对于硬件的要求不高，同时还具有迭代效率高、计算准确度较高的优点。

Description

CT散射校正方法及系统
技术领域
本公开属于辐射成像技术领域，涉及一种CT散射校正方法及系统。
背景技术
计算机断层成像(CT，Computed Tomography)是使用X射线束和探 测器扫描物体的某一截面(断层)，利用人体器官或组织对X射线的吸收 系数不同的特性获得反映此截面物理或化学特性的投影数据集合，通过计 算机运算获得截面上任意位置的参数值，并由此得到断层图像。
现有CT系统的主要部件包括：X射线源、探测器、以及旋转装置。 在CT成像过程中，X射线源与探测器围绕物体进行相对旋转运动，以此 获得不同旋转角度下的CT数据，也称CT投影值。在CT成像中，X射线 的衰减服从指数规律，CT投影值是间接获取的，通常需要“负对数”预 处理。
光子散射是CT自诞生以来就存在的影响CT图像质量的基础物理挑 战，会造成图像伪影和CT值不准确等问题。提高CT成像性能的最核心 问题之一是去除或减少射线散射。去除散射的方法大致可以分为：基于硬 件的直接散射测量和基于算法的散射估计两大类，基于硬件的直接散射测 量的方法中的硬件可以利用散射遮挡块/散射遮挡条等，基于算法的散射 估计的方法中的算法包括：基于物理的解析模型玻尔兹曼输运方程/蒙特卡 洛计算、投影域卷积滤波等。一般而言，直接测量类散射校正方法精度高， 但对硬件有额外要求，且往往需要二次扫描，可能增加剂量；而算法估计 类方法对硬件没有额外要求，也不需要二次扫描，但校正效果可能差一些， 或者计算复杂度高。
因此，亟需提出一种对于硬件要求不高且计算较为简单的散射校正方 法。
发明内容
(一)要解决的技术问题
本公开提供了一种CT散射校正方法及系统，以至少部分解决以上所 提出的技术问题。
(二)技术方案
根据本公开的一个方面，提供了一种CT散射校正方法，包括：基于 预先建立的投影数据值与散射值映射关系查找表和测量的投影值得到散 射估计值，映射关系中的所述散射值为经过能谱校正后的散射值；或者基 于能谱物理对散射值的分布函数进行校正后迭代求解得到散射估计值；
其中，所述基于能谱物理对散射值的分布函数进行校正后迭代求解得 到散射估计值，包括：将被扫描物体放于多光源CT成像系统中，获取至 少两个不同位置的光源在探测器同一像素点处的投影数据；对至少两个不 同位置的光源在探测器同一像素点处的投影数据进行能谱校正，获得能谱 校正后探测器像素探测得到的总光子数的表达式，该表达式中包括：经过 能谱校正后散射值的分布函数，其中所述经过能谱校正后散射值的分布函数基于能谱硬化对于散射分布的影响确定；以及根据能谱校正后散射值的 分布函数求解出能谱校正后散射值的表达式，并根据能谱校正后散射值的 表达式进行迭代求解得到散射估计值。
在本公开的一实施例中，所述能谱校正后探测器像素探测得到的总光 子数的表达式为：能谱校正后探测器像素探测得到的总光子数为能谱校正 后的透射光子数与能谱校正后散射值分布函数的叠加。
在本公开的一实施例中，对应于两个不同位置的光源的情况，所述能 谱校正后探测器像素探测得到的总光子数的表达式满足：
Figure BDA0002395479180000021
其中，
Figure BDA0002395479180000022
表示当光源Pi闪烁时，经过能谱校正后的探测器像素探测得 到的总光子数，i表示不同位置光源的序号，个数为至少两个，这里i＝1，2；
Figure BDA0002395479180000023
表示光源Pi发光时穿过滤波栅格的光子数；fi(·)为光源Pi闪烁时对应的 能谱校正函数；
Figure BDA0002395479180000024
表示光源Pi闪烁时实际测算得到的投影值；
Figure BDA0002395479180000025
表示光 源Pi闪烁时对应的透射投影值；Gi(Is)表示光源Pi闪烁时经过能谱校正之后 关于散射值IS的分布函数，该分布函数的表达式有多种，根据能谱硬化对 于散射分布的影响确定；
对应于个数大于两个的不同位置的光源的情况，按照上述总光子数的 表达式类推增加方程的个数。
在本公开的一实施例中，所述基于预先建立的投影数据值与散射值映 射关系查找表和测量的投影值得到散射估计值，包括：
预先建立投影数据值与散射值之间的映射关系查找表：将被扫描物体 放于多光源CT成像系统中，获取至少两个不同位置的光源在不设置散射 遮挡硬件时，不放置和放置被扫描物体对应在探测器同一像素点处的投影 数据值，不放置被扫描物体对应的投影数据值为光源Pi发光时穿过滤波栅 格的光子数
Figure BDA0002395479180000031
放置被扫描物体对应的投影数据值为探测器像素探测得到 的总光子数
Figure BDA0002395479180000032
i表示不同位置光源的序号，个数为至少两个；获取在设置 散射遮挡硬件时放置被扫描物体对应在探测器同一像素点的透射光子数
Figure BDA0002395479180000033
根据不设置散射遮挡硬件时光源发光时穿过滤波栅格的光子数
Figure BDA0002395479180000034
不设置散射遮挡硬件时探测器像素探测得到的总光子数
Figure BDA0002395479180000035
和设置散射遮 挡硬件时同一探测器像素探测得到的透射光子数
Figure BDA0002395479180000036
得到经过能谱校正后 的探测器像素探测得到的总光子数
Figure BDA0002395479180000037
和透射光子数
Figure BDA0002395479180000038
利用经过能谱校正 后的探测器像素探测得到的总光子数
Figure BDA0002395479180000039
减去经过能谱校正后的探测器像 素探测得到的透射光子数
Figure BDA00023954791800000310
求得散射值，该散射值为经过能谱校正后的散 射值；建立探测器像素探测得到的总光子数
Figure BDA00023954791800000311
和所述散射值之间的映射关 系查找表；以及
将被扫描物体放置于多光源CT成像系统中，实际测量至少两个不同 位置的光源对应探测器同一像素点的投影值，并基于所述查找表根据测量 得到的投影值得到该投影值对应的散射估计值。
在本公开的一实施例中，所述经过能谱校正后散射值的分布函数基于 能谱硬化对于散射分布的影响确定，包括：在能谱硬化对于散射分布具有 影响的情况下，则假设经过能谱校正后散射值的分布函数具有变化，将实 际测量投影值表示为透射投影值减去散射因素引起的偏差的情况下，将能 谱校正后实际测量投影值采用一阶泰勒展开，建立能谱校正后实际测量投 影值与能谱校正后透射投影值间的方程，将能谱校正后实际测量总光子数 的表达式与所述方程联立进而得到能谱校正后的散射分布函数。
在本公开的一实施例中，所述经过能谱校正后散射值的分布函数基于 能谱硬化对于散射分布的影响确定，包括：在能谱硬化对于散射分布无影 响的情况下，则假设经过能谱校正后散射值的分布函数无变化，即得到如 下关系：能谱校正后散射值的分布函数为初始散射值。
在本公开的一实施例中，所述根据能谱校正后散射值的表达式进行迭 代求解得到散射估计值，包括：在每个迭代循环内，基于能谱校正后散射 值的表达式获取上一步得到的累积散射，对于第一步迭代，上一步得到的 累积散射为散射值初值，每一步迭代利用上一步得到的累积散射计算散射 去除后的投影值，并通过该投影值计算当前步残留的散射，最后更新当前 步的累积散射，直到到达给定的迭代次数或者当前步计算的残留散射小于 给定值时停止迭代。
在本公开的一实施例中，所述对至少两个不同位置的光源在探测器同 一像素点处的投影数据进行能谱校正的步骤中，基于如下基本假设：处于 不同位置的光源分别发光时，对于同一个探测器像素而言，射线经过被扫 描物体的路径长度相等，同时处于不同位置的光源对应的散射光子数相 等。
在本公开的一实施例中，该CT散射校正方法还包括：利用求解得到 的散射估计值对投影数据进行校正，得到经过散射校正后的投影值。
在本公开的一实施例中，所述将被扫描物体放于多光源CT成像系统 中，获取至少两个不同位置的光源在探测器同一像素点处的投影数据，包 括：将被扫描物体放于多光源CT成像系统中，不同位置光源发射的X射 线经过非均匀滤波栅格进行能谱调制，能谱调制后的X射线经过被扫描物 体吸收之后由探测器进行探测，分别获取至少两个不同位置光源对应在探 测器同一像素点处的投影数据。
根据本公开的另一个方面，提供了一种CT散射校正系统，用于执 行本公开提及的任一种CT散射校正方法，所述CT散射校正系统包括： 多光源CT成像系统，用于放置被扫描物体，以获取至少两个不同位置的 光源在探测器同一像素点处的投影数据；能谱校正模块，用于对至少两个 不同位置的光源在探测器同一像素点处的投影数据进行能谱校正，获得能 谱校正后探测器像素探测得到的总光子数的表达式，该表达式中包括：经 过能谱校正后散射值的分布函数，其中所述经过能谱校正后散射值的分布 函数基于能谱硬化对于散射分布的影响确定；以及散射估计值求解模块， 用于根据能谱校正后散射值的分布函数求解出能谱校正后散射值的表达 式，并根据能谱校正后散射值的表达式进行迭代求解得到散射估计值；和 /或，所述CT散射校正系统包括：存储模块，用于存储预先建立的投影数据值与散射值之间的映射关系查找表，映射关系中的所述散射值为经过能 谱校正后的散射值；多光源CT成像系统，用于放置被扫描物体，以获取 至少两个不同位置的光源对应探测器同一像素点的投影值；以及散射估计 值计算模块，用于基于所述查找表根据测量得到的投影值得到该投影值对 应的散射估计值。
在本公开的一实施例中，所述CT散射校正系统还包括：投影数据校 正模块，用于利用求解得到的散射估计值对投影数据进行校正，得到经过 散射校正后的投影值。
在本公开的一实施例中，所述多光源CT成像系统包括：多个光源、 非均匀滤波栅格以及探测器，将被扫描物体放置于非均匀滤波栅格与探测 器之间；其中，所述多个光源为分布式射线源，该分布式射线源中各个射 线源按照设定依次闪烁；或者所述多个光源由单个射线源的焦点快速移动 形成。
(三)有益效果
从上述技术方案可以看出，本公开提供的CT散射校正方法及系统， 具有以下有益效果：
提出两种方式对散射值进行估计，一种方式是基于预先建立的投影数 据值与散射值映射关系查找表和测量的投影值得到散射估计值，映射关系 中的所述散射值为经过能谱校正后的散射值；另一种方式是基于能谱物理 对散射值的分布函数进行校正后迭代求解得到散射估计值，通过获取至少 两个不同位置的光源在探测器同一像素点处的投影数据，对至少两个不同 位置的光源在探测器同一像素点处的投影数据进行能谱校正，校正后的表 达式中涵盖了能谱硬化对于散射分布的影响，根据能谱硬化对于散射分布 的影响可以确定经过能谱校正后散射值的分布函数，进而根据能谱校正后 散射值的分布函数求出能谱校正后散射值的表达式，对能谱校正后散射值 的表达式进行迭代求解可以得到散射估计值。两种方法均通过在散射估计 过程中加入能谱校正的因素，给出了能谱硬化对于散射分布影响的计算表 达形式或者测量得到方式，基于不同的能谱硬化对于散射分布影响的假 设，按照迭代求解的方式或者基于散射遮挡硬件测量的方式可以得到能谱 校正后的散射估计值，对于硬件的要求不高，同时还具有迭代效率高，计 算准确度较高的优点。
附图说明
图1为现有技术中常用CT系统结构及散射分布示意图。
图2为根据本公开一实施例所示的多光源CT成像系统的结构示意图。
图3为根据本公开一实施例所示的CT散射校正方法流程图。
图4为根据本公开一实施例所示的CT散射校正方法还包括预先获取 能谱校正函数的步骤的流程图。
图5为根据本公开一实施例所示的CT散射校正系统的结构框图。
【符号说明】
2-CT校正系统；
21-多光源CT成像系统；
22a-能谱校正模块；
22b-存储模块；
23-散射估计值求解模块；
24-投影数据校正模块。
具体实施方式
图1为现有技术中常用CT系统结构及散射分布示意图，参照图1所 示，进行CT成像时，由射线源发出X射线，X射线经过被扫描物体的衰 减之后被探测器探测得到投影数据。X射线在经过被扫描物体时有部分光 子发生透射，部分光子发生散射，散射分布相比于通常穿过物体进入探测 器的投影数据分布，变化比较缓慢，因此以高频主束分布来表示透射信息， 以低频散射分布来表示散射信息，由探测器测试得到的信息中包括高频主 束分布和低频散射分布。其中，发生透射的透射光子为有用信息，透射光 子连同散射光子一同穿过待测物体进入探测器，散热光子会造成图像伪影 和CT值不准确等问题。
目前X射线作为CT成像的射线源，由X射线管发射的X射线是分 布在频谱上的量，并非严格的单值。当X射线通过吸收物质时，X射线频 谱的不同部分被吸收(或称为衰减)的量并不相同，较低能量的X射线优 先被吸收，对应X射线频谱分布的峰值会向高能方向移动，这一过程称之 为能谱硬化。射线束硬化使得X射线对应的光子吸收不均衡，这种非线性衰减不作修正，会产生环状或条状伪影，由此导致对CT成像的质量造成 很大影响。
本公开提供了一种CT散射校正方法及系统，考虑到能谱硬化和X射 线散射引起的图像伪影和CT值不准确的问题，将被扫描物体放于多光源 CT成像系统中，不同位置光源发射的X射线经过非均匀滤波栅格进行能 谱调制，能谱调制后的X射线经过被扫描物体吸收之后由探测器进行探 测，一种途径为采用查找表的方法，预先建立投影数据值与散射值之间的 映射关系查找表；以及将被扫描物体放置于多光源CT成像系统中，实际 测量至少两个不同位置的光源对应探测器同一像素点的投影值，并基于所 述查找表根据测量得到的投影值得到该投影值对应的散射估计值。预先建 立投影数据值与散射值之间的映射关系查找表的方法包括：将被扫描物体 放于多光源CT成像系统中，获取至少两个不同位置的光源在不设置散射 遮挡硬件时，不放置和放置被扫描物体对应在探测器同一像素点处的投影 数据值，不放置被扫描物体对应的投影数据值为光源发光时穿过滤波栅格 的光子数，放置被扫描物体对应的投影数据值为探测器像素探测得到的总 光子数；获取在设置散射遮挡硬件时放置被扫描物体对应在探测器同一像 素点的透射光子数，根据不设置散射遮挡硬件时光源发光时穿过滤波栅格 的光子数、不设置散射遮挡硬件时探测器像素探测得到的总光子数和设置 散射遮挡硬件时同一探测器像素探测得到的透射光子数得到经过能谱校 正后的探测器像素探测得到的总光子数和经过能谱校正后的探测器像素 探测到的透射光子数，利用经过能谱校正后的探测器像素探测得到的总光 子数减去经过能谱校正后的探测器像素探测得到的透射光子数求得散射 值，该散射值为经过能谱校正后的散射值；建立探测器像素探测得到的总 光子数和所述散射值之间的映射关系查找表。
另外一种途径为采用基于能谱物理的散射估计方法，将被扫描物体放 于多光源CT成像系统中，获取至少两个不同位置的光源在探测器同一像 素点处的投影数据；对至少两个不同位置的光源在探测器同一像素点处的 投影数据进行能谱校正，获得能谱校正后探测器像素探测得到的总光子数 的表达式，该表达式中包括：经过能谱校正后散射值的分布函数，其中所 述经过能谱校正后散射值的分布函数基于能谱硬化对于散射分布的影响确定；以及根据能谱校正后散射值的分布函数求解出能谱校正后散射值的 表达式，并根据能谱校正后散射值的表达式进行迭代求解得到散射估计 值。上述CT散射校正方法中，不局限于两个不同位置光源对应的情况， 还可以拓展至N个不同位置光源的情况进行散射校正。此种散射估计方法 对应的迭代次数较少，迭代效率很高。当然，在下面的实施例中，对于能 谱校正后散射值的分布函数仅提供了两种示例性的实例，本公开的能谱校 正后散射值的分布函数Gi(Is)的表达式不局限于列举的实施例，在本公开 的其他实施方式中，可以对能谱校正后散射值的分布函数表示为其他表达 式，进行散射值估计的过程与示例的方式一致，仅需要更换能谱校正后散 射值的分布函数Gi(Is)以及散射值IS的表达式即可，这里不再穷举。
为使本公开的目的、技术方案和优点更加清楚明白，以下结合具体实 施例，并参照附图，对本公开进一步详细说明。
第一实施例
在本公开的第一个示例性实施例中，提供了一种CT散射校正方法。
图2为根据本公开一实施例所示的多光源CT成像系统的结构示意图。 图3为根据本公开一实施例所示的CT散射校正方法流程图。
本实施例中，基于多光源CT成像系统来实施上述CT散射校正方法。 参照图2所示，该多光源CT成像系统包括：多个光源、非均匀滤波栅格 以及探测器。在CT成像过程中，X射线源与探测器围绕物体进行相对旋 转运动，以此获得不同旋转角度下的CT数据，即为CT投影值。不同位 置光源发射的X射线经过非均匀滤波栅格进行能谱调制，能谱调制后的X 射线经过被扫描物体吸收之后由探测器进行探测。
本实施例中，将被扫描物体放于多光源CT成像系统中，例如图2中 所示意的将被扫描物体放置于非均匀滤波栅格与平板探测器之间。多光源 扫描的一种实现方式为：多个光源为分布式射线源，该分布式射线源中各 个射线源依次闪烁，射线源可为X光管、碳纳米管或加速器，光源分布形 式及分布间距可以根据需求进行定义。多光源扫描的另一种实现方式为： 单个射线源的焦点快速移动，以此改变射线源发射射线的初始位置，焦点 可以根据需求自定义移动方向。
非均匀滤波栅格为非均匀排布的滤波栅格，滤波栅格由能够吸收部分 射线的材料制造，能够对入射的射线能谱进行调制，根据与射线源焦点的 相对位置产生不同能谱的射线，实现对于入射射线的能谱调制。滤波栅格 的材料包括但不限于为如下材料的一种或几种：铜、银、锡以及钼等。
参照图3所示，本实施例的CT散射校正方法，包括：步骤S111、步 骤S112、步骤S113和步骤S114。
在步骤S111，将被扫描物体放于多光源CT成像系统中，获取至少两 个不同位置的光源在探测器同一像素点处的投影数据。下面以两个不同位 置的光源对应的情况作为示例，需要说明的是，不同位置的光源的个数可 以拓展至多个(大于2个)。
参照图2所示，将被扫描物体放于多光源CT成像系统中，不同位置 光源发射的X射线经过非均匀滤波栅格进行能谱调制，能谱调制后的X 射线经过被扫描物体吸收之后由探测器进行探测，分别获取至少两个不同 位置光源对应在探测器同一像素点处的投影数据。
本实施例中，在同一个探测器像素点可以获取至少两个不同位置的光 源下的X光信号，该探测器像素获取的X光信号即为投影数据。该投影 数据中包含透射光子的贡献和散射光子的贡献，在后续的步骤中需要去除 散射光子的部分，以校正投影数据。
参照图2所示，该多光源CT成像系统中，示意了相邻的光源Pn和光 源Pn+1经过非均匀滤波栅格和被扫描物体后入射到探测器的路径，n为正 整数，其他任意两个不同位置的光源Pa和Pb的情况与上述光源Pn和光源 Pn+1的情况类似，这里仅以光源Pn和光源Pn+1这两个光源作为示例。当两 个任意位置的光源分别闪烁时，例如光源P1和光源P2分别闪烁，对于同一个探测器像素，采集得到的X光信号包括透射光子的贡献以及散射光子 的贡献，可以简单地表示为如下表达式：
Figure BDA0002395479180000101
其中，
Figure BDA0002395479180000102
表示当光源Pi闪烁时，探测器像素探测得到的总光子数，其 中i＝1，2；
Figure BDA0002395479180000103
表示当光源Pi闪烁时对应的散射光子数；
Figure BDA0002395479180000104
表示当光源Pi闪烁 时对应的透射光子数。
透射光子数可以进一步表示为如下表达式：
Figure BDA0002395479180000105
其中，Si(E)表示光源Pi照向探测器的入射能谱，um(E)表示滤波栅格 在能量E下的线性衰减系数；
Figure RE-GDA0002581429730000106
表示射线路径上滤波栅格的厚度；uo(E)表 示被扫描物体在能量E下的线性衰减系数；Li表示射线穿过被扫描物体的 路径长度。考虑到CT成像系统中不同方向上的出射能谱差异性以及不同 探测器之间的不一致性，作为优选，公式(2)中的入射能谱可以是等效 入射能谱，通过对每个探测器像素探测得到的初始能谱进行优化而得到该 等效入射能谱。另外，光源在不同位置下射线所经过滤波栅格的厚度可以 通过事先标定和优化得到。
上述公式(2)中的透射光子数的表达式中包含滤波栅格和被扫描物 体的信息，后续步骤旨在将散射光子数进行能谱校正，从而对应得到去除 了散射光子贡献的投影数据，即得到相对准确的透射光子数，从而可以根 据公式(2)来进行CT成像的过程。
术语“入射能谱”是指相对于探测器为入射的方向，由探测器探测到 的入射能谱；相对于光源来说，从光源发出。
在步骤S112，对至少两个不同位置的光源在探测器同一像素点处的投 影数据进行能谱校正，获得能谱校正后探测器像素探测得到的总光子数的 表达式，该表达式中包括：经过能谱校正后散射值的分布函数。
本实施例中，进行能谱校正后的表达式将探测器像素探测得到的总光 子数近似表示为两部分的加和，一部分为光源发光时穿过滤波栅格的光子 数与能谱校正因子的乘积；另一部分为经过能谱校正后散射值的分布函 数。
对于同一个探测器像素而言，当来自不同位置光源的经过扫描物体的 射线路径在空间上比较接近时，散射基本相同：
L1＝L2＝L， (3)
IS1＝IS2＝IS， (4)
其中，L为射线穿过被扫描物体的路径长度值；IS为散射值。
以上可以通过合理配置不同光源的间距实现。此时可以认为当光源P1和光源P2分别发光时，对于同一个探测器像素而言，射线经过被扫描物体 的路径长度基本相等，不同光源闪烁时对应的散射光子数也相等。
将公式(4)带入至公式(1)中，得到如下表达式：
Figure BDA0002395479180000111
考虑到能谱硬化对于散射分布的影响，利用能谱校正的方法对探测得 到的光子数进行能谱校正。能谱校正后探测器像素探测得到的总光子数的 表达式为：能谱校正后探测器像素探测得到的总光子数为能谱校正后的透 射光子数与能谱校正后散射值分布函数的叠加。
在一实施例中，对公式(5)进行能谱校正之后，得到如下表达式：
Figure BDA0002395479180000112
其中，
Figure BDA0002395479180000113
表示当光源Pi闪烁时，经过能谱校正后的探测器像素探测得 到的总光子数，i＝1，2；
Figure BDA0002395479180000114
表示光源Pi发光时穿过滤波栅格的光子数；fi(·) 为光源Pi闪烁时对应的能谱校正函数，可以由预先模拟或实验标定得到；
Figure BDA0002395479180000115
表示光源Pi闪烁时实际测算得到的投影值；
Figure BDA0002395479180000116
表示光源Pi闪烁时对应的 透射投影值；Gi(Is)表示光源Pi闪烁时经过能谱校正之后关于散射值IS的分 布函数。上述约等号可以视为等号，表示进行能谱校正的过程。
图4为根据本公开一实施例所示的CT散射校正方法还包括预先获取 能谱校正函数的步骤的流程图。
需要说明的是，上述能谱校正函数fi(·)可以由预先模拟或实验标定得 到，例如参照图4所示，可以在步骤S101中，基于入射能谱(优选等效 入射能谱)和滤波栅格厚度预先得到能谱校正函数。
其中，
Figure BDA0002395479180000117
和
Figure BDA0002395479180000118
分别满足如下表达式：
Figure BDA0002395479180000121
Figure BDA0002395479180000122
Figure BDA0002395479180000123
在步骤S113，对经过能谱校正后散射值的分布函数的变化情况进行假 设，根据对应的假设推导得出能谱校正后散射值的表达式，并根据能谱校 正后散射值的表达式进行迭代算法求解得出散射估计值。该步骤S113中， 由于散射值的分布函数已经经过了能谱校正，则根据对应假设推导得出的 散射值已经是经过能谱校正后的散射值，由该散射值最终迭代求解得到的 散射估计值已经涵盖了能谱硬化的影响。
基于Gi(Is)取值的不同，可以进行不同的迭代求解过程。
下面采用两个实施例来示例性说明两种迭代方法对散射值进行估计 的过程。其他形式的Gi(Is)的表达式对应的求解过程与下面的示例类似， 只需要代入相应Gi(Is)和Is的表达式中即可。该步骤S113中，对经过能谱 校正后散射值的分布函数的变化情况进行假设，一种假设是能谱硬化对于 散射分布没有影响或者影响很小，即经过能谱校正后散射值的分布函数不 发生变化。另一种假设是能谱硬化对于散射分布具有影响，即经过能谱校正后散射值的分布函数发生变化。
下面针对不同的假设来详细介绍如何根据对应的假设推导得出散射 值的表达式，并根据散射值的表达式对散射值进行迭代算法求解得出散射 值估计的过程。根据前述能谱硬化对于散射分布的影响得出的能谱校正后 散射值的分布函数Gi(Is)，可以求解得出由投影数据计算散射的公式，即 得到散射值的表达式，该散射值为经过能谱校正后的散射值，然后基于散 射值的表达式可以通过迭代算法求解得出散射估计值。
第一种迭代方法示例如下，所述经过能谱校正后散射值的分布函数基 于能谱硬化对于散射分布的影响确定，包括：在能谱硬化对于散射分布无 影响的情况下，则假设经过能谱校正后散射值的分布函数无变化，即得到 如下关系：能谱校正后散射值的分布函数不发生变化，该能谱校正后散射 值的分布函数表示为初始散射值，这里的初始散射值含义为并未发生变化 的散射值Is，在公式(10)中以Is表示初始散射值。
该实施例中，可以简单地认为经过能谱校正之后散射的分布没有改 变，即存在如下关系式：
Gi(Is)＝Is， (10)
则将公式(10)代入至公式(6)中，得到如下表达式：
Figure BDA0002395479180000131
对于同一个探测器像素而言，两个不同位置的光源出射的射线经过扫 描物体的路径长度基本一致，那么经过能谱校正函数会映射到相同的投影 值，即存在如下关系式：
Figure BDA0002395479180000132
那么根据公式(6)、(11)和(12)可以求解得到能谱校正后Is的表达 式如下：
Figure BDA0002395479180000133
针对公式(13)所示的能谱校正后Is的表达式，采用迭代算法进行求 解。
在本公开的一实施例中，根据能谱校正后散射值的表达式进行迭代求 解得到散射估计值，包括：在每个迭代循环内，基于能谱校正后散射值的 表达式获取上一步得到的累积散射，对于第一步迭代，上一步得到的累积 散射为散射值初值，每一步迭代利用上一步得到的累积散射计算散射去除 后的投影值，并通过该投影值计算当前步残留的散射，最后更新当前步的 累积散射，直到到达给定的迭代次数或者当前步计算的残留散射小于给定 值时停止迭代。
下面以具体迭代的实例进行说明。在一实施例中，根据散射值的表达 式对散射值进行迭代算法求解得出散射估计值，包括：
给散射值赋予初值Is 0，初始化设定停止准则stop_c以及最大迭代步数 Niter；迭代步数从1开始，每次循环增加1步，从1增加至Niter，即满足： 迭代步数iter＝1：Niter，这里的等号表示赋值；
在每个迭代循环内，将特定迭代步数下的散射值代入至
Figure BDA0002395479180000141
的表达式中 得到该特定迭代步数下
Figure BDA0002395479180000142
的值，将
Figure BDA0002395479180000143
的值代入至所述散射值的表达式中 得到特定迭代步数下的散射值残差ΔIs iter；则根据得到的特定迭代步数下的 散射值残差ΔIs iter与上一步迭代得到的散射值
Figure BDA0002395479180000144
计算得到本次迭代循环 的散射值
Figure BDA0002395479180000145
Figure BDA0002395479180000146
在特定迭代步数下的表达式如下：
Figure BDA0002395479180000147
本次迭代循环的散射值
Figure BDA0002395479180000148
满足：
Figure BDA0002395479180000149
如果在某步迭代后，|ΔIs|小于停止准则stop_c，则停止迭代，得到最 终的散射估计值为
Figure BDA00023954791800001410
否则，一直迭代至最大迭代步数，返回最终的散 射估计值为
Figure BDA00023954791800001411
下面采用伪代码的形式进行示例性说明。在一实例中，迭代算法求解 散射值Is的伪代码(Pseudocode)如下：
(a1)初始化Is 0、停止准则stop_c以及最大迭代步数Niter；
(b1)for iter＝1：Niter：
(c1)
Figure BDA00023954791800001412
(d1)
Figure BDA00023954791800001413
(e1)
Figure BDA00023954791800001414
(f1)
Figure BDA00023954791800001415
(g1)如果|ΔIs|小于stop_c，则停止迭代，最后的散射估计值 为
Figure BDA00023954791800001416
(h1)end for.返回最后的散射估计值为
Figure BDA00023954791800001417
第二种迭代方法示例如下，在能谱硬化对于散射分布具有影响的情况 下，则假设经过能谱校正后散射值的分布函数具有变化，将实际测量投影 值表示为透射投影值减去散射因素引起的偏差的情况下，将能谱校正后实 际测量投影值采用一阶泰勒展开，建立能谱校正后实际测量投影值与能谱 校正后透射投影值间的方程，将能谱校正后实际测量总光子数的表达式与 所述方程联立进而得到能谱校正后的散射分布函数。
具体而言，该实施例中，考虑到能谱硬化对于散射分布的影响，并对 散射进行相应的补偿。透射投影值
Figure BDA0002395479180000151
和实际测量投影值
Figure BDA0002395479180000152
之间存在散射 因素引起的偏差Δp，利用公式表示如下：
Figure BDA0002395479180000153
Figure BDA0002395479180000154
对于实际测量投影值
Figure BDA0002395479180000155
进行能谱校正，并应用一阶泰勒展开进行近 似，能谱校正后实际测量投影值采用一阶泰勒展开的表达式如下：
Figure BDA0002395479180000156
其中，f′表示f的一阶导数。
在上述(16)表达式的基础上，可以建立能谱校正后实际测量投影值 与能谱校正后透射投影值间的方程。具体建立过程为：根据公式(16)可 以进行推导，得到能谱校正之后的探测器像素探测得到的总光子数
Figure BDA0002395479180000159
的表 达式如下：
Figure BDA0002395479180000157
基于公式(17)中的泰勒近似，对应公式(6)可以得到如下表达式：
Figure BDA0002395479180000158
方程(17)示例性描述了本实施例中能谱校正后实际测量投影值与能 谱校正后透射投影值间的方程，其中本领域技术人员公知，实际测量投影 值与实际测量时探测器像素探测得到的总光子数之间存在负对数关系。将 能谱校正后实际测量总光子数对应的表达式(6)与所述方程(17)联立 进而得到能谱校正后散射值的分布函数(散射值的分布函数也可以称为散 射分布函数)，即得到公式(18)。
然后基于能谱校正后散射值的分布函数求解出散射值的表达式。在本 实施例中，可以将不同像素点对应的能谱校正后散射值的分布函数代入能 谱校正后实际测量时探测器像素探测得到的总光子数的表达式中，从而求 解得出能谱校正后散射值的表达式。具体而言，将公式(18)代入公式(6) 中，可以得到如下表达式：
Figure BDA0002395479180000161
根据公式(19)列出的方程可以解出Is的表达式。下面列出求解出的 能谱校正后Is的表达式如下：
Figure BDA0002395479180000162
针对公式(20)所示的能谱校正后Is的表达式，采用迭代算法进行求 解。
根据能谱校正后散射值的表达式进行迭代求解得到散射估计值，包 括：在每个迭代循环内，基于能谱校正后散射值的表达式获取上一步得到 的累积散射，对于第一步迭代，上一步得到的累积散射为散射值初值，每 一步迭代利用上一步得到的累积散射计算散射去除后的投影值，并通过该 投影值计算当前步残留的散射，最后更新当前步的累积散射，直到到达给 定的迭代次数或者当前步计算的残留散射小于给定值时停止迭代。
下面以具体实例进行介绍。在一实例中，迭代算法求解散射值Is的伪 代码(Pseudocode)如下：
(a2)初始化Is 0、停止准则stop_c以及最大迭代步数Niter；
(b2)for iter＝1：Niter：
(c2)
Figure BDA0002395479180000163
(d2)
Figure BDA0002395479180000164
(e2)
Figure BDA0002395479180000165
(f2)
Figure BDA0002395479180000171
(g2)如果|ΔIs|小于stop_c，则停止迭代，最后的散射估计值 为
Figure BDA0002395479180000172
(h2)end for.返回最后的散射估计值为
Figure BDA0002395479180000173
上述实施例中，第二种迭代方法具有更高的迭代效率，可以很快逼近 散射估计值。
在步骤S114，利用求解得到的散射估计值对投影数据进行校正，得到 经过散射校正后的投影值。
本实施例中，步骤S114中的投影数据中包含透射光子的贡献和散射 光子的贡献，经过散射校正后的投影值仅包含透射光子的贡献。其中利用 求解得到的散射估计值对投影数据校正，包括：在投影数据中减去散射估 计值，以实现对于投影数据的校正。
当然，上述实施例中，以两个不同位置的光源在探测器同一像素点处 的投影数据列出两个方程进行示例，本公开中，i表示不同位置光源的序 号，个数为至少两个，这里仅以i＝1，2作为示例，实际上，对于个数大于 两个的情况，i＝1，2，…，k，可以按照公式(6)中的形式进行类推得到包含k 个方程的方程组，k＞2；
Figure BDA0002395479180000174
在不同位置的光源的个数大于2的情况下，不仅可以按照上述方法求 解得到散射估计值，还可以得出其他参数。
第二实施例
本公开的第二个示例性实施例中，提供了一种CT散射校正方法，本 实施例的CT散射校正方法与第一实施例不同，是利用查找表法来对散射 进行校正。
本实施例中，基于预先建立的投影数据值与散射值映射关系查找表和 测量的投影值得到散射估计值，包括：预先建立投影数据值与散射值之间 的映射关系查找表；以及实际测量被扫描物体放置于多光源CT成像系统 中，至少两个不同位置的光源对应探测器同一像素点的投影值，并基于所 述查找表根据测量得到的投影值得到该投影值对应的散射估计值。
在一实例中，预先建立投影数据值与散射值之间的映射关系查找表， 包括：将被扫描物体放于多光源CT成像系统中，获取至少两个不同位置 的光源在不设置散射遮挡硬件时，不放置和放置被扫描物体对应在探测器 同一像素点处的投影数据值。其中不放置被扫描物体对应的投影数据值为 光源Pi发光时穿过滤波栅格的光子数
Figure BDA0002395479180000181
放置被扫描物体对应的投影数据 值为探测器像素探测得到的总光子数
Figure BDA0002395479180000182
i表示不同位置光源的序号，个数 为至少两个；获取在设置散射遮挡硬件时放置被扫描物体对应在探测器同 一像素点的透射光子数
Figure BDA0002395479180000183
根据
Figure BDA0002395479180000184
和
Figure BDA0002395479180000185
得到经过能谱校正后的探测 器像素探测得到的总光子数
Figure BDA0002395479180000186
和经过能谱校正后的探测器像素探测得到 的透射光子数
Figure BDA0002395479180000187
可以参照前述公式(6)的求得。利用经过能谱校正后 的探测器像素探测得到的总光子数
Figure BDA0002395479180000188
减去经过能谱校正后的探测器像素 探测得到的透射光子数
Figure BDA0002395479180000189
求得散射值，该散射值为经过能谱校正后的散射 值。基于上述，可以建立探测器像素探测得到的总光子数
Figure BDA00023954791800001810
和能谱校正后 散射值之间的映射关系查找表。也就是说，建立了不设置散射遮挡硬件、 放置被扫描物体时对应的实测投影数据值和能谱校正后散射值之间的查 找表，直接可以利用实际测量的投影值和查找表获取散射估计值，该散射 估计值为能谱校正后的散射估计值，无需重复进行测量或者其他复杂计算 过程。
然后，将被扫描物体放置于多光源CT成像系统中，实际测量至少两 个不同位置的光源对应探测器同一像素点的投影值，并基于所述查找表根 据测量得到的投影值得到该投影值对应的散射估计值。
第三实施例
本公开的第三个示例性实施例中，提供了一种CT散射校正系统，用 于执行本公开提及的任一种CT散射校正方法。
图5为根据本公开一实施例所示的CT散射校正系统的结构框图。
参照图5所示，本实施例的CT散射校正系统2包括：多光源CT成 像系统21、能谱校正模块22a和散射估计值求解模块23。
其中，多光源CT成像系统21用于放置被扫描物体，以获取至少两个 不同位置的光源在探测器同一像素点处的投影数据。能谱校正模块22a用 于对至少两个不同位置的光源在探测器同一像素点处的投影数据进行能 谱校正，获得能谱校正后探测器像素探测得到的总光子数的表达式。该表 达式中包括：经过能谱校正后散射值的分布函数。其中，经过能谱校正后 散射值的分布函数基于能谱硬化对于散射分布的影响确定。散射估计值求 解模块23用于根据能谱校正后散射值的分布函数求解出能谱校正后散射 值的表达式，并根据能谱校正后散射值的表达式进行迭代求解得到散射估 计值。
在另一实施例中，CT散射校正系统2可以用于执行查找表法来实现 散射值的估计校正。该CT散射校正系统2还可以包括：存储模块22b， 该存储模块22b用于存储预先建立的投影数据值与散射值之间的映射关系 查找表，映射关系中的所述散射值为经过能谱校正后的散射值。上述散射 估计值计算模块23还用于基于所述查找表根据测量得到的投影值得到该 投影值对应的散射估计值。当然，该CT散射校正系统可以同时包括：多 光源CT成像系统21、能谱校正模块22a、散射估计值求解模块23和存储 模块22b，或者，可以包括：多光源CT成像系统21、能谱校正模块22a 和散射估计值求解模块23或者包括：多光源CT成像系统21、存储模块 22b和散射估计值求解模块23。
在本公开的一实施例中，参照图5所示，CT散射校正系统2还包括： 投影数据校正模块24。该投影数据校正模块24用于利用求解得到的散射 估计值对投影数据进行校正，得到经过散射校正后的投影值。
在本公开的一实施例中，参照图1所示，多光源CT成像系统21包括： 多个光源、非均匀滤波栅格以及探测器，将被扫描物体放置于非均匀滤波 栅格与探测器之间。本实施例中，多个光源可以为分布式射线源，该分布 式射线源中各个射线源按照设定依次闪烁；或者所述多个光源可以由单个 射线源的焦点快速移动形成。
综上所述，本公开提供了一种CT散射校正方法及系统，提供了两种 途径来进行散射分布的估计，一种是通过查找表的方法，另一种是通过基 于能谱物理的散射估计方法，在基于能谱物理的散射估计方法中，通过获 取至少两个不同位置的光源在探测器同一像素点处的投影数据，对至少两 个不同位置的光源在探测器同一像素点处的投影数据进行能谱校正，校正 后的表达式中涵盖了能谱硬化对于散射的影响，根据能谱硬化对于散射分布的影响可以确定经过能谱校正后散射值的分布函数，进而根据能谱校正 后散射值的分布函数求解出能谱校正后散射值的表达式，并根据能谱校正 后散射值的表达式进行迭代求解可以得到散射估计值，两种方法均通过在 散射估计过程中加入能谱校正的因素，给出了能谱硬化对于散射分布影响 的计算表达形式或者测量得到方式，基于不同的能谱硬化对于散射分布影 响的假设，按照迭代求解的方式或者基于散射遮挡硬件测量的方式可以得 到能谱校正后的散射估计值，对于硬件的要求不高，同时还具有迭代效率 高，计算准确度较高的优点。
本领域那些技术人员可以理解，可以对实施例中的系统中的模块进行 自适应性地改变并且把它们设置在与该实施例不同的一个或多个模块中。 可以把实施例中的模块组合成一个模块，以及此外可以把它们分成多个子 模块。除了这样的特征和/或过程或者单元中的至少一些是相互排斥之外， 可以采用任何组合对本说明书(包括伴随的权利要求、摘要和附图)中公 开的所有特征以及如此公开的任何方法或者装置的所有过程或单元进行 组合。除非另外明确陈述，本说明书(包括伴随的权利要求、摘要和附图) 中公开的每个特征可以由提供相同、等同或相似目的的替代特征来代替。 根据本公开实施例的模块中的任意一个或多个可以至少被部分地实现为 硬件电路，例如现场可编程门阵列(FPGA)、可编程逻辑阵列(PLA)、 片上系统、基板上的系统、封装上的系统、专用集成电路(ASIC)，或可以通过对电路进行集成或封装的任何其他的合理方式的硬件或固件来实 现，或以软件、硬件以及固件三种实现方式中任意一种或以其中任意几种 的适当组合来实现。或者，根据本公开实施例的模块中的一个或多个可以 至少被部分地实现为计算机程序模块，当该计算机程序模块被运行时，可 以执行相应的功能。
再者，单词“包含”或“包括”不排除存在未列在权利要求中的元件 或步骤。位于元件之前的单词“一”或“一个”不排除存在多个这样的元 件。
除非存在技术障碍或矛盾，本发明的上述各种实施方式可以自由组合 以形成另外的实施例，这些另外的实施例均在本发明的保护范围中。
以上所述的具体实施例，对本公开的目的、技术方案和有益效果进行 了进一步详细说明，所应理解的是，以上所述仅为本公开的具体实施例而 已，并不用于限制本公开，凡在本公开的精神和原则之内，所做的任何修 改、等同替换、改进等，均应包含在本公开的保护范围之内。

Claims (10)

1.一种CT散射校正方法，其特征在于，包括：
基于预先建立的投影数据值与散射值映射关系查找表和测量的投影值得到散射估计值，映射关系中的所述散射值为经过能谱校正后的散射值；或者，
基于能谱物理对散射值的分布函数进行校正后迭代求解得到散射估计值；
其中，所述基于能谱物理对散射值的分布函数进行校正后迭代求解得到散射估计值，包括：
将被扫描物体放于多光源CT成像系统中，获取至少两个不同位置的光源在探测器同一像素点处的投影数据；
对至少两个不同位置的光源在探测器同一像素点处的投影数据进行能谱校正，获得能谱校正后探测器像素探测得到的总光子数的表达式，该表达式中包括：经过能谱校正后散射值的分布函数，其中所述经过能谱校正后散射值的分布函数基于能谱硬化对于散射分布的影响确定；以及
根据能谱校正后散射值的分布函数求解出能谱校正后散射值的表达式，并根据能谱校正后散射值的表达式进行迭代求解得到散射估计值。
2.根据权利要求1所述的CT散射校正方法，其特征在于，所述能谱校正后探测器像素探测得到的总光子数的表达式为：能谱校正后探测器像素探测得到的总光子数为能谱校正后的透射光子数与能谱校正后散射值分布函数的叠加。
3.根据权利要求1或2所述的CT散射校正方法，其特征在于，所述基于预先建立的投影数据值与散射值映射关系查找表和测量的投影值得到散射估计值，包括：
预先建立投影数据值与散射值之间的映射关系查找表：将被扫描物体放于多光源CT成像系统中，获取至少两个不同位置的光源在不设置散射遮挡硬件时，不放置和放置被扫描物体对应在探测器同一像素点处的投影数据值，不放置被扫描物体对应的投影数据值为光源发光时穿过滤波栅格的光子数，放置被扫描物体对应的投影数据值为探测器像素探测得到的总光子数；获取在设置散射遮挡硬件时放置被扫描物体对应在探测器同一像素点探测得到的透射光子数；根据不设置散射遮挡硬件时光源发光时穿过滤波栅格的光子数、不设置散射遮挡硬件时探测器像素探测得到的总光子数和设置散射遮挡硬件时同一探测器像素探测得到的透射光子数得到经过能谱校正后的探测器像素探测得到的总光子数和透射光子数；利用能谱校正后的探测器像素探测得到的总光子数减去能谱校正后的探测器像素探测得到的透射光子数求得散射值，该散射值为经过能谱校正后的散射值；建立探测器像素探测得到的总光子数和所述散射值之间的映射关系查找表；以及
将被扫描物体放置于多光源CT成像系统中，实际测量至少两个不同位置的光源对应探测器同一像素点的投影值，并基于所述查找表根据测量得到的投影值得到该投影值对应的散射估计值。
4.根据权利要求1或2所述的CT校正方法，其特征在于，
所述经过能谱校正后散射值的分布函数基于能谱硬化对于散射分布的影响确定，包括：在能谱硬化对于散射分布具有影响的情况下，则假设经过能谱校正后散射值的分布函数具有变化，将实际测量投影值表示为透射投影值减去散射因素引起的偏差的情况下，将能谱校正后实际测量投影值采用一阶泰勒展开，建立能谱校正后实际测量投影值与能谱校正后透射投影值间的方程，将能谱校正后实际测量总光子数的表达式与所述方程联立进而得到能谱校正后的散射分布函数。
5.根据权利要求1或2所述的CT散射校正方法，其特征在于，所述经过能谱校正后散射值的分布函数基于能谱硬化对于散射分布的影响确定，包括：在能谱硬化对于散射分布无影响的情况下，则假设经过能谱校正后散射值的分布函数无变化，即得到如下关系：能谱校正后散射值的分布函数为初始散射值。
6.根据权利要求1或2所述的CT散射校正方法，其特征在于，所述根据能谱校正后散射值的表达式进行迭代求解得到散射估计值，包括：在每个迭代循环内，基于能谱校正后散射值的表达式获取上一步得到的累积散射，对于第一步迭代，上一步得到的累积散射为散射值初值，每一步迭代利用上一步得到的累积散射计算散射去除后的投影值，并通过该投影值计算当前步残留的散射，最后更新当前步的累积散射，直到到达给定的迭代次数或者当前步计算的残留散射小于给定值时停止迭代。
7.根据权利要求1所述的CT散射校正方法，其特征在于，还包括：
利用求解得到的散射估计值对投影数据进行校正，得到经过散射校正后的投影值。
8.根据权利要求1所述的CT散射校正方法，其特征在于，所述将被扫描物体放于多光源CT成像系统中，获取至少两个不同位置的光源在探测器同一像素点处的投影数据，包括：
将被扫描物体放于多光源CT成像系统中，不同位置光源发射的X射线经过非均匀滤波栅格进行能谱调制，能谱调制后的X射线经过被扫描物体吸收之后由探测器进行探测，分别获取至少两个不同位置光源对应在探测器同一像素点处的投影数据。
9.一种CT散射校正系统，其特征在于，用于执行权利要求1-8中任一项所述的CT散射校正方法，所述CT散射校正系统包括：
多光源CT成像系统，用于放置被扫描物体，以获取至少两个不同位置的光源在探测器同一像素点处的投影数据；
能谱校正模块，用于对至少两个不同位置的光源在探测器同一像素点处的投影数据进行能谱校正，获得能谱校正后探测器像素探测得到的总光子数的表达式，该表达式中包括：经过能谱校正后散射值的分布函数，其中所述经过能谱校正后散射值的分布函数基于能谱硬化对于散射分布的影响确定；以及
散射估计值求解模块，用于根据能谱校正后散射值的分布函数求解出能谱校正后散射值的表达式，并根据能谱校正后散射值的表达式进行迭代求解得到散射估计值；
和/或，所述CT散射校正系统包括：
存储模块，用于存储预先建立的投影数据值与散射值之间的映射关系查找表，映射关系中的所述散射值为经过能谱校正后的散射值；
多光源CT成像系统，用于放置被扫描物体，以获取至少两个不同位置的光源对应探测器同一像素点的投影值；以及
散射估计值计算模块，用于基于所述查找表根据测量得到的投影值得到该投影值对应的散射估计值。
10.根据权利要求9所述的CT散射校正系统，其特征在于，
所述多光源CT成像系统包括：多个光源、非均匀滤波栅格以及探测器，将被扫描物体放置于非均匀滤波栅格与探测器之间；其中，所述多个光源为分布式射线源，该分布式射线源中各个射线源按照设定依次闪烁；或者所述多个光源由单个射线源的焦点快速移动形成。
CN202010133896.6A 2020-02-28 2020-02-28 Ct散射校正方法及系统 Active CN111643104B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010133896.6A CN111643104B (zh) 2020-02-28 2020-02-28 Ct散射校正方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010133896.6A CN111643104B (zh) 2020-02-28 2020-02-28 Ct散射校正方法及系统

Publications (2)

Publication Number Publication Date
CN111643104A true CN111643104A (zh) 2020-09-11
CN111643104B CN111643104B (zh) 2022-03-18

Family

ID=72346270

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010133896.6A Active CN111643104B (zh) 2020-02-28 2020-02-28 Ct散射校正方法及系统

Country Status (1)

Country Link
CN (1) CN111643104B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112415569A (zh) * 2020-10-30 2021-02-26 湖北锐世数字医学影像科技有限公司 能量校正方法、装置、设备、系统及计算机可读存储介质
CN113237903A (zh) * 2021-06-15 2021-08-10 清华大学 基于双层平板探测器锥束ct的散射校正方法及装置
CN113409414A (zh) * 2021-06-08 2021-09-17 江苏一影医疗设备有限公司 X线图像的散射修正方法、装置、电子设备、存储介质
CN115270075A (zh) * 2022-08-03 2022-11-01 北京朗视仪器股份有限公司 基于钨靶x射线球管的双能ct系统的能谱数据校正方法
CN116413290A (zh) * 2023-02-22 2023-07-11 奥影检测科技(上海)有限公司 一种工业ct散射校正方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1934590A (zh) * 2004-03-17 2007-03-21 皇家飞利浦电子股份有限公司 用于相干散射ct的射束硬化和衰减校正
CN101226642A (zh) * 2008-01-25 2008-07-23 西安交通大学 基于ct数据一致性的投影射束硬化校正方法
US20180235562A1 (en) * 2017-02-17 2018-08-23 Toshiba Medical Systems Corporation Combined sinogram- and image-domain material decomposition for spectral computed tomography (ct)
CN109363703A (zh) * 2018-10-18 2019-02-22 清华大学 Ct系统能谱不一致性的校正方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN1934590A (zh) * 2004-03-17 2007-03-21 皇家飞利浦电子股份有限公司 用于相干散射ct的射束硬化和衰减校正
CN101226642A (zh) * 2008-01-25 2008-07-23 西安交通大学 基于ct数据一致性的投影射束硬化校正方法
US20180235562A1 (en) * 2017-02-17 2018-08-23 Toshiba Medical Systems Corporation Combined sinogram- and image-domain material decomposition for spectral computed tomography (ct)
CN109363703A (zh) * 2018-10-18 2019-02-22 清华大学 Ct系统能谱不一致性的校正方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
HEWEI GAO等: "Modulator design for x-ray scatter correction using primary modulation material selection", 《MEDICAL PHYSICS》 *
HEWEI GAO等: "Physics-based spectral compensation algorithm for x-ray CT with primary modulator", 《PHYSICS IN MEDICINE AND BIOLOGY》 *
TAO ZHANG等: "An Analysis of Scatter Characteristics in X-ray CT Spectral", 《ARXIV》 *

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112415569A (zh) * 2020-10-30 2021-02-26 湖北锐世数字医学影像科技有限公司 能量校正方法、装置、设备、系统及计算机可读存储介质
CN112415569B (zh) * 2020-10-30 2022-06-17 湖北锐世数字医学影像科技有限公司 能量校正方法、装置、设备、系统及计算机可读存储介质
CN113409414A (zh) * 2021-06-08 2021-09-17 江苏一影医疗设备有限公司 X线图像的散射修正方法、装置、电子设备、存储介质
CN113409414B (zh) * 2021-06-08 2024-03-26 江苏一影医疗设备有限公司 X线图像的散射修正方法、装置、电子设备、存储介质
CN113237903A (zh) * 2021-06-15 2021-08-10 清华大学 基于双层平板探测器锥束ct的散射校正方法及装置
CN113237903B (zh) * 2021-06-15 2022-06-10 清华大学 基于双层平板探测器锥束ct的散射校正方法及装置
CN115270075A (zh) * 2022-08-03 2022-11-01 北京朗视仪器股份有限公司 基于钨靶x射线球管的双能ct系统的能谱数据校正方法
CN116413290A (zh) * 2023-02-22 2023-07-11 奥影检测科技(上海)有限公司 一种工业ct散射校正方法
CN116413290B (zh) * 2023-02-22 2023-11-28 奥影检测科技(上海)有限公司 一种工业ct散射校正方法

Also Published As

Publication number Publication date
CN111643104B (zh) 2022-03-18

Similar Documents

Publication Publication Date Title
CN111643104B (zh) Ct散射校正方法及系统
GB2551029B (en) Method for performing material decoposition using a dual-energy X-Ray CT and corresponding Dual-energy X-Ray CT Apparatus
US11331064B2 (en) Method and apparatus for computed tomography (CT) and material decomposition with pile-up correction calibrated using real pulse pileup effect and detector response
JP7139119B2 (ja) 医用画像生成装置及び医用画像生成方法
US9476993B2 (en) Apparatus and method for computing detector response of a photon-counting detector
US9687207B2 (en) Pre-reconstruction calibration, data correction, and material decomposition method and apparatus for photon-counting spectrally-resolving X-ray detectors and X-ray imaging
US8705827B2 (en) Scatter correction methods
US9020092B2 (en) Apparatus and method for angular response calibration of photon-counting detectors in sparse spectral computed tomography imaging
US8958524B2 (en) Correction of projection data in radiation system
US9285326B2 (en) Sparse and energy discriminating collimated detector elements to assist scatter evaluation in CT imaging
JP6640498B2 (ja) Ｘ線コンピュータ断層撮影装置、画像再構成方法及び画像再構成プログラム
WO2014171539A1 (ja) Ｘ線コンピュータ断層撮影装置及び補正方法
US10679385B1 (en) System and method for statistical iterative reconstruction and material decomposition
JPWO2013089155A1 (ja) Ｘ線ｃｔ装置および散乱ｘ線補正方法
US20160217594A1 (en) Cost-function based method and apparatus for projection-domain basis decomposition in spectral computed tomography
US10398405B2 (en) Gain calibration and correction in radiation system
JPWO2021108715A5 (zh)
JP2009082616A (ja) 多管球ｘ線ｃｔにおける散乱線強度分布のスケーリング方法および多管球ｘ線ｃｔ装置
US9230348B2 (en) Imaging system for imaging a region of interest
US10646186B2 (en) X-ray CT apparatus, information processing device and information processing method
JP2016104125A (ja) Ｘ線ｃｔ装置、画像処理装置およびプログラム
US10497153B2 (en) Heel effect correction in computed tomography
US11013487B2 (en) Method and apparatus for computed tomography (CT) and material decomposition with count-rate dependent pileup correction
CN112562025A (zh) X射线成像中的数据校正
CN111134709B (zh) 一种多能量ct基材料分解方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant