CN102819826B - 降噪方法 - Google Patents

降噪方法 Download PDF

Info

Publication number
CN102819826B
CN102819826B CN201210185139.9A CN201210185139A CN102819826B CN 102819826 B CN102819826 B CN 102819826B CN 201210185139 A CN201210185139 A CN 201210185139A CN 102819826 B CN102819826 B CN 102819826B
Authority
CN
China
Prior art keywords
mentioned
noise
wave filter
variance
determination signal
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
CN201210185139.9A
Other languages
English (en)
Other versions
CN102819826A (zh
Inventor
杨智
A·察莫耶廷
邹宇
M·D·西尔弗
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.)
Canon Medical Systems Corp
Original Assignee
Toshiba Medical Systems Corp
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 Toshiba Medical Systems Corp filed Critical Toshiba Medical Systems Corp
Publication of CN102819826A publication Critical patent/CN102819826A/zh
Application granted granted Critical
Publication of CN102819826B publication Critical patent/CN102819826B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/70Denoising; Smoothing
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00Image enhancement or restoration
    • G06T5/20Image enhancement or restoration using local operators
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10072Tomographic images
    • G06T2207/10081Computed x-ray tomography [CT]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/20Special algorithmic details
    • G06T2207/20004Adaptive image processing
    • G06T2207/20008Globally adaptive

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Image Processing (AREA)
  • Noise Elimination (AREA)

Abstract

在降低计算机断层摄影(CT)等的测定信号的噪声的方法中,抑制滤波器效果对于滤波器参数的依赖性。降噪方法包含:a)根据规定的噪声模型,决定测定信号各自的噪声方差的步骤;b)根据噪声方差的P次方,产生针对各测定信号的滤波器的离散核的步骤,该针对各测定信号的滤波器的离散核表示滤波器对于各测定信号的频率响应;以及c)应用分别对应的上述离散核,对测定信号进行滤波的步骤。

Description

