CN109919868B - 一种锥束ct射束硬化曲线侦测及投影加权校正方法 - Google Patents

一种锥束ct射束硬化曲线侦测及投影加权校正方法 Download PDF

Info

Publication number
CN109919868B
CN109919868B CN201910143911.2A CN201910143911A CN109919868B CN 109919868 B CN109919868 B CN 109919868B CN 201910143911 A CN201910143911 A CN 201910143911A CN 109919868 B CN109919868 B CN 109919868B
Authority
CN
China
Prior art keywords
projection
ray
image
hardening
correction
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.)
Active
Application number
CN201910143911.2A
Other languages
English (en)
Other versions
CN109919868A (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.)
Northwestern Polytechnical University
Original Assignee
Northwestern Polytechnical 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 Northwestern Polytechnical University filed Critical Northwestern Polytechnical University
Priority to CN201910143911.2A priority Critical patent/CN109919868B/zh
Publication of CN109919868A publication Critical patent/CN109919868A/zh
Application granted granted Critical
Publication of CN109919868B publication Critical patent/CN109919868B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Apparatus For Radiation Diagnosis (AREA)
  • Analysing Materials By The Use Of Radiation (AREA)
  • Image Processing (AREA)

Abstract

本发明公开了一种锥束CT射束硬化曲线侦测及投影加权校正方法,利用射线采样矩阵信息先验知识和预先重建的切片图像,以射线与二值化图像进行求取交线长度,侦测得到射束硬化曲线。通过对拟合函数进行投影加权修正,实现多能投影的硬化校正,改善因多能投影重建引起的杯状伪影,提高图像质量。本发明提供的锥束CT射束硬化曲线侦测及投影加权校正方法,适用于任意复杂度的被测物体的硬化曲线侦测及杯状伪影校正,方法的可靠性、稳定性、抗噪性好,可在很大程度上减少锥束CT硬化投影杯状伪影对图像的干扰与影响,明显改善锥束CT图像质量。

Description

