图像特征点|SIFT特征点位置与方向
本次主要为小伙伴们讲解，如何求取关键点的位置和方向。
空间极值点（即关键点）检测
关键点是由DOG空间的局部极值点组成的，关键点的初步探查是通过同一组内各DoG相邻两层图像之间比较完成的。为了寻找DoG函数的极值点，每一个像素点要和它所有的相邻点比较，看其是否比它的图像域和尺度域的相邻点大或者小。如图下图所示，中间的检测点和它同尺度的8个相邻点和上下相邻尺度对应的9×2个点共26个点比较，以确保在尺度空间和二维图像空间都检测到极值点。
2.1、极值点检测过程
2.1.1、极值点检测示意
2.1.2、极值点检测源码分析
if (val > 0)//极大值检测{
for (i = -1; i <= 1; i++)
for (j = -1; j <= 1; j++)
for (k = -1; k <= 1; k++)
if (val < pixval32f(dog_pyr[octv][intvl + i], r + j, c + k))//pixval32f为提取图像像素位置上的灰度值
return 0;}
else /* check for minimum */
{
for (i = -1; i <= 1; i++)
for (j = -1; j <= 1; j++)
for (k = -1; k <= 1; k++)
if (val > pixval32f(dog_pyr[octv][intvl + i], r + j, c + k))//r c为图像的行数和列数，dog_pyr为高斯差分图
return 0;2.2、关键点定位
以上方法检测到的极值点是离散空间的极值点，以下通过拟合三维二次函数来精确确定关键点的位置和尺度，同时去除低对比度的关键点和不稳定的边缘响应点(因为DoG算子会产生较强的边缘响应)，以增强匹配稳定性、提高抗噪声能力。
2.2.1、关键点精确定位
离散空间的极值点并不是真正的极值点，下图显示了二维函数离散空间得到的极值点与连续空间极值点的差别。利用已知的离散空间点插值得到的连续空间极值点的方法叫做子像素插值。
为了提高关键点的稳定性，需要对尺度空间DoG函数进行曲线插值。利用DoG函数在尺度空间的Taylor展开式(插值函数)为：
上面算式的矩阵表示如下：
其中，X求导并让方程等于零，可以得到极值点的偏移量为：
对应极值点，方程的值为：
其中, X^代表相对插值中心的偏移量，当它在任一维度上的偏移量大于0.5时（即x或y或 σ），意味着插值中心已经偏移到它的邻近点上，所以必须改变当前关键点的位置。同时在新的位置上反复插值直到收敛；也有可能超出所设定的迭代次数或者超出图像边界的范围，此时这样的点应该删除，在Lowe中进行了5次迭代。另外，过小的点易受噪声的干扰而变得不稳定，所以将 小于某个经验值(Lowe论文中使用0.03，Rob Hess等人实现时使用0.04/S)的极值点删除。同时，在此过程中获取特征点的精确位置(原位置加上拟合的偏移量)以及尺度(σ)。
2.2.2、消除边缘响应
一个定义不好的高斯差分算子的极值在横跨边缘的地方有较大的主曲率，而在垂直边缘的方向有较小的主曲率。DOG算子会产生较强的边缘响应，需要剔除不稳定的边缘响应点。获取特征点处的Hessian矩阵，主曲率通过一个2x2 的Hessian矩阵H求出（D的主曲率和H的特征值成正比）：
假设H的特征值为α和β（α、β代表x和y方向的梯度）且α>β。令α=rβ则有：
其中Tr(H)求取H的对角元素和；Det(H)为求H的行列式值。
则公式(r+1)^2/r的值在两个特征值相等时最小，随着的增大而增大。值越大，说明两个特征值的比值越大，即在某一个方向的梯度值越大，而在另一个方向的梯度值越小，而边缘恰恰就是这种情况。所以为了剔除边缘响应点，需要让该比值小于一定的阈值，因此，为了检测主曲率是否在某域值r下，只需检测:
论文建议r=10,OpenCv也采用r=10
2.2.3、精确定位中的泰勒插值源码分析
while (i < SIFT_MAX_INTERP_STEPS)//SIFT_MAX_INTERP_STEPS=5为最大迭代次数，避免长时迭代
{
interp_step(dog_pyr, octv, intvl, r, c, &xi, &xr, &xc);// 泰勒展开拟合，xi,xr,xc依次为x、y、σ方向偏移量,
if (ABS(xi) < 0.5 && ABS(xr) < 0.5 && ABS(xc) < 0.5)//如果当前偏移量绝对值中的每个值均小于0.5，退出迭代
break;
c += cvRound(xc);//计算行坐标，cvRound 为四舍五入。
r += cvRound(xr);
intvl += cvRound(xi);
if (intvl < 1 ||//不在计算的图像层中
intvl > intvls ||//高斯差分每组的层数为intvls
c < SIFT_IMG_BORDER ||//靠近图像边缘5个像素的区域不做检测，SIFT_IMG_BORDER=5,
r < SIFT_IMG_BORDER ||
c >= dog_pyr[octv][0]->width - SIFT_IMG_BORDER ||//靠近图像边缘5个像素的区域不做检测
r >= dog_pyr[octv][0]->height - SIFT_IMG_BORDER)
{
return NULL;
}
i++;//迭代计数
}static void interp_step(IplImage*** dog_pyr, int octv, int intvl, int r, int c,double* xi, double* xr, double* xc)
{
CvMat* dD, *H, *H_inv, X;
double x[3] = { 0 };
dD = deriv_3D(dog_pyr, octv, intvl, r, c);//一阶偏导数
H = hessian_3D(dog_pyr, octv, intvl, r, c);//Hessian 矩阵即二阶导数组成的矩阵
H_inv = cvCreateMat(3, 3, CV_64FC1);
cvInvert(H, H_inv, CV_SVD);//求Hessian矩阵的逆矩阵
cvInitMatHeader(&X, 3, 1, CV_64FC1, x, CV_AUTOSTEP);
cvGEMM(H_inv, dD, -1, NULL, 0, &X, 0); //cvGEMM为矩阵乘法，//第一个矩阵的系数；//H_inv、dD第一二个矩阵//-1矩阵前的常数//X结果矩阵
cvReleaseMat(&dD);
cvReleaseMat(&H);
cvReleaseMat(&H_inv);
*xi = x[2];
*xr = x[1];
*xc = x[0];
}3、关键点方向分配
为了使描述符具有旋转不变性，需要利用图像的局部特征为给每一个关键点分配一个基准方向。使用图像梯度的方法求取局部结构的稳定方向。
3.1、特征点的梯度
3.1.1、梯度的计算
对于在DOG金字塔中检测出的关键点点，采集其所在高斯金字塔图像3σ领域窗口内像素的梯度和方向分布特征。梯度的模值和方向如下：
L为关键点所在的尺度空间值，按Lowe的建议，梯度的模值m(x,y)按 σ=1.5σ_oct 的高斯分布加成，按尺度采样的3σ原则，领域窗口半径为 3x1.5σ_oct。
3.1.1、梯度直方图
在完成关键点的梯度计算后，使用直方图统计领域内像素的梯度和方向。梯度直方图将0~360度的方向范围分为36个柱(bins)，其中每柱10度。如图5.1所示，直方图的峰值方向代表了关键点的主方向，(为简化，图中只画了八个方向的直方图)。
3.2、特征点主方向的确定
方向直方图的峰值则代表了该特征点处邻域梯度的方向，以直方图中最大值作为该关键点的主方向。为了增强匹配的鲁棒性，只保留峰值大于主方向峰值80％的方向作为该关键点的辅方向。因此，对于同一梯度值的多个峰值的关键点位置，在相同位置和尺度将会有多个关键点被创建但方向不同。仅有15％的关键点被赋予多个方向，但可以明显的提高关键点匹配的稳定性。实际编程实现中，就是把该关键点复制成多份关键点，并将方向值分别赋给这些复制后的关键点，并且，离散的梯度方向直方图要进行插值拟合处理，来求得更精确的方向角度值。
3.2.1、梯度图像的平滑处理
为了防止某个梯度方向角度因受到噪声的干扰而突变，我们还需要对梯度方向直方图进行平滑处理。Opencv 所使用的平滑公式为：
其中i∈[0,35],h 和H 分别表示平滑前和平滑后的直方图。由于角度是循环的，即00=3600，如果出现h(j)，j超出了(0,…,35)的范围，那么可以通过圆周循环的方法找到它所对应的、在00=3600之间的值，如h(-1) = h(35)。
3.2.2、梯度直方图抛物线插值
假设我们在第i个小柱子要找一个精确的方向，那么由上面分析知道：
设插值抛物线方程为h(t)=at^2+bt+c,其中a、b、c为抛物线的系数，t为自变量，t∈[-1,1],此抛物线求导并令它等于0。即h(t)´=0 得tmax=-b/(2a)
现在把这三个插值点带入方程可得：
3.2.3、抛物线插值源码分析
#define interp_hist_peak( l, c, r ) ( 0.5 * ((l)-(r)) / ((l) - 2.0*(c) + (r)) )//插值计算式，l为左侧柱子值，r为左侧柱子值
static void add_good_ori_features(CvSeq* features, double* hist, int n,
double mag_thr, struct feature* feat)//精确主方向及辅方向
{
struct feature* new_feat;
double bin, PI2 = CV_PI * 2.0;//CV_PI=pi
int l, r, i;
for (i = 0; i < n; i++)// 直方图有n=36个小柱子
{
l = (i == 0) ? n - 1 : i - 1;//把小柱子看成是循环的，角度的取值为0-360即一个圆周
r = (i + 1) % n;
//只对小柱子的值大于等于主峰80%且此小柱子比左右两边小柱子都高的柱子进行抛物线插值
if (hist[i] > hist[l] && hist[i] > hist[r] && hist[i] >= mag_thr)// mag_thr为>=80%的最高峰值
{
bin = i + interp_hist_peak(hist[l], hist[i], hist[r]);//interp_hist_peak 插值函数
bin = (bin < 0) ? n + bin : (bin >= n) ? bin - n : bin;//角度取值约束在0-360之间，且是连续循环的
new_feat = clone_feature(feat);//幅值特征点
new_feat->ori = ((PI2 * bin) / n) - CV_PI;//？
cvSeqPush(features, new_feat);
free(new_feat);
}至此，图像的关键点已检测完毕，每个关键点有三个信息：位置、所处尺度、方向。由此可以确定一个SIFT特征区域。
本文由作者首发在CSDN，授权转载，禁止二次转载 https://blog.csdn.net/lingyunxianhe/article/details/79063547