Понарамное фото с использованием OpenCV

Используя понятие однородных координат и матрицы гомографии попробуем написать не сложное приложение используя библиотеку компьютерного зрения OpenCV, которое из ряда последовательных фотографий будет строить одну большую фотографию:

С начала будет дан необходимый теоретический минимум, а затем некоторые примеры кода на С++.
Для того что бы склеить две фотографии необходимо выполнить следующие действия:

1. Для каждого из двух изображений найти особые точки (о том что это такое ниже) и их дескрипторы.
2. Сопоставление особых точек двух изображений.
3. Отсеивание ложных совпадений точек.
4. Построение матрицы гомографии. 
5. Проецирование с помощью полученной матрицы гомографии одного изображения на другое.

Теперь подробнее по каждому пункту.

1. Для склейки двух фотографий, нам вообще говоря, надо как-то определить, что за фото перед нами и из этой информации сделать вывод, можно ли их склеить (вдруг перед нами фото молекулы водорода и боинга 777 и склеивать нечего) и если можно, то в каких местах их склеивать. Для этой цели используются некоторые признаки изображений. Одним из таких признаков являются особые точки изображения. Иначе говоря, особые точки, это точки которые каким-то образом характеризующие данное изображение. Каждая особая точка имеет некоторый параметр характеризующей ее (по мимо координат). Этот параметр называется дескриптором особой точки. 

Есть несколько алгоритмов решающих задачу нахождения особых точек изображения. Те что я знаю, это алгоритм Speeded Up Robust Features (SURF) и Scale Invariant Feature Transform (SIFT). Описанием этих алгоритмов я заниматься не буду. Достаточно хорошо оба алгоритма описаны на хабре (SURF, SIFT). Для выявления особых точек в нашей программе будем использовать SURF, хотя можно было бы использовать и SIFT.

2. Из первого шага мы имеем особые точки первого изображения и особые точки второго изображения. Т.к. особые точки являются характеристикой изображения, мы как раз и будем делать вывод о возможностях склеивания изображений исходя из сравнения этих характеристик. Т.е. на данном шаге будем сравнивать дескрипторы особых точек.

3. Из предыдущего шага мы получили набор особых точек которые сравнивая по дескрипторам являются общими для двух изображений. В этом наборе точек могут оказаться ложные совпадения. Т.е. может оказаться что дескрипторы точек совпадают, но на самом деле эти точки принадлежат совершенно разным частям изображений. Этих точек будет совсем не много в сравнении с "правильно" совпавшим точками, но нам их надо каким-то образом отсеять. Для этого мы будем использовать еще один алгоритм который называется RANdom SAmple Consensus (RANSAC), который поможет отсеять ложно совпавшие точки. 

4. На этом шаге мы построим матрицу гомографии для двух изображений.

5. И наконец, умножая каждую точку с первого изображения на матрицу гомографии, мы получим получим это изображение в плоскости второго изображения с совпадением общих особых точек двух изображений.

Давайте теперь перейдем к практике. В программе будет использоваться библиотека компьютерного зрения OpenCV, которая уже содержит в себе все необходимые для программы алгоритмы и структуры для работы с изображениями. (Как пользоваться OpenCV можно ознакомится на RoboCraft). И еще, программа будет рассчитана на любое количество изображений, склеивая каждую следующую фотографию с результатом склейки двух предыдущих. Начнем с функции main:

1. int main(int argc, char *argv[])
2. {
3.         if (argc < 3) return -1;
4.  
5.         for (int i = 1; i < argc; ++i) {
6.                 cvNamedWindow(argv[i]);
7.                 images.push_back(cvLoadImage(argv[i]));
8.                 cvShowImage(argv[i], images.at(i-1));
9.         }
10.  
11.         IplImage *result;
12.  
13.         std::cout << "Make a panorama from images:n";
14.         for (int i = 0; i < images.size() - 1; ++i)
15.         {
16.                 if (i > 0) images.at(i) = result;
17.                 remapping(images.at(i), images.at(i+1), &result);
18.         }
19.  
20.         cvNamedWindow("win");
21.         cvShowImage("win", result);
22.  
23.         cvWaitKey();
24.  
25.         for (int i = 0; i < images.size(); ++i)
26.                 cvReleaseImage(&images.at(i));
27.  
28.         cvReleaseImage(&result);
29.         cvDestroyAllWindows();
30.  
31.         return 0;
32. }

В качестве входных параметров программа принимает имена файлов изображений (при условии что файлы находятся в той же папке что и  исполняемый файл) в том порядке котором они будут склеиваться. Каждое изображение будет перед склейкой показано в отельном окне (6 строка), а так же загружено в вектор (7). Вывод изображений в 8й строке. В 11й строке объявляется структура которая будет содержать результат склейки. Функция remapping в строке 17 склеивает два изображения. Результат склейки всех изображений будет выведен на экран в строке 21. В конце освобождается память.
Рассмотрим функцию remapping:

