susan角点检测
suan(smallest univalue segment assimilating nucleus 最小核值相似区)
基本原理
SUSAN简单而言就是是利用一个固定尺寸的圆形窗口,但是其没有滑动,直接比较圆形窗口内像素与中心像素(Nucleus)的灰度差异,统计整个圆形窗口内灰度差异情况,利用阈值判断进行角点检测。
USAN(核值相似区)
对于图像中非纹理区域(不具备明显重复或规律性结构的区域)的任一点,在以它为中心的模板窗中存在一 块亮度与其相同的区域,这块区域即为 USAN 区域。
USAN区域
USAN的典型区域如图所示。模板在图像上移动时,
- 当圆形模板完全在背景或者目标区域时,其USAN区域最大,如图(a);
- 当核心在边缘时,USAN区域减少一半,如图(c);
- 当核心在角点时, USAN区域最小,如图(d)。基于这一原理, Smith提出了最小核值相似区角点检测算法。
模板窗口
实现步骤
-
用高斯滤波器对图像进行平滑处理
-
在图像上放置一个36个像素的圆形模板,模板在图像上滑动,依次比较模板内各个像素点的灰度与模板核的灰度,判断是否属于USAN区域。判别函数如下:
-
统计圆形模板中和核心点有相似亮度值的像素个数n(r0)。
其中,D(r0)是以r0为中心的圆形模板区域。 -
使用如下角点响应函数。若某个像素点的USAN值小于某一特定阈值,则该点被认为是初始角点,其中,g可以设定为USAN的最大面积的一半。
-
对初始角点进行非极值抑制来求得最后的角点。
代码实现
void susanCornerDetection(cv::Mat& src, cv::Mat& dst, int g = 18)
{
if (src.empty())
{
std::cout << "Failed to open file!!!" << std::endl;
return;
}
//使用高斯滤波器对图像进行处理
cv::Mat dst1;
cv::GaussianBlur(src, dst1,cv::Size(5, 5), 0, 0);
cv::Mat grayImg;
dst = src.clone();
cv::cvtColor(dst1, grayImg, cv::COLOR_BGR2GRAY);
//在x,y方向的模板
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 };
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 };
//对图像边缘进行填充
cv::Mat grayImg_padded;
copyMakeBorder(grayImg, grayImg_padded, 3, 3, 3, 3, cv::BORDER_REFLECT);
//寻找阈值
int max = 0, min = 255;
for (int i = 0; i < grayImg.rows; i++)
{
for (int j = 0; j < grayImg.cols; j++)
{
if (grayImg.at<uchar>(i, j) > max)
{
max = grayImg.at<uchar>(i, j);
}
if (grayImg.at<uchar>(i, j) < min)
{
min = grayImg.at<uchar>(i, j);
}
}
}
int t = (max - min) / 10;
//统计圆形模板中和核心点有相似亮度值的像素个数
for (int i = 3; i < grayImg_padded.rows - 3; i++)
{
for (int j = 3; j < grayImg_padded.cols - 3; j++)
{
int same = 0;//same表示近似像素的数量
for (int k = 0; k < 37; k++)
{
if (abs(grayImg_padded.at<uchar>(i + OffSetX[k], j + OffSetY[k]) - grayImg_padded.at<uchar>(i, j)) <= t)
same++;
}
if (same < g)
{
grayImg_padded.at<uchar>(i, j) = g - same;
}
else
{
grayImg_padded.at<uchar>(i, j) = 0;
}
}
}
//非极大值抑制
int x[8] = { -1, -1, -1, 0, 0, 1, 1, 1 };
int y[8] = { -1, 0, 1, -1, 1, -1, 0, 1 };
for (int i = 3; i < grayImg_padded.rows - 3; i++)
{
for (int j = 3; j < grayImg_padded.cols - 3; j++)
{
int flag = 0;
for (int k = 0; k < 8; k++)
{
if (grayImg_padded.at<uchar>(i, j) <= grayImg_padded.at<uchar>(i + x[k], j + y[k]))
{
flag = 1;
break;
}
}
if (flag == 0)
{
line(dst, cv::Point(j - 3, i - 3), cv::Point(j + 3, i - 3), cv::Scalar(255, 0, 255), 2, 8);
line(dst, cv::Point(j - 3, i - 3), cv::Point(j - 3, i + 3), cv::Scalar(255, 0, 255), 2, 8);
}
}
}
}
实验结果分析
不同尺寸高斯滤波器下的实验结果
图像:rectangle.bmp
15* 15
5 *5
7 *7
不同图像变换下的检测效果
rectangle
5*5 15 23
rectangle_smooth1
rectangle_smooth2
rectangle_smooth3
结果分析:USAN算法是基于像素灰度值的局部统计特征来判断角点的,当图像模糊时,像素灰度值的变化会变得不明显,导致角点的响应值减小。因此,图像模糊会降低SUSAN角点检测算法的准确性和稳定性。
harris角点检测
基本原理
实现步骤
1.计算x y 方向的梯度值Mat_x,Mat_y
2.计算Mat_xx,Mat_yy,Mat_xy
3.利用高斯函数对Mat_xx,Mat_yy,Mat_xy进行滤波
4.计算局部特征结果矩阵M的特征值和响应函数
C(i,j) = Det(M) - k(trace(M)^2) k∈(0.04,0.06]
5.将计算完的响应函数的值C进行非极大值抑制,滤除边缘与非角点的点,保留满足大于设定的阈值的区域
6.找出角点
代码实现
void myDetecHarrisCornerAlgorithm(const Mat& src, vector<cv::Point> &points, double k)
{
Mat gray;
if (src.channels() == 3)
{
cvtColor(src, gray, COLOR_BGR2GRAY);
}
else if (src.channels() == 1)
{
gray = src.clone();
}
else
{
cout << "Image channnels is Error! " << endl;
return;
}
gray.convertTo(gray, CV_64F);
//Step1
//1.1 创建x与y方向内核
Mat xKernel = (Mat_<double>(1, 3) << -1, 0, 1);
//[-1,0,1] x方向
Mat yKernel = xKernel.t();//反转矩阵。 该方法通过矩阵表达式进行矩阵求逆。
//Mat yKernel = (Mat_<double>(3, 1) << -1, 0, 1);
//[-1
// 0
// 1] y方向
//1.2卷积获取x与y方向的梯度值
Mat Ix, Iy;
filter2D(gray, Ix, CV_64F, xKernel);
filter2D(gray, Iy, CV_64F, yKernel);
//Step2
//计算Mat_xx,Mat_yy,Mat_xy
Mat Ix2, Iy2, Ixy;
Ix2 = Ix.mul(Ix);// 执行两个矩阵按元素相乘 获取Mat_xx。
Iy2 = Iy.mul(Iy);// 执行两个矩阵按元素相乘 获取Mat_yy。
Ixy = Ix.mul(Iy);// 执行两个矩阵按元素相乘 获取Mat_xy。
//Step3
//3.1获取高斯滤波内核
Mat gaussKernel = getGaussianKernel(7, 1);
//3.2利用高斯函数对Mat_xx,Mat_yy,Mat_xy进行滤波
filter2D(Ix2, Ix2, CV_64F, gaussKernel);
filter2D(Iy2, Iy2, CV_64F, gaussKernel);
filter2D(Ixy, Ixy, CV_64F, gaussKernel);
//Step4
//计算局部特征结果矩阵M的特征值和响应函数
//C(i,j) = Det(M) - k(trace(M)^2) k∈(0.04,0.06]
Mat cornerStrength(gray.size(), CV_64F);
int width = gray.size().width;
int height = gray.size().height;
for (int h = 0; h < height; h++)
{
for (int w = 0; w < width; w++)
{
//M = [Ix2,Ixy
// Ixy,Iy2]
//det = Ix2 * Ix2 - Ixy^2
//trace = Ix2 + Iy2
//C = det - k * trace
double det_m = Ix2.at<double>(h, w) * Iy2.at<double>(h, w) - pow(Ixy.at<double>(h, w), 2);
double trace_m = Ix2.at<double>(h, w) + Iy2.at<double>(h, w);
cornerStrength.at<double>(h, w) = det_m - k * trace_m * trace_m;
}
}
//Step5
//5.1寻找最大值
double maxStrength;
minMaxLoc(cornerStrength, NULL, &maxStrength, NULL, NULL);
//5.2非极大值抑制
Mat dilated;
dilate(cornerStrength, dilated, Mat());
Mat localMax;
compare(cornerStrength, dilated, localMax, CMP_EQ);
//5.3保留满足大于设定的阈值
Mat cornerMap;
double qualityLevel = 0.01;
double thresh = qualityLevel * maxStrength;
cornerMap = cornerStrength > thresh;// 大于标识符重载函数
// 等同于threshold(cornerStrength,cornerMap,thresh,255,THRESH_BINARY)
bitwise_and(cornerMap, localMax, cornerMap);
//Step6
// Iterate over the pixels to obtain all feature points 迭代像素以获得所有特征点
for (int y = 0; y < cornerMap.rows; y++) {
const uchar* rowPtr = cornerMap.ptr<uchar>(y); //行指针
for (int x = 0; x < cornerMap.cols; x++) {
// if it is a feature point 如果是特征点(像素值非0值为特征点)
if (rowPtr[x]) {
points.push_back(cv::Point(x, y));
}
}
}
}
实验结果分析
不同图像变换下的实验结果
rectangle
rectangle_smooth1
rectangle_smooth2
rectangle_smooth3
结果分析:图像模糊会导致角点的局部结构信息模糊化,从而使得角点的响应值减小。Harris算法是基于图像局部结构的二阶导数来判断角点的,当图像模糊时,局部结构的变化会变得不明显,导致角点的响应值减小。因此,图像模糊会降低Harris角点检测算法的准确性和稳定性。