Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4ChargeExchangeProcess Class Reference

#include <G4ChargeExchangeProcess.hh>

Inheritance diagram for G4ChargeExchangeProcess:
Collaboration 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=nullptr)
 
G4double GetMicroscopicCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
void DumpPhysicsTable (const G4ParticleDefinition &p)
 
void AddDataSet (G4VCrossSectionDataSet *aDataSet)
 
std::vector
< G4HadronicInteraction * > & 
GetHadronicInteractionList ()
 
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *)
 
const G4NucleusGetTargetNucleus () const
 
const G4IsotopeGetTargetIsotope ()
 
virtual void ProcessDescription (std::ostream &outFile) const
 
void BiasCrossSectionByFactor (G4double aScale)
 
void SetIntegral (G4bool val)
 
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 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 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 &)
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Member Functions inherited from G4HadronicProcess
G4HadronicInteractionChooseHadronicInteraction (const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, G4Material *aMaterial, G4Element *anElement)
 
G4NucleusGetTargetNucleusPointer ()
 
void DumpState (const G4Track &, const G4String &, G4ExceptionDescription &)
 
G4HadronicInteractionGetHadronicInteraction () const
 
G4double GetLastCrossSection ()
 
void FillResult (G4HadFinalState *aR, const G4Track &aT)
 
G4HadFinalStateCheckResult (const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
 
void CheckEnergyMomentumConservation (const G4Track &, const G4Nucleus &)
 
- 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 ( const G4String procName = "chargeExchange")

Definition at line 58 of file G4ChargeExchangeProcess.cc.

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

Here is the call graph for this function:

G4ChargeExchangeProcess::~G4ChargeExchangeProcess ( )
virtual

Definition at line 96 of file G4ChargeExchangeProcess.cc.

97 {
98  if (factors) delete factors;
99 }

Member Function Documentation

void G4ChargeExchangeProcess::BuildPhysicsTable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 102 of file G4ChargeExchangeProcess.cc.

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

Here is the call graph for this function:

void G4ChargeExchangeProcess::DumpPhysicsTable ( const G4ParticleDefinition aParticleType)
virtual

Definition at line 196 of file G4ChargeExchangeProcess.cc.

197 {
198  store->DumpPhysicsTable(aParticleType);
199 }
void DumpPhysicsTable(const G4ParticleDefinition &)

Here is the call graph for this function:

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

Definition at line 136 of file G4ChargeExchangeProcess.cc.

140 {
141  // gives the microscopic cross section in GEANT4 internal units
142  G4double Z = elm->GetZ();
143  G4int iz = G4int(Z);
144  G4double x = 0.0;
145 
146  // The process is effective only above the threshold
147  if(iz == 1 || dp->GetKineticEnergy() < thEnergy) return x;
148 
149  if(verboseLevel>1)
150  G4cout << "G4ChargeExchangeProcess compute GHAD CS for element "
151  << elm->GetName()
152  << G4endl;
153  x = store->GetCrossSection(dp, elm, mat);
154 
155  if(verboseLevel>1)
156  G4cout << "G4ChargeExchangeProcess cross(mb)= " << x/millibarn
157  << " E(MeV)= " << dp->GetKineticEnergy()
158  << " " << theParticle->GetParticleName()
159  << " in Z= " << iz
160  << G4endl;
161  G4bool b;
162  G4double A = elm->GetN();
163  G4double ptot = dp->GetTotalMomentum();
164  x *= factors->GetValue(ptot, b)/G4Pow::GetInstance()->powA(A, 0.42);
165  if(theParticle == thePiPlus || theParticle == theProton ||
166  theParticle == theKPlus || theParticle == theANeutron)
167  { x *= (1.0 - Z/A); }
168 
169  else if(theParticle == thePiMinus || theParticle == theNeutron ||
170  theParticle == theKMinus || theParticle == theAProton)
171  { x *= Z/A; }
172 
173  if(theParticle->GetPDGMass() < GeV) {
174  if(ptot > 2.*GeV) x *= 4.0*GeV*GeV/(ptot*ptot);
175  }
176 
177  if(verboseLevel>1)
178  G4cout << "Corrected cross(mb)= " << x/millibarn << G4endl;
179 
180  return x;
181 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
G4double GetValue(G4double theEnergy, G4bool &isOutRange) const
G4int verboseLevel
Definition: G4VProcess.hh:368
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
bool G4bool
Definition: G4Types.hh:79
G4double GetCrossSection(const G4DynamicParticle *, const G4Material *)
G4double GetPDGMass() const
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static constexpr double millibarn
Definition: G4SIunits.hh:106

Here is the call graph for this function:

G4bool G4ChargeExchangeProcess::IsApplicable ( const G4ParticleDefinition aParticleType)
virtual

Reimplemented from G4VProcess.

Definition at line 184 of file G4ChargeExchangeProcess.cc.

185 {
186  const G4ParticleDefinition* p = &aParticleType;
187  return (p == thePiPlus || p == thePiMinus ||
188  p == theProton || p == theNeutron ||
189  p == theAProton|| p == theANeutron||
190  p == theKPlus || p == theKMinus ||
191  p == theK0S || p == theK0L ||
192  p == theL);
193 }
const char * p
Definition: xmltok.h:285

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