Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4KDTree.hh
이 파일의 문서화 페이지로 가기
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 110873 2018-06-22 13:11:22Z 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> class G4KDNode;
56 
57 //__________________________________
58 // Methods to act on kdnode
59 // Methods defined in G4KDNode.cc :
61 void Free(G4KDNode_Base*&);
62 //void* GetData(G4KDNode*);
63 //const double* GetNodePosition(G4KDNode_Base*);
64 //__________________________________
65 
71 class G4KDTree
72 {
73  friend class G4KDNode_Base;
74 public:
75  G4KDTree(size_t dim = 3);
76  ~G4KDTree();
77  void Clear();
78 
79  void Print(std::ostream& out = G4cout) const;
80  void Build();
82  {
84  if (fNbActiveNodes <= 0) Clear();
85  }
86 
87  size_t GetDim() const
88  {
89  return fDim;
90  }
91  int GetNbNodes() const
92  {
93  return fNbNodes;
94  }
96  {
97  return fRoot;
98  }
99 
100  template<typename PointT>
101  G4KDNode_Base* InsertMap(PointT* pos);
102 
103  // Insert and attache the data to a node at the specified position
104  // In return, it gives you the corresponding node
105  template<typename PointT> G4KDNode_Base* Insert(PointT* pos); // 3D
106 
107  template<typename PointT> G4KDNode_Base* Insert(const PointT& pos); // 3D
108 
109  /* Find one of the nearest nodes from the specified point.
110  *
111  * This function returns a pointer to a result set with at most one element.
112  */
113  template<typename Position> G4KDTreeResultHandle Nearest(const Position& pos);
115 
116  /* Find any nearest nodes from the specified point within a range.
117  *
118  * This function returns a pointer to a result set, which can be manipulated
119  * by the G4KDTreeResult.
120  * The returned pointer can be null as an indication of an error. Otherwise
121  * a valid result set is always returned which may contain 0 or more elements.
122  */
123  template<typename Position>
124  G4KDTreeResultHandle NearestInRange(const Position& pos,
125  const double& range);
127 
128  void *operator new(size_t);
129  void operator delete(void *);
130 
131 protected:
132 
133  //______________________________________________________________________
134  class HyperRect
135  {
136  public:
137  HyperRect(size_t dim)
138  {
139  fDim = dim;
140  fMin = new double[fDim];
141  fMax = new double[fDim];
142  }
143 
144  template<typename Position>
145  void SetMinMax(const Position& min, const Position& max)
146  {
147  for (size_t i = 0; i < fDim; i++)
148  {
149  fMin[i] = min[i];
150  fMax[i] = max[i];
151  }
152  }
153 
155  {
156  delete[] fMin;
157  delete[] fMax;
158  }
159 
160  HyperRect(const HyperRect& rect)
161  {
162  fDim = rect.fDim;
163  fMin = new double[fDim];
164  fMax = new double[fDim];
165 
166  for (size_t i = 0; i < fDim; i++)
167  {
168  fMin[i] = rect.fMin[i];
169  fMax[i] = rect.fMax[i];
170  }
171  }
172 
173  template<typename Position>
174  void Extend(const Position& pos)
175  {
176  for (size_t i = 0; i < fDim; i++)
177  {
178  if (pos[i] < fMin[i])
179  {
180  fMin[i] = pos[i];
181  }
182  if (pos[i] > fMax[i])
183  {
184  fMax[i] = pos[i];
185  }
186  }
187  }
188 
189  template<typename Position>
190  bool CompareDistSqr(const Position& pos, const double* bestmatch)
191  {
192  double result = 0;
193 
194  for (size_t i = 0; i < fDim; i++)
195  {
196  if (pos[i] < fMin[i])
197  {
198  result += sqr(fMin[i] - pos[i]);
199  }
200  else if (pos[i] > fMax[i])
201  {
202  result += sqr(fMax[i] - pos[i]);
203  }
204 
205  if (result >= *bestmatch) return false;
206  }
207 
208  return true;
209  }
210 
211  size_t GetDim()
212  {
213  return fDim;
214  }
215  double* GetMin()
216  {
217  return fMin;
218  }
219  double* GetMax()
220  {
221  return fMax;
222  }
223 
224  protected:
225  size_t fDim;
226  double *fMin, *fMax; /* minimum/maximum coords */
227 
228  private:
229  // should not be used
231  {
232  if (this == &rhs) return *this;
233  return *this;
234  }
235  };
236 
237 protected:
238  void __InsertMap(G4KDNode_Base *node);
239  void __Clear_Rec(G4KDNode_Base *node);
240 
241  template<typename Position>
243  const Position& pos,
244  const double& range_sq,
245  const double& range,
246  G4KDTreeResult& list,
247  int ordered,
248  G4KDNode_Base *source_node = 0);
249 
250  template<typename Position>
252  const Position& pos,
254  double *result_dist_sq,
255  HyperRect* fRect);
256 
257  template<typename Position>
258  void __NearestToNode(G4KDNode_Base *source_node,
259  G4KDNode_Base *node,
260  const Position& pos,
261  std::vector<G4KDNode_Base*>& result,
262  double *result_dist_sq,
263  HyperRect* fRect,
264  int& nbresult);
265 
266 protected:
269  size_t fDim;
270  int fNbNodes;
273 
275 };
276 
277 #include "G4KDTree.icc"
278 
279 #endif // G4KDTREE_HH
void Free(G4KDNode_Base *&)
Definition: G4KDNode.cc:65
~G4KDTree()
Definition: G4KDTree.cc:64
T max(const T t1, const T t2)
brief Return the largest of the two arguments
static G4Allocator< G4KDTree > *& fgAllocator()
Definition: G4KDTree.cc:46
void __NearestToPosition(G4KDNode_Base *node, const Position &pos, G4KDNode_Base *&result, double *result_dist_sq, HyperRect *fRect)
static const G4double pos
HyperRect * fRect
Definition: G4KDTree.hh:267
G4KDTree(size_t dim=3)
Definition: G4KDTree.cc:54
G4KDNode_Base * fRoot
Definition: G4KDTree.hh:268
G4KDTreeResultHandle NearestInRange(const Position &pos, const double &range)
int __NearestInRange(G4KDNode_Base *node, const Position &pos, const double &range_sq, const double &range, G4KDTreeResult &list, int ordered, G4KDNode_Base *source_node=0)
G4KDNode_Base * Insert(PointT *pos)
HyperRect & operator=(const HyperRect &rhs)
Definition: G4KDTree.hh:230
double * GetMax()
Definition: G4KDTree.hh:219
G4KDMap * fKDMap
Definition: G4KDTree.hh:272
HyperRect(const HyperRect &rect)
Definition: G4KDTree.hh:160
void Print(std::ostream &out=G4cout) const
Definition: G4KDTree.cc:89
void __InsertMap(G4KDNode_Base *node)
Definition: G4KDTree.cc:117
void NoticeNodeDeactivation()
Definition: G4KDTree.hh:81
int fNbNodes
Definition: G4KDTree.hh:270
const G4ParticleDefinition const G4Material *G4double range
size_t GetDim() const
Definition: G4KDTree.hh:87
G4KDNode_Base * InsertMap(PointT *pos)
G4KDTreeResultHandle Nearest(const Position &pos)
void SetMinMax(const Position &min, const Position &max)
Definition: G4KDTree.hh:145
void Extend(const Position &pos)
Definition: G4KDTree.hh:174
bool CompareDistSqr(const Position &pos, const double *bestmatch)
Definition: G4KDTree.hh:190
G4double G4ParticleHPJENDLHEData::G4double result
void Build()
Definition: G4KDTree.cc:122
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)
G4GLOB_DLL std::ostream G4cout
G4KDNode_Base * GetRoot()
Definition: G4KDTree.hh:95
HyperRect(size_t dim)
Definition: G4KDTree.hh:137
int GetNbNodes() const
Definition: G4KDTree.hh:91
void InactiveNode(G4KDNode_Base *)
Definition: G4KDNode.cc:58
T sqr(const T &x)
Definition: templates.hh:145
size_t fDim
Definition: G4KDTree.hh:269
void Clear()
Definition: G4KDTree.cc:94
int fNbActiveNodes
Definition: G4KDTree.hh:271
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
double * GetMin()
Definition: G4KDTree.hh:215
void __Clear_Rec(G4KDNode_Base *node)
Definition: G4KDTree.cc:107