一种锥束CT射束硬化曲线侦测及投影加权校正方法
技术领域
本发明属于锥束CT应用相关的医学成像和工业无损检测领域,涉及一种锥束CT射束硬化曲线侦测及投影加权校正方法。
背景技术
锥束CT(Cone Beam Computed Tomography,CBCT)作为一种先进的医学成像和工业无损检测技术,在不破坏物体情况下,以二维或三维断层图像的形式清晰、准确、直观地展现被检测物体的内部结构,定量地提供物体内部缺陷位置和尺寸。
实际工业射线检测中使用的射线不是单色射线,而是宽束、连续谱的多色射线。且工业 CT扫描能量较高,由于不同能量射线衰减程度不同,射线穿越物体后能谱发生变化。射束硬化伪影是影响重建图像质量的重要因素。在多色X射线谱中,低能光子更容易被吸收,随着射线穿透深度的增加,高能光子所占比例增加,射线平均能量升高,重建时依然以单能等效能量进行计算,导致重建图像出现杯状伪影,各点灰度值同实际物体的吸收系数之间出现偏差,降低了图像质量。
目前,X射线的射束硬化校正方法主要分为单能法和双能法。多项式拟合法是一种单能射束硬化校正方法,虽然实现起来简单,但灵活性差。Zou X等人在Journal ofMedical Imaging &Health Informatics(2016,6(7):1701-1707)的文章“TV-BasedCorrection for Beam Hardening in Computed Tomography”中提到的后处理法是基于重建-滤波-投影的迭代重建校正方法。一般先用前投影法校正,然后对重建出的切片进行门限值滤波,然后再进行投影、重建。该方法对各部分的材料吸收系数差别较大的物体有一定的效果,但实施起来复杂、费时。Abbema J等人在Physica Medica(2012,28(1):25-32)的文章“Feasibility and accuracy of tissue characterization with dual sourcecomputed tomography”中提出的双能量法是理论上消除射束硬化的最佳方法,它是对引起射束硬化现象的原因进行建模,实现物体位置函数和能量函数的分离,但该方法需要知道X射线的能谱分布。
上述不同方法虽然取得了一定的射束硬化伪影抑制效果,但对于锥束CT医学成像要求和工业无损检测需求,已有方法的通用性较差,实际应用中往往受到诸多限制。
发明内容
针对多能X射线衰减系数随着射线穿越路径的变长而减小,出现投影灰度随射线路径长度的非线性变化,重建时产生的杯状伪影问题,本发明提供一种锥束CT射束硬化曲线侦测及投影加权校正方法。利用射线采样矩阵信息先验知识和预先重建的切片图像,以射线与二值化图像进行求取交线长度,侦测得到射束硬化曲线。通过对拟合函数进行投影加权修正,实现多能投影的硬化校正,改善因多能投影重建引起的杯状伪影,提高图像质量。
本发明解决其技术问题所采用的技术方案包括以下步骤:
(1)设置投影采样几何参数;
(2)对多色投影Pp进行预重建,获取中心切片图像S对应的二值化图像B;
(3)获取射线与二值化图像的交线长度,拟合得到侦测硬化曲线;
(4)对侦测硬化曲线进行投影加权及修正;
(5)以修正后的函数模型对多能射线进行校正,三维重建获取切片。
在上述步骤(1)中,投影采样几何参数包括:射线源到探测器中心距离、射线源到旋转中心距离、重建分辨率及体素大小,采样幅数,采样的起始位置。
在上述步骤(1)中,获取射线采样的矩阵信息的具体步骤包括:
1)将待重建图像离散化,以M×N表示离散网格图像,称之为D,M和N分别表示沿 X方向和Y方向的网格数目;
2)构建平面坐标系,设置射线源起始坐标及探测器探元坐标,建立射线源与探测器像元之间的射线方程;
3)计算射线与二维离散网格图像D中每个像素的交线长度,以此交线长度作为对应像素位置ID的权值系数,不相交的像素位置,记为0;
4)遍历每条射线,以每条射线与二维离散网格图像D中每个像素的交线长度信息作为先验知识,存入采样矩阵。
在上述步骤(2)中,对多色投影预重建之前,需要先对投影进行散射校正,散射校正完成后需要对校正投影进行对数变换然后重建,重建时只需重建中心层切片即可,重建方法可以选择滤波反投影算法,也可以选择迭代类算法,然后对中心切片图像S根据OTSU大律法得到二值化图像B。
在上述步骤(3)中,获取射线与二值化图像的交线长度的具体步骤包括:
1)构建平面坐标系,设置射线源起始坐标及探测器探元坐标,建立射线源与探测器像元之间的射线方程,根据平面坐标系,将步骤2中得到的二值化图像B离散化,以M×N表示离散网格图像,称之为B1,M和N分别表示沿X方向和Y方向的网格数目;
2)选择不同角度下的射线,射线索引,记为R,求射线方程与二值化图像的离散网格图像B1的交点,将单个网格内的线段长度作为射线与该网格的交线长度,并记录该网格索引为
Figure GDA0003615298940000031
以此类推,获取同一射线下获得的全部网格索引
Figure GDA0003615298940000032
并对
Figure GDA0003615298940000033
对应的交线长度求和,作为该条射线穿越路径长度。
在上述步骤(3)中,拟合得到侦测硬化曲线的具体步骤包括:
1)首先根据射线索引R,提取射线索引R在探测器上投影的灰度值
Figure GDA0003615298940000034
获取射线索引R 对应的射线穿越路径长度Ls;
2)建立射线穿越路径长度Ls和投影灰度
Figure GDA0003615298940000035
的对应关系,即硬化离散数据点
Figure GDA0003615298940000036
3)选择函数L=AtB+Ct,以最小二乘法拟合离散数据点
Figure GDA0003615298940000037
得到射线穿越长度Ls 关于投影灰度
Figure GDA0003615298940000038
变化的射束硬化侦测曲线,以A、B、C作为拟合系数,L表示射线穿越物体长度,t为多色投影灰度值。
在上述步骤(4)中,对侦测硬化曲线进行修正及投影加权的具体步骤包括:
1)以拟合系数C对函数模型就行修正,得到更新模型
Figure GDA0003615298940000039
2)对更新模型增加投影加权补偿项λt,得到最终模型
Figure GDA00036152989400000310
其中λ表示加权系数,计算方法为:
Figure GDA0003615298940000041
在上述方法中,如果检测对象的扫描方式及扫描参数与上次检测相同或批量检测时,则步骤(3)中侦测的硬化曲线只需要获取一次,后面的不同检测对象可以直接应用,无需再次获取。
本发明的有益效果是:本发明提供的锥束CT射束硬化曲线侦测及投影加权校正方法,适用于任意复杂度的被测物体的硬化曲线侦测及杯状伪影校正,方法的可靠性、稳定性、抗噪性好,可在很大程度上减少锥束CT硬化投影杯状伪影对图像的干扰与影响,明显改善锥束CT图像质量。
下面结合附图和具体实施方式对本发明进行详细说明。
附图说明
图1为本发明算法流程图。
图2为重建图像杯状伪影校正前后相同位置的线性灰度比较。
具体实施方式
通过现有的工业锥束CT设备(X射线源为Comet的MXR-451HP/11,平板探测器为PerkinElmer的XRD 1621AN15 ES,并具备扫描机构、系统控制及计算用计算机),对两个零件进行投影采样,应用本发明方法对锥束CT射束硬化曲线侦测及投影加权校正方法,执行以下步骤:
(1)通过多能谱射线源工业锥束CT设备,选择射线源电压420kV和电流0.75mA,扫描几何参数为:射线源到探测器距离1241.832886mm,射线源到旋转中心距离965.323242mm;重建分辨率为512×512,圆周扫描得到检测对象锥束CT投影360幅。
1)图像网格离散化成512×512,图像网格与投影采样中的分辨率设置相同;
2)构建平面坐标系,设置射线源起始坐标及探测器探元坐标,建立射线源与探测器像元之间的射线方程;
3)计算射线与二维离散网格图像D中每个像素的交线长度,以此交线长度作为对应像素位置ID的权值系数,不相交的像素位置,记为0;
4)遍历每条射线,以每条射线与二维离散网格图像D中每个像素的交线长度信息作为先验知识,构建采样矩阵M,大小为360×262144(512×512=262144)。
(2)对多色投影预重建之前,需要先对投影进行散射校正,散射校正完成后需要对校正投影进行对数变换然后重建,重建时只需重建中心层切片即可,重建方法可以选择滤波反投影算法,也可以选择迭代类算法,然后对中心切片图像S根据OTSU大律法得到二值化图像B。
(3)获取射线与二值化图像的交线长度,拟合得到侦测硬化曲线,具体步骤包括:
1)构建平面坐标系,设置射线源起始坐标及探测器探元坐标,建立射线源与探测器像元之间的射线方程,根据平面坐标系,将步骤2中得到的二值化图像B离散化成512×512的 B1图像;
2)选择不同角度下的射线,射线索引,记为R,求射线方程与二值化图像的离散网格图像B1的交点,将单个网格内的线段长度作为射线与该网格的交线长度,并记录该网格索引为
Figure GDA0003615298940000051
3)以此类推,获取同一射线下获得的全部网格索引
Figure GDA0003615298940000052
并对
Figure GDA0003615298940000053
对应的交线长度求和,作为该条射线穿越路径长度;
4)根据射线索引R,提取射线索引R在探测器上投影的灰度值
Figure GDA0003615298940000054
获取射线索引R对应的射线穿越路径长度Ls,建立射线穿越长度和对应投影灰度关系,得到硬化离散数据点;
5)选择函数L=AtB+Ct,以最小二乘法拟合得到射线穿越长度关于多色投影灰度变化的射束硬化侦测曲线,以A、B、C作为拟合系数,L表示射线穿越物体长度,t为多色投影灰度值。
(4)对侦测硬化曲线进行投影加权及修正,具体步骤包括:
1)以拟合系数C对函数模型就行修正,得到更新模型
Figure GDA0003615298940000055
2)对更新模型增加投影加权补偿项λt,得到最终模型
Figure GDA0003615298940000056
其中λ表示加权系数,计算方法为:
Figure GDA0003615298940000057
(5)以修正后的函数模型对多能射线进行校正,三维重建获取切片。
本实施例中,基于多能谱射线硬化投影的锥束CT射束硬化曲线侦测及投影加权校正方法特点在于:
(1)根据检测对象的投影信息结合先验知识侦测得到射束硬化曲线,描述并表示了具体对象的射束硬化行为;
(2)通过对硬化曲线模型拟合完成了模型修正,经过投影加权补偿项,使得硬化校正模型更精准的建立多能投影到单能投影之间的映射;
(3)以最终模型完成多能谱硬化投影校正,抑制杯状伪影,获得高质量图像。
图2为重建图像杯状伪影校正前后相同位置的线性灰度比较,可见本发明方法可使硬化投影重建图像的杯状伪影进行抑制,使得图像轮廓对比度和清晰度得到显著提高。

