CN1505804A - 图像处理方法和图像处理设备 - Google Patents

Abstract

一种图像处理方法，在该图像处理方法，通过进行统计测试，对构成图像的第一像素与第二像素之间的相似性进行数字化，而且在数字化的相似性高时，平均第一像素值和第二像素值，而在所确定的相似性低时，不平均它们，从而在不损害空间分辨率和时间分辨率的情况下，降低噪声。

Description

图像处理方法和图像处理设备

技术领域

本发明涉及一种图像处理方法和图像处理设备，更具体地说，本发明涉及一种可以降低图像上存在的噪声的方法和设备。

背景技术

当前，图像处理技术已经开始应用于各种领域。

为了解决例如利用磁带录像机、数码相机等获得的图像的退化和改进，进行图像处理，或者为了例如清楚地掌握结构对象本身的图形或结构，进行图像处理以检查该结构是否是根据设计制造的。

此外，在诸如X射线CT设备、SPECT设备以及MRI设备的各种医学图像诊断设备中，执行多方面的图像处理。在执行例如描绘血流或造影剂流，采样病变部分，或者绘制内部器官的轮廓方面，在多方面受益。

图像处理技术包括各种构成技术，例如噪声抑制技术、特征提取技术以及图形识别技术，而且可以单独使用各种技术，或者以正确组合方式使用各种技术。顺便提一句，在这些构成技术中，特别是用于降低图像中的随机噪声的技术是一种必不可少的技术，它可以更清楚地再现被成像、重构等的对象。

然而，要求进一步改进现有技术的图像处理技术，特别是降噪技术。作为例子，所谓“平滑”通常被称为降噪技术。以这样的方式进行平滑，即在存在特定像素(i，j)的输入值f(i，j)时，提供该像素(i，j)附近的平均密度作为该像素(i，j)的输出值(i，j)。具体地说，假定使用像素(i，j)附近的n×n个像素，其输出值被表示为：

$g(i,j)=\underset{k=a}{\overset{b}{\mathrm{\Σ}}}\underset{j=c}{\overset{d}{\mathrm{\Σ}}}\frac{1}{(b-a+1)(d-c+1)}\•f(i+k,j+1)---\left(1\right)$

在此，上面等式(1)中的字母a、b、c和d表示整数。此外，等式(1)中的1/(b-a+1)(d-c+1)被称为“权重”。顺便提一句，图13示出a、b、c和d＝-1、1、-1和1的情况。

同时，通常已知，在计算从其方差为σ²的总体分布中单独采集的n个采样的平均值时，平均值的方差变成σ²/n。因此，根据等式(1)，上述“总体”(population)和“其方差σ²”分别对应于其随机变量是起因于包含在每个像素(i，j)的值内噪声的分量的概率分布及其方差，因此可以降低噪声产生的每个像素的值f(i，j)的分量。

然而，仅简单应用平滑过程，出现所谓“边缘不可视性”(edge obscurity)，降低图像的空间分辨率，而且会给人整个图像不清晰的感觉。例如，对于上述医学图像，即使在以最低可能噪声绘制微细血管结构的情况下，根据等式(1)的噪声抑制处理过程仍对原本不绘制血管结构的像素执行平均(平滑)过程。因此，虽然抑制了噪声，但通过进行平滑，也降低了表示血管结构的对比度，因此有时难以绘制微细血管结构。

根据上述情况提出本发明，而且本发明的目的在于提供一种在不会引起图像的不清晰的情况下，可以满意地抑制噪声，而且可以有效促进其他图像处理技术，例如图形识别技术的图像处理方法和图像处理设备，

发明内容

(本发明的基本组成)

本发明的要点是一种图像处理方法，其特征在于，在其中对构成特定图像的像素建立是矢量值或标量值的像素值、其中对该像素值与另一个根据前者像素值单独建立的像素值之间的拟和优度进行量化以及其中利用另一个像素值建立该像素的新像素值的情况下，通过在拟和优度高时，使另一个像素值的贡献大，而在拟和优度低时，使另一个像素值的贡献小，建立新像素值。

此外，在本发明中，可以根据构成该图像的其他像素很好地建立“单独建立的另一个像素值”。即，在这种情况下，对于一个像素和其他像素，像素值被分别表示为分别是矢量值或标量值的“一个像素值”和“其他像素值”。

此外，在本发明中，特别是在以将作为拟和优度的函数的权重函数应用于对各其他像素值计算的拟和优度值，以确定各其他像素值的权重的方式，以及以利用权重顺序计算其他像素值的加权平均值的方式来建立新像素值的情况下，通过在拟和优度高时，使权重大，以使相应其他像素值在其他像素值的加权平均值中的贡献大，而在拟和优度低时使权重小，以使相应其他像素值在其他像素值的加权平均值中的贡献小，建立新像素值。

根据该处理过程，利用对另一个被判定与一个像素“相似”的像素附加的重要性，建立新像素值。因此，不破坏空间分辨率，与在现有技术中相同，或者在图像是动态图像情况下，不破坏时间分辨率。

(权重函数和权重)

在本发明中，权重函数最好应该是关于拟和优度的非负单调递增函数。

此外，利用x和y分别表示上述一个像素和另一个像素时，利用v(x)＝(v₁(x)，v₂(x)，...，v_k(x))和v(y)＝(v₁(y)，v₂(y)，...，v_k(y))表示它们的像素值，而利用v’(x)＝(v’₁(x)，v’₂(x)，...，v’_k(x))表示要建立的新像素值。此外，ρ(x，y)表示拟和优度，w表示权重函数以及w(ρ(x，y))表示通过对拟和优度应用权重函数w获得的权重。然后，最好应该以下式表示的形式进行取加权平均值的处理过程：

${v\′}_k\left(x\right)=\frac{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}{v}_k\left(y\right)\•w\left(\mathrm{\ρ}(x,y)\right)}{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}w\left(\mathrm{\ρ}(x,y)\right)}---\left(2\right)$

(其中N(x)表示包括其他像素的范围)。

(以下将类似于用于建立新像素值v’(x)的上式的右侧的形式称为“根据本发明的相干滤波器”。)

(可以应用本发明的图像处理技术的例子)

本发明可以应用于各种应用中，例如以上描述的处理过程，其中通过例如图像的整个区域建立新像素值，从而降低图像的噪声，或者在图像内进行图形识别。如下所述，本发明在这些图像处理技术的分支技术中均显示良好效果。

(拟和优度)

在本发明中，根据对分别构成上述一个像素值和另一个像素值的标量值应用统计测试方法获得的临界概率，可以很好地量化拟和优度。顺便提一句，作为在此所称的“统计测试方法”，可以考虑例如“x平方测试方法”，将在以下进行说明。此外，通常，在此引入的临界概率与拟和优度的关系尤其包括其中在它们中的一个增加时，而另一个降低的情况。即，在建立新像素值的过程中，在临界概率升高(拟和优度降低)时，另一个像素的权重变小(另一个像素对要建立的新像素值的贡献变小)，而在临界概率降低(拟和优度升高)时，另一个像素的权重变大(对要建立的新像素值的贡献变大)。

(“另一个像素”的选择范围)

此外，本发明可以很好地从围绕上述一个像素的预定区域内选择另一个像素，特别是，预定区域是以上述一个像素为中心的邻近区域。利用这种结构，被判定与一个像素相似的另一个像素通常很可能存在与上述一个像素的周围，因此，在选择上述区域时，可以省略不经济的计算过程，以有效建立新像素值。

(建立像素值的方法)

在本发明中，可以选择建立(construct)像素值(一个像素和另一个像素)的各种方法。更具体地说，首先，在该图像包括多个图像的情况下，可以将像素值设置为其中分别排列了位于所有多个图像中的同一个点上的各像素的标量值的矢量值。

这时，多个图像可以很好地是构成动态图像的多个静态图像。在这种情况下，在建立新像素值的过程中，既不破坏空间分辨率，又不破坏时间分辨率。顺便提一句，这种情况下的良好应用就是其中动态图像是例如利用医学图像诊断设备获取的动态CT图像的情况。

此外，建立像素值的第二种方法对应于该图像包括同一个受治疗者的多种图像的情况。在此，如上所述，可以将像素值设置为其中分别排列了位于所有多种图像中的同一个点上的各像素的标量值的矢量值。顺便提一句，这种情况下的良好应用是例如其中该图像是利用SPECT设备或PET设备进行多窗***线照相获得的图像，而且其中多种图像是根据不同放射性同位素发出的γ射线互相不同的图像的情况。另一种良好应用是例如其中该图像是彩色图像，而且其中多种图像是通过将光分解为3原色获得的互相不同的图像。

此外，作为建立像素值的第三种方法，该图像是一个图像，而且可以根据一个图像建立像素值。更具体地说，可以将像素值表示为例如其中分别排列了一个图像上的特定像素的标量值和存在于该特定像素附近的像素的标量值的矢量值。顺便提一句，根据建立像素值的该方法，本发明可以应用于任何数字图像。

(用于实现本发明的图像处理过程的优选设备结构)

最后，用于实现该图像处理方法的设备结构最好应该包括：拟和优度量化装置，用于对构成特定图像的像素，判定被表示为矢量值或标量值的像素值与单独建立的另一个像素值之间的拟和优度；以及像素值计算装置，在利用另一个像素值，对该像素建立新像素值的过程中，该像素值计算装置运行，用于通过在拟和优度高时，使另一个像素值的贡献大，而在拟和优度低时，使另一个像素值的贡献小，建立新像素值。

附图说明

图1是示出临界概率p(x，y)、权重(x，y)以及虚假设滞后的必然关系的表。

图2是示出通过根据像素值v(x)的建立方法的差别，对它们进行分类，根据本发明的相干滤波器实施例应用中的各种例子的表。

图3是示出其中引入了相干滤波器的X射线CT设备的原理配置例子的示意图。

图4是用于从理论上解释基于根据本发明的相干滤波器的图像处理过程的说明图。

图5示出其中基于根据本发明的相干滤波器对构成动态CT图像的多个图像进行噪声抑制处理的例子，其中(a)示出还未被处理的原始图像，而(b)和(c)分别示出被处理的图像。

图6是示出在实施例中基于相干滤波器的噪声抑制处理流程的流程图。

图7示出其中基于根据本发明的相干滤波器对一个CT图像进行噪声抑制处理的例子；其中(a)示出还未被处理的原始图像，(b)示出处理的图像，而(c)示出通过平均(a)所示图像和(b)所示图像获得的图像。

图8(a)、(b)、(c)和(d)是分别示出3D采样血管例子的增强图像，而图8(e)是作为另一个采样例子，通过采样人脑的灰质获得的图像。

图9是示出根据一个图像建立像素值v(x)的方法的一个例子的说明图。

图10是示出根据一个图像建立像素值v(x)的方法的另一个例子的说明图。

图11是示出在将相干回归方法与AUC方法组合在一起的情况下获得的时间密度曲线，以及将现有技术方法与AUC方法组合在一起的情况下获得的时间密度曲线的曲线图。

图12涉及根据说明书主体内容中的等式(4)等估计方差σ2的另一种估计方法，而且是用于解释根据其他估计方法计算方差σ2++的计算方法的说明图。

图13是用于说明作为现有技术中的噪声抑制方法的“平滑处理”的说明图。

具体实施方式

现在，将参考附图说明本发明实施例。顺便提一句，以下将首先更一般地说明(以下称为项目“I”)根据本发明“相干滤波器”(coherent filter)的概况，此后将顺序说明对各图像处理方面应用相干滤波器的例子(项目“II”和“III”-“VIII”)。

在说明了各种应用的例子之后说明根据实施例的最通用形式的相干滤波器及其应用(项目“IX”)的例子，最后概述该实施例的辅助项目(项目“X”)。

(I根据本发明相干滤波器的一般说明)

首先，为了说明根据本发明的“相干滤波器”，做预备说明。

(I-1像素值v(x))

通常，利用多个像素构造使用诸如相机之类的成像装置获得的数字图像(或者，可以将该图像看作像素集)。在随后的说明中，将每个像素的位置表示为矢量x(即，坐标值的矢量)，将该像素X所具有的数值(例如，密度的数值表示)表示为k维矢量。对于二维图像，像素x是二维矢量，其表示表示在图像上位置的坐标值(x，y)。对特定像素x定义的“像素值v(x)”被表示为：

v(x)＝(v₁(x)，v₂(x)，...，v_k(x))              ...(3)

顺便提一句，将位于该等式(3)的右侧的v₁(x)，v₂(x)，...和v_k(x)分别称为该像素x的“标量值”(scalar value)，如下所述。

作为例子，在图像是“彩色图像”时，每个像素具有3原色(红、绿和蓝)的亮度(标量值)，因此可以将每个像素的像素值v(x)看作其维数为k＝3的矢量(假定其中在上述等式(3)右侧的各项中，下标是例如“红”、“绿”和“蓝”，请参考后面说明的等式(3.1))。此外，在例如其中图像是包括k个静态图像的动态图像，而且其中第n个图像的每个像素具有标量值v_n(x)的情况下，通过排列位于k个静态图像上的各相同公用点(相同坐标)上的各像素的像素值(标量值)建立的k维矢量值v_n(x)＝(v₁(x)，v₂(x)，...，v_k(x))是作为矢量值的像素值，如下所述。

(I-2拟和优度(goodness of fit)概率或临界概率p(x，y)以及权重w(p(x，y)))

对于像素x，要考虑到适当像素集N(x)(像素集N(x)可包括像素x)。随后，要考虑到像素x与每个作为像素集N(x)的元素的像素y之间的权重w(p(x，y))。权重w(p(x，y))具有以下说明的属性。

(I-2-1拟和优度概率或临界概率p(x，y))

首先，将说明用于控制权重w(p(x，y))的值的函数p(x，y)的重要性。该函数p(x，y)是用于量化本发明所称“拟和优度”的工具。一般的说，它给出用于指出在任何意义上像素x与像素y∈N(x)的相似程度的具体数值(例如，在像素x和像素y的像素值v(x)与v(y)之间注意到的统计偏差程度)。

更具体地说，作为例子，在函数p(x，y)给出小值时，判定在像素值v(x)与v(y)之间“不存在统计上的有效偏差(拟和优度高)”，而且判定像素x和y很可能相似，而在函数p(x，y)给出大值时，判定“存在统计上的有效偏差＝拟和优度低”。

同时，必须考虑到噪声不明显包含在像素值v(x)和v(y)中(或者标量值v₁(x)，...，v_k(x)以及v₁(y)，...，v_k(y))。作为例子，在考虑到利用CCD成像器件获取图像的情况时，对于构成图像的各像素，存在因为该器件内的暗电流和从外部进入的光量的不规则波动产生的噪声。

通常，对于所有像素，这种噪声取互异值。因此，即使像素x和y反映同一个对象(在外部)，在实际观测的图像上，它们有时仍不能具有同样的值。反过来说，假定在从反映同一个对象的像素x和y中去除了噪声的情况下，在表示同一个对象时，在图像上显示(识别)它们，而且这两个像素应该具有同样(或非常接近)的像素值。

