Geant4  10.02.p01
G4KDTree.cc
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 /*
27  * G4KDTree.cc
28  *
29  * Created on: 22 oct. 2013
30  * Author: kara
31  */
32 
33 #include "globals.hh"
34 #include <stdio.h>
35 #include <stdlib.h>
36 #include <cmath>
37 #include "G4KDTree.hh"
38 #include "G4KDMap.hh"
39 #include "G4KDNode.hh"
40 #include "G4KDTreeResult.hh"
41 #include <list>
42 #include <iostream>
43 
44 using namespace std;
45 
47 
48 //______________________________________________________________________
49 // KDTree methods
50 G4KDTree::G4KDTree(size_t k) :
51  fKDMap(new G4KDMap(k))
52 {
53  fDim = k;
54  fRoot = 0;
55  fRect = 0;
56  fNbNodes = 0;
57  fNbActiveNodes = 0;
58 }
59 
61 {
62  if (fRoot) __Clear_Rec(fRoot);
63  fRoot = 0;
64 
65  if (fRect)
66  {
67  delete fRect;
68  fRect = 0;
69  }
70 
71  if (fKDMap) delete fKDMap;
72 }
73 
74 void* G4KDTree::operator new(size_t)
75 {
76  if (!fgAllocator) fgAllocator = new G4Allocator<G4KDTree>;
77  return (void *) fgAllocator->MallocSingle();
78 }
79 
80 void G4KDTree::operator delete(void *aNode)
81 {
82  fgAllocator->FreeSingle((G4KDTree*) aNode);
83 }
84 
85 void G4KDTree::Print(std::ostream& out) const
86 {
87  if (fRoot) fRoot->Print(out);
88 }
89 
91 {
93  fRoot = 0;
94  fNbNodes = 0;
95 
96  if (fRect)
97  {
98  delete fRect;
99  fRect = 0;
100  }
101 }
102 
104 {
105  if (!node) return;
106 
107  if (node->GetLeft()) __Clear_Rec(node->GetLeft());
108  if (node->GetRight()) __Clear_Rec(node->GetRight());
109 
110  delete node;
111 }
112 
114 {
115  fKDMap->Insert(node);
116 }
117 
119 {
120  size_t Nnodes = fKDMap->GetSize();
121 
122  G4cout << "********************" << G4endl;
123  G4cout << "template<typename PointT> G4KDTree<PointT>::Build" << G4endl;
124  G4cout << "Map size = " << Nnodes << G4endl;
125 
126  G4KDNode_Base* root = fKDMap->PopOutMiddle(0);
127 
128  if(root == 0) return;
129 
130  fRoot = root;
131  fNbActiveNodes++;
132  fRect = new HyperRect(fDim);
133  fRect->SetMinMax(*fRoot, *fRoot);
134 
135  Nnodes--;
136 
137  G4KDNode_Base* parent = fRoot;
138 
139  for (size_t n = 0; n < Nnodes; n += fDim)
140  {
141  for (size_t dim = 0; dim < fDim; dim++)
142  {
143  G4KDNode_Base* node = fKDMap->PopOutMiddle(dim);
144  if (node)
145  {
146  parent->Insert(node);
147  fNbActiveNodes++;
148  fRect->Extend(*node);
149  parent = node;
150  }
151  }
152  }
153 }
154 
156 {
157  // G4cout << "Nearest(node)" << G4endl ;
158  if (!fRect)
159  {
160  G4cout << "Tree empty" << G4endl;
161  return 0;
162  }
163 
164  std::vector<G4KDNode_Base* > result;
165  double dist_sq = DBL_MAX;
166 
167  /* Duplicate the bounding hyperrectangle, we will work on the copy */
168  HyperRect *newrect = new HyperRect(*fRect);
169 
170  /* Search for the nearest neighbour recursively */
171  int nbresult = 0;
172 
173  __NearestToNode(node, fRoot, *node, result, &dist_sq, newrect, nbresult);
174 
175  /* Free the copy of the hyperrect */
176  delete newrect;
177 
178  /* Store the result */
179  if (!result.empty())
180  {
181  G4KDTreeResultHandle rset(new G4KDTreeResult(this));
182  int j = 0;
183  while (j<nbresult)
184  {
185  rset->Insert(dist_sq, result[j]);
186  j++;
187  }
188  rset->Rewind();
189 
190  return rset;
191  }
192  else
193  {
194  return 0;
195  }
196 }
197 
199  const double& range)
200 {
201  if (!node) return 0;
202  int ret(-1);
203 
204  G4KDTreeResult *rset = new G4KDTreeResult(this);
205 
206  const double range_sq = sqr(range);
207 
208  if ((ret = __NearestInRange(fRoot, *node, range_sq, range, *rset, 0, node)) == -1)
209  {
210  delete rset;
211  return 0;
212  }
213  rset->Sort();
214  rset->Rewind();
215  return rset;
216 }
Type * MallocSingle()
Definition: G4Allocator.hh:202
HyperRect * fRect
Definition: G4KDTree.hh:268
void Clear()
Definition: G4KDTree.cc:90
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
G4KDNode_Base * Insert(PointT *point)
G4KDNode_Base * GetRight()
Definition: G4KDNode.hh:84
size_t GetSize()
Definition: G4KDMap.hh:110
G4KDTree is used by the ITManager to locate the neareast neighbours.
Definition: G4KDTree.hh:72
int fNbNodes
Definition: G4KDTree.hh:271
G4KDNode_Base * PopOutMiddle(size_t dimension)
Definition: G4KDMap.cc:138
G4KDTree(size_t dim=3)
Definition: G4KDTree.cc:50
#define G4ThreadLocal
Definition: tls.hh:89
void Print(std::ostream &out, int level=0) const
Definition: G4KDNode.cc:178
G4KDTreeResultHandle Nearest(const Position &pos)
int fNbActiveNodes
Definition: G4KDTree.hh:272
G4GLOB_DLL std::ostream G4cout
G4KDNode_Base * fRoot
Definition: G4KDTree.hh:269
G4KDMap * fKDMap
Definition: G4KDTree.hh:273
void Extend(const Position &pos)
Definition: G4KDTree.hh:175
G4KDTreeResultHandle NearestInRange(const Position &pos, const double &range)
const G4int n
void Build()
Definition: G4KDTree.cc:118
void __InsertMap(G4KDNode_Base *node)
Definition: G4KDTree.cc:113
#define G4endl
Definition: G4ios.hh:61
~G4KDTree()
Definition: G4KDTree.cc:60
void Insert(G4KDNode_Base *pos)
Definition: G4KDMap.cc:90
void __Clear_Rec(G4KDNode_Base *node)
Definition: G4KDTree.cc:103
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)
#define DBL_MAX
Definition: templates.hh:83
G4ThreadLocalStatic G4Allocator< G4KDTree > * fgAllocator
Definition: G4KDTree.hh:275
G4KDNode_Base * GetLeft()
Definition: G4KDNode.hh:83
G4KDTreeResult enables to go through the nearest entities found by G4KDTree.