CN102376084A - 使用各向异性噪声模型对ct图像的迭代图像滤波 - Google Patents

Abstract

使用各向异性噪声模型对CT图像的迭代图像滤波。本发明涉及一种用于处理计算机断层造影图像数据(f)的方法，其中，通过对图像数据(f)的减少噪声的处理来获得新的图像数据(f)。在此，对所述图像数据(f)进行加权的高通滤波，其中，所述加权逐个像素地考虑在不同方向

中的相应像素上的图像噪声，使得增大的噪声导致加强的高通效果。在使用加权的高通滤波的情况下，对所述图像数据(f)进行减少噪声的平滑。

Description

使用各向异性噪声模型对CT图像的迭代图像滤波

技术领域

本发明涉及计算机断层造影图像数据的用于减少噪声的处理。

背景技术

断层造影成像方法的特征在于，可以检查检查对象的内部结构而无需在此过程中对该检查对象进行手术干预。断层造影图像产生的一种可能方式在于，从不同角度拍摄待检查对象的一定数量的投影。从这些投影中可以计算出检查对象的二维截面图像或者三维立体图像。

这样的断层造影成像方法的一个例子是计算机断层造影。利用CT系统对检查对象进行扫描的方法是普遍公知的。在此，例如使用圆形扫描、具有进动的顺序圆形扫描、或者螺旋扫描。不基于圆形移动的其它类型的扫描也可以，例如具有线性片段的扫描。借助至少一个X射线源和至少一个与该至少一个X射线源对置的检测器从不同的拍摄角度拍摄检查对象的吸收数据，并且借助相应的重建方法将这样收集的吸收数据或投影结算为所述检查对象的截面图像。

为了由计算机断层造影设备(CT设备)的X射线CT数据组、也就是由所采集的投影来重建出计算机断层造影图像，目前采用所谓的滤波反投影方法(Filtered Back Projection；FBP)作为标准方法。在采集数据之后，通常执行所谓的“重排(rebinning)”步骤，在该步骤中对利用从源扇形扩展的射线所产生的数据进行重排，使得该数据就如同检测器被平行地进入该检测器的X射线击中的那样的方式存在。然后这些数据被变换到频域。在频域中进行滤波，接着对经过滤波的数据进行反变换。由此借助这样经过重排和滤波的数据，反投影到感兴趣体积内的各个体素。

在采集CT测量数据期间，通常是患者的检查对象被暴露在X射线辐射剂量中。由于该辐射对检查对象来说一般并不是无害的，因此力求设法采用尽可能小的辐射负担。但是，所使用的剂量直接与由CT测量数据重建的图像数据中的图像噪声相关：剂量的减小导致噪声的增加。因此，为了尽可能好地利用特定的辐射剂量，值得追求的是采用有效减少CT图像中的噪声的图像重建或处理方法。

发明内容

本发明要解决的技术问题是，展示一种用于处理计算机断层造影图像数据的方法，该方法带来噪声的减少。此外，还应当展示一种相应的控制和计算单元、一种CT系统、一种计算机程序以及一种计算机程序产品。

该技术问题通过具有本发明的特征的方法以及通过具有本发明并列的特征的控制和计算单元、CT系统、计算机程序以及计算机程序产品来解决。

在本发明的用于处理计算机断层造影图像数据的方法中，通过对图像数据的减少噪声的处理来获得新的图像数据，在该处理中对图像数据进行加权的高通滤波，其中，所述加权逐个像素地(bildpunktweise)考虑在不同方向中的相应像素上的图像噪声，使得增大的噪声导致加强的高通效果。在使用加权的高通滤波的情况下，对图像数据进行减少噪声的平滑。

所述处理涉及计算机断层造影图像数据，也就是由测量数据重建出的检查对象的图像数据，该测量数据是利用计算机断层造影系统采集的。通常，在测量数据采集时在计算机断层造影系统的辐射源与检查对象之间进行相对的旋转移动。为了由该测量数据计算出图像数据，可以采用本身公知的重建方法，尤其是滤波反投影方法。所得到的图像可以是二维的或三维的。该图像现在被以减少图像噪声为目的地处理。

减少噪声的处理包括至少两个步骤。一方面，对图像数据进行加权的高通滤波。该加权在此引起与例如在使用拉普拉斯滤波器的情况下普通高通滤波的偏差。与普通高通滤波的偏差至少取决于图像噪声及其它。因此，在加权时，对每个像素考虑该像素上的噪声。该噪声不必是各向同性的大小，各向同性的大小也就是与所考察的方向无关地具有相同的值。因此，不仅考察每个像素的单个噪声值，而且还考察不同方向中的噪声值。也就是说，基于相应的像素考察不同的方向，并且对于不同方向中的每个方向都使用一个噪声值。因此在高通滤波的加权中，噪声被用作与方向有关的参数。噪声值对不同的方向可以是不同的；但是还可能的是噪声值对一些或者所有方向是相同的。

