Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 100802 2016-11-02 14:55:27Z gcosmo $
27 //
28 // Author: Mathieu Karamitros
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  template<typename PointT> G4KDNode_Base* Insert(const PointT& pos); // 3D
109 
110  /* Find one of the nearest nodes from the specified point.
111  *
112  * This function returns a pointer to a result set with at most one element.
113  */
114  template<typename Position> G4KDTreeResultHandle Nearest(const Position& pos);
116 
117  /* Find any nearest nodes from the specified point within a range.
118  *
119  * This function returns a pointer to a result set, which can be manipulated
120  * by the G4KDTreeResult.
121  * The returned pointer can be null as an indication of an error. Otherwise
122  * a valid result set is always returned which may contain 0 or more elements.
123  */
124  template<typename Position>
125  G4KDTreeResultHandle NearestInRange(const Position& pos,
126  const double& range);
128 
129  void *operator new(size_t);
130  void operator delete(void *);
131 
132 protected:
133 
134  //______________________________________________________________________
135  class HyperRect
136  {
137  public:
138  HyperRect(size_t dim)
139  {
140  fDim = dim;
141  fMin = new double[fDim];
142  fMax = new double[fDim];
143  }
144 
145  template<typename Position>
146  void SetMinMax(const Position& min, const Position& max)
147  {
148  for (size_t i = 0; i < fDim; i++)
149  {
150  fMin[i] = min[i];
151  fMax[i] = max[i];
152  }
153  }
154 
156  {
157  delete[] fMin;
158  delete[] fMax;
159  }
160 
161  HyperRect(const HyperRect& rect)
162  {
163  fDim = rect.fDim;
164  fMin = new double[fDim];
165  fMax = new double[fDim];
166 
167  for (size_t i = 0; i < fDim; i++)
168  {
169  fMin[i] = rect.fMin[i];
170  fMax[i] = rect.fMax[i];
171  }
172  }
173 
174  template<typename Position>
175  void Extend(const Position& pos)
176  {
177  for (size_t i = 0; i < fDim; i++)
178  {
179  if (pos[i] < fMin[i])
180  {
181  fMin[i] = pos[i];
182  }
183  if (pos[i] > fMax[i])
184  {
185  fMax[i] = pos[i];
186  }
187  }
188  }
189 
190  template<typename Position>
191  bool CompareDistSqr(const Position& pos, const double* bestmatch)
192  {
193  double result = 0;
194 
195  for (size_t i = 0; i < fDim; i++)
196  {
197  if (pos[i] < fMin[i])
198  {
199  result += sqr(fMin[i] - pos[i]);
200  }
201  else if (pos[i] > fMax[i])
202  {
203  result += sqr(fMax[i] - pos[i]);
204  }
205 
206  if (result >= *bestmatch) return false;
207  }
208 
209  return true;
210  }
211 
212  size_t GetDim()
213  {
214  return fDim;
215  }
216  double* GetMin()
217  {
218  return fMin;
219  }
220  double* GetMax()
221  {
222  return fMax;
223  }
224 
225  protected:
226  size_t fDim;
227  double *fMin, *fMax; /* minimum/maximum coords */
228 
229  private:
230  // should not be used
231  HyperRect& operator=(const HyperRect& rhs)
232  {
233  if (this == &rhs) return *this;
234  return *this;
235  }
236  };
237 
238 protected:
239  void __InsertMap(G4KDNode_Base *node);
240  void __Clear_Rec(G4KDNode_Base *node);
241 
242  template<typename Position>
244  const Position& pos,
245  const double& range_sq,
246  const double& range,
247  G4KDTreeResult& list,
248  int ordered,
249  G4KDNode_Base *source_node = 0);
250 
251  template<typename Position>
253  const Position& pos,
255  double *result_dist_sq,
256  HyperRect* fRect);
257 
258  template<typename Position>
259  void __NearestToNode(G4KDNode_Base *source_node,
260  G4KDNode_Base *node,
261  const Position& pos,
262  std::vector<G4KDNode_Base*>& result,
263  double *result_dist_sq,
264  HyperRect* fRect,
265  int& nbresult);
266 
267 protected:
270  size_t fDim;
271  int fNbNodes;
274 
276 };
277 
278 #include "G4KDTree.icc"
279 
280 #endif // G4KDTREE_HH
G4double G4ParticleHPJENDLHEData::G4double result
HyperRect * fRect
Definition: G4KDTree.hh:268
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:146
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:270
int GetNbNodes() const
Definition: G4KDTree.hh:92
int fNbNodes
Definition: G4KDTree.hh:271
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:161
int fNbActiveNodes
Definition: G4KDTree.hh:272
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition const G4Material *G4double range
G4KDNode_Base * fRoot
Definition: G4KDTree.hh:269
bool CompareDistSqr(const Position &pos, const double *bestmatch)
Definition: G4KDTree.hh:191
G4KDMap * fKDMap
Definition: G4KDTree.hh:273
void InactiveNode(G4KDNode_Base *)
Definition: G4KDNode.cc:58
double * GetMin()
Definition: G4KDTree.hh:216
void Extend(const Position &pos)
Definition: G4KDTree.hh:175
void Free(G4KDNode_Base *&)
Definition: G4KDNode.cc:65
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:220
~G4KDTree()
Definition: G4KDTree.cc:60
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:138
G4ThreadLocalStatic G4Allocator< G4KDTree > * fgAllocator
Definition: G4KDTree.hh:275
void NoticeNodeDeactivation()
Definition: G4KDTree.hh:82
static const G4double pos
G4KDNode_Base * InsertMap(PointT *pos)