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