在计算了加权的高通滤波之后，呈现经过高通滤波的图像。该图像被用于执行图像数据的减少噪声的平滑并且由此获得新的图像数据。

在本发明的扩展中，所述加权逐个像素地考虑相应的像素和其它像素的像素值之间的差异，使得增大的差异导致减弱的高通滤波效果。除了与方向有关的噪声之外，在加权时还考察不同像素的像素值之间的差异：取决于两个像素的像素值的差异有多大，这些像素之一或多或少地用于计算另一个像素的经过高通滤波的值。在此，至少在像素值差异的特定值域中，增大的差异会导致减弱的高通效果。

根据本发明的一种扩展，通过从图像数据中去掉经过高通滤波的图像数据，在使用加权的高通滤波的情况下进行减少噪声的平滑。在这种形成差的情况下，必要时可以采用加权因子，也就是说经过高通滤波的图像数据可以与必要时对不同的像素可以不同的因子相乘，以便这样相乘地从图像数据中除去。从图像数据中除去经过高通滤波的图像数据相当于低通滤波效果。按照这种方式实现了平滑。

特别有利的是，所述不同的方向是至相应像素的相邻像素的方向。在此，可以考察在相应像素的图像平面内的近邻，也就是在二维中或者是在三维中的近邻。优选地，相邻像素是像素的直接近邻；但是还可以涉及相隔更远的近邻。所考察的近邻优选是那些在高通滤波中被使用的像素；也就是说，在高通滤波中只考虑像素的直接近邻，从而由此在所述加权时也只采用这些方向中的噪声。

在本发明的扩展中，按照以下方式确定不同方向中的在相应像素上的图像噪声，即，通过在测量数据空间中的考察来获得穿过相应像素的属于不同方向的衰减值。通过考察原始的测量数据(由这些原始的测量数据重建出了待平滑的图像数据)，可以获得该测量数据空间。作为对此的替换，还可以将图像数据前向投影到该测量数据空间中。在计算机断层造影的情况下，通过测量穿过检查对象的衰减积分。在测量数据空间中可以显示该衰减积分或者由该衰减积分导出的参数，例如负对数。因此在测量数据空间中可以考察那些涉及特定像素的衰减值，也就是穿过该像素分布的衰减积分。

在本发明的实施方式中，确定最大和最小的衰减值，以便由此确定在不同方向中的相应像素上的图像噪声。因此，对每个像素都确定衰减值中的最大值和最小值。该最大值和最小值可以在大值和小值的范围内通过形成平均值来确定。然后可以从这两个值中确定不同方向中的衰减值。不同方向的数量优选地大于2。

根据本发明的另一实施方式，对于所述不同方向中的每一个方向都确定一个衰减值，以便由此确定相应方向中的在相应像素上的图像噪声。因此与上面的措施相反，不是确定两个衰减值，而是确定与不同方向的数量相应数量的衰减值。由于方向是已知的，因此可以在测量数据空间中简单地确定所属的衰减值。

根据本发明的一种扩展，图像噪声作为线性因子影响高通滤波效果的强度。这意味着，较大的噪声按照线性程度引起较强的高通滤波。除了该线性影响之外，图像噪声还可以通过其它方式影响高通滤波效果。

根据本发明的一种特别优选的扩展，借助一种函数来进行所述加权，该函数首先线性增加，并且在自变量变大的情况下比线性弱地增加。如果在该函数的自变量增大的方向上来考察该函数，则该函数在自变量的值为0的附近首先是是线性的。稍后，也就是在自变量变大的情况下，函数的斜率相对于线性变化过程来说减小。在此可能的是，该函数首先线性增加，并且在自变量变大时比线性弱地增加以及在自变量变得更大时减小；因此在这种情况下，该函数的斜率改变了其符号。此外还可能的是，该函数首先线性增加，并且在自变量变大时比线性弱地增加，在自变量变得更大时减小以及在自变量变得还要大时改变其符号。所述函数的这些不同的实施方式使得在该函数的自变量包含图像噪声时可以通过不同的方式让该噪声进入到高通滤波的加权中。优选地，所述函数的自变量包含两个像素的像素值之差与图像噪声之间的比例。按照这种方式，高通滤波效果取决于在相应像素上的各向异性的对比度与噪声之比。

特别有利的是，通过处理图像数据来获得新图像数据而不采用测量数据。这与迭代的重建算法相反，在迭代的重建算法中在基于所述图像的图像计算之后计算投影数据并且将该投影数据与测量数据进行比较，以便在使用计算出的投影数据与测量数据之间存在的偏差的情况下计算新的图像数据。与此相反，在这里，图像数据只需要用于由此计算改善的图像数据，而无需重新考虑测量数据。

