CN102509266B - 一种图像快速保边滤波方法 - Google Patents

Abstract

本发明公开一种图像处理技术领域的图像快速保边滤波方法，首先对滤波器采用设定的核函数；然后用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)，最后得到滤波结果；所述设定的核函数具有四点性质：在滤波模板内积分和为1、具有保边滤波的性质、可用FastBoxFilter加速计算使得计算复杂度为O(1)、可调节参数获得不同滤波效果。整个方法的时间复杂为O(1)，即与滤波核大小无关，且能达到与双边滤波相当的效果，并且结果证明其不会出现边缘人造物和光晕现象。

Description

一种图像快速保边滤波方法

技术领域

本发明属于图像处理领域，具体涉及的是一种图像保边滤波方法。

背景技术

图像滤波是图像处理以及计算机视觉邻域最基本的操作，其能滤除图像中的噪声而使得图像更加光滑和清晰。常用的滤波方法有线性滤波器和非线性滤波器。线性滤波器有均值滤波器，带权均值滤波器，以及高斯滤波器，而非线性滤波器的典型代表是中值滤波器。对于线性滤波器，如均值滤波和高斯滤波，其在去除噪声的同时使得整个图像变得模糊，包括图像的边界，在有些需要滤除噪声并保持区域边界的地方运用受到限制，如轮廓提取，纹理提取，区域分割等。故针对此种情况，C.Tomasi和R.Manduchi在1998年提出了能保持图像边缘的滤波器——双边滤波器(Bilateral Filtering)，其内在想法是在图像的空域和值域上都进行滤波，其表达式为：

$y_h\left(p\right)=\frac{1}{K_b}{\mathrm{\Σ}}_{k\∈S}f\left(k\right)g(I\left(p\right)-I(p+k))I(p+k)$

其中p，p+k表示图像位置向量，k为位置增量，S为以原点为中心，半径为r的正方形邻域，I(p)表示原图像在p位置处的像素值；f，g为空域和值域的滤波核函数，在C.Tomasi等人的方法中取g为高斯核函数：

g(I(p)-I(p+k))＝exp(-x[I(p)-I(p+k)]²)

而f可取高斯核函数，也可取任意线性滤波器的核函数。

双边滤波能够在保持强边缘的同时有效的对图像的细小变化进行平滑，所以其被用于噪声去除，细节分解，HDR压缩和图像抽象。但是其也存在着一些不足，其中一点是当它被用于细节分解和HDR压缩时，会引起梯度反转人造物，常表现为光晕，这是因为当一个像素点周围的相似像素很少时，其高斯加权平均是不稳定的。

另外双边滤波的效率偏低，其计算时间复杂度与核半径r有关，当调节滤波器的核半径r时会明显感受到相应的计算时间上的变化，这是双边滤波器技术的主要特征差异。现有的双边滤波是结合了空间域和像素值域的滤波，所以它的计算量要求非常大。一个标准的双边滤波器的实现计算整个图像输出像素的复杂O(Nr²)，当固定图像大小N时，复杂度为O(r²)，意味着当r过大时滤波器的计算将会非常的缓慢。

近年来有不少关于双边滤波器的加速办法，众所周知的Photoshop CS2的表面模糊功能所用到的双边滤波器时间复杂度为O(r)，Pham和Vliet在05年提出用1维的双边滤波器分别对各个空间方向进行滤波，其是对2维双边滤波器的近似，计算速度很快，但在图像边缘会产生人工痕迹，B.Weiss在06年提出用一个基于3D直方图的方法，并用此方法来加速正方形核空间上的二维双边滤波，其思想是改变图像直方图的表达方式来达到加速图像处理的过程，同时达到在任意位宽的图像，任意大小的滤波核的情况下其时间复杂度都为O(logr)。在08年，F.Porikli同样介绍了基于三维直方图的一种低计算复杂度的方法，F.Porikli提出常数时间双边滤波方法，即其计算时间复杂度为常数O(1)，与核大小r无关，其文章针对不同特征的双边滤波(即不同类型核函数)采用了不同的证明方法，结果均能使得双边滤波达到O(1)，但是其计算时间却与直方图个数有关，且是对双边滤波的一种近似，并且在HDR压缩中可能产生人造物(表现为光晕)。

发明内容

