| 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