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 // 28 // 29 30 #include "G4ParallelWorldProcess.hh" 31 32 #include "G4FieldTrackUpdator.hh" 33 #include "G4Material.hh" 34 #include "G4Navigator.hh" 35 #include "G4ParallelWorldProcessStore.hh" 36 #include "G4ParticleChange.hh" 37 #include "G4PathFinder.hh" 38 #include "G4ProductionCuts.hh" 39 #include "G4ProductionCutsTable.hh" 40 #include "G4SDManager.hh" 41 #include "G4Step.hh" 42 #include "G4StepPoint.hh" 43 #include "G4TransportationManager.hh" 44 #include "G4TransportationProcessType.hh" 45 #include "G4VPhysicalVolume.hh" 46 #include "G4VSensitiveDetector.hh" 47 #include "G4VTouchable.hh" 48 #include "G4ios.hh" 49 50 G4ThreadLocal G4Step* G4ParallelWorldProcess::fpHyperStep = nullptr; 51 G4ThreadLocal G4int G4ParallelWorldProcess::nParallelWorlds = 0; 52 G4ThreadLocal G4int G4ParallelWorldProcess::fNavIDHyp = 0; 53 54 const G4Step* G4ParallelWorldProcess::GetHyperStep() 55 { 56 return fpHyperStep; 57 } 58 59 G4int G4ParallelWorldProcess::GetHypNavigatorID() 60 { 61 return fNavIDHyp; 62 } 63 64 G4ParallelWorldProcess::G4ParallelWorldProcess(const G4String& processName, G4ProcessType theType) 65 : G4VProcess(processName, theType), fFieldTrack('0') 66 { 67 SetProcessSubType(PARALLEL_WORLD_PROCESS); 68 if (fpHyperStep == nullptr) fpHyperStep = new G4Step(); 69 iParallelWorld = ++nParallelWorlds; 70 71 pParticleChange = &aDummyParticleChange; 72 73 fGhostStep = new G4Step(); 74 fGhostPreStepPoint = fGhostStep->GetPreStepPoint(); 75 fGhostPostStepPoint = fGhostStep->GetPostStepPoint(); 76 77 fTransportationManager = G4TransportationManager::GetTransportationManager(); 78 fTransportationManager->GetNavigatorForTracking()->SetPushVerbosity(false); 79 fPathFinder = G4PathFinder::GetInstance(); 80 81 fGhostWorldName = "** NotDefined **"; 82 G4ParallelWorldProcessStore::GetInstance()->SetParallelWorld(this, processName); 83 84 if (verboseLevel > 0) { 85 G4cout << GetProcessName() << " is created " << G4endl; 86 } 87 } 88 89 G4ParallelWorldProcess::~G4ParallelWorldProcess() 90 { 91 delete fGhostStep; 92 nParallelWorlds--; 93 if (nParallelWorlds == 0) { 94 delete fpHyperStep; 95 fpHyperStep = nullptr; 96 } 97 } 98 99 void G4ParallelWorldProcess::SetParallelWorld(G4String parallelWorldName) 100 { 101 fGhostWorldName = parallelWorldName; 102 fGhostWorld = fTransportationManager->GetParallelWorld(fGhostWorldName); 103 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld); 104 fGhostNavigator->SetPushVerbosity(false); 105 } 106 107 void G4ParallelWorldProcess::SetParallelWorld(G4VPhysicalVolume* parallelWorld) 108 { 109 fGhostWorldName = parallelWorld->GetName(); 110 fGhostWorld = parallelWorld; 111 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld); 112 fGhostNavigator->SetPushVerbosity(false); 113 } 114 115 void G4ParallelWorldProcess::StartTracking(G4Track* trk) 116 { 117 if (fGhostNavigator != nullptr) { 118 fNavigatorID = fTransportationManager->ActivateNavigator(fGhostNavigator); 119 } 120 else { 121 G4Exception( 122 "G4ParallelWorldProcess::StartTracking", "ProcParaWorld000", FatalException, 123 "G4ParallelWorldProcess is used for tracking without having a parallel world assigned"); 124 } 125 fPathFinder->PrepareNewTrack(trk->GetPosition(), trk->GetMomentumDirection()); 126 127 fOldGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID); 128 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable); 129 fNewGhostTouchable = fOldGhostTouchable; 130 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable); 131 132 fGhostSafety = -1.; 133 fOnBoundary = false; 134 fGhostPreStepPoint->SetStepStatus(fUndefined); 135 fGhostPostStepPoint->SetStepStatus(fUndefined); 136 137 *(fpHyperStep->GetPostStepPoint()) = *(trk->GetStep()->GetPostStepPoint()); 138 if (layeredMaterialFlag) { 139 G4StepPoint* realWorldPostStepPoint = trk->GetStep()->GetPostStepPoint(); 140 SwitchMaterial(realWorldPostStepPoint); 141 G4StepPoint* realWorldPreStepPoint = trk->GetStep()->GetPreStepPoint(); 142 SwitchMaterial(realWorldPreStepPoint); 143 G4double velocity = trk->CalculateVelocity(); 144 realWorldPostStepPoint->SetVelocity(velocity); 145 realWorldPreStepPoint->SetVelocity(velocity); 146 trk->SetVelocity(velocity); 147 } 148 *(fpHyperStep->GetPreStepPoint()) = *(fpHyperStep->GetPostStepPoint()); 149 } 150 151 G4double G4ParallelWorldProcess::AtRestGetPhysicalInteractionLength(const G4Track& /*track*/, 152 G4ForceCondition* condition) 153 { 154 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 155 // At Rest must be registered ONLY for the particle which has other At Rest 156 // process(es). 157 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 158 *condition = Forced; 159 return DBL_MAX; 160 } 161 162 G4VParticleChange* G4ParallelWorldProcess::AtRestDoIt(const G4Track& track, const G4Step& step) 163 { 164 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 165 // At Rest must be registered ONLY for the particle which has other At Rest 166 // process(es). 167 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 168 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle(); 169 G4VSensitiveDetector* aSD = nullptr; 170 if (fOldGhostTouchable->GetVolume() != nullptr) { 171 aSD = fOldGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector(); 172 } 173 fOnBoundary = false; 174 if (aSD != nullptr) { 175 CopyStep(step); 176 fGhostPreStepPoint->SetSensitiveDetector(aSD); 177 178 fNewGhostTouchable = fOldGhostTouchable; 179 180 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable); 181 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable); 182 if (fNewGhostTouchable->GetVolume() != nullptr) { 183 fGhostPostStepPoint->SetSensitiveDetector( 184 fNewGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector()); 185 } 186 else { 187 fGhostPostStepPoint->SetSensitiveDetector(nullptr); 188 } 189 190 aSD->Hit(fGhostStep); 191 } 192 193 pParticleChange->Initialize(track); 194 return pParticleChange; 195 } 196 197 G4double G4ParallelWorldProcess::PostStepGetPhysicalInteractionLength(const G4Track& /*track*/, 198 G4double /*previousStepSize*/, 199 G4ForceCondition* condition) 200 { 201 *condition = StronglyForced; 202 return DBL_MAX; 203 } 204 205 G4VParticleChange* G4ParallelWorldProcess::PostStepDoIt(const G4Track& track, const G4Step& step) 206 { 207 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle(); 208 G4VSensitiveDetector* aSD = nullptr; 209 if (fOldGhostTouchable->GetVolume() != nullptr) { 210 aSD = fOldGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector(); 211 } 212 CopyStep(step); 213 fGhostPreStepPoint->SetSensitiveDetector(aSD); 214 215 if (fOnBoundary) { 216 fNewGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID); 217 } 218 else { 219 fNewGhostTouchable = fOldGhostTouchable; 220 } 221 222 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable); 223 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable); 224 225 if (fNewGhostTouchable->GetVolume() != nullptr) { 226 fGhostPostStepPoint->SetSensitiveDetector( 227 fNewGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector()); 228 } 229 else { 230 fGhostPostStepPoint->SetSensitiveDetector(nullptr); 231 } 232 233 G4VSensitiveDetector* sd = fGhostPreStepPoint->GetSensitiveDetector(); 234 if (sd != nullptr) { 235 sd->Hit(fGhostStep); 236 } 237 238 pParticleChange->Initialize(track); 239 if (layeredMaterialFlag) { 240 G4StepPoint* realWorldPostStepPoint = ((G4Step*)(track.GetStep()))->GetPostStepPoint(); 241 SwitchMaterial(realWorldPostStepPoint); 242 } 243 return pParticleChange; 244 } 245 246 G4double G4ParallelWorldProcess::AlongStepGetPhysicalInteractionLength(const G4Track& track, 247 G4double previousStepSize, 248 G4double currentMinimumStep, 249 G4double& proposedSafety, 250 G4GPILSelection* selection) 251 { 252 static G4ThreadLocal G4FieldTrack* endTrack_G4MT_TLS_ = nullptr; 253 if (endTrack_G4MT_TLS_ == nullptr) endTrack_G4MT_TLS_ = new G4FieldTrack('0'); 254 G4FieldTrack& endTrack = *endTrack_G4MT_TLS_; 255 // static ELimited eLimited; 256 ELimited eLimited; 257 ELimited eLim = kUndefLimited; 258 259 *selection = NotCandidateForSelection; 260 G4double returnedStep = DBL_MAX; 261 262 if (previousStepSize > 0.) { 263 fGhostSafety -= previousStepSize; 264 } 265 if (fGhostSafety < 0.) fGhostSafety = 0.0; 266 267 if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.) { 268 // I have no chance to limit 269 returnedStep = currentMinimumStep; 270 fOnBoundary = false; 271 proposedSafety = fGhostSafety - currentMinimumStep; 272 eLim = kDoNot; 273 } 274 else { 275 G4FieldTrackUpdator::Update(&fFieldTrack, &track); 276 277 #ifdef G4DEBUG_PARALLEL_WORLD_PROCESS 278 if (verboseLevel > 0) { 279 int localVerb = verboseLevel - 1; 280 281 if (localVerb == 1) { 282 G4cout << " Pll Wrl proc::AlongStepGPIL " << this->GetProcessName() << G4endl; 283 } 284 else if (localVerb > 1) { 285 G4cout << "----------------------------------------------" << G4endl; 286 G4cout << " ParallelWorldProcess: field Track set to : " << G4endl; 287 G4cout << "----------------------------------------------" << G4endl; 288 G4cout << fFieldTrack << G4endl; 289 G4cout << "----------------------------------------------" << G4endl; 290 } 291 } 292 #endif 293 294 returnedStep = fPathFinder->ComputeStep(fFieldTrack, currentMinimumStep, fNavigatorID, 295 track.GetCurrentStepNumber(), fGhostSafety, eLimited, 296 endTrack, track.GetVolume()); 297 if (eLimited == kDoNot) { 298 fOnBoundary = false; 299 fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition()); 300 } 301 else { 302 fOnBoundary = true; 303 // fGhostSafetyEnd = 0.0; // At end-point of expected step only 304 } 305 proposedSafety = fGhostSafety; 306 if (eLimited == kUnique || eLimited == kSharedOther) { 307 *selection = CandidateForSelection; 308 } 309 else if (eLimited == kSharedTransport) { 310 returnedStep *= (1.0 + 1.0e-9); 311 } 312 eLim = eLimited; 313 } 314 315 if (iParallelWorld == nParallelWorlds) fNavIDHyp = 0; 316 if (eLim == kUnique || eLim == kSharedOther) fNavIDHyp = fNavigatorID; 317 return returnedStep; 318 } 319 320 G4VParticleChange* G4ParallelWorldProcess::AlongStepDoIt(const G4Track& track, const G4Step&) 321 { 322 pParticleChange->Initialize(track); 323 return pParticleChange; 324 } 325 326 void G4ParallelWorldProcess::CopyStep(const G4Step& step) 327 { 328 G4StepStatus prevStat = fGhostPostStepPoint->GetStepStatus(); 329 330 fGhostStep->SetTrack(step.GetTrack()); 331 fGhostStep->SetStepLength(step.GetStepLength()); 332 fGhostStep->SetTotalEnergyDeposit(step.GetTotalEnergyDeposit()); 333 fGhostStep->SetNonIonizingEnergyDeposit(step.GetNonIonizingEnergyDeposit()); 334 fGhostStep->SetControlFlag(step.GetControlFlag()); 335 fGhostStep->SetSecondary((const_cast<G4Step&>(step)).GetfSecondary()); 336 337 *fGhostPreStepPoint = *(step.GetPreStepPoint()); 338 *fGhostPostStepPoint = *(step.GetPostStepPoint()); 339 340 fGhostPreStepPoint->SetStepStatus(prevStat); 341 if (fOnBoundary) { 342 fGhostPostStepPoint->SetStepStatus(fGeomBoundary); 343 } 344 else if (fGhostPostStepPoint->GetStepStatus() == fGeomBoundary) { 345 fGhostPostStepPoint->SetStepStatus(fPostStepDoItProc); 346 } 347 348 if (iParallelWorld == 1) { 349 G4StepStatus prevStatHyp = fpHyperStep->GetPostStepPoint()->GetStepStatus(); 350 351 fpHyperStep->SetTrack(step.GetTrack()); 352 fpHyperStep->SetStepLength(step.GetStepLength()); 353 fpHyperStep->SetTotalEnergyDeposit(step.GetTotalEnergyDeposit()); 354 fpHyperStep->SetNonIonizingEnergyDeposit(step.GetNonIonizingEnergyDeposit()); 355 fpHyperStep->SetControlFlag(step.GetControlFlag()); 356 357 *(fpHyperStep->GetPreStepPoint()) = *(fpHyperStep->GetPostStepPoint()); 358 *(fpHyperStep->GetPostStepPoint()) = *(step.GetPostStepPoint()); 359 360 fpHyperStep->GetPreStepPoint()->SetStepStatus(prevStatHyp); 361 } 362 363 if (fOnBoundary) { 364 fpHyperStep->GetPostStepPoint()->SetStepStatus(fGeomBoundary); 365 } 366 } 367 368 void G4ParallelWorldProcess::SwitchMaterial(G4StepPoint* realWorldStepPoint) 369 { 370 if (realWorldStepPoint->GetStepStatus() == fWorldBoundary) return; 371 G4VPhysicalVolume* thePhys = fNewGhostTouchable->GetVolume(); 372 if (thePhys != nullptr) { 373 G4Material* ghostMaterial = thePhys->GetLogicalVolume()->GetMaterial(); 374 if (ghostMaterial != nullptr) { 375 G4Region* ghostRegion = thePhys->GetLogicalVolume()->GetRegion(); 376 G4ProductionCuts* prodCuts = realWorldStepPoint->GetMaterialCutsCouple()->GetProductionCuts(); 377 if (ghostRegion != nullptr) { 378 G4ProductionCuts* ghostProdCuts = ghostRegion->GetProductionCuts(); 379 if (ghostProdCuts != nullptr) prodCuts = ghostProdCuts; 380 } 381 const G4MaterialCutsCouple* ghostMCCouple = 382 G4ProductionCutsTable::GetProductionCutsTable()->GetMaterialCutsCouple(ghostMaterial, 383 prodCuts); 384 if (ghostMCCouple != nullptr) { 385 realWorldStepPoint->SetMaterial(ghostMaterial); 386 realWorldStepPoint->SetMaterialCutsCouple(ghostMCCouple); 387 *(fpHyperStep->GetPostStepPoint()) = *(fGhostPostStepPoint); 388 fpHyperStep->GetPostStepPoint()->SetMaterial(ghostMaterial); 389 fpHyperStep->GetPostStepPoint()->SetMaterialCutsCouple(ghostMCCouple); 390 } 391 else { 392 G4cout << "!!! MaterialCutsCouple is not found for " << ghostMaterial->GetName() << "." 393 << G4endl << " Material in real world (" 394 << realWorldStepPoint->GetMaterial()->GetName() << ") is used." << G4endl; 395 } 396 } 397 } 398 } 399 400 G4bool G4ParallelWorldProcess::IsAtRestRequired(G4ParticleDefinition* partDef) 401 { 402 G4int pdgCode = partDef->GetPDGEncoding(); 403 if (pdgCode == 0) { 404 G4String partName = partDef->GetParticleName(); 405 if (partName == "geantino") return false; 406 if (partName == "chargedgeantino") return false; 407 } 408 else { 409 if (pdgCode == 11 || pdgCode == 2212) return false; // electrons and proton 410 pdgCode = std::abs(pdgCode); 411 if (pdgCode == 22) return false; // gamma and optical photons 412 if (pdgCode == 12 || pdgCode == 14 || pdgCode == 16) return false; // all neutronos 413 } 414 return true; 415 } 416