| 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 <cstdio>
35 #include <cmath> 35 #include <cmath>
36 #include "G4KDTree.hh" 36 #include "G4KDTree.hh"
37 #include "G4KDMap.hh" 37 #include "G4KDMap.hh"
38 #include "G4KDNode.hh" 38 #include "G4KDNode.hh"
39 #include "G4KDTreeResult.hh" 39 #include "G4KDTreeResult.hh"
40 #include <list> 40 #include <list>
41 #include <iostream> 41 #include <iostream>
42 42
43 using namespace std; 43 using namespace std;
44 44
45 G4Allocator<G4KDTree>*& G4KDTree::fgAllocator( 45 G4Allocator<G4KDTree>*& G4KDTree::fgAllocator()
46 { 46 {
47 G4ThreadLocalStatic G4Allocator<G4KDTree>* _ 47 G4ThreadLocalStatic G4Allocator<G4KDTree>* _instance = nullptr;
48 return _instance; 48 return _instance;
49 } 49 }
50 50
51 //____________________________________________ 51 //______________________________________________________________________
52 // KDTree methods 52 // KDTree methods
53 G4KDTree::G4KDTree(size_t k) 53 G4KDTree::G4KDTree(size_t k)
54 : fDim(k) 54 : fDim(k)
55 ,fKDMap(new G4KDMap(k)) 55 ,fKDMap(new G4KDMap(k))
56 {} 56 {}
57 57
58 G4KDTree::~G4KDTree() 58 G4KDTree::~G4KDTree()
59 { 59 {
60 if(fRoot != nullptr){ 60 if(fRoot != nullptr){
61 __Clear_Rec(fRoot); 61 __Clear_Rec(fRoot);
62 fRoot = nullptr; 62 fRoot = nullptr;
63 } 63 }
64 64
65 if(fRect != nullptr){ 65 if(fRect != nullptr){
66 delete fRect; 66 delete fRect;
67 fRect = nullptr; 67 fRect = nullptr;
68 } 68 }
69 69
70 if(fKDMap != nullptr){ 70 if(fKDMap != nullptr){
71 delete fKDMap; 71 delete fKDMap;
72 fKDMap = nullptr; 72 fKDMap = nullptr;
73 } 73 }
74 } 74 }
75 75
76 void* G4KDTree::operator new(size_t) 76 void* G4KDTree::operator new(size_t)
77 { 77 {
78 if(fgAllocator() == nullptr){ 78 if(fgAllocator() == nullptr){
79 fgAllocator() = new G4Allocator<G4KDTree>; 79 fgAllocator() = new G4Allocator<G4KDTree>;
80 } 80 }
81 return (void*) fgAllocator()->MallocSingle() 81 return (void*) fgAllocator()->MallocSingle();
82 } 82 }
83 83
84 void G4KDTree::operator delete(void* aNode) 84 void G4KDTree::operator delete(void* aNode)
85 { 85 {
86 fgAllocator()->FreeSingle((G4KDTree*) aNode) 86 fgAllocator()->FreeSingle((G4KDTree*) aNode);
87 } 87 }
88 88
89 void G4KDTree::Print(std::ostream& out) const 89 void G4KDTree::Print(std::ostream& out) const
90 { 90 {
91 if(fRoot != nullptr){ 91 if(fRoot != nullptr){
92 fRoot->Print(out); 92 fRoot->Print(out);
93 } 93 }
94 } 94 }
95 95
96 void G4KDTree::Clear() 96 void G4KDTree::Clear()
97 { 97 {
98 __Clear_Rec(fRoot); 98 __Clear_Rec(fRoot);
99 fRoot = nullptr; 99 fRoot = nullptr;
100 fNbNodes = 0; 100 fNbNodes = 0;
101 101
102 if(fRect != nullptr) 102 if(fRect != nullptr)
103 { 103 {
104 delete fRect; 104 delete fRect;
105 fRect = nullptr; 105 fRect = nullptr;
106 } 106 }
107 } 107 }
108 108
109 void G4KDTree::__Clear_Rec(G4KDNode_Base* node 109 void G4KDTree::__Clear_Rec(G4KDNode_Base* node)
110 { 110 {
111 if(node == nullptr) 111 if(node == nullptr)
112 { 112 {
113 return; 113 return;
114 } 114 }
115 115
116 if(node->GetLeft() != nullptr) 116 if(node->GetLeft() != nullptr)
117 { 117 {
118 __Clear_Rec(node->GetLeft()); 118 __Clear_Rec(node->GetLeft());
119 } 119 }
120 if(node->GetRight() != nullptr) 120 if(node->GetRight() != nullptr)
121 { 121 {
122 __Clear_Rec(node->GetRight()); 122 __Clear_Rec(node->GetRight());
123 } 123 }
124 124
125 delete node; 125 delete node;
126 } 126 }
127 127
128 void G4KDTree::__InsertMap(G4KDNode_Base* node 128 void G4KDTree::__InsertMap(G4KDNode_Base* node) { fKDMap->Insert(node); }
129 129
130 void G4KDTree::Build() 130 void G4KDTree::Build()
131 { 131 {
132 size_t Nnodes = fKDMap->GetSize(); 132 size_t Nnodes = fKDMap->GetSize();
133 133
134 G4cout << "********************" << G4endl; 134 G4cout << "********************" << G4endl;
135 G4cout << "template<typename PointT> G4KDTre 135 G4cout << "template<typename PointT> G4KDTree<PointT>::Build" << G4endl;
136 G4cout << "Map size = " << Nnodes << G4endl; 136 G4cout << "Map size = " << Nnodes << G4endl;
137 137
138 G4KDNode_Base* root = fKDMap->PopOutMiddle(0 138 G4KDNode_Base* root = fKDMap->PopOutMiddle(0);
139 139
140 if(root == nullptr) 140 if(root == nullptr)
141 { 141 {
142 return; 142 return;
143 } 143 }
144 144
145 fRoot = root; 145 fRoot = root;
146 fNbActiveNodes++; 146 fNbActiveNodes++;
147 fRect = new HyperRect(fDim); 147 fRect = new HyperRect(fDim);
148 fRect->SetMinMax(*fRoot, *fRoot); 148 fRect->SetMinMax(*fRoot, *fRoot);
149 149
150 Nnodes--; 150 Nnodes--;
151 151
152 G4KDNode_Base* parent = fRoot; 152 G4KDNode_Base* parent = fRoot;
153 153
154 for(size_t n = 0; n < Nnodes; n += fDim) 154 for(size_t n = 0; n < Nnodes; n += fDim)
155 { 155 {
156 for(size_t dim = 0; dim < fDim; dim++) 156 for(size_t dim = 0; dim < fDim; dim++)
157 { 157 {
158 G4KDNode_Base* node = fKDMap->PopOutMidd 158 G4KDNode_Base* node = fKDMap->PopOutMiddle(dim);
159 if(node != nullptr) 159 if(node != nullptr)
160 { 160 {
161 parent->Insert(node); 161 parent->Insert(node);
162 fNbActiveNodes++; 162 fNbActiveNodes++;
163 fRect->Extend(*node); 163 fRect->Extend(*node);
164 parent = node; 164 parent = node;
165 } 165 }
166 } 166 }
167 } 167 }
168 } 168 }
169 169
170 G4KDTreeResultHandle G4KDTree::Nearest(G4KDNod 170 G4KDTreeResultHandle G4KDTree::Nearest(G4KDNode_Base* node)
171 { 171 {
172 if(fRect == nullptr) 172 if(fRect == nullptr)
173 { 173 {
174 return nullptr; 174 return nullptr;
175 } 175 }
176 176
177 std::vector<G4KDNode_Base*> result; 177 std::vector<G4KDNode_Base*> result;
178 G4double dist_sq = DBL_MAX; 178 G4double dist_sq = DBL_MAX;
179 179
180 /* Duplicate the bounding hyperrectangle, we 180 /* Duplicate the bounding hyperrectangle, we will work on the copy */
181 auto newrect = new HyperRect(*fRect); 181 auto newrect = new HyperRect(*fRect);
182 182
183 /* Search for the nearest neighbour recursiv 183 /* Search for the nearest neighbour recursively */
184 G4int nbresult = 0; 184 G4int nbresult = 0;
185 185
186 __NearestToNode(node, fRoot, *node, result, 186 __NearestToNode(node, fRoot, *node, result, &dist_sq, newrect, nbresult);
187 187
188 /* Free the copy of the hyperrect */ 188 /* Free the copy of the hyperrect */
189 delete newrect; 189 delete newrect;
190 190
191 /* Store the result */ 191 /* Store the result */
192 if(!result.empty()) 192 if(!result.empty())
193 { 193 {
194 G4KDTreeResultHandle rset(new G4KDTreeResu 194 G4KDTreeResultHandle rset(new G4KDTreeResult(this));
195 G4int j = 0; 195 G4int j = 0;
196 while(j < nbresult) 196 while(j < nbresult)
197 { 197 {
198 rset->Insert(dist_sq, result[j]); 198 rset->Insert(dist_sq, result[j]);
199 j++; 199 j++;
200 } 200 }
201 rset->Rewind(); 201 rset->Rewind();
202 202
203 return rset; 203 return rset;
204 } 204 }
205 205
206 return nullptr; 206 return nullptr;
207 } 207 }
208 208
209 G4KDTreeResultHandle G4KDTree::NearestInRange( 209 G4KDTreeResultHandle G4KDTree::NearestInRange(G4KDNode_Base* node,
210 210 const G4double& range)
211 { 211 {
212 if(node == nullptr) 212 if(node == nullptr)
213 { 213 {
214 return nullptr; 214 return nullptr;
215 } 215 }
216 G4int ret(-1); 216 G4int ret(-1);
217 217
218 auto* rset = new G4KDTreeResult(this); 218 auto* rset = new G4KDTreeResult(this);
219 219
220 const G4double range_sq = sqr(range); 220 const G4double range_sq = sqr(range);
221 221
222 if((ret = __NearestInRange(fRoot, *node, ran 222 if((ret = __NearestInRange(fRoot, *node, range_sq, range, *rset, 0, node)) ==
223 -1) 223 -1)
224 { 224 {
225 delete rset; 225 delete rset;
226 return nullptr; 226 return nullptr;
227 } 227 }
228 rset->Sort(); 228 rset->Sort();
229 rset->Rewind(); 229 rset->Rewind();
230 return rset; 230 return rset;
231 } 231 }
232 232