降噪方法
技术领域
本发明的实施方式涉及一种降噪方法。
背景技术
受到弱光引起的噪声及条纹(streak)的影响,X射线CT图像的质量有时候会严重变差。这个问题虽然可以通过增加X射线的剂量得到缓解,但从患者安全性方面考虑,这种方法在临床上是不被允许的。为了在安全辐射剂量等级的前提下达到在诊断上有用的图像质量,在以往寻求用于大幅度地降低噪声及条纹的优选的解决策略,近数十年来进行了各种各样的尝试。受到关于低辐射剂量的意识提升的影响,所提出的上述研究在最近重要性显著增加,备受关注。
为了改良结果,在以往的尝试中,利用自适应滤波器来代替固定滤波器。作为固定滤波器的若干示例有三角滤波器及双向滤波器。作为自适应滤波器的若干示例有自适应型高斯滤波器、自适应型裁剪均值滤波器。
在以往技术的自适应滤波器中,其滤波器参数需要适当地进行调整。依照经验,对每个图像分别实现此调整的最佳化。在其他的例示的滤波器中,高斯核的方差和数据的噪声方差相同。同样地,由于是对图像分别进行调整,因此利用上述以往技术的自适应滤波器执行的噪声和/或条纹的降低并不充分。
依然期望有一种有效的滤波器,该滤波器虽然使噪声及条纹的降低达到最大程度但其滤波器参数却不用复杂或特定的方法来进行调整。
发明内容
目的是在降低计算机断层摄影(CT)等的测定信号的噪声的方法中,抑制滤波器效果对于滤波器参数的依赖性。
本实施方式涉及的降噪方法包含:a)根据规定的噪声模型,决定测定信号各自的噪声方差的步骤;b)根据上述噪声方差的P次方,产生针对各上述测定信号的滤波器的离散核的步骤,该针对各上述测定信号的滤波器的离散核表示针对各上述测定信号的滤波器的频率响应;以及c)应用分别对应的上述离散核,对上述测定信号进行滤波的步骤。
附图说明
图1是表示基于本实施方式的多切片X射线CT装置或扫描器的结构的图。
图2是表示利用噪声模型预测的方差(y轴)和实测方差(x轴)的关系的图表。
图3是表示在图2的计数范围、及利用噪声模型预测的方差(y轴)和实测方差(x轴)的关系变得格外非线性的计数较少的延伸范围内的上述关系的图表。
图4是表示使用规定高斯滤波器的与维N、VarScale及VarPower相关的高斯滤波器的有效滤波器大小和计数的关系的图表。
图5是表示使用基于本实施方式的具有第1组参数值的规定滤波器的噪声和/或条纹减少处理部的一实施方式中的滤波后的平均计数和对数转换后方差的关联的图表。
图6是表示使用基于本实施方式的具有第2组参数值的规定滤波器的噪声和/或条纹减少处理部的一实施方式中的滤波后的平均计数和对数转换后方差的关联的图表。
图7是表示使用基于本实施方式的具有第3组参数值的规定滤波器的噪声和/或条纹减少处理部的一实施方式中的滤波后的平均计数和对数转换后方差的关联的图表。
图8是表示与使利用本实施方式重建CT图像前测定的数据的噪声和/或条纹大幅度地降低的一个例示性处理相关的步骤的流程图。
图9A是作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善的肩部图像。
图9B是作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善的肩部图像。
图9C是作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善的肩部图像。
图10A是作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善的其他肩部图像。
图10B是作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善的其他肩部图像。
图10C是作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善的其他肩部图像。
图11A是基于图10A、10B及10C的同一数据而以肺窗(lung window)表示的缩放后的重建的图像。
图11B是基于图10A、10B、及10C的同一数据而以肺窗表示的缩放后的重建的图像。
图11C是基于图10A、10B、及10C的同一数据而以肺窗表示的缩放后的重建的图像。
具体实施方式
以下,一边参照附图,一边说明本实施方式。
在此,参照附图,相同参照编号表示所有附图中对应的构造,尤其参照图1,图中表示包含台架100和其他处理部或单元的、基于本实施方式的多切片X射线CT装置或扫描器的一实施方式。如前视图所示,台架100还包含X射线管101、环状框架102、多列或二维排列型X射线检测器103。X射线管101及X射线检测器103是在直径方向上隔着被摄体(subject)S而安装于环状框架102,环状框架102绕着轴RA旋转。被摄体S沿着轴RA在图示的页面内外移动期间,旋转单元107使框架102以0.4秒/旋转等高速旋转。
多切片X射线CT装置还包含电流调整器113和高电压产生装置109,该高电压产生装置109对X射线管101施加管电压,使X射线管101生成X射线。在一实施方式中,高电压产生装置109被安装于框架102。X射线向被摄体S辐射,用圆表示被摄体S的截面积。X射线检测器103为了检测透过被摄体S的辐射出的X射线,隔着被摄体S而位于X射线管101的相反侧。
另外,参照图1,X射线CT装置或扫描器还包含用来处理来自X射线检测器103的被检测出的信号的其他处理部。数据收集电路或数据收集系统(DAS)104针对各通道将X射线检测器103输出的信号转换成电压信号,对其进行放大后再进一步转换成数字信号。X射线检测器103及DAS104构成为对每1旋转的规定的所有投影数(TPPR:total number ofprojections per rotation)进行操作。
上述数据通过非接触式数据发送器105,被发送给台架100外部的控制中心内所容纳的前处理处理部106。前处理处理部106执行对原始数据的灵敏度校正等特定校正。接着,存储处理部112在即将重建处理前的阶段,存储还被称为投影数据的作为结果而生成的数据。存储处理部112和重建处理部114、显示处理部116、输入处理部115、及扫描计划帮助装置200一起,经由数据/控制总线而连接于系统控制器110。扫描计划帮助装置200包含为制定扫描计划而帮助摄像技师的功能。
本实施方式的一实施方式还包含各种软件模块和硬件构成要素,使用规定滤波器使计算机断层摄影图像的条纹和/或噪声大幅度地降低。根据本实施方式的一方式,CT装置的降噪处理部117有利地执行噪声和/或条纹的减少。在一实施方式中,降噪处理部117经由数据/控制总线,动作性地连接于其他软件模块、和/或存储处理部112、重建处理部114、显示处理部116、及输入处理部115等系统构成要素。关于这一点,在本实施方式的其他实施方式中,未必由降噪处理部117单独地执行降噪功能和/或与此关联的任务。而且,降噪处理部117在本实施方式的代替实施方式中,通过任意选择而变成重建处理部114等其他处理部的一部分。
一般来说,CT内的投影数据可以在规定的对数转换处理之后利用。在该对数转换处理中,将由于扫描对象而衰减的被测定后的X射线强度信号转换成线积分数据。然后,根据该线积分数据,利用已知的数学反演(mathematical inversion)方法来重建CT图像。在基于本实施方式的噪声/条纹减少系统的例示性一实施方式中,降噪处理部117转换投影数据而使其返回至原始的X射线强度数据或光子计数测定值。这种情况下,降噪处理部117在转换步骤中需要与系统校正处理相关的一些信息。或者,降噪处理部117对测定出的X射线强度信号直接存取。
降噪处理部117基于X射线强度信号或光子计数,决定对数转换后的数据的噪声方差(V)。在本实施方式中,以对数转换处理后使噪声均匀化的方式来决定噪声滤波器参数。
为了理解对数转换对于测定出的数据的效果,关于对数转换前及后的方差调用噪声模型。对数转换前(before-log)的噪声方差VarBL是利用由式(1)定义的对数转换前噪声模型来推测的。
VarBL=Ve+WI (1)
VarBL是对数转换前的全部噪声方差,Ve是电子噪声的方差(electronic noisevariance),I是平均计数。W是检测器的增益,这是通道、片段、数据收集系统(DAS)、和/或校准的函数。另一方面,对数转换后的噪声方差VarAL是利用由式(2)定义的对数转换后的噪声模型来推测的。
上述两个式子均公开于「Adaptive streak artifact reduction in CTresulting from excessive X-ray photon noise」、Jiang Hsieh(GE)、Med.Phys.25(11)、2139~47、1998年。
噪声方差被决定之后,自适应型高斯滤波器被应用于X射线强度数据。自适应型高斯滤波器(G)的例示性一方式是由式(3)定义。
式中,标准偏差(SD)σG决定有效的核大小,x是核内的任意位置和中央位置的距离。自适应型高斯滤波的良好应用的关键在于,决定噪声SDσG的值而作为数据或噪声的局部特性的函数。在基于本实施方式的噪声和/或条纹减少处理部的实施方式中,核的方差VG(或σ2 G)是由式(4)定义的,自适应地决定为噪声方差VarAL的函数f。
VG=f(VarAL) (4)
关于上述一般性的式子,进一步说明基于本实施方式的噪声和/或条纹减少处理部的特定的实现方式。式(4)的例示性一函数在以下由式(5)定义,但此例并不用来限定本实施方式的其他实现方式。
VG=K(VarAL)P (5)
式中,K是用于控制滤波强度的参数。在一实施方式中,K设定为1。在其他实施方式中,参数K根据应用种类而决定,参数K的值是任意的绝对值。而且,其他参数P是其典型范围包含为0至2之间的乘方(power)或指数参数。在一实施方式中,P的值为0.5。在用于适当的应用例的其他实施方式中,P=1等这样的其他P值也可以任意选择。
图2表示与式(2)的对数转换后噪声模型相关的光子计数和测定出的方差的关联。由线表示的对数转换后噪声模型在计数为44以上的情况下,密切关联真的对数转换后方差值而表示(应留意Ve约为30~40)。另一方面,如图3所示,若计数减少为低于44(即,减少至电子噪声的下限),则测定出的方差关于对数转换后噪声模型而变得格外的非线性。也就是说,对数转换后噪声模型在计数较少时并不准确。较少计数时方差会变大,因此在利用对数转换后数据决定噪声方差的情况下,较少计数时需要减少方差。方差用表示为基于任意选择的较大值的分贝(dB)的相对方差值来表示。
在基于本实施方式的噪声和/或条纹减少处理部的实施方式中,噪声方差是基于对数转换后的投影数据而推测的,滤波器基于对数转换后方差而作成的。之后,上述滤波器应用于对数转换前测定出的数据或X射线强度数据。最后,转换上述滤波后的X射线强度数据,返回至投影数据(即,线积分)区域,然后通过本实施方式重建CT图像。因此,以下的例示性滤波器作成方法是使用高斯滤波器来表示。
在第1滤波器作成方法中,计数I0时的噪声为能够允许的,所需的方差等级Var0假定为由式(6)定义。
Var0=(Ve0+W0I0)/I0 2 (6)
式中,Ve0及W0分别是电子噪声的方差的平均值及检测器的增益。
然后,利用低通滤波来降低方差。也就是说,对数转换后的方差VarAL通过一定的比率VRR而降低为由式(7)定义的滤波后的方差VarF
VarF=VarAL VRR (7)
式中,在滤波后的方差VarF为所需的方差等级Var0以下(VarF≤Var0)时,所需的缩小率VRR由式(8)决定。
VRR=Var0/VarAL=Var0I2/VarBL (8)
式中,VarBL是对数转换前全部噪声方差,VarAL是对数转换后全部噪声方差,I是平均计数。
上述低通滤波器假定为通过规定的一组系数{ck}所提供。在此,k=1…N、ck>0、及∑kck=1。关于上述系数,所需的缩小率VRR是通过式(9)定义的。
VRR=∑kck 2. (9)
然后,系数{ck}任意选择,是N维高斯滤波器,其滤波器方差VG是通过式(10)定义。
kck 2=1/(4πVG)N/2. (10)
根据式(8)~(10),高斯滤波器的方差VG在维数为N的情况下通过式(11A)及(11B)推测。
1/(4πVG)N/2=Var0I2/VarBL. (11A)
VG=1/(4π)(VarBL/(Var0Mean2))2/N (11B)
高斯滤波器的方差VG是基于式(11B),简化表示为利用了2个变量VarScale及VarPower的式(11C)。实际上,式(11C)等同于式(5),这里式(5)的参数K及P分别相当于式(11C)的VarScale及VarPower。在本申请中,式(5)的参数K及P和式(11C)的VarScale及VarPower可以交换。根据第1滤波器作成方法,式(11C)定义高斯滤波器的方差VG
VG=VarScale VarAL VarPower (11C)
式中,VarAL已经由式(2)定义。最后,高斯滤波器的方差VG是利用2个变量VarScale(尺度方差)及VarPower而根据对数转换后噪声方差VarAL来决定。而且,变量VarScale是通过与以前由式(6)定义的Var0相关的式(12)来定义。根据第1滤波器作成方法,式(12)定义变量VarScale。
VarScale=1/(4π)(1/Var0)2/N (12)
另一方面,另一个变量VarPower是由式(13)定义,其值在基于本实施方式的噪声和/或条纹减少处理部的实施方式中相对于2D滤波器、3D滤波器、及4D滤波器容易被决定。根据第1滤波器作成方法,式(13)定义变量VarPower。
VarPower=2/N (13)
=1用于2D高斯滤波器
=2/3=0.67用于3D高斯滤波器
=1/2用于4D高斯滤波器
接着,参照图4,图表表示使用规定高斯滤波器的与维数为N、VarScale、及VarPower相关的高斯滤波器的有效滤波器大小和计数的关系。在维数为N的1至4的范围内的特定值下,VarPower的值为2至0.5的范围内。另一方面,VarScale的值具有406.0075至0.6726这样更广的范围。如图所示,在N=4时,方差在计数最广的范围中是固定的。另一方面,在N=1时,方差在计数最窄的范围中是固定的。为了实施本实施方式,在其他实施方式中,VarPower及VarScale的范围并不限定为上述公开的值的范围。实际上,根据本实施方式的特定的其他方式,VarPower在1至10000的范围任意选择,VarPower在0.2至1的范围任意选择。
如上上述,在基于本实施方式的噪声和/或条纹减少处理部的实施方式中,噪声方差是基于对数转换后的投影数据而推测的,滤波器是基于对数转换后方差而作成的。然后,上述滤波器应用于对数转换前测定出的数据或X射线强度数据。最后,转换上述滤波后的X射线强度数据,返回至投影数据(即,线积分)区域,然后通过本实施方式重建CT图像。上述滤波器应用的结果为,对数转换后的投影数据实质上具有均匀的方差。
一般来说,测定后的数据的全部方差受泊松噪声及高斯噪声影响。来自检测器的计数数据理想情况下是依照泊松分布,实际数据和数据收集系统(DAS)电路引起的高斯分布的电子噪声混合。在重建后的图像中,利用较少计数数据生成条纹和无法容忍的噪声,电子噪声在较少的计数下已经无法忽视。因此,准确的噪声模型需要考虑泊松噪声和高斯噪声这两者。
在式(11A)、(11B)、及(11C)中,假定无限的高斯核。实际上,高斯核被限定为5×7等有限的掩膜。在降噪处理部117等的噪声和/或条纹减少处理部的另一实现方式中,式(14)定义用于决定针对测定出的各数据值的离散性高斯滤波的一个装置。也就是说,对于特定检测器元件(i)的滤波器的离散性的核是通过下式定义。
式中,Δxi是一维(1D)检测器中的第i个像素离规定的基准像素(referencepixel)i0的距离,Vi0是基准像素i0的滤波器的频率响应。VR是滤波器的参数。
通过相同标记,降噪处理部117决定针对二维(2D)检测器测定出的各数据值的离散性高斯滤波。也就是说,滤波器对于特定的检测器元件(i,j)的滤波器的离散性的核通过式(15)定义。
式中,Δxi,j是2D检测器中的(第i个,第j个)像素离规定的基准像素(i0,j0)的距离,Vi0,j0是基准像素(i0,j0)的滤波器的频率响应。VR是滤波器的参数。
和以往技术的尝试相反,在基于本实施方式的噪声和/或条纹减少处理部的实施方式中,基于对数转换后的投影数据的噪声特性而推测噪声方差,滤波器基于对数转换推测后的后方差而作成。上述对数转换推测后的后方差是
上述I是上述测定信号,上述Ve是已知的电子噪声的值。上述滤波器是高通型高斯滤波器。然后,上述滤波器应用于对数转换前原始的测定出的数据或X射线强度数据。最后,转换上述滤波后的X射线强度数据,返回至投影数据(即,线积分)区域,通过本实施方式而重建CT图像。
在特定的实施方式中,为了实现期望的降噪效果,使噪声和/或条纹实质上最小的上述步骤以任意选择的方式而重复执行多次。在这些实施方式中,参数p及K为任意选择,为了所需的解决策略而变化。而且,参数p及K为任意选择,针对每次重复被执行的情况而发生变化。
为了决定上述噪声和/或条纹减少处理中的式(5)的参数K及P或式(11C)的VarScale及VarPower的临床上有用的值,最佳的噪声滤波器作为逐次近似重建(IR:iterative reconstruction)噪声模型滤波被作成。一般来说,纯粹的原始数据滤波的作用是除去电子噪声或弱光的影响。对数转换的结果为,值较低的原始数据变得无法信赖,因此对这些数据分配了较低的统计上的权重。因此,纯粹的原始数据滤波和统计上的权重相互补充。由于这些理由,在基于本实施方式的噪声和/或条纹减少处理部的实施方式中,不需要基于total variation(TV:全变分)或自适应型加权各向异性扩散法(AWAD:AdaptiveWeighted Anisotropic Diffusion)的图像的规格化(regularization)等强力的数据滤波。总而言之,IR噪声模型滤波通过使用较弱的滤波参数而得到改良,维持低辐射剂量数据的空间分辨率。
通过以下示例说明参数的其他变动。接着,参照图5,图表表示基于本实施方式的使用规定滤波器的噪声和/或条纹减少处理部的一实施方式中的平均计数和对数转换后方差的关联。假定无限的高斯核,将式(5)的参数K或式(11C)的VarScale调整为1、10、及100,式(5)的参数P或式(11C)的VarPower固定为1。在没有滤波器的情况下,好像在平均计数下的对数转换后方差的噪声没有均匀化。k=1时,在小于30的较小的低平均计数区域内的对数转换后方差中,好像有一些使噪声均匀化的效果。k=10时,在小于160的较广的平均计数区域内的对数转换后方差中,好像有较大地使噪声均匀化的效果。最后,k=100时,在小于200的所有平均计数区域内的对数转换后方差中,好像有使噪声较大地均匀化的效果。方差是用表示为基于任意选择的较大值的分贝(dB)的相对方差值来表示的。
接着,参照图6,图表表示基于本实施方式的使用规定滤波器的噪声和/或条纹减少处理部的一实施方式中的平均计数和对数转换后方差的关联。假定无限的高斯核,将式(5)的参数P或式(11C)的VarPower调整为1、0.75、及0.5,式(5)的参数K或式(11C)的VarScale固定为1。在没有滤波器的情况下,好像在平均计数下的对数转换后方差的噪声没有均匀化。P=1时,在小于30的较小的低平均计数区域内的对数转换后方差中,好像有一些使噪声均匀化的效果。P=0.75时,小于50的略广的平均计数区域内的对数转换后方差中,也好像有一些使噪声均匀化的效果,但在此均匀化效果中并未实现完全平坦且均匀化的区域。最后,P=0.5时,小于200的相当广的平均计数区域内的对数转换后方差中,好像具有较大地使噪声均匀化的效果。在60至200之间的平均计数时,好像使噪声均匀化的效果显著。方差是以表示为基于任意选择的较大值的分贝(dB)的相对方差值来表示。
接着,参照图7,图表表示基于本实施方式的使用规定的滤波器的噪声和/或条纹减少处理部的一实施方式中的平均计数和对数转换后方差的关联。假定无限的高斯核,将式(5)的参数K或式(11C)的VarScale调整为1、1.5及2,式(5)的参数P或式(11C)的VarPower固定为0.5。在没有滤波器的情况下,好像在平均计数下的对数转换后方差的噪声没有均匀化。K=1时,在小于180的较小的低平均计数区域内的对数转换后方差中,好像有一些使噪声均匀化的效果。K=1.5时,小于380的略广的平均计数区域内的对数转换后方差中,好像有一些使噪声均匀化的效果。最后,K=2.0时,小于800的较广的平均计数区域内的对数转换后方差中,好像有使噪声均匀化的效果。方差是以表示为基于任意选择出的较大值的分贝(dB)的相对方差值来表示。
接下来,参照图8,表示与使通过本实施方式重建CT图像前测定出的数据的噪声和/或条纹大幅度地降低的一个例示性处理相关的步骤。一般来说,CT内的投影数据可以在规定的对数转换处理之后利用。在该对数转换处理中,将由于扫描的对象而衰减的测定出的X射线强度信号转换为线积分数据。然后,根据该线积分数据,通过已知的数学反演方法重建CT图像。
在基于本实施方式的噪声/条纹减少处理的例示性一实施方式中,以下的步骤是通过软件和硬件的规定组合而执行。该处理的实现方式并不限定于哪种特定的软件模块或硬件模块。
在步骤S20中,通过降噪处理,基于X射线强度信号或光子计数而决定测定数据的噪声方差(V)。对该噪声方差进行计算,以使得在对数转换处理后使噪声均匀化。对数转换后噪声方差VarAL是由上述式(2)定义的噪声模型来决定的,该式(2)是根据式1将测定出的X射线强度或光子计数I和数据收集系统等的电子噪声Ve加入考虑而得到的。
在步骤S20中决定对数转换后噪声方差VarAL之后,在步骤S30中自适应型高斯滤波器应用于来自步骤S10的X射线强度数据。自适应型高斯滤波器(G)的例示性一方式是通过上述式(3)定义。自适应型高斯滤波器依赖于决定有效的核大小的标准偏差σ、以及核内的任意位置和中央位置之间的距离x。自适应型高斯滤波作为数据或噪声的局部特性的函数而正常应用于噪声SDσ的值。在通过本实施方式使噪声和/或条纹大幅度降低的例示性处理中,核的方差VG(或σ2 G)是通过上述式(4)定义的作为噪声方差VarAL的函数f而自适应地被决定。噪声和/或条纹减少处理的特定的实现方式需要用于控制滤波强度的参数K。在一实施方式中,K设定为1。在其他实施方式中,参数K是由应用种类决定,参数K的值为任意选择的位数。而且,在噪声和/或条纹减少处理的特定的实现方式中,还需要作为包含其典型范围为0至2之间的乘方或指数参数的第2参数P。在一处理中,指数P的值为0.5。在用于适当的应用例的其他处理中,P=1等的其他p值也是能够任意选择地利用。上述参数K及P分别在通过本实施方式而使噪声和条纹大幅度地降低的例示性处理中,可以替换为上述式(11C)的VarScale及VarPower。
在步骤S30的另一实现方式的详细内容中,用于决定针对测定出的各数据值的离散性高斯滤波的一个装置,是由上述式(14)对于特定的一维(1D)检测器元件(i)来定义的。同样地,用于决定针对测定出的各数据值的离散性高斯滤波的一个装置,是由上述式(15)对于特定的二维(2D)检测器元件(i,j)来定义的。步骤S30中,离散性高斯核是针对测定出的各数据而作成,并应用于这些数据。和以往技术的尝试相反,在通过本实施方式而使噪声和/或条纹大幅度地降低的上述例示性处理中,步骤S20中基于对数转换后的投影数据的噪声特性而推测噪声方差,基于对数转换后推测出的方差而作成滤波器。然后,在步骤S30中,上述滤波器被应用于对数转换前原始的测定出的数据或X射线强度数据。
以任意选择的方式,在基于本实施方式的例示性处理中,重复地反复步骤S20及S30。为了以任意选择的方式进行重复,在步骤S40中决定例示性处理是否返回到步骤S20。步骤S40中,在决定为滤波器应重复地被应用的情况下,例示性处理从步骤S20开始反复。另一方面,步骤S40中,在决定为滤波器不重复应用或重复结束的情况下,例示性处理进入到步骤S50中。最后,在步骤S50中,转换上述滤波后的X射线强度数据,返回至投影数据(即,线积分)区域内,通过本实施方式重建CT图像。
接下来,参照图9A、9B及9C,肩部图像作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善。DAS电子噪声的方差Ve是由不旋转X射线管而取得的实际扫描数据所决定。唯一的滤波器参数为K,在图9A、9B、及9C的所有测试案例中均设定为1。而且,图9A、9B、及9C所示的图像是在120kv、195mAs下以+69.5mm/床旋转速度(rotation couch speed)利用160列检测器所取得。图9A表示未利用基于本实施方式的自适应型高斯滤波器对数据区域进行滤波的重建图像。图9B表示利用基于本实施方式的自适应型高斯滤波器对数据区域进行滤波后的重建图像。图9B表示条纹及噪声通过基于本实施方式的自适应型高斯滤波器而得到较大抑制的情况。换句话说,图9A中,由于条纹及噪声导致临床上关联的信息被遮盖和/或不明确。另一方面,图9B中,描绘出了上述临床上关联的信息。最后,图9C表示构造信息基本上没有损失地示出的差分图像。
接下来,参照图10A、10B、及10C,其他肩部图像作为基于本实施方式的自适应型高斯滤波器的应用结果而示出在临床上有意义的改善。质量改良在所有案例中均大体一致。DAS电子噪声的方差Ve是根据不旋转X射线管所取得的实际扫描数据而决定。唯一的滤波器参数为K,在10A、10B、及10C的所有测试案例中均设定为1。而且,这些图像是在120kv、30mAs下以+47.5mm/床旋转速度而利用64列检测器来取得。图10A表示未通过基于本实施方式的自适应型高斯滤波器对数据区域进行滤波的重建图像。图10B表示通过基于本实施方式的自适应型高斯滤波器对数据区域进行滤波的重建图像。图10B表示条纹及噪声通过基于本实施方式的自适应型高斯滤波器得到大幅度地抑制的情况。换句话说,图10A中,由于条纹及噪声导致临床上关联的信息被遮盖和/或不明确。另一方面,图10B中,描绘出上述临床上关联的信息。最后,图10C表示构造信息基本上没有损失地被示出的差分图像。
接下来,参照图11A、11B、及11C,图像是基于图10A、10B、及10C的同一数据表示但是以肺窗表示的缩放后的重建图像。可明确地知道,图11B中条纹在肺组织中可不使分辨率明显变差地很好地被除去。图11B中在下部肺区域清晰地识别出云状污染物,而该污染物在图11A中并不明了。在肋骨处识别出若干模糊效果。这是由作为骨强度增加的最下端的极低的观测窗等级(-700HU)所引起的。在作为差分图像的图11B中,除了中间略偏上的部分处的血管较暗的征兆,并无明确的构造信息。
虽然对本发明的若干实施方式进行了说明,但这些实施方式是作为示例而提示的,并不试图限定本发明的范围。这些实施方式可以通过其他各种方式来实施,在不脱离发明主旨的范围内可以进行各种省略、置换、变更。这些实施方式及其变形包含于发明范围及主旨内,同样包含于权利要求记载的发明及其均等范围内。

