00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038 #include "G4CrossSectionPairGG.hh"
00039
00040 #include "globals.hh"
00041 #include "G4PhysicalConstants.hh"
00042 #include "G4SystemOfUnits.hh"
00043 #include "G4HadTmpUtil.hh"
00044 #include "G4NistManager.hh"
00045 #include "G4ThreeVector.hh"
00046 #include "G4NistManager.hh"
00047
00048 G4CrossSectionPairGG::G4CrossSectionPairGG(G4VCrossSectionDataSet* low,
00049 G4double Etransit) :
00050 G4VCrossSectionDataSet("G4CrossSectionPairGG"), theLowX(low), ETransition(
00051 Etransit) {
00052 NistMan = G4NistManager::Instance();
00053 theHighX = new G4GlauberGribovCrossSection();
00054 verboseLevel = 0;
00055 }
00056
00057 G4CrossSectionPairGG::~G4CrossSectionPairGG() {
00058
00059
00060
00061 }
00062
00063 void G4CrossSectionPairGG::CrossSectionDescription(
00064 std::ostream& outFile) const {
00065 outFile << "G4CrossSectionPairGG is used to add the relativistic rise to\n"
00066 << "hadronic cross section data sets above a given energy. In this\n"
00067 << "case, the Glauber-Gribov cross section is used above 91 GeV.\n"
00068 << "At this energy the low energy cross section is smoothly joined\n"
00069 << "to the high energy cross section. Below 91 GeV the Barashenkov\n"
00070 << "cross section is used for pions (G4PiNuclearCrossSection), the\n"
00071 << "Axen-Wellisch cross section is used for protons\n"
00072 << "(G4ProtonInelasticCrossSection), and the Wellisch-Laidlaw cross\n"
00073 << "section is used for neutrons (G4NeutronInelasticCrossSection).\n";
00074 }
00075
00076 G4bool G4CrossSectionPairGG::IsElementApplicable(
00077 const G4DynamicParticle* aParticle, G4int Z, const G4Material* mat) {
00078 G4bool isApplicable(false);
00079 G4double Ekin = aParticle->GetKineticEnergy();
00080 if (Ekin <= ETransition) {
00081 isApplicable = theLowX->IsElementApplicable(aParticle, Z, mat);
00082 } else if (Z > 1) {
00083 isApplicable = true;
00084 }
00085 return isApplicable;
00086 }
00087
00088 G4double G4CrossSectionPairGG::GetElementCrossSection(
00089 const G4DynamicParticle* aParticle, G4int ZZ, const G4Material* mat)
00090 {
00091 G4double Xsec(0.);
00092
00093 if (aParticle->GetKineticEnergy() < ETransition)
00094 {
00095 Xsec = theLowX->GetElementCrossSection(aParticle, ZZ, mat);
00096 } else {
00097
00098 std::vector<ParticleXScale>::iterator iter = scale_factors.begin();
00099 G4ParticleDefinition * pDef = aParticle->GetDefinition();
00100 while (iter != scale_factors.end() && (*iter).first != pDef)
00101 {
00102 ++iter;
00103 }
00104 if (iter != scale_factors.end() )
00105 {
00106 G4int AA = G4lrint(NistMan->GetAtomicMassAmu(ZZ));
00107 Xsec = theHighX->GetInelasticGlauberGribov(aParticle, ZZ, AA)
00108 * (*iter).second[ZZ];
00109 if (verboseLevel > 2)
00110 {
00111 G4cout << " scaling .." << ZZ << " " << AA << " "
00112 << (*iter).second[ZZ] << " "
00113 << theHighX->GetInelasticGlauberGribov(aParticle, ZZ, AA)
00114 << " " << Xsec << G4endl;
00115 }
00116 } else {
00117
00118
00119 BuildPhysicsTable(*pDef);
00120 Xsec=GetElementCrossSection(aParticle, ZZ, mat);
00121 }
00122 }
00123
00124 return Xsec;
00125 }
00126
00127 void G4CrossSectionPairGG::BuildPhysicsTable(const G4ParticleDefinition& pDef) {
00128 theLowX->BuildPhysicsTable(pDef);
00129 theHighX->BuildPhysicsTable(pDef);
00130
00131 if (verboseLevel > 0) {
00132 G4cout << "G4CrossSectionPairGG::BuildPhysicsTable "
00133 << theLowX->GetName() << " " << theHighX->GetName() << G4endl;
00134 }
00135
00136 G4ParticleDefinition * myDef = const_cast<G4ParticleDefinition*>(&pDef);
00137 std::vector<ParticleXScale>::iterator iter;
00138 iter = scale_factors.begin();
00139 while (iter != scale_factors.end() && (*iter).first != myDef) {
00140 ++iter;
00141 }
00142
00143
00144
00145 G4Material* mat = 0;
00146
00147 if (iter == scale_factors.end()) {
00148 XS_factors factors(93);
00149 G4ThreeVector mom(0.0, 0.0, 1.0);
00150 G4DynamicParticle DynPart(myDef, mom, ETransition);
00151
00152 if (verboseLevel > 0) {
00153 G4cout << "G4CrossSectionPairGG::BuildPhysicsTable for particle "
00154 << pDef.GetParticleName() << G4endl;
00155 }
00156 for (G4int aZ = 1; aZ < 93; ++aZ) {
00157 factors[aZ] = 1.;
00158 G4int AA = G4lrint(NistMan->GetAtomicMassAmu(aZ));
00159 G4bool isApplicable = theLowX->IsElementApplicable(&DynPart, aZ,
00160 mat) && (aZ > 1);
00161
00162 if (isApplicable) {
00163 factors[aZ] = theLowX->GetElementCrossSection(&DynPart, aZ, mat)
00164 / theHighX->GetInelasticGlauberGribov(&DynPart, aZ, AA);
00165
00166 }
00167 if (verboseLevel > 0) {
00168 G4cout << "Z=" << aZ << ", A=" << AA << ", scale="
00169 << factors[aZ];
00170 if (verboseLevel == 1) {
00171 G4cout << G4endl;
00172 } else {
00173 if (isApplicable) {
00174 G4cout << ", low / high "
00175 << theLowX->GetElementCrossSection(&DynPart, aZ,
00176 mat) << " "
00177 << theHighX->GetInelasticGlauberGribov(&DynPart,
00178 aZ, AA) << G4endl;
00179 } else {
00180 G4cout << ", N/A" << G4endl;
00181 }
00182 }
00183 }
00184 }
00185 ParticleXScale forPart(myDef, factors);
00186 scale_factors.push_back(forPart);
00187 }
00188 }
00189
00190
00191
00192
00193
00194
00195
00196
00197
00198
00199
00200 void G4CrossSectionPairGG::DumpPhysicsTable(const G4ParticleDefinition&) {
00201 G4cout << std::setw(24) << " " << " G4CrossSectionPairGG: "
00202 << theLowX->GetName() << " cross sections " << G4endl;
00203 G4cout << std::setw(27) << " " << "below " << ETransition / GeV
00204 << " GeV, Glauber-Gribov above " << G4endl;
00205 }