【OpenCV】邻域滤波：方框、高斯、中值、双边滤波

1、邻域滤波（卷积） 邻域算子值利用给定像素周围像素的值决定此像素的最终输出。如图左边图像与中间图像 卷积禅城右边图像。目标图像中绿色的像素由原图像中蓝色标记的像素计算得到。 通用线性邻域滤波是一种常用的邻域算子，输入像素加权得到输出像素： 斡J） = F f （F+砒+级化。 其中权重核 城V）为“滤波系数”。上面的式子可以简记为： g = 【方框滤波】 最简单的线性滤波是移动平均或方框滤波，用KkK窗口中的像素值平均后输出，核函 数为： 1 L 1 - 1 h-1 1 L 1 1 L hsize.yvith * hisze.height n 1 L 1 . 1 L1 其实等价于图像与全部元素值

2、为 1的核函数进行卷积再进行尺度缩放。 代码 OpenCV中的blur函数是进行标准方框滤波: cpp view pla in copy 1. void cv:blur( InputArray src, OutputArray dst, 2. Size ksize, Point anchor, int borderType ) 3. 4. boxFilter( src, dst, -1, ksize, anchor, true , borderType ); 5. 而boxFilter函数源码如下: cpp view pla in copy 1. cv:Ptr cv:createBoxFilt

