该博客介绍了如何利用opencv2.4.10和Visual Studio 2012结合RobHess的SIFT特征算法来实现全景图像的拼接。详细步骤包括opencv配置、SIFT源码下载与测试、创建新项目、添加源文件、取消预编译头、C/C++混合编程以及关键步骤的代码实现，如特征提取、匹配查找、RANSAC筛选和图像融合。
利用opencv2.4.10+VS2012和RobHess的sift特征算法实现全景图像拼接
一、保证opencv配置成功，可以参考以下网址来配置。
openCV的配置
二、下载RobHess的SIFT源码并测试使用其源代码。
1、下载
下载“sift-latest_win.zip”这个文件（下载链接）。
2、测试（自信opencv配置成功的话可以跳过这一步）
测试使用其源码可参考此网址。
在测试过程可能会遇到如下错误：
“error LNK1104: 无法打开文件“cv200d.lib” D:\sift-latest_win\sift-1.1.2_20101207_win\match\LINK match”
在菜单中点击”项目”->”属性”->”配置属性”->”链接器”->”输入”，把”附加依赖项”的其中四个lib文件删除，文件名如图
如果能成功运行程序，证明环境配置成功，可以使用RobHess的源码。
三、实现全景拼接的具体操作。
1、新建一个控制台应用程序并给它命名，例如“RobHess_SIFT”。
2、添加所需项。
从下载的文件夹当中找到所需的头文件（.h）和源文件（.c）拷贝到项目的默认路径，然后将其添加到”解决方案资源管理器“，如图
（图中显示的即为所需文件）
3、取消.c文件的预编译。
“项目”->”属性”->”配置属性”->”C/C++”->”预编译头”，在”预编译头”选项中选择”不适用预编译头”。
4、C语法设定，实现C语言在C++文件里面混合编程。
方法：分别打开xform.h、sift.h、kdtree.h、utils.h、imgfeatures.h（以后的代码中有用到这些头文件里面声明的函数，所以在主函数所在文件要把这些头文件都include进来），把所有函数包含在
#ifdef __cplusplus
extern "C" {
#endif和
#ifdef __cplusplus
}
#endif之间。
5、（如果跳过了测试步骤的要做这一步）
在util.c中添加“#include < stdargh.h>,不然编译时候会报错。
6、利用源代码实现全景拼接。
主要分为以下几个步骤：
(1) 读入两张图片并分别提取SIFT特征
(2) 利用k-d tree和BBF算法进行特征匹配查找
(3) 利用RANSAC算法筛选匹配点并计算变换矩阵
(4) 图像融合
(1) 读入两张图片并分别提取SIFT特征
调用RobHess源码中的sift_features()函数进行默认参数的SIFT特征提取，主要代码如下：
img1_Feat = cvCloneImage(img1);//复制图1，深拷贝，用来画特征点
img2_Feat = cvCloneImage(img2);//复制图2，深拷贝，用来画特征点
//默认提取的是LOWE格式的SIFT特征点
//提取并显示第1幅图片上的特征点
n1 = sift_features( img1, &feat1 );//检测图1中的SIFT特征点,n1是图1的特征点个数
export_features("feature1.txt",feat1,n1);//将特征向量数据写入到文件
draw_features( img1_Feat, feat1, n1 );//画出特征点
cvNamedWindow(IMG1_FEAT);//创建窗口
cvShowImage(IMG1_FEAT,img1_Feat);//显示
//提取并显示第2幅图片上的特征点
n2 = sift_features( img2, &feat2 );//检测图2中的SIFT特征点，n2是图2的特征点个数
export_features("feature2.txt",feat2,n2);//将特征向量数据写入到文件
draw_features( img2_Feat, feat2, n2 );//画出特征点
cvNamedWindow(IMG2_FEAT);//创建窗口
cvShowImage(IMG2_FEAT,img2_Feat);//显示 以上要添加定义的变量（定义为全局变量，以后定义的也是）有
IplImage *img1, *img2, *img1_Feat, *img2_Feat;
struct feature *feat1, *feat2;
int n1, n2;其中img1，img2保存载入的两张待拼接的图像。
检测出的SIFT特征点如下：
(2) 利用k-d tree和BBF算法进行特征匹配查找，并根据最近邻和次近邻距离比值进行初步筛选
也是调用RobHess源码中的函数，加上之后的一些筛选处理，主要代码如下：
//根据图1的特征点集feat1建立k-d树，返回k-d树根给kd_root
kd_root = kdtree_build( feat1, n1 );
Point pt1,pt2;//连线的两个端点
double d0,d1;//feat2中每个特征点到最近邻和次近邻的距离
int matchNum = 0;//经距离比值法筛选后的匹配点对的个数
//遍历特征点集feat2，针对feat2中每个特征点feat，选取符合距离比值条件的匹配点，放到feat的fwd_match域中
for(int i = 0; i < n2; i++ )
{
feat = feat2+i;//第i个特征点的指针
//在kd_root中搜索目标点feat的2个最近邻点，存放在nbrs中，返回实际找到的近邻点个数
int k = kdtree_bbf_knn( kd_root, feat, 2, &nbrs, KDTREE_BBF_MAX_NN_CHKS );
if( k == 2 )
{
d0 = descr_dist_sq( feat, nbrs[0] );//feat与最近邻点的距离的平方
d1 = descr_dist_sq( feat, nbrs[1] );//feat与次近邻点的距离的平方
//若d0和d1的比值小于阈值NN_SQ_DIST_RATIO_THR，则接受此匹配，否则剔除
if( d0 < d1 * NN_SQ_DIST_RATIO_THR )
{ //将目标点feat和最近邻点作为匹配点对
pt2 = Point( cvRound( feat->x ), cvRound( feat->y ) );//图2中点的坐标
pt1 = Point( cvRound( nbrs[0]->x ), cvRound( nbrs[0]->y ) );//图1中点的坐标(feat的最近邻点)
pt2.x += img1->width;//由于两幅图是左右排列的，pt2的横坐标加上图1的宽度，作为连线的终点
cvLine( stacked, pt1, pt2, CV_RGB(255,0,255), 1, 8, 0 );//画出连线
matchNum++;//统计匹配点对的个数
feat2[i].fwd_match = nbrs[0];//使点feat的fwd_match域指向其对应的匹配点
}
}
free( nbrs );//释放近邻数组
}
//显示并保存经距离比值法筛选后的匹配图
cvNamedWindow(IMG_MATCH1);//创建窗口
cvShowImage(IMG_MATCH1,stacked);//显示 以上要添加定义的变量有
struct kd_node *kd_root;
struct feature *feat;
struct feature** nbrs;
IplImage *stacked = stack_imgs1( img1, img2 );其中stacked变量需要用到的stack_imgs1函数（用途：把两张图直接左右拼在一起）在sift源文件里并没有，需要模仿utils文件里stack_imgs函数（用途：把两张图直接上下拼在一起）来定义。而如果使用的图片用于上下拼接，则只需调用原有的stack_imgs函数即可
定义方法：
1、打开utils.h文件，找到“extern IplImage* stack_imgs( IplImage* img1, IplImage* img2 );”，在其下方插入”extern IplImage* stack_imgs1( IplImage* img1, IplImage* img2 );“。
2、打开utils.c文件，在定义方法的地方添加以下代码。
extern IplImage* stack_imgs1( IplImage* img1, IplImage* img2 )
{
IplImage* expanded = cvCreateImage( cvSize(
img1->width + img2->width,MAX(img1->height, img2->height)),
IPL_DEPTH_8U, 3 );
cvZero( expanded );
cvSetImageROI( expanded, cvRect( 0, 0, img1->width, img1->height ) );
cvAdd( img1, expanded, expanded, NULL );
cvSetImageROI( expanded, cvRect(img1->width, 0, img2->width, img2->height) );
cvAdd( img2, expanded, expanded, NULL );
cvResetImageROI( expanded );
return expanded;
}此外，还需宏定义添加两个标识符。
#define KDTREE_BBF_MAX_NN_CHKS 200 //用BBF算法搜素最近邻近点的最大搜索次数
#define NN_SQ_DIST_RATIO_THR 0.49 //为最近邻点和次近邻点间距离平方所设的阀值匹配结果如下：
(3) 利用RANSAC算法筛选匹配点并计算变换矩阵
此部分最主要的是RobHess源码中的ransac_xform()函数，此函数实现了用RANSAC算法筛选匹配点，返回结果是计算好的变换矩阵。
此部分中利用匹配点的坐标关系，对输入的两幅图像的左右关系进行了判断，并根据结果选择使用矩阵H或H的逆阵进行变换,所以读入的两幅要拼接的图像的左右位置关系可以随意，程序中可自动调整。
主要代码如下：
//利用RANSAC算法筛选匹配点,计算变换矩阵H，
//无论img1和img2的左右顺序，计算出的H永远是将feat2中的特征点变换为其匹配点，即将img2中的点变换为img1中的对应点
H = ransac_xform(feat2,n2,FEATURE_FWD_MATCH,lsq_homog,4,0.01,homog_xfer_err,3.0,&inliers,&n_inliers);
//若能成功计算出变换矩阵，即两幅图中有共同区域
if( H )
{
cout<<"经RANSAC算法筛选后的匹配点对个数："<<n_inliers<<endl; //输出筛选后的匹配点对个数
int invertNum = 0;//统计pt2.x > pt1.x的匹配点对的个数，来判断img1中是否右图
//遍历经RANSAC算法筛选后的特征点集合inliers，找到每个特征点的匹配点，画出连线
for(int i=0; i<n_inliers; i++)
{
feat = inliers[i];//第i个特征点
pt2 = Point(cvRound(feat->x), cvRound(feat->y));//图2中点的坐标
pt1 = Point(cvRound(feat->fwd_match->x), cvRound(feat->fwd_match->y));//图1中点的坐标(feat的匹配点)
//统计匹配点的左右位置关系，来判断图1和图2的左右位置关系
if(pt2.x > pt1.x)
invertNum++;
pt2.x += img1->width;//由于两幅图是左右排列的，pt2的横坐标加上图1的宽度，作为连线的终点
cvLine(stacked_ransac,pt1,pt2,CV_RGB(255,0,255),1,8,0);//在匹配图上画出连线
}
cvNamedWindow(IMG_MATCH2);//创建窗口
cvShowImage(IMG_MATCH2,stacked_ransac);//显示经RANSAC算法筛选后的匹配图
/*程序中计算出的变换矩阵H用来将img2中的点变换为img1中的点，正常情况下img1应该是左图，img2应该是右图。
此时img2中的点pt2和img1中的对应点pt1的x坐标的关系基本都是：pt2.x < pt1.x
若用户打开的img1是右图，img2是左图，则img2中的点pt2和img1中的对应点pt1的x坐标的关系基本都是：pt2.x > pt1.x
所以通过统计对应点变换前后x坐标大小关系，可以知道img1是不是右图。
如果img1是右图，将img1中的匹配点经H的逆阵H_IVT变换后可得到img2中的匹配点*/
//若pt2.x > pt1.x的点的个数大于内点个数的80%，则认定img1中是右图
if(invertNum > n_inliers * 0.8)
{
CvMat * H_IVT = cvCreateMat(3, 3, CV_64FC1);//变换矩阵的逆矩阵
//求H的逆阵H_IVT时，若成功求出，返回非零值
if( cvInvert(H,H_IVT) )
{
cvReleaseMat(&H);//释放变换矩阵H，因为用不到了
H = cvCloneMat(H_IVT);//将H的逆阵H_IVT中的数据拷贝到H中
cvReleaseMat(&H_IVT);//释放逆阵H_IVT
//将img1和img2对调
IplImage * temp = img2;
img2 = img1;
img1 = temp;
}
else//H不可逆时，返回0
{
cvReleaseMat(&H_IVT);//释放逆阵H_IVT
cout<<"变换矩阵H不可逆"<<endl;
}
}
}
else //无法计算出变换矩阵，即两幅图中没有重合区域
{
cout<<"Warning: 两图中无公共区域";
} 此处要添加定义的变量：
struct feature **inliers;
int n_inliers;
CvMat *H;
IplImage *stacked_ransac = stack_imgs1( img1, img2 );经RANSAC筛选后的匹配结果如下图：
(4) 图像融合
这里有两种拼接方法：
① 简易拼接方法的过程是：首先将右图img2经变换矩阵H变换到一个新图像中，然后直接将左图img1加到新图像中，这样拼接出来会有明显的拼接缝，但也是一个初步的成品了。
② 另一种方法首先也是将右图img2经变换矩阵H变换到一个新图像中，然后图像的融合过程将目标图像分为三部分，最左边完全取自img1中的数据，中间的重合部分是两幅图像的加权平均，重合区域右边的部分完全取自img2经变换后的图像。加权平均的权重选择也有好多方法，比如可以使用最基本的取两张图像的平均值，但这样会有明显的拼接缝。这里首先计算出拼接区域的宽度，设d1，d2分别是重叠区域中的点到重叠区域左边界和右边界的距离，则使用如下公式计算重叠区域的像素值：
主要代码如下：
//若能成功计算出变换矩阵，即两幅图中有共同区域，才可以进行全景拼接
if(H)
{
//拼接图像，img1是左图，img2是右图
CalcFourCorner();//计算图2的四个角经变换后的坐标
//为拼接结果图xformed分配空间,高度为图1图2高度的较小者，根据图2右上角和右下角变换后的点的位置决定拼接图的宽度
xformed = cvCreateImage(cvSize(MIN(rightTop.x,rightBottom.x),MIN(img1->height,img2->height)),IPL_DEPTH_8U,3);
//用变换矩阵H对右图img2做投影变换(变换后会有坐标右移)，结果放到xformed中
cvWarpPerspective(img2,xformed,H,CV_INTER_LINEAR + CV_WARP_FILL_OUTLIERS,cvScalarAll(0));
cvNamedWindow(IMG_MOSAIC_TEMP); //显示临时图,即只将图2变换后的图
cvShowImage(IMG_MOSAIC_TEMP,xformed);
//简易拼接法：直接将将左图img1叠加到xformed的左边
xformed_simple = cvCloneImage(xformed);//简易拼接图，克隆自xformed
cvSetImageROI(xformed_simple,cvRect(0,0,img1->width,img1->height));
cvAddWeighted(img1,1,xformed_simple,0,0,xformed_simple);
cvResetImageROI(xformed_simple);
cvNamedWindow(IMG_MOSAIC_SIMPLE);//创建窗口
cvShowImage(IMG_MOSAIC_SIMPLE,xformed_simple);//显示简易拼接图
//处理后的拼接图，克隆自xformed
xformed_proc = cvCloneImage(xformed);
//重叠区域左边的部分完全取自图1
cvSetImageROI(img1,cvRect(0,0,MIN(leftTop.x,leftBottom.x),xformed_proc->height));
cvSetImageROI(xformed,cvRect(0,0,MIN(leftTop.x,leftBottom.x),xformed_proc->height));
cvSetImageROI(xformed_proc,cvRect(0,0,MIN(leftTop.x,leftBottom.x),xformed_proc->height));
cvAddWeighted(img1,1,xformed,0,0,xformed_proc);
cvResetImageROI(img1);
cvResetImageROI(xformed);
cvResetImageROI(xformed_proc);
cvNamedWindow(IMG_MOSAIC_BEFORE_FUSION);
cvShowImage(IMG_MOSAIC_BEFORE_FUSION,xformed_proc);//显示融合之前的拼接图
//采用加权平均的方法融合重叠区域
int start = MIN(leftTop.x,leftBottom.x) ;//开始位置，即重叠区域的左边界
double processWidth = img1->width - start;//重叠区域的宽度
double alpha = 1;//img1中像素的权重
for(int i=0; i<xformed_proc->height; i++)//遍历行
{
const uchar * pixel_img1 = ((uchar *)(img1->imageData + img1->widthStep * i));//img1中第i行数据的指针
const uchar * pixel_xformed = ((uchar *)(xformed->imageData + xformed->widthStep * i));//xformed中第i行数据的指针
uchar * pixel_xformed_proc = ((uchar *)(xformed_proc->imageData + xformed_proc->widthStep * i));//xformed_proc中第i行数据的指针
for(int j=start; j<img1->width; j++)//遍历重叠区域的列
{
//如果遇到图像xformed中无像素的黑点，则完全拷贝图1中的数据
if(pixel_xformed[j*3] < 50 && pixel_xformed[j*3+1] < 50 && pixel_xformed[j*3+2] < 50 )
{
alpha = 1;
}
else
{ //img1中像素的权重，与当前处理点距重叠区域左边界的距离成正比
alpha = (processWidth-(j-start)) / processWidth ;
}
pixel_xformed_proc[j*3] = pixel_img1[j*3] * alpha + pixel_xformed[j*3] * (1-alpha);//B通道
pixel_xformed_proc[j*3+1] = pixel_img1[j*3+1] * alpha + pixel_xformed[j*3+1] * (1-alpha);//G通道
pixel_xformed_proc[j*3+2] = pixel_img1[j*3+2] * alpha + pixel_xformed[j*3+2] * (1-alpha);//R通道
}
}
cvNamedWindow(IMG_MOSAIC_PROC);//创建窗口
cvShowImage(IMG_MOSAIC_PROC,xformed_proc);//显示处理后的拼接图
/*重叠区域取两幅图像的平均值，效果不好
//设置ROI，是包含重叠区域的矩形
cvSetImageROI(xformed_proc,cvRect(MIN(leftTop.x,leftBottom.x),0,img1->width-MIN(leftTop.x,leftBottom.x),xformed_proc->height));
cvSetImageROI(img1,cvRect(MIN(leftTop.x,leftBottom.x),0,img1->width-MIN(leftTop.x,leftBottom.x),xformed_proc->height));
cvSetImageROI(xformed,cvRect(MIN(leftTop.x,leftBottom.x),0,img1->width-MIN(leftTop.x,leftBottom.x),xformed_proc->height));
cvAddWeighted(img1,0.5,xformed,0.5,0,xformed_proc);
cvResetImageROI(xformed_proc);
cvResetImageROI(img1);
cvResetImageROI(xformed); //*/
/*对拼接缝周围区域进行滤波来消除拼接缝，效果不好
//在处理前后的图上分别设置横跨拼接缝的矩形ROI
cvSetImageROI(xformed_proc,cvRect(img1->width-10,0,img1->width+10,xformed->height));
cvSetImageROI(xformed,cvRect(img1->width-10,0,img1->width+10,xformed->height));
cvSmooth(xformed,xformed_proc,CV_MEDIAN,5);//对拼接缝周围区域进行中值滤波
cvResetImageROI(xformed);
cvResetImageROI(xformed_proc);
cvShowImage(IMG_MOSAIC_PROC,xformed_proc);//显示处理后的拼接图 */
/*想通过锐化解决变换后的图像失真的问题，对于扭曲过大的图像，效果不好
double a[]={ 0, -1, 0, -1, 5, -1, 0, -1, 0 };//拉普拉斯滤波核的数据
CvMat kernel = cvMat(3,3,CV_64FC1,a);//拉普拉斯滤波核
cvFilter2D(xformed_proc,xformed_proc,&kernel);//滤波
cvShowImage(IMG_MOSAIC_PROC,xformed_proc);//显示处理后的拼接图*/
}此处要定义的函数及变量
IplImage *xformed, *xformed_simple, *xformed_proc;
Point leftTop,leftBottom,rightTop,rightBottom;
void CalcFourCorner();其中CalcFourCorner函数主体为：
void CalcFourCorner()
{
//计算图2的四个角经矩阵H变换后的坐标
double v2[]={0,0,1};//左上角
double v1[3];//变换后的坐标值
CvMat V2 = cvMat(3,1,CV_64FC1,v2);
CvMat V1 = cvMat(3,1,CV_64FC1,v1);
cvGEMM(H,&V2,1,0,1,&V1);//矩阵乘法
leftTop.x = cvRound(v1[0]/v1[2]);
leftTop.y = cvRound(v1[1]/v1[2]);
//cvCircle(xformed,leftTop,7,CV_RGB(255,0,0),2);
//将v2中数据设为左下角坐标
v2[0] = 0;
v2[1] = img2->height;
V2 = cvMat(3,1,CV_64FC1,v2);
V1 = cvMat(3,1,CV_64FC1,v1);
cvGEMM(H,&V2,1,0,1,&V1);
leftBottom.x = cvRound(v1[0]/v1[2]);
leftBottom.y = cvRound(v1[1]/v1[2]);
//cvCircle(xformed,leftBottom,7,CV_RGB(255,0,0),2);
//将v2中数据设为右上角坐标
v2[0] = img2->width;
v2[1] = 0;
V2 = cvMat(3,1,CV_64FC1,v2);
V1 = cvMat(3,1,CV_64FC1,v1);
cvGEMM(H,&V2,1,0,1,&V1);
rightTop.x = cvRound(v1[0]/v1[2]);
rightTop.y = cvRound(v1[1]/v1[2]);
//cvCircle(xformed,rightTop,7,CV_RGB(255,0,0),2);
//将v2中数据设为右下角坐标
v2[0] = img2->width;
v2[1] = img2->height;
V2 = cvMat(3,1,CV_64FC1,v2);
V1 = cvMat(3,1,CV_64FC1,v1);
cvGEMM(H,&V2,1,0,1,&V1);
rightBottom.x = cvRound(v1[0]/v1[2]);
rightBottom.y = cvRound(v1[1]/v1[2]);
//cvCircle(xformed,rightBottom,7,CV_RGB(255,0,0),2);
}右图经变换后的结果如下图：
简易拼接结果如下图：
使用第二种方法时，重合区域融合之前如下图：
加权平均融合之后如下图：
本文参考博文：http://blog.youkuaiyun.com/masibuaa/article/details/9246493
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
