双边滤波(bilateralFilter)原理及C++实现

写在前面
双边滤波是一种非线性滤波，能够达到去噪保边的效果。相比高斯滤波，双边滤波多了一种掩膜，也就是还考虑了灰度相似性，所以双边滤波是结合图像的空间邻近度和像素值相似度的一种折衷处理。

先看看对比效果：wsize=23*23, sigma(space)=10，sigma(color)=35

在同样的掩膜尺寸、sigma(space)的情况下，双边滤波能去除噪声污染并能保持边缘。

应用：图像增强、图像去雾、图像去噪、图像恢复、美颜等等

opencv函数:

void bilateralFilter(InputArray src, OutputArray dst, int d, double sigmaColor, 
                     double sigmaSpace, int borderType=BORDER_DEFAULT )

Github: 源码

σ的意义及选取
空间域sigma(space）选取：

和高斯滤波一样，sigma（space）越大，图像越平滑，趋于无穷大时，每个权重都一样，类似均值滤波。

sigma（space）越小，中心点权重越大，周围点权重越小，对图像的滤波作用越小，趋于零时，输出等同于原图。

值域sigma(color）选取：

Sigma(color）越大，边缘越模糊，极限情况为simga无穷大，值域系数近似相等（忽略常数时，将近为exp（0）= 1），与高斯模板（空间域模板）相乘后可认为等效于高斯滤波。

Sigma(color）越小，边缘越清晰，极限情况为simga无限接近0，值域系数除了中心位置，其他近似为0（接近exp（-∞） = 0），与高斯模板（空间域模板）相乘进行滤波的结果等效于源图像。

基于OpenCV的c++实现
获取值域高斯核模板
优化：可用查表的方式计算，因为灰度差的范围是固定的。

 
//获取色彩模板（值域模板）
///
void getColorMask(std::vector<double> &colorMask , double colorSigma){
 
	for (int i = 0; i < 256; ++i){
		double colordiff = exp(-(i*i) / (2 * colorSigma * colorSigma));
		colorMask.push_back(colordiff);
	}
 
}

Note:注意啦！！！(x3)，我用的表的范围是0_{255，大家是不是也看到有博客用的是0}255channels呢，其实OpenCV中用的也是0~255channels。一开始，我还没弄明白，灰度值差异最大不就是255么，为啥要乘通道数？其实我发现是我处理彩色图像实现的方式和OpenCV的方式不同，我是每个通道分别计算其值域高斯核模板，所以灰度差异最大是255，但是opencv实现是这样的，我贴一下opencv源码:（源码在modules/imgproc/src/smooth.cpp 从2828行开始就是双边滤波）

        float w = space_weight[k]*color_weight[std::abs(b - b0) 
        + std::abs(g - g0) + std::abs(r - r0)];

OpenCV是三个通道的像素值先加起来再相减的，所以差异最大的当然是255*channels了。

不过我觉得两种方式都差不多，只要原理没错就行，后面我会做一个效果对比。

获取空间域高斯核模板（高斯滤波模板）

 
//获取高斯模板（空间模板）
///
void getGausssianMask(cv::Mat& Mask, cv::Size wsize, double spaceSigma){
	Mask.create(wsize, CV_64F);
	int h = wsize.height;
	int w = wsize.width;
	int center_h = (h - 1) / 2;
	int center_w = (w - 1) / 2;
	double sum = 0.0;
	double x, y;
 
	for (int i = 0; i < h; ++i){
		y = pow(i - center_h, 2);
		double* Maskdate = Mask.ptr<double>(i);
		for (int j = 0; j < w; ++j){
			x = pow(j - center_w, 2);
			double g = exp(-(x + y) / (2 * spaceSigma * spaceSigma));
			Maskdate[j] = g;
			sum += g;
		}
	}
}

双边滤波

注意：两个模板权重先相乘，然后还要求和，最后再归一化！（我一开始就忘了归一化了…）；求灰度差异时要用绝对值（有时候忘记了…）。

 
//双边滤波
///
void bilateralfiter(cv::Mat& src, cv::Mat& dst, cv::Size wsize, double spaceSigma, double colorSigma){
	cv::Mat spaceMask;
	std::vector<double> colorMask;
	cv::Mat Mask0 = cv::Mat::zeros(wsize, CV_64F);
	cv::Mat Mask1 = cv::Mat::zeros(wsize, CV_64F);
	cv::Mat Mask2 = cv::Mat::zeros(wsize, CV_64F);
 
	getGausssianMask(spaceMask, wsize, spaceSigma);//空间模板
	getColorMask(colorMask, colorSigma);//值域模板
	int hh = (wsize.height - 1) / 2;
	int ww = (wsize.width - 1) / 2;
	dst.create(src.size(), src.type());
	//边界填充
	cv::Mat Newsrc;
	cv::copyMakeBorder(src, Newsrc, hh, hh, ww, ww, cv::BORDER_REPLICATE);//边界复制;
 
	for (int i = hh; i < src.rows + hh; ++i){
		for (int j = ww; j < src.cols + ww; ++j){
			double sum[3] = { 0 };
			int graydiff[3] = { 0 };
			double space_color_sum[3] = { 0.0 };
 
			for (int r = -hh; r <= hh; ++r){
				for (int c = -ww; c <= ww; ++c){
					if (src.channels() == 1){
						int centerPix = Newsrc.at<uchar>(i, j);
						int pix = Newsrc.at<uchar>(i + r, j + c);
						graydiff[0] = abs(pix - centerPix);
						double colorWeight = colorMask[graydiff[0]];
						Mask0.at<double>(r + hh, c + ww) = colorWeight * spaceMask.at<double>(r + hh, c + ww);//滤波模板
						space_color_sum[0] = space_color_sum[0] + Mask0.at<double>(r + hh, c + ww);
 
					}
					else if (src.channels() == 3){
						cv::Vec3b centerPix = Newsrc.at<cv::Vec3b>(i, j);
						cv::Vec3b bgr = Newsrc.at<cv::Vec3b>(i + r, j + c);
						graydiff[0] = abs(bgr[0] - centerPix[0]); graydiff[1] = abs(bgr[1] - centerPix[1]); graydiff[2] = abs(bgr[2] - centerPix[2]);
						double colorWeight0 = colorMask[graydiff[0]];
						double colorWeight1 = colorMask[graydiff[1]];
						double colorWeight2 = colorMask[graydiff[2]];
						Mask0.at<double>(r + hh, c + ww) = colorWeight0 * spaceMask.at<double>(r + hh, c + ww);//滤波模板
						Mask1.at<double>(r + hh, c + ww) = colorWeight1 * spaceMask.at<double>(r + hh, c + ww);
						Mask2.at<double>(r + hh, c + ww) = colorWeight2 * spaceMask.at<double>(r + hh, c + ww);
						space_color_sum[0] = space_color_sum[0] + Mask0.at<double>(r + hh, c + ww);
						space_color_sum[1] = space_color_sum[1] + Mask1.at<double>(r + hh, c + ww);
						space_color_sum[2] = space_color_sum[2] + Mask2.at<double>(r + hh, c + ww);
 
					}
				}
			}
 
			//滤波模板归一化
			if(src.channels() == 1)
			Mask0 = Mask0 / space_color_sum[0]; 
			else{
				Mask0 = Mask0 / space_color_sum[0];
				Mask1 = Mask1 / space_color_sum[1];
				Mask2 = Mask2 / space_color_sum[2];
			}
				
 
			for (int r = -hh; r <= hh; ++r){
				for (int c = -ww; c <= ww; ++c){
 
					if (src.channels() == 1){
						sum[0] = sum[0] + Newsrc.at<uchar>(i + r, j + c) * Mask0.at<double>(r + hh, c + ww); //滤波
					}
					else if(src.channels() == 3){
						cv::Vec3b bgr = Newsrc.at<cv::Vec3b>(i + r, j + c); //滤波
						sum[0] = sum[0] + bgr[0] * Mask0.at<double>(r + hh, c + ww);//B
						sum[1] = sum[1] + bgr[1] * Mask1.at<double>(r + hh, c + ww);//G
						sum[2] = sum[2] + bgr[2] * Mask2.at<double>(r + hh, c + ww);//R
					}
				}
			}
 
			for (int k = 0; k < src.channels(); ++k){
				if (sum[k] < 0)
					sum[k] = 0;
				else if (sum[k]>255)
					sum[k] = 255;
			}
			if (src.channels() == 1)
			{
				dst.at<uchar>(i - hh, j - ww) = static_cast<uchar>(sum[0]);
			}
			else if (src.channels() == 3)
			{
				cv::Vec3b bgr = { static_cast<uchar>(sum[0]), static_cast<uchar>(sum[1]), static_cast<uchar>(sum[2]) };
				dst.at<cv::Vec3b>(i - hh, j - ww) = bgr;
			}
 
		}
	}
 
}

效果
wsize=23*23， sigma(space)=10，sigma(color)=35

灰度图：效果差别不大，因为灰度差异最大都是255，所以计算值域高斯核模板的结果都一样。

彩色图：有区别，Opencv自带的函数平滑程度低一些，例如脸上的黑点还存在。

参考文章：
https://www.cnblogs.com/wangguchangqing/p/6416401.html（有理论、有实验、讲的很好）

https://blog.csdn.net/Jfuck/article/details/8932978（图文并茂）