SIFT+RANSAC+反距离加权图像拼接融合（代码记录）

单纯记录一下之前的代码

#include <opencv2/imgproc/imgproc.hpp> 
#include <opencv2/opencv.hpp>  
#include<opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <cv.h>
#include <iostream>  
#include<fstream>
#include <stdio.h>  
//#include<opencv2/legacy/legacy.hpp>
#include <string>
#include <vector>
#include<opencv2/xfeatures2d.hpp>
#include<algorithm>
using namespace cv;
using namespace std;
using namespace cv::xfeatures2d;
using namespace cv::ml;

int main()
{

	int64 begintime = getTickCount();

	ofstream outfile1, outfile2;
	outfile1.open("matrix.txt");
	//outfile2.open("mmright.txt");



	//读入图像并进行灰度处理 
	//Mat img1 = imread("D:/SU/Tao13_LF_Depth_Light/DEPTH_FROM_CORR_DEFOC_JPEGONLY/Desk/9168/AllFocus.tif");
	//Mat img2 = imread("D:/SU/Tao13_LF_Depth_Light/DEPTH_FROM_CORR_DEFOC_JPEGONLY/Desk/9169/AllFocus.tif");
	Mat img1 = imread("panoLeft.tif");
	Mat img2 = imread("panoRight.tif");
	
	
	//第一步，用SIFT算子检测关键点；

	Ptr<SiftFeatureDetector>detector=SiftFeatureDetector::create();
	std::vector<KeyPoint> m_LeftKey, m_RightKey;//构造2个专门由点组成的点向量用来存储特征点
	detector->detect(img1, m_LeftKey);//将img1图像中检测到的特征点存储起来放在m_LeftKey中
	detector->detect(img2, m_RightKey);//同理
	cout << "图像1特征点的个数：" << m_LeftKey.size() << endl;
	cout << "图像2特征点的个数：" << m_RightKey.size() << endl;
	//计算特征向量

	
	cv::Mat descriptors1, descriptors2;//存放特征向量的矩阵
	detector->compute(img1, m_LeftKey, descriptors1);
	detector->compute(img2, m_RightKey, descriptors2);
	
	//画出特征点  
	Mat img_m_LeftKey, img_m_RightKey;
	drawKeypoints(img1, m_LeftKey, img_m_LeftKey, Scalar::all(-1), 0);  //cvScalar(255,0,0)画的圈圈是蓝色，对应于特征点的颜色,DrawMatchesFlags::DRAW_RICH_KEYPOINTS表示关键点上圆圈的尺寸与特征的尺度成正比，对应于0，是“标志位”的意思
	drawKeypoints(img2, m_RightKey, img_m_RightKey, Scalar::all(-1), 0);
	imshow("Keysrc1", img_m_LeftKey);
	imshow("Keysrc2", img_m_RightKey);
	imwrite("图像2的特征点.tif", img_m_LeftKey);
	imwrite("图像3的特征点.tif", img_m_RightKey);
	
	BFMatcher matcher(NORM_L2, false);
	vector<vector<DMatch>> matches2;
	vector<DMatch>matches;
	matcher.match(descriptors1, descriptors2, matches);
	const float ratio = 0.7;
	matches.clear();
	matcher.knnMatch(descriptors1, descriptors2, matches2, 2);
	for (int n = 0; n < matches2.size(); n++)
	{
		DMatch& bestmatch = matches2[n][0];
		DMatch& bettermatch = matches2[n][1];
		if (bestmatch.distance < ratio*bettermatch.distance)
		{
			matches.push_back(bestmatch);
		}
	}
	cout << "match个数:" << matches.size() << endl;
	//画出直接匹配结果
	Mat key_matches;
	drawMatches(img1, m_LeftKey, img2, m_RightKey, matches, key_matches, CV_RGB(255, 0, 255), CV_RGB(0, 255, 0), Mat(), 2);
	imshow("ratio_matches", key_matches);
	imwrite("ratio_matches.tif", key_matches);

	//RANSAC匹配过程
	vector<DMatch> m_Matches =matches;

	//cout<<"m_Matches="<<m_Matches.size()<<endl;
	// 分配空间
	int ptCount = (int)m_Matches.size();
	//cout<<"m_Matches="<<ptCount<<endl;
	Mat p1(ptCount, 2, CV_32F);
	Mat p2(ptCount, 2, CV_32F);
	//cout<<"p1="<<p1<<endl;
	// 把Keypoint转换为Mat
	Point2f pt;
	for (int i = 0; i<ptCount; i++)
	{
		pt = m_LeftKey[m_Matches[i].queryIdx].pt;
		p1.at<float>(i, 0) = pt.x;
		p1.at<float>(i, 1) = pt.y;
		pt = m_RightKey[m_Matches[i].trainIdx].pt;
		p2.at<float>(i, 0) = pt.x;
		p2.at<float>(i, 1) = pt.y;
	}
	//cout<<"p1="<<p1<<endl;//图1的匹配点坐标
	//cout<<"p2="<<p2<<endl;//图2的匹配点坐标
	// 用RANSAC方法计算F(基础矩阵)
	Mat m_Fundamental;
	vector<uchar> m_RANSACStatus;       // 这个变量用于存储RANSAC后每个点的状态
	findFundamentalMat(p1, p2, m_RANSACStatus, FM_RANSAC);
	// 计算内点个数
	int OutlinerCount = 0;
	for (int i = 0; i<ptCount; i++)
	{
		if (m_RANSACStatus[i] == 0)    // 状态为0表示外点
		{
			OutlinerCount++;
		}
	}
	int InlinerCount = ptCount - OutlinerCount;   // 计算内点	
	cout << "内点数为：" << InlinerCount << endl;
	cout << "外点数为：" << OutlinerCount << endl;
	// 这三个变量用于保存内点和匹配关系
	vector<Point2f> m_LeftInlier;
	vector<Point2f> m_RightInlier;
	vector<DMatch> m_InlierMatches;
	m_InlierMatches.resize(InlinerCount);
	m_LeftInlier.resize(InlinerCount);
	m_RightInlier.resize(InlinerCount);
	InlinerCount = 0;
	float inlier_minRx = img1.cols;        //用于存储内点中右图最小横坐标，以便后续融合
										   //cout<<"inlier="<<inlier_minRx<<endl;
	for (int i = 0; i<ptCount; i++)
	{
		if (m_RANSACStatus[i] != 0)
		{
			m_LeftInlier[InlinerCount].x = p1.at<float>(i, 0);
			m_LeftInlier[InlinerCount].y = p1.at<float>(i, 1);
			m_RightInlier[InlinerCount].x = p2.at<float>(i, 0);
			m_RightInlier[InlinerCount].y = p2.at<float>(i, 1);
			m_InlierMatches[InlinerCount].queryIdx = InlinerCount;
			m_InlierMatches[InlinerCount].trainIdx = InlinerCount;

		



			if (m_RightInlier[InlinerCount].x<inlier_minRx) inlier_minRx = m_RightInlier[InlinerCount].x;   //存储内点中右图最小横坐标
			InlinerCount++;
		}
	}
	cout<<"inlier="<<inlier_minRx<<endl;//最小横坐标	

	// 把内点转换为drawMatches可以使用的格式
	vector<KeyPoint> key1(InlinerCount);
	vector<KeyPoint> key2(InlinerCount);
	KeyPoint::convert(m_LeftInlier, key1);
	KeyPoint::convert(m_RightInlier, key2);
	
	// 显示计算F过后的内点匹配
	Mat OutImage;
	drawMatches(img1, key1, img2, key2, m_InlierMatches, OutImage, CV_RGB(255, 0, 255), CV_RGB(0, 255, 0));
	imshow("RANSAC match features", OutImage);
	imwrite("RansacMatchb.tif", OutImage);

	//矩阵H用以存储RANSAC得到的单应矩阵
	Mat H = findHomography(m_LeftInlier, m_RightInlier, RANSAC);
	outfile1 << H << endl;
	cout << H << endl;
	//存储左图四角，及其变换到右图位置
	std::vector<Point2f> obj_corners(4);
	obj_corners[0] = Point(0, 0); obj_corners[1] = Point(img1.cols, 0);
	obj_corners[2] = Point(img1.cols, img1.rows); obj_corners[3] = Point(0, img1.rows);
	std::vector<Point2f> scene_corners(4);
	perspectiveTransform(obj_corners, scene_corners, H);
/*	//画出变换后图像位置
	Point2f offset((float)img1.cols, 0);
	line(OutImage, scene_corners[0] + offset, scene_corners[1] + offset, Scalar(0, 255, 0), 4);
	line(OutImage, scene_corners[1] + offset, scene_corners[2] + offset, Scalar(0, 255, 0), 4);
	line(OutImage, scene_corners[2] + offset, scene_corners[3] + offset, Scalar(0, 255, 0), 4);
	line(OutImage, scene_corners[3] + offset, scene_corners[0] + offset, Scalar(0, 255, 0), 4);
	imshow("Good Matches & Object detection", OutImage);
	imwrite("transform.tif", OutImage);*/
	int drift = scene_corners[1].x;
	//储存偏移量
	cout << "scene=" << scene_corners << endl;
	cout << "scene0=" << scene_corners[0].x << endl;
	cout << "scene1=" << scene_corners[1].x << endl;
	cout << "scene2=" << scene_corners[2].x << endl;
	cout << "scene3=" << scene_corners[3].x << endl;
	//新建一个矩阵存储配准后四角的位置
	int width = int(max(abs(scene_corners[1].x), abs(scene_corners[2].x)));//等于左图右边两个顶点中横坐标较大的那个
	int height = img1.rows;                                                                  //或者：int height = int(max(abs(scene_corners[2].y), abs(scene_corners[3].y)));
	float origin_x = 0, origin_y = 0;
	if (scene_corners[0].x<0) {
		if (scene_corners[3].x<0) origin_x += min(scene_corners[0].x, scene_corners[3].x);
		else origin_x += scene_corners[0].x;
	}
	width -= int(origin_x);//图像的宽度变为左图变换之后左边的定点坐标小的数值+右边的横坐标较大的数值
	if (scene_corners[0].y<0) {
		if (scene_corners[1].y) origin_y += min(scene_corners[0].y, scene_corners[1].y);
		else origin_y += scene_corners[0].y;
	}
	//height-=int(origin_y);

	//imageturn为变换只有的左图的大小
	Mat imageturn = Mat::zeros(width, height, img1.type());
	cout << "PerspectiveWidth: " << width << endl;
	cout << "PerspectiveHeight: " << height << endl;
	cout << "img1.type(): " << img1.type() << endl;
	
	//获取新的变换矩阵，使图像完整显示
	for (int i = 0; i<4; i++) { scene_corners[i].x -= origin_x; }//相当于把整个图像往右移动至右图上，移动距离为初次变换之后最小的横坐标	
	Mat H1 = getPerspectiveTransform(obj_corners, scene_corners);//然后移动之后的与右图再求一次矩阵
	//进行图像变换，显示效果
	warpPerspective(img1, imageturn, H1, Size(width, height));
	imshow("image_Perspective", imageturn);
	imwrite("perpective.tif", imageturn);

	cout << "origin_x=" << origin_x << endl;
	cout << "origin_y=" << origin_y << endl;
	cout << "img1.width=" << img1.cols << endl;
	cout << "inlier_minRx=" << inlier_minRx << endl;
	cout << "scene_corners= " << scene_corners << endl;
	cout << "单应矩阵H=" << H << endl;
	cout << "变换矩阵H1=" << H1 << endl;
	
	//图像融合
	int width_ol = width - int(inlier_minRx - origin_x);//重叠区域的宽
	int start_x = int(inlier_minRx - origin_x);//imageturn上的重叠区域的起始坐标  右侧内点的最小点-origin_x
	cout << "start_X:" << start_x << endl;
	/*//imageturn为左图图转化之后的 重叠区域 开始的横坐标为右图上的最小横坐标的内点在imageturn上的坐标，右侧为左图的右边界
	Mat subValue_image = Mat::zeros(imageturn.rows, width_ol, imageturn.type());//创建一个图用于保存差值图像

	Mat sub1 = imageturn(Rect(start_x, 0, width_ol, imageturn.rows));//start_x  重叠区域的宽 高不变
	Mat sub2 = img2(Rect(inlier_minRx, 0, width_ol, img2.rows));
	
	//求两个图像插值
	absdiff(sub1, sub2, subValue_image);
	imshow("差值图", subValue_image);

	//创建一个mask
	Mat mask = Mat::zeros(height, width_ol, img1.type());//重叠区域大小
	Mat imgL = imageturn(Rect(start_x, 0, width_ol, height));
	Mat imgR = img2(Rect(0, 0, width_ol, height));
	mask(Rect(0, 0, width_ol*0.5, height)).setTo(255);
	imshow("mask", mask);

	Point Center(imgL.cols*0.5, imgL.rows*0.5);

	Mat resultImage;
	resultImage.create(imgL.size(), imgL.type());
	seamlessClone(imgL, imgR, mask, Center, resultImage, NORMAL_CLONE);

	imshow("blending", resultImage);*/


	uchar* ptr = imageturn.data;//perspective的图像数据
	double alpha = 0, beta = 1;//定义权重
	for (int row = 0; row<height; row++) //图像转换之后的高度
	{
		//step可以理解为Mat矩阵中每一行的“步长”，以字节为基本单位，
		//每一行中所有元素的字节总量，是累计了一行中所有元素、所有通道、所有通道的elemSize1之后的值；	
		ptr = imageturn.data + row*imageturn.step + (start_x)*imageturn.elemSize(); //data   uchar类型的指针，指向Mat数据矩阵的首地址。可以理解为标示一个房屋的门牌号；
		for (int col = 0; col<width_ol; col++)//图像转换之后宽度			
		{
			uchar* ptr_c1 = ptr + imageturn.elemSize1();
			uchar*  ptr_c2 = ptr_c1 + imageturn.elemSize1();
			uchar* ptr2 = img2.data + row*img2.step + (col + int(inlier_minRx))*img2.elemSize();
			uchar* ptr2_c1 = ptr2 + img2.elemSize1();  
			uchar* ptr2_c2 = ptr2_c1 + img2.elemSize1();
			alpha = double(col) / double(width_ol); 
			beta = 1 - alpha;
			if (*ptr == 0 && *ptr_c1 == 0 && *ptr_c2 == 0) {//左图中转换了之后的无像素值的黑点，则完全拷贝右侧图的像素值
				*ptr = (*ptr2);//该像素第一通道
				*ptr_c1 = (*ptr2_c1);//该像素第二通道
				*ptr_c2 = (*ptr2_c2);//该像素第三通道
			}
			//如不为0
			*ptr = (*ptr)*beta + (*ptr2)*alpha;//左图通道1的像素值乘以权重+右侧通道1像素值乘以权重
			*ptr_c1 = (*ptr_c1)*beta + (*ptr2_c1)*alpha;
			*ptr_c2 = (*ptr_c2)*beta + (*ptr2_c2)*alpha;
			ptr += imageturn.elemSize();//指针指向下一个元素
		}
	}
	imshow("image_overlap", imageturn);

	Mat img_result = Mat::zeros(height, width + img2.cols - drift, img1.type());//int drift = scene_corners[1].x;
	cout << "宽：" << width + img2.cols - drift << endl;
	cout << "高：" << height << endl;
	uchar* ptr_r = imageturn.data;
	for (int row = 0; row<height; row++) {
		ptr_r = img_result.data + row*img_result.step;
		for (int col = 0; col<imageturn.cols; col++)
		{
			uchar* ptr_rc1 = ptr_r + imageturn.elemSize1();  uchar*  ptr_rc2 = ptr_rc1 + imageturn.elemSize1();
			uchar* ptr = imageturn.data + row*imageturn.step + col*imageturn.elemSize();
			uchar* ptr_c1 = ptr + imageturn.elemSize1();  uchar*  ptr_c2 = ptr_c1 + imageturn.elemSize1();
			*ptr_r = *ptr;
			*ptr_rc1 = *ptr_c1;
			*ptr_rc2 = *ptr_c2;
			ptr_r += img_result.elemSize();
		}
		ptr_r = img_result.data + row*img_result.step + imageturn.cols*img_result.elemSize();
		for (int col = imageturn.cols; col<img_result.cols; col++)
		{
			uchar* ptr_rc1 = ptr_r + imageturn.elemSize1();  uchar*  ptr_rc2 = ptr_rc1 + imageturn.elemSize1();
			uchar* ptr2 = img2.data + row*img2.step + (col - imageturn.cols + drift)*img2.elemSize();
			uchar* ptr2_c1 = ptr2 + img2.elemSize1();  uchar* ptr2_c2 = ptr2_c1 + img2.elemSize1();
			*ptr_r = *ptr2;
			*ptr_rc1 = *ptr2_c1;
			*ptr_rc2 = *ptr2_c2;
			ptr_r += img_result.elemSize();
		}
	}
	imshow("image_result", img_result);
	imwrite("panoAll.tif", img_result);
	//imwrite("14.jpg",img_result);
	int64 endtime = getTickCount();
	cout << "时间：" << (endtime - begintime) / getTickFrequency() << "s" << endl;
	waitKey(0);
	return 0;
}