Geant4 Cross Reference |
1 // 1 2 // ******************************************* 3 // * License and Disclaimer 4 // * 5 // * The Geant4 software is copyright of th 6 // * the Geant4 Collaboration. It is provided 7 // * conditions of the Geant4 Software License 8 // * LICENSE and available at http://cern.ch/ 9 // * include a list of copyright holders. 10 // * 11 // * Neither the authors of this software syst 12 // * institutes,nor the agencies providing fin 13 // * work make any representation or warran 14 // * regarding this software system or assum 15 // * use. Please see the license in the file 16 // * for the full disclaimer and the limitatio 17 // * 18 // * This code implementation is the result 19 // * technical work of the GEANT4 collaboratio 20 // * By using, copying, modifying or distri 21 // * any work based on the software) you ag 22 // * use in resulting scientific publicati 23 // * acceptance of all terms of the Geant4 Sof 24 // ******************************************* 25 // 26 #include "globals.hh" 27 #include "G4ios.hh" 28 #include "G4Exp.hh" 29 #include "G4Log.hh" 30 #include "G4SystemOfUnits.hh" 31 #include "G4XNNElasticLowE.hh" 32 #include "G4KineticTrack.hh" 33 #include "G4ParticleDefinition.hh" 34 #include "G4PhysicsVector.hh" 35 #include "G4PhysicsLogVector.hh" 36 #include "G4Proton.hh" 37 #include "G4Neutron.hh" 38 39 const G4double G4XNNElasticLowE::_lowLimit = 0 40 const G4double G4XNNElasticLowE::_highLimit = 41 42 // Low energy limit of the cross-section table 43 // Units are assigned when filling the Physics 44 const G4double G4XNNElasticLowE::_eMinTable = 45 const G4double G4XNNElasticLowE::_eStepLog = 0 46 47 // Cross-sections in mb 48 // Units are assigned when filling the Physics 49 50 const G4int G4XNNElasticLowE::tableSize = 101; 51 52 const G4double G4XNNElasticLowE::ppTable[101] 53 { 54 60.00, //was 0. 55 33.48, 26.76, 25.26, 24.55, 23.94, 23.77, 23 56 25.48, 27.52, 27.72, 27.21, 25.80, 26.00, 24 57 24.37, 24.36, 23.13, 22.43, 21.71, 21.01, 20 58 20.25, 20.10, 20.59, 20.04, 20.83, 20.84, 21 59 20.79, 21.88, 21.15, 20.92, 19.00, 18.60, 17 60 16.70, 16.50, 16.20, 15.80, 15.57, 15.20, 15 61 14.20, 14.00, 13.80, 13.60, 13.40, 13.20, 13 62 12.70, 12.60, 12.50, 12.40, 12.30, 12.20, 12 63 11.90, 11.80, 11.75, 11.70, 11.64, 11.53, 11 64 11.22, 11.13, 11.05, 10.97, 10.89, 10.82, 10 65 10.61, 10.54, 10.48, 10.41, 10.35, 10.28, 10 66 10.13, 10.10, 10.08, 10.05, 10.02, 9.99, 9 67 9.90, 9.87, 9.84, 9.80 68 }; 69 70 const G4double G4XNNElasticLowE::npTable[101] 71 { 72 1500.00, // was 0. 73 248.20, 93.38, 55.26, 44.50, 41.33, 38.48, 3 74 35.02, 34.47, 32.48, 30.76, 29.46, 28.53, 27 75 26.53, 25.95, 25.59, 25.46, 25.00, 24.49, 24 76 23.17, 22.70, 21.88, 21.48, 20.22, 19.75, 18 77 17.98, 17.63, 17.21, 16.72, 16.68, 16.58, 16 78 15.98, 15.71, 15.42, 15.14, 14.87, 14.65, 14 79 14.10, 13.95, 13.80, 13.64, 13.47, 13.29, 13 80 12.68, 12.47, 12.27, 12.06, 11.84, 11.76, 11 81 11.50, 11.41, 11.29, 11.17, 11.06, 10.93, 10 82 10.56, 10.44, 10.33, 10.21, 10.12, 10.03, 9 83 9.83, 9.80, 9.77, 9.75, 9.74, 9.74, 9. 84 9.73, 9.70, 9.68, 9.65, 9.63, 9.60, 9. 85 9.52, 9.49, 9.46, 9.43 86 }; 87 88 89 G4XNNElasticLowE::G4XNNElasticLowE() 90 { 91 // Cross-sections are available in the range 92 93 _eMin = _eMinTable * GeV; 94 _eMax = G4Exp(G4Log(_eMinTable) + tableSize 95 if (_eMin < _lowLimit) 96 throw G4HadronicException(__FILE__, __LINE 97 if (_highLimit > _eMax) 98 throw G4HadronicException(__FILE__, __LINE 99 G4PhysicsVector* pp = new G4PhysicsLogVector 100 101 _eMin = G4Exp(G4Log(_eMinTable)-_eStepLog)*G 102 if (_eMin < _lowLimit) 103 throw G4HadronicException(__FILE__, __LINE 104 G4PhysicsVector* np = new G4PhysicsLogVector 105 106 G4int i; 107 for (i=0; i<tableSize; i++) 108 { 109 G4double value = ppTable[i] * millibarn; 110 pp->PutValue(i,value); 111 value = npTable[i] * millibarn; 112 np->PutValue(i,value); 113 } 114 xMap[G4Proton::ProtonDefinition()] = pp; 115 xMap[G4Neutron::NeutronDefinition()] = np; 116 } 117 118 119 G4XNNElasticLowE::~G4XNNElasticLowE() 120 { 121 delete xMap[G4Proton::ProtonDefinition()]; 122 delete xMap[G4Neutron::NeutronDefinition()]; 123 } 124 125 126 G4bool G4XNNElasticLowE::operator==(const G4XN 127 { 128 return (this == (G4XNNElasticLowE *) &right) 129 } 130 131 132 G4bool G4XNNElasticLowE::operator!=(const G4XN 133 { 134 135 return (this != (G4XNNElasticLowE *) &right) 136 } 137 138 139 G4double G4XNNElasticLowE::CrossSection(const 140 { 141 G4double sigma = 0.; 142 G4double sqrtS = (trk1.Get4Momentum() + trk2 143 G4bool dummy = false; 144 145 const G4ParticleDefinition * key = FindKeyPa 146 147 typedef std::map <const G4ParticleDefinition 148 149 if (xMap.find(key)!= xMap.end()) 150 { 151 152 StringPhysMap::const_iterator iter; 153 for (iter = xMap.begin(); iter != xMap.e 154 { 155 const G4ParticleDefinition * str = (*iter) 156 if (str == key) 157 { 158 G4PhysicsVector* physVector = (*iter). 159 // G4PhysicsVector* physVector = x 160 if (sqrtS >= _eMin && sqrtS <= _eMax) 161 { 162 sigma = physVector->GetValue(sqrtS,dummy 163 } else if ( sqrtS < _eMin ) 164 { 165 sigma = physVector->GetValue 166 } 167 //G4cout << " sqrtS / sigma " << sqrtS/GeV 168 // sigma/millibarn << G4endl; 169 } 170 } 171 } 172 return sigma; 173 } 174 175 176 void G4XNNElasticLowE::Print() const 177 { 178 // Dump the pp cross-section table 179 180 G4cout << Name() << ", pp cross-section: " < 181 182 G4bool dummy = false; 183 G4int i; 184 const G4ParticleDefinition * key = G4Proton: 185 G4PhysicsVector* pp = 0; 186 187 typedef std::map <const G4ParticleDefinition 188 StringPhysMap::const_iterator iter; 189 190 for (iter = xMap.begin(); iter != xMap.end() 191 { 192 const G4ParticleDefinition * str = (*ite 193 if (str == key) 194 { 195 pp = (*iter).second; 196 } 197 } 198 199 if (pp != 0) 200 { 201 for (i=0; i<tableSize; i++) 202 { 203 G4double e = pp->GetLowEdgeEnergy(i); 204 G4double sigma = pp->GetValue(e,dummy) / m 205 G4cout << i << ") e = " << e / GeV << " Ge 206 } 207 } 208 209 // Dump the np cross-section table 210 211 G4cout << Name() << ", np cross-section: " < 212 213 key = G4Neutron::NeutronDefinition(); 214 G4PhysicsVector* np = 0; 215 for (iter = xMap.begin(); iter != xMap.end() 216 { 217 const G4ParticleDefinition * str = (*ite 218 if (str == key) 219 { 220 np = (*iter).second; 221 } 222 } 223 224 // G4PhysicsVector* np = xMap[G4Neutron::Ne 225 226 if (np != 0) 227 { 228 for (i=0; i<tableSize; i++) 229 { 230 G4double e = np->GetLowEdgeEnergy(i); 231 G4double sigma = np->GetValue(e,dummy) / m 232 G4cout << i << ") e = " << e / GeV << " Ge 233 } 234 } 235 G4VCrossSectionSource::Print(); 236 } 237 238 239 G4String G4XNNElasticLowE::Name() const 240 { 241 G4String name("NNElasticLowE"); 242 return name; 243 } 244 245 246 247 G4bool G4XNNElasticLowE::IsValid(G4double e) c 248 { 249 G4bool answer = InLimits(e,_lowLimit,_highLi 250 251 return answer; 252 } 253 254 255