Claims (3)

1.一种锥束CT射束硬化曲线侦测及投影加权校正方法,其特征在于包括下述步骤:
步骤1:设置投影采样几何参数;
步骤2:对多色投影Pp进行预重建,获取中心切片图像S对应的二值化图像B;
步骤3:获取射线与二值化图像的交线长度,拟合得到侦测硬化曲线;
步骤4:对侦测硬化曲线进行投影加权及修正;
步骤5:以修正后的函数模型对多能射线进行校正,三维重建获取切片;
在所述步骤1中,投影采样几何参数包括:射线源到探测器中心距离、射线源到旋转中心距离、重建分辨率及体素大小,采样幅数,采样的起始位置;
在所述步骤2中,对多色投影预重建之前,需要先对投影进行散射校正,散射校正完成后需要对校正投影进行对数变换然后重建,重建时只需重建中心层切片即可,重建方法可以选择滤波反投影算法,也可以选择迭代类算法,然后对中心切片图像S根据OTSU大律法得到二值化图像B;
在所述步骤3中,获取射线与二值化图像的交线长度的具体步骤包括:
1)构建平面坐标系,设置射线源起始坐标及探测器探元坐标,建立射线源与探测器像元之间的射线方程;
2)根据平面坐标系,将步骤2中得到的二值化图像B离散化,以M×N表示离散网格图像,称之为B1,M和N分别表示沿X方向和Y方向的网格数目;
3)选择不同角度下的射线,射线索引,记为R,求射线方程与二值化图像的离散网格图像B1的交点,将单个网格内的线段长度作为射线与该网格的交线长度,并记录该网格索引为
Figure FDA0003615298930000011
4)以此类推,获取同一射线下获得的全部网格索引
Figure FDA0003615298930000012
并对
Figure FDA0003615298930000013
对应的交线长度求和,作为该条射线穿越路径长度;
在所述步骤3中,拟合得到侦测硬化曲线的具体步骤包括:
1)首先根据射线索引R,提取射线索引R在探测器上投影的灰度值
Figure FDA0003615298930000014
2)获取射线索引R对应的射线穿越路径长度Ls;
3)建立射线穿越路径长度Ls和投影灰度
Figure FDA0003615298930000015
的对应关系,即硬化离散数据点
Figure FDA0003615298930000016
4)选择函数L=AtB+Ct,以最小二乘法拟合离散数据点
Figure FDA0003615298930000017
得到射线穿越长度Ls关于投影灰度
Figure FDA0003615298930000018
变化的射束硬化侦测曲线,以A、B、C作为拟合系数,L表示射线穿越物体长度,t为多色投影灰度值;
在所述步骤4中,对侦测硬化曲线进行修正及投影加权的具体步骤包括:
1)以拟合系数C对函数模型就行修正,得到更新模型
Figure FDA0003615298930000021
2)对更新模型增加投影加权补偿项λt,得到最终模型
Figure FDA0003615298930000022
其中λ表示加权系数,计算方法为:
Figure FDA0003615298930000023
2.根据权利要求1所述的一种锥束CT射束硬化曲线侦测及投影加权校正方法,其特征在于:在所述方法中,如果检测对象的扫描方式及扫描参数与上次检测相同或批量检测时,则步骤3中侦测的硬化曲线只需要获取一次,后面的不同检测对象可以直接应用,无需再次获取。
3.根据权利要求1所述的一种锥束CT射束硬化曲线侦测及投影加权校正方法,其特征在于:基于多能谱射线硬化投影的锥束CT射束硬化曲线侦测及投影加权校正方法中:
(1)根据检测对象的投影信息结合先验知识侦测得到射束硬化曲线,描述并表示了具体对象的射束硬化行为;
(2)通过对硬化曲线模型拟合完成了模型修正,经过投影加权补偿项,使得硬化校正模型更精准的建立多能投影到单能投影之间的映射;
(3)以最终模型完成多能谱硬化投影校正,抑制杯状伪影,获得高质量图像。
CN201910143911.2A 2019-02-27 2019-02-27 一种锥束ct射束硬化曲线侦测及投影加权校正方法 Active CN109919868B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910143911.2A CN109919868B (zh) 2019-02-27 2019-02-27 一种锥束ct射束硬化曲线侦测及投影加权校正方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910143911.2A CN109919868B (zh) 2019-02-27 2019-02-27 一种锥束ct射束硬化曲线侦测及投影加权校正方法

