G4VCrossSectionDataSet.cc

Go to the documentation of this file.
00001 //
00002 // ********************************************************************
00003 // * License and Disclaimer                                           *
00004 // *                                                                  *
00005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
00006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
00007 // * conditions of the Geant4 Software License,  included in the file *
00008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
00009 // * include a list of copyright holders.                             *
00010 // *                                                                  *
00011 // * Neither the authors of this software system, nor their employing *
00012 // * institutes,nor the agencies providing financial support for this *
00013 // * work  make  any representation or  warranty, express or implied, *
00014 // * regarding  this  software system or assume any liability for its *
00015 // * use.  Please see the license in the file  LICENSE  and URL above *
00016 // * for the full disclaimer and the limitation of liability.         *
00017 // *                                                                  *
00018 // * This  code  implementation is the result of  the  scientific and *
00019 // * technical work of the GEANT4 collaboration.                      *
00020 // * By using,  copying,  modifying or  distributing the software (or *
00021 // * any work based  on the software)  you  agree  to acknowledge its *
00022 // * use  in  resulting  scientific  publications,  and indicate your *
00023 // * acceptance of all terms of the Geant4 Software license.          *
00024 // ********************************************************************
00025 //
00026 // $Id$
00027 //
00028 // -------------------------------------------------------------------
00029 //
00030 // GEANT4 Class file
00031 //
00032 //
00033 // File name:    G4VCrossSectionDataSet
00034 //
00035 // Author  F.W. Jones, TRIUMF, 20-JAN-97
00036 //
00037 // Modifications:
00038 // 23.01.2009 V.Ivanchenko move constructor and destructor to source
00039 // 12.08.2011 G.Folger, V.Ivanchenko, T.Koi, D.Wright redesign the class
00040 //
00041 
00042 #include "G4VCrossSectionDataSet.hh"
00043 #include "G4SystemOfUnits.hh"
00044 #include "G4CrossSectionDataSetRegistry.hh"
00045 #include "G4DynamicParticle.hh"
00046 #include "G4Material.hh"
00047 #include "G4Element.hh"
00048 #include "G4Isotope.hh"
00049 #include "G4NistManager.hh"
00050 #include "G4HadronicException.hh"
00051 #include "G4HadTmpUtil.hh"
00052 #include "Randomize.hh"
00053 
00054 G4VCrossSectionDataSet::G4VCrossSectionDataSet(const G4String& nam) :
00055   verboseLevel(0),minKinEnergy(0.0),maxKinEnergy(100*TeV),name(nam) 
00056 {
00057   G4CrossSectionDataSetRegistry::Instance()->Register(this);
00058 }
00059 
00060 G4VCrossSectionDataSet::~G4VCrossSectionDataSet()
00061 {
00062   G4CrossSectionDataSetRegistry::Instance()->DeRegister(this);
00063 }
00064 
00065 G4bool 
00066 G4VCrossSectionDataSet::IsElementApplicable(const G4DynamicParticle*, 
00067                                             G4int,
00068                                             const G4Material*)
00069 {
00070   return false;
00071 }
00072 
00073 G4bool 
00074 G4VCrossSectionDataSet::IsIsoApplicable(const G4DynamicParticle*, 
00075                                         G4int, G4int,
00076                                         const G4Element*,  
00077                                         const G4Material*)
00078 {
00079   return false;
00080 }
00081 
00082 G4double 
00083 G4VCrossSectionDataSet::ComputeCrossSection(const G4DynamicParticle* part, 
00084                                             const G4Element* elm,
00085                                             const G4Material* mat)
00086 {
00087   G4int Z = G4lrint(elm->GetZ());
00088 
00089   if (IsElementApplicable(part, Z, mat)) { 
00090     return GetElementCrossSection(part, Z, mat);
00091   }
00092 
00093   // isotope-wise cross section making sum over available
00094   // isotope cross sections, which may be incomplete, so
00095   // the result is corrected 
00096   G4int nIso = elm->GetNumberOfIsotopes();    
00097   G4double fact = 0.0;
00098   G4double xsec = 0.0;
00099   G4Isotope* iso = 0;
00100 
00101   if (0 < nIso) { 
00102 
00103     // user-defined isotope abundances        
00104     G4IsotopeVector* isoVector = elm->GetIsotopeVector();
00105     G4double* abundVector = elm->GetRelativeAbundanceVector();
00106 
00107     for (G4int j = 0; j<nIso; ++j) {
00108       iso = (*isoVector)[j];
00109       G4int A = iso->GetN();
00110       if(abundVector[j] > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
00111         fact += abundVector[j];
00112         xsec += abundVector[j]*GetIsoCrossSection(part, Z, A, iso, elm, mat);
00113       }
00114     }
00115 
00116   } else {
00117 
00118     // natural isotope abundances
00119     G4NistManager* nist = G4NistManager::Instance();
00120     G4int n0 = nist->GetNistFirstIsotopeN(Z);
00121     G4int nn = nist->GetNumberOfNistIsotopes(Z);
00122     for (G4int A = n0; A < n0+nn; ++A) {
00123       G4double abund = nist->GetIsotopeAbundance(Z, A);
00124       if(abund > 0.0 && IsIsoApplicable(part, Z, A, elm, mat)) {
00125         fact += abund;
00126         xsec += abund*GetIsoCrossSection(part, Z, A, iso, elm, mat);
00127       }
00128     }
00129   }
00130   if(fact > 0.0) { xsec /= fact; }
00131   return xsec;
00132 }
00133 
00134 G4double 
00135 G4VCrossSectionDataSet::GetElementCrossSection(const G4DynamicParticle* dynPart,
00136                                                G4int Z,
00137                                                const G4Material* mat)
00138 {
00139   G4cout << "G4VCrossSectionDataSet::GetCrossSection per element ERROR: "
00140          << " there is no cross section for "
00141          << dynPart->GetDefinition()->GetParticleName()
00142          << "  E(MeV)= "  << dynPart->GetKineticEnergy()/MeV;
00143   if(mat) { G4cout << "  inside " << mat->GetName(); }
00144   G4cout << " for Z= " << Z << G4endl;
00145   throw G4HadronicException(__FILE__, __LINE__,
00146         "G4VCrossSectionDataSet::GetElementCrossSection is absent");
00147   return 0.0;
00148 }
00149 
00150 G4double 
00151 G4VCrossSectionDataSet::GetIsoCrossSection(const G4DynamicParticle* dynPart,
00152                                            G4int Z, G4int A,
00153                                            const G4Isotope*,
00154                                            const G4Element* elm,
00155                                            const G4Material* mat)
00156 {
00157   G4cout << "G4VCrossSectionDataSet::GetCrossSection per isotope ERROR: "
00158          << " there is no cross section for "
00159          << dynPart->GetDefinition()->GetParticleName()
00160          << "  E(MeV)= "  << dynPart->GetKineticEnergy()/MeV;
00161   if(mat) { G4cout << "  inside " << mat->GetName(); }
00162   if(elm) { G4cout << " for " << elm->GetName(); }
00163   G4cout << "  Z= " << Z << " A= " << A << G4endl;
00164   throw G4HadronicException(__FILE__, __LINE__,
00165         "G4VCrossSectionDataSet::GetIsoCrossSection is absent");
00166   return 0.0;
00167 }
00168 
00169 G4Isotope* 
00170 G4VCrossSectionDataSet::SelectIsotope(const G4Element* anElement, G4double)
00171 {
00172   G4int nIso = anElement->GetNumberOfIsotopes();
00173   G4IsotopeVector* isoVector = anElement->GetIsotopeVector();
00174   G4Isotope* iso = (*isoVector)[0];
00175 
00176   // more than 1 isotope
00177   if(1 < nIso) {
00178     G4double* abundVector = anElement->GetRelativeAbundanceVector();
00179     G4double sum = 0.0;
00180     G4double q = G4UniformRand();
00181     for (G4int j = 0; j<nIso; ++j) {
00182       sum += abundVector[j];
00183       if(q <= sum) {
00184         iso = (*isoVector)[j];
00185         break;
00186       }
00187     }
00188   }
00189   return iso;
00190 }
00191 
00192 void G4VCrossSectionDataSet::BuildPhysicsTable(const G4ParticleDefinition&)
00193 {}
00194 
00195 void G4VCrossSectionDataSet::DumpPhysicsTable(const G4ParticleDefinition&)
00196 {}
00197 
00198 void G4VCrossSectionDataSet::CrossSectionDescription(std::ostream& outFile) const
00199 {
00200   outFile << "The description for this cross section data set has not been written yet.\n";
00201 }

Generated on Mon May 27 17:50:10 2013 for Geant4 by  doxygen 1.4.7