OpenCV图像角点检测(19)

亚像素级角点检测

数字图像通常是离散化成像素；每个像素对应一个整数坐标位置；整数坐标位置对于很多应用并不精确，比如跟踪、相机标定、图像配准、图像拼接以及三维重构；为达到有些应用的精确性，需要精确到浮点坐标位置；所以亚像素定位问题。亚像素定位就是计算特征所在图像中的真实位置，而真实位置有时候并不在像素所在整数坐标位置上，而是在像素的内部。

计算过程

• 插值
• 基于图像矩阵计算
• 曲线拟合

左图中点p附近的图像是均匀的，其梯度为0，右图中边缘的梯度与沿边缘方向的q-p向量正交，p点梯度与q-p向量的点积均为0

我们假设起始角点q在实际亚像素角点的附近。检测所有的q-p向量。若点p位于一个均匀区域，则点p的梯度为0。若q-p向量的方向与边缘的方向一致，则此边缘上p点处的梯度与q-p向量正交，在这两种情况下，p点处的梯度与q-p向量的点积为0.我们可以在p点周围找到很多组梯度以及相关的向量q-p，令其点集为0，然后可以通过求解方程组，方程组的解即为角点q的亚像素级精度的位置，即精确的角点位置。

API介绍

void cornerSubPix( InputArray image, InputOutputArray corners,Size winSize, Size zeroZone,TermCriteria criteria );
/*******************************************************************
*			image: 				输入图				
*			corners:	        角点信息
*			winSize:			计算亚像素角点区域大小
*			zeroZone:  			搜索区域
*			criteria:			点精准化迭代过程的终止条件			
*********************************************************************/
//TermCriteria
class CV_EXPORTS TermCriteria
{
public:
    enum Type
    {
        COUNT=1, 			//计算元素或者迭代次数最小值
        MAX_ITER=COUNT,		//最大迭代次数
        EPS=2 				//当满足该精确度时，迭代算法停止
    };
    //TermCriteria::COUNT + TermCriteria::EPS 作为判定条件
    TermCriteria();
    TermCriteria(int type, int maxCount, double epsilon);
    inline bool isValid() const
    {
        const bool isCount = (type & COUNT) && maxCount > 0;
        const bool isEps = (type & EPS) && !cvIsNaN(epsilon);
        return isCount || isEps;
    }
    int type; 			//终止条件类型
    int maxCount; 		//计算的迭代数或者最大元素数
    double epsilon; 	//当达到要求的精确度或参数的变化范围时，迭代算法停止
};

综合代码

#include <iostream>
#include <map>
#include <iomanip>
#include <vector>
#include <opencv2/opencv.hpp>
using namespace std;
using namespace cv;
class AngularPoint
{
public:
	AngularPoint() :img(imread("test.jpg")) 
	{
		result["Harris"] = img;
	}
	void AngularByHarris() 
	{
		Mat norm, harris;
		//图像灰度处理
		cvtColor(img, result["Gray"], COLOR_BGR2GRAY);
		cornerHarris(result["Gray"], harris, 5, 3, 0.04);

		//归一化处理
		normalize(harris, norm, 0, 255, NORM_MINMAX, CV_32FC1);
		convertScaleAbs(norm, result["8U"]);
		//绘制角点
		for (int i = 0; i < norm.rows; i++) 
		{
			for (int j = 0; j < norm.cols; j++) 
			{
				if (norm.at<float>(i, j) > 40) 
				{
					circle(result["Harris"], Point(j, i), 1, Scalar(0, 0, 255));
					circle(result["8U"], Point(j, i), 1, Scalar(0, 0, 255));
				}
			}
		}
	}
	void AngTomasi() 
	{
		result["Tomasi"] = img;
		vector<Point2f> corner;
		//角点检测采用的是灰度图
		goodFeaturesToTrack(result["Gray"], corner, 40, 0.01, 10);
		for (int i = 0; i < corner.size(); i++) 
		{
			circle(result["Tomasi"], corner[i],1, Scalar(0, 255, 0));
		}
	}
	void AngularSubPix() 
	{
		result["SubPix"] = img;
		vector<Point2f> corner;
		TermCriteria tc = TermCriteria(TermCriteria::EPS + TermCriteria::MAX_ITER, 40, 0.001);
		goodFeaturesToTrack(result["Gray"], corner, 40, 0.01, 10);				//先获取角点信息
		cornerSubPix(result["Gray"], corner, Size(5, 5), Size(-1, -1), tc);		//矫正过程
		for (int i = 0; i < corner.size(); i++) 
		{
			circle(result["SubPix"], corner[i], 1, Scalar(255, 0, 0));
		}

	}

