lsh_table.h 19 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513
  1. /***********************************************************************
  2. * Software License Agreement (BSD License)
  3. *
  4. * Copyright 2008-2009 Marius Muja (mariusm@cs.ubc.ca). All rights reserved.
  5. * Copyright 2008-2009 David G. Lowe (lowe@cs.ubc.ca). All rights reserved.
  6. *
  7. * THE BSD LICENSE
  8. *
  9. * Redistribution and use in source and binary forms, with or without
  10. * modification, are permitted provided that the following conditions
  11. * are met:
  12. *
  13. * 1. Redistributions of source code must retain the above copyright
  14. * notice, this list of conditions and the following disclaimer.
  15. * 2. Redistributions in binary form must reproduce the above copyright
  16. * notice, this list of conditions and the following disclaimer in the
  17. * documentation and/or other materials provided with the distribution.
  18. *
  19. * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
  20. * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
  21. * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
  22. * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
  23. * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
  24. * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
  25. * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
  26. * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
  27. * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
  28. * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  29. *************************************************************************/
  30. /***********************************************************************
  31. * Author: Vincent Rabaud
  32. *************************************************************************/
  33. #ifndef OPENCV_FLANN_LSH_TABLE_H_
  34. #define OPENCV_FLANN_LSH_TABLE_H_
  35. #include <algorithm>
  36. #include <iostream>
  37. #include <iomanip>
  38. #include <limits.h>
  39. // TODO as soon as we use C++0x, use the code in USE_UNORDERED_MAP
  40. #ifdef __GXX_EXPERIMENTAL_CXX0X__
  41. # define USE_UNORDERED_MAP 1
  42. #else
  43. # define USE_UNORDERED_MAP 0
  44. #endif
  45. #if USE_UNORDERED_MAP
  46. #include <unordered_map>
  47. #else
  48. #include <map>
  49. #endif
  50. #include <math.h>
  51. #include <stddef.h>
  52. #include "dynamic_bitset.h"
  53. #include "matrix.h"
  54. namespace cvflann
  55. {
  56. namespace lsh
  57. {
  58. ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  59. /** What is stored in an LSH bucket
  60. */
  61. typedef uint32_t FeatureIndex;
  62. /** The id from which we can get a bucket back in an LSH table
  63. */
  64. typedef unsigned int BucketKey;
  65. /** A bucket in an LSH table
  66. */
  67. typedef std::vector<FeatureIndex> Bucket;
  68. ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  69. /** POD for stats about an LSH table
  70. */
  71. struct LshStats
  72. {
  73. std::vector<unsigned int> bucket_sizes_;
  74. size_t n_buckets_;
  75. size_t bucket_size_mean_;
  76. size_t bucket_size_median_;
  77. size_t bucket_size_min_;
  78. size_t bucket_size_max_;
  79. size_t bucket_size_std_dev;
  80. /** Each contained vector contains three value: beginning/end for interval, number of elements in the bin
  81. */
  82. std::vector<std::vector<unsigned int> > size_histogram_;
  83. };
  84. /** Overload the << operator for LshStats
  85. * @param out the streams
  86. * @param stats the stats to display
  87. * @return the streams
  88. */
  89. inline std::ostream& operator <<(std::ostream& out, const LshStats& stats)
  90. {
  91. int w = 20;
  92. out << "Lsh Table Stats:\n" << std::setw(w) << std::setiosflags(std::ios::right) << "N buckets : "
  93. << stats.n_buckets_ << "\n" << std::setw(w) << std::setiosflags(std::ios::right) << "mean size : "
  94. << std::setiosflags(std::ios::left) << stats.bucket_size_mean_ << "\n" << std::setw(w)
  95. << std::setiosflags(std::ios::right) << "median size : " << stats.bucket_size_median_ << "\n" << std::setw(w)
  96. << std::setiosflags(std::ios::right) << "min size : " << std::setiosflags(std::ios::left)
  97. << stats.bucket_size_min_ << "\n" << std::setw(w) << std::setiosflags(std::ios::right) << "max size : "
  98. << std::setiosflags(std::ios::left) << stats.bucket_size_max_;
  99. // Display the histogram
  100. out << std::endl << std::setw(w) << std::setiosflags(std::ios::right) << "histogram : "
  101. << std::setiosflags(std::ios::left);
  102. for (std::vector<std::vector<unsigned int> >::const_iterator iterator = stats.size_histogram_.begin(), end =
  103. stats.size_histogram_.end(); iterator != end; ++iterator) out << (*iterator)[0] << "-" << (*iterator)[1] << ": " << (*iterator)[2] << ", ";
  104. return out;
  105. }
  106. ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  107. /** Lsh hash table. As its key is a sub-feature, and as usually
  108. * the size of it is pretty small, we keep it as a continuous memory array.
  109. * The value is an index in the corpus of features (we keep it as an unsigned
  110. * int for pure memory reasons, it could be a size_t)
  111. */
  112. template<typename ElementType>
  113. class LshTable
  114. {
  115. public:
  116. /** A container of all the feature indices. Optimized for space
  117. */
  118. #if USE_UNORDERED_MAP
  119. typedef std::unordered_map<BucketKey, Bucket> BucketsSpace;
  120. #else
  121. typedef std::map<BucketKey, Bucket> BucketsSpace;
  122. #endif
  123. /** A container of all the feature indices. Optimized for speed
  124. */
  125. typedef std::vector<Bucket> BucketsSpeed;
  126. /** Default constructor
  127. */
  128. LshTable()
  129. {
  130. key_size_ = 0;
  131. feature_size_ = 0;
  132. speed_level_ = kArray;
  133. }
  134. /** Default constructor
  135. * Create the mask and allocate the memory
  136. * @param feature_size is the size of the feature (considered as a ElementType[])
  137. * @param key_size is the number of bits that are turned on in the feature
  138. */
  139. LshTable(unsigned int feature_size, unsigned int key_size)
  140. {
  141. feature_size_ = feature_size;
  142. CV_UNUSED(key_size);
  143. std::cerr << "LSH is not implemented for that type" << std::endl;
  144. assert(0);
  145. }
  146. /** Add a feature to the table
  147. * @param value the value to store for that feature
  148. * @param feature the feature itself
  149. */
  150. void add(unsigned int value, const ElementType* feature)
  151. {
  152. // Add the value to the corresponding bucket
  153. BucketKey key = (lsh::BucketKey)getKey(feature);
  154. switch (speed_level_) {
  155. case kArray:
  156. // That means we get the buckets from an array
  157. buckets_speed_[key].push_back(value);
  158. break;
  159. case kBitsetHash:
  160. // That means we can check the bitset for the presence of a key
  161. key_bitset_.set(key);
  162. buckets_space_[key].push_back(value);
  163. break;
  164. case kHash:
  165. {
  166. // That means we have to check for the hash table for the presence of a key
  167. buckets_space_[key].push_back(value);
  168. break;
  169. }
  170. }
  171. }
  172. /** Add a set of features to the table
  173. * @param dataset the values to store
  174. */
  175. void add(Matrix<ElementType> dataset)
  176. {
  177. #if USE_UNORDERED_MAP
  178. buckets_space_.rehash((buckets_space_.size() + dataset.rows) * 1.2);
  179. #endif
  180. // Add the features to the table
  181. for (unsigned int i = 0; i < dataset.rows; ++i) add(i, dataset[i]);
  182. // Now that the table is full, optimize it for speed/space
  183. optimize();
  184. }
  185. /** Get a bucket given the key
  186. * @param key
  187. * @return
  188. */
  189. inline const Bucket* getBucketFromKey(BucketKey key) const
  190. {
  191. // Generate other buckets
  192. switch (speed_level_) {
  193. case kArray:
  194. // That means we get the buckets from an array
  195. return &buckets_speed_[key];
  196. break;
  197. case kBitsetHash:
  198. // That means we can check the bitset for the presence of a key
  199. if (key_bitset_.test(key)) return &buckets_space_.find(key)->second;
  200. else return 0;
  201. break;
  202. case kHash:
  203. {
  204. // That means we have to check for the hash table for the presence of a key
  205. BucketsSpace::const_iterator bucket_it, bucket_end = buckets_space_.end();
  206. bucket_it = buckets_space_.find(key);
  207. // Stop here if that bucket does not exist
  208. if (bucket_it == bucket_end) return 0;
  209. else return &bucket_it->second;
  210. break;
  211. }
  212. }
  213. return 0;
  214. }
  215. /** Compute the sub-signature of a feature
  216. */
  217. size_t getKey(const ElementType* /*feature*/) const
  218. {
  219. std::cerr << "LSH is not implemented for that type" << std::endl;
  220. assert(0);
  221. return 1;
  222. }
  223. /** Get statistics about the table
  224. * @return
  225. */
  226. LshStats getStats() const;
  227. private:
  228. /** defines the speed fo the implementation
  229. * kArray uses a vector for storing data
  230. * kBitsetHash uses a hash map but checks for the validity of a key with a bitset
  231. * kHash uses a hash map only
  232. */
  233. enum SpeedLevel
  234. {
  235. kArray, kBitsetHash, kHash
  236. };
  237. /** Initialize some variables
  238. */
  239. void initialize(size_t key_size)
  240. {
  241. const size_t key_size_lower_bound = 1;
  242. //a value (size_t(1) << key_size) must fit the size_t type so key_size has to be strictly less than size of size_t
  243. const size_t key_size_upper_bound = (std::min)(sizeof(BucketKey) * CHAR_BIT + 1, sizeof(size_t) * CHAR_BIT);
  244. if (key_size < key_size_lower_bound || key_size >= key_size_upper_bound)
  245. {
  246. CV_Error(cv::Error::StsBadArg, cv::format("Invalid key_size (=%d). Valid values for your system are %d <= key_size < %d.", (int)key_size, (int)key_size_lower_bound, (int)key_size_upper_bound));
  247. }
  248. speed_level_ = kHash;
  249. key_size_ = (unsigned)key_size;
  250. }
  251. /** Optimize the table for speed/space
  252. */
  253. void optimize()
  254. {
  255. // If we are already using the fast storage, no need to do anything
  256. if (speed_level_ == kArray) return;
  257. // Use an array if it will be more than half full
  258. if (buckets_space_.size() > ((size_t(1) << key_size_) / 2)) {
  259. speed_level_ = kArray;
  260. // Fill the array version of it
  261. buckets_speed_.resize(size_t(1) << key_size_);
  262. for (BucketsSpace::const_iterator key_bucket = buckets_space_.begin(); key_bucket != buckets_space_.end(); ++key_bucket) buckets_speed_[key_bucket->first] = key_bucket->second;
  263. // Empty the hash table
  264. buckets_space_.clear();
  265. return;
  266. }
  267. // If the bitset is going to use less than 10% of the RAM of the hash map (at least 1 size_t for the key and two
  268. // for the vector) or less than 512MB (key_size_ <= 30)
  269. if (((std::max(buckets_space_.size(), buckets_speed_.size()) * CHAR_BIT * 3 * sizeof(BucketKey)) / 10
  270. >= (size_t(1) << key_size_)) || (key_size_ <= 32)) {
  271. speed_level_ = kBitsetHash;
  272. key_bitset_.resize(size_t(1) << key_size_);
  273. key_bitset_.reset();
  274. // Try with the BucketsSpace
  275. for (BucketsSpace::const_iterator key_bucket = buckets_space_.begin(); key_bucket != buckets_space_.end(); ++key_bucket) key_bitset_.set(key_bucket->first);
  276. }
  277. else {
  278. speed_level_ = kHash;
  279. key_bitset_.clear();
  280. }
  281. }
  282. /** The vector of all the buckets if they are held for speed
  283. */
  284. BucketsSpeed buckets_speed_;
  285. /** The hash table of all the buckets in case we cannot use the speed version
  286. */
  287. BucketsSpace buckets_space_;
  288. /** What is used to store the data */
  289. SpeedLevel speed_level_;
  290. /** If the subkey is small enough, it will keep track of which subkeys are set through that bitset
  291. * That is just a speedup so that we don't look in the hash table (which can be mush slower that checking a bitset)
  292. */
  293. DynamicBitset key_bitset_;
  294. /** The size of the sub-signature in bits
  295. */
  296. unsigned int key_size_;
  297. unsigned int feature_size_;
  298. // Members only used for the unsigned char specialization
  299. /** The mask to apply to a feature to get the hash key
  300. * Only used in the unsigned char case
  301. */
  302. std::vector<size_t> mask_;
  303. };
  304. ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  305. // Specialization for unsigned char
  306. template<>
  307. inline LshTable<unsigned char>::LshTable(unsigned int feature_size, unsigned int subsignature_size)
  308. {
  309. feature_size_ = feature_size;
  310. initialize(subsignature_size);
  311. // Allocate the mask
  312. mask_ = std::vector<size_t>((feature_size * sizeof(char) + sizeof(size_t) - 1) / sizeof(size_t), 0);
  313. // A bit brutal but fast to code
  314. std::vector<int> indices(feature_size * CHAR_BIT);
  315. for (size_t i = 0; i < feature_size * CHAR_BIT; ++i) indices[i] = (int)i;
  316. #ifndef OPENCV_FLANN_USE_STD_RAND
  317. cv::randShuffle(indices);
  318. #else
  319. std::random_shuffle(indices.begin(), indices.end());
  320. #endif
  321. // Generate a random set of order of subsignature_size_ bits
  322. for (unsigned int i = 0; i < key_size_; ++i) {
  323. size_t index = indices[i];
  324. // Set that bit in the mask
  325. size_t divisor = CHAR_BIT * sizeof(size_t);
  326. size_t idx = index / divisor; //pick the right size_t index
  327. mask_[idx] |= size_t(1) << (index % divisor); //use modulo to find the bit offset
  328. }
  329. // Set to 1 if you want to display the mask for debug
  330. #if 0
  331. {
  332. size_t bcount = 0;
  333. BOOST_FOREACH(size_t mask_block, mask_){
  334. out << std::setw(sizeof(size_t) * CHAR_BIT / 4) << std::setfill('0') << std::hex << mask_block
  335. << std::endl;
  336. bcount += __builtin_popcountll(mask_block);
  337. }
  338. out << "bit count : " << std::dec << bcount << std::endl;
  339. out << "mask size : " << mask_.size() << std::endl;
  340. return out;
  341. }
  342. #endif
  343. }
  344. /** Return the Subsignature of a feature
  345. * @param feature the feature to analyze
  346. */
  347. template<>
  348. inline size_t LshTable<unsigned char>::getKey(const unsigned char* feature) const
  349. {
  350. // no need to check if T is dividable by sizeof(size_t) like in the Hamming
  351. // distance computation as we have a mask
  352. // FIXIT: This is bad assumption, because we reading tail bytes after of the allocated features buffer
  353. const size_t* feature_block_ptr = reinterpret_cast<const size_t*> ((const void*)feature);
  354. // Figure out the subsignature of the feature
  355. // Given the feature ABCDEF, and the mask 001011, the output will be
  356. // 000CEF
  357. size_t subsignature = 0;
  358. size_t bit_index = 1;
  359. for (unsigned i = 0; i < feature_size_; i += sizeof(size_t)) {
  360. // get the mask and signature blocks
  361. size_t feature_block;
  362. if (i <= feature_size_ - sizeof(size_t))
  363. {
  364. feature_block = *feature_block_ptr;
  365. }
  366. else
  367. {
  368. size_t tmp = 0;
  369. memcpy(&tmp, feature_block_ptr, feature_size_ - i); // preserve bytes order
  370. feature_block = tmp;
  371. }
  372. size_t mask_block = mask_[i / sizeof(size_t)];
  373. while (mask_block) {
  374. // Get the lowest set bit in the mask block
  375. size_t lowest_bit = mask_block & (-(ptrdiff_t)mask_block);
  376. // Add it to the current subsignature if necessary
  377. subsignature += (feature_block & lowest_bit) ? bit_index : 0;
  378. // Reset the bit in the mask block
  379. mask_block ^= lowest_bit;
  380. // increment the bit index for the subsignature
  381. bit_index <<= 1;
  382. }
  383. // Check the next feature block
  384. ++feature_block_ptr;
  385. }
  386. return subsignature;
  387. }
  388. template<>
  389. inline LshStats LshTable<unsigned char>::getStats() const
  390. {
  391. LshStats stats;
  392. stats.bucket_size_mean_ = 0;
  393. if ((buckets_speed_.empty()) && (buckets_space_.empty())) {
  394. stats.n_buckets_ = 0;
  395. stats.bucket_size_median_ = 0;
  396. stats.bucket_size_min_ = 0;
  397. stats.bucket_size_max_ = 0;
  398. return stats;
  399. }
  400. if (!buckets_speed_.empty()) {
  401. for (BucketsSpeed::const_iterator pbucket = buckets_speed_.begin(); pbucket != buckets_speed_.end(); ++pbucket) {
  402. stats.bucket_sizes_.push_back((lsh::FeatureIndex)pbucket->size());
  403. stats.bucket_size_mean_ += pbucket->size();
  404. }
  405. stats.bucket_size_mean_ /= buckets_speed_.size();
  406. stats.n_buckets_ = buckets_speed_.size();
  407. }
  408. else {
  409. for (BucketsSpace::const_iterator x = buckets_space_.begin(); x != buckets_space_.end(); ++x) {
  410. stats.bucket_sizes_.push_back((lsh::FeatureIndex)x->second.size());
  411. stats.bucket_size_mean_ += x->second.size();
  412. }
  413. stats.bucket_size_mean_ /= buckets_space_.size();
  414. stats.n_buckets_ = buckets_space_.size();
  415. }
  416. std::sort(stats.bucket_sizes_.begin(), stats.bucket_sizes_.end());
  417. // BOOST_FOREACH(int size, stats.bucket_sizes_)
  418. // std::cout << size << " ";
  419. // std::cout << std::endl;
  420. stats.bucket_size_median_ = stats.bucket_sizes_[stats.bucket_sizes_.size() / 2];
  421. stats.bucket_size_min_ = stats.bucket_sizes_.front();
  422. stats.bucket_size_max_ = stats.bucket_sizes_.back();
  423. // TODO compute mean and std
  424. /*float mean, stddev;
  425. stats.bucket_size_mean_ = mean;
  426. stats.bucket_size_std_dev = stddev;*/
  427. // Include a histogram of the buckets
  428. unsigned int bin_start = 0;
  429. unsigned int bin_end = 20;
  430. bool is_new_bin = true;
  431. for (std::vector<unsigned int>::iterator iterator = stats.bucket_sizes_.begin(), end = stats.bucket_sizes_.end(); iterator
  432. != end; )
  433. if (*iterator < bin_end) {
  434. if (is_new_bin) {
  435. stats.size_histogram_.push_back(std::vector<unsigned int>(3, 0));
  436. stats.size_histogram_.back()[0] = bin_start;
  437. stats.size_histogram_.back()[1] = bin_end - 1;
  438. is_new_bin = false;
  439. }
  440. ++stats.size_histogram_.back()[2];
  441. ++iterator;
  442. }
  443. else {
  444. bin_start += 20;
  445. bin_end += 20;
  446. is_new_bin = true;
  447. }
  448. return stats;
  449. }
  450. // End the two namespaces
  451. }
  452. }
  453. ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  454. #endif /* OPENCV_FLANN_LSH_TABLE_H_ */