Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4UnknownDecay Class Reference

#include <G4UnknownDecay.hh>

Inheritance diagram for G4UnknownDecay:
Collaboration diagram for G4UnknownDecay:

Public Member Functions

 G4UnknownDecay (const G4String &processName="UnknownDecay")
 
virtual ~G4UnknownDecay ()
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

virtual G4VParticleChangeDecayIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

Definition at line 44 of file G4UnknownDecay.hh.

Constructor & Destructor Documentation

G4UnknownDecay::G4UnknownDecay ( const G4String processName = "UnknownDecay")

Definition at line 47 of file G4UnknownDecay.cc.

48  :G4VDiscreteProcess(processName, fDecay),
49  verboseLevel(1),
50  HighestValue(20.0)
51 {
52  // set Process Sub Type
53  SetProcessSubType(static_cast<int>(DECAY_Unknown));
54 
55 #ifdef G4VERBOSE
56  if (GetVerboseLevel()>1) {
57  G4cout << "G4UnknownDecay constructor " << " Name:" << processName << G4endl;
58  }
59 #endif
60  pParticleChange = &fParticleChangeForDecay;
61 }
G4int GetVerboseLevel() const
G4GLOB_DLL std::ostream G4cout
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:432
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

G4UnknownDecay::~G4UnknownDecay ( )
virtual

Definition at line 63 of file G4UnknownDecay.cc.

64 {
65 }

Member Function Documentation

void G4UnknownDecay::BuildPhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 78 of file G4UnknownDecay.cc.

79 {
80  return;
81 }
G4VParticleChange * G4UnknownDecay::DecayIt ( const G4Track aTrack,
const G4Step aStep 
)
protectedvirtual

Definition at line 83 of file G4UnknownDecay.cc.

84 {
85  // The DecayIt() method returns by pointer a particle-change object.
86  // Units are expressed in GEANT4 internal units.
87 
88  // Initialize ParticleChange
89  // all members of G4VParticleChange are set to equal to
90  // corresponding member in G4Track
91  fParticleChangeForDecay.Initialize(aTrack);
92 
93  // get particle
94  const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
95 
96  //check if thePreAssignedDecayProducts exists
97  const G4DecayProducts* o_products = (aParticle->GetPreAssignedDecayProducts());
98  G4bool isPreAssigned = (o_products != 0);
99  G4DecayProducts* products = 0;
100 
101  if (!isPreAssigned ){
102  fParticleChangeForDecay.SetNumberOfSecondaries(0);
103  // Kill the parent particle
104  fParticleChangeForDecay.ProposeTrackStatus( fStopAndKill ) ;
105  fParticleChangeForDecay.ProposeLocalEnergyDeposit(0.0);
106 
108  return &fParticleChangeForDecay ;
109  }
110 
111  // copy decay products
112  products = new G4DecayProducts(*o_products);
113 
114  // get parent particle information ...................................
115  G4double ParentEnergy = aParticle->GetTotalEnergy();
116  G4double ParentMass = aParticle->GetMass();
117  if (ParentEnergy < ParentMass) {
118  ParentEnergy = ParentMass;
119 #ifdef G4VERBOSE
120  if (GetVerboseLevel()>1) {
121  G4cout << "G4UnknownDecay::DoIt : Total Energy is less than its mass" << G4endl;
122  G4cout << " Particle: " << aParticle->GetDefinition()->GetParticleName();
123  G4cout << " Energy:" << ParentEnergy/MeV << "[MeV]";
124  G4cout << " Mass:" << ParentMass/MeV << "[MeV]";
125  G4cout << G4endl;
126  }
127 #endif
128  }
129  G4ThreeVector ParentDirection(aParticle->GetMomentumDirection());
130 
131  G4double energyDeposit = 0.0;
132  G4double finalGlobalTime = aTrack.GetGlobalTime();
133  //boost all decay products to laboratory frame
134  //if the particle has traveled
135  if(aParticle->GetPreAssignedDecayProperTime()>0.) {
136  products->Boost( ParentEnergy, ParentDirection);
137  }
138 
139  //add products in fParticleChangeForDecay
140  G4int numberOfSecondaries = products->entries();
141  fParticleChangeForDecay.SetNumberOfSecondaries(numberOfSecondaries);
142 #ifdef G4VERBOSE
143  if (GetVerboseLevel()>1) {
144  G4cout << "G4UnknownDecay::DoIt : Decay vertex :";
145  G4cout << " Time: " << finalGlobalTime/ns << "[ns]";
146  G4cout << " X:" << (aTrack.GetPosition()).x() /cm << "[cm]";
147  G4cout << " Y:" << (aTrack.GetPosition()).y() /cm << "[cm]";
148  G4cout << " Z:" << (aTrack.GetPosition()).z() /cm << "[cm]";
149  G4cout << G4endl;
150  G4cout << "G4UnknownDecay::DoIt : decay products in Lab. Frame" << G4endl;
151  products->DumpInfo();
152  }
153 #endif
154  G4int index;
155  G4ThreeVector currentPosition;
156  const G4TouchableHandle thand = aTrack.GetTouchableHandle();
157  for (index=0; index < numberOfSecondaries; index++)
158  {
159  // get current position of the track
160  currentPosition = aTrack.GetPosition();
161  // create a new track object
162  G4Track* secondary = new G4Track( products->PopProducts(),
163  finalGlobalTime ,
164  currentPosition );
165  // switch on good for tracking flag
166  secondary->SetGoodForTrackingFlag();
167  secondary->SetTouchableHandle(thand);
168  // add the secondary track in the List
169  fParticleChangeForDecay.AddSecondary(secondary);
170  }
171  delete products;
172 
173  // Kill the parent particle
174  fParticleChangeForDecay.ProposeTrackStatus( fStopAndKill ) ;
175  fParticleChangeForDecay.ProposeLocalEnergyDeposit(energyDeposit);
176  fParticleChangeForDecay.ProposeGlobalTime( finalGlobalTime );
177  // reset NumberOfInteractionLengthLeft
179 
180  return &fParticleChangeForDecay ;
181 }
G4int GetVerboseLevel() const
G4double GetTotalEnergy() const
const G4DynamicParticle * GetDynamicParticle() const
const G4ThreeVector & GetPosition() const
const G4DecayProducts * GetPreAssignedDecayProducts() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4ParticleDefinition * GetDefinition() const
tuple x
Definition: test.py:50
void Boost(G4double totalEnergy, const G4ThreeVector &momentumDirection)
void ClearNumberOfInteractionLengthLeft()
Definition: G4VProcess.hh:447
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
virtual void Initialize(const G4Track &)
G4GLOB_DLL std::ostream G4cout
G4double GetMass() const
bool G4bool
Definition: G4Types.hh:79
void DumpInfo() const
const G4ThreeVector & GetMomentumDirection() const
static constexpr double cm
Definition: G4SIunits.hh:119
G4double GetGlobalTime() const
void AddSecondary(G4Track *aSecondary)
const G4TouchableHandle & GetTouchableHandle() const
void SetNumberOfSecondaries(G4int totSecondaries)
G4DynamicParticle * PopProducts()
tuple z
Definition: test.py:28
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
G4int entries() const
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
G4double GetPreAssignedDecayProperTime() const
#define ns
Definition: xmlparse.cc:614
void SetGoodForTrackingFlag(G4bool value=true)

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4UnknownDecay::GetMeanFreePath ( const G4Track aTrack,
G4double  previousStepSize,
G4ForceCondition condition 
)
protectedvirtual

Implements G4VDiscreteProcess.

Definition at line 73 of file G4UnknownDecay.cc.

74 {
75  return DBL_MIN;
76 }
#define DBL_MIN
Definition: templates.hh:75
G4int G4UnknownDecay::GetVerboseLevel ( ) const
inline

Definition at line 154 of file G4UnknownDecay.hh.

154 { return verboseLevel; }

Here is the caller graph for this function:

G4bool G4UnknownDecay::IsApplicable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4VProcess.

Definition at line 67 of file G4UnknownDecay.cc.

68 {
69  if(aParticleType.GetParticleName()=="unknown") return true;
70  return false;
71 }
const G4String & GetParticleName() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4VParticleChange * G4UnknownDecay::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
inlinevirtual

Reimplemented from G4VDiscreteProcess.

Definition at line 157 of file G4UnknownDecay.hh.

161 {
162  return DecayIt(aTrack, aStep);
163 }
virtual G4VParticleChange * DecayIt(const G4Track &aTrack, const G4Step &aStep)

Here is the call graph for this function:

G4double G4UnknownDecay::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
inlinevirtual

Reimplemented from G4VDiscreteProcess.

Definition at line 126 of file G4UnknownDecay.hh.

131 {
132  // pre-assigned UnknownDecay time
134 
135  if (pTime < 0.) pTime = DBL_MIN;
136 
137  // condition is set to "Not Forced"
138  *condition = NotForced;
139 
140  // reminder proper time
141  G4double remainder = pTime - track.GetProperTime();
142  if (remainder <= 0.0) remainder = DBL_MIN;
143 
144  // use pre-assigned Decay time to determine PIL
145  //return GetMeanFreePath(track, previousStepSize, condition);
146  return remainder*CLHEP::c_light;
147 
148 }
G4double condition(const G4ErrorSymMatrix &m)
G4double GetProperTime() const
const G4DynamicParticle * GetDynamicParticle() const
static constexpr double c_light
#define DBL_MIN
Definition: templates.hh:75
double G4double
Definition: G4Types.hh:76
G4double GetPreAssignedDecayProperTime() const

Here is the call graph for this function:

void G4UnknownDecay::SetVerboseLevel ( G4int  value)
inline

Definition at line 151 of file G4UnknownDecay.hh.

151 { verboseLevel = value; }
const XML_Char int const XML_Char * value
Definition: expat.h:331

The documentation for this class was generated from the following files: