#include "opencv2/opencv.hpp"

#include <opencv2/core/core.hpp>

#include <opencv2/features2d/features2d.hpp>

#include <opencv2/highgui/highgui.hpp>

#include <opencv2/nonfree/features2d.hpp>

#include <opencv2/legacy/legacy.hpp>

#include <iostream>

#include <vector>

#define __BEGIN__ 50





using namespace std;

using namespace cv;





char img1_file[] = "27.jpg";			//读取原图像

char img2_file[] = "28.jpg";





IplImage* stack_imgs1( IplImage* img1, IplImage* img2 )

{

	IplImage* stacked = cvCreateImage( cvSize( img1->width + img2->width,

										MAX (img1->height , img2->height) ),

										IPL_DEPTH_8U, 3 );



	cvZero( stacked );						//清空数组

	cvSetImageROI( stacked, cvRect( 0, 0, img1->width, img1->height ) );			//基于给定的矩形设置感兴趣区域

	cvAdd( img1, stacked, stacked, NULL );					//计算两个数组中每个元素的和

	cvSetImageROI( stacked, cvRect(img1->width, 0, img2->width, img2->height) );

	cvAdd( img2, stacked, stacked, NULL );

	cvResetImageROI( stacked );							//释放图像的ROI 



	return stacked;

}



int main(int argc, char* argv[])

{

	IplImage* img1_cvLoadImage, * img2_cvLoadImage, * stacked;    //表示图像的结构体  声明结构体指针，指向图像





		     //图像融合时用到的参数

			unsigned char* data1;   //temp_xformed_2

			int step1;

			int channels1;



			unsigned char* data2 ;   //img2

			int step2;

			int channels2;



			unsigned char* data3 ;   //temp_xformed_1

			int step3;

			int channels3;



			double w;   //使用加权平均法进行图像融合时的权值





	img1_cvLoadImage = cvLoadImage( img1_file, 1 );   //从文件中读取图像  从指定文件读入图像，返回读入图像的指针

	if( ! img1_cvLoadImage )

		printf( "unable to load image from %s", img1_file );

	img2_cvLoadImage = cvLoadImage( img2_file, 1 );

	if( ! img2_cvLoadImage )

		printf( "unable to load image from %s", img2_file );



    cvNamedWindow( "原图1 img1_cvLoadImage    ", 1 );   

	cvShowImage( "原图1 img1_cvLoadImage ", img1_cvLoadImage );



	cvNamedWindow( "原图2 img2_cvLoadImage  ", 1 );

	cvShowImage( "原图2 img2_cvLoadImage ", img2_cvLoadImage );



	//stacked = stack_imgs( img2_cvLoadImage, img1_cvLoadImage );			//上下组合两幅图像

	stacked = stack_imgs1( img2_cvLoadImage, img1_cvLoadImage );			//左右组合两幅图像



	

	vector<KeyPoint> keys1, keys2;   //存储关键点

	Mat descriptors1, descriptors2;

	

	SurfFeatureDetector detector_Freak(2000);

	BruteForceMatcher<L2<float>> matcher;

     std::vector< DMatch > matches;

	 

	FREAK freak;

	int64 st, et;

	st = cvGetTickCount();

	

	double t = (double)getTickCount();

	// detect

	t = (double)getTickCount();

	detector_Freak.detect(img1_cvLoadImage,keys1);

	detector_Freak.detect(img2_cvLoadImage, keys2);

	t = ((double)getTickCount() - t)/getTickFrequency();

	cout << "FREAK detection time [s]: " << t/1.0 << endl;



	// extract

	t = (double)getTickCount();

	freak.compute(img1_cvLoadImage, keys1, descriptors1);

	freak.compute(img2_cvLoadImage, keys2, descriptors2);

	t = ((double)getTickCount() - t)/getTickFrequency();

	cout << "FREAK extraction time [s]: " << t << endl;

	matcher.match(descriptors1, descriptors2, matches);

	

	double max_dist = 0;

	double min_dist = 100;

	

	for (int i=0; i<descriptors1.rows; i++)

	{ 

		double dist = matches[i].distance;

		if (dist < min_dist) min_dist = dist;

		if(dist > max_dist) max_dist = dist;

	}

	printf("-- Max dist : %f \n", max_dist);

	printf("-- Min dist : %f \n", min_dist);

	

	vector<DMatch> good_matches_Freak;

	for (int i=0; i<descriptors1.rows; i++)

	{ 

		if(matches[i].distance < 0.8*max_dist)

		{ 

			good_matches_Freak.push_back(matches[i]); 

		}

	}

	

	Mat imgMatch_Freak;

	drawMatches(img1_cvLoadImage, keys1, img2_cvLoadImage, keys2, good_matches_Freak, imgMatch_Freak,

		Scalar::all(-1), Scalar::all(-1),

		vector<char>(), DrawMatchesFlags::NOT_DRAW_SINGLE_POINTS);



	// display

	imshow("matchFREAK", imgMatch_Freak);



	

	//waitKey();

	

	vector<Point2f>	pt1,pt2;   //定义点容器，存储匹配的那些特征点

	pt1.resize(good_matches_Freak.size());

	pt2.resize(good_matches_Freak.size());//设置断点于此，前面运行正确

	for(int i = 0; i < (int)good_matches_Freak.size(); i++){

		pt1.push_back(Point2f(keys2[good_matches_Freak[i].queryIdx].pt.x, keys2[good_matches_Freak[i].queryIdx].pt.y));    



		pt2.push_back(Point2f(keys1[good_matches_Freak[i].trainIdx].pt.x, keys1[good_matches_Freak[i].trainIdx].pt.y));

	}





	Mat homo;



	//st = cvGetTickCount();

	homo = findHomography(pt1, pt2,CV_RANSAC);   //寻找单应变换矩阵

	//et = cvGetTickCount();

	//printf("ransac time: %f\n", (et-st)/(double)cvGetTickFrequency()/1000.);   //ransac所花费的时间



	printf("homo\n"     //是一个3*3的变换矩阵

		"%f %f %f\n"

		"%f %f %f\n"

		"%f %f %f\n",

		homo.at<double>(0,0), homo.at<double>(0,1), homo.at<double>(0,2),

		homo.at<double>(1,0), homo.at<double>(1,1), homo.at<double>(1,2),

		homo.at<double>(2,0), homo.at<double>(2,1), homo.at<double>(2,2));