| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 /* 26 /*
27 * G4KDTree.cc 27 * G4KDTree.cc
28 * 28 *
29 * Created on: 22 oct. 2013 29 * Created on: 22 oct. 2013
30 * Author: kara 30 * Author: kara
31 */ 31 */
32 32
33 #include "globals.hh" 33 #include "globals.hh"
34 #include <cstdio> << 34 #include <stdio.h>
>> 35 #include <stdlib.h>
35 #include <cmath> 36 #include <cmath>
36 #include "G4KDTree.hh" 37 #include "G4KDTree.hh"
37 #include "G4KDMap.hh" 38 #include "G4KDMap.hh"
38 #include "G4KDNode.hh" 39 #include "G4KDNode.hh"
39 #include "G4KDTreeResult.hh" 40 #include "G4KDTreeResult.hh"
40 #include <list> 41 #include <list>
41 #include <iostream> 42 #include <iostream>
42 43
43 using namespace std; 44 using namespace std;
44 45
45 G4Allocator<G4KDTree>*& G4KDTree::fgAllocator( << 46 G4ThreadLocal G4Allocator<G4KDTree>* G4KDTree::fgAllocator(0);
46 { <<
47 G4ThreadLocalStatic G4Allocator<G4KDTree>* _ <<
48 return _instance; <<
49 } <<
50 47
51 //____________________________________________ 48 //______________________________________________________________________
52 // KDTree methods 49 // KDTree methods
53 G4KDTree::G4KDTree(size_t k) << 50 G4KDTree::G4KDTree(size_t k) :
54 : fDim(k) << 51 fKDMap(new G4KDMap(k))
55 ,fKDMap(new G4KDMap(k)) << 52 {
56 {} << 53 fDim = k;
>> 54 fRoot = 0;
>> 55 fRect = 0;
>> 56 fNbNodes = 0;
>> 57 fNbActiveNodes = 0;
>> 58 }
57 59
58 G4KDTree::~G4KDTree() 60 G4KDTree::~G4KDTree()
59 { 61 {
60 if(fRoot != nullptr){ << 62 if (fRoot) __Clear_Rec(fRoot);
61 __Clear_Rec(fRoot); << 63 fRoot = 0;
62 fRoot = nullptr; <<
63 } <<
64 64
65 if(fRect != nullptr){ << 65 if (fRect)
>> 66 {
66 delete fRect; 67 delete fRect;
67 fRect = nullptr; << 68 fRect = 0;
68 } 69 }
69 70
70 if(fKDMap != nullptr){ << 71 if (fKDMap) delete fKDMap;
71 delete fKDMap; <<
72 fKDMap = nullptr; <<
73 } <<
74 } 72 }
75 73
76 void* G4KDTree::operator new(size_t) 74 void* G4KDTree::operator new(size_t)
77 { 75 {
78 if(fgAllocator() == nullptr){ << 76 if (!fgAllocator) fgAllocator = new G4Allocator<G4KDTree>;
79 fgAllocator() = new G4Allocator<G4KDTree>; << 77 return (void *) fgAllocator->MallocSingle();
80 } <<
81 return (void*) fgAllocator()->MallocSingle() <<
82 } 78 }
83 79
84 void G4KDTree::operator delete(void* aNode) << 80 void G4KDTree::operator delete(void *aNode)
85 { 81 {
86 fgAllocator()->FreeSingle((G4KDTree*) aNode) << 82 fgAllocator->FreeSingle((G4KDTree*) aNode);
87 } 83 }
88 84
89 void G4KDTree::Print(std::ostream& out) const 85 void G4KDTree::Print(std::ostream& out) const
90 { 86 {
91 if(fRoot != nullptr){ << 87 if (fRoot) fRoot->Print(out);
92 fRoot->Print(out); <<
93 } <<
94 } 88 }
95 89
96 void G4KDTree::Clear() 90 void G4KDTree::Clear()
97 { 91 {
98 __Clear_Rec(fRoot); 92 __Clear_Rec(fRoot);
99 fRoot = nullptr; << 93 fRoot = 0;
100 fNbNodes = 0; 94 fNbNodes = 0;
101 95
102 if(fRect != nullptr) << 96 if (fRect)
103 { 97 {
104 delete fRect; 98 delete fRect;
105 fRect = nullptr; << 99 fRect = 0;
106 } 100 }
107 } 101 }
108 102
109 void G4KDTree::__Clear_Rec(G4KDNode_Base* node 103 void G4KDTree::__Clear_Rec(G4KDNode_Base* node)
110 { 104 {
111 if(node == nullptr) << 105 if (!node) return;
112 { <<
113 return; <<
114 } <<
115 106
116 if(node->GetLeft() != nullptr) << 107 if (node->GetLeft()) __Clear_Rec(node->GetLeft());
117 { << 108 if (node->GetRight()) __Clear_Rec(node->GetRight());
118 __Clear_Rec(node->GetLeft()); <<
119 } <<
120 if(node->GetRight() != nullptr) <<
121 { <<
122 __Clear_Rec(node->GetRight()); <<
123 } <<
124 109
125 delete node; 110 delete node;
126 } 111 }
127 112
128 void G4KDTree::__InsertMap(G4KDNode_Base* node << 113 void G4KDTree::__InsertMap(G4KDNode_Base *node)
>> 114 {
>> 115 fKDMap->Insert(node);
>> 116 }
129 117
130 void G4KDTree::Build() 118 void G4KDTree::Build()
131 { 119 {
132 size_t Nnodes = fKDMap->GetSize(); 120 size_t Nnodes = fKDMap->GetSize();
133 121
134 G4cout << "********************" << G4endl; 122 G4cout << "********************" << G4endl;
135 G4cout << "template<typename PointT> G4KDTre 123 G4cout << "template<typename PointT> G4KDTree<PointT>::Build" << G4endl;
136 G4cout << "Map size = " << Nnodes << G4endl; 124 G4cout << "Map size = " << Nnodes << G4endl;
137 125
138 G4KDNode_Base* root = fKDMap->PopOutMiddle(0 126 G4KDNode_Base* root = fKDMap->PopOutMiddle(0);
139 127
140 if(root == nullptr) << 128 if(root == 0) return;
141 { <<
142 return; <<
143 } <<
144 129
145 fRoot = root; 130 fRoot = root;
146 fNbActiveNodes++; 131 fNbActiveNodes++;
147 fRect = new HyperRect(fDim); 132 fRect = new HyperRect(fDim);
148 fRect->SetMinMax(*fRoot, *fRoot); 133 fRect->SetMinMax(*fRoot, *fRoot);
149 134
150 Nnodes--; 135 Nnodes--;
151 136
152 G4KDNode_Base* parent = fRoot; 137 G4KDNode_Base* parent = fRoot;
153 138
154 for(size_t n = 0; n < Nnodes; n += fDim) << 139 for (size_t n = 0; n < Nnodes; n += fDim)
155 { 140 {
156 for(size_t dim = 0; dim < fDim; dim++) << 141 for (size_t dim = 0; dim < fDim; dim++)
157 { 142 {
158 G4KDNode_Base* node = fKDMap->PopOutMidd 143 G4KDNode_Base* node = fKDMap->PopOutMiddle(dim);
159 if(node != nullptr) << 144 if (node)
160 { 145 {
161 parent->Insert(node); 146 parent->Insert(node);
162 fNbActiveNodes++; 147 fNbActiveNodes++;
163 fRect->Extend(*node); 148 fRect->Extend(*node);
164 parent = node; 149 parent = node;
165 } 150 }
166 } 151 }
167 } 152 }
168 } 153 }
169 154
170 G4KDTreeResultHandle G4KDTree::Nearest(G4KDNod 155 G4KDTreeResultHandle G4KDTree::Nearest(G4KDNode_Base* node)
171 { 156 {
172 if(fRect == nullptr) << 157 // G4cout << "Nearest(node)" << G4endl ;
>> 158 if (!fRect)
173 { 159 {
174 return nullptr; << 160 G4cout << "Tree empty" << G4endl;
>> 161 return 0;
175 } 162 }
176 163
177 std::vector<G4KDNode_Base*> result; << 164 std::vector<G4KDNode_Base* > result;
178 G4double dist_sq = DBL_MAX; << 165 double dist_sq = DBL_MAX;
179 166
180 /* Duplicate the bounding hyperrectangle, we 167 /* Duplicate the bounding hyperrectangle, we will work on the copy */
181 auto newrect = new HyperRect(*fRect); << 168 HyperRect *newrect = new HyperRect(*fRect);
182 169
183 /* Search for the nearest neighbour recursiv 170 /* Search for the nearest neighbour recursively */
184 G4int nbresult = 0; << 171 int nbresult = 0;
185 172
186 __NearestToNode(node, fRoot, *node, result, 173 __NearestToNode(node, fRoot, *node, result, &dist_sq, newrect, nbresult);
187 174
188 /* Free the copy of the hyperrect */ 175 /* Free the copy of the hyperrect */
189 delete newrect; 176 delete newrect;
190 177
191 /* Store the result */ 178 /* Store the result */
192 if(!result.empty()) << 179 if (!result.empty())
193 { 180 {
194 G4KDTreeResultHandle rset(new G4KDTreeResu 181 G4KDTreeResultHandle rset(new G4KDTreeResult(this));
195 G4int j = 0; << 182 int j = 0;
196 while(j < nbresult) << 183 while (j<nbresult)
197 { 184 {
198 rset->Insert(dist_sq, result[j]); 185 rset->Insert(dist_sq, result[j]);
199 j++; 186 j++;
200 } 187 }
201 rset->Rewind(); 188 rset->Rewind();
202 189
203 return rset; 190 return rset;
204 } 191 }
205 << 192 else
206 return nullptr; << 193 {
>> 194 return 0;
>> 195 }
207 } 196 }
208 197
209 G4KDTreeResultHandle G4KDTree::NearestInRange( 198 G4KDTreeResultHandle G4KDTree::NearestInRange(G4KDNode_Base* node,
210 << 199 const double& range)
211 { 200 {
212 if(node == nullptr) << 201 if (!node) return 0;
213 { << 202 int ret(-1);
214 return nullptr; <<
215 } <<
216 G4int ret(-1); <<
217 203
218 auto* rset = new G4KDTreeResult(this); << 204 G4KDTreeResult *rset = new G4KDTreeResult(this);
219 205
220 const G4double range_sq = sqr(range); << 206 const double range_sq = sqr(range);
221 207
222 if((ret = __NearestInRange(fRoot, *node, ran << 208 if ((ret = __NearestInRange(fRoot, *node, range_sq, range, *rset, 0, node)) == -1)
223 -1) <<
224 { 209 {
225 delete rset; 210 delete rset;
226 return nullptr; << 211 return 0;
227 } 212 }
228 rset->Sort(); 213 rset->Sort();
229 rset->Rewind(); 214 rset->Rewind();
230 return rset; 215 return rset;
231 } 216 }
232 217