1/*M///////////////////////////////////////////////////////////////////////////////////////
2//
3//  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4//
5//  By downloading, copying, installing or using the software you agree to this license.
6//  If you do not agree to this license, do not download, install,
7//  copy or use the software.
8//
9//
10//                           License Agreement
11//                For Open Source Computer Vision Library
12//
13// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
14// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
15// Third party copyrights are property of their respective owners.
16//
17// Redistribution and use in source and binary forms, with or without modification,
18// are permitted provided that the following conditions are met:
19//
20//   * Redistribution's of source code must retain the above copyright notice,
21//     this list of conditions and the following disclaimer.
22//
23//   * Redistribution's in binary form must reproduce the above copyright notice,
24//     this list of conditions and the following disclaimer in the documentation
25//     and/or other materials provided with the distribution.
26//
27//   * The name of the copyright holders may not be used to endorse or promote products
28//     derived from this software without specific prior written permission.
29//
30// This software is provided by the copyright holders and contributors "as is" and
31// any express or implied warranties, including, but not limited to, the implied
32// warranties of merchantability and fitness for a particular purpose are disclaimed.
33// In no event shall the Intel Corporation or contributors be liable for any direct,
34// indirect, incidental, special, exemplary, or consequential damages
35// (including, but not limited to, procurement of substitute goods or services;
36// loss of use, data, or profits; or business interruption) however caused
37// and on any theory of liability, whether in contract, strict liability,
38// or tort (including negligence or otherwise) arising in any way out of
39// the use of this software, even if advised of the possibility of such damage.
40//
41//M*/
42
43#include "test_precomp.hpp"
44#include <time.h>
45
46#define CALIB3D_HOMOGRAPHY_ERROR_MATRIX_SIZE 1
47#define CALIB3D_HOMOGRAPHY_ERROR_MATRIX_DIFF 2
48#define CALIB3D_HOMOGRAPHY_ERROR_REPROJ_DIFF 3
49#define CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK 4
50#define CALIB3D_HOMOGRAPHY_ERROR_RANSAC_DIFF 5
51
52#define MESSAGE_MATRIX_SIZE "Homography matrix must have 3*3 sizes."
53#define MESSAGE_MATRIX_DIFF "Accuracy of homography transformation matrix less than required."
54#define MESSAGE_REPROJ_DIFF_1 "Reprojection error for current pair of points more than required."
55#define MESSAGE_REPROJ_DIFF_2 "Reprojection error is not optimal."
56#define MESSAGE_RANSAC_MASK_1 "Sizes of inliers/outliers mask are incorrect."
57#define MESSAGE_RANSAC_MASK_2 "Mask mustn't have any outliers."
58#define MESSAGE_RANSAC_MASK_3 "All values of mask must be 1 (true) or 0 (false)."
59#define MESSAGE_RANSAC_MASK_4 "Mask of inliers/outliers is incorrect."
60#define MESSAGE_RANSAC_MASK_5 "Inlier in original mask shouldn't be outlier in found mask."
61#define MESSAGE_RANSAC_DIFF "Reprojection error for current pair of points more than required."
62
63#define MAX_COUNT_OF_POINTS 303
64#define COUNT_NORM_TYPES 3
65#define METHODS_COUNT 4
66
67int NORM_TYPE[COUNT_NORM_TYPES] = {cv::NORM_L1, cv::NORM_L2, cv::NORM_INF};
68int METHOD[METHODS_COUNT] = {0, cv::RANSAC, cv::LMEDS, cv::RHO};
69
70using namespace cv;
71using namespace std;
72
73class CV_HomographyTest: public cvtest::ArrayTest
74{
75public:
76    CV_HomographyTest();
77    ~CV_HomographyTest();
78
79    void run (int);
80
81protected:
82
83    int method;
84    int image_size;
85    double reproj_threshold;
86    double sigma;
87
88private:
89    float max_diff, max_2diff;
90    bool check_matrix_size(const cv::Mat& H);
91    bool check_matrix_diff(const cv::Mat& original, const cv::Mat& found, const int norm_type, double &diff);
92    int check_ransac_mask_1(const Mat& src, const Mat& mask);
93    int check_ransac_mask_2(const Mat& original_mask, const Mat& found_mask);
94
95    void print_information_1(int j, int N, int method, const Mat& H);
96    void print_information_2(int j, int N, int method, const Mat& H, const Mat& H_res, int k, double diff);
97    void print_information_3(int method, int j, int N, const Mat& mask);
98    void print_information_4(int method, int j, int N, int k, int l, double diff);
99    void print_information_5(int method, int j, int N, int l, double diff);
100    void print_information_6(int method, int j, int N, int k, double diff, bool value);
101    void print_information_7(int method, int j, int N, int k, double diff, bool original_value, bool found_value);
102    void print_information_8(int method, int j, int N, int k, int l, double diff);
103};
104
105CV_HomographyTest::CV_HomographyTest() : max_diff(1e-2f), max_2diff(2e-2f)
106{
107    method = 0;
108    image_size = 100;
109    reproj_threshold = 3.0;
110    sigma = 0.01;
111}
112
113CV_HomographyTest::~CV_HomographyTest() {}
114
115bool CV_HomographyTest::check_matrix_size(const cv::Mat& H)
116{
117    return (H.rows == 3) && (H.cols == 3);
118}
119
120bool CV_HomographyTest::check_matrix_diff(const cv::Mat& original, const cv::Mat& found, const int norm_type, double &diff)
121{
122    diff = cvtest::norm(original, found, norm_type);
123    return diff <= max_diff;
124}
125
126int CV_HomographyTest::check_ransac_mask_1(const Mat& src, const Mat& mask)
127{
128    if (!(mask.cols == 1) && (mask.rows == src.cols)) return 1;
129    if (countNonZero(mask) < mask.rows) return 2;
130    for (int i = 0; i < mask.rows; ++i) if (mask.at<uchar>(i, 0) > 1) return 3;
131    return 0;
132}
133
134int CV_HomographyTest::check_ransac_mask_2(const Mat& original_mask, const Mat& found_mask)
135{
136    if (!(found_mask.cols == 1) && (found_mask.rows == original_mask.rows)) return 1;
137    for (int i = 0; i < found_mask.rows; ++i) if (found_mask.at<uchar>(i, 0) > 1) return 2;
138    return 0;
139}
140
141void CV_HomographyTest::print_information_1(int j, int N, int _method, const Mat& H)
142{
143    cout << endl; cout << "Checking for homography matrix sizes..." << endl; cout << endl;
144    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
145    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
146    cout << "Count of points: " << N << endl; cout << endl;
147    cout << "Method: "; if (_method == 0) cout << 0; else if (_method == 8) cout << "RANSAC"; else if (_method == cv::RHO) cout << "RHO"; else cout << "LMEDS"; cout << endl;
148    cout << "Homography matrix:" << endl; cout << endl;
149    cout << H << endl; cout << endl;
150    cout << "Number of rows: " << H.rows << "   Number of cols: " << H.cols << endl; cout << endl;
151}
152
153void CV_HomographyTest::print_information_2(int j, int N, int _method, const Mat& H, const Mat& H_res, int k, double diff)
154{
155    cout << endl; cout << "Checking for accuracy of homography matrix computing..." << endl; cout << endl;
156    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
157    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
158    cout << "Count of points: " << N << endl; cout << endl;
159    cout << "Method: "; if (_method == 0) cout << 0; else if (_method == 8) cout << "RANSAC"; else if (_method == cv::RHO) cout << "RHO"; else cout << "LMEDS"; cout << endl;
160    cout << "Original matrix:" << endl; cout << endl;
161    cout << H << endl; cout << endl;
162    cout << "Found matrix:" << endl; cout << endl;
163    cout << H_res << endl; cout << endl;
164    cout << "Norm type using in criteria: "; if (NORM_TYPE[k] == 1) cout << "INF"; else if (NORM_TYPE[k] == 2) cout << "L1"; else cout << "L2"; cout << endl;
165    cout << "Difference between matrices: " << diff << endl;
166    cout << "Maximum allowed difference: " << max_diff << endl; cout << endl;
167}
168
169void CV_HomographyTest::print_information_3(int _method, int j, int N, const Mat& mask)
170{
171    cout << endl; cout << "Checking for inliers/outliers mask..." << endl; cout << endl;
172    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
173    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
174    cout << "Count of points: " << N << endl; cout << endl;
175    cout << "Method: "; if (_method == RANSAC) cout << "RANSAC" << endl; else if (_method == cv::RHO) cout << "RHO" << endl; else cout << _method << endl;
176    cout << "Found mask:" << endl; cout << endl;
177    cout << mask << endl; cout << endl;
178    cout << "Number of rows: " << mask.rows << "   Number of cols: " << mask.cols << endl; cout << endl;
179}
180
181void CV_HomographyTest::print_information_4(int _method, int j, int N, int k, int l, double diff)
182{
183    cout << endl; cout << "Checking for accuracy of reprojection error computing..." << endl; cout << endl;
184    cout << "Method: "; if (_method == 0) cout << 0 << endl; else cout << "CV_LMEDS" << endl;
185    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
186    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
187    cout << "Sigma of normal noise: " << sigma << endl;
188    cout << "Count of points: " << N << endl;
189    cout << "Number of point: " << k << endl;
190    cout << "Norm type using in criteria: "; if (NORM_TYPE[l] == 1) cout << "INF"; else if (NORM_TYPE[l] == 2) cout << "L1"; else cout << "L2"; cout << endl;
191    cout << "Difference with noise of point: " << diff << endl;
192    cout << "Maxumum allowed difference: " << max_2diff << endl; cout << endl;
193}
194
195void CV_HomographyTest::print_information_5(int _method, int j, int N, int l, double diff)
196{
197    cout << endl; cout << "Checking for accuracy of reprojection error computing..." << endl; cout << endl;
198    cout << "Method: "; if (_method == 0) cout << 0 << endl; else cout << "CV_LMEDS" << endl;
199    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
200    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
201    cout << "Sigma of normal noise: " << sigma << endl;
202    cout << "Count of points: " << N << endl;
203    cout << "Norm type using in criteria: "; if (NORM_TYPE[l] == 1) cout << "INF"; else if (NORM_TYPE[l] == 2) cout << "L1"; else cout << "L2"; cout << endl;
204    cout << "Difference with noise of points: " << diff << endl;
205    cout << "Maxumum allowed difference: " << max_diff << endl; cout << endl;
206}
207
208void CV_HomographyTest::print_information_6(int _method, int j, int N, int k, double diff, bool value)
209{
210    cout << endl; cout << "Checking for inliers/outliers mask..." << endl; cout << endl;
211    cout << "Method: "; if (_method == RANSAC) cout << "RANSAC" << endl; else if (_method == cv::RHO) cout << "RHO" << endl; else cout << _method << endl;
212    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
213    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
214    cout << "Count of points: " << N << "   " << endl;
215    cout << "Number of point: " << k << "   " << endl;
216    cout << "Reprojection error for this point: " << diff << "   " << endl;
217    cout << "Reprojection error threshold: " << reproj_threshold << "   " << endl;
218    cout << "Value of found mask: "<< value << endl; cout << endl;
219}
220
221void CV_HomographyTest::print_information_7(int _method, int j, int N, int k, double diff, bool original_value, bool found_value)
222{
223    cout << endl; cout << "Checking for inliers/outliers mask..." << endl; cout << endl;
224    cout << "Method: "; if (_method == RANSAC) cout << "RANSAC" << endl; else if (_method == cv::RHO) cout << "RHO" << endl; else cout << _method << endl;
225    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
226    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
227    cout << "Count of points: " << N << "   " << endl;
228    cout << "Number of point: " << k << "   " << endl;
229    cout << "Reprojection error for this point: " << diff << "   " << endl;
230    cout << "Reprojection error threshold: " << reproj_threshold << "   " << endl;
231    cout << "Value of original mask: "<< original_value << "   Value of found mask: " << found_value << endl; cout << endl;
232}
233
234void CV_HomographyTest::print_information_8(int _method, int j, int N, int k, int l, double diff)
235{
236    cout << endl; cout << "Checking for reprojection error of inlier..." << endl; cout << endl;
237    cout << "Method: "; if (_method == RANSAC) cout << "RANSAC" << endl; else if (_method == cv::RHO) cout << "RHO" << endl; else cout << _method << endl;
238    cout << "Sigma of normal noise: " << sigma << endl;
239    cout << "Type of srcPoints: "; if ((j>-1) && (j<2)) cout << "Mat of CV_32FC2"; else  cout << "vector <Point2f>";
240    cout << "   Type of dstPoints: "; if (j % 2 == 0) cout << "Mat of CV_32FC2"; else cout << "vector <Point2f>"; cout << endl;
241    cout << "Count of points: " << N << "   " << endl;
242    cout << "Number of point: " << k << "   " << endl;
243    cout << "Norm type using in criteria: "; if (NORM_TYPE[l] == 1) cout << "INF"; else if (NORM_TYPE[l] == 2) cout << "L1"; else cout << "L2"; cout << endl;
244    cout << "Difference with noise of point: " << diff << endl;
245    cout << "Maxumum allowed difference: " << max_2diff << endl; cout << endl;
246}
247
248void CV_HomographyTest::run(int)
249{
250    for (int N = 4; N <= MAX_COUNT_OF_POINTS; ++N)
251    {
252        RNG& rng = ts->get_rng();
253
254        float *src_data = new float [2*N];
255
256        for (int i = 0; i < N; ++i)
257        {
258            src_data[2*i] = (float)cvtest::randReal(rng)*image_size;
259            src_data[2*i+1] = (float)cvtest::randReal(rng)*image_size;
260        }
261
262        cv::Mat src_mat_2f(1, N, CV_32FC2, src_data),
263        src_mat_2d(2, N, CV_32F, src_data),
264        src_mat_3d(3, N, CV_32F);
265        cv::Mat dst_mat_2f, dst_mat_2d, dst_mat_3d;
266
267        vector <Point2f> src_vec, dst_vec;
268
269        for (int i = 0; i < N; ++i)
270        {
271            float *tmp = src_mat_2d.ptr<float>()+2*i;
272            src_mat_3d.at<float>(0, i) = tmp[0];
273            src_mat_3d.at<float>(1, i) = tmp[1];
274            src_mat_3d.at<float>(2, i) = 1.0f;
275
276            src_vec.push_back(Point2f(tmp[0], tmp[1]));
277        }
278
279        double fi = cvtest::randReal(rng)*2*CV_PI;
280
281        double t_x = cvtest::randReal(rng)*sqrt(image_size*1.0),
282        t_y = cvtest::randReal(rng)*sqrt(image_size*1.0);
283
284        double Hdata[9] = { cos(fi), -sin(fi), t_x,
285                            sin(fi),  cos(fi), t_y,
286                            0.0f,     0.0f, 1.0f };
287
288        cv::Mat H_64(3, 3, CV_64F, Hdata), H_32;
289
290        H_64.convertTo(H_32, CV_32F);
291
292        dst_mat_3d = H_32*src_mat_3d;
293
294        dst_mat_2d.create(2, N, CV_32F); dst_mat_2f.create(1, N, CV_32FC2);
295
296        for (int i = 0; i < N; ++i)
297        {
298            float *tmp_2f = dst_mat_2f.ptr<float>()+2*i;
299            tmp_2f[0] = dst_mat_2d.at<float>(0, i) = dst_mat_3d.at<float>(0, i) /= dst_mat_3d.at<float>(2, i);
300            tmp_2f[1] = dst_mat_2d.at<float>(1, i) = dst_mat_3d.at<float>(1, i) /= dst_mat_3d.at<float>(2, i);
301            dst_mat_3d.at<float>(2, i) = 1.0f;
302
303            dst_vec.push_back(Point2f(tmp_2f[0], tmp_2f[1]));
304        }
305
306        for (int i = 0; i < METHODS_COUNT; ++i)
307        {
308            method = METHOD[i];
309            switch (method)
310            {
311            case 0:
312            case LMEDS:
313                {
314                    Mat H_res_64 [4] = { cv::findHomography(src_mat_2f, dst_mat_2f, method),
315                                         cv::findHomography(src_mat_2f, dst_vec, method),
316                                         cv::findHomography(src_vec, dst_mat_2f, method),
317                                         cv::findHomography(src_vec, dst_vec, method) };
318
319                    for (int j = 0; j < 4; ++j)
320                    {
321
322                        if (!check_matrix_size(H_res_64[j]))
323                        {
324                            print_information_1(j, N, method, H_res_64[j]);
325                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_MATRIX_SIZE, MESSAGE_MATRIX_SIZE);
326                            return;
327                        }
328
329                        double diff;
330
331                        for (int k = 0; k < COUNT_NORM_TYPES; ++k)
332                            if (!check_matrix_diff(H_64, H_res_64[j], NORM_TYPE[k], diff))
333                            {
334                            print_information_2(j, N, method, H_64, H_res_64[j], k, diff);
335                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_MATRIX_DIFF, MESSAGE_MATRIX_DIFF);
336                            return;
337                        }
338                    }
339
340                    continue;
341                }
342            case cv::RHO:
343            case RANSAC:
344                {
345                    cv::Mat mask [4]; double diff;
346
347                    Mat H_res_64 [4] = { cv::findHomography(src_mat_2f, dst_mat_2f, method, reproj_threshold, mask[0]),
348                                         cv::findHomography(src_mat_2f, dst_vec, method, reproj_threshold, mask[1]),
349                                         cv::findHomography(src_vec, dst_mat_2f, method, reproj_threshold, mask[2]),
350                                         cv::findHomography(src_vec, dst_vec, method, reproj_threshold, mask[3]) };
351
352                    for (int j = 0; j < 4; ++j)
353                    {
354
355                        if (!check_matrix_size(H_res_64[j]))
356                        {
357                            print_information_1(j, N, method, H_res_64[j]);
358                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_MATRIX_SIZE, MESSAGE_MATRIX_SIZE);
359                            return;
360                        }
361
362                        for (int k = 0; k < COUNT_NORM_TYPES; ++k)
363                            if (!check_matrix_diff(H_64, H_res_64[j], NORM_TYPE[k], diff))
364                            {
365                            print_information_2(j, N, method, H_64, H_res_64[j], k, diff);
366                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_MATRIX_DIFF, MESSAGE_MATRIX_DIFF);
367                            return;
368                        }
369
370                        int code = check_ransac_mask_1(src_mat_2f, mask[j]);
371
372                        if (code)
373                        {
374                            print_information_3(method, j, N, mask[j]);
375
376                            switch (code)
377                            {
378                            case 1: { CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_1); break; }
379                            case 2: { CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_2); break; }
380                            case 3: { CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_3); break; }
381
382                            default: break;
383                            }
384
385                            return;
386                        }
387
388                    }
389
390                    continue;
391                }
392
393            default: continue;
394            }
395        }
396
397        Mat noise_2f(1, N, CV_32FC2);
398        rng.fill(noise_2f, RNG::NORMAL, Scalar::all(0), Scalar::all(sigma));
399
400        cv::Mat mask(N, 1, CV_8UC1);
401
402        for (int i = 0; i < N; ++i)
403        {
404            float *a = noise_2f.ptr<float>()+2*i, *_2f = dst_mat_2f.ptr<float>()+2*i;
405            _2f[0] += a[0]; _2f[1] += a[1];
406            mask.at<bool>(i, 0) = !(sqrt(a[0]*a[0]+a[1]*a[1]) > reproj_threshold);
407        }
408
409        for (int i = 0; i < METHODS_COUNT; ++i)
410        {
411            method = METHOD[i];
412            switch (method)
413            {
414            case 0:
415            case LMEDS:
416                {
417                    Mat H_res_64 [4] = { cv::findHomography(src_mat_2f, dst_mat_2f),
418                                         cv::findHomography(src_mat_2f, dst_vec),
419                                         cv::findHomography(src_vec, dst_mat_2f),
420                                         cv::findHomography(src_vec, dst_vec) };
421
422                    for (int j = 0; j < 4; ++j)
423                    {
424
425                        if (!check_matrix_size(H_res_64[j]))
426                        {
427                            print_information_1(j, N, method, H_res_64[j]);
428                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_MATRIX_SIZE, MESSAGE_MATRIX_SIZE);
429                            return;
430                        }
431
432                        Mat H_res_32; H_res_64[j].convertTo(H_res_32, CV_32F);
433
434                        cv::Mat dst_res_3d(3, N, CV_32F), noise_2d(2, N, CV_32F);
435
436                        for (int k = 0; k < N; ++k)
437                        {
438
439                            Mat tmp_mat_3d = H_res_32*src_mat_3d.col(k);
440
441                            dst_res_3d.at<float>(0, k) = tmp_mat_3d.at<float>(0, 0) /= tmp_mat_3d.at<float>(2, 0);
442                            dst_res_3d.at<float>(1, k) = tmp_mat_3d.at<float>(1, 0) /= tmp_mat_3d.at<float>(2, 0);
443                            dst_res_3d.at<float>(2, k) = tmp_mat_3d.at<float>(2, 0) = 1.0f;
444
445                            float *a = noise_2f.ptr<float>()+2*k;
446                            noise_2d.at<float>(0, k) = a[0]; noise_2d.at<float>(1, k) = a[1];
447
448                            for (int l = 0; l < COUNT_NORM_TYPES; ++l)
449                                if (cv::norm(tmp_mat_3d, dst_mat_3d.col(k), NORM_TYPE[l]) - cv::norm(noise_2d.col(k), NORM_TYPE[l]) > max_2diff)
450                                {
451                                print_information_4(method, j, N, k, l, cv::norm(tmp_mat_3d, dst_mat_3d.col(k), NORM_TYPE[l]) - cv::norm(noise_2d.col(k), NORM_TYPE[l]));
452                                CV_Error(CALIB3D_HOMOGRAPHY_ERROR_REPROJ_DIFF, MESSAGE_REPROJ_DIFF_1);
453                                return;
454                            }
455
456                        }
457
458                        for (int l = 0; l < COUNT_NORM_TYPES; ++l)
459                            if (cv::norm(dst_res_3d, dst_mat_3d, NORM_TYPE[l]) - cv::norm(noise_2d, NORM_TYPE[l]) > max_diff)
460                            {
461                            print_information_5(method, j, N, l, cv::norm(dst_res_3d, dst_mat_3d, NORM_TYPE[l]) - cv::norm(noise_2d, NORM_TYPE[l]));
462                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_REPROJ_DIFF, MESSAGE_REPROJ_DIFF_2);
463                            return;
464                        }
465
466                    }
467
468                    continue;
469                }
470            case cv::RHO:
471            case RANSAC:
472                {
473                    cv::Mat mask_res [4];
474
475                    Mat H_res_64 [4] = { cv::findHomography(src_mat_2f, dst_mat_2f, method, reproj_threshold, mask_res[0]),
476                                         cv::findHomography(src_mat_2f, dst_vec, method, reproj_threshold, mask_res[1]),
477                                         cv::findHomography(src_vec, dst_mat_2f, method, reproj_threshold, mask_res[2]),
478                                         cv::findHomography(src_vec, dst_vec, method, reproj_threshold, mask_res[3]) };
479
480                    for (int j = 0; j < 4; ++j)
481                    {
482                        if (!check_matrix_size(H_res_64[j]))
483                        {
484                            print_information_1(j, N, method, H_res_64[j]);
485                            CV_Error(CALIB3D_HOMOGRAPHY_ERROR_MATRIX_SIZE, MESSAGE_MATRIX_SIZE);
486                            return;
487                        }
488
489                        int code = check_ransac_mask_2(mask, mask_res[j]);
490
491                        if (code)
492                        {
493                            print_information_3(method, j, N, mask_res[j]);
494
495                            switch (code)
496                            {
497                            case 1: { CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_1); break; }
498                            case 2: { CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_3); break; }
499
500                            default: break;
501                            }
502
503                            return;
504                        }
505
506                        cv::Mat H_res_32; H_res_64[j].convertTo(H_res_32, CV_32F);
507
508                        cv::Mat dst_res_3d = H_res_32*src_mat_3d;
509
510                        for (int k = 0; k < N; ++k)
511                        {
512                            dst_res_3d.at<float>(0, k) /= dst_res_3d.at<float>(2, k);
513                            dst_res_3d.at<float>(1, k) /= dst_res_3d.at<float>(2, k);
514                            dst_res_3d.at<float>(2, k) = 1.0f;
515
516                            float *p = dst_mat_2f.ptr<float>()+2*k;
517
518                            dst_mat_3d.at<float>(0, k) = p[0];
519                            dst_mat_3d.at<float>(1, k) = p[1];
520
521                            double diff = cv::norm(dst_res_3d.col(k), dst_mat_3d.col(k), NORM_L2);
522
523                            if (mask_res[j].at<bool>(k, 0) != (diff <= reproj_threshold))
524                            {
525                                print_information_6(method, j, N, k, diff, mask_res[j].at<bool>(k, 0));
526                                CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_4);
527                                return;
528                            }
529
530                            if (mask.at<bool>(k, 0) && !mask_res[j].at<bool>(k, 0))
531                            {
532                                print_information_7(method, j, N, k, diff, mask.at<bool>(k, 0), mask_res[j].at<bool>(k, 0));
533                                CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_MASK, MESSAGE_RANSAC_MASK_5);
534                                return;
535                            }
536
537                            if (mask_res[j].at<bool>(k, 0))
538                            {
539                                float *a = noise_2f.ptr<float>()+2*k;
540                                dst_mat_3d.at<float>(0, k) -= a[0];
541                                dst_mat_3d.at<float>(1, k) -= a[1];
542
543                                cv::Mat noise_2d(2, 1, CV_32F);
544                                noise_2d.at<float>(0, 0) = a[0]; noise_2d.at<float>(1, 0) = a[1];
545
546                                for (int l = 0; l < COUNT_NORM_TYPES; ++l)
547                                {
548                                    diff = cv::norm(dst_res_3d.col(k), dst_mat_3d.col(k), NORM_TYPE[l]);
549
550                                    if (diff - cv::norm(noise_2d, NORM_TYPE[l]) > max_2diff)
551                                    {
552                                        print_information_8(method, j, N, k, l, diff - cv::norm(noise_2d, NORM_TYPE[l]));
553                                        CV_Error(CALIB3D_HOMOGRAPHY_ERROR_RANSAC_DIFF, MESSAGE_RANSAC_DIFF);
554                                        return;
555                                    }
556                                }
557                            }
558                        }
559                    }
560
561                    continue;
562                }
563
564            default: continue;
565            }
566        }
567    }
568}
569
570TEST(Calib3d_Homography, accuracy) { CV_HomographyTest test; test.safe_run(); }
571
572TEST(Calib3d_Homography, EKcase)
573{
574    float pt1data[] =
575    {
576        2.80073029e+002f, 2.39591217e+002f, 2.21912201e+002f, 2.59783997e+002f,
577        2.16053192e+002f, 2.78826569e+002f, 2.22782532e+002f, 2.82330383e+002f,
578        2.09924820e+002f, 2.89122559e+002f, 2.11077698e+002f, 2.89384674e+002f,
579        2.25287689e+002f, 2.88795532e+002f, 2.11180801e+002f, 2.89653503e+002f,
580        2.24126404e+002f, 2.90466064e+002f, 2.10914429e+002f, 2.90886963e+002f,
581        2.23439362e+002f, 2.91657715e+002f, 2.24809387e+002f, 2.91891602e+002f,
582        2.09809082e+002f, 2.92891113e+002f, 2.08771164e+002f, 2.93093231e+002f,
583        2.23160095e+002f, 2.93259460e+002f, 2.07874023e+002f, 2.93989990e+002f,
584        2.08963638e+002f, 2.94209839e+002f, 2.23963165e+002f, 2.94479645e+002f,
585        2.23241791e+002f, 2.94887817e+002f, 2.09438782e+002f, 2.95233337e+002f,
586        2.08901886e+002f, 2.95762878e+002f, 2.21867981e+002f, 2.95747711e+002f,
587        2.24195511e+002f, 2.98270905e+002f, 2.09331345e+002f, 3.05958191e+002f,
588        2.24727875e+002f, 3.07186035e+002f, 2.26718842e+002f, 3.08095795e+002f,
589        2.25363953e+002f, 3.08200226e+002f, 2.19897797e+002f, 3.13845093e+002f,
590        2.25013474e+002f, 3.15558777e+002f
591    };
592
593    float pt2data[] =
594    {
595        1.84072723e+002f, 1.43591202e+002f, 1.25912483e+002f, 1.63783859e+002f,
596        2.06439407e+002f, 2.20573929e+002f, 1.43801437e+002f, 1.80703903e+002f,
597        9.77904129e+000f, 2.49660202e+002f, 1.38458405e+001f, 2.14502701e+002f,
598        1.50636337e+002f, 2.15597183e+002f, 6.43103180e+001f, 2.51667648e+002f,
599        1.54952499e+002f, 2.20780014e+002f, 1.26638412e+002f, 2.43040924e+002f,
600        3.67568909e+002f, 1.83624954e+001f, 1.60657944e+002f, 2.21794052e+002f,
601        -1.29507828e+000f, 3.32472443e+002f, 8.51442242e+000f, 4.15561554e+002f,
602        1.27161377e+002f, 1.97260361e+002f, 5.40714645e+000f, 4.90978302e+002f,
603        2.25571690e+001f, 3.96912415e+002f, 2.95664978e+002f, 7.36064959e+000f,
604        1.27241104e+002f, 1.98887573e+002f, -1.25569367e+000f, 3.87713226e+002f,
605        1.04194012e+001f, 4.31495758e+002f, 1.25868874e+002f, 1.99751617e+002f,
606        1.28195480e+002f, 2.02270355e+002f, 2.23436356e+002f, 1.80489182e+002f,
607        1.28727692e+002f, 2.11185410e+002f, 2.03336639e+002f, 2.52182083e+002f,
608        1.29366486e+002f, 2.12201904e+002f, 1.23897598e+002f, 2.17847351e+002f,
609        1.29015259e+002f, 2.19560623e+002f
610    };
611
612    int npoints = (int)(sizeof(pt1data)/sizeof(pt1data[0])/2);
613
614    Mat p1(1, npoints, CV_32FC2, pt1data);
615    Mat p2(1, npoints, CV_32FC2, pt2data);
616    Mat mask;
617
618    Mat h = findHomography(p1, p2, RANSAC, 0.01, mask);
619    ASSERT_TRUE(!h.empty());
620
621    transpose(mask, mask);
622    Mat p3, mask2;
623    int ninliers = countNonZero(mask);
624    Mat nmask[] = { mask, mask };
625    merge(nmask, 2, mask2);
626    perspectiveTransform(p1, p3, h);
627    mask2 = mask2.reshape(1);
628    p2 = p2.reshape(1);
629    p3 = p3.reshape(1);
630    double err = norm(p2, p3, NORM_INF, mask2);
631
632    printf("ninliers: %d, inliers err: %.2g\n", ninliers, err);
633    ASSERT_GE(ninliers, 10);
634    ASSERT_LE(err, 0.01);
635}
636
637TEST(Calib3d_Homography, fromImages)
638{
639    Mat img_1 = imread(cvtest::TS::ptr()->get_data_path() + "cv/optflow/image1.png", 0);
640    Mat img_2 = imread(cvtest::TS::ptr()->get_data_path() + "cv/optflow/image2.png", 0);
641    Ptr<ORB> orb = ORB::create();
642    vector<KeyPoint> keypoints_1, keypoints_2;
643    Mat descriptors_1, descriptors_2;
644    orb->detectAndCompute( img_1, Mat(), keypoints_1, descriptors_1, false );
645    orb->detectAndCompute( img_2, Mat(), keypoints_2, descriptors_2, false );
646
647    //-- Step 3: Matching descriptor vectors using Brute Force matcher
648    BFMatcher  matcher(NORM_HAMMING,false);
649    std::vector< DMatch > matches;
650    matcher.match( descriptors_1, descriptors_2, matches );
651
652    double max_dist = 0; double min_dist = 100;
653    //-- Quick calculation of max and min distances between keypoints
654    for( int i = 0; i < descriptors_1.rows; i++ )
655    {
656        double dist = matches[i].distance;
657        if( dist < min_dist ) min_dist = dist;
658        if( dist > max_dist ) max_dist = dist;
659    }
660
661    //-- Draw only "good" matches (i.e. whose distance is less than 3*min_dist )
662    std::vector< DMatch > good_matches;
663    for( int i = 0; i < descriptors_1.rows; i++ )
664    {
665        if( matches[i].distance <= 100 )
666            good_matches.push_back( matches[i]);
667    }
668
669    //-- Localize the model
670    std::vector<Point2f> pointframe1;
671    std::vector<Point2f> pointframe2;
672    for( int i = 0; i < (int)good_matches.size(); i++ )
673    {
674        //-- Get the keypoints from the good matches
675        pointframe1.push_back( keypoints_1[ good_matches[i].queryIdx ].pt );
676        pointframe2.push_back( keypoints_2[ good_matches[i].trainIdx ].pt );
677    }
678
679    Mat H0, H1, inliers0, inliers1;
680    double min_t0 = DBL_MAX, min_t1 = DBL_MAX;
681    for( int i = 0; i < 10; i++ )
682    {
683        double t = (double)getTickCount();
684        H0 = findHomography( pointframe1, pointframe2, RANSAC, 3.0, inliers0 );
685        t = (double)getTickCount() - t;
686        min_t0 = std::min(min_t0, t);
687    }
688    int ninliers0 = countNonZero(inliers0);
689    for( int i = 0; i < 10; i++ )
690    {
691        double t = (double)getTickCount();
692        H1 = findHomography( pointframe1, pointframe2, RHO, 3.0, inliers1 );
693        t = (double)getTickCount() - t;
694        min_t1 = std::min(min_t1, t);
695    }
696    int ninliers1 = countNonZero(inliers1);
697    double freq = getTickFrequency();
698    printf("nfeatures1 = %d, nfeatures2=%d, matches=%d, ninliers(RANSAC)=%d, "
699           "time(RANSAC)=%.2fmsec, ninliers(RHO)=%d, time(RHO)=%.2fmsec\n",
700           (int)keypoints_1.size(), (int)keypoints_2.size(),
701           (int)good_matches.size(), ninliers0, min_t0*1000./freq, ninliers1, min_t1*1000./freq);
702
703    ASSERT_TRUE(!H0.empty());
704    ASSERT_GE(ninliers0, 80);
705    ASSERT_TRUE(!H1.empty());
706    ASSERT_GE(ninliers1, 80);
707}
708