数字图像处理基本算法实现（3）--section3.5-6空间滤波处理

         这次是《DIP》的section3.5-3.6，也是第3章的最后一个部分了，最后还有一个使用模糊技术的，看着比较麻烦（主要是没大看懂）就不去实现了，下次直接到频域滤波。

        这次的空间滤波主要就是平滑滤波以及一二阶微分的处理，其中线性平滑滤波只实现了简单的3*3滑动平均和高斯平均滤波，还有3*3和5*5的中值滤波，微分处理上，实现了一下二阶的拉普拉斯算子，一阶的Sobel算子，最后将3.6.3的高提升滤波也顺便实现了一下看了一下效果。其实以前学习机器视觉的时候，各种检测边缘的算子：拉普拉斯，Sobel，Prewitt，Robert，Canny等也都学过，也都懒得去实现，现在把简单的实现了一下，就当练练手了，工作了可能会用到，可没有现成的库让咱用。

       先是头文件

class section_3_5_6
{
public:
	enum FilterType{Avg3, Gauss3, Median,Laplace,Laplace8,Sobel};
	Mat& filterImg(Mat &src, FilterType ftype, int fsize=3, bool addsrc=false);
	Mat& promoteFilter(Mat &src, float k, Mat &temp=Mat());
private:
	Mat& linearFilter(Mat &src, FilterType ftype);
	Mat& medianFilter(Mat &src, int fsize);
	Mat& laplaceFilter(Mat &src, FilterType ftype,bool addsrc=false);
	Mat& sobelFilter(Mat &src);
	Mat m_destMat;
};

 

然后是实现，对于边缘的像素点，采用了与前面的几部分一样的复制边缘的方式来扩充。

Mat §ion_3_5_6::filterImg(Mat &src, FilterType ftype, int fsize, bool addsrc)
{
	if(src.channels() != 1)
		return m_destMat;

	switch(ftype)
	{
	case Avg3:
	case Gauss3:
		return linearFilter(src, ftype);
		break;
	case Median:
		return medianFilter(src, fsize);
		break;
	case Laplace:
	case Laplace8:
		return laplaceFilter(src, ftype, addsrc);
		break;
	case Sobel:
		return sobelFilter(src);
		break;
	default:
		break;
	}
	return m_destMat;
}

Mat& section_3_5_6::linearFilter(Mat &src, FilterType ftype)
{
	m_destMat.create(src.rows, src.cols, CV_8UC1);
	
	//offset box;
	const int fsize = 3;
	int box[fsize] = {-1,0,1};
	
	//weight numberator and denominator, default is guassian3;
	int w_num[fsize*fsize]={1,2,1,
							2,4,2,
							1,2,1};
	int w_deno = 16;
	if(ftype==Avg3)
	{
		for(int i=0; i<sizeof(w_num)/sizeof(w_num[0]); ++i)
			w_num[i] = 1;
		w_deno = 9;
	}

	//scanning image and process
	int step = src.step[0];
	int offset_i=0,offset_j=0;
	int count=0, inc = 0;;
	for (int i=0; i<src.rows; ++i)
	{
		for (int j=0; j<src.cols; ++j)
		{
			count = 0;
			inc = 0;
			for (int a=0; a<fsize; ++a)
			{
				for (int b=0;b<fsize; ++b)
				{
					offset_i = std::min(std::max(0,i+box[a]), src.rows-1);
					offset_j = std::min(std::max(0,j+box[b]), src.cols-1);

					count += *(src.data+step*offset_i+offset_j)*w_num[inc++];
				}
			}
			*(m_destMat.data+step*i+j) = saturate_cast<uchar>(count/w_deno);
		}
	}
		
	return m_destMat;
}
//local function used for qsort;
static int compare_uchar(const void *elem1, const void *elem2)
{
	return *(uchar*)elem1 - *(uchar*)elem2;
}

Mat& section_3_5_6::medianFilter(Mat &src, int fsize)
{
	if(fsize!=3 && fsize!=5)
		return m_destMat;

	m_destMat.create(src.rows, src.cols, CV_8UC1);
	
	//array save the pixels
	int rgsize = fsize*fsize;
	uchar *rg = new uchar[rgsize];
	memset(rg, 0, rgsize);
	//offset box
	int *box = new int[fsize];
	for (int i=0; i<fsize; ++i)
		box[i] = i-fsize/2;
	
	int median = fsize==3 ? 4 : 12;

	int step = src.step[0];
	int offset_i=0,offset_j=0;
	int cnt = 0;	
	for (int i=0; i<src.rows; ++i)
	{
		for (int j=0; j<src.cols; ++j)
		{
			cnt = 0;	
			for (int a=0; a<fsize; ++a)
			{
				for (int b=0;b<fsize; ++b)
				{
					offset_i = std::min(std::max(0,i+box[a]), src.rows-1);
					offset_j = std::min(std::max(0,j+box[b]), src.cols-1);

					rg[cnt++] = *(src.data+step*offset_i+offset_j);
				}
			}
			
			qsort((void*)rg, rgsize, sizeof(uchar), compare_uchar);

			*(m_destMat.data+i*step+j) = rg[median];
		}
	}

	return m_destMat;
}