因此，在考虑到上述噪声的属性，对函数p(x，y)应用统计测试方法中众所周知的“虚假设”(null hypothesis)概念时，以下具体解释该函数p(x，y)。如果建立虚假设H，即“在去除了像素x和y的噪声时，像素x和y具有同样的像素值”，  换句话说，即“如果起因于这两个像素的噪声的偏差除外，则v(x)＝v(y)”(即，在肯定该命题时，可以认为“像素x与y相似(＝拟和优度高)”)，在放弃假设H的情况下(在这种情况下，p(x，y)被定义为其范围为[0，1]的函数(p(x，y)∈[0，1]))，可以将函数p(x，y)称为临界概率(或有效度)。

因此，在临界概率p(x，y)高的情况下，即在该放弃很可能是错误的情况下，可以说实现假设H的概率高，相反，在临界概率p(x，y)低的情况下，即，在该放弃几乎不可能是错误的情况下，可以说不能实现假设H的概率高。(顺便提一句，在统计测试方法中众所周知，即使在不“放弃”假设H时，也不意味着该假设为“真”。在这种情况下，仅意味着可以否定利用假设H表示的命题)。

(I-2-2权重w(p(x，y)))

同时，从其表达式可以看出，权重w(p(x，y))是以上说明的临界概率p(x，y)的函数(在更一般情况下，可以将它构造为w(ρ(x，y))，作为拟和优度的函数(其中ρ(x，y)表示拟和优度))。此外，为了建立权重w(p(x，y))而对像素x和y的每个组合获得的临界概率p(x，y)作用的权重函数w，通常具有体现上述“放弃”的作用。具体地说，以这样的方式调整权重函数w，即，在临界概率p(x，y)高的情况下，权重函数w，即权重w(p(x，y))的值取大的正值，而在相反的情况下，它取小的正值(或“0”)(以下将说明权重函数w的具体形式)。也就是说，在像素x和y可能实现假设H表示的命题的情况下，权重w(p(x，y))取大值，而在相反情况下，它取小值。作为例子，函数w很可能是特殊构造的函数，以使它可以仅取两个值，即“0”和非“0”的常量值。

顺便提一句，图1对以上说明的假设H、临界概率p(x，y)以及权重w(p(x，y))之间的关系进行了概括。此外，更一般的说，还可以将权重函数w(t)称为“对t∈[0，1]定义的非负单调递增函数”，而且这至少可能是该函数w(t)要满足的属性。

(I-3相干滤波器)

由于在以上进行了预备说明，所以可以得出根据本发明的“相干滤波器”，如下所述。首先，对于构成图像的特定像素x，对作为像素集N(x)的元素的所有像素y计算以上说明的权重w(p(x，y))。随后，利用多个权重w(p(x，y))，根据下面的等式计算构成像素x的新标量值v’_k(x)：

${v^{\′}}_k\left(x\right)=\frac{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}{v}_k\left(y\right)w\left(p(x,y)\right)}{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}w\left(p(x,y)\right)}---\left(4\right)$

在此，k＝1，2，...，K。此外，利用根据该等式获得的v’_k(x)，对像素x进行变换后获得的像素值(新像素值)v’(x)被表示为：

v’(x)＝(v’₁(x)，v’₂(x)，...，v’_k(x))    ...(5)

其中，在该实施例中，用于如上述等式(4)将像素值v(y)＝(v₁(y)，v₂(y)，...，v_k(y))(包括其中y＝x成立的情况)变换为v’(x)＝(v’₁(x)，v’₂(x)，...，v’_k(x))的滤波器是根据本发明的“相干滤波器”。从其表达式可以看出，新像素值表示构成该像素值的标量值v_k(y)的加权平均值。

该处理过程产生以下结果。像素值v’(x)表示由加权平均值v’_k(x)构成的矢量，该加权平均值v’_k(x)具有对去除噪声后可能与像素x的像素值具有同样像素值的像素y附加的重要性(＝非常可能实现假设H的命题)。此外，如果存在足够数量的这种像素y，则像素值v’(x)将具有仅利用上述平均作用抑制了噪声的值，而不偏离像素x应该具有的真值。

顺便提一句，即使在临界概率p(x，y)低的情况下，仍相应“放弃”虚假设H，因此权重w(p(x，y))变小，不始终完全“放弃”该标题，正如根据上述说明或等式(3)的表达式所理解的那样。即使在临界概率p(x，y)接近“0”(＝0％)时，仍可以设置权重w(p(x，y))≠0(小于概率p(x，y)接近“1”情况下的正值)，但是这依赖于以下说明的权重函数w的具体形式。(顺便提一句，p(x，y)＝1对应于如下所述的保持v(x)＝v(y)的情况。)也就是说，与完全放弃不同，很可能允许小贡献(contribution)。(顺便提一句，如果在这种情况下设置权重w(p(x，y))＝0，则这与完全放弃具有同样的意义。请参考后面的等式(14)。)

以下对该处理过程做一般性说明。可以说在图像处理方法中，在存在构成特定图像的(多个)像素x时，对每个像素x与要求的特定像素y(在上述说明中被设置为y∈N(x))之间的拟和优度进行量化(根据以上说明的p(x，y))，而且在拟和优度高时，在利用像素值v(y)进行加权平均处理的过程中，允许像素y具有大贡献，而在拟和优度低时，仅允许小贡献，因此可以有效抑制像素x的噪声。换句话说，在像素x和y“互相相似”时，使得像素y对平均处理起更大作用，而在它们“互相不相似”时，几乎或者完全忽略像素y。

通过对整个图像执行该处理过程，显示非常高效抑制噪声，而且几乎不会使图像不清晰。此外，不仅在应用于噪声抑制过程中，而且在应用于例如图形识别领域时，通过以优选具体形式实现权重函数或相干滤波器，可以显示出良好的处理效果。顺便提一句，在以下对各种图像处理方面所做的说明中，将说明这些效果具体是怎样的。

(II各种应用例子；根据构建像素值v(x)的方法的差别进行分类)

同时，根据上述通用形式的相干滤波器，本发明可以应用于各种领域。作为例子，图2所示的表概括了各种应用。顺便提一句，在图2所示的表中，根据建立像素值v(x)的方法的差别对各领域内的各种应用例子进行分类。根据图2，可以明白，像素值v(x)是以上说明的矢量值，可以利用各种解决方案，“根据多个图像”(利用例如一次射线照相获得的)、“根据多种图像”或者“根据一个图像”建立该像素值v(x)，而且可以明白，涉及该实施例的各种应用例子被归入各种解决方案。

现在，将按照图2所示分类表，顺序说明其中具有上述通用形式的相干滤波器应用于更实际情况的各种应用例子，以及如何具体实现建立图2所示像素值v(x)的3种方法(与图2所示内容不同)。

(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)

首先，说明在X射线CT设备中将相干滤波器应用于所谓“动态CT”射线照相的情况。在此，如图3所示，X射线CT设备100包括X射线管101以及X射线检测器102、数据采集***103(DAS)、预处理单元104、存储单元105、重构单元106和图像显示单元107。此外还包括用于控制这些组成部分的控制单元108和可以对X辐照条件、射线照相模式等进行各种设置和输入的输入单元109。

根据X射线CT设备100，在X射线管101和X射线检测器102绕病人P旋转时，对受治疗者或病人P照射X射线，因此可以将病人P的内部状态显示为CT图像，例如层析X射线照片。

这时，上述组成部分工作，如下所述。通过施加高压发生器101a产生的电压，由X射线管101发出的通过病人P的X射线被X射线检测器102变换为模拟电信号。此后，数据采集***103对该模拟电信号进行数字变换处理，预处理单元104对其进行各种校正处理，然后将获得的数字信号存储到存储单元105内作为投影数据。此外，根据以例如使X射线管101和X射线检测器102绕病人旋转一圈(＝360°)的方式获得的投影数据，重构单元106重构层析X射线照片或任何其他CT图像，然后，图像显示单元107显示该CT图像。顺便提一句，还可以将重构的CT图像存储到存储装置10M。

在此，以上说明的“动态CT”射线照相表示了一种射线照相解决方案，其中X射线管101和X射线检测器102反复对病人P的同一个部分进行射线照相(通过反复扫描，在连续旋转CT设备时，通常执行基于连续旋转的反复射线照相)，以连续采集投影数据，然后，根据该投影数据连续执行重构过程，以按时间顺序获得一系列图像。  (在这种情况下，例如在未示出的计数器的控制下，图像显示在图像显示单元107上，以便在与采集相关图像的源投影数据相关的扫描开始点或扫描结束点之后的预定时间进行显示。)

因此，这样获取和显示的图像变成象运动图像等一样包括多个时间顺序(time-series)静态图像的所谓“动态图像”。顺便提一句，为了通过将造影剂注入病人P并观察和探查图像随时间的变化，来分析病变部分的病变，例如血管的收缩或阻塞，通常采用这种射线照相解决方案。此外，在广义上，可以将其中仅在注入造影剂之前和之后，对同一个部分进行两次CT射线照相的解决方案看作动态CT射线照相。

同时，在现有技术中，如果例如在进行K次射线照相时，在上述“动态CT”射线照相模式下，病人P经受任何变化(例如，通常考虑到造影剂的密度变化，或呼吸运动)，则为了抑制图像噪声而不破坏空间分辨率，必然要在时间方向进行平滑。因此，已经不可能避免破坏时间分辨率(temporal resolution)的缺陷。

在这方面，利用动态CT射线照相获得的图像是动态图像，而且意在仔细观察时间变化，如上所述，因此不能将破坏时间分辨率称为良好情况。

在采用根据本发明的相干滤波器时，可以执行以下说明的、可以对所有K个静态图像(多个图像)抑制噪声而不破坏时间分辨率的处理(以下称为“动态相干滤波”)。

首先，关于对以上述方式获得的、构成动态图像的K个静态图像定义的像素x，可以建立已经说明的像素值v(x)：

v(x)＝(v₁(x)，v₂(x)，...，v_k(x))    ...(上面说明的3)

在此，右侧各项的下标1，2，...和k是对K个静态图像分配的编号(参考以上对等式(3)所做的说明)。

随后，在这种情况下，利用例如下面的等式给出具体形式的权重函数w1：

$w1\left(p(x,y)\right)=\exp \[-\{\sqrt{\underset{k=1}{\overset{k}{\mathrm{\Σ}}}\frac{\{v_k\left(x\right)-v_k\left(y\right){\}}^2}{k{\left(2{\mathrm{\σ}}_k\right)}^2}}{\}}^c\]---\left(6\right)$

在此，y∈N(x)成立，而对于像素x，可以任意设置像素集N(x)(＝可以根据任意判据进行设置)。然而，实际上，可以认为像素x和其位置远离像素x的像素y通常以低概率实现假设“除了因为这两个像素的噪声产生的差别，v(x)＝v(y)”，因此，它具有实际意义，例如提高计算速度，这样利用作为接近像素x的像素集判据限制像素集N(x)。

因此，在此，作为例子，假定像素集N(x)是其各像素包括在以像素x为中心的矩形区域内的一组像素。更具体地说，在例如构成现在注意到的一个静态图像的所有像素为128×128像素情况下，可以选择像素集N(x)作为围绕像素x的3×3像素区域。此外，在所有像素为512×512像素情况下，可以选择像素集N(x)为围绕像素x的13×13像素区域。

此外，在上面的等式(6)中，σ_K表示假定第k个静态图像的各像素具有某种一般程度的噪声时估计的标准偏差，而C表示在将权重函数w1(p(x，y))代入上述等式(4)情况下，可以确定并调节作用程度的参数。

现在，将顺序说明参数σ_K和C。

首先，解释等式(6)中的参数σ_K(以下将它作为方差σ_K ²说明)。如上所述，方差σ_K ²是第k个静态图像上各像素的标量值具有的噪声分量的方差。此外，已假定第k个图像上各像素的标量值含有作为常量值的方差σ_K ²的噪声，估计了上面的等式(6)中的方差σ_K ²。通常，根据以下说明的背景技术，该假定具有令人满意的理由。

首先，在病人P的尺寸以及X射线管101和X射线检测器102、重构单元106等的结构不变，而且辐照X射线的能力保持不变的情况下，利用X射线辐照量，即X射线管101内的管电流与辐照时间周期的乘积(被称为“电流时间积(mA·s)”)确定CT图像的噪声，它与X射线辐照量成正比。另一方面，还知道CT图像的噪声均是加性的(additive)，而且与高斯分布大致一致。更具体地说，对于构成某一像素x的像素值v(x)的任意标量值v_n(x)(n＝1，2，...，K)，在设v_n ⁰(x)表示真值(通过消除噪声的贡献分量获得的值)时，它们之间的差值，即v_n(x)-v_n ⁰(x)与平均值0和方差σ_K ²的高斯分布大致一致(顺便提一句，X射线辐照量或电流时间积mA·s与该噪声的方差σ_K ²大致成反比关系)。

此外，尽管方差σ_K ²还取决于像素x的位置(例如，上述坐标值x＝(x，y))，但是在普通X射线CT设备100中，可以忽略它，因为将用于调节X射线辐照量的物理X射线滤波器(被称为利用例如铜线圈或金属块构造的“***物”(wedge)或“X射线滤波器”)(未示出)插在X射线管101与X射线检测器102之间。以下说明其原因。***物具有调节辐照的X射线的一部分或全部X射线的功能，因此根据病人P由其密度大致等于水的物质构成这个事实，在任意X射线检测器102内可以检测到近似同样的X射线辐照量。因此，这种***物产生效果是，使噪声的方差σ_K ²变成几乎与像素x的位置无关的常量值(顺便提一句，通常利用有效使用X射线检测器102的动态范围的源对象设置***物)。

根据以上事实，在通过动态CT射线照相获得的K个静止图像中，对第k个静态图像上的所有像素估计方差σ_K ²为常量值是合理的。当然，对于每个像素的方差不相同的情况，可以容易地扩展该实施例(将在后面的X-3小节进行说明)。

随后，该问题变成，利用什么数值作为方差σ_K ²具体计算上面的等式(5)。该问题起因于，尽管可以假设噪声分布的形状(在以上说明中，为高斯分布)，但是通常并不清楚方差σ_K ²的具体数值。

此外，通常，通过在每次进行射线照相时改变射线辐照量(X射线管电流×辐照时间周期(mA·s))，可以很好地进行射线照相。

现在，设σ_K ²表示第k个图像中的各像素的标量值具有的噪声方差，而R_K表示在射线照相第k个图像时使用的射线辐照量，σ_K ²与R_K成正比。因此，在可以对至少一个图像k＝k₀指定σK₀ ²时，根据下面的等式，也可以精确估计与其他图像k有关的方差σ_k ²：

${{\mathrm{\σ}}_k}^2={{\mathrm{\σk}}_0}^2\frac{R_k}{R_{k0}}$

在该实施例(对其应用这种情况)中，对于至少一个图像k，利用如下所述的方法可以估计方差σ_K ²的具体数值。

