Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-18 09:59:19

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // -------------------------------------------------------------------
0027 //
0028 // GEANT4 Class header file
0029 //
0030 //
0031 // File name:    G4VCrossSectionDataSet
0032 //
0033 // Author  F.W. Jones, TRIUMF, 20-JAN-97
0034 //
0035 // Modifications:
0036 // 23.01.2009 V.Ivanchenko move constructor and destructor to source
0037 // 05.07.2010 V.Ivanchenko added name, min and max energy limit and
0038 //            corresponding access methods
0039 // 12.08.2011 G.Folger, V.Ivanchenko, T.Koi, D.Wright redesign the class
0040 // 
0041 //
0042 // Class Description
0043 // This is a base class for hadronic cross section data sets.  Users may 
0044 // derive specialized cross section classes and register them with the
0045 // appropriate process, or use provided data sets.
0046 //
0047 // Each cross section should have unique name
0048 // Minimal and maximal energy for the cross section will be used in run
0049 // time before IsApplicable method is called
0050 // 
0051 // Both the name and the energy interval will be used for documentation 
0052 //
0053 // Class Description - End
0054 
0055 #ifndef G4VCrossSectionDataSet_h
0056 #define G4VCrossSectionDataSet_h 1
0057 
0058 #include "globals.hh"
0059 #include "G4ParticleDefinition.hh"
0060 #include "G4DynamicParticle.hh"
0061 #include "G4Element.hh"
0062 #include <iostream>
0063 
0064 class G4DynamicParticle;
0065 class G4Isotope;
0066 class G4Material;
0067 class G4CrossSectionDataSetRegistry;
0068 
0069 class G4VCrossSectionDataSet
0070 {
0071 public: //with description
0072 
0073   G4VCrossSectionDataSet(const G4String& nam = "");
0074 
0075   virtual ~G4VCrossSectionDataSet();
0076 
0077   //============== Is Applicable methods ===============================
0078   // The following three methods have default implementations returning
0079   // "false".  Derived classes should implement only needed methods.
0080  
0081   // Element-wise cross section
0082   virtual
0083   G4bool IsElementApplicable(const G4DynamicParticle*, G4int Z, 
0084                  const G4Material* mat = nullptr);
0085 
0086   // Derived classes should implement this method if they provide isotope-wise 
0087   // cross sections.  Default arguments G4Element and G4Material are needed to 
0088   // access low-energy neutron cross sections, but are not required for others.
0089   virtual
0090   G4bool IsIsoApplicable(const G4DynamicParticle*, G4int Z, G4int A,    
0091              const G4Element* elm = nullptr,
0092              const G4Material* mat = nullptr);
0093 
0094   //============== GetCrossSection methods ===============================
0095 
0096   // This is a generic method to access cross section per element
0097   // This method should not be overwritten in a derived class
0098   inline G4double GetCrossSection(const G4DynamicParticle*, const G4Element*,
0099                                   const G4Material* mat = nullptr);
0100 
0101   // This is a generic method to compute cross section per element
0102   // If the DataSet is not applicable the method returns zero
0103   // This method should not be overwritten in a derived class
0104   G4double ComputeCrossSection(const G4DynamicParticle*, 
0105                                const G4Element*,
0106                                const G4Material* mat = nullptr);
0107 
0108   // Implement element cross section, IsApplicable does not checked.
0109   // In the default implementation a sum of isotope cross sections is computed
0110   virtual
0111   G4double ComputeCrossSectionPerElement(G4double kinEnergy, G4double loge,
0112                                          const G4ParticleDefinition*, 
0113                                          const G4Element*,
0114                                          const G4Material* mat = nullptr);
0115 
0116   // Implement these methods for element-wise cross section 
0117   virtual
0118   G4double GetElementCrossSection(const G4DynamicParticle*, G4int Z,
0119                   const G4Material* mat = nullptr);
0120 
0121   // Derived classes should implement these methods if they provide isotope-wise
0122   // cross sections. Extra arguments G4Isotope, G4Element, and G4Material are 
0123   // needed to access low-energy neutron cross sections, but not in other cases. 
0124   virtual
0125   G4double GetIsoCrossSection(const G4DynamicParticle*, G4int Z, G4int A,  
0126                   const G4Isotope* iso = nullptr,
0127                   const G4Element* elm = nullptr,
0128                   const G4Material* mat = nullptr);
0129 
0130   virtual
0131   G4double ComputeIsoCrossSection(G4double kinEnergy, G4double loge,
0132                                   const G4ParticleDefinition*, G4int Z, G4int A,  
0133                       const G4Isotope* iso = nullptr,
0134                       const G4Element* elm = nullptr,
0135                       const G4Material* mat = nullptr);
0136 
0137   //=====================================================================
0138 
0139   // Implement this method if needed
0140   // This method is called for element-wise cross section
0141   // Default implementation assumes equal cross sections for all isotopes 
0142   virtual const G4Isotope* SelectIsotope(const G4Element*, G4double kinEnergy, 
0143                                          G4double logE);
0144 
0145   // Implement this method if needed
0146   virtual
0147   void BuildPhysicsTable(const G4ParticleDefinition&);
0148 
0149   // Implement this method if needed
0150   // Default implementation will provide a dump of the cross section 
0151   // in logarithmic scale in the interval of applicability 
0152   virtual
0153   void DumpPhysicsTable(const G4ParticleDefinition&);
0154 
0155   virtual void CrossSectionDescription(std::ostream&) const;
0156 
0157   virtual void SetVerboseLevel(G4int value);
0158 
0159   inline G4double GetMinKinEnergy() const;
0160 
0161   inline void SetMinKinEnergy(G4double value);
0162 
0163   inline G4double GetMaxKinEnergy() const;
0164 
0165   inline void SetMaxKinEnergy(G4double value);
0166 
0167   inline bool ForAllAtomsAndEnergies() const;
0168 
0169   inline void SetForAllAtomsAndEnergies(G4bool val);
0170 
0171   inline const G4String& GetName() const;
0172 
0173   inline void SetName(const G4String& nam);
0174 
0175   G4VCrossSectionDataSet & operator=
0176   (const G4VCrossSectionDataSet &right) = delete;
0177   G4VCrossSectionDataSet(const G4VCrossSectionDataSet&) = delete;
0178 
0179 protected:
0180 
0181   G4int verboseLevel;
0182 
0183   G4String name;
0184 
0185 private:
0186 
0187   G4CrossSectionDataSetRegistry* registry;
0188 
0189   G4double minKinEnergy;
0190   G4double maxKinEnergy;
0191 
0192   G4bool isForAllAtomsAndEnergies;
0193 
0194 };
0195 
0196 inline G4double 
0197 G4VCrossSectionDataSet::GetCrossSection(const G4DynamicParticle* dp, 
0198                                         const G4Element* elm,
0199                                         const G4Material* mat)
0200 {
0201   return ComputeCrossSection(dp, elm, mat);
0202 }
0203 
0204 inline void G4VCrossSectionDataSet::SetVerboseLevel(G4int value)
0205 {
0206   verboseLevel = value;
0207 }
0208 
0209 inline void G4VCrossSectionDataSet::SetMinKinEnergy(G4double value)
0210 {
0211   minKinEnergy = value;
0212 }
0213 
0214 inline G4double G4VCrossSectionDataSet::GetMinKinEnergy() const
0215 {
0216   return minKinEnergy;
0217 }
0218 
0219 inline void G4VCrossSectionDataSet::SetMaxKinEnergy(G4double value)
0220 {
0221   maxKinEnergy = value;
0222 }
0223 
0224 inline G4double G4VCrossSectionDataSet::GetMaxKinEnergy() const
0225 {
0226   return maxKinEnergy;
0227 }
0228 
0229 inline const G4String& G4VCrossSectionDataSet::GetName() const
0230 {
0231   return name;
0232 }
0233 
0234 inline bool G4VCrossSectionDataSet::ForAllAtomsAndEnergies() const
0235 {
0236   return isForAllAtomsAndEnergies;
0237 }
0238 
0239 inline void G4VCrossSectionDataSet::SetForAllAtomsAndEnergies(G4bool val)
0240 {
0241   isForAllAtomsAndEnergies = val;
0242 }
0243 
0244 inline void G4VCrossSectionDataSet::SetName(const G4String& nam)
0245 {
0246   name = nam;
0247 }
0248 
0249 #endif