| 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 // Author: HoangTRAN, 20/2/2019 27 // Author: HoangTRAN, 20/2/2019
28 28
29 #define OCTREE G4Octree<Iterator,Extractor,Poi 29 #define OCTREE G4Octree<Iterator,Extractor,Point>
30 #define OCTREE_TEMPLATE typename Iterator, cla 30 #define OCTREE_TEMPLATE typename Iterator, class Extractor,typename Point
31 31
32 //....oooOO0OOooo........oooOO0OOooo........oo 32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
33 33
34 template <OCTREE_TEMPLATE> 34 template <OCTREE_TEMPLATE>
35 G4ThreadLocal G4Allocator<OCTREE>* OCTREE::fgA 35 G4ThreadLocal G4Allocator<OCTREE>* OCTREE::fgAllocator = nullptr;
36 36
37 //....oooOO0OOooo........oooOO0OOooo........oo 37 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
38 38
39 template <OCTREE_TEMPLATE> 39 template <OCTREE_TEMPLATE>
40 OCTREE::G4Octree() 40 OCTREE::G4Octree()
41 : functor_(Extractor()) 41 : functor_(Extractor())
42 , head_(nullptr) 42 , head_(nullptr)
43 , size_(0) 43 , size_(0)
44 {} 44 {}
45 45
46 //....oooOO0OOooo........oooOO0OOooo........oo 46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47 47
48 template <OCTREE_TEMPLATE> 48 template <OCTREE_TEMPLATE>
49 OCTREE::G4Octree(Iterator begin, Iterator end) 49 OCTREE::G4Octree(Iterator begin, Iterator end)
50 : G4Octree(begin,end,Extractor()) 50 : G4Octree(begin,end,Extractor())
51 { 51 {
52 head_ = nullptr; 52 head_ = nullptr;
53 size_ = 0; 53 size_ = 0;
54 } 54 }
55 55
56 //....oooOO0OOooo........oooOO0OOooo........oo 56 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
57 57
58 template <OCTREE_TEMPLATE> 58 template <OCTREE_TEMPLATE>
59 OCTREE::G4Octree(Iterator begin, Iterator end, 59 OCTREE::G4Octree(Iterator begin, Iterator end, Extractor f)
60 : functor_(std::move(f)),head_(nullptr),si << 60 : functor_(f),head_(nullptr),size_(0)
61 { 61 {
62 std::vector<std::pair<Iterator,Point>> v; 62 std::vector<std::pair<Iterator,Point>> v;
63 for(auto it=begin;it!=end;++it) 63 for(auto it=begin;it!=end;++it)
64 { 64 {
65 v.push_back(std::pair<Iterator,Point>( 65 v.push_back(std::pair<Iterator,Point>(it,functor_(it)));
66 } 66 }
67 size_ = v.size(); 67 size_ = v.size();
68 head_ = new Node(v); 68 head_ = new Node(v);
69 } 69 }
70 70
71 //....oooOO0OOooo........oooOO0OOooo........oo 71 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
72 72
73 template <OCTREE_TEMPLATE> 73 template <OCTREE_TEMPLATE>
74 OCTREE::G4Octree(OCTREE::tree_type&& rhs) 74 OCTREE::G4Octree(OCTREE::tree_type&& rhs)
75 : functor_(rhs.functor_) 75 : functor_(rhs.functor_)
76 , head_(rhs.head_) 76 , head_(rhs.head_)
77 , size_(rhs.size_) 77 , size_(rhs.size_)
78 {} 78 {}
79 79
80 //....oooOO0OOooo........oooOO0OOooo........oo 80 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
81 81
82 template <OCTREE_TEMPLATE> 82 template <OCTREE_TEMPLATE>
83 void OCTREE::swap(OCTREE::tree_type& rhs) 83 void OCTREE::swap(OCTREE::tree_type& rhs)
84 { 84 {
85 std::swap(head_, rhs.head_); 85 std::swap(head_, rhs.head_);
86 std::swap(functor_, rhs.functor_); 86 std::swap(functor_, rhs.functor_);
87 std::swap(size_, rhs.size_); 87 std::swap(size_, rhs.size_);
88 } 88 }
89 89
90 //....oooOO0OOooo........oooOO0OOooo........oo 90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
91 91
92 template <OCTREE_TEMPLATE> 92 template <OCTREE_TEMPLATE>
93 typename OCTREE::tree_type& OCTREE::operator=( 93 typename OCTREE::tree_type& OCTREE::operator=(typename OCTREE::tree_type rhs)
94 { 94 {
95 swap(rhs); 95 swap(rhs);
96 return *this; 96 return *this;
97 } 97 }
98 98
99 //....oooOO0OOooo........oooOO0OOooo........oo 99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
100 100
101 template <OCTREE_TEMPLATE> 101 template <OCTREE_TEMPLATE>
102 typename OCTREE::tree_type& OCTREE::operator=( 102 typename OCTREE::tree_type& OCTREE::operator=(typename OCTREE::tree_type&& rhs)
103 { 103 {
104 swap(rhs); 104 swap(rhs);
105 return *this; 105 return *this;
106 } 106 }
107 //....oooOO0OOooo........oooOO0OOooo........oo 107 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
108 108
109 template <OCTREE_TEMPLATE> 109 template <OCTREE_TEMPLATE>
110 OCTREE::~G4Octree() 110 OCTREE::~G4Octree()
111 { 111 {
112 delete head_; 112 delete head_;
113 } 113 }
114 114
115 template <OCTREE_TEMPLATE> 115 template <OCTREE_TEMPLATE>
116 size_t OCTREE::size() const 116 size_t OCTREE::size() const
117 { 117 {
118 return size_; 118 return size_;
119 } 119 }
120 120
121 template <OCTREE_TEMPLATE> 121 template <OCTREE_TEMPLATE>
122 OCTREE::Node::Node(const NodeVector& input_val 122 OCTREE::Node::Node(const NodeVector& input_values)
123 : Node(input_values, 123 : Node(input_values,
124 G4DNABoundingBox(InnerIterator(input_ 124 G4DNABoundingBox(InnerIterator(input_values.begin()),
125 InnerIterator(input_va 125 InnerIterator(input_values.end())),
126 0) 126 0)
127 {} 127 {}
128 128
129 //....oooOO0OOooo........oooOO0OOooo........oo 129 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
130 130
131 template <OCTREE_TEMPLATE> 131 template <OCTREE_TEMPLATE>
132 OCTREE::Node::Node( 132 OCTREE::Node::Node(
133 const NodeVector& input_values, 133 const NodeVector& input_values,
134 const G4DNABoundingBox& box, 134 const G4DNABoundingBox& box,
135 size_t current_depth): 135 size_t current_depth):
136 fpValue(nullptr), 136 fpValue(nullptr),
137 fBigVolume(box), 137 fBigVolume(box),
138 fNodeType(DEFAULT) 138 fNodeType(DEFAULT)
139 { 139 {
140 if (current_depth > max_depth) 140 if (current_depth > max_depth)
141 { 141 {
142 init_max_depth_leaf(input_values); 142 init_max_depth_leaf(input_values);
143 } 143 }
144 else if (input_values.size() <= max_per_no 144 else if (input_values.size() <= max_per_node)
145 { 145 {
146 init_leaf(input_values); 146 init_leaf(input_values);
147 } 147 }
148 else 148 else
149 { 149 {
150 init_internal(input_values, current_de 150 init_internal(input_values, current_depth);
151 } 151 }
152 } 152 }
153 153
154 //....oooOO0OOooo........oooOO0OOooo........oo 154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
155 template <OCTREE_TEMPLATE> 155 template <OCTREE_TEMPLATE>
156 OCTREE::Node::~Node() 156 OCTREE::Node::~Node()
157 { 157 {
158 if (fNodeType == NodeTypes::INTERNAL) 158 if (fNodeType == NodeTypes::INTERNAL)
159 { 159 {
160 childNodeArray& children = *static_cas 160 childNodeArray& children = *static_cast<childNodeArray*>(fpValue);
161 for (size_t i = 0; i < 8; ++i) 161 for (size_t i = 0; i < 8; ++i)
162 { 162 {
163 if (children[i] != nullptr) 163 if (children[i] != nullptr)
164 { 164 {
165 delete children[i]; 165 delete children[i];
166 children[i] = nullptr; 166 children[i] = nullptr;
167 } 167 }
168 } 168 }
169 delete &children; 169 delete &children;
170 } 170 }
171 else if (fNodeType == NodeTypes::LEAF) 171 else if (fNodeType == NodeTypes::LEAF)
172 { 172 {
173 auto toDelete = static_cast<LeafValues 173 auto toDelete = static_cast<LeafValues*>(fpValue);
174 toDelete->size_ = 0; 174 toDelete->size_ = 0;
175 delete static_cast<LeafValues*>(fpValu 175 delete static_cast<LeafValues*>(fpValue);
176 } 176 }
177 else if (fNodeType == NodeTypes::MAX_DEPTH 177 else if (fNodeType == NodeTypes::MAX_DEPTH_LEAF)
178 { 178 {
179 delete static_cast<NodeVector*>(fpValu 179 delete static_cast<NodeVector*>(fpValue);
180 } 180 }
181 fpValue = nullptr; 181 fpValue = nullptr;
182 } 182 }
183 183
184 //....oooOO0OOooo........oooOO0OOooo........oo 184 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
185 185
186 template <OCTREE_TEMPLATE> 186 template <OCTREE_TEMPLATE>
187 void OCTREE::Node::init_max_depth_leaf( 187 void OCTREE::Node::init_max_depth_leaf(
188 const NodeVector& input_values) 188 const NodeVector& input_values)
189 { 189 {
190 fpValue = new NodeVector(input_values); 190 fpValue = new NodeVector(input_values);
191 fNodeType = NodeTypes::MAX_DEPTH_LEAF; 191 fNodeType = NodeTypes::MAX_DEPTH_LEAF;
192 } 192 }
193 //....oooOO0OOooo........oooOO0OOooo........oo 193 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
194 template <OCTREE_TEMPLATE> 194 template <OCTREE_TEMPLATE>
195 void OCTREE::Node::init_leaf(const NodeVector& 195 void OCTREE::Node::init_leaf(const NodeVector& input_values)
196 { 196 {
197 std::array<std::pair<Iterator, Point>, max 197 std::array<std::pair<Iterator, Point>, max_per_node> a;
198 std::copy(input_values.begin(), input_valu 198 std::copy(input_values.begin(), input_values.end(), a.begin());
199 fpValue = new LeafValues{a, input_values.s 199 fpValue = new LeafValues{a, input_values.size()};
200 fNodeType = NodeTypes::LEAF; 200 fNodeType = NodeTypes::LEAF;
201 } 201 }
202 //....oooOO0OOooo........oooOO0OOooo........oo 202 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
203 template <OCTREE_TEMPLATE> 203 template <OCTREE_TEMPLATE>
204 void OCTREE::Node::init_internal( 204 void OCTREE::Node::init_internal(
205 const NodeVector& input_values, 205 const NodeVector& input_values,
206 size_t current_depth) 206 size_t current_depth)
207 { 207 {
208 std::array<NodeVector, 8> childVectors; 208 std::array<NodeVector, 8> childVectors;
209 std::array<G4DNABoundingBox, 8> boxes = fB 209 std::array<G4DNABoundingBox, 8> boxes = fBigVolume.partition();
210 std::array<Node*, 8> children{}; 210 std::array<Node*, 8> children{};
211 for (size_t child = 0; child < 8; ++child) 211 for (size_t child = 0; child < 8; ++child)
212 { 212 {
213 NodeVector& childVector = childVectors 213 NodeVector& childVector = childVectors[child];
214 childVector.reserve(input_values.size( 214 childVector.reserve(input_values.size()/8);
215 std::copy_if(input_values.begin(),inpu 215 std::copy_if(input_values.begin(),input_values.end(),
216 std::back_inserter(childVe 216 std::back_inserter(childVector),
217 [&boxes, child](const std::pai 217 [&boxes, child](const std::pair<Iterator, Point>& element)
218 -> G4bool 218 -> G4bool
219 { 219 {
220 const Point& p = element.second; 220 const Point& p = element.second;
221 return boxes[child].contains(p); 221 return boxes[child].contains(p);
222 } 222 }
223 ); 223 );
224 children[child] = childVector.empty() 224 children[child] = childVector.empty()
225 ? nullptr 225 ? nullptr
226 : new Node(childVector, boxes[chil 226 : new Node(childVector, boxes[child], ++current_depth);
227 } 227 }
228 fpValue = new std::array<Node*, 8>(childre 228 fpValue = new std::array<Node*, 8>(children);
229 fNodeType = NodeTypes::INTERNAL; 229 fNodeType = NodeTypes::INTERNAL;
230 } 230 }
231 //....oooOO0OOooo........oooOO0OOooo........oo 231 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
232 232
233 template <OCTREE_TEMPLATE> 233 template <OCTREE_TEMPLATE>
234 template <typename OutPutContainer> 234 template <typename OutPutContainer>
235 G4bool OCTREE::Node::radiusNeighbors(const Poi 235 G4bool OCTREE::Node::radiusNeighbors(const Point& query,
236 G4double r 236 G4double radius,
237 OutPutContaine 237 OutPutContainer& resultIndices) const
238 { 238 {
239 G4bool success = false; 239 G4bool success = false;
240 G4double distance = 0; 240 G4double distance = 0;
241 if (fNodeType == NodeTypes::INTERNAL) 241 if (fNodeType == NodeTypes::INTERNAL)
242 { 242 {
243 childNodeArray& children = *static_cas 243 childNodeArray& children = *static_cast<childNodeArray*>(fpValue);
244 for (auto eachChild : children) 244 for (auto eachChild : children)
245 { 245 {
246 if (eachChild == nullptr) 246 if (eachChild == nullptr)
247 { 247 {
248 continue; 248 continue;
249 } 249 }
250 if(!eachChild->fBigVolume.overlap( 250 if(!eachChild->fBigVolume.overlap(query,radius))
251 { 251 {
252 continue; 252 continue;
253 } 253 }
254 success = eachChild->radiusNeighbo 254 success = eachChild->radiusNeighbors(query, radius, resultIndices) || success;
255 } 255 }
256 } 256 }
257 else if (fNodeType == NodeTypes::LEAF) 257 else if (fNodeType == NodeTypes::LEAF)
258 { 258 {
259 if(fpValue != nullptr) 259 if(fpValue != nullptr)
260 { 260 {
261 LeafValues& children = *static_cas 261 LeafValues& children = *static_cast<LeafValues*>(fpValue);
262 for (size_t i = 0; i < children.si 262 for (size_t i = 0; i < children.size_; ++i)
263 { 263 {
264 distance = (query - std::get<1 264 distance = (query - std::get<1>(children.values_[i])).mag();
265 265
266 if(distance != 0) 266 if(distance != 0)
267 { 267 {
268 if( distance < radius )//T 268 if( distance < radius )//TODO: find another solution for this using boundingbox
269 { 269 {
270 resultIndices.push_bac 270 resultIndices.push_back(std::make_pair(std::get<0>(children.values_[i]),distance));
271 success = true; 271 success = true;
272 } 272 }
273 } 273 }
274 } 274 }
275 } 275 }
276 } 276 }
277 else if (fNodeType == NodeTypes::MAX_DEPTH 277 else if (fNodeType == NodeTypes::MAX_DEPTH_LEAF)
278 { 278 {
279 NodeVector& children = *static_cast<No 279 NodeVector& children = *static_cast<NodeVector*>(fpValue);
280 for (auto & child : children) 280 for (auto & child : children)
281 { 281 {
282 const Point& point = std::get<1>(c 282 const Point& point = std::get<1>(child);
283 //if (this->fBigVolume.contains(qu 283 //if (this->fBigVolume.contains(query, point, radius))
284 distance = (query - point).mag(); 284 distance = (query - point).mag();
285 if( distance == 0. ) 285 if( distance == 0. )
286 { 286 {
287 continue; 287 continue;
288 } 288 }
289 if( distance < radius ) 289 if( distance < radius )
290 { 290 {
291 if(distance == 0) 291 if(distance == 0)
292 { 292 {
293 throw std::runtime_error(" 293 throw std::runtime_error("distance == 0 => OCTREE::Node::radiusNeighbors : find itself");
294 } 294 }
295 Iterator resultIndex = std::ge 295 Iterator resultIndex = std::get<0>(child);
296 resultIndices.push_back(std::m 296 resultIndices.push_back(std::make_pair(resultIndex,distance));
297 success = true; 297 success = true;
298 } 298 }
299 } 299 }
300 } 300 }
301 else 301 else
302 { 302 {
303 throw std::runtime_error("fNodeType is 303 throw std::runtime_error("fNodeType is not set : find itself");
304 } 304 }
305 return success; 305 return success;
306 } 306 }
307 307
308 //....oooOO0OOooo........oooOO0OOooo........oo 308 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
309 309
310 template <OCTREE_TEMPLATE> 310 template <OCTREE_TEMPLATE>
311 template <typename OutPutContainer> 311 template <typename OutPutContainer>
312 void OCTREE::radiusNeighbors(const Point& quer 312 void OCTREE::radiusNeighbors(const Point& query,
313 const G4double& r 313 const G4double& radius,
314 OutPutContainer& 314 OutPutContainer& resultIndices) const
315 { 315 {
316 resultIndices.clear(); 316 resultIndices.clear();
317 if (head_ == nullptr) 317 if (head_ == nullptr)
318 { 318 {
319 return; 319 return;
320 } 320 }
321 head_->radiusNeighbors(query, radius, resu 321 head_->radiusNeighbors(query, radius, resultIndices);
322 } 322 }
323 323
324 //....oooOO0OOooo........oooOO0OOooo........oo 324 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
325 325
326 template <OCTREE_TEMPLATE> 326 template <OCTREE_TEMPLATE>
327 void* OCTREE::operator new(size_t) 327 void* OCTREE::operator new(size_t)
328 { 328 {
329 if (!fgAllocator) 329 if (!fgAllocator)
330 { 330 {
331 fgAllocator = new G4Allocator<OCTREE>; 331 fgAllocator = new G4Allocator<OCTREE>;
332 } 332 }
333 return (void *) fgAllocator->MallocSingle( 333 return (void *) fgAllocator->MallocSingle();
334 } 334 }
335 335
336 template <OCTREE_TEMPLATE> 336 template <OCTREE_TEMPLATE>
337 void OCTREE::operator delete(void *a) 337 void OCTREE::operator delete(void *a)
338 { 338 {
339 fgAllocator->FreeSingle((OCTREE*)a); 339 fgAllocator->FreeSingle((OCTREE*)a);
340 } 340 }
341 341