#include “stdafx.h”#include#include #include #include #include #include #include using namespace std;IplImage *image = 0;double compareSURFDescriptors( const float* d1, const float* d2, double best, int length ){double total_cost = 0;assert( length % 4 == 0 );for( int i = 0; i < length; i += 4 ){double t0 = d1[i] – d2[i];double t1 = d1[i+1] – d2[i+1];double t2 = d1[i+2] – d2[i+2];double t3 = d1[i+3] – d2[i+3];total_cost += t0*t0 + t1*t1 + t2*t2 + t3*t3;if( total_cost > best )break;}return total_cost;}int naiveNearestNeighbor( const float* vec, int laplacian,const CvSeq* model_keypoints,const CvSeq* model_descriptors ){int length = (int)(model_descriptors->elem_size/sizeof(float));int i, neighbor = -1;double d, dist1 = 1e6, dist2 = 1e6;CvSeqReader reader, kreader;cvStartReadSeq( model_keypoints, &kreader, 0 );cvStartReadSeq( model_descriptors, &reader, 0 );for( i = 0; i < model_descriptors->total; i++ ){const CvSURFPoint* kp = (const CvSURFPoint*)kreader.ptr;const float* mvec = (const float*)reader.ptr;CV_NEXT_SEQ_ELEM( kreader.seq->elem_size, kreader );CV_NEXT_SEQ_ELEM( reader.seq->elem_size, reader );if( laplacian != kp->laplacian )continue;d = compareSURFDescriptors( vec, mvec, dist2, length );if( d < dist1 ){dist2 = dist1;dist1 = d;neighbor = i;}else if ( d < dist2 )dist2 = d;}if ( dist1 < 0.6*dist2 )return neighbor;return -1;}//用于找到两幅图像之间匹配的点对，并把匹配的点对存储在 ptpairs 向量中，其中物体(object)图像的特征点//及其相应的描述器（局部特征）分别存储在 objectKeypoints 和 objectDescriptors，场景(image)图像的特//征点及其相应的描述器（局部特征）分别存储在 imageKeypoints和 imageDescriptorsvoid findPairs( const CvSeq* objectKeypoints, const CvSeq* objectDescriptors,const CvSeq* imageKeypoints, const CvSeq* imageDescriptors, vector & ptpairs ){int i;CvSeqReader reader, kreader;cvStartReadSeq( objectKeypoints, &kreader );cvStartReadSeq( objectDescriptors, &reader );ptpairs.clear();for( i = 0; i < objectDescriptors->total; i++ ){const CvSURFPoint* kp = (const CvSURFPoint*)kreader.ptr;const float* descriptor = (const float*)reader.ptr;CV_NEXT_SEQ_ELEM( kreader.seq->elem_size, kreader );CV_NEXT_SEQ_ELEM( reader.seq->elem_size, reader );int nearest_neighbor = naiveNearestNeighbor( descriptor, kp->laplacian, imageKeypoints, imageDescriptors );if( nearest_neighbor >= 0 ){ptpairs.push_back(i);ptpairs.push_back(nearest_neighbor);}}}//用于寻找物体(object)在场景(image)中的位置,位置信息保存在参数dst_corners中，参数src_corners由物//体(object的width几height等决定，其他部分参数如上findPairsint locatePlanarObject( const CvSeq* objectKeypoints, const CvSeq* objectDescriptors,const CvSeq* imageKeypoints, const CvSeq* imageDescriptors,const CvPoint src_corners[4], CvPoint dst_corners[4] ){double h[9];CvMat _h = cvMat(3, 3, CV_64F, h);vector ptpairs;vector pt1, pt2;CvMat _pt1, _pt2;int i, n;findPairs( objectKeypoints, objectDescriptors, imageKeypoints, imageDescriptors, ptpairs );n = ptpairs.size()/2;if( n < 4 )return 0;pt1.resize(n);pt2.resize(n);for( i = 0; i < n; i++ ){pt1[i] = ((CvSURFPoint*)cvGetSeqElem(objectKeypoints,ptpairs[i*2]))->pt;pt2[i] = ((CvSURFPoint*)cvGetSeqElem(imageKeypoints,ptpairs[i*2+1]))->pt;}_pt1 = cvMat(1, n, CV_32FC2, &pt1[0] );_pt2 = cvMat(1, n, CV_32FC2, &pt2[0] );if( !cvFindHomography( &_pt1, &_pt2, &_h, CV_RANSAC, 5 ))return 0;for( i = 0; i < 4; i++ ){double x = src_corners[i].x, y = src_corners[i].y;double Z = 1./(h[6]*x + h[7]*y + h[8]);double X = (h[0]*x + h[1]*y + h[2])*Z;double Y = (h[3]*x + h[4]*y + h[5])*Z;dst_corners[i] = cvPoint(cvRound(X), cvRound(Y));}return 1;}int main(int argc, char** argv){//物体(object)和场景(scene)的图像向来源const char* object_filename = argc == 3 ? argv[1] : “box.png”;const char* scene_filename = argc == 3 ? argv[2] : “box_in_scene.png”;//内存存储器CvMemStorage* storage = cvCreateMemStorage(0);cvNamedWindow(“Object”, 1);cvNamedWindow(“Object Correspond”, 1);//颜色值static CvScalar colors[] ={{ {0,0,255}},{ {0,128,255}},{ {0,255,255}},{ {0,255,0}},{ {255,128,0}},{ {255,255,0}},{ {255,0,0}},{ {255,0,255}},{ {255,255,255}}};IplImage* object = cvLoadImage( object_filename, CV_LOAD_IMAGE_GRAYSCALE );IplImage* image = cvLoadImage( scene_filename, CV_LOAD_IMAGE_GRAYSCALE );if( !object || !image ){fprintf( stderr, “Can not load %s and/or %s\n”“Usage: find_obj [ ]\n”,object_filename, scene_filename );exit(-1);}IplImage* object_color = cvCreateImage(cvGetSize(object), 8, 3);cvCvtColor( object, object_color, CV_GRAY2BGR );//物体(object)和场景(scene)的图像的特征点CvSeq *objectKeypoints = 0, *objectDescriptors = 0;CvSeq *imageKeypoints = 0, *imageDescriptors = 0;int i;//定义Surf算法要用的参数分别为 threshold 和 extendedCvSURFParams params = cvSURFParams(500, 1);double tt = (double)cvGetTickCount();//提取物体(object)和场景(scene)的图像的特征点及其描述器cvExtractSURF( object, 0, &objectKeypoints, &objectDescriptors, storage, params );printf(“Object Descriptors: %d\n”, objectDescriptors->total);cvExtractSURF( image, 0, &imageKeypoints, &imageDescriptors, storage, params );printf(“Image Descriptors: %d\n”, imageDescriptors->total);//计算所消耗的时间tt = (double)cvGetTickCount() – tt;printf( “Extraction time = %gms\n”, tt/(cvGetTickFrequency()*1000.));CvPoint src_corners[4] = { {0,0}, {object->width,0}, {object->width, object->height}, {0, object->height}};//定义感兴趣的区域CvPoint dst_corners[4];IplImage* correspond = cvCreateImage( cvSize(image->width, object->height+image->height), 8, 1 );//设置感兴趣区域cvSetImageROI( correspond, cvRect( 0, 0, object->width, object->height ) );cvCopy( object, correspond );cvSetImageROI( correspond, cvRect( 0, object->height, correspond->width, correspond->height ) );cvCopy( image, correspond );cvResetImageROI( correspond );//寻找物体(object)在场景(image)中的位置，并将信息保存if( locatePlanarObject( objectKeypoints, objectDescriptors, imageKeypoints,imageDescriptors, src_corners, dst_corners )){for( i = 0; i < 4; i++ ){CvPoint r1 = dst_corners[i%4];CvPoint r2 = dst_corners[(i+1)%4];cvLine( correspond, cvPoint(r1.x, r1.y+object->height ),cvPoint(r2.x, r2.y+object->height ), colors[8] );}}//定义并保存物体(object)在场景(image)图形之间的匹配点对，并将其存储在向量 ptpairs 中，之后可以对//ptpairs 进行操作vector ptpairs;findPairs( objectKeypoints, objectDescriptors, imageKeypoints, imageDescriptors, ptpairs );//显示匹配结果for( i = 0; i < (int)ptpairs.size(); i += 2 ){CvSURFPoint* r1 = (CvSURFPoint*)cvGetSeqElem( objectKeypoints, ptpairs[i] );CvSURFPoint* r2 = (CvSURFPoint*)cvGetSeqElem( imageKeypoints, ptpairs[i+1] );cvLine( correspond, cvPointFrom32f(r1->pt),cvPoint(cvRound(r2->pt.x), cvRound(r2->pt.y+object->height)), colors[8] );}cvShowImage( “Object Correspond”, correspond );//显示物体(object)的所有特征点for( i = 0; i < objectKeypoints->total; i++ ){CvSURFPoint* r = (CvSURFPoint*)cvGetSeqElem( objectKeypoints, i );CvPoint center;int radius;center.x = cvRound(r->pt.x);center.y = cvRound(r->pt.y);radius = cvRound(r->size*1.2/9.*2);cvCircle( object_color, center, radius, colors[0], 1, 8, 0 );}cvShowImage( “Object”, object_color );cvWaitKey(0);//释放窗口所占用的内存cvDestroyWindow(“Object”);cvDestroyWindow(“Object Correspond”);return 0;}
surf OpenCV Demo分析（find_obj.cpp）
OpenCV2.1中有关于Surf算法的简单示例(1.1以上的版本都添加了这个算法)，在路径：C:\Program Files\OpenCV2.1\samples\c下，名为find_obj.cpp，运行它可以直接观察到相应结果。为了便于介绍这个示例，简单做了如下修改（只是删掉一些代码，但是对于如何使用Surf算法没有影响）。
修改后的代码及其注释如下：（主要是介绍这个main函数）
#include “stdafx.h”
#include <cv.h>
#include <highgui.h>
#include <ctype.h>
#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include <vector>
using namespace std;
IplImage *image = 0;
double compareSURFDescriptors( const float* d1, const float* d2, double best, int length )
{
double total_cost = 0;
assert( length % 4 == 0 );
for( int i = 0; i < length; i += 4 )
{
double t0 = d1[i] – d2[i];
double t1 = d1[i+1] – d2[i+1];
double t2 = d1[i+2] – d2[i+2];
double t3 = d1[i+3] – d2[i+3];
total_cost += t0*t0 + t1*t1 + t2*t2 + t3*t3;
if( total_cost > best )
break;
}
return total_cost;
}
int naiveNearestNeighbor( const float* vec, int laplacian,
const CvSeq* model_keypoints,
const CvSeq* model_descriptors )
{
int length = (int)(model_descriptors->elem_size/sizeof(float));
int i, neighbor = -1;
double d, dist1 = 1e6, dist2 = 1e6;
CvSeqReader reader, kreader;
cvStartReadSeq( model_keypoints, &kreader, 0 );
cvStartReadSeq( model_descriptors, &reader, 0 );
for( i = 0; i < model_descriptors->total; i++ )
{
const CvSURFPoint* kp = (const CvSURFPoint*)kreader.ptr;
const float* mvec = (const float*)reader.ptr;
CV_NEXT_SEQ_ELEM( kreader.seq->elem_size, kreader );
CV_NEXT_SEQ_ELEM( reader.seq->elem_size, reader );
if( laplacian != kp->laplacian )
continue;
d = compareSURFDescriptors( vec, mvec, dist2, length );
if( d < dist1 )
{
dist2 = dist1;
dist1 = d;
neighbor = i;
}
else if ( d < dist2 )
dist2 = d;
}
if ( dist1 < 0.6*dist2 )
return neighbor;
return -1;
}
//用于找到两幅图像之间匹配的点对，并把匹配的点对存储在 ptpairs 向量中，其中物体(object)图像的特征点
//及其相应的描述器（局部特征）分别存储在 objectKeypoints 和 objectDescriptors，场景(image)图像的特
//征点及其相应的描述器（局部特征）分别存储在 imageKeypoints和 imageDescriptors
void findPairs( const CvSeq* objectKeypoints, const CvSeq* objectDescriptors,
const CvSeq* imageKeypoints, const CvSeq* imageDescriptors, vector<int>& ptpairs )
{
int i;
CvSeqReader reader, kreader;
cvStartReadSeq( objectKeypoints, &kreader );
cvStartReadSeq( objectDescriptors, &reader );
ptpairs.clear();
for( i = 0; i < objectDescriptors->total; i++ )
{
const CvSURFPoint* kp = (const CvSURFPoint*)kreader.ptr;
const float* descriptor = (const float*)reader.ptr;
CV_NEXT_SEQ_ELEM( kreader.seq->elem_size, kreader );
CV_NEXT_SEQ_ELEM( reader.seq->elem_size, reader );
int nearest_neighbor = naiveNearestNeighbor( descriptor, kp->laplacian, imageKeypoints, imageDescriptors );
if( nearest_neighbor >= 0 )
{
ptpairs.push_back(i);
ptpairs.push_back(nearest_neighbor);
}
}
}
//用于寻找物体(object)在场景(image)中的位置,位置信息保存在参数dst_corners中，参数src_corners由物
//体(object的width几height等决定，其他部分参数如上findPairs
int locatePlanarObject( const CvSeq* objectKeypoints, const CvSeq* objectDescriptors,
const CvSeq* imageKeypoints, const CvSeq* imageDescriptors,
const CvPoint src_corners[4], CvPoint dst_corners[4] )
{
double h[9];
CvMat _h = cvMat(3, 3, CV_64F, h);
vector<int> ptpairs;
vector<CvPoint2D32f> pt1, pt2;
CvMat _pt1, _pt2;
int i, n;
findPairs( objectKeypoints, objectDescriptors, imageKeypoints, imageDescriptors, ptpairs );
n = ptpairs.size()/2;
if( n < 4 )
return 0;
pt1.resize(n);
pt2.resize(n);
for( i = 0; i < n; i++ )
{
pt1[i] = ((CvSURFPoint*)cvGetSeqElem(objectKeypoints,ptpairs[i*2]))->pt;
pt2[i] = ((CvSURFPoint*)cvGetSeqElem(imageKeypoints,ptpairs[i*2+1]))->pt;
}
_pt1 = cvMat(1, n, CV_32FC2, &pt1[0] );
_pt2 = cvMat(1, n, CV_32FC2, &pt2[0] );
if( !cvFindHomography( &_pt1, &_pt2, &_h, CV_RANSAC, 5 ))
return 0;
for( i = 0; i < 4; i++ )
{
double x = src_corners[i].x, y = src_corners[i].y;
double Z = 1./(h[6]*x + h[7]*y + h[8]);
double X = (h[0]*x + h[1]*y + h[2])*Z;
double Y = (h[3]*x + h[4]*y + h[5])*Z;
dst_corners[i] = cvPoint(cvRound(X), cvRound(Y));
}
return 1;
}
int main(int argc, char** argv)
{
//物体(object)和场景(scene)的图像向来源
const char* object_filename = argc == 3 ? argv[1] : “box.png”;
const char* scene_filename = argc == 3 ? argv[2] : “box_in_scene.png”;
//内存存储器
CvMemStorage* storage = cvCreateMemStorage(0);
cvNamedWindow(“Object”, 1);
cvNamedWindow(“Object Correspond”, 1);
//颜色值
static CvScalar colors[] =
{
{
{0,0,255}},{
{0,128,255}},{
{0,255,255}},{
{0,255,0}},{
{255,128,0}},{
{255,255,0}},{
{255,0,0}},{
{255,0,255}},{
{255,255,255}}};
IplImage* object = cvLoadImage( object_filename, CV_LOAD_IMAGE_GRAYSCALE );
IplImage* image = cvLoadImage( scene_filename, CV_LOAD_IMAGE_GRAYSCALE );
if( !object || !image )
{
fprintf( stderr, “Can not load %s and/or %s\n”
“Usage: find_obj [<object_filename> <scene_filename>]\n”,
object_filename, scene_filename );
exit(-1);
}
IplImage* object_color = cvCreateImage(cvGetSize(object), 8, 3);
cvCvtColor( object, object_color, CV_GRAY2BGR );
//物体(object)和场景(scene)的图像的特征点
CvSeq *objectKeypoints = 0, *objectDescriptors = 0;
CvSeq *imageKeypoints = 0, *imageDescriptors = 0;
int i;
//定义Surf算法要用的参数分别为 threshold 和 extended
CvSURFParams params = cvSURFParams(500, 1);
double tt = (double)cvGetTickCount();
//提取物体(object)和场景(scene)的图像的特征点及其描述器
cvExtractSURF( object, 0, &objectKeypoints, &objectDescriptors, storage, params );
printf(“Object Descriptors: %d\n”, objectDescriptors->total);
cvExtractSURF( image, 0, &imageKeypoints, &imageDescriptors, storage, params );
printf(“Image Descriptors: %d\n”, imageDescriptors->total);
//计算所消耗的时间
tt = (double)cvGetTickCount() – tt;
printf( “Extraction time = %gms\n”, tt/(cvGetTickFrequency()*1000.));
CvPoint src_corners[4] = {
{0,0}, {object->width,0}, {object->width, object->height}, {0, object->height}};
//定义感兴趣的区域
CvPoint dst_corners[4];
IplImage* correspond = cvCreateImage( cvSize(image->width, object->height+image->height), 8, 1 );
//设置感兴趣区域
cvSetImageROI( correspond, cvRect( 0, 0, object->width, object->height ) );
cvCopy( object, correspond );
cvSetImageROI( correspond, cvRect( 0, object->height, correspond->width, correspond->height ) );
cvCopy( image, correspond );
cvResetImageROI( correspond );
//寻找物体(object)在场景(image)中的位置，并将信息保存
if( locatePlanarObject( objectKeypoints, objectDescriptors, imageKeypoints,
imageDescriptors, src_corners, dst_corners ))
{
for( i = 0; i < 4; i++ )
{
CvPoint r1 = dst_corners[i%4];
CvPoint r2 = dst_corners[(i+1)%4];
cvLine( correspond, cvPoint(r1.x, r1.y+object->height ),
cvPoint(r2.x, r2.y+object->height ), colors[8] );
}
}
//定义并保存物体(object)在场景(image)图形之间的匹配点对，并将其存储在向量 ptpairs 中，之后可以对
//ptpairs 进行操作
vector<int> ptpairs;
findPairs( objectKeypoints, objectDescriptors, imageKeypoints, imageDescriptors, ptpairs );
//显示匹配结果
for( i = 0; i < (int)ptpairs.size(); i += 2 )
{
CvSURFPoint* r1 = (CvSURFPoint*)cvGetSeqElem( objectKeypoints, ptpairs[i] );
CvSURFPoint* r2 = (CvSURFPoint*)cvGetSeqElem( imageKeypoints, ptpairs[i+1] );
cvLine( correspond, cvPointFrom32f(r1->pt),
cvPoint(cvRound(r2->pt.x), cvRound(r2->pt.y+object->height)), colors[8] );
}
cvShowImage( “Object Correspond”, correspond );
//显示物体(object)的所有特征点
for( i = 0; i < objectKeypoints->total; i++ )
{
CvSURFPoint* r = (CvSURFPoint*)cvGetSeqElem( objectKeypoints, i );
CvPoint center;
int radius;
center.x = cvRound(r->pt.x);
center.y = cvRound(r->pt.y);
radius = cvRound(r->size*1.2/9.*2);
cvCircle( object_color, center, radius, colors[0], 1, 8, 0 );
}
cvShowImage( “Object”, object_color );
cvWaitKey(0);
//释放窗口所占用的内存
cvDestroyWindow(“Object”);
cvDestroyWindow(“Object Correspond”);
return 0;
}
通过调试运行，可以得到dst_corners中的数据如下：
ptpairs 中的数据如下：
也就是：
[78]
(
29 484
77 134
82 274
206 797
228 210
243 203
244 203
249 404
295 105
347 451
360 142
417 190
427 191
436 198
445 204
452 211
466 218
473 105
486 684
502 133
521 169
522 178
527 190
530 190
532 450
533 198
535 197
539 205
542 202
544 208
547 483
558 412
559 412
583 623
586 624
587 624
594 748
595 654
597 657
)
总共有39对匹配点，第一列表示物体(object)图像中匹配上的点，第一列表示场景(image)图像中匹配的点，其实也就是物体(object)图像中第28个特征点和场景(image)图像中第484个特征点相匹配。
通过这种索引（可以这么说ptpairs中存储的是索引）可以求的那个特征点的坐标，如下：
//取得图像中第i个特征点
CvSURFPoint* r = (CvSURFPoint*)cvGetSeqElem( objectKeypoints, i );
或
//通过ptpairs取得图像中第ptpairs[i] 个特征点，这个特征点是匹配上的点
CvSURFPoint* r1 = (CvSURFPoint*)cvGetSeqElem( objectKeypoints, ptpairs[i] );
运行示意图如下：
关于find_obj.cpp中的
#ifdef USE_FLANN
flannFindPairs( objectKeypoints, objectDescriptors, imageKeypoints, imageDescriptors, ptpairs );
#else
findPairs( objectKeypoints, objectDescriptors, imageKeypoints, imageDescriptors, ptpairs );
#endif
不外乎就是在与说明是否采用 approximate nearest-neighbor 方法来寻找匹配点对。