Geant4_10
G4CrossSectionPairGG.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 // $Id: G4CrossSectionPairGG.cc 66241 2012-12-13 18:34:42Z gunter $
27 // $ GEANT4 tag $Name: not supported by cvs2svn $
28 //
29 // Class G4CrossSectionPairGG
30 //
31 // smoothly join two cross section sets by scaling the second at a given
32 // transition energy to match the first.
33 //
34 // Author: Gunter Folger
35 // November 2009
36 //
37 
38 #include "G4CrossSectionPairGG.hh"
39 
40 #include "globals.hh"
41 #include "G4PhysicalConstants.hh"
42 #include "G4SystemOfUnits.hh"
43 #include "G4HadTmpUtil.hh"
44 #include "G4NistManager.hh"
45 #include "G4ThreeVector.hh"
46 #include "G4NistManager.hh"
47 
48 G4CrossSectionPairGG::G4CrossSectionPairGG(G4VCrossSectionDataSet* low,
49  G4double Etransit) :
50  G4VCrossSectionDataSet("G4CrossSectionPairGG"), theLowX(low), ETransition(
51  Etransit) {
52  NistMan = G4NistManager::Instance();
53  theHighX = new G4GlauberGribovCrossSection();
54  verboseLevel = 0;
55 }
56 
58  // The cross section registry wil delete these
59  // delete theLowX;
60  // delete theHighX;
61 }
62 
64  std::ostream& outFile) const {
65  outFile << "G4CrossSectionPairGG is used to add the relativistic rise to\n"
66  << "hadronic cross section data sets above a given energy. In this\n"
67  << "case, the Glauber-Gribov cross section is used above 91 GeV.\n"
68  << "At this energy the low energy cross section is smoothly joined\n"
69  << "to the high energy cross section. Below 91 GeV the Barashenkov\n"
70  << "cross section is used for pions (G4PiNuclearCrossSection), the\n"
71  << "Axen-Wellisch cross section is used for protons\n"
72  << "(G4ProtonInelasticCrossSection), and the Wellisch-Laidlaw cross\n"
73  << "section is used for neutrons (G4NeutronInelasticCrossSection).\n";
74 }
75 
77  const G4DynamicParticle* aParticle, G4int Z, const G4Material* mat) {
78  G4bool isApplicable(false);
79  G4double Ekin = aParticle->GetKineticEnergy();
80  if (Ekin <= ETransition) {
81  isApplicable = theLowX->IsElementApplicable(aParticle, Z, mat);
82  } else if (Z > 1) {
83  isApplicable = true;
84  }
85  return isApplicable;
86 }
87 
89  const G4DynamicParticle* aParticle, G4int ZZ, const G4Material* mat)
90 {
91  G4double Xsec(0.);
92 
93  if (aParticle->GetKineticEnergy() < ETransition)
94  {
95  Xsec = theLowX->GetElementCrossSection(aParticle, ZZ, mat);
96  } else {
97 
98  std::vector<ParticleXScale>::iterator iter = scale_factors.begin();
99  G4ParticleDefinition * pDef = aParticle->GetDefinition();
100  while (iter != scale_factors.end() && (*iter).first != pDef)
101  {
102  ++iter;
103  }
104  if (iter != scale_factors.end() )
105  {
106  G4int AA = G4lrint(NistMan->GetAtomicMassAmu(ZZ));
107  Xsec = theHighX->GetInelasticGlauberGribov(aParticle, ZZ, AA)
108  * (*iter).second[ZZ];
109  if (verboseLevel > 2)
110  {
111  G4cout << " scaling .." << ZZ << " " << AA << " "
112  << (*iter).second[ZZ] << " "
113  << theHighX->GetInelasticGlauberGribov(aParticle, ZZ, AA)
114  << " " << Xsec << G4endl;
115  }
116  } else {
117  // BuildPhysicsTable not done for pDef=aParticle->GetDefinition
118  // build table, and recurse
119  BuildPhysicsTable(*pDef);
120  Xsec=GetElementCrossSection(aParticle, ZZ, mat);
121  }
122  }
123 
124  return Xsec;
125 }
126 
128  theLowX->BuildPhysicsTable(pDef);
129  theHighX->BuildPhysicsTable(pDef);
130 
131  if (verboseLevel > 0) {
132  G4cout << "G4CrossSectionPairGG::BuildPhysicsTable "
133  << theLowX->GetName() << " " << theHighX->GetName() << G4endl;
134  }
135 
136  G4ParticleDefinition * myDef = const_cast<G4ParticleDefinition*>(&pDef);
137  std::vector<ParticleXScale>::iterator iter;
138  iter = scale_factors.begin();
139  while (iter != scale_factors.end() && (*iter).first != myDef) {
140  ++iter;
141  }
142 
143  // new particle, initialise
144 
145  G4Material* mat = 0;
146 
147  if (iter == scale_factors.end()) {
148  XS_factors factors(93);
149  G4ThreeVector mom(0.0, 0.0, 1.0);
150  G4DynamicParticle DynPart(myDef, mom, ETransition); // last is kinetic Energy
151 
152  if (verboseLevel > 0) {
153  G4cout << "G4CrossSectionPairGG::BuildPhysicsTable for particle "
154  << pDef.GetParticleName() << G4endl;
155  }
156  for (G4int aZ = 1; aZ < 93; ++aZ) {
157  factors[aZ] = 1.; // default, to give reasonable value if only high is applicable
158  G4int AA = G4lrint(NistMan->GetAtomicMassAmu(aZ));
159  G4bool isApplicable = theLowX->IsElementApplicable(&DynPart, aZ,
160  mat) && (aZ > 1);
161 
162  if (isApplicable) {
163  factors[aZ] = theLowX->GetElementCrossSection(&DynPart, aZ, mat)
164  / theHighX->GetInelasticGlauberGribov(&DynPart, aZ, AA);
165 
166  }
167  if (verboseLevel > 0) {
168  G4cout << "Z=" << aZ << ", A=" << AA << ", scale="
169  << factors[aZ];
170  if (verboseLevel == 1) {
171  G4cout << G4endl;
172  } else {
173  if (isApplicable) {
174  G4cout << ", low / high "
175  << theLowX->GetElementCrossSection(&DynPart, aZ,
176  mat) << " "
177  << theHighX->GetInelasticGlauberGribov(&DynPart,
178  aZ, AA) << G4endl;
179  } else {
180  G4cout << ", N/A" << G4endl;
181  }
182  }
183  }
184  }
185  ParticleXScale forPart(myDef, factors);
186  scale_factors.push_back(forPart);
187  }
188 }
189 
190 /*
191  void G4CrossSectionPairGG::DumpHtml(const G4ParticleDefinition&,
192  std::ofstream outFile)
193  {
194  outFile << " <li><b>"
195  << " G4CrossSectionPairGG: " << theLowX->GetName() << " cross sections \n";
196  outFile << "below " << ETransition/GeV << " GeV, Glauber-Gribov above \n";
197  }
198  */
199 
201  G4cout << std::setw(24) << " " << " G4CrossSectionPairGG: "
202  << theLowX->GetName() << " cross sections " << G4endl;
203  G4cout << std::setw(27) << " " << "below " << ETransition / GeV
204  << " GeV, Glauber-Gribov above " << G4endl;
205 }
G4double GetKineticEnergy() const
std::ofstream outFile
Definition: GammaRayTel.cc:68
virtual void CrossSectionDescription(std::ostream &) const
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
virtual G4double GetElementCrossSection(const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
G4ParticleDefinition * GetDefinition() const
const G4String & GetName() const
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
const G4String & GetParticleName() const
Float_t mat
Definition: plot.C:40
virtual G4bool IsElementApplicable(const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
G4GLOB_DLL std::ostream G4cout
Float_t Z
Definition: plot.C:39
bool G4bool
Definition: G4Types.hh:79
virtual G4bool IsElementApplicable(const G4DynamicParticle *, G4int Z, const G4Material *mat=0)
G4double GetInelasticGlauberGribov(const G4DynamicParticle *, G4int Z, G4int A)
int G4lrint(double ad)
Definition: templates.hh:163
G4double GetAtomicMassAmu(const G4String &symb) const
virtual void DumpPhysicsTable(const G4ParticleDefinition &)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
virtual G4double GetElementCrossSection(const G4DynamicParticle *, G4int Z, const G4Material *mat=0)