Claims (21)

1.一种降噪方法,包含:
a)根据规定的噪声模型,决定针对测定信号各自的对数转换处理后的噪声方差的步骤;
b)根据上述噪声方差的P次方,产生针对各上述测定信号的滤波器的离散核的步骤,该针对各上述测定信号的滤波器的离散核表示滤波器针对各上述测定信号的频率响应;以及
c)应用分别对应的上述离散核,对上述测定信号进行滤波的步骤。
2.根据权利要求1所述的降噪方法,其中
上述规定的噪声模型是通过VarAL=(Ve+WI)/I2定义的,
上述I是上述测定信号,上述Ve是已知的电子噪声的值,上述W是信号增益。
3.根据权利要求2所述的降噪方法,其中
针对各上述测定信号的上述滤波器的上述离散核是通过对上述噪声方差的P次方乘以常数K而生成的。
4.根据权利要求3所述的降噪方法,其中
上述P是根据包含线性扩散式的尺度空间的导出而被决定的。
5.根据权利要求3所述的降噪方法,其中
上述P具有0.2至1的范围内的值,上述常数K具有1至10000的范围内的值。
6.根据权利要求5所述的降噪方法,其中
相对于上述测定信号,上述P为0.5,上述K为1。
7.根据权利要求1所述的降噪方法,其中
针对用于产生上述测定信号的检测器内的各检测器元件(i)的上述滤波器的上述离散核是由
定义的,上述Δxi是1维图像中的第i个像素和规定的基准像素i0之间的距离,上述是上述基准像素i0中的上述滤波器的上述频率响应,上述VR是上述滤波器的参数。
8.根据权利要求7所述的降噪方法,其中
上述参数VR等于由VarAL=(Ve+WI)/I2决定的对数转换后的方差,上述I是上述测定信号,上述Ve是已知的电子噪声的值。
9.根据权利要求1所述的降噪方法,其中
针对用于产生上述测定信号的检测器内的各检测器元件(i,j)的上述滤波器的上述离散核是由
定义的,上述Δxi,j是2维图像中的(i,j)像素和规定的基准像素(i0,j0)之间的距离,其中,i和j是表示像素中的坐标,上述是上述基准像素(i0,j0)中的上述滤波器的频率响应,上述VR是上述滤波器的参数。
10.根据权利要求9所述的降噪方法,其中
上述参数VR等于由VarAL=(Ve+WI)/I2决定的对数转换后的方差,上述I是上述测定信号,上述Ve是已知的电子噪声的值。
11.根据权利要求1所述的降噪方法,其中
重复执行上述步骤a)、b)及c)。
12.根据权利要求1所述的降噪方法,其中
上述滤波器是低通型的高斯滤波器。
13.根据权利要求1所述的降噪方法,其中
上述滤波器是为了使上述测定信号中的对数转换处理后方差下的噪声均匀化而设置的。
14.根据权利要求2所述的降噪方法,其中
上述滤波器是高斯滤波器,
针对各上述测定信号的上述高斯滤波器的方差VG是根据对数转换后的噪声方差VarAL的规定的VarPower,通过乘以作为各上述测定信号中的上述滤波器的上述频率响应而使用的规定常数VarScale而生成,所述VarPower是利用2/N来定义的,N是上述高斯滤波器的维数。
15.根据权利要求14所述的降噪方法,其中
上述VarScale是由VarScale=1/(4π)(1/Var0)2/N定义的,Var0为(Ve0+W0I0)/I0 2,式中Ve0及W0分别是电子噪声的方差的平均值及检测器的增益,I0是能够容许的噪声计数。
16.一种降噪方法,包含:
a)根据规定的噪声模型,决定多个测定信号各自中的对数转换推测后的后方差的步骤;
b)根据上述对数转换推测后的后方差,决定用于使各上述测定信号中的对数转换推测后的后方差下的噪声大体均匀化的滤波器的离散核的步骤;以及
c)应用分别对应的上述离散核,对上述测定信号进行滤波的步骤。
17.根据权利要求16所述的降噪方法,其中
上述对数转换推测后的后方差是
上述I是上述测定信号,上述Ve是已知的电子噪声的值。
18.根据权利要求17所述的降噪方法,其中
针对用于产生上述测定信号的检测器内的各检测器元件(i)的上述滤波器的上述离散核是由
定义的,上述Δxi是1维图像中的第i个像素和规定的基准像素i0之间的距离,上述是上述基准像素i0中的上述滤波器的频率响应,上述VR是上述滤波器的参数。
19.根据权利要求17所述的降噪方法,其中
针对用于产生上述测定信号的检测器内的各检测器元件(i,j)的上述滤波器的上述离散核是由
定义的,上述Δxi,j是2维图像中的(i,j)像素和规定的基准像素(i0,j0)之间的距离,其中,i和j是表示像素中的坐标,上述 是上述基准像素(i0,j0)中的上述滤波器的频率响应,上述VR是上述滤波器的参数。
20.根据权利要求16所述的降噪方法,其中
重复执行上述步骤a)、b)及c)。
21.根据权利要求16所述的降噪方法,其中
上述滤波器是高通型高斯滤波器。
CN201210185139.9A 2011-06-06 2012-06-06 降噪方法 Active CN102819826B (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
US13/154,124 US8538114B2 (en) 2011-06-06 2011-06-06 Method and system utilizing parameter-less filter for substantially reducing streak and or noise in computer tomography (CT) images
US13/154,124 2011-06-06

Publications (2)

Publication Number Publication Date
CN102819826A CN102819826A (zh) 2012-12-12
CN102819826B true CN102819826B (zh) 2017-07-11

Family

ID=46851793

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210185139.9A Active CN102819826B (zh) 2011-06-06 2012-06-06 降噪方法

Country Status (4)

Country Link
US (2) US8538114B2 (zh)
EP (1) EP2533197B1 (zh)
JP (1) JP6021448B2 (zh)
CN (1) CN102819826B (zh)

Families Citing this family (23)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP2583241A1 (en) * 2010-06-21 2013-04-24 Koninklijke Philips Electronics N.V. Method and system for noise reduction in low dose computed tomography
EA017302B1 (ru) * 2011-10-07 2012-11-30 Закрытое Акционерное Общество "Импульс" Способ подавления шума серий цифровых рентгенограмм
US9031297B2 (en) * 2013-02-01 2015-05-12 Kabushiki Kaisha Toshiba Alternative noise map estimation methods for CT images
CA2906973C (en) * 2013-04-04 2020-10-27 Illinois Tool Works Inc. Helical computed tomography
JP6492005B2 (ja) * 2013-04-08 2019-03-27 株式会社日立製作所 X線ct装置、再構成演算装置、及び再構成演算方法
EP3039650B1 (en) * 2013-08-30 2020-07-08 Koninklijke Philips N.V. Spectral projection data de-noising with anti-correlation filter
CN105682559B (zh) * 2013-09-25 2019-07-09 瓦里安医疗系统公司 用于估计散射的方法和系统
KR102144994B1 (ko) 2013-09-30 2020-08-14 삼성전자주식회사 영상의 노이즈를 저감하는 방법 및 이를 이용한 영상 처리 장치
US10192328B2 (en) * 2013-10-24 2019-01-29 Toshiba Medical Systems Corporation Method for statistical weights design in iterative reconstruction algorithms
JP6386762B2 (ja) * 2014-03-24 2018-09-05 キヤノンメディカルシステムズ株式会社 放射線診断装置、放射線検出装置、及び放射線検出データ処理方法
CN105374012B (zh) * 2014-08-27 2018-11-27 通用电气公司 用于消除由性能差异的探测器单元所致的条状伪影的方法
JPWO2016132880A1 (ja) * 2015-02-16 2017-11-30 株式会社日立製作所 演算装置、x線ct装置、及び画像再構成方法
US9805448B2 (en) * 2015-02-26 2017-10-31 Toshiba Medical Systems Corporation Method and apparatus for computed tomography using asymetric filter for volume half reconstruction
CN106157275A (zh) * 2015-04-21 2016-11-23 杭州与盟医疗技术有限公司 基于热扩散方程的多尺度ct图像去噪算法
CN105049846B (zh) * 2015-08-14 2019-05-21 广东中星微电子有限公司 图像和视频编解码的方法和设备
GB2542774B (en) 2015-09-25 2021-08-18 Smiths Heimann Sas Denoising and/or zooming of inspection images
US10282831B2 (en) * 2015-12-28 2019-05-07 Novatek Microelectronics Corp. Method and apparatus for motion compensated noise reduction
US10213176B2 (en) 2016-04-27 2019-02-26 Toshiba Medical Systems Corporation Apparatus and method for hybrid pre-log and post-log iterative image reconstruction for computed tomography
CN113873237B (zh) * 2016-12-01 2023-12-29 谷歌有限责任公司 用于恢复由重构产生的劣化帧的劣化图块的方法和装置
US11861765B2 (en) 2018-07-09 2024-01-02 Koninklijke Philips N.V. Imaging system detector clipping-induced bias correction
US11937009B2 (en) * 2019-02-07 2024-03-19 Mitsubishi Electric Corporation Infrared imaging device and non-transitory computer-readable storage medium
CN109887051B (zh) * 2019-02-28 2023-07-25 沈阳开普医疗影像技术有限公司 Ct图像重建反投影过程中的线性插值方法
CN117370737B (zh) * 2023-12-08 2024-02-06 成都信息工程大学 一种基于自适应高斯滤波器的非稳态非高斯噪声去除方法

Family Cites Families (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5276614A (en) * 1989-11-17 1994-01-04 Picker International, Inc. Dynamic bandwidth reconstruction
DE69331719T2 (de) * 1992-06-19 2002-10-24 Agfa Gevaert Nv Verfahren und Vorrichtung zur Geräuschunterdrückung
US5764307A (en) * 1995-07-24 1998-06-09 Motorola, Inc. Method and apparatus for spatially adaptive filtering for video encoding
US6963670B2 (en) 2001-11-21 2005-11-08 Ge Medical Systems Global Technology Company, Llc CT dose reduction filter with a computationally efficient implementation
US6493416B1 (en) * 2001-11-21 2002-12-10 Ge Medical Systems Global Technology Company, Llc Method and apparatus for noise reduction in computed tomographic systems
US7599579B2 (en) 2002-07-11 2009-10-06 Ge Medical Systems Global Technology Company, Llc Interpolated image filtering method and apparatus
US20050201605A1 (en) 2004-03-11 2005-09-15 Jianying Li Methods and apparatus for CT smoothing to reduce artifacts
CN100563570C (zh) * 2004-07-07 2009-12-02 皇家飞利浦电子股份有限公司 心脏锥面光束ct重建中条纹伪影的减少
JP2006279162A (ja) * 2005-03-28 2006-10-12 Fuji Xerox Co Ltd 画像入力方法及びこれを用いた画像入力装置
US7672421B2 (en) 2005-10-12 2010-03-02 Siemens Medical Solutions Usa, Inc. Reduction of streak artifacts in low dose CT imaging through multi image compounding
EP1892666A1 (en) 2006-08-16 2008-02-27 Toyota Motor Europe NV A method, an apparatus and a computer-readable medium for processing an image dataset
DE102006046191B4 (de) * 2006-09-29 2017-02-02 Siemens Healthcare Gmbh Streustrahlungskorrektur in der Radiographie und Computertomographie mit Flächendetektoren
JP4480760B2 (ja) * 2007-12-29 2010-06-16 株式会社モルフォ 画像データ処理方法および画像処理装置
JP5028528B2 (ja) * 2008-09-30 2012-09-19 株式会社日立メディコ X線ct装置
US8965078B2 (en) * 2009-02-20 2015-02-24 Mayo Foundation For Medical Education And Research Projection-space denoising with bilateral filtering in computed tomography
RU2586968C2 (ru) 2010-10-27 2016-06-10 Конинклейке Филипс Электроникс Н.В. Уменьшение уровня шума в низкодозной компьютерной томографии

Also Published As

Publication number Publication date
EP2533197B1 (en) 2020-10-28
US8538114B2 (en) 2013-09-17
US20120308104A1 (en) 2012-12-06
EP2533197A2 (en) 2012-12-12
US8965144B2 (en) 2015-02-24
JP6021448B2 (ja) 2016-11-09
CN102819826A (zh) 2012-12-12
EP2533197A3 (en) 2013-05-29
JP2012250043A (ja) 2012-12-20
US20130243349A1 (en) 2013-09-19

Similar Documents

Publication Publication Date Title
CN102819826B (zh) 降噪方法
Manduca et al. Projection space denoising with bilateral filtering and CT noise modeling for dose reduction in CT
Solomon et al. Quantitative comparison of noise texture across CT scanners from different manufacturers
US8965078B2 (en) Projection-space denoising with bilateral filtering in computed tomography
US10643319B2 (en) Apparatus and method for context-oriented blending of reconstructed images
Zhang et al. Statistical image reconstruction for low-dose CT using nonlocal means-based regularization. Part II: An adaptive approach
CN112651885A (zh) 用于减少图像记录噪声的方法和装置
US6778692B1 (en) Image processing method and apparatus including image improving circuit
CN108292430A (zh) 用于对功能医学成像中的定量图生成进行自动优化的方法
Prabhat et al. Deep neural networks-based denoising models for CT imaging and their efficacy
Li et al. Incorporation of residual attention modules into two neural networks for low‐dose CT denoising
US20080118128A1 (en) Methods and systems for enhanced accuracy image noise addition
Hayes et al. Low‐dose cone‐beam CT via raw counts domain low‐signal correction schemes: Performance assessment and task‐based parameter optimization (Part I: Assessment of spatial resolution and noise performance)
WO2022045210A1 (ja) 画像処理装置、画像処理方法、学習装置、学習方法、及びプログラム
Hao et al. A wavelet transform-based photon starvation artifacts suppression algorithm in CT imaging
Zheng et al. Improving spatial adaptivity of nonlocal means in low‐dosed CT imaging using pointwise fractal dimension
Al-Hinnawi et al. Assessment of bilateral filter on 1/2-dose chest-pelvis CT views
Wang et al. Locally linear transform based three‐dimensional gradient‐norm minimization for spectral CT reconstruction
Us et al. Combining dual-tree complex wavelets and multiresolution in iterative CT reconstruction with application to metal artifact reduction
Gomez-Cardona et al. Low signal correction scheme for low dose CBCT: the good, the bad, and the ugly
Huang et al. An iterative reconstruction method for sparse-projection data for low-dose CT
US20220414832A1 (en) X-ray imaging restoration using deep learning algorithms
Jiang et al. A novel method to improve the visual quality of X-ray CR Images
Zhang et al. Nonlocal means-based regularizations for statistical CT reconstruction
JP4707471B2 (ja) 画像処理プログラム、装置及び方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C41 Transfer of patent application or patent right or utility model
TA01 Transfer of patent application right

Effective date of registration: 20160718

Address after: Japan Tochigi

Applicant after: Toshiba Medical System Co., Ltd.

Address before: Tokyo, Japan, Japan

Applicant before: Toshiba Corp

Applicant before: Toshiba Medical System Co., Ltd.

GR01 Patent grant
GR01 Patent grant