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

#include <G4ChargeExchangeProcess.hh>

+ Inheritance diagram for G4ChargeExchangeProcess:

Public Member Functions

 G4ChargeExchangeProcess (const G4String &procName="chargeExchange")
 
virtual ~G4ChargeExchangeProcess ()
 
virtual G4bool IsApplicable (const G4ParticleDefinition &aParticleType)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &aParticleType)
 
virtual void DumpPhysicsTable (const G4ParticleDefinition &aParticleType)
 
virtual G4double GetElementCrossSection (const G4DynamicParticle *aParticle, const G4Element *anElement, const G4Material *mat=0)
 
- Public Member Functions inherited from G4HadronicProcess
 G4HadronicProcess (const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)
 
 G4HadronicProcess (const G4String &processName, G4HadronicProcessType subType)
 
virtual ~G4HadronicProcess ()
 
void RegisterMe (G4HadronicInteraction *a)
 
G4double GetElementCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=0)
 
G4double GetMicroscopicCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=0)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
void DumpPhysicsTable (const G4ParticleDefinition &p)
 
void AddDataSet (G4VCrossSectionDataSet *aDataSet)
 
G4EnergyRangeManagerGetManagerPointer ()
 
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *)
 
const G4NucleusGetTargetNucleus () const
 
const G4IsotopeGetTargetIsotope ()
 
virtual void ProcessDescription (std::ostream &outFile) const
 
void BiasCrossSectionByFactor (G4double aScale)
 
void SetEpReportLevel (G4int level)
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
std::pair< G4double, G4doubleGetEnergyMomentumCheckLevels () const
 
G4CrossSectionDataStoreGetCrossSectionDataStore ()
 
void MultiplyCrossSectionBy (G4double factor)
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
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
 
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)
 
- Protected Member Functions inherited from G4HadronicProcess
G4HadronicInteractionChooseHadronicInteraction (G4double kineticEnergy, G4Material *aMaterial, G4Element *anElement)
 
G4NucleusGetTargetNucleusPointer ()
 
void DumpState (const G4Track &, const G4String &, G4ExceptionDescription &)
 
const G4EnergyRangeManagerGetEnergyRangeManager () const
 
void SetEnergyRangeManager (const G4EnergyRangeManager &value)
 
G4HadronicInteractionGetHadronicInteraction () const
 
G4double GetLastCrossSection ()
 
void FillResult (G4HadFinalState *aR, const G4Track &aT)
 
G4HadFinalStateCheckResult (const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result) const
 
void CheckEnergyMomentumConservation (const G4Track &, const G4Nucleus &)
 
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4HadronicProcess
G4HadProjectile thePro
 
G4ParticleChangetheTotalResult
 
G4int epReportLevel
 
- 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 53 of file G4ChargeExchangeProcess.hh.

Constructor & Destructor Documentation

◆ G4ChargeExchangeProcess()

G4ChargeExchangeProcess::G4ChargeExchangeProcess ( const G4String procName = "chargeExchange")

Definition at line 54 of file G4ChargeExchangeProcess.cc.

55 : G4HadronicProcess(procName,fChargeExchange), first(true)
56{
57 thEnergy = 20.*MeV;
58 pPDG = 0;
59 verboseLevel= 1;
61 theProton = G4Proton::Proton();
62 theNeutron = G4Neutron::Neutron();
63 theAProton = G4AntiProton::AntiProton();
64 theANeutron = G4AntiNeutron::AntiNeutron();
65 thePiPlus = G4PionPlus::PionPlus();
66 thePiMinus = G4PionMinus::PionMinus();
67 thePiZero = G4PionZero::PionZero();
68 theKPlus = G4KaonPlus::KaonPlus();
69 theKMinus = G4KaonMinus::KaonMinus();
72 theL = G4Lambda::Lambda();
73 theAntiL = G4AntiLambda::AntiLambda();
74 theSPlus = G4SigmaPlus::SigmaPlus();
76 theSMinus = G4SigmaMinus::SigmaMinus();
78 theS0 = G4SigmaZero::SigmaZero();
80 theXiMinus = G4XiMinus::XiMinus();
81 theXi0 = G4XiZero::XiZero();
82 theAXiMinus = G4AntiXiMinus::AntiXiMinus();
83 theAXi0 = G4AntiXiZero::AntiXiZero();
84 theOmega = G4OmegaMinus::OmegaMinus();
86 theD = G4Deuteron::Deuteron();
87 theT = G4Triton::Triton();
88 theA = G4Alpha::Alpha();
89 theHe3 = G4He3::He3();
90}
@ fChargeExchange
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
static G4AntiLambda * AntiLambda()
static G4AntiNeutron * AntiNeutron()
static G4AntiOmegaMinus * AntiOmegaMinus()
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
static G4AntiSigmaMinus * AntiSigmaMinus()
static G4AntiSigmaPlus * AntiSigmaPlus()
static G4AntiSigmaZero * AntiSigmaZero()
static G4AntiXiMinus * AntiXiMinus()
static G4AntiXiZero * AntiXiZero()
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
static G4He3 * He3()
Definition: G4He3.cc:94
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:113
static G4KaonPlus * KaonPlus()
Definition: G4KaonPlus.cc:113
static G4KaonZeroLong * KaonZeroLong()
static G4KaonZeroShort * KaonZeroShort()
static G4Lambda * Lambda()
Definition: G4Lambda.cc:108
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4OmegaMinus * OmegaMinus()
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:98
static G4PionPlus * PionPlus()
Definition: G4PionPlus.cc:98
static G4PionZero * PionZero()
Definition: G4PionZero.cc:104
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4SigmaMinus * SigmaMinus()
static G4SigmaPlus * SigmaPlus()
Definition: G4SigmaPlus.cc:108
static G4SigmaZero * SigmaZero()
Definition: G4SigmaZero.cc:99
static G4Triton * Triton()
Definition: G4Triton.cc:95
G4int verboseLevel
Definition: G4VProcess.hh:368
static G4XiMinus * XiMinus()
Definition: G4XiMinus.cc:106
static G4XiZero * XiZero()
Definition: G4XiZero.cc:106

