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 // GEANT4 Class file 28 // GEANT4 Class file 29 // 29 // 30 // File name: G4VCrossSectionDataSet 30 // File name: G4VCrossSectionDataSet 31 // 31 // 32 // Author F.W. Jones, TRIUMF, 20-JAN-97 32 // Author F.W. Jones, TRIUMF, 20-JAN-97 33 // 33 // 34 // Modifications: 34 // Modifications: 35 // 23.01.2009 V.Ivanchenko move constructor an 35 // 23.01.2009 V.Ivanchenko move constructor and destructor to source 36 // 12.08.2011 G.Folger, V.Ivanchenko, T.Koi, D 36 // 12.08.2011 G.Folger, V.Ivanchenko, T.Koi, D.Wright redesign the class 37 // 37 // 38 38 39 #include "G4VCrossSectionDataSet.hh" 39 #include "G4VCrossSectionDataSet.hh" 40 #include "G4SystemOfUnits.hh" 40 #include "G4SystemOfUnits.hh" 41 #include "G4CrossSectionDataSetRegistry.hh" 41 #include "G4CrossSectionDataSetRegistry.hh" 42 #include "G4Material.hh" 42 #include "G4Material.hh" 43 #include "G4Element.hh" 43 #include "G4Element.hh" 44 #include "G4Isotope.hh" 44 #include "G4Isotope.hh" 45 #include "G4NistManager.hh" 45 #include "G4NistManager.hh" 46 #include "Randomize.hh" 46 #include "Randomize.hh" 47 #include "G4HadronicParameters.hh" 47 #include "G4HadronicParameters.hh" 48 48 49 G4VCrossSectionDataSet::G4VCrossSectionDataSet 49 G4VCrossSectionDataSet::G4VCrossSectionDataSet(const G4String& nam) : 50 name(nam), << 50 verboseLevel(0),name(nam),minKinEnergy(0.0), 51 maxKinEnergy(G4HadronicParameters::Instance( << 51 maxKinEnergy(G4HadronicParameters::Instance()->GetMaxEnergy()), >> 52 isForAllAtomsAndEnergies(false) 52 { 53 { 53 registry = G4CrossSectionDataSetRegistry::In 54 registry = G4CrossSectionDataSetRegistry::Instance(); 54 registry->Register(this); 55 registry->Register(this); 55 } 56 } 56 57 57 G4VCrossSectionDataSet::~G4VCrossSectionDataSe 58 G4VCrossSectionDataSet::~G4VCrossSectionDataSet() 58 { 59 { 59 registry->DeRegister(this); 60 registry->DeRegister(this); 60 } 61 } 61 62 62 G4bool 63 G4bool 63 G4VCrossSectionDataSet::IsElementApplicable(co 64 G4VCrossSectionDataSet::IsElementApplicable(const G4DynamicParticle*, 64 G4int, 65 G4int, 65 const G4Material*) 66 const G4Material*) 66 { 67 { 67 return false; 68 return false; 68 } 69 } 69 70 70 G4bool 71 G4bool 71 G4VCrossSectionDataSet::IsIsoApplicable(const 72 G4VCrossSectionDataSet::IsIsoApplicable(const G4DynamicParticle*, 72 G4int, 73 G4int, G4int, 73 const 74 const G4Element*, 74 const 75 const G4Material*) 75 { 76 { 76 return false; 77 return false; 77 } 78 } 78 79 79 G4double 80 G4double 80 G4VCrossSectionDataSet::ComputeCrossSection(co 81 G4VCrossSectionDataSet::ComputeCrossSection(const G4DynamicParticle* part, 81 co 82 const G4Element* elm, 82 co 83 const G4Material* mat) 83 { 84 { 84 G4int Z = elm->GetZasInt(); 85 G4int Z = elm->GetZasInt(); 85 86 86 if (IsElementApplicable(part, Z, mat)) { 87 if (IsElementApplicable(part, Z, mat)) { 87 return GetElementCrossSection(part, Z, mat 88 return GetElementCrossSection(part, Z, mat); 88 } 89 } 89 90 90 // isotope-wise cross section making sum ove 91 // isotope-wise cross section making sum over available 91 // isotope cross sections, which may be inco 92 // isotope cross sections, which may be incomplete, so 92 // the result is corrected 93 // the result is corrected 93 std::size_t nIso = elm->GetNumberOfIsotopes( 94 std::size_t nIso = elm->GetNumberOfIsotopes(); 94 G4double fact = 0.0; 95 G4double fact = 0.0; 95 G4double xsec = 0.0; 96 G4double xsec = 0.0; 96 97 97 // user-defined isotope abundances 98 // user-defined isotope abundances 98 const G4IsotopeVector* isoVector = elm->GetI 99 const G4IsotopeVector* isoVector = elm->GetIsotopeVector(); 99 const G4double* abundVector = elm->GetRelati 100 const G4double* abundVector = elm->GetRelativeAbundanceVector(); 100 for (std::size_t j=0; j<nIso; ++j) { 101 for (std::size_t j=0; j<nIso; ++j) { 101 const G4Isotope* iso = (*isoVector)[j]; 102 const G4Isotope* iso = (*isoVector)[j]; 102 G4int A = iso->GetN(); 103 G4int A = iso->GetN(); 103 if(abundVector[j] > 0.0 && IsIsoApplicable 104 if(abundVector[j] > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) { 104 fact += abundVector[j]; 105 fact += abundVector[j]; 105 xsec += abundVector[j]*GetIsoCrossSectio 106 xsec += abundVector[j]*GetIsoCrossSection(part, Z, A, iso, elm, mat); 106 } 107 } 107 } 108 } 108 return (fact > 0.0) ? xsec/fact : 0.0; 109 return (fact > 0.0) ? xsec/fact : 0.0; 109 } 110 } 110 111 111 G4double 112 G4double 112 G4VCrossSectionDataSet::ComputeCrossSectionPer 113 G4VCrossSectionDataSet::ComputeCrossSectionPerElement( 113 G4double kinEnergy, G4 114 G4double kinEnergy, G4double loge, 114 const G4ParticleDefini 115 const G4ParticleDefinition* pd, 115 const G4Element* elm, 116 const G4Element* elm, const G4Material* mat) 116 { 117 { 117 G4int Z = elm->GetZasInt(); 118 G4int Z = elm->GetZasInt(); 118 std::size_t nIso = elm->GetNumberOfIsotopes( 119 std::size_t nIso = elm->GetNumberOfIsotopes(); 119 G4double xsec = 0.0; 120 G4double xsec = 0.0; 120 const G4IsotopeVector* isoVector = elm->GetI 121 const G4IsotopeVector* isoVector = elm->GetIsotopeVector(); 121 const G4double* abundVector = elm->GetRelati 122 const G4double* abundVector = elm->GetRelativeAbundanceVector(); 122 for (std::size_t j=0; j<nIso; ++j) { 123 for (std::size_t j=0; j<nIso; ++j) { 123 const G4Isotope* iso = (*isoVector)[j]; 124 const G4Isotope* iso = (*isoVector)[j]; 124 G4int A = iso->GetN(); 125 G4int A = iso->GetN(); 125 xsec += abundVector[j]* 126 xsec += abundVector[j]* 126 ComputeIsoCrossSection(kinEnergy, loge, 127 ComputeIsoCrossSection(kinEnergy, loge, pd, Z, A, iso, elm, mat); 127 } 128 } 128 return xsec; 129 return xsec; 129 } 130 } 130 131 131 G4double 132 G4double 132 G4VCrossSectionDataSet::GetElementCrossSection 133 G4VCrossSectionDataSet::GetElementCrossSection(const G4DynamicParticle* dynPart, 133 G4int Z, 134 G4int Z, 134 const G4Material* mat) 135 const G4Material* mat) 135 { 136 { 136 G4ExceptionDescription ed; 137 G4ExceptionDescription ed; 137 ed << "GetElementCrossSection is not impleme 138 ed << "GetElementCrossSection is not implemented in <" << name << ">\n" 138 << "Particle: " << dynPart->GetDefinition 139 << "Particle: " << dynPart->GetDefinition()->GetParticleName() 139 << " Ekin(MeV)= " << dynPart->GetKineti 140 << " Ekin(MeV)= " << dynPart->GetKineticEnergy()/MeV; 140 if(nullptr != mat) { ed << " material: " << 141 if(nullptr != mat) { ed << " material: " << mat->GetName(); } 141 ed << " target Z= " << Z << G4endl; 142 ed << " target Z= " << Z << G4endl; 142 G4Exception("G4VCrossSectionDataSet::GetElem 143 G4Exception("G4VCrossSectionDataSet::GetElementCrossSection", "had001", 143 FatalException, ed); 144 FatalException, ed); 144 return 0.0; 145 return 0.0; 145 } 146 } 146 147 147 G4double 148 G4double 148 G4VCrossSectionDataSet::GetIsoCrossSection(con 149 G4VCrossSectionDataSet::GetIsoCrossSection(const G4DynamicParticle* dynPart, 149 G4int Z, G4int A, 150 G4int Z, G4int A, 150 const G4Isotope*, 151 const G4Isotope*, 151 const G4Element* elm, 152 const G4Element* elm, 152 const G4Material* mat) 153 const G4Material* mat) 153 { 154 { 154 G4ExceptionDescription ed; 155 G4ExceptionDescription ed; 155 ed << "GetIsoCrossSection is not implemented 156 ed << "GetIsoCrossSection is not implemented in <" << name << ">\n" 156 << "Particle: " << dynPart->GetDefinition 157 << "Particle: " << dynPart->GetDefinition()->GetParticleName() 157 << " Ekin(MeV)= " << dynPart->GetKineti 158 << " Ekin(MeV)= " << dynPart->GetKineticEnergy()/MeV; 158 if(nullptr != mat) { ed << " material: " << 159 if(nullptr != mat) { ed << " material: " << mat->GetName(); } 159 if(nullptr != elm) { ed << " element: " << e 160 if(nullptr != elm) { ed << " element: " << elm->GetName(); } 160 ed << " target Z= " << Z << " A= " << A << G 161 ed << " target Z= " << Z << " A= " << A << G4endl; 161 G4Exception("G4VCrossSectionDataSet::GetIsoC 162 G4Exception("G4VCrossSectionDataSet::GetIsoCrossSection", "had001", 162 FatalException, ed); 163 FatalException, ed); 163 return 0.0; 164 return 0.0; 164 } 165 } 165 166 166 G4double 167 G4double 167 G4VCrossSectionDataSet::ComputeIsoCrossSection 168 G4VCrossSectionDataSet::ComputeIsoCrossSection(G4double kinEnergy, G4double, 168 const G4ParticleDefinition* p 169 const G4ParticleDefinition* pd, 169 G4int Z, G4int A, 170 G4int Z, G4int A, 170 const G4Isotope*, 171 const G4Isotope*, 171 const G4Element* elm, 172 const G4Element* elm, 172 const G4Material* mat) 173 const G4Material* mat) 173 { 174 { 174 G4ExceptionDescription ed; 175 G4ExceptionDescription ed; 175 ed << "GetIsoCrossSection is not implemented 176 ed << "GetIsoCrossSection is not implemented in <" << name << ">\n" 176 << "Particle: " << pd->GetParticleName() 177 << "Particle: " << pd->GetParticleName() 177 << " Ekin(MeV)= " << kinEnergy/CLHEP::M 178 << " Ekin(MeV)= " << kinEnergy/CLHEP::MeV; 178 if(nullptr != mat) { ed << " material: " << 179 if(nullptr != mat) { ed << " material: " << mat->GetName(); } 179 if(nullptr != elm) { ed << " element: " << e 180 if(nullptr != elm) { ed << " element: " << elm->GetName(); } 180 ed << " target Z= " << Z << " A= " << A << G 181 ed << " target Z= " << Z << " A= " << A << G4endl; 181 G4Exception("G4VCrossSectionDataSet::GetIsoC 182 G4Exception("G4VCrossSectionDataSet::GetIsoCrossSection", "had001", 182 FatalException, ed); 183 FatalException, ed); 183 return 0.0; 184 return 0.0; 184 } 185 } 185 186 186 const G4Isotope* 187 const G4Isotope* 187 G4VCrossSectionDataSet::SelectIsotope(const G4 188 G4VCrossSectionDataSet::SelectIsotope(const G4Element* anElement, 188 G4double 189 G4double, G4double) 189 { 190 { 190 G4int nIso = (G4int)anElement->GetNumberOfIs 191 G4int nIso = (G4int)anElement->GetNumberOfIsotopes(); 191 const G4Isotope* iso = anElement->GetIsotope 192 const G4Isotope* iso = anElement->GetIsotope(0); 192 193 193 // more than 1 isotope 194 // more than 1 isotope 194 if(1 < nIso) { 195 if(1 < nIso) { 195 const G4double* abundVector = anElement->G 196 const G4double* abundVector = anElement->GetRelativeAbundanceVector(); 196 G4double sum = 0.0; 197 G4double sum = 0.0; 197 G4double q = G4UniformRand(); 198 G4double q = G4UniformRand(); 198 for (G4int j=0; j<nIso; ++j) { 199 for (G4int j=0; j<nIso; ++j) { 199 sum += abundVector[j]; 200 sum += abundVector[j]; 200 if(q <= sum) { 201 if(q <= sum) { 201 iso = anElement->GetIsotope(j); 202 iso = anElement->GetIsotope(j); 202 break; 203 break; 203 } 204 } 204 } 205 } 205 } 206 } 206 return iso; 207 return iso; 207 } 208 } 208 209 209 void G4VCrossSectionDataSet::BuildPhysicsTable 210 void G4VCrossSectionDataSet::BuildPhysicsTable(const G4ParticleDefinition&) 210 {} 211 {} 211 212 212 void G4VCrossSectionDataSet::DumpPhysicsTable( 213 void G4VCrossSectionDataSet::DumpPhysicsTable(const G4ParticleDefinition&) 213 {} 214 {} 214 215 215 void G4VCrossSectionDataSet::CrossSectionDescr 216 void G4VCrossSectionDataSet::CrossSectionDescription(std::ostream& outFile) const 216 { 217 { 217 outFile << "The description for this cross s 218 outFile << "The description for this cross section data set has not been written yet.\n"; 218 } 219 } 219 220