1. int remapping(IplImage *img1, IplImage *img2, IplImage **result)
2. {
3.         // features points
4.         CvSeq *keyPoints1, *keyPoints2;
5.         CvSeq *desPoints1, *desPoints2;
6.         // for cvExtractSURF
7.         CvMemStorage *storage = cvCreateMemStorage();
8.         CvSURFParams params = cvSURFParams(500,1);
9.  
10.         std::vector<int> ptPairs;
11.  
12.         std::cout << "Remapping...n";
13.  
14.  
15.         cvExtractSURF(img1, 0, &keyPoints1, &desPoints1, storage, params);
16.         cvExtractSURF(img2, 0, &keyPoints2, &desPoints2, storage, params);
17.  
18.         flannFindPairs(keyPoints1, desPoints1, keyPoints2, desPoints2, ptPairs);
19.  
20.         int n = (int)(ptPairs.size()/2);
21.         if (n < 4) return 0;
22.  
23.         std::vector<CvPoint2D32f> pt1, pt2;
24.         pt1.resize(n);  pt2.resize(n);
25.  
26.         for(int j = 0; j < n; ++j) {
27.                 pt1[j] = ((CvSURFPoint*)cvGetSeqElem(keyPoints1,ptPairs[j*2]))->pt;
28.                 pt2[j] = ((CvSURFPoint*)cvGetSeqElem(keyPoints2,ptPairs[j*2+1]))->pt;
29.         }
30.  
31.         srcPoints = cv::Mat(1, n, CV_32FC2, &pt1[0]);
32.         disPoints = cv::Mat(1, n, CV_32FC2, &pt2[0]);
33.         // Calculating matrix of homography
34.         H = cv::findHomography(srcPoints, disPoints, CV_RANSAC, 5);
35.  
36.         std::cout << "Perspective transformation H between two images:n"; 
37.         std::cout << "H = n" << H << "n";
38.  
39.         int maxX, maxY, minX, minY;
40.         int maxX2, maxY2, minX2, minY2;
41.         // Define image's size
42.         ImageSizePoint(img1, &H, &minX, &minY, &maxX, &maxY);
43.         ImageSizePoint(img2, 0, &minX2, &minY2, &maxX2, &maxY2);
44.         // Select the most minimal
45.         minX = std::min(minX, minX2); minY = std::min(minY, minY2);
46.         // if the image is out of the window border             // the move images
47.         if (minX < 0) minX = std::abs(minX);
48.         else    minX = 0;
49.         if (minY < 0) minY = std::abs(minY);
50.         else    minY = 0;
51.         globalMaxX = std::max(maxX, maxX2);
52.         globalMaxY = std::max(maxY, maxY2);
53.         // Create final image
54.         *result = cvCreateImage(cvSize(  std::max(maxX,maxX2) + minX,
55.                                         std::max(maxY,maxY2) + minY)
56.                                         ,img2->depth,img2->nChannels);
57.  
58.         // H^(-1)
59.         H=H.inv();
60.  
61.         CvScalar s;
62.         // Display original 2nd image
63.         for (int i = 0; i < img2->height; ++i)
64.                 for (int j = 0; j < img2->width; ++j)
65.                 {
66.                         s = cvGet2D(img2, i, j);
67.                         cvSet2D(*result, i+minY, j+minX, s);
68.                 }
69.  
70.         vector<double> x;
71.         // Remapping 1st image
72.         for (int i = -minY; i < ((*result)->height - minY); ++i)
73.                 for (int j = -minX; j < ((*result)->width - minX); ++j)
74.                 {
75.                         // Compute point's coordinates at remapping image
76.                         // x = H^(-1)*x'
77.                         x = matrixMultiplication(&H, j, i);
78.  
79.                         // If does computed point exists at 1st image
80.                         // note: x(0)=j' x(1)=i'
81.                         if (    x(1)>=0 && x(1)<img1->height &&
82.                                 x(0)>=0 && x(0)<img1->width)
83.                         {
84.                                 s = cvGet2D(img1, x(1), x(0));
85.                                 cvSet2D(*result, i+minY, j+minX, s);
86.                         }
87.                 }
88.         return 1;
89. }

Все волшебство выполняется именно в этой функции. В качестве параметров принимается два изображения для которых в строках 15 и 16 с помощью библиотечных функций OpenCV ищутся особые точки и их дескрипторы по алгоритму SURF. В строке 18 с помощью функции flannFindPairs сравниваются полученных ключевые точки. Если общих точек оказалось совсем мало, то в строке 21 мы завершаем работу функции, т.е. данные два изображения разные и склеить их не удастся.  В строке 34, мы находим матрицу гомографии с помощью библиотечной функции findHomography. Обратите внимание, что в параметрах указанно значение CV_RANSAC, т.е. для поиска матрицы гомографии будет использоваться алгоритм RANSAC. В строках с 39 по 56 мы вычисляем размеры конечного изображения и создаем соответствующего размера окно.
По скольку мы собираемся проецировать первое изображение в плоскость второго нам необходимо найти обратную матрицу к матрице гомографии, что и выполняется в строке 59.
Второе изображение оставляем исходным, т.е. мы его просто копируем в результирующее изображение (62-68), а вот первое копируем, но уже умножая на обратную матрицу гомографии (72-87). При этом надо отметить, что копирования и первого и второго изображения происходит с некоторым сдвигом. Это выполняется для этого что бы все изображение поместилось в результирующее окно. Вот собственно и все.
Так же в программе используются вспомогательные функции:

1. void flannFindPairs(...)
2. void ImageSizePoint(...)
3. vector<double> matrixMultiplication(...)

Первая функция сравнивает наборы особых точек двух изображений по дескрипторам. Вторая функция используется для вычисления размеров результирующего изображения. И наконец, третья функция представляет из себя простое матричное произведение. Результат работы программы можно посмотреть в начале текста, а исходный код программы можно скачать тут.