CN101739670B - 非局部均值空域时变图像滤波方法 - Google Patents

Abstract

本发明公开了一种非局部均值空域时变图像滤波方法，主要解决现有非局部均值空域时变滤波器的计算量大和应用范围受限的问题。其滤波过程是：初始化原始图像所有点处的权值归一化系数和非归一化滤波后的值为0；对于搜索域内每个坐标偏置，分别将原始图像所有的像素点进行统一的预处理，然后快速计算原始图像所有像素点在该搜索域坐标处的加权核权值；根据该加权核权值，更新权值归一化系数和非归一化滤波后的值；根据权值归一化系数和非归一化滤波后的值，计算得到滤波后的图像。本发明极大降低了现有非局部均值空域时变滤波方法的计算复杂度，可应用于图像恢复、图像去噪、图像超分辨率重建领域。

Description

非局部均值空域时变图像滤波方法

技术领域

本发明属于图像处理领域，涉及图像滤波方法，具体地说是一种非局部均值空域时变图像滤波方法，用于对数字图像和视频图像的处理。

技术背景

随着信息与多媒体技术的发展，数字图像处理技术的应用越来越广泛。空域滤波技术是数字图像处理的基本方法之一，在数字图像领域有着广泛的应用，比如图像恢复、图像去噪、图像超分辨率重建等。空域滤波的基本思想是利用临近像素点之间的关系来得到相应的处理目的。临近的像素点组成的形状被称为搜索域。对于图像中的任意像素点，该点滤波后的值等于以该点为中心的搜索域的所有点进行归一化加权并求和后的值。如果用I，

分别表示原始图像和空域滤波后的图像，则有空域滤波器的基本表达式为：

$\hat{I}(x,y)=\frac{1}{\mathrm{\α}(x,y)}\underset{(u,v)\∈V_R}{\mathrm{\Σ}}w(x,y,u,v)I(x+u,y+v),$

其中，V_R是临近像素点组成的搜索域，u、v为V_R内行、列坐标，x、y分别为像素点的行、列坐标，w(x，y，u，v)是加权核，

是归一化系数。当加权核w(x，y，u，v)为一个与x，y无关的固定函数时，称上述过程为空域时不变滤波；当加权核w(x，y，u，v)为一个与x，y相关的变化函数时，称上述过程为空域时变滤波。

非局部均值滤波器是由Buades等人在文章“A review of imagedenoising algorithms，with a new one”中提出的空域时变滤波器。非局部均值滤波器搜索域是半径为R的矩形窗，即

非局部均值滤波器的加权核其中，S(x，y，u，v)称为图像偏置相似函数， $S(x,y,u,v)={\left|\right|B_r(x,y)-B_r(x+u,y+v)\left|\right|}_p^p,$ h是加权核宽度，

表示l_p范数值的p次方，p≥1，exp(·)为指数函数，B_r(x，y)为表示以(x，y)为中心，以r为半径的矩形块的匹配窗。

参见图1，现有非局部均值空域时变图像滤波流程为：

[1]获取原始图像信息，设置搜索域形状、匹配窗形状、加权核宽度h和相似性比较范数参数p；

[2]对原始图像中每个坐标点，分别执行步骤[3]～[6]进行滤波；

[3]对原始图像中以当前滤波坐标为中心的搜索域内的每个像素点，分别计算对应的加权核权值：

$w(x,y,u,v)=\exp (-\frac{S(x,y,u,v)}{{2h}^p}),$

其中，(u，v)为搜索域坐标， $S(x,y,u,v)={\left|\right|B_r(x,y)-B_r(x+u,y+v)\left|\right|}_p^p;$

[4]计算原始图像中当前滤波坐标为中心的搜素域范围内所有加权核权值的和，得到归一化系数；

[5]计算原始图像中以当前滤波坐标为中心的搜素域内，每个像素点对应的权

值乘以对应点信息值的和，作为当前滤波坐标处的非归一化滤波后的值；

[6]将原始图像中在当前滤波坐标处的非归一化滤波后的值除以归一化系数，获得当前滤波坐标处的滤波后的值；