在本发明的实施方式中，将新的图像数据作为结果图像数据输出。因此，在进行了减少噪声的处理之后就已经给出了不再基于用于减小图像噪声的进一步计算的图像数据。作为对此的替换，接着还可以对新的图像数据进行减少噪声的处理。后者意味着，现在可以执行先前基于图像数据被执行来用于计算新的图像数据的相同方法步骤，以便进一步处理新的图像数据。因此这是迭代的图像处理。

在迭代的图像处理中，可能的是对于不同的减少噪声的处理改变在不同方向中的相应像素上的图像噪声。这可以通过以下方式进行，即，对图像噪声的重新计算如在第一次迭代中所见的那样执行。作为对此的替换，可以基于关于由于减少噪声的处理而引起的噪声变化的特定假设来改变噪声值。

根据本发明的另一扩展，减少噪声的处理引起图像数据的取决于对比度的噪声减小。因此，不是均匀地在整个图像上进行平滑，由此可能失去清晰度；而是特别是去除在图像数据的对比度差的位置上的噪声，而在对比度强的位置上则注意保持清晰度。

本发明的控制和计算单元用于由CT系统的测量数据来重建检查对象的图像数据。该控制和计算单元包括用于存储程序代码的程序存储器，其中在该程序存储器中存在程序代码(必要时还有其它)，该程序代码适于执行以上描述类型的方法或者引起或控制所示实施方式。本发明的CT系统包括这样的控制和计算单元。此外，该CT系统可以包含其它例如用于采集测量数据所需要的部件。

本发明的计算机程序具有程序代码，该程序代码适于当该计算机程序在计算机上执行时执行上述类型的方法。

本发明的计算机程序产品包括存储在计算机可读数据载体上的程序代码，该程序代码适于当该计算机程序在计算机上执行时执行上述类型的方法。

附图说明

下面借助实施例更详细地解释本发明。在此：

图1示出了具有图像重建组件的计算机断层造影系统的实施例的第一示意图示，

图2示出了具有图像重建组件的计算机断层造影系统的实施例的第二示意图示，

图3示出了第一影响函数(Influenz-Funktion)，

图4示出了第二影响函数，

图5示出了在一个平面内的相邻像素，

图6示出了CT图像。

具体实施方式

在图1中首先示意性地显示了具有图像重建装置C21的计算机断层造影系统C1。该计算机断层造影系统是所谓的第三代CT设备，但是，本发明不限于该第三代CT设备。在支架外壳C6中具有在此未示出的封闭的支架，在该支架上设置了第一X射线管C2以及与该第一X射线管对置的检测器C3。可选地，在这里所示的CT系统中设置了第二X射线管C4以及与该第二X射线管对置的检测器C5，从而可以通过附加提供的辐射器/检测器组合实现更大的时间分辨率，或者在辐射器/检测器系统中使用不同的X射线能量谱的情况下还可以执行“双能量”检查。

此外，CT系统C1具有患者卧榻C8，患者在检查时可以在该患者卧榻上沿着系统轴C9(也称为z轴)被移动到测量场中，其中，扫描本身可以作为没有患者进动的纯圆形扫描仅在感兴趣检查区域中进行。患者卧榻C8相对于支架的移动通过合适的机械化来引起。在该移动期间，X射线源C2或C4分别围绕患者旋转。在此，在X射线源C2或C4对面检测器C3或C5并行地一起运转，以便采集投影测量数据，然后这些投影测量数据被用于重建截面图像。当然，替换其中患者在各次扫描之间逐步地移动穿过检查场的顺序扫描，还给出了螺旋扫描的可能性，其中患者在利用X射线进行旋转扫描期间连续地沿着系统轴C9移动穿过在X射线管C2或C4与检测器C3或C5之间的检查场。通过患者沿着轴C9的移动以及X射线源C2或C4的同时旋转，在测量期间在X射线源C2或C4相对于患者的螺旋扫描中产生了螺旋线轨道。该轨道也可以通过以下方式实现，即支架在患者未移动的情况下沿着轴C9移动。此外还可能的是，患者连续地以及周期地在两个点之间来回移动。

通过控制和计算单元C10来控制CT系统10，该控制和计算单元具有存在于存储器中的计算机程序代码Prg₁至Prg_n。要指出的是，这些计算机程序代码Prg₁至Prg_n当然还可以包含在外部存储介质上并且在需要时加载到控制和计算单元C10中。

从控制和计算单元C10可以通过控制接口24传输采集控制信号AS，以根据特定的测量协议控制CT系统C1。采集控制信号AS在此例如涉及X射线管C2和C4，其中可以预先给定它们的功率以及它们的接通时刻和关闭时刻，并且涉及支架，其中可以预先给定支架的旋转速度，以及涉及卧榻进动。

由于控制和计算单元C10具有输入控制台，因此可以由CT设备C1的使用者或操作员输入测量参数，这些测量参数然后以采集控制信号AS的形式控制数据采集。关于当前使用的测量参数的信息可以显示在控制和计算单元C10的显示屏上；附加地，还可以显示其它对操作员来说相关的信息。

将由检测器C3或C5所采集的投影测量数据p或者说原始数据通过原始数据接口C23传送给控制和计算单元C10。然后，该原始数据p(必要时在合适的预处理之后)在图像重建组件C21中被进一步处理。图像重建组件C21在该实施例中在控制和计算单元C10中以软件的形式实现在处理器上，例如按照所述计算机程序代码Prg₁至Prg_n中的一个或多个的形式实现。对于图像重建，正如已经针对测量过程的控制所解释的那样，这些计算机程序代码Prg₁至Prg_n还可以包含在外部存储介质上并且在需要时加载到控制和计算单元C10中。此外可能的是，由不同的计算单元来执行对测量过程的控制和图像重建。

然后，将由图像重建组件C21重建的图像数据f存储在控制和计算单元C10的存储器C22中，和/或按照通常的方式在控制和计算单元C10的显示屏上输出。该图像数据还可以通过在图1中未示出的接口馈入连接到计算机断层造影系统C1的网络，例如放射信息系统(RIS)，并且存储在该网络中可访问的海量存储器中或者作为图像输出。

控制和计算单元C10附加地还可以执行EKG的功能，其中，将导线C12用于在患者与控制和计算单元C10之间传导EKG电势。附加地，图1中所示的CT系统C1还具有造影剂注射器C11，通过该造影剂注射器可以附加地将造影剂注入到患者的血液循环中，从而例如患者的血管，尤其是跳动的心脏的心室能够被更好地显示。此外，由此还存在执行灌注测量的可能性，对于灌注测量来说所建议的方法同样适用。

图2示出C形臂系统，其中与图1的CT系统相反，外壳C6承载C形臂C7，在该C形臂上一侧安装了X射线管C2，而另一侧则安装了与该X射线管对置的检测器C3。C形臂C7对于扫描来说同样围绕系统轴C9摆动，从而可以从多个扫描角度进行扫描并且可以从多个投影角度获得相应的投影数据p。图2的C形臂系统C1与图1的CT系统一样也具有针对图1所描述类型的控制和计算单元C10。

本发明可以用在图1和图2所示的两种系统中。此外，本发明原则上还可以用于其它CT系统，例如用于具有形成一个完整环的检测器的CT系统。

由于在由计算单元C10重建的CT图像中包含临床上的相关信息，因此特别重要的是这些图像是具有说服力的。例如，还应当可以在图像中识别出小的肿瘤，也就是说，能很好地与周围组织区分开来并且就肿瘤的大小和位置来说是可识别的。因此，致力于在保持或甚至提高细节信息的可见度的同时减少CT图像中的噪声。由此可以在降低辐射剂量的情况下获得相同的图像质量或在相同的剂量的情况下获得更高的图像质量。

原则上可以通过以下方式减少CT图像中的噪声，即按照线性低通滤波器的形式来应用进行平滑的图像滤波器。但是在此有缺陷的是同时图像清晰度降低，由此从该CT图像中去除了细节信息。为了估计图像清晰度，可以考察CT图像中的与实际检查对象内的理想边缘突变相应的边缘陡度。CT图像内的边缘越陡，CT图像就越清晰。基于噪声的减少而进行的平滑导致该边缘的模糊，从而该边缘的陡度减小以及由此图像清晰度降低。

下面描述一种特别有效地减少CT图像的噪声的方法。该方法在此获得CT图像的细节信息，也就是说，清晰度不会或者几乎不被降低。同时，边缘周围的区域在减少噪声时不被排除在外；而是噪声减少也在这些图像区域中起作用。此外，该噪声减少方法不改变CT典型的噪声结构。因为CT图像具有典型的噪声功率频谱，分析CT图像的人员，尤其是放射科医师已经习惯于这种噪声功率频谱。由于这种习惯效应或训练效应，不期望图像噪声的统计特性发生根本性改变。

使用迭代的、非线性的图像滤波器来减少噪声。为此，首先由测量数据重建出CT图像V₀，下面称为第0次迭代的图像。为了重建图像可以使用本身公知的重建算法，例如FBP或迭代的重建方法。图像V₀可以是检查对象的二维截面图像或三维立体图像。相应的也适用于如下面描述的那样由图像V₀计算出的图像V_k。

对第0次迭代的图像以及接着对由该第0次迭代的图像计算出的图像V_k应用根据以下更新等式的迭代图像滤波器：

$V_{k+1}=V_k-{\mathrm{\γ}}_k\·V_k^E$ 等式(1)

在此，V_k+1是第(k+1)次迭代的图像，该图像是从第k次迭代的图像V_k中计算出来的。

γ_k是在迭代k中的滤波器强度。该滤波器强度是对不同的迭代可能改变的数字；但是该滤波器强度也可以保持恒定，从而γ_k＝γ。如果γ_k依据迭代地改变，则γ值例如可以在迭代运行时减小，使得正则化图像(Regularisierungs-bild)

对下次迭代的图像的作用越来越小。

