双边滤波_双边滤波matlab - CSDN
精华内容
参与话题
  • 双边滤波算法介绍与实现

    千次阅读 2017-01-16 22:30:59
    Bilateral filter&Adaptive bilateral filtering双边滤波的思想是抑制与中心像素差别太大的像素,简单的来说:起到滤波保边的效果。目录Bilateral filterAdaptive bilateral filtering 双边滤波算法 思想 数学表达式...

    Bilateral filter&Adaptive bilateral filtering

    双边滤波的思想是抑制与中心像素差别太大的像素,简单的来说:起到滤波保边的效果。

    目录

    双边滤波算法

    思想:

    双边滤波算法考虑到了空间域和值域两个方面,实现双边滤波的基本思路是:

    对于模板的各个点来说:
    
        1、从空间域出发,计算出模板的点与目标点的距离(利用勾股定理),然后计算得出空域权重
        2、从值域出发,计算出模板的点与目标点值的差值(取绝对值),然后计算得到值域权重
        4、同时让模板各个点的像素值乘以该点的w,并加起来得到sum_pixel
        5、最后让sum_pixel除以sumw就得到目标点最终的像素值了
    

    数学表达式:

    这里写图片描述

    其中i,j是变化的,(i, j)对应模板的各个点,(k,l)是你要求的像素点,也就是模板的中心点

    代码实现

    /*
    sigma_color参数是上式中的σd, sigma_space参数是上式中的σr
    */
    #pragma once
    #include "core/core.hpp"  //3个opencv库的头文件
    #include "highgui/highgui.hpp"  
    #include "imgproc/imgproc.hpp"  
    #include <iostream>  
    
    using namespace cv;
    static void bilateralFilter_8u(const Mat& src, Mat& dst, int d, double sigma_color, double sigma_space, int borderType)
    {
        int cn = src.channels();//获得图片的通道数
        int i, j, k, maxk, radius;
        Size size = src.size();
        CV_Assert((src.type() == CV_8UC1 || src.type() == CV_8UC3) &&
            src.type() == dst.type() && src.size() == dst.size() &&
            src.data != dst.data);
    
    
        if (sigma_color <= 0)
            sigma_color = 1;
        if (sigma_space <= 0)
            sigma_space = 1;
    
    
        double gauss_color_coeff = -0.5 / (sigma_color*sigma_color);//分母值
        double gauss_space_coeff = -0.5 / (sigma_space*sigma_space);//分母值
    
    
        if (d <= 0)
            radius = cvRound(sigma_space*1.5);//进行四舍五入
        else
            radius = d / 2;
        radius = MAX(radius, 1);
        d = radius * 2 + 1;
    
    
        Mat temp;
        copyMakeBorder(src, temp, radius, radius, radius, radius, borderType);
        /*复制图像并且制作边界。(处理边界卷积)
        目的是为了放大图像好做图像的边界处理*/
    
    
        vector<float> _color_weight(cn * 256);//用来存放值域差值对应的权重
        vector<float> _space_weight(d*d);//用来存放空间距离对应的权重
        vector<int> _space_ofs(d*d);//用来存放模板各点与锚点(中心点)的偏移量
        float* color_weight = &_color_weight[0];
        float* space_weight = &_space_weight[0];
        int* space_ofs = &_space_ofs[0];
    
    
        // initialize color-related bilateral filter coefficients  函数1
        //由于sigma_color已经给定,所以可以先算出差值为0-255时,对应的高斯相似度权重。 
        for (i = 0; i < 256 * cn; i++)
            color_weight[i] = (float)std::exp(i*i*gauss_color_coeff);
    
    
        // initialize space-related bilateral filter coefficients  函数2
        //由于sigma_space已经给定,所以一旦选定好模板就可计算出高斯距离权重了
        for (i = -radius, maxk = 0; i <= radius; i++)
            for (j = -radius; j <= radius; j++)
            {
                double r = std::sqrt((double)i*i + (double)j*j);
                /*if (r > radius)
                    continue;*/
                space_weight[maxk] = (float)std::exp(r*r*gauss_space_coeff);
                space_ofs[maxk++] = (int)(i*temp.step + j*cn);
            }
    
    
        for (i = 0; i < size.height; i++)
        {
            const uchar* sptr = temp.data + (i + radius)*temp.step + radius*cn;
            uchar* dptr = dst.data + i*dst.step;
        ///灰度图
            if (cn == 1)
            {
                for (j = 0; j < size.width; j++)
                {
                    float sum = 0, wsum = 0;
                    int val0 = sptr[j];
                    for (k = 0; k < maxk; k++)
                    {
                        int val = sptr[j + space_ofs[k]];
                        float w = space_weight[k] * color_weight[std::abs(val - val0)];
                        sum += val*w;
                        wsum += w;
                    }
                    // overflow is not possible here => there is no need to use CV_CAST_8U  
                    dptr[j] = (uchar)cvRound(sum / wsum);
                }
            }
            else
            {
            //彩色图
                assert(cn == 3);
                for (j = 0; j < size.width * 3; j += 3)
                {
                    float sum_b = 0, sum_g = 0, sum_r = 0, wsum = 0;
                    int b0 = sptr[j], g0 = sptr[j + 1], r0 = sptr[j + 2];
                    for (k = 0; k < maxk; k++)
                    {
                        const uchar* sptr_k = sptr + j + space_ofs[k];
                        int b = sptr_k[0], g = sptr_k[1], r = sptr_k[2];
                        float w = space_weight[k] * color_weight[std::abs(b - b0) +
                            std::abs(g - g0) + std::abs(r - r0)];
                        sum_b += b*w; sum_g += g*w; sum_r += r*w;
                        wsum += w;
                    }
                    wsum = 1.f / wsum;
                    b0 = cvRound(sum_b*wsum);
                    g0 = cvRound(sum_g*wsum);
                    r0 = cvRound(sum_r*wsum);
                    dptr[j] = (uchar)b0; dptr[j + 1] = (uchar)g0; dptr[j + 2] = (uchar)r0;
                }
            }
        }
    } 

    自适应双边滤波

    介绍

    高斯距离权重值还受到模块的的大小影响( space_weight[maxk] = (float)std::exp(r*r*gauss_space_coeff);),因此高斯距离权值受到一定的约束,但是高斯相似度权值则没有这样的约束。
    因此为了给高斯相似度权值增加约束,让其根据图片进行自适应,而不人为进行指定,就得根据式子:

    sigma_color =((模板内所有的像素值的平方的和)*(模板的点的数量)-(模板所有的像素值的和)^2 ) / (模板的点的数量)^2

    算出sigma_color的值:

    double sigma_color= ((sumsqr*n) - sum*sum) / ((double)(n*n));

    然后再计算gauss_color_cofee的值

    double gauss_color_coeff = -0.5 / (double)(sigma_color*sigma_color);

    之后根据模板的点与运算点的差值,算出相似度权值:

    colorr_Weight = (double)std::exp((std::abs(r - r0))*(std::abs(r - r0))*gauss_colorr_coeff);

    代码实现:

    (与opencv源码的实现方式有所不同,当算法的思想是一致的)

    //获得自适应值
    double adaptive_gausscolor(float sum, float sumsqr, int n, double maxSigma_color)
    {
        double sigma_color= ((sumsqr*n) - sum*sum) / ((double)(n*n));
        if (sigma_color < 0.01)sigma_color = 0.01;
        else
            if (sigma_color >(double)maxSigma_color)
                sigma_color = (double)maxSigma_color;
        double gauss_color_coeff = -0.5 / (double)(sigma_color*sigma_color);
        return gauss_color_coeff;
    }
    自适应双边滤波算法
    void adaptivebilateralFilter(Mat& src, Mat& dst, int d, double sigma_space, double sigmacolor_max, int borderType)
    {
        int cn = src.channels();//获得图片的通道数
        int i, j, k, maxk, radius;
        Size size = src.size();
        CV_Assert((src.type() == CV_8UC1 || src.type() == CV_8UC3) &&
            src.type() == dst.type() && src.size() == dst.size() &&
            src.data != dst.data);
    
        if (sigma_space <= 0)
            sigma_space = 1;
        CV_Assert(d & 1);//确保为奇数
        double gauss_space_coeff = -0.5 / (sigma_space*sigma_space);
    
        if (d <= 0)
            radius = cvRound(sigma_space*1.5);//进行四舍五入
        else
            radius = d / 2;
        radius = MAX(radius, 1);
        d = radius * 2 + 1;
    
        Mat temp;
        copyMakeBorder(src, temp, radius, radius, radius, radius, borderType);
        /*复制图像并且制作边界。(处理边界卷积)
        目的是为了放大图像好做图像的边界处理*/
    
        vector<float> _space_weight(d*d);
    
    
        vector<int> _space_ofs(d*d);
        float* space_weight = &_space_weight[0];
        int* space_ofs = &_space_ofs[0];
         // initialize color-related bilateral filter coefficients  函数1
        //由于sigma_color没有给定,所以无法算出差值为0-255时,对应的高斯相似度权重。 
        /*for (i = 0; i < 256 * cn; i++)
            color_weight[i] = (float)std::exp(i*i*gauss_color_coeff);*/
    
     // initialize space-related bilateral filter coefficients 
        //由于sigma_space已经给定,所以一旦选定好模板就可计算出高斯距离权重了
    
        for (i = -radius, maxk = 0; i <= radius; i++)
            for (j = -radius; j <= radius; j++)
            {
                double r = std::sqrt((double)i*i + (double)j*j);
                space_weight[maxk] = (float)std::exp(r*r*gauss_space_coeff);
                space_ofs[maxk++] = (int)(i*temp.step + j*cn);
            }
    
        for (i = 0; i < size.height; i++)
        {
            const uchar* sptr = temp.data + (i + radius)*temp.step + radius*cn;
            uchar* dptr = dst.data + i*dst.step;
    
            if (cn == 1)
            {
                for (j = 0; j < size.width; j++)
                {
                    float sum = 0, wsum = 0, sumValsqr = 0;
                    int val0 = sptr[j];
                    //获得自适应的高斯相似度的sigma值并计算相似度
                    for (k = 0; k < maxk; k++)
                    {
                        int val = sptr[j + space_ofs[k]];
                        sum += val;
                        sumValsqr += (val*val);
                    }
                    double gauss_color_coeff=adaptive_gausscolor(sum, sumValsqr, maxk, sigmacolor_max);
    
                    sum = 0;
                    for (k = 0; k < maxk; k++)
                    {
                        int val = sptr[j + space_ofs[k]];
                        int temp = std::abs(val - val0);
                        float color_Weight =(float)std::exp((float)temp*temp*gauss_color_coeff);
                        float w = space_weight[k] * color_Weight;
                        sum += val*w;
                        wsum += w;
                    }
                    // overflow is not possible here => there is no need to use CV_CAST_8U
                    dptr[j] = (uchar)cvRound(sum / wsum);
                }
            }
            else
            {
                assert(cn == 3);
                for (j = 0; j < size.width * 3; j += 3)
                {
                    float sum_b = 0, sum_g = 0, sum_r = 0, wbsum = 0, wgsum = 0, wrsum = 0, sum_bsqr = 0, sum_gsqr = 0, sum_rsqr = 0;
                    int b0 = sptr[j], g0 = sptr[j + 1], r0 = sptr[j + 2];
                    //获得自适应的高斯相似度的sigma值并计算相似度
                    for (k = 0; k < maxk; k++)
                    {
                        const uchar* sptr_k = sptr + j + space_ofs[k];
                        int b = sptr_k[0], g = sptr_k[1], r = sptr_k[2];
                        sum_b += b; sum_g += g; sum_r += r;
                        sum_bsqr += b*b; sum_gsqr += g*g; sum_rsqr += r*r;
                    }
                    double gauss_colorb_coeff=adaptive_gausscolor(sum_b, sum_bsqr, maxk, sigmacolor_max);
                    double gauss_colorg_coeff=adaptive_gausscolor(sum_g, sum_gsqr, maxk, sigmacolor_max);
                    double gauss_colorr_coeff=adaptive_gausscolor(sum_r, sum_rsqr, maxk, sigmacolor_max);
    
                    sum_b = 0; sum_g = 0; sum_r = 0;
                    for (k = 0; k < maxk; k++)
                    {
                        const uchar* sptr_k = sptr + j + space_ofs[k];
                        int b = sptr_k[0], g = sptr_k[1], r = sptr_k[2];
    
                        double colorb_Weight = (double)std::exp((std::abs(b - b0))*(std::abs(b - b0))*gauss_colorb_coeff);
                        double colorg_Weight = (double)std::exp((std::abs(g - g0))*(std::abs(g - g0))*gauss_colorg_coeff);
                        double colorr_Weight = (double)std::exp((std::abs(r - r0))*(std::abs(r - r0))*gauss_colorr_coeff);
    
    
                        float wb = space_weight[k] * colorb_Weight;
                        float wg= space_weight[k] * colorg_Weight;
                        float wr=space_weight[k] * colorr_Weight;
                        sum_b += b*wb; sum_g += g*wg; sum_r += r*wr;
    
                        wbsum += wb;
                        wgsum += wg;
                        wrsum += wr;
    
                    }
                    wbsum = 1.f / wbsum;
                    wgsum = 1.f / wgsum;
                    wrsum = 1.f / wrsum;
                    b0 = cvRound(sum_b*wbsum);
                    g0 = cvRound(sum_g*wgsum);
                    r0 = cvRound(sum_r*wrsum);
    
                    dptr[j] = (uchar)b0; dptr[j + 1] = (uchar)g0; dptr[j + 2] = (uchar)r0;
                }
            }
        }
    }
    
    展开全文
  • 双边滤波

    千次阅读 2018-01-23 17:06:20
    双边滤波(Bilateral filter)是一种可以保边去噪的滤波器。之所以可以达到此去噪效果,是因为滤波器是由两个函数构成。一个函数是由几何空间距离决定滤波器系数。另一个由像素差值决定滤波器系数。可以与其相比较的...

    http://blog.csdn.net/abcjennifer/article/details/7616663

    双边滤波器是什么?

    双边滤波(Bilateral filter)是一种可以保边去噪的滤波器。之所以可以达到此去噪效果,是因为滤波器是由两个函数构成。一个函数是由几何空间距离决定滤波器系数。另一个由像素差值决定滤波器系数。可以与其相比较的两个filter:高斯低通滤波器(http://en.wikipedia.org/wiki/Gaussian_filter)和α-截尾均值滤波器(去掉百分率为α的最小值和最大之后剩下像素的均值作为滤波器),后文中将结合公式做详细介绍。


    双边滤波器中,输出像素的值依赖于邻域像素的值的加权组合,


    权重系数w(i,j,k,l)取决于定义域核


    和值域核

    的乘积


    同时考虑了空间域与值域的差别,而Gaussian Filter和α均值滤波分别只考虑了空间域和值域差别。


    =======================================================================

    双边滤波器的实现(MATLAB):function B = bfilter2(A,w,sigma)

    CopyRight:

    % Douglas R. Lanman, Brown University, September 2006.
    % dlanman@brown.edu, http://mesh.brown.edu/dlanman


    具体请见function B = bfltGray(A,w,sigma_d,sigma_r)函数说明。


    [cpp] view plain copy
    1. %简单地说:  
    2. %A为给定图像,归一化到[0,1]的矩阵  
    3. %W为双边滤波器(核)的边长/2  
    4. %定义域方差σd记为SIGMA(1),值域方差σr记为SIGMA(2)  
    5.   
    6. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
    7. % Pre-process input and select appropriate filter.  
    8. function B = bfilter2(A,w,sigma)  
    9.   
    10. % Verify that the input image exists and is valid.  
    11. if ~exist('A','var') || isempty(A)  
    12.    error('Input image A is undefined or invalid.');  
    13. end  
    14. if ~isfloat(A) || ~sum([1,3] == size(A,3)) || ...  
    15.       min(A(:)) < 0 || max(A(:)) > 1  
    16.    error(['Input image A must be a double precision ',...  
    17.           'matrix of size NxMx1 or NxMx3 on the closed ',...  
    18.           'interval [0,1].']);        
    19. end  
    20.   
    21. % Verify bilateral filter window size.  
    22. if ~exist('w','var') || isempty(w) || ...  
    23.       numel(w) ~= 1 || w < 1  
    24.    w = 5;  
    25. end  
    26. w = ceil(w);  
    27.   
    28. % Verify bilateral filter standard deviations.  
    29. if ~exist('sigma','var') || isempty(sigma) || ...  
    30.       numel(sigma) ~= 2 || sigma(1) <= 0 || sigma(2) <= 0  
    31.    sigma = [3 0.1];  
    32. end  
    33.   
    34. % Apply either grayscale or color bilateral filtering.  
    35. if size(A,3) == 1  
    36.    B = bfltGray(A,w,sigma(1),sigma(2));  
    37. else  
    38.    B = bfltColor(A,w,sigma(1),sigma(2));  
    39. end  
    40.   
    41.   
    42. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
    43. % Implements bilateral filtering for grayscale images.  
    44. function B = bfltGray(A,w,sigma_d,sigma_r)  
    45.   
    46. % Pre-compute Gaussian distance weights.  
    47. [X,Y] = meshgrid(-w:w,-w:w);  
    48. %创建核距离矩阵,e.g.  
    49. %  [x,y]=meshgrid(-1:1,-1:1)  
    50. %   
    51. % x =  
    52. %   
    53. %     -1     0     1  
    54. %     -1     0     1  
    55. %     -1     0     1  
    56. %   
    57. %   
    58. % y =  
    59. %   
    60. %     -1    -1    -1  
    61. %      0     0     0  
    62. %      1     1     1  
    63. %计算定义域核  
    64. G = exp(-(X.^2+Y.^2)/(2*sigma_d^2));  
    65.   
    66. % Create waitbar.  
    67. h = waitbar(0,'Applying bilateral filter...');  
    68. set(h,'Name','Bilateral Filter Progress');  
    69.   
    70. % Apply bilateral filter.  
    71. %计算值域核H 并与定义域核G 乘积得到双边权重函数F  
    72. dim = size(A);  
    73. B = zeros(dim);  
    74. for i = 1:dim(1)  
    75.    for j = 1:dim(2)  
    76.         
    77.          % Extract local region.  
    78.          iMin = max(i-w,1);  
    79.          iMax = min(i+w,dim(1));  
    80.          jMin = max(j-w,1);  
    81.          jMax = min(j+w,dim(2));  
    82.          %定义当前核所作用的区域为(iMin:iMax,jMin:jMax)  
    83.          I = A(iMin:iMax,jMin:jMax);%提取该区域的源图像值赋给I  
    84.         
    85.          % Compute Gaussian intensity weights.  
    86.          H = exp(-(I-A(i,j)).^2/(2*sigma_r^2));  
    87.         
    88.          % Calculate bilateral filter response.  
    89.          F = H.*G((iMin:iMax)-i+w+1,(jMin:jMax)-j+w+1);  
    90.          B(i,j) = sum(F(:).*I(:))/sum(F(:));  
    91.                  
    92.    end  
    93.    waitbar(i/dim(1));  
    94. end  
    95.   
    96. % Close waitbar.  
    97. close(h);  
    98.   
    99.   
    100. %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
    101. % Implements bilateral filter for color images.  
    102. function B = bfltColor(A,w,sigma_d,sigma_r)  
    103.   
    104. % Convert input sRGB image to CIELab color space.  
    105. if exist('applycform','file')  
    106.    A = applycform(A,makecform('srgb2lab'));  
    107. else  
    108.    A = colorspace('Lab<-RGB',A);  
    109. end  
    110.   
    111. % Pre-compute Gaussian domain weights.  
    112. [X,Y] = meshgrid(-w:w,-w:w);  
    113. G = exp(-(X.^2+Y.^2)/(2*sigma_d^2));  
    114.   
    115. % Rescale range variance (using maximum luminance).  
    116. sigma_r = 100*sigma_r;  
    117.   
    118. % Create waitbar.  
    119. h = waitbar(0,'Applying bilateral filter...');  
    120. set(h,'Name','Bilateral Filter Progress');  
    121.   
    122. % Apply bilateral filter.  
    123. dim = size(A);  
    124. B = zeros(dim);  
    125. for i = 1:dim(1)  
    126.    for j = 1:dim(2)  
    127.         
    128.          % Extract local region.  
    129.          iMin = max(i-w,1);  
    130.          iMax = min(i+w,dim(1));  
    131.          jMin = max(j-w,1);  
    132.          jMax = min(j+w,dim(2));  
    133.          I = A(iMin:iMax,jMin:jMax,:);  
    134.         
    135.          % Compute Gaussian range weights.  
    136.          dL = I(:,:,1)-A(i,j,1);  
    137.          da = I(:,:,2)-A(i,j,2);  
    138.          db = I(:,:,3)-A(i,j,3);  
    139.          H = exp(-(dL.^2+da.^2+db.^2)/(2*sigma_r^2));  
    140.         
    141.          % Calculate bilateral filter response.  
    142.          F = H.*G((iMin:iMax)-i+w+1,(jMin:jMax)-j+w+1);  
    143.          norm_F = sum(F(:));  
    144.          B(i,j,1) = sum(sum(F.*I(:,:,1)))/norm_F;  
    145.          B(i,j,2) = sum(sum(F.*I(:,:,2)))/norm_F;  
    146.          B(i,j,3) = sum(sum(F.*I(:,:,3)))/norm_F;  
    147.                   
    148.    end  
    149.    waitbar(i/dim(1));  
    150. end  
    151.   
    152. % Convert filtered image back to sRGB color space.  
    153. if exist('applycform','file')  
    154.    B = applycform(B,makecform('lab2srgb'));  
    155. else    
    156.    B = colorspace('RGB<-Lab',B);  
    157. end  
    158.   
    159. % Close waitbar.  
    160. close(h);  


    调用方法:

    [cpp] view plain copy
    1. I=imread('einstein.jpg');  
    2. I=double(I)/255;  
    3.   
    4. w     = 5;       % bilateral filter half-width  
    5. sigma = [3 0.1]; % bilateral filter standard deviations  
    6.   
    7. I1=bfilter2(I,w,sigma);  
    8.   
    9. subplot(1,2,1);  
    10. imshow(I);  
    11. subplot(1,2,2);  
    12. imshow(I1)  

    实验结果:


    参考资料:

    1.《Computer Vision Algorithms and Applications》

    2. http://de.wikipedia.org/wiki/Bilaterale_Filterung

    3.http://www.cs.duke.edu/~tomasi/papers/tomasi/tomasiIccv98.pdf

    4. http://homepages.inf.ed.ac.uk/rbf/CVonline/LOCAL_COPIES/MANDUCHI1/Bilateral_Filtering.html

    5. http://mesh.brown.edu/dlanman


    展开全文
  • 双边滤波算法原理

    万次阅读 多人点赞 2016-03-10 22:15:40
     双边滤波在图像处理领域中有着广泛的应用,比如去噪、去马赛克、光流估计等等,最近,比较流行的Non-Local算法也可以看成是双边滤波的一种扩展。自从Tomasi et al等人提出该算法那一天起,如何快速的实现他,一直...

    1. 简介

    图像平滑是一个重要的操作,而且有多种成熟的算法。这里主要简单介绍一下Bilateral方法(双边滤波),这主要是由于前段时间做了SSAO,需要用bilateral blur 算法进行降噪。Bilateral blur相对于传统的高斯blur来说很重要的一个特性即可可以保持边缘(Edge Perseving),这个特点对于一些图像模糊来说很有用。一般的高斯模糊在进行采样时主要考虑了像素间的空间距离关系,但是却并没有考虑像素值之间的相似程度,因此这样我们得到的模糊结果通常是整张图片一团模糊。Bilateral blur的改进就在于在采样时不仅考虑像素在空间距离上的关系,同时加入了像素间的相似程度考虑,因而可以保持原始图像的大体分块进而保持边缘。在于游戏引擎的post blur算法中,bilateral blur常常被用到,比如对SSAO的降噪。

    2. 原理

    滤波算法中,目标点上的像素值通常是由其所在位置上的周围的一个小局部邻居像素的值所决定。在2D高斯滤波中的具体实现就是对周围的一定范围内的像素值分别赋以不同的高斯权重值,并在加权平均后得到当前点的最终结果。而这里的高斯权重因子是利用两个像素之间的空间距离(在图像中为2D)关系来生成。通过高斯分布的曲线可以发现,离目标像素越近的点对最终结果的贡献越大,反之则越小。其公式化的描述一般如下所述:

    其中的c即为基于空间距离的高斯权重,而用来对结果进行单位化。

    高斯滤波在低通滤波算法中有不错的表现,但是其却有另外一个问题,那就是只考虑了像素间的空间位置上的关系,因此滤波的结果会丢失边缘的信息。这里的边缘主要是指图像中主要的不同颜色区域(比如蓝色的天空,黑色的头发等),而Bilateral就是在Gaussian blur中加入了另外的一个权重分部来解决这一问题。Bilateral滤波中对于边缘的保持通过下述表达式来实现:

    其中的s为基于像素间相似程度的高斯权重,同样用来对结果进行单位化。对两者进行结合即可以得到基于空间距离、相似程度综合考量的Bilateral滤波:

    上式中的单位化分部综合了两种高斯权重于一起而得到,其中的cs计算可以详细描述如下:

    且有

    且有

    上述给出的表达式均是在空间上的无限积分,而在像素化的图像中当然无法这么做,而且也没必要如此做,因而在使用前需要对其进行离散化。而且也不需要对于每个局部像素从整张图像上进行加权操作,距离超过一定程度的像素实际上对当前的目标像素影响很小,可以忽略的。限定局部子区域后的离散化公就可以简化为如下形式:

    上述理论公式就构成了Bilateral滤波实现的基础。为了直观地了解高斯滤波与双边滤波的区别,我们可以从下列图示中看出依据。假设目标源图像为下述左右区域分明的带有噪声的图像(由程序自动生成),蓝色框的中心即为目标像素所在的位置,那么当前像素处所对应的高斯权重与双边权重因子3D可视化后的形状如后边两图所示:            

    左图为原始的噪声图像;中间为高斯采样的权重;右图为Bilateral采样的权重。从图中可以看出Bilateral加入了相似程度分部以后可以将源图像左侧那些跟当前像素差值过大的点给滤去,这样就很好地保持了边缘。为了更加形象地观察两者间的区别,使用Matlab将该图在两种不同方式下的高度图3D绘制出来,如下:

      

    上述三图从左到右依次为:双边滤波,原始图像,高斯滤波。从高度图中可以明显看出Bilateral和Gaussian两种方法的区别,前者较好地保持了边缘处的梯度,而在高斯滤波中,由于其在边缘处的变化是线性的,因而就使用连累的梯度呈现出渐变的状态,而这表现在图像中的话就是边界的丢失(图像的示例可见于后述)。                                              

    3. 代码实现

    有了上述理论以后实现Bilateral Filter就比较简单了,其实它也与普通的Gaussian Blur没有太大的区别。这里主要包括3部分的操作: 基于空间距离的权重因子生成;基于相似度的权重因子的生成;最终filter颜色的计算。

    3.1 Spatial Weight

    这就是通常的Gaussian Blur中使用的计算高斯权重的方法,其主要通过两个pixel之间的距离并使用如下公式计算而来:

    其中的就表示两个像素间的距离,比如当前像素与其右边紧邻的一个像素之间的距离我们就可以用来计算,也即两个二维向量{0 , 0}以及{0 , 1}之间的欧氏距离。直接计算一个区域上的高斯权重并单位化后就可以进行高斯模糊了。

    3.2 Similarity Weight

    与基于距离的高斯权重计算类似,只不过此处不再根据两个pixel之间的空间距离,而是根据其相似程度(或者两个pixel的值之间的距离)。

    其中的表示两个像素值之间的距离,可以直接使用其灰度值之间的差值或者RGB向量之间的欧氏距离。

    3.3 Color Filtering

    有了上述两部分所必需的权重因子之后,那么具体的双边滤波的实现即与普通的高斯滤波无异。主要部分代码如下述:

    [cpp] view plain copy
    1. UCHAR3 BBColor(int posX , int posY)  
    2. {  
    3.     int centerItemIndex = posY * picWidth4 + posX * 3 , neighbourItemIndex;  
    4.     int weightIndex;  
    5.     double gsAccumWeight = 0;  
    6.     double accumColor = 0;  
    7.   
    8.     // 计算各个采样点处的Gaussian权重,包括closeness,similarity  
    9.     for(int i = -number ; i <= number ; ++i)  
    10.     {  
    11.         for(int j = -number ; j <= number ; ++j)  
    12.         {  
    13.             weightIndex = (i + number) * (number * 2 + 1) + (j + number);  
    14.             neighbourItemIndex = min(noiseImageHeight - 1 , max(0 , posY + j * radius)) * picWidth4 +  
    15.                              min(noiseImageWidth - 1  , max(0 , posX + i * radius)) * 3;  
    16.               
    17.             pCSWeight[weightIndex] = LookupGSWeightTable(pSrcDataBuffer[neighbourItemIndex] , pSrcDataBuffer[centerItemIndex]);  
    18.             pCSWeight[weightIndex] = pGSWeight[weightIndex] * pGCWeight[weightIndex];  
    19.             gsAccumWeight += pCSWeight[weightIndex];  
    20.         }  
    21.     }  
    22.       
    23.     // 单位化权重因子  
    24.     gsAccumWeight = 1 / gsAccumWeight;  
    25.     for(int i = -number ; i <= number ; ++i)  
    26.     {  
    27.         for(int j = -number ; j <= number ; ++j)  
    28.         {  
    29.             weightIndex = (i + number) * (number * 2 + 1) + (j + number);  
    30.             pCSWeight[weightIndex] *= gsAccumWeight;  
    31.         }  
    32.     }  
    33.       
    34.     // 计算最终的颜色并返回  
    35.     for(int i = -number ; i <= number ; ++i)  
    36.     {  
    37.         for(int j = -number ; j <= number ; ++j)  
    38.         {  
    39.             weightIndex = (i + number) * (number * 2 + 1) + (j + number);  
    40.             neighbourItemIndex = min(noiseImageHeight - 1 , max(0 , posY + j * radius)) * picWidth4 +  
    41.                                  min(noiseImageWidth - 1  , max(0 , posX + i * radius)) * 3;  
    42.             accumColor += pSrcDataBuffer[neighbourItemIndex + 0] * pCSWeight[weightIndex];  
    43.         }  
    44.     }  
    45.   
    46.     return UCHAR3(accumColor , accumColor , accumColor);  
    47. }  

     其中的相似度分部的权重s主要根据两个Pixel之间的颜色差值计算面来。对于灰度图而言,这个差值的范围是可以预知的,即[-255, 255],因而为了提高计算的效率我们可以将该部分权重因子预计算生成并存表,在使用时快速查询即可。使用上述实现的算法对几张带有噪声的图像进行滤波后的结果如下所示:

       

        

    上图从左到右分别为:双边滤波;原始图像;高斯滤波。从图片中可以较为明显地看出两种算法的区别,最直观的感受差别就是使用高斯算法后整张图片都是一团模糊的状态;而双边滤波则可以较好地保持原始图像中的区域信息,看起来仍然嘴是嘴、眼是眼(特别是在第一张美女图像上的效果!看来PS是灰常重要啊~~^o^)。

    4. 在SSAO中的使用

    在上述实现中的边缘判定函数主要是通过两个像素值之间的差异来决定,这也是我们观察普通图片的一种普遍感知方式。当然,也可以根据使用的需求情况来使用其它的方式判断其它定义下的边缘,比如使用场景的depth或是normal。比如在对SSAO进行滤波时可以直接使用Depth值来行边缘判断。首先,设置一个深度的阈值,在作边缘检测时比较两点间的depth差值,如果差值大于阈值,则认为其属于不同的区域,则此处就应为边界。使用此方法得到的效果可见于下图所示:

    高斯滤波

    双边滤波 

    引用:

    【1】http://blog.csdn.net/bugrunner/article/details/7170471

    【2】http://blog.csdn.net/bugrunner/article/details/7170471

    展开全文
  • 双边滤波——原理及matlab实现

    万次阅读 2018-08-18 21:05:42
    1、双边滤波简介:  双边滤波(Bilateral filter)是一种非线性滤波方法(空间权值+相似权值)——空间权值:模糊去噪;相似权值:保护边缘。 2、双边滤波原理  双边滤波具有两个权重:空间权重与相似权重  ...

       思维闭塞时可外出采采风。

    1、双边滤波简介:

        双边滤波(Bilateral filter)是一种非线性滤波方法(空间权值+相似权值)——空间权值:模糊去噪;相似权值:保护边缘。


    2、双边滤波原理

        双边滤波具有两个权重:空间权重与相似权重

        1)空间权重:与像素位置有关,为像素之间的距离(欧式距离,空间度量),故可定义为全局变量放在循环外,通常定义为

    c(\xi ,x)=e^{-\tfrac{1}{2}(\tfrac{d(\xi ,x)}{\sigma_{d}})^{2}}

    d(\xi,x)=d(\xi-x)=||\xi -x||

        其中d(\xi,x)表示两个像素间的距离(欧式距离)。该过程滤波如下:

    h(x)=k_{d}^{-1}(x)\int_{-\infty }^{\infty }\int_{-\infty }^{\infty }f(\xi )c(\xi ,x)d\xi

        权值为:

    k_{d}(x)=\int_{-\infty }^{\infty }\int_{-\infty }^{\infty }c(\xi ,x)d\xi

      2)相似权重:与像素值大小有关,为像素值之间的距离(辐射距离,相似性度量),根据像素值不同而不同,需要放在循环内,通常定义为

    s(f(\xi) ,f(x))=e^{-\frac{1}{2}(\frac{\sigma (f(\xi ),f(x))}{\sigma _{\tau }})^{2}}

    \sigma(\xi ,x)=\sigma(\xi-x)=||\xi -x||

    其中\sigma (f(\xi ),f(x))表示两个像素值之间的距离。该过程滤波如下:

    h(x)=k_{\tau }^{-1}(x)\int_{-\infty }^{\infty }\int_{-\infty }^{\infty }f(\xi )s(f(\xi) ,f(x))d\xi

        权值为:

    k_{\tau }(x)=\int_{-\infty }^{\infty }\int_{-\infty }^{\infty }s(f(\xi) ,f(x))d\xi

        3)两者结合,得到基于空间距离、相似程度整体考虑的双边滤波如下:

    h(x)=k^{-1}(x)\int_{-\infty }^{\infty }\int_{-\infty }^{\infty }f(\xi )c(\xi ,x)s(f(\xi) ,f(x))d\xi

        权值为:

    k(x)=\int_{-\infty }^{\infty }\int_{-\infty }^{\infty }c(\xi,x)s(f(\xi) ,f(x))d\xi


    3、双边滤波实现:

        实际应用时,根据需要对积分采用离散形式表示,滤波半径根据需要自行设置。

        在进行滤波前,需将数据转换为浮点型等。


    4、matlab源代码

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%主函数%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    % 读取
    f = imread(filename);
    
    % 设置参数
    r = 5;% 滤波半径
    a = 3;% 全局方差
    b = 0.1;% 局部方差
    
    % 判断二维图还是三维图
    if ismatrix(f)
        g = bfilt_gray(f,r,a,b);
    else
        g = bfilt_rgb(f,r,a,b);
    end
    
    % 显示
    subplot(121)
    imshow(f)
    subplot(122)
    imshow(g)
    
    
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%灰度图双边滤波%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    function g = bfilt_gray(f,r,a,b)
    % f灰度图;r滤波半径;a全局方差;b局部方差
    [x,y] = meshgrid(-r:r);
    w1 = exp(-(x.^2+y.^2)/(2*a^2));
    f = tofloat(f);%f = im2double(f);
    
    h = waitbar(0,'Applying bilateral filter...');
    set(h,'Name','Bilateral Filter Progress');
    
    [m,n] = size(f);
    f_temp = padarray(f,[r r],'symmetric');
    g = zeros(m,n);
    for i = r+1:m+r
        for j = r+1:n+r
            temp = f_temp(i-r:i+r,j-r:j+r);
            w2 = exp(-(temp-f(i-r,j-r)).^2/(2*b^2));
            w = w1.*w2;
            s = temp.*w;
            g(i-r,j-r) = sum(s(:))/sum(w(:));
        end
        waitbar((i-r)/m);
    end
    % g = revertclass(g);
    
    
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%彩色图双边滤波%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    function g = bfilt_rgb(f,r,a,b)
    % f灰度图;r滤波半径;a全局方差;b局部方差
    [x,y] = meshgrid(-r:r);
    w1 = exp(-(x.^2+y.^2)/(2*a^2));
    f = tofloat(f);%f = im2double(f);
    
    h = waitbar(0,'Applying bilateral filter...');
    set(h,'Name','Bilateral Filter Progress');
    
    fr = f(:,:,1);
    fg = f(:,:,2);
    fb = f(:,:,3);
    [m,n] = size(fr);
    fr_temp = padarray(fr,[r r],'symmetric');
    fg_temp = padarray(fg,[r r],'symmetric');
    fb_temp = padarray(fb,[r r],'symmetric');
    [gr,gg,gb] = deal(zeros(size(fr)));
    
    
    for i = r+1:m+r
        for j = r+1:n+r
            temp1 = fr_temp(i-r:i+r,j-r:j+r);
            temp2 = fg_temp(i-r:i+r,j-r:j+r);
            temp3 = fb_temp(i-r:i+r,j-r:j+r);
            dr = temp1 - fr_temp(i,j);
            dg = temp2 - fg_temp(i,j);
            db = temp3 - fb_temp(i,j);
            w2 = exp(-(dr.^2+dg.^2+db.^2)/(2*b^2));
            w = w1.*w2;
            gr(i-r,j-r) = sum(sum(temp1.*w))/sum(w(:));
            gg(i-r,j-r) = sum(sum(temp2.*w))/sum(w(:));
            gb(i-r,j-r) = sum(sum(temp3.*w))/sum(w(:));
        end
        waitbar((i-r)/n);
    end
    g = cat(3,gr,gg,gb);
    % g = revertclass(g);
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%可以用到的函数%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    function [out,revertclass]=tofloat(in)
    identity=@(x) x;
    tosingle=@im2single;
    
    table={'uint8',tosingle,@im2uint8
    'uint16',tosingle,@im2uint16
    'int16',tosingle,@im2int16
    'logical',tosingle,@logical
    'double',identity,identity
    'single',identity,identity};
    
    classIndex=find(strcmp(class(in),table(:,1)));
    
    if isempty(classIndex)
        error('unsupported input immage class.');
    end
    
    out=table{classIndex,2}(in);
    revertclass=table{classIndex,3}; 
    
    
    

     5、实验结果

     

    6、参考资料

    https://blog.csdn.net/abcjennifer/article/details/7616663

    https://blog.csdn.net/bugrunner/article/details/7170471

    https://www.cnblogs.com/qiqibaby/p/5296681.html

    展开全文
  • Bilateral Filters(双边滤波算法)原理及实现

    万次阅读 多人点赞 2017-10-21 23:13:26
    双边滤波算法原理: 双边滤波算法实现: 双边滤波算法实例: 参考: http://people.csail.mit.edu/sparis/bf/ http://blog.csdn.net/fightingforcv/article/details/52723376 ...
  • Bilateral Filters(双边滤波算法)原理及实现(一)

    万次阅读 多人点赞 2019-02-25 16:27:17
    双边滤波算法原理 双边滤波是一种非线性滤波器,它可以达到保持边缘、降噪平滑的效果。和其他滤波原理一样,双边滤波也是采用加权平均的方法,用周边像素亮度值的加权平均代表某个像素的强度,所用的加权平均基于...
  • OpenCV双边滤波详解及实代码实现

    万次阅读 多人点赞 2018-10-07 21:30:55
    双边滤波(Bilateral Filter)是非线性滤波中的一种。这是一种结合图像的空间邻近度与像素值相似度的处理办法。在滤波时,该滤波方法同时考虑空间临近信息与颜色相似信息,在滤除噪声、平滑图像的同时,又做到边缘...
  • 双边滤波与引导滤波

    万次阅读 多人点赞 2014-05-22 15:41:20
    双边滤波与引导滤波 分类: AI and Computer Vision2014-03-07 17:04344人阅读评论(0)收藏举报 图像处理滤波 双边滤波 双边滤波很有名,使用广泛,简单的说就是一种同时考虑了像素空间差异与强度差异的...
  • 【图像处理】——双边滤波

    千次阅读 热门讨论 2019-09-05 10:39:52
    双边滤波    由于高斯噪声在信号采集系统中往往是无法避免的,所以在信号处理中高斯滤波器是最常用的滤波器之一;数字图像也是一种信号,高斯滤波也是最常用的图像去噪方法。但是在之前的博客中提到过,高斯滤波是...
  • http://blog.csdn.net/pi9nc/article/details/26592377
  • python 双边滤波与高斯滤波

    千次阅读 2018-04-07 17:58:17
    python 双边滤波与高斯滤波高斯滤波就是对整幅图像进行加权平均的过程。每个像素点的值,都由其本身和邻域内的其它像素值经过加权平均后得到。高斯滤波的详细操作是:用一个模板(或称卷积、掩模)扫描图像中的每个...
  • C/C++ OpenCV中值滤波&双边滤波

    千次阅读 2017-01-09 21:48:21
    C/C++ OpenCV中值滤波&双边滤波
  • 通过对局部线性关系求导可以判断哪些边缘需要保留,给出的a,b两个参数则是以像素k为窗口的周围的权重的均值2、引导滤波的优点引导滤波相对于双边滤波最大的优点在于算法的复杂度与窗口的大小无关,对于处理较为大型...
  • 图像处理-双边滤波 Bilateral Filtering

    千次阅读 2018-03-23 11:01:42
    在Shader中实现双边滤波的时候,总感觉理解的不太透彻,这里写博客记录一下。 参考资料: 算法原理 GPUImage中Bilateral Filtering的实现 高斯滤波 空间域 在理解双边滤波之前,先来理解上面是高斯滤波。 ...
  • 一维信号分别经过高斯滤波和双边滤波后,得到的曲线和原曲线之间的差异,图片如![图片说明](https://img-ask.csdn.net/upload/201512/17/1450320427_527179.jpg),下面蓝色的是滤波后与滤波前的差值,matlab程序如何...
  • python&opencv 图像的双边滤波

    千次阅读 2018-05-28 18:43:15
    双边滤波的操作主要是ccv2.bilateralFilter()函数来操作,它能够保持边界清晰的情况下有效的去除噪声,但是这种操作比较慢。它拥有着美颜的效果: 下面是代码演示: import cv2 def bi_demo(image):#高斯双边滤波 ...
  • 滤波算法主要包括均值滤波,高斯滤波,中值滤波和双边滤波。 每种算法都有自己的特点,建议从原理上了解每种算法的优缺点。上图给出简洁版的总结。 以下是代码: import numpy as np import cv2 import ...
  • opencv 中的双边滤波用法总结(10)

    千次阅读 2018-08-14 09:19:33
    双边滤波】结合空间临近度和像素值相似度的一种折中处理 原型:void bilateralFilter( InputArray src,  OutputArray dst, int d, double sigmaColor, double sigmaSpace, int borderType=BORDER_DEFAULT ); ...
  • 双边滤波参数选择

    千次阅读 2018-01-09 15:39:50
    双边滤波部分总结 参数ksizesigmaCsigmaS Ksize的选择:双边滤波的作用式图像去噪,因此该窗口的大小要盖住要除去的噪声大小,根据图中噪声的大小,设置双边的窗口大小。窗口要是过小的话,去噪效果不明显;窗口过...
  • 双边滤波原理(Bilateral Filtering)

    万次阅读 2017-12-18 22:49:15
    双边滤波原理(Bilateral Filtering)
1 2 3 4 5 ... 20
收藏数 5,048
精华内容 2,019
关键字:

双边滤波