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 46using namespace cv; 47using namespace std; 48 49#define sign(a) a > 0 ? 1 : a == 0 ? 0 : -1 50 51#define CORE_EIGEN_ERROR_COUNT 1 52#define CORE_EIGEN_ERROR_SIZE 2 53#define CORE_EIGEN_ERROR_DIFF 3 54#define CORE_EIGEN_ERROR_ORTHO 4 55#define CORE_EIGEN_ERROR_ORDER 5 56 57#define MESSAGE_ERROR_COUNT "Matrix of eigen values must have the same rows as source matrix and 1 column." 58#define MESSAGE_ERROR_SIZE "Source matrix and matrix of eigen vectors must have the same sizes." 59#define MESSAGE_ERROR_DIFF_1 "Accurasy of eigen values computing less than required." 60#define MESSAGE_ERROR_DIFF_2 "Accuracy of eigen vectors computing less than required." 61#define MESSAGE_ERROR_ORTHO "Matrix of eigen vectors is not orthogonal." 62#define MESSAGE_ERROR_ORDER "Eigen values are not sorted in ascending order." 63 64const int COUNT_NORM_TYPES = 3; 65const int NORM_TYPE[COUNT_NORM_TYPES] = {cv::NORM_L1, cv::NORM_L2, cv::NORM_INF}; 66 67enum TASK_TYPE_EIGEN {VALUES, VECTORS}; 68 69class Core_EigenTest: public cvtest::BaseTest 70{ 71public: 72 73 Core_EigenTest(); 74 ~Core_EigenTest(); 75 76protected: 77 78 bool test_values(const cv::Mat& src); // complex test for eigen without vectors 79 bool check_full(int type); // compex test for symmetric matrix 80 virtual void run (int) = 0; // main testing method 81 82protected: 83 84 float eps_val_32, eps_vec_32; 85 float eps_val_64, eps_vec_64; 86 int ntests; 87 88 bool check_pair_count(const cv::Mat& src, const cv::Mat& evalues, int low_index = -1, int high_index = -1); 89 bool check_pair_count(const cv::Mat& src, const cv::Mat& evalues, const cv::Mat& evectors, int low_index = -1, int high_index = -1); 90 bool check_pairs_order(const cv::Mat& eigen_values); // checking order of eigen values & vectors (it should be none up) 91 bool check_orthogonality(const cv::Mat& U); // checking is matrix of eigen vectors orthogonal 92 bool test_pairs(const cv::Mat& src); // complex test for eigen with vectors 93 94 void print_information(const size_t norm_idx, const cv::Mat& src, double diff, double max_diff); 95}; 96 97class Core_EigenTest_Scalar : public Core_EigenTest 98{ 99public: 100 Core_EigenTest_Scalar() : Core_EigenTest() {} 101 ~Core_EigenTest_Scalar(); 102 103 virtual void run(int) = 0; 104}; 105 106class Core_EigenTest_Scalar_32 : public Core_EigenTest_Scalar 107{ 108public: 109 Core_EigenTest_Scalar_32() : Core_EigenTest_Scalar() {} 110 ~Core_EigenTest_Scalar_32(); 111 112 void run(int); 113}; 114 115class Core_EigenTest_Scalar_64 : public Core_EigenTest_Scalar 116{ 117public: 118 Core_EigenTest_Scalar_64() : Core_EigenTest_Scalar() {} 119 ~Core_EigenTest_Scalar_64(); 120 void run(int); 121}; 122 123class Core_EigenTest_32 : public Core_EigenTest 124{ 125public: 126 Core_EigenTest_32(): Core_EigenTest() {} 127 ~Core_EigenTest_32() {} 128 void run(int); 129}; 130 131class Core_EigenTest_64 : public Core_EigenTest 132{ 133public: 134 Core_EigenTest_64(): Core_EigenTest() {} 135 ~Core_EigenTest_64() {} 136 void run(int); 137}; 138 139Core_EigenTest_Scalar::~Core_EigenTest_Scalar() {} 140Core_EigenTest_Scalar_32::~Core_EigenTest_Scalar_32() {} 141Core_EigenTest_Scalar_64::~Core_EigenTest_Scalar_64() {} 142 143void Core_EigenTest_Scalar_32::run(int) 144{ 145 for (int i = 0; i < ntests; ++i) 146 { 147 float value = cv::randu<float>(); 148 cv::Mat src(1, 1, CV_32FC1, Scalar::all((float)value)); 149 test_values(src); 150 } 151} 152 153void Core_EigenTest_Scalar_64::run(int) 154{ 155 for (int i = 0; i < ntests; ++i) 156 { 157 float value = cv::randu<float>(); 158 cv::Mat src(1, 1, CV_64FC1, Scalar::all((double)value)); 159 test_values(src); 160 } 161} 162 163void Core_EigenTest_32::run(int) { check_full(CV_32FC1); } 164void Core_EigenTest_64::run(int) { check_full(CV_64FC1); } 165 166Core_EigenTest::Core_EigenTest() 167: eps_val_32(1e-3f), eps_vec_32(12e-3f), 168 eps_val_64(1e-4f), eps_vec_64(1e-3f), ntests(100) {} 169Core_EigenTest::~Core_EigenTest() {} 170 171bool Core_EigenTest::check_pair_count(const cv::Mat& src, const cv::Mat& evalues, int low_index, int high_index) 172{ 173 int n = src.rows, s = sign(high_index); 174 if (!( (evalues.rows == n - max<int>(0, low_index) - ((int)((n/2.0)*(s*s-s)) + (1+s-s*s)*(n - (high_index+1)))) && (evalues.cols == 1))) 175 { 176 std::cout << endl; std::cout << "Checking sizes of eigen values matrix " << evalues << "..." << endl; 177 std::cout << "Number of rows: " << evalues.rows << " Number of cols: " << evalues.cols << endl; 178 std::cout << "Size of src symmetric matrix: " << src.rows << " * " << src.cols << endl; std::cout << endl; 179 CV_Error(CORE_EIGEN_ERROR_COUNT, MESSAGE_ERROR_COUNT); 180 return false; 181 } 182 return true; 183} 184 185bool Core_EigenTest::check_pair_count(const cv::Mat& src, const cv::Mat& evalues, const cv::Mat& evectors, int low_index, int high_index) 186{ 187 int n = src.rows, s = sign(high_index); 188 int right_eigen_pair_count = n - max<int>(0, low_index) - ((int)((n/2.0)*(s*s-s)) + (1+s-s*s)*(n - (high_index+1))); 189 190 if (!(evectors.rows == right_eigen_pair_count && evectors.cols == right_eigen_pair_count)) 191 { 192 std::cout << endl; std::cout << "Checking sizes of eigen vectors matrix " << evectors << "..." << endl; 193 std::cout << "Number of rows: " << evectors.rows << " Number of cols: " << evectors.cols << endl; 194 std:: cout << "Size of src symmetric matrix: " << src.rows << " * " << src.cols << endl; std::cout << endl; 195 CV_Error (CORE_EIGEN_ERROR_SIZE, MESSAGE_ERROR_SIZE); 196 return false; 197 } 198 199 if (!(evalues.rows == right_eigen_pair_count && evalues.cols == 1)) 200 { 201 std::cout << endl; std::cout << "Checking sizes of eigen values matrix " << evalues << "..." << endl; 202 std::cout << "Number of rows: " << evalues.rows << " Number of cols: " << evalues.cols << endl; 203 std:: cout << "Size of src symmetric matrix: " << src.rows << " * " << src.cols << endl; std::cout << endl; 204 CV_Error (CORE_EIGEN_ERROR_COUNT, MESSAGE_ERROR_COUNT); 205 return false; 206 } 207 208 return true; 209} 210 211void Core_EigenTest::print_information(const size_t norm_idx, const cv::Mat& src, double diff, double max_diff) 212{ 213 switch (NORM_TYPE[norm_idx]) 214 { 215 case cv::NORM_L1: std::cout << "L1"; break; 216 case cv::NORM_L2: std::cout << "L2"; break; 217 case cv::NORM_INF: std::cout << "INF"; break; 218 default: break; 219 } 220 221 cout << "-criteria... " << endl; 222 cout << "Source size: " << src.rows << " * " << src.cols << endl; 223 cout << "Difference between original eigen vectors matrix and result: " << diff << endl; 224 cout << "Maximum allowed difference: " << max_diff << endl; cout << endl; 225} 226 227bool Core_EigenTest::check_orthogonality(const cv::Mat& U) 228{ 229 int type = U.type(); 230 double eps_vec = type == CV_32FC1 ? eps_vec_32 : eps_vec_64; 231 cv::Mat UUt; cv::mulTransposed(U, UUt, false); 232 233 cv::Mat E = Mat::eye(U.rows, U.cols, type); 234 235 for (int i = 0; i < COUNT_NORM_TYPES; ++i) 236 { 237 double diff = cvtest::norm(UUt, E, NORM_TYPE[i]); 238 if (diff > eps_vec) 239 { 240 std::cout << endl; std::cout << "Checking orthogonality of matrix " << U << ": "; 241 print_information(i, U, diff, eps_vec); 242 CV_Error(CORE_EIGEN_ERROR_ORTHO, MESSAGE_ERROR_ORTHO); 243 return false; 244 } 245 } 246 247 return true; 248} 249 250bool Core_EigenTest::check_pairs_order(const cv::Mat& eigen_values) 251{ 252 switch (eigen_values.type()) 253 { 254 case CV_32FC1: 255 { 256 for (int i = 0; i < (int)(eigen_values.total() - 1); ++i) 257 if (!(eigen_values.at<float>(i, 0) > eigen_values.at<float>(i+1, 0))) 258 { 259 std::cout << endl; std::cout << "Checking order of eigen values vector " << eigen_values << "..." << endl; 260 std::cout << "Pair of indexes with non ascending of eigen values: (" << i << ", " << i+1 << ")." << endl; 261 std::cout << endl; 262 CV_Error(CORE_EIGEN_ERROR_ORDER, MESSAGE_ERROR_ORDER); 263 return false; 264 } 265 266 break; 267 } 268 269 case CV_64FC1: 270 { 271 for (int i = 0; i < (int)(eigen_values.total() - 1); ++i) 272 if (!(eigen_values.at<double>(i, 0) > eigen_values.at<double>(i+1, 0))) 273 { 274 std::cout << endl; std::cout << "Checking order of eigen values vector " << eigen_values << "..." << endl; 275 std::cout << "Pair of indexes with non ascending of eigen values: (" << i << ", " << i+1 << ")." << endl; 276 std::cout << endl; 277 CV_Error(CORE_EIGEN_ERROR_ORDER, "Eigen values are not sorted in ascending order."); 278 return false; 279 } 280 281 break; 282 } 283 284 default:; 285 } 286 287 return true; 288} 289 290bool Core_EigenTest::test_pairs(const cv::Mat& src) 291{ 292 int type = src.type(); 293 double eps_vec = type == CV_32FC1 ? eps_vec_32 : eps_vec_64; 294 295 cv::Mat eigen_values, eigen_vectors; 296 297 cv::eigen(src, eigen_values, eigen_vectors); 298 299 if (!check_pair_count(src, eigen_values, eigen_vectors)) 300 return false; 301 302 if (!check_orthogonality (eigen_vectors)) 303 return false; 304 305 if (!check_pairs_order(eigen_values)) 306 return false; 307 308 cv::Mat eigen_vectors_t; cv::transpose(eigen_vectors, eigen_vectors_t); 309 310 cv::Mat src_evec(src.rows, src.cols, type); 311 src_evec = src*eigen_vectors_t; 312 313 cv::Mat eval_evec(src.rows, src.cols, type); 314 315 switch (type) 316 { 317 case CV_32FC1: 318 { 319 for (int i = 0; i < src.cols; ++i) 320 { 321 cv::Mat tmp = eigen_values.at<float>(i, 0) * eigen_vectors_t.col(i); 322 for (int j = 0; j < src.rows; ++j) eval_evec.at<float>(j, i) = tmp.at<float>(j, 0); 323 } 324 325 break; 326 } 327 328 case CV_64FC1: 329 { 330 for (int i = 0; i < src.cols; ++i) 331 { 332 cv::Mat tmp = eigen_values.at<double>(i, 0) * eigen_vectors_t.col(i); 333 for (int j = 0; j < src.rows; ++j) eval_evec.at<double>(j, i) = tmp.at<double>(j, 0); 334 } 335 336 break; 337 } 338 339 default:; 340 } 341 342 cv::Mat disparity = src_evec - eval_evec; 343 344 for (int i = 0; i < COUNT_NORM_TYPES; ++i) 345 { 346 double diff = cvtest::norm(disparity, NORM_TYPE[i]); 347 if (diff > eps_vec) 348 { 349 std::cout << endl; std::cout << "Checking accuracy of eigen vectors computing for matrix " << src << ": "; 350 print_information(i, src, diff, eps_vec); 351 CV_Error(CORE_EIGEN_ERROR_DIFF, MESSAGE_ERROR_DIFF_2); 352 return false; 353 } 354 } 355 356 return true; 357} 358 359bool Core_EigenTest::test_values(const cv::Mat& src) 360{ 361 int type = src.type(); 362 double eps_val = type == CV_32FC1 ? eps_val_32 : eps_val_64; 363 364 cv::Mat eigen_values_1, eigen_values_2, eigen_vectors; 365 366 if (!test_pairs(src)) return false; 367 368 cv::eigen(src, eigen_values_1, eigen_vectors); 369 cv::eigen(src, eigen_values_2); 370 371 if (!check_pair_count(src, eigen_values_2)) return false; 372 373 for (int i = 0; i < COUNT_NORM_TYPES; ++i) 374 { 375 double diff = cvtest::norm(eigen_values_1, eigen_values_2, NORM_TYPE[i]); 376 if (diff > eps_val) 377 { 378 std::cout << endl; std::cout << "Checking accuracy of eigen values computing for matrix " << src << ": "; 379 print_information(i, src, diff, eps_val); 380 CV_Error(CORE_EIGEN_ERROR_DIFF, MESSAGE_ERROR_DIFF_1); 381 return false; 382 } 383 } 384 385 return true; 386} 387 388bool Core_EigenTest::check_full(int type) 389{ 390 const int MAX_DEGREE = 7; 391 392 srand((unsigned int)time(0)); 393 394 for (int i = 0; i < ntests; ++i) 395 { 396 int src_size = (int)(std::pow(2.0, (rand()%MAX_DEGREE)+1.)); 397 398 cv::Mat src(src_size, src_size, type); 399 400 for (int j = 0; j < src.rows; ++j) 401 for (int k = j; k < src.cols; ++k) 402 if (type == CV_32FC1) src.at<float>(k, j) = src.at<float>(j, k) = cv::randu<float>(); 403 else src.at<double>(k, j) = src.at<double>(j, k) = cv::randu<double>(); 404 405 if (!test_values(src)) return false; 406 } 407 408 return true; 409} 410 411TEST(Core_Eigen, scalar_32) {Core_EigenTest_Scalar_32 test; test.safe_run(); } 412TEST(Core_Eigen, scalar_64) {Core_EigenTest_Scalar_64 test; test.safe_run(); } 413TEST(Core_Eigen, vector_32) { Core_EigenTest_32 test; test.safe_run(); } 414TEST(Core_Eigen, vector_64) { Core_EigenTest_64 test; test.safe_run(); } 415