是第k次迭代的图像的经过高通滤波的版本，也称为正则化图像。以下等式成立：

$V_k^E\left(i\right)=\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{d}_{\mathrm{ij}}\·{\mathrm{df}}_{\mathrm{ij}}\·H({\mathrm{df}}_{\mathrm{ij}}/{\mathrm{\σ}}^j\left(i\right))$ 等式(2)

i在此是正则化图像

的像素的标号，其中每个图像V_k和

的像素都从1数到N。对所有像素j进行求和。还可以代替对全部像素求和而对与像素i相邻的所有像素j求和。例如可以使用像素i周围的面积，例如3x3个像素的大小，或体积，例如3x3x3个像素的大小，并且仅在图像的该有限部分内进行求和。

图像

的高通特征通过域滤波器(Domaine-Filter)d_ij来产生，该域滤波器例如可以通过像素相互之间的反比距离(inverse Abstand)来给定。还可以代替该反比距离而使用具有高通特性的其它类型的滤波函数。

df_ij是像素i的像素值与像素j的像素值之差，也就是关于两个像素i和j的对比度。

σ^j(i)是在像素j方向中的像素i上的噪声。在CT图像的情况下，可以假定至少在一些图像区域中不存在各向同性的噪声，而存在各向异性的并因此是定向的噪声。这尤其适用于检查对象的具有各向异性衰减的部位。具有高噪声的方向与具有高统计不确定性的方向相应，该统计不确定性通过X射线的强烈衰减产生。

因此，在像素i中确定针对相邻的或者必要时还更远距离的像素j的局部对比度-噪声之比df_ij/σ^j(i)。用特征线H对该对比度-噪声之比df_ij/σ^j(i)加权。在此，该特征线H是影响函数的组成部分：

影响函数G＝df_ij·H(df_ij/σ^j(i))依据每个相邻像素i的图像值与直接相邻或更远相邻的像素j的图像值之差对高通滤波d_ij加权。该影响函数有利地满足以下特性：

G(-x)＝-G(x)，也就是该函数是非对称的，以及

G(ε)＞0，也就是对于小的正值来说G为正。

通过应用d_ij，将对图像进行“正常的”高通滤波。影响函数引起与“正常的”高通滤波之间的偏差，而且是依赖于CNR(对比度-噪声之比)。因为G的自变量一方面是图像值的差df_ij。另一方面，σ^j(i)也被用作G的自变量；有利的是图像值之差涉及局部噪声值，以便按照这种方式与局部噪声无关地执行正则化。为此定义：

$G(t;{\mathrm{\σ}}^j)={\mathrm{\σ}}^j\left(i\right)\·\tilde{G}\left(\frac{t}{{\mathrm{\σ}}^j\left(i\right)}\right)$ 等式(3)

针对影响函数

的一个具体示例是：

$\tilde{G}\left(t\right)=\frac{t}{1+{\left(\frac{\left|t\right|}{c_0}\right)}^p}$ 等式(4)

在图3中示出该函数，其中选择了c₀＝4和p＝2。可以识别出，影响函数

在开始时，也就是在值0的附近是线性的。也就是说，具有与各考察的像素类似的图像值的像素线性地被进行高通滤波；这相当于高通滤波的“正常”应用。随着自变量的变大(这相当于像素的图像值与现场各考察的像素的图像值之间的值偏差增大)，所述影响函数偏离线性：该影响函数首先比线性弱地增加，以便最后甚至下降。G值下降意味着相应的像素在高通滤波计算中被较少地考虑。由于σ^j(i)也进入G中，实现了随着CNR增加而下降的高通滤波效果。

高通滤波效果的下降在图3中通过曲线

表征。该曲线与值1偏离越多，影响函数

就偏离线性越多。在等式(4)中，c₀＝4是其中该影响函数对值的效果已经下降1/2的CNR，也就是相对于“正常的”高通滤波来说，像素仅还贡献大约一半。

因此，来自图3的影响函数的效果是，根据高通滤波在经过高通滤波的图像中包含小的和中等的边缘。相反，较大的边缘在高通滤波中仅被很少考虑，从而较大的边缘在经过高通滤波的图像中几乎看不见。所述影响函数越靠近值0，相应的边缘在经过高通滤波的图像中就越不重要。

甚至可以选择影响函数，使得该影响函数在自变量较大时改变其符号。一种示例是：

$\tilde{G}\left(t\right)=\frac{t\·(1-\frac{t}{c_1})}{1+{\left(\frac{t}{c_0}\right)}^p}$ 等式(5)

在图4中针对c₀＝4，c₁＝10和p＝1.5示出了该函数。如果该影响函数改变了符号，这会引起高通滤波的滤波器系数的符号改变。因此，针对所述较大的边缘进行具有负号的高通滤波。

原则上应当这样来选择影响函数，即，该影响函数首先线性增加，以便稍后比线性弱地增加。与线性增加之间的偏差可以表现得很小，使得该影响函数的增加即使对于大的自变量来说也还是正的。但是该增加也可以改变其符号(参见图3)，并且该影响函数在这种情况下也可以改变其符号(参见图4)。

所描述的影响函数只是示例性的。影响函数的其它合适的类别例如是GGMRF-先验(Generalized Gaussian Markov Random Field priors，通用高斯马可夫随机场先验)。

通过将利用经过修改的高通滤波进行滤波的图像根据等式(1)从图像V_k中除去，其效果如同低通滤波，并且由此如同对图像V_k进行平滑。在此，大的噪声引起更强的高通滤波效果并且相应地引起更强的平滑；对于较大的对比度值来说则相反。

等式(2)的前提是，各向异性的噪声σ^j(i)是已知的；下面描述如何能获得各向异性的噪声σ^j(i)。在此介绍两种可替换的措施，但是它们的一些部分也可被添加或者可被更换。

对于两种措施都成立的是，噪声σ^j(i)仅在二维，也就是在像素i的图像平面中确定。这在图像V_k是二维截面图像的情况下当然如此。相反，如果V_k是由多幅二维截面图像组成的三维立体图像，则为了确定噪声σ^j(i)要考察与像素i的截面图像相应的平面。对于位于该截面图像之外的像素j来说，该像素j与z轴平行地投影到像素i的图像平面上。然后作为针对该像素j的σ^j(i)，使用以下像素的噪声，即像素j的投影落到这些像素上。

对于第一种措施下式成立：

等式(6)

也就是说，由方位向量

获得在像素j方向中的像素方差，其中

表示从像素i到像素j的方向。

在图5中示出了方向

在该图中对位于中心的像素i示出8个相邻的像素。在三维中存在26个相邻像素。

σ(i)在此是这样一个参数，该参数是对假设为像素i上各向同性的噪声的度量。可以采用不同的方法来进行噪声估计，例如以下方法：对每个像素都计算沿着若干穿过该像素的直线的方差。因此沿着一条直线的图像值被视为统计的系综(Ensemble)并且计算该系综的方差。由此计算在不同的空间方向中的多个一维方差。这些方差中的最小方差将被作为针对相应像素i的结果值σ(i)输出。使用该最小值的原因在于，对于穿过边缘或已有结构的直线获得大的方差值。但是，通过σ(i)应当不是在相应像素的周围环境中识别出结构的存在，而是提供一种对噪声的度量。在最小方差值的情况下，可以假定该最小方差值典型地通过噪声而不是通过结构来受到影响。

替换地，还可以针对σ(i)使用与像素i无关的参数。

方位向量(Aspektvektor)

再现了投影到笛卡尔坐标轴

和(参见图5)上的、噪声椭圆的主轴。

通过以下方式获得该方位向量：

图像V_k被前向投影到数据空间中，从而达到正弦图空间。正弦图空间对每个探测器行来说都是一个二维空间，该二维空间一方面通过投影角度(也就是X射线源相对于检查对象的角度位置)以及另一方面通过X射线内的扇形角度(也就是通过检测器像素在通道方向中的位置)来展开。因此，所述正弦图空间是测量数据的域，而图像空间是图像数据的域。

在正弦图中，穿过像素i的所有投影的测量值都位于一条正弦形线上。通过沿着该正弦形线移动，查找最大衰减积分和最小衰减积分。其原因是，强衰减相应于高噪声值，小衰减相应于低噪声值。由此获得：

$k_{\mathrm{low}}=\frac{1}{p^{\mathrm{low}}+p^{\mathrm{high}}}\·p^{\mathrm{low}}$ 和 $k_{\mathrm{high}}=\frac{1}{p^{\mathrm{low}}+p^{\mathrm{high}}}\·p^{\mathrm{high}}$ 等式(7)

在此，p^low是沿着正弦形线的最小衰减积分，而p^high是沿着正弦形线的最大衰减积分。代替p^low和p^high的单个值，还可以使用在该正弦形线上的相邻测量值的区域上形成的平均值。在这种情况下下式成立：

$k_{\mathrm{low}}=\frac{1}{\mathrm{\Σ}{p}^{\mathrm{low}}+p^{\mathrm{high}}}\·\mathrm{\Σ}{p}^{\mathrm{low}}$ 和 $k_{\mathrm{high}}=\frac{1}{\mathrm{\Σ}{p}^{\mathrm{low}}+p^{\mathrm{high}}}\·{\mathrm{\Σp}}^{\mathrm{high}}$ 等式(8)

向量

给出噪声椭圆的主轴长度。在此，k_low表示最小衰减值，k_high表示最大衰减值。

在确定了主轴长度之后，从正弦图中获得最小衰减的所属方向

和最强衰减的所属方向

为此，只需要从正弦图中读取出p_low和p_high的y值，该y值与X射线内的相应扇形角度相应并由此与检测器像素在通道方向中的相应位置(也就是投影角度)对应。正如也在主轴长度

情况下那样，为了确定