该方法是有效的，其中利用可以假定在K次射线照相中病人P几乎不经受变化的N(1＜N≤K)个图像，通过实际进行测量可以发现方差σ_K ²的期望值E[σ_K ²]。为了简化后面的说明，假定N个图像基于同样的射线辐照量，因此对于k＝1，2...，N，方差σ_K ²为常数(被书写为σ²)。可以推测，包含在构成N个图像中的某个像素xf的像素值v(xf)的各标量值v₁(xf)、v₂(xf)、...以及v_K(xf)内噪声与平均值0和方差σ²的高斯分布一致，如上所述。因此，如下利用标量值的平均值：

$v^{*}\left(\mathrm{xf}\right)=\frac{1}{N}\underset{k=1}{\overset{N}{\mathrm{\Σ}}}{v}_k\left(\mathrm{xf}\right)---\left(7\right)$

可以利用下式计算真值方差σ²的期望值E[σ²]：

$E\left[{\mathrm{\σ}}^2\right]$ $=$ $\frac{1}{N-1}\underset{N-1}{\overset{N}{\mathrm{\Σ}}}\{v_k\left(\mathrm{xf}\right)-v^{*}\left(\mathrm{xf}\right){\}}^2---\left(8\right)$

在此，对于如上所述所有K个静态图像上的所有像素x，可以认为该方差的期望值E[σ²]是合理的，而且它是在利用其代替真值方差σ²时其可信度确保在预定程度或在预定程度之上的值。因此，在实际计算上面的等式(6)时，可以将期望值E[σ²]代入等式(6)的项σ²。

更具体地说，利用根据例如K个静态图像中的第一和第二静态图像的实际测量值可以很好地计算该期望值E[σ²](在上面的等式(7)和(8)中，这对应于设置N＝2)。此外，关于实际计算等式(7)和(8)中使用的像素xf，可以很好地设想，仅选择例如成像空气或骨骼的各部分之外的适当像素xf(在已经选择了多个像素的情况下，对获得的所有期望值E[σ²]计算平均值)。此外，通常，最好设想一种方法抑制因为病人P的移动而产生的干扰。(顺便提一句，估计等式(6)中噪声的方差σ²可以重新参考最后说明的“实施例的辅助项”。)

可以容易地建议，即使在N个图像的射线照相步骤，射线辐照量不是常数的情况下，利用其对射线辐照量R_K的比值，可以精确估计方差θ_R ²。

接着，将解释等式(6)中的参数C。首先，在等式(6)中，包括了在通用形式的说明中说明的临界概率p(x，y)的思想，如下所述。等式(6)的右侧的分子中的根与x²值一致，x²值与所谓“x-平方分布”一致。以利用x²除以(2σ)²，然后将大括号内的全部内容放置到e的肩部获得的值就是临界概率p1(x，y)本身。即：

$p1(x,y)=\mathrm{Aexp}\{-\underset{k=1}{\overset{k}{\mathrm{\Σ}}}\frac{\{v_k\left(x\right)-v_k\left(y\right){\}}^2}{K{\left({2\mathrm{\σ}}_k\right)}^2}\}---\left(9\right)$

此外，对于利用该等式(9)表示的p1(x，y)，上面的等式(6)可仅表示为下式：

$w1\left(p(x,y)\right)=\{p1(x,y){\}}^{-\frac{c}{2}}---\left(10\right)$

顺便提一句，字母A表示用于归一化临界概率p1以变成(0至1)的值的常数。

总之，等式(6)未确实指出通用形式的上述临界概率p(x，y)，但是可以认为权重w1(p(x，y))的实际状态是如上所述的临界概率p1(x，y)的极函数(very function)(等式(10))，即，它是“拟和优度的函数”。(如上所述，临界概率与拟和优度之间的关系是，在一个升高时，另一个也升高。)

此外，从上面的等式(10)可以看出，参数C具有确定权重w1(p(x，y))如何对临界概率p1(x，y)敏感地起反应的作用。即，在参数C增大时，仅通过少许降低临界概率p1(x，y)，权重w1(p(x，y))接近0。相反，在参数C变小时，可以抑制过于敏感的反应。顺便提一句，可以将参数C具体设置为约1至10，而且应该最好设置为C＝3。

在该实施例中，根据上述描述可以明白，利用所谓“x平方测试方法”(统计测试方法)，根据临界概率p1(x，y)，对两个像素x和y进行相似性判定，换句话说，对这两个像素x和y放弃上述虚假设H进行判定。

此外，从上面的等式(6)的表达式可以看出，在本发明中，不需要始终跟踪对像素x、y的各组合计算临界概率p(x，y)，然后得出权重w(p(x，y))的各步骤，但是可以直接计算复合函数(wOp)，而无需具体计算临界概率p(x，y)。

因此，如上所述，在等式(8)中估计(例如，E[σ²])方差σ²，然后正确确定参数C(例如，C＝3)，从而利用等式(6)，对包括在对某一像素x定义的像素集N(x)内的所有像素y，得出具体权重w1(p(x，y))(例如，如上所述对应于围绕像素x的3×3像素的区域)。此后，利用权重w1(p(x，y))代替上面的等式(4)内的w(p(x，y))，因此可以执行相干滤波器的具体数值计算过程。因此，具有噪声的像素值v’(x)＝(v’₁(x)，v’₂(x)，...，v’_k(x))被强烈抑制，换句话说，可以获得该k个静态图像或者该动态图像，同时不仅不破坏时间分辨率，而且不破坏空间分辨率。

为了容易从理论上掌握它，图4示出了该图像处理过程。首先，参考图4(a)，在第一、第二，...以及第K个图像中，对某一像素x假定，以该像素x为中心并对应于3×3个像素的矩形区域N^{3 ×3}(x)。设y₁表示矩形区域N^3×3(x)的左上角的像素，该像素y₁具有像素值v(y1)，如图4所示。

其中，利用上面的等式(6)，根据构成像素值v(y₁)的标量值v₁(y₁)，v₂(y₁)，...，v_k(y₁)以及像素值v(x)的标量值v₁(x)，v₂(x)，...，v_k(x)，计算权重w1(p(x，y))(图4(b))。此外，对矩形区域N^3×3(x)内的剩余像素y₂，...和y₈进行类似处理，  因此最后获得了权重w1(p(x，y₁))、...、w1(p(x，y₈))和w1(p(x，x))，如图4(b)所示。(在这种情况下，根据等式(9)，临界概率p(x，x)为“1”，因此，根据等式(10)，权重w1(p(x，x))也是“1”(＝实现了最大权重)。)

随后，将这样获得的权重w1(p(x，y₁))、...、w1(p(x，y₈))和w1(p(x，x))分别乘以第k个图像中相应像素的标量值v_k(y₂)、...v_k(y₈)以及v_k(x)，对所获得的乘积求和(对应于上面的等式(4)中的分子)，将所获得的和数除以矩形区域N^3×3(x)的权重w1的和数(对应于同一个等式(4)中的分母)。然后，可以计算抑制了第k个图像的像素x的噪声的标量值v’_k(x)(图4(c))。此外，利用同样的权重w1(p(x，y₁))、...、w1(p(x，y₈))和w1(p(x，x))，可以对k＝1、2、...、K的所有图像计算抑制了噪声的标量值v’_k(x)，从而获得抑制了像素x内的噪声的像素值v’_k(x)＝(v’₁(x)，v’₂(x)，...，v’_k(x))。在对所有像素x重复进行上述计算时，可以获得抑制了噪声的K个图像。

利用这样计算的像素值v(x)构造的图像变成例如图5所示的图像。图5(a)示出全部k＝23个静态图像中，应用相干滤波器之前的原始图像(k＝11)，而(b)和(c)分别示出对k＝11和k＝23应用相干滤波器之后的图像。图5(b)和(c)所示的图像充分抑制了在图5(a)所示的原始图像中可以看到的随机噪声。

顺便提一句，根据例如图6所示的流程图，可以执行以上解释的各处理步骤。此外，为了根据处理步骤在实际X射线CT设备100上进行计算、图像显示等，可以通过设置例如图像处理单元110进行计算、图像显示等，该图像处理单元110配置了方差估计部分111、权重计算部分112以及像素值计算部分113，如图3所示。

其中，这样构造权重计算部分112，即如以上说明的过程，可以利用像素值v(x)和v(y)直接计算权重w1(p(x，y))。因此，计算部分112是不具体计算临界概率p1(x，y)的值，而直接计算权重的装置(即，通过包括临界概率计算部分(在本发明中称为“拟和优度量化部分”))。顺便提一句，与上面说明的结构不同，可以采用这样的结构，它可以对“临界概率计算部分(拟和优度量化部分)”和“权重计算部分”的两级进行跟踪，“临界概率计算部分(拟和优度量化部分)”具体计算临界概率p1(x，y)的值，而“权重计算部分”根据前者部分的输出计算权重w1(p(x，y))。不论怎样，利用方差估计部分111估计的方差σ²以及像素值v(x)和v(y)，权重计算部分112可以计算权重w1(p(x，y))。

此外，利用像素值v(x)和v(y)以及权重计算部分112数值计算的权重w1(p(x，y))，像素值计算部分113计算像素值v’(x)。换句话说，计算部分113实际执行处理抑制原始图像的噪声的过程，即，应用相干滤波器(以下将该过程表达为“应用相干滤波器”)。

在如上所述进行动态相干滤波的情况下，对由K个静态图像构成的动态图像应用相干滤波器，可以在图像处理单元110内很好地进行该处理过程，以致可以一次性重构所有动态图像，此后将它们存储到存储装置10M作以后的后处理。然而，本发明并不局限于该方面，在上述的连续扫描、连续投影数据采集、以及连续显示的流程中，可以很好地实时执行应用相干滤波器的处理过程(以下将该处理过程称为“实时相干滤波”)。

在实时相干滤波的优选实施例中，在每次获得并重构新图像时，执行如下所述的处理。在具有第M、M-1、...以及M-K+1个图像的K个静态图像上，在通过最新图像(第M个图像)首先获得的图像(第1个图像)中，对同一个公共点(相同坐标)上的各像素x所具有的像素值(标量值)进行排列以建立K维矢量值v(x)＝(v_M(x)，v_M-1(x)，...，v_M-K+1(x))。这样，可以以与上述“动态相干滤波处理”非常类似的方式应用相干滤波器。在此，像素值计算部分113不实际计算像素值v’(x)的所有元素，而是仅计算对应于最新图像(第M个图像)的标量值v_M’(x)。因此，提高了计算速度，可以实时显示抑制了噪声的最新图像。

可以很好地这样建立“实时相干滤波”的另一个优选实施例，即在首先获得K个图像时，以与上述非常类似的方式应用相干滤波器，以求得v₁’(x)，...，以及v_K’(x)，此后，利用具有第M、M-1、...以及M-K+1个图像的K个静态图像，由v(x)＝(v_M(x)，v_{M- 1}’(x)，...，v_M-K+1’(x))建立K维矢量值，然后，对K维矢量值应用实时相干滤波。顺便提一句，其中在实时相干滤波过程中，适于有时通过进行手动设置或自动设置来改变K维像素值矢量v(x)的结构是非常方便的。

(IV对X射线CT射线照相应用本发明(以获得一个被降噪的图像)的情况2)

接着，将说明一个切实可行的例子，在该例子中，与在动态CT内不同，获得被降噪的图像。在此之前，通常，在获得噪声小的一个CT图像情况下，已经考虑到实现两种解决方案，在一个解决方案中，通过重复K次利用普通X射线管电流进行的射线照相获得的K个CT图像具有沿时间方向取的简单平均值(simple average)(作为例子，在采用本实施例的像素值v(x)概念时，简单平均值是{v₁(x)+...+v_K(x)}/K)，结果获得了该单个图像，在另一个解决方案中，利用是普通X射线管电流的K倍的X射线管电流进行一次射线照相获得该单个图像。在这两种解决方案中，噪声方差是利用普通X射线管电流进行一次射线照相获得的单个图像的噪声方差的1/K倍(因为该原因，在前者解决方案中，对噪声进行平均，而在后者解决方案中，X射线辐照量与噪声方差成反比关系，如上所述。)。因此，在此之前，考虑到工作效率，只要射线管的性能容许，通常采用进行一次射线照相的后者解决方案获得单个CT图像，或者在X射线管的性能不能承受大量X射线管电流的情况下，通常采用重复进行K次射线照相的前者解决方案获得单个CT图像。

然而，利用这两个解决方案，噪声抑制得越多，X射线辐照量增加得越多。即，对病人P的照射量增加得越多。

在此，在采用根据本发明的相干滤波器时，可以抑制该单个CT图像内的噪声，而不增加对病人P的照射量。

首先，在该实施例中，在对于病人P完全可以认为没有变化发生的时间周期内，利用比普通X射线辐照量小的X射线辐照获得总共K个CT图像。因此，在这种情况下，对于对K个CT图像确定的像素x，可以建立像素值v(x)：

v(x)＝(v₁(x)，v₂(x)，...，v_k(x))          ...(上述3)

与上面描述的动态相干滤波过程相同。顺便提一句，在这种情况下，X射线辐照量小，因此，包含在每个单个CT图像内的噪声的方差σ²较大。

随后，在这种情况下，对权重函数给与与上面的等式(6)类似的具体形式。

然而，由于在这种情况下将在最后获得的图像是单个CT图像，考虑采用如下所述的解决方案。

例如，可以采用的解决方案是：其中与动态CT射线照相相同，对所有K个图像应用相干滤波器，而且对所获得的K个图像计算加权平均值或简单平均值，以建立一个要求的CT图像的解决方案；或者其中仅对从K个图像中任意选择的单个图像应用相干滤波器的解决方案(仅对选择的单个图像执行图4(c)所示的计算)。

顺便提一句，根据前者解决方案，利用例如下式具体进行计算，以利用通过对所有K个图像应用相干滤波器获得的K个图像，获得单个平均图像：

$v^{*\′}\left(x\right)=\frac{1}{k}\underset{k=1}{\overset{K}{\mathrm{\Σ}}}{v_k}^{\′}\left(x\right)---\left(11\right)$

其中，位于等式(11)右侧的v^*′(x)表示通过进行简单平均计算获得的一个图像上的像素x的像素值(标量值)。此时，包含在v^*′(x)内的噪声的方差(σ²)大致变成：

${\left({\mathrm{\σ}}^2\right)}^{*}=\frac{{\mathrm{\σ}}^2}{K}\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}\[\frac{w\left(p(x,y)\right)}{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}w\left(p(x,y)\right)}{\]}^2---\left(12\right)$

此外，由于这里w(p(x，y))≥0成立，所以等式(12)右侧的求和项始终小于1。因此，根据该实施例的噪声抑制效果优于利用现有技术解决方案获得的噪声抑制效果(σ²/K)。反过来说，如果采用本实施例实现的噪声抑制效果等效于采用现有技术实现的噪声抑制效果，则所需的X射线辐照量可以小于现有技术解决方案所需的X射线辐照量，因此，可以减少对病人P的照射量。

例如图7示出采用该实施例实现的效果。图7(a)示出利用其中在时间方向对根据重复K次射线照相获得的图像简单平均的现有技术解决方案获得的图像，而图7(b)示出根据该实施例，对所有K个图像应用相干滤波器并取所获得的K个图像的简单平均值建立的图像。显然，在图7(b)中充分抑制了在图7(a)中可以看到的随机噪声。(顺便提一句，图7示出在进行3次(K＝3)0.5秒扫描的射线照相、对应于13×13像素的矩形区域被用作像素集N(x)以及设置参数C＝3条件下获得的图像。)