[7]所有坐标处滤波后的值组成滤波后的图像。

上述的非局部均值滤波方法在计算相邻的像素点和搜索域偏置组合的图像偏置函数的时候，会出现图像匹配窗覆盖的现象，从而会导致覆盖的匹配窗内像素点和加搜索域偏置像素组合之间，多次执行比较计算的计算冗余问题。而上述非局部均值滤波器的在逐点滤波过程中，每个像素点的每个搜索域偏置点的加权核权值都需要重新计算，这样会造成大量的重复计算。若图像大小为M×N，则非局部均值滤波需要的计算复杂度为O[(2r+1)²(2R+1)²MN]。这种极大的计算复杂度，造成计算时间过长，极大的限制了非局部均值滤波器的应用，比如在实时嵌入式系统中就不能应用这种计算复杂度的非局部均值滤波方法。

发明内容

本发明的目的在于克服上述已有技术的问题，提出一种非局部均值空域时变图像滤波方法，以减小非局部均值空域时变滤波器的计算复杂度，扩大非局部均值滤波器的应用范围。

实现本发明的技术思路为：对于搜索域内每个坐标，分别将原始图像所有的像素点进行统一的预处理，然后计算原始图像所有像素点在该搜索域坐标处的加权核权值；根据搜索域所有坐标处对应的原始图像每个点的加权核权值，计算出滤波后的图像信息。其具体实现方案如下：

技术方案1，一种非局部均值空域时变图像滤波方法，其包括如下步骤：

(1)获取原始图像信息，设置搜索域形状、匹配窗大小、加权核宽度h和相似性比较范数参数p；

(2)初始化非归一化滤波后的图像

和权值归一化系数α(x，y)均为0，其中(x，y)为图像坐标；

(3)以搜索域内的点作为偏置，对所有点执行步骤(4)～(9)的操作；

(4)以当前偏置(u，v)为偏置量，计算原始图像与偏置后图像之间所有对应点的数据差的绝对值的p次方之值，得到p范数差值图像Δ^(u，v)(x，y)；

(5)计算当前偏置(u，v)下差值图像Δ^(u，v)(x，y)的积分图像S^(u，v)(x，y)，计算公式为：

$S^{(u,v)}(x,y)=\underset{y^{\′}\≤y}{\mathrm{\Σ}}\underset{x^{\′}\≤x}{\mathrm{\Σ}}{\mathrm{\Δ}}^{(u,v)}(x^{\′},y^{\′});$

(6)通过下式计算在当前偏置坐标(u，v)下的原始图像每个坐标处的图像偏置相似函数值：

D(x，y，u，v)＝S^(u，v)(x_l，y_t)+S^(u，v)(x_r，y_b)-S^(u，v)(x_l，y_b)-S^(u，v)(x_r，y_t)，

其中，(x，y)为图像坐标，(x_l，y_t)，(x_r，y_b)分别表示以(x，y)中心的匹配窗的左上角坐标和右下角坐标；

(7)通过下式计算在当前偏置(u，v)下，原始图像每个坐标处的加权核权值：

$w(x,y,u,v)=\exp (-\frac{D(x,y,u,v)}{{2h}^p}),$

其中，w(x，y，u，v)为加权核权值，h为加权核宽度，exp(·)为指数函数；

(8)在当前偏置坐标(u，v)下，通过下式更新原始图像每个坐标处的权值归一化系数：

α(x，y)＝α(x，y)+w(x，y，u，v)，

其中，α(x，y)为归一化系数；

(9)在当前偏置坐标下，通过下式更新原始图像每个坐标处的非归一化滤波后的值：

$\stackrel{\‾}{I}(x,y)=\stackrel{\‾}{I}(x,y)+w(x,y,u,v)\×I(x+u,y+v),$

其中，

为非归一化滤波后的值，I为原始图像信息；

(10)通过下式对原始图像每个坐标进行归一化处理，得到滤波后图像

