Loading...
Searching...
No Matches
Kd_tree_search.h
/* This file is part of the Gudhi Library - https://gudhi.inria.fr/ - which is released under MIT.
* See file LICENSE or go to https://gudhi.inria.fr/licensing/ for full license details.
* Author(s): Clement Jamin
*
* Copyright (C) 2016 Inria
*
* Modification(s):
* - 2019/08 Vincent Rouvreau: Fix issue #10 for CGAL and Eigen3
* - YYYY/MM Author: Description of the modification
*/
#ifndef KD_TREE_SEARCH_H_
#define KD_TREE_SEARCH_H_
#include <gudhi/Debug_utils.h>
#include <CGAL/Orthogonal_k_neighbor_search.h>
#include <CGAL/Orthogonal_incremental_neighbor_search.h>
#include <CGAL/Search_traits.h>
#include <CGAL/Search_traits_adapter.h>
#include <CGAL/property_map.h>
#include <CGAL/version.h> // for CGAL_VERSION_NR
#include <Eigen/src/Core/util/Macros.h> // for EIGEN_VERSION_AT_LEAST
#include <boost/property_map/property_map.hpp>
#include <boost/iterator/counting_iterator.hpp>
#include <cstddef>
#include <vector>
// Make compilation fail - required for external projects - https://github.com/GUDHI/gudhi-devel/issues/10
#if CGAL_VERSION_NR < 1041101000
# error Kd_tree_search is only available for CGAL >= 4.11
#endif
#if !EIGEN_VERSION_AT_LEAST(3,1,0)
# error Kd_tree_search is only available for Eigen3 >= 3.1.0 installed with CGAL
#endif
namespace Gudhi {
namespace spatial_searching {
template <typename Search_traits, typename Point_range>
class Kd_tree_search {
typedef boost::iterator_property_map<
typename Point_range::const_iterator,
CGAL::Identity_property_map<std::ptrdiff_t> > Point_property_map;
public:
typedef Search_traits Traits;
typedef typename Traits::FT FT;
typedef typename Point_range::value_type Point;
typedef CGAL::Search_traits<
FT, Point,
typename Traits::Cartesian_const_iterator_d,
typename Traits::Construct_cartesian_const_iterator_d,
typename Traits::Dimension> Traits_base;
typedef CGAL::Search_traits_adapter<
std::ptrdiff_t,
Point_property_map,
Traits_base> STraits;
typedef CGAL::Distance_adapter<
std::ptrdiff_t,
Point_property_map,
CGAL::Euclidean_distance<Traits_base> > Orthogonal_distance;
typedef CGAL::Orthogonal_k_neighbor_search<STraits> K_neighbor_search;
typedef typename K_neighbor_search::Tree Tree;
typedef typename K_neighbor_search::Distance Distance;
typedef K_neighbor_search KNS_range;
typedef CGAL::Orthogonal_incremental_neighbor_search<
STraits, Distance, CGAL::Sliding_midpoint<STraits>, Tree>
Incremental_neighbor_search;
typedef Incremental_neighbor_search INS_range;
// Because CGAL::Fuzzy_sphere takes the radius and not its square
struct Sphere_for_kdtree_search
{
typedef typename Traits::Point_d Point_d;
typedef typename Traits::FT FT;
typedef typename Traits::Dimension D;
typedef D Dimension;
private:
STraits traits;
Point_d c;
FT sqradmin, sqradmax;
bool use_max;
public:
// `prefer_max` means that we prefer outputting more points at squared distance between r2min and r2max,
// while `!prefer_max` means we prefer fewer.
Sphere_for_kdtree_search(Point_d const& c_, FT const& r2min, FT const& r2max, bool prefer_max=true, STraits const& traits_ = {})
: traits(traits_), c(c_), sqradmin(r2min), sqradmax(r2max), use_max(prefer_max)
{ GUDHI_CHECK(r2min >= 0 && r2max >= r2min, "0 <= r2min <= r2max"); }
bool contains(std::ptrdiff_t i) const {
const Point_d& p = get(traits.point_property_map(), i);
auto ccci = traits.construct_cartesian_const_iterator_d_object();
return contains_point_given_as_coordinates(ccci(p), ccci(p, 0));
}
template <typename Coord_iterator>
bool contains_point_given_as_coordinates(Coord_iterator pi, Coord_iterator) const {
FT distance = 0;
auto ccci = traits.construct_cartesian_const_iterator_d_object();
auto ci = ccci(c);
auto ce = ccci(c, 0);
FT const& limit = use_max ? sqradmax : sqradmin;
while (ci != ce) {
distance += CGAL::square(*pi++ - *ci++);
// I think Clément advised to check the distance at every step instead of
// just at the end, especially when the dimension becomes large. Distance
// isn't part of the concept anyway.
if (distance > limit) return false;
}
return true;
}
bool inner_range_intersects(CGAL::Kd_tree_rectangle<FT, D> const& rect) const {
auto ccci = traits.construct_cartesian_const_iterator_d_object();
FT distance = 0;
auto ci = ccci(c);
auto ce = ccci(c, 0);
for (int i = 0; ci != ce; ++i, ++ci) {
distance += CGAL::square(CGAL::max<FT>(CGAL::max<FT>(*ci - rect.max_coord(i), rect.min_coord(i) - *ci), 0 ));
if (distance > sqradmin) return false;
}
return true;
}
bool outer_range_contains(CGAL::Kd_tree_rectangle<FT, D> const& rect) const {
auto ccci = traits.construct_cartesian_const_iterator_d_object();
FT distance = 0;
auto ci = ccci(c);
auto ce = ccci(c, 0);
for (int i = 0; ci != ce; ++i, ++ci) {
distance += CGAL::square(CGAL::max<FT>(*ci - rect.min_coord(i), rect.max_coord(i) - *ci));
if (distance > sqradmax) return false;
}
return true;
}
};
Kd_tree_search(Point_range const& points)
: m_points(points),
m_tree(boost::counting_iterator<std::ptrdiff_t>(0),
boost::counting_iterator<std::ptrdiff_t>(points.size()),
typename Tree::Splitter(),
STraits(std::begin(points))) {
// Build the tree now (we don't want to wait for the first query)
m_tree.build();
}
template <typename Point_indices_range>
Kd_tree_search(
Point_range const& points,
Point_indices_range const& only_these_points)
: m_points(points),
m_tree(
only_these_points.begin(), only_these_points.end(),
typename Tree::Splitter(),
STraits(std::begin(points))) {
// Build the tree now (we don't want to wait for the first query)
m_tree.build();
}
Kd_tree_search(
Point_range const& points,
std::size_t begin_idx, std::size_t past_the_end_idx)
: m_points(points),
m_tree(
boost::counting_iterator<std::ptrdiff_t>(begin_idx),
boost::counting_iterator<std::ptrdiff_t>(past_the_end_idx),
typename Tree::Splitter(),
STraits(std::begin(points))) {
// Build the tree now (we don't want to wait for the first query)
m_tree.build();
}
// Be careful, this function invalidates the tree,
// which will be recomputed at the next query
void insert(std::ptrdiff_t point_idx) {
m_tree.insert(point_idx);
}
KNS_range k_nearest_neighbors(
Point const& p,
unsigned int k,
bool sorted = true,
FT eps = FT(0)) const {
// Initialize the search structure, and search all N points
// Note that we need to pass the Distance explicitly since it needs to
// know the property map
K_neighbor_search search(
m_tree,
p,
k,
eps,
true,
Orthogonal_distance(std::begin(m_points)), sorted);
return search;
}
INS_range incremental_nearest_neighbors(Point const& p, FT eps = FT(0)) const {
// Initialize the search structure, and search all N points
// Note that we need to pass the Distance explicitly since it needs to
// know the property map
Incremental_neighbor_search search(
m_tree,
p,
eps,
true,
Orthogonal_distance(std::begin(m_points)) );
return search;
}
KNS_range k_furthest_neighbors(
Point const& p,
unsigned int k,
bool sorted = true,
FT eps = FT(0)) const {
// Initialize the search structure, and search all N points
// Note that we need to pass the Distance explicitly since it needs to
// know the property map
K_neighbor_search search(
m_tree,
p,
k,
eps,
false,
Orthogonal_distance(std::begin(m_points)), sorted);
return search;
}
INS_range incremental_furthest_neighbors(Point const& p, FT eps = FT(0)) const {
// Initialize the search structure, and search all N points
// Note that we need to pass the Distance explicitly since it needs to
// know the property map
Incremental_neighbor_search search(
m_tree,
p,
eps,
false,
Orthogonal_distance(std::begin(m_points)) );
return search;
}
template <typename OutputIterator>
void all_near_neighbors(Point const& p,
FT const& radius,
OutputIterator it,
FT eps = FT(0)) const {
all_near_neighbors2(p, CGAL::square(radius - eps), CGAL::square(radius + eps), it);
}
template <typename OutputIterator>
void all_near_neighbors2(Point const& p,
FT const& sq_radius_min,
FT const& sq_radius_max,
OutputIterator it) const {
m_tree.search(it, Sphere_for_kdtree_search(p, sq_radius_min, sq_radius_max, true, m_tree.traits()));
}
int tree_depth() const {
return m_tree.root()->depth();
}
private:
Point_range const& m_points;
Tree m_tree;
};
} // namespace spatial_searching
} // namespace Gudhi
#endif // KD_TREE_SEARCH_H_
Gudhi::spatial_searching::Kd_tree_search< Kernel_, Point_range >::Point
Point_range::value_type Point
Definition Kd_tree_search.h:81
Gudhi::spatial_searching::Kd_tree_search::k_furthest_neighbors
KNS_range k_furthest_neighbors(Point const &p, unsigned int k, bool sorted=true, FT eps=FT(0)) const
Search for the k-furthest points from a query point.
Definition Kd_tree_search.h:291
Gudhi::spatial_searching::Kd_tree_search< Kernel_, Point_range >::KNS_range
K_neighbor_search KNS_range
Definition Kd_tree_search.h:104
Gudhi::spatial_searching::Kd_tree_search::all_near_neighbors
void all_near_neighbors(Point const &p, FT const &radius, OutputIterator it, FT eps=FT(0)) const
Search for all the neighbors in a ball.
Definition Kd_tree_search.h:338
Gudhi::spatial_searching::Kd_tree_search< Kernel_, Point_range >::FT
Traits::FT FT
Definition Kd_tree_search.h:79
Gudhi::spatial_searching::Kd_tree_search::incremental_nearest_neighbors
INS_range incremental_nearest_neighbors(Point const &p, FT eps=FT(0)) const
Search incrementally for the nearest neighbors from a query point.
Definition Kd_tree_search.h:271
Gudhi::spatial_searching::Kd_tree_search::Kd_tree_search
Kd_tree_search(Point_range const &points)
Constructor.
Definition Kd_tree_search.h:187
Gudhi::spatial_searching::Kd_tree_search< Kernel_, Point_range >::Traits
Kernel_ Traits
Definition Kd_tree_search.h:77
Gudhi::spatial_searching::Kd_tree_search::Kd_tree_search
Kd_tree_search(Point_range const &points, Point_indices_range const &only_these_points)
Constructor.
Definition Kd_tree_search.h:203
Gudhi::spatial_searching::Kd_tree_search::k_nearest_neighbors
KNS_range k_nearest_neighbors(Point const &p, unsigned int k, bool sorted=true, FT eps=FT(0)) const
Search for the k-nearest neighbors from a query point.
Definition Kd_tree_search.h:245
Gudhi::spatial_searching::Kd_tree_search::Kd_tree_search
Kd_tree_search(Point_range const &points, std::size_t begin_idx, std::size_t past_the_end_idx)
Constructor.
Definition Kd_tree_search.h:220
Gudhi::spatial_searching::Kd_tree_search::all_near_neighbors2
void all_near_neighbors2(Point const &p, FT const &sq_radius_min, FT const &sq_radius_max, OutputIterator it) const
Search for all the neighbors in a ball. This is similar to all_near_neighbors but takes directly the ...
Definition Kd_tree_search.h:354
Gudhi::spatial_searching::Kd_tree_search< Kernel_, Point_range >::INS_range
Incremental_neighbor_search INS_range
Definition Kd_tree_search.h:112
Gudhi::spatial_searching::Kd_tree_search::incremental_furthest_neighbors
INS_range incremental_furthest_neighbors(Point const &p, FT eps=FT(0)) const
Search incrementally for the furthest neighbors from a query point.
Definition Kd_tree_search.h:317
Gudhi
Gudhi namespace.
Definition SimplicialComplexForAlpha.h:14
std
STL namespace.