opencv学习_13 (trajkovic 角点检测)

来源:http://blog.csdn.net/songzitea/article/details/13614977

背景引言

本节主要内容来源于是由 Miroslav Trajkovic和Mark Hedley[1]在1998年提出Trajkovic算子，其论文为FastCorner Detection.和Trajkovic Operator (4-Neighbours)[3](注：本节图片主要是来源于此)。Trajkovic算子角点提取方法存在价值是在角点提取的效果上，它优于同时期的其他角点(如Moravec角点，Harris角点)提取方法，同时，从算法的运行速度角度而言，它比同是期的其他角点提取方法要快很多。

基本理论

作者对角点的定义，与Moravec角点、 Harris 角点的定义是一样的：图像灰度值在各个方向变化都比较大的点，即认为是角点。角点量的定义也类似于Moravec,即角点量的是在各个方向上灰度变化的最小值（当在某个方向上，灰度值的变化最小，并且这个最小值也大于某某个设定的阈值，那么认为这个点就是一个角点）。与Moravec 相比，Trajkovic 的性能更好，因为Trajkovic通过使用像素插值的方法比较每一个方向上的灰度值变化程度，而moravec只有限个方向的比较，与harris相比，两都效果差不多，但是Trajkovic的速度远快于haris.

Trajkovic在速度方面的表现于决定它能够应用对实时要求比较高的系统中。但是Trajkovic也有缺点是它不具备旋转不变性，对噪声比较敏感，对对角点的响应值比较大等。对于这些缺点也会有相应的方法，但不能完全的解决这些缺点。Trajkovic算子计算角点量是在一个小的圆形窗口内，并且考虑了所有通过圆心的直线，即计算了每一个可能的方向上的灰度值变化程度）。假设用C表示圆的圆心，那么通过圆心C的任意一条直线会与圆在两个交点，如图1所示：

Figure 1: Notation for Trajkovic Operator

那么，Trajkovic对角点量的定义是

C(x,y)表示的是图像I上，任意一点(x,y)处的角点量。这个公式能断出是否是角点。点(x,y)在图像目的位置有四种情况。如下图所示：

第一情况：点(x,y)在内部区域：此时经点(x,y)位中心的圆形窗中完全在一个比较平坦的区域，或者大部分的圆形窗口是在这个平坦区域的（如上图A所示）。这种情况下，至少存大一条直线过中心点C的直线使得Ic约等于Ip和Ic约等于Ip‘,根据角点量计算公式，计算出的C(x,y)一定是一个很小的值，那么就可以排除点(x,y)不是角点。并且，大多数的情况下，会有多个这样的点P存大，使得对噪声一定的承受能力。

第二情况：点(x,y)在边缘上：此时圆形窗口的中心点(x,y)刚好落在边缘上（如上图B所示）。仅存在一个点对P,P'使得Ic约等于Ip和Ic约等于Ip‘。在理想状况下，根据角点量的公式计算出来的角点量C(x,y)是比较小的。但是由于只存在一个条这样直线，所以容易受到噪声的干扰。

第三情况：点(x,y)在角点上：当圆窗口的中心(x, y)位于角点上时（如上图C 所示）。任意一条通过点(x,y)与圆相交的交点P、P', P、P' 中至少有一个点的灰度值与中心的点的灰度值差别比较大。所以计算出来的角点量也比较大（大于设定的阈值，所以会被认为是角点）。

第四情况：点(x,y)在一个孤立点上：如果在孤立点上（如上图D所示），对于每一个直线上的Ip, Ip’,都分与Ic的值相差比较大，此时计算出的角点量也比较大。此时这个点是噪声,但是由于角点量比较大，会被认为是角点。为了削除这种噪声的影响，可以先使用高斯平滑进行去噪。

根据上面的分析，Trajkovic 受噪声的影响比较大，所以，可以先实现高斯平滑去噪再使用Trajkovic算子进行角点检测。

角点量的计算方法

目前，关于Trajkovic计算角点量的并未具体介绍如何时计算。首先，如下图2所示，水平方向的灰度值变化量rA,垂直方向的灰度值变化量rB很容易计算。

如何计算任意一个方向上灰度值的变化程度呢？首先，我们角点量就取其最小值即可：

Figure 2: Interpixel positions Figure3: Interpixel approximation for a 3x3 window using 4-neighbours

现在，我们可将问题可以转换为，任意一条通过点C 的直线与圆的交点是P，P',Q,Q'即可：

把公式(2)代入公式(1) 中即可：

那么，最小值是即可结论是

多格算法（Multigrid Algorithm）

