Geant4_10
G4ChargeExchangeProcess.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4ChargeExchangeProcess.cc 66892 2013-01-17 10:57:59Z gunter $
28 //
29 //
30 // Geant4 Hadron Charge Exchange Process -- source file
31 //
32 // Created 21 April 2006 V.Ivanchenko
33 //
34 // Modified:
35 // 24-Apr-06 V.Ivanchenko add neutron scattering on hydrogen from CHIPS
36 // 07-Jun-06 V.Ivanchenko fix problem of rotation of final state
37 // 25-Jul-06 V.Ivanchenko add 19 MeV low energy for CHIPS
38 // 23-Jan-07 V.Ivanchenko add cross section interfaces with Z and A
39 // and do not use CHIPS for cross sections
40 // 14-Sep-12 M.Kelsey -- Pass subType code to base ctor
41 
43 #include "globals.hh"
44 #include "G4SystemOfUnits.hh"
47 #include "G4Element.hh"
48 #include "G4ElementVector.hh"
49 #include "G4IsotopeVector.hh"
50 #include "G4Neutron.hh"
51 #include "G4Proton.hh"
52 #include "G4PhysicsLinearVector.hh"
53 
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();
75  theASPlus = G4AntiSigmaPlus::AntiSigmaPlus();
76  theSMinus = G4SigmaMinus::SigmaMinus();
77  theASMinus = G4AntiSigmaMinus::AntiSigmaMinus();
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 }
91 
93 {
94  delete factors;
95 }
96 
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 }
131 
133  const G4DynamicParticle* dp,
134  const G4Element* elm,
135  const G4Material* mat)
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 }
178 
180 IsApplicable(const G4ParticleDefinition& aParticleType)
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 }
190 
193 {
194  store->DumpPhysicsTable(aParticleType);
195 }
virtual G4bool IsApplicable(const G4ParticleDefinition &aParticleType)
G4double GetValue(G4double theEnergy, G4bool &isOutRange) const
G4int verboseLevel
Definition: G4VProcess.hh:368
static G4AntiOmegaMinus * AntiOmegaMinus()
G4double GetKineticEnergy() const
G4double GetN() const
Definition: G4Element.hh:134
const char * p
Definition: xmltok.h:285
static G4OmegaMinus * OmegaMinus()
G4double GetZ() const
Definition: G4Element.hh:131
static G4KaonZeroLong * KaonZeroLong()
tuple x
Definition: test.py:50
static G4AntiSigmaPlus * AntiSigmaPlus()
int G4int
Definition: G4Types.hh:78
int millibarn
Definition: hepunit.py:40
const G4String & GetParticleName() const
static G4SigmaZero * SigmaZero()
Definition: G4SigmaZero.cc:99
static G4KaonMinus * KaonMinus()
Definition: G4KaonMinus.cc:113
Float_t mat
Definition: plot.C:40
G4double GetTotalMomentum() const
static G4AntiSigmaMinus * AntiSigmaMinus()
G4ChargeExchangeProcess(const G4String &procName="chargeExchange")
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
static G4XiZero * XiZero()
Definition: G4XiZero.cc:106
tuple b
Definition: test.py:12
Char_t n[5]
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
G4GLOB_DLL std::ostream G4cout
virtual G4double GetElementCrossSection(const G4DynamicParticle *aParticle, const G4Element *anElement, const G4Material *mat=0)
Float_t Z
Definition: plot.C:39
static G4KaonZeroShort * KaonZeroShort()
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
bool G4bool
Definition: G4Types.hh:79
G4double iz
Definition: TRTMaterials.hh:39
static G4XiMinus * XiMinus()
Definition: G4XiMinus.cc:106
G4CrossSectionDataStore * GetCrossSectionDataStore()
void PutValue(size_t index, G4double theValue)
static G4AntiXiMinus * AntiXiMinus()
static G4Triton * Triton()
Definition: G4Triton.cc:95
G4double GetCrossSection(const G4DynamicParticle *, const G4Material *)
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:104
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
static G4SigmaMinus * SigmaMinus()
void DumpPhysicsTable(const G4ParticleDefinition &)
G4double GetPDGMass() const
static G4AntiLambda * AntiLambda()
static G4PionMinus * PionMinus()
Definition: G4PionMinus.cc:98
static G4AntiSigmaZero * AntiSigmaZero()
G4int first
static G4AntiXiZero * AntiXiZero()
#define G4endl
Definition: G4ios.hh:61
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
static G4SigmaPlus * SigmaPlus()
Definition: G4SigmaPlus.cc:108
static G4Lambda * Lambda()
Definition: G4Lambda.cc:108
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Element.hh:127
static G4KaonPlus * KaonPlus()
Definition: G4KaonPlus.cc:113
static G4He3 * He3()
Definition: G4He3.cc:94
virtual void BuildPhysicsTable(const G4ParticleDefinition &aParticleType)
static G4AntiNeutron * AntiNeutron()
virtual void DumpPhysicsTable(const G4ParticleDefinition &aParticleType)