$\hat{I}(x,y)=\frac{\stackrel{\‾}{I}(x,y)}{\mathrm{\α}(x,y)}.$

技术方案2，一种非局部均值空域时变图像滤波方法，包括如下步骤：

1)获取原始图像信息，设置搜索域形状、匹配窗大小、加权核宽度h和相似性比较范数参数p；

2)以搜索域内的点作为偏置，对所有点执行步骤3)～6)的操作；

3)以当前偏置(u，v)为偏置量，计算原始图像与偏置后图像之间所有对应点的数据差的绝对值的p次方之值，得到p范数差值图像Δ^(u，v)(x，y)，其中(x，y)为图像坐标；

4)计算当前偏置(u，v)下差值图像Δ^(u，v)(x，y)的积分图像S^(u，v)(x，y)，计算公式为：

$S^{(u,v)}(x,y)=\underset{y^{\′}\≤y}{\mathrm{\Σ}}\underset{x^{\′}\≤x}{\mathrm{\Σ}}{\mathrm{\Δ}}^{(u,v)}(x^{\′},y^{\′});$

5)通过下式计算在当前偏置坐标(u，v)下的原始图像每个坐标处的图像偏置相似函数值：

D(x，y，u，v)＝S^(u，v)(x_l，y_t)+S^(u，v)(x_r，y_b)-S^(u，v)(x_l，y_b)-S^(u，v)(x_r，y_t)，

其中，(x，y)为图像坐标，(x_l，y_t)，(x_r，y_b)分别表示以(x，y)中心的匹配窗的左上角坐标和右下角坐标；

6)通过下式计算在当前偏置(u，v)下，原始图像每个坐标处的加权核权值：

$w(x,y,u,v)=\exp (-\frac{D(x,y,u,v)}{2h^p}),$

其中，w(x，y，u，v)为加权核权值，h为加权核宽度，exp(·)为指数函数；

7)通过下式计算原始图像每个坐标处的权值归一化系数α(x，y)：

$\mathrm{\α}(x,y)=\underset{(u,v)\∈V_R}{\mathrm{\Σ}}w(x,y,u,v),$

其中，(u，v)为搜索域内的坐标，V_R为搜索域形状，w(x，y，u，v，j)为加权核权值；

8)通过下式计算计算原始图像每个坐标处的非归一化滤波后的值

$\stackrel{\‾}{I}(x,y)=\underset{(u,v)\∈V_R}{\mathrm{\Σ}}[w(x,y,u,v)\×I(x+u,y+v)],$

其中，I为原始图像信息；

9)通过下式对原始图像每个坐标进行归一化处理，得到滤波后图像

$\hat{I}(x,y)=\frac{\stackrel{\‾}{I}(x,y)}{\mathrm{\α}(x,y)}.$

本发明具有如下优点：

(1)现有非局部均值空域时变滤波方法在计算相邻的像素点和搜索域偏置组合的加权核权值的时候，会出现图像匹配窗覆盖的现象，从而会导致覆盖的匹配窗内像素点和加搜索域偏置像素组合之间，多次执行比较计算的计算冗余问题。而本发明由于对于搜索域内每个偏置坐标，采用技术方案1的步骤(4)～(7)或者技术方案2的步骤3)～6)，对于每个像素点和加搜索域偏置像素组合，只计算一次比较计算操作，消除了现有非局部均值空域时变滤波方法存在的由匹配窗覆盖导致的计算冗余的问题。所以本发明提出的非局部均值空域时变图像滤波方法，在不影响滤波性能的情况下，极大的降低了现有非局部均值空域时变滤波方法的计算量，加快了滤波速度，拓宽了非局部均值滤波器的应用范围。

附图说明

图1为现有非局部均值空域时变图像滤波方法流程图；

图2为本发明技术方案1的图像滤波方法流程图；

图3为本发明技术方案2的图像滤波方法流程图；

图4为图2和图3中，计算原始图像在当前偏置坐标下所有权值的子流程图；

图5为原始图像生成偏置图像原理图；

图6为差值图像中矩形区域描述图；