此外，图7(c)示出考虑到由于图7(b)所示的噪声抑制效果非常良好，所以与现有技术的图像相比，而认不出该图像的问题建立的图像，它是通过平均图7(a)和7(b)具体获得的图像。更具体地说，设图7(c)所示图像的像素值为v^**′(x)，则图7(a)所示像素值为：

v(x)＝({v₁(x)+...+v_K(x)}/K)

通过在下面的等式中设u＝1/2，可以获得图7(c)所示的图像：

v^**′(x)＝u·v(x)+(1-u)·v^*′(x)           ...(13)

顺便提一句，在本发明中，该实施例本身(可以获得图7(b)所示的图像)不是限制性的，而且采用可以根据诊断者的偏爱调节例如上面的等式(13)中的u的值(满足0≤u≤1)的建立过程也是更可取的。

(V对X射线CT射线照相应用本发明；进行多级能量射线照相(以降低多种图像的噪声)的情况3)

现在，将说明将本发明应用于所谓“多级能量射线照相”的情况。在此，多级能量射线照相表示了一种其中以改变对X射线管101施加的电压的方式，或者以将由金属薄板等构成的X射线滤波器插在X射线管101与病人P之间的方式，从备方面改变照射X射线的能量分布，而且其中在这样获得的多种能量分布下对同一个病人P进行射线照相的解决方案。在这种情况下，尽管病人P未发生变化，但是图像仍根据X射线的能量分布发生变化，结果，获得了多种图像。

在此，假定在改变能量分布的情况下，对病人P进行K次射线照相，通过对进行K次射线照相获得的K种图像的像素x的标量值进行排列，可以建立像素值v(x)＝(v₁(x)，v_x(x)，...，v_K(x))。因此，利用与在小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的进行动态相干滤波使用的装置相同的装置，可以对像素值v(x)应用相干滤波器，结果，获得了抑制了随机噪声的K种图像。

(VI对X射线CT设备之外的医学图像诊断设备应用本发明的情况)

现在，将说明其中将根据本发明的相干滤波器应用于X射线CT设备之外的医学图像诊断设备，例如，核磁共振成像设备(所谓“MRI设备”)、医疗放射学诊断设备(“闪烁照相机(scintillation camera)、SPECT设备或PET设备”)以及荧光镜设备。

(VI-1对MRI设备应用本发明的情况)

首先说明对MRI设备应用本发明的例子。关于MRI设备，本发明可以应用于至今被称为“平均”、“快速射线照相”以及“三维射线照相”的射线照相解决方案，而且还可以应用于“采用多个脉冲序列的射线照相解决方案”等。

首先，平均(averaging)是一种射线照相解决方案，在该射线照相解决方案中，正如它的名字所表示的那样，目的是通过平均抑制图像噪声，在病人P不发生变化时，重复进行射线照相，而且所获得的多个(K)个图像取其简单平均值，因此可以获得单个图像。在排列K个图像的像素x的标量值时，可以建立像素值v(x)＝(v₁(x)，v_x(x)，...，v_K(x))。此外，关于权重函数w，作为例子可以很好地采用上面的等式(5)。因此，利用与在小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的进行动态相干滤波使用的装置完全相同的装置，可以对像素值v(x)应用相干滤波器，结果，获得了多个抑制了噪声、优于根据平均过程获得的图像的图像。作为一种选择，利用与在小节(IV对X射线CT射线照相应用本发明(以获得一个被降噪的图像)的情况2)描述的处理过程完全相同的处理过程，可以获得一个抑制了噪声的图像。

此外，快速射线照相是一种解决方案，在该解决方案中，在病人P发生变化时，重复对病人P进行射线照相，从而获得表示随时间发生变化的动态图像，而且还将它称为“动态MRI”。此外，在这种情况下，同样可以获得多个(K)静态图像，因此，以与在小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的进行动态相干滤波的方式完全相同的方式，应用相干滤波器，可以抑制噪声，正如平均过程的情况那样。此外，在快速射线照相中，还容许应用与在小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的进行实时相干滤波的处理过程完全相同的处理过程。

三维射线照相是一种一次性获得三维图像的射线照相方法。此外，在这种射线照相方法中，有时执行平均过程(即，其中对于抑制了噪声的对象，重复进行射线照相以取所获得的图像的简单平均值的处理过程)或快速射线照相过程(其中在病人P发生变化时，重复对病人P进行射线照相)。在这种情况下，除了将像素x表示为三维矢量x＝(x，y，z)这点之外，处理过程与上述对平均过程或快速射线照相过程应用相干滤波器的方法完全相同。此外，为了在仅进行一次三维射线照相时应用相干滤波器，可以应用将在后面的小节(VIII根据一个图像建立像素值v(x))描述的解决方案。

同时，在MRI设备中，通过改变脉冲序列，可以使不同的物理参数对图像起作用。“脉冲序列”表示将应用RF＝脉冲的序列与梯度磁场、应用的强度和时间周期以及静止期间组合在一起，通过改变脉冲序列可以获得各种图像。作为通常使用的非常少的脉冲序列的例子，已知有被称为“T2^*-加权图像”、“T1-加权图像”以及“质子密度图像”的脉冲序列，在“T2^*-加权图像”、“T1-加权图像”以及“质子密度图像”的前面添加对图像的对比度影响最大的物理参数。尽管病人未发生变化，但是根据该脉冲序列，可以获得不同图像，结果，实际获得了多种图像。因此，可以以与小节(V对X射线CT射线照相应用本发明；以降低多种图像的噪声的情况3)描述的方式完全相同的方式，应用相干滤波器。

更具体地说，假定在改变脉冲序列时，已经对病人进行了K次射线照相，则通过排列进行K次射线照相获得的图像的像素x的标量值，可以建立矢量v(x)＝(v₁(x)，v_x(x)，...，v_K(x))。因此，在这种情况下，还可以以这样的方式抑制噪声，即利用与动态相干滤波过程中使用的装置完全相同的装置，应用相干滤波器，如平均过程的情况那样。

(VI-2将本发明应用于闪烁照相机、SPECT设备或PET设备的情况)

其次，将说明对闪烁照相机、SPECT设备或PET设备应用本发明的例子。通过将含有放射性同位素(以下称为“RI”)的药剂注入受治疗者体内，然后检测因为RI的自然衰变发出的γ射线，这种设备使在受治疗者体内药剂(放射性药剂)的分布转化为图像。特别是在受治疗者是活体的情况下，使用该设备的目的是获得表示活体的生化代谢程度和血液循环程度的所谓“功能映像”。闪烁照相机获得药剂的分布作为二维图像，好象利用荧光镜对受治疗者进行检查。此外，利用计算的层析X射线照相法，SPECT设备和PET设备获得表示药剂的三维分布的三维图像。

在闪烁照相机、SPECT设备或PET设备中，检测到的γ射线量通常非常少，因此在图像上出现因为RI的量子力学属性导致的许多光子的扩散而产生的所谓“光子噪声”。

至此，为了抑制光子噪声又不降低图像的分辨率，唯一的方法就是增加待检测的γ射线量。具体地说，存在增加注入的RI数量的方法，和延长射线照相时间周期的方法。

在可以增加注入的RI数量的情况下(容许对例如进行医学研究的实验动物增加剂量)，通过增加注入的RI数量，可以降低光子噪声，但是这种方法导致增加对操作员的照射量。此外，在诊断人体的情况下(＝人体就是受治疗者)，首先不容许不受限制地如上所述增加RI数量(受法律制约)。此外，在利用包含在自然物体中的RI发出的放射线进行射线照相的情况下，不能增加放射线数量。

另一方面，如果采用延长射线照相时间周期的方法，则为了实现足够的降噪效果，要求受治疗者非常长时间不能动(例如，几个小时)。

在活体(特别是人体)是受治疗者的情况下，这种要求难以实现。另一个问题是，因为射线照相的时间周期长，而降低了射线照相设备的处理能力(或性能价格比)。因此，实际上，只能设置受容许时间周期限制的适当射线照相时间周期(例如，在一个小时之内)，而强噪声不可避免地出现在该图像上。因此，仅通过平滑该图像，抑制噪声。因此，可以认为破坏图像的分辨率(伴随平滑产生的缺点)是不可避免的。

在这种情况下，根据本发明，可以执行如上所述的处理。在等效于以上说明的适当射线照相时间周期的时间周期内，可以在短时间周期内进行多次射线照相，从而获得多个含有较多噪声的图像。随后，应用与小节(IV对X射线CT射线照相应用本发明(以获得一个被降噪的图像)的情况2)说明的处理过程完全相同的处理过程，因此可以在不破坏空间分辨率的情况下，获得被更多抑制了噪声的图像。此外，可以在不增加RI数量的情况下，进行该处理。

此外，还是在闪烁照相机、SPECT设备或PET设备中，存在一种射线照相解决方案，在该射线照相解决方案中，与在X射线CT等中相同，对同一个受治疗者或病人重复进行射线照相，对病人体内的药剂分布随时间发生的变化进行跟踪，即，被称为  “动态闪烁图”、“动态SPECT”或“动态PET”。此外，在这种情况下，可以应用与小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的动态相干滤波或实时相干滤波的处理过程完全相同的处理过程，而且可以有效抑制图像噪声。

此外，在PET设备或SPECT设备中，有时采用所谓“多窗***线照相”，在所谓“多窗***线照相”中，利用RI发出其能量为其种类所特有的γ射线这个事实，将多种RI同时注入受治疗者或病人，以便一次性对多种RI获得相应分布图像。多窗***线照相基于一种其中对所检测的γ射线的能量进行测量，而且根据能级是否在人为设置的数值范围(能量窗口)内，识别发出γ射线的多种RI的技术。

更具体地说，设置适当数量的能量窗口(例如，窗口W1和W2)，分别形成各图像，例如基于对应于其中一个窗口(例如，窗口W1)的γ射线的图像以及基于对应于另一个窗口(例如，窗口W2)的γ射线的图像，从而通过进行一次射线照相获得如上所述的多个RI分布图像。

也就是说，利用成像设备对同一个受治疗者进行射线照相，而且在不同的处理条件下，获得多种图像。更一般地说，多窗***线照相是一种“其中利用K种RI，通过进行一次射线照相，获得K种互相不同图像的射线照相解决方案”。因此，可以以与小节(V对X射线CT射线照相应用本发明；以降低多种图像的噪声的情况3)描述的方式完全相同方式，应用相干滤波器。

因此，还是在这种情况下，对于这K种图像上的像素x，可以建立其像素值v(x)＝(v₁(x)，v₂(x)，...，v_k(x))。此外，尽管包含在像素值v(x)内的噪声方差随像素x的位置和图像的种类的不同而不同，但是可以认为它与对各种图像检测的γ射线的总量大致成正比，因此可以根据该属性估计它。因此，如上所述，通过应用本发明，相干滤波器还可以应用于“多窗***线照相”，其结果是可以获得抑制了随机噪声的K种图像。

顺便提一句，在医疗放射学诊断设备中，存在有效提高分辨率的重构方法(扇形射束的MU-EN重构方法、扇形射束的超高分辨率重构方法或平行孔光管(collimator)等)。从原理上说，并不完全局限于将这些方法与根据本发明的相干滤波器组合。因此，可以产生具有高分辨率、而且充分抑制了其随机噪声的图像。

(VI-3将本发明应用于荧光镜设备的情况)

下面，将说明对荧光镜设备应用本发明的例子。假定在此所称“荧光镜设备”是一种其中在以低剂量对病人连续照射X射线时，每个预定时间周期(例如，1/24秒)重复进行一次射线照相，而且其中各“帧”(每“帧”(frame)分别表示一次射线照相)所获得X射线图像被连续显示，从而观察所谓“荧光镜图像(动态图像)”的设备。因此，包括被称为“C臂”或“Ω臂”结构的设备通常相当于荧光镜设备。此外，如下所述的方面自然可以应用于上述X射线CT设备100，只要安装了上述功能部件即可。

在这种荧光镜设备中，荧光镜图像含有大量因为上述低剂量产生的光子噪声而导致的随机噪声。在为了解决该问题而充分增加X射线辐照量时，不仅对病人，而且对检查者显著增加照射量，其原因是连续进行X射线辐照(irradiation)。

在应用根据本发明的相干滤波器时，可以解决这些问题，而且可以提供满意抑制了其噪声的荧光镜图像。

具体地说，为了便于说明，对各帧获得的图像连续指定为第1，2，...和M个图像。将最新的K个图像，即第M，M-1，...，和M-K+1各图像的图像存储到存储器中，而且对这K个图像执行与小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的动态相干滤波或实时相干滤波的处理过程完全相同的处理过程。

顺便提一句，除了医学图像诊断设备，上述方法还可以应用于用于进行所谓“无损检验”的“荧光镜设备”。

此外，本发明同样可以应用于利用超声波诊断设备获得的层析X射线照片(tomogram)。还是在这种情况下，可以确保在对受治疗者连续进行射线照相时，连续产生图像。因此，可以利用与荧光镜设备的结构相同的结构应用相干滤波器，而且可以抑制噪声。

(VII将本发明更一般地应用于彩色图像等的情况)

本发明还可以应用于不是上述医学图像诊断设备、用于处理图像的“通用设备”或“一般图像”。现在，将说明对这种通用设备或图像应用本发明的例子。顺便提一句，处理一般“彩色图像”的情况属于特别重要的领域，所以在此主要对其进行说明。

在彩色图像中，如上所述可以对构成该图像的每个像素x，建立其元素是光的3原色亮度的三维矢量，即：

v(x)＝(v_R(x)，v_G(x)，v_B(x))    ...(3.1)

(其中V_R、V_G和V_B表示作为光的3原色的R(红)、G(绿)和B(蓝))。从图2中可以看出，这相当于“根据多种图像建立像素值v(x)”的情况。

此外，在利用例如由CCD成像器件、CMOS成像器件等构成的数码相机拍摄彩色图像的情况下，利用进入像素x的光量确定包含在每个像素的像素值中的噪声方差。具体地说，标量值v_R(x)，v_G(x)和v_B(x)与包含在标量值内的噪声方差σ² _R(x)、σ² _G(x)和σ² _B(x)大致分别成正比关系。利用该属性，可以求得噪声方差的估计值。因此，给出与上面的等式(6)中的权重函数具体形式相同的权重函数w的具体形式，而且该权重函数w的具体形式用于建立上述等式(4)的相干滤波器，从而获得抑制了噪声的像素值v’(x)＝(v’_R(x)，v’_G(x)，v’_B(x))。