角点可以分为两类：几何角点和纹理角点。几何角点是由图像中物体边缘的相交。而纹理角点有由物体表面的纹理产生的角点（如：草地，衣服的纹理等）。通常情况下，一幅图像中的几何角点的数量要远少于纹理角点。多格算法的目的是希望能多检测到几何角点，少检测到纹理点。因Trajkovic 认为几何角点比纹理角点更加的稳定。所以，从实践使用来看，减少纹理角点是合理的。

经过观察发现，纹理角点一般都是非常的密集，并且是在一个很小区域内灰度值发生变化，所以，采用把原图缩小以后，再来提取角点，通过缩小图像，可以消除区域内灰度值的变换。缩小图像的方法不是使用插值法，而是使用平均法。如下图所示：

Trajkovic角点检测，首先使用多格法对原如图像进行缩小，这样不仅能够减少检测到纹理角点数目，而且加快了Trajkovic角点检测的速度，在Trajkovic角点检测缩小版的图像上，使用角为简单的角点量计算公式初步判断是不是侯选角点。如果是，再到原始图像上进一步判定。

Figure 4: Corner points detected at different image resolutions

如图4所示，原始的大小是256*256,图像中有一片草地，而草地上的角点，对于我们进行角点匹配时，是没有什么意义。当把图像缩小到128*128时，草地上的角点已经消除了。当然，图像也不能缩小太多，否则，几何角点也会消失，当缩小到64*64时，部分几何角点也消失了。

算法流程

代码：

[cpp] view plain copy print ?

