diff options
Diffstat (limited to 'src/Spatial_searching')
6 files changed, 443 insertions, 0 deletions
diff --git a/src/Spatial_searching/doc/Intro_spatial_searching.h b/src/Spatial_searching/doc/Intro_spatial_searching.h new file mode 100644 index 00000000..0a857221 --- /dev/null +++ b/src/Spatial_searching/doc/Intro_spatial_searching.h @@ -0,0 +1,58 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Clement Jamin + * + * Copyright (C) 2016 INRIA + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + */ + +#ifndef DOC_SPATIAL_SEARCHING_INTRO_SPATIAL_SEARCHING_H_ +#define DOC_SPATIAL_SEARCHING_INTRO_SPATIAL_SEARCHING_H_ + +// needs namespaces for Doxygen to link on classes +namespace Gudhi { +namespace spatial_searching { + +/** \defgroup spatial_searching Spatial_searching + * + * \author Clément Jamin + * + * @{ + * + * \section introduction Introduction + * + * This Gudhi component is a wrapper around + * <a target="_blank" href="http://doc.cgal.org/latest/Spatial_searching/index.html">CGAL dD spatial searching algorithms</a>. + * It provides a simplified API to perform (approximate) nearest neighbor searches. Contrary to CGAL default behavior, the tree + * does not store the points themselves, but stores indices. + * + * \section spatial_searching_examples Example + * + * This example generates 500 random points, then performs queries for nearest points using different methods. + * + * \include Spatial_searching/example_spatial_searching.cpp + * + * \copyright GNU General Public License v3. + * \verbatim Contact: gudhi-users@lists.gforge.inria.fr \endverbatim + */ +/** @} */ // end defgroup spatial_searching + +} // namespace spatial_searching + +} // namespace Gudhi + +#endif // DOC_SPATIAL_SEARCHING_INTRO_SPATIAL_SEARCHING_H_ diff --git a/src/Spatial_searching/example/CMakeLists.txt b/src/Spatial_searching/example/CMakeLists.txt new file mode 100644 index 00000000..0885c24c --- /dev/null +++ b/src/Spatial_searching/example/CMakeLists.txt @@ -0,0 +1,24 @@ +cmake_minimum_required(VERSION 2.6) +project(Spatial_searching_examples) + +if(CGAL_FOUND) + if (NOT CGAL_VERSION VERSION_LESS 4.8.0) + message(STATUS "CGAL version: ${CGAL_VERSION}.") + + find_package(Eigen3 3.1.0) + if (EIGEN3_FOUND) + message(STATUS "Eigen3 version: ${EIGEN3_VERSION}.") + include( ${EIGEN3_USE_FILE} ) + include_directories (BEFORE "../../include") + + add_executable( Spatial_searching_example_spatial_searching example_spatial_searching.cpp ) + target_link_libraries(Spatial_searching_example_spatial_searching ${CGAL_LIBRARY}) + else() + message(WARNING "Eigen3 not found. Version 3.1.0 is required for the Tangential_complex examples.") + endif() + else() + message(WARNING "CGAL version: ${CGAL_VERSION} is too old to compile Spatial_searching examples. Version 4.8.0 is required.") + endif () +else() + message(WARNING "CGAL not found. It is required for the Spatial_searching examples.") +endif() diff --git a/src/Spatial_searching/example/example_spatial_searching.cpp b/src/Spatial_searching/example/example_spatial_searching.cpp new file mode 100644 index 00000000..1a807b90 --- /dev/null +++ b/src/Spatial_searching/example/example_spatial_searching.cpp @@ -0,0 +1,42 @@ +#include <gudhi/Spatial_tree_data_structure.h> + +#include <CGAL/Epick_d.h> +#include <CGAL/Random.h> + +#include <array> +#include <vector> + +namespace gss = Gudhi::spatial_searching; + +int main (void) +{ + typedef CGAL::Epick_d<CGAL::Dimension_tag<4> > K; + typedef typename K::FT FT; + typedef typename K::Point_d Point; + typedef std::vector<Point> Points; + + typedef gss::Spatial_tree_data_structure<K, Points> Points_ds; + + CGAL::Random rd; + + Points points; + for (int i = 0; i < 500; ++i) + points.push_back(Point(std::array<FT, 4>({ rd.get_double(-1.,1),rd.get_double(-1.,1),rd.get_double(-1.,1),rd.get_double(-1.,1) }))); + + Points_ds points_ds(points); + + // 20-nearest neighbor query + std::cout << "20 nearest neighbors:\n"; + auto kns_range = points_ds.query_ANN(points[20], 10, true); + for (auto const& nghb : kns_range) + std::cout << nghb.first << " (sq. dist. = " << nghb.second << ")\n"; + + // Incremental nearest neighbor query + std::cout << "Incremental nearest neighbors:\n"; + auto ins_range = points_ds.query_incremental_ANN(points[45]); + // Get all the neighbors that are closer than 0.5 + for (auto ins_iterator = ins_range.begin(); ins_iterator->second < 0.5*0.5 ; ++ins_iterator) + std::cout << ins_iterator->first << " (sq. dist. = " << ins_iterator->second << ")\n"; + + return 0; +} diff --git a/src/Spatial_searching/include/gudhi/Spatial_tree_data_structure.h b/src/Spatial_searching/include/gudhi/Spatial_tree_data_structure.h new file mode 100644 index 00000000..5a29b153 --- /dev/null +++ b/src/Spatial_searching/include/gudhi/Spatial_tree_data_structure.h @@ -0,0 +1,226 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Clement Jamin + * + * Copyright (C) 2016 INRIA + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + */ + +#ifndef GUDHI_SPATIAL_TREE_DS_H_ +#define GUDHI_SPATIAL_TREE_DS_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 <boost/iterator/counting_iterator.hpp> + +#include <cstddef> +#include <vector> + +namespace Gudhi { +namespace spatial_searching { + + + /** + * \class Spatial_tree_data_structure Spatial_tree_data_structure.h gudhi/Spatial_tree_data_structure.h + * \brief Spatial tree data structure to perform (approximate) nearest neighbor search. + * + * \ingroup spatial_searching + * + * \details + * The class Spatial_tree_data_structure is a tree-based data structure, based on + * <a target="_blank" href="http://doc.cgal.org/latest/Spatial_searching/index.html">CGAL dD spatial searching data structures</a>. + * It provides a simplified API to perform (approximate) nearest neighbor searches. Contrary to CGAL default behavior, the tree + * does not store the points themselves, but stores indices. + * + * There are two types of queries: the <i>k-nearest neighbor query</i>, where <i>k</i> is fixed and the <i>k</i> nearest points are + * computed right away, + * and the <i>incremental nearest neighbor query</i>, where no number of neighbors is provided during the call, as the + * neighbors will be computed incrementally when the iterator on the range is incremented. + * + * \tparam K must be a model of the <a target="_blank" + * href="http://doc.cgal.org/latest/Spatial_searching/classSearchTraits.html">SearchTraits</a> + * concept, such as the <a target="_blank" + * href="http://doc.cgal.org/latest/Kernel_d/classCGAL_1_1Epick__d.html">CGAL::Epick_d</a> class, which + * can be static if you know the ambiant dimension at compile-time, or dynamic if you don't. + * \tparam Point_container_ is the type of the container where points are stored (on the user side). + * It must provide random-access via `operator[]` and the points should be stored contiguously in memory. + * `std::vector` is a good candidate. + */ +template <typename K, typename Point_container_> +class Spatial_tree_data_structure +{ +public: + typedef typename Point_container_::value_type Point; + typedef K Kernel; + /// Number type used for distances. + typedef typename Kernel::FT FT; + + typedef CGAL::Search_traits< + FT, Point, + typename Kernel::Cartesian_const_iterator_d, + typename Kernel::Construct_cartesian_const_iterator_d> Traits_base; + // using a pointer as a special property map type + typedef CGAL::Search_traits_adapter< + std::ptrdiff_t, Point*, Traits_base> STraits; + + typedef CGAL::Orthogonal_k_neighbor_search<STraits> K_neighbor_search; + typedef typename K_neighbor_search::Tree Tree; + typedef typename K_neighbor_search::Distance Distance; + /// \brief The range returned by a k-nearest neighbor search. + /// Its value type is `std::pair<std::size_t, FT>` where `first` is the index + /// of a point P and `second` is the squared distance between P and the query point. + typedef K_neighbor_search KNS_range; + + typedef CGAL::Orthogonal_incremental_neighbor_search< + STraits, Distance, CGAL::Sliding_midpoint<STraits>, Tree> + Incremental_neighbor_search; + /// \brief The range returned by an incremental nearest neighbor search. + /// Its value type is `std::pair<std::size_t, FT>` where `first` is the index + /// of a point P and `second` is the squared distance between P and the query point. + typedef Incremental_neighbor_search INS_range; + + /// \brief Constructor + /// @param[in] points Const reference to the point container. This container + /// is not copied, so it should not be destroyed or modified afterwards. + Spatial_tree_data_structure(Point_container_ 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((Point*)&(points[0])) ) + { + // Build the tree now (we don't want to wait for the first query) + m_tree.build(); + } + + /// \brief Constructor + /// @param[in] points Const reference to the point container. This container + /// is not copied, so it should not be destroyed or modified afterwards. + /// @param[in] Only_these_points Specifies the indices of the points that + /// should be actually inserted into the tree. The other points are ignored. + template <typename Point_indices_range> + Spatial_tree_data_structure( + Point_container_ 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((Point*)&(points[0]))) + { + // Build the tree now (we don't want to wait for the first query) + m_tree.build(); + } + + /// \brief Constructor + /// @param[in] points Const reference to the point container. This container + /// is not copied, so it should not be destroyed or modified afterwards. + /// @param[in] begin_idx, past_the_end_idx Define the subset of the points that + /// should be actually inserted into the tree. The other points are ignored. + Spatial_tree_data_structure( + Point_container_ 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((Point*)&(points[0])) ) + { + // Build the tree now (we don't want to wait for the first query) + m_tree.build(); + } + + /*Point_container_ &points() + { + return m_points; + } + + const Point_container_ &points() const + { + return m_points; + }*/ + + // 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); + } + + /// \brief Search for the k-nearest neighbors from a query point. + /// @param[in] p The query point. + /// @param[in] k Number of nearest points to search. + /// @param[in] sorted Indicates if the computed sequence of k-nearest neighbors needs to be sorted. + /// @param[in] eps Approximation factor. + /// @return A range containing the k-nearest neighbors. + KNS_range query_ANN(const + Point &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, + CGAL::Distance_adapter<std::ptrdiff_t,Point*,CGAL::Euclidean_distance<Traits_base> >( + (Point*)&(m_points[0])), + sorted); + + return search; + } + + /// \brief Search incrementally for the nearest neighbors from a query point. + /// @param[in] p The query point. + /// @param[in] eps Approximation factor. + /// @return A range containing the neighbors sorted by their distance to p. + /// All the neighbors are not computed by this function, but they will be + /// computed incrementally when the iterator on the range is incremented. + INS_range query_incremental_ANN(const Point &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, + CGAL::Distance_adapter<std::ptrdiff_t, Point*, CGAL::Euclidean_distance<Traits_base> >( + (Point*)&(m_points[0])) ); + + return search; + } + +protected: + Point_container_ const& m_points; + Tree m_tree; +}; + +} // namespace spatial_searching +} // namespace Gudhi + +#endif // GUDHI_SPATIAL_TREE_DS_H_ diff --git a/src/Spatial_searching/test/CMakeLists.txt b/src/Spatial_searching/test/CMakeLists.txt new file mode 100644 index 00000000..41059286 --- /dev/null +++ b/src/Spatial_searching/test/CMakeLists.txt @@ -0,0 +1,28 @@ +cmake_minimum_required(VERSION 2.6) +project(Spatial_searching_tests) + +if (GCOVR_PATH) + # for gcovr to make coverage reports - Corbera Jenkins plugin + set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fprofile-arcs -ftest-coverage") +endif() +if (GPROF_PATH) + # for gprof to make coverage reports - Jenkins + set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -pg") +endif() + +if(CGAL_FOUND) + find_package(Eigen3 3.1.0) + if (EIGEN3_FOUND) + message(STATUS "Eigen3 version: ${EIGEN3_VERSION}.") + include( ${EIGEN3_USE_FILE} ) + include_directories (BEFORE "../../include") + + add_executable( Spatial_searching_test_Spatial_tree_data_structure test_Spatial_tree_data_structure.cpp ) + target_link_libraries(Spatial_searching_test_Spatial_tree_data_structure + ${CGAL_LIBRARY} ${Boost_SYSTEM_LIBRARY} ${Boost_UNIT_TEST_FRAMEWORK_LIBRARY}) + else() + message(WARNING "Eigen3 not found. Version 3.1.0 is required for the Spatial_searching tests.") + endif() +else() + message(WARNING "CGAL not found. It is required for the Spatial_searching tests.") +endif() diff --git a/src/Spatial_searching/test/test_Spatial_tree_data_structure.cpp b/src/Spatial_searching/test/test_Spatial_tree_data_structure.cpp new file mode 100644 index 00000000..916710ef --- /dev/null +++ b/src/Spatial_searching/test/test_Spatial_tree_data_structure.cpp @@ -0,0 +1,65 @@ +/* This file is part of the Gudhi Library. The Gudhi library + * (Geometric Understanding in Higher Dimensions) is a generic C++ + * library for computational topology. + * + * Author(s): Clement Jamin + * + * Copyright (C) 2016 INRIA Sophia-Antipolis (France) + * + * This program is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * This program is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see <http://www.gnu.org/licenses/>. + */ + +#define BOOST_TEST_DYN_LINK +#define BOOST_TEST_MODULE Spatial_searching - test Spatial_tree_data_structure +#include <boost/test/unit_test.hpp> + +#include <gudhi/Spatial_tree_data_structure.h> + +#include <CGAL/Epick_d.h> +#include <CGAL/Random.h> + +#include <array> +#include <vector> + +BOOST_AUTO_TEST_CASE(test_Spatial_tree_data_structure) +{ + typedef CGAL::Epick_d<CGAL::Dimension_tag<4> > K; + typedef K::FT FT; + typedef K::Point_d Point; + typedef std::vector<Point> Points; + + typedef Gudhi::spatial_searching::Spatial_tree_data_structure< + K, Points> Points_ds; + + CGAL::Random rd; + + Points points; + for (int i = 0 ; i < 500 ; ++i) + points.push_back(Point(std::array<FT,4>({rd.get_double(-1.,1),rd.get_double(-1.,1),rd.get_double(-1.,1),rd.get_double(-1.,1)}))); + + Points_ds points_ds(points); + + std::size_t closest_pt_index = + points_ds.query_ANN(points[10], 1, false).begin()->first; + BOOST_CHECK(closest_pt_index == 10); + + auto kns_range = points_ds.query_ANN(points[20], 10, true); + + FT last_dist = -1.; + for (auto const& nghb : kns_range) + { + BOOST_CHECK(nghb.second > last_dist); + last_dist = nghb.second; + } +} |