Publications (2)

Publication Number Publication Date
CN109919868A CN109919868A (zh) 2019-06-21
CN109919868B true CN109919868B (zh) 2022-10-04

Family

ID=66962533

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910143911.2A Active CN109919868B (zh) 2019-02-27 2019-02-27 一种锥束ct射束硬化曲线侦测及投影加权校正方法

Country Status (1)

Country Link
CN (1) CN109919868B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111487265B (zh) * 2020-05-25 2023-06-16 中国人民解放军战略支援部队信息工程大学 一种结合投影一致性的锥束ct硬化伪影校正方法
CN113034636B (zh) * 2021-03-09 2022-04-29 浙江大学 基于跨尺度多能谱ct标签的锥束ct图像质量改善方法和装置
CN116543088B (zh) * 2023-07-07 2023-09-19 有方(合肥)医疗科技有限公司 Cbct图像重建方法及装置

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101126722A (zh) * 2007-09-30 2008-02-20 西北工业大学 基于配准模型仿真的锥束ct射束硬化校正方法
CN101509879A (zh) * 2009-03-17 2009-08-19 西北工业大学 一种ct快速批量扫描与校正方法
CN101533098A (zh) * 2009-04-08 2009-09-16 西北工业大学 Ct射束硬化校正中的噪声抑制方法
CN102609908A (zh) * 2012-01-13 2012-07-25 中国人民解放军信息工程大学 基于基图像tv模型的ct射束硬化校正方法
CN105528771A (zh) * 2016-01-19 2016-04-27 南京邮电大学 一种使用能量函数方法的锥束ct中杯状伪影的校正方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7778392B1 (en) * 2004-11-02 2010-08-17 Pme Ip Australia Pty Ltd Method of reconstructing computed tomography (CT) volumes suitable for execution on commodity central processing units (CPUs) and graphics processors, and apparatus operating in accord with those methods (rotational X-ray on GPUs)
US10605933B2 (en) * 2016-10-21 2020-03-31 Carestream Health, Inc. X-ray spectral calibration technique for cone-beam CT

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101126722A (zh) * 2007-09-30 2008-02-20 西北工业大学 基于配准模型仿真的锥束ct射束硬化校正方法
CN101509879A (zh) * 2009-03-17 2009-08-19 西北工业大学 一种ct快速批量扫描与校正方法
CN101533098A (zh) * 2009-04-08 2009-09-16 西北工业大学 Ct射束硬化校正中的噪声抑制方法
CN102609908A (zh) * 2012-01-13 2012-07-25 中国人民解放军信息工程大学 基于基图像tv模型的ct射束硬化校正方法
CN105528771A (zh) * 2016-01-19 2016-04-27 南京邮电大学 一种使用能量函数方法的锥束ct中杯状伪影的校正方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
Geometric Calibration Method for Multiple-Head Cone-Beam SPECT System;Ph. Rizo 等;《IEEE TRANSACTIONS ON NUCLEAR SCIENCE》;19941231;第41卷(第6期);第2748-2757页 *
基于K-N模型的锥束CT散射伪影校正方法;刘建邦 等;《光学学报》;20181130;第38卷(第11期);第1-9页 *
基于基图像TV模型的CT射束硬化校正方法;李庆亮 等;《核电子学与探测技术》;20120331;第32卷(第3期);第287-291页 *
基于线框模型的锥束CT几何参数校正方法;王珏 等;《仪器仪表学报》;20180228;第39卷(第2期);第177-184页 *

