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

#include <G4ChargeExchange.hh>

Inheritance diagram for G4ChargeExchange:
Collaboration diagram for G4ChargeExchange:

Public Member Functions

 G4ChargeExchange ()
 
virtual ~G4ChargeExchange ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
void SetLowestEnergyLimit (G4double value)
 
void SetRecoilKinEnergyLimit (G4double value)
 
G4double SampleT (G4double p, G4double A)
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair
< G4double, G4double
GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double,
G4double
GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 52 of file G4ChargeExchange.hh.

Constructor & Destructor Documentation

G4ChargeExchange::G4ChargeExchange ( )

Definition at line 54 of file G4ChargeExchange.cc.

54  : G4HadronicInteraction("Charge Exchange")
55 {
56  SetMinEnergy( 0.0*GeV );
57  SetMaxEnergy( 100.*TeV );
58 
59  lowEnergyRecoilLimit = 100.*keV;
60  lowestEnergyLimit = 1.*MeV;
61 
62  theProton = G4Proton::Proton();
63  theNeutron = G4Neutron::Neutron();
64  theAProton = G4AntiProton::AntiProton();
65  theANeutron = G4AntiNeutron::AntiNeutron();
66  thePiPlus = G4PionPlus::PionPlus();
67  thePiMinus = G4PionMinus::PionMinus();
68  thePiZero = G4PionZero::PionZero();
69  theKPlus = G4KaonPlus::KaonPlus();
70  theKMinus = G4KaonMinus::KaonMinus();
73  theL = G4Lambda::Lambda();
74  theAntiL = G4AntiLambda::AntiLambda();
75  theSPlus = G4SigmaPlus::SigmaPlus();
76  theASPlus = G4AntiSigmaPlus::AntiSigmaPlus();
77  theSMinus = G4SigmaMinus::SigmaMinus();
78  theASMinus = G4AntiSigmaMinus::AntiSigmaMinus();
79  theS0 = G4SigmaZero::SigmaZero();
81  theXiMinus = G4XiMinus::XiMinus();
82  theXi0 = G4XiZero::XiZero();
83  theAXiMinus = G4AntiXiMinus::AntiXiMinus();
84  theAXi0 = G4AntiXiZero::AntiXiZero();
85  theOmega = G4OmegaMinus::OmegaMinus();
87  theD = G4Deuteron::Deuteron();
88  theT = G4Triton::Triton();
89  theA = G4Alpha::Alpha();
90  theHe3 = G4He3::He3();
91 }
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 constexpr double TeV
Definition: G4SIunits.hh:218
static G4AntiSigmaMinus * AntiSigmaMinus()
static G4XiZero * XiZero()
Definition: G4XiZero.cc:106
void SetMinEnergy(G4double anEnergy)
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
G4HadronicInteraction(const G4String &modelName="HadronicModel")
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 constexpr double GeV
Definition: G4SIunits.hh:217
void SetMaxEnergy(const G4double anEnergy)
static G4AntiXiZero * AntiXiZero()
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 constexpr double keV
Definition: G4SIunits.hh:216
static G4AntiNeutron * AntiNeutron()

Here is the call graph for this function:

G4ChargeExchange::~G4ChargeExchange ( )
virtual

Definition at line 93 of file G4ChargeExchange.cc.

94 {}

Member Function Documentation

G4HadFinalState * G4ChargeExchange::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtual

Implements G4HadronicInteraction.

Definition at line 96 of file G4ChargeExchange.cc.

98 {
100  const G4HadProjectile* aParticle = &aTrack;
101  G4double ekin = aParticle->GetKineticEnergy();
102 
103  G4int A = targetNucleus.GetA_asInt();
104  G4int Z = targetNucleus.GetZ_asInt();
105 
106  if(ekin <= lowestEnergyLimit || A < 3) {
109  return &theParticleChange;
110  }
111 
112  G4double plab = aParticle->GetTotalMomentum();
113 
114  if (verboseLevel > 1)
115  G4cout << "G4ChargeExchange::DoIt: Incident particle plab="
116  << plab/GeV << " GeV/c "
117  << " ekin(MeV) = " << ekin/MeV << " "
118  << aParticle->GetDefinition()->GetParticleName() << G4endl;
119 
120  // Scattered particle referred to axis of incident particle
121  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
122 
123  G4int N = A - Z;
124  G4int projPDG = theParticle->GetPDGEncoding();
125  if (verboseLevel > 1)
126  G4cout << "G4ChargeExchange for " << theParticle->GetParticleName()
127  << " PDGcode= " << projPDG << " on nucleus Z= " << Z
128  << " A= " << A << " N= " << N
129  << G4endl;
130 
131  G4ParticleDefinition * theDef = 0;
132 
134  G4LorentzVector lv1 = aParticle->Get4Momentum();
135  G4LorentzVector lv0(0.0,0.0,0.0,mass2);
136 
137  G4LorentzVector lv = lv0 + lv1;
138  G4ThreeVector bst = lv.boostVector();
139  lv1.boost(-bst);
140  lv0.boost(-bst);
141 
142  // Sample final particles
143  G4bool theHyperon = false;
144  G4ParticleDefinition* theRecoil = 0;
145  G4ParticleDefinition* theSecondary = 0;
146 
147  if(theParticle == theProton) {
148  theSecondary = theNeutron;
149  Z++;
150  } else if(theParticle == theNeutron) {
151  theSecondary = theProton;
152  Z--;
153  } else if(theParticle == thePiPlus) {
154  theSecondary = thePiZero;
155  Z++;
156  } else if(theParticle == thePiMinus) {
157  theSecondary = thePiZero;
158  Z--;
159  } else if(theParticle == theKPlus) {
160  if(G4UniformRand()<0.5) theSecondary = theK0S;
161  else theSecondary = theK0L;
162  Z++;
163  } else if(theParticle == theKMinus) {
164  if(G4UniformRand()<0.5) theSecondary = theK0S;
165  else theSecondary = theK0L;
166  Z--;
167  } else if(theParticle == theK0S || theParticle == theK0L) {
168  if(G4UniformRand()*A < G4double(Z)) {
169  theSecondary = theKPlus;
170  Z--;
171  } else {
172  theSecondary = theKMinus;
173  Z++;
174  }
175  } else if(theParticle == theANeutron) {
176  theSecondary = theAProton;
177  Z++;
178  } else if(theParticle == theAProton) {
179  theSecondary = theANeutron;
180  Z--;
181  } else if(theParticle == theL) {
183  if(G4UniformRand()*A < G4double(Z)) {
184  if(x < 0.2) {
185  theSecondary = theS0;
186  } else if (x < 0.4) {
187  theSecondary = theSPlus;
188  Z--;
189  } else if (x < 0.6) {
190  theSecondary = theProton;
191  theRecoil = theL;
192  theHyperon = true;
193  A--;
194  } else if (x < 0.8) {
195  theSecondary = theProton;
196  theRecoil = theS0;
197  theHyperon = true;
198  A--;
199  } else {
200  theSecondary = theNeutron;
201  theRecoil = theSPlus;
202  theHyperon = true;
203  A--;
204  }
205  } else {
206  if(x < 0.2) {
207  theSecondary = theS0;
208  } else if (x < 0.4) {
209  theSecondary = theSMinus;
210  Z++;
211  } else if (x < 0.6) {
212  theSecondary = theNeutron;
213  theRecoil = theL;
214  A--;
215  theHyperon = true;
216  } else if (x < 0.8) {
217  theSecondary = theNeutron;
218  theRecoil = theS0;
219  theHyperon = true;
220  A--;
221  } else {
222  theSecondary = theProton;
223  theRecoil = theSMinus;
224  theHyperon = true;
225  A--;
226  }
227  }
228  }
229 
230  if (Z == 1 && A == 2) theDef = theD;
231  else if (Z == 1 && A == 3) theDef = theT;
232  else if (Z == 2 && A == 3) theDef = theHe3;
233  else if (Z == 2 && A == 4) theDef = theA;
234  else {
235  theDef =
237  }
238  if(!theSecondary) { return &theParticleChange; }
239 
240  G4double m11 = theSecondary->GetPDGMass();
241  G4double m21 = theDef->GetPDGMass();
242  if(theRecoil) { m21 += theRecoil->GetPDGMass(); }
243  else { theRecoil = theDef; }
244 
245  G4double etot = lv0.e() + lv1.e();
246 
247  // kinematiacally impossible
248  if(etot < m11 + m21) {
251  return &theParticleChange;
252  }
253 
254  G4ThreeVector p1 = lv1.vect();
255  G4double e1 = 0.5*etot*(1.0 - (m21*m21 - m11*m11)/(etot*etot));
256  // G4double e2 = etot - e1;
257  G4double ptot = std::sqrt(e1*e1 - m11*m11);
258 
259  G4double tmax = 4.0*ptot*ptot;
260  G4double g2 = GeV*GeV;
261 
262  G4double t = g2*SampleT(tmax/g2, A);
263 
264  if(verboseLevel>1)
265  G4cout <<"## G4ChargeExchange t= " << t << " tmax= " << tmax
266  << " ptot= " << ptot << G4endl;
267 
268  // Sampling in CM system
269  G4double phi = G4UniformRand()*twopi;
270  G4double cost = 1. - 2.0*t/tmax;
271  if(std::abs(cost) > 1.0) cost = 1.0;
272  G4double sint = std::sqrt((1.0-cost)*(1.0+cost));
273 
274  //if (verboseLevel > 1)
275  // G4cout << "cos(t)=" << cost << " std::sin(t)=" << sint << G4endl;
276 
277  G4ThreeVector v1(sint*std::cos(phi),sint*std::sin(phi),cost);
278  v1 *= ptot;
279  G4LorentzVector nlv1(v1.x(),v1.y(),v1.z(),e1);
280  G4LorentzVector nlv0 = lv0 + lv1 - nlv1;
281 
282  nlv0.boost(bst);
283  nlv1.boost(bst);
284 
287  G4DynamicParticle * aSec = new G4DynamicParticle(theSecondary, nlv1);
289 
290  G4double erec = nlv0.e() - m21;
291 
292  //G4cout << "erec= " <<erec << " Esec= " << aSec->GetKineticEnergy() << G4endl;
293 
294  if(theHyperon) {
296  aSec = new G4DynamicParticle();
297  aSec->SetDefinition(theRecoil);
298  aSec->SetKineticEnergy(0.0);
299  } else if(erec > lowEnergyRecoilLimit) {
300  aSec = new G4DynamicParticle(theRecoil, nlv0);
302  } else {
303  if(erec < 0.0) erec = 0.0;
305  }
306  return &theParticleChange;
307 }
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
Hep3Vector boostVector() const
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
tuple x
Definition: test.py:50
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
static constexpr double twopi
Definition: G4SIunits.hh:76
void SetStatusChange(G4HadFinalStateStatus aS)
G4IonTable * GetIonTable() const
Hep3Vector vect() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
const G4ParticleDefinition * GetDefinition() const
bool G4bool
Definition: G4Types.hh:79
HepLorentzVector & boost(double, double, double)
G4double GetKineticEnergy() const
const G4LorentzVector & Get4Momentum() const
void SetKineticEnergy(G4double aEnergy)
void SetEnergyChange(G4double anEnergy)
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
G4double SampleT(G4double p, G4double A)
Hep3Vector unit() const
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
static constexpr double GeV
Definition: G4SIunits.hh:217
void SetLocalEnergyDeposit(G4double aE)
#define G4endl
Definition: G4ios.hh:61
**D E S C R I P T I O N
Definition: HEPEvtcom.cc:77
static constexpr double MeV
Definition: G4SIunits.hh:214
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
void SetDefinition(const G4ParticleDefinition *aParticleDefinition)
void SetMomentumChange(const G4ThreeVector &aV)
G4double GetTotalMomentum() const

Here is the call graph for this function:

G4double G4ChargeExchange::SampleT ( G4double  p,
G4double  A 
)

Definition at line 309 of file G4ChargeExchange.cc.

310 {
311  G4double aa, bb, cc, dd;
312  if (A <= 62.) {
313  aa = G4Pow::GetInstance()->powA(A, 1.63);
314  bb = 14.5*G4Pow::GetInstance()->powA(A, 0.66);
315  cc = 1.4*G4Pow::GetInstance()->powA(A, 0.33);
316  dd = 10.;
317  } else {
318  aa = G4Pow::GetInstance()->powA(A, 1.33);
319  bb = 60.*G4Pow::GetInstance()->powA(A, 0.33);
320  cc = 0.4*G4Pow::GetInstance()->powA(A, 0.40);
321  dd = 10.;
322  }
323  G4double x1 = (1.0 - G4Exp(-tmax*bb))*aa/bb;
324  G4double x2 = (1.0 - G4Exp(-tmax*dd))*cc/dd;
325 
326  G4double t;
327  G4double y = bb;
328  if(G4UniformRand()*(x1 + x2) < x2) y = dd;
329 
330  const G4int maxNumberOfLoops = 10000;
331  G4int loopCounter = 0;
332  do {
333  t = -G4Log(G4UniformRand())/y;
334  } while ( (t > tmax) &&
335  ++loopCounter < maxNumberOfLoops ); /* Loop checking, 10.08.2015, A.Ribon */
336  if ( loopCounter >= maxNumberOfLoops ) {
337  t = 0.0;
338  }
339 
340  return t;
341 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
G4double powA(G4double A, G4double y) const
Definition: G4Pow.hh:259
int G4int
Definition: G4Types.hh:78
static ulg bb
Definition: csz_inflate.cc:344
#define G4UniformRand()
Definition: Randomize.hh:97
double A(double temperature)
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4ChargeExchange::SetLowestEnergyLimit ( G4double  value)
inline

Definition at line 112 of file G4ChargeExchange.hh.

113 {
114  lowestEnergyLimit = value;
115 }
const XML_Char int const XML_Char * value
Definition: expat.h:331
void G4ChargeExchange::SetRecoilKinEnergyLimit ( G4double  value)
inline

Definition at line 107 of file G4ChargeExchange.hh.

108 {
109  lowEnergyRecoilLimit = value;
110 }
const XML_Char int const XML_Char * value
Definition: expat.h:331

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