本发明要解决的技术问题是提出一种图像快速保边滤波方法，能达到双边滤波器的滤波效果，对图像可以使用任意大的r获得不同的保边滤波效果，而其计算时间不变。

本发明是通过以下技术方案来实现的，本发明首先对滤波器采用设定的核函数；然后用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)，最后得到滤波结果。

所述的设定的核函数是指：

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}]$

其中I_i为模板W中心处的像素，I_j为以I_i为中心、r为半径的模板W里的任一像素，μ_i为以I_i为中心、r为半径的模板W中像素的均值，

为模板W中像素的方差，|W|为模板W的面积，ε为大于0的调整参数。本函数可以做各类简易的变换；如将函数式进行通分，多项式展开等操作；又如将ε取为关于

的非负减函数，函数值随着

的增加而减小，减小的幅度则决定了不同区域滤波后界限的分明程度，减小得越快，则界限越分明；又如将函数扩充到RGB彩色空间，对应核函数为

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+{(I_i-{\mathrm{\μ}}_i)}^T{(H_i+\mathrm{\ϵ}*E)}^{-1}(I_j-{\mathrm{\μ}}_i)]$

这里I_i，I_j是3×1颜色向量，μ_i是模板中3×1的均值向量，H_i是模板中3×3的颜色协方差矩阵，E是3×3的单位矩阵。

所述的用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)具体是指：

首先给出对图像采用半径为r的FastBoxFilter滤波方法：

1)计算图像的累积和矩阵T

2)设输出图像为Q，则将输出图像坐标为(i，j)的元素Q(i，j)用T(i+r，j+r)，T(i-r，j+r)，T(i+r，j-r)，T(i-r，j-r)来计算，计算公式为：

Q(i，j)＝T(i+r，j+r)-T(i-r，j+r)-T(i+r，j-r)+T(i-r，j-r)

设对图像I(I为矩阵)运用半径为r的FastBoxFilter得到的结果矩阵表示为FastBoxFilter(I，r)，则对本发明提出的滤波器，设输入图像为I，计算滤波后图像的步骤如下：

1)计算FastBoxFilter(I，r)，FastBoxFilter(I.*I，r)的结果，其中I.*I是指图像I每一个对应像素值的平方所得到的矩阵。

2)计算均值矩阵和方差矩阵

μ＝FastBoxFilter(I，r)/(2r+1)²

σ²＝FastBoxFilter(I.*I，r)/(2r+1)²-μ.*μ

3)将计算到的各个参数带入∑_jw_ij*I_j中，∑_jw_ij*I_j如下式所示：

${\mathrm{\Σ}}_j{w}_{\mathrm{ij}}*I_j={\mathrm{\μ}}_i+\frac{{\mathrm{\σ}}_i^2*(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}$

根据下式计算输出图像Q：

Q＝μ+(σ².*(I-μ))./(σ²+ε)

其中A.*B表示A，B矩阵对应位置元素相乘得到的矩阵，A./B表示A，B矩阵对应位置元素相除得到的矩阵。

本发明通过使用所提出的滤波器核函数，该核函数具有四点性质：在滤波模板内积分和为1、具有保边滤波的性质、可用FastBoxFilter加速计算使得计算复杂度为O(1)、可调节参数获得不同滤波效果。

采用上述技术方案，本发明的有益效果是：整个方法的时间复杂为O(1)，与滤波核大小无关，因此对图像可以使用任意大的r获得不同的保边滤波效果，而其计算时间不变；且能达到与双边滤波相当的效果，并且结果证明其不会出现人造物(表现为光晕)现象。

附图说明

图1为本发明实施例中核函数模板示意图。

图2为本发明实施例中核函数保边滤波示意图。

图3为采用计算累积图像的方法的原理图。

图4为本发明实施例结果示意图一。

图5为本发明实施例结果示意图二。

具体实施方式

下面对本发明的实施例作详细说明，本实施例以本发明技术方案为前提进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

本发明首先对滤波器采用设定的核函数；然后用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)，最后得到滤波结果。

如图1所示，为本实施例核函数模板示意图，所述的核函数为：

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}]$

其中I_i为模板W中心处的像素，I_j为以I_i为中心、r为半径的模板W里的任一像素，μ_i为以I_i为中心、r为半径的模板W中像素的均值，

为模板W中像素的方差，|W|为模板W的面积，ε为大于0的调整参数。