3、er(int srcType, int dstType, Size ksize, 2. Point anchor, bool normalize, int borderType ) 3. 4. int sdepth = CV_MAT_DEPTH(srcType); 5. int cn = CV_MAT_CN(srcType), sumType = CV_64F; 6. if ( sdepth = CV_32S 10. sumType = CV_MAKETYPE( sumType, cn ); 11. 11. Ptr rowFilter = getRowSumFilter(srcType, su

4、mType, ksize.width, anc hor.x ); 12. Ptr columnFilter = getColumnSumFilter(sumType, 13. dstType, ksize.height, anchor.y, normalize ? 1./(ksize.width*ksize.height) : 1); 15. 14. return Ptr( new FilterEngine(Ptr(0), rowFilter, columnFilter, 15. srcType, dstType, sumType, borderType ); 16. 这里blur和boxFi

5、lter的区别是，blur是标准化后的boxFilter，即boxFilter的核函 数: normii = true cr = hsizeAvidth * hsize.height 其中, 11af/ipnvis cpp view pla in copy 1. blur( src,dst,Size(1,1 ),Point(-1,-1); 2. blur( 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); 实

6、验结果 F图是对一幅图像分别用 标准方框滤波后的图像: 1*1,4*4,8*8,16*16 【高斯滤波】 高斯滤波器是一类根据高斯函数的形状来选择权值的线性平滑滤波器。它对去除服从正态 分布的噪声很有效。 常用的零均值离散高斯滤波器函数： g(x) = exp(- 7) 2 CT 2D图像中表示为： 代码 cpp view pla in copy 1. 2. 3. 4. Gaussian Blur */ / 5. cv:Mat cv:getGaussianKernel(int n, double sigma, int ktype ) 6. 7. const int SMALL_GAUSSIAN

7、_SIZE = 7; 8. static const float small_gaussian_tabSMALL_GAUSSIAN_SIZE = 9. 10. 1.f, 11. 0.25f, 0.5f, 0.25f, 12. 0.0625f, 0.25f, 0.375f, 0.25f, 0.0625f, 13. 0.03125f, 0.109375f, 0.21875f, 0.28125f, 0.21875f, 0.109375f, 0.03125f 14. ； 15. 16. const float * fixed_kernel = n % 2 = 1 18. 19. CV_Assert(

8、ktype = CV_32F | ktype = CV_64F ); 20. Mat kernel(n, 1, ktype); 21. float * cf = ( float *)kernel.data; 22. double * cd = ( double *)kernel.data; 23. 24. 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; 27. 28. int i; 29. for

9、( i = 0; i n; i+ ) 30. 31. double x = i - (n-1)*0.5; 32. double t = fixed_kernel ? ( double )fixed_kerneli : std:exp(scale2X*x*x); 33. if ( ktype = CV_32F ) 34. 35. cfi = ( float )t; 36. sum += cfi; 37. 38. else 39. 40. cdi = t; 41. sum += cdi; 42. 43. 44. 45. 46. 47. 48. 49. 50. 51. 52. 53. 54. 55.

10、 56. 57. 58. 59. 60. 61. 62. 63. 64. 65. 66. 67. 68. 69. 70. 71. 72. 73. 74. 75. 76. 77. 78. 79. 80. 81. 82. 83. 84. 85. 86. 87. sum = 1./sum; for ( i = 0; i n; i+ ) if ( ktype = CV_32F ) cfi = ( float )(cfi*sum); else cdi *= sum; return kernel; cv:Ptr cv:createGaussianFilter(int type, Size ksize, d

11、ouble sigma1, double sigma2, int borderType ) int depth = CV_MAT_DEPTH(type); if ( sigma2 = 0 ) sigma2 = sigma1; / automatic detection of kernel size from sigma if ( ksize.width 0 ) ksize.width = cvRound(sigma1*(depth = CV_8U ? 3 : 4)*2 + 1)|1; if ( ksize.height 0 ) ksize.height = cvRound(sigma2*(de

12、pth = CV_8U ? 3 : 4)*2 + 1)|1; CV_Assert( ksize.width 0 sigma1 = std:max( sigma1,0.); sigma2 = std:max( sigma2, 0.); Mat kx = getGaussianKernel( ksize.width, sigma1, std:max(depth, CV_32F); Mat ky; if ( ksize.height = ksize.width else ky = getGaussianKernel( ksize.height, sigma2, std:max(depth, CV_3

13、2F); return createSeparableLinearFilter( type, type, kx, ky, Point(-1,-1), 0, borderType ); 88. 89. void cv:GaussianBlur( InputArray _src, OutputArray _dst. Size ksize, 90. double sigma1, double sigma2, 91. int borderType ) 92. 93. Mat src = _src.getMat(); 94. _dst.create( src.size(), src.type(); 95

14、. Mat dst = _dst.getMat(); 96. 96. if ( borderType != BORDER_CONSTANT ) 97. 98. if ( src.rows = 1 ) 99. ksize.height = 1; 100. if ( src.cols = 1 ) 101. ksize.width = 1; 102. 104. 103. if ( ksize.width = 1 106. return ; 107. 110. 108. #ifdef HAVE_TEGRA_OPTIMIZATION 109. if (sigma1 = 0 111. #endif 115

15、. 112. Ptr f = createGaussianFilter( src.type(), ksize, sigma1, sigma2, borderT ype ); 113. f-apply( src, dst ); 114. 实验结果 F图是对一幅图像分别用1*1,3*3,5*5,9*9标准方框滤波后的图像: 非线性滤波 线性滤波易于构造，且易于从频率响应的角度分析，但如果噪声是散粒噪声而非高斯噪声 时线性滤波不能去除噪声。如图像突然出现很大的值，线性滤波只是转换为柔和但仍可见 的散粒。这时需要非线性滤波。 简单的非线性滤波有 中值滤波，立-截尾均值滤波，定义域滤波和值域滤波。 中值

16、滤波选择每个邻域像素的中值输出；口-截尾均值滤波是指去掉百分率为圧的最小 值和最大值；定义域滤波中沿着边界的数字是像素的距离；值域就是去掉值域外的像素 值。 中值滤波代码 cpp view pla in copy 1. medianBlur ( src, dst, i ); 中值滤波实验 下图是对一幅图像分别用3*3,5*5,7*7,9*9(这里必须是奇数)标准方框滤波后的图像： 【双边滤波】 双边滤波的思想是抑制与中心像素值差别太大的像素，输出像素值依赖于邻域像素值的加 权合： 工仗)诫丿和) 权重系数 取决于定义域核 犯丿肿)= exp( ) 和依赖于数据的值域核 rQ.jJ,/) = e

17、xp(- 2刁 的乘积。相乘后会产生依赖于数据的双边权重函数: 亠匸珂*十(丿厅 双边滤波源码 cpp view pla in copy Bilateral Filtering namespace cv static void bilateralFilter_8u( const Mat int i, j, k, maxk, radius; Size size = src.size(); CV_Assert( (src.type() = CV_8UC1 | src.type() = CV_8UC3) if ( sigma_color = 0 ) sigma_color = 1; if ( sig

18、ma_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, rad

19、ius, radius, radius, radius, borderType ); vector _color_weight(cn*256); vector _space_weight(d*d); vector _space_ofs(d*d); 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19. 20. 21. 22. 23. 24. 25. 26. 27. 28. 29. 30. 31. 32. 33. 34. 35. 36. 37. 38. 39. 40. 41. 42. 43. 44. 45. float *

20、color_weight = float * space_weight = int * space_ofs = 46. 47. 48. 49. 50. 51. 52. 53. 54. 55. 56. 57. 58. 59. 60. 61. 62. 63. 64. 65. 66. 67. 68. 69. 70. 71. 72. 73. 74. 75. 76. 77. 78. 79. 80. 81. 82. 83. 84. 85. 86. 87. 88. 89. /initialize color-related bilateral filter coefficients for ( i = 0;

21、 i 256*cn; i+ ) color_weighti = ( float )std:exp(i*i*gauss_color_coeff); / initialize space-related bilateral filter coefficients for ( i = -radius, maxk = 0; i = radius; i+ ) for ( j = -radius; j radius ) continue ; space_weightmaxk = ( float )std:exp(r*r*gauss_space_coeff); space_ofsmaxk+ = ( int

22、)(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 = sptrj; for ( k = 0; k there is no need to use CV_CAST_8U dp

23、trj = (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 = sptrj, g0 = sptrj+1, r0 = sptrj+2; 90. for ( k = 0; k maxk; k+ ) 91. 92. const uchar* sptr_k = sptr + j + space_ofsk; 93. int b = sptr_k0, g = sptr_k1, r = sptr_k2; 94. float w = space_weightk*color_weightstd:abs(b - b0) + 95. std:abs(g - g0) + std:abs(r - r0); 96. sum_b += b*w; sum_g += g*w; sum_r += r*w; 97. wsum += w;