#include <iostream>
#include "cv.h"
#include "highgui.h"
#include "cxcore.h"
using namespace std;
void getTrajkovic4NCorner(IplImage* src,IplImage* srcResize , float T1 , float T2 ,CvSeq* corners,int maximumSize)
{
int x,y,maxChar=255,scaleX = src->width/srcResize->width,scaleY = src->height/srcResize->height;
IplImage* srcResizeMap ,*srcMap;
srcResizeMap = cvCreateImage(cvGetSize(srcResize),8,1); // 用来保存缩小图像的角点量
srcMap = cvCreateImage(cvGetSize(src),32,1); // 用来保存源图像计算得到的角点量
cvZero(srcMap);
cvZero(srcResizeMap);
for( y=1;y<srcResize->height-1;y++)
{
uchar* preRow = (uchar*)(srcResize->imageData + (y-1)*srcResize->widthStep);
uchar* curRow = (uchar*)(srcResize->imageData + y*srcResize->widthStep);
uchar* nextRow = (uchar*)(srcResize->imageData + (y+1)*srcResize->widthStep);
uchar* MapData = (uchar*)(srcResizeMap->imageData + y*srcResizeMap->widthStep);
for(x=1;x<srcResize->width-1;x++)
{
int IC,IA,IB,IAA,IBB;
int rA,rB,C_Simple;
// curRow[x]是uchar型,而IC是int型,如果直接赋值int类型的其它三个字节可能会产生随机值
IC = curRow[x]&maxChar;
IA = curRow[x+1]&maxChar;
IAA = curRow[x-1]&maxChar;
IB = preRow[x]&maxChar;
IBB = nextRow[x]&maxChar;
rA = (IA-IC)*(IA-IC) + (IAA-IC)*(IAA-IC);
rB = (IB-IC)*(IB-IC) + (IBB-IC)*(IBB-IC);
C_Simple = rA < rB ? rA : rB;
if(C_Simple > T1)
{
MapData[x]=1;
}
}
}
for( y=1;y<srcResize->height-1;y++)
{
uchar* srcResizeMapData = (uchar*)(srcResizeMap->imageData + y*srcResizeMap->widthStep);
//uchar* srcMapData = (uchar*)(srcMap->imageData + y*srcMap->widthStep);
for(x=1;x<srcResize->width-1;x++)
{
if(srcResizeMapData[x]==0)
{
continue;
}
int originX,originY,IC,IA,IB,IAA,IBB;
int rA,rB,C_Simple;
originX = x*scaleX;
originY = y*scaleY;
uchar* srcPreRow = (uchar*)(src->imageData + (originY-1)*src->widthStep);
uchar* srcCurRow = (uchar*)(src->imageData + originY*src->widthStep);
uchar* srcNextRow = (uchar*)(src->imageData + (originY+1)*src->widthStep);
IC = srcCurRow[originX]&maxChar;
IA = srcCurRow[originX+1]&maxChar;
IAA = srcCurRow[originX-1]&maxChar;
IB = srcPreRow[x]&maxChar;
IBB = srcNextRow[x]&maxChar;
rA = (IA-IC)*(IA-IC) + (IAA-IC)*(IAA-IC);
rB = (IB-IC)*(IB-IC) + (IBB-IC)*(IBB-IC);
C_Simple = rA < rB ? rA : rB;
if(C_Simple>T2)
{
float B1,B2,C,A,B,C_InterPixel;
B1 = (IB-IA)*(IA-IC)+(IBB-IAA)*(IAA-IC);
B2 = (IB-IAA)*(IAA-IC)+(IBB-IA)*(IA-IC);
C = rA;
B = B1<B2 ? B1 : B2;
A = rB-rA-2*B;
if(B<0 &&(A+B)>0)
{
C_InterPixel = C-(B*B)/A;
}
else
{
C_InterPixel=C_Simple;
}
if(C_InterPixel>T2)
{
float* srcMapData = (float*)(srcMap->imageData + originY*srcMap->widthStep);
srcMapData[originX]=C_InterPixel;
//cvSetReal2D(srcMap,originY,originX,C_InterPixel);
}
}
}
}
//计算局部极大值及极大值是否大于阈值
int beginY,endY,beginX,endX;
int halfWinSize = maximumSize/2;
beginY = halfWinSize;
endY = srcMap->height - halfWinSize;
beginX = halfWinSize;
endX = srcMap->width - halfWinSize;
for(y=beginY;y<endY;)
{
for(x=beginX;x<endX;)
{
//寻找局部极大值及其位置信息
float maxValue=0;
int flag = 0 ;
CvPoint maxLoc;
maxLoc.x = -1;
maxLoc.y = -1;
//首先计算以点(x,y)位中心的maximumSize*maximumSize的窗口内部的局部极大值
for(int winy=-halfWinSize;winy<=halfWinSize;winy++)
{
for(int winx=-halfWinSize;winx<=halfWinSize;winx++)
{
float value ;
value = cvGetReal2D(srcMap,y+winy,x+winx);
//计算该窗口内最大值保存到max 并保存其坐标到maxLoc
if(value>maxValue)
{
maxValue = value;
maxLoc.x = x+winx;
maxLoc.y = y+winy;
flag = 1;
}
}
}
//如果找到局部极大值并且该值大于预先设定的阈值则认为是角点
if(flag==1 && maxValue>T2)
{
cvSeqPush(corners,&maxLoc);
}
x = x+halfWinSize;
}
y = y + halfWinSize;
}
cvReleaseImage(&srcResizeMap);
cvReleaseImage(&srcMap);
}
int main(int argc,char* argv[])
{
//相关变量
int scale = 2;
IplImage* src,*srcGray,*srcGrayResize;
CvMemStorage* mem = cvCreateMemStorage(0);
CvSeq* TrajkovicPoints;
src = cvLoadImage("E:\\study_opencv_video\\lesson17_2\\2.jpg");//源图像
srcGray = cvCreateImage(cvGetSize(src),8,1);
if(!src)
{
cout << "src is null" << endl;
return 0;
}
cvCvtColor(src,srcGray,CV_BGR2GRAY);
srcGrayResize = cvCreateImage(cvSize(srcGray->width/scale,srcGray->height/scale),8,1);
cvResize(srcGray,srcGrayResize);//请将resize修改为多格算法
//Trajkovic 4 邻域角点角点保存的空间角点坐标保存在一个序列中
TrajkovicPoints = cvCreateSeq(0,sizeof(CvSeq),sizeof(CvPoint),mem);
//进行Trajkovic 4 邻域角点检测
int T1=120,T2=150,localArea = 8; // localArea 是局部极大值抑制时窗口的大小
getTrajkovic4NCorner(srcGray,srcGrayResize,T1,T2,TrajkovicPoints,localArea);
//获取每一个角点的坐标
for(int x=0;x<TrajkovicPoints->total;x++)
{
//获取第x个角点的坐标
CvPoint* pt = (CvPoint*)cvGetSeqElem(TrajkovicPoints,x);
//以角点坐标为中心绘制一个半径为5的圆
cvCircle(src,*pt,2,cvScalar(255,0,255,0));
}
//显示图像
cvNamedWindow("dst");
cvShowImage("dst",src);
cvWaitKey(0);
//释放资源
cvReleaseImage(&src);
cvReleaseImage(&srcGray);
cvReleaseImage(&srcGrayResize);
cvReleaseMemStorage(&mem);
return 0;
}

英文链接： http://kiwi.cs.dal.ca/~dparks/CornerDetection/trajkovic.htm