Geant4  10.01.p02
G4KDTree.hh
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4KDTree.hh 85244 2014-10-27 08:24:13Z gcosmo $
27 //
28 // Author: Mathieu Karamitros, kara@cenbg.in2p3.fr
29 
30 // The code is developed in the framework of the ESA AO7146
31 //
32 // We would be very happy hearing from you, send us your feedback! :)
33 //
34 // In order for Geant4-DNA to be maintained and still open-source,
35 // article citations are crucial.
36 // If you use Geant4-DNA chemistry and you publish papers about your software,
37 // in addition to the general paper on Geant4-DNA:
38 //
39 // Int. J. Model. Simul. Sci. Comput. 1 (2010) 157–178
40 //
41 // we would be very happy if you could please also cite the following
42 // reference papers on chemistry:
43 //
44 // J. Comput. Phys. 274 (2014) 841-882
45 // Prog. Nucl. Sci. Tec. 2 (2011) 503-508
46 
47 #ifndef G4KDTREE_HH
48 #define G4KDTREE_HH 1
49 
50 #include <vector>
51 #include "G4KDNode.hh"
52 #include "G4KDTreeResult.hh"
53 
54 class G4KDMap;
55 template<typename PointT>
56  class G4KDNode;
57 
58 //__________________________________
59 // Methods to act on kdnode
60 // Methods defined in G4KDNode.cc :
62 void Free(G4KDNode_Base*&);
63 //void* GetData(G4KDNode*);
64 //const double* GetNodePosition(G4KDNode_Base*);
65 //__________________________________
66 
72 class G4KDTree
73 {
74  friend class G4KDNode_Base;
75 public:
76  G4KDTree(size_t dim = 3);
77  ~G4KDTree();
78  void Clear();
79 
80  void Print(std::ostream& out = G4cout) const;
81  void Build();
83  {
85  if (fNbActiveNodes <= 0) Clear();
86  }
87 
88  size_t GetDim() const
89  {
90  return fDim;
91  }
92  int GetNbNodes() const
93  {
94  return fNbNodes;
95  }
97  {
98  return fRoot;
99  }
100 
101  template<typename PointT>
102  G4KDNode_Base* InsertMap(PointT* pos);
103 
104  // Insert and attache the data to a node at the specified position
105  // In return, it gives you the corresponding node
106  template<typename PointT> G4KDNode_Base* Insert(PointT* pos); // 3D
107 
108  /* Find one of the nearest nodes from the specified point.
109  *
110  * This function returns a pointer to a result set with at most one element.
111  */
112  template<typename Position> G4KDTreeResultHandle Nearest(const Position& pos);
114 
115  /* Find any nearest nodes from the specified point within a range.
116  *
117  * This function returns a pointer to a result set, which can be manipulated
118  * by the G4KDTreeResult.
119  * The returned pointer can be null as an indication of an error. Otherwise
120  * a valid result set is always returned which may contain 0 or more elements.
121  */
122  template<typename Position>
123  G4KDTreeResultHandle NearestInRange(const Position& pos,
124  const double& range);
125  G4KDTreeResultHandle NearestInRange(G4KDNode_Base* node, const double& range);
126 
127  void *operator new(size_t);
128  void operator delete(void *);
129 
130 protected:
131 
132  //______________________________________________________________________
133  class HyperRect
134  {
135  public:
136  HyperRect(size_t dim)
137  {
138  fDim = dim;
139  fMin = new double[fDim];
140  fMax = new double[fDim];
141  }
142 
143  template<typename Position>
144  void SetMinMax(const Position& min, const Position& max)
145  {
146  for (size_t i = 0; i < fDim; i++)
147  {
148  fMin[i] = min[i];
149  fMax[i] = max[i];
150  }
151  }
152 
154  {
155  delete[] fMin;
156  delete[] fMax;
157  }
158 
159  HyperRect(const HyperRect& rect)
160  {
161  fDim = rect.fDim;
162  fMin = new double[fDim];
163  fMax = new double[fDim];
164 
165  for (size_t i = 0; i < fDim; i++)
166  {
167  fMin[i] = rect.fMin[i];
168  fMax[i] = rect.fMax[i];
169  }
170  }
171 
172  template<typename Position>
173  void Extend(const Position& pos)
174  {
175  for (size_t i = 0; i < fDim; i++)
176  {
177  if (pos[i] < fMin[i])
178  {
179  fMin[i] = pos[i];
180  }
181  if (pos[i] > fMax[i])
182  {
183  fMax[i] = pos[i];
184  }
185  }
186  }
187 
188  template<typename Position>
189  bool CompareDistSqr(const Position& pos, const double* bestmatch)
190  {
191  double result = 0;
192 
193  for (size_t i = 0; i < fDim; i++)
194  {
195  if (pos[i] < fMin[i])
196  {
197  result += sqr(fMin[i] - pos[i]);
198  }
199  else if (pos[i] > fMax[i])
200  {
201  result += sqr(fMax[i] - pos[i]);
202  }
203 
204  if (result >= *bestmatch) return false;
205  }
206 
207  return true;
208  }
209 
210  size_t GetDim()
211  {
212  return fDim;
213  }
214  double* GetMin()
215  {
216  return fMin;
217  }
218  double* GetMax()
219  {
220  return fMax;
221  }
222 
223  protected:
224  size_t fDim;
225  double *fMin, *fMax; /* minimum/maximum coords */
226 
227  private:
228  // should not be used
230  {
231  if (this == &rhs) return *this;
232  return *this;
233  }
234  };
235 
236 protected:
237  void __InsertMap(G4KDNode_Base *node);
238  void __Clear_Rec(G4KDNode_Base *node);
239 
240  template<typename Position>
242  const Position& pos,
243  const double& range_sq,
244  const double& range,
245  G4KDTreeResult& list,
246  int ordered,
247  G4KDNode_Base *source_node = 0);
248 
249  template<typename Position>
251  const Position& pos,
252  G4KDNode_Base *&result,
253  double *result_dist_sq,
254  HyperRect* fRect);
255 
256  template<typename Position>
257  void __NearestToNode(G4KDNode_Base *source_node,
258  G4KDNode_Base *node,
259  const Position& pos,
260  std::vector<G4KDNode_Base*>& result,
261  double *result_dist_sq,
262  HyperRect* fRect,
263  int& nbresult);
264 
265 protected:
268  size_t fDim;
269  int fNbNodes;
272 
274 };
275 
276 #include "G4KDTree.icc"
277 
278 #endif // G4KDTREE_HH
HyperRect * fRect
Definition: G4KDTree.hh:266
void Clear()
Definition: G4KDTree.cc:90
G4KDNode_Base * GetRoot()
Definition: G4KDTree.hh:96
void Print(std::ostream &out=G4cout) const
Definition: G4KDTree.cc:85
void SetMinMax(const Position &min, const Position &max)
Definition: G4KDTree.hh:144
void __NearestToNode(G4KDNode_Base *source_node, G4KDNode_Base *node, const Position &pos, std::vector< G4KDNode_Base * > &result, double *result_dist_sq, HyperRect *fRect, int &nbresult)
size_t fDim
Definition: G4KDTree.hh:268
int GetNbNodes() const
Definition: G4KDTree.hh:92
G4KDTree is used by the ITManager to locate the neareast neighbours.
Definition: G4KDTree.hh:72
int fNbNodes
Definition: G4KDTree.hh:269
G4KDTree(size_t dim=3)
Definition: G4KDTree.cc:50
#define G4ThreadLocalStatic
Definition: tls.hh:88
G4KDTreeResultHandle Nearest(const Position &pos)
HyperRect(const HyperRect &rect)
Definition: G4KDTree.hh:159
int fNbActiveNodes
Definition: G4KDTree.hh:270
G4GLOB_DLL std::ostream G4cout
G4KDNode_Base * fRoot
Definition: G4KDTree.hh:267
bool CompareDistSqr(const Position &pos, const double *bestmatch)
Definition: G4KDTree.hh:189
G4KDMap * fKDMap
Definition: G4KDTree.hh:271
void InactiveNode(G4KDNode_Base *)
Definition: G4KDNode.cc:58
double * GetMin()
Definition: G4KDTree.hh:214
G4KDNode stores one entity in G4KDTree This class is for internal use only.
Definition: G4KDNode.hh:128
void Extend(const Position &pos)
Definition: G4KDTree.hh:173
void Free(G4KDNode_Base *&)
Definition: G4KDNode.cc:64
G4KDTreeResultHandle NearestInRange(const Position &pos, const double &range)
void Build()
Definition: G4KDTree.cc:118
void __InsertMap(G4KDNode_Base *node)
Definition: G4KDTree.cc:113
T max(const T t1, const T t2)
brief Return the largest of the two arguments
size_t GetDim() const
Definition: G4KDTree.hh:88
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double * GetMax()
Definition: G4KDTree.hh:218
~G4KDTree()
Definition: G4KDTree.cc:60
HyperRect & operator=(const HyperRect &rhs)
Definition: G4KDTree.hh:229
void __Clear_Rec(G4KDNode_Base *node)
Definition: G4KDTree.cc:103
void __NearestToPosition(G4KDNode_Base *node, const Position &pos, G4KDNode_Base *&result, double *result_dist_sq, HyperRect *fRect)
G4KDNode_Base * Insert(PointT *pos)
T sqr(const T &x)
Definition: templates.hh:145
int __NearestInRange(G4KDNode_Base *node, const Position &pos, const double &range_sq, const double &range, G4KDTreeResult &list, int ordered, G4KDNode_Base *source_node=0)
HyperRect(size_t dim)
Definition: G4KDTree.hh:136
G4ThreadLocalStatic G4Allocator< G4KDTree > * fgAllocator
Definition: G4KDTree.hh:273
void NoticeNodeDeactivation()
Definition: G4KDTree.hh:82
static const G4double pos
G4KDNode_Base * InsertMap(PointT *pos)
G4KDTreeResult enables to go through the nearest entities found by G4KDTree.