Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/dna/management/include/G4Octree.icc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

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