Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/cross_sections/src/G4BGGPionInelasticXS.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 //
 28 // GEANT4 Class file
 29 //
 30 //
 31 // File name:     G4BGGPionInelasticXS
 32 //
 33 // Author:        Vladimir Ivanchenko
 34 //
 35 // Creation date: 01.10.2003
 36 // Modifications:
 37 //
 38 // -------------------------------------------------------------------
 39 //
 40 
 41 #include "G4BGGPionInelasticXS.hh"
 42 #include "G4SystemOfUnits.hh"
 43 #include "G4ComponentGGHadronNucleusXsc.hh"
 44 #include "G4UPiNuclearCrossSection.hh"
 45 #include "G4HadronNucleonXsc.hh"
 46 #include "G4NuclearRadii.hh"
 47 
 48 #include "G4Proton.hh"
 49 #include "G4PionPlus.hh"
 50 #include "G4PionMinus.hh"
 51 #include "G4NistManager.hh"
 52 #include "G4Pow.hh"
 53 #include "G4HadronicParameters.hh"
 54 
 55 G4double G4BGGPionInelasticXS::theGlauberFacPiPlus[93] = {0.0};
 56 G4double G4BGGPionInelasticXS::theGlauberFacPiMinus[93] = {0.0};
 57 G4double G4BGGPionInelasticXS::theLowEPiPlus[93] = {0.0};
 58 G4double G4BGGPionInelasticXS::theLowEPiMinus[93] = {0.0};
 59 G4int    G4BGGPionInelasticXS::theA[93] = {0};
 60 
 61 G4BGGPionInelasticXS::G4BGGPionInelasticXS(const G4ParticleDefinition* p) 
 62  : G4VCrossSectionDataSet("BarashenkovGlauberGribov")
 63 {
 64   verboseLevel = 0;
 65   fGlauberEnergy = 91.*CLHEP::GeV;
 66   fLowEnergy = 20.*CLHEP::MeV;
 67   fLowestEnergy = 1.*CLHEP::MeV;
 68   SetMinKinEnergy(0.0);
 69   SetMaxKinEnergy( G4HadronicParameters::Instance()->GetMaxEnergy() );
 70 
 71   fPion    = new G4UPiNuclearCrossSection();
 72   fGlauber = new G4ComponentGGHadronNucleusXsc();
 73   fHadron  = new G4HadronNucleonXsc();
 74 
 75   fG4pow   = G4Pow::GetInstance();
 76 
 77   theProton = G4Proton::Proton();
 78   thePiPlus = G4PionPlus::PionPlus();
 79   isPiplus = (p == thePiPlus);
 80   SetForAllAtomsAndEnergies(true);
 81 
 82   if (0 == theA[0]) { Initialise(); } 
 83 }
 84 
 85 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 86 
 87 G4BGGPionInelasticXS::~G4BGGPionInelasticXS()
 88 {
 89   delete fHadron;
 90 }
 91 
 92 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 93 
 94 G4bool 
 95 G4BGGPionInelasticXS::IsElementApplicable(const G4DynamicParticle*, G4int,
 96                                           const G4Material*)
 97 {
 98   return true;
 99 }
100 
101 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
102 
103 G4bool G4BGGPionInelasticXS::IsIsoApplicable(const G4DynamicParticle*, 
104                                              G4int Z, G4int,  
105                                              const G4Element*,
106                                              const G4Material*)
107 {
108   return (1 == Z);
109 }
110 
111 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
112 
113 G4double
114 G4BGGPionInelasticXS::GetElementCrossSection(const G4DynamicParticle* dp,
115                                              G4int ZZ, const G4Material*)
116 {
117   // this method should be called only for Z > 1
118 
119   G4double cross = 0.0;
120   G4double ekin = std::max(dp->GetKineticEnergy(), fLowestEnergy);
121   G4int Z = std::min(ZZ, 92);
122 
123   if(1 == Z) {
124     cross = 1.0115*GetIsoCrossSection(dp,1,1);
125   } else if(ekin < fLowEnergy) {
126     cross = (isPiplus) ? theLowEPiPlus[Z]*CoulombFactorPiPlus(ekin, Z)
127       : theLowEPiMinus[Z]*FactorPiMinus(ekin);
128   } else if(ekin > fGlauberEnergy) {
129     cross = (isPiplus) ? theGlauberFacPiPlus[Z] : theGlauberFacPiMinus[Z];
130     cross *= fGlauber->GetInelasticGlauberGribov(dp, Z, theA[Z]);
131   } else {
132     cross = fPion->GetInelasticCrossSection(dp, Z, theA[Z]);
133   }
134 #ifdef G4VERBOSE
135   if(verboseLevel > 1) {
136     G4cout << "G4BGGPionInelasticXS::GetCrossSection  for "
137            << dp->GetDefinition()->GetParticleName()
138            << "  Ekin(GeV)= " << dp->GetKineticEnergy()
139            << " in nucleus Z= " << Z << "  A= " << theA[Z]
140            << " XS(b)= " << cross/barn 
141            << G4endl;
142   }
143 #endif
144   return cross;
145 }
146 
147 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
148 
149 G4double
150 G4BGGPionInelasticXS::GetIsoCrossSection(const G4DynamicParticle* dp, 
151                                          G4int, G4int A, 
152                                          const G4Isotope*,
153                                          const G4Element*,
154                                          const G4Material*)
155 {
156   // this method should be called only for Z = 1
157   fHadron->HadronNucleonXscNS(dp->GetDefinition(), theProton, 
158                               dp->GetKineticEnergy());
159   G4double cross = A*fHadron->GetInelasticHadronNucleonXsc();
160 
161 #ifdef G4VERBOSE
162   if(verboseLevel > 1) {
163     G4cout << "G4BGGPionInelasticXS::GetCrossSection  for "
164            << dp->GetDefinition()->GetParticleName()
165            << "  Ekin(GeV)= " << dp->GetKineticEnergy()
166            << " in nucleus Z=1  A=" << A
167            << " XS(b)= " << cross/barn 
168            << G4endl;
169   }
170 #endif
171   return cross;
172 }
173 
174 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
175 
176 void G4BGGPionInelasticXS::BuildPhysicsTable(const G4ParticleDefinition& p)
177 {
178   if(verboseLevel > 1) {
179     G4cout << "G4BGGPionInelasticXS::BuildPhysicsTable for " 
180            << p.GetParticleName() << G4endl;
181   } 
182   if(&p == G4PionPlus::PionPlus() || &p == G4PionMinus::PionMinus()) {
183     isPiplus = (&p == G4PionPlus::PionPlus());
184   } else {
185     G4ExceptionDescription ed;
186     ed << "This BGG cross section is applicable only to pions and not to " 
187        << p.GetParticleName() << G4endl; 
188     G4Exception("G4BGGPionInelasticXS::BuildPhysicsTable", "had001", 
189                 FatalException, ed);
190   }
191 }
192 
193 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
194 
195 void G4BGGPionInelasticXS::Initialise()
196 {
197   theA[0] = theA[1] = 1;
198   G4ThreeVector mom(0.0,0.0,1.0);
199   G4DynamicParticle dp(thePiPlus, mom, fGlauberEnergy);
200 
201   G4NistManager* nist = G4NistManager::Instance();
202   G4double csup, csdn;
203 
204   for (G4int iz=2; iz<93; ++iz) {
205     G4int A = G4lrint(nist->GetAtomicMassAmu(iz));
206     theA[iz] = A;
207 
208     csup = fGlauber->GetInelasticGlauberGribov(&dp, iz, A);
209     csdn = fPion->GetInelasticCrossSection(&dp, iz, A);
210     theGlauberFacPiPlus[iz] = csdn/csup;
211   }
212 
213   dp.SetDefinition(G4PionMinus::PionMinus());
214   for (G4int iz=2; iz<93; ++iz) {
215     csup = fGlauber->GetInelasticGlauberGribov(&dp, iz, theA[iz]);
216     csdn = fPion->GetInelasticCrossSection(&dp, iz, theA[iz]);
217     theGlauberFacPiMinus[iz] = csdn/csup;
218 
219     if(verboseLevel > 1) {
220       G4cout << "Z= " << iz <<  "  A= " << theA[iz] 
221        << " factorPiPlus= " << theGlauberFacPiPlus[iz] 
222        << " factorPiMinus= " << theGlauberFacPiMinus[iz] 
223        << G4endl;
224     }
225   }
226 
227   theLowEPiPlus[1] = theLowEPiMinus[1]= 1.0;
228   dp.SetDefinition(thePiPlus);
229   dp.SetKineticEnergy(fLowEnergy);
230   for (G4int iz=2; iz<93; ++iz) {
231     theLowEPiPlus[iz] = fPion->GetInelasticCrossSection(&dp, iz, theA[iz])
232       /CoulombFactorPiPlus(fLowEnergy, iz);
233   }
234 
235   dp.SetDefinition(G4PionMinus::PionMinus());
236   for (G4int iz=2; iz<93; ++iz) {
237     theLowEPiMinus[iz] = fPion->GetInelasticCrossSection(&dp, iz, theA[iz])
238       /FactorPiMinus(fLowEnergy);
239     
240     if (verboseLevel > 1) {
241       G4cout << "Z= " << iz <<  "  A= " << theA[iz] 
242        << " LowEtorPiPlus= " << theLowEPiPlus[iz] 
243        << " LowEtorPiMinus= " << theLowEPiMinus[iz] 
244        << G4endl;
245     }
246   }
247 }
248 
249 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
250 
251 G4double G4BGGPionInelasticXS::CoulombFactorPiPlus(G4double kinEnergy, G4int Z)
252 {
253   return (kinEnergy > 0.0) ? 
254     G4NuclearRadii::CoulombFactor(Z, theA[Z], thePiPlus, kinEnergy) : 0.0;
255 }
256 
257 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
258 
259 G4double G4BGGPionInelasticXS::FactorPiMinus(G4double kinEnergy)
260 {
261   return 1.0/std::sqrt(kinEnergy);
262 }
263 
264 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
265 
266 void
267 G4BGGPionInelasticXS::CrossSectionDescription(std::ostream& outFile) const 
268 {
269   outFile << "The Barashenkov-Glauber-Gribov cross section handles inelastic\n"
270           << "pion scattering from nuclei at all energies.  The Barashenkov\n"
271           << "parameterization is used below 91 GeV and the Glauber-Gribov\n"
272           << "parameterization is used above 91 GeV.\n";
273 }
274 
275 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
276