例如，如果为了全面理解上述红光、绿光和蓝光之外的红外线和紫外线而从最广意义上理解单词“彩色”，还可以将等式(3.1)表示为更高阶的矢量或表示为与彩色R、G和B无关的矢量。作为例子，已知利用矢量值的像素值表示彩色图像的方法并不始终仅局限于3原色，而且本发明可以很好地应用于采用这种矢量表示而不使用3原色的情况。此外，可以考虑包括紫外线VU的4原色、具有分离的紫色VP的5原色(顺便提一句，这相当于蝴蝶图像(vision ofbutterflies))等的情况。对于例如5原色情况，将像素x表示为5维矢量值V(x)＝(V_R(x)，V_G(x)，V_B(x)，V_U(x)，V_P(x))，因此可以象在上述结构中那样，应用相干滤波器。

此外，在其中为了测量例如受治疗者的反射光谱，利用单色照相机(monochroimatic camera)，在用于照射受治疗者的照射光的波长是逐步或连续转换时，以各照射波长对受治疗者重复进行射线照相，或者通过逐步转换滤色镜对受治疗者重复进行射线照相，从而以多种波长进行射线照相的射线照相设备中，可以根据大量图像(每个图像分别是一个波长或一个颜色的图像)建立像素值v(x)。即，利用与小节(V对X射线CT射线照相应用本发明以降低多种图像的噪声的情况3)中的方法完全相同的方法，可以应用本发明，结果是获得了降低了随机噪声的多种图像。

此外，许多这种成像设备通常进行一次射线照相以获得一个图像。然而，在应用本发明时，通过在预定时间周期内重复进行短时间周期的射线照相可以更好地获得多种(K)图像，如在小节(IV对X射线CT射线照相应用本发明(以获得一个被降噪的图像)的情况2)对“X射线CT射线照相”所做的说明。其原因是，(尽管因为短时间射线照相导致光量变少，使得K个图像含有大量随机噪声)通过对这K个图像应用相干滤波器，可以建立更强烈抑制了随机噪声的一个图像。此外，在其中K个图像是例如彩色图像的情况下，可以将它们表示为“多个其每个分别包括对应于光的3原色获得的多种(3种)图像的组合(K个组合)中的”图像。更具体地说，可以对构成每个图像的各像素x建立其各元素是光的3原色亮度的三维矢量，因此可以对所有K个图像建立3K维矢量的矢量值v(x)。相干滤波器可以很好地应用于3K维矢量的像素值。

此外，在用于拍摄动态图像的、诸如摄像机的成像设备中，可以应用在小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的动态相干滤波或实时相干滤波。此外，通过将动态图像作为“多个其每个分别包括对应于光的3原色获得的多种(3种)图像的组合(K个组合)中的”图像进行处理，如上所述，可以应用相干滤波器。

顺便提一句，对于不是以上说明的数码相机的摄像机、红外线摄像机、快速照相摄像机(camera)等(在以下通常称为“成像设备”)拍摄的图像，可以以上述同样的方式应用本发明。特别是，红外线摄像机或快速照相摄像机是在每个帧的光量少的情况下，首先使用的摄像机，而且图像中含有大量随机噪声。然而，它无助于破坏空间分辨率，而且无助于破坏时间分辨率，特别是在快速照相摄像机中。因此，红外线摄像机或快速照相摄像机是一种最适于应用本发明的设备。

因此，还是在这种情况下，利用这种照相方法可以获得抑制了更多噪声的图像。

同时，在不进行“成像”，但是可以根据记录在任何记录介质上的信号或发送的信号再现图像的设备，例如磁带录像机或DVD播放器中，在对要再现的图像应用根据本发明的相干滤波器时，可以获得抑制了随机噪声的再现图像。

此外，在在这种设备中产生并显示静态图像(所谓“暂停图像”或用于打印的图像)的情况下，通过应用相干滤波器可以获得抑制了随机噪声的静态图像，而且该图像还适于打印等。为此，可以采用小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的方法、小节(V对X射线CT射线照相应用本发明，进行多级能量射线照相，以降低多种图像的噪声的情况3)描述的方法、将在下面的小节(VIII根据一个图像建立像素值v(x))描述的方法或这些方法的组合。

(VIII根据一个图像建立像素值v(x))

这样到此为止描述的应用例子对应于图2所示分类中的那些例子，在这些例子中，“根据多个图像”以及“根据多种图像”分别建立的像素值v(x)。在本发明中，除了建立像素值v(x)的两种方法，图2还示出“根据一个图像建立像素值v(x)”的方法，而且可以如下所述实现该建立方法。

在一个二维图像，例如利用X射线CT设备100、MRI设备等获得的一个CT图像中，或者在利用成像设备等进行单色照相获得的一个图像中，将像素x的值仅由一个标量值构成。换句话说，使作为矢量值的像素值v(x)变成单维矢量。

在这种情况下，应用根据本发明的相干滤波器的一种方法是将像素值v(x)本身看作“单维矢量”。即，在上面的等式(3)中，设置K＝1以获得：

v(x)＝v(x)                             ...(3.2)

此后，可以以上述过程中的同样方式建立相干滤波器或抑制图像噪声。本发明还包括作为“矢量值”的“像素值”v(x)是单维矢量，即标量值的情况。

此外，可以利用例如下面的等式计算这种情况下的权重w3：

$w3\left(p(x,y)\right)=$

在此，D的属性类似于等式(6)中的参数C的属性。即，在D小时，权重w3(p(x，y))变小，而在D大时，权重也变大。

作为例子，基于根据该权重w3(p(x，y))和上面的等式(4)计算的像素值v’(x)的图像变成图8所示的图像。图8(a)、(b)、(c)和(d)是分别示出血管3D采样的例子的增强图像，而图8(e)是作为另一个采样例子通过采样人脑的灰质(gray matter)获得的图像。图8(a)示出将被相干滤波器处理的CT图像(原始图像)。图8(b)和(d)示出利用3×3像素的矩形区域作为等式(6)中的像素集N(x)，应用相干滤波器获得的结果。图8(b)所示的结果对应于其中在上面的等式(14)中设置D＝9的情况，而图8(d)所示的结果是对通过应用相干滤波器获得的(c)所示图像进一步进行阈值处理获得的图像。此外，图8(c)所示的图像是对(a)所示的原始图像进行阈值处理获得的图像。

在将图8(a)与8(b)所示的图像进行比较时，可以明白，由于应用了相干滤波器，所以与图8(a)所示的图像相比，图8(b)所示的图像被适当抑制了噪声，因此可以确保满足诊断的图像质量。

在将图8(c)与8(d)所示的图像进行比较时，可以明白，通过仅对原始图像进行阈值处理获得的图8(c)所示图像难以根据阈值进行血管采样，而在应用相干滤波器后进行阈值处理获得的图8(d)所示图像允许根据阈值进行血管采样，而不会忽略非常微细的血管。

此外，根据图8(e)，可以明白，已经采样了灰质区域。

根据该解决方案，可以真实获得图8所示的效果。然而，由于单维矢量值v(x)＝v(x)用于测量像素x与y之间的拟和优度，所以拟和优度的识别性能较低，换句话说，仅获得统计数字不确定的临界概率p(x，y)。因此，不能说不破坏分辨率而抑制随机噪声的效果非常好。因此，在这种情况下，应该更优先采用如下所述的另一种解决方案。其中利用包括像素x附近的一些像素的像素集Z(x)建立多维矢量，从而将该多维矢量用作像素x的像素值v(x)的解决方案是有效的。作为例子，如图9所示，我们研究其中利用包括4个围绕像素x＝(x，y)并与其邻接的像素(x+1，y)、(x-1，y)、(x，y+1)和(x，y-1)的像素集Z(x)(以下将作为围绕像素x的像素集z(x)的各元素的像素表示为“z₁(x)，z₂(x)，z₃(x)和z₄(x))的5维矢量，这4个像素的标量值被加到像素x本身的标量值，即v(x)＝v(x，y)。即：

v(x)＝(v(x，y)，v(x+1，y)，v(x-1，y)，v(x，y+1)，v(x，y-1))

                                                            ...(3.3)

对于构成图9所示图像G的所有像素，同样可以建立这种矢量(顺便提一句，图9仅以放大比例示出整个图像的一部分)。

这样，利用根据一个图像建立的多维矢量的像素值v(x)，可以建立根据本发明的相干滤波器，因此可以以比采用单维矢量像素值的情况的精度高的精度，判定拟和优度。结果，可以更强烈抑制随机噪声，而不破坏图像的分辨率。

顺便提一句，在这样建立像素值的情况下，上述用于量化拟和优度的虚假设H、“在消除了像素x和y的噪声后，像素x和y具有同样的像素值”仅表示“除去噪声，像素x和周围像素z₁(x)至z₄(x)总共5个像素以及像素y和周围像素z₁(y)至z₄(y)总共5个像素与相应像素分别具有同样的像素值”的虚假设H’。

更具体地说，在如图9所示，存在于像素集z(x)(即，例如，如上所述围绕像素x的区域)内的特定像素x和特定像素y均存在于例如该图所示图像G上的对象E1的区域E1内的情况下，利用例如等式(9)计算的临界概率p(x，y)的值变小，而利用例如等式(6)计算的权重w(p(x，y))变大。其原因是，在这种情况下，像素z₂(x)和z₄(x)以及像素z₂(y)和z₄(y)位于特定对象的图像E1区域内，而像素z₁(x)和z₃(x)以及像素z₁(y)和z₃(y)位于对象E2的图像内，因此利用等式(3.3)表示的矢量值V(x)与V(y)相似。

现在，在上面的实施例中，在特别注意位于该对象的图像E1的边缘上的像素x时，可以推测通常不存在非常大量的对于其在放弃虚假设H’的情况下临界概率p(x，y)低的像素y(作为例子，在图9中，即使对于均存在于对象E1的图像边缘E1E内的像素x和像素y_p，在放弃虚假设H’情况下的临界概率p(x，y)仍比该图所示的像素x和y的临界概率高)。

更一般地说，对于存在于特定对象的图像边缘上的像素x和大多数y∈N(x)，利用例如等式(9)计算的临界概率p(x，y)变大，而利用例如等式(6)计算的权重w(p(x，y))变小。因此，即使进行加权平均，都可以使平滑效果变弱，因此基于利用例如等式(4)计算的相干滤波器的噪声抑制效果不够好。

因此，以下将说明进一步改进上述解决方案获得的另一种解决方案。在该解决方案中，像素值v(x)被构建为：

v(x)＝(v(x，y)，v₁，v₂，v₃，v₄)           ...(3.4)

其中v₁至v₄是根据预定规则(例如，以大小值顺序)排列上面的等式(3.3)中的v(x+1，y)、v(x-1，y)、v(x，y+1)以及v(x，y-1)获得的。

这样，虚假设H’转变为更适度的虚假设H”，即“除去噪声，像素x和周围像素z₁(x)至z₄(x)总共5个像素与像素y和周围像素z₁(y)至z₄(y)总共5个像素分别具有同样的像素值[但是容许改变像素z₁(x)至z₄(x)以及像素z₁(y)至z₄(y)的顺序]”(已经增加了括号[]内的内容)。

此外，在这种情况下，在存在于例如图10所示对象的图像E2的边缘E2E上的像素x与y之间的关系中，在放弃虚假设H”情况下的临界概率p(x，y)变低(因为像素z₂(x)、z₃(x)和z₄(x)以及像素z₁(y)、z₂(y)和z₄(y)与上面的等式(3.4)中的v₁至v₃一致，而像素z₁(x)和像素z₃(x)同样与V₄一致)，结果是，其噪声抑制效果比基于根据上面的等式(3.3)建立的像素值v(x)的解决方案中的噪声抑制效果要好，而且不破坏空间分辨率。

如上所述，在本发明中，可以根据一个图像建立作为矢量值的像素值v(x)，因此利用相干滤波器，可以抑制噪声，而不破坏空间分辨率。由于该解决方案只需要一个图像，因此，通常，它基本上可以应用于任何图像。即，该解决方案不仅可以应用于上述各种应用例子(X射线CT设备，医学图像诊断设备、彩色图像等)，而且可以应用于所有图像。此外，当然，即使在图像是动态图像(可以看作多个静态图像，以上说明了几次)的情况下，该解决方案(其中“根据一个图像建立像素值v(x)”)仍可以应用于每个静态图像。此外，在对三维图像或多维图像采用该解决方案的情况下，由例如像素x₁本身和位于接触像素x₁位置的6个像素可以构成像素集z(x)。可以非常容易地建议这些扩展方法，而且这些扩展方法是显而易见的。

顺便提一句，在本发明中，利用基于如图2所示的“根据多个图像”、“根据多种图像”以及“根据一个图像”的每种解决方案建立的像素值v(x)，当然可以建立相干滤波器，但是还可以根据通过适当组合各解决方案建立的、被表示为高维矢量值的像素值v(x)建立相干滤波器。

正如在小节VII中所述的那样，作为例子，假定存在关于同一个受治疗者的两个彩色图像，对于特定像素x，因为像素值v(x)是6维矢量，所以利用两个图像获得位于上面的等式(2.2)的右侧的标量值v_R(x)、v_G(x)和v_B(x)，作为像素x的像素值v(x)：

v(x)＝(v_R1(x)，v_G1(x)，v_B1(x)，v_R2(x)，v_G2(x)，v_B2(x))

                                                   ...(3.5)

在以这种方式使用利用高维矢量表示的像素值v(x)时，可以对临界概率p(x，y)计算更可靠的值，它是利用例如等式(9)和(6)计算的。因此，可以更强烈抑制随机噪声，而不破坏图像的分辨率。

(IX最通用形式相干滤波器的实施例)

现在，除了顺着图2所做的以上说明之外，将说明在上面的等式(4)中引入的更通用形式的相干滤波器。

在通用形式(等式(4))的相干滤波器中，或者在各种应用例子中，采用了在开头说明的“在已经去除像素x和y的噪声时像素x和y具有同样的值”的，或者“除了因为两个像素的噪声产生的差别v(x)＝v(y)”的虚假设H，但是可以给出该虚假设的更通用表述，如下所述。

引入包括未知的M维参数a＝(a₁，a₂，...，a_M)(其中M＜K)的模型函数f(a，v(y))，因此，对于表示像素位置的两个特定矢量x和y，建立虚假设H₀“v(x)＝f(a，v(y))+ξ(其中ξ与(已知的)概率分布一致)”。此时，利用适当的拟和方法(例如，如下所述的最小二乘方法)，根据v(x)和v(y)确定参数a的估计值a~，因此，换句话说，可以将上述虚假设表示为其等效命题H₀“v(x)＝f(a~，v(y))+ξ(其中ξ符合(已知的)概率分布)”。在此，上述f(a~，v(y))是这样的，即通过最佳调节函数f描述的模型容许的自由度(即，参数a)，使像素y的像素值v(y)进行变换以与像素x的像素值v(x)具有最高拟和优度。此外，这种最佳参数是a~。

因此，可以利用已知的概率分布(假定如果虚假设正确，则ξ与其一致)具体计算临界概率p(x，y)，因此可以计算权重w(p(x，y))。因此，利用下面的等式可以计算加权平均值v’_K(x)：

$v^{\′}k\left(x\right)\frac{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}w\left(p(x,y)\right)f(a^{~},v\left(y\right))}{\underset{y\∈N\left(x\right)}{\mathrm{\Σ}}w\left(p(x,y)\right)}---\left(15\right)$