图7为在p＝1时，本发明与现有非局部均值空域时变滤波方法的滤波对比图；

图8为在p＝2时，本发明与现有非局部均值空域时变滤波方法的滤波对比图。

具体实施方式

参见图2，本发明的技术方案1提出的一种非局部均值空域时变图像滤波方法，包括如下步骤：

步骤一，获取原始图像信息，设置滤波相关参数，初始化非归一化滤波后的图像

和权值归一化系数α(x，y)。

设置的滤波相关参数包括搜索域形状、匹配窗大小、加权核宽度h和相似性比较范数参数p；该搜索域为矩形、三角形、圆形或任意多边形，该匹配窗形状为矩形，设置在根据图像偏置相似函数值计算加权核权值时的加权核宽度h和相似性比较范数参数p；该加权核宽度h为任意正实数，该相似性比较范数参数p为任意非负实数；将图像每个坐标处非归一化滤波后的图像的值

和权值归一化系数α(x，y)，均初始为0，其中(x，y)为图像坐标。

步骤二，以搜索域内的点作为偏置，对所有点执行步骤三至步骤五的操作；如果搜索域内所有点均已执行步骤三至步骤五的操作，跳转到步骤六。

步骤三，在当前偏置坐标下，计算原始图像所有坐标点处的加权核权值。

参见图4，本步骤的具体实现方式如下：

3.1)如图5所示，白色框表示为原始图像信息I，灰色框为加了(u，v)偏置的图像白色框和灰色框对应区域值相等，即偏置图像

坐标(x，y)处的值为：

$I_s^{(u,v)}(x,y)=I(x+u,y+v);$

差值图像Δ^(u，v)(x，y)通过下式计算得到：

Δ^(u，v)(x，y)＝|I(x，y)-I(x+u，y+v)|^p，

其中I为原始图像信息，|·|^p为求绝对值的p次方操作；

3.2)计算当前偏置(u，v)下差值图像Δ^(u，v)(x，y)的积分图像S^(u，v)(x，y)，计算公式为：

$S^{(u,v)}(x,y)=\underset{y^{\′}\≤y}{\mathrm{\Σ}}\underset{x^{\′}\≤x}{\mathrm{\Σ}}{\mathrm{\Δ}}^{(u,v)}(x^{\′},y^{\′});$

3.3)通过下式计算在当前偏置坐标下的原始图像每个坐标处的图像偏置相似函数值：

D(x，y，u，v)＝S^(u，v)(x_l，y_t)+S^(u，v)(x_r，y_b)-S^(u，v)(x_l，y_b)-S^(u，v)(x_r，y_t)，

其中(x，y)、(u，v)分别表示图像坐标和当前偏置坐标，D(x，y，u，v)为图像偏置相似函数值，S^(u，v)表示在当前偏置坐标下的积分图像，(x_l，y_t)，(x_r，y_b)分别表示以(x，y)中心的匹配窗形状的左上角坐标和右下角坐标，(x_r，y_t)和(x_l，y_b)分别表示以(x，y)中心的匹配窗形状的右上角坐标和左下角坐标，匹配窗形状在步骤一中被设置；

3.4)通过下式计算在当前偏置坐标下，原始图像每个坐标处的加权核权值：

$w(x,y,u,v)=\exp (-\frac{D(x,y,u,v)}{2h^p}),$

其中(x，y)、(u，v)分别表示图像坐标和当前偏置坐标，w(x，y，u，v)为加权核权值，h为加权核宽度，exp(·)为指数函数。

步骤四，在当前偏置坐标(u，v)下，通过下式更新原始图像每个坐标处的权值归一化系数：

α(x，y)＝α(x，y)+w(x，y，u，v)，

其中，α(x，y)为归一化系数。

步骤五，在当前偏置坐标下，通过下式更新原始图像每个坐标处的非归一化滤波后的值：

$\stackrel{\‾}{I}(x,y)=\stackrel{\‾}{I}(x,y)+w(x,y,u,v)\×I(x+u,y+v),$