	void Show()
	{
		for (auto& v : result)
		{
			imshow(v.first, v.second);
		}
		waitKey(0);
	}

protected:
	Mat img;
	map<string, Mat> result;
};

int main() 
{
	AngularPoint* p = new AngularPoint;
	p->AngularByHarris();
	p->AngTomasi();
	p->AngularSubPix();
	p->Show();
	return 0;
}

opencv图像纹理

纹理是一种反映图像中同质现象的视觉特征，它体现了物体表面的具有缓慢变化或者周期性变化的表面结构组织排列属性，纹理特征不是基于像素点的特征，它需要在包含多个像素点的区域中进行统计计算，通常图像纹理具有以特性:

• 某种局部序列性不断重复
• 非随机排列
• 纹理区域内大致为均匀的统一体
• 具有旋转不变性，且对噪声有较强的抵抗能力

适用于检索具有粗细、疏密等方面较大差别的纹理图像

GLCM描述

灰度共生矩阵法(GLCM， Gray-level co-occurrence matrix)，就是通过计算灰度图像得到它的共生矩阵，然后透过计算该共生矩阵得到矩阵的部分特征值，来分别代表图像的某些纹理特征。灰度共生矩阵能反映图像灰度关于方向、相邻间隔、变化幅度等综合信息，它是分析图像的局部模式和它们排列规则的基础。

计算过程

灰度共生矩阵就是从N×N的图像f(x,y)的灰度为i的像素计算与i距离为δ，灰度为 j的像素同时出现的概率P(i,j,δ,θ)

图示图下:

上述表述可能会比较抽象，接下来我们举一个例子表述 一下： 这是原始的图片像素矩阵，可以看到其最大像素为3，最小像素为0，因此可知最后得出的共生矩阵应为4x4的矩阵(行为像素i,列为像素j,共生矩阵中坐标（i,j）表示像素为i的点走到距离为δ处灰度为 j的像素出现的次数)，注意：不记录重复的情况

0度方向，只能左右走一个单位长度。此时dx=1或dx=-1,dy=0，共生矩阵[0][0]表示像素值为 0的点左或者右移动一格像素值仍然为1的次数(矩阵中[0][0]向右是[0][1]像素为 0记一次，但是 [0][1]向左到[0][0]像素为0就不计数了)

常用的还有45度方向，90度方向和135度方向检索方式，这里不做例举。

GLCM特征

GLCM纹理特征主要有以下四点，数学计算方式如下图:

• 熵：灰度级分布随机性

• 纹理少，熵值接近为零
• 细小纹理多，熵值最大

• 能量：均匀性和平滑性
• 对比度：图像点对中前后点间灰度差的度量
• 均匀度：图像的均匀程度

GLCM实现

#include <iostream>
#include <map>
#include <iomanip>
#include <cmath>
#include <algorithm>
#include <numeric>
#include <opencv2/opencv.hpp>
using namespace std;
using namespace cv;

