]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Minor code updates to the kdtree implementation. 4937/head
authorWolfgang Bangerth <bangerth@colostate.edu>
Wed, 23 Aug 2017 02:41:38 +0000 (20:41 -0600)
committerWolfgang Bangerth <bangerth@colostate.edu>
Thu, 24 Aug 2017 13:41:45 +0000 (07:41 -0600)
include/deal.II/numerics/kdtree.h
source/numerics/kdtree.cc

index d87cde2156ca1daf11d432a04eb1e6486c55a82a..ea8a15be66cc403c36d920ad469973f9a4ab3211 100644 (file)
@@ -34,30 +34,31 @@ DEAL_II_NAMESPACE_OPEN
 /**
  * A wrapper for the nanoflann library, used to compute the distance from a
  * collection of points, and to efficiently return nearest neighbors to a
- * target point. This function uses nanoflann to efficiently partition the
- * space in a k-dimensional tree. The cost of each query is then roughly of
- * order log(n), where n is the number of points stored in this class.
+ * target point. This class uses nanoflann to efficiently partition the
+ * space in a $k$-dimensional tree. The cost of each query is then roughly of
+ * order $\log(n)$, where $n$ is the number of points stored in this class.
  *
  * The wrapper provides methods that give access to some of the functionalities
- * of the nanoflann library, like searching the n nearest neighbors, or
+ * of the nanoflann library, like searching the $p$ nearest neighbors of
+ * a given point, or
  * searching the points that fall within a radius of a target point.
  *
  * @quotation
  * From wikipedia (https://en.wikipedia.org/wiki/K-d_tree):
  *
- * A k-d tree is a binary tree in which every node is a k-dimensional point.
+ * A k-d tree is a binary tree in which every node is a $k$-dimensional point.
  * Every non-leaf node can be thought of as implicitly generating a splitting
  * hyperplane that divides the space into two parts, known as half-spaces.
  * Points to the left of this hyperplane are represented by the left subtree of
  * that node and points right of the hyperplane are represented by the right
  * subtree. The hyperplane direction is chosen in the following way: every node
- * in the tree is associated with one of the k-dimensions, with the hyperplane
+ * in the tree is associated with one of the $k$-dimensions, with the hyperplane
  * perpendicular to that dimension's axis. So, for example, if for a particular
  * split the "x" axis is chosen, all points in the subtree with a smaller "x"
  * value than the node will appear in the left subtree and all points with
  * larger "x" value will be in the right subtree. In such a case, the
- * hyperplane would be set by the x-value of the point, and its normal would be
- * the unit x-axis.
+ * hyperplane would be set by the $x$-value of the point, and its normal would be
+ * the unit $x$-axis.
  * @endquotation
  *
  * @author Luca Heltai, 2017.
@@ -67,28 +68,33 @@ class KDTree
 {
 public:
   /**
-   * The max leaf parameter is used to decide how many points per leaf
+   * Constructor.
+   *
+   * @param[in] max_leaf_size A number denoting how many points per leaf
    * are used in the kdtree algorithm.
    *
-   * If the points are not passed to this constructor, then you have
+   * @param[in] pts A vector of points that are to be represented by
+   * the current object. If no points are passed to this constructor
+   * (or if the default value of the argument is used), then you have
    * to pass them later to this object by calling the set_points()
    * method.
    *
    * Access to any of the methods without first passing a reference to
    * a vector of points will result in an exception. Only a reference
    * to the points is stored, so you should make sure that the life of
-   * the the vector you pass is longer than the life of this class, or
-   * you'll get undefined behaviour.
+   * the vector you pass is longer than the life of this class, or
+   * you will get undefined behaviour.
    *
-   * If you update the vector of points in someway, remember to call
-   * again the set_points() method. The tree and the index are
-   * constructed only once, when you pass the points (either at
+   * @warning If you change the contents of the vector of points that you
+   * passed either to the constructor or to set_points(), remember to call
+   * the set_points() method again. The tree and the index are
+   * constructed only once when you pass the points (either at
    * construction time, or when you call set_points()). If you update
-   * your points, and do not call again set_points(), your results
+   * your points, and do not call set_points() again, then all following results
    * will likely be wrong.
    */
-  KDTree(const unsigned int &max_leaf_size=10,
-         const std::vector<Point<dim> > &pts=std::vector<Point<dim> >());
+  KDTree(const unsigned int             &max_leaf_size = 10,
+         const std::vector<Point<dim> > &pts           = std::vector<Point<dim> >());
 
 
   /**
@@ -115,7 +121,7 @@ public:
      * The constructor needs the vector of points from which we want to build
      * the tree.
      */
-    PointCloudAdaptor(const std::vector<Point<dim> > &_points);
+    PointCloudAdaptor (const std::vector<Point<dim> > &_points);
 
 
     /**
@@ -127,13 +133,16 @@ public:
     /**
      * Return the L2 distance between points
      */