Mat& section_3_5_6::laplaceFilter(Mat &src, FilterType ftype, bool addsrc)
{
	m_destMat.create(src.rows, src.cols, CV_8UC1);
	//init maskbox and offsetbox;
	const int masksize = 9;
	const int fsize = 3;
	int maskbox[masksize] = {0, -1, 0,
							-1,  4, -1,
							 0,	-1,	0};
	if(ftype == Laplace8)
	{
		for(int i=0; i<masksize; ++i)
			maskbox[i] = -1;
		maskbox[4] = 8;
	}
	int offsetbox[fsize] = {-1, 0, 1};

	//scanning image
	int step = src.step[0];
	int offset_i=0,offset_j=0;
	int cnt = 0, result=0;	
	for (int i=0; i<src.rows; ++i)
	{
		for (int j=0; j<src.cols; ++j)
		{
			cnt = 0;	
			result = 0;

			for (int a=0; a<fsize; ++a)
			{
				for (int b=0;b<fsize; ++b)
				{
					offset_i = std::min(std::max(0,i+offsetbox[a]), src.rows-1);
					offset_j = std::min(std::max(0,j+offsetbox[b]), src.cols-1);

					result += *(src.data+step*offset_i+offset_j)*maskbox[cnt++];
				}
			}
			if(addsrc)
				*(m_destMat.data+i*step+j) = saturate_cast<uchar>(*(src.data+i*step+j)+result);
			else
				*(m_destMat.data+i*step+j) = result;
		}
	}

	return m_destMat;
}

Mat& section_3_5_6::sobelFilter(Mat &src)
{
	m_destMat.create(src.rows, src.cols, CV_8UC1);

	//init maskbox and offsetbox;
	const int masksize = 9;
	const int fsize = 3;
	int mask_x[masksize] = {-1, -2, -1,
							 0,  0,  0,
							 1,	 2,	 1};

	int mask_y[masksize] = {-1,  0, 1,
							-2,  0, 2,
							-1,	 0,	1};

	int offsetbox[fsize] = {-1, 0, 1};

	//scanning image
	int step = src.step[0];
	int offset_i=0,offset_j=0;
	int cnt=0;
	int gx=0, gy=0;	
	for (int i=0; i<src.rows; ++i)
	{
		for (int j=0; j<src.cols; ++j)
		{
			cnt = 0;	
			gx = gy = 0;

			for (int a=0; a<fsize; ++a)
			{
				for (int b=0;b<fsize; ++b)
				{
					offset_i = std::min(std::max(0,i+offsetbox[a]), src.rows-1);
					offset_j = std::min(std::max(0,j+offsetbox[b]), src.cols-1);
					
					gx += *(src.data+step*offset_i+offset_j)*mask_x[cnt];
					gy += *(src.data+step*offset_i+offset_j)*mask_y[cnt];
					++cnt;
				}
			}
			//m=|gx|+|gy|;
			*(m_destMat.data+i*step+j) = saturate_cast<uchar>(abs(gx)+abs(gy));
			//or m=sqrt(gx*gx+gy*gy);
//			*(m_destMat.data+i*step+j) = saturate_cast<uchar>(sqrt(float(gx*gx + gy*gy)));
		}
	}

	return m_destMat;
}

Mat& section_3_5_6::promoteFilter(Mat &src, float k, Mat &temp)
{
	//first blur the image using gaussian 3*3;
	linearFilter(src, Gauss3);
	if(m_destMat.empty())
		return m_destMat;

	Mat blured;
	m_destMat.copyTo(blured);

	//temp is used to save the temporary difference mat
	temp.create(src.rows, src.cols, CV_8UC1);
	
	for (int i=0; i<src.rows; ++i)
	{
		uchar *src_ptr = src.ptr<uchar>(i);
		uchar *dest_ptr = m_destMat.ptr<uchar>(i);
		uchar *blur_ptr = blured.ptr<uchar>(i);
		uchar *temp_ptr = temp.ptr<uchar>(i);
		for (int j=0; j<src.cols; ++j)
		{
			dest_ptr[j] = MIN(255, MAX(0, src_ptr[j]+k*(src_ptr[j]-blur_ptr[j])));	
			temp_ptr[j] = MIN(255, MAX(0, k*(src_ptr[j]-blur_ptr[j])));
		}
	}

	return m_destMat;
}

下面是以dog为例的几张处理后的图像：

滑动平均滤波，高斯滤波，中值滤波5*5，以及模板中间为8的Laplace+原图

最后是sobel算子和k=2时的高提升滤波

相比原始图像进行高提升滤波后的图像增强了细节边缘，跟上面的Laplace+原图的效果类似。