可证明该核函数具有下述4个性质：

1.其在模板内的积分为1

证明如下：

${\mathrm{\Σ}}_j{w}_{\mathrm{ij}}={\mathrm{\Σ}}_j\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_k^2+\mathrm{\ϵ}}]=\frac{1}{\left|W\right|}[{\mathrm{\Σ}}_{j}1+{\mathrm{\Σ}}_j\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}]$

$=\frac{1}{\left|W\right|}\left[\right|W|+\frac{(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}{\mathrm{\Σ}}_j(i_j-{\mathrm{\μ}}_i)]=\frac{1}{\left|W\right|}*\left|W\right|=1$

在模板内的积分为1，则不需要额外的对模板中的核函数进行归一化，简化了计算，这是本发明用到核函数较其它核函数的一个优点。

2.其具有保边滤波的性质，达到双边滤波器相当的效果

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}]$

如图2所示，取1维信号为例子，当I_i恰好取在边界附近时，如下图所示，当I_j与I_i在边界的同侧时I_i-μ_i与I_j-μ_i同号，式

接近1，此时I_j的权值w_ij较大；当I_j与I_i在边界的异侧时，I_i-μ_i与I_j-μ_i反号，

接近-1，此时I_j的权值w_ij较小(接近于0)。

综上其具有保边滤波的性质，即与模板中心点在边界同侧的像素系数较大，而在异侧的像素系数较小。

3.其可采用计算累积图像(Integral Image)的方法使得时间复杂度为O(1)

如图3所示，为一幅图像的累积和矩阵T，要计算原图像中对应白色区域内的像素和Sum(Area)，则在累积和矩阵中只需要用到区域四个角上的值即，T(Xr，Yt)，T(Xr，Yb)，T(Xl，Yt)，T(Xl，Yb)，则白色区域的像素和由下式给出：

Sum(Area)＝T(Xr，Yt)-T(Xr，Yb)-T(Xl，Yt)+T(Xl，Yb)

由于一幅图像的累积和矩阵可固定计算，白色区域的面积再大，也只需要用到四个角上的值，故可知某一区域的像素和计算的时间复杂度与区域大小r无关。将用上述方法计算像素和的方法称FastBoxFilter，则用FastBoxFilter对图像滤波的计算时间复杂度与模板半径r无关，可达到O(1)。

下面说明上述的滤波核函数可用FastBoxFilter计算使其时间复杂达到O(1)，对于输出像素Q_i＝∑_jw_ij*I_j，代入w_ij的定义得到

$Q_i={\mathrm{\Σ}}_j{w}_{\mathrm{ij}}*I_j={\mathrm{\Σ}}_j\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_k^2+\mathrm{\ϵ}}]*I_j$

