Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VITDiscreteProcess.cc
Go to the documentation of this file.
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 // $Id: G4VITDiscreteProcess.cc 71231 2013-06-12 13:06:28Z gcosmo $
28 //
29 //
30 // --------------------------------------------------------------
31 // GEANT 4 class implementation file
32 //
33 // History: first implementation, based on object model of
34 // 2nd December 1995, G.Cosmo
35 // --------------------------------------------------------------
36 // New Physics scheme 8 Jan. 1997 H.Kurahige
37 // ------------------------------------------------------------
38 
39 #include "G4VITDiscreteProcess.hh"
40 #include "G4SystemOfUnits.hh"
41 
42 #include "G4Step.hh"
43 #include "G4Track.hh"
44 #include "G4MaterialTable.hh"
45 #include "G4VParticleChange.hh"
46 
48  G4VITProcess("No Name Discrete Process")
49 {
50  G4Exception("G4VDiscreteProcess::G4VDiscreteProcess()",
51  "ProcMan102",
53  "Default constructor is called");
54 }
55 
56 //------------------------------------------------------------------------------
57 
59  G4ProcessType aType) :
60  G4VITProcess(aName, aType)
61 {
62  enableAtRestDoIt = false;
63  enableAlongStepDoIt = false;
64 }
65 
66 //------------------------------------------------------------------------------
67 
69 {
70 }
71 
72 //------------------------------------------------------------------------------
73 
75  G4VITProcess(right)
76 {
77 }
78 
79 //------------------------------------------------------------------------------
80 
84  G4double previousStepSize,
86 {
87  if((previousStepSize < 0.0)
88  || (fpState->theNumberOfInteractionLengthLeft <= 0.0))
89  {
90  // beginning of tracking (or just after DoIt of this process)
92  } else if(previousStepSize > 0.0)
93  {
94  // subtract NumberOfInteractionLengthLeft
95  SubtractNumberOfInteractionLengthLeft(previousStepSize);
96  } else
97  {
98  // zero step
99  // DO NOTHING
100  }
101 
102  // condition is set to "Not Forced"
103  *condition = NotForced;
104 
105  // get mean free path
106  fpState->currentInteractionLength = GetMeanFreePath(track,
107  previousStepSize,
108  condition);
109 
110  G4double value;
111  if(fpState->currentInteractionLength < DBL_MAX)
112  {
113  value = fpState->theNumberOfInteractionLengthLeft
114  * fpState->currentInteractionLength;
115  } else
116  {
117  value = DBL_MAX;
118  }
119 #ifdef G4VERBOSE
120  if(verboseLevel > 1)
121  {
122  G4cout << "G4VDiscreteProcess::PostStepGetPhysicalInteractionLength ";
123  G4cout << "[ " << GetProcessName() << "]" << G4endl;
124  track.GetDynamicParticle()->DumpInfo();
125  G4cout << " in Material " << track.GetMaterial()->GetName() << G4endl;
126  G4cout << "InteractionLength= " << value / cm << "[cm] " << G4endl;
127  }
128 #endif
129  return value;
130 }
131 
132 //------------------------------------------------------------------------------
133 
135  const G4Step&)
136 {
137 // clear NumberOfInteractionLengthLeft
139 
140  return pParticleChange;
141 }
G4double condition(const G4ErrorSymMatrix &m)
virtual void ClearNumberOfInteractionLengthLeft()
G4int verboseLevel
Definition: G4VProcess.hh:368
const G4DynamicParticle * GetDynamicParticle() const
virtual void ResetNumberOfInteractionLengthLeft()
const G4String & GetName() const
Definition: G4Material.hh:178
void DumpInfo(G4int mode=0) const
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
virtual void SubtractNumberOfInteractionLengthLeft(G4double previousStepSize)
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)
G4GLOB_DLL std::ostream G4cout
const XML_Char int const XML_Char * value
Definition: expat.h:331
static constexpr double cm
Definition: G4SIunits.hh:119
G4shared_ptr< G4ProcessState > fpState
Definition: G4Step.hh:76
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
G4Material * GetMaterial() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
virtual G4double GetMeanFreePath(const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
#define G4endl
Definition: G4ios.hh:61
G4bool enableAlongStepDoIt
Definition: G4VProcess.hh:351
double G4double
Definition: G4Types.hh:76
G4VITDiscreteProcess(const G4String &, G4ProcessType aType=fNotDefined)
G4ForceCondition
#define DBL_MAX
Definition: templates.hh:83
G4ProcessType