其中，

为非归一化滤波后的值，I为原始图像信息。

通过步骤二至五，计算得到了非归一化滤波后的图像

和权值归一化系数α(x，y)。

步骤六，通过下式对原始图像每个坐标进行归一化处理，得到滤波后图像 $\hat{I}(x,y):$

$\hat{I}(x,y)=\frac{\stackrel{\‾}{I}(x,y)}{\mathrm{\α}(x,y)}.$

上述技术方案1所述的图像滤波方法，其中所述的步骤四与步骤五的顺序能够互换。

上述技术方案1所述的图像滤波方法，其中步骤三至五与步骤六，使用串行方式、并行方式或者串行与并行混合的方式进行计算。

上述技术方案1方法处理的原始图像为直接采集到的数字图像，或者是经过一定处理后的数字图像。

上述技术方案1可以应用于类似于非局部均值滤波器加权核权值计算方法的空域滤波方法中。步骤3.1)中，差值图像可以由原始图像和偏置图像之间所有对应点之间的任意初等数学变换获得。步骤3.4)中，加权核权值可以由对图像偏置相似函数值作任意初等数学变换获得。

上述技术方案1的滤波过程，也可以用下述技术方案2的步骤实现。

参见图3，本发明技术方案2提出的一种非局部均值空域时变图像滤波方法，包括如下步骤：

步骤1，获取原始图像信息，设置滤波相关参数。

设置的滤波相关参数包括搜索域形状、匹配窗大小、加权核宽度h和相似性比较范数参数p；该搜索域为矩形、三角形、圆形或任意多边形，该匹配窗形状为矩形；设置在根据图像偏置相似函数值计算加权核权值时的加权核宽度h和相似性比较范数参数p；该加权核宽度h为任意正实数，该相似性比较范数参数p为任意非负实数。

步骤2，以搜索域内的点作为偏置，对所有点执行步骤3操作；如果搜索域内所有点均已执行步骤3操作，跳转到步骤4。

步骤3，在当前偏置坐标下，计算原始图像所有坐标点处的加权核权值。

参见图4，本步骤的具体实现方式如下：

3a)如图5所示，白色框表示为原始图像信息I，灰色框为加了(u，v)偏置的图像

白色框和灰色框对应区域值相等，即偏置图像

坐标(x，y)处的值为：

$I_s^{(u,v)}(x,y)=I(x+u,y+v);$

差值图像Δ^(u，v)(x，y)通过下式计算得到：

Δ^(u，v)(x，y)＝|I(x，y)-I(x+u，y+v)|^p，

其中I为原始图像信息，|·|^p为求绝对值的p次方操作；

3b)计算当前偏置(u，v)下差值图像Δ^(u，v)(x，y)的积分图像S^(u，v)(x，y)，计算公式为：

$S^{(u,v)}(x,y)=\underset{y^{\′}\≤y}{\mathrm{\Σ}}\underset{x^{\′}\≤x}{\mathrm{\Σ}}{\mathrm{\Δ}}^{(u,v)}(x^{\′},y^{\′});$

3c)通过下式计算在当前偏置坐标下的原始图像每个坐标处的图像偏置相似函数值：

D(x，y，u，v)＝S^(u，v)(x_l，y_t)+S^(u，v)(x_r，y_b)-S^(u，v)(x_l，y_b)-S^(u，v)(x_r，y_t)，

其中(x，y)、(u，v)分别表示图像坐标和当前偏置坐标，D(x，y，u，v)为图像偏置相似函数值，S^(u，v)表示在当前偏置坐标下的积分图像，(x_l，y_t)，(x_r，y_b)分别表示以(x，y)中心的匹配窗形状的左上角坐标和右下角坐标，(x_r，y_t)和(x_l，y_b)分别表示以(x，y)中心的匹配窗形状的右上角坐标和左下角坐标，匹配窗形状在步骤一中被设置；

3d)通过下式计算在当前偏置坐标下，原始图像每个坐标处的加权核权值：

