| 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 // Author: Mathieu Karamitros (kara (AT) cenbg
28 //
29 // History:
30 // -----------
31 // 10 Oct 2011 M.Karamitros created
32 //
33 // -------------------------------------------
34
35 #include "G4ITStepProcessor.hh"
36
37 #include "G4ForceCondition.hh"
38 #include "G4GPILSelection.hh"
39 #include "G4GeometryTolerance.hh"
40 #include "G4IT.hh"
41 #include "G4ITNavigator.hh" // Inc
42 #include "G4ITSteppingVerbose.hh"
43 #include "G4ITTrackHolder.hh"
44 #include "G4ITTrackingInteractivity.hh"
45 #include "G4ITTrackingManager.hh"
46 #include "G4ITTransportation.hh"
47 #include "G4ITTransportationManager.hh"
48 #include "G4ParticleTable.hh"
49 #include "G4TrackingInformation.hh"
50 #include "G4UImanager.hh"
51 #include "G4VITProcess.hh"
52 #include "G4VITSteppingVerbose.hh"
53 #include "G4VProcess.hh"
54
55 #include <iomanip> // Include fro
56 #include <vector> // Include fro
57
58 using namespace std;
59
60 static const std::size_t SizeOfSelectedDoItVec
61
62 template<typename T>
63 inline bool IsInf(T value)
64 {
65 return std::numeric_limits<T>::has_infinit
66 && value == std::numeric_limits<T>::in
67 }
68
69 //____________________________________________
70
71 G4ITStepProcessor::G4ITStepProcessor()
72 {
73 fpVerbose = nullptr;
74 // fpUserSteppingAction = 0 ;
75 fStoreTrajectory = 0;
76 fpTrackingManager = nullptr;
77 fpNavigator = nullptr;
78 kCarTolerance = -1.;
79 fInitialized = false;
80 fPreviousTimeStep = DBL_MAX;
81 fILTimeStep = DBL_MAX;
82 fpTrackContainer = nullptr;
83
84 CleanProcessor();
85 ResetSecondaries();
86 }
87
88 //____________________________________________
89
90 //G4ITStepProcessor::
91 G4ITStepProcessorState::G4ITStepProcessorState
92 G4ITStepProcessorState_Lock(),
93 fSelectedAtRestDoItVector(G4VITProcess::Ge
94 fSelectedPostStepDoItVector(G4VITProcess::
95 {
96 fPhysicalStep = -1.;
97 fPreviousStepSize = -1.;
98
99 fSafety = -1.;
100 fProposedSafety = -1.;
101 fEndpointSafety = -1;
102
103 fStepStatus = fUndefined;
104
105 fTouchableHandle = nullptr;
106 }
107
108 //____________________________________________
109
110 //G4ITStepProcessor::
111 G4ITStepProcessorState::G4ITStepProcessorState
112 G4ITStepProcessorState_Lock(),
113 fSelectedAtRestDoItVector(right.fSelectedA
114 fSelectedPostStepDoItVector(right.fSelecte
115 {
116 fPhysicalStep = right.fPhysicalStep;
117 fPreviousStepSize = right.fPreviousStepSize;
118
119 fSafety = right.fSafety;
120 fProposedSafety = right.fProposedSafety;
121 fEndpointSafety = right.fEndpointSafety;
122
123 fStepStatus = right.fStepStatus;
124
125 fTouchableHandle = right.fTouchableHandle;
126 }
127
128 //____________________________________________
129
130 //G4ITStepProcessor::
131 G4ITStepProcessorState&
132 //G4ITStepProcessor::
133 G4ITStepProcessorState::operator=(const G4ITSt
134 {
135 if(this == &right) return *this;
136
137 fSelectedAtRestDoItVector.clear();
138 fSelectedAtRestDoItVector = right.fSelectedA
139 fSelectedPostStepDoItVector.clear();
140 fSelectedPostStepDoItVector = right.fSelecte
141
142 fPhysicalStep = right.fPhysicalStep;
143 fPreviousStepSize = right.fPreviousStepSize;
144
145 fSafety = right.fSafety;
146 fProposedSafety = right.fProposedSafety;
147 fEndpointSafety = right.fEndpointSafety;
148
149 fStepStatus = right.fStepStatus;
150
151 fTouchableHandle = right.fTouchableHandle;
152 return *this;
153 }
154
155 //____________________________________________
156
157 //G4ITStepProcessor::
158 G4ITStepProcessorState::~G4ITStepProcessorStat
159 {
160 ;
161 }
162
163 //____________________________________________
164
165 void G4ITStepProcessor::ClearProcessInfo()
166 {
167 std::map<const G4ParticleDefinition*, Proces
168
169 for(it = fProcessGeneralInfoMap.begin(); it
170 it++)
171 {
172 if(it->second != nullptr)
173 {
174 delete it->second;
175 it->second = 0;
176 }
177 }
178
179 fProcessGeneralInfoMap.clear();
180 }
181
182 //____________________________________________
183
184 void G4ITStepProcessor::ForceReInitialization(
185 {
186 fInitialized = false;
187 ClearProcessInfo();
188 Initialize();
189 }
190
191 //____________________________________________
192
193 void G4ITStepProcessor::Initialize()
194 {
195 CleanProcessor();
196 if(fInitialized) return;
197 // ActiveOnlyITProcess();
198
199 SetNavigator(G4ITTransportationManager::GetT
200 ->GetNavigatorForTracking());
201
202 fPhysIntLength = DBL_MAX;
203 kCarTolerance = 0.5
204 * G4GeometryTolerance::GetInstance()->Ge
205
206 if(fpVerbose == nullptr)
207 {
208 G4ITTrackingInteractivity* interactivity =
209
210 if(interactivity != nullptr)
211 {
212 fpVerbose = interactivity->GetSteppingVe
213 fpVerbose->SetStepProcessor(this);
214 }
215 }
216
217 fpTrackContainer = G4ITTrackHolder::Instance
218
219 fInitialized = true;
220 }
221 //____________________________________________
222
223 G4ITStepProcessor::~G4ITStepProcessor()
224 {
225 if(fpStep != nullptr)
226 {
227 fpStep->DeleteSecondaryVector();
228 delete fpStep;
229 }
230
231 delete fpSecondary;
232 ClearProcessInfo();
233 //G4ITTransportationManager::DeleteInstance(
234
235 // if(fpUserSteppingAction) d
236 }
237 //____________________________________________
238 // should not be used
239 G4ITStepProcessor::G4ITStepProcessor(const G4I
240 {
241 fpVerbose = rhs.fpVerbose;
242 fStoreTrajectory = rhs.fStoreTrajectory;
243
244 // fpUserSteppingAction = 0 ;
245 fpTrackingManager = nullptr;
246 fpNavigator = nullptr;
247 fInitialized = false;
248
249 kCarTolerance = rhs.kCarTolerance;
250 fInitialized = false;
251 fPreviousTimeStep = DBL_MAX;
252
253 CleanProcessor();
254 ResetSecondaries();
255 fpTrackContainer = nullptr;
256 fILTimeStep = DBL_MAX;
257 }
258
259 //____________________________________________
260
261 void G4ITStepProcessor::ResetLeadingTracks()
262 {
263 fLeadingTracks.Reset();
264 }
265
266 //____________________________________________
267
268 void G4ITStepProcessor::PrepareLeadingTracks()
269 {
270 fLeadingTracks.PrepareLeadingTracks();
271 }
272
273 //____________________________________________
274
275 G4ITStepProcessor& G4ITStepProcessor::operator
276 {
277 if(this == &rhs) return *this; // handle sel
278 //assignment operator
279 return *this;
280 }
281
282 //____________________________________________
283
284 void G4ITStepProcessor::ActiveOnlyITProcess()
285 {
286 // Method not used for the time being
287 #ifdef debug
288 G4cout<<"G4ITStepProcessor::CloneProcesses:
289 #endif
290
291 G4ParticleTable* theParticleTable = G4Partic
292 G4ParticleTable::G4PTblDicIterator* theParti
293 ->GetIterator();
294
295 theParticleIterator->reset();
296 // TODO : Ne faire la boucle que sur les IT
297 while((*theParticleIterator)())
298 {
299 G4ParticleDefinition* particle = thePartic
300 G4ProcessManager* pm = particle->GetProces
301
302 if(pm == nullptr)
303 {
304 G4cerr << "ERROR - G4ITStepProcessor::Ge
305 << particle->GetParticleName() << ", PDG
306 << particle->GetPDGEncoding() << G4endl;
307 G4Exception("G4ITStepProcessor::GetProce
308 FatalException, "Process Manager is
309 return;
310 }
311
312 ActiveOnlyITProcess(pm);
313 }
314 }
315
316 //____________________________________________
317
318 void G4ITStepProcessor::ActiveOnlyITProcess(G4
319 {
320 // Method not used for the time being
321 G4ProcessVector* processVector = processMana
322
323 G4VITProcess* itProcess = nullptr;
324 for(G4int i = 0; i < (G4int)processVector->s
325 {
326 G4VProcess* base_process = (*processVector
327 itProcess = dynamic_cast<G4VITProcess*>(ba
328
329 if(itProcess == nullptr)
330 {
331 processManager->SetProcessActivation(bas
332 }
333 }
334 }
335
336 //____________________________________________
337
338 void G4ITStepProcessor::SetupGeneralProcessInf
339
340 {
341
342 #ifdef debug
343 G4cout<<"G4ITStepProcessor::GetProcessNumber
344 #endif
345 if(pm == nullptr)
346 {
347 G4cerr << "ERROR - G4SteppingManager::GetP
348 << particle->GetParticleName() << ", PDG_c
349 << particle->GetPDGEncoding() << G4endl;
350 G4Exception("G4SteppingManager::GetProcess
351 FatalException, "Process Manager is no
352 return;
353 }
354
355 auto it =
356 fProcessGeneralInfoMap.find(particle);
357 if(it != fProcessGeneralInfoMap.end())
358 {
359 G4Exception("G4SteppingManager::SetupGener
360 "ITStepProcessor0003",
361 FatalException, "Process info already
362 return;
363 }
364
365 // here used as temporary
366 fpProcessInfo = new ProcessGeneralInfo();
367
368 // AtRestDoits
369 fpProcessInfo->MAXofAtRestLoops = pm->GetAtR
370 fpProcessInfo->fpAtRestDoItVector = pm->GetA
371 fpProcessInfo->fpAtRestGetPhysIntVector =
372 pm->GetAtRestProcessVector(typeGPIL);
373 #ifdef debug
374 G4cout << "G4ITStepProcessor::GetProcessNumb
375 << fpProcessInfo->MAXofAtRestLoops << G4endl
376 #endif
377
378 // AlongStepDoits
379 fpProcessInfo->MAXofAlongStepLoops =
380 pm->GetAlongStepProcessVector()->entries();
381 fpProcessInfo->fpAlongStepDoItVector =
382 pm->GetAlongStepProcessVector(typeDoIt);
383 fpProcessInfo->fpAlongStepGetPhysIntVector =
384 pm->GetAlongStepProcessVector(typeGPIL);
385 #ifdef debug
386 G4cout << "G4ITStepProcessor::GetProcessNumb
387 << fpProcessInfo->MAXofAlongStepLoops << G4e
388 #endif
389
390 // PostStepDoits
391 fpProcessInfo->MAXofPostStepLoops =
392 pm->GetPostStepProcessVector()->entries();
393 fpProcessInfo->fpPostStepDoItVector = pm->Ge
394 fpProcessInfo->fpPostStepGetPhysIntVector =
395 pm->GetPostStepProcessVector(typeGPIL);
396 #ifdef debug
397 G4cout << "G4ITStepProcessor::GetProcessNumb
398 << fpProcessInfo->MAXofPostStepLoops << G4en
399 #endif
400
401 if (SizeOfSelectedDoItVector<fpProcessInfo->
402 SizeOfSelectedDoItVector<fpProcessInfo->
403 SizeOfSelectedDoItVector<fpProcessInfo->
404 {
405 G4cerr << "ERROR - G4ITStepProcessor::GetP
406 << " SizeOfSelectedDoItVector= " <<
407 << " ; is smaller then one of MAXofAtRestL
408 << fpProcessInfo->MAXofAtRestLoops << G4en
409 << " or MAXofAlongStepLoops= " << f
410 << " or MAXofPostStepLoops= " << fpProcess
411 G4Exception("G4ITStepProcessor::GetProcess
412 "ITStepProcessor0004", FatalException,
413 "The array size is smaller than the ac
414 }
415
416 if((fpProcessInfo->fpAtRestDoItVector == nul
417 (fpProcessInfo->fpAlongStepDoItVector ==
418 (fpProcessInfo->fpPostStepDoItVector ==
419 {
420 G4ExceptionDescription exceptionDescriptio
421 exceptionDescription << "No DoIt process f
422 G4Exception("G4ITStepProcessor::DoStepping
423 FatalErrorInArgument,exceptionDescript
424 return;
425 }
426
427 if((fpProcessInfo->fpAlongStepGetPhysIntVect
428 && fpProcessInfo->MAXofAlongStepLoops>0)
429 {
430 fpProcessInfo->fpTransportation = dynamic_
431 ((*fpProcessInfo->fpAlongStepGetPhysIntVec
432 [G4int(fpProcessInfo->MAXofAlongStepLo
433
434 if(fpProcessInfo->fpTransportation == null
435 {
436 G4ExceptionDescription exceptionDescript
437 exceptionDescription << "No transportati
438 G4Exception("G4ITStepProcessor::SetupGen
439 "ITStepProcessor0006",
440 FatalErrorInArgument,exceptionDescri
441 }
442 }
443 fProcessGeneralInfoMap[particle] = fpProcess
444 // fpProcessInfo = 0;
445 }
446
447 //____________________________________________
448
449 void G4ITStepProcessor::SetTrack(G4Track* trac
450 {
451 fpTrack = track;
452 if(fpTrack != nullptr)
453 {
454 fpITrack = GetIT(fpTrack);
455 fpStep = const_cast<G4Step*>(fpTrack->GetS
456
457 if(fpITrack != nullptr)
458 {
459 fpTrackingInfo = fpITrack->GetTrackingIn
460 }
461 else
462 {
463 fpTrackingInfo = nullptr;
464 G4cerr << "Track ID : " << fpTrack->GetT
465
466 G4ExceptionDescription errMsg;
467 errMsg << "No IT pointer was attached to
468 G4Exception("G4ITStepProcessor::SetTrack
469 "ITStepProcessor0007",
470 FatalErrorInArgument,
471 errMsg);
472 }
473 }
474 else
475 {
476 fpITrack = nullptr;
477 fpStep = nullptr;
478 }
479 }
480 //____________________________________________
481
482 void G4ITStepProcessor::GetProcessInfo()
483 {
484 G4ParticleDefinition* particle = fpTrack->Ge
485 auto it =
486 fProcessGeneralInfoMap.find(particle);
487
488 if(it == fProcessGeneralInfoMap.end())
489 {
490 SetupGeneralProcessInfo(particle,
491 fpTrack->GetDefini
492 if(fpProcessInfo == nullptr)
493 {
494 G4ExceptionDescription exceptionDescript
495 G4Exception("G4ITStepProcessor::GetProce
496 "ITStepProcessor0008",
497 FatalErrorInArgument,
498 exceptionDescription);
499 return;
500 }
501 }
502 else
503 {
504 fpProcessInfo = it->second;
505 }
506 }
507
508 //____________________________________________
509
510 void G4ITStepProcessor::SetupMembers()
511 {
512 fpSecondary = fpStep->GetfSecondary();
513 fpPreStepPoint = fpStep->GetPreStepPoint();
514 fpPostStepPoint = fpStep->GetPostStepPoint()
515
516 fpState = (G4ITStepProcessorState*) fpITrack
517 ->GetStepProcessorState();
518
519 GetProcessInfo();
520 ResetSecondaries();
521 }
522
523 //____________________________________________
524
525 void G4ITStepProcessor::ResetSecondaries()
526 {
527 // Reset the secondary particles
528 fN2ndariesAtRestDoIt = 0;
529 fN2ndariesAlongStepDoIt = 0;
530 fN2ndariesPostStepDoIt = 0;
531 }
532
533 //____________________________________________
534
535 void G4ITStepProcessor::GetAtRestIL()
536 {
537 // Select the rest process which has the sho
538 // it is invoked. In rest processes, GPIL()
539 // returns the time before a process occurs.
540 G4double lifeTime(DBL_MAX), shortestLifeTime
541
542 fAtRestDoItProcTriggered = 0;
543 shortestLifeTime = DBL_MAX;
544
545 unsigned int NofInactiveProc=0;
546
547 for( G4int ri=0; ri < (G4int)fpProcessInfo->
548 {
549 fpCurrentProcess = dynamic_cast<G4VITProce
550 if (fpCurrentProcess== nullptr)
551 {
552 (fpState->fSelectedAtRestDoItVector)[ri]
553 NofInactiveProc++;
554 continue;
555 } // NULL means the process is inactivated
556
557 fCondition=NotForced;
558 fpCurrentProcess->SetProcessState(
559 fpTrackingInfo->GetProcessState(fpCurr
560
561 lifeTime = fpCurrentProcess->AtRestGPIL( *
562 fpCurrentProcess->ResetProcessState();
563
564 if(fCondition==Forced)
565 {
566 (fpState->fSelectedAtRestDoItVector)[ri]
567 }
568 else
569 {
570 (fpState->fSelectedAtRestDoItVector)[ri]
571 if(lifeTime < shortestLifeTime )
572 {
573 shortestLifeTime = lifeTime;
574 fAtRestDoItProcTriggered = G4int(ri);
575 }
576 }
577 }
578
579 (fpState->fSelectedAtRestDoItVector)[fAtRest
580
581 // G4cout << " --> Selected at rest process :
582 // << (*fpProcessInfo->fpAtRestGetPhys
583 // ->GetProcessName()
584 // << G4endl;
585
586 fTimeStep = shortestLifeTime;
587
588 // at least one process is necessary to dest
589 // exit with warning
590 if(NofInactiveProc==fpProcessInfo->MAXofAtRe
591 {
592 G4cerr << "ERROR - G4ITStepProcessor::Invo
593 << " No AtRestDoIt process is activ
594 }
595 }
596
597 //____________________________________________
598 G4double G4ITStepProcessor::ComputeInteraction
599 {
600 G4TrackManyList* mainList = fpTrackContainer
601 G4TrackManyList::iterator it = mainList ->be
602 G4TrackManyList::iterator end = mainList ->e
603
604 SetPreviousStepTime(previousTimeStep);
605
606 fILTimeStep = DBL_MAX;
607
608 for (; it != end; )
609 {
610 G4Track * track = *it;
611
612 #ifdef DEBUG
613 G4cout << "*CIL* " << GetIT(track)->GetNam
614 << " ID: " << track->GetTrackID()
615 << " at time : " << track->GetGlobalT
616 << G4endl;
617 #endif
618
619 ++it;
620 DefinePhysicalStepLength(track);
621
622 ExtractILData();
623 }
624
625 return fILTimeStep;
626 }
627
628 //____________________________________________
629
630 void G4ITStepProcessor::ExtractILData()
631 {
632 assert(fpTrack != 0);
633 if (fpTrack == nullptr)
634 {
635 CleanProcessor();
636 return;
637 }
638
639 // assert(fpTrack->GetTrackStatus() != fStop
640
641 if (fpTrack->GetTrackStatus() == fStopAndKil
642 {
643 // trackContainer->GetMainList()->pop(fpTra
644 fpTrackingManager->EndTracking(fpTrack);
645 CleanProcessor();
646 return;
647 }
648
649 if (IsInf(fTimeStep))
650 {
651 // G4cout << "!!!!!!!!!!!! IS INF " << tra
652 CleanProcessor();
653 return;
654 }
655 if (fTimeStep < fILTimeStep - DBL_EPSILON)
656 {
657 // G4cout << "!!!!!!!!!!!! TEMPS DIFFERENT
658 // << track->GetTrackID() << G4endl;
659
660 fLeadingTracks.Reset();
661
662 fILTimeStep = GetInteractionTime();
663
664 // G4cout << "Will set leading step to true
665 // << SP -> GetInteractionTime() <<
666 // << G4BestUnit(fILTimeStep, "Time"
667 // << G4endl;
668
669 // GetIT(fpTrack)->GetTrackingInfo()->SetLe
670 fLeadingTracks.Push(fpTrack);
671 }
672 else if(fabs(fILTimeStep - fTimeStep) < DBL_
673 {
674
675 // G4cout << "!!!!!!!!!!!! MEME TEMPS " <<
676 // G4cout << "Will set leading step to tru
677 // << SP -> GetInteractionTime() <<
678 // << fTimeStep << " the trackID is
679 // GetIT(fpTrack)->GetTrackingInfo()->SetLe
680 fLeadingTracks.Push(fpTrack);
681 }
682 // else
683 // {
684 // G4cout << "!!!! Bigger time : " << "c
685 // << " proposedTime : " << SP -> GetIntera
686 // }
687
688 CleanProcessor();
689 }
690
691 //____________________________________________
692
693 void G4ITStepProcessor::DefinePhysicalStepLeng
694 {
695 SetTrack(track);
696 DoDefinePhysicalStepLength();
697 }
698
699 //____________________________________________
700
701 void G4ITStepProcessor::SetInitialStep()
702 {
703 // DEBUG
704 // G4cout << "SetInitialStep for : " << f
705 //__________________________________________
706 // Initialize geometry
707
708 if(!fpTrack->GetTouchableHandle())
709 {
710 //========================================
711 // Create navigator state and Locate parti
712 //========================================
713 /*
714 fpNavigator->NewNavigatorStateAndLocate(f
715 fpTrack->GetMomentumDirection());
716
717 fpITrack->GetTrackingInfo()->
718 SetNavigatorState(fpNavigator->GetNavigat
719 */
720 fpNavigator->NewNavigatorState();
721 fpITrack->GetTrackingInfo()->SetNavigatorS
722 ->GetNavigatorState());
723
724 G4ThreeVector direction = fpTrack->GetMome
725 fpNavigator->LocateGlobalPointAndSetup(fpT
726 &di
727 fal
728 fal
729
730 fpState->fTouchableHandle = fpNavigator->C
731
732 fpTrack->SetTouchableHandle(fpState->fTouc
733 fpTrack->SetNextTouchableHandle(fpState->f
734 }
735 else
736 {
737 fpState->fTouchableHandle = fpTrack->GetTo
738 fpTrack->SetNextTouchableHandle(fpState->f
739
740 //========================================
741 // Create OR set navigator state
742 //========================================
743
744 if(fpITrack->GetTrackingInfo()->GetNavigat
745 {
746 fpNavigator->SetNavigatorState(fpITrack-
747 ->GetNavigatorState());
748 fpITrack->GetTrackingInfo()->SetNavigato
749 ->GetNavigatorState());
750 }
751 else
752 {
753 fpNavigator->NewNavigatorState(*((G4Touc
754 ->fTouchableHandle()));
755 fpITrack->GetTrackingInfo()->SetNavigato
756 ->GetNavigatorState());
757 }
758
759 G4VPhysicalVolume* oldTopVolume =
760 fpTrack->GetTouchableHandle()->GetVolu
761
762 //========================================
763 // Locate particle in geometry
764 //========================================
765
766 // G4VPhysicalVolume* newTopVolume =
767 // fpNavigator->LocateGlobalPointAndSet
768 // fpTrack->GetPosition(),
769 // &fpTrack->GetMomentumDirection()
770 // true, false);
771
772 G4VPhysicalVolume* newTopVolume =
773 fpNavigator->ResetHierarchyAndLocate(f
774 f
775 *
776
777
778 if(newTopVolume != oldTopVolume || oldTopV
779 == 1)
780 {
781 fpState->fTouchableHandle = fpNavigator-
782 fpTrack->SetTouchableHandle(fpState->fTo
783 fpTrack->SetNextTouchableHandle(fpState-
784 }
785 }
786
787 fpCurrentVolume = fpState->fTouchableHandle-
788
789 //__________________________________________
790 // If the primary track has 'Suspend' or 'Po
791 // set the track state to 'Alive'.
792 if((fpTrack->GetTrackStatus() == fSuspend) |
793 == fPostponeToNextEvent))
794 {
795 fpTrack->SetTrackStatus(fAlive);
796 }
797
798 //HoangTRAN: it's better to check the status
799 if(fpTrack->GetTrackStatus() == fStopAndKill
800
801 // If the primary track has 'zero' kinetic e
802 // state to 'StopButAlive'.
803 if(fpTrack->GetKineticEnergy() <= 0.0)
804 {
805 fpTrack->SetTrackStatus(fStopButAlive);
806 }
807 //__________________________________________
808 // Set vertex information of G4Track at here
809 if(fpTrack->GetCurrentStepNumber() == 0)
810 {
811 fpTrack->SetVertexPosition(fpTrack->GetPos
812 fpTrack->SetVertexMomentumDirection(fpTrac
813 fpTrack->SetVertexKineticEnergy(fpTrack->G
814 fpTrack->SetLogicalVolumeAtVertex(fpTrack-
815 }
816 //__________________________________________
817 // If track is already outside the world bou
818 if(fpCurrentVolume == nullptr)
819 {
820 // If the track is a primary, stop process
821 if(fpTrack->GetParentID() == 0)
822 {
823 G4cerr << "ERROR - G4ITStepProcessor::Se
824 << fpTrack->GetPosition()
825 << " - is outside of the world volume."
826 G4Exception("G4ITStepProcessor::SetIniti
827 FatalException, "Primary vertex outs
828 }
829
830 fpTrack->SetTrackStatus( fStopAndKill );
831 G4cout << "WARNING - G4ITStepProcessor::Se
832 << " Initial track position is ou
833 << fpTrack->GetPosition() << G4endl;
834 }
835 else
836 {
837 // Initial set up for attribues of 'Step'
838 fpStep->InitializeStep( fpTrack );
839 }
840
841 fpState->fStepStatus = fUndefined;
842 }
843 //____________________________________________
844
845 void G4ITStepProcessor::InitDefineStep()
846 {
847
848 if(fpStep == nullptr)
849 {
850 // Create new Step and give it to the trac
851 fpStep = new G4Step();
852 fpTrack->SetStep(fpStep);
853 fpSecondary = fpStep->NewSecondaryVector()
854
855 // Create new state and set it in the trac
856 fpState = new G4ITStepProcessorState();
857 fpITrack->GetTrackingInfo()->SetStepProces
858
859 SetupMembers();
860 SetInitialStep();
861
862 fpTrackingManager->StartTracking(fpTrack);
863 }
864 else
865 {
866 SetupMembers();
867
868 fpState->fPreviousStepSize = fpTrack->GetS
869 /***
870 // Send G4Step information to Hit/Dig if
871 fpCurrentVolume = fpStep->GetPreStepPoint
872 StepControlFlag = fpStep->GetControlFlag
873 if( fpCurrentVolume != 0 && StepControlFl
874 {
875 fpSensitive = fpStep->GetPreStepPoint()->
876 GetSensitiveDetector();
877
878 // if( fSensitive != 0 ) {
879 // fSensitive->Hit(fStep);
880 // }
881 }
882 ***/
883 // Store last PostStepPoint to PreStepPoin
884 // volume information of G4Track. Reset to
885 fpStep->CopyPostToPreStepPoint();
886 fpStep->ResetTotalEnergyDeposit();
887
888 //JA Set the volume before it is used (in
889 fpCurrentVolume = fpStep->GetPreStepPoint(
890 /*
891 G4cout << G4endl;
892 G4cout << "!!!!!!!!!!!!!!!!!!!!!!!!!!!" <
893 G4cout << "PreStepPoint Volume : "
894 << fpCurrentVolume->GetName() << G4endl;
895 G4cout << "Track Touchable : "
896 << fpTrack->GetTouchableHandle()->GetVolu
897 G4cout << "Track NextTouchable : "
898 << fpTrack->GetNextTouchableHandle()->Get
899 << G4endl;
900 */
901 // Reset the step's auxiliary points vecto
902 fpStep->SetPointerToVectorOfAuxiliaryPoint
903
904 // Switch next touchable in track to curre
905 fpTrack->SetTouchableHandle(fpTrack->GetNe
906 fpState->fTouchableHandle = fpTrack->GetTo
907 fpTrack->SetNextTouchableHandle(fpState->f
908
909 //! ADDED BACK
910 /*
911 G4VPhysicalVolume* oldTopVolume =
912 fpTrack->GetTouchableHandle()->GetVolume(
913 fpNavigator->SetNavigatorState(
914 fpITrack->GetTrackingInfo()->GetNavigator
915
916 G4VPhysicalVolume* newTopVolume = fpNavig
917 fpTrack->GetPosition(), fpTrack->GetMomen
918 *((G4TouchableHistory*) fpTrack->GetTouch
919
920 // G4VPhysicalVolume* newTopVolume=
921 // fpNavigator->LocateGlobalPointAndS
922 //
923 //
924
925 // G4cout << "New Top Volume : " <
926
927 if (newTopVolume != oldTopVolume || oldTo
928 == 1)
929 {
930 fpState->fTouchableHandle = fpNavigator->
931 fpTrack->SetTouchableHandle(fpState->fTou
932 fpTrack->SetNextTouchableHandle(fpState->
933 }
934
935 */
936 //! ADDED BACK
937 //========================================
938 // Only reset navigator state + reset volu
939 // No need to relocate
940 //========================================
941 fpNavigator->SetNavigatorState(fpITrack->G
942 ->GetNavigatorState());
943 }
944 }
945
946 //____________________________________________
947
948 // -------------------------------------------
949 // Compute Interaction Length
950 // -------------------------------------------
951 void G4ITStepProcessor::DoDefinePhysicalStepLe
952 {
953
954 InitDefineStep();
955
956 #ifdef G4VERBOSE
957 // !!!!! Verbose
958 if(fpVerbose != nullptr) fpVerbose->DPSLStar
959 #endif
960
961 G4TrackStatus trackStatus = fpTrack->GetTrac
962
963 if(trackStatus == fStopAndKill)
964 {
965 return;
966 }
967
968 if(trackStatus == fStopButAlive)
969 {
970 fpITrack->GetTrackingInfo()->SetNavigatorS
971 ->GetNavigatorState());
972 fpNavigator->ResetNavigatorState();
973 return GetAtRestIL();
974 }
975
976 // Find minimum Step length and correspondin
977 // demanded by active disc./cont. processes
978
979 // ReSet the counter etc.
980 fpState->fPhysicalStep = DBL_MAX; // Initial
981 fPhysIntLength = DBL_MAX; // Initialize by a
982
983 G4double proposedTimeStep = DBL_MAX;
984 G4VProcess* processWithPostStepGivenByTimeSt
985
986 // GPIL for PostStep
987 fPostStepDoItProcTriggered = fpProcessInfo->
988 fPostStepAtTimeDoItProcTriggered = fpProcess
989
990 // G4cout << "fpProcessInfo->MAXofPostStepL
991 // << fpProcessInfo->MAXofPostStepLo
992 // << " mol : " << fpITrack -> GetNa
993 // << " id : " << fpTrack->GetTrackI
994 // << G4endl;
995
996 for(std::size_t np = 0; np < fpProcessInfo->
997 {
998 fpCurrentProcess = dynamic_cast<G4VITProce
999 ->fpPostStepGetPhysIntVector)[(G4int)n
1000 if(fpCurrentProcess == nullptr)
1001 {
1002 (fpState->fSelectedPostStepDoItVector)[
1003 continue;
1004 } // NULL means the process is inactivate
1005
1006 fCondition = NotForced;
1007 fpCurrentProcess->SetProcessState(fpTrack
1008 ->GetProcessID()));
1009
1010 // G4cout << "Is going to call : "
1011 // << fpCurrentProcess -> GetProce
1012 // << G4endl;
1013 fPhysIntLength = fpCurrentProcess->PostSt
1014
1015
1016
1017 #ifdef G4VERBOSE
1018 // !!!!! Verbose
1019 if(fpVerbose != nullptr) fpVerbose->DPSLP
1020 #endif
1021
1022 fpCurrentProcess->ResetProcessState();
1023 //fpCurrentProcess->SetProcessState(0);
1024
1025 switch(fCondition)
1026 {
1027 case ExclusivelyForced: // Will need sp
1028 (fpState->fSelectedPostStepDoItVector
1029 fpState->fStepStatus = fExclusivelyFo
1030 fpStep->GetPostStepPoint()->SetProces
1031 break;
1032
1033 case Conditionally:
1034 // (fpState->fSelectedPostStepD
1035 G4Exception("G4ITStepProcessor::Defin
1036 "ITStepProcessor0008",
1037 FatalException,
1038 "This feature is no more
1039 break;
1040
1041 case Forced:
1042 (fpState->fSelectedPostStepDoItVector
1043 break;
1044
1045 case StronglyForced:
1046 (fpState->fSelectedPostStepDoItVector
1047 break;
1048
1049 default:
1050 (fpState->fSelectedPostStepDoItVector
1051 break;
1052 }
1053
1054 if(fCondition == ExclusivelyForced)
1055 {
1056 for(std::size_t nrest = np + 1; nrest <
1057 ++nrest)
1058 {
1059 (fpState->fSelectedPostStepDoItVector
1060 }
1061 return; // Please note the 'return' at
1062 }
1063
1064 if(fPhysIntLength < fpState->fPhysicalSte
1065 {
1066 // To avoid checking whether the proces
1067 // proposing a time step, the returned
1068 // negative (just for tagging)
1069 if(fpCurrentProcess->ProposesTimeStep()
1070 {
1071 fPhysIntLength *= -1;
1072 if(fPhysIntLength < proposedTimeStep)
1073 {
1074 proposedTimeStep = fPhysIntLength;
1075 fPostStepAtTimeDoItProcTriggered =
1076 processWithPostStepGivenByTimeStep
1077 }
1078 }
1079 else
1080 {
1081 fpState->fPhysicalStep = fPhysIntLeng
1082 fpState->fStepStatus = fPostStepDoItP
1083 fPostStepDoItProcTriggered = G4int(np
1084 fpStep->GetPostStepPoint()->SetProces
1085 }
1086 }
1087 }
1088
1089 // GPIL for AlongStep
1090 fpState->fProposedSafety = DBL_MAX;
1091 G4double safetyProposedToAndByProcess = fpS
1092
1093 for(std::size_t kp = 0; kp < fpProcessInfo-
1094 {
1095 fpCurrentProcess = dynamic_cast<G4VITProc
1096 ->fpAlongStepGetPhysIntVector)[(G4int
1097 if(fpCurrentProcess == nullptr) continue;
1098 // NULL means the process is inactivated
1099
1100 fpCurrentProcess->SetProcessState(
1101 fpTrackingInfo->GetProcessState(fpCur
1102 fPhysIntLength =
1103 fpCurrentProcess->AlongStepGPIL(*fpTr
1104 fpSta
1105 fpSta
1106 safet
1107 &fGPI
1108
1109 #ifdef G4VERBOSE
1110 // !!!!! Verbose
1111 if(fpVerbose != nullptr) fpVerbose->DPSLA
1112 #endif
1113
1114 if(fPhysIntLength < fpState->fPhysicalSte
1115 {
1116 fpState->fPhysicalStep = fPhysIntLength
1117 // Should save PS and TS in IT
1118
1119 // Check if the process wants to be the
1120 // multi-scattering proposes Step limit
1121 if(fGPILSelection == CandidateForSelect
1122 {
1123 fpState->fStepStatus = fAlongStepDoIt
1124 fpStep->GetPostStepPoint()->SetProces
1125 }
1126
1127 // Transportation is assumed to be the
1128 if(kp == fpProcessInfo->MAXofAlongStepL
1129 {
1130 fpTransportation = dynamic_cast<G4ITT
1131
1132 if(fpTransportation == nullptr)
1133 {
1134 G4ExceptionDescription exceptionDes
1135 exceptionDescription << "No transpo
1136 G4Exception("G4ITStepProcessor::DoD
1137 "ITStepProcessor0009",
1138 FatalErrorInArgument,
1139 exceptionDescription);
1140 }
1141
1142 fTimeStep = fpTransportation->GetInte
1143
1144 if(fpTrack->GetNextVolume() != nullpt
1145 else fpState->fStepStatus = fWorldBou
1146 }
1147 }
1148 else
1149 {
1150 if(kp == fpProcessInfo->MAXofAlongStepL
1151 {
1152 fpTransportation = dynamic_cast<G4ITT
1153
1154 if(fpTransportation == nullptr)
1155 {
1156 G4ExceptionDescription exceptionDes
1157 exceptionDescription << "No transpo
1158 G4Exception("G4ITStepProcessor::DoD
1159 "ITStepProcessor0010",
1160 FatalErrorInArgument,
1161 exceptionDescription);
1162 }
1163
1164 fTimeStep = fpTransportation->GetInte
1165 }
1166 }
1167
1168 // Handle PostStep processes sending back
1169 if(proposedTimeStep < fTimeStep)
1170 {
1171 if(fPostStepAtTimeDoItProcTriggered < f
1172 {
1173 if((fpState->fSelectedPostStepDoItVec
1174 {
1175 (fpState->fSelectedPostStepDoItVect
1176 NotForced;
1177 // (fpState->fSelectedPostStepDoItV
1178
1179 fpState->fStepStatus = fPostStepDoI
1180 fpStep->GetPostStepPoint()->SetProc
1181
1182 fTimeStep = proposedTimeStep;
1183
1184 fpTransportation->ComputeStep(*fpTr
1185 *fpSt
1186 fTime
1187 fpSta
1188 }
1189 }
1190 }
1191 else
1192 {
1193 if(fPostStepDoItProcTriggered < fpProce
1194 {
1195 if((fpState->fSelectedPostStepDoItVec
1196 {
1197 (fpState->fSelectedPostStepDoItVect
1198 NotForced;
1199 }
1200 }
1201 }
1202
1203 // fpCurrentProcess->SetProcessState(0);
1204 fpCurrentProcess->ResetProcessState();
1205
1206 // Make sure to check the safety, even if
1207 // by this process.
1208 //
1209 if(safetyProposedToAndByProcess < fpState
1210 // proposedSafety keeps the smallest valu
1211 fpState->fProposedSafety = safetyProposed
1212 else
1213 // safetyProposedToAndByProcess always pr
1214 safetyProposedToAndByProcess = fpState->f
1215
1216 }
1217
1218 fpITrack->GetTrackingInfo()->SetNavigatorSt
1219 fpNavigator->ResetNavigatorState();
1220 }
1221
1222 //___________________________________________
1223