Geant4 Cross Reference |
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: Mathieu Karamitros (kara (AT) cenbg . in2p3 . fr) 28 // 29 // WARNING : This class is released as a prototype. 30 // It might strongly evolve or even disapear in the next releases. 31 // 32 // History: 33 // ----------- 34 // 10 Oct 2011 M.Karamitros created 35 // 36 // ------------------------------------------------------------------- 37 38 TEMPLATE 39 G4ThreadLocal G4ITMANAGER* G4ITMANAGER::fInstance = nullptr; 40 41 TEMPLATE 42 G4ITMANAGER* G4ITMANAGER::Instance() 43 { 44 if(!fInstance) 45 { 46 fInstance = new G4ITFinder(); 47 } 48 return fInstance; 49 } 50 51 TEMPLATE G4ITMANAGER::G4ITFinder() 52 : G4VITFinder() 53 { 54 G4AllITFinder::Instance()->RegisterManager(this); 55 } 56 57 TEMPLATE G4ITMANAGER::~G4ITFinder() 58 { 59 for(auto & it : fTree) 60 { 61 delete it.second; 62 } 63 fTree.clear(); 64 fInstance = nullptr; 65 } 66 67 TEMPLATE 68 void G4ITMANAGER::Clear() 69 { 70 for(auto & it : fTree) 71 { 72 delete it.second; 73 } 74 fTree.clear(); 75 } 76 77 TEMPLATE 78 void G4ITMANAGER::Push(G4Track* track) 79 { 80 /* 81 * If you want to use this class with another type than G4Molecule 82 * inheriting as well from G4IT, replace T::GetMolecule by GetIT 83 * and aIT->GetMoleculeID by GetSubITType 84 */ 85 T* aIT = T::GetMolecule(track); 86 aIT->RecordCurrentPositionNTime(); 87 88 G4int key = aIT->GetMoleculeID(); 89 90 if(!(aIT->GetNode())) 91 { 92 G4KDNode_Base* node; 93 94 auto it_fTree = fTree.find(key); 95 96 if(it_fTree != fTree.end()) 97 { 98 node = it_fTree->second->Insert(aIT); 99 } 100 else 101 { 102 auto aTree = new G4KDTree(); 103 fTree.insert(std::make_pair(key, aTree)); 104 node = aTree->Insert(aIT); 105 } 106 aIT->SetNode(node); 107 } 108 } 109 110 TEMPLATE 111 G4KDTreeResultHandle G4ITMANAGER::FindNearest(const G4ThreeVector& position, 112 G4int key) 113 { 114 auto it = fTree.find(key); 115 if(it != fTree.end()) 116 { 117 return it->second->Nearest(position); 118 } 119 120 return nullptr; 121 } 122 123 TEMPLATE 124 G4KDTreeResultHandle G4ITMANAGER::FindNearest(const T* point0, G4int key) 125 { 126 if(G4int(*point0) == key) 127 { 128 G4KDNode_Base* node0 = point0->GetNode(); 129 if(node0 == nullptr) 130 { 131 G4ExceptionDescription exceptionDescription( 132 "Bad request : no node found in the IT you are searching " 133 "closest neighbourg for"); 134 G4Exception("G4ITManager::FindNearest", "ITManager002", 135 FatalErrorInArgument, exceptionDescription); 136 return nullptr; // coverity 137 } 138 auto it = fTree.find(key); 139 if(it != fTree.end()) 140 { 141 G4KDTreeResultHandle output(it->second->Nearest(node0)); 142 if(!output) 143 { 144 return nullptr; 145 } 146 return output; 147 } 148 149 return nullptr; 150 } 151 152 auto it = fTree.find(key); 153 if(it != fTree.end()) 154 { 155 G4KDTreeResultHandle output(it->second->Nearest(*point0)); 156 if(!output) 157 { 158 return nullptr; 159 } 160 return output; 161 } 162 163 return nullptr; 164 } 165 166 TEMPLATE 167 G4KDTreeResultHandle G4ITMANAGER::FindNearest(const T* source, const T* type) 168 { 169 return FindNearest(source, G4int(*type)); 170 } 171 172 TEMPLATE 173 G4KDTreeResultHandle G4ITMANAGER::FindNearestInRange( 174 const G4ThreeVector& position, G4int key, G4double R) 175 { 176 auto it = fTree.find(key); 177 if(it != fTree.end()) 178 { 179 return it->second->NearestInRange(position, R); 180 } 181 182 return nullptr; 183 } 184 185 TEMPLATE 186 G4KDTreeResultHandle G4ITMANAGER::FindNearestInRange(const T* point0, G4int key, 187 G4double R) 188 { 189 if(point0->GetMoleculeID() == key) 190 { 191 G4KDNode_Base* node0 = point0->GetNode(); 192 auto it = fTree.find(key); 193 if(it != fTree.end()) 194 return it->second->NearestInRange(node0, R); 195 196 return nullptr; 197 } 198 199 auto it = fTree.find(key); 200 if(it != fTree.end()) 201 return it->second->NearestInRange(*point0, R); 202 203 return nullptr; 204 } 205 206 //#define DEBUG_MEM 207 208 #ifdef DEBUG_MEM 209 # include "G4MemStat.hh" 210 using namespace G4MemStat; 211 #endif 212 TEMPLATE 213 void G4ITMANAGER::UpdatePositionMap() 214 { 215 #if defined(DEBUG_MEM) 216 MemStat mem_first, mem_second, mem_diff; 217 #endif 218 219 #if defined(DEBUG_MEM) 220 mem_first = MemoryUsage(); 221 #endif 222 Clear(); 223 const std::map<G4int, PriorityList*>& listMap = 224 G4ITTrackHolder::Instance()->GetLists(); 225 auto it = listMap.begin(); 226 auto end = listMap.end(); 227 228 for(; it != end; it++) 229 { 230 G4int key = it->first; 231 PriorityList* listUnion = it->second; 232 if(listUnion == nullptr || listUnion->GetMainList() == nullptr || 233 listUnion->GetMainList()->empty()) 234 { 235 #if defined(DEBUG_MEM) 236 mem_second = MemoryUsage(); 237 mem_diff = mem_second - mem_first; 238 G4cout << "\t || MEM || G4ITMANAGER::UpdatePositionMap || " 239 "In if(currentBox->Empty()), diff is : " 240 << mem_diff << G4endl; 241 #endif 242 243 continue; 244 } 245 246 auto currentTree = new G4KDTree(); 247 fTree[key] = currentTree; 248 #if defined(DEBUG_MEM) 249 mem_second = MemoryUsage(); 250 mem_diff = mem_second - mem_first; 251 G4cout << "\t || MEM || G4ITMANAGER::UpdatePositionMap || " 252 "after creating tree, diff is : " 253 << mem_diff << G4endl; 254 #endif 255 256 G4TrackList* trackList = listUnion->GetMainList(); 257 G4TrackList::iterator __it = trackList->begin(); 258 G4TrackList::iterator __end = trackList->end(); 259 260 for(; __it != __end; __it++) 261 { 262 G4IT* currentIT = GetIT(*__it); 263 G4KDNode_Base* currentNode = currentTree->Insert(currentIT); 264 currentIT->SetNode(currentNode); 265 266 #if defined(DEBUG_MEM) 267 mem_second = MemoryUsage(); 268 mem_diff = mem_second - mem_first; 269 G4cout << "\t || MEM || G4ITMANAGER::UpdatePositionMap || " 270 "after currentIT->SetNode(currentNode), diff is : " 271 << mem_diff << G4endl; 272 #endif 273 } 274 275 #if defined(DEBUG_MEM) 276 mem_second = MemoryUsage(); 277 mem_diff = mem_second - mem_first; 278 G4cout << "\t || MEM || G4ITMANAGER::UpdatePositionMap || " 279 "In else{...}, diff is : " 280 << mem_diff << G4endl; 281 #endif 282 283 } 284 } 285