Also Published As

Publication number Publication date
CN109919868A (zh) 2019-06-21

Similar Documents

Publication Publication Date Title
US10274439B2 (en) System and method for spectral x-ray imaging
JP7202302B2 (ja) 断層撮影再構成に使用するためのデータのディープラーニングに基づく推定
JP7139119B2 (ja) 医用画像生成装置及び医用画像生成方法
US8000435B2 (en) Method and system for error compensation
JP2018140165A (ja) 医用画像生成装置
US10111638B2 (en) Apparatus and method for registration and reprojection-based material decomposition for spectrally resolved computed tomography
CN109919868B (zh) 一种锥束ct射束硬化曲线侦测及投影加权校正方法
US8666137B2 (en) Apparatus and method for processing projection data
WO2007148263A1 (en) Method and system for error compensation
CN109920020B (zh) 一种锥束ct病态投影重建伪影抑制方法
CN107106109B (zh) 计算机断层扫描系统
JP2019111346A (ja) 医用処理装置及び放射線診断装置
Yang et al. Cupping artifacts correction for polychromatic X-ray cone-beam computed tomography based on projection compensation and hardening behavior
Li et al. Multienergy cone-beam computed tomography reconstruction with a spatial spectral nonlocal means algorithm
Xia et al. Patient‐bounded extrapolation using low‐dose priors for volume‐of‐interest imaging in C‐arm CT
US9916669B2 (en) Projection data correction and computed tomography value computation
CN106228584B (zh) 锥束ct圆加直线轨迹反投影滤波重建方法
CN109870471B (zh) 一种单光栅侦测的锥束ct角度序列散射获取方法
CN111899312B (zh) 一种迭代补偿的有限角度ct投影重建方法
EP3920144A1 (en) Methods and systems for multi-material decomposition
Guo et al. Improved iterative image reconstruction algorithm for the exterior problem of computed tomography
CN116172594B (zh) 用于生成患者的结果图像数据集的方法及装置
US20230083935A1 (en) Method and apparatus for partial volume identification from photon-counting macro-pixel measurements
EP4047557A1 (en) Projection-domain material decomposition for spectral imaging

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