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

#include <G4HadronElastic.hh>

Inheritance diagram for G4HadronElastic:
Collaboration diagram for G4HadronElastic:

Public Member Functions

 G4HadronElastic (const G4String &name="hElasticLHEP")
 
virtual ~G4HadronElastic ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
void SetLowestEnergyLimit (G4double value)
 
G4double LowestEnergyLimit () const
 
G4double ComputeMomentumCMS (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual void ModelDescription (std::ostream &) const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
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 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 50 of file G4HadronElastic.hh.

Constructor & Destructor Documentation

G4HadronElastic::G4HadronElastic ( const G4String name = "hElasticLHEP")

Definition at line 48 of file G4HadronElastic.cc.

49  : G4HadronicInteraction(name)
50 {
51  SetMinEnergy( 0.0*GeV );
52  SetMaxEnergy( 100.*TeV );
53  lowestEnergyLimit= 1.e-6*eV;
54 
55  theProton = G4Proton::Proton();
56  theNeutron = G4Neutron::Neutron();
57  theDeuteron = G4Deuteron::Deuteron();
58  theAlpha = G4Alpha::Alpha();
59  //Description();
60 }
static constexpr double TeV
Definition: G4SIunits.hh:218
void SetMinEnergy(G4double anEnergy)
static G4Proton * Proton()
Definition: G4Proton.cc:93
static constexpr double eV
Definition: G4SIunits.hh:215
G4HadronicInteraction(const G4String &modelName="HadronicModel")
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static constexpr double GeV
Definition: G4SIunits.hh:217
void SetMaxEnergy(const G4double anEnergy)
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89

Here is the call graph for this function:

G4HadronElastic::~G4HadronElastic ( )
virtual

Definition at line 63 of file G4HadronElastic.cc.

64 {}

Member Function Documentation

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

Implements G4HadronicInteraction.

Reimplemented in G4LEnp, G4LEHadronProtonElastic, and G4LEpp.

Definition at line 80 of file G4HadronElastic.cc.

82 {
84 
85  const G4HadProjectile* aParticle = &aTrack;
86  G4double ekin = aParticle->GetKineticEnergy();
87  if(ekin <= lowestEnergyLimit) {
90  return &theParticleChange;
91  }
92 
93  G4int A = targetNucleus.GetA_asInt();
94  G4int Z = targetNucleus.GetZ_asInt();
95 
96  G4double plab = aParticle->GetTotalMomentum();
97 
98  // Scattered particle referred to axis of incident particle
99  const G4ParticleDefinition* theParticle = aParticle->GetDefinition();
100  G4double m1 = theParticle->GetPDGMass();
101 
102  if (verboseLevel>1) {
103  G4cout << "G4HadronElastic: "
104  << aParticle->GetDefinition()->GetParticleName()
105  << " Plab(GeV/c)= " << plab/GeV
106  << " Ekin(MeV) = " << ekin/MeV
107  << " scattered off Z= " << Z
108  << " A= " << A
109  << G4endl;
110  }
111 
113  G4LorentzVector lv1 = aParticle->Get4Momentum();
114  G4LorentzVector lv(0.0,0.0,0.0,mass2);
115  lv += lv1;
116 
117  G4ThreeVector bst = lv.boostVector();
118  lv1.boost(-bst);
119 
120  G4ThreeVector p1 = lv1.vect();
121  G4double momentumCMS = p1.mag();
122  G4double tmax = 4.0*momentumCMS*momentumCMS;
123 
124  // Sampling in CM system
125  G4double t = SampleInvariantT(theParticle, plab, Z, A);
127  G4double cost = 1. - 2.0*t/tmax;
128  G4double sint;
129 
130  // problem in sampling
131  if(cost > 1.0 || cost < -1.0) {
132  //if(verboseLevel > 0) {
133  G4cout << "G4HadronElastic WARNING (1 - cost)= " << 1 - cost
134  << " after scattering of "
135  << aParticle->GetDefinition()->GetParticleName()
136  << " p(GeV/c)= " << plab/GeV
137  << " on an ion Z= " << Z << " A= " << A
138  << G4endl;
139  //}
140  cost = 1.0;
141  sint = 0.0;
142 
143  // normal situation
144  } else {
145  sint = std::sqrt((1.0-cost)*(1.0+cost));
146  }
147  if (verboseLevel>1) {
148  G4cout << " t= " << t << " tmax(GeV^2)= " << tmax/(GeV*GeV)
149  << " Pcms(GeV)= " << momentumCMS/GeV << " cos(t)=" << cost
150  << " sin(t)=" << sint << G4endl;
151  }
152  G4ThreeVector v1(sint*std::cos(phi),sint*std::sin(phi),cost);
153  v1 *= momentumCMS;
154  G4LorentzVector nlv1(v1.x(),v1.y(),v1.z(),
155  std::sqrt(momentumCMS*momentumCMS + m1*m1));
156 
157  nlv1.boost(bst);
158 
159  G4double eFinal = nlv1.e() - m1;
160  if (verboseLevel > 1) {
161  G4cout <<" m= " << m1 << " Efin(MeV)= " << eFinal
162  << " Proj: 4-mom " << lv1 << " Final: " << nlv1
163  << G4endl;
164  }
165  if(eFinal <= lowestEnergyLimit) {
166  if(eFinal < 0.0 && verboseLevel > 0) {
167  G4cout << "G4HadronElastic WARNING Efinal= " << eFinal
168  << " after scattering of "
169  << aParticle->GetDefinition()->GetParticleName()
170  << " p(GeV/c)= " << plab/GeV
171  << " on an ion Z= " << Z << " A= " << A
172  << G4endl;
173  }
175  nlv1 = G4LorentzVector(0.0,0.0,0.0,m1);
176 
177  } else {
178  theParticleChange.SetMomentumChange(nlv1.vect().unit());
180  }
181 
182  lv -= nlv1;
183  G4double erec = lv.e() - mass2;
184  if (verboseLevel > 1) {
185  G4cout << "Recoil: " <<" m= " << mass2 << " Erec(MeV)= " << erec
186  << " 4-mom: " << lv
187  << G4endl;
188  }
189 
190  if(erec > GetRecoilEnergyThreshold()) {
191  G4ParticleDefinition * theDef = 0;
192  if(Z == 1 && A == 1) { theDef = theProton; }
193  else if (Z == 1 && A == 2) { theDef = theDeuteron; }
194  else if (Z == 1 && A == 3) { theDef = G4Triton::Triton(); }
195  else if (Z == 2 && A == 3) { theDef = G4He3::He3(); }
196  else if (Z == 2 && A == 4) { theDef = theAlpha; }
197  else {
198  theDef =
200  }
201  G4DynamicParticle * aSec = new G4DynamicParticle(theDef, lv);
203  } else if(erec > 0.0) {
205  }
206 
207  return &theParticleChange;
208 }
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
static G4double GetNuclearMass(const G4double A, const G4double Z)
virtual G4double SampleInvariantT(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
int G4int
Definition: G4Types.hh:78
G4double GetRecoilEnergyThreshold() const
const G4String & GetParticleName() const
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
HepLorentzVector & boost(double, double, double)
G4double GetKineticEnergy() const
static G4Triton * Triton()
Definition: G4Triton.cc:95
const G4LorentzVector & Get4Momentum() const
void SetEnergyChange(G4double anEnergy)
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
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
static constexpr double MeV
Definition: G4SIunits.hh:214
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
static G4He3 * He3()
Definition: G4He3.cc:94
double mag() const
void SetMomentumChange(const G4ThreeVector &aV)
static constexpr double twopi
Definition: SystemOfUnits.h:55
G4double GetTotalMomentum() const
CLHEP::HepLorentzVector G4LorentzVector

Here is the call graph for this function:

G4double G4HadronElastic::ComputeMomentumCMS ( const G4ParticleDefinition p,
G4double  plab,
G4int  Z,
G4int  A 
)
inline

Definition at line 98 of file G4HadronElastic.hh.

100 {
101  G4double m1 = p->GetPDGMass();
102  G4double m12= m1*m1;
104  return plab*mass2/std::sqrt(m12 + mass2*mass2 + 2.*mass2*std::sqrt(m12 + plab*plab));
105 }
static G4double GetNuclearMass(const G4double A, const G4double Z)
double A(double temperature)
G4double GetPDGMass() const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4HadronElastic::LowestEnergyLimit ( ) const
inline

Definition at line 92 of file G4HadronElastic.hh.

93 {
94  return lowestEnergyLimit;
95 }
void G4HadronElastic::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronicInteraction.

Reimplemented in G4ElasticHadrNucleusHE, and G4ChipsElasticModel.

Definition at line 67 of file G4HadronElastic.cc.

68 {
69 
70  outFile << "G4HadronElastic is a hadron-nucleus elastic scattering\n"
71  << "model which uses the Gheisha two-exponential momentum\n"
72  << "transfer parameterization. The model is fully relativistic\n"
73  << "as opposed to the original Gheisha model which was not.\n"
74  << "This model may be used for all long-lived hadrons at all\n"
75  << "incident energies.\n";
76 
77 }
G4double G4HadronElastic::SampleInvariantT ( const G4ParticleDefinition p,
G4double  plab,
G4int  Z,
G4int  A 
)
virtual

Reimplemented from G4HadronicInteraction.

Reimplemented in G4ElasticHadrNucleusHE, G4hhElastic, G4DiffuseElastic, G4NuclNuclDiffuseElastic, G4LEnp, G4LEHadronProtonElastic, G4LEpp, G4ChipsElasticModel, and G4AntiNuclElastic.

Definition at line 212 of file G4HadronElastic.cc.

215 {
216  static const G4double GeV2 = GeV*GeV;
217  G4double momentumCMS = ComputeMomentumCMS(p,plab,Z,A);
218  G4double tmax = 4.0*momentumCMS*momentumCMS/GeV2;
219  G4double aa, bb, cc;
220  G4double dd = 10.;
221  G4Pow* g4pow = G4Pow::GetInstance();
222  if (A <= 62) {
223  bb = 14.5*g4pow->Z23(A);
224  aa = g4pow->powZ(A, 1.63)/bb;
225  cc = 1.4*g4pow->Z13(A)/dd;
226  } else {
227  bb = 60.*g4pow->Z13(A);
228  aa = g4pow->powZ(A, 1.33)/bb;
229  cc = 0.4*g4pow->powZ(A, 0.4)/dd;
230  }
231  G4double q1 = 1.0 - G4Exp(-bb*tmax);
232  G4double q2 = 1.0 - G4Exp(-dd*tmax);
233  G4double s1 = q1*aa;
234  G4double s2 = q2*cc;
235  if((s1 + s2)*G4UniformRand() < s2) {
236  q1 = q2;
237  bb = dd;
238  }
239  return -GeV2*G4Log(1.0 - G4UniformRand()*q1)/bb;
240 }
static G4Pow * GetInstance()
Definition: G4Pow.cc:55
Definition: G4Pow.hh:56
static ulg bb
Definition: csz_inflate.cc:344
#define G4UniformRand()
Definition: Randomize.hh:97
G4double Z13(G4int Z) const
Definition: G4Pow.hh:127
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
G4double ComputeMomentumCMS(const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double Z23(G4int Z) const
Definition: G4Pow.hh:154
G4double powZ(G4int Z, G4double y) const
Definition: G4Pow.hh:254
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4HadronElastic::SetLowestEnergyLimit ( G4double  value)
inline

Definition at line 87 of file G4HadronElastic.hh.

88 {
89  lowestEnergyLimit = value;
90 }
const XML_Char int const XML_Char * value
Definition: expat.h:331

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