这是该实施例中的最通用形式相干滤波器。顺便提一句，在将上面的等式(4)与该等式(15)进行比较时，可以明白，前者仅是将下式实际代入后者获得的更具体形式：

f(a~，v(y))＝v(y)             ...(16)

(IX-1根据上面的等式(15)的应用例子，根据动态图像建立时间-密度曲线)

将说明采用等式(15)所示形式的相干滤波器的本发明的另一个应用例子。在此，该应用例子是一个其中在量化测量用于处理利用医学图像诊断设备获得的动态CT图像的时间-密度曲线的过程中，试图通过应用本发明获得有效抑制了随机噪声的时间-密度曲线的实施例。

首先，在此所称的“动态图像”是对同一个受治疗者重复进行射线照相获得的一系列图像(第1、2、...K个图像)与各图像的射线照相时间(t₁，t₂，...，t_K)的组合。它是例如利用医学图像诊断设备，即X射线CT设备100获得的动态CT图像，或者利用荧光镜设备、闪烁照相机、MRI设备、SPECT设备或PET设备、超声波诊断设备等进行射线照相。此外，利用通用照相设备，即数码相机、摄像机、红外线摄像机、快速照相摄像机等重复进行照相获得的一系列图像可以很好地构成动态图像。

此外，在此所称的“时间-密度曲线”是表示动态图像中指定部分的图像密度值随时间发生变化的曲线。特别是在医学图像诊断设备中，为了检查例如血流流动、人体组织内的代谢功能等的详情，测量人体特定组织内的造影剂密度随时间的变化作为时间-密度曲线。此外，在进行天文观测等的过程中，为了分析例如指定天体的亮度变化，利用时间-密度曲线。在进行更正式的说明时，将“时间-密度曲线”表示为{<t_k，d_k>(k＝1，2，...，K)}对的序列，其中d_k表示在时间t_k时某一部分的密度值。此外，在时间-密度曲线的许多应用中，不是总需要密度值d_k的绝对值，但是它足以仅获得参考第一图像1计算的增量(d_k-d₁)。此外，在许多这种应用中，它足以仅获得与增量(d_k-d₁)成正比的数据A(d_k-d₁)(其中A表示未知的比例系数)。因此，在这种情况下，{<t_k，A(d_k-d₁)>(k＝1，2，...，K)}对的序列就是要获得的时间-密度曲线。

为了获得这种时间-密度曲线，从理论上说，在构成动态图像的每个图像k(k＝1，2，...，K)中，利用包括在要对其测量时间-密度曲线的部分内的像素x的标量值v_k(x)，可以建立{<t_k，v_k(x)>}对或{<t_k，A(v_k(x)-v₁(x))>}对的序列。

然而，在实际使用中，存在的问题是，由于利用医学图像诊断设备等获得的动态图像中含有随机噪声，所以不能精确获得最初要测量的时间-密度曲线。

此外，在实际使用中，在动态图像中产生所谓“部分体效应(partial volume effect)”。部分体效应是一种受治疗者或病人内的微小对象的图像被表示为动态图像上的少量像素，而且这些少量像素受受治疗者体内邻近对象的干扰，因此这些少量像素的像素值仅显示微小波动的现象(但是它们与最初要测量的密度值的波动成正比)。换句话说，这些少量像素的像素值仅含有微小信号。因此，在出现部分体效应的情况下，{<t_k，v_k(x)>}对序列在任何像素x均具有非常低的信号电平，而且它受不是最初要测量的组织内的密度值的变化的干扰。此外，还存在随机噪声。这样产生的问题是，不能精确获得最初要测量的时间-密度曲线{<t_k，d_k>}。

因此，为了抑制随机噪声，至今采用基于时间平均的方法，如下所述。首先，对包含在要测量部分内的像素x建立{<t_k，v_k(x)>(k＝1，2，...，K)}对序列。随后，利用例如下式计算时间平均值v_k ^#(x)：

$v_k\#\left(x\right)=\frac{\underset{j=m}{\overset{n}{\mathrm{\&Sigma;}}}{v}_{k+j}+\left(x\right)G_j}{\underset{j=m}{\overset{n}{\mathrm{\&Sigma;}}}{G}_j}---\left(17\right)$

此外，利用下式计算对应于上述时间平均值的备时间的时间平均值：

$t_k\#\left(x\right)=\frac{\underset{j=m}{\overset{n}{\mathrm{\&Sigma;}}}{t}_{k+j}\left(x\right)G_j}{\underset{j=m}{\overset{n}{\mathrm{\&Sigma;}}}{G}_j}---\left(18\right)$