$w(x,y,u,v)=\exp (-\frac{D(x,y,u,v)}{2h^p}),$

其中(x，y)、(u，v)分别表示图像坐标和当前偏置坐标，w(x，y，u，v)为加权核权值，h为加权核宽度，exp(·)为指数函数。

通过步骤2～3，计算得到了在每个搜索域偏置下的所有原始图像坐标处的加权核权值，也就是获得了每个原始图像坐标处在所有搜索域偏置下的加权核权值。

步骤4，通过下式计算原始图像每个坐标处的权值归一化系数α(x，y)：

$\mathrm{\α}(x,y)=\underset{(u,v)\∈V_R}{\mathrm{\Σ}}w(x,y,u,v),$

其中，(u，v)为搜索域内的坐标，V_R为搜索域形状，w(x，y，u，v)为加权核权值。

步骤5，通过下式计算计算原始图像每个坐标处的非归一化滤波后的值

$\stackrel{\‾}{I}(x,y)=\underset{(u,v)\∈V_R}{\mathrm{\Σ}}[w(x,y,u,v)\×I(x+u,y+v)],$

步骤6，通过下式对原始图像每个坐标进行归一化处理，得到滤波后图像 $\hat{I}(x,y):$

$\hat{I}(x,y)=\frac{\stackrel{\‾}{I}(x,y)}{\mathrm{\α}(x,y)}.$

上述技术方案2所述的图像滤波方法，其中所述的步骤4与步骤5的顺序能够互换。

上述技术方案2所述的图像滤波方法，其中步骤3与步骤4～6，使用串行方式、并行方式或者串行与并行混合的方式进行计算。

上述技术方案2方法处理的原始图像为直接采集到的数字图像，或者是经过一定处理后的数字图像。

上述技术方案2的可以应用于类似于非局部均值滤波器加权核权值计算方法的空域滤波方法中。步骤3a)中，差值图像可以由原始图像和偏置图像之间所有对应点之间的任意初等数学变换获得。步骤3d)中，加权核权值可以由对图像偏置相似函数值通过任何初等数学变换获得。

下面以技术方案1为例，说明本发明的滤波方法能够实现非局部均值空域时变滤波器的滤波效果：

在技术方案1步骤一里设置的矩形区域形状为半径为r的正方形。

对于步骤3.1)获得的差值图像Δ^(u，v)(x，y)，如图6所示，若要计算差值图像Δ^(u，v)(x，y)中矩形块D内所有坐标点处数值之和，可以通过下面的式子获得：

D_s＝A_s+(A+B+C+D)_s-(A+B)_s-(A+C)_s，

其中D_s，A_s，(A+B+C+B)_s，(A+B)_s，(A+C)_s分别为相应矩形范围内所有坐标处数值之和。

显然，A_s，(A+B+C+B)_s，(A+B)_s，(A+C)_s的值分别为积分图像S^(u，v)中S^(u，v)(x_l，y_t)，S^(u，v)(x_r，y_b)，S^(u，v)(x_l，y_b)，S^(u，v)(x_r，y_t)的值。对比步骤3.3)可知，D(x，y，u，v)的值为差值图像Δ^(u，v)(x，y)中以(x，y)为中心的矩形区域的数据之和。结合步骤3.1)差值图像的计算过程，可知步骤3.3)得到的图像偏置相似函数值也可以表示为：

$D(x,y,u,v)={\left|\right|B_r(x,y)-B_r(x+u,y+v)\left|\right|}_p^p,$

其中B_r(x，y)为表示以(x，y)为中心，以r为半径的矩形块的匹配窗，

表示l_p范数值的p次方。

上面的结果说明，本发明的方法能够获得与现有非局部均值滤波空域滤波器相同的图像偏置相似函数值，据此也就获得了相同的加权核权值。

显然，在相同的加权核权值下，本发明的滤波方法与现有非局部均值空域时变滤波器的图像滤波方法有相同的滤波效果。

