图像处理知识点4：高斯滤波之两个一维卷积实现（C++）

       根据高斯函数的可分离特性，二维高斯卷积可以变成两个一维的实现，而且现实中我们不一定都用两维实现，一个一维就达到了预期平滑效果为什么还用二维的呢？但是高斯滤波不具备保护边缘的特性，如果需要保护边缘的滤波方式可以用双边滤波，因为双边滤波再高斯滤波的基础上增加了像素差值域的判断，具备边缘保护效果。双边滤波的文章留着下一次写吧。其实我之前的Canny算子里具有高斯滤波，只不过那个是直接二维实现的。

1.第一部分：图像读写

#define _CRT_SECURE_NO_WARNINGS
#include <windows.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <time.h>
#include<iostream>
#include<assert.h>
#include<vector>
#define WIDTHBYTES(bits) (((bits)+31)/32*4)
 
using namespace std;
 
BYTE* Read8BitBmpFile2Img(const char* filename, int* width, int* height)
{
	FILE* BinFile;
	BITMAPFILEHEADER FileHeader;
	BITMAPINFOHEADER BmpHeader;
	BYTE* plmg;
	unsigned int size;
	int Suc = 1, w, h;
 
	//open file
	*width = *height = 0;
	if ((BinFile = fopen(filename, "rb")) == NULL) return NULL;
	//read struct info
	if (fread((void*)&FileHeader, 1, sizeof(FileHeader), BinFile) != sizeof(FileHeader)) Suc = -1;
	if (fread((void*)&BmpHeader, 1, sizeof(BmpHeader), BinFile) != sizeof(BmpHeader)) Suc = -1;
	if ((Suc == -1) || (FileHeader.bfOffBits < sizeof(FileHeader) + sizeof(BmpHeader)))
	{
		fclose(BinFile);
		return NULL;
	}
	//read Image data
	*width = w = BmpHeader.biWidth;
	*height = h = BmpHeader.biHeight;
	size = w * h;
	fseek(BinFile, FileHeader.bfOffBits, SEEK_SET);
	if ((plmg = new BYTE[size]) != NULL)
	{
		for (int i = 0; i < h; i++)  //0,1,2,3,4(5);400-499
		{
			if (fread(plmg + (h - 1 - i) * w, sizeof(BYTE), w, BinFile) != w)
			{
				fclose(BinFile);
				delete plmg;
				plmg = NULL;
				return NULL;
			}
			fseek(BinFile, (w + 3) / 4 * 4 - w, SEEK_CUR);
		}
	}
	fclose(BinFile);
	return plmg;
}
 
bool Write8BitImg2BmpFile(BYTE* pImg, int width, int height, const char* filename)
//当宽度不是4的倍数时自动添加成4的倍数
{
	FILE* BinFile;
	BITMAPFILEHEADER FileHeader;
	BITMAPINFOHEADER BmpHeader;
	int i, extend;
	bool Suc = true;
	BYTE p[4], * pCur;
 
	//Open File
	if ((BinFile = fopen(filename, "w+b")) == NULL) { return false; }
	//Fill the FileHeader
	FileHeader.bfType = ((WORD)('M' << 8) | 'B');
	FileHeader.bfOffBits = sizeof(BITMAPFILEHEADER) + sizeof(BITMAPINFOHEADER) + 256 * 4L;
	FileHeader.bfSize = FileHeader.bfOffBits + width * height;
	FileHeader.bfReserved1 = 0;
	FileHeader.bfReserved2 = 0;
	if (fwrite((void*)&FileHeader, 1, sizeof(FileHeader), BinFile) != sizeof(FileHeader)) Suc = false;
	//Fill the ImgHeader
	BmpHeader.biSize = 40;
	BmpHeader.biWidth = width;
	BmpHeader.biHeight = height;
	BmpHeader.biPlanes = 1;
	BmpHeader.biBitCount = 8;
	BmpHeader.biCompression = 0;
	BmpHeader.biSizeImage = 0;
	BmpHeader.biXPelsPerMeter = 0;
	BmpHeader.biYPelsPerMeter = 0;
	BmpHeader.biClrUsed = 0;
	BmpHeader.biClrImportant = 0;
	if (fwrite((void*)&BmpHeader, 1, sizeof(BmpHeader), BinFile) != sizeof(BmpHeader)) Suc = false;
	//write Pallete
	for (i = 0, p[3] = 0; i < 256; i++)
	{
		p[3] = 0;
		p[0] = p[1] = p[2] = i;//blue,green,red
		if (fwrite((void*)p, 1, 4, BinFile) != 4) { Suc = false; break; }
	}
	//write image data
	extend = (width + 3) / 4 * 4 - width;
	if (extend == 0)
	{
		for (pCur = pImg + (height - 1) * width; pCur >= pImg; pCur -= width)
		{
			if (fwrite((void*)pCur, 1, width, BinFile) != (unsigned int)width) Suc = false;//真实的数据
		}
	}
	else
	{
		for (pCur = pImg + (height - 1) * width; pCur >= pImg; pCur -= width)
		{
			if (fwrite((void*)pCur, 1, width+1, BinFile) != (unsigned int)width) Suc = false;//真实的数据
			for (i = 0; i < extend; i++)//扩充的数据
				if (fwrite((void*)(pCur + width - 1), 1, 1, BinFile) != 1) Suc = false;
		}
	}
	//return
	fclose(BinFile);
	return Suc;
}

2.第二部分：图像处理高斯滤波函数