和

在衰减最大值和衰减最小值周围的测量值中形成平均值。在这种情况下下式成立：

${\stackrel{\→}{e}}_{\mathrm{low}}\left(i\right)=\frac{1}{\mathrm{\Σ}{p}^{\mathrm{low}}}\·\mathrm{\Σ}{p}^{\mathrm{low}}\·{\stackrel{\→}{e}}^{\mathrm{low}}\left(i\right)$ 和 ${\stackrel{\→}{e}}_{\mathrm{high}}\left(i\right)=\frac{1}{\mathrm{\Σ}{p}^{\mathrm{high}}}\·\mathrm{\Σ}{p}^{\mathrm{high}}\·{\stackrel{\→}{e}}^{\mathrm{high}}\left(i\right)$ 等式(9)

在此，

和

分别表示在最小衰减方向和最大衰减方向中的单位向量。

方位向量

现在通过

的坐标变换给出：

等式(10)

等式(10)现在可以被用在等式(6)中以便获得σ^j(i)。

迄今为止假定在像素i上的噪声具有椭圆形状。这相当于只确定了两个衰减值以及属于该两个衰减值的方向。与此不同的是，可以考察更多数量的方向，从而一般来说不存在椭圆形状。

在用于确定σ^j(i)的第二种措施中，对针对像素j的每个方向都用角度

确定方向θ_j的标称噪声参数μ_i(θ_j)：

${\mathrm{\μ}}_i\left({\mathrm{\θ}}_j\right)=\frac{{\left({\mathrm{\σ}}^j\left(i\right)\right)}^2}{{\mathrm{\σ}}_0^2}=\frac{\exp \left({-p}_0\right)}{\exp \left({-p}_i\left({\mathrm{\θ}}_j\right)\right)}$ 等式(11)

参数exp(-p)是强度并因此是实际的测量参数，其中p是线积分并因此是衰减值。通常在正弦图中记录的不是强度，而是强度的负对数，也就是衰减值p。将在等式(11)中使用的p₀和p_i(θ_j)的值正如在第一种措施中那样从正弦图中获得。

因此，噪声(σ^j(i))²可以由测得的强度exp(-p_i(θ_j))与标称强度exp(-p₀)之比来给出。这样，在强衰减的方向上，投影数据具有较小的统计可靠性并因此具有增大的噪声。

p₀是标称衰减，例如在具有直径为30cm的水模体(Phantom)中的衰减，对于该衰减在给定核的情况下达到噪声σ₀。

p_i(θ_j)是在正弦图空间中在像素i的正弦线上的衰减值，其y值与方向θ_j相应。图6示出了一幅CT图像，在该CT图像中绘制了属于图5的相邻像素的投影。从所述正弦图中获得与这些投影相应的衰减值p_i(θ_j)。

因此，代替其中从正弦图中仅获得两个值(也就是最强衰减和最小衰减)的第一种措施，现在对每个像素j使用一个衰减值。

由此等式(2)变成：

${V_k}^E\left(i\right)=\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{d}_{\mathrm{ij}}\·{\mathrm{df}}_{\mathrm{ij}}\·H\left(\frac{{\mathrm{df}}_{\mathrm{ij}}}{{\mathrm{\σ}}_0\·{\mathrm{\μ}}_i\left({\mathrm{\θ}}_j\right)}\right)$ 等式(12)

在此，如上已经解释的，表达式

引起平滑优选在噪声较大的方向上进行。

在使用以下正则化项的情况下也获得在质量上相同的效果：

${V_k}^E\left(i\right)=\underset{j=1}{\overset{N}{\mathrm{\Σ}}}{d}_{\mathrm{ij}}\·{\mathrm{df}}_{\mathrm{ij}}\·{\mathrm{\μ}}_i\left({\mathrm{\θ}}_j\right)\·H\left(\frac{{\mathrm{df}}_{\mathrm{ij}}}{{\mathrm{\σ}}_0\·{\mathrm{\μ}}_i\left({\mathrm{\θ}}_j\right)}\right)$ 等式(13)

根据等式(13)，μ_i(θ_j)上的正则化强度也取决于方向。这意味着，正则化图像

的混合的程度(这与平滑的程度相应)取决于在特定方向中的噪声大小。一方面在特征线函数H中使用μ_i(θ_j)以及另一方面在正则化强度中使用μ_i(θ_j)在质量上是相同的，从而它们两个互相加强。

在所介绍的确定σ^j(i)的两种变型方案中，从原始数据获得关于噪声方向的信息。

在确定了各向异性的噪声σ^j(i)之后，可以借助等式(1)和(2)计算新的迭代图像V_k。在等式(1)中的表达式

的效果是对上一幅迭代图像V_k-1的各向异性的平滑。通过该平滑改变了图像V_k相对于图像V_k-1的噪声特性。强衰减的投影导致噪声在其相应的方向上增大。通过在低通滤波中包含噪声的各向异性(参见)，在这些方向上得到了更强的平滑。

