GetFEM  5.4.2
bgeot_kdtree.h
Go to the documentation of this file.
1 /* -*- c++ -*- (enables emacs c++ mode) */
2 /*===========================================================================
3 
4  Copyright (C) 2004-2020 Julien Pommier
5 
6  This file is a part of GetFEM
7 
8  GetFEM is free software; you can redistribute it and/or modify it
9  under the terms of the GNU Lesser General Public License as published
10  by the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version along with the GCC Runtime Library
12  Exception either version 3.1 or (at your option) any later version.
13  This program is distributed in the hope that it will be useful, but
14  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  License and GCC Runtime Library Exception for more details.
17  You should have received a copy of the GNU Lesser General Public License
18  along with this program; if not, write to the Free Software Foundation,
19  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
20 
21  As a special exception, you may use this file as it is a part of a free
22  software library without restriction. Specifically, if other files
23  instantiate templates or use macros or inline functions from this file,
24  or you compile this file and link it with other files to produce an
25  executable, this file does not by itself cause the resulting executable
26  to be covered by the GNU Lesser General Public License. This exception
27  does not however invalidate any other reasons why the executable file
28  might be covered by the GNU Lesser General Public License.
29 
30 ===========================================================================*/
31 
32 #ifndef BGEOT_KDTREE_H
33 #define BGEOT_KDTREE_H
34 
35 /** @file bgeot_kdtree.h
36  @author Julien Pommier <Julien.Pommier@insa-toulouse.fr>
37  @date January 2004.
38  @brief Simple implementation of a KD-tree.
39 
40  Basically, a KD-tree is a balanced N-dimensional tree.
41 */
42 #include "bgeot_small_vector.h"
43 
44 namespace bgeot {
45 
46  /* generic node for the kdtree */
47  struct kdtree_elt_base {
48  enum { PTS_PER_LEAF=8 };
49  unsigned n; /* 0 => is a tree node, != 0 => tree leaf storing n points */
50  bool isleaf() const { return (n != 0); }
51  kdtree_elt_base(unsigned n_) : n(n_) {}
52  virtual ~kdtree_elt_base() {}
53  };
54 
55  /// store a point and the associated index for the kdtree.
56  /* std::pair<size_type,base_node> is not ok since it does not
57  have a suitable overloaded swap function ...
58  */
59  struct index_node_pair {
60  size_type i;
61  base_node n;
62  index_node_pair() {}
63  index_node_pair(size_type i_, base_node n_) : i(i_), n(n_) {}
64  void swap(index_node_pair& other) { std::swap(i,other.i); n.swap(other.n);}
65  };
66 
67  /// store a set of points with associated indexes.
68  typedef std::vector<index_node_pair> kdtree_tab_type;
69 
70  /** Balanced tree over a set of points.
71 
72  Once the tree have been built, it is possible to query very
73  quickly for the list of points lying in a given box. Note that
74  this is not a dynamic structure: once you start to call
75  kdtree::points_in_box, you should not use anymore kdtree::add_point.
76 
77  Here is an example of use (which tries to find the mapping between
78  the dof of the mesh_fem and the node numbers of its mesh):
79 
80  @code
81  void dof_to_nodes(const getfem::mesh_fem &mf) {
82  const getfem::getfem_mesh &m = mf.linked_mesh();
83  bgeot::kdtree tree;
84  for (dal::bv_visitor ip(m.points().index()); !ip.finished(); ++ip) {
85  tree.add_point_with_id(m.points()[ip], ip);
86  }
87  bgeot::kdtree_tab_type t;
88  for (size_type d = 0; d < mf.nb_dof(); ++d) {
89  getfem::base_node P(mf.point_of_dof(d)), P2(P);
90  for (unsigned i=0; i < P.size(); ++i) {
91  P[i] -= 1e-12; P2[i]+= 1e-12;
92  }
93  tree.points_in_box(t, P, P2);
94  if (t.size()) {
95  assert(t.size() == 1);
96  cout << " dof " << d << " maps to mesh node " << t[0].i << "\n";
97  }
98  }
99  }
100  @endcode
101  */
102  class kdtree {
103  dim_type N; /* dimension of points */
104  std::unique_ptr<kdtree_elt_base> tree;
105  kdtree_tab_type pts;
106  public:
107  kdtree() : N(0) {}
108 
109  kdtree(const kdtree&) = delete;
110  kdtree &operator = (const kdtree&) = delete;
111 
112  /// reset the tree, remove all points
113  void clear() { clear_tree(); pts = kdtree_tab_type(); N = 0; }
114  void reserve(size_type n) { pts.reserve(n); }
115  /// insert a new point
117  size_type i = pts.size(); add_point_with_id(n,i); return i;
118  }
119  /// insert a new point, with an associated number.
121  if (pts.size() == 0)
122  N = n.size();
123  else
124  GMM_ASSERT2(N == n.size(), "invalid dimension");
125  if (tree) clear_tree();
126  pts.push_back(index_node_pair(i, n));
127  }
128  size_type nb_points() const { return pts.size(); }
129  const kdtree_tab_type &points() const { return pts; }
130  /* fills ipts with the indexes of points in the box
131  [min,max] */
132  void points_in_box(kdtree_tab_type &ipts,
133  const base_node &min,
134  const base_node &max);
135  /* assigns at ipt the index of the nearest neighbor at location
136  pos and returns the square of the distance to this point*/
137  scalar_type nearest_neighbor(index_node_pair &ipt,
138  const base_node &pos);
139  private:
140  typedef std::vector<size_type>::const_iterator ITER;
141  void clear_tree();
142  };
143 }
144 
145 #endif
bgeot::kdtree_tab_type
std::vector< index_node_pair > kdtree_tab_type
store a set of points with associated indexes.
Definition: bgeot_kdtree.h:68
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
bgeot::kdtree::add_point
size_type add_point(const base_node &n)
insert a new point
Definition: bgeot_kdtree.h:116
bgeot::kdtree
Balanced tree over a set of points.
Definition: bgeot_kdtree.h:102
bgeot::kdtree::clear
void clear()
reset the tree, remove all points
Definition: bgeot_kdtree.h:113
bgeot::small_vector< scalar_type >
bgeot_small_vector.h
Small (dim < 8) vectors.
bgeot
Basic Geometric Tools.
Definition: bgeot_convex_ref.cc:27
bgeot::index_node_pair
store a point and the associated index for the kdtree.
Definition: bgeot_kdtree.h:59
bgeot::kdtree::add_point_with_id
void add_point_with_id(const base_node &n, size_type i)
insert a new point, with an associated number.
Definition: bgeot_kdtree.h:120