Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4KaonMinusAbsorption Class Reference

#include <G4KaonMinusAbsorption.hh>

+ Inheritance diagram for G4KaonMinusAbsorption:

Public Member Functions

 G4KaonMinusAbsorption (const G4String &processName="KaonMinusAbsorption", G4ProcessType aType=fHadronic)
 
 ~G4KaonMinusAbsorption ()
 
G4bool IsApplicable (const G4ParticleDefinition &)
 
void PreparePhysicsTable (const G4ParticleDefinition &)
 
void BuildPhysicsTable (const G4ParticleDefinition &)
 
G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
G4double GetMeanLifeTime (const G4Track &, G4ForceCondition *)
 
G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
G4int GetNumberOfSecondaries ()
 
G4GHEKinematicsVectorGetSecondaryKinematics ()
 
- Public Member Functions inherited from G4VRestProcess
 G4VRestProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VRestProcess (G4VRestProcess &)
 
virtual ~G4VRestProcess ()
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &, G4double, G4ForceCondition *)
 
virtual G4VParticleChangePostStepDoIt (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
 
virtual G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)=0
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)=0
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)=0
 
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 G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
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
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
virtual G4double GetMeanLifeTime (const G4Track &aTrack, G4ForceCondition *condition)=0
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- 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 47 of file G4KaonMinusAbsorption.hh.

Constructor & Destructor Documentation

◆ G4KaonMinusAbsorption()

G4KaonMinusAbsorption::G4KaonMinusAbsorption ( const G4String processName = "KaonMinusAbsorption",
G4ProcessType  aType = fHadronic 
)

Definition at line 46 of file G4KaonMinusAbsorption.cc.

47 :
48 G4VRestProcess (processName, aType), // initialization
49 massKaonMinus(G4KaonMinus::KaonMinus()->GetPDGMass()/GeV),
50 massGamma(G4Gamma::Gamma()->GetPDGMass()/GeV),
51 massPionZero(G4PionZero::PionZero()->GetPDGMass()/GeV),
52 massProton(G4Proton::Proton()->GetPDGMass()/GeV),
53 massLambda(G4Lambda::Lambda()->GetPDGMass()/GeV),
54 pdefKaonMinus(G4KaonMinus::KaonMinus()),
55 pdefGamma(G4Gamma::Gamma()),
56 pdefPionZero(G4PionZero::PionZero()),
57 pdefProton(G4Proton::Proton()),
58 pdefNeutron(G4Neutron::Neutron()),
59 pdefLambda(G4Lambda::Lambda()),
60 pdefDeuteron(G4Deuteron::Deuteron()),
61 pdefTriton(G4Triton::Triton()),
62 pdefAlpha(G4Alpha::Alpha())
63{
64 G4HadronicDeprecate("G4KaonMinusAbsorption");
65 if (verboseLevel>0) {
66 G4cout << GetProcessName() << " is created "<< G4endl;
67 }
72
74}
#define MAX_SECONDARIES
#define G4HadronicDeprecate(name)
@ fHadronAtRest
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
void RegisterExtraProcess(G4VProcess *)
static G4HadronicProcessStore * Instance()
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:113
static G4Lambda * Lambda()
Definition: G4Lambda.cc:108
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4PionZero * PionZero()
Definition: G4PionZero.cc:104
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Triton * Triton()
Definition: G4Triton.cc:95
G4int verboseLevel
Definition: G4VProcess.hh:368
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:403
const G4String & GetProcessName() const
Definition: G4VProcess.hh:379

◆ ~G4KaonMinusAbsorption()

G4KaonMinusAbsorption::~G4KaonMinusAbsorption ( )

Definition at line 78 of file G4KaonMinusAbsorption.cc.

79{
81 delete [] pv;
82 delete [] eve;
83 delete [] gkin;
84}
void DeRegisterExtraProcess(G4VProcess *)

Member Function Documentation

◆ AtRestDoIt()

G4VParticleChange * G4KaonMinusAbsorption::AtRestDoIt ( const G4Track track,
const G4Step  
)
virtual

Reimplemented from G4VRestProcess.

Definition at line 146 of file G4KaonMinusAbsorption.cc.

155{
156
157// Initialize ParticleChange
158// all members of G4VParticleChange are set to equal to
159// corresponding member in G4Track
160
162
163// Store some global quantities that depend on current material and particle
164
165 globalTime = track.GetGlobalTime()/s;
166 G4Material * aMaterial = track.GetMaterial();
167 const G4int numberOfElements = aMaterial->GetNumberOfElements();
168 const G4ElementVector* theElementVector = aMaterial->GetElementVector();
169
170 const G4double* theAtomicNumberDensity = aMaterial->GetAtomicNumDensityVector();
171 G4double normalization = 0;
172 for ( G4int i1=0; i1 < numberOfElements; i1++ )
173 {
174 normalization += theAtomicNumberDensity[i1] ; // change when nucleon specific
175 // probabilities are included.
176 }
177 G4double runningSum= 0.;
178 G4double random = G4UniformRand()*normalization;
179 for ( G4int i2=0; i2 < numberOfElements; i2++ )
180 {
181 runningSum += theAtomicNumberDensity[i2]; // change when nucleon specific
182 // probabilities are included.
183 if (random<=runningSum)
184 {
185 targetCharge = G4double( ((*theElementVector)[i2])->GetZ());
186 targetAtomicMass = (*theElementVector)[i2]->GetN();
187 }
188 }
189 if (random>runningSum)
190 {
191 targetCharge = G4double((*theElementVector)[numberOfElements-1]->GetZ());
192 targetAtomicMass = (*theElementVector)[numberOfElements-1]->GetN();
193
194 }
195
196 if (verboseLevel>1) {
197 G4cout << "G4KaonMinusAbsorption::AtRestDoIt is invoked " <<G4endl;
198 }
199
200 G4ParticleMomentum momentum;
201 G4float localtime;
202
204
205 GenerateSecondaries(); // Generate secondaries
206
208
209 for ( G4int isec = 0; isec < ngkine; isec++ ) {
210 G4DynamicParticle* aNewParticle = new G4DynamicParticle;
211 aNewParticle->SetDefinition( gkin[isec].GetParticleDef() );
212 aNewParticle->SetMomentum( gkin[isec].GetMomentum() * GeV );
213
214 localtime = globalTime + gkin[isec].GetTOF();
215
216 G4Track* aNewTrack = new G4Track( aNewParticle, localtime*s, position );
217 aNewTrack->SetTouchableHandle(track.GetTouchableHandle());
218 aParticleChange.AddSecondary( aNewTrack );
219
220 }
221
223
224 aParticleChange.ProposeTrackStatus(fStopAndKill); // Kill the incident KaonMinus
225
226// clear InteractionLengthLeft
227
229
230 return &aParticleChange;
231
232}
std::vector< G4Element * > G4ElementVector
@ fStopAndKill
double G4double
Definition: G4Types.hh:64
float G4float
Definition: G4Types.hh:65
int G4int
Definition: G4Types.hh:66
#define G4UniformRand()
Definition: Randomize.hh:53
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
void SetMomentum(const G4ThreeVector &momentum)
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:189
size_t GetNumberOfElements() const
Definition: G4Material.hh:185
const G4double * GetAtomicNumDensityVector() const
Definition: G4Material.hh:215
void AddSecondary(G4Track *aSecondary)
virtual void Initialize(const G4Track &)
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
G4double GetGlobalTime() const
G4Material * GetMaterial() const
const G4TouchableHandle & GetTouchableHandle() const
void ProposeTrackStatus(G4TrackStatus status)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetNumberOfSecondaries(G4int totSecondaries)
virtual void ResetNumberOfInteractionLengthLeft()
Definition: G4VProcess.cc:92
G4ParticleChange aParticleChange
Definition: G4VProcess.hh:289

◆ AtRestGetPhysicalInteractionLength()

G4double G4KaonMinusAbsorption::AtRestGetPhysicalInteractionLength ( const G4Track track,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VRestProcess.

Definition at line 120 of file G4KaonMinusAbsorption.cc.

124{
125 // beggining of tracking
127
128 // condition is set to "Not Forced"
130
131 // get mean life time
133
134 if ((currentInteractionLength <0.0) || (verboseLevel>2)){
135 G4cout << "G4KaonMinusAbsorptionProcess::AtRestGetPhysicalInteractionLength ";
136 G4cout << "[ " << GetProcessName() << "]" <<G4endl;
137 track.GetDynamicParticle()->DumpInfo();
138 G4cout << " in Material " << track.GetMaterial()->GetName() <<G4endl;
139 G4cout << "MeanLifeTime = " << currentInteractionLength/ns << "[ns]" <<G4endl;
140 }
141
143
144}
G4double condition(const G4ErrorSymMatrix &m)
@ NotForced
void DumpInfo(G4int mode=0) const
G4double GetMeanLifeTime(const G4Track &, G4ForceCondition *)
const G4String & GetName() const
Definition: G4Material.hh:177
const G4DynamicParticle * GetDynamicParticle() const
G4double currentInteractionLength
Definition: G4VProcess.hh:297
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
#define ns
Definition: xmlparse.cc:597

◆ BuildPhysicsTable()

void G4KaonMinusAbsorption::BuildPhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4VProcess.

Definition at line 91 of file G4KaonMinusAbsorption.cc.

92{
94}
void PrintInfo(const G4ParticleDefinition *)

◆ GetMeanLifeTime()

G4double G4KaonMinusAbsorption::GetMeanLifeTime ( const G4Track ,
G4ForceCondition  
)
inlinevirtual

Implements G4VRestProcess.

Definition at line 72 of file G4KaonMinusAbsorption.hh.

73 {return 0.0;}

Referenced by AtRestGetPhysicalInteractionLength().

◆ GetNumberOfSecondaries()

G4int G4KaonMinusAbsorption::GetNumberOfSecondaries ( )

Definition at line 107 of file G4KaonMinusAbsorption.cc.

108{
109 return ( ngkine );
110
111}

◆ GetSecondaryKinematics()

G4GHEKinematicsVector * G4KaonMinusAbsorption::GetSecondaryKinematics ( )

Definition at line 114 of file G4KaonMinusAbsorption.cc.

115{
116 return ( &gkin[0] );
117
118}

◆ IsApplicable()

G4bool G4KaonMinusAbsorption::IsApplicable ( const G4ParticleDefinition particle)
virtual

Reimplemented from G4VProcess.

Definition at line 98 of file G4KaonMinusAbsorption.cc.

101{
102 return ( &particle == pdefKaonMinus );
103
104}

◆ PreparePhysicsTable()

void G4KaonMinusAbsorption::PreparePhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4VProcess.

Definition at line 86 of file G4KaonMinusAbsorption.cc.

87{
89}
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)

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