$=\frac{1}{\left|W\right|}[{\mathrm{\Σ}}_j{I}_j+\frac{(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}{\mathrm{\Σ}}_j(I_j-{\mathrm{\μ}}_i)*I_j]$

$=\frac{1}{\left|W\right|}[{\mathrm{\Σ}}_j{I}_j+\frac{(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}({\mathrm{\Σ}}_j{I}_j*I_j-{\mathrm{\μ}}_i{\mathrm{\Σ}}_j{I}_j)]$

$=\frac{1}{\left|W\right|}[{\mathrm{\Σ}}_j{i}_j+\frac{(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}*|W\left|{\mathrm{\σ}}_i^2\right]={\mathrm{\μ}}_i+\frac{{\mathrm{\σ}}_i^2*(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}$

其中∑_jI_j与∑_jI_j*I_j均可以用FastBoxFilter计算，故参数μ_i与

也均可以用FastBoxFilter计算，所以对整个输出图像Q_i的计算与核大小r无关，实际上该核函数w_ij是关于I_j的1次多项式核函数，也就是关于I_j-I_i的1次多项式核函数，可以证明，关于I_j的任一多项式核函数均可用FastBoxFilter优化，即其时间复杂度与核大小r无关。

4.调整参数ε(ε＞0)用于调节核函数大小以取得不同的滤波效果

根据下式

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}]$

$Q_i={\mathrm{\μ}}_i+\frac{{\mathrm{\σ}}_i^2*(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}$

固定滤波核模板半径r，对于w_ij，当ε增大时，若I_i与I_j在边界的同侧，则(I_i-μ_i)(I_j-μ_i)为正，w_ij减小；若I_i与I_j在边界的异侧，则(I_i-μ_i)(I_j-μ_i)为负，w_ij增大；换言之ε增大使得边界两边像素的系数差距减小，使得滤波效果变得更平滑，当ε趋于+∞时，本发明提出滤波器则退化为均值滤波器。当ε＝0时，Q_i＝I_i，即输出图像等于输入图像，则不具有滤波效果，且当模板内像素值均相同时，除数

所以一般取ε＞0，若在模板中

远大于ε，则Q_i接近于I_i，在此模板中几乎不具有滤波效果；若在模板中

远小于ε，则Q_i接近于μ_i，在此模板中具有近似均值滤波效果。

所述的用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)具体是指：

首先给出对图像采用半径为r的FastBoxFilter滤波方法

1.计算图像的累积和矩阵T

2.设输出图像为Q，则将输出图像坐标为(i，j)的元素Q(i，j)用T(i+r，j+r)，T(i-r，j+r)，T(i+r，j-r)，T(i-r，j-r)来计算，计算公式为：

Q(i，j)＝T(i+r，j+r)-T(i-r，j+r)-T(i+r，j-r)+T(i-r，j-r)

设对图像I(I为矩阵)运用半径为r的FastBoxFilter得到的结果矩阵表示为FastBoxFilter(I，r)，则对本发明提出的滤波器，设输入图像为I，计算滤波后图像的步骤如下：

1.计算FastBoxFilter(I，r)，FastBoxFilter(I.*I，r)的结果，其中I.*I是指图像I每一个对应像素值的平方所得到的矩阵。

2.计算均值矩阵和方差矩阵

μ＝FastBoxFilter(I，r)/(2r+1)²

σ²＝FastBoxFilter(I.*I，r)/(2r+1)²-μ.*μ

3.将计算到的各个参数带入∑_jw_ij*I_j中，∑_jw_ij*I_j如下式所示：

${\mathrm{\Σ}}_j{w}_{\mathrm{ij}}*I_j={\mathrm{\μ}}_i+\frac{{\mathrm{\σ}}_i^2*(I_i-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}$

根据下式计算输出图像Q：

Q＝μ+(σ².*(I-μ))./(σ²+ε)

其中A.*B表示A，B矩阵对应位置元素相乘得到的矩阵，A./B表示A，B矩阵对应位置元素相除得到的矩阵。

本实施例采用上述核函数对输入图像进行保边滤波，设对图像I(I为矩阵)运用半径为r的FastBoxFilter得到的结果矩阵表示为FastBoxFilter(I，r)，具体实施步骤如下：

1.以矩阵的形式读入一幅图像，将其像素值归一化后存入矩阵I。

2.计算FastBoxFilter(I，r)，FastBoxFilter(I.*I，r)的结果，其中I.*I是指图像I每一个对应像素值的平方所得到的矩阵。

3.计算均值矩阵和方差矩阵

3.μ＝FastBoxFilter(I，r)/(2r+1)²

4.σ²＝FastBoxFilter(I.*I，r)/(2r+1)²-μ.*μ

4.根据下式计算输出图像Q

Q＝μ+(σ².*(I-μ))./(σ²+ε)

并将Q的像素值范围重新调整至0-255。

图4与图5给出了实验结果，如图4所示：原图(左图)，经过本实施例提出的保边滤波(r＝10ε-0.04)(中图)，与均值滤波(r＝10)(右图)。如图5所示：原图(左上图)，经过本实施例提出的保边滤波(r＝10ε-0.01)(右上图)，双边滤波(左下图)与均值滤波(r＝10)(右下图)结果。

以上实施例仅对灰度图像进行，如果输入是彩色图像，则只需将本发明提出的保边滤波核函数扩展到RGB空间，w_ij与Q_i对应表达式为：

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+{(I_i-{\mathrm{\μ}}_i)}^T{(H_i+\mathrm{\ϵ}*E)}^{-1}(I_j-{\mathrm{\μ}}_i)]$

Q_i＝μ_i+H_i(H_i+ε*E)^-1(I_i-μ_i)

这里I_i，I_j是3×1颜色向量，μ_i是模板中3×1的均值向量，H_i是模板中3×3的颜色协方差矩阵，E是3×3的单位矩阵。

通过以上实施例可以看出，本发明使用上述滤波器核函数，再使用FastBoxFilter对滤波过程进行加速，能达到双边滤波器的滤波效果，方法简单高效，不需要额外的核函数归一化，并且计算时间复杂度与滤波器的核半径r无关，可证明其为O(N)，当图像大小N固定时，其复杂度为O(1)，因此对图像可以使用任意大的r获得不同的保边滤波效果，而其计算时间不变，并且结果证明其不会出现人造物(表现为光晕)现象，所得图像比现有技术处理后的图像更加光滑和清晰。

当然，本发明提出上述滤波器核函数还应当包括该函数所做的各类简易的变换：如将函数式进行通分，多项式展开等操作；又如将ε取为关于

的非负减函数，函数值随着

的增加而减小，减小的幅度则决定了不同区域滤波后界限的分明程度，减小得越快，则界限越分明。这些对于本领域的技术人员来说是很容易做到的。

尽管本发明的内容已经通过上述优选实施例作了详细介绍，但应当认识到上述的描述不应被认为是对本发明的限制。在本领域技术人员阅读了上述内容后，对于本发明的多种修改和替代都将是显而易见的。因此，本发明的保护范围应由所附的权利要求来限定。

Claims (3)

1.一种图像快速保边滤波方法，其特征在于具体为：

首先对滤波器采用设定的核函数；

然后用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)，最后得到滤波结果；

所述设定的核函数为：

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+\frac{(I_i-{\mathrm{\μ}}_i)(I_j-{\mathrm{\μ}}_i)}{{\mathrm{\σ}}_i^2+\mathrm{\ϵ}}]$

