Merge pull request #6043 from danston/Kd_tree-nice_print-danston

Kd Tree Nice Print

Installation/CHANGES.md

@@ -1,6 +1,12 @@
 Release History
 ===============
 
+
+### [dD Spatial Searching](https://doc.cgal.org/5.5/Manual/packages.html#PkgSpatialSearchingD)
+
+-   Added the member function `write_graphviz()` to the class The Kd_tree` that writes the tree in a stream in the [Graphviz](https://graphviz.org/) format.
+
+
 [Release 5.4](https://github.com/CGAL/cgal/releases/tag/v5.4)
 -----------
 

Spatial_searching/benchmark/Spatial_searching/CMakeLists.txt

@@ -8,29 +8,28 @@ find_package(CGAL REQUIRED COMPONENTS Core)
 
 include(${CGAL_USE_FILE})
 
-find_package(Eigen3 3.1.91) #(requires 3.2.0 or greater)
+find_package(Eigen3 3.1.91) # (requires 3.2.0 or greater)
 include(CGAL_Eigen3_support)
 
 include_directories(BEFORE "include")
 
-create_single_source_cgal_program("Compare_ANN_STANN_CGAL.cpp")
+# create_single_source_cgal_program("Compare_ANN_STANN_CGAL.cpp") # does not compile, missing dependency
 create_single_source_cgal_program("nanoflan.cpp")
 create_single_source_cgal_program("binary.cpp")
 create_single_source_cgal_program("nearest_neighbor_searching_50.cpp")
 create_single_source_cgal_program("nearest_neighbor_searching_inplace_50.cpp")
 create_single_source_cgal_program("Nearest_neighbor_searching.cpp")
 create_single_source_cgal_program("Nearest_neighbor_searching_2D.cpp")
-create_single_source_cgal_program(
-  "Nearest_neighbor_searching_2D_user_defined.cpp")
+create_single_source_cgal_program("Nearest_neighbor_searching_2D_user_defined.cpp")
 create_single_source_cgal_program("Split_data.cpp")
 create_single_source_cgal_program("nn3cgal.cpp")
-create_single_source_cgal_program("nn4cgal.cpp")
+# create_single_source_cgal_program("nn4cgal.cpp") # this file does not exist for some reason
 create_single_source_cgal_program("nn3nanoflan.cpp")
 create_single_source_cgal_program("sizeof.cpp")
-create_single_source_cgal_program("deque.cpp")
+# create_single_source_cgal_program("deque.cpp") # does not compile, lots of errors
 foreach(
   target
-  Compare_ANN_STANN_CGAL
+  # Compare_ANN_STANN_CGAL # see above
   nanoflan
   binary
   nearest_neighbor_searching_50
@@ -40,9 +39,10 @@ foreach(
   Nearest_neighbor_searching_2D_user_defined
   Split_data
   nn3cgal
-  nn4cgal
+  # nn4cgal # see above
   nn3nanoflan
   sizeof
-  deque)
+  # deque # see above
+  )
   target_link_libraries(${target} PUBLIC CGAL::Eigen3_support)
 endforeach()

Spatial_searching/benchmark/Spatial_searching/Distance.h

@@ -34,7 +34,7 @@ struct Distance {
   }
 
   double new_distance(double& dist, double old_off, double new_off,
-                      int cutting_dimension)  const {
+                      int /* cutting_dimension */)  const {
     return dist + new_off*new_off - old_off*old_off;
   }
 

Spatial_searching/benchmark/Spatial_searching/binary.cpp

@@ -9,7 +9,7 @@
 
 typedef CGAL::Simple_cartesian<double>::Point_3 Point_3;
 
-int main(int argc, char* argv[])
+int main(int /* argc */, char* argv[])
 {
 #if 0
   std::cerr << "ASCII to binary\n";

Spatial_searching/benchmark/Spatial_searching/include/nanoflann.hpp

@@ -485,6 +485,7 @@ namespace nanoflann
                  */
                 PooledAllocator(const size_t blocksize_ = BLOCKSIZE) : blocksize(blocksize_) {
                         internal_init();
+                        blocksize = blocksize_;
                 }
 
                 /**
@@ -963,7 +964,7 @@ namespace nanoflann
                  *  \sa radiusSearch, findNeighbors
                  * \note nChecks_IGNORED is ignored but kept for compatibility with the original FLANN interface.
                  */
-                inline void knnSearch(const ElementType *query_point, const size_t num_closest, IndexType *out_indices, DistanceType *out_distances_sq,const int nChecks_IGNORED = 10) const
+                inline void knnSearch(const ElementType *query_point, const size_t num_closest, IndexType *out_indices, DistanceType *out_distances_sq,const int /* nChecks_IGNORED = 10 */) const
                 {
                         nanoflann::KNNResultSet<DistanceType,IndexType> resultSet(num_closest);
                         resultSet.init(out_indices, out_distances_sq);
@@ -1412,7 +1413,7 @@ namespace nanoflann
                 index_t* index; //! The kd-tree index for the user to call its methods as usual with any other FLANN index.
 
                 /// Constructor: takes a const ref to the matrix object with the data points
-                KDTreeEigenMatrixAdaptor(const int dimensionality, const MatrixType &mat, const int leaf_max_size = 10) : m_data_matrix(mat)
+                KDTreeEigenMatrixAdaptor(const int /* dimensionality */, const MatrixType &mat, const int leaf_max_size = 10) : m_data_matrix(mat)
                 {
                         const size_t dims = mat.cols();
                         if (DIM>0 && static_cast<int>(dims)!=DIM)
@@ -1436,7 +1437,7 @@ namespace nanoflann
                   *  The user can also call index->... methods as desired.
                   * \note nChecks_IGNORED is ignored but kept for compatibility with the original FLANN interface.
                   */
-                inline void query(const num_t *query_point, const size_t num_closest, IndexType *out_indices, num_t *out_distances_sq, const int nChecks_IGNORED = 10) const
+                inline void query(const num_t *query_point, const size_t num_closest, IndexType *out_indices, num_t *out_distances_sq, const int /* nChecks_IGNORED = 10 */) const
                 {
                         nanoflann::KNNResultSet<typename MatrixType::Scalar,IndexType> resultSet(num_closest);
                         resultSet.init(out_indices, out_distances_sq);
@@ -1478,7 +1479,7 @@ namespace nanoflann
                 //   Return true if the BBOX was already computed by the class and returned in "bb" so it can be avoided to redo it again.
                 //   Look at bb.size() to find out the expected dimensionality (e.g. 2 or 3 for point clouds)
                 template <class BBOX>
-                bool kdtree_get_bbox(BBOX &bb) const {
+                bool kdtree_get_bbox(BBOX & /* bb */) const {
                         return false;
                 }
 

Spatial_searching/benchmark/Spatial_searching/nanoflan.cpp

@@ -39,7 +39,7 @@ struct PointCloud
         inline size_t kdtree_get_point_count() const { return pts.size(); }
 
         // Returns the distance between the vector "p1[0:size-1]" and the data point with index "idx_p2" stored in the class:
-        inline T kdtree_distance(const T *p1, const size_t idx_p2,size_t size) const
+        inline T kdtree_distance(const T *p1, const size_t idx_p2,size_t /* size */) const
         {
                 const T d0=p1[0]-pts[idx_p2].x;
                 const T d1=p1[1]-pts[idx_p2].y;
@@ -61,7 +61,7 @@ struct PointCloud
         //   Return true if the BBOX was already computed by the class and returned in "bb" so it can be avoided to redo it again.
         //   Look at bb.size() to find out the expected dimensionality (e.g. 2 or 3 for point clouds)
         template <class BBOX>
-        bool kdtree_get_bbox(BBOX &bb) const { return false; }
+        bool kdtree_get_bbox(BBOX & /* bb */) const { return false; }
 
 };
 
@@ -76,7 +76,7 @@ void generateRandomPointCloud(PointCloud<T> &point, istream& is)
 }
 
 template <typename num_t>
-void kdtree_demo(const size_t N)
+void kdtree_demo(const size_t /* N */)
 {
         PointCloud<num_t> cloud;
 
@@ -131,14 +131,14 @@ void kdtree_demo(const size_t N)
           //std::cout << "knnSearch(nn="<<num_results<<"): \n";
           //std::cout << "ret_index=" << ret_index << " out_dist_sqr=" << out_dist_sqr << endl;
           //std::cout << cloud.pts[ret_index] << std::endl;
-          for (int k=0; k<num_results; ++k)
+          for (size_t k=0; k<num_results; ++k)
             sum += cloud.pts[ret_index[k]].x;
         }
         timer.stop();
         std::cout << sum << " done in " << timer.time() << " sec."<< std::endl;
 }
 
-int main(int argc, char** argv)
+int main(int /* argc */, char** /* argv */)
 {
   //kdtree_demo<float>(100000);
         kdtree_demo<double>(100000);

Spatial_searching/benchmark/Spatial_searching/nearest_neighbor_searching_50.cpp

@@ -15,7 +15,7 @@ typedef Neighbor_search::Tree Tree;
 
 
 int main() {
-  const unsigned int N = 50;
+  // const unsigned int N = 50;
   CGAL::Timer timer;
   timer.start();
 
@@ -51,7 +51,7 @@ int main() {
   double d = 0;
   timer.reset();
   timer.start();
-  for(int i = 0; i < queries.size(); i++){
+  for(std::size_t i = 0; i < queries.size(); i++){
     Neighbor_search search(tree, queries[i], 50, 0);
 
     // report the N nearest neighbors and their distance

Spatial_searching/benchmark/Spatial_searching/nearest_neighbor_searching_inplace_50.cpp

@@ -41,7 +41,7 @@ typedef K_neighbor_search::Distance Distance;
 
 
 int main() {
-  const unsigned int N = 50;
+  // const unsigned int N = 50;
   CGAL::Timer timer;
   timer.start();
 
@@ -82,7 +82,7 @@ int main() {
   double d = 0;
   timer.reset();
   timer.start();
-  for(int i = 0; i < queries.size(); i++){
+  for(std::size_t i = 0; i < queries.size(); i++){
     K_neighbor_search search(tree, queries[i], 50, 0, true, tr_dist);
 
     // report the N nearest neighbors and their distance

Spatial_searching/benchmark/Spatial_searching/nn3nanoflan.cpp

@@ -10,6 +10,7 @@
 #include <fstream>
 #include <vector>
 #include <CGAL/Memory_sizer.h>
+#include <CGAL/IO/io.h>
 
 using namespace std;
 using namespace nanoflann;
@@ -42,7 +43,7 @@ struct PointCloud
         inline size_t kdtree_get_point_count() const { return pts.size(); }
 
         // Returns the distance between the vector "p1[0:size-1]" and the data point with index "idx_p2" stored in the class:
-        inline T kdtree_distance(const T *p1, const size_t idx_p2,size_t size) const
+        inline T kdtree_distance(const T *p1, const size_t idx_p2,size_t /* size */) const
         {
                 const T d0=p1[0]-pts[idx_p2].x;
                 const T d1=p1[1]-pts[idx_p2].y;
@@ -64,7 +65,7 @@ struct PointCloud
         //   Return true if the BBOX was already computed by the class and returned in "bb" so it can be avoided to redo it again.
         //   Look at bb.size() to find out the expected dimensionality (e.g. 2 or 3 for point clouds)
         template <class BBOX>
-        bool kdtree_get_bbox(BBOX &bb) const { return false; }
+        bool kdtree_get_bbox(BBOX & /* bb */) const { return false; }
 
 };
 
@@ -162,7 +163,7 @@ void kdtree_demo(int argc, char** argv)
             index.findNeighbors(resultSet, &query_pt[0], nanoflann::SearchParams(10,0));
             timer.stop();
 
-            for (int k=0; k<num_results; ++k){
+            for (size_t k=0; k<num_results; ++k){
               if(dump)
                 std::cerr <<cloud.pts[ret_index[k]] << std::endl;
             }

Spatial_searching/doc/Spatial_searching/CGAL/Kd_tree.h

@@ -223,7 +223,7 @@ Returns the number of points that are stored in the tree.
 size_type size() const;
 
 /*!
-return the instance of the traits used to construct the tree.
+Returns the instance of the traits used to construct the tree.
 */
 Traits traits() const;
 
@@ -248,6 +248,12 @@ Inserts statistics of the tree into the output stream `s`.
 */
 std::ostream& statistics(std::ostream& s) const;
 
+/*!
+Inserts the tree in the <a href="https://graphviz.org/">`Graphviz`</a> format
+into the output stream `s`.
+*/
+std::ostream& write_graphviz(std::ostream& s) const;
+
 /// @}
 
 }; /* end Kd_tree */

Spatial_searching/include/CGAL/Kd_tree.h

@@ -21,6 +21,9 @@
 #include <CGAL/basic.h>
 #include <CGAL/assertions.h>
 #include <vector>
+#include <string>
+#include <unordered_map>
+#include <ostream>
 
 #include <CGAL/algorithm.h>
 #include <CGAL/Kd_tree_node.h>
@@ -377,6 +380,39 @@ public:
     return dim_;
   }
 
+  std::ostream&
+  write_graphviz(std::ostream& s) const
+  {
+    int counter = -1;
+    std::unordered_map<const Node*, int> node_to_index;
+    tree_root->get_indices(counter, node_to_index);
+
+    const auto node_name = [&](const Node* node) {
+      const int index = node_to_index.at(node);
+      std::string node_name = "default_name";
+      if (node->is_leaf()) { // leaf node
+        node_name = "L" + std::to_string(index);
+      } else {
+        if (index == 0) { // root node
+          node_name = "R" + std::to_string(index);
+        } else { // internal node
+          node_name = "N" + std::to_string(index);
+        }
+      }
+      CGAL_assertion(node_name != "default_name");
+      return node_name;
+    };
+
+    s << "graph G" << std::endl;
+    s << "{" << std::endl << std::endl;
+    s << "label=\"Graph G. Num leaves: " << tree_root->num_nodes() << ". ";
+    s << "Num items: " << tree_root->num_items() << ".\"" << std::endl;
+    s << node_name(tree_root) + " ;";
+    tree_root->print(s, node_name);
+    s << std::endl << "}" << std::endl << std::endl;
+    return s;
+  }
+
 private:
   //any call to this function is for the moment not threadsafe
   void const_build() const {

Spatial_searching/include/CGAL/Kd_tree_node.h

@@ -15,7 +15,8 @@
 
 #include <CGAL/license/Spatial_searching.h>
 
-
+#include <unordered_map>
+#include <ostream>
 
 #include <CGAL/Splitters.h>
 #include <CGAL/Compact_container.h>
@@ -58,6 +59,55 @@ namespace CGAL {
 
     bool is_leaf() const { return leaf; }
 
+    void
+    get_indices(int& index, std::unordered_map<const Kd_tree_node*, int>& node_to_index) const
+    {
+      if (is_leaf()) {
+        Leaf_node_const_handle node =
+          static_cast<Leaf_node_const_handle>(this);
+        ++index;
+        node_to_index[node] = index;
+      } else {
+        Internal_node_const_handle node =
+          static_cast<Internal_node_const_handle>(this);
+        ++index;
+        node_to_index[node] = index;
+        node->lower()->get_indices(index, node_to_index);
+        node->upper()->get_indices(index, node_to_index);
+      }
+    }
+
+    template<typename Node_name>
+    void
+    print(std::ostream& s, const Node_name& node_name) const
+    {
+      if (is_leaf()) { // draw leaf nodes
+
+        Leaf_node_const_handle node =
+          static_cast<Leaf_node_const_handle>(this);
+
+        s << std::endl;
+        if (node->size() > 0) {
+          s << node_name(node) << " [label=\"" << node_name(node) << ", Size: "
+          << node->size() << "\"] ;" << std::endl;
+        } else {
+          CGAL_assertion_msg(false, "ERROR: NODE SIZE IS ZERO!");
+        }
+
+      } else { // draw internal nodes
+
+        Internal_node_const_handle node =
+          static_cast<Internal_node_const_handle>(this);
+
+        s << std::endl;
+        s << node_name(node) << " [label=\"" << node_name(node) << "\"] ;" << std::endl;
+        s << node_name(node) << " -- " << node_name(node->lower()) << " ;";
+        node->lower()->print(s, node_name);
+        s << node_name(node) << " -- " << node_name(node->upper()) << " ;";
+        node->upper()->print(s, node_name);
+      }
+    }
+
     std::size_t
     num_items() const
     {

Spatial_searching/test/Spatial_searching/data/balanced.xyz

@@ -0,0 +1,25 @@
+2 3 3
+5 4 2
+9 6 7
+4 7 9
+8 1 5
+7 2 6
+9 4 1
+8 4 2
+9 7 8
+6 3 1
+3 4 5
+1 6 8
+9 5 3
+2 1 3
+8 7 6
+5 4 2
+6 3 1
+8 7 6
+9 6 7
+2 1 3
+7 2 6
+4 7 9
+1 6 8
+3 4 5
+9 4 1

Spatial_searching/test/Spatial_searching/test_print.cpp

@@ -0,0 +1,67 @@
+#include <cmath>
+#include <vector>
+#include <fstream>
+#include <iostream>
+
+#include <CGAL/Real_timer.h>
+#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+#include <CGAL/Search_traits_3.h>
+#include <CGAL/Splitters.h>
+#include <CGAL/Kd_tree.h>
+
+using Kernel  = CGAL::Exact_predicates_inexact_constructions_kernel;
+using Point_3 = typename Kernel::Point_3;
+
+using Traits   = CGAL::Search_traits_3<Kernel>;
+using Splitter = CGAL::Sliding_midpoint<Traits>;
+using Kd_tree  = CGAL::Kd_tree<Traits, Splitter>;
+
+using Timer = CGAL::Real_timer;
+
+void test_print(const std::string filename) {
+
+  std::cout << std::endl;
+  std::cout << "- testing " << filename << " ... " << std::endl;
+
+  std::vector<Point_3> points;
+  std::ifstream in(filename);
+  CGAL::IO::set_ascii_mode(in);
+
+  Point_3 p; std::size_t count = 0;
+  while (in >> p) {
+    points.push_back(p);
+    ++count;
+  }
+  assert(points.size() > 0);
+
+  std::cout << "* num points: " << points.size() << std::endl;
+  std::cout << "* building the tree ... " << std::endl;
+
+  // Building the tree.
+  Splitter splitter(5); // bucket size = 5
+  Kd_tree tree(points.begin(), points.end(), splitter);
+
+  Timer timer;
+  timer.start();
+  tree.build();
+  timer.stop();
+  std::cout << "* building done in " << timer.time() << " sec." << std::endl;
+
+  tree.statistics(std::cout);
+  // Use this command to print in png:
+  // dot -Tpng tree.graphviz > tree.png
+  std::ofstream outfile("tree.graphviz");
+  tree.write_graphviz(outfile);
+
+  assert(tree.root()->num_items() == points.size());
+  assert(tree.root()->num_nodes() == 8);
+  assert(tree.root()->depth() == 5);
+  std::cout << std::endl;
+}
+
+int main() {
+
+  std::cout.precision(20);
+  test_print("data/balanced.xyz");
+  return EXIT_SUCCESS;
+}