OpenCV鄰域濾波方框高斯中值雙邊濾波

1、鄰域濾波（卷積）鄰域算子值利用給定像素周圍像素的值決定此像素的最終輸出。如圖左邊圖像與中間圖像卷積禪城右邊圖像。目標圖像中綠色的像素由原圖像中藍色標記的像素計算得到。通用線性鄰域濾波是一種常用的鄰域算子，輸入像素加權(quán)得到輸出像素：其中權(quán)重核   為“濾波系數(shù)”。上面的式子可以簡記為：【方框濾波】最簡單的線性濾波是移動平均或方框濾波，用 窗口中的像素值平均后輸出，核函數(shù)為：其實等價于圖像與全部元素值為1的核函數(shù)進行卷積再進行尺度縮放。代碼OpenCV中的 blur函數(shù)是進行標準方框濾波：cpp view plain copy1. voi

2、d cv:blur( InputArray src, OutputArray dst,  2.            Size ksize, Point anchor, int borderType )  3.   4.     boxFilter( src, 

3、;dst, -1, ksize, anchor, true, borderType );  5.   而boxFilter函數(shù)源碼如下：cpp view plain copy1. cv:Ptr<cv:FilterEngine> cv:createBoxFilter( int srcType, int dstType, Size ksize,  2.   &

4、#160;                 Point anchor, bool normalize, int borderType )  3.   4.     int sdepth = CV_MAT_DEPTH(srcType);  5

5、.     int cn = CV_MAT_CN(srcType), sumType = CV_64F;  6.     if( sdepth <= CV_32S && (!normalize |  7.         ksize.width*ksize.

6、height <= (sdepth = CV_8U ? (1<<23) :  8.             sdepth = CV_16U ? (1 << 15) : (1 << 16) )  9.  

7、60;      sumType = CV_32S;  10.     sumType = CV_MAKETYPE( sumType, cn );  11.   12.     Ptr<BaseRowFilter> rowFilter = getRowSumFilter(srcType,

8、 sumType, ksize.width, anchor.x );  13.     Ptr<BaseColumnFilter> columnFilter = getColumnSumFilter(sumType,  14.         dstType, ksize.height, anchor.y, normalize

9、 ? 1./(ksize.width*ksize.height) : 1);  15.   16.     return Ptr<FilterEngine>(new FilterEngine(Ptr<BaseFilter>(0), rowFilter, columnFilter,  17.         

10、0;  srcType, dstType, sumType, borderType );  18.   這里 blur 和 boxFilter 的區(qū)別是，blur是標準化后的 boxFilter，即boxFilter的核函數(shù)：其中，cpp view plain copy1. blur( src, dst, Size( 1, 1 ), Point(-1,-1);  2. blur(

11、0;src, dst, Size( 4, 4 ), Point(-1,-1);  3. blur( src, dst, Size( 8, 8 ), Point(-1,-1);  4. blur( src, dst, Size( 16, 16 ), Point(-1,-1);  實驗結(jié)果下圖是對一幅圖像分別用1*1,4*4,8*8,16*1

12、6標準方框濾波后的圖像：      【高斯濾波】高斯濾波器是一類根據(jù)高斯函數(shù)的形狀來選擇權(quán)值的線性平滑濾波器。它對去除服從正態(tài)分布的噪聲很有效。常用的零均值離散高斯濾波器函數(shù)：2D圖像中表示為：代碼cpp view plain copy1. /* 2.                      &#

13、160;               Gaussian Blur 3. */  4.   5. cv:Mat cv:getGaussianKernel( int n, double sigma, int ktype )  6.   7.   

14、0; const int SMALL_GAUSSIAN_SIZE = 7;  8.     static const float small_gaussian_tabSMALL_GAUSSIAN_SIZE =  9.       10.         1.f,  11

15、.         0.25f, 0.5f, 0.25f,  12.         0.0625f, 0.25f, 0.375f, 0.25f, 0.0625f,  13.         0.03125f, 0.109375f, 0.2

16、1875f, 0.28125f, 0.21875f, 0.109375f, 0.03125f  14.       15.   16.     const float* fixed_kernel = n % 2 = 1 && n <= SMALL_GAUSSIAN_SIZE&#

17、160;&& sigma <= 0 ?  17.         small_gaussian_tabn>>1 : 0;  18.   19.     CV_Assert( ktype = CV_32F | ktype = CV_64F );&#

18、160; 20.     Mat kernel(n, 1, ktype);  21.     float* cf = (float*)kernel.data;  22.     double* cd = (double*)kernel.data;  23.   24.   

19、60; double sigmaX = sigma > 0 ? sigma : (n-1)*0.5 - 1)*0.3 + 0.8;  25.     double scale2X = -0.5/(sigmaX*sigmaX);  26.     double sum = 0

20、;  27.   28.     int i;  29.     for( i = 0; i < n; i+ )  30.       31.         double x =&#

21、160;i - (n-1)*0.5;  32.         double t = fixed_kernel ? (double)fixed_kerneli : std:exp(scale2X*x*x);  33.         if( ktype = CV_32F 

22、;)  34.           35.             cfi = (float)t;  36.             sum += cfi;  

23、37.           38.         else  39.           40.             cdi = t; 

24、0;41.             sum += cdi;  42.           43.       44.   45.     sum = 1./sum;  46. &#

25、160;   for( i = 0; i < n; i+ )  47.       48.         if( ktype = CV_32F )  49.         

26、0;   cfi = (float)(cfi*sum);  50.         else  51.             cdi *= sum;  52.       53.   54

27、.     return kernel;  55.   56.   57.   58. cv:Ptr<cv:FilterEngine> cv:createGaussianFilter( int type, Size ksize,  59.             &#

28、160;                           double sigma1, double sigma2,  60.             &#

29、160;                           int borderType )  61.   62.     int depth = CV_MAT_DEPTH(type); 

30、 63.     if( sigma2 <= 0 )  64.         sigma2 = sigma1;  65.   66.     / automatic detection of kernel size from sig

31、ma  67.     if( ksize.width <= 0 && sigma1 > 0 )  68.         ksize.width = cvRound(sigma1*(depth = CV_8U ? 3 : 4)*2 +

32、60;1)|1;  69.     if( ksize.height <= 0 && sigma2 > 0 )  70.         ksize.height = cvRound(sigma2*(depth = CV_8U ? 3 : 4)*2&#

33、160;+ 1)|1;  71.   72.     CV_Assert( ksize.width > 0 && ksize.width % 2 = 1 &&  73.         ksize.height > 0 &&am

34、p; ksize.height % 2 = 1 );  74.   75.     sigma1 = std:max( sigma1, 0. );  76.     sigma2 = std:max( sigma2, 0. );  77.   78. 

35、0;   Mat kx = getGaussianKernel( ksize.width, sigma1, std:max(depth, CV_32F) );  79.     Mat ky;  80.     if( ksize.height = ksize.width && std:abs

36、(sigma1 - sigma2) < DBL_EPSILON )  81.         ky = kx;  82.     else  83.         ky = getGaussianKernel( ksize.heig

37、ht, sigma2, std:max(depth, CV_32F) );  84.   85.     return createSeparableLinearFilter( type, type, kx, ky, Point(-1,-1), 0, borderType );  86.   87.   88.  &#

38、160;89. void cv:GaussianBlur( InputArray _src, OutputArray _dst, Size ksize,  90.                    double sigma1, double sigma2,  91.

39、                    int borderType )  92.   93.     Mat src = _src.getMat();  94.     _dst.create( src.si

40、ze(), src.type() );  95.     Mat dst = _dst.getMat();  96.   97.     if( borderType != BORDER_CONSTANT )  98.       99.     &

41、#160;   if( src.rows = 1 )  100.             ksize.height = 1;  101.         if( src.cols = 1 )  102.

42、60;           ksize.width = 1;  103.       104.   105.     if( ksize.width = 1 && ksize.height = 1 )  106.

43、      107.         src.copyTo(dst);  108.         return;  109.       110.   111. #ifdef HAVE_TEGRA_OPTIMIZATION  1

44、12.     if(sigma1 = 0 && sigma2 = 0 && tegra:gaussian(src, dst, ksize, borderType)  113.         return;  114. #endif  115.   116

45、.     Ptr<FilterEngine> f = createGaussianFilter( src.type(), ksize, sigma1, sigma2, borderType );  117.     f->apply( src, dst );  118.   實驗結(jié)果下圖是對一幅圖像分別用1*1,3*3

46、,5*5,9*9標準方框濾波后的圖像：      非線性濾波線性濾波易于構(gòu)造，且易于從頻率響應的角度分析，但如果噪聲是散粒噪聲而非高斯噪聲時線性濾波不能去除噪聲。如圖像突然出現(xiàn)很大的值，線性濾波只是轉(zhuǎn)換為柔和但仍可見的散粒。這時需要非線性濾波。簡單的非線性濾波有 中值濾波， -截尾均值濾波，定義域濾波 和值域濾波 。中值濾波選擇每個鄰域像素的中值輸出； -截尾均值濾波是指去掉百分率為 的最小值和最大值；定義域濾波中沿著邊界的數(shù)字是像素的距離；值域就是去掉值域外的像素值。中值濾

47、波代碼cpp view plain copy1. medianBlur ( src, dst, i );  中值濾波實驗下圖是對一幅圖像分別用3*3,5*5,7*7,9*9（這里必須是奇數(shù)）標準方框濾波后的圖像：      【雙邊濾波】雙邊濾波的思想是抑制與中心像素值差別太大的像素，輸出像素值依賴于鄰域像素值的加權(quán)合：權(quán)重系數(shù) 取決于定義域核 和依賴于數(shù)據(jù)的值域核 的乘積。相乘后會產(chǎn)生依賴于數(shù)據(jù)的雙邊權(quán)重函數(shù)：雙邊濾波源碼cp

48、p view plain copy1. /* 2.                                    Bilateral Filtering 3. */  4.

49、   5. namespace cv  6.   7.   8. static void  9. bilateralFilter_8u( const Mat& src, Mat& dst, int d,  10.              &#

50、160;      double sigma_color, double sigma_space,  11.                     int borderType )  12.   13.   

51、  int cn = src.channels();  14.     int i, j, k, maxk, radius;  15.     Size size = src.size();  16.   17.     CV_Assert( (src.typ

52、e() = CV_8UC1 | src.type() = CV_8UC3) &&  18.         src.type() = dst.type() && src.size() = dst.size() &&  19.      

53、60;  src.data != dst.data );  20.   21.     if( sigma_color <= 0 )  22.         sigma_color = 1;  23.     if( sigma_s

54、pace <= 0 )  24.         sigma_space = 1;  25.   26.     double gauss_color_coeff = -0.5/(sigma_color*sigma_color);  27.     double 

55、gauss_space_coeff = -0.5/(sigma_space*sigma_space);  28.   29.     if( d <= 0 )  30.         radius = cvRound(sigma_space*1.5);  31.     

56、else  32.         radius = d/2;  33.     radius = MAX(radius, 1);  34.     d = radius*2 + 1;  35.   36.    

57、60;Mat temp;  37.     copyMakeBorder( src, temp, radius, radius, radius, radius, borderType );  38.   39.     vector<float> _color_weight(cn*256);  40.   

58、  vector<float> _space_weight(d*d);  41.     vector<int> _space_ofs(d*d);  42.     float* color_weight = &_color_weight0;  43.     float* space_weight

59、0;= &_space_weight0;  44.     int* space_ofs = &_space_ofs0;  45.   46.     / initialize color-related bilateral filter coefficients  47.     for(&

60、#160;i = 0; i < 256*cn; i+ )  48.         color_weighti = (float)std:exp(i*i*gauss_color_coeff);  49.   50.     / initialize space-related bilatera

61、l filter coefficients  51.     for( i = -radius, maxk = 0; i <= radius; i+ )  52.         for( j = -radius; j <= radius

62、; j+ )  53.           54.             double r = std:sqrt(double)i*i + (double)j*j);  55.        

63、60;    if( r > radius )  56.                 continue;  57.             space_weightmaxk =&#

64、160;(float)std:exp(r*r*gauss_space_coeff);  58.             space_ofsmaxk+ = (int)(i*temp.step + j*cn);  59.           60.   61.   

65、;  for( i = 0; i < size.height; i+ )  62.       63.         const uchar* sptr = temp.data + (i+radius)*temp.step + radius*cn;&#

66、160; 64.         uchar* dptr = dst.data + i*dst.step;  65.   66.         if( cn = 1 )  67.         &

67、#160; 68.             for( j = 0; j < size.width; j+ )  69.               70.      

68、0;          float sum = 0, wsum = 0;  71.                 int val0 = sptrj;  72.     

69、            for( k = 0; k < maxk; k+ )  73.                   74.      &

70、#160;              int val = sptrj + space_ofsk;  75.                     float w =&#

71、160;space_weightk*color_weightstd:abs(val - val0);  76.                     sum += val*w;  77.            

72、;         wsum += w;  78.                   79.                 /

73、 overflow is not possible here => there is no need to use CV_CAST_8U  80.                 dptrj = (uchar)cvRound(sum/wsum); 

74、60;81.               82.           83.         else  84.           85.  

75、60;          assert( cn = 3 );  86.             for( j = 0; j < size.width*3; j += 3 )  87.  

76、;             88.                 float sum_b = 0, sum_g = 0, sum_r = 0, wsum = 0;  8

77、9.                 int b0 = sptrj, g0 = sptrj+1, r0 = sptrj+2;  90.                 f

78、or( k = 0; k < maxk; k+ )  91.                   92.                  &#

79、160;  const uchar* sptr_k = sptr + j + space_ofsk;  93.                     int b = sptr_k0, g = sptr_k1, r = sptr_k2;  94.    &