nanoflann
C++ header-only ANN library
Loading...
Searching...
No Matches
nanoflann.hpp
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 * Copyright 2011-2025 Jose Luis Blanco (joseluisblancoc@gmail.com).
7 * All rights reserved.
8 *
9 * THE BSD LICENSE
10 *
11 * Redistribution and use in source and binary forms, with or without
12 * modification, are permitted provided that the following conditions
13 * are met:
14 *
15 * 1. Redistributions of source code must retain the above copyright
16 * notice, this list of conditions and the following disclaimer.
17 * 2. Redistributions in binary form must reproduce the above copyright
18 * notice, this list of conditions and the following disclaimer in the
19 * documentation and/or other materials provided with the distribution.
20 *
21 * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
22 * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
23 * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
24 * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
25 * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
26 * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
27 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
28 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
29 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
30 * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
31 *************************************************************************/
32
45#pragma once
46
47#include <algorithm>
48#include <array>
49#include <atomic>
50#include <cassert>
51#include <cmath> // for abs()
52#include <cstdint>
53#include <cstdlib> // for abs()
54#include <functional> // std::reference_wrapper
55#include <future>
56#include <istream>
57#include <limits> // std::numeric_limits
58#include <ostream>
59#include <stdexcept>
60#include <unordered_set>
61#include <vector>
62
64#define NANOFLANN_VERSION 0x171
65
66// Avoid conflicting declaration of min/max macros in Windows headers
67#if !defined(NOMINMAX) && \
68 (defined(_WIN32) || defined(_WIN32_) || defined(WIN32) || defined(_WIN64))
69#define NOMINMAX
70#ifdef max
71#undef max
72#undef min
73#endif
74#endif
75// Avoid conflicts with X11 headers
76#ifdef None
77#undef None
78#endif
79
80namespace nanoflann
81{
86template <typename T>
88{
89 return static_cast<T>(3.14159265358979323846);
90}
91
96template <typename T, typename = int>
97struct has_resize : std::false_type
98{
99};
100
101template <typename T>
102struct has_resize<T, decltype((void)std::declval<T>().resize(1), 0)>
103 : std::true_type
104{
105};
106
107template <typename T, typename = int>
108struct has_assign : std::false_type
109{
110};
111
112template <typename T>
113struct has_assign<T, decltype((void)std::declval<T>().assign(1, 0), 0)>
114 : std::true_type
115{
116};
117
121template <typename Container>
122inline typename std::enable_if<has_resize<Container>::value, void>::type resize(
123 Container& c, const size_t nElements)
124{
125 c.resize(nElements);
126}
127
132template <typename Container>
133inline typename std::enable_if<!has_resize<Container>::value, void>::type
134 resize(Container& c, const size_t nElements)
135{
136 if (nElements != c.size())
137 throw std::logic_error("Try to change the size of a std::array.");
138}
139
143template <typename Container, typename T>
144inline typename std::enable_if<has_assign<Container>::value, void>::type assign(
145 Container& c, const size_t nElements, const T& value)
146{
147 c.assign(nElements, value);
148}
149
153template <typename Container, typename T>
154inline typename std::enable_if<!has_assign<Container>::value, void>::type
155 assign(Container& c, const size_t nElements, const T& value)
156{
157 for (size_t i = 0; i < nElements; i++) c[i] = value;
158}
159
162{
164 template <typename PairType>
165 bool operator()(const PairType& p1, const PairType& p2) const
166 {
167 return p1.second < p2.second;
168 }
169};
170
179template <typename IndexType = size_t, typename DistanceType = double>
181{
182 ResultItem() = default;
183 ResultItem(const IndexType index, const DistanceType distance)
184 : first(index), second(distance)
185 {
186 }
187
188 IndexType first;
189 DistanceType second;
190};
191
196template <
197 typename _DistanceType, typename _IndexType = size_t,
198 typename _CountType = size_t>
200{
201 public:
202 using DistanceType = _DistanceType;
203 using IndexType = _IndexType;
204 using CountType = _CountType;
205
206 private:
207 IndexType* indices;
208 DistanceType* dists;
209 CountType capacity;
210 CountType count;
211
212 public:
213 explicit KNNResultSet(CountType capacity_)
214 : indices(nullptr), dists(nullptr), capacity(capacity_), count(0)
215 {
216 }
217
218 void init(IndexType* indices_, DistanceType* dists_)
219 {
220 indices = indices_;
221 dists = dists_;
222 count = 0;
223 }
224
225 CountType size() const { return count; }
226 bool empty() const { return count == 0; }
227 bool full() const { return count == capacity; }
228
234 bool addPoint(DistanceType dist, IndexType index)
235 {
236 CountType i;
237 for (i = count; i > 0; --i)
238 {
241#ifdef NANOFLANN_FIRST_MATCH
242 if ((dists[i - 1] > dist) ||
243 ((dist == dists[i - 1]) && (indices[i - 1] > index)))
244 {
245#else
246 if (dists[i - 1] > dist)
247 {
248#endif
249 if (i < capacity)
250 {
251 dists[i] = dists[i - 1];
252 indices[i] = indices[i - 1];
253 }
254 }
255 else
256 break;
257 }
258 if (i < capacity)
259 {
260 dists[i] = dist;
261 indices[i] = index;
262 }
263 if (count < capacity) count++;
264
265 // tell caller that the search shall continue
266 return true;
267 }
268
271 DistanceType worstDist() const
272 {
273 return (count < capacity || !count)
274 ? std::numeric_limits<DistanceType>::max()
275 : dists[count - 1];
276 }
277
278 void sort()
279 {
280 // already sorted
281 }
282};
283
285template <
286 typename _DistanceType, typename _IndexType = size_t,
287 typename _CountType = size_t>
289{
290 public:
291 using DistanceType = _DistanceType;
292 using IndexType = _IndexType;
293 using CountType = _CountType;
294
295 private:
296 IndexType* indices;
297 DistanceType* dists;
298 CountType capacity;
299 CountType count;
300 DistanceType maximumSearchDistanceSquared;
301
302 public:
303 explicit RKNNResultSet(
304 CountType capacity_, DistanceType maximumSearchDistanceSquared_)
305 : indices(nullptr),
306 dists(nullptr),
307 capacity(capacity_),
308 count(0),
309 maximumSearchDistanceSquared(maximumSearchDistanceSquared_)
310 {
311 }
312
313 void init(IndexType* indices_, DistanceType* dists_)
314 {
315 indices = indices_;
316 dists = dists_;
317 count = 0;
318 if (capacity) dists[capacity - 1] = maximumSearchDistanceSquared;
319 }
320
321 CountType size() const { return count; }
322 bool empty() const { return count == 0; }
323 bool full() const { return count == capacity; }
324
330 bool addPoint(DistanceType dist, IndexType index)
331 {
332 CountType i;
333 for (i = count; i > 0; --i)
334 {
337#ifdef NANOFLANN_FIRST_MATCH
338 if ((dists[i - 1] > dist) ||
339 ((dist == dists[i - 1]) && (indices[i - 1] > index)))
340 {
341#else
342 if (dists[i - 1] > dist)
343 {
344#endif
345 if (i < capacity)
346 {
347 dists[i] = dists[i - 1];
348 indices[i] = indices[i - 1];
349 }
350 }
351 else
352 break;
353 }
354 if (i < capacity)
355 {
356 dists[i] = dist;
357 indices[i] = index;
358 }
359 if (count < capacity) count++;
360
361 // tell caller that the search shall continue
362 return true;
363 }
364
367 DistanceType worstDist() const
368 {
369 return (count < capacity || !count) ? maximumSearchDistanceSquared
370 : dists[count - 1];
371 }
372
373 void sort()
374 {
375 // already sorted
376 }
377};
378
382template <typename _DistanceType, typename _IndexType = size_t>
384{
385 public:
386 using DistanceType = _DistanceType;
387 using IndexType = _IndexType;
388
389 public:
390 const DistanceType radius;
391
392 std::vector<ResultItem<IndexType, DistanceType>>& m_indices_dists;
393
394 explicit RadiusResultSet(
395 DistanceType radius_,
396 std::vector<ResultItem<IndexType, DistanceType>>& indices_dists)
397 : radius(radius_), m_indices_dists(indices_dists)
398 {
399 init();
400 }
401
402 void init() { clear(); }
403 void clear() { m_indices_dists.clear(); }
404
405 size_t size() const { return m_indices_dists.size(); }
406 size_t empty() const { return m_indices_dists.empty(); }
407
408 bool full() const { return true; }
409
415 bool addPoint(DistanceType dist, IndexType index)
416 {
417 if (dist < radius) m_indices_dists.emplace_back(index, dist);
418 return true;
419 }
420
421 DistanceType worstDist() const { return radius; }
422
428 {
429 if (m_indices_dists.empty())
430 throw std::runtime_error(
431 "Cannot invoke RadiusResultSet::worst_item() on "
432 "an empty list of results.");
433 auto it = std::max_element(
434 m_indices_dists.begin(), m_indices_dists.end(), IndexDist_Sorter());
435 return *it;
436 }
437
438 void sort()
439 {
440 std::sort(
441 m_indices_dists.begin(), m_indices_dists.end(), IndexDist_Sorter());
442 }
443};
444
449template <typename T>
450void save_value(std::ostream& stream, const T& value)
451{
452 stream.write(reinterpret_cast<const char*>(&value), sizeof(T));
453}
454
455template <typename T>
456void save_value(std::ostream& stream, const std::vector<T>& value)
457{
458 size_t size = value.size();
459 stream.write(reinterpret_cast<const char*>(&size), sizeof(size_t));
460 stream.write(reinterpret_cast<const char*>(value.data()), sizeof(T) * size);
461}
462
463template <typename T>
464void load_value(std::istream& stream, T& value)
465{
466 stream.read(reinterpret_cast<char*>(&value), sizeof(T));
467}
468
469template <typename T>
470void load_value(std::istream& stream, std::vector<T>& value)
471{
472 size_t size;
473 stream.read(reinterpret_cast<char*>(&size), sizeof(size_t));
474 value.resize(size);
475 stream.read(reinterpret_cast<char*>(value.data()), sizeof(T) * size);
476}
482struct Metric
483{
484};
485
496template <
497 class T, class DataSource, typename _DistanceType = T,
498 typename IndexType = uint32_t>
500{
501 using ElementType = T;
502 using DistanceType = _DistanceType;
503
504 const DataSource& data_source;
505
506 L1_Adaptor(const DataSource& _data_source) : data_source(_data_source) {}
507
508 DistanceType evalMetric(
509 const T* a, const IndexType b_idx, size_t size,
510 DistanceType worst_dist = -1) const
511 {
512 DistanceType result = DistanceType();
513 const T* last = a + size;
514 const T* lastgroup = last - 3;
515 size_t d = 0;
516
517 /* Process 4 items with each loop for efficiency. */
518 while (a < lastgroup)
519 {
520 const DistanceType diff0 =
521 std::abs(a[0] - data_source.kdtree_get_pt(b_idx, d++));
522 const DistanceType diff1 =
523 std::abs(a[1] - data_source.kdtree_get_pt(b_idx, d++));
524 const DistanceType diff2 =
525 std::abs(a[2] - data_source.kdtree_get_pt(b_idx, d++));
526 const DistanceType diff3 =
527 std::abs(a[3] - data_source.kdtree_get_pt(b_idx, d++));
528 result += diff0 + diff1 + diff2 + diff3;
529 a += 4;
530 if ((worst_dist > 0) && (result > worst_dist)) { return result; }
531 }
532 /* Process last 0-3 components. Not needed for standard vector lengths.
533 */
534 while (a < last)
535 {
536 result += std::abs(*a++ - data_source.kdtree_get_pt(b_idx, d++));
537 }
538 return result;
539 }
540
541 template <typename U, typename V>
542 DistanceType accum_dist(const U a, const V b, const size_t) const
543 {
544 return std::abs(a - b);
545 }
546};
547
558template <
559 class T, class DataSource, typename _DistanceType = T,
560 typename IndexType = uint32_t>
562{
563 using ElementType = T;
564 using DistanceType = _DistanceType;
565
566 const DataSource& data_source;
567
568 L2_Adaptor(const DataSource& _data_source) : data_source(_data_source) {}
569
570 DistanceType evalMetric(
571 const T* a, const IndexType b_idx, size_t size,
572 DistanceType worst_dist = -1) const
573 {
574 DistanceType result = DistanceType();
575 const T* last = a + size;
576 const T* lastgroup = last - 3;
577 size_t d = 0;
578
579 /* Process 4 items with each loop for efficiency. */
580 while (a < lastgroup)
581 {
582 const DistanceType diff0 =
583 a[0] - data_source.kdtree_get_pt(b_idx, d++);
584 const DistanceType diff1 =
585 a[1] - data_source.kdtree_get_pt(b_idx, d++);
586 const DistanceType diff2 =
587 a[2] - data_source.kdtree_get_pt(b_idx, d++);
588 const DistanceType diff3 =
589 a[3] - data_source.kdtree_get_pt(b_idx, d++);
590 result +=
591 diff0 * diff0 + diff1 * diff1 + diff2 * diff2 + diff3 * diff3;
592 a += 4;
593 if ((worst_dist > 0) && (result > worst_dist)) { return result; }
594 }
595 /* Process last 0-3 components. Not needed for standard vector lengths.
596 */
597 while (a < last)
598 {
599 const DistanceType diff0 =
600 *a++ - data_source.kdtree_get_pt(b_idx, d++);
601 result += diff0 * diff0;
602 }
603 return result;
604 }
605
606 template <typename U, typename V>
607 DistanceType accum_dist(const U a, const V b, const size_t) const
608 {
609 return (a - b) * (a - b);
610 }
611};
612
623template <
624 class T, class DataSource, typename _DistanceType = T,
625 typename IndexType = uint32_t>
627{
628 using ElementType = T;
629 using DistanceType = _DistanceType;
630
631 const DataSource& data_source;
632
633 L2_Simple_Adaptor(const DataSource& _data_source)
634 : data_source(_data_source)
635 {
636 }
637
638 DistanceType evalMetric(
639 const T* a, const IndexType b_idx, size_t size) const
640 {
641 DistanceType result = DistanceType();
642 for (size_t i = 0; i < size; ++i)
643 {
644 const DistanceType diff =
645 a[i] - data_source.kdtree_get_pt(b_idx, i);
646 result += diff * diff;
647 }
648 return result;
649 }
650
651 template <typename U, typename V>
652 DistanceType accum_dist(const U a, const V b, const size_t) const
653 {
654 return (a - b) * (a - b);
655 }
656};
657
668template <
669 class T, class DataSource, typename _DistanceType = T,
670 typename IndexType = uint32_t>
672{
673 using ElementType = T;
674 using DistanceType = _DistanceType;
675
676 const DataSource& data_source;
677
678 SO2_Adaptor(const DataSource& _data_source) : data_source(_data_source) {}
679
680 DistanceType evalMetric(
681 const T* a, const IndexType b_idx, size_t size) const
682 {
683 return accum_dist(
684 a[size - 1], data_source.kdtree_get_pt(b_idx, size - 1), size - 1);
685 }
686
689 template <typename U, typename V>
690 DistanceType accum_dist(const U a, const V b, const size_t) const
691 {
692 DistanceType result = DistanceType();
693 DistanceType PI = pi_const<DistanceType>();
694 result = b - a;
695 if (result > PI)
696 result -= 2 * PI;
697 else if (result < -PI)
698 result += 2 * PI;
699 return result;
700 }
701};
702
713template <
714 class T, class DataSource, typename _DistanceType = T,
715 typename IndexType = uint32_t>
717{
718 using ElementType = T;
719 using DistanceType = _DistanceType;
720
722 distance_L2_Simple;
723
724 SO3_Adaptor(const DataSource& _data_source)
725 : distance_L2_Simple(_data_source)
726 {
727 }
728
729 DistanceType evalMetric(
730 const T* a, const IndexType b_idx, size_t size) const
731 {
732 return distance_L2_Simple.evalMetric(a, b_idx, size);
733 }
734
735 template <typename U, typename V>
736 DistanceType accum_dist(const U a, const V b, const size_t idx) const
737 {
738 return distance_L2_Simple.accum_dist(a, b, idx);
739 }
740};
741
743struct metric_L1 : public Metric
744{
745 template <class T, class DataSource, typename IndexType = uint32_t>
750};
753struct metric_L2 : public Metric
754{
755 template <class T, class DataSource, typename IndexType = uint32_t>
760};
764{
765 template <class T, class DataSource, typename IndexType = uint32_t>
770};
772struct metric_SO2 : public Metric
773{
774 template <class T, class DataSource, typename IndexType = uint32_t>
779};
781struct metric_SO3 : public Metric
782{
783 template <class T, class DataSource, typename IndexType = uint32_t>
788};
789
795enum class KDTreeSingleIndexAdaptorFlags
796{
797 None = 0,
798 SkipInitialBuildIndex = 1
799};
800
801inline std::underlying_type<KDTreeSingleIndexAdaptorFlags>::type operator&(
802 KDTreeSingleIndexAdaptorFlags lhs, KDTreeSingleIndexAdaptorFlags rhs)
803{
804 using underlying =
805 typename std::underlying_type<KDTreeSingleIndexAdaptorFlags>::type;
806 return static_cast<underlying>(lhs) & static_cast<underlying>(rhs);
807}
808
811{
813 size_t _leaf_max_size = 10,
814 KDTreeSingleIndexAdaptorFlags _flags =
815 KDTreeSingleIndexAdaptorFlags::None,
816 unsigned int _n_thread_build = 1)
817 : leaf_max_size(_leaf_max_size),
818 flags(_flags),
819 n_thread_build(_n_thread_build)
820 {
821 }
822
823 size_t leaf_max_size;
824 KDTreeSingleIndexAdaptorFlags flags;
825 unsigned int n_thread_build;
826};
827
830{
831 SearchParameters(float eps_ = 0, bool sorted_ = true)
832 : eps(eps_), sorted(sorted_)
833 {
834 }
835
836 float eps;
837 bool sorted;
839};
860{
861 static constexpr size_t WORDSIZE = 16; // WORDSIZE must >= 8
862 static constexpr size_t BLOCKSIZE = 8192;
863
864 /* We maintain memory alignment to word boundaries by requiring that all
865 allocations be in multiples of the machine wordsize. */
866 /* Size of machine word in bytes. Must be power of 2. */
867 /* Minimum number of bytes requested at a time from the system. Must be
868 * multiple of WORDSIZE. */
869
870 using Size = size_t;
871
872 Size remaining_ = 0;
873 void* base_ = nullptr;
874 void* loc_ = nullptr;
875
876 void internal_init()
877 {
878 remaining_ = 0;
879 base_ = nullptr;
880 usedMemory = 0;
881 wastedMemory = 0;
882 }
883
884 public:
885 Size usedMemory = 0;
886 Size wastedMemory = 0;
887
891 PooledAllocator() { internal_init(); }
892
896 ~PooledAllocator() { free_all(); }
897
899 void free_all()
900 {
901 while (base_ != nullptr)
902 {
903 // Get pointer to prev block
904 void* prev = *(static_cast<void**>(base_));
905 ::free(base_);
906 base_ = prev;
907 }
908 internal_init();
909 }
910
915 void* malloc(const size_t req_size)
916 {
917 /* Round size up to a multiple of wordsize. The following expression
918 only works for WORDSIZE that is a power of 2, by masking last bits
919 of incremented size to zero.
920 */
921 const Size size = (req_size + (WORDSIZE - 1)) & ~(WORDSIZE - 1);
922
923 /* Check whether a new block must be allocated. Note that the first
924 word of a block is reserved for a pointer to the previous block.
925 */
926 if (size > remaining_)
927 {
928 wastedMemory += remaining_;
929
930 /* Allocate new storage. */
931 const Size blocksize =
932 size > BLOCKSIZE ? size + WORDSIZE : BLOCKSIZE + WORDSIZE;
933
934 // use the standard C malloc to allocate memory
935 void* m = ::malloc(blocksize);
936 if (!m)
937 {
938 fprintf(stderr, "Failed to allocate memory.\n");
939 throw std::bad_alloc();
940 }
941
942 /* Fill first word of new block with pointer to previous block. */
943 static_cast<void**>(m)[0] = base_;
944 base_ = m;
945
946 remaining_ = blocksize - WORDSIZE;
947 loc_ = static_cast<char*>(m) + WORDSIZE;
948 }
949 void* rloc = loc_;
950 loc_ = static_cast<char*>(loc_) + size;
951 remaining_ -= size;
952
953 usedMemory += size;
954
955 return rloc;
956 }
957
965 template <typename T>
966 T* allocate(const size_t count = 1)
967 {
968 T* mem = static_cast<T*>(this->malloc(sizeof(T) * count));
969 return mem;
970 }
971};
980template <int32_t DIM, typename T>
982{
983 using type = std::array<T, DIM>;
984};
986template <typename T>
987struct array_or_vector<-1, T>
988{
989 using type = std::vector<T>;
990};
991
1008template <
1009 class Derived, typename Distance, class DatasetAdaptor, int32_t DIM = -1,
1010 typename index_t = uint32_t>
1012{
1013 public:
1016 void freeIndex(Derived& obj)
1017 {
1018 obj.pool_.free_all();
1019 obj.root_node_ = nullptr;
1020 obj.size_at_index_build_ = 0;
1021 }
1022
1023 using ElementType = typename Distance::ElementType;
1024 using DistanceType = typename Distance::DistanceType;
1025 using IndexType = index_t;
1026
1030 std::vector<IndexType> vAcc_;
1031
1032 using Offset = typename decltype(vAcc_)::size_type;
1033 using Size = typename decltype(vAcc_)::size_type;
1034 using Dimension = int32_t;
1035
1036 /*---------------------------
1037 * Internal Data Structures
1038 * --------------------------*/
1039 struct Node
1040 {
1043 union
1044 {
1045 struct leaf
1046 {
1047 Offset left, right;
1048 } lr;
1049 struct nonleaf
1050 {
1051 Dimension divfeat;
1053 DistanceType divlow, divhigh;
1054 } sub;
1055 } node_type;
1056
1058 Node *child1 = nullptr, *child2 = nullptr;
1059 };
1060
1061 using NodePtr = Node*;
1062 using NodeConstPtr = const Node*;
1063
1065 {
1066 ElementType low, high;
1067 };
1068
1069 NodePtr root_node_ = nullptr;
1070
1071 Size leaf_max_size_ = 0;
1072
1074 Size n_thread_build_ = 1;
1076 Size size_ = 0;
1078 Size size_at_index_build_ = 0;
1079 Dimension dim_ = 0;
1080
1083 using BoundingBox = typename array_or_vector<DIM, Interval>::type;
1084
1087 using distance_vector_t = typename array_or_vector<DIM, DistanceType>::type;
1088
1091
1100
1102 Size size(const Derived& obj) const { return obj.size_; }
1103
1105 Size veclen(const Derived& obj) { return DIM > 0 ? DIM : obj.dim; }
1106
1108 ElementType dataset_get(
1109 const Derived& obj, IndexType element, Dimension component) const
1110 {
1111 return obj.dataset_.kdtree_get_pt(element, component);
1112 }
1113
1118 Size usedMemory(Derived& obj)
1119 {
1120 return obj.pool_.usedMemory + obj.pool_.wastedMemory +
1121 obj.dataset_.kdtree_get_point_count() *
1122 sizeof(IndexType); // pool memory and vind array memory
1123 }
1124
1125 void computeMinMax(
1126 const Derived& obj, Offset ind, Size count, Dimension element,
1127 ElementType& min_elem, ElementType& max_elem)
1128 {
1129 min_elem = dataset_get(obj, vAcc_[ind], element);
1130 max_elem = min_elem;
1131 for (Offset i = 1; i < count; ++i)
1132 {
1133 ElementType val = dataset_get(obj, vAcc_[ind + i], element);
1134 if (val < min_elem) min_elem = val;
1135 if (val > max_elem) max_elem = val;
1136 }
1137 }
1138
1147 Derived& obj, const Offset left, const Offset right, BoundingBox& bbox)
1148 {
1149 assert(left < obj.dataset_.kdtree_get_point_count());
1150
1151 NodePtr node = obj.pool_.template allocate<Node>(); // allocate memory
1152 const auto dims = (DIM > 0 ? DIM : obj.dim_);
1153
1154 /* If too few exemplars remain, then make this a leaf node. */
1155 if ((right - left) <= static_cast<Offset>(obj.leaf_max_size_))
1156 {
1157 node->child1 = node->child2 = nullptr; /* Mark as leaf node. */
1158 node->node_type.lr.left = left;
1159 node->node_type.lr.right = right;
1160
1161 // compute bounding-box of leaf points
1162 for (Dimension i = 0; i < dims; ++i)
1163 {
1164 bbox[i].low = dataset_get(obj, obj.vAcc_[left], i);
1165 bbox[i].high = dataset_get(obj, obj.vAcc_[left], i);
1166 }
1167 for (Offset k = left + 1; k < right; ++k)
1168 {
1169 for (Dimension i = 0; i < dims; ++i)
1170 {
1171 const auto val = dataset_get(obj, obj.vAcc_[k], i);
1172 if (bbox[i].low > val) bbox[i].low = val;
1173 if (bbox[i].high < val) bbox[i].high = val;
1174 }
1175 }
1176 }
1177 else
1178 {
1179 Offset idx;
1180 Dimension cutfeat;
1181 DistanceType cutval;
1182 middleSplit_(obj, left, right - left, idx, cutfeat, cutval, bbox);
1183
1184 node->node_type.sub.divfeat = cutfeat;
1185
1186 BoundingBox left_bbox(bbox);
1187 left_bbox[cutfeat].high = cutval;
1188 node->child1 = this->divideTree(obj, left, left + idx, left_bbox);
1189
1190 BoundingBox right_bbox(bbox);
1191 right_bbox[cutfeat].low = cutval;
1192 node->child2 = this->divideTree(obj, left + idx, right, right_bbox);
1193
1194 node->node_type.sub.divlow = left_bbox[cutfeat].high;
1195 node->node_type.sub.divhigh = right_bbox[cutfeat].low;
1196
1197 for (Dimension i = 0; i < dims; ++i)
1198 {
1199 bbox[i].low = std::min(left_bbox[i].low, right_bbox[i].low);
1200 bbox[i].high = std::max(left_bbox[i].high, right_bbox[i].high);
1201 }
1202 }
1203
1204 return node;
1205 }
1206
1218 Derived& obj, const Offset left, const Offset right, BoundingBox& bbox,
1219 std::atomic<unsigned int>& thread_count, std::mutex& mutex)
1220 {
1221 std::unique_lock<std::mutex> lock(mutex);
1222 NodePtr node = obj.pool_.template allocate<Node>(); // allocate memory
1223 lock.unlock();
1224
1225 const auto dims = (DIM > 0 ? DIM : obj.dim_);
1226
1227 /* If too few exemplars remain, then make this a leaf node. */
1228 if ((right - left) <= static_cast<Offset>(obj.leaf_max_size_))
1229 {
1230 node->child1 = node->child2 = nullptr; /* Mark as leaf node. */
1231 node->node_type.lr.left = left;
1232 node->node_type.lr.right = right;
1233
1234 // compute bounding-box of leaf points
1235 for (Dimension i = 0; i < dims; ++i)
1236 {
1237 bbox[i].low = dataset_get(obj, obj.vAcc_[left], i);
1238 bbox[i].high = dataset_get(obj, obj.vAcc_[left], i);
1239 }
1240 for (Offset k = left + 1; k < right; ++k)
1241 {
1242 for (Dimension i = 0; i < dims; ++i)
1243 {
1244 const auto val = dataset_get(obj, obj.vAcc_[k], i);
1245 if (bbox[i].low > val) bbox[i].low = val;
1246 if (bbox[i].high < val) bbox[i].high = val;
1247 }
1248 }
1249 }
1250 else
1251 {
1252 Offset idx;
1253 Dimension cutfeat;
1254 DistanceType cutval;
1255 middleSplit_(obj, left, right - left, idx, cutfeat, cutval, bbox);
1256
1257 node->node_type.sub.divfeat = cutfeat;
1258
1259 std::future<NodePtr> right_future;
1260
1261 BoundingBox right_bbox(bbox);
1262 right_bbox[cutfeat].low = cutval;
1263 if (++thread_count < n_thread_build_)
1264 {
1265 // Concurrent right sub-tree
1266 right_future = std::async(
1267 std::launch::async, &KDTreeBaseClass::divideTreeConcurrent,
1268 this, std::ref(obj), left + idx, right,
1269 std::ref(right_bbox), std::ref(thread_count),
1270 std::ref(mutex));
1271 }
1272 else { --thread_count; }
1273
1274 BoundingBox left_bbox(bbox);
1275 left_bbox[cutfeat].high = cutval;
1276 node->child1 = this->divideTreeConcurrent(
1277 obj, left, left + idx, left_bbox, thread_count, mutex);
1278
1279 if (right_future.valid())
1280 {
1281 // Block and wait for concurrent right sub-tree
1282 node->child2 = right_future.get();
1283 --thread_count;
1284 }
1285 else
1286 {
1287 node->child2 = this->divideTreeConcurrent(
1288 obj, left + idx, right, right_bbox, thread_count, mutex);
1289 }
1290
1291 node->node_type.sub.divlow = left_bbox[cutfeat].high;
1292 node->node_type.sub.divhigh = right_bbox[cutfeat].low;
1293
1294 for (Dimension i = 0; i < dims; ++i)
1295 {
1296 bbox[i].low = std::min(left_bbox[i].low, right_bbox[i].low);
1297 bbox[i].high = std::max(left_bbox[i].high, right_bbox[i].high);
1298 }
1299 }
1300
1301 return node;
1302 }
1303
1304 void middleSplit_(
1305 const Derived& obj, const Offset ind, const Size count, Offset& index,
1306 Dimension& cutfeat, DistanceType& cutval, const BoundingBox& bbox)
1307 {
1308 const auto dims = (DIM > 0 ? DIM : obj.dim_);
1309 const auto EPS = static_cast<DistanceType>(0.00001);
1310 ElementType max_span = bbox[0].high - bbox[0].low;
1311 for (Dimension i = 1; i < dims; ++i)
1312 {
1313 ElementType span = bbox[i].high - bbox[i].low;
1314 if (span > max_span) { max_span = span; }
1315 }
1316 ElementType max_spread = -1;
1317 cutfeat = 0;
1318 ElementType min_elem = 0, max_elem = 0;
1319 for (Dimension i = 0; i < dims; ++i)
1320 {
1321 ElementType span = bbox[i].high - bbox[i].low;
1322 if (span >= (1 - EPS) * max_span)
1323 {
1324 ElementType min_elem_, max_elem_;
1325 computeMinMax(obj, ind, count, i, min_elem_, max_elem_);
1326 ElementType spread = max_elem_ - min_elem_;
1327 if (spread > max_spread)
1328 {
1329 cutfeat = i;
1330 max_spread = spread;
1331 min_elem = min_elem_;
1332 max_elem = max_elem_;
1333 }
1334 }
1335 }
1336 // split in the middle
1337 DistanceType split_val = (bbox[cutfeat].low + bbox[cutfeat].high) / 2;
1338
1339 if (split_val < min_elem)
1340 cutval = min_elem;
1341 else if (split_val > max_elem)
1342 cutval = max_elem;
1343 else
1344 cutval = split_val;
1345
1346 Offset lim1, lim2;
1347 planeSplit(obj, ind, count, cutfeat, cutval, lim1, lim2);
1348
1349 if (lim1 > count / 2)
1350 index = lim1;
1351 else if (lim2 < count / 2)
1352 index = lim2;
1353 else
1354 index = count / 2;
1355 }
1356
1367 const Derived& obj, const Offset ind, const Size count,
1368 const Dimension cutfeat, const DistanceType& cutval, Offset& lim1,
1369 Offset& lim2)
1370 {
1371 /* Move vector indices for left subtree to front of list. */
1372 Offset left = 0;
1373 Offset right = count - 1;
1374 for (;;)
1375 {
1376 while (left <= right &&
1377 dataset_get(obj, vAcc_[ind + left], cutfeat) < cutval)
1378 ++left;
1379 while (right && left <= right &&
1380 dataset_get(obj, vAcc_[ind + right], cutfeat) >= cutval)
1381 --right;
1382 if (left > right || !right)
1383 break; // "!right" was added to support unsigned Index types
1384 std::swap(vAcc_[ind + left], vAcc_[ind + right]);
1385 ++left;
1386 --right;
1387 }
1388 /* If either list is empty, it means that all remaining features
1389 * are identical. Split in the middle to maintain a balanced tree.
1390 */
1391 lim1 = left;
1392 right = count - 1;
1393 for (;;)
1394 {
1395 while (left <= right &&
1396 dataset_get(obj, vAcc_[ind + left], cutfeat) <= cutval)
1397 ++left;
1398 while (right && left <= right &&
1399 dataset_get(obj, vAcc_[ind + right], cutfeat) > cutval)
1400 --right;
1401 if (left > right || !right)
1402 break; // "!right" was added to support unsigned Index types
1403 std::swap(vAcc_[ind + left], vAcc_[ind + right]);
1404 ++left;
1405 --right;
1406 }
1407 lim2 = left;
1408 }
1409
1410 DistanceType computeInitialDistances(
1411 const Derived& obj, const ElementType* vec,
1412 distance_vector_t& dists) const
1413 {
1414 assert(vec);
1415 DistanceType dist = DistanceType();
1416
1417 for (Dimension i = 0; i < (DIM > 0 ? DIM : obj.dim_); ++i)
1418 {
1419 if (vec[i] < obj.root_bbox_[i].low)
1420 {
1421 dists[i] =
1422 obj.distance_.accum_dist(vec[i], obj.root_bbox_[i].low, i);
1423 dist += dists[i];
1424 }
1425 if (vec[i] > obj.root_bbox_[i].high)
1426 {
1427 dists[i] =
1428 obj.distance_.accum_dist(vec[i], obj.root_bbox_[i].high, i);
1429 dist += dists[i];
1430 }
1431 }
1432 return dist;
1433 }
1434
1435 static void save_tree(
1436 const Derived& obj, std::ostream& stream, const NodeConstPtr tree)
1437 {
1438 save_value(stream, *tree);
1439 if (tree->child1 != nullptr) { save_tree(obj, stream, tree->child1); }
1440 if (tree->child2 != nullptr) { save_tree(obj, stream, tree->child2); }
1441 }
1442
1443 static void load_tree(Derived& obj, std::istream& stream, NodePtr& tree)
1444 {
1445 tree = obj.pool_.template allocate<Node>();
1446 load_value(stream, *tree);
1447 if (tree->child1 != nullptr) { load_tree(obj, stream, tree->child1); }
1448 if (tree->child2 != nullptr) { load_tree(obj, stream, tree->child2); }
1449 }
1450
1456 void saveIndex(const Derived& obj, std::ostream& stream) const
1457 {
1458 save_value(stream, obj.size_);
1459 save_value(stream, obj.dim_);
1460 save_value(stream, obj.root_bbox_);
1461 save_value(stream, obj.leaf_max_size_);
1462 save_value(stream, obj.vAcc_);
1463 if (obj.root_node_) save_tree(obj, stream, obj.root_node_);
1464 }
1465
1471 void loadIndex(Derived& obj, std::istream& stream)
1472 {
1473 load_value(stream, obj.size_);
1474 load_value(stream, obj.dim_);
1475 load_value(stream, obj.root_bbox_);
1476 load_value(stream, obj.leaf_max_size_);
1477 load_value(stream, obj.vAcc_);
1478 load_tree(obj, stream, obj.root_node_);
1479 }
1480};
1481
1523template <
1524 typename Distance, class DatasetAdaptor, int32_t DIM = -1,
1525 typename index_t = uint32_t>
1527 : public KDTreeBaseClass<
1528 KDTreeSingleIndexAdaptor<Distance, DatasetAdaptor, DIM, index_t>,
1529 Distance, DatasetAdaptor, DIM, index_t>
1530{
1531 public:
1535 Distance, DatasetAdaptor, DIM, index_t>&) = delete;
1536
1538 const DatasetAdaptor& dataset_;
1539
1540 const KDTreeSingleIndexAdaptorParams indexParams;
1541
1542 Distance distance_;
1543
1544 using Base = typename nanoflann::KDTreeBaseClass<
1546 Distance, DatasetAdaptor, DIM, index_t>,
1547 Distance, DatasetAdaptor, DIM, index_t>;
1548
1549 using Offset = typename Base::Offset;
1550 using Size = typename Base::Size;
1551 using Dimension = typename Base::Dimension;
1552
1553 using ElementType = typename Base::ElementType;
1554 using DistanceType = typename Base::DistanceType;
1555 using IndexType = typename Base::IndexType;
1556
1557 using Node = typename Base::Node;
1558 using NodePtr = Node*;
1559
1560 using Interval = typename Base::Interval;
1561
1564 using BoundingBox = typename Base::BoundingBox;
1565
1568 using distance_vector_t = typename Base::distance_vector_t;
1569
1590 template <class... Args>
1592 const Dimension dimensionality, const DatasetAdaptor& inputData,
1593 const KDTreeSingleIndexAdaptorParams& params, Args&&... args)
1594 : dataset_(inputData),
1595 indexParams(params),
1596 distance_(inputData, std::forward<Args>(args)...)
1597 {
1598 init(dimensionality, params);
1599 }
1600
1601 explicit KDTreeSingleIndexAdaptor(
1602 const Dimension dimensionality, const DatasetAdaptor& inputData,
1603 const KDTreeSingleIndexAdaptorParams& params = {})
1604 : dataset_(inputData), indexParams(params), distance_(inputData)
1605 {
1606 init(dimensionality, params);
1607 }
1608
1609 private:
1610 void init(
1611 const Dimension dimensionality,
1612 const KDTreeSingleIndexAdaptorParams& params)
1613 {
1614 Base::size_ = dataset_.kdtree_get_point_count();
1615 Base::size_at_index_build_ = Base::size_;
1616 Base::dim_ = dimensionality;
1617 if (DIM > 0) Base::dim_ = DIM;
1618 Base::leaf_max_size_ = params.leaf_max_size;
1619 if (params.n_thread_build > 0)
1620 {
1621 Base::n_thread_build_ = params.n_thread_build;
1622 }
1623 else
1624 {
1625 Base::n_thread_build_ =
1626 std::max(std::thread::hardware_concurrency(), 1u);
1627 }
1628
1629 if (!(params.flags &
1630 KDTreeSingleIndexAdaptorFlags::SkipInitialBuildIndex))
1631 {
1632 // Build KD-tree:
1633 buildIndex();
1634 }
1635 }
1636
1637 public:
1642 {
1643 Base::size_ = dataset_.kdtree_get_point_count();
1644 Base::size_at_index_build_ = Base::size_;
1645 init_vind();
1646 this->freeIndex(*this);
1647 Base::size_at_index_build_ = Base::size_;
1648 if (Base::size_ == 0) return;
1649 computeBoundingBox(Base::root_bbox_);
1650 // construct the tree
1651 if (Base::n_thread_build_ == 1)
1652 {
1653 Base::root_node_ =
1654 this->divideTree(*this, 0, Base::size_, Base::root_bbox_);
1655 }
1656 else
1657 {
1658#ifndef NANOFLANN_NO_THREADS
1659 std::atomic<unsigned int> thread_count(0u);
1660 std::mutex mutex;
1661 Base::root_node_ = this->divideTreeConcurrent(
1662 *this, 0, Base::size_, Base::root_bbox_, thread_count, mutex);
1663#else /* NANOFLANN_NO_THREADS */
1664 throw std::runtime_error("Multithreading is disabled");
1665#endif /* NANOFLANN_NO_THREADS */
1666 }
1667 }
1668
1688 template <typename RESULTSET>
1690 RESULTSET& result, const ElementType* vec,
1691 const SearchParameters& searchParams = {}) const
1692 {
1693 assert(vec);
1694 if (this->size(*this) == 0) return false;
1695 if (!Base::root_node_)
1696 throw std::runtime_error(
1697 "[nanoflann] findNeighbors() called before building the "
1698 "index.");
1699 float epsError = 1 + searchParams.eps;
1700
1701 // fixed or variable-sized container (depending on DIM)
1702 distance_vector_t dists;
1703 // Fill it with zeros.
1704 auto zero = static_cast<typename RESULTSET::DistanceType>(0);
1705 assign(dists, (DIM > 0 ? DIM : Base::dim_), zero);
1706 DistanceType dist = this->computeInitialDistances(*this, vec, dists);
1707 searchLevel(result, vec, Base::root_node_, dist, dists, epsError);
1708
1709 if (searchParams.sorted) result.sort();
1710
1711 return result.full();
1712 }
1713
1730 const ElementType* query_point, const Size num_closest,
1731 IndexType* out_indices, DistanceType* out_distances) const
1732 {
1734 resultSet.init(out_indices, out_distances);
1735 findNeighbors(resultSet, query_point);
1736 return resultSet.size();
1737 }
1738
1759 const ElementType* query_point, const DistanceType& radius,
1760 std::vector<ResultItem<IndexType, DistanceType>>& IndicesDists,
1761 const SearchParameters& searchParams = {}) const
1762 {
1764 radius, IndicesDists);
1765 const Size nFound =
1766 radiusSearchCustomCallback(query_point, resultSet, searchParams);
1767 return nFound;
1768 }
1769
1775 template <class SEARCH_CALLBACK>
1777 const ElementType* query_point, SEARCH_CALLBACK& resultSet,
1778 const SearchParameters& searchParams = {}) const
1779 {
1780 findNeighbors(resultSet, query_point, searchParams);
1781 return resultSet.size();
1782 }
1783
1801 const ElementType* query_point, const Size num_closest,
1802 IndexType* out_indices, DistanceType* out_distances,
1803 const DistanceType& radius) const
1804 {
1806 num_closest, radius);
1807 resultSet.init(out_indices, out_distances);
1808 findNeighbors(resultSet, query_point);
1809 return resultSet.size();
1810 }
1811
1814 public:
1818 {
1819 // Create a permutable array of indices to the input vectors.
1820 Base::size_ = dataset_.kdtree_get_point_count();
1821 if (Base::vAcc_.size() != Base::size_) Base::vAcc_.resize(Base::size_);
1822 for (IndexType i = 0; i < static_cast<IndexType>(Base::size_); i++)
1823 Base::vAcc_[i] = i;
1824 }
1825
1826 void computeBoundingBox(BoundingBox& bbox)
1827 {
1828 const auto dims = (DIM > 0 ? DIM : Base::dim_);
1829 resize(bbox, dims);
1830 if (dataset_.kdtree_get_bbox(bbox))
1831 {
1832 // Done! It was implemented in derived class
1833 }
1834 else
1835 {
1836 const Size N = dataset_.kdtree_get_point_count();
1837 if (!N)
1838 throw std::runtime_error(
1839 "[nanoflann] computeBoundingBox() called but "
1840 "no data points found.");
1841 for (Dimension i = 0; i < dims; ++i)
1842 {
1843 bbox[i].low = bbox[i].high =
1844 this->dataset_get(*this, Base::vAcc_[0], i);
1845 }
1846 for (Offset k = 1; k < N; ++k)
1847 {
1848 for (Dimension i = 0; i < dims; ++i)
1849 {
1850 const auto val =
1851 this->dataset_get(*this, Base::vAcc_[k], i);
1852 if (val < bbox[i].low) bbox[i].low = val;
1853 if (val > bbox[i].high) bbox[i].high = val;
1854 }
1855 }
1856 }
1857 }
1858
1865 template <class RESULTSET>
1867 RESULTSET& result_set, const ElementType* vec, const NodePtr node,
1868 DistanceType mindist, distance_vector_t& dists,
1869 const float epsError) const
1870 {
1871 /* If this is a leaf node, then do check and return. */
1872 if ((node->child1 == nullptr) && (node->child2 == nullptr))
1873 {
1874 DistanceType worst_dist = result_set.worstDist();
1875 for (Offset i = node->node_type.lr.left;
1876 i < node->node_type.lr.right; ++i)
1877 {
1878 const IndexType accessor = Base::vAcc_[i]; // reorder... : i;
1879 DistanceType dist = distance_.evalMetric(
1880 vec, accessor, (DIM > 0 ? DIM : Base::dim_));
1881 if (dist < worst_dist)
1882 {
1883 if (!result_set.addPoint(dist, Base::vAcc_[i]))
1884 {
1885 // the resultset doesn't want to receive any more
1886 // points, we're done searching!
1887 return false;
1888 }
1889 }
1890 }
1891 return true;
1892 }
1893
1894 /* Which child branch should be taken first? */
1895 Dimension idx = node->node_type.sub.divfeat;
1896 ElementType val = vec[idx];
1897 DistanceType diff1 = val - node->node_type.sub.divlow;
1898 DistanceType diff2 = val - node->node_type.sub.divhigh;
1899
1900 NodePtr bestChild;
1901 NodePtr otherChild;
1902 DistanceType cut_dist;
1903 if ((diff1 + diff2) < 0)
1904 {
1905 bestChild = node->child1;
1906 otherChild = node->child2;
1907 cut_dist =
1908 distance_.accum_dist(val, node->node_type.sub.divhigh, idx);
1909 }
1910 else
1911 {
1912 bestChild = node->child2;
1913 otherChild = node->child1;
1914 cut_dist =
1915 distance_.accum_dist(val, node->node_type.sub.divlow, idx);
1916 }
1917
1918 /* Call recursively to search next level down. */
1919 if (!searchLevel(result_set, vec, bestChild, mindist, dists, epsError))
1920 {
1921 // the resultset doesn't want to receive any more points, we're done
1922 // searching!
1923 return false;
1924 }
1925
1926 DistanceType dst = dists[idx];
1927 mindist = mindist + cut_dist - dst;
1928 dists[idx] = cut_dist;
1929 if (mindist * epsError <= result_set.worstDist())
1930 {
1931 if (!searchLevel(
1932 result_set, vec, otherChild, mindist, dists, epsError))
1933 {
1934 // the resultset doesn't want to receive any more points, we're
1935 // done searching!
1936 return false;
1937 }
1938 }
1939 dists[idx] = dst;
1940 return true;
1941 }
1942
1943 public:
1949 void saveIndex(std::ostream& stream) const
1950 {
1951 Base::saveIndex(*this, stream);
1952 }
1953
1959 void loadIndex(std::istream& stream) { Base::loadIndex(*this, stream); }
1960
1961}; // class KDTree
1962
2000template <
2001 typename Distance, class DatasetAdaptor, int32_t DIM = -1,
2002 typename IndexType = uint32_t>
2004 : public KDTreeBaseClass<
2005 KDTreeSingleIndexDynamicAdaptor_<
2006 Distance, DatasetAdaptor, DIM, IndexType>,
2007 Distance, DatasetAdaptor, DIM, IndexType>
2008{
2009 public:
2013 const DatasetAdaptor& dataset_;
2014
2015 KDTreeSingleIndexAdaptorParams index_params_;
2016
2017 std::vector<int>& treeIndex_;
2018
2019 Distance distance_;
2020
2021 using Base = typename nanoflann::KDTreeBaseClass<
2023 Distance, DatasetAdaptor, DIM, IndexType>,
2024 Distance, DatasetAdaptor, DIM, IndexType>;
2025
2026 using ElementType = typename Base::ElementType;
2027 using DistanceType = typename Base::DistanceType;
2028
2029 using Offset = typename Base::Offset;
2030 using Size = typename Base::Size;
2031 using Dimension = typename Base::Dimension;
2032
2033 using Node = typename Base::Node;
2034 using NodePtr = Node*;
2035
2036 using Interval = typename Base::Interval;
2039 using BoundingBox = typename Base::BoundingBox;
2040
2043 using distance_vector_t = typename Base::distance_vector_t;
2044
2061 const Dimension dimensionality, const DatasetAdaptor& inputData,
2062 std::vector<int>& treeIndex,
2063 const KDTreeSingleIndexAdaptorParams& params =
2065 : dataset_(inputData),
2066 index_params_(params),
2067 treeIndex_(treeIndex),
2068 distance_(inputData)
2069 {
2070 Base::size_ = 0;
2071 Base::size_at_index_build_ = 0;
2072 for (auto& v : Base::root_bbox_) v = {};
2073 Base::dim_ = dimensionality;
2074 if (DIM > 0) Base::dim_ = DIM;
2075 Base::leaf_max_size_ = params.leaf_max_size;
2076 if (params.n_thread_build > 0)
2077 {
2078 Base::n_thread_build_ = params.n_thread_build;
2079 }
2080 else
2081 {
2082 Base::n_thread_build_ =
2083 std::max(std::thread::hardware_concurrency(), 1u);
2084 }
2085 }
2086
2089 const KDTreeSingleIndexDynamicAdaptor_& rhs) = default;
2090
2094 {
2096 std::swap(Base::vAcc_, tmp.Base::vAcc_);
2097 std::swap(Base::leaf_max_size_, tmp.Base::leaf_max_size_);
2098 std::swap(index_params_, tmp.index_params_);
2099 std::swap(treeIndex_, tmp.treeIndex_);
2100 std::swap(Base::size_, tmp.Base::size_);
2101 std::swap(Base::size_at_index_build_, tmp.Base::size_at_index_build_);
2102 std::swap(Base::root_node_, tmp.Base::root_node_);
2103 std::swap(Base::root_bbox_, tmp.Base::root_bbox_);
2104 std::swap(Base::pool_, tmp.Base::pool_);
2105 return *this;
2106 }
2107
2112 {
2113 Base::size_ = Base::vAcc_.size();
2114 this->freeIndex(*this);
2115 Base::size_at_index_build_ = Base::size_;
2116 if (Base::size_ == 0) return;
2117 computeBoundingBox(Base::root_bbox_);
2118 // construct the tree
2119 if (Base::n_thread_build_ == 1)
2120 {
2121 Base::root_node_ =
2122 this->divideTree(*this, 0, Base::size_, Base::root_bbox_);
2123 }
2124 else
2125 {
2126#ifndef NANOFLANN_NO_THREADS
2127 std::atomic<unsigned int> thread_count(0u);
2128 std::mutex mutex;
2129 Base::root_node_ = this->divideTreeConcurrent(
2130 *this, 0, Base::size_, Base::root_bbox_, thread_count, mutex);
2131#else /* NANOFLANN_NO_THREADS */
2132 throw std::runtime_error("Multithreading is disabled");
2133#endif /* NANOFLANN_NO_THREADS */
2134 }
2135 }
2136
2160 template <typename RESULTSET>
2162 RESULTSET& result, const ElementType* vec,
2163 const SearchParameters& searchParams = {}) const
2164 {
2165 assert(vec);
2166 if (this->size(*this) == 0) return false;
2167 if (!Base::root_node_) return false;
2168 float epsError = 1 + searchParams.eps;
2169
2170 // fixed or variable-sized container (depending on DIM)
2171 distance_vector_t dists;
2172 // Fill it with zeros.
2173 assign(
2174 dists, (DIM > 0 ? DIM : Base::dim_),
2175 static_cast<typename distance_vector_t::value_type>(0));
2176 DistanceType dist = this->computeInitialDistances(*this, vec, dists);
2177 searchLevel(result, vec, Base::root_node_, dist, dists, epsError);
2178 return result.full();
2179 }
2180
2196 const ElementType* query_point, const Size num_closest,
2197 IndexType* out_indices, DistanceType* out_distances,
2198 const SearchParameters& searchParams = {}) const
2199 {
2201 resultSet.init(out_indices, out_distances);
2202 findNeighbors(resultSet, query_point, searchParams);
2203 return resultSet.size();
2204 }
2205
2226 const ElementType* query_point, const DistanceType& radius,
2227 std::vector<ResultItem<IndexType, DistanceType>>& IndicesDists,
2228 const SearchParameters& searchParams = {}) const
2229 {
2231 radius, IndicesDists);
2232 const size_t nFound =
2233 radiusSearchCustomCallback(query_point, resultSet, searchParams);
2234 return nFound;
2235 }
2236
2242 template <class SEARCH_CALLBACK>
2244 const ElementType* query_point, SEARCH_CALLBACK& resultSet,
2245 const SearchParameters& searchParams = {}) const
2246 {
2247 findNeighbors(resultSet, query_point, searchParams);
2248 return resultSet.size();
2249 }
2250
2253 public:
2254 void computeBoundingBox(BoundingBox& bbox)
2255 {
2256 const auto dims = (DIM > 0 ? DIM : Base::dim_);
2257 resize(bbox, dims);
2258
2259 if (dataset_.kdtree_get_bbox(bbox))
2260 {
2261 // Done! It was implemented in derived class
2262 }
2263 else
2264 {
2265 const Size N = Base::size_;
2266 if (!N)
2267 throw std::runtime_error(
2268 "[nanoflann] computeBoundingBox() called but "
2269 "no data points found.");
2270 for (Dimension i = 0; i < dims; ++i)
2271 {
2272 bbox[i].low = bbox[i].high =
2273 this->dataset_get(*this, Base::vAcc_[0], i);
2274 }
2275 for (Offset k = 1; k < N; ++k)
2276 {
2277 for (Dimension i = 0; i < dims; ++i)
2278 {
2279 const auto val =
2280 this->dataset_get(*this, Base::vAcc_[k], i);
2281 if (val < bbox[i].low) bbox[i].low = val;
2282 if (val > bbox[i].high) bbox[i].high = val;
2283 }
2284 }
2285 }
2286 }
2287
2292 template <class RESULTSET>
2294 RESULTSET& result_set, const ElementType* vec, const NodePtr node,
2295 DistanceType mindist, distance_vector_t& dists,
2296 const float epsError) const
2297 {
2298 /* If this is a leaf node, then do check and return. */
2299 if ((node->child1 == nullptr) && (node->child2 == nullptr))
2300 {
2301 DistanceType worst_dist = result_set.worstDist();
2302 for (Offset i = node->node_type.lr.left;
2303 i < node->node_type.lr.right; ++i)
2304 {
2305 const IndexType index = Base::vAcc_[i]; // reorder... : i;
2306 if (treeIndex_[index] == -1) continue;
2307 DistanceType dist = distance_.evalMetric(
2308 vec, index, (DIM > 0 ? DIM : Base::dim_));
2309 if (dist < worst_dist)
2310 {
2311 if (!result_set.addPoint(
2312 static_cast<typename RESULTSET::DistanceType>(dist),
2313 static_cast<typename RESULTSET::IndexType>(
2314 Base::vAcc_[i])))
2315 {
2316 // the resultset doesn't want to receive any more
2317 // points, we're done searching!
2318 return; // false;
2319 }
2320 }
2321 }
2322 return;
2323 }
2324
2325 /* Which child branch should be taken first? */
2326 Dimension idx = node->node_type.sub.divfeat;
2327 ElementType val = vec[idx];
2328 DistanceType diff1 = val - node->node_type.sub.divlow;
2329 DistanceType diff2 = val - node->node_type.sub.divhigh;
2330
2331 NodePtr bestChild;
2332 NodePtr otherChild;
2333 DistanceType cut_dist;
2334 if ((diff1 + diff2) < 0)
2335 {
2336 bestChild = node->child1;
2337 otherChild = node->child2;
2338 cut_dist =
2339 distance_.accum_dist(val, node->node_type.sub.divhigh, idx);
2340 }
2341 else
2342 {
2343 bestChild = node->child2;
2344 otherChild = node->child1;
2345 cut_dist =
2346 distance_.accum_dist(val, node->node_type.sub.divlow, idx);
2347 }
2348
2349 /* Call recursively to search next level down. */
2350 searchLevel(result_set, vec, bestChild, mindist, dists, epsError);
2351
2352 DistanceType dst = dists[idx];
2353 mindist = mindist + cut_dist - dst;
2354 dists[idx] = cut_dist;
2355 if (mindist * epsError <= result_set.worstDist())
2356 {
2357 searchLevel(result_set, vec, otherChild, mindist, dists, epsError);
2358 }
2359 dists[idx] = dst;
2360 }
2361
2362 public:
2368 void saveIndex(std::ostream& stream) { saveIndex(*this, stream); }
2369
2375 void loadIndex(std::istream& stream) { loadIndex(*this, stream); }
2376};
2377
2392template <
2393 typename Distance, class DatasetAdaptor, int32_t DIM = -1,
2394 typename IndexType = uint32_t>
2396{
2397 public:
2398 using ElementType = typename Distance::ElementType;
2399 using DistanceType = typename Distance::DistanceType;
2400
2401 using Offset = typename KDTreeSingleIndexDynamicAdaptor_<
2402 Distance, DatasetAdaptor, DIM>::Offset;
2403 using Size = typename KDTreeSingleIndexDynamicAdaptor_<
2404 Distance, DatasetAdaptor, DIM>::Size;
2405 using Dimension = typename KDTreeSingleIndexDynamicAdaptor_<
2406 Distance, DatasetAdaptor, DIM>::Dimension;
2407
2408 protected:
2409 Size leaf_max_size_;
2410 Size treeCount_;
2411 Size pointCount_;
2412
2416 const DatasetAdaptor& dataset_;
2417
2420 std::vector<int> treeIndex_;
2421 std::unordered_set<int> removedPoints_;
2422
2423 KDTreeSingleIndexAdaptorParams index_params_;
2424
2425 Dimension dim_;
2426
2428 Distance, DatasetAdaptor, DIM, IndexType>;
2429 std::vector<index_container_t> index_;
2430
2431 public:
2434 const std::vector<index_container_t>& getAllIndices() const
2435 {
2436 return index_;
2437 }
2438
2439 private:
2441 int First0Bit(IndexType num)
2442 {
2443 int pos = 0;
2444 while (num & 1)
2445 {
2446 num = num >> 1;
2447 pos++;
2448 }
2449 return pos;
2450 }
2451
2453 void init()
2454 {
2455 using my_kd_tree_t = KDTreeSingleIndexDynamicAdaptor_<
2456 Distance, DatasetAdaptor, DIM, IndexType>;
2457 std::vector<my_kd_tree_t> index(
2458 treeCount_,
2459 my_kd_tree_t(dim_ /*dim*/, dataset_, treeIndex_, index_params_));
2460 index_ = index;
2461 }
2462
2463 public:
2464 Distance distance_;
2465
2482 const int dimensionality, const DatasetAdaptor& inputData,
2483 const KDTreeSingleIndexAdaptorParams& params =
2485 const size_t maximumPointCount = 1000000000U)
2486 : dataset_(inputData), index_params_(params), distance_(inputData)
2487 {
2488 treeCount_ = static_cast<size_t>(std::log2(maximumPointCount)) + 1;
2489 pointCount_ = 0U;
2490 dim_ = dimensionality;
2491 treeIndex_.clear();
2492 if (DIM > 0) dim_ = DIM;
2493 leaf_max_size_ = params.leaf_max_size;
2494 init();
2495 const size_t num_initial_points = dataset_.kdtree_get_point_count();
2496 if (num_initial_points > 0) { addPoints(0, num_initial_points - 1); }
2497 }
2498
2502 Distance, DatasetAdaptor, DIM, IndexType>&) = delete;
2503
2505 void addPoints(IndexType start, IndexType end)
2506 {
2507 const Size count = end - start + 1;
2508 int maxIndex = 0;
2509 treeIndex_.resize(treeIndex_.size() + count);
2510 for (IndexType idx = start; idx <= end; idx++)
2511 {
2512 const int pos = First0Bit(pointCount_);
2513 maxIndex = std::max(pos, maxIndex);
2514 treeIndex_[pointCount_] = pos;
2515
2516 const auto it = removedPoints_.find(idx);
2517 if (it != removedPoints_.end())
2518 {
2519 removedPoints_.erase(it);
2520 treeIndex_[idx] = pos;
2521 }
2522
2523 for (int i = 0; i < pos; i++)
2524 {
2525 for (int j = 0; j < static_cast<int>(index_[i].vAcc_.size());
2526 j++)
2527 {
2528 index_[pos].vAcc_.push_back(index_[i].vAcc_[j]);
2529 if (treeIndex_[index_[i].vAcc_[j]] != -1)
2530 treeIndex_[index_[i].vAcc_[j]] = pos;
2531 }
2532 index_[i].vAcc_.clear();
2533 }
2534 index_[pos].vAcc_.push_back(idx);
2535 pointCount_++;
2536 }
2537
2538 for (int i = 0; i <= maxIndex; ++i)
2539 {
2540 index_[i].freeIndex(index_[i]);
2541 if (!index_[i].vAcc_.empty()) index_[i].buildIndex();
2542 }
2543 }
2544
2546 void removePoint(size_t idx)
2547 {
2548 if (idx >= pointCount_) return;
2549 removedPoints_.insert(idx);
2550 treeIndex_[idx] = -1;
2551 }
2552
2569 template <typename RESULTSET>
2571 RESULTSET& result, const ElementType* vec,
2572 const SearchParameters& searchParams = {}) const
2573 {
2574 for (size_t i = 0; i < treeCount_; i++)
2575 {
2576 index_[i].findNeighbors(result, &vec[0], searchParams);
2577 }
2578 return result.full();
2579 }
2580};
2581
2607template <
2608 class MatrixType, int32_t DIM = -1, class Distance = nanoflann::metric_L2,
2609 bool row_major = true>
2611{
2612 using self_t =
2614 using num_t = typename MatrixType::Scalar;
2615 using IndexType = typename MatrixType::Index;
2616 using metric_t = typename Distance::template traits<
2617 num_t, self_t, IndexType>::distance_t;
2618
2620 metric_t, self_t,
2621 row_major ? MatrixType::ColsAtCompileTime
2622 : MatrixType::RowsAtCompileTime,
2623 IndexType>;
2624
2625 index_t* index_;
2627
2628 using Offset = typename index_t::Offset;
2629 using Size = typename index_t::Size;
2630 using Dimension = typename index_t::Dimension;
2631
2634 const Dimension dimensionality,
2635 const std::reference_wrapper<const MatrixType>& mat,
2636 const int leaf_max_size = 10, const unsigned int n_thread_build = 1)
2637 : m_data_matrix(mat)
2638 {
2639 const auto dims = row_major ? mat.get().cols() : mat.get().rows();
2640 if (static_cast<Dimension>(dims) != dimensionality)
2641 throw std::runtime_error(
2642 "Error: 'dimensionality' must match column count in data "
2643 "matrix");
2644 if (DIM > 0 && static_cast<int32_t>(dims) != DIM)
2645 throw std::runtime_error(
2646 "Data set dimensionality does not match the 'DIM' template "
2647 "argument");
2648 index_ = new index_t(
2649 dims, *this /* adaptor */,
2651 leaf_max_size, nanoflann::KDTreeSingleIndexAdaptorFlags::None,
2652 n_thread_build));
2653 }
2654
2655 public:
2658
2659 ~KDTreeEigenMatrixAdaptor() { delete index_; }
2660
2661 const std::reference_wrapper<const MatrixType> m_data_matrix;
2662
2671 void query(
2672 const num_t* query_point, const Size num_closest,
2673 IndexType* out_indices, num_t* out_distances) const
2674 {
2675 nanoflann::KNNResultSet<num_t, IndexType> resultSet(num_closest);
2676 resultSet.init(out_indices, out_distances);
2677 index_->findNeighbors(resultSet, query_point);
2678 }
2679
2683 const self_t& derived() const { return *this; }
2684 self_t& derived() { return *this; }
2685
2686 // Must return the number of data points
2687 Size kdtree_get_point_count() const
2688 {
2689 if (row_major)
2690 return m_data_matrix.get().rows();
2691 else
2692 return m_data_matrix.get().cols();
2693 }
2694
2695 // Returns the dim'th component of the idx'th point in the class:
2696 num_t kdtree_get_pt(const IndexType idx, size_t dim) const
2697 {
2698 if (row_major)
2699 return m_data_matrix.get().coeff(idx, IndexType(dim));
2700 else
2701 return m_data_matrix.get().coeff(IndexType(dim), idx);
2702 }
2703
2704 // Optional bounding-box computation: return false to default to a standard
2705 // bbox computation loop.
2706 // Return true if the BBOX was already computed by the class and returned
2707 // in "bb" so it can be avoided to redo it again. Look at bb.size() to
2708 // find out the expected dimensionality (e.g. 2 or 3 for point clouds)
2709 template <class BBOX>
2710 bool kdtree_get_bbox(BBOX& /*bb*/) const
2711 {
2712 return false;
2713 }
2714
2717}; // end of KDTreeEigenMatrixAdaptor
// end of grouping
2721} // namespace nanoflann
Definition nanoflann.hpp:1012
void freeIndex(Derived &obj)
Definition nanoflann.hpp:1016
BoundingBox root_bbox_
Definition nanoflann.hpp:1090
Size veclen(const Derived &obj)
Definition nanoflann.hpp:1105
void saveIndex(const Derived &obj, std::ostream &stream) const
Definition nanoflann.hpp:1456
Size usedMemory(Derived &obj)
Definition nanoflann.hpp:1118
typename array_or_vector< DIM, DistanceType >::type distance_vector_t
Definition nanoflann.hpp:1087
void planeSplit(const Derived &obj, const Offset ind, const Size count, const Dimension cutfeat, const DistanceType &cutval, Offset &lim1, Offset &lim2)
Definition nanoflann.hpp:1366
NodePtr divideTree(Derived &obj, const Offset left, const Offset right, BoundingBox &bbox)
Definition nanoflann.hpp:1146
std::vector< IndexType > vAcc_
Definition nanoflann.hpp:1030
Size size(const Derived &obj) const
Definition nanoflann.hpp:1102
NodePtr divideTreeConcurrent(Derived &obj, const Offset left, const Offset right, BoundingBox &bbox, std::atomic< unsigned int > &thread_count, std::mutex &mutex)
Definition nanoflann.hpp:1217
void loadIndex(Derived &obj, std::istream &stream)
Definition nanoflann.hpp:1471
PooledAllocator pool_
Definition nanoflann.hpp:1099
ElementType dataset_get(const Derived &obj, IndexType element, Dimension component) const
Helper accessor to the dataset points:
Definition nanoflann.hpp:1108
typename array_or_vector< DIM, Interval >::type BoundingBox
Definition nanoflann.hpp:1083
Definition nanoflann.hpp:1530
bool searchLevel(RESULTSET &result_set, const ElementType *vec, const NodePtr node, DistanceType mindist, distance_vector_t &dists, const float epsError) const
Definition nanoflann.hpp:1866
void saveIndex(std::ostream &stream) const
Definition nanoflann.hpp:1949
Size radiusSearch(const ElementType *query_point, const DistanceType &radius, std::vector< ResultItem< IndexType, DistanceType > > &IndicesDists, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:1758
void init_vind()
Definition nanoflann.hpp:1817
void buildIndex()
Definition nanoflann.hpp:1641
const DatasetAdaptor & dataset_
Definition nanoflann.hpp:1538
KDTreeSingleIndexAdaptor(const KDTreeSingleIndexAdaptor< Distance, DatasetAdaptor, DIM, index_t > &)=delete
bool findNeighbors(RESULTSET &result, const ElementType *vec, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:1689
Size rknnSearch(const ElementType *query_point, const Size num_closest, IndexType *out_indices, DistanceType *out_distances, const DistanceType &radius) const
Definition nanoflann.hpp:1800
Size radiusSearchCustomCallback(const ElementType *query_point, SEARCH_CALLBACK &resultSet, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:1776
typename Base::distance_vector_t distance_vector_t
Definition nanoflann.hpp:1568
void loadIndex(std::istream &stream)
Definition nanoflann.hpp:1959
typename Base::BoundingBox BoundingBox
Definition nanoflann.hpp:1564
KDTreeSingleIndexAdaptor(const Dimension dimensionality, const DatasetAdaptor &inputData, const KDTreeSingleIndexAdaptorParams &params, Args &&... args)
Definition nanoflann.hpp:1591
Size knnSearch(const ElementType *query_point, const Size num_closest, IndexType *out_indices, DistanceType *out_distances) const
Definition nanoflann.hpp:1729
Definition nanoflann.hpp:2008
Size radiusSearch(const ElementType *query_point, const DistanceType &radius, std::vector< ResultItem< IndexType, DistanceType > > &IndicesDists, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:2225
KDTreeSingleIndexDynamicAdaptor_(const Dimension dimensionality, const DatasetAdaptor &inputData, std::vector< int > &treeIndex, const KDTreeSingleIndexAdaptorParams &params=KDTreeSingleIndexAdaptorParams())
Definition nanoflann.hpp:2060
typename Base::BoundingBox BoundingBox
Definition nanoflann.hpp:2039
const DatasetAdaptor & dataset_
The source of our data.
Definition nanoflann.hpp:2013
Size radiusSearchCustomCallback(const ElementType *query_point, SEARCH_CALLBACK &resultSet, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:2243
KDTreeSingleIndexDynamicAdaptor_(const KDTreeSingleIndexDynamicAdaptor_ &rhs)=default
void buildIndex()
Definition nanoflann.hpp:2111
void saveIndex(std::ostream &stream)
Definition nanoflann.hpp:2368
typename Base::distance_vector_t distance_vector_t
Definition nanoflann.hpp:2043
void loadIndex(std::istream &stream)
Definition nanoflann.hpp:2375
Size knnSearch(const ElementType *query_point, const Size num_closest, IndexType *out_indices, DistanceType *out_distances, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:2195
void searchLevel(RESULTSET &result_set, const ElementType *vec, const NodePtr node, DistanceType mindist, distance_vector_t &dists, const float epsError) const
Definition nanoflann.hpp:2293
bool findNeighbors(RESULTSET &result, const ElementType *vec, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:2161
KDTreeSingleIndexDynamicAdaptor_ operator=(const KDTreeSingleIndexDynamicAdaptor_ &rhs)
Definition nanoflann.hpp:2092
Definition nanoflann.hpp:2396
bool findNeighbors(RESULTSET &result, const ElementType *vec, const SearchParameters &searchParams={}) const
Definition nanoflann.hpp:2570
const DatasetAdaptor & dataset_
The source of our data.
Definition nanoflann.hpp:2416
void removePoint(size_t idx)
Definition nanoflann.hpp:2546
void addPoints(IndexType start, IndexType end)
Definition nanoflann.hpp:2505
KDTreeSingleIndexDynamicAdaptor(const int dimensionality, const DatasetAdaptor &inputData, const KDTreeSingleIndexAdaptorParams &params=KDTreeSingleIndexAdaptorParams(), const size_t maximumPointCount=1000000000U)
Definition nanoflann.hpp:2481
std::vector< int > treeIndex_
Definition nanoflann.hpp:2420
const std::vector< index_container_t > & getAllIndices() const
Definition nanoflann.hpp:2434
Dimension dim_
Dimensionality of each data point.
Definition nanoflann.hpp:2425
KDTreeSingleIndexDynamicAdaptor(const KDTreeSingleIndexDynamicAdaptor< Distance, DatasetAdaptor, DIM, IndexType > &)=delete
Definition nanoflann.hpp:200
bool addPoint(DistanceType dist, IndexType index)
Definition nanoflann.hpp:234
DistanceType worstDist() const
Definition nanoflann.hpp:271
Definition nanoflann.hpp:860
~PooledAllocator()
Definition nanoflann.hpp:896
void free_all()
Definition nanoflann.hpp:899
void * malloc(const size_t req_size)
Definition nanoflann.hpp:915
T * allocate(const size_t count=1)
Definition nanoflann.hpp:966
PooledAllocator()
Definition nanoflann.hpp:891
Definition nanoflann.hpp:289
bool addPoint(DistanceType dist, IndexType index)
Definition nanoflann.hpp:330
DistanceType worstDist() const
Definition nanoflann.hpp:367
Definition nanoflann.hpp:384
ResultItem< IndexType, DistanceType > worst_item() const
Definition nanoflann.hpp:427
bool addPoint(DistanceType dist, IndexType index)
Definition nanoflann.hpp:415
std::enable_if< has_assign< Container >::value, void >::type assign(Container &c, const size_t nElements, const T &value)
Definition nanoflann.hpp:144
T pi_const()
Definition nanoflann.hpp:87
std::enable_if< has_resize< Container >::value, void >::type resize(Container &c, const size_t nElements)
Definition nanoflann.hpp:122
Definition nanoflann.hpp:162
bool operator()(const PairType &p1, const PairType &p2) const
Definition nanoflann.hpp:165
Definition nanoflann.hpp:1065
Definition nanoflann.hpp:1040
DistanceType divlow
The values used for subdivision.
Definition nanoflann.hpp:1053
Offset right
Indices of points in leaf node.
Definition nanoflann.hpp:1047
union nanoflann::KDTreeBaseClass::Node::@0 node_type
Dimension divfeat
Definition nanoflann.hpp:1051
Node * child1
Definition nanoflann.hpp:1058
Definition nanoflann.hpp:2611
void query(const num_t *query_point, const Size num_closest, IndexType *out_indices, num_t *out_distances) const
Definition nanoflann.hpp:2671
KDTreeEigenMatrixAdaptor(const self_t &)=delete
typename index_t::Offset Offset
Definition nanoflann.hpp:2628
KDTreeEigenMatrixAdaptor(const Dimension dimensionality, const std::reference_wrapper< const MatrixType > &mat, const int leaf_max_size=10, const unsigned int n_thread_build=1)
Constructor: takes a const ref to the matrix object with the data points.
Definition nanoflann.hpp:2633
Definition nanoflann.hpp:811
Definition nanoflann.hpp:500
Definition nanoflann.hpp:562
Definition nanoflann.hpp:627
Definition nanoflann.hpp:483
Definition nanoflann.hpp:181
DistanceType second
Distance from sample to query point.
Definition nanoflann.hpp:189
IndexType first
Index of the sample in the dataset.
Definition nanoflann.hpp:188
Definition nanoflann.hpp:672
DistanceType accum_dist(const U a, const V b, const size_t) const
Definition nanoflann.hpp:690
Definition nanoflann.hpp:717
Definition nanoflann.hpp:830
bool sorted
distance (default: true)
Definition nanoflann.hpp:837
float eps
search for eps-approximate neighbours (default: 0)
Definition nanoflann.hpp:836
Definition nanoflann.hpp:982
Definition nanoflann.hpp:109
Definition nanoflann.hpp:98
Definition nanoflann.hpp:747
Definition nanoflann.hpp:744
Definition nanoflann.hpp:757
Definition nanoflann.hpp:767
Definition nanoflann.hpp:764
Definition nanoflann.hpp:754
Definition nanoflann.hpp:776
Definition nanoflann.hpp:773
Definition nanoflann.hpp:785
Definition nanoflann.hpp:782