图像处理知识点2：均值滤波算法实现（C++）

今天还有些时间，就把我之前的实现均值滤波之后仿照他人进行优化的思想进行实现，下面是均值滤波的实现及优化，优化采用积分图的优化方式，最近在学CPU指令集SSE和AVX优化，等我学的差不多再说，均值滤波的原理就不用说了，网上的资料不少。好了废话不多说开始。

第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;
}

第二部分均值滤波算法

void meanFilter(unsigned char* corrupted, unsigned char* smooth, int width, int height)
{

	memcpy(smooth, corrupted, width * height * sizeof(unsigned char));  //复制了一个原图

	for (int j = 1; j < height - 1; j++)       //边界不处理，因为起始是（1，1）
	{
		for (int i = 1; i < width - 1; i++)
		{
			smooth[j * width + i] = (corrupted[(j - 1) * width + (i - 1)] + corrupted[(j - 1) * width + i] + corrupted[(j - 1) * width + (i + 1)] +
				corrupted[j * width + (i - 1)] + corrupted[j * width + i] + corrupted[j * width + (i + 1)] +
				corrupted[(j + 1) * width + (i - 1)] + corrupted[(j + 1) * width + i] + corrupted[(j + 1) * width + (i + 1)]) / 9;
		}
	}
}

第三部分均值滤波算法2

void MeanWave(unsigned char* img, unsigned char* out, int size,int height,int width)
{
	vector<vector<int> > a(height, vector<int>(width));

		for (int i = 0; i < height; i++)
		{
			for (int j = 0; j < width; j++)
			{
				a[i][j] = 0;
			}
		}
		for (int i = 0; i <height; i++)
		{
			int sum = 0;
			for (int j = 0; j < width; j++)
			{
				sum += img[i*width+j];
				if (i - 1 >= 0)
					a[i][j] = a[i - 1][j] + sum;
				else
					a[i][j] = sum;
			}
		}
		for (int i = 0; i < height; i++)
		{
			for (int j = 0; j < width; j++)
			{
				int w = size;

				if (i - w - 1 >= 0 && i + w < height && j - w - 1 >= 0 && j + w < width)
				{
					int aa = a[i - w - 1][j - w - 1], bb = a[i - w - 1][j + w];
					int cc = a[i + w][j - w - 1], dd = a[i + w][j + w];
					out[i*width+j] = (dd + aa - bb - cc) / pow(2 * size + 1, 2);
				}
				else
					out[i*width+j] = img[i*width+j];
			}
		}
}

第三部分均值滤波优化算法3

/*********************************************求图像的积分图******************************************/
void GetGrayIntegralImage(unsigned char* Src, int* Integral, int Width, int Height, int Stride)
{
	int* ColSum = (int*)calloc(Width, sizeof(int));        //    用的calloc函数哦，自动内存清0
	memset(Integral, 0, (Width + 1) * sizeof(int));
	for (int Y = 0; Y < Height; Y++)
	{
		unsigned char* LinePS = Src + Y * Stride;
		int* LinePL = Integral + Y * (Width + 1) + 1;
		int* LinePD = Integral + (Y + 1) * (Width + 1) + 1;
		LinePD[-1] = 0;
		for (int X = 0; X < Width; X++)
		{
			ColSum[X] += LinePS[X];
			LinePD[X] = LinePD[X - 1] + ColSum[X];
		}
	}
	free(ColSum);
}
/**********************************3积分图进行快速均值滤波算法（在上一步的基础上再进行优化优化2）（Stride为一行像素所占字节数，Radius为滤波核半径）*******************************************/
void BoxBlur(unsigned char* Src, unsigned char* Dest, int Width, int Height, int Stride, int Radius)
{
	int* Integral = (int*)malloc((Width + 1) * (Height + 1) * sizeof(int));
	GetGrayIntegralImage(Src, Integral, Width, Height, Stride);
	for (int Y = 0; Y < Height; Y++)
	{
		int Y1 = max(Y - Radius, 0);
		int Y2 = min(Y + Radius + 1, Height - 1);
		int* LineP1 = Integral + Y1 * (Width + 1);
		int* LineP2 = Integral + Y2 * (Width + 1);
		unsigned char* LinePD = Dest + Y * Stride;
		for (int X = 0; X < Width; X++)
		{
			int X1 = max(X - Radius, 0);
			int X2 = min(X + Radius + 1, Width);
			int Sum = LineP2[X2] - LineP1[X2] - LineP2[X1] + LineP1[X1];
			int PixelCount = (X2 - X1) * (Y2 - Y1);
			LinePD[X] = (Sum + (PixelCount >> 1)) / PixelCount;
		}

		Dest= LinePD- Y * Stride;
	}
	free(Integral);
}

第四部分主函数

int main()
{
	int width, height, func = 0, threshold = 0;
	clock_t startTime, endTime;
    startTime = clock();       //计时开始
	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);

	//meanFilter(pGryImg, dstImg, width, height);

	//积分图实现不受滤波核大小影响
	//MeanWave(pGryImg, dstImg, 9, height, width);

	BoxBlur(pGryImg, dstImg, width, height, width, 4);      

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

	endTime = clock();    //计时结束
    cout << "The run time is: " << (float)(endTime - startTime) / CLOCKS_PER_SEC << "s" << endl;
	delete[] dstImg;
	delete[] pGryImg;
}

实现结果：通过HALCON21.05进行观察就好了，灰度剖面图就可以看是否实现了均值滤波效果。通过halcon可以观察你的算法实现的时候准确。下图就是灰度分布曲线例子。

积分图实现均值滤波看了博客上的大佬作者名（imageshop）。搞算法优化可以看他。有时间再继续更新其他图像处理算法吧。最近看引导滤波实现呢。