Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** 25 // 25 // 26 // 26 // 27 // 27 // 28 // 28 // 29 // Geant4 Hadron Charge Exchange Process -- so 29 // Geant4 Hadron Charge Exchange Process -- source file 30 // 30 // 31 // Created 21 April 2006 V.Ivanchenko 31 // Created 21 April 2006 V.Ivanchenko 32 // 32 // 33 // Modified: 33 // Modified: 34 // 24-Apr-06 V.Ivanchenko add neutron scatteri 34 // 24-Apr-06 V.Ivanchenko add neutron scattering on hydrogen from CHIPS 35 // 07-Jun-06 V.Ivanchenko fix problem of rotat 35 // 07-Jun-06 V.Ivanchenko fix problem of rotation of final state 36 // 25-Jul-06 V.Ivanchenko add 19 MeV low energ 36 // 25-Jul-06 V.Ivanchenko add 19 MeV low energy for CHIPS 37 // 23-Jan-07 V.Ivanchenko add cross section in 37 // 23-Jan-07 V.Ivanchenko add cross section interfaces with Z and A 38 // and do not use CHIPS 38 // and do not use CHIPS for cross sections 39 // 14-Sep-12 M.Kelsey -- Pass subType code to 39 // 14-Sep-12 M.Kelsey -- Pass subType code to base ctor 40 // 06-Aug-15 A.Ribon migrating to G4Pow 40 // 06-Aug-15 A.Ribon migrating to G4Pow 41 41 42 #include "G4ChargeExchangeProcess.hh" 42 #include "G4ChargeExchangeProcess.hh" 43 #include "globals.hh" 43 #include "globals.hh" 44 #include "G4SystemOfUnits.hh" 44 #include "G4SystemOfUnits.hh" 45 #include "G4CrossSectionDataStore.hh" 45 #include "G4CrossSectionDataStore.hh" 46 #include "G4CrossSectionElastic.hh" << 46 #include "G4HadronElasticDataSet.hh" 47 #include "G4ComponentGGHadronNucleusXsc.hh" << 48 #include "G4Element.hh" 47 #include "G4Element.hh" 49 #include "G4ElementVector.hh" 48 #include "G4ElementVector.hh" 50 #include "G4IsotopeVector.hh" 49 #include "G4IsotopeVector.hh" 51 #include "G4Neutron.hh" 50 #include "G4Neutron.hh" 52 #include "G4Proton.hh" 51 #include "G4Proton.hh" 53 #include "G4PhysicsLinearVector.hh" 52 #include "G4PhysicsLinearVector.hh" 54 53 55 #include "G4Pow.hh" 54 #include "G4Pow.hh" 56 55 57 56 58 G4ChargeExchangeProcess::G4ChargeExchangeProce 57 G4ChargeExchangeProcess::G4ChargeExchangeProcess(const G4String& procName) 59 : G4HadronicProcess(procName,fChargeExchange 58 : G4HadronicProcess(procName,fChargeExchange), first(true) 60 { 59 { 61 G4cout << "###=== The class G4ChargeExchange << 62 G4cout << "###=== It will be removed at the << 63 thEnergy = 20.*MeV; 60 thEnergy = 20.*MeV; 64 pPDG = 0; 61 pPDG = 0; 65 verboseLevel= 1; 62 verboseLevel= 1; 66 AddDataSet( new G4CrossSectionElastic( new G << 63 AddDataSet(new G4HadronElasticDataSet); 67 theProton = G4Proton::Proton(); 64 theProton = G4Proton::Proton(); 68 theNeutron = G4Neutron::Neutron(); 65 theNeutron = G4Neutron::Neutron(); 69 theAProton = G4AntiProton::AntiProton(); 66 theAProton = G4AntiProton::AntiProton(); 70 theANeutron = G4AntiNeutron::AntiNeutron(); 67 theANeutron = G4AntiNeutron::AntiNeutron(); 71 thePiPlus = G4PionPlus::PionPlus(); 68 thePiPlus = G4PionPlus::PionPlus(); 72 thePiMinus = G4PionMinus::PionMinus(); 69 thePiMinus = G4PionMinus::PionMinus(); 73 thePiZero = G4PionZero::PionZero(); 70 thePiZero = G4PionZero::PionZero(); 74 theKPlus = G4KaonPlus::KaonPlus(); 71 theKPlus = G4KaonPlus::KaonPlus(); 75 theKMinus = G4KaonMinus::KaonMinus(); 72 theKMinus = G4KaonMinus::KaonMinus(); 76 theK0S = G4KaonZeroShort::KaonZeroShort 73 theK0S = G4KaonZeroShort::KaonZeroShort(); 77 theK0L = G4KaonZeroLong::KaonZeroLong() 74 theK0L = G4KaonZeroLong::KaonZeroLong(); 78 theL = G4Lambda::Lambda(); 75 theL = G4Lambda::Lambda(); 79 theAntiL = G4AntiLambda::AntiLambda(); 76 theAntiL = G4AntiLambda::AntiLambda(); 80 theSPlus = G4SigmaPlus::SigmaPlus(); 77 theSPlus = G4SigmaPlus::SigmaPlus(); 81 theASPlus = G4AntiSigmaPlus::AntiSigmaPlus 78 theASPlus = G4AntiSigmaPlus::AntiSigmaPlus(); 82 theSMinus = G4SigmaMinus::SigmaMinus(); 79 theSMinus = G4SigmaMinus::SigmaMinus(); 83 theASMinus = G4AntiSigmaMinus::AntiSigmaMin 80 theASMinus = G4AntiSigmaMinus::AntiSigmaMinus(); 84 theS0 = G4SigmaZero::SigmaZero(); 81 theS0 = G4SigmaZero::SigmaZero(); 85 theAS0 = G4AntiSigmaZero::AntiSigmaZero 82 theAS0 = G4AntiSigmaZero::AntiSigmaZero(); 86 theXiMinus = G4XiMinus::XiMinus(); 83 theXiMinus = G4XiMinus::XiMinus(); 87 theXi0 = G4XiZero::XiZero(); 84 theXi0 = G4XiZero::XiZero(); 88 theAXiMinus = G4AntiXiMinus::AntiXiMinus(); 85 theAXiMinus = G4AntiXiMinus::AntiXiMinus(); 89 theAXi0 = G4AntiXiZero::AntiXiZero(); 86 theAXi0 = G4AntiXiZero::AntiXiZero(); 90 theOmega = G4OmegaMinus::OmegaMinus(); 87 theOmega = G4OmegaMinus::OmegaMinus(); 91 theAOmega = G4AntiOmegaMinus::AntiOmegaMin 88 theAOmega = G4AntiOmegaMinus::AntiOmegaMinus(); 92 theD = G4Deuteron::Deuteron(); 89 theD = G4Deuteron::Deuteron(); 93 theT = G4Triton::Triton(); 90 theT = G4Triton::Triton(); 94 theA = G4Alpha::Alpha(); 91 theA = G4Alpha::Alpha(); 95 theHe3 = G4He3::He3(); 92 theHe3 = G4He3::He3(); 96 } 93 } 97 94 98 G4ChargeExchangeProcess::~G4ChargeExchangeProc 95 G4ChargeExchangeProcess::~G4ChargeExchangeProcess() 99 { 96 { 100 if (factors) delete factors; 97 if (factors) delete factors; 101 } 98 } 102 99 103 void G4ChargeExchangeProcess:: 100 void G4ChargeExchangeProcess:: 104 BuildPhysicsTable(const G4ParticleDefinition& 101 BuildPhysicsTable(const G4ParticleDefinition& aParticleType) 105 { 102 { 106 if(first) { 103 if(first) { 107 first = false; 104 first = false; 108 theParticle = &aParticleType; 105 theParticle = &aParticleType; 109 pPDG = theParticle->GetPDGEncoding(); 106 pPDG = theParticle->GetPDGEncoding(); 110 107 111 store = G4HadronicProcess::GetCrossSection 108 store = G4HadronicProcess::GetCrossSectionDataStore(); 112 109 113 const size_t n = 10; 110 const size_t n = 10; 114 if(theParticle == thePiPlus || theParticle 111 if(theParticle == thePiPlus || theParticle == thePiMinus || 115 theParticle == theKPlus || theParticle 112 theParticle == theKPlus || theParticle == theKMinus || 116 theParticle == theK0S || theParticle 113 theParticle == theK0S || theParticle == theK0L) { 117 114 118 G4double F[n] = {0.33,0.27,0.29,0.31,0.2 115 G4double F[n] = {0.33,0.27,0.29,0.31,0.27,0.18,0.13,0.1,0.09,0.07}; 119 factors = new G4PhysicsLinearVector(0.0, 116 factors = new G4PhysicsLinearVector(0.0,2.0*GeV,n); 120 for(size_t i=0; i<n; i++) {factors->PutV 117 for(size_t i=0; i<n; i++) {factors->PutValue(i,F[i]);} 121 118 122 } else { 119 } else { 123 120 124 G4double F[n] = {0.50,0.45,0.40,0.35,0.3 121 G4double F[n] = {0.50,0.45,0.40,0.35,0.30,0.25,0.06,0.04,0.005,0.0}; 125 factors = new G4PhysicsLinearVector(0.0, 122 factors = new G4PhysicsLinearVector(0.0,4.0*GeV,n); 126 for(size_t i=0; i<n; i++) {factors->PutV 123 for(size_t i=0; i<n; i++) {factors->PutValue(i,F[i]);} 127 } 124 } >> 125 //factors->SetSpline(true); 128 126 129 if(verboseLevel>1) 127 if(verboseLevel>1) 130 G4cout << "G4ChargeExchangeProcess for " 128 G4cout << "G4ChargeExchangeProcess for " 131 << theParticle->GetParticleName() 129 << theParticle->GetParticleName() 132 << G4endl; 130 << G4endl; 133 } 131 } 134 G4HadronicProcess::BuildPhysicsTable(aPartic 132 G4HadronicProcess::BuildPhysicsTable(aParticleType); 135 } 133 } 136 134 137 G4double G4ChargeExchangeProcess::GetElementCr 135 G4double G4ChargeExchangeProcess::GetElementCrossSection( 138 const G4Dyna 136 const G4DynamicParticle* dp, 139 const G4Element* elm, 137 const G4Element* elm, 140 const G4Material* mat) 138 const G4Material* mat) 141 { 139 { 142 // gives the microscopic cross section in GE 140 // gives the microscopic cross section in GEANT4 internal units 143 G4double Z = elm->GetZ(); 141 G4double Z = elm->GetZ(); 144 G4int iz = G4int(Z); 142 G4int iz = G4int(Z); 145 G4double x = 0.0; 143 G4double x = 0.0; 146 144 147 // The process is effective only above the t 145 // The process is effective only above the threshold 148 if(iz == 1 || dp->GetKineticEnergy() < thEne 146 if(iz == 1 || dp->GetKineticEnergy() < thEnergy) return x; 149 147 150 if(verboseLevel>1) 148 if(verboseLevel>1) 151 G4cout << "G4ChargeExchangeProcess compute 149 G4cout << "G4ChargeExchangeProcess compute GHAD CS for element " 152 << elm->GetName() 150 << elm->GetName() 153 << G4endl; 151 << G4endl; 154 x = store->GetCrossSection(dp, elm, mat); 152 x = store->GetCrossSection(dp, elm, mat); 155 153 156 if(verboseLevel>1) 154 if(verboseLevel>1) 157 G4cout << "G4ChargeExchangeProcess cross(m 155 G4cout << "G4ChargeExchangeProcess cross(mb)= " << x/millibarn 158 << " E(MeV)= " << dp->GetKineticEn 156 << " E(MeV)= " << dp->GetKineticEnergy() 159 << " " << theParticle->GetParticleName() 157 << " " << theParticle->GetParticleName() 160 << " in Z= " << iz 158 << " in Z= " << iz 161 << G4endl; 159 << G4endl; 162 G4bool b; 160 G4bool b; 163 G4double A = elm->GetN(); 161 G4double A = elm->GetN(); 164 G4double ptot = dp->GetTotalMomentum(); 162 G4double ptot = dp->GetTotalMomentum(); 165 x *= factors->GetValue(ptot, b)/G4Pow::GetIn 163 x *= factors->GetValue(ptot, b)/G4Pow::GetInstance()->powA(A, 0.42); 166 if(theParticle == thePiPlus || theParticle = 164 if(theParticle == thePiPlus || theParticle == theProton || 167 theParticle == theKPlus || theParticle = 165 theParticle == theKPlus || theParticle == theANeutron) 168 { x *= (1.0 - Z/A); } 166 { x *= (1.0 - Z/A); } 169 167 170 else if(theParticle == thePiMinus || thePart 168 else if(theParticle == thePiMinus || theParticle == theNeutron || 171 theParticle == theKMinus || thePart 169 theParticle == theKMinus || theParticle == theAProton) 172 { x *= Z/A; } 170 { x *= Z/A; } 173 171 174 if(theParticle->GetPDGMass() < GeV) { 172 if(theParticle->GetPDGMass() < GeV) { 175 if(ptot > 2.*GeV) x *= 4.0*GeV*GeV/(ptot*p 173 if(ptot > 2.*GeV) x *= 4.0*GeV*GeV/(ptot*ptot); 176 } 174 } 177 175 178 if(verboseLevel>1) 176 if(verboseLevel>1) 179 G4cout << "Corrected cross(mb)= " << x/mil 177 G4cout << "Corrected cross(mb)= " << x/millibarn << G4endl; 180 178 181 return x; 179 return x; 182 } 180 } 183 181 184 G4bool G4ChargeExchangeProcess:: 182 G4bool G4ChargeExchangeProcess:: 185 IsApplicable(const G4ParticleDefinition& aPart 183 IsApplicable(const G4ParticleDefinition& aParticleType) 186 { 184 { 187 const G4ParticleDefinition* p = &aParticleTy 185 const G4ParticleDefinition* p = &aParticleType; 188 return (p == thePiPlus || p == thePiMinus || 186 return (p == thePiPlus || p == thePiMinus || 189 p == theProton || p == theNeutron || 187 p == theProton || p == theNeutron || 190 p == theAProton|| p == theANeutron|| 188 p == theAProton|| p == theANeutron|| 191 p == theKPlus || p == theKMinus || 189 p == theKPlus || p == theKMinus || 192 p == theK0S || p == theK0L || 190 p == theK0S || p == theK0L || 193 p == theL); 191 p == theL); 194 } 192 } 195 193 196 void G4ChargeExchangeProcess:: 194 void G4ChargeExchangeProcess:: 197 DumpPhysicsTable(const G4ParticleDefinition& a 195 DumpPhysicsTable(const G4ParticleDefinition& aParticleType) 198 { 196 { 199 store->DumpPhysicsTable(aParticleType); 197 store->DumpPhysicsTable(aParticleType); 200 } 198 } 201 199