susan角点检测算法

 SUSAN算法是1997年牛津大学的Smith等人提出的一种处理灰度图像的方法,它主要是用来计算图像中的角点特征。SUSAN算法选用圆形模板(如图1所示)。将位于圆形窗口模板中心等待检测的像素点称为核心点。核心点的邻域被划分为两个区域:亮度值相似于核心点亮度的区域即核值相似区(Univalue SegmentAs-similatingNueleus,USAN)和亮度值不相似于核心点亮度的区域。

图1 圆形模板

       USAN的典型区域如图2所示。模板在图像上移动时,当圆形模板完全在背景或者目标区域时，其USAN区域最大,如图2(a)；当核心在边缘时,USAN区域减少一半,如图2(c)；当核心在角点时, USAN区域最小，如图2(d)。基于这一原理, Smith提出了最小核值相似区角点检测算法。

图2 典型区域

SUSAN角点检测算法的具体步骤如下:
(1)在图像上放置一个37个像素的圆形模板，模板在图像上滑动，依次比较模板内各个像素点的灰度与模板核的灰度，判断是否属于USAN区域。判别函数如下：

 

(2)统计圆形模板中和核心点有相似亮度值的像素个数n(r0)。

其中，D(r0)是以r0为中心的圆形模板区域

(3)使用如下角点响应函数。若某个像素点的USAN值小于某一特定阈值，则该点被认为是初始角点，其中，g可以设定为USAN的最大面积的一半。

 

 

(4)对初始角点进行非极值抑制来求得最后的角点。

 

算法实现如下:

IplImage* SUSAN(IplImage* src)
{
	uchar *data0, *data1, *data2;
	int same, max, min, thresh;
	IplImage* img = cvCreateImage(cvGetSize(src), 8, 1); //中间图像
	IplImage* dst = cvCreateImage(cvGetSize(src), 8, 1); //结果图像
	int height = src->height;
	int width = src->width;
	int step = src->widthStep / sizeof(uchar);
	int channels = src->nChannels;
	data0 = (uchar*)src->imageData;
	data1 = (uchar*)img->imageData;
	data2 = (uchar*)dst->imageData;
	int  g = 18; //核值相似区域中点个数的阈值

	//模版 x 和 y的坐标的偏移量
	int OffSetX[37] =
	{ 
			-1, 0, 1,
		-2, -1, 0, 1, 2,
	-3, -2, -1, 0, 1, 2, 3,
	-3, -2, -1, 0, 1, 2, 3,
	-3, -2, -1, 0, 1, 2, 3,
		-2, -1, 0, 1, 2,
			-1, 0, 1
	};
	int OffSetY[37] =
	{ 
			-3, -3, -3,
		-2, -2, -2, -2, -2,
	-1, -1, -1, -1, -1, -1, -1,
		 0, 0, 0, 0, 0, 0, 0,
		 1, 1, 1, 1, 1, 1, 1,
			2, 2, 2, 2, 2,
			   3, 3, 3
	};

	//求阈值
	max = min = data0[0];
	for (int i = 0; i < height; i++)
	{
		for (int j = 0; j<width; j++)
		{
			if (data0[i*step + j]>max) max = data0[i*step + j];
			if (data0[i*step + j]<min) min = data0[i*step + j];
		}
	}
	thresh = (max - min) / 10;//可选取其他方法

	for (int i = 3; i < height - 3; i++)
	{
		for (int j = 3; j < width - 3; j++)
		{
			same = 0;
			for (int k = 0; k < 37; k++)
			{
				if (abs(data0[(i + OffSetY[k])*step + (j + OffSetX[k])] - data0[i*step + j]) < thresh)
					same++;
			}
			 
			if (same < g)//g值可改 
				data1[i*step + j] = g-same;
			else
				data1[i*step + j] = 0;
		}
	}

	//非极大值抑制
	int i_s[8] = { -1, -1, -1, 0, 0, 1, 1, 1 };
	int j_s[8] = { -1, 0, 1, -1, 1, -1, 0, 1 };
	int flag;
	for (int i = 4; i < height - 4; i++)
	{
		for (int j = 4; j < width -4; j++)
		{
			flag = 0;
			for (int k = 0; k < 8; k++)
			{
				if (data1[i*step + j] <= data1[ (i + i_s[k]) * step + (j + j_s[k]) ] )
				{
					flag = 1;
					break;
				}
			}
			if (flag == 0)
			{
				data2[i*step + j] = 255;
			}
			else
			{
				data2[i*step + j] = 0;
			}
		}
	}

	return dst;
}