opencv/src/cvaux/cvoneway.cpp

   1 /*
   2  *  cvoneway.cpp
   3  *  one_way_sample
   4  *
   5  *  Created by Victor  Eruhimov on 3/23/10.
   6  *  Copyright 2010 Argus Corp. All rights reserved.
   7  *
   8  */
   9
  10 #include <highgui.h>
  11 #include <cxcore.h>
  12 #include <cv.h>
  13 #include <cvaux.h>
  14
  15 namespace cv{
  16
  17     inline int round(float value)
  18     {
  19         if(value > 0)
  20         {
  21             return int(value + 0.5f);
  22         }
  23         else
  24         {
  25             return int(value - 0.5f);
  26         }
  27     }
  28
  29     inline CvRect resize_rect(CvRect rect, float alpha)
  30     {
  31         return cvRect(rect.x + round((float)(0.5*(1 - alpha)*rect.width)), rect.y + round((float)(0.5*(1 - alpha)*rect.height)),
  32                       round(rect.width*alpha), round(rect.height*alpha));
  33     }
  34
  35     CvMat* ConvertImageToMatrix(IplImage* patch);
  36
  37     class CvCameraPose
  38         {
  39         public:
  40             CvCameraPose()
  41             {
  42                 m_rotation = cvCreateMat(1, 3, CV_32FC1);
  43                 m_translation = cvCreateMat(1, 3, CV_32FC1);
  44             };
  45
  46             ~CvCameraPose()
  47             {
  48                 cvReleaseMat(&m_rotation);
  49                 cvReleaseMat(&m_translation);
  50             };
  51
  52             void SetPose(CvMat* rotation, CvMat* translation)
  53             {
  54                 cvCopy(rotation, m_rotation);
  55                 cvCopy(translation, m_translation);
  56             };
  57
  58             CvMat* GetRotation() {return m_rotation;};
  59             CvMat* GetTranslation() {return m_translation;};
  60
  61         protected:
  62             CvMat* m_rotation;
  63             CvMat* m_translation;
  64         };
  65
  66     // AffineTransformPatch: generates an affine transformed image patch.
  67     // - src: source image (roi is supported)
  68     // - dst: output image. ROI of dst image should be 2 times smaller than ROI of src.
  69     // - pose: parameters of an affine transformation
  70     void AffineTransformPatch(IplImage* src, IplImage* dst, CvAffinePose pose);
  71
  72     // GenerateAffineTransformFromPose: generates an affine transformation matrix from CvAffinePose instance
  73     // - size: the size of image patch
  74     // - pose: affine transformation
  75     // - transform: 2x3 transformation matrix
  76     void GenerateAffineTransformFromPose(CvSize size, CvAffinePose pose, CvMat* transform);
  77
  78     // Generates a random affine pose
  79     CvAffinePose GenRandomAffinePose();
  80
  81
  82     const static int num_mean_components = 500;
  83     const static float noise_intensity = 0.15f;
  84
  85
  86     static inline CvPoint rect_center(CvRect rect)
  87     {
  88         return cvPoint(rect.x + rect.width/2, rect.y + rect.height/2);
  89     }
  90
  91     void homography_transform(IplImage* frontal, IplImage* result, CvMat* homography)
  92     {
  93         cvWarpPerspective(frontal, result, homography);
  94     }
  95
  96     CvAffinePose perturbate_pose(CvAffinePose pose, float noise)
  97     {
  98         // perturbate the matrix
  99         float noise_mult_factor = 1 + (0.5f - float(rand())/RAND_MAX)*noise;
 100         float noise_add_factor = noise_mult_factor - 1;
 101
 102         CvAffinePose pose_pert = pose;
 103         pose_pert.phi += noise_add_factor;
 104         pose_pert.theta += noise_mult_factor;
 105         pose_pert.lambda1 *= noise_mult_factor;
 106         pose_pert.lambda2 *= noise_mult_factor;
 107
 108         return pose_pert;
 109     }
 110
 111     void generate_mean_patch(IplImage* frontal, IplImage* result, CvAffinePose pose, int pose_count, float noise)
 112     {
 113         IplImage* sum = cvCreateImage(cvSize(result->width, result->height), IPL_DEPTH_32F, 1);
 114         IplImage* workspace = cvCloneImage(result);
 115         IplImage* workspace_float = cvCloneImage(sum);
 116
 117         cvSetZero(sum);
 118         for(int i = 0; i < pose_count; i++)
 119         {
 120             CvAffinePose pose_pert = perturbate_pose(pose, noise);
 121
 122             AffineTransformPatch(frontal, workspace, pose_pert);
 123             cvConvertScale(workspace, workspace_float);
 124             cvAdd(sum, workspace_float, sum);
 125         }
 126
 127         cvConvertScale(sum, result, 1.0f/pose_count);
 128
 129         cvReleaseImage(&workspace);
 130         cvReleaseImage(&sum);
 131         cvReleaseImage(&workspace_float);
 132     }
 133
 134     void generate_mean_patch_fast(IplImage* frontal, IplImage* result, CvAffinePose pose,
 135                                   CvMat* pca_hr_avg, CvMat* pca_hr_eigenvectors, const OneWayDescriptor* pca_descriptors)
 136     {
 137         for(int i = 0; i < pca_hr_eigenvectors->cols; i++)
 138         {
 139
 140         }
 141     }
 142
 143     void readPCAFeatures(const char* filename, CvMat** avg, CvMat** eigenvectors);
 144     void eigenvector2image(CvMat* eigenvector, IplImage* img);
 145
 146     void FindOneWayDescriptor(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch, int& desc_idx, int& pose_idx, float& distance,
 147                               CvMat* avg = 0, CvMat* eigenvalues = 0);
 148
 149     void FindOneWayDescriptor(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch, int n,
 150                               std::vector<int>& desc_idxs, std::vector<int>&  pose_idxs, std::vector<float>& distances,
 151                               CvMat* avg = 0, CvMat* eigenvalues = 0);
 152
 153     void FindOneWayDescriptor(cv::flann::Index* m_pca_descriptors_tree, CvSize patch_size, int m_pca_dim_low, int m_pose_count, IplImage* patch, int& desc_idx, int& pose_idx, float& distance,
 154                               CvMat* avg = 0, CvMat* eigenvalues = 0);
 155
 156     void FindOneWayDescriptorEx(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch,
 157                                 float scale_min, float scale_max, float scale_step,
 158                                 int& desc_idx, int& pose_idx, float& distance, float& scale,
 159                                 CvMat* avg, CvMat* eigenvectors);
 160
 161     void FindOneWayDescriptorEx(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch,
 162                                 float scale_min, float scale_max, float scale_step,
 163                                 int n, std::vector<int>& desc_idxs, std::vector<int>& pose_idxs,
 164                                 std::vector<float>& distances, std::vector<float>& scales,
 165                                 CvMat* avg, CvMat* eigenvectors);
 166
 167     void FindOneWayDescriptorEx(cv::flann::Index* m_pca_descriptors_tree, CvSize patch_size, int m_pca_dim_low, int m_pose_count, IplImage* patch,
 168                                 float scale_min, float scale_max, float scale_step,
 169                                 int& desc_idx, int& pose_idx, float& distance, float& scale,
 170                                 CvMat* avg, CvMat* eigenvectors);
 171
 172     inline CvRect fit_rect_roi_fixedsize(CvRect rect, CvRect roi)
 173     {
 174         CvRect fit = rect;
 175         fit.x = MAX(fit.x, roi.x);
 176         fit.y = MAX(fit.y, roi.y);
 177         fit.x = MIN(fit.x, roi.x + roi.width - fit.width - 1);
 178         fit.y = MIN(fit.y, roi.y + roi.height - fit.height - 1);
 179         return(fit);
 180     }
 181
 182     inline CvRect fit_rect_fixedsize(CvRect rect, IplImage* img)
 183     {
 184         CvRect roi = cvGetImageROI(img);
 185         return fit_rect_roi_fixedsize(rect, roi);
 186     }
 187
 188     OneWayDescriptor::OneWayDescriptor()
 189     {
 190         m_pose_count = 0;
 191         m_samples = 0;
 192         m_input_patch = 0;
 193         m_train_patch = 0;
 194         m_pca_coeffs = 0;
 195         m_affine_poses = 0;
 196         m_transforms = 0;
 197         m_pca_dim_low = 100;
 198         m_pca_dim_high = 100;
 199     }
 200
 201     OneWayDescriptor::~OneWayDescriptor()
 202     {
 203         if(m_pose_count)
 204         {
 205             for(int i = 0; i < m_pose_count; i++)
 206             {
 207                 cvReleaseImage(&m_samples[i]);
 208                 cvReleaseMat(&m_pca_coeffs[i]);
 209             }
 210             cvReleaseImage(&m_input_patch);
 211             cvReleaseImage(&m_train_patch);
 212             delete []m_samples;
 213             delete []m_pca_coeffs;
 214
 215             if(!m_transforms)
 216             {
 217                 delete []m_affine_poses;
 218             }
 219         }
 220     }
 221
 222     void OneWayDescriptor::Allocate(int pose_count, CvSize size, int nChannels)
 223     {
 224         m_pose_count = pose_count;
 225         m_samples = new IplImage* [m_pose_count];
 226         m_pca_coeffs = new CvMat* [m_pose_count];
 227         m_patch_size = cvSize(size.width/2, size.height/2);
 228
 229         if(!m_transforms)
 230         {
 231             m_affine_poses = new CvAffinePose[m_pose_count];
 232         }
 233
 234         int length = m_pca_dim_low;//roi.width*roi.height;
 235         for(int i = 0; i < m_pose_count; i++)
 236         {
 237             m_samples[i] = cvCreateImage(cvSize(size.width/2, size.height/2), IPL_DEPTH_32F, nChannels);
 238             m_pca_coeffs[i] = cvCreateMat(1, length, CV_32FC1);
 239         }
 240
 241         m_input_patch = cvCreateImage(GetPatchSize(), IPL_DEPTH_8U, 1);
 242         m_train_patch = cvCreateImage(GetInputPatchSize(), IPL_DEPTH_8U, 1);
 243     }
 244
 245     void cvmSet2DPoint(CvMat* matrix, int row, int col, CvPoint2D32f point)
 246     {
 247         cvmSet(matrix, row, col, point.x);
 248         cvmSet(matrix, row, col + 1, point.y);
 249     }
 250
 251     void cvmSet3DPoint(CvMat* matrix, int row, int col, CvPoint3D32f point)
 252     {
 253         cvmSet(matrix, row, col, point.x);
 254         cvmSet(matrix, row, col + 1, point.y);
 255         cvmSet(matrix, row, col + 2, point.z);
 256     }
 257
 258     CvAffinePose GenRandomAffinePose()
 259     {
 260         const float scale_min = 0.8f;
 261         const float scale_max = 1.2f;
 262         CvAffinePose pose;
 263         pose.theta = float(rand())/RAND_MAX*120 - 60;
 264         pose.phi = float(rand())/RAND_MAX*360;
 265         pose.lambda1 = scale_min + float(rand())/RAND_MAX*(scale_max - scale_min);
 266         pose.lambda2 = scale_min + float(rand())/RAND_MAX*(scale_max - scale_min);
 267
 268         return pose;
 269     }
 270
 271     void GenerateAffineTransformFromPose(CvSize size, CvAffinePose pose, CvMat* transform)
 272     {
 273         CvMat* temp = cvCreateMat(3, 3, CV_32FC1);
 274         CvMat* final = cvCreateMat(3, 3, CV_32FC1);
 275         cvmSet(temp, 2, 0, 0.0f);
 276         cvmSet(temp, 2, 1, 0.0f);
 277         cvmSet(temp, 2, 2, 1.0f);
 278
 279         CvMat rotation;
 280         cvGetSubRect(temp, &rotation, cvRect(0, 0, 3, 2));
 281
 282         cv2DRotationMatrix(cvPoint2D32f(size.width/2, size.height/2), pose.phi, 1.0, &rotation);
 283         cvCopy(temp, final);
 284
 285         cvmSet(temp, 0, 0, pose.lambda1);
 286         cvmSet(temp, 0, 1, 0.0f);
 287         cvmSet(temp, 1, 0, 0.0f);
 288         cvmSet(temp, 1, 1, pose.lambda2);
 289         cvmSet(temp, 0, 2, size.width/2*(1 - pose.lambda1));
 290         cvmSet(temp, 1, 2, size.height/2*(1 - pose.lambda2));
 291         cvMatMul(temp, final, final);
 292
 293         cv2DRotationMatrix(cvPoint2D32f(size.width/2, size.height/2), pose.theta - pose.phi, 1.0, &rotation);
 294         cvMatMul(temp, final, final);
 295
 296         cvGetSubRect(final, &rotation, cvRect(0, 0, 3, 2));
 297         cvCopy(&rotation, transform);
 298
 299         cvReleaseMat(&temp);
 300         cvReleaseMat(&final);
 301     }
 302
 303     void AffineTransformPatch(IplImage* src, IplImage* dst, CvAffinePose pose)
 304     {
 305         CvRect src_large_roi = cvGetImageROI(src);
 306
 307         IplImage* temp = cvCreateImage(cvSize(src_large_roi.width, src_large_roi.height), IPL_DEPTH_32F, src->nChannels);
 308         cvSetZero(temp);
 309         IplImage* temp2 = cvCloneImage(temp);
 310         CvMat* rotation_phi = cvCreateMat(2, 3, CV_32FC1);
 311
 312         CvSize new_size = cvSize(temp->width*pose.lambda1, temp->height*pose.lambda2);
 313         IplImage* temp3 = cvCreateImage(new_size, IPL_DEPTH_32F, src->nChannels);
 314
 315         cvConvertScale(src, temp);
 316         cvResetImageROI(temp);
 317
 318
 319         cv2DRotationMatrix(cvPoint2D32f(temp->width/2, temp->height/2), pose.phi, 1.0, rotation_phi);
 320         cvWarpAffine(temp, temp2, rotation_phi);
 321
 322         cvSetZero(temp);
 323
 324         cvResize(temp2, temp3);
 325
 326         cv2DRotationMatrix(cvPoint2D32f(temp3->width/2, temp3->height/2), pose.theta - pose.phi, 1.0, rotation_phi);
 327         cvWarpAffine(temp3, temp, rotation_phi);
 328
 329         cvSetImageROI(temp, cvRect(temp->width/2 - src_large_roi.width/4, temp->height/2 - src_large_roi.height/4,
 330                                    src_large_roi.width/2, src_large_roi.height/2));
 331         cvConvertScale(temp, dst);
 332         cvReleaseMat(&rotation_phi);
 333
 334         cvReleaseImage(&temp3);
 335         cvReleaseImage(&temp2);
 336         cvReleaseImage(&temp);
 337     }
 338
 339     void OneWayDescriptor::GenerateSamples(int pose_count, IplImage* frontal, int norm)
 340     {
 341         /*    if(m_transforms)
 342          {
 343          GenerateSamplesWithTransforms(pose_count, frontal);
 344          return;
 345          }
 346          */
 347         CvRect roi = cvGetImageROI(frontal);
 348         IplImage* patch_8u = cvCreateImage(cvSize(roi.width/2, roi.height/2), frontal->depth, frontal->nChannels);
 349         for(int i = 0; i < pose_count; i++)
 350         {
 351             if(!m_transforms)
 352             {
 353                 m_affine_poses[i] = GenRandomAffinePose();
 354             }
 355             //AffineTransformPatch(frontal, patch_8u, m_affine_poses[i]);
 356             generate_mean_patch(frontal, patch_8u, m_affine_poses[i], num_mean_components, noise_intensity);
 357
 358             float scale = 1.0f;
 359             if(norm)
 360             {
 361                 float sum = cvSum(patch_8u).val[0];
 362                 scale = 1/sum;
 363             }
 364             cvConvertScale(patch_8u, m_samples[i], scale);
 365
 366 #if 0
 367             double maxval;
 368             cvMinMaxLoc(m_samples[i], 0, &maxval);
 369             IplImage* test = cvCreateImage(cvSize(roi.width/2, roi.height/2), IPL_DEPTH_8U, 1);
 370             cvConvertScale(m_samples[i], test, 255.0/maxval);
 371             cvNamedWindow("1", 1);
 372             cvShowImage("1", test);
 373             cvWaitKey(0);
 374 #endif
 375         }
 376         cvReleaseImage(&patch_8u);
 377     }
 378
 379     void OneWayDescriptor::GenerateSamplesFast(IplImage* frontal, CvMat* pca_hr_avg,
 380                                                CvMat* pca_hr_eigenvectors, OneWayDescriptor* pca_descriptors)
 381     {
 382         CvRect roi = cvGetImageROI(frontal);
 383         if(roi.width != GetInputPatchSize().width || roi.height != GetInputPatchSize().height)
 384         {
 385             cvResize(frontal, m_train_patch);
 386             frontal = m_train_patch;
 387         }
 388
 389         CvMat* pca_coeffs = cvCreateMat(1, pca_hr_eigenvectors->cols, CV_32FC1);
 390         double maxval;
 391         cvMinMaxLoc(frontal, 0, &maxval);
 392         CvMat* frontal_data = ConvertImageToMatrix(frontal);
 393
 394         float sum = cvSum(frontal_data).val[0];
 395         cvConvertScale(frontal_data, frontal_data, 1.0f/sum);
 396         cvProjectPCA(frontal_data, pca_hr_avg, pca_hr_eigenvectors, pca_coeffs);
 397         for(int i = 0; i < m_pose_count; i++)
 398         {
 399             cvSetZero(m_samples[i]);
 400             for(int j = 0; j < m_pca_dim_high; j++)
 401             {
 402                 float coeff = cvmGet(pca_coeffs, 0, j);
 403                 IplImage* patch = pca_descriptors[j + 1].GetPatch(i);
 404                 cvAddWeighted(m_samples[i], 1.0, patch, coeff, 0, m_samples[i]);
 405
 406 #if 0
 407                 printf("coeff%d = %f\n", j, coeff);
 408                 IplImage* test = cvCreateImage(cvSize(12, 12), IPL_DEPTH_8U, 1);
 409                 double maxval;
 410                 cvMinMaxLoc(patch, 0, &maxval);
 411                 cvConvertScale(patch, test, 255.0/maxval);
 412                 cvNamedWindow("1", 1);
 413                 cvShowImage("1", test);
 414                 cvWaitKey(0);
 415 #endif
 416             }
 417
 418             cvAdd(pca_descriptors[0].GetPatch(i), m_samples[i], m_samples[i]);
 419             float sum = cvSum(m_samples[i]).val[0];
 420             cvConvertScale(m_samples[i], m_samples[i], 1.0/sum);
 421
 422 #if 0
 423             IplImage* test = cvCreateImage(cvSize(12, 12), IPL_DEPTH_8U, 1);
 424             /*        IplImage* temp1 = cvCreateImage(cvSize(12, 12), IPL_DEPTH_32F, 1);
 425              eigenvector2image(pca_hr_avg, temp1);
 426              IplImage* test = cvCreateImage(cvSize(12, 12), IPL_DEPTH_8U, 1);
 427              cvAdd(m_samples[i], temp1, temp1);
 428              cvMinMaxLoc(temp1, 0, &maxval);
 429              cvConvertScale(temp1, test, 255.0/maxval);*/
 430             cvMinMaxLoc(m_samples[i], 0, &maxval);
 431             cvConvertScale(m_samples[i], test, 255.0/maxval);
 432
 433             cvNamedWindow("1", 1);
 434             cvShowImage("1", frontal);
 435             cvNamedWindow("2", 1);
 436             cvShowImage("2", test);
 437             cvWaitKey(0);
 438 #endif
 439         }
 440
 441         cvReleaseMat(&pca_coeffs);
 442         cvReleaseMat(&frontal_data);
 443     }
 444
 445     void OneWayDescriptor::SetTransforms(CvAffinePose* poses, CvMat** transforms)
 446     {
 447         if(m_affine_poses)
 448         {
 449             delete []m_affine_poses;
 450         }
 451
 452         m_affine_poses = poses;
 453         m_transforms = transforms;
 454     }
 455
 456     void OneWayDescriptor::Initialize(int pose_count, IplImage* frontal, const char* feature_name, int norm)
 457     {
 458         m_feature_name = std::string(feature_name);
 459         CvRect roi = cvGetImageROI(frontal);
 460         m_center = rect_center(roi);
 461
 462         Allocate(pose_count, cvSize(roi.width, roi.height), frontal->nChannels);
 463
 464         GenerateSamples(pose_count, frontal, norm);
 465     }
 466
 467     void OneWayDescriptor::InitializeFast(int pose_count, IplImage* frontal, const char* feature_name,
 468                                           CvMat* pca_hr_avg, CvMat* pca_hr_eigenvectors, OneWayDescriptor* pca_descriptors)
 469     {
 470         if(pca_hr_avg == 0)
 471         {
 472             Initialize(pose_count, frontal, feature_name, 1);
 473             return;
 474         }
 475         m_feature_name = std::string(feature_name);
 476         CvRect roi = cvGetImageROI(frontal);
 477         m_center = rect_center(roi);
 478
 479         Allocate(pose_count, cvSize(roi.width, roi.height), frontal->nChannels);
 480
 481         GenerateSamplesFast(frontal, pca_hr_avg, pca_hr_eigenvectors, pca_descriptors);
 482     }
 483
 484     void OneWayDescriptor::InitializePCACoeffs(CvMat* avg, CvMat* eigenvectors)
 485     {
 486         for(int i = 0; i < m_pose_count; i++)
 487         {
 488             ProjectPCASample(m_samples[i], avg, eigenvectors, m_pca_coeffs[i]);
 489         }
 490     }
 491
 492     void OneWayDescriptor::ProjectPCASample(IplImage* patch, CvMat* avg, CvMat* eigenvectors, CvMat* pca_coeffs) const
 493     {
 494         CvMat* patch_mat = ConvertImageToMatrix(patch);
 495         //    CvMat eigenvectorsr;
 496         //    cvGetSubRect(eigenvectors, &eigenvectorsr, cvRect(0, 0, eigenvectors->cols, pca_coeffs->cols));
 497         CvMat* temp = cvCreateMat(1, eigenvectors->cols, CV_32FC1);
 498         cvProjectPCA(patch_mat, avg, eigenvectors, temp);
 499         CvMat temp1;
 500         cvGetSubRect(temp, &temp1, cvRect(0, 0, pca_coeffs->cols, 1));
 501         cvCopy(&temp1, pca_coeffs);
 502
 503         cvReleaseMat(&temp);
 504         cvReleaseMat(&patch_mat);
 505     }
 506
 507     void OneWayDescriptor::EstimatePosePCA(CvArr* patch, int& pose_idx, float& distance, CvMat* avg, CvMat* eigenvectors) const
 508     {
 509         if(avg == 0)
 510         {
 511             // do not use pca
 512             if (!CV_IS_MAT(patch))
 513             {
 514                 EstimatePose((IplImage*)patch, pose_idx, distance);
 515             }
 516             else
 517             {
 518
 519             }
 520             return;
 521         }
 522         CvRect roi;
 523         if (!CV_IS_MAT(patch))
 524         {
 525             roi = cvGetImageROI((IplImage*)patch);
 526             if(roi.width != GetPatchSize().width || roi.height != GetPatchSize().height)
 527             {
 528                 cvResize(patch, m_input_patch);
 529                 patch = m_input_patch;
 530                 roi = cvGetImageROI((IplImage*)patch);
 531             }
 532         }
 533
 534         CvMat* pca_coeffs = cvCreateMat(1, m_pca_dim_low, CV_32FC1);
 535
 536         if (CV_IS_MAT(patch))
 537         {
 538             cvCopy((CvMat*)patch, pca_coeffs);
 539         }
 540         else
 541         {
 542             IplImage* patch_32f = cvCreateImage(cvSize(roi.width, roi.height), IPL_DEPTH_32F, 1);
 543             float sum = cvSum(patch).val[0];
 544             cvConvertScale(patch, patch_32f, 1.0f/sum);
 545             ProjectPCASample(patch_32f, avg, eigenvectors, pca_coeffs);
 546             cvReleaseImage(&patch_32f);
 547         }
 548
 549
 550         distance = 1e10;
 551         pose_idx = -1;
 552
 553         for(int i = 0; i < m_pose_count; i++)
 554         {
 555             float dist = cvNorm(m_pca_coeffs[i], pca_coeffs);
 556             //          float dist = 0;
 557             //          float data1, data2;
 558             //          //CvMat* pose_pca_coeffs = m_pca_coeffs[i];
 559             //          for (int x=0; x < pca_coeffs->width; x++)
 560             //                  for (int y =0 ; y < pca_coeffs->height; y++)
 561             //                  {
 562             //                          data1 = ((float*)(pca_coeffs->data.ptr + pca_coeffs->step*x))[y];
 563             //                          data2 = ((float*)(m_pca_coeffs[i]->data.ptr + m_pca_coeffs[i]->step*x))[y];
 564             //                          dist+=(data1-data2)*(data1-data2);
 565             //                  }
 566             ////#if 1
 567             //          for (int j = 0; j < m_pca_dim_low; j++)
 568             //          {
 569             //                  dist += (pose_pca_coeffs->data.fl[j]- pca_coeffs->data.fl[j])*(pose_pca_coeffs->data.fl[j]- pca_coeffs->data.fl[j]);
 570             //          }
 571             //#else
 572             //          for (int j = 0; j <= m_pca_dim_low - 4; j += 4)
 573             //          {
 574             //                  dist += (pose_pca_coeffs->data.fl[j]- pca_coeffs->data.fl[j])*
 575             //                          (pose_pca_coeffs->data.fl[j]- pca_coeffs->data.fl[j]);
 576             //                  dist += (pose_pca_coeffs->data.fl[j+1]- pca_coeffs->data.fl[j+1])*
 577             //                          (pose_pca_coeffs->data.fl[j+1]- pca_coeffs->data.fl[j+1]);
 578             //                  dist += (pose_pca_coeffs->data.fl[j+2]- pca_coeffs->data.fl[j+2])*
 579             //                          (pose_pca_coeffs->data.fl[j+2]- pca_coeffs->data.fl[j+2]);
 580             //                  dist += (pose_pca_coeffs->data.fl[j+3]- pca_coeffs->data.fl[j+3])*
 581             //                          (pose_pca_coeffs->data.fl[j+3]- pca_coeffs->data.fl[j+3]);
 582             //          }
 583             //#endif
 584             if(dist < distance)
 585             {
 586                 distance = dist;
 587                 pose_idx = i;
 588             }
 589         }
 590
 591         cvReleaseMat(&pca_coeffs);
 592     }
 593
 594     void OneWayDescriptor::EstimatePose(IplImage* patch, int& pose_idx, float& distance) const
 595     {
 596         distance = 1e10;
 597         pose_idx = -1;
 598
 599         CvRect roi = cvGetImageROI(patch);
 600         IplImage* patch_32f = cvCreateImage(cvSize(roi.width, roi.height), IPL_DEPTH_32F, patch->nChannels);
 601         float sum = cvSum(patch).val[0];
 602         cvConvertScale(patch, patch_32f, 1/sum);
 603
 604         for(int i = 0; i < m_pose_count; i++)
 605         {
 606             if(m_samples[i]->width != patch_32f->width || m_samples[i]->height != patch_32f->height)
 607             {
 608                 continue;
 609             }
 610             float dist = cvNorm(m_samples[i], patch_32f);
 611             //float dist = 0.0f;
 612             //float i1,i2;
 613
 614             //for (int y = 0; y<patch_32f->height; y++)
 615             //  for (int x = 0; x< patch_32f->width; x++)
 616             //  {
 617             //          i1 = ((float*)(m_samples[i]->imageData + m_samples[i]->widthStep*y))[x];
 618             //          i2 = ((float*)(patch_32f->imageData + patch_32f->widthStep*y))[x];
 619             //          dist+= (i1-i2)*(i1-i2);
 620             //  }
 621
 622             if(dist < distance)
 623             {
 624                 distance = dist;
 625                 pose_idx = i;
 626             }
 627
 628 #if 0
 629             IplImage* img1 = cvCreateImage(cvSize(roi.width, roi.height), IPL_DEPTH_8U, 1);
 630             IplImage* img2 = cvCreateImage(cvSize(roi.width, roi.height), IPL_DEPTH_8U, 1);
 631             double maxval;
 632             cvMinMaxLoc(m_samples[i], 0, &maxval);
 633             cvConvertScale(m_samples[i], img1, 255.0/maxval);
 634             cvMinMaxLoc(patch_32f, 0, &maxval);
 635             cvConvertScale(patch_32f, img2, 255.0/maxval);
 636
 637             cvNamedWindow("1", 1);
 638             cvShowImage("1", img1);
 639             cvNamedWindow("2", 1);
 640             cvShowImage("2", img2);
 641             printf("Distance = %f\n", dist);
 642             cvWaitKey(0);
 643 #endif
 644         }
 645
 646         cvReleaseImage(&patch_32f);
 647     }
 648
 649     void OneWayDescriptor::Save(const char* path)
 650     {
 651         for(int i = 0; i < m_pose_count; i++)
 652         {
 653             char buf[1024];
 654             sprintf(buf, "%s/patch_%04d.jpg", path, i);
 655             IplImage* patch = cvCreateImage(cvSize(m_samples[i]->width, m_samples[i]->height), IPL_DEPTH_8U, m_samples[i]->nChannels);
 656
 657             double maxval;
 658             cvMinMaxLoc(m_samples[i], 0, &maxval);
 659             cvConvertScale(m_samples[i], patch, 255/maxval);
 660
 661             cvSaveImage(buf, patch);
 662
 663             cvReleaseImage(&patch);
 664         }
 665     }
 666
 667     void OneWayDescriptor::Write(CvFileStorage* fs, const char* name)
 668     {
 669         CvMat* mat = cvCreateMat(m_pose_count, m_samples[0]->width*m_samples[0]->height, CV_32FC1);
 670
 671         // prepare data to write as a single matrix
 672         for(int i = 0; i < m_pose_count; i++)
 673         {
 674             for(int y = 0; y < m_samples[i]->height; y++)
 675             {
 676                 for(int x = 0; x < m_samples[i]->width; x++)
 677                 {
 678                     float val = *((float*)(m_samples[i]->imageData + m_samples[i]->widthStep*y) + x);
 679                     cvmSet(mat, i, y*m_samples[i]->width + x, val);
 680                 }
 681             }
 682         }
 683
 684         cvWrite(fs, name, mat);
 685
 686         cvReleaseMat(&mat);
 687     }
 688
 689     int OneWayDescriptor::ReadByName(CvFileStorage* fs, CvFileNode* parent, const char* name)
 690     {
 691         CvMat* mat = (CvMat*)cvReadByName(fs, parent, name);
 692         if(!mat)
 693         {
 694             return 0;
 695         }
 696
 697
 698         for(int i = 0; i < m_pose_count; i++)
 699         {
 700             for(int y = 0; y < m_samples[i]->height; y++)
 701             {
 702                 for(int x = 0; x < m_samples[i]->width; x++)
 703                 {
 704                     float val = cvmGet(mat, i, y*m_samples[i]->width + x);
 705                     *((float*)(m_samples[i]->imageData + y*m_samples[i]->widthStep) + x) = val;
 706                 }
 707             }
 708         }
 709
 710         cvReleaseMat(&mat);
 711         return 1;
 712     }
 713
 714     IplImage* OneWayDescriptor::GetPatch(int index)
 715     {
 716         return m_samples[index];
 717     }
 718
 719     CvAffinePose OneWayDescriptor::GetPose(int index) const
 720     {
 721         return m_affine_poses[index];
 722     }
 723
 724     void FindOneWayDescriptor(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch, int& desc_idx, int& pose_idx, float& distance,
 725                               CvMat* avg, CvMat* eigenvectors)
 726     {
 727         desc_idx = -1;
 728         pose_idx = -1;
 729         distance = 1e10;
 730         //--------
 731         //PCA_coeffs precalculating
 732         int m_pca_dim_low = descriptors[0].GetPCADimLow();
 733         CvMat* pca_coeffs = cvCreateMat(1, m_pca_dim_low, CV_32FC1);
 734         int patch_width = descriptors[0].GetPatchSize().width;
 735         int patch_height = descriptors[0].GetPatchSize().height;
 736         if (avg)
 737         {
 738             CvRect _roi = cvGetImageROI((IplImage*)patch);
 739             IplImage* test_img = cvCreateImage(cvSize(patch_width,patch_height), IPL_DEPTH_8U, 1);
 740             if(_roi.width != patch_width|| _roi.height != patch_height)
 741             {
 742
 743                 cvResize(patch, test_img);
 744                 _roi = cvGetImageROI(test_img);
 745             }
 746             else
 747             {
 748                 cvCopy(patch,test_img);
 749             }
 750             IplImage* patch_32f = cvCreateImage(cvSize(_roi.width, _roi.height), IPL_DEPTH_32F, 1);
 751             float sum = cvSum(test_img).val[0];
 752             cvConvertScale(test_img, patch_32f, 1.0f/sum);
 753
 754             //ProjectPCASample(patch_32f, avg, eigenvectors, pca_coeffs);
 755             //Projecting PCA
 756             CvMat* patch_mat = ConvertImageToMatrix(patch_32f);
 757             CvMat* temp = cvCreateMat(1, eigenvectors->cols, CV_32FC1);
 758             cvProjectPCA(patch_mat, avg, eigenvectors, temp);
 759             CvMat temp1;
 760             cvGetSubRect(temp, &temp1, cvRect(0, 0, pca_coeffs->cols, 1));
 761             cvCopy(&temp1, pca_coeffs);
 762             cvReleaseMat(&temp);
 763             cvReleaseMat(&patch_mat);
 764             //End of projecting
 765
 766             cvReleaseImage(&patch_32f);
 767             cvReleaseImage(&test_img);
 768         }
 769
 770         //--------
 771
 772
 773
 774         for(int i = 0; i < desc_count; i++)
 775         {
 776             int _pose_idx = -1;
 777             float _distance = 0;
 778
 779 #if 0
 780             descriptors[i].EstimatePose(patch, _pose_idx, _distance);
 781 #else
 782             if (!avg)
 783             {
 784                 descriptors[i].EstimatePosePCA(patch, _pose_idx, _distance, avg, eigenvectors);
 785             }
 786             else
 787             {
 788                 descriptors[i].EstimatePosePCA(pca_coeffs, _pose_idx, _distance, avg, eigenvectors);
 789             }
 790 #endif
 791
 792             if(_distance < distance)
 793             {
 794                 desc_idx = i;
 795                 pose_idx = _pose_idx;
 796                 distance = _distance;
 797             }
 798         }
 799         cvReleaseMat(&pca_coeffs);
 800     }
 801
 802 #if defined(_KDTREE)
 803
 804     void FindOneWayDescriptor(cv::flann::Index* m_pca_descriptors_tree, CvSize patch_size, int m_pca_dim_low, int m_pose_count, IplImage* patch, int& desc_idx, int& pose_idx, float& distance,
 805                               CvMat* avg, CvMat* eigenvectors)
 806     {
 807         desc_idx = -1;
 808         pose_idx = -1;
 809         distance = 1e10;
 810         //--------
 811         //PCA_coeffs precalculating
 812         CvMat* pca_coeffs = cvCreateMat(1, m_pca_dim_low, CV_32FC1);
 813         int patch_width = patch_size.width;
 814         int patch_height = patch_size.height;
 815         //if (avg)
 816         //{
 817                 CvRect _roi = cvGetImageROI((IplImage*)patch);
 818                 IplImage* test_img = cvCreateImage(cvSize(patch_width,patch_height), IPL_DEPTH_8U, 1);
 819                 if(_roi.width != patch_width|| _roi.height != patch_height)
 820                 {
 821
 822                         cvResize(patch, test_img);
 823                         _roi = cvGetImageROI(test_img);
 824                 }
 825                 else
 826                 {
 827                         cvCopy(patch,test_img);
 828                 }
 829                 IplImage* patch_32f = cvCreateImage(cvSize(_roi.width, _roi.height), IPL_DEPTH_32F, 1);
 830                 float sum = cvSum(test_img).val[0];
 831                 cvConvertScale(test_img, patch_32f, 1.0f/sum);
 832
 833                 //ProjectPCASample(patch_32f, avg, eigenvectors, pca_coeffs);
 834                 //Projecting PCA
 835                 CvMat* patch_mat = ConvertImageToMatrix(patch_32f);
 836                 CvMat* temp = cvCreateMat(1, eigenvectors->cols, CV_32FC1);
 837                 cvProjectPCA(patch_mat, avg, eigenvectors, temp);
 838                 CvMat temp1;
 839                 cvGetSubRect(temp, &temp1, cvRect(0, 0, pca_coeffs->cols, 1));
 840                 cvCopy(&temp1, pca_coeffs);
 841                 cvReleaseMat(&temp);
 842                 cvReleaseMat(&patch_mat);
 843                 //End of projecting
 844
 845                 cvReleaseImage(&patch_32f);
 846                 cvReleaseImage(&test_img);
 847         //      }
 848
 849         //--------
 850
 851                 //float* target = new float[m_pca_dim_low];
 852                 //::flann::KNNResultSet res(1,pca_coeffs->data.fl,m_pca_dim_low);
 853                 //::flann::SearchParams params;
 854                 //params.checks = -1;
 855
 856                 //int maxDepth = 1000000;
 857                 //int neighbors_count = 1;
 858                 //int* neighborsIdx = new int[neighbors_count];
 859                 //float* distances = new float[neighbors_count];
 860                 //if (m_pca_descriptors_tree->findNearest(pca_coeffs->data.fl,neighbors_count,maxDepth,neighborsIdx,0,distances) > 0)
 861                 //{
 862                 //      desc_idx = neighborsIdx[0] / m_pose_count;
 863                 //      pose_idx = neighborsIdx[0] % m_pose_count;
 864                 //      distance = distances[0];
 865                 //}
 866                 //delete[] neighborsIdx;
 867                 //delete[] distances;
 868
 869                 cv::Mat m_object(1, m_pca_dim_low, CV_32F);
 870                 cv::Mat m_indices(1, 1, CV_32S);
 871                 cv::Mat m_dists(1, 1, CV_32F);
 872
 873                 float* object_ptr = m_object.ptr<float>(0);
 874                 for (int i=0;i<m_pca_dim_low;i++)
 875                 {
 876                         object_ptr[i] = pca_coeffs->data.fl[i];
 877                 }
 878
 879                 m_pca_descriptors_tree->knnSearch(m_object, m_indices, m_dists, 1, cv::flann::SearchParams(-1) );
 880
 881                 desc_idx = ((int*)(m_indices.ptr<int>(0)))[0] / m_pose_count;
 882                 pose_idx = ((int*)(m_indices.ptr<int>(0)))[0] % m_pose_count;
 883                 distance = ((float*)(m_dists.ptr<float>(0)))[0];
 884
 885         //      delete[] target;
 886
 887
 888         //    for(int i = 0; i < desc_count; i++)
 889         //    {
 890         //        int _pose_idx = -1;
 891         //        float _distance = 0;
 892         //
 893         //#if 0
 894         //        descriptors[i].EstimatePose(patch, _pose_idx, _distance);
 895         //#else
 896         //              if (!avg)
 897         //              {
 898         //                      descriptors[i].EstimatePosePCA(patch, _pose_idx, _distance, avg, eigenvectors);
 899         //              }
 900         //              else
 901         //              {
 902         //                      descriptors[i].EstimatePosePCA(pca_coeffs, _pose_idx, _distance, avg, eigenvectors);
 903         //              }
 904         //#endif
 905         //
 906         //        if(_distance < distance)
 907         //        {
 908         //            desc_idx = i;
 909         //            pose_idx = _pose_idx;
 910         //            distance = _distance;
 911         //        }
 912         //    }
 913         cvReleaseMat(&pca_coeffs);
 914     }
 915 #endif
 916     //**
 917     void FindOneWayDescriptor(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch, int n,
 918                               std::vector<int>& desc_idxs, std::vector<int>&  pose_idxs, std::vector<float>& distances,
 919                               CvMat* avg, CvMat* eigenvectors)
 920     {
 921         for (int i=0;i<n;i++)
 922         {
 923             desc_idxs[i] = -1;
 924             pose_idxs[i] = -1;
 925             distances[i] = 1e10;
 926         }
 927         //--------
 928         //PCA_coeffs precalculating
 929         int m_pca_dim_low = descriptors[0].GetPCADimLow();
 930         CvMat* pca_coeffs = cvCreateMat(1, m_pca_dim_low, CV_32FC1);
 931         int patch_width = descriptors[0].GetPatchSize().width;
 932         int patch_height = descriptors[0].GetPatchSize().height;
 933         if (avg)
 934         {
 935             CvRect _roi = cvGetImageROI((IplImage*)patch);
 936             IplImage* test_img = cvCreateImage(cvSize(patch_width,patch_height), IPL_DEPTH_8U, 1);
 937             if(_roi.width != patch_width|| _roi.height != patch_height)
 938             {
 939
 940                 cvResize(patch, test_img);
 941                 _roi = cvGetImageROI(test_img);
 942             }
 943             else
 944             {
 945                 cvCopy(patch,test_img);
 946             }
 947             IplImage* patch_32f = cvCreateImage(cvSize(_roi.width, _roi.height), IPL_DEPTH_32F, 1);
 948             float sum = cvSum(test_img).val[0];
 949             cvConvertScale(test_img, patch_32f, 1.0f/sum);
 950
 951             //ProjectPCASample(patch_32f, avg, eigenvectors, pca_coeffs);
 952             //Projecting PCA
 953             CvMat* patch_mat = ConvertImageToMatrix(patch_32f);
 954             CvMat* temp = cvCreateMat(1, eigenvectors->cols, CV_32FC1);
 955             cvProjectPCA(patch_mat, avg, eigenvectors, temp);
 956             CvMat temp1;
 957             cvGetSubRect(temp, &temp1, cvRect(0, 0, pca_coeffs->cols, 1));
 958             cvCopy(&temp1, pca_coeffs);
 959             cvReleaseMat(&temp);
 960             cvReleaseMat(&patch_mat);
 961             //End of projecting
 962
 963             cvReleaseImage(&patch_32f);
 964             cvReleaseImage(&test_img);
 965         }
 966         //--------
 967
 968
 969
 970         for(int i = 0; i < desc_count; i++)
 971         {
 972             int _pose_idx = -1;
 973             float _distance = 0;
 974
 975 #if 0
 976             descriptors[i].EstimatePose(patch, _pose_idx, _distance);
 977 #else
 978             if (!avg)
 979             {
 980                 descriptors[i].EstimatePosePCA(patch, _pose_idx, _distance, avg, eigenvectors);
 981             }
 982             else
 983             {
 984                 descriptors[i].EstimatePosePCA(pca_coeffs, _pose_idx, _distance, avg, eigenvectors);
 985             }
 986 #endif
 987
 988             for (int j=0;j<n;j++)
 989             {
 990                 if(_distance < distances[j])
 991                 {
 992                     for (int k=(n-1);k > j;k--)
 993                     {
 994                         desc_idxs[k] = desc_idxs[k-1];
 995                         pose_idxs[k] = pose_idxs[k-1];
 996                         distances[k] = distances[k-1];
 997                     }
 998                     desc_idxs[j] = i;
 999                     pose_idxs[j] = _pose_idx;
1000                     distances[j] = _distance;
1001                     break;
1002                 }
1003             }
1004         }
1005         cvReleaseMat(&pca_coeffs);
1006     }
1007
1008     void FindOneWayDescriptorEx(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch,
1009                                 float scale_min, float scale_max, float scale_step,
1010                                 int& desc_idx, int& pose_idx, float& distance, float& scale,
1011                                 CvMat* avg, CvMat* eigenvectors)
1012     {
1013         CvSize patch_size = descriptors[0].GetPatchSize();
1014         IplImage* input_patch;
1015         CvRect roi;
1016
1017         input_patch= cvCreateImage(patch_size, IPL_DEPTH_8U, 1);
1018         roi = cvGetImageROI((IplImage*)patch);
1019
1020         int _desc_idx, _pose_idx;
1021         float _distance;
1022         distance = 1e10;
1023         for(float cur_scale = scale_min; cur_scale < scale_max; cur_scale *= scale_step)
1024         {
1025             //        printf("Scale = %f\n", cur_scale);
1026
1027             CvRect roi_scaled = resize_rect(roi, cur_scale);
1028             cvSetImageROI(patch, roi_scaled);
1029             cvResize(patch, input_patch);
1030
1031
1032 #if 0
1033             if(roi.x > 244 && roi.y < 200)
1034             {
1035                 cvNamedWindow("1", 1);
1036                 cvShowImage("1", input_patch);
1037                 cvWaitKey(0);
1038             }
1039 #endif
1040
1041             FindOneWayDescriptor(desc_count, descriptors, input_patch, _desc_idx, _pose_idx, _distance, avg, eigenvectors);
1042             if(_distance < distance)
1043             {
1044                 distance = _distance;
1045                 desc_idx = _desc_idx;
1046                 pose_idx = _pose_idx;
1047                 scale = cur_scale;
1048             }
1049         }
1050
1051
1052         cvSetImageROI((IplImage*)patch, roi);
1053         cvReleaseImage(&input_patch);
1054
1055     }
1056
1057     void FindOneWayDescriptorEx(int desc_count, const OneWayDescriptor* descriptors, IplImage* patch,
1058                                 float scale_min, float scale_max, float scale_step,
1059                                 int n, std::vector<int>& desc_idxs, std::vector<int>& pose_idxs,
1060                                 std::vector<float>& distances, std::vector<float>& scales,
1061                                 CvMat* avg, CvMat* eigenvectors)
1062     {
1063         CvSize patch_size = descriptors[0].GetPatchSize();
1064         IplImage* input_patch;
1065         CvRect roi;
1066
1067         input_patch= cvCreateImage(patch_size, IPL_DEPTH_8U, 1);
1068         roi = cvGetImageROI((IplImage*)patch);
1069
1070         //  float min_distance = 1e10;
1071         std::vector<int> _desc_idxs;
1072         _desc_idxs.resize(n);
1073         std::vector<int> _pose_idxs;
1074         _pose_idxs.resize(n);
1075         std::vector<float> _distances;
1076         _distances.resize(n);
1077
1078
1079         for (int i=0;i<n;i++)
1080         {
1081             distances[i] = 1e10;
1082         }
1083
1084         for(float cur_scale = scale_min; cur_scale < scale_max; cur_scale *= scale_step)
1085         {
1086
1087             CvRect roi_scaled = resize_rect(roi, cur_scale);
1088             cvSetImageROI(patch, roi_scaled);
1089             cvResize(patch, input_patch);
1090
1091
1092
1093             FindOneWayDescriptor(desc_count, descriptors, input_patch, n,_desc_idxs, _pose_idxs, _distances, avg, eigenvectors);
1094             for (int i=0;i<n;i++)
1095             {
1096                 if(_distances[i] < distances[i])
1097                 {
1098                     distances[i] = _distances[i];
1099                     desc_idxs[i] = _desc_idxs[i];
1100                     pose_idxs[i] = _pose_idxs[i];
1101                     scales[i] = cur_scale;
1102                 }
1103             }
1104         }
1105
1106
1107
1108         cvSetImageROI((IplImage*)patch, roi);
1109         cvReleaseImage(&input_patch);
1110     }
1111
1112 #if defined(_KDTREE)
1113     void FindOneWayDescriptorEx(cv::flann::Index* m_pca_descriptors_tree, CvSize patch_size, int m_pca_dim_low,
1114                                 int m_pose_count, IplImage* patch,
1115                                 float scale_min, float scale_max, float scale_step,
1116                                 int& desc_idx, int& pose_idx, float& distance, float& scale,
1117                                 CvMat* avg, CvMat* eigenvectors)
1118     {
1119         IplImage* input_patch;
1120         CvRect roi;
1121
1122         input_patch= cvCreateImage(patch_size, IPL_DEPTH_8U, 1);
1123         roi = cvGetImageROI((IplImage*)patch);
1124
1125         int _desc_idx, _pose_idx;
1126         float _distance;
1127         distance = 1e10;
1128         for(float cur_scale = scale_min; cur_scale < scale_max; cur_scale *= scale_step)
1129         {
1130             //        printf("Scale = %f\n", cur_scale);
1131
1132             CvRect roi_scaled = resize_rect(roi, cur_scale);
1133             cvSetImageROI(patch, roi_scaled);
1134             cvResize(patch, input_patch);
1135
1136             FindOneWayDescriptor(m_pca_descriptors_tree, patch_size, m_pca_dim_low, m_pose_count, input_patch, _desc_idx, _pose_idx, _distance, avg, eigenvectors);
1137             if(_distance < distance)
1138             {
1139                 distance = _distance;
1140                 desc_idx = _desc_idx;
1141                 pose_idx = _pose_idx;
1142                 scale = cur_scale;
1143             }
1144         }
1145
1146
1147         cvSetImageROI((IplImage*)patch, roi);
1148         cvReleaseImage(&input_patch);
1149
1150     }
1151 #endif
1152
1153     const char* OneWayDescriptor::GetFeatureName() const
1154     {
1155         return m_feature_name.c_str();
1156     }
1157
1158     CvPoint OneWayDescriptor::GetCenter() const
1159     {
1160         return m_center;
1161     }
1162
1163     int OneWayDescriptor::GetPCADimLow() const
1164     {
1165         return m_pca_dim_low;
1166     }
1167
1168     int OneWayDescriptor::GetPCADimHigh() const
1169     {
1170         return m_pca_dim_high;
1171     }
1172
1173     CvMat* ConvertImageToMatrix(IplImage* patch)
1174     {
1175         CvRect roi = cvGetImageROI(patch);
1176         CvMat* mat = cvCreateMat(1, roi.width*roi.height, CV_32FC1);
1177
1178         if(patch->depth == 32)
1179         {
1180             for(int y = 0; y < roi.height; y++)
1181             {
1182                 for(int x = 0; x < roi.width; x++)
1183                 {
1184                     mat->data.fl[y*roi.width + x] = *((float*)(patch->imageData + (y + roi.y)*patch->widthStep) + x + roi.x);
1185                 }
1186             }
1187         }
1188         else if(patch->depth == 8)
1189         {
1190             for(int y = 0; y < roi.height; y++)
1191             {
1192                 for(int x = 0; x < roi.width; x++)
1193                 {
1194                     mat->data.fl[y*roi.width + x] = (float)(unsigned char)patch->imageData[(y + roi.y)*patch->widthStep + x + roi.x];
1195                 }
1196             }
1197         }
1198         else
1199         {
1200             printf("Image depth %d is not supported\n", patch->depth);
1201             return 0;
1202         }
1203
1204         return mat;
1205     }
1206
1207     OneWayDescriptorBase::OneWayDescriptorBase(CvSize patch_size, int pose_count, const char* train_path,
1208                                                const char* pca_config, const char* pca_hr_config,
1209                                                const char* pca_desc_config, int pyr_levels,
1210                                                int pca_dim_high, int pca_dim_low) : m_pca_dim_high(pca_dim_high), m_pca_dim_low(pca_dim_low)
1211     {
1212         //      m_pca_descriptors_matrix = 0;
1213         m_patch_size = patch_size;
1214         m_pose_count = pose_count;
1215         m_pyr_levels = pyr_levels;
1216         m_poses = 0;
1217         m_transforms = 0;
1218
1219         m_pca_avg = 0;
1220         m_pca_eigenvectors = 0;
1221         m_pca_hr_avg = 0;
1222         m_pca_hr_eigenvectors = 0;
1223         m_pca_descriptors = 0;
1224
1225         m_descriptors = 0;
1226
1227         if(train_path == 0 || strlen(train_path) == 0)
1228         {
1229             // skip pca loading
1230             return;
1231         }
1232         char pca_config_filename[1024];
1233         sprintf(pca_config_filename, "%s/%s", train_path, pca_config);
1234         readPCAFeatures(pca_config_filename, &m_pca_avg, &m_pca_eigenvectors);
1235         if(pca_hr_config && strlen(pca_hr_config) > 0)
1236         {
1237             char pca_hr_config_filename[1024];
1238             sprintf(pca_hr_config_filename, "%s/%s", train_path, pca_hr_config);
1239             readPCAFeatures(pca_hr_config_filename, &m_pca_hr_avg, &m_pca_hr_eigenvectors);
1240         }
1241
1242         m_pca_descriptors = new OneWayDescriptor[m_pca_dim_high + 1];
1243
1244 #if !defined(_GH_REGIONS)
1245         if(pca_desc_config && strlen(pca_desc_config) > 0)
1246             //    if(0)
1247         {
1248             //printf("Loading the descriptors...");
1249             char pca_desc_config_filename[1024];
1250             sprintf(pca_desc_config_filename, "%s/%s", train_path, pca_desc_config);
1251             LoadPCADescriptors(pca_desc_config_filename);
1252             //printf("done.\n");
1253         }
1254         else
1255         {
1256             printf("Initializing the descriptors...\n");
1257             InitializePoseTransforms();
1258             CreatePCADescriptors();
1259             SavePCADescriptors("pca_descriptors.yml");
1260         }
1261 #endif //_GH_REGIONS
1262         //    SavePCADescriptors("./pca_descriptors.yml");
1263
1264     }
1265
1266     OneWayDescriptorBase::~OneWayDescriptorBase()
1267     {
1268         cvReleaseMat(&m_pca_avg);
1269         cvReleaseMat(&m_pca_eigenvectors);
1270
1271         if(m_pca_hr_eigenvectors)
1272         {
1273             delete[] m_pca_descriptors;
1274             cvReleaseMat(&m_pca_hr_avg);
1275             cvReleaseMat(&m_pca_hr_eigenvectors);
1276         }
1277
1278
1279         delete []m_descriptors;
1280         delete []m_poses;
1281
1282         for(int i = 0; i < m_pose_count; i++)
1283         {
1284             cvReleaseMat(&m_transforms[i]);
1285         }
1286         delete []m_transforms;
1287
1288 #if defined(_KDTREE)
1289         //      if (m_pca_descriptors_matrix)
1290         //              delete m_pca_descriptors_matrix;
1291         cvReleaseMat(&m_pca_descriptors_matrix);
1292         delete m_pca_descriptors_tree;
1293 #endif
1294     }
1295
1296     void OneWayDescriptorBase::InitializePoses()
1297     {
1298         m_poses = new CvAffinePose[m_pose_count];
1299         for(int i = 0; i < m_pose_count; i++)
1300         {
1301             m_poses[i] = GenRandomAffinePose();
1302         }
1303     }
1304
1305     void OneWayDescriptorBase::InitializeTransformsFromPoses()
1306     {
1307         m_transforms = new CvMat*[m_pose_count];
1308         for(int i = 0; i < m_pose_count; i++)
1309         {
1310             m_transforms[i] = cvCreateMat(2, 3, CV_32FC1);
1311             GenerateAffineTransformFromPose(cvSize(m_patch_size.width*2, m_patch_size.height*2), m_poses[i], m_transforms[i]);
1312         }
1313     }
1314
1315     void OneWayDescriptorBase::InitializePoseTransforms()
1316     {
1317         InitializePoses();
1318         InitializeTransformsFromPoses();
1319     }
1320
1321     void OneWayDescriptorBase::InitializeDescriptor(int desc_idx, IplImage* train_image, const KeyPoint& keypoint, const char* feature_label)
1322     {
1323
1324         // TBD add support for octave != 0
1325         CvPoint center = keypoint.pt;
1326
1327         CvRect roi = cvRect(center.x - m_patch_size.width/2, center.y - m_patch_size.height/2, m_patch_size.width, m_patch_size.height);
1328         cvResetImageROI(train_image);
1329         roi = fit_rect_fixedsize(roi, train_image);
1330         cvSetImageROI(train_image, roi);
1331         if(roi.width != m_patch_size.width || roi.height != m_patch_size.height)
1332         {
1333             return;
1334         }
1335
1336         InitializeDescriptor(desc_idx, train_image, feature_label);
1337         cvResetImageROI(train_image);
1338     }
1339
1340     void OneWayDescriptorBase::InitializeDescriptor(int desc_idx, IplImage* train_image, const char* feature_label)
1341     {
1342         m_descriptors[desc_idx].SetPCADimHigh(m_pca_dim_high);
1343         m_descriptors[desc_idx].SetPCADimLow(m_pca_dim_low);
1344         m_descriptors[desc_idx].SetTransforms(m_poses, m_transforms);
1345
1346         if(!m_pca_hr_eigenvectors)
1347         {
1348             m_descriptors[desc_idx].Initialize(m_pose_count, train_image, feature_label);
1349         }
1350         else
1351         {
1352             m_descriptors[desc_idx].InitializeFast(m_pose_count, train_image, feature_label,
1353                                                    m_pca_hr_avg, m_pca_hr_eigenvectors, m_pca_descriptors);
1354         }
1355
1356         if(m_pca_avg)
1357         {
1358             m_descriptors[desc_idx].InitializePCACoeffs(m_pca_avg, m_pca_eigenvectors);
1359         }
1360     }
1361
1362     void OneWayDescriptorBase::FindDescriptor(IplImage* src, cv::Point2f pt, int& desc_idx, int& pose_idx, float& distance) const
1363     {
1364
1365         CvRect roi = cvRect(pt.x - m_patch_size.width/4, pt.y - m_patch_size.height/4, m_patch_size.width/2, m_patch_size.height/2);
1366         cvSetImageROI(src, roi);
1367
1368
1369         FindDescriptor(src, desc_idx, pose_idx, distance);
1370
1371         cvResetImageROI(src);
1372
1373     }
1374
1375     void OneWayDescriptorBase::FindDescriptor(IplImage* patch, int& desc_idx, int& pose_idx, float& distance, float* _scale, float* scale_ranges) const
1376     {
1377 #if 0
1378         ::FindOneWayDescriptor(m_train_feature_count, m_descriptors, patch, desc_idx, pose_idx, distance, m_pca_avg, m_pca_eigenvectors);
1379 #else
1380         float scale_min = 0.7f;
1381         float scale_max = 2.0f;
1382         float scale_step = 1.2f;
1383
1384         if (scale_ranges)
1385         {
1386             scale_min = scale_ranges[0];
1387             scale_max = scale_ranges[1];
1388         }
1389
1390         float scale = 1.0f;
1391
1392 #if !defined(_KDTREE)
1393         cv::FindOneWayDescriptorEx(m_train_feature_count, m_descriptors, patch,
1394                                    scale_min, scale_max, scale_step, desc_idx, pose_idx, distance, scale,
1395                                    m_pca_avg, m_pca_eigenvectors);
1396 #else
1397         cv::FindOneWayDescriptorEx(m_pca_descriptors_tree, m_descriptors[0].GetPatchSize(), m_descriptors[0].GetPCADimLow(), m_pose_count, patch,
1398                                    scale_min, scale_max, scale_step, desc_idx, pose_idx, distance, scale,
1399                                    m_pca_avg, m_pca_eigenvectors);
1400 #endif
1401
1402         if (_scale)
1403             *_scale = scale;
1404
1405 #endif
1406     }
1407
1408     void OneWayDescriptorBase::FindDescriptor(IplImage* patch, int n, std::vector<int>& desc_idxs, std::vector<int>& pose_idxs,
1409                                               std::vector<float>& distances, std::vector<float>& _scales, float* scale_ranges) const
1410     {
1411         float scale_min = 0.7f;
1412         float scale_max = 2.5f;
1413         float scale_step = 1.2f;
1414
1415         if (scale_ranges)
1416         {
1417             scale_min = scale_ranges[0];
1418             scale_max = scale_ranges[1];
1419         }
1420
1421         distances.resize(n);
1422         _scales.resize(n);
1423         desc_idxs.resize(n);
1424         pose_idxs.resize(n);
1425         /*float scales = 1.0f;*/
1426
1427         cv::FindOneWayDescriptorEx(m_train_feature_count, m_descriptors, patch,
1428                                    scale_min, scale_max, scale_step ,n, desc_idxs, pose_idxs, distances, _scales,
1429                                    m_pca_avg, m_pca_eigenvectors);
1430
1431     }
1432
1433     void OneWayDescriptorBase::SetPCAHigh(CvMat* avg, CvMat* eigenvectors)
1434     {
1435         m_pca_hr_avg = cvCloneMat(avg);
1436         m_pca_hr_eigenvectors = cvCloneMat(eigenvectors);
1437     }
1438
1439     void OneWayDescriptorBase::SetPCALow(CvMat* avg, CvMat* eigenvectors)
1440     {
1441         m_pca_avg = cvCloneMat(avg);
1442         m_pca_eigenvectors = cvCloneMat(eigenvectors);
1443     }
1444
1445     void OneWayDescriptorBase::AllocatePCADescriptors()
1446     {
1447         m_pca_descriptors = new OneWayDescriptor[m_pca_dim_high + 1];
1448         for(int i = 0; i < m_pca_dim_high + 1; i++)
1449         {
1450             m_pca_descriptors[i].SetPCADimHigh(m_pca_dim_high);
1451             m_pca_descriptors[i].SetPCADimLow(m_pca_dim_low);
1452         }
1453     }
1454
1455     void OneWayDescriptorBase::CreatePCADescriptors()
1456     {
1457         if(m_pca_descriptors == 0)
1458         {
1459             AllocatePCADescriptors();
1460         }
1461         IplImage* frontal = cvCreateImage(m_patch_size, IPL_DEPTH_32F, 1);
1462
1463         eigenvector2image(m_pca_hr_avg, frontal);
1464         m_pca_descriptors[0].SetTransforms(m_poses, m_transforms);
1465         m_pca_descriptors[0].Initialize(m_pose_count, frontal, "", 0);
1466
1467         for(int j = 0; j < m_pca_dim_high; j++)
1468         {
1469             CvMat eigenvector;
1470             cvGetSubRect(m_pca_hr_eigenvectors, &eigenvector, cvRect(0, j, m_pca_hr_eigenvectors->cols, 1));
1471             eigenvector2image(&eigenvector, frontal);
1472
1473             m_pca_descriptors[j + 1].SetTransforms(m_poses, m_transforms);
1474             m_pca_descriptors[j + 1].Initialize(m_pose_count, frontal, "", 0);
1475
1476             printf("Created descriptor for PCA component %d\n", j);
1477         }
1478
1479         cvReleaseImage(&frontal);
1480     }
1481
1482
1483     int OneWayDescriptorBase::LoadPCADescriptors(const char* filename)
1484     {
1485         CvMemStorage* storage = cvCreateMemStorage();
1486         CvFileStorage* fs = cvOpenFileStorage(filename, storage, CV_STORAGE_READ);
1487         if(!fs)
1488         {
1489             cvReleaseMemStorage(&storage);
1490             printf("File %s not found...\n", filename);
1491             return 0;
1492         }
1493
1494         // read affine poses
1495         CvFileNode* node = cvGetFileNodeByName(fs, 0, "affine poses");
1496         if(node != 0)
1497         {
1498             CvMat* poses = (CvMat*)cvRead(fs, node);
1499             //if(poses->rows != m_pose_count)
1500             //{
1501             //    printf("Inconsistency in the number of poses between the class instance and the file! Exiting...\n");
1502             //    cvReleaseMat(&poses);
1503             //    cvReleaseFileStorage(&fs);
1504             //    cvReleaseMemStorage(&storage);
1505             //}
1506
1507             if(m_poses)
1508             {
1509                 delete m_poses;
1510             }
1511             m_poses = new CvAffinePose[m_pose_count];
1512             for(int i = 0; i < m_pose_count; i++)
1513             {
1514                 m_poses[i].phi = cvmGet(poses, i, 0);
1515                 m_poses[i].theta = cvmGet(poses, i, 1);
1516                 m_poses[i].lambda1 = cvmGet(poses, i, 2);
1517                 m_poses[i].lambda2 = cvmGet(poses, i, 3);
1518             }
1519             cvReleaseMat(&poses);
1520
1521             // now initialize pose transforms
1522             InitializeTransformsFromPoses();
1523         }
1524         else
1525         {
1526             printf("Node \"affine poses\" not found...\n");
1527         }
1528
1529         node = cvGetFileNodeByName(fs, 0, "pca components number");
1530         if(node != 0)
1531         {
1532
1533             m_pca_dim_high = cvReadInt(node);
1534             if(m_pca_descriptors)
1535             {
1536                 delete []m_pca_descriptors;
1537             }
1538             AllocatePCADescriptors();
1539             for(int i = 0; i < m_pca_dim_high + 1; i++)
1540             {
1541                 m_pca_descriptors[i].Allocate(m_pose_count, m_patch_size, 1);
1542                 m_pca_descriptors[i].SetTransforms(m_poses, m_transforms);
1543                 char buf[1024];
1544                 sprintf(buf, "descriptor for pca component %d", i);
1545                 m_pca_descriptors[i].ReadByName(fs, 0, buf);
1546             }
1547         }
1548         else
1549         {
1550             printf("Node \"pca components number\" not found...\n");
1551         }
1552
1553         cvReleaseFileStorage(&fs);
1554         cvReleaseMemStorage(&storage);
1555
1556         printf("Successfully read %d pca components\n", m_pca_dim_high);
1557
1558         return 1;
1559     }
1560
1561     void OneWayDescriptorBase::SavePCADescriptors(const char* filename)
1562     {
1563         CvMemStorage* storage = cvCreateMemStorage();
1564         CvFileStorage* fs = cvOpenFileStorage(filename, storage, CV_STORAGE_WRITE);
1565
1566         cvWriteInt(fs, "pca components number", m_pca_dim_high);
1567         cvWriteComment(fs, "The first component is the average Vector, so the total number of components is <pca components number> + 1", 0);
1568         cvWriteInt(fs, "patch width", m_patch_size.width);
1569         cvWriteInt(fs, "patch height", m_patch_size.height);
1570
1571         // pack the affine transforms into a single CvMat and write them
1572         CvMat* poses = cvCreateMat(m_pose_count, 4, CV_32FC1);
1573         for(int i = 0; i < m_pose_count; i++)
1574         {
1575             cvmSet(poses, i, 0, m_poses[i].phi);
1576             cvmSet(poses, i, 1, m_poses[i].theta);
1577             cvmSet(poses, i, 2, m_poses[i].lambda1);
1578             cvmSet(poses, i, 3, m_poses[i].lambda2);
1579         }
1580         cvWrite(fs, "affine poses", poses);
1581         cvReleaseMat(&poses);
1582
1583         for(int i = 0; i < m_pca_dim_high + 1; i++)
1584         {
1585             char buf[1024];
1586             sprintf(buf, "descriptor for pca component %d", i);
1587             m_pca_descriptors[i].Write(fs, buf);
1588         }
1589
1590         cvReleaseMemStorage(&storage);
1591         cvReleaseFileStorage(&fs);
1592     }
1593
1594     void OneWayDescriptorBase::Allocate(int train_feature_count)
1595     {
1596         m_train_feature_count = train_feature_count;
1597         m_descriptors = new OneWayDescriptor[m_train_feature_count];
1598         for(int i = 0; i < m_train_feature_count; i++)
1599         {
1600             m_descriptors[i].SetPCADimHigh(m_pca_dim_high);
1601             m_descriptors[i].SetPCADimLow(m_pca_dim_low);
1602         }
1603     }
1604
1605     void OneWayDescriptorBase::InitializeDescriptors(IplImage* train_image, const vector<KeyPoint>& features,
1606                                                      const char* feature_label, int desc_start_idx)
1607     {
1608         for(int i = 0; i < (int)features.size(); i++)
1609         {
1610             InitializeDescriptor(desc_start_idx + i, train_image, features[i], feature_label);
1611
1612         }
1613         cvResetImageROI(train_image);
1614
1615 #if defined(_KDTREE)
1616         ConvertDescriptorsArrayToTree();
1617 #endif
1618     }
1619
1620     void OneWayDescriptorBase::CreateDescriptorsFromImage(IplImage* src, const std::vector<KeyPoint>& features)
1621     {
1622         m_train_feature_count = (int)features.size();
1623
1624         m_descriptors = new OneWayDescriptor[m_train_feature_count];
1625
1626         InitializeDescriptors(src, features);
1627
1628     }
1629
1630 #if defined(_KDTREE)
1631     void OneWayDescriptorBase::ConvertDescriptorsArrayToTree()
1632     {
1633         int n = this->GetDescriptorCount();
1634         if (n<1)
1635             return;
1636         int pca_dim_low = this->GetDescriptor(0)->GetPCADimLow();
1637
1638         //if (!m_pca_descriptors_matrix)
1639         //      m_pca_descriptors_matrix = new ::flann::Matrix<float>(n*m_pose_count,pca_dim_low);
1640         //else
1641         //{
1642         //      if ((m_pca_descriptors_matrix->cols != pca_dim_low)&&(m_pca_descriptors_matrix->rows != n*m_pose_count))
1643         //      {
1644         //              delete m_pca_descriptors_matrix;
1645         //              m_pca_descriptors_matrix = new ::flann::Matrix<float>(n*m_pose_count,pca_dim_low);
1646         //      }
1647         //}
1648
1649         m_pca_descriptors_matrix = cvCreateMat(n*m_pose_count,pca_dim_low,CV_32FC1);
1650         for (int i=0;i<n;i++)
1651         {
1652             CvMat** pca_coeffs = m_descriptors[i].GetPCACoeffs();
1653             for (int j = 0;j<m_pose_count;j++)
1654             {
1655                 for (int k=0;k<pca_dim_low;k++)
1656                 {
1657                     m_pca_descriptors_matrix->data.fl[(i*m_pose_count+j)*m_pca_dim_low + k] = pca_coeffs[j]->data.fl[k];
1658                 }
1659             }
1660         }
1661         cv::Mat pca_descriptors_mat(m_pca_descriptors_matrix,false);
1662
1663         //::flann::KDTreeIndexParams params;
1664         //params.trees = 1;
1665         //m_pca_descriptors_tree = new KDTree(pca_descriptors_mat);
1666         m_pca_descriptors_tree = new cv::flann::Index(pca_descriptors_mat,cv::flann::KDTreeIndexParams(1));
1667         //cvReleaseMat(&m_pca_descriptors_matrix);
1668         //m_pca_descriptors_tree->buildIndex();
1669     }
1670 #endif
1671
1672     void OneWayDescriptorObject::Allocate(int train_feature_count, int object_feature_count)
1673     {
1674         OneWayDescriptorBase::Allocate(train_feature_count);
1675         m_object_feature_count = object_feature_count;
1676
1677         m_part_id = new int[m_object_feature_count];
1678     }
1679
1680
1681     void OneWayDescriptorObject::InitializeObjectDescriptors(IplImage* train_image, const vector<KeyPoint>& features,
1682                                                              const char* feature_label, int desc_start_idx, float scale, int is_background)
1683     {
1684         InitializeDescriptors(train_image, features, feature_label, desc_start_idx);
1685
1686         for(int i = 0; i < (int)features.size(); i++)
1687         {
1688             CvPoint center = features[i].pt;
1689
1690             if(!is_background)
1691             {
1692                 // remember descriptor part id
1693                 CvPoint center_scaled = cvPoint(round(center.x*scale), round(center.y*scale));
1694                 m_part_id[i + desc_start_idx] = MatchPointToPart(center_scaled);
1695             }
1696         }
1697         cvResetImageROI(train_image);
1698     }
1699
1700     int OneWayDescriptorObject::IsDescriptorObject(int desc_idx) const
1701     {
1702         return desc_idx < m_object_feature_count ? 1 : 0;
1703     }
1704
1705     int OneWayDescriptorObject::MatchPointToPart(CvPoint pt) const
1706     {
1707         int idx = -1;
1708         const int max_dist = 10;
1709         for(int i = 0; i < (int)m_train_features.size(); i++)
1710         {
1711             if(norm(Point2f(pt) - m_train_features[i].pt) < max_dist)
1712             {
1713                 idx = i;
1714                 break;
1715             }
1716         }
1717
1718         return idx;
1719     }
1720
1721     int OneWayDescriptorObject::GetDescriptorPart(int desc_idx) const
1722     {
1723         //    return MatchPointToPart(GetDescriptor(desc_idx)->GetCenter());
1724         return desc_idx < m_object_feature_count ? m_part_id[desc_idx] : -1;
1725     }
1726
1727     OneWayDescriptorObject::OneWayDescriptorObject(CvSize patch_size, int pose_count, const char* train_path,
1728                                                    const char* pca_config, const char* pca_hr_config, const char* pca_desc_config, int pyr_levels) :
1729     OneWayDescriptorBase(patch_size, pose_count, train_path, pca_config, pca_hr_config, pca_desc_config, pyr_levels)
1730     {
1731         m_part_id = 0;
1732     }
1733
1734     OneWayDescriptorObject::~OneWayDescriptorObject()
1735     {
1736         delete []m_part_id;
1737     }
1738
1739     vector<KeyPoint> OneWayDescriptorObject::_GetLabeledFeatures() const
1740     {
1741         vector<KeyPoint> features;
1742         for(size_t i = 0; i < m_train_features.size(); i++)
1743         {
1744             features.push_back(m_train_features[i]);
1745         }
1746
1747         return features;
1748     }
1749
1750     void eigenvector2image(CvMat* eigenvector, IplImage* img)
1751     {
1752         CvRect roi = cvGetImageROI(img);
1753         if(img->depth == 32)
1754         {
1755             for(int y = 0; y < roi.height; y++)
1756             {
1757                 for(int x = 0; x < roi.width; x++)
1758                 {
1759                     float val = cvmGet(eigenvector, 0, roi.width*y + x);
1760                     *((float*)(img->imageData + (roi.y + y)*img->widthStep) + roi.x + x) = val;
1761                 }
1762             }
1763         }
1764         else
1765         {
1766             for(int y = 0; y < roi.height; y++)
1767             {
1768                 for(int x = 0; x < roi.width; x++)
1769                 {
1770                     float val = cvmGet(eigenvector, 0, roi.width*y + x);
1771                     img->imageData[(roi.y + y)*img->widthStep + roi.x + x] = (unsigned char)val;
1772                 }
1773             }
1774         }
1775     }
1776
1777     void readPCAFeatures(const char* filename, CvMat** avg, CvMat** eigenvectors)
1778     {
1779         CvMemStorage* storage = cvCreateMemStorage();
1780         CvFileStorage* fs = cvOpenFileStorage(filename, storage, CV_STORAGE_READ);
1781         if(!fs)
1782         {
1783             printf("Cannot open file %s! Exiting!", filename);
1784             cvReleaseMemStorage(&storage);
1785         }
1786
1787         CvFileNode* node = cvGetFileNodeByName(fs, 0, "avg");
1788         CvMat* _avg = (CvMat*)cvRead(fs, node);
1789         node = cvGetFileNodeByName(fs, 0, "eigenvectors");
1790         CvMat* _eigenvectors = (CvMat*)cvRead(fs, node);
1791
1792         *avg = cvCloneMat(_avg);
1793         *eigenvectors = cvCloneMat(_eigenvectors);
1794
1795         cvReleaseMat(&_avg);
1796         cvReleaseMat(&_eigenvectors);
1797         cvReleaseFileStorage(&fs);
1798         cvReleaseMemStorage(&storage);
1799     }
1800 }