Main MRPT website > C++ reference for MRPT 1.9.9
CFeatureExtraction_SIFT.cpp
Go to the documentation of this file.
1 /* +------------------------------------------------------------------------+
2  | Mobile Robot Programming Toolkit (MRPT) |
3  | http://www.mrpt.org/ |
4  | |
5  | Copyright (c) 2005-2017, Individual contributors, see AUTHORS file |
6  | See: http://www.mrpt.org/Authors - All rights reserved. |
7  | Released under BSD License. See details in http://www.mrpt.org/License |
8  +------------------------------------------------------------------------+ */
9 
10 #include "vision-precomp.h" // Precompiled headers
11 
12 #include <mrpt/system/os.h>
14 
15 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
16 #include "sift-hess/sift.h"
17 #include "sift-hess/imgfeatures.h"
18 #include "sift-hess/utils.h"
19 #endif
20 
21 // Universal include for all versions of OpenCV
22 #include <mrpt/otherlibs/do_opencv_includes.h>
23 #ifdef HAVE_OPENCV_NONFREE // MRPT_HAS_OPENCV_NONFREE
24 #include <opencv2/nonfree/nonfree.hpp>
25 #endif
26 #ifdef HAVE_OPENCV_FEATURES2D
27 #include <opencv2/features2d/features2d.hpp>
28 #endif
29 #ifdef HAVE_OPENCV_XFEATURES2D
30 #include <opencv2/xfeatures2d.hpp>
31 #endif
32 
33 // TODO: Remove, it's just for GetTempPathA
34 #ifdef MRPT_OS_WINDOWS
35 #include <windows.h>
36 #endif
37 
38 using namespace mrpt;
39 using namespace mrpt::vision;
40 using namespace mrpt::system;
41 using namespace mrpt::utils;
42 using namespace mrpt::math;
43 using namespace std;
44 
45 #if MRPT_HAS_OPENCV
46 #if MRPT_OPENCV_VERSION_NUM >= 0x211
47 using namespace cv;
48 #endif
49 #endif
50 
51 /************************* Local Function Prototypes for Hess' SIFT
52  * *************************/
53 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
54 extern "C" // This is mandatory since the implementations are in ".c" files
55 {
56 IplImage* create_init_img(const IplImage*, int, double);
57 IplImage*** build_gauss_pyr(IplImage*, int, int, double);
58 IplImage*** build_dog_pyr(IplImage***, int, int);
59 CvSeq* scale_space_extrema(IplImage***, int, int, double, int, CvMemStorage*);
60 int is_extremum(IplImage***, int, int, int, int);
61 struct feature* new_feature(void);
62 void calc_feature_scales(CvSeq*, double, int);
63 void adjust_for_img_dbl(CvSeq*);
64 void calc_feature_oris(CvSeq*, IplImage***);
65 void compute_descriptors(CvSeq*, IplImage***, int, int);
66 int feature_cmp(void*, void*, void*);
67 void release_pyr(IplImage****, int, int);
68 }
69 #endif // MRPT_HAS_OPENCV
70 
71 /************************************************************************************************
72 * extractFeaturesSIFT *
73 ************************************************************************************************/
75  const CImage& img, CFeatureList& feats, unsigned int init_ID,
76  unsigned int nDesiredFeatures, const TImageROI& ROI) const
77 {
78  bool usingROI = false;
79  if (ROI.xMin != 0 || ROI.xMax != 0 || ROI.yMin != 0 || ROI.yMax != 0)
80  usingROI = true; // A ROI has been defined
81 
82  // ROI can not be managed properly (yet) with these method, so we extract a
83  // subimage
84 
85  // use a smart pointer so we just copy the pointer if the image is
86  // grayscale, or we'll create a new one if it was RGB:
87  CImage img_grayscale(
88  img, FAST_REF_OR_CONVERT_TO_GRAY); // Was: auxImg::Ptr;
89  if (usingROI)
90  {
91  ASSERT_(
92  ROI.xMin >= 0 && ROI.xMin < ROI.xMax && ROI.xMax < img.getWidth() &&
93  ROI.yMin >= 0 && ROI.yMax < img.getHeight() && ROI.yMin < ROI.yMax);
94  CImage auximg;
95  img_grayscale.extract_patch(
96  auximg, ROI.xMin, ROI.yMin, ROI.xMax - ROI.xMin + 1,
97  ROI.yMax - ROI.yMin + 1); // Subimage in "auxImg"
98  img_grayscale.swap(auximg);
99  }
100 
101  switch (options.SIFTOptions.implementation)
102  {
103  // --------------------------------------------------------------------------------------
104  // Binary in C# -> OPTIONAL: Feature position already computed
105  // --------------------------------------------------------------------------------------
106 
107  case CSBinary:
108  {
109 #ifdef MRPT_OS_WINDOWS
110 
111  char filImg[2000], filOut[2000], filFeat[2000];
112  char paramImg[2000];
113 
114  GetTempPathA(1000, filOut);
115  os::strcat(filOut, 1000, "temp_out.txt"); // OUTPUT FILE
116  GetTempPathA(1000, filImg);
117  os::strcat(
118  filImg, 1000,
119  "temp_img.bmp"); // INPUT IMAGE (BMP) FOR BINARY IN (C#)
120 
121  bool onlyDesc = feats.size() > 0 ? true : false;
122 
123  if (onlyDesc)
124  {
125  GetTempPathA(1000, filFeat);
126  os::strcat(
127  filFeat, 1000, "temp_feats.txt"); // KEYPOINTS INPUT FILE
128  CMatrix listPoints(feats.size(), 2);
129  for (size_t i = 0; i < feats.size(); i++)
130  {
131  listPoints(i, 0) = feats[i]->x;
132  listPoints(i, 1) = feats[i]->y;
133  }
134  listPoints.saveToTextFile(
135  filFeat, MATRIX_FORMAT_FIXED /*Float format*/);
136  } // end if
137  // -------------------------------------------
138  // CALL TO "extractSIFT.exe"
139  // -------------------------------------------
140  img_grayscale.saveToFile(filImg);
141 
142  // ------------------------------------
143  // Version with "CreateProcess":
144  // ------------------------------------
145  os::strcpy(paramImg, 1000, "extractSIFT.exe -i");
146  os::strcat(paramImg, 1000, filImg);
147  os::strcat(paramImg, 1000, " -f");
148  os::strcat(paramImg, 1000, filOut);
149  os::strcat(paramImg, 1000, " -l");
150  os::strcat(paramImg, 1000, filFeat);
151 
152  // ------------------------------------
153  // Launch process
154  // ------------------------------------
155  bool ret = mrpt::system::launchProcess(paramImg);
156 
157  if (!ret)
159  "[extractFeaturesSIFT] Could not launch external "
160  "process... (extractSIFT.exe)")
161 
162  // Process Results
163  CFeatureList::iterator itFeat = feats.begin();
164  size_t nFeats;
165 
166  CMatrix aux;
167  aux.loadFromTextFile(filOut);
168  std::cout << "[computeSiftFeatures] " << aux.getRowCount()
169  << " features." << std::endl;
170 
171  if (onlyDesc)
172  nFeats = feats.size();
173  else
174  {
175  nFeats = aux.getRowCount();
176  feats.resize(nFeats);
177  }
178 
179  for (size_t i = 0; itFeat != feats.end(); i++, itFeat++)
180  {
181  (*itFeat)->type = featSIFT;
182  (*itFeat)->x = usingROI ? aux(i, 0) + ROI.xMin : aux(i, 0);
183  (*itFeat)->y = usingROI ? aux(i, 1) + ROI.yMin : aux(i, 1);
184  (*itFeat)->orientation = aux(i, 2);
185  (*itFeat)->scale = aux(i, 3);
186  (*itFeat)->ID = init_ID + i;
187 
188  // The descriptor:
189  aux.extractRow(i, (*itFeat)->descriptors.SIFT, 4);
190  }
191  remove(filImg);
192  remove(filOut);
193 #else
195  "Unfortunately, this SIFT Implementation only runs in Windows "
196  "OS, try Hess implementation");
197 #endif
198  break;
199  } // end case Binary in C#
200  case VedaldiBinary:
201  {
202 // --------------------------------------------------------------------------------------
203 // Binary by Vedaldi: NOT IMPLEMENTED YET. Input in PGM format
204 // --------------------------------------------------------------------------------------
205 #ifdef MRPT_OS_WINDOWS
207  "Usage of Vedaldi Binary not implemented yet, please, try "
208  "another one");
209 #else
211  "Unfortunately, this SIFT Implementation only runs in Windows "
212  "OS, try Hess implementation");
213 #endif
214  break;
215  } // end case Binary by Vedaldi
216  // --------------------------------------------------------------------------------------
217  // Binary by David Lowe
218  // --------------------------------------------------------------------------------------
219  case LoweBinary: // Binary by Lowe
220  {
221 #ifdef MRPT_OS_WINDOWS
222  char filImg[2000], filOut[2000];
223  char paramImg[2000];
224 
225  feats.clear();
226 
227  GetTempPathA(1000, filOut);
228  os::strcat(filOut, 1000, "temp_out.txt"); // OUTPUT FILE
229  GetTempPathA(1000, filImg);
230  os::strcat(filImg, 1000, "temp_img.pgm"); // INPUT IMAGE (PGM) FOR
231  // ORIGINAL BINARY BY
232  // LOWE
233 
234  bool valid = img_grayscale.saveToFile(filImg);
235  if (!valid)
237  "An error occurred when saving input image into a .pgm "
238  "file");
239 
240  // CONVERT TO UNCOMPRESSED RAW PGM (TODO: Solve in a better way this
241  // issue)
242  os::strcpy(
243  paramImg, 1000,
244  format("cmd /C gmic.exe %s -o %s -quiet", filImg, filImg)
245  .c_str());
246 
247  bool ret = mrpt::system::launchProcess(paramImg);
248 
249  if (!ret)
251  "[extractFeaturesSIFT] Could not launch external "
252  "process... (gmic.exe)");
253 
254  // ------------------------------------
255  // Version with "CreateProcess":
256  // ------------------------------------
257  os::strcpy(paramImg, 1000, "cmd /C siftWin32.exe <");
258  os::strcat(paramImg, 1000, filImg);
259  os::strcat(paramImg, 1000, " >");
260  os::strcat(paramImg, 1000, filOut);
261 
262  ret = mrpt::system::launchProcess(paramImg);
263 
264  if (!ret)
266  "[extractFeaturesSIFT] Could not launch external "
267  "process... (siftWin32.exe)");
268 
269  // ------------------------------------
270  // Process Results
271  // ------------------------------------
272  unsigned int dLen, nFeats;
273  FILE* f = os::fopen(filOut, "rt");
274  if (!f)
276  "Error in extract SIFT with Lowe binary, output file not "
277  "found!");
278  fscanf(
279  f, "%u %u", &nFeats,
280  &dLen); // Number of feats and length of the descriptor
281 
282  for (size_t i = 0; i < nFeats; i++)
283  {
284  CFeature::Ptr feat = mrpt::make_aligned_shared<CFeature>();
285 
286  feat->type = featSIFT; // Type
287  feat->ID = init_ID + i; // Identifier
288 
289  // Position, orientation and scale
290  // IMPORTANTE NOTE: Lowe format stores first the 'y' coordinate
291  // and then the 'x' one
292  float fx, fy, fo, fs;
293  fscanf(f, "%f %f %f %f", &fy, &fx, &fo, &fs);
294 
295  feat->x = usingROI ? fx + ROI.xMin : fx;
296  feat->y = usingROI ? fy + ROI.yMin : fy;
297  feat->orientation = fo;
298  feat->scale = fs;
299 
300  // The descriptor
301  feat->descriptors.SIFT.resize(dLen);
302  unsigned int c;
303  for (unsigned int k = 0; k < dLen; k++)
304  {
305  fscanf(f, "%u", &c);
306  feat->descriptors.SIFT[k] = (unsigned char)c;
307  }
308 
309  feats.push_back(feat);
310  } // end for
311  os::fclose(f);
312  remove(filImg);
313  remove(filOut);
314 #else
316  "Unfortunately, this SIFT Implementation only runs in Windows "
317  "OS, try Hess implementation");
318 #endif
319  break;
320  } // end case Binary by Lowe
321  // --------------------------------------------------------------------------------------
322  // Hess implementation
323  // --------------------------------------------------------------------------------------
324  case Hess: // Implementation by Robert Hess
325  {
326 #if !MRPT_HAS_SIFT_HESS
328  "Method not available since MRPT has been compiled without "
329  "Hess' SIFT library")
330 #elif MRPT_HAS_OPENCV // OK, we have Hess' sift:
331  IplImage* init_img;
332  IplImage ***gauss_pyr, ***dog_pyr;
333  CvMemStorage* storage;
334  CvSeq* features;
335  int octvs;
336  /* check arguments */
337  ASSERT_(
338  img_grayscale.getWidth() != 0 &&
339  img_grayscale.getHeight() != 0);
340  /* build scale space pyramid; smallest dimension of top level is ~4
341  * pixels */
342  const IplImage* ipl_im = img_grayscale.getAs<IplImage>();
343  init_img = create_init_img(ipl_im, SIFT_IMG_DBL, SIFT_SIGMA);
344  octvs = log((float)(MIN(init_img->width, init_img->height))) /
345  log((float)2) -
346  2;
347  gauss_pyr =
348  build_gauss_pyr(init_img, octvs, SIFT_INTVLS, SIFT_SIGMA);
349  dog_pyr = build_dog_pyr(gauss_pyr, octvs, SIFT_INTVLS);
350  storage = cvCreateMemStorage(0);
351  features = scale_space_extrema(
352  dog_pyr, octvs, SIFT_INTVLS,
353  options.SIFTOptions.threshold, // SIFT_CONTR_THR,
354  options.SIFTOptions.edgeThreshold, // SIFT_CURV_THR
355  storage);
356  calc_feature_scales(features, SIFT_SIGMA, SIFT_INTVLS);
357  if (SIFT_IMG_DBL) adjust_for_img_dbl(features);
358  calc_feature_oris(features, gauss_pyr);
359  compute_descriptors(
360  features, gauss_pyr, SIFT_DESCR_WIDTH, SIFT_DESCR_HIST_BINS);
361 
362  /* sort features by decreasing scale and move from CvSeq to array */
363  cvSeqSort(features, (CvCmpFunc)feature_cmp, nullptr);
364 
365  /* get only the desired features */
366  if (nDesiredFeatures > 0)
367  {
368  if (nDesiredFeatures < (unsigned int)features->total)
369  cvSeqPopMulti(
370  features, nullptr, features->total - nDesiredFeatures);
371  else
372  cout << "[Warning] Detected less features than the "
373  "requested "
374  << features->total << " vs " << nDesiredFeatures
375  << endl;
376  } // end if
377 
378  /* convert CvSeq into a FeatureList */
379  convertCvSeqInCFeatureList(features, feats, init_ID, ROI);
380 
381  // clear Hess-features
382  cvClearSeq(features);
383  cvReleaseMemStorage(&storage);
384  cvReleaseImage(&init_img);
385  release_pyr(&gauss_pyr, octvs, SIFT_INTVLS + 3);
386  release_pyr(&dog_pyr, octvs, SIFT_INTVLS + 2);
387 #else
389  "Method not available since MRPT has been compiled without "
390  "OpenCV")
391 #endif // MRPT_HAS_OPENCV
392  break;
393  } // end case Hess
394  //***********************************************************************************************
395  // USING OPENCV
396  //***********************************************************************************************
397  case OpenCV:
398  {
399 #if defined(HAVE_OPENCV_NONFREE) || \
400  defined(HAVE_OPENCV_XFEATURES2D) // MRPT_HAS_OPENCV_NONFREE
401 
402 #if MRPT_OPENCV_VERSION_NUM >= 0x211 && MRPT_OPENCV_VERSION_NUM < 0x300
403 
404  SiftFeatureDetector SIFTDetector(
405  options.SIFTOptions
406  .threshold, // SIFT::DetectorParams::GET_DEFAULT_THRESHOLD(),
407  options.SIFTOptions
408  .edgeThreshold // SIFT::DetectorParams::GET_DEFAULT_EDGE_THRESHOLD()
409  // );
410  );
411 
412  SiftDescriptorExtractor SIFTDescriptor;
413 
414  vector<KeyPoint> cv_feats; // The OpenCV output feature list
415 
416  const IplImage* cGrey = img_grayscale.getAs<IplImage>();
417 
418  Mat theImg = cvarrToMat(cGrey);
419  SIFTDetector.detect(theImg, cv_feats);
420 
421  Mat desc;
422  SIFTDescriptor.compute(theImg, cv_feats, desc);
423 
424  // fromOpenCVToMRPT( theImg, cv_feats, desc, nDesiredFeatures,
425  // outList );
426  const size_t N = cv_feats.size();
427  unsigned int nMax = nDesiredFeatures != 0 && N > nDesiredFeatures
428  ? nDesiredFeatures
429  : N;
430  const int offset = (int)this->options.patchSize / 2 + 1;
431  const size_t size_2 = options.patchSize / 2;
432  const size_t imgH = img.getHeight();
433  const size_t imgW = img.getWidth();
434  unsigned int i = 0;
435  unsigned int cont = 0;
436  TFeatureID nextID = init_ID;
437  feats.clear();
438  while (cont != nMax && i != N)
439  {
440  const int xBorderInf = (int)floor(cv_feats[i].pt.x - size_2);
441  const int xBorderSup = (int)floor(cv_feats[i].pt.x + size_2);
442  const int yBorderInf = (int)floor(cv_feats[i].pt.y - size_2);
443  const int yBorderSup = (int)floor(cv_feats[i].pt.y + size_2);
444 
445  if (options.patchSize == 0 ||
446  ((xBorderSup < (int)imgW) && (xBorderInf > 0) &&
447  (yBorderSup < (int)imgH) && (yBorderInf > 0)))
448  {
449  CFeature::Ptr ft = mrpt::make_aligned_shared<CFeature>();
450  ft->type = featSIFT;
451  ft->ID = nextID++;
452  ft->x = cv_feats[i].pt.x;
453  ft->y = cv_feats[i].pt.y;
454  ft->response = cv_feats[i].response;
455  ft->orientation = cv_feats[i].angle;
456  ft->scale = cv_feats[i].size;
457  ft->patchSize =
458  options.patchSize; // The size of the feature patch
459  ft->descriptors.SIFT.resize(128);
460  memcpy(
461  &(ft->descriptors.SIFT[0]), &desc.data[128 * i],
462  128 *
463  sizeof(ft->descriptors.SIFT[0])); // The descriptor
464 
465  if (options.patchSize > 0)
466  {
467  img.extract_patch(
468  ft->patch, round(ft->x) - offset,
469  round(ft->y) - offset, options.patchSize,
470  options.patchSize); // Image patch surronding the
471  // feature
472  }
473  feats.push_back(ft);
474  ++cont;
475  }
476  ++i;
477  }
478  feats.resize(cont);
479 #else
480  // MRPT_OPENCV_VERSION_NUM >= 0x300
481  using namespace cv;
482  vector<KeyPoint> cv_feats;
483 
484  cv::Ptr<cv::xfeatures2d::SIFT> sift = cv::xfeatures2d::SIFT::create(
485  nDesiredFeatures, 3, options.SIFTOptions.threshold,
486  options.SIFTOptions.edgeThreshold, 1.6); // gb
487  const IplImage* cGrey = img_grayscale.getAs<IplImage>();
488  Mat theImg = cvarrToMat(cGrey);
489  // SIFTDetector.detect(theImg, cv_feats);
490  sift->detect(theImg, cv_feats); // gb
491  Mat desc;
492  // SIFTDescriptor.compute(theImg, cv_feats, desc);
493  sift->compute(theImg, cv_feats, desc);
494 
495  // fromOpenCVToMRPT( theImg, cv_feats, desc, nDesiredFeatures,
496  // outList );
497  const size_t N = cv_feats.size();
498  unsigned int nMax = nDesiredFeatures != 0 && N > nDesiredFeatures
499  ? nDesiredFeatures
500  : N;
501  const int offset = (int)this->options.patchSize / 2 + 1;
502  const size_t size_2 = options.patchSize / 2;
503  const size_t imgH = img.getHeight();
504  const size_t imgW = img.getWidth();
505  unsigned int i = 0;
506  unsigned int cont = 0;
507  TFeatureID nextID = init_ID;
508  feats.clear();
509 
510  while (cont != nMax && i != N)
511  {
512  const int xBorderInf = (int)floor(cv_feats[i].pt.x - size_2);
513  const int xBorderSup = (int)floor(cv_feats[i].pt.x + size_2);
514  const int yBorderInf = (int)floor(cv_feats[i].pt.y - size_2);
515  const int yBorderSup = (int)floor(cv_feats[i].pt.y + size_2);
516 
517  if (options.patchSize == 0 ||
518  ((xBorderSup < (int)imgW) && (xBorderInf > 0) &&
519  (yBorderSup < (int)imgH) && (yBorderInf > 0)))
520  {
521  CFeature::Ptr ft = mrpt::make_aligned_shared<CFeature>();
522  ft->type = featSIFT;
523  ft->ID = nextID++;
524  ft->x = cv_feats[i].pt.x;
525  ft->y = cv_feats[i].pt.y;
526  ft->response = cv_feats[i].response;
527  ft->orientation = cv_feats[i].angle;
528  ft->scale = cv_feats[i].size;
529  ft->patchSize =
530  options.patchSize; // The size of the feature patch
531  ft->descriptors.SIFT.resize(128);
532  memcpy(
533  &(ft->descriptors.SIFT[0]), &desc.data[128 * i],
534  128 *
535  sizeof(ft->descriptors.SIFT[0])); // The descriptor
536 
537  if (options.patchSize > 0)
538  {
539  img.extract_patch(
540  ft->patch, round(ft->x) - offset,
541  round(ft->y) - offset, options.patchSize,
542  options.patchSize); // Image patch surronding the
543  // feature
544  }
545  feats.push_back(ft);
546  ++cont;
547  }
548  ++i;
549  }
550  feats.resize(cont);
551 #endif
552 #else
554  "This method requires OpenCV >= 2.1.1 with nonfree module")
555 #endif
556  break;
557  } // end case OpenCV
558  return;
559  default:
560  {
561  break;
562  } // end default
563  } // end switch
564 } // end extractFeaturesSIFT
565 
566 /************************************************************************************************
567 * computeSiftDescriptors *
568 ************************************************************************************************/
569 // Compute SIFT descriptors on a set of already localized points
571  const CImage& in_img, CFeatureList& in_features) const
572 {
573  ASSERT_(in_features.size() > 0);
574  switch (options.SIFTOptions.implementation)
575  {
576  case CSBinary:
577  {
578 #ifdef MRPT_OS_WINDOWS
579  char filImg[2000], filOut[2000], filFeat[2000];
580  char paramImg[2000];
581 
583 
584  // Save to temporary file
585  GetTempPathA(1000, filImg);
586  os::strcat(filImg, 1000, "temp_img.bmp");
587  GetTempPathA(1000, filOut);
588  os::strcat(filOut, 1000, "temp_feats.txt");
589  GetTempPathA(1000, filFeat);
590  os::strcat(filFeat, 1000, "temp_KLT_feats.txt");
591 
592  // Fill the input file
593  FILE* fout = os::fopen(filFeat, "wt");
594  for (feat = in_features.begin(); feat != in_features.end(); feat++)
595  os::fprintf(fout, "%.6f %.6f\n", (*feat)->x, (*feat)->y);
596 
597  os::fclose(fout);
598 
599  // -------------------------------------------
600  // CALL TO "extractSIFT.exe"
601  // -------------------------------------------
602  in_img.saveToFile(filImg);
603 
604  // Version with "CreateProcess":
605  // ------------------------------------
606  os::strcpy(paramImg, 1000, "extractSIFT.exe -i");
607  os::strcat(paramImg, 1000, filImg);
608  os::strcat(paramImg, 1000, " -f");
609  os::strcat(paramImg, 1000, filOut);
610  os::strcat(paramImg, 1000, " -l");
611  os::strcat(paramImg, 1000, filFeat);
612 
613  bool ret = mrpt::system::launchProcess(paramImg);
614 
615  if (!ret)
617  "[extractFeaturesSIFT] Could not launch external "
618  "process... (extractSIFT.exe)");
619 
620  MRPT_START
621  // Load the results:
622  CMatrix aux;
623  aux.loadFromTextFile(filOut);
624  size_t nRows = aux.getRowCount();
625 
626  std::cout << "[computeSiftFeatures1] " << nRows << " features.\n";
627 
628  unsigned int i;
629  float lx, ly;
630  lx = 0.0;
631  ly = 0.0;
632  feat = in_features.begin();
633 
634  for (i = 0; i < nRows; i++)
635  {
636  if (aux(i, 0) != lx ||
637  aux(i, 1) != ly) // Only one descriptor for feature
638  {
639  (*feat)->type = featSIFT;
640  (*feat)->orientation = aux(i, 2);
641  (*feat)->scale = aux(i, 3);
642 
643  // The descriptor:
644  aux.extractRow(i, (*feat)->descriptors.SIFT, 4);
645 
646  lx = aux(i, 0);
647  ly = aux(i, 1);
648 
649  feat++;
650  } // end if
651  } // end for
652  MRPT_END
653 
654  break;
655 #else
656  THROW_EXCEPTION("TO DO @ linux OS!");
657 #endif
658  } // end case CSBinary
659  case Hess: // Implementation by Hess
660  {
661 #if !MRPT_HAS_SIFT_HESS
663  "Method not available since MRPT has been compiled without "
664  "Hess' SIFT library")
665 #elif MRPT_HAS_OPENCV // OK, we have Hess' sift:
666  IplImage* init_img;
667  IplImage ***gauss_pyr, ***dog_pyr;
668  CvMemStorage* storage;
669  CvSeq* features;
670  int octvs; //, n = 0;
671 
672  /* check arguments */
673  ASSERT_(in_img.getWidth() != 0 && in_img.getHeight() != 0);
674 
675  /* build scale space pyramid; smallest dimension of top level is ~4
676  * pixels */
677  const CImage img_grayscale(in_img, FAST_REF_OR_CONVERT_TO_GRAY);
678  const IplImage* ipl_im = img_grayscale.getAs<IplImage>();
679 
680  init_img = create_init_img(ipl_im, SIFT_IMG_DBL, SIFT_SIGMA);
681  octvs = log((float)(MIN(init_img->width, init_img->height))) /
682  log((float)2) -
683  2;
684  gauss_pyr =
685  build_gauss_pyr(init_img, octvs, SIFT_INTVLS, SIFT_SIGMA);
686  dog_pyr = build_dog_pyr(gauss_pyr, octvs, SIFT_INTVLS);
687 
688  storage = cvCreateMemStorage(0);
689  features = static_cast<CvSeq*>(
690  my_scale_space_extrema(
691  in_features, dog_pyr, octvs, SIFT_INTVLS,
692  options.SIFTOptions.threshold, // SIFT_CONTR_THR,
693  options.SIFTOptions.edgeThreshold, // SIFT_CURV_THR
694  storage));
695  calc_feature_scales(features, SIFT_SIGMA, SIFT_INTVLS);
696  if (SIFT_IMG_DBL) my_adjust_for_img_dbl(features);
697  calc_feature_oris(features, gauss_pyr);
698  compute_descriptors(
699  features, gauss_pyr, SIFT_DESCR_WIDTH, SIFT_DESCR_HIST_BINS);
700 
701  // merge Hess-features and MRPT-features
702  insertCvSeqInCFeatureList(features, in_features);
703 
704  // clear Hess-features
705  cvClearSeq(features);
706 
707  // Free memory
708  cvReleaseMemStorage(&storage);
709  cvReleaseImage(&init_img);
710  release_pyr(&gauss_pyr, octvs, SIFT_INTVLS + 3);
711  release_pyr(&dog_pyr, octvs, SIFT_INTVLS + 2);
712 
713 /* sort features by decreasing scale and move from CvSeq to array */
714 // cvSeqSort( features, (CvCmpFunc)feature_cmp, nullptr );
715 // n = features->total;
716 
717 // struct feature** feat;
718 //*feat = (feature*)calloc( n, sizeof(struct feature) );
719 //*feat = (feature*)cvCvtSeqToArray( features, *feat, CV_WHOLE_SEQ );
720 // for( i = 0; i < n; i++ )
721 //{
722 // free( (*feat)[i].feature_data );
723 // (*feat)[i].feature_data = nullptr;
724 //}
725 
726 //// The number of features must be the same
727 // cout << "NUMBER OF INPUT FEATURES: " << in_features.size() << endl;
728 // cout << "NUMBER OF OUTPUT FEATURES: " << n << endl;
729 
730 //// Transform to CFeatureList:
731 // CFeatureList::iterator itFeat;
732 // int m;
733 // for( itFeat = in_features.begin(), m = 0; itFeat != in_features.end();
734 // itFeat++, m++)
735 //{
736 // (*itFeat)->type = featSIFT;
737 // (*itFeat)->x = (*feat)[m].x;
738 // (*itFeat)->y = (*feat)[m].y;
739 // (*itFeat)->orientation = (*feat)[m].ori;
740 // (*itFeat)->scale = (*feat)[m].scl;
741 // for( unsigned int k = 0; k < FEATURE_MAX_D; k++)
742 // (*itFeat)->descriptors.SIFT[k] = (*feat)[m].descr[k];
743 // (*itFeat)->hasDescriptor = true;
744 //} // end for features
745 #else
747  "Method not available since MRPT has been compiled without "
748  "OpenCV")
749 #endif // MRPT_HAS_OPENCV
750  break;
751  } // end case Hess
752  default:
753  {
754  cout << "SIFT Extraction method not supported for features with "
755  "already known image coordinates"
756  << endl;
757  break;
758  }
759  } // end switch
760 
761 } // end computeSiftDescriptors
762 
763 // ------------------------------------------------------------------------------------
764 // convertCvSeqInCFeatureList
765 // ------------------------------------------------------------------------------------
767  void* features_, CFeatureList& list, unsigned int init_ID,
768  const TImageROI& ROI) const
769 {
770 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
771  CvSeq* features = reinterpret_cast<CvSeq*>(features_);
772 
773  // Is there a defined ROI?
774  bool usingROI = false;
775  if (ROI.xMin != 0 || ROI.xMax != 0 || ROI.yMin != 0 || ROI.yMax != 0)
776  usingROI = true;
777 
778  int n = features->total;
779  ASSERT_(n > 0);
780 
781  list.clear();
782  struct feature* thisFeat;
783  for (int k = 0; k < n; k++)
784  {
785  thisFeat = (feature*)cvGetSeqElem(features, k);
786  CFeature::Ptr feat = mrpt::make_aligned_shared<CFeature>();
787  feat->ID = (TFeatureID)(k + init_ID);
788  feat->x = usingROI ? thisFeat->x + ROI.xMin : thisFeat->x;
789  feat->y = usingROI ? thisFeat->y + ROI.yMin : thisFeat->y;
790  feat->type = featSIFT;
791  feat->orientation = thisFeat->ori;
792  feat->scale = thisFeat->scl;
793  feat->descriptors.SIFT.resize(thisFeat->d);
794  for (int i = 0; i < thisFeat->d; i++)
795  feat->descriptors.SIFT[i] = (unsigned char)thisFeat->descr[i];
796 
797  list.push_back(feat);
798  } // end for
799 #else
801  "Method not available since MRPT has been compiled without OpenCV")
802 #endif // MRPT_HAS_OPENCV
803 }
804 // ------------------------------------------------------------------------------------
805 // insertCvSeqInCFeatureList
806 // ------------------------------------------------------------------------------------
808  void* features_, CFeatureList& list, unsigned int init_ID) const
809 {
810 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
811  CvSeq* features = reinterpret_cast<CvSeq*>(features_);
812 
813  int n = features->total;
814  ASSERT_(n > 0);
815 
816  CFeatureList::iterator itFeat;
817  struct feature* thisFeat;
818  int k;
819  for (itFeat = list.begin(), k = 0; itFeat != list.end() && k < n; k++)
820  {
821  thisFeat = (feature*)cvGetSeqElem(features, k);
822  if ((*itFeat)->x == thisFeat->x && (*itFeat)->y == thisFeat->y)
823  {
824  (*itFeat)->ID = (TFeatureID)(k + init_ID);
825  (*itFeat)->orientation = thisFeat->ori;
826  (*itFeat)->scale = thisFeat->scl;
827  (*itFeat)->descriptors.SIFT.resize(thisFeat->d);
828  for (int i = 0; i < thisFeat->d; i++)
829  (*itFeat)->descriptors.SIFT[i] =
830  (unsigned char)thisFeat->descr[i];
831 
832  itFeat++;
833  } // end if
834  } // end for
835 #else
837  "Method not available since MRPT has been compiled without OpenCV")
838 #endif // MRPT_HAS_OPENCV
839 }
840 
841 /*
842 Halves feature coordinates and scale in case the input image was doubled
843 prior to scale space construction.
844 
845 @param features array of features
846 */
847 void CFeatureExtraction::my_adjust_for_img_dbl(void* features_) const
848 {
849 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
850  CvSeq* features = reinterpret_cast<CvSeq*>(features_);
851  struct feature* feat;
852  int i, n;
853 
854  n = features->total;
855  for (i = 0; i < n; i++)
856  {
857  feat = CV_GET_SEQ_ELEM(struct feature, features, i);
858  feat->scl /= 2.0;
859  }
860 #else
862  "Method not available since MRPT has been compiled without OpenCV")
863 #endif // MRPT_HAS_OPENCV
864 }
865 
866 // ------------------------------------------------------------------------------------
867 // my_scale_space_extrema
868 // ------------------------------------------------------------------------------------
870  CFeatureList& featList, void* dog_pyr_, int octvs, int intvls,
871  double contr_thr, int curv_thr, void* storage_) const
872 {
873  MRPT_UNUSED_PARAM(contr_thr);
874  MRPT_UNUSED_PARAM(curv_thr);
875 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
876  CvMemStorage* storage = reinterpret_cast<CvMemStorage*>(storage_);
877  IplImage*** dog_pyr = reinterpret_cast<IplImage***>(dog_pyr_);
878 
879  CvSeq* features;
880  // double prelim_contr_thr = 0.5 * contr_thr / intvls;
881  struct feature* feat;
882  struct detection_data* ddata;
883  int o, i;
884 
885  CFeatureList::iterator itFeats;
886 
887  features = cvCreateSeq(0, sizeof(CvSeq), sizeof(struct feature), storage);
888  for (itFeats = featList.begin(); itFeats != featList.end(); itFeats++)
889  {
890  // unsigned int nMax, nMin;
891  double gExt = -1e5;
892  int go = 0, gi = 1;
893  float gr = 0.0, gc = 0.0;
894 
895  if ((*itFeats)->y < 0 || (*itFeats)->x < 0)
896  continue; // Do not process this feature
897 
898  // Find the best octave and interval where the feature is the 'best'
899  // extrema
900  for (o = 0; o < octvs; o++)
901  {
902  float r = (*itFeats)->y / pow(2.0, o); // Dog pyramid coordinates
903  float c = (*itFeats)->x / pow(2.0, o);
904 
905  for (i = 1; i <= intvls; i++)
906  {
907  double s = SIFT_SIGMA * pow(2.0, o + i / intvls);
908  double val = s * getLaplacianValue(dog_pyr, o, i, r, c);
909  if (val > gExt)
910  {
911  gExt = val;
912  go = o;
913  gi = i;
914  gr = r;
915  gc = c;
916  }
917  // getTimesExtrema( dog_pyr, o, i, r, c, nMax, nMin ); //
918  // Get
919  // number of times the point is bigger or lower than the
920  // surroundings
921  // if( nMax > gExt || nMin > gExt )
922  //{
923  // gExt = max(nMax,nMin);
924  // go = o;
925  // gi = i;
926  // gr = r;
927  // gc = c;
928  //} // end if
929  } // end for intervals
930  } // end for octaves
931 
932  feat = new_feature(); // Create the new feature
933  ddata = feat_detection_data(feat); // Feature data
934  feat->img_pt.x = feat->x =
935  (*itFeats)->x; // Feature Coordinates in the original image
936  feat->img_pt.y = feat->y = (*itFeats)->y;
937  ddata->r = (int)gr; // Feature Coordinates in the proper octave
938  ddata->c = (int)gc;
939  ddata->octv = go; // Octave
940  ddata->intvl = gi; // Interval
941  ddata->subintvl = 0.0; // Subinterval (not used)
942 
943  cvSeqPush(features, feat); // Insert into the feature sequence
944  } // end for features in the list
945 
946  return features;
947 #else
949  "Method not available since MRPT has been compiled without OpenCV")
950 #endif // MRPT_HAS_OPENCV
951 }
952 
953 /************************************************************************************************
954 * getLaplaceValue *
955 ************************************************************************************************/
957  void* dog_pyr_, int octv, int intvl, float row, float col) const
958 {
959 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
960  IplImage*** dog_pyr = reinterpret_cast<IplImage***>(dog_pyr_);
961  double value = 0.0;
962  int j, k;
963 
964  for (j = -1; j <= 1; j++)
965  for (k = -1; k <= 1; k++)
966  {
967  // Apply convolution mask for computing the Laplacian:
968  // -1 -1 -1
969  // -1 +8 -1
970  // -1 -1 -1
971  if (k == 0 && j == 0)
972  value += 8 * pixval32f(
973  dog_pyr[octv][intvl], (int)(row + j),
974  (int)(col + k));
975  else
976  value -= pixval32f(
977  dog_pyr[octv][intvl], (int)(row + j), (int)(col + k));
978  }
979  return value;
980 #else
982  "Method not available since MRPT has been compiled without OpenCV")
983 #endif // MRPT_HAS_OPENCV
984 } // end getLaplacianValue
985 
986 /************************************************************************************************
987 * getTimesExtrema *
988 ************************************************************************************************/
989 // void CFeatureExtraction::getTimesExtrema( IplImage*** dog_pyr, int octv, int
990 // intvl, float row, float col, unsigned int &nMin, unsigned int &nMax )
992  void* dog_pyr_, int octv, int intvl, float row, float col,
993  unsigned int& nMin, unsigned int& nMax) const
994 {
995 #if MRPT_HAS_OPENCV && MRPT_HAS_SIFT_HESS
996  IplImage*** dog_pyr = reinterpret_cast<IplImage***>(dog_pyr_);
997 
998  float val = pixval32f(dog_pyr[octv][intvl], (int)row, (int)col);
999  int i, j, k;
1000 
1001  nMax = 0;
1002  nMin = 0;
1003  /* check for extrema */
1004  // cout << "VAL: " << pixval32f( dog_pyr[octv][intvl], row, col );
1005  for (i = -1; i <= 1; i++)
1006  for (j = -1; j <= 1; j++)
1007  for (k = -1; k <= 1; k++)
1008  {
1009  // cout << pixval32f( dog_pyr[octv][intvl+i], row + j, col + k )
1010  // << endl;
1011  if (val > pixval32f(
1012  dog_pyr[octv][intvl + i], (int)(row + j),
1013  (int)(col + k)))
1014  nMax++;
1015  if (val < pixval32f(
1016  dog_pyr[octv][intvl + i], (int)(row + j),
1017  (int)(col + k)))
1018  nMin++;
1019  }
1020 #else
1022  "Method not available since MRPT has been compiled without OpenCV")
1023 #endif // MRPT_HAS_OPENCV
1024 } // end getTimesExtrema
std::shared_ptr< CFeature > Ptr
Definition: CFeature.h:58
double scl
scale of a Lowe-style feature
Definition: imgfeatures.h:56
Classes for serialization, sockets, ini-file manipulation, streams, list of properties-values, timewatch, extensions to STL.
void internal_computeSiftDescriptors(const mrpt::utils::CImage &in_img, CFeatureList &in_features) const
Compute the SIFT descriptor of the provided features into the input image.
int octv
Definition: sift.h:28
double subintvl
Definition: sift.h:30
void getTimesExtrema(void *dog_pyr, int octvs, int intvls, float row, float col, unsigned int &nMin, unsigned int &nMax) const
Gets the number of times that a point in the image is higher or lower than the surroundings in the im...
This namespace provides a OS-independent interface to many useful functions: filenames manipulation...
Definition: math_frwds.h:30
int void fclose(FILE *f)
An OS-independent version of fclose.
Definition: os.cpp:272
A class for storing images as grayscale or RGB bitmaps.
Definition: CImage.h:118
static __inline float pixval32f(IplImage *img, int r, int c)
A function to get a pixel value from a 32-bit floating-point image.
double x
x coord
Definition: imgfeatures.h:51
size_t size() const
Definition: CFeature.h:387
#define THROW_EXCEPTION(msg)
#define SIFT_DESCR_HIST_BINS
default number of bins per histogram in descriptor array
Definition: sift.h:57
GLenum GLsizei n
Definition: glext.h:5074
GLintptr offset
Definition: glext.h:3925
CvPoint2D64f img_pt
location in image
Definition: imgfeatures.h:65
double ori
orientation of a Lowe-style feature
Definition: imgfeatures.h:57
#define MIN(a, b)
Definition: jpegint.h:282
STL namespace.
GLdouble s
Definition: glext.h:3676
const T * getAs() const
Returns a pointer to a const T* containing the image - the idea is to call like "img.getAs<IplImage>()" so we can avoid here including OpenCV&#39;s headers.
Definition: CImage.h:587
float yMin
Y coordinate limits [0,imageHeight)
int intvl
Definition: sift.h:29
A structure for defining a ROI within an image.
void extractFeaturesSIFT(const mrpt::utils::CImage &img, CFeatureList &feats, unsigned int init_ID=0, unsigned int nDesiredFeatures=0, const TImageROI &ROI=TImageROI()) const
Extract features from the image based on the SIFT method.
int d
descriptor length
Definition: imgfeatures.h:58
This base provides a set of functions for maths stuff.
Definition: CArrayNumeric.h:19
size_t getWidth() const override
Returns the width of the image in pixels.
Definition: CImage.cpp:869
void resize(size_t N)
Definition: CFeature.h:393
#define MRPT_END
#define MRPT_UNUSED_PARAM(a)
Can be used to avoid "not used parameters" warnings from the compiler.
double getLaplacianValue(void *dog_pyr, int octvs, int intvls, float row, float col) const
Computes the Laplacian value of the feature in the corresponing image in the pyramid.
const GLubyte * c
Definition: glext.h:6313
GLint GLvoid * img
Definition: glext.h:3763
Scale Invariant Feature Transform [LOWE&#39;04].
size_t getHeight() const override
Returns the height of the image in pixels.
Definition: CImage.cpp:897
int val
Definition: mrpt_jpeglib.h:955
std::string format(const char *fmt,...) MRPT_printf_format_check(1
A std::string version of C sprintf.
Definition: format.cpp:19
Classes for computer vision, detectors, features, etc.
#define SIFT_DESCR_WIDTH
default width of descriptor histogram array
Definition: sift.h:54
void insertCvSeqInCFeatureList(void *features, CFeatureList &list, unsigned int init_ID=0) const
Append a sequence of openCV features into an MRPT feature list.
uint64_t TFeatureID
Definition of a feature ID.
holds feature data relevant to detection
Definition: sift.h:24
void swap(CImage &o)
Very efficient swap of two images (just swap the internal pointers)
Definition: CImage.cpp:138
A list of visual features, to be used as output by detectors, as input/output by trackers, etc.
Definition: CFeature.h:305
#define feat_detection_data(f)
Definition: sift.h:93
double descr[FEATURE_MAX_D]
descriptor
Definition: imgfeatures.h:59
float xMin
X coordinate limits [0,imageWidth)
int fprintf(FILE *fil, const char *format,...) noexcept MRPT_printf_format_check(2
An OS-independent version of fprintf.
Definition: os.cpp:405
void my_adjust_for_img_dbl(void *features) const
Adjust scale if the image was initially doubled.
void saveToTextFile(const std::string &fileName, bool APPEND=false)
Save feature list to a text file.
Definition: CFeature.cpp:1035
#define MRPT_START
void * my_scale_space_extrema(CFeatureList &featList, void *dog_pyr, int octvs, int intvls, double contr_thr, int curv_thr, void *storage) const
Computes extrema in the scale space.
This is the global namespace for all Mobile Robot Programming Toolkit (MRPT) libraries.
GLdouble GLdouble GLdouble r
Definition: glext.h:3705
char * strcat(char *dest, size_t destSize, const char *source) noexcept
An OS-independent version of strcat.
Definition: os.cpp:281
Structure to represent an affine invariant image feature.
Definition: imgfeatures.h:49
#define SIFT_INTVLS
default number of sampled intervals per octave
Definition: sift.h:39
GLenum GLenum GLvoid * row
Definition: glext.h:3576
void extract_patch(CImage &patch, const unsigned int col=0, const unsigned int row=0, const unsigned int width=1, const unsigned int height=1) const
Extract a patch from this image, saveing it into "patch" (its previous contents will be overwritten)...
Definition: CImage.cpp:1367
#define ASSERT_(f)
int round(const T value)
Returns the closer integer (int) to x.
Definition: round.h:25
bool saveToFile(const std::string &fileName, int jpeg_quality=95) const
Save the image to a file, whose format is determined from the extension (internally uses OpenCV)...
Definition: CImage.cpp:301
TInternalFeatList::iterator iterator
Definition: CFeature.h:366
FILE * fopen(const char *fileName, const char *mode) noexcept
An OS-independent version of fopen.
Definition: os.cpp:254
GLsizei const GLfloat * value
Definition: glext.h:4117
#define SIFT_SIGMA
default sigma for initial gaussian smoothing
Definition: sift.h:42
#define SIFT_IMG_DBL
default threshold on keypoint contrast |D(x)|
Definition: sift.h:51
double y
y coord
Definition: imgfeatures.h:52
char * strcpy(char *dest, size_t destSize, const char *source) noexcept
An OS-independent version of strcpy.
Definition: os.cpp:296
This class is a "CSerializable" wrapper for "CMatrixFloat".
Definition: CMatrix.h:25
bool launchProcess(const std::string &command)
Executes the given command (which may contain a program + arguments), and waits until it finishes...
Definition: os.cpp:559
void convertCvSeqInCFeatureList(void *features, CFeatureList &list, unsigned int init_ID=0, const TImageROI &ROI=TImageROI()) const
Converts a sequence of openCV features into an MRPT feature list.
fixed floating point &#39;f&#39;
Definition: math_frwds.h:76
void push_back(const CFeature::Ptr &f)
Definition: CFeature.h:399
void memcpy(void *dest, size_t destSize, const void *src, size_t copyCount) noexcept
An OS and compiler independent version of "memcpy".
Definition: os.cpp:355



Page generated by Doxygen 1.8.14 for MRPT 1.9.9 Git: ae4571287 Thu Nov 23 00:06:53 2017 +0100 at dom oct 27 23:51:55 CET 2019