深入剖析高斯滤波

1、高斯滤波（高斯平滑）是数字图像处理和计算机视觉里面最常见的操作。平时，我们都是用Matlab 或者OpenCV 的函数调用：imfilter 或者cvSmooth ,并不关心底层的实现。然 而当开发者要自己实现高斯滤波的时候，往往就会很迷惘，经常会被下面几个问题所困扰：1、给定sigma，即标准偏差，怎么确定离散化后的滤波器窗口大小？2、 给定窗口大小，怎么计算高斯核的sigma ,即标准偏差？3、怎么实现可分离滤波器？在网上搜了一下，还真没几个人把实现的细节讲清楚了。这里，我尝试结合三份源码，做个小小的总结。三份源码分别是：1、OpenCV 中的 cvfilter.cpp2、autopano

2、-sift-c 中的 GaussianConvolution.c3、GIMP 中的 blur-gauss.c 和 unsharp-mask.c在图像处理中，高斯滤波一般有两种实现方式，一种是用离散化的滑动窗口进行卷积运算，另一种则是通过傅里叶变换来实现。最常见的就是第一种卷积核实现，只有当离散化的窗口非常大，计算非常耗时（这时可使用可分离滤波器）的情况下，可能会考虑基于傅里叶变换 的实现方法。这里我们只讨论第一种方法。二维高斯函数的形式是这样的：I (J-Jo)2 , (v-yp2 y) = Ae ' 气 2try有着如下的形状:基本上，离散化的主旨就是保留高斯函数中心能量最集中的中间

3、部分，忽略四周能量很小的平坦区域。这只是个很感性的描述， 具体实现起来，就会出现千奇百怪的版本。下面结合三 份源码，看看现实世界里的高斯平滑到底长的什么样子。首先是第一个问题：给定 sigma ,怎么计算窗口大小？直接上OpenCV 的源码，在cvFilter函数中：paraml = cvRound(sigma1*(depth = CV_8U ? 3 : 4)*2 + 1)|1;OpenCV认为半径为3*sigma的窗口就是高斯函数的能量最集中的区域。(在图像位深度不是8的时候，使用4*sigma 半径的窗口？)autopan0-sift-c是图像拼接软件 hugin里面的sift实现，在实现

4、DoG的时候需要做不同尺度的高斯平滑，在 GaussianConvolution_new1函数中实现如下：dim = 1 + 2 * (int) (3.0 * sigma);可见autopano也是实现的 3*sigma 半径的窗口。在GIMP里，实现比较奇特，在 blur_gauss.c 的make_rle_curve函数里面，const gdouble sigma2 = 2 * sigma * sigma;const gdouble l = sqrt (-sigma2 * log (1.0 / 255.0);int n = ceil (l) * 2;if (n % 2) = 0)n +=

5、1;从效果来看，这个实现的窗口半径是约等于2.2*sigma 。然后是第二个问题：给定窗口大小，怎么计算 sigma ?OpenCV 的实现，在 cvFilter.cpp 的 init_gaussian_kernel函数中:sigmaX = sigma > 0 ? sigma : (n/2 1)*0.3 + 0.8;sigma的大小约为窗口半径的0.3倍再加上0.8 。autopano没有实现这个特性。GIMP的实现：/* we want to generate a matrix that goes out a certain radius* from the center, so we

6、 have to go out ceil(rad-0.5) pixels,* inlcuding the center pixel. Of course, thatn 6 nes otiniction,* so we have to go the same amount in the other direction, but not count* the center pixel again. So we double the previous result and subtract* one.* The radius parameter that is passed to this func

7、tion is used ass a little confusing.* the standard deviation, and the radius of effect is the* standard deviation * 2. It*/ radius = fabs (radius) + 1.0;std_dev = radius;radius = std_dev * 2;/* go out ' radius ' in each direction */matrix_length = 2 * ceil (radius 0.5) + 1;注释讲的很清楚了，基本上就是认为si

8、gma应该等于窗口半径的一半。看完这三份源码，结论就是，关于 sigma和窗口的大小，你爱怎么算都可以，这个由你的 实际需要决定，别太离谱就行。(根据概率理论，高斯分布的样本绝大部分集中在3倍标准差以内的区域。)第三个问题是可分离滤波器：首先说明为什么要使用可分离滤波器。实际上，模板运算(滑动窗口卷积)在数字图像处理中是一项非常耗时的运算。1 X 16。2 1'2 42121I以上图中的3*3高斯模板为例，每个像素完成一次模板操作要用9个乘法、8个加法和1个除法。对于一幅n*n的图像，大约就是 9n2个乘法，8n 2个加法和n2个除法，这对于比 较大的图像来说，是非常可怕的。而且随着模

9、板大小的增加，计算量是呈指数增长的。 那么有没有一种办法能够减少计算量呢？答案是肯定的。由于高斯函数可以写成可分离的形式，因此可以采用可分离滤波器实现来加速。所谓的可分离滤波器，就是可以把一个多维的卷积化成多个一维的卷积。具体到二维的高斯滤波，就是指先对行做一维卷积，再对列做一维卷 积。这样就可以将计算复杂度从O(M*M*N*N) 降到O(2*M*M*N), M, N分别是图像和滤波器的窗口大小。问题是实现的时候怎么计算一维的卷积核呢？其实很简单，按照前面计算出来的窗口大小，将二维的高斯模板合并成一维，计算所有离散点上一维高斯函数的权值，最后将权值之和归一化到1。下面是来自OpenCV 的源码：for( i = 0; i <= n/2; i+ )(double t = fixed_kernel ? (double)fixed_kerneli : exp(scale2X*i*i);if( type = CV_32FC1 )(cf(n/2+i)*step = (float)t;sum += cf(n/2+i)*step*2;)else(cd(n/2+i)*step = t;sum += cd(n/2+i)*step*2;)sum = 1./sum;for( i = 0; i &l