现有非局部均值空域时变滤波方法在计算相邻的像素点和搜索域偏置组合的加权核权值的时候，会出现图像匹配窗覆盖的现象，从而会导致覆盖的匹配窗内像素点和加搜索域偏置像素组合之间，多次执行比较计算的计算冗余问题。而本发明由于对于搜索域内每个偏置坐标，采用技术方案1的步骤3.1)～3.4)或者技术方案2的步骤3a)～3d)，对于每个像素点和加搜索域偏置像素组合，只计算一次比较计算操作，消除了现有非局部均值空域时变滤波方法存在的由匹配窗覆盖导致的计算冗余的问题。所以本发明提出的非局部均值空域时变图像滤波方法，在不影响滤波性能的情况下，极大的降低了现有非局部均值空域时变滤波方法的计算量，加快了滤波速度，拓宽了非局部均值滤波器的应用范围。

本发明的快速滤波效果可以通过图7和图8的仿真结果进一步说明。

在运行于3.0GHz Intel E8400GPU上的64位GNU/Linux操作系统上，分别用现有的非局部均值空域时变图像滤波方法和本发明提出的图像滤波方法，对同一组图像序列进行滤波处理的仿真实验。仿真程序采用的C/C++程序设计语言实现。

在仿真试验中，搜索域V_R形状设置为半径为R的正方形，匹配窗形状设置为半径为r的正方形。分别在相似性比较范数参数p＝1和p＝2的情况下，统计在R值为4、6、8、10和12，r为2、3、4、5、6和7的情况下，现有滤波方法和本发明提出的滤波方法的计算时间。p＝1时的结果见图7，p＝2时的结果见图8。从图7和图8的结果可以看出，无论在p＝1还是p＝2的情况下，本发明提出的方法比现有滤波方法的计算时间都有极大的降低；在固定的R值下，现有滤波方法的计算时间随着r的增加而迅速增加，而本发明的滤波方法的计算时间几乎不随r的增加而增加。特别是在R和r取较大值的情况下，本发明提出的方法快速滤波效果更为明显。

Claims (8)

1.一种非局部均值空域时变图像滤波方法，包括如下步骤：

(1)获取原始图像信息，设置搜索域形状、匹配窗大小、加权核宽度h和相似性比较范数参数p；

(2)初始化非归一化滤波后的图像

和权值归一化系数α(x，y)均为0，其中(x，y)为图像坐标；

(3)以搜索域内的点作为偏置，对所有点执行步骤(4)～(9)的操作；

(4)以当前偏置(u，v)为偏置量，计算原始图像与偏置后图像之间所有对应点的数据差的绝对值的p次方之值，得到p范数差值图像Δ^(u，v)(x，y)；

(5)计算当前偏置(u，v)下差值图像Δ^(u，v)(x，y)的积分图像S^(u，v)(x，y)，计算公式为：

$S^{(u,v)}(x,y)=\underset{y^{\′}\≤y}{\mathrm{\Σ}}\underset{x^{\′}\≤x}{\mathrm{\Σ}}{\mathrm{\Δ}}^{(u,v)}(x^{\′},y^{\′});$

(6)通过下式计算在当前偏置坐标(u，v)下的原始图像每个坐标处的图像偏置相似函数值：

D(x，y，u，v)＝S^(u，v)(x_l，y_t)+S^(u，v)(x_r，y_b)-S^(u，v)(x_l，y_b)-S^(u，v)(x_r，y_t)，

其中，(x，y)为图像坐标，(x_l，y_t)，(x_r，y_b)分别表示以(x，y)中心的匹配窗的左上角坐标和右下角坐标；

(7)通过下式计算在当前偏置(u，v)下，原始图像每个坐标处的加权核权值：

$w(x,y,u,v)=\exp (-\frac{D(x,y,u,v)}{2h^p}),$

其中，w(x，y，u，v)为加权核权值，h为加权核宽度，exp(·)为指数函数；

(8)在当前偏置坐标(u，v)下，通过下式更新原始图像每个坐标处的权值归一化系数：

α(x，y)＝α(x，y)+w(x，y，u，v)，

其中，α(x，y)为归一化系数；