-    coord_t kdtree_distance(const coord_t *p1, const size_t idx_p2,size_t size) const;
+    coord_t kdtree_distance (const coord_t *p1,
+                             const size_t idx_p2,
+                             const size_t size) const;
 
 
     /**
      * Return the d-th component of the idx-th point in the class.
      */
-    coord_t kdtree_get_pt(const size_t idx, int d) const;
+    coord_t kdtree_get_pt (const size_t idx,
+                           const int d) const;
 
 
     /**
@@ -144,15 +153,17 @@ public:
      * expected dimensionality (e.g. 2 or 3 for point clouds).
      */
     template <class BBOX>
-    bool kdtree_get_bbox(BBOX &) const;
+    bool kdtree_get_bbox (BBOX &) const;
   };
 
 
   /**
    * A typedef for the actual KDTree object.
    */
-  typedef typename nanoflann::KDTreeSingleIndexAdaptor<nanoflann::L2_Simple_Adaptor<double, PointCloudAdaptor> ,
-          PointCloudAdaptor, dim, unsigned int>  NanoFlannKDTree;
+  typedef
+  typename nanoflann::KDTreeSingleIndexAdaptor<nanoflann::L2_Simple_Adaptor<double, PointCloudAdaptor> ,
+           PointCloudAdaptor, dim, unsigned int>
+           NanoFlannKDTree;
 
 
   /**
@@ -163,7 +174,7 @@ public:
    * the get_closest_points() methods.
    *
    * Notice that the constructor calls this method internally if you
-   * pass it a non empty vector of points.
+   * pass it a non-empty vector of points.
    *
    * Whenever your points change, you should call this method again,
    * since this is the method responsible for building the index and
@@ -171,56 +182,57 @@ public:
    * don't call again this method, any function you call later will
    * happily return wrong values without you noticing.
    *
-   * @param pts: a collection of points
+   * @param[in] pts A collection of points
    */
-  void set_points(const std::vector<Point<dim> > &pts);
+  void set_points (const std::vector<Point<dim> > &pts);
 
 
   /**
-   * A const accessor to the underlying points.
+   * A const accessor to the @p i'th one among the underlying points.
    */
-  const Point<dim> &operator[](unsigned int i) const;
+  const Point<dim> &operator[] (const unsigned int i) const;
 
 
   /**
-   * The size of the vector stored by this class.
+   * The number of points currently stored by this class.
    */
   unsigned int size() const;
 
 
   /**
-   * Fill a vector with the indices and the distance of the points
+   * Fill and return a vector with the indices and the distance of the points
    * that are at distance less than or equal to the given radius from
-   * the target point. Consider preallocating the size of the return
-   * vector if you have a wild guess of how many should be there.
+   * the target point.
    *
-   * @param[in] point: the target point
-   * @param[in] radius: the radius of the ball
-   * @param[in] sorted: sort the output results in ascending order with respect to distances
+   * @param[in] point The target point
+   * @param[in] radius The radius of the ball
+   * @param[in] sorted Sort the output results in ascending order with respect to distances
    *
-   * @return vector of indices and distances of the matching points
+   * @return vector of indices and distances of the matching points
    */
-  std::vector<std::pair<unsigned int, double> > get_points_within_ball(const Point<dim> &target,
-      const double &radius,
-      bool sorted=false) const;
+  std::vector<std::pair<unsigned int, double> >
+  get_points_within_ball (const Point<dim> &target,
+                          const double     &radius,
+                          const bool        sorted=false) const;
 
   /**
-   * Fill a vectors with the indices and distances of the closest @p n_points
+   * Fill and return a vector with the indices and distances of the closest @p n_points
    * points to the given target point.
    *
-   * @param[in] target: the target point
-   * @param[in] n_points: the number of requested points
+   * @param[in] target The target point
+   * @param[in] n_points The number of requested points
    *
-   * @return a vector of pairs of indices and distances of the matching points
+   * @return A vector of pairs of indices and distances of the matching points
    */
-  std::vector<std::pair<unsigned int, double> >  get_closest_points(const Point<dim> &target,
-      const unsigned int n_points) const;
+  std::vector<std::pair<unsigned int, double> >
+  get_closest_points (const Point<dim>  &target,
+                      const unsigned int n_points) const;
 
 private:
   /**
-   * Max number of points per leaf.
+   * Max number of points per leaf as set in the constructor.
    */
