计算机视觉
图像处理

基于SIFT特征的全景图像拼接

主要分为以下几个步骤:

(1) 读入两张图片并分别提取SIFT特征

(2) 利用k-d tree和BBF算法进行特征匹配查找

(3) 利用RANSAC算法筛选匹配点并计算变换矩阵

(3) 图像融合

SIFT算法以及RANSAC算法都是利用的RobHess的SIFT源码,前三个步骤RobHess的源码中都有自带的示例。

(1) SIFT特征提取

直接调用RobHess源码(RobHess的SIFT源码分析:综述)中的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);//显示

检测出的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);//显示

匹配结果如下:


(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 )
    {
        qDebug()<<tr("经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的匹配点)
            //qDebug()<<"pt2:("<<pt2.x<<","<<pt2.y<<")--->pt1:("<<pt1.x<<","<<pt1.y<<")";//输出对应点对 //统计匹配点的左右位置关系,来判断图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;
                ui->mosaicButton->setEnabled(true);//激活全景拼接按钮
            }
            else//H不可逆时,返回0
            {
                cvReleaseMat(&H_IVT);//释放逆阵H_IVT
                QMessageBox::warning(this,tr("警告"),tr("变换矩阵H不可逆"));
            }
        }
        else
            ui->mosaicButton->setEnabled(true);//激活全景拼接按钮
    }
    else //无法计算出变换矩阵,即两幅图中没有重合区域
    {
        QMessageBox::warning(this,tr("警告"),tr("两图中无公共区域"));
    }

经RANSAC筛选后的匹配结果如下图:


(3) 图像融合

这里有两种拼接方法:

① 简易拼接方法的过程是:首先将右图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; iheight; 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; jwidth; 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);//显示处理后的拼接图*/

    }

右图经变换后的结果如下图:

简易拼接结果如下图:

使用第二种方法时,重合区域融合之前如下图:

加权平均融合之后如下图:

用Qt做了个简单的界面,如下:

还有很多不足,经常有黑边无法去除,望大家多多指正。

实验环境:

IDE是 Qt Creator 2.4.1Qt版本是4.8.4OpenCV版本是2.4.4,SIFT源码来自RobHess给出的源码

源码下载:基于SIFT特征的全景图像拼接,Qt工程:

!!温馨提示: 此处需要评论本文后才能查看.

转载注明来源:CV视觉网 » 基于SIFT特征的全景图像拼接

分享到:更多 ()
扫描二维码,给作者 打赏
pay_weixinpay_weixin

请选择你看完该文章的感受:

3不错 3超赞 无聊 扯淡 不解 路过

评论 5

评论前必须登录!

 

  1. #4

    好厉害!

    xiemomo3天前
  2. #3

    好厉害!

    manu_98765个月前 (05-17)
  3. #2

    可以添加 微信ID:readsense 一起探讨计算机视觉

    宁传鹏2年前 (2016-09-26)
  4. #1

    呵呵,原来是周涛教授写的啊

    谭卓2年前 (2016-09-12)