//1.描述GLCM特征
struct GLCMFeatures 
{
	double energy = 0.0;			//能量
	double entropy = 0.0;			//熵
	double contract = 0.0;			//对比度
	double idMoment = 0.0;			//均匀度
	void PrintFeatures() 
	{
		cout << "energy:" << energy << endl;
		cout << "entropy:" << entropy << endl;
		cout << "contract:" << contract << endl;
		cout << "idMoment:" << idMoment << endl;

	}
};
class GLCM
{
public:
	enum class Type
	{
		GLCM_ANGLE_0 = 0,				//水平
		GLCM_ANGLE_45 = 1,			//斜角
		GLCM_ANGLE_90 = 2,			//垂直
		GLCM_ANGLE_135 = 3
	};
	GLCM(string url = "mm.jpg", int grade = 256) :img(imread(url, IMREAD_GRAYSCALE)), grade(grade)
	{
		//便于采用下标法填充矩阵
		vecGLCM.resize(grade);
		for (int i = 0; i < grade; i++)
		{
			vecGLCM[i].resize(grade);
			for (int j = 0; j < grade; j++)
			{
				vecGLCM[i][j] = 0;
			}
		}
	}
	void GLCM_0() 
	{
		uchar* p = nullptr;
		int row = 0;
		int cols = 0;
		for (int i = 0; i < img.rows; i++)
		{
			p = img.ptr<uchar>(i);
			for (int j = 0; j < img.cols-1; j++) 
			{
				row = p[j];					//第一个点
				cols = p[j + 1];			//第二个点
				vecGLCM[row][cols]++;
			}
		}
	}
	void GLCM_45()
	{
		uchar* p = nullptr;
		int row = 0;
		int cols = 0;
		for (int i = 0; i < img.rows-1; i++)
		{
			p = img.ptr<uchar>(i);
			for (int j = 0; j < img.cols - 1; j++)
			{
				row = p[j];					//第一个点
				cols = (p+1)[j + 1];			//第二个点
				vecGLCM[row][cols]++;
			}
		}
	}
	void GLCM_90()
	{
		uchar* p = nullptr;
		int row = 0;
		int cols = 0;
		for (int i = 0; i < img.rows - 1; i++)
		{
			p = img.ptr<uchar>(i);
			for (int j = 0; j < img.cols; j++)
			{
				row = p[j];					//第一个点
				cols = (p + 1)[j];			//第二个点
				vecGLCM[row][cols]++;
			}
		}
	}
	void GLCM_135()
	{
		uchar* p = nullptr;
		int row = 0;
		int cols = 0;
		for (int i = 0; i < img.rows - 1; i++)
		{
			p = img.ptr<uchar>(i);
			for (int j = 1; j < img.cols; j++)
			{
				row = p[j];						//第一个点
				cols = (p + 1)[j-1];			//第二个点
				vecGLCM[row][cols]++;
			}
		}
	}
	//求四个数学表达式
	void  GetGLCMFeature()
	{
		int total = 0;
		for (auto v : vecGLCM) 
		{
			total += accumulate(v.begin(), v.end(), 0);
		}
		vector<vector<double>> temp(vecGLCM.size());
		for (auto& v : temp) 
		{
			v.resize(vecGLCM.size());
		}
		//归一化处理
		for (int i = 0; i < vecGLCM.size(); i++) 
		{
			for (int j = 0; j < vecGLCM[i].size(); j++) 
			{
				temp[i][j] = double(vecGLCM[i][j]) / double(total);
			}
		}
		for (int i = 0; i < temp.size(); i++) 
		{
			for (int j = 0; j < temp[i].size(); j++) 
			{
				features.energy += temp[i][j] * temp[i][j];   
				if (temp[i][j] > 0) 
				{
					features.entropy -= temp[i][j] * log(temp[i][j]);
				}
				features.contract += double(i - j * 1.0) * double(i - j * 1.0) * temp[i][j];
				features.idMoment += temp[i][j] / (1 + double(i - j * 1.0)* double(i - j * 1.0));
			}
		}
	}
	//封装打印
	void  GetGLCM(GLCM::Type type = GLCM::Type::GLCM_ANGLE_0) 
	{
		switch (type) 
		{
		case Type::GLCM_ANGLE_0:
			GLCM_0();
			break;
		case Type::GLCM_ANGLE_45:
			GLCM_45();
			break;
		case Type::GLCM_ANGLE_90:
			GLCM_90();
			break;
		case Type::GLCM_ANGLE_135:
			GLCM_135();
			break;
		}
		GetGLCMFeature();
		features.PrintFeatures();
	}

protected:
	vector<vector<uchar>> vecGLCM;
	GLCMFeatures features;
	Mat img;
	int grade;
};
int main() 
{
	unique_ptr<GLCM> p(new GLCM);
	p->GetGLCM(GLCM::Type::GLCM_ANGLE_0);
	return 0;
}