这意味着，σ^j(i)也随着迭代不同而改变。为了考虑这一点，可以在每次迭代中重新计算σ^j(i)。但是由于这在计算上很费事，因此简化地建议使用取决于迭代的延伸因子δ_k。

由此对于第一种措施，等式(6)在这种情况下变成：

等式(14)

第二种措施的等式(11)变成：

${\mathrm{\μ}}_i\left({\mathrm{\θ}}_j\right)=\frac{{({\mathrm{\σ}}^j\left(i\right)\·{\mathrm{\δ}}_k)}^2}{{\mathrm{\σ}}_0^2}$ 等式(15)

δ_k可以随迭代不同而改变，以考虑迭代图像的噪声特性的变化。在此要考虑到，噪声的方向性随迭代不同而下降，也就是随着迭代的增加噪声与各向同性的特性匹配。

上面针对实施例描述了本发明。应当理解，可以进行众多的变化和修改而不脱离本发明的范围。

Claims (22)

1.一种用于处理计算机断层造影图像数据(f)的方法，

其中，通过对图像数据(f)的减少噪声的处理来获得新的图像数据(f)，在该处理中

-对所述图像数据(f)进行加权的高通滤波，其中，所述加权逐个像素地考虑在不同方向

中的相应像素上的图像噪声，使得增大的噪声导致加强的高通效果，

-在使用加权的高通滤波的情况下，对所述图像数据(f)进行减少噪声的平滑。

2.根据权利要求1所述的方法，其中，所述加权逐个像素地考虑相应的像素和其它像素的像素值之间的差异，使得增大的差异导致减弱的高通滤波效果。

3.根据权利要求1或2所述的方法，其中，通过从所述图像数据(f)中去掉经过高通滤波的图像数据，在使用加权的高通滤波的情况下进行减少噪声的平滑。

4.根据权利要求1至3之一所述的方法，其中，所述不同的方向

是至相应像素的相邻像素的方向。

5.根据权利要求1至4之一所述的方法，其中，按照以下方式确定不同方向

中的在相应像素上的图像噪声：即，通过在测量数据空间中的考察来确定穿过相应像素的属于不同方向的衰减值。

6.根据权利要求5所述的方法，其中，确定最大和最小的衰减值，以便由此确定在不同方向

中的相应像素上的图像噪声。

7.根据权利要求5所述的方法，其中，对于所述不同方向中的每一个方向都确定一个衰减值，以便由此确定相应方向

中的在相应像素上的图像噪声。

8.根据权利要求1至7之一所述的方法，其中，所述图像噪声作为线性因子影响高通滤波效果的强度。

9.根据权利要求1至8之一所述的方法，其中，所述加权借助一种函数

来进行，该函数首先线性增加，并且在自变量变大的情况下比线性弱地增加。

10.根据权利要求9所述的方法，其中，所述函数

首先线性增加，并且在自变量变大时比线性弱地增加以及在自变量变得更大时减小。

11.根据权利要求9所述的方法，其中，所述函数

首先线性增加，并且在自变量变大时比线性弱地增加，在自变量变得更大时减小以及在自变量变得还要大时改变其符号。

12.根据权利要求9至11之一所述的方法，其中，所述函数的自变量包含图像噪声。

13.根据权利要求9至12之一所述的方法，其中，所述函数

的自变量包含两个像素的像素值之差与图像噪声之间的比例。

14.根据权利要求1至13之一所述的方法，其中，通过处理所述图像数据(f)来获得新图像数据(f)而不采用原始的测量数据(p)。

15.根据权利要求1至14之一所述的方法，其中，将新的图像数据(f)作为结果图像数据(f)输出。

16.根据权利要求1至15之一所述的方法，其中，接着还对新的图像数据(f)进行减少噪声的处理。

17.根据权利要求16所述的方法，其中，对于不同的减少噪声的处理改变在不同方向

中的相应像素上的图像噪声。

18.根据权利要求1至17之一所述的方法，其中，减少噪声的处理引起所述图像数据(f)的取决于对比度的噪声减小。

19.一种控制和计算单元(C10)，用于从CT系统(C1)的测量数据(p)中来重建检查对象的图像数据(f)，

该控制和计算单元包括用于存储程序代码(Prg₁-Prg_n)的程序存储器，

其中，在该程序存储器中存在程序代码(Prg₁-Prg_n)，该程序代码执行根据权利要求1至18之一所述的方法。

20.一种具有根据权利要求19所述的控制和计算单元(C10)的CT系统。

21.一种具有程序代码(Prg₁-Prg_n)的计算机程序，用于当该计算机程序在计算机上执行时执行根据权利要求1至18之一所述的方法。

22.一种计算机程序产品，包括存储在计算机可读数据载体上的程序代码(Prg₁-Prg_n)，用于当该计算机程序在计算机上执行时执行根据权利要求1至18之一所述的方法。

Images