(9)在当前偏置坐标下，通过下式更新原始图像每个坐标处的非归一化滤波后的值：

$\stackrel{\‾}{I}(x,y)=\stackrel{\‾}{I}(x,y)+w(x,y,u,v)\×I(x+u,y+v),$

其中，为非归一化滤波后的值，I为原始图像信息；

(10)通过下式对原始图像每个坐标进行归一化处理，得到滤波后图像

$\hat{I}(x,y)=\frac{\stackrel{\‾}{I}(x,y)}{\mathrm{\α}(x,y)}.$

2.根据权利要求1所述的图像滤波方法，其中的步骤(8)与步骤(9)的顺序能够互换。

3.根据权利要求1所述的图像滤波方法，其中步骤(1)所述的搜索域的形状设置为圆形或任意多边形。

4.根据权利要求1所述的图像滤波方法，其中步骤(4)～(9)与步骤(10)，使用串行方式、并行方式或者串行与并行混合的方式进行计算。

5.一种非局部均值空域时变图像滤波方法，包括如下步骤：

1)获取原始图像信息，设置搜索域形状、匹配窗大小、加权核宽度h和相似性比较范数参数p；

2)以搜索域内的点作为偏置，对所有点执行步骤3)～6)的操作；

3)以当前偏置(u，v)为偏置量，计算原始图像与偏置后图像之间所有对应点的数据差的绝对值的p次方之值，得到p范数差值图像Δ^(u，v)(x，y)，其中(x，y)为图像坐标；

4)计算当前偏置(u，v)下差值图像Δ^(u，v)(x，y)的积分图像S^(u，v)(x，y)，计算公式为：

$S^{(u,v)}(x,y)=\underset{y^{\′}\≤y}{\mathrm{\Σ}}\underset{x^{\′}\≤x}{\mathrm{\Σ}}{\mathrm{\Δ}}^{(u,v)}(x^{\′},y^{\′});$

5)通过下式计算在当前偏置坐标(u，v)下的原始图像每个坐标处的图像偏置相似函数值：

D(x，y，u，v)＝S^(u，v)(x_l，y_t)+S^(u，v)(x_r，y_b)-S^(u，v)(x_l，y_b)-S^(u，v)(x_r，y_t)，

其中，(x，y)为图像坐标，(x_l，y_t)，(x_r，y_b)分别表示以(x，y)中心的匹配窗的左上角坐标和右下角坐标；

6)通过下式计算在当前偏置(u，v)下，原始图像每个坐标处的加权核权值：

$w(x,y,u,v)=\exp (-\frac{D(x,y,u,v)}{2h^p}),$

其中，w(x，y，u，v)为加权核权值，h为加权核宽度，exp(·)为指数函数；

7)通过下式计算原始图像每个坐标处的权值归一化系数α(x，y)：

$\mathrm{\α}(x,y)=\underset{(u,v)\∈V_R}{\mathrm{\Σ}}w(x,y,u,v),$

其中，(u，v)为搜索域内的坐标，V_R为搜索域形状，w(x，y，u，v)为加权核权值；

8)通过下式计算计算原始图像每个坐标处的非归一化滤波后的值

$\stackrel{\‾}{I}(x,y)=\underset{(u,v)\∈V_R}{\mathrm{\Σ}}[w(x,y,u,v)\×I(x+u,y+v)],$

其中，I为原始图像信息；

9)通过下式对原始图像每个坐标进行归一化处理，得到滤波后图像

$\hat{I}(x,y)=\frac{\stackrel{\‾}{I}(x,y)}{\mathrm{\α}(x,y)}.$

6.根据权利要求5所述的图像滤波方法，其中的步骤7)与步骤8)的顺序能够互换。

7.根据权利要求5所述的图像滤波方法，其中步骤1)所述的搜索域的形状设置为圆形或任意多边形。

8.根据权利要求5所述的图像滤波方法，其中步骤3)～6)与步骤7)～9)，使用串行方式、并行方式或者串行与并行混合的方式进行计算。

Images