Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/dna/management/include/G4ITStepProcessor.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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
 28 
 29 // The code is developed in the framework of the ESA AO7146
 30 //
 31 // We would be very happy hearing from you, send us your feedback! :)
 32 //
 33 // In order for Geant4-DNA to be maintained and still open-source,
 34 // article citations are crucial. 
 35 // If you use Geant4-DNA chemistry and you publish papers about your software, 
 36 // in addition to the general paper on Geant4-DNA:
 37 //
 38 // Int. J. Model. Simul. Sci. Comput. 1 (2010) 157–178
 39 //
 40 // we would be very happy if you could please also cite the following
 41 // reference papers on chemistry:
 42 //
 43 // J. Comput. Phys. 274 (2014) 841-882
 44 // Prog. Nucl. Sci. Tec. 2 (2011) 503-508 
 45 
 46 #ifndef G4ITSTEPPROCESSOR_H
 47 #define G4ITSTEPPROCESSOR_H
 48 
 49 #include "G4ios.hh"                   // Include from 'system'
 50 #include "globals.hh"                 // Include from 'global'
 51 #include "Randomize.hh"               // Include from 'global'
 52 
 53 #include "G4LogicalVolume.hh"         // Include from 'geometry'
 54 #include "G4VPhysicalVolume.hh"       // Include from 'geometry'
 55 #include "G4ProcessManager.hh"        // Include from 'piim'
 56 
 57 #include "G4Track.hh"                 // Include from 'track'
 58 #include "G4TrackVector.hh"           // Include from 'track'
 59 #include "G4TrackStatus.hh"           // Include from 'track'
 60 #include "G4StepStatus.hh"            // Include from 'track'
 61 #include "G4Step.hh"                  // Include from 'track'
 62 #include "G4StepPoint.hh"             // Include from 'track'
 63 #include "G4TouchableHandle.hh"       // Include from 'geometry'
 64 
 65 #include "G4ITStepProcessorState_Lock.hh"
 66 #include "G4ITLeadingTracks.hh"
 67 
 68 #include <vector>
 69 
 70 class G4ITNavigator;
 71 class G4ParticleDefinition;
 72 class G4ITTrackingManager;
 73 class G4IT;
 74 class G4TrackingInformation;
 75 class G4ITTransportation;
 76 class G4VITProcess;
 77 class G4VITSteppingVerbose;
 78 class G4ITTrackHolder;
 79 using G4SelectedAtRestDoItVector = std::vector<int, std::allocator<int>>;
 80 using G4SelectedAlongStepDoItVector = std::vector<int, std::allocator<int>>;
 81 using G4SelectedPostStepDoItVector = std::vector<int, std::allocator<int>>;
 82 
 83 //________________________________________________
 84 //
 85 // Members related to ParticleDefinition and not
 86 // proper to a track
 87 //________________________________________________
 88 struct ProcessGeneralInfo
 89 {
 90   G4ProcessVector* fpAtRestDoItVector;
 91   G4ProcessVector* fpAlongStepDoItVector;
 92   G4ProcessVector* fpPostStepDoItVector;
 93 
 94   G4ProcessVector* fpAtRestGetPhysIntVector;
 95   G4ProcessVector* fpAlongStepGetPhysIntVector;
 96   G4ProcessVector* fpPostStepGetPhysIntVector;
 97   //
 98   // Note: DoItVector has inverse order against GetPhysIntVector
 99   //       and SelectedPostStepDoItVector.
100   //
101   // * Max number of processes
102   std::size_t MAXofAtRestLoops;
103   std::size_t MAXofAlongStepLoops;
104   std::size_t MAXofPostStepLoops;
105   // Maximum number of processes for each type of process
106   // These depend on the G4ParticleDefinition, so on the track
107 
108   // * Transportation process
109   G4ITTransportation* fpTransportation;
110 };
111 
112 //________________________________________________
113 //
114 //          Members proper to a track
115 //________________________________________________
116 class G4ITStepProcessorState : public G4ITStepProcessorState_Lock
117 {
118 public:
119   G4ITStepProcessorState();
120   ~G4ITStepProcessorState() override;
121 
122   G4ITStepProcessorState(const G4ITStepProcessorState&);
123   G4ITStepProcessorState& operator=(const G4ITStepProcessorState&);
124 
125   // * Max Number of Process
126   G4SelectedAtRestDoItVector fSelectedAtRestDoItVector;
127   G4SelectedPostStepDoItVector fSelectedPostStepDoItVector;
128 
129   G4double fPhysicalStep;
130   G4double fPreviousStepSize;
131   G4double fSafety;
132 
133   G4StepStatus fStepStatus;
134 
135   // * Safety
136   G4double fProposedSafety;
137   // This keeps the minimum safety value proposed by AlongStepGPILs.
138   G4ThreeVector fEndpointSafOrigin;
139   G4double fEndpointSafety;
140   // To get the true safety value at the PostStepPoint, you have
141   // to subtract the distance to 'endpointSafOrigin' from this value.
142 
143   G4TouchableHandle fTouchableHandle;
144 };
145 
146 /**
147  * Its role is the same as G4StepManager :
148  * - Find the minimum physical length and corresponding time step
149  * - Step one track BUT on a given time step.
150  */
151 
152 class G4ITStepProcessor
153 {
154 friend class G4Scheduler;
155 public:
156   G4ITStepProcessor();
157   virtual ~G4ITStepProcessor();
158 
159   inline void SetPreviousStepTime(G4double);
160 
161   inline G4Track* GetTrack()
162   {
163     return fpTrack;
164   }
165   inline G4Step* GetStep()
166   {
167     return fpStep;
168   }
169   inline const G4Step* GetStep() const
170   {
171     return fpStep;
172   }
173   inline void SetStep(G4Step* val)
174   {
175     fpStep = val;
176   }
177 
178   inline G4TrackVector* GetSecondaries() const
179   {
180     return fpSecondary;
181   }
182   inline void SetTrackingManager(G4ITTrackingManager* trackMan)
183   {
184     fpTrackingManager = trackMan;
185   }
186   inline G4ITTrackingManager* GetTrackingManager()
187   {
188     return fpTrackingManager;
189   }
190 
191   //___________________________________
192 
193   virtual void Initialize();
194   void ForceReInitialization();
195 
196   void ResetLeadingTracks();
197   void PrepareLeadingTracks();
198 
199   //___________________________________
200   G4double ComputeInteractionLength(double previousTimeStep);
201   void DefinePhysicalStepLength(G4Track*);
202   G4double GetILTimeStep()
203   {
204     return fILTimeStep;
205   }
206 
207   //___________________________________
208   // DoIt
209   void DoIt(double timeStep);
210   void ExtractDoItData();
211   void Stepping(G4Track*, const double&);
212   void FindTransportationStep();
213   //___________________________________
214 
215   inline double GetInteractionTime();
216   inline const G4Track* GetTrack() const;
217   inline void CleanProcessor();
218 
219   std::size_t GetAtRestDoItProcTriggered() const
220   {
221     return fAtRestDoItProcTriggered;
222   }
223 
224   G4GPILSelection GetGPILSelection() const
225   {
226     return fGPILSelection;
227   }
228 
229   G4int GetN2ndariesAlongStepDoIt() const
230   {
231     return fN2ndariesAlongStepDoIt;
232   }
233 
234   G4int GetN2ndariesAtRestDoIt() const
235   {
236     return fN2ndariesAtRestDoIt;
237   }
238 
239   G4int GetN2ndariesPostStepDoIt() const
240   {
241     return fN2ndariesPostStepDoIt;
242   }
243 
244   const G4VITProcess* GetCurrentProcess() const
245   {
246     return fpCurrentProcess;
247   }
248 
249   G4double GetPhysIntLength() const
250   {
251     return fPhysIntLength;
252   }
253 
254   std::size_t GetPostStepAtTimeDoItProcTriggered() const
255   {
256     return fPostStepAtTimeDoItProcTriggered;
257   }
258 
259   std::size_t GetPostStepDoItProcTriggered() const
260   {
261     return fPostStepDoItProcTriggered;
262   }
263 
264   const ProcessGeneralInfo* GetCurrentProcessInfo() const
265   {
266     return fpProcessInfo;
267   }
268 
269   const G4ITStepProcessorState* GetProcessorState() const
270   {
271     return fpState;
272   }
273 
274   const G4VParticleChange* GetParticleChange() const
275   {
276     return fpParticleChange;
277   }
278 
279   const G4VPhysicalVolume* GetCurrentVolume() const
280   {
281     return fpCurrentVolume;
282   }
283 
284   G4ForceCondition GetCondition() const
285   {
286     return fCondition;
287   }
288 
289 protected:
290 
291   void ExtractILData();
292 
293   void SetupGeneralProcessInfo(G4ParticleDefinition*, G4ProcessManager*);
294   void ClearProcessInfo();
295   void SetTrack(G4Track*);
296 
297   void GetProcessInfo();
298 
299   void SetupMembers();
300   void ResetSecondaries();
301   void InitDefineStep();
302 
303   void SetInitialStep();
304 
305   void GetAtRestIL();
306   void DoDefinePhysicalStepLength();
307   void DoStepping();
308   void PushSecondaries();
309 
310   // void CalculateStep();
311   // void DoCalculateStep();
312 
313   // void CloneProcesses();
314   void ActiveOnlyITProcess();
315   void ActiveOnlyITProcess(G4ProcessManager*);
316 
317   void DealWithSecondaries(G4int&);
318   void InvokeAtRestDoItProcs();
319   void InvokeAlongStepDoItProcs();
320   void InvokePostStepDoItProcs();
321   void InvokePSDIP(std::size_t); //
322   void InvokeTransportationProc();
323   void SetNavigator(G4ITNavigator *value);
324   G4double CalculateSafety();
325 
326   // Return the estimated safety value at the PostStepPoint
327   void ApplyProductionCut(G4Track*);
328 
329   G4ITStepProcessor(const G4ITStepProcessor& other);
330   G4ITStepProcessor& operator=(const G4ITStepProcessor& other);
331 
332 private:
333   //________________________________________________
334   //
335   //              General members
336   //________________________________________________
337 
338   G4bool fInitialized;
339 
340   G4ITTrackingManager* fpTrackingManager;
341 
342   G4double kCarTolerance;
343   // Cached geometrical tolerance on surface
344 
345   G4ITNavigator* fpNavigator;
346   G4int fStoreTrajectory;
347   G4VITSteppingVerbose* fpVerbose;
348 
349   G4ITTrackHolder* fpTrackContainer;
350   G4ITLeadingTracks fLeadingTracks;
351 
352   //________________________________________________
353   //
354   // Members used as temporaries (= not proper to a track)
355   //________________________________________________
356 
357   G4double fTimeStep; // not proper to a track
358   G4double fILTimeStep; // proper to a track ensemble
359 
360   G4double fPreviousTimeStep;
361   G4TrackVector* fpSecondary; // get from fpStep at every configuration setup
362   G4VParticleChange* fpParticleChange;
363 
364   G4VITProcess* fpCurrentProcess;
365   // The pointer to the process of which DoIt or
366   // GetPhysicalInteractionLength has been just executed
367 
368   // * Secondaries
369   G4int fN2ndariesAtRestDoIt;
370   G4int fN2ndariesAlongStepDoIt;
371   G4int fN2ndariesPostStepDoIt;
372   // These are the numbers of secondaries generated by the process
373   // just executed.
374 
375   // * Process selection
376   std::size_t fAtRestDoItProcTriggered;
377   std::size_t fPostStepDoItProcTriggered;
378   std::size_t fPostStepAtTimeDoItProcTriggered;
379   // Record the selected process
380 
381   G4ForceCondition fCondition;
382   G4GPILSelection fGPILSelection;
383   // Above three variables are for the method
384   // DefinePhysicalStepLength(). To pass these information to
385   // the method Verbose, they are kept at here. Need a more
386   // elegant mechanism.
387 
388   G4double fPhysIntLength;
389   // The minimum physical interaction length over all possible processes
390 
391   // * Sensitive detector
392   //    G4SteppingControl StepControlFlag;
393   //    G4VSensitiveDetector*   fpSensitive;
394 
395   G4VPhysicalVolume* fpCurrentVolume;
396   // Get from fpStep or touchable, keep as member for user interface
397 
398   //________________________________________________
399   //
400   // Members related to ParticleDefinition and not
401   // proper to a track
402   //________________________________________________
403 
404   std::map<const G4ParticleDefinition*, ProcessGeneralInfo*> fProcessGeneralInfoMap;
405   ProcessGeneralInfo* fpProcessInfo;
406   G4ITTransportation* fpTransportation;
407 
408   //________________________________________________
409   //
410   // Members used for setting up the processor
411   //________________________________________________
412 
413   G4Track* fpTrack; // Set track
414   G4IT* fpITrack; // Set track
415   G4TrackingInformation* fpTrackingInfo; // Set track
416 
417   G4ITStepProcessorState* fpState; // SetupMembers or InitDefineStep
418   G4Step* fpStep; // Set track or InitDefineStep
419 
420   G4StepPoint* fpPreStepPoint; // SetupMembers
421   G4StepPoint* fpPostStepPoint; // SetupMembers
422 };
423 
424 //______________________________________________________________________________
425 
426 inline void G4ITStepProcessor::SetPreviousStepTime(G4double previousTimeStep)
427 {
428   fPreviousTimeStep = previousTimeStep;
429 }
430 
431 //______________________________________________________________________________
432 
433 inline const G4Track* G4ITStepProcessor::GetTrack() const
434 {
435   return fpTrack;
436 }
437 
438 //______________________________________________________________________________
439 
440 inline G4double G4ITStepProcessor::CalculateSafety()
441 {
442   return std::max(fpState->fEndpointSafety - (fpState->fEndpointSafOrigin
443                       - fpPostStepPoint->GetPosition()).mag(),
444                   kCarTolerance);
445 }
446 
447 //______________________________________________________________________________
448 
449 inline void G4ITStepProcessor::SetNavigator(G4ITNavigator *value)
450 {
451   fpNavigator = value;
452 }
453 
454 //______________________________________________________________________________
455 
456 inline void G4ITStepProcessor::CleanProcessor()
457 {
458   fTimeStep = DBL_MAX;
459   fPhysIntLength = DBL_MAX;
460 
461   fpState = nullptr;
462   fpTrack = nullptr;
463   fpTrackingInfo = nullptr;
464   fpITrack = nullptr;
465   fpStep = nullptr;
466   fpPreStepPoint = nullptr;
467   fpPostStepPoint = nullptr;
468 
469   fpParticleChange = nullptr;
470 
471   fpCurrentVolume = nullptr;
472   //    fpSensitive = 0;
473 
474   fpSecondary = nullptr;
475 
476   fpTransportation = nullptr;
477 
478   fpCurrentProcess= nullptr;
479   fpProcessInfo = nullptr;
480 
481   fAtRestDoItProcTriggered = INT_MAX;
482   fPostStepDoItProcTriggered = INT_MAX;
483   fPostStepAtTimeDoItProcTriggered = INT_MAX;
484   fGPILSelection = NotCandidateForSelection;
485   fCondition = NotForced;
486 }
487 
488 //______________________________________________________________________________
489 
490 inline double G4ITStepProcessor::GetInteractionTime()
491 {
492   return fTimeStep;
493 }
494 
495 #endif // G4ITSTEPPROCESSOR_H
496