(其中m＜n表示适当整数，而G_j表示适当加权系数)。因此，可以建立并采用时间-密度曲线{<t_k ^#，v_k ^#(x)>(k＝1，2，...，K)}。

作为一种选择，为了抑制随机噪声，至今采用基于空间平均的方法，如下所述。建立像素集R，像素集R与要测量的部分(即，动态图像上的ROI(感兴趣区域))大致对应。在第k个图像上，利用下面给出的等式(19)，计算包括在该像素集内的所有像素x∈R的标量值v_k(x)的空间平均值。此后，建立{<t_k，v_k(R)>(k＝1，2，...，K)}对序列，并将它用作时间-密度曲线。

$v_k\left(R\right)=\frac{\underset{x\&Element;R}{\mathrm{\&Sigma;}}{v}_k\left(x\right)}{\underset{x\&Element;R}{\mathrm{\&Sigma;}}1}---\left(19\right)$

此外，还采用其中利用时间平均和空间平均更多抑制随机噪声的解决方案。

然而，这些现有技术的解决方案存在的问题是，在执行时间平均时，破坏时间分辨率，而在执行空间平均时，最初要测量的部分之外的各部分的密度随时间的变化混入测量值。

为了解决这些问题，而且为了获得更精确的时间-密度曲线，可以根据小节(IX最通用形式相干滤波器的实施例)应用根据本发明的相干滤波器。

首先，将说明用于该实施例的相干滤波器的虚假设。假定要测量部分的真实时间-密度曲线是{<t_k，d_k>(k＝l，2，...，K)}，则在意在测量作为真实时间-密度曲线的线性变换的{<t_k，A(d_k-d₁)>(k＝1，2，...，K)}(其中A表示未知系数)的情况下，建立大致对应于要测量部分的像素集R。对于作为像素集R的元素的任何像素x∈R，如果满足条件Q：“像素x很好地反映真实时间-密度曲线，而且几乎不受其他部分内的密度随时间变化的干扰”，则考虑到部分体效应和随机噪声对像素值(是矢量值)v(x)＝(v₁(x)，v₂(x)，...，v_k(x))的干扰，可以认为下面给出的等式(20)成立：

        v_k(x)＝p(x)d_k+q(x)+γ_k(x)

               (k＝1，2，...，K)                ...(20)

在此，p(x)和q(x)表示对各像素x不同、但是又不根据第k个图像(即，射线照相时间t_k)发生变化的未知系数，通过对部分体效应建模可以获得它们。此外，通过对随机噪声建模获得γ_k(x)。尽管γ_k(x)的值对各像素x和各第k个图像不同，但是其期望值为0(零)，而且其静态分布既不取决于像素x又不取决于第k个图像。(顺便提一句，在上面的描述中，为了简洁起见，举例说明了静态分布是平均值0和方差σ²的正态分布的情况，但是仅要求统计分布大致接近任何已知的分布，而且如何修改下面的实施例以适应该情况是显而易见的。)

根据上面的假定，如果对于作为像素集R的任何两个元素x和y，命题“像素x和y均满足条件Q(如上所述)”成立，则它可以证实下式等式关系成立：

      v_k(x)＝a₁v_k(y)+a₂+ξ_k

             (k＝1，2，...，K)                 ...(21)

在此，a₁和a₂表示对于像素的每个组合x、y不同，但是不根据第k个图像(即，射线照相时间t_k)发生变化的未知系数。此外，ξ_k表示随机噪声，而且对于各像素组合x、y和各第k个图像，其值不同，但是其期望值为0。

如以下所述推导出等式(21)。在通过将y代入等式(19)内的x获得的等式：

     v_k(y)＝p(y)d_k+q(y)+γ_k(y)             ...(22)

进行转换时，获得下面的等式：

$d_k=\frac{v_k\left(y\right)-q\left(y\right)-{\mathrm{\&gamma;}}_k\left(y\right)}{p\left(y\right)}---(22\text{'})$

通过将该等式代入等式(20)，获得下面的等式：

$v_k\left(x\right)=\frac{p\left(x\right)}{p\left(y\right)}{v}_k\left(y\right)+(q\left(x\right)-\frac{p\left(x\right)}{p\left(x\right)}q\left(y\right))+({\mathrm{\&gamma;}}_k\left(x\right)-\frac{p\left(x\right)}{p\left(x\right)}{\mathrm{\&gamma;}}_k\left(y\right))---\left(23\right)$

因此，通过设：

推导出等式(21)。其中，等式(24)中的a₁和a₂代表表示部分体效应的参数，而ξ_k同样表示随机噪声。

以上表明命题“像素x和y均满足条件Q”等效于虚假设H₀’“v_k(x)＝a₁v_k(y)+a₂+ξ_k(k＝1，2，...，K)成立”。

  接着，将说明一种将虚假设H₀’“v_k(x)＝a₁v_k(y)+a₂+ξ_k(k＝1，2，...，K)成立”变换为与其实际等效而且可以实际测试的命题形式。在利用严格的数学表达式重新表示该虚假设时，它就变成虚假设H₀’“存在特定常数a₁和a₂，而且ξ_k＝v_k(x)-a₁v_k(y)-a₂(k＝1，2，...，K)符合平均值0和方差(σ²h(a₁)的正态分布)”。在此，系数h(a₁)是：

     h(a₁)＝1+a₁ ²                            ...(25)

(根据给出a₁和ξ_k的定义的等式(24)以及随机变量的方差所具有的一般属性，可以立即推导出等式(25)。)此外，在实际使用过程中，利用例如在上面的小节(III将本发明应用于X射线CT射线照相、动态CT射线照相(以消除多个图像中的噪声)的情况1)描述的方法，可以以令人满意的精度、方便地估计上面的方差σ²。

根据上述说明可以明白，如果可以确定常数a₁和a₂，则可以测试虚假设H₀’。在此，实际上足以获得这些常数的最佳估计值a₁~和a₂~。

可以直接利用已知的拟和方法计算常数a₁和a₂的这种最佳估计值。因此，在下面的说明中，作为这种拟和方法的具体典型例子，将概括说明采用线性最小二乘方法的情况。在对该实施例应用线性最小二乘方法的过程中，以下仅将上述虚假设中的ξ_k的平方和定义为S(a)：

$\left(a\right)=\underset{k=1}{\overset{K}{\mathrm{\&Sigma;}}}\{v_k\left(x\right)-a_1{v}_k\left(y\right)-a_2{\}}^2---\left(26\right)$

平方和S(a)的值取决于常矢量a＝(a₁，a₂)，即常数a₁和a₂的值。在计算根据其平方和S(a)取最小值的常矢量a时，可以获得无偏估计常数a₁和a₂意义上的最佳估计值a₁~和a₂~。顺便提一向，可以将各种已知的方法用作线性最小二乘方法的具体计算方法，此外，已知任何的计算方法均非常简单，而且所需计算时间非常短。

这样，就计算了常数a₁和a₂的最佳估计值a₁~和a₂~，因此，可以具体计算下面的等式定义的残数(residual)：

     r_k~(x，y)＝v_k(x)-a₁~v_k(y)-a₂~               ...(27)

因此，换句话说，利用残数r_k~，可以将上述虚假设H₀’表示为与其实际等效的虚假设H₀”“r_k~(x，y)(k＝1，...，K)符合平均值0和方差(1+(a₁~)²)σ²的正态分布”。这是可以进行要实际执行的测量计算的具体命题。

此外，换句话说，通过引入矢量表达式表示虚假设H₀’：

(其中矢量a和ξ取决于像素的组合x，y)并采用利用下面的等式确定的矢量值函数f：

     f(a~，v(y))＝a₁~v(y)+a₂~                    ...(29)

然后，虚假设H₀”变成“v(x)＝f(a~，v(y))+ξ)(其中ξ符合平均值0和方差(1+(a₁~)²)σ²的正态分布)”，而且它具有与在小节(IX最通用形式相干滤波器的实施例)描述的虚假设H₀的形式完全相同的形式。因此，显然，该实施例是一种根据本发明的相干滤波器的变换实施例。顺便提一句，在此，函数f(a~，v(y))表示对于像素y的像素值v(y)，已经对其进行最佳调节和变换从而与像素x的像素值v(x)具有最高拟和优度的、用于表示部分体效应的参数a。

接着，将说明一种其中在该实施例中，利用具有上述虚假设H₀”的相干滤波器获取时间-密度曲线的方法。对于大致对应于要测量的部分的像素的像素集R，对包括在该像素集R内的、与作为像素集R的元素的所有像素y∈R相关的特定像素x∈R进行如下所述的计算。利用上述方法实际计算残数r_k~(x，y)(k＝1，...，K)。随后，根据等式(6)具体计算放弃上述虚假设H₀”“r_k~(x，y)(k＝1，...，K)符合平均值0和方差(1+(a₁~)²)σ²的正态分布”情况下的临界概率p(x，y)或权重w(p(x，y))。此外，根据下面的等式(15’)计算加权平均值v’_k(x)，以建立像素x的时间-密度曲线{<t_k，v’_k(x)-v’₁(x)>(k＝1，2，...，K)}。

$\left(x\right)=\frac{\underset{y\&Element;R}{\mathrm{\&Sigma;}}w(p(x,y)f(a^{~},-v\left(y\right))}{\underset{y\&Element;R}{\mathrm{\&Sigma;}}w\left(p(x,y)\right)}---(15\text{'})$

这样获得的时间-密度曲线是近似于作为像素x的真实时间-密度曲线{<t_k，d_k>}的线性变换的{<t_k，A(d_k-d₁)>}(其中A表示未知系数)的测量值，而且利用基于等式(15)的加权平均值的作用，抑制了其随机噪声。此外，根据等式(15)可以明白，在根据等式(15)的计算过程中使用的其他像素y的像素值矢量受到被补偿的部分体效应的干扰。此外，该实施例具有属性“根本不使用时间平均值，而利用基于与像素x的拟和优度的权重计算空间平均值”，这是相干滤波器的共同特征。因此，利用该实施例可以获得不破坏时间分辨率、其受部分体效应的干扰受到抑制以及其随机噪声受到抑制的时间-密度曲线。顺便提一句，将这样获取时间-密度曲线的解决方案特别称为“相干回归方法”。

接着，将具体说明在进行医学处理的X射线CT中，在例如动态CT射线照相过程获得的动态图像中临床应用时间-密度曲线的例子。在该应用例子中，在将造影剂迅速注入血管时，进行动态CT射线照相等，测量存在于人体组织内的动脉图像的密度变化作为时间-密度曲线，从而对该组织内的血流流动作出诊断。

在该应用例子中，人体组织内的动脉在许多情况下通常非常细微，出现在基于CT的层析X射线照片上的动脉图像受到部分体效应的影响。此外，随机噪声当然地存在于该图像中。因此，利用现有技术，难以获得精度令人满意的、动脉的时间-密度曲线。在冒险测量时间-密度曲线时，只能获得在某种程度上近似于该动脉的真实时间-密度曲线<t_k，D_k>的线性变换<t_k，A(D_k-D₁)>的测量值<t_k，(v_k(x)-v₁(x))>(其中D_k表示在时间t_k对应于动脉图像的一组像素的像素值(是标量值)，而且k＝1，2，...，K保持)。该测量值含有随机噪声。此外，系数A仍为未知，因为受到部分体效应的干扰。

因此，在应用根据本发明的上述解决方案时，可以获得测量值<t_k，(v’_k(x)-v’₁(x))>(k＝1，2，...，K)，它是线性变换<t_k，A(D_k-D₁)>的满意近似。另一方面，可以在同一个层析X射线照片中观察到的静脉显著包括厚静脉。因此，对于这种静脉，利用现有技术方法可以获得时间-密度曲线的令人满意的良好近似值<t_k，(J_k-J₁)>(k＝1，2，...，K)。其中J_K表示在时间t_k对应于静脉图像的一组像素的像素值。

同时，关于血液循环的时间-密度曲线，已知命题S的固有属性：“如果在时间t₁，血液中造影剂的密度为0(零)，则曲线下面的面积(AUC)与任何血管d的时间-密度曲线<t_k，(d_k-d₁)>一致”。在此提到的“曲线下的面积”表示时间-密度曲线<t_k，(d_k-d₁)>对时间t的积分。

因此，利用例如下面的等式，可以近似计算特定血管d的时间-密度曲线<t_k，(d_k-d₁)>的曲线下的面积(area-under-curve)AUC(d)：

$\mathrm{AUC}\left(d\right)\&cong;\underset{k=1}{\overset{K-1}{\mathrm{\&Sigma;}}}\frac{\{(d_{k+1}-d_1)+(d_k-d_1)\}(t_{k+1}-t_k)}{2}---\left(30\right)$

因此，利用等式(30)(J可以代入d)可以计算利用现有技术方法对静脉获得的时间-密度曲线{<t_k，(J_k-J₁)>}的曲线下的面积AUC(J)。此外，关于动脉，如果时间-密度曲线{<t_k，(D_k-D₁)>}已知，则同样可以利用等式(30)计算曲线下的面积AUC(D)。此外，根据上述命题S，下式应当成立：

$\mathrm{AUC}\left(d\right)\&cong;\mathrm{AUC}\left(J\right)---\left(31\right)$

然而，实际上，不知道该时间-密度曲线{<t_k，(d_k-d₁)>}，因此不能计算曲线下的面积AUC(D)。

相反，利用根据本发明的解决方案获得的时间-密度曲线<t_k，(v’_k(x)-v’₁(x))>近似于<t_k，A(D_k-D₁)>，而后者包括未知系数A。因此，可以利用等式(30)、根据曲线{<t_k，(v’_k(x)-v’₁(x))>}具体计算的曲线下的面积AUC(v’)一定是面积AUC(D)的A倍。即：

$\mathrm{AUC}(v\&prime;)\&cong;\mathrm{AAUC}\left(D\right)---\left(32\right)$

根据等式(31)和(32)，以下关系成立：

$A\&cong;\mathrm{AUC}(V\text{'})/\mathrm{AUC}\left(J\right)---\left(33\right)$

由于可以利用等式(30)具体计算等式(33)的右侧，所以可以具体确定未知系数A的值。因此，在利用系数A的值建立时间-密度曲线<t_k，(v’_k(x)-v’₁(x))/A>时，只能近似计算动脉的时间-密度曲线<t_k，(D_k-D₁)>。将以这种方式利用曲线下的面积，确定其未知的比例系数A的值来建立时间-密度曲线的方法称为“AUC方法”。

因此，如上所述，在临床使用利用动态CT射线照相等获得的动态图像中的时间-密度曲线时，还可以进一步将AUC方法与相干回归方法组合在一起使用，这样，即使对于利用现有技术方法难以或不可能测量的微小动脉的时间-密度曲线，仍可以获得不包括未知比例系数A、其排除了部分体效应和随机噪声的干扰的测量值。

顺便提一句，AUC方法当然可以单独应用于利用现有技术方法测量的动脉的时间-密度曲线<t_k，(v’_k(x)-v’₁(x))>，而且可以建立确定了其未知比例系数A的值的时间-密度曲线(但是不能排除随机噪声和部分体效应的干扰)。

对于特定像素x，进一步将AUC方法与相干回归方法组合在一起获得的时间-密度曲线如例如图11中的符号(A)所示。顺便提一句，该图中的(B)还示出利用现有技术方法建立的、仅对其应用AUC方法获得的时间-密度曲线。显然，在曲线(B)中，可以清楚地看到随机噪声的干扰，而在曲线(A)中，在根本不破坏时间分辨率的情况下，该噪声被充分抑制。

顺便提一句，尽管在上述说明中，将标量值d_k用作时间-密度曲线的密度值，但是该实施例并不局限于这种情况。显然，可以容易地扩展该实施例，而且该实施例可以容易地应用于例如构成动态图像的各图像是彩色图像，更一般地是多种图像的组合、并将时间-密度曲线的密度值表示为矢量值的情况。

(X各实施例的辅助项目)

现在，将说明上述实施例的辅助项目。

(X-1通用辅助项目)

首先，在上述实施例中，假定放弃虚假设H情况下的“临界概率”p(x，y)是用于量化在本发明中被称为“拟和优度”的手段，但是本发明并不局限于该方面。正如在开头说明的那样，“拟和优度”可以是用于表示像素x与y在各种意义上是否相似的数值“指标”(index)。

此外，在这方面，为了在放弃几个虚假设情况下计算临界概率(请参考等式(6)和(9))，采用x平方测试方法，但是本发明并不局限于此。换句话说，本发明不仅仅局限于采用“统计测试方法”或“x平方测试方法”的方面，该方面是本发明的一种具体方面。

此外，在上述实施例中，作为具体形式的权重函数w(p(x，y))，所描述的权重函数总共有两种，w1(等式(6)或(9)代入其内的等式(10))和w3(等式(14))，但是本发明并不局限于此。如上所述，作为临界概率p(x，y)∈[0，1]的函数的权重函数w满足的属性仅是权重函数w(t)是域t∈[0，1]定义的非负单调递增函数这个事实。此外，更一般地说，权重函数w可以是“关于该拟和优度的非负单调递增函数”。

(X-2关于拟和优度的最通用方面)

在本发明中所称的“拟和优度”并不局限于以上说明的其中在放弃虚假设的情况下根据临界概率p(x，y)给出其数值的解决方案。因此，“权重”也不局限于其中通过使权重函数w作用于临界概率p(x，y)计算它的解决方案。在下面的说明中，将结合实际例子说明这种拟和优度等的更通用方面。

首先，一般地说，本发明通常采用其中对于根据图像处理的用途而正确设置的特定命题，采用符合该用途的正确解决方案，利用数值给出在要处理的图像内该命题成立的确实性，然后利用该数值进行图像处理的建立过程。作为例子，为了利用数值给出像素x与y之间的拟和优度，在本发明中可以很好地建立更适度的比例。作为其实际例子，考虑到其中采用模糊逻辑等中的成员函数，利用数值给出拟和优度的建立过程。

(X-2-1将本发明应用于图形识别的例子)

在此，作为切实可行的例子，关于拟和优度等的更通用方面，将对一个在将本发明应用于用于识别并提取图像上的特征图形的所谓“图形识别”的情况下的例子进行说明。

在这种情况下，根据要识别的图形的特征，事先确定用于表示命题“像素x是构成该图形的像素”的确实性的函数m(x)。在此，m(x)∈[0，1]成立。更详细地说，函数m(x)是用于根据作为像素x的像素值的矢量值v(x)、位于像素x附近包括在像素集N(x)内的每个像素y的像素值v(y)、表示该矢量的位置的矢量x自身等计算表示该命题的确实性的数值的函数。因此，需要根据要提取的图形与不要提取的图形之间的差别，正确设置函数m(x)。此外，利用下式定义权重函数：

此外，在建立对于其函数m(x)大于特定“阈值”的像素集X时，获得对应于该图形的像素x的区域。

现在，在给出要处理的图像P时，对该图像的所有像素x，计算其中将具体阈值T和函数m(x)代入上面的权重函数中获得的w(T，m(x))。(因此，确定其成员函数是w(T，m(x))(即，对于其w(T，m(x))＝1成立的像素x的像素集)的像素集X，而且该像素集X仅是该图形所占据的图像P的区域。)此外，还产生新图像P’，如下所述。在函数w(T，m(x))＝1成立时，图像P’的像素x具有像素值v(x)，而在所有其他情况下，其像素值为0(零矢量)。这样产生的图像P’就是利用0擦除图像P上不对应于要识别的图形的各像素的像素值获得的图像。

当然，可以很好地顺序执行计算函数m(x)、计算权重w(T，m(x))以及建立图像P’的处理步骤，也可以利用集成了这些处理步骤的软件执行这些处理步骤。

在该例子中，通过检查该图形的特征、作为像素x的像素值的矢量值v(x)以及矢量x本身的值，该函数m(x)利用数值给出图形与像素x之间的“拟和优度”。此外，利用等式(34)定义的权重函数计算权重w(T，m(x))，而利用该权重和图像P建立作为处理结果的新图像P’。

作为该处理过程又一个切实可行的例子，将说明为了提取航摄照片、金属表面的显微照片等中的细线而进行图像处理的建立过程。

首先，将说明用于建立利用数值表示命题“像素x是构成该细线的像素”的确实性的指标(index)的方法。

事先确定位于像素x附近的各像素的像素集N(x)。通常，可以将该像素集N(x)定义为以像素x为中心的矩形区域。此外，还要定义其中排列了像素集N(x)内的各像素的像素值(标量值或矢量值)的多维矢量值v(x)。K表示矢量的维数。通常，最好将数字K设置在约25至169之间。

随后，如下所述，定义作为K维矢量函数、在后面的计算中辅助使用的归一化函数nrm。在v≠0的条件下，对任何一个K维矢量v给出归一化函数nrm：

${\mathrm{nrm}}_k\left(v\right)=\frac{v_k-\frac{1}{K}\underset{i=1}{\overset{K}{\mathrm{\&Sigma;}}}{v}_i}{\sqrt{\underset{j=1}{\overset{K}{\mathrm{\&Sigma;}}}{(v_j-\frac{1}{K}\underset{i=1}{\overset{K}{\mathrm{\&Sigma;}}}{v}_i)}^2}}---\left(35\right)$

其中k＝1，2，...，K成立。因为该归一化函数，可以确保，对于任何一个K维矢量v，下式成立(只要矢量v不是零矢量)：

$\underset{i=1}{\overset{K}{\mathrm{\&Sigma;}}}\mathrm{nr}{m}_k\left(v\right)=0$

$\underset{i=1}{\overset{K}{\mathrm{\&Sigma;}}}\{\mathrm{nr}{m}_k\left(v\right){\}}^2=1$

随后，随意确定一个像素x，以适当编号J采集每个分别具有0(零)背景而且仅包括一条通过像素x的细直线的各图像的典型例子，并对其分配图像号1，2，...，J。计算第J个图像中的v(x)的值，从而根据下式具体建立矢量r(j)：

      r(j)＝nmrj(v(x))

对于j＝1，2，...，J执行该计算过程。(以下将这样建立的J矢量r(j)称为“图形矢量”。顺便提一句，该图形矢量对应于在本发明中所称的“分别建立的像素值”。)

此外，利用下面的等式确定函数m(x)：

      m(x)＝max{r(j)·nrm(v(x))|j＝1，2，...，J}

其中“·”表示K维矢量之间的内积，而“max”表示用于获得像素集内各元素的最大值的函数。在像素x与任何一个图形矢量r(j)(j＝1，2，...，J)匹配或相似的情况下，这样建立的函数m(x)取大值，而在其他情况下，它取小值。因此，该函数m(x)是利用数值表示命题“像素x是构成细线的像素”的确实性的指标，即，表示“拟和优度”的函数。此外。利用下式确定权重函数：

(再一次提到的34)

其中T表示用作阈值的常数，并对其设置适当值。

现在，在给出要处理的图像P时，对该图像的所有像素x，计算其中将具体阈值T和函数m(x)代入上面的权重函数获得的w(T，m(x))。(因此，确定其成员函数是w(T，m(x))(即，对于其w(T，m(x))＝1成立的像素x的像素集)的像素集X，而且该像素集X仅是该图形所占据的图像P的区域。)

此外，还产生新图像P’，如下所述。在函数w(T，m(x))＝1成立时，图像P’的像素x具有像素值v(x)，而在所有其他情况下，其像素值为0(零矢量)。这样产生的图像P’就是利用0擦除图像P上不对应于要识别的图形的各像素的像素值获得的图像。

此外，可以很好地将w(T，t)构建为模糊识别函数，例如：

(其中C表示正常数)，可以很好地建立图像P’的像素x以使w(T，m(x)V(x))作为其像素值。这样就降低了因为不与“细线”显著对应而利用模糊逻辑调节的像素的对比度，从而产生了建立图像P’的效果，因此可以使“细线”部分特别醒目。

(X-3估计方差σ²的估计方法的辅助项目)

接着，将辅助说明用于估计在上面的等式(6)或(14)中使用的方差σ²的方法。

在上述“将本发明应用于X射线CT射线照相的情况”中，根据对于所有K个图像上的所有像素x，利用上面的等式(6)和(7)获得的方差的期望值E[σ²]是正确的背景，该期望值E[σ²]被不加区别地用于所有像素x。然而，还存在不满足这种假设条件的情况，即，对于各像素x(在此，它们是x₁，...，和x_F)，通常假定特定方差σ²(x₁)，...，和σ²(x_F)的情况(即使分布形状相同(例如，高斯分布))。

在这种情况下，对于各像素x，可以根据上面的等式(6)和(7)，分别计算各方差的期望值E[σ²(x₁)]，...和E[σ²(x_F)]。此后，根据各种情况(＝根据注意的像素x₁，...，和x_F)，可以将这样分别计算的期望值E[σ²(x₁)]，...和E[σ²(x_F)]代入等式(3)，以获得v’(x)＝(v’₁(x)，v’₂(x)，...，v’_K(x))。

此外，在假定在所有像素的像素值内包含的噪声具有共同方差σ²的情况下，可以采用如下所述的解决方案计算方差σ²。首先，对于所有像素x(＝{x₁,...，x_F})的像素集Ω，获得利用上面的等式(7)和(8)计算的各方差的期望值E[σ²(x)]的平均值(σ²)⁺，并将它用作等式(6)中的σ²的第一近似估计值。即：

${\left({\mathrm{\&sigma;}}^2\right)}^{+}=\frac{1}{\left|\mathrm{\&Omega;}\right|}\underset{x\&Element;\mathrm{\&Omega;}}{\mathrm{\&Sigma;}}E\left[{\mathrm{\&sigma;}}^2\left(x\right)\right]---\left(35\right)$

在上面的等式中，|Ω|表示像素集Ω内的元素数量(其中|Ω|＝F成立)。随后，设M表示对于其期望值E[σ²(x)]在例如上述平均值(σ²)⁺的几倍范围内的像素x的像素集。在设(σ²)⁺⁺表示仅对使用像素集M的元素的像素x取的新平均值时，它变成：

${\left({\mathrm{\&sigma;}}^2\right)}^{++}=\frac{1}{\left|M\right|}\underset{x\&Element;M}{\mathrm{\&Sigma;}}E\left[{\mathrm{\&sigma;}}^2\left(x\right)\right]---\left(36\right)$

可以将该平均值用作σ的更适当估计值。

在例如以下具体说明的情况下，即在如图12(a)所示存在构成动态图像的2个静态图像G1和G2，而且获取图像G1与G2之间的差别如图12(b)所示，图像Z1和Z2在某种程度上互相偏离的情况(换句话说，移动图像Z1作为图Z2)，该处理过程是有效的。

上面的等式(35)的方差(σ²)⁺对应于根据构成图像G3的所有像素x(＝Ω)计算的方差，如图12(b)所示。然而，在这种情况下，使因为图像G1和G2不互相重叠的部分ζ₁和ζ₂产生的两个峰值ζ₁’和ζ₂’出现在图G3的概率分布图内，如图12(c)所示。因此，对方差(6²)⁺的估计过高。因此，对于所有静态图像，不能原封不动地正确使用方差(σ²)⁺。另一方面，对于图12(b)所示的部分ζ_M，可以认为图12(c)所示的中心峰值ζ_M’与其对应。因此，在通过去除两个峰值ζ₁’和ζ₂’而保留该峰值ζ_M’的情况下，更优先计算概率分布图的方差。该方差确实与上面的等式(32)中的(σ²)⁺⁺对应。

换句话说，仅对“对于其期望值E[σ²(x)]在平均值(σ²)⁺”几倍范围内的像素x”的像素集M，取平均值的上述处理过程仅是通过仅提取图12(c)所示的峰值ζ_M’获取平均值的处理过程。因此，方差(σ²)⁺⁺的计算过程等效于在去除图12(c)所示峰值ζ₁’和ζ₂’情况下计算概率分布图的方差的计算过程。

因此，在这种情况下，将平均值(σ²)⁺⁺代入上面的等式(5)等，  因此可以获得各像素x₁，...，和x_F的变换像素值v’(x)＝(v’₁(x)，v’₂(x)，...，v’_K(x))的更适当值。

顺便提一句，关于噪声估计过程，通常，其中将均匀亮度的对象置于位于例如图像中心的位置，然后测量位于相应部分的像素的像素集的值的分布，从而估计噪声分布的方法也非常有效。

在本发明中，用于估计噪声分布或噪声方差的方法基本上可以基于任何一种解决方案，而且并不局限于以上说明的几个解决方案。实际上，在可用的一实施例中，应当正确采用其中可以最优先利用实验信息、成像处理原理、测量值或实际测量值等计算分布或方差的解决方案。顺便提一句，这时，最好，根据每个实施例要求的实际精度，简化设备的结构并提高尽可能容易地建立估计方法的处理速度。

(X本发明的另一个实施例)

最后，将辅助说明在以上说明中未说明的本发明的另一个实施例。

(X-4-1数据压缩预处理的实施例)

通常，在建立其中对图像Q附加了随机噪声的图像Q’时，图像Q’具有比图像Q大得多的信息内容(例外情况是，原始图像Q含有非常大量的噪声)。因此，如果要传送该图像Q’或者要将它记录到记录介质上，在通信或记录过程中，用于描述没有实质意义的噪声的信息内容占用位数是一种浪费。反过来说，在根据图像Q’建立抑制了噪声的图像Q”，并为了进行通信或记录而对其进行处理时，可以利用小信息内容表示具有同样有效内容(而且鲜艳的)的图像。

通常，在在通信过程中或记录过程中处理图像的情况下，通常要对图像进行所谓“压缩处理”。该处理过程意在最大限度地减少要传送或记录的信息内容，从而实现高速通信或大量记录(对一个记录介质)。在实际使用中，正是在“动态图像”情况下，尤其要减少信息内容。这些情况在不可逆图像压缩技术中众所周知。

根据这些项目，通过压缩抑制了噪声的图像Q”，而非通过压缩含有噪声的图像Q’，可以减少进行通信或记录所需的位数(信息内容)。

因此，在通过应用根据本发明的相干滤波器(具体地说，利用硬件或软件建立“预处理设备”，这同样适用于以下内容)，对要通信或要记录的图像Q’进行预处理以抑制噪声时，获得抑制了噪声的图像Q”，而不破坏时间分辨率或空间分辨率，而且可以以高压缩比压缩该图像Q”，因此可以显著节省进行通信或记录所需的时间周期和存储容量。

(X-4-2二维或三维局部提取处理的预处理过程的实施例)

通常执行从航摄照片等的二维图像或从基于MRI设备等的三维分布图像中提取指定区域的处理过程。作为例子，对于前者的航摄照片，仅提取公路区域，而对于后者基于MRI设备等获得的人头部的三维分布图像，仅提取血管区域，然后，对提取的区域进行再现(rendering)以将它们显示为计算机图形。特别是，作为例子，可以相对旋转或从各方向观察后者中的血管，其射线照片被认为是诊断脑血管疾病等的重要技术。

在进行这种区域提取的情况下，最简单、最基本的处理就是所谓“阈值处理”。阈值处理是其中根据其标量值是否大于某个特定阈值，对构成二维图像或三维分布图像等的各像素x进行分类的处理。此外，在对仅包括属于规定类的像素x的像素集进行再现时，可以显示区域提取后的图像，例如“仅包括血管的三维图像”(顺便提一句，作为现有技术，还使用非“阈值处理”的其他解决方案)。

同时，在该区域提取处理过程中，原始图像(上述二维图像或三维分布图像等)中最好几乎不含有噪声。这是因为，在图像中存在噪声时，出现这样的情况，即，位于不应该是例如血管的地方的像素具有就好象它是血管一样的像素值，而位于应该是血管的地方的像素具有就好象它不是血管一样的像素值。因此，在这种情况下，需要利用复杂的后处理过程执行去除被错误分类的像素的步骤。此外，在实际情况中，为了消除原始图像中的噪声，通常执行现有技术小节中描述的“平滑处理”。然而，在这种情况下，存在的问题是破坏了空间分辨率，因此要提取的区域上的细微结构(包括例如微细血管)被丢失。

因此，在应用根据本发明的相干滤波器，对原始图像，即二维图像或三维分布图像进行预处理以抑制噪声时，可以阈值图像的噪声，而又不破坏空间分辨率，因此可以实现精确区域提取处理。

顺便提一句，在上面的小节“根据一个图像建立像素值v(x)的情况(参考图8(c))”中已经对与这种提取处理类似的实施例进行了说明。

Claims (32)

1、一种图像处理设备，其特征在于包括：

判定装置，用于通过进行测试，判定构成图像的第一像素与第二像素之间的相似性；以及

平均装置，用于根据判定装置的判定结果，对第一像素和第二像素进行加权平均。

2、一种图像处理设备，其特征在于包括：

数值装置，用于以数值给出利用构成图像的第一像素与第二像素之间的相似性确定的权重；以及

平均装置，用于利用以数值给出的权重，对第一像素和第二像素的像素值进行加权平均。

3、根据权利要求1或2所述的图像处理设备，其特征在于，平均装置获得第一像素的新像素值作为加权平均值。

4、根据权利要求1所述的图像处理设备，其特征在于，平均装置包括：确定装置，用于根据判定结果，确定加权系数；以及乘法装置，用于将第一像素和第二像素的值乘以加权系数。

5、根据权利要求1、2、3和4之任一所述的图像处理设备，其特征在于：

在第一像素与第二像素之间的相似性高时，平均装置利用乘法装置使第二像素值乘以大加权系数，而在相似性低时，使第二像素值乘以小加权系数；以及

对被加权的第一像素值和第二像素值进行加权平均。

6、根据权利要求1所述的图像处理设备，其特征在于，根据以利用成像设备在不同时间对同一个受治疗者进行照相的方式获得的多个图像，判定装置判定第一像素与第二像素之间的相似性。

7、根据权利要求1或6所述的图像处理设备，其特征在于，根据以利用成像设备在不同照相条件下对同一个受治疗者进行照相的方式获得的多个图像，判定装置判定第一像素与第二像素之间的相似性。

8、根据权利要求1、6和7之任一所述的图像处理设备，其特征在于，根据利用成像设备在不同处理条件下对同一个受治疗者进行照相获得的多个图像，判定装置判定第一像素与第二像素之间的相似性。

9、根据权利要求1所述的图像处理设备，其特征在于，利用对通过排列利用成像设备对受治疗者进行照相获得的单个图像中的第一像素的标量值建立的矢量值与通过排列第二像素的标量值建立的矢量值进行比较获得的比较结果，判定装置判定第一像素与第二像素之间的相似性。

10、根据权利要求1、6、7和8之任一所述的图像处理设备，其特征在于，利用对通过将多个图像的相应像素值排列到第一像素坐标中建立的第一矢量与通过将多个图像的相应像素值排列到第二像素坐标中建立的第二矢量进行比较获得的比较结果，判定装置判定多个图像中的一个预定图像的第一像素与第二像素之间的相似性。

11、根据权利要求1或2所述的图像处理设备，其特征在于，该图像处理设备进一步包括图像平均装置，用于平均利用平均装置对图像进行平均处理获得的第一图像和对该图像进行不同于该平均处理的处理获得的第二图像。

12、根据权利要求1或2所述的图像处理设备，其特征在于，该图像处理设备进一步包括图像平均装置，用于平均平均装置进行平均处理的第一图像和对第一图像进行平均处理之前的第二图像。

13、一种图像处理设备，其特征在于包括：

确定装置，用于通过进行统计测试，确定并以数值给出构成图像的第一像素与第二像素之间的相似性；以及

平均装置，用于根据该数值，对第一像素和第二像素进行加权平均。

14、一种图像处理设备，其特征在于包括：

数值装置，用于通过进行统计测试，以数值给出构成图像的第一像素与第二像素之间的相似性；以及

平均装置，用于在数值装置中以数值表示的相似性高时，对第一像素和第二像素的值进行平均，而在确定的相似性低时，不平均第一像素值和第二像素值。

15、一种图像处理方法，其特征在于：

对构成图像的第一像素和第二像素，计算分别通过排列标量值建立的各矢量值；以及

在量化分别根据第一像素值建立的第一像素的值和第二像素的值之间的相似性，然后利用第二像素值，对第一像素建立新像素值的过程中，通过在相似性高时，使第二像素的贡献大，而在相似性低时，使第二像素的贡献小，建立新像素值。

16、一种图像处理方法，其特征在于：

对构成图像的第一像素和第二像素，建立分别是通过排列标量值建立的矢量值的第一像素值和第二像素值；

量化第一像素值与第二像素值之间的相似性；以及

利用第二像素值建立第一像素的新像素值，

通过在相似性高时，使第二像素的贡献大，而在相似性低时，使第二像素的贡献小，建立该新像素值。

17、一种图像处理方法，其特征在于：

对构成图像的第一像素和第二像素，建立分别是通过排列标量值建立的矢量值的第一像素值和第二像素值；

量化第一像素值与第二像素值之间的相似性；以及

在通过对各第二像素值计算的相似性应用作为相似性的函数的权重函数，从而确定各其他像素值的权重，然后根据权重顺序计算其他像素值的加权平均值，来建立第一像素的新像素值过程中，通过在相似性高时，使权重大来使其他像素在加权平均值中的贡献大，而在相似性低时，使权重小来使其他像素在加权平均值中的贡献小，建立新像素值。

18、根据权利要求17所述的图像处理方法，其特征在于，权重函数是关于相似性的非负单调递增函数。

19、根据权利要求17所述的图像处理方法，其特征在于，

在利用x和y分别表示第一像素和第二像素时，分别利用v(x)＝(v₁(x)，v₂(x)，...，v_k(x))和v(y)＝(v₁(y)，v₂(y)，...，v_k(y))表示第一像素值和第二像素值，而利用v’(x)＝(v’₁(x)，v’₂(x)，...，v’_k(x))表示要建立的新像素值，而且在ρ(x，y)表示相似性，w表示权重函数以及w(ρ(x，y))表示通过对相似性应用权重函数w获得的权重时，将进行加权平均的处理过程表示为：

${v\&prime;}_k\left(x\right)=\frac{\underset{y\&Element;N\left(x\right)}{\mathrm{\&Sigma;}}{v}_k\left(y\right)w\left(\mathrm{\&rho;}(x,y)\right)}{\underset{y\&Element;N\left(x\right)}{\mathrm{\&Sigma;}}w\left(\mathrm{\&rho;}(x,y)\right)}$

(其中N(x)表示包括其他像素的范围)。

20、根据权利要求15至17之任一所述的图像处理方法，其特征在于，通过建立新像素值，抑制图像的噪声。

21、根据权利要求15至17之任一所述的图像处理方法，其特征在于，通过建立新像素值，在图像内进行图形识别。

22、根据权利要求15至17之任一所述的图像处理方法，其特征在于，根据通过对分别构成第一像素值和第二像素值的标量值应用统计测试方法获得的临界概率，量化相似性。

23、根据权利要求16或17所述的图像处理方法，其特征在于，在围绕第一像素的预定范围内选择第二像素。

24、根据权利要求23所述的图像处理方法，其特征在于，预定范围是以第一像素为中心的矩形区域。

25、根据权利要求16或17所述的图像处理方法，其特征在于，所述图像包括多个图像，而且其特征还在于，像素值是通过排列多个图像的同一个点上的各像素的标量值获得的矢量值。

26、根据权利要求25所述的图像处理方法，其特征在于，多个图像是构成动态图像的多个静态图像。

27、根据权利要求26所述的图像处理方法，其特征在于，动态图像是医学图像诊断设备获取的动态CT图像。

28、根据权利要求16或17所述的图像处理方法，其特征在于，该图像是同一个受治疗者的多种图像，而且其特征还在于，像素值是通过排列多种图像的同一个点上的各像素的标量值获得的矢量值。

29、根据权利要求28所述的图像处理方法，其特征在于，该图像是通过在SPECT设备或PET设备内进行多窗***线照相获取的图像，而且其特征还在于，多种图像分别是根据不同放射性同位素分别发出的γ射线而不同的图像。

30、根据权利要求28所述的图像处理方法，其特征在于，图像是彩色图像，而且其特征还在于，多种图像分别是通过将光分解为3原色获取的不同的图像。

31、根据权利要求15至17之任一所述的图像处理方法，其特征在于，所述图像是单个图像，而且其特征还在于，根据所述单个图像建立像素值。

32、根据权利要求31所述的图像处理方法，其特征在于，像素值是通过所述单个图像上的特定像素的标量值以及位于该特定像素附近的各像素的标量值获得的矢量值。

Images