/*************************生成高斯滤波核***********************************/
void MakeGauss(double sigma, double** pdKernel, int* pnWindowSize)
{
	int i;
	int nCenter;
	// 数组的某一点到中心点的距离
	double dDis;
	double PI = 3.14159;
	double dValue;
	double dSum;
	dSum = 0;
	// 数组长度，根据概率论的知识，选取[-3*sigma, 3*sigma]以内的数据。
	// 这些数据会覆盖绝大部分的滤波系数
	*pnWindowSize = 1 + 2 * ceil(3 * sigma);
	nCenter = (*pnWindowSize) / 2;
	*pdKernel = new double[*pnWindowSize];
	for (i = 0; i < (*pnWindowSize); i++)
	{
		dDis = (double)(i - nCenter);
		dValue = exp(-(1 / 2) * dDis * dDis / (sigma * sigma)) / (sqrt(2 * PI) * sigma);
		(*pdKernel)[i] = dValue;
		dSum += dValue;
	}
	//归一化
	for (i = 0; i < (*pnWindowSize); i++)
	{
		(*pdKernel)[i] /= dSum;
	}
}

/*************************************************************************
*
*\函数名称：
*   GaussianSmooth()
*
* \输入参数:
*   unsigned char * pUnchImg    - 指向图象数据的指针
*   int nWidth           - 图象数据宽度
*   int nHeight          - 图象数据高度
*   double dSigma         - 高斯函数的标准差
*   unsigned char * pUnchSmthdImg - 指向经过平滑之后的图象数据
* \说明:
*   为了抑止噪声，采用高斯滤波对图象进行滤波，滤波先对x方向进行，然后对
*   y方向进行。
*
*************************************************************************
*/
void GaussianSmooth(BYTE* pUnchImg, int nWidth, int nHeight, double sigma, BYTE* pUnchSmthdImg)
{
	// 循环控制变量
	int y;
	int x;
	int i;
	// 高斯滤波器的数组长度
	int nWindowSize;
	// 窗口长度的1/2
	int nHalfLen;
	// 一维高斯数据滤波器
	double* pdKernel;
	// 高斯系数与图象数据的点乘
	double dDotMul;
	// 高斯滤波系数的总和
	double dWeightSum;
	// 中间变量
	double* pdTmp;
	// 分配内存
	pdTmp = new double[nWidth * nHeight];
	// 产生一维高斯数据滤波器
	MakeGauss(sigma, &pdKernel, &nWindowSize);
	// MakeGauss返回窗口的长度，利用此变量计算窗口的半长
	nHalfLen = nWindowSize / 2;
	// x方向进行滤波
	for (y = 0; y < nHeight; y++)
	{
		for (x = 0; x < nWidth; x++)
		{
			dDotMul = 0;
			dWeightSum = 0;
			for (i = (-nHalfLen); i <= nHalfLen; i++)
			{
				// 判断是否在图象内部
				if ((i + x) >= 0 && (i + x) < nWidth)
				{
					dDotMul += (double)pUnchImg[y * nWidth + (i + x)] * pdKernel[nHalfLen + i];
					dWeightSum += pdKernel[nHalfLen + i];
				}
			}
			pdTmp[y * nWidth + x] = dDotMul / dWeightSum;
		}
	}
	// y方向进行滤波
	for (x = 0; x < nWidth; x++)
	{
		for (y = 0; y < nHeight; y++)
		{
			dDotMul = 0;
			dWeightSum = 0;
			for (i = (-nHalfLen); i <= nHalfLen; i++)
			{
				// 判断是否在图象内部
				if ((i + y) >= 0 && (i + y) < nHeight)
				{
					dDotMul += (double)pdTmp[(y + i) * nWidth + x] * pdKernel[nHalfLen + i];
					dWeightSum += pdKernel[nHalfLen + i];
				}
			}
			pUnchSmthdImg[y * nWidth + x] = (unsigned char)(int)dDotMul / dWeightSum;
			//	pUnchSmthdImg[y*nWidth + x] = (int)dDotMul/dWeightSum ;
		}
	}

	// 释放内存
	delete[]pdKernel;
	pdKernel = NULL;
	delete[]pdTmp;
	pdTmp = NULL;
}

3.第三部分：主函数

int main()
{
	int width, height, func = 0, threshold = 0;

	char readPath[] = "D:\\机器视觉学习工程\\lianxi\\me.bmp";
	BYTE* pGryImg = Read8BitBmpFile2Img(readPath, &width, &height);        //注意：下面对应位置也要该改名。

	BYTE* dstImg = new BYTE[width * height];
	memset(dstImg, 0, width * height);
	//printf("%d,%d\n", width, height);

	int size = 7; //定义卷积核大小
	double** gaus = new double* [size];  //卷积核数组
	for (int i = 0; i < size; i++)
	{
		gaus[i] = new double[size];  //动态生成矩阵
	}

	MakeGauss(1, gaus, 7);
	GaussianSmooth(pGryImg, width, height, gaus, 7, dstImg);

	cout << "操作完成" << endl;
	char writePath[] = "D:\\机器视觉学习工程\\lianxi\\megao.bmp";
	Write8BitImg2BmpFile(dstImg, width, height, writePath);
	ShellExecuteA(nullptr, "open", writePath, "", "", SW_SHOW);     //打开写入的位图结果

	delete[] dstImg;
	delete[] pGryImg;
}

完事！一次性写两章，下次有时间再写