其中I_i为模板W中心处的像素，I_j为以I_i为中心、r为半径的模板W里的任一像素，μ_i为以I_i为中心、r为半径的模板W中像素的均值，

为模板W中像素的方差，|W|为模板W的面积，ε为大于0的调整参数；

所述用FastBoxFilter对滤波过程进行加速，使得其计算时间复杂度为O(1)，是指：

(1)计算FastBoxFilter(I,r)，FastBoxFilter(I.*I,r)的结果，其中I.*I是指图像I每一个对应像素值的平方所得到的矩阵；

(2)计算均值矩阵和方差矩阵

1.μ=FastBoxFilter(I,r)/(2r+1)²

2.σ²=FastBoxFilter(I.*I,r)/(2r+1)²-μ.*μ

(3)根据下式计算输出图像Q

Q=μ+(σ².*(I-μ))./(σ²+ε)

其中A.*B表示A，B矩阵对应位置元素相乘得到的矩阵，A./B表示A，B矩阵对应位置元素相除得到的矩阵；

对图像采用半径为r的FastBoxFilter滤波，其中：

1)计算图像的累积和矩阵T；

2)设输出图像为Q，则将输出图像坐标为(i,j)的元素Q(i,j)用T(i+r,j+r),T(i-r,j+r),T(i+r,j-r),T(i-r,j-r)来计算，计算公式为：

Q(i,j)=T(i+r,j+r)-T(i-r,j+r)-T(i+r,j-r)+T(i-r,j-r)

设对图像I，I为矩阵，运用半径为r的FastBoxFilter得到的结果矩阵表示为FastBoxFilter(I,r)。

2.根据权利要求1所述的图像快速保边滤波方法，其特征在于：所述设定的核函数包括对该函数式进行通分，多项式展开操作，或者将ε取为关于

的非负减函数，函数值随着

的增加而减小，减小的幅度则决定了不同区域滤波后界限的分明程度，减小得越快，则界限越分明；或者将函数扩充到RGB彩色空间，对应核函数为

$w_{\mathrm{ij}}=\frac{1}{\left|W\right|}[1+{(I_i-{\mathrm{\μ}}_i)}^T{(H_i+\mathrm{\ϵ}*E)}^{-1}(I_j-{\mathrm{\μ}}_i)]$

这里I_i，I_j是3×1颜色向量,μ_i是模板中3×1的均值向量，H_i是模板中3×3的颜色协方差矩阵,E是3×3的单位矩阵。

3.根据权利要求1或2所述的图像快速保边滤波方法，其特征在于：所述设定的核函数具有四点性质：在滤波模板内积分和为1，具有保边滤波的性质，可用FastBoxFilter加速计算使得计算复杂度为O(1)，可调节参数获得不同滤波效果。

Images