#include <G4ParallelWorldScoringProcess.hh>
Inheritance diagram for G4ParallelWorldScoringProcess:
Public Member Functions | |
G4ParallelWorldScoringProcess (const G4String &processName="ParaWorldScore", G4ProcessType theType=fParameterisation) | |
virtual | ~G4ParallelWorldScoringProcess () |
void | SetParallelWorld (G4String parallelWorldName) |
void | SetParallelWorld (G4VPhysicalVolume *parallelWorld) |
G4bool | IsAtRestRequired (G4ParticleDefinition *partDef) |
void | StartTracking (G4Track *) |
G4double | AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *) |
G4VParticleChange * | AtRestDoIt (const G4Track &, const G4Step &) |
G4double | AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *) |
G4VParticleChange * | AlongStepDoIt (const G4Track &, const G4Step &) |
G4double | PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition) |
G4VParticleChange * | PostStepDoIt (const G4Track &, const G4Step &) |
void | Verbose (const G4Step &) const |
Definition at line 68 of file G4ParallelWorldScoringProcess.hh.
G4ParallelWorldScoringProcess::G4ParallelWorldScoringProcess | ( | const G4String & | processName = "ParaWorldScore" , |
|
G4ProcessType | theType = fParameterisation | |||
) |
Definition at line 54 of file G4ParallelWorldScoringProcess.cc.
References G4cout, G4endl, G4PathFinder::GetInstance(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4VProcess::GetProcessName(), G4TransportationManager::GetTransportationManager(), G4VProcess::pParticleChange, and G4VProcess::verboseLevel.
00055 :G4VProcess(processName,theType), fGhostNavigator(0), fNavigatorID(-1), fFieldTrack('0') 00056 { 00057 pParticleChange = &aDummyParticleChange; 00058 00059 fGhostStep = new G4Step(); 00060 fGhostPreStepPoint = fGhostStep->GetPreStepPoint(); 00061 fGhostPostStepPoint = fGhostStep->GetPostStepPoint(); 00062 00063 fTransportationManager = G4TransportationManager::GetTransportationManager(); 00064 fPathFinder = G4PathFinder::GetInstance(); 00065 00066 if (verboseLevel>0) 00067 { 00068 G4cout << GetProcessName() << " is created " << G4endl; 00069 } 00070 }
G4ParallelWorldScoringProcess::~G4ParallelWorldScoringProcess | ( | ) | [virtual] |
G4VParticleChange * G4ParallelWorldScoringProcess::AlongStepDoIt | ( | const G4Track & | , | |
const G4Step & | ||||
) | [virtual] |
Implements G4VProcess.
Definition at line 390 of file G4ParallelWorldScoringProcess.cc.
References G4VParticleChange::Initialize(), and G4VProcess::pParticleChange.
00392 { 00393 // Dummy ParticleChange ie: does nothing 00394 // Expecting G4Transportation to move the track 00395 pParticleChange->Initialize(track); 00396 return pParticleChange; 00397 }
G4double G4ParallelWorldScoringProcess::AlongStepGetPhysicalInteractionLength | ( | const G4Track & | , | |
G4double | , | |||
G4double | , | |||
G4double & | , | |||
G4GPILSelection * | ||||
) | [virtual] |
Implements G4VProcess.
Definition at line 324 of file G4ParallelWorldScoringProcess.cc.
References CandidateForSelection, G4Navigator::ComputeSafety(), G4PathFinder::ComputeStep(), DBL_MAX, G4Track::GetCurrentStepNumber(), G4Track::GetVolume(), kDoNot, kSharedOther, kSharedTransport, kUnique, NotCandidateForSelection, and G4FieldTrackUpdator::Update().
00327 { 00328 static G4FieldTrack endTrack('0'); 00329 static ELimited eLimited; 00330 00331 *selection = NotCandidateForSelection; 00332 G4double returnedStep = DBL_MAX; 00333 00334 if (previousStepSize > 0.) 00335 { fGhostSafety -= previousStepSize; } 00336 // else 00337 // { fGhostSafety = -1.; } 00338 if (fGhostSafety < 0.) fGhostSafety = 0.0; 00339 00340 // ------------------------------------------ 00341 // Determination of the proposed STEP LENGTH: 00342 // ------------------------------------------ 00343 if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.) 00344 { 00345 // I have no chance to limit 00346 returnedStep = currentMinimumStep; 00347 fOnBoundary = false; 00348 proposedSafety = fGhostSafety - currentMinimumStep; 00349 } 00350 else // (currentMinimumStep > fGhostSafety: I may limit the Step) 00351 { 00352 G4FieldTrackUpdator::Update(&fFieldTrack,&track); 00353 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00354 // ComputeStep 00355 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00356 returnedStep 00357 = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID, 00358 track.GetCurrentStepNumber(),fGhostSafety,eLimited, 00359 endTrack,track.GetVolume()); 00360 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00361 if(eLimited == kDoNot) 00362 { 00363 // Track is not on the boundary 00364 fOnBoundary = false; 00365 fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition()); 00366 } 00367 else 00368 { 00369 // Track is on the boundary 00370 fOnBoundary = true; 00371 // proposedSafety = fGhostSafety; 00372 } 00373 proposedSafety = fGhostSafety; 00374 if(eLimited == kUnique || eLimited == kSharedOther) { 00375 *selection = CandidateForSelection; 00376 }else if (eLimited == kSharedTransport) { 00377 returnedStep *= (1.0 + 1.0e-9); 00378 // Expand to disable its selection in Step Manager comparison 00379 } 00380 } 00381 00382 // ---------------------------------------------- 00383 // Returns the fGhostSafety as the proposedSafety 00384 // The SteppingManager will take care of keeping 00385 // the smallest one. 00386 // ---------------------------------------------- 00387 return returnedStep; 00388 }
G4VParticleChange * G4ParallelWorldScoringProcess::AtRestDoIt | ( | const G4Track & | , | |
const G4Step & | ||||
) | [virtual] |
Implements G4VProcess.
Definition at line 198 of file G4ParallelWorldScoringProcess.cc.
References G4StepPoint::GetSensitiveDetector(), G4StepPoint::GetTouchableHandle(), G4VSensitiveDetector::Hit(), G4VParticleChange::Initialize(), G4VProcess::pParticleChange, G4StepPoint::SetSensitiveDetector(), G4StepPoint::SetTouchableHandle(), Verbose(), and G4VProcess::verboseLevel.
00201 { 00202 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle(); 00203 G4VSensitiveDetector* aSD = 0; 00204 if(fOldGhostTouchable->GetVolume()) 00205 { aSD = fOldGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector(); } 00206 fOnBoundary = false; 00207 CopyStep(step); 00208 fGhostPreStepPoint->SetSensitiveDetector(aSD); 00209 00210 fNewGhostTouchable = fOldGhostTouchable; 00211 00212 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable); 00213 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable); 00214 if(fNewGhostTouchable->GetVolume()) 00215 { 00216 fGhostPostStepPoint->SetSensitiveDetector( 00217 fNewGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector()); 00218 } 00219 else 00220 { fGhostPostStepPoint->SetSensitiveDetector(0); } 00221 00222 if (verboseLevel>1) Verbose(step); 00223 00224 G4VSensitiveDetector* sd = fGhostPreStepPoint->GetSensitiveDetector(); 00225 if(sd) 00226 { 00227 sd->Hit(fGhostStep); 00228 } 00229 00230 pParticleChange->Initialize(track); 00231 return pParticleChange; 00232 }
G4double G4ParallelWorldScoringProcess::AtRestGetPhysicalInteractionLength | ( | const G4Track & | , | |
G4ForceCondition * | ||||
) | [virtual] |
Implements G4VProcess.
Definition at line 185 of file G4ParallelWorldScoringProcess.cc.
References DBL_MAX, and Forced.
G4bool G4ParallelWorldScoringProcess::IsAtRestRequired | ( | G4ParticleDefinition * | partDef | ) |
Definition at line 110 of file G4ParallelWorldScoringProcess.cc.
References G4ParticleDefinition::GetParticleName(), and G4ParticleDefinition::GetPDGEncoding().
00111 { 00112 G4int pdgCode = partDef->GetPDGEncoding(); 00113 if(pdgCode==0) 00114 { 00115 G4String partName = partDef->GetParticleName(); 00116 if(partName=="opticalphoton") return false; 00117 if(partName=="geantino") return false; 00118 if(partName=="chargedgeantino") return false; 00119 } 00120 else 00121 { 00122 if(pdgCode==22) return false; // gamma 00123 if(pdgCode==11) return false; // electron 00124 if(pdgCode==2212) return false; // proton 00125 if(pdgCode==-12) return false; // anti_nu_e 00126 if(pdgCode==12) return false; // nu_e 00127 if(pdgCode==-14) return false; // anti_nu_mu 00128 if(pdgCode==14) return false; // nu_mu 00129 if(pdgCode==-16) return false; // anti_nu_tau 00130 if(pdgCode==16) return false; // nu_tau 00131 } 00132 return true; 00133 }
G4VParticleChange * G4ParallelWorldScoringProcess::PostStepDoIt | ( | const G4Track & | , | |
const G4Step & | ||||
) | [virtual] |
Implements G4VProcess.
Definition at line 255 of file G4ParallelWorldScoringProcess.cc.
References G4PathFinder::CreateTouchableHandle(), G4StepPoint::GetSensitiveDetector(), G4StepPoint::GetTouchableHandle(), G4VSensitiveDetector::Hit(), G4VParticleChange::Initialize(), G4VProcess::pParticleChange, G4StepPoint::SetSensitiveDetector(), G4StepPoint::SetTouchableHandle(), Verbose(), and G4VProcess::verboseLevel.
00258 { 00259 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle(); 00260 G4VSensitiveDetector* aSD = 0; 00261 if(fOldGhostTouchable->GetVolume()) 00262 { aSD = fOldGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector(); } 00263 CopyStep(step); 00264 fGhostPreStepPoint->SetSensitiveDetector(aSD); 00265 00266 // fPathFinder->Locate( track.GetPosition(), 00267 // track.GetMomentumDirection(), 00268 // true); 00269 00270 // fPathFinder->Locate(step.GetPostStepPoint()->GetPosition(), 00271 // step.GetPostStepPoint()->GetMomentumDirection()); 00272 00273 if(fOnBoundary) 00274 { 00275 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00276 // Locate the point and get new touchable 00277 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00278 //?? fPathFinder->Locate(step.GetPostStepPoint()->GetPosition(), 00279 //?? step.GetPostStepPoint()->GetMomentumDirection()); 00280 fNewGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID); 00281 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00282 } 00283 else 00284 { 00285 // Do I need this ?????????????????????????????????????????????????????????? 00286 // fGhostNavigator->LocateGlobalPointWithinVolume(track.GetPosition()); 00287 // ????????????????????????????????????????????????????????????????????????? 00288 00289 // fPathFinder->ReLocate(track.GetPosition()); 00290 00291 // reuse the touchable 00292 fNewGhostTouchable = fOldGhostTouchable; 00293 } 00294 00295 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable); 00296 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable); 00297 00298 if(fNewGhostTouchable->GetVolume()) 00299 { 00300 fGhostPostStepPoint->SetSensitiveDetector( 00301 fNewGhostTouchable->GetVolume()->GetLogicalVolume()->GetSensitiveDetector()); 00302 } 00303 else 00304 { fGhostPostStepPoint->SetSensitiveDetector(0); } 00305 00306 if (verboseLevel>1) Verbose(step); 00307 00308 G4VSensitiveDetector* sd = fGhostPreStepPoint->GetSensitiveDetector(); 00309 if(sd) 00310 { 00311 sd->Hit(fGhostStep); 00312 } 00313 00314 pParticleChange->Initialize(track); // Does not change the track properties 00315 return pParticleChange; 00316 }
G4double G4ParallelWorldScoringProcess::PostStepGetPhysicalInteractionLength | ( | const G4Track & | track, | |
G4double | previousStepSize, | |||
G4ForceCondition * | condition | |||
) | [virtual] |
Implements G4VProcess.
Definition at line 240 of file G4ParallelWorldScoringProcess.cc.
References DBL_MAX, and StronglyForced.
00244 { 00245 // I must be invoked anyway to score the hit. 00246 *condition = StronglyForced; 00247 return DBL_MAX; 00248 }
void G4ParallelWorldScoringProcess::SetParallelWorld | ( | G4VPhysicalVolume * | parallelWorld | ) |
Definition at line 98 of file G4ParallelWorldScoringProcess.cc.
References G4VPhysicalVolume::GetName(), and G4TransportationManager::GetNavigator().
00099 { 00100 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00101 // Get pointer of navigator 00102 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00103 fGhostWorldName = parallelWorld->GetName(); 00104 fGhostWorld = parallelWorld; 00105 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld); 00106 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00107 }
void G4ParallelWorldScoringProcess::SetParallelWorld | ( | G4String | parallelWorldName | ) |
Definition at line 86 of file G4ParallelWorldScoringProcess.cc.
References G4TransportationManager::GetNavigator(), and G4TransportationManager::GetParallelWorld().
Referenced by G4RunManager::ConstructScoringWorlds().
00087 { 00088 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00089 // Get pointers of the parallel world and its navigator 00090 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00091 fGhostWorldName = parallelWorldName; 00092 fGhostWorld = fTransportationManager->GetParallelWorld(fGhostWorldName); 00093 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld); 00094 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00095 }
void G4ParallelWorldScoringProcess::StartTracking | ( | G4Track * | ) | [virtual] |
Reimplemented from G4VProcess.
Definition at line 141 of file G4ParallelWorldScoringProcess.cc.
References G4TransportationManager::ActivateNavigator(), G4PathFinder::CreateTouchableHandle(), FatalException, fUndefined, G4Exception(), G4Track::GetMomentumDirection(), G4Track::GetPosition(), G4PathFinder::PrepareNewTrack(), G4StepPoint::SetStepStatus(), and G4StepPoint::SetTouchableHandle().
00142 { 00143 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00144 // Activate navigator and get the navigator ID 00145 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00146 // G4cout << " G4ParallelWorldScoringProcess::StartTracking" << G4endl; 00147 if(fGhostNavigator) 00148 { fNavigatorID = fTransportationManager->ActivateNavigator(fGhostNavigator); } 00149 else 00150 { 00151 G4Exception("G4ParallelWorldScoringProcess::StartTracking", 00152 "ProcParaWorld000",FatalException, 00153 "G4ParallelWorldScoringProcess is used for tracking without having a parallel world assigned"); 00154 } 00155 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00156 00157 // G4cout << "G4ParallelWorldScoringProcess::StartTracking <<<<<<<<<<<<<<<<<< " << G4endl; 00158 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00159 // Let PathFinder initialize 00160 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00161 fPathFinder->PrepareNewTrack(trk->GetPosition(),trk->GetMomentumDirection()); 00162 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00163 00164 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00165 // Setup initial touchables for the first step 00166 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 00167 fOldGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID); 00168 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable); 00169 fNewGhostTouchable = fOldGhostTouchable; 00170 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable); 00171 00172 // Initialize 00173 fGhostSafety = -1.; 00174 fOnBoundary = false; 00175 fGhostPreStepPoint->SetStepStatus(fUndefined); 00176 fGhostPostStepPoint->SetStepStatus(fUndefined); 00177 }
void G4ParallelWorldScoringProcess::Verbose | ( | const G4Step & | ) | const |
Definition at line 424 of file G4ParallelWorldScoringProcess.cc.
References G4cout, G4endl, G4Track::GetMomentumDirection(), G4VPhysicalVolume::GetName(), G4StepPoint::GetPhysicalVolume(), G4StepPoint::GetPosition(), G4Step::GetPostStepPoint(), G4Step::GetPreStepPoint(), G4StepPoint::GetProcessDefinedStep(), G4VProcess::GetProcessName(), G4VTouchable::GetReplicaNumber(), G4Step::GetStepLength(), G4Step::GetTotalEnergyDeposit(), G4StepPoint::GetTouchable(), and G4Step::GetTrack().
Referenced by AtRestDoIt(), and PostStepDoIt().
00425 { 00426 G4cout << "In mass geometry ------------------------------------------------" << G4endl; 00427 G4cout << " StepLength : " << step.GetStepLength()/mm << " TotalEnergyDeposit : " 00428 << step.GetTotalEnergyDeposit()/MeV << G4endl; 00429 G4cout << " PreStepPoint : " 00430 << step.GetPreStepPoint()->GetPhysicalVolume()->GetName() << " - "; 00431 if(step.GetPreStepPoint()->GetProcessDefinedStep()) 00432 { G4cout << step.GetPreStepPoint()->GetProcessDefinedStep()->GetProcessName(); } 00433 else 00434 { G4cout << "NoProcessAssigned"; } 00435 G4cout << G4endl; 00436 G4cout << " " << step.GetPreStepPoint()->GetPosition() << G4endl; 00437 G4cout << " PostStepPoint : "; 00438 if(step.GetPostStepPoint()->GetPhysicalVolume()) 00439 { G4cout << step.GetPostStepPoint()->GetPhysicalVolume()->GetName(); } 00440 else 00441 { G4cout << "OutOfWorld"; } 00442 G4cout << " - "; 00443 if(step.GetPostStepPoint()->GetProcessDefinedStep()) 00444 { G4cout << step.GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName(); } 00445 else 00446 { G4cout << "NoProcessAssigned"; } 00447 G4cout << G4endl; 00448 G4cout << " " << step.GetPostStepPoint()->GetPosition() << G4endl; 00449 00450 G4cout << "In ghost geometry ------------------------------------------------" << G4endl; 00451 G4cout << " StepLength : " << fGhostStep->GetStepLength()/mm 00452 << " TotalEnergyDeposit : " 00453 << fGhostStep->GetTotalEnergyDeposit()/MeV << G4endl; 00454 G4cout << " PreStepPoint : " 00455 << fGhostStep->GetPreStepPoint()->GetPhysicalVolume()->GetName() << " [" 00456 << fGhostStep->GetPreStepPoint()->GetTouchable()->GetReplicaNumber() 00457 << " ]" << " - "; 00458 if(fGhostStep->GetPreStepPoint()->GetProcessDefinedStep()) 00459 { G4cout << fGhostStep->GetPreStepPoint()->GetProcessDefinedStep()->GetProcessName(); } 00460 else 00461 { G4cout << "NoProcessAssigned"; } 00462 G4cout << G4endl; 00463 G4cout << " " << fGhostStep->GetPreStepPoint()->GetPosition() << G4endl; 00464 G4cout << " PostStepPoint : "; 00465 if(fGhostStep->GetPostStepPoint()->GetPhysicalVolume()) 00466 { 00467 G4cout << fGhostStep->GetPostStepPoint()->GetPhysicalVolume()->GetName() << " [" 00468 << fGhostStep->GetPostStepPoint()->GetTouchable()->GetReplicaNumber() 00469 << " ]"; 00470 } 00471 else 00472 { G4cout << "OutOfWorld"; } 00473 G4cout << " - "; 00474 if(fGhostStep->GetPostStepPoint()->GetProcessDefinedStep()) 00475 { G4cout << fGhostStep->GetPostStepPoint()->GetProcessDefinedStep()->GetProcessName(); } 00476 else 00477 { G4cout << "NoProcessAssigned"; } 00478 G4cout << G4endl; 00479 G4cout << " " << fGhostStep->GetPostStepPoint()->GetPosition() << " == " 00480 << fGhostStep->GetTrack()->GetMomentumDirection() 00481 << G4endl; 00482 00483 }