◆ ~G4ChargeExchangeProcess()

G4ChargeExchangeProcess::~G4ChargeExchangeProcess ( )
virtual

Definition at line 92 of file G4ChargeExchangeProcess.cc.

93{
94 delete factors;
95}

Member Function Documentation

◆ BuildPhysicsTable()

void G4ChargeExchangeProcess::BuildPhysicsTable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 97 of file G4ChargeExchangeProcess.cc.

99{
100 if(first) {
101 first = false;
102 theParticle = &aParticleType;
103 pPDG = theParticle->GetPDGEncoding();
104
106
107 const size_t n = 10;
108 if(theParticle == thePiPlus || theParticle == thePiMinus ||
109 theParticle == theKPlus || theParticle == theKMinus ||
110 theParticle == theK0S || theParticle == theK0L) {
111
112 G4double F[n] = {0.33,0.27,0.29,0.31,0.27,0.18,0.13,0.1,0.09,0.07};
113 factors = new G4PhysicsLinearVector(0.0,2.0*GeV,n);
114 for(size_t i=0; i<n; i++) {factors->PutValue(i,F[i]);}
115
116 } else {
117
118 G4double F[n] = {0.50,0.45,0.40,0.35,0.30,0.25,0.06,0.04,0.005,0.0};
119 factors = new G4PhysicsLinearVector(0.0,4.0*GeV,n);
120 for(size_t i=0; i<n; i++) {factors->PutValue(i,F[i]);}
121 }
122 //factors->SetSpline(true);
123
124 if(verboseLevel>1)
125 G4cout << "G4ChargeExchangeProcess for "
126 << theParticle->GetParticleName()
127 << G4endl;
128 }
130}
double G4double
Definition: G4Types.hh:64
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
G4CrossSectionDataStore * GetCrossSectionDataStore()
const G4String & GetParticleName() const
void PutValue(size_t index, G4double theValue)

◆ DumpPhysicsTable()

void G4ChargeExchangeProcess::DumpPhysicsTable ( const G4ParticleDefinition aParticleType)
virtual

Definition at line 191 of file G4ChargeExchangeProcess.cc.

193{
194 store->DumpPhysicsTable(aParticleType);
195}
void DumpPhysicsTable(const G4ParticleDefinition &)

◆ GetElementCrossSection()

G4double G4ChargeExchangeProcess::GetElementCrossSection ( const G4DynamicParticle aParticle,
const G4Element anElement,
const G4Material mat = 0 
)
virtual

Definition at line 132 of file G4ChargeExchangeProcess.cc.

136{
137 // gives the microscopic cross section in GEANT4 internal units
138 G4double Z = elm->GetZ();
139 G4int iz = G4int(Z);
140 G4double x = 0.0;
141
142 // The process is effective only above the threshold
143 if(iz == 1 || dp->GetKineticEnergy() < thEnergy) return x;
144
145 if(verboseLevel>1)
146 G4cout << "G4ChargeExchangeProcess compute GHAD CS for element "
147 << elm->GetName()
148 << G4endl;
149 x = store->GetCrossSection(dp, elm, mat);
150
151 if(verboseLevel>1)
152 G4cout << "G4ChargeExchangeProcess cross(mb)= " << x/millibarn
153 << " E(MeV)= " << dp->GetKineticEnergy()
154 << " " << theParticle->GetParticleName()
155 << " in Z= " << iz
156 << G4endl;
157 G4bool b;
158 G4double A = elm->GetN();
159 G4double ptot = dp->GetTotalMomentum();
160 x *= factors->GetValue(ptot, b)/std::pow(A, 0.42);
161 if(theParticle == thePiPlus || theParticle == theProton ||
162 theParticle == theKPlus || theParticle == theANeutron)
163 { x *= (1.0 - Z/A); }
164
165 else if(theParticle == thePiMinus || theParticle == theNeutron ||
166 theParticle == theKMinus || theParticle == theAProton)
167 { x *= Z/A; }
168
169 if(theParticle->GetPDGMass() < GeV) {
170 if(ptot > 2.*GeV) x *= 4.0*GeV*GeV/(ptot*ptot);
171 }
172
173 if(verboseLevel>1)
174 G4cout << "Corrected cross(mb)= " << x/millibarn << G4endl;
175
176 return x;
177}
int G4int
Definition: G4Types.hh:66
bool G4bool
Definition: G4Types.hh:67
G4double GetCrossSection(const G4DynamicParticle *, const G4Material *)
G4double GetValue(G4double theEnergy, G4bool &isOutRange)

◆ IsApplicable()

G4bool G4ChargeExchangeProcess::IsApplicable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4VProcess.

Definition at line 179 of file G4ChargeExchangeProcess.cc.

181{
182 const G4ParticleDefinition* p = &aParticleType;
183 return (p == thePiPlus || p == thePiMinus ||
184 p == theProton || p == theNeutron ||
185 p == theAProton|| p == theANeutron||
186 p == theKPlus || p == theKMinus ||
187 p == theK0S || p == theK0L ||
188 p == theL);
189}

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