-  unsigned int max_leaf_size;
+  const unsigned int max_leaf_size;
 
 
   /**
@@ -253,7 +265,7 @@ unsigned int KDTree<dim>::size() const
 
 template <int dim>
 inline const Point<dim> &
-KDTree<dim>::operator[](unsigned int i) const
+KDTree<dim>::operator[] (const unsigned int i) const
 {
   AssertIndexRange(i, size());
   return adaptor->points[i];
@@ -263,7 +275,8 @@ KDTree<dim>::operator[](unsigned int i) const
 
 template <int dim>
 KDTree<dim>::PointCloudAdaptor::PointCloudAdaptor(const std::vector<Point<dim> > &_points)
-  : points(_points)
+  :
+  points(_points)
 {}
 
 
@@ -290,7 +303,7 @@ KDTree<dim>::PointCloudAdaptor::kdtree_get_pt(const size_t idx, int d) const
 template <int dim>
 template <class BBOX>
 inline bool
-KDTree<dim>::PointCloudAdaptor::kdtree_get_bbox(BBOX &) const
+KDTree<dim>::PointCloudAdaptor::kdtree_get_bbox (BBOX &) const
 {
   return false;
 }
@@ -299,12 +312,14 @@ KDTree<dim>::PointCloudAdaptor::kdtree_get_bbox(BBOX &) const
 
 template <int dim>
 inline double
-KDTree<dim>::PointCloudAdaptor::kdtree_distance(const double *p1, const size_t idx_p2,size_t size) const
+KDTree<dim>::PointCloudAdaptor::kdtree_distance (const double *p1,
+                                                 const size_t idx_p2,
+                                                 const size_t size) const
 {
   AssertDimension(size, dim);
   double res=0.0;
   for (size_t d=0; d<size; ++d)
-    res += (p1[d]-points[idx_p2][d])*(p1[d]-points[idx_p2][d]);
+    res += (p1[d]-points[idx_p2][d]) * (p1[d]-points[idx_p2][d]);
   return std::sqrt(res);
 }
 #endif
index ee6228bc720c61d881e6ca5b9741564cb965d776..8cdab9311b1b8904d292da253e22d882e51ff09c 100644 (file)
@@ -8,21 +8,23 @@ DEAL_II_NAMESPACE_OPEN
 
 
 template <int dim>
-KDTree<dim>::KDTree(const unsigned int &max_leaf_size,
+KDTree<dim>::KDTree(const unsigned int             &max_leaf_size,
                     const std::vector<Point<dim> > &pts)
-  : max_leaf_size(max_leaf_size)
+  :
+  max_leaf_size(max_leaf_size)
 {
   if (pts.size() > 0)
     set_points(pts);
 }
 
 
+
 template <int dim>
-std::vector<std::pair<unsigned int, double> > KDTree<dim>::get_points_within_ball(const Point<dim> &center,
-    const double &radius,
-    bool sorted) const
+std::vector<std::pair<unsigned int, double> >
+KDTree<dim>::get_points_within_ball(const Point<dim> &center,
+                                    const double &radius,
+                                    bool sorted) const
 {
-  std::vector<std::pair<unsigned int, double> > matches;
   Assert(adaptor, ExcNotInitialized());
   Assert(kdtree, ExcInternalError());
 
@@ -31,32 +33,47 @@ std::vector<std::pair<unsigned int, double> > KDTree<dim>::get_points_within_bal
 
   nanoflann::SearchParams params;
   params.sorted = sorted;
+
+  std::vector<std::pair<unsigned int, double> > matches;
   kdtree->radiusSearch(&center[0], radius, matches, params);
+
   return matches;
 }
 
+
+
 template <int dim>
-std::vector<std::pair<unsigned int, double> > KDTree<dim>::get_closest_points(const Point<dim> &target,
-    const unsigned int n_points) const
+std::vector<std::pair<unsigned int, double> >
+KDTree<dim>::get_closest_points(const Point<dim> &target,
+                                const unsigned int n_points) const
 {
   Assert(adaptor, ExcNotInitialized());
   Assert(kdtree, ExcInternalError());
+
+  // get the information out of nanoflann
   std::vector<unsigned int> indices(n_points);
   std::vector<double> distances(n_points);
-  std::vector<std::pair<unsigned int, double> > matches(n_points);
 
   kdtree->knnSearch(&target[0], n_points, &indices[0], &distances[0]);
+
+  // convert it to the format we want to return
+  std::vector<std::pair<unsigned int, double> > matches(n_points);
   for (unsigned int i=0; i<n_points; ++i)
     matches[i] = std::make_pair(indices[i], distances[i]);
+
   return matches;
 }
 
+
+
 template <int dim>
 void KDTree<dim>::set_points(const std::vector<Point<dim> > &pts)
 {
   Assert(pts.size() > 0, ExcMessage("Expecting a non zero set of points."));
   adaptor = std_cxx14::make_unique<PointCloudAdaptor>(pts);
-  kdtree = std_cxx14::make_unique<NanoFlannKDTree>(dim, *adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(max_leaf_size));
+  kdtree = std_cxx14::make_unique<NanoFlannKDTree>(dim,
+                                                   *adaptor,
+                                                   nanoflann::KDTreeSingleIndexAdaptorParams(max_leaf_size));
   kdtree->buildIndex();
 }
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.