CN106796733B - 用于谱成像的噪声降低 - Google Patents

用于谱成像的噪声降低 Download PDF

Info

Publication number
CN106796733B
CN106796733B CN201580055482.3A CN201580055482A CN106796733B CN 106796733 B CN106796733 B CN 106796733B CN 201580055482 A CN201580055482 A CN 201580055482A CN 106796733 B CN106796733 B CN 106796733B
Authority
CN
China
Prior art keywords
projection data
spectral
pairs
noise
polarity
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
CN201580055482.3A
Other languages
English (en)
Other versions
CN106796733A (zh
Inventor
K·M·布朗
R·莱文森
G·谢克特
M·P·弗赖曼
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.)
Koninklijke Philips NV
Original Assignee
Koninklijke Philips NV
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 Koninklijke Philips NV filed Critical Koninklijke Philips NV
Publication of CN106796733A publication Critical patent/CN106796733A/zh
Application granted granted Critical
Publication of CN106796733B publication Critical patent/CN106796733B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/002D [Two Dimensional] image generation
    • G06T11/003Reconstruction from projections, e.g. tomography
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T11/002D [Two Dimensional] image generation
    • G06T11/003Reconstruction from projections, e.g. tomography
    • G06T11/005Specific pre-processing for tomographic reconstruction, e.g. calibration, source positioning, rebinning, scatter correction, retrospective gating
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N23/00Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00
    • G01N23/02Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material
    • G01N23/04Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material
    • G01N23/046Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material using tomography, e.g. computed tomography [CT]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2211/00Image generation
    • G06T2211/40Computed tomography
    • G06T2211/408Dual energy

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Pulmonology (AREA)
  • Radiology & Medical Imaging (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Biochemistry (AREA)
  • General Health & Medical Sciences (AREA)
  • Immunology (AREA)
  • Pathology (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Analysing Materials By The Use Of Radiation (AREA)

Abstract

一种系统,包括:存储器(420),其具有用于以下中的至少一项的指令:处理谱CT投影数据以减轻所述谱CT投影数据的噪声或所述谱CT投影数据的噪声引发的偏差中的至少一个,或生成减轻所述谱CT投影数据的所述噪声引发的偏差的分解算法。所述系统还包括:处理器(418),其运行所述指令并且进行以下中的至少一项:处理所述谱CT投影数据,或者生成所述分解算法并将所述谱CT投影数据分解为基材料。所述系统还包括:重建器(434),其重建所述基材料,从而生成谱图像。

Description

用于谱成像的噪声降低
技术领域
以下总体上涉及谱成像,并且利用对计算机断层摄影(CT)的特定应用来进行描述。
背景技术
谱(或多能量)CT利用在多个不同光子能量处采集的多个衰减值来求解光电效应、康普顿散射和(一个或多个)其他分量(例如,一个或多个K边缘)对材料的质量衰减系数的贡献。这允许对虚拟单色图像、碘浓度图、虚拟非对比图像等以及常规CT图像的重建。存在用于采集这样的衰减值的若干方法,包括使用多个X射线管、kVp切换、多层探测器,以及光子计数探测器。
利用双能量CT,文献已经指示材料的能量依赖性衰减(μi(E))能够被近似为如在公式1中示出的在水中的光电效应
Figure GDA0002918977520000011
贡献与康普顿散射
Figure GDA0002918977520000012
贡献的衰减分布图的线性组合:
公式1:
Figure GDA0002918977520000013
能够根据对应的线积分(lp和ls)来重建具有表示光电效应和康普顿效应的系数Cp和Cp的像素值的图像的集合。能够通过公式2从两个不同能量谱FH和FL(L=低能量,并且H=高能量)来获得线积分:
公式2:
Figure GDA0002918977520000021
这里,FH/L(E)项是高谱和低谱的通量,所述高谱和低谱的通量是通过对从上游射束朝向探测器像素发出的通量以及探测器响应进行建模来获得的。因管调制或因两个探测器层而造成两个谱之间的差。
lp/s项是水中的光电和散射等效路径,并且
Figure GDA0002918977520000022
是这两种机制的能量依赖性衰减系数。能够通过转换公式2来求解lp和ls,其得到从ps和pp的域到pL和pH的域的映射函数:
Figure GDA0002918977520000023
这里,ps和pp是通过将lp/s乘以固定的标量衰减系数从lp/s获得的,以便稍后为了方便在无维度预备值上进行处理。已经通过多项式函数或查找表求解出
Figure GDA0002918977520000024
遗憾的是,分解函数
Figure GDA0002918977520000025
是非线性的。这样,在存在有噪声的数据的情况下,分解的平均值并不等于输入的平均值的分解,或者
Figure GDA0002918977520000026
Figure GDA0002918977520000027
结果,分解将包括偏差。该噪声引发的偏差在正弦图内变化并传播到重建的光电图像和散射图像中的人为偏差。此外,这导致不准确的碘浓度估计结果,并且导致在虚拟单色图像中的可见的图像伪影。
低能量投影和高能量投影中的每个都好具有因X射线光子的偶然性(光子射击噪声)产生的统计分布,其能够由高斯函数在良好的近似内进行描述,作为光子泊松行为的结果。在单色近似内,并且忽略电子噪声,投影值将具有正态分布,所述正态分布具有等于光子计数值(μ)均方根的倒数的标准偏差(σ)。低能量投影和高能量投影对的分布并不相关且具有随机极性,每个极性组合(+,+)、(-,+)、(-,-)和(+,-)是相等的。根据图1,低能量投影和高能量投影对被均匀地分布在相同的极性象限Ⅰ和Ⅲ以及相反的极性象限Ⅱ和Ⅳ上。
利用双能量CT,分解算法生成两个数据集,例如,光电数据集(PE)和康普顿散射数据集(Sc)。注意,PE/Sc分解能够被认为等于基材料分解(即,水/碘)到一“角度”常量内。针对每个投影“p”和视图“v”,两个数据集被成对产生:(低能量投影,高能量投影)p,v=>(PE,Sc)p,v。分解创建光电效应数据集与康普顿散射数据集之间的负相关。相关因子是光电分量和康普顿散射分量的双能量衰减系数的函数。图2示出了PE、Sc对)。
图3示出了经由分解映射308被映射到Sc和PE的低能量投影数据和高能量投影数据的值306的低能量投影数据和高能量投影数据对302,它们具有距平均值304相等的变化。遗憾的是,分解对低能量投影数据和高能量投影数据对的相对极性敏感,并且针对相反极性的低能量投影数据和高能量投影数据对噪声被放大。这通过象限Ⅱ和Ⅳ中的大的变化310和312能够看出。归因于光子射击噪声,时间的大致50%发生“相反极性”变化。
US 2012/0134561 A1公开了用于对相关联的噪声进行抑制的方法和系统。所述相关联的噪声抑制技术估计分别对应于第一基材料和第二基材料的第一MD图像与第二MD图像中的噪声值之间的关联方向。使用所估计的关联方向来扩散这两幅MD图像,以生成第一扩散图像和第二扩散图像。另外,通过从对应的MD图像减去扩散图像来生成第一噪声掩模和第二噪声掩模。分别利用第一噪声掩模和第二噪声掩模来处理第一MD图像和第二MD图像中的边缘,以生成最终的第一噪声掩模和最终的第二噪声掩模。然后利用最终的第二噪声掩模来处理第一MD图像以生成最终的第一MD图像,并且利用最终的第一噪声掩模来处理第二MD图像以生成最终的第二MD图像。
James T.Dobbins等人的“Dose reduction in CT with correlated-polaritynoise reduction:comparable image quality at half the dose with projectionspace processing”(Proceedings of Spie,第8668卷,第8668页)公开了一种相关联的极性噪声降低(CPNR),其是使用统计学方法来降低噪声同时保持优秀的解析和“正常”噪声外观的一种新颖的噪声降低技术。在该公开中,首次展现了在降低CT图像中的噪声中使用CPNR作为降低CT中的剂量的手段。CPNR在投影图像空间中被应用到半剂量图像,并且然后利用Feldkamp方法使用滤波反投影来重建图像。
James T.Dobbins的“Correlated polarity noise reduction for dual-energyimaging”(Radiological Society of North America 82nd scientific assembly andannual meeting的摘要,12月1-6日,1996年)公开了通过其极性匹配图像噪声的极性的高斯噪声在双能量放射摄影中降低噪声的研究报告。还研究了一种新技术,其中,从每个像素建成区随机变量,假设只有噪声分量的极性而不是幅度是已知的。
本文中描述的各方面解决了上述问题和其他问题。
发明内容
在一个方面中,一种系统,包括:存储器,其具有用于以下的至少一项的指令:通过以下来处理谱CT投影数据以减轻所述谱CT投影数据的噪声或所述谱CT投影数据的噪声引发的偏差中的至少一个:将所述谱CT投影数据的低能量投影数据和高能量投影数据对分类为相同极性对或相反极性对,将所述相反极性对变换为相同极性对,以创建经修改的谱CT投影数据,并且分解所述经修改的谱CT投影数据,以生成基材料;或者生成减轻所述谱CT投影数据的所述噪声引发的偏差的分解算法。所述系统还包括:处理器,其运行所述指令并且进行以下中的至少一项:处理所述谱CT投影数据,或生成所述分解算法并将所述谱CT投影数据或所述经修改的谱CT投影数据分解到基材料。所述系统还包括:重建器,其重建所述基材料,从而生成谱图像。
在另一方面中,一种方法,包括:进行以下中的至少一项:通过以下来处理谱CT投影数据以减少噪声或噪声引发的偏差:将所述谱CT投影数据的低能量投影数据和高能量投影数据对分类为相同极性对或相反极性对,将所述相反极性对变换为相同极性对,并且创建经修改的谱CT投影数据;或生成减小针对所述谱CT投影数据的所述噪声引发的偏差的分解算法。所述方法还包括:分解所述谱CT投影数据或所述经修改的谱CT投影数据以生成基材料。所述方法还包括:重建所述基材料以生成谱图像。
在另一方面中,一种计算机可读存储介质,其被编码有计算机可读指令。所述计算机可读指令当由处理器运行时令所述处理器:进行以下中的至少一项:通过以下来处理双能量CT投影数据:将所述谱CT投影数据的低能量投影数据和高能量投影数据对分类为相同极性对或相反极性对,将所述相反极性对变换为相同极性对,并且创建经修改的谱CT投影数据,或生成用于所述双能量CT投影数据的分解算法;分解所述双能量CT投影数据或所述经修改的谱CT投影数据以生成基材料;并且通过重建所述基材料来生成谱图像。
附图说明
本发明可以采取各种部件和各部件的布置,以及各个步骤和各步骤的安排的形式。附图仅出于图示优选实施例的目的,并且不应被解释为限制本发明。
图1示出了被均匀分布在相同的极性象限Ⅰ和Ⅲ以及相反的极性象限Ⅱ和Ⅳ中的低能量投影与高能量投影对。
图2示出了光电分量与康普顿散射分量对。
图3示出了被映射到低能量投影数据和高能量投影数据的值的低能量投影数据与高能量投影数据对,它们具有距平均值相等的变化。
图4示意性地图示了具有控制台的范例成像系统,所述控制台具有材料基分解器和投影数据处理器。
图5示意性地图示了投影数据处理器的范例。
图6示出了根据空间相邻项的集合来计算局部平均值的范例。
图7示出了根据时间相邻项的集合来计算局部平均值的范例。
图8示意性地图示了投影数据处理器的另一范例。
图9示出了分解输入范围和它的列中的一个的位置。
图10示出了表示沿着它的列中的一个的去卷积的解的曲线。
图11示出了针对光电效应LUT的解。
图12示出了针对光电效应LUT的范例值。
图13示出了使用并不包括去偏差的分解算法的现有技术方法。
图14示出了对应于图13的散点图。
图15示出了使用本文讨论的包括去偏差的分解算法的光电效应图像。
图16示出了对应于图15的散点图。
图17图示了范例方法。
图18图示了用于模拟投影域中的减小的造影剂图像的范例方法。
图19图示了用于模拟投影域中的减小的造影剂图像的范例方法。
具体实施方式
首先参考图4,图示了诸如计算机断层摄影(CT)扫描器的成像系统400。成像系统400包括固定机架402和旋转机架404,所述旋转机架404由固定机架402可旋转地支撑并关于z轴围绕检查区域406旋转。诸如卧榻的对象支撑物407在检查区域406中支撑对象或目标。对象支撑物407能协调扫描而移动,以便相对于检查区域406引导对象或目标以用于扫描对象或目标。
诸如X射线管的辐射源408由旋转机架404可旋转地支撑,与旋转机架404一起旋转,并发出贯穿检查区域406的辐射。在一个实例中,控制器控制辐射源408的均值发射电压或峰值发射电压。这包括在一积分时间段内在两个或更多个发射电压(例如,80kVp与140kVp、100kVp与120kVp等)之间和/或以其他方式切换发射电压。在变型中,成像系统400包括以两个不同的发射电压发射辐射的至少两个辐射源408。在另一变型中,辐射源408包括单宽谱X射线管。
探测器阵列412相对于辐射源408在检查区域406对面对向一角度弧。探测器阵列412探测贯穿检查区域406的辐射并生成指示该辐射的投影数据。当在至少两个发射电压之间切换辐射源电压和/或两个或更多个X射线管以两个不同的发射电压发射辐射的情况下,探测器阵列412生成针对辐射源电压中的每个的投影数据。针对单宽谱X射线管,探测器阵列412包括产生谱投影数据的能量解析探测器(例如,多层、光子计数等)。
计算系统服务于操作者控制台416并且包括至少一个处理器418(例如,微处理器、中央处理单元等),所述至少一个处理器418运行被存储在计算机可读存储介质(“存储器”)420中的至少一条计算机可读指令,所述计算机可读存储介质排除瞬态介质并包括物理存储器和/或其他非瞬态介质。微处理器418也可以运行由载波、信号或其他瞬态介质承载的一条或多条计算机可读指令。计算系统416还包括(一个或多个)输出设备422(例如,显示监视器、影片放映机等)和(一个或多个)输入设备424(例如,鼠标、键盘等)。
在所图示的范例中,至少一个计算机可读指令实施分解器426和投影数据处理器430,所述分解器426基于(一个或多个)分解算法428来分解谱投影数据,所述投影数据处理器430基于(一个或多个)投影数据处理算法432来处理谱投影数据。如在下文中更加详细地描述的,在一个实例中,分解器426使用(一个或多个)分解算法将谱投影数据分解成基材料,例如,光电分量和康普顿散射分量和/或其他基材料,所述(一个或多个)分解算法减少噪声引发的偏差。也如在下文中更加详细地描述的,在一个实例中,投影数据处理器430使用(一个或多个)处理算法432处理谱投影数据,所述(一个或多个)处理算法432减小噪声引发的偏差和/或减小投影数据噪声。
(一个或多个)分解算法428减小噪声引发的偏差,或另外(一个或多个)分解算法428能够与经处理的谱投影数据一起使用。在变型中,分解器426或投影数据处理器430中的一个或多个被定位在与控制台416分离且不同的计算装置中。此外,在另一变型中,投影数据处理器430被省略。利用该变型,减小噪声引发的偏差的(一个或多个)分解算法428被利用。能够基于默认算法、指示用户选择的感兴趣算法的输入、成像协议等来采用特定的分解算法。同样地,能够基于默认算法、指示用户选择的感兴趣算法的输入、成像协议等来采用特定的投影数据处理算法。
重建器434重建经分解的谱投影数据。在一个实例中,重建器434重建能量依赖性偏差材料中的一个或多个,生成对应于一个或多个不同能量的一幅或多幅图像。额外地或备选地,重建器434组合经分解的投影数据并在整个能量谱上重建非谱(或常规的)图像数据。在又一实例中,对对应于一个或多个不同能量的一幅或多幅图像进行组合以产生非谱(或常规的)图像数据。经分解的谱投影数据也能够用于重建虚拟单色图像、碘浓度图、虚拟非对比图像等。
所图示的控制台416处理从成像系统400获得的谱投影数据。在变型中,谱投影数据是从不同的成像系统和/或数据存储库获得的,所述数据存储库例如为图片存档与通信系统(PACS)、放射科信息系统(RIS)、医院信息系统(HIS)、电子病历(EMR)、数据库、服务器、成像系统、计算机和/或其他数据存储库。谱投影数据能够经由医学数字成像与通信(DICOM)、健康水平7(HL7)和/或其他协议来传输。
图5示意性地图示了投影数据处理器430的范例。所图示的投影数据处理器430包括投影数据极性修改器506和校正因子确定器508。出于简明和清晰的目的,下文结合双能量投影数据来进行描述。然而,应当理解,投影数据处理器430能够投影针对三个或更多个不同的能量谱投影数据。
投影数据极性分类器502将低能量投影数据与高能量投影数据对分类为属于相同极性象限(图2和图3中的象限Ⅰ和Ⅲ)或属于相反极性象限(图2和图3中的象限Ⅱ和Ⅳ)。分类能够通过将测得的低能量投影值和高能量投影值与低能量投影数据和高能量投影数据的局部平均值进行比较来完成。能够根据空间相邻项和/或时间相邻项的集合来计算局部平均值。
图6示出了根据空间相邻项的集合来计算局部平均值的范例。在该范例中,针对像素(n)600的局部平均值被确定为8个空间相邻项(n1(602)、n2(604)、n3(606)、n4(608)、n5(610)、n6(612)、n7(614)、n8(616))的平均值(ni),其中,i=8。图7示出了根据时间相邻项的集合计算局部平均值的范例。在该范例中,针对像素(n)700和视图“j”的局部平均值被确定为8个时间相邻项(vj-4(702)、vj-3(704)、vj-2(706)、vj-1(708)、vj+1(710)、vj+2(712)、vj+3(714)、vj+4(716))的平均值(vj),其中,j=-4,4。本文也预期其他方法。
转到图5,局部低能量投影平均值与局部高能量投影平均值并不是真的低能量投影数据值和高能量投影数据值。然而,它们在统计上比个体的低能量投影数据值和高能量投影数据值更接近真值。分类过程不必具有100%的准确率。对于泊松统计行为(光子射击噪声),针对低频率数据,基于8个最接近的空间相邻项和8个最接近的时间相邻项的成功率高于90%。
当代替相反极性而计算相同极性时,误分类发生,或者反之亦然。对于前者,不采取行动,使用原始投影数据,并且光电数据与康普顿散射数据对的变化不减小。对于后者,对高能量投影数据进行增量,低能量投影数据与高能量投影数据对现在处于错误的象限中,并且光电与康普顿散射对(或光电表面密度和康普顿散射表面密度)的变化增加。误分类发生的时间小于10%并且引起总体噪声减小的小的减小。
投影数据相反极性识别器504基于分类来识别属于相反极性的低能量投影数据与高能量投影数据对。针对具有相同极性的低能量投影数据与高能量投影数据对(图3中的象限Ⅰ和Ⅲ),投影数据极性修改器506并不修改低能量投影数据与高能量投影数据对极性。然而,针对具有相反极性的低能量投影数据与高能量投影数据对(图3中的象限Ⅱ和Ⅳ),极性修改器506改变极性。
例如,如果高能量投影值大于局部平均值(其指示相反极性象限Ⅱ中的高能量投影值),则投影数据极性修改器506将高能量投影值修改为相同极性象限Ⅲ中的高能量投影值。在一个实例中,投影数据极性修改器506修改高能量投影值,如公式3所示:
公式3:
高_能量_投影_值=高_能量_投影_值-2σ
在另一范例中,如果高能量投影值小于局部平均值(其指示相反极性象限Ⅳ中的高能量投影值),投影数据极性修改器506将高能量投影值修改为相同极性象限Ⅰ中的高能量投影值。在一个实例中,投影数据极性修改器506修改高能量投影值,如公式4所示:
公式4:
高_能量_投影_值=高_能量_投影_值+2σ
2σ校正因子表示“最佳猜想”。归因于泊松行为,高能量投影值距“真值”的变化将在统计上进行分布,并且平均变化为σ。修改创建具有相同极性变化的低能量投影数据与高能量投影数据对。这样,为了校正高能量投影值变化的极性,使用两倍的平均变化值(2σ)。2σ校正值将具有0到2σ的高能量投影值变化的低能量投影数据与高能量投影数据对从象限Ⅱ移动到象限Ⅲ或从象限Ⅳ移动到象限Ⅰ。
投影数据校正因子确定器508确定校正因子。在一个非限制性实例中,投影数据校正因子确定器508确定校正因子,如公式5所示:
公式5:
2*σ=2*SQRT[(kV/mA因子)*((exp(高_能量_投影)-1)]
其中,“kV/mA因子”项是针对每个CT扫描器测得的经验值,并且通常是在探测器表面处的(空气中的)X射线管通量的量度。
材料基分解器426(图4)分解极性调节的投影数据。该投影数据处理器430并不要求对分解算法或处理链的其他部分的任何改变。
图8示意性地图示了投影数据处理器430的另一范例。在该范例中,投影数据处理器430包括投影数据滤波器802、噪声确定器804,以及投影数据值选择器806。
投影数据滤波器802滤除低能量投影数据和高能量投影数据。合适的滤波器的范例包括,但不限于,高斯滤波器(例如,3×3等)、自适应滤波器、中值滤波器等。这样的滤波器对低能量投影数据和高能量投影数据进行平滑。
噪声确定器804确定原始投影数据与经滤波的投影数据之间的噪声值。这能够通过确定原始投影数据与经滤波的投影数据之间的差(例如,通过减去、加上负值,等)来实现。
投影数据值选择器806选择投影数据以基于噪声和预定噪声阈值的集合808来进行分解。在噪声满足预定噪声阈值的情况下,投影数据值选择器806选择经滤波的投影。否则,投影数据值选择器806选择原始(或未经滤波的)投影。
在一个实例中,针对低能量投影数据和高能量投影数据设定不同阈值808,并且当低能量投影数据噪声为正且高能量投影数据噪声为负时选择经滤波的数据。在存在噪声的情况下,信号pL和pH上的梯度将几乎总是处于相同的方向。如果在扫描状况下低信号和高信号在相反方向上移动,则很高的可能性是这是归因于噪声而并不是归因于真实的目标特征。
在另一实例中,只有在低信号误差大于第一预定阈值且高信号小于第二预定阈值时,低值才被设定到经滤波的低值,并且只有在高信号误差小于第三阈值(其在一个实例中为负的第一阈值)且低信号误差大于第四阈值(其在一个实例中为负的第二阈值)时,高值才被设定到滤波器高值。
在又一实例中,在低误差与高误差之间的差的绝对值大于预定阈值时低值被设定到经滤波的低值。本文也预期其他选择标准。材料基分解器426(图4)分解所选择的投影数据。该投影数据处理器430并不要求对分解算法或处理链的其他部分的任何改变。
图9图示了其中材料基分解器426使用减小噪声引发的偏差的(一个或多个)分解算法428来分解谱投影数据的范例。利用该配置,投影数据处理器430和(一个或多个)投影数据处理算法432能够被省略。材料基分解器426包括分解查找表(LUT)生成器902、卷积核确定器904以及卷积器906。
通常,(一个或多个)分解算法减轻在估计的光电和散射等效路径中的人为偏差。适合用于用谱投影数据替换无噪声输入谱投影数据的分解函数不仅考虑了低谱投影和高谱投影,而且还考虑了它们的差的估计的噪声。在一个实例中,通过对适当的1D去卷积进行求解来执行抑制偏差。
简洁地转到图10,示出了分解的输入范围。在图10中,第一轴1002和第二轴1004分别表示高能量投影数据和低能量投影数据。区域1006表示对应于朝向给定的探测器像素发射的射线路径的低投影数据值和高投影数据值所有实际对。区域1006归因于具有扇形角度和锥形角度的谱的变化而从一个探测器像素改变到另一个探测器像素。区域1008表示针对所有探测器像素取单位的预计算的查找表(LUT)边界。LUT输入范围被对准到p1轴和p2轴。线1010表示一个特定的LUT列。
当在两个谱的投影之间不存在噪声相关的情况下,p1和p2是不相关的并且包含相似的噪声水平。计算被施加到适合用于关于p1和p2无噪声输入LUT值的离散的二阶导数矩阵的本征向量和本征值示出LUT局部非线性主轴非常接近p2。因此,只有p2的噪声分量负责人们想要避免的噪声引发的偏差。沿着p2的非线性表达针对与高原子量材料相交的射线路径的谱分离的损失。其在p2的负值和p1的高正值处变得主导。
返回图9,分解LUT生成器902生成针对散射分量和光电分量两者的从
Figure GDA0002918977520000117
Figure GDA0002918977520000116
的LUT。除了常规的分解输入,即,低能量投影值和高能量投影值,LUT采取针对p2的噪声标准偏差(STD)的估计值作为输入,其由
Figure GDA0002918977520000111
来指代。针对
Figure GDA0002918977520000112
的每个可能值,逐列计算LUT,即,通过针对p1的每个固定值计算针对所有p2的LUT。1D函数(由
Figure GDA0002918977520000113
指代)通过满足公式6来构建:
公式6:
Figure GDA0002918977520000114
在公式6中,p2被表达为由s2指代的无噪声信号分量与由n2指代的噪声分量的总和:p2=s2+n2。左手侧上的上标代表其下面的统计预期值。右手侧的变量
Figure GDA0002918977520000121
是LUT列与区域1006的边缘的相交,如在图10中能够看见的。右手侧的函数L(p1,s2)包含适合无噪声输入值的LUT列值。
最终的值是通过对如公式2所示给定的一个有关的倒数进行数值求解来获得。公式7的左手侧能够被写为无噪声解L(p1,s)与噪声分布核之间的卷积。排除探测器处的低信号的情况和/或在分解之前不对预备值施加任何去噪声算法,卷积核确定器904可以将卷积核
Figure GDA0002918977520000122
确定为公式7中示出的近似的高斯项。
公式7:
Figure GDA0002918977520000123
将公式6与公式7进行组合,
Figure GDA0002918977520000124
能够被表达为针对在公式8中示出的去卷积的解,这由去卷积器906来执行:
公式8:
Figure GDA0002918977520000125
为了避免断点,关于p2限制
Figure GDA0002918977520000126
的全局变化。为了不增加去分解的正弦图内的噪声,获得针对
Figure GDA0002918977520000127
的平滑且非振荡解。
在这些指导之后,
Figure GDA0002918977520000128
是在公式9中示出的二次三项代价函数的最小化器:
公式9:
Figure GDA0002918977520000129
C(f)=F(f)+S(f)+R(f);
其中,F(f)是保真度项,S(f)是有限输入范围的辅助函数,并且R(f)是平滑性正则化项。考虑到沿着p2的分解LUT的采样点的离散性,这三项能够如公式10、11和12来表达:
公式10:
Figure GDA0002918977520000131
公式11:
Figure GDA0002918977520000132
以及
公式12:
Figure GDA0002918977520000133
其中,mf和ml
Figure GDA0002918977520000134
附近的采样点的指数,βS和βR是相似度和平滑性正则化参数,M是一列中的LUT点的数量,并且A是等于点
Figure GDA0002918977520000135
与点
Figure GDA0002918977520000136
之间的L(p1,s2)的辅助函数且从该段推被外推。
针对
Figure GDA0002918977520000137
的解和无噪声情况解L,以及辅助函数A被展示在图11中以用于光电LUT。LUT边界处的A的导数是特意地约为零。该外推的选择得到在接近边界处表现相似的
Figure GDA0002918977520000138
的解。这允许我们通过在LUT边界以外给其指定恒定值以用于计算保真度项F(f)来以平滑的方式沿着p2从LUT外推出
Figure GDA0002918977520000139
能够防止A在p2的低值处超过10000的典型值。
对针对LUTL(p1,s2)的无噪声解利用噪声核
Figure GDA00029189775200001310
进行卷积将其从其本身移位,参见黑色虚线曲线。相比之下,利用
Figure GDA00029189775200001311
Figure GDA00029189775200001312
的进行卷积得到蓝色虚线曲线,其与L(p1,s2)的原始分布图匹配地相当好。这展示了对获得的噪声引发的偏差的抑制。归因于代价函数C(f)的二次性质,公式10能够由标准加权的最小二乘拟合来执行。为了改善数值稳定性,能够使用共轭梯度方法。
在图12中示出了针对光电LUT的βR的范例值。βS的值被选取为比βR的值小大约两个数量级。利用公式10计算LUT值以用于
Figure GDA0002918977520000141
值的离散集合。针对这样的集合的范例具有0.02无维度自然预备单元的最小值、0.02的增量,以及0.7的最大值。采取
Figure GDA0002918977520000142
以及从p和p直接获得的p1和p2作为我们的三个输入参数,能够针对每个探测器像素构建3D LUT。然而,为了减小计算时间和限额,能够针对包含在探测器像素的2D子集中的探测器像素的子集生成3D LUT。多个这些预计算的3D LUT然后能够被在线内插以用于所有剩余的探测器像素。
然而,对于由这些LUT分解任何给定的双能量正弦图所剩下的是在线计算针对
Figure GDA0002918977520000143
的估计的正弦图。能够通过并入谱模型的噪声传播方法基于p和p的正弦图以及已知的瞬时管电流来完成对该正弦图的估计。
图13和图15示出了针对圆形扫描的包含碘、空气和CaCl2管的固态水体模的光电效应图像,其中,mAs=133。窗口水平/宽度为50/-947HU。
图13是在没有去偏差的情况下从沿着p2剪裁分解LUT而获得的。图像经受源自于水中的光等效路径的断点的条纹。图像也经受人工白化。图14示出了对应的散点图。在图14中,y轴1402表示光电效应(以HU为单位),x轴1404表示康普顿散射(以HU为单位),并且点1406表示mAs。
图15是利用本文讨论的分解方法获得的并且缺少图13的伪影。在由紫色圆圈环绕的感兴趣区域内的平均HU等于-942(针对a)和-955(针对b)。570mg/cc CaCl2的杆由黄色虚线圆圈环绕。图16示出了对应的散点图。在图16中,y轴1402表示光电效应(以HU为单位),x轴1404表示康普顿散射(以HU为单位),并且点1602表示mAs。图14和图16的散点图示出了本文描述的方法能够如何抑制在该杆内的mAs依赖性噪声引发的HU偏差。
图17图示了用于减小在光电分量和康普顿散射分量中的谱投影数据噪声的范例方法。
应当理解,动作的排序并非限制。这样,本文预期其他排序。另外,可以省略一个或多个动作和/或可以包括一个或多个额外的动作。
在1702处,获得谱CT投影数据。谱投影数据包括至少两个不同的能量谱。
在1704处,谱投影数据被分类成在相同极性象限和相反极性象限上的低能量投影数据和高能量投影数据对。
在1706处,相反极性对被变换到相同极性对。
在1708处,经调节的投影数据被分解为光电分量和康普顿散射分量。
在1710处,重建光电分量和康普顿散射分量,产生谱图像。
图18图示了用于模拟投影域中的减小的造影剂图像的范例方法。
应当理解,动作的排序并非限制。这样,本文预期其他排序。另外,可以省略一个或多个动作和/或可以包括一个或多个额外的动作。
在1802处,获得谱CT投影数据。谱投影数据包括至少两个不同的能量谱。
在1804处,通过向数据施加平滑滤波器来对谱投影数据进行平滑。
在1806处,确定原始数据与经平滑的数据之间的差。
在1808处,选择差满足预定阈值的经滤波的数据,并且以其他方式选择原始数据。
在1810处,将经选择的投影数据分解为光电分量和康普顿散射分量。
在1812处,重建光电分量和康普顿散射分量,产生谱图像。
图19图示了用于模拟投影域中的减小的造影剂图像的范例方法。
应当理解,动作的排序并非限制。正因如此,本文预期其他排序。另外,可以省略一个或多个动作和/或可以包括一个或多个额外的动作。
在1902处,获得谱CT投影数据。谱投影数据包括至少两个不同的能量谱。
在1904处,估计谱投影数据的噪声。
在1906处,基于噪声生成卷积核。
在1908处,通过对卷积核进行去卷积来确定分解结果。
在1910处,分解结果被用于分解谱投影数据并生成光电分量和康普顿散射分量。
在1912处,重建光电分量和康普顿散射分量,产生谱图像。
本文的方法可以借助于计算机可读指令来实施,所述计算机指令被编码或被嵌入在计算机可读存储介质上,所述计算机可读指令当由(一个或多个)计算机处理器运行时令(一个或多个)处理器执行所描述的动作。额外地或备选地,计算机可读指令中的至少一个有信号、载波或其他瞬态介质来承载。
应当理解,本文描述的实施例应用于诸如光电分量和康普顿散射分量的基材料和/或其他基材料对。
已经参考优选实施例描述了本发明。他人在阅读和理解前面的具体描述的情况下可以想到修改和替代。本文旨在将本发明解释为包括所有这样的修改和替代,只要它们落入权利要求书及其等价方案的范围内。

Claims (10)

1.一种成像系统,包括:
存储器(420),其包括用于以下的指令:处理谱CT投影数据以减轻所述谱CT投影数据的噪声或所述谱CT投影数据的噪声引发的偏差中的至少一个;
处理器(418),其被配置为运行所述指令,并且被配置为通过以下来处理所述谱CT投影数据:将所述谱CT投影数据的低能量投影数据和高能量投影数据对分类为相同极性对或相反极性对,将所述相反极性对变换为相同极性对,以创建经修改的谱CT投影数据,并且分解所述经修改的谱CT投影数据,以生成基材料;以及
重建器(434),其被配置为重建所述基材料,从而生成谱图像。
2.根据权利要求1所述的系统,其中,所述处理器被配置为通过将所述低能量和所述高能量对的值与所述低能量和所述高能量对的局部平均值进行比较来对所述低能量和所述高能量对进行分类。
3.根据权利要求2所述的系统,其中,所述处理器被配置为基于关于感兴趣的低能量和高能量对的预定空间邻域或预定时间邻域中的至少一项来确定所述局部平均值。
4.根据权利要求2至3中的任一项所述的系统,其中,所述处理器被配置为通过响应于高能量值大于所述局部平均值而从所述高能量值减去所述高能量值的标准偏差的两倍来将相反极性对变换为相同极性对。
5.根据权利要求2至3中的任一项所述的系统,其中,所述处理器被配置为通过响应于高能量值小于所述局部平均值而将所述高能量值加上所述高能量值的标准偏差的两倍来将相反极性对变换为相同极性对。
6.根据权利要求1-3中的任一项所述的系统,其中,所述处理器被配置为通过考虑所述谱CT投影数据的估计的噪声来修改适合用于无噪声谱投影数据的分解结果,以创建针对所述谱CT投影数据的经修改的分解结果。
7.根据权利要求6所述的系统,其中,所述处理器被配置为基于所述估计的噪声来确定卷积核。
8.根据权利要求1-3中的任一项所述的系统,其中,所述处理器被配置为通过将适当的代价函数最小化来对分解函数进行去卷积。
9.一种成像方法,包括:
通过以下来处理谱CT投影数据以降低噪声或噪声引发的偏差:将所述谱CT投影数据的低能量和高能量对分类为相同极性对或相反极性对,将所述相反极性对变换为相同极性对,创建经修改的谱CT投影数据;
分解所述经修改的谱CT投影数据,以生成基材料;并且
重建所述基材料,以生成谱图像。
10.一种编码有计算机可读指令的计算机可读存储介质,所述计算机可读指令当由处理器运行时令所述处理器:
通过以下处理谱CT投影数据以降低噪声或噪声引发的偏差:将所述谱CT投影数据的低能量和高能量对分类为相同极性对或相反极性对,将所述相反极性对变换为相同极性对,创建经修改的谱CT投影数据;分解所述经修改的谱CT投影数据以生成基材料;并且通过重建所述基材料来生成谱图像。
CN201580055482.3A 2014-10-13 2015-10-10 用于谱成像的噪声降低 Active CN106796733B (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
US201462062985P 2014-10-13 2014-10-13
US62/062,985 2014-10-13
PCT/IB2015/057755 WO2016059527A2 (en) 2014-10-13 2015-10-10 Spectral imaging

Publications (2)

Publication Number Publication Date
CN106796733A CN106796733A (zh) 2017-05-31
CN106796733B true CN106796733B (zh) 2021-05-25

Family

ID=54360497

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201580055482.3A Active CN106796733B (zh) 2014-10-13 2015-10-10 用于谱成像的噪声降低

Country Status (4)

Country Link
US (1) US10282870B2 (zh)
EP (1) EP3207524B1 (zh)
CN (1) CN106796733B (zh)
WO (1) WO2016059527A2 (zh)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10820865B2 (en) 2016-07-13 2020-11-03 Koninklijke Philiops N.V. Spectral computed tomography fingerprinting
WO2019232650A1 (en) * 2018-06-08 2019-12-12 Ka Imaging Inc. Method and system for determining virtual outputs for a multi-energy x-ray imaging apparatus
JP2020103571A (ja) * 2018-12-27 2020-07-09 キヤノンメディカルシステムズ株式会社 医用処理装置及びx線診断システム
CN112446930A (zh) * 2019-08-30 2021-03-05 通用电气精准医疗有限责任公司 用于图像重建的方法、成像系统以及存储有对应程序的介质
US11350895B2 (en) * 2019-11-29 2022-06-07 Wisconsin Alumni Research Foundation System and method for spectral computed tomography using single polychromatic x-ray spectrum acquisition
CN111134709B (zh) * 2020-01-17 2021-09-14 清华大学 一种多能量ct基材料分解方法
CN111896566B (zh) * 2020-07-20 2023-07-18 上海交通大学 一种增加同步辐射光源成像范围的装置及方法
EP3992912A1 (en) * 2020-11-03 2022-05-04 Koninklijke Philips N.V. Methods and systems for generating a spectral computed tomography image

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101647706A (zh) * 2008-08-13 2010-02-17 清华大学 高能双能ct系统的图象重建方法
CN102682434A (zh) * 2012-05-14 2012-09-19 西安电子科技大学 基于边缘先验和nsct域gsm的图像去噪方法
CN103559699A (zh) * 2013-11-18 2014-02-05 首都师范大学 一种基于投影估计的多能谱ct图像重建方法
CN103559729A (zh) * 2013-11-18 2014-02-05 首都师范大学 一种双能谱ct图像迭代重建方法

Family Cites Families (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4789930A (en) 1985-11-15 1988-12-06 Picker International, Inc. Energy dependent gain correction for radiation detection
US6052433A (en) 1995-12-29 2000-04-18 Advanced Optical Technologies, Inc. Apparatus and method for dual-energy x-ray imaging
WO2002028154A1 (en) 2000-09-29 2002-04-04 Analogic Corporation Method of and system for improving the signal to noise characteristics of images from a digital
US8263938B2 (en) 2004-03-01 2012-09-11 Varian Medical Systems, Inc. Dual energy radiation scanning of objects
US8199874B2 (en) 2009-12-11 2012-06-12 General Electric Company System and method of mitigating low signal data for dual energy CT
US8199875B2 (en) 2009-12-11 2012-06-12 General Electric Company System and method of acquiring multi-energy CT imaging data
US8401266B2 (en) * 2010-11-29 2013-03-19 General Electric Company Method and system for correlated noise suppression in dual energy imaging
JP6200428B2 (ja) * 2011-12-19 2017-09-20 コーニンクレッカ フィリップス エヌ ヴェKoninklijke Philips N.V. X線検出器、x線装置、プロセッサ、コンピュータプログラム、及びその関連方法
US20140169520A1 (en) * 2012-12-19 2014-06-19 Morpho Detection, Inc. Systems and methods for dual energy imaging
US8885910B2 (en) * 2012-12-28 2014-11-11 General Electric Company Systems and methods for X-ray imaging
US9700264B2 (en) * 2013-10-25 2017-07-11 The Johns Hopkins University Joint estimation of tissue types and linear attenuation coefficients for computed tomography

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101647706A (zh) * 2008-08-13 2010-02-17 清华大学 高能双能ct系统的图象重建方法
CN102682434A (zh) * 2012-05-14 2012-09-19 西安电子科技大学 基于边缘先验和nsct域gsm的图像去噪方法
CN103559699A (zh) * 2013-11-18 2014-02-05 首都师范大学 一种基于投影估计的多能谱ct图像重建方法
CN103559729A (zh) * 2013-11-18 2014-02-05 首都师范大学 一种双能谱ct图像迭代重建方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
双能CT图像重建算法研究;高洋;《中国优秀硕士学位论文全文数据库 信息科技辑》;20130315;正文第2-4章 *

Also Published As

Publication number Publication date
WO2016059527A3 (en) 2016-06-09
EP3207524A2 (en) 2017-08-23
WO2016059527A2 (en) 2016-04-21
US20170278278A1 (en) 2017-09-28
US10282870B2 (en) 2019-05-07
CN106796733A (zh) 2017-05-31
EP3207524B1 (en) 2019-10-02

Similar Documents

Publication Publication Date Title
CN106796733B (zh) 用于谱成像的噪声降低
Maier et al. Real‐time scatter estimation for medical CT using the deep scatter estimation: method and robustness analysis with respect to different anatomies, dose levels, tube voltages, and data truncation
Rührnschopf et al. A general framework and review of scatter correction methods in x‐ray cone‐beam computerized tomography. Part 1: scatter compensation approaches
US10147168B2 (en) Spectral CT
Schüller et al. Segmentation‐free empirical beam hardening correction for CT
Ritschl et al. Robust primary modulation‐based scatter estimation for cone‐beam CT
EP2702563B1 (en) Multi-energy imaging
JP7346429B2 (ja) スペクトルコンピュータ断層撮影(ct)スキャナによって生成される画質が向上したバーチャル非造影画像
JP6431068B2 (ja) 反相関フィルタによるスペクトルプロジェクションデータノイズ除去
CN107430779B (zh) 多能量(谱)图像数据处理
CN107106105B (zh) 用于物体的x射线成像装置
EP2671069A1 (en) Detection values processing apparatus
CN110574073B (zh) 能谱计算机断层摄影(ct)成像中的残余碘伪影的探测和/或校正
Humphries et al. Superiorized method for metal artifact reduction
Lee et al. Ultra-low-dose spectral CT based on a multi-level wavelet convolutional neural network
Mouton et al. A distance driven method for metal artefact reduction in computed tomography
CN111344742A (zh) 针对多种不同类型的投影数据每个体素具有一次几何计算的单个ct反投影器
Burton et al. Theoretical and experimental comparison of image signal and noise for dual-energy subtraction angiography and conventional x-ray angiography
Chukalina et al. A Way To Reduce The Artifacts Caused By Intensely Absorbing Areas In Computed Tomography.
Park et al. Pseudo-monochromatic Imaging in Industrial X-Ray Computed Tomography
Lu Multi-dimensional extension of the alternating minimization algorithm in x-ray computed tomography

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