Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/management/src/G4EnergyRangeManager.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  // Hadronic Process: Energy Range Manager
 28  // original by H.P. Wellisch
 29  // modified by J.L. Chuma, TRIUMF, 22-Nov-1996
 30  // Last modified: 24-Mar-1997
 31  // fix in the counter-hndling: H.P. Wellisch 04-Apr-97
 32  // throw an exception if no model found:  J.L. Chuma  04-Apr-97
 33  
 34 #include "G4EnergyRangeManager.hh"
 35 #include "Randomize.hh"
 36 #include "G4HadronicException.hh"
 37 #include "G4SystemOfUnits.hh"
 38 
 39 G4EnergyRangeManager::G4EnergyRangeManager()
 40  : theHadronicInteractionCounter(0)
 41 {}
 42 
 43 G4EnergyRangeManager::~G4EnergyRangeManager()
 44 {}
 45 
 46 void G4EnergyRangeManager::RegisterMe(G4HadronicInteraction* a)
 47 {
 48   if(!a) { return; }
 49   if(0 < theHadronicInteractionCounter) {
 50     for(G4int i=0; i<theHadronicInteractionCounter; ++i) {
 51       if(a == theHadronicInteraction[i]) { return; }
 52     }
 53   }
 54   theHadronicInteraction.push_back(a);
 55   ++theHadronicInteractionCounter;
 56 }
 57 
 58 G4HadronicInteraction*
 59 G4EnergyRangeManager::GetHadronicInteraction(const G4HadProjectile & aHadProjectile, 
 60                                              G4Nucleus & aTargetNucleus,
 61                                              const G4Material* aMaterial,
 62                                              const G4Element* anElement) const
 63 {
 64   // VI shortcut: if only one interaction is registered skip all checks
 65   if(1 == theHadronicInteractionCounter) { return theHadronicInteraction[0]; }
 66   else if(0 == theHadronicInteractionCounter) {
 67     G4cout << "G4EnergyRangeManager::GetHadronicInteraction: "
 68      << "no models defined for a process" << G4endl;
 69     return nullptr;
 70   }
 71 
 72   G4double kineticEnergy = aHadProjectile.GetKineticEnergy();
 73   // For ions, get kinetic energy per nucleon
 74   if ( std::abs( aHadProjectile.GetDefinition()->GetBaryonNumber() ) > 1 ) {
 75     kineticEnergy /= static_cast< G4double >( std::abs( aHadProjectile.GetDefinition()->GetBaryonNumber() ) );
 76   }
 77 
 78   G4int cou = 0, memory = 0, memor2 = 0;
 79   G4double emi1 = 0.0, ema1 = 0.0, emi2 = 0.0, ema2 = 0.0;
 80 
 81   for (G4int i = 0; i<theHadronicInteractionCounter; ++i) {
 82     if ( theHadronicInteraction[i]->IsApplicable( aHadProjectile, aTargetNucleus ) ) {
 83       G4double low  = theHadronicInteraction[i]->GetMinEnergy( aMaterial, anElement );
 84       G4double high = theHadronicInteraction[i]->GetMaxEnergy( aMaterial, anElement );
 85       if (low <= kineticEnergy && high >= kineticEnergy) {
 86         ++cou;
 87         emi2 = emi1;
 88         ema2 = ema1;
 89         emi1 = low;
 90         ema1 = high;
 91         memor2 = memory;
 92         memory = i;
 93       }
 94     }
 95   }
 96 
 97   G4HadronicInteraction* hi = nullptr;
 98   switch (cou) {
 99     case 0:
100       G4cout << "No model found out of " << theHadronicInteractionCounter << G4endl;
101       for( G4int j=0; j<theHadronicInteractionCounter; ++j) {
102   G4HadronicInteraction* hint=theHadronicInteraction[j];
103   G4cout << "   "<< j << ".  Elow= " << hint->GetMinEnergy(aMaterial,anElement)
104          <<", Ehigh= " << hint->GetMaxEnergy(aMaterial,anElement)
105          <<"   " << hint->GetModelName() << G4endl;
106       }
107       break;
108 
109     case 1:
110       hi = theHadronicInteraction[memory];
111       break;
112 
113     case 2:
114       if( (emi2<=emi1 && ema2>=ema1) || (emi2>=emi1 && ema2<=ema1) ) {
115   G4cout << "Energy ranges of two models fully overlapping " << G4endl;
116   for( G4int j=0; j<theHadronicInteractionCounter; ++j) {
117     G4HadronicInteraction* hint=theHadronicInteraction[j];
118     G4cout << "   "<< j << ".  Elow= " << hint->GetMinEnergy(aMaterial,anElement)
119      <<", Ehigh= " << hint->GetMaxEnergy(aMaterial,anElement)
120      <<"   " << hint->GetModelName() << G4endl;
121   }
122       } else {
123   G4double rand = G4UniformRand();
124   G4int mem;
125   if( emi1 < emi2 ) {
126     if( (ema1-kineticEnergy) < rand*(ema1-emi2) ) {
127       mem = memor2;
128     } else {
129       mem = memory;
130     }
131   } else {
132     if( (ema2-kineticEnergy) < rand*(ema2-emi1) ) {
133       mem = memory;
134     } else {
135       mem = memor2;
136     }
137   }
138   hi = theHadronicInteraction[mem];
139       }
140       break;
141 
142     default:
143       G4cout << "More than two competing models for this energy" << G4endl;
144       for( G4int j=0; j<theHadronicInteractionCounter; ++j) {
145   G4HadronicInteraction* hint=theHadronicInteraction[j];
146   G4cout << "   "<< j << ".  Elow= " << hint->GetMinEnergy(aMaterial,anElement)
147          <<", Ehigh= " << hint->GetMaxEnergy(aMaterial,anElement)
148          <<"   " << hint->GetModelName() << G4endl;
149       }
150       break;
151   }
152   return hi;
153 } 
154 
155 std::vector<G4HadronicInteraction*>& 
156 G4EnergyRangeManager::GetHadronicInteractionList()
157 {
158   return theHadronicInteraction;
159 }
160 
161 void G4EnergyRangeManager::Dump( G4int verbose )
162 {
163   G4cout << "G4EnergyRangeManager " << this << G4endl;
164   for (G4int i = 0 ; i < theHadronicInteractionCounter; i++) {
165     G4cout << "   HadronicModel " << i <<":"
166            << theHadronicInteraction[i]->GetModelName() << G4endl;
167     if (verbose > 0) {
168       G4cout << "      Minimum Energy " 
169        << theHadronicInteraction[i]->GetMinEnergy()/GeV << " [GeV], "
170              << "Maximum Energy " 
171        << theHadronicInteraction[i]->GetMaxEnergy()/GeV << " [GeV]"
172              << G4endl;
173     }
174   }
175 }
176 
177 void
178 G4EnergyRangeManager::BuildPhysicsTable(const G4ParticleDefinition& aParticleType)
179 {
180   for (auto & hadi : theHadronicInteraction) {
181     hadi->BuildPhysicsTable( aParticleType );
182   }
183 }
184 
185  
186