G4ComponentGGNuclNuclXsc.hh

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 // Calculation of the nucleus-nucleus total, inelastic, production, 
00027 // elastic and quasi-elastic  cross-sections
00028 // based on parametrisations of nucleon-nucleon
00029 // cross-sections  in 
00030 // the framework of simplified Glauber-Gribov approach
00031 //
00032 //
00033 // 24.11.08 V. Grichine - first implementation based on G4GlauberGribovCrossSection
00034 //
00035 //
00036 
00037 #ifndef G4ComponentGGNuclNuclXsc_h
00038 #define G4ComponentGGNuclNuclXsc_h
00039 
00040 #include "globals.hh"
00041 #include "G4Proton.hh"
00042 #include "G4Nucleus.hh"
00043 #include "G4NistManager.hh"
00044 
00045 #include "G4VComponentCrossSection.hh"
00046 
00047 class G4ParticleDefinition;
00048 class G4HadronNucleonXsc;
00049 
00050 class G4ComponentGGNuclNuclXsc : public G4VComponentCrossSection
00051 {
00052 public:
00053 
00054   G4ComponentGGNuclNuclXsc ();
00055   virtual ~G4ComponentGGNuclNuclXsc ();
00056 
00057 
00058   // virtual interface methods
00059 
00060   virtual
00061   G4double GetTotalIsotopeCrossSection(const G4ParticleDefinition* aParticle,
00062                                        G4double kinEnergy,
00063                                        G4int Z, G4int A);
00064 
00065   virtual
00066   G4double GetTotalElementCrossSection(const G4ParticleDefinition* aParticle,
00067                                        G4double kinEnergy, 
00068                                        G4int Z, G4double A);
00069 
00070   virtual
00071   G4double GetInelasticIsotopeCrossSection(const G4ParticleDefinition* aParticle,
00072                                            G4double kinEnergy, 
00073                                            G4int Z, G4int A);
00074 
00075   virtual
00076   G4double GetInelasticElementCrossSection(const G4ParticleDefinition* aParticle,
00077                                            G4double kinEnergy, 
00078                                            G4int Z, G4double A);
00079 
00080   virtual
00081   G4double GetElasticElementCrossSection(const G4ParticleDefinition* aParticle,
00082                                          G4double kinEnergy, 
00083                                          G4int Z, G4double A);
00084 
00085   virtual
00086   G4double GetElasticIsotopeCrossSection(const G4ParticleDefinition* aParticle,
00087                                          G4double kinEnergy, 
00088                                          G4int Z, G4int A);
00089  
00090   virtual
00091   G4double ComputeQuasiElasticRatio(const G4ParticleDefinition* aParticle,
00092                                          G4double kinEnergy, 
00093                                          G4int Z, G4int A);
00094  
00095    
00096   //  virtual
00097   G4bool IsElementApplicable(const G4DynamicParticle*, 
00098                              G4int Z, const G4Material*);
00099 
00100   //  virtual
00101   G4double GetElementCrossSection(const G4DynamicParticle*, 
00102                                   G4int Z, const G4Material*);
00103 
00104   G4double GetZandACrossSection(const G4DynamicParticle*, 
00105                                 G4int Z, G4int A);
00106 
00107   G4double GetCoulombBarier(const G4DynamicParticle*, 
00108                             G4double Z, G4double A, G4double pR, G4double tR);
00109 
00110   virtual
00111   void BuildPhysicsTable(const G4ParticleDefinition&)
00112   {}
00113 
00114   virtual
00115   void DumpPhysicsTable(const G4ParticleDefinition&) 
00116   {G4cout << "G4NuclNuclCrossSection: uses Glauber-Gribov formula"<<G4endl;}
00117 
00118   virtual void CrossSectionDescription(std::ostream&) const;
00119 
00120   G4double GetRatioSD(const G4DynamicParticle*, G4double At, G4double Zt);
00121   G4double GetRatioQE(const G4DynamicParticle*, G4double At, G4double Zt);
00122 
00123   G4double GetHadronNucleonXsc(const G4DynamicParticle*, const G4Element*);
00124   G4double GetHadronNucleonXsc(const G4DynamicParticle*, G4int At, G4int Zt);
00125 
00126  
00127   G4double GetHadronNucleonXscPDG(G4ParticleDefinition*,G4double sMand, G4ParticleDefinition*);
00128   G4double GetHadronNucleonXscNS(G4ParticleDefinition*,G4double pTkin, G4ParticleDefinition*);
00129 
00130   G4double GetHNinelasticXscVU(const G4DynamicParticle*, G4int At, G4int Zt);
00131   G4double CalculateEcmValue(const G4double, const G4double, const G4double); 
00132   G4double CalcMandelstamS( const G4double , const G4double , const G4double );
00133 
00134   G4double GetElasticGlauberGribov(const G4DynamicParticle*,G4int Z, G4int A);
00135   G4double GetInelasticGlauberGribov(const G4DynamicParticle*,G4int Z, G4int A);
00136 
00137   G4double GetTotalGlauberGribovXsc()    { return fTotalXsc;     }; 
00138   G4double GetElasticGlauberGribovXsc()  { return fElasticXsc;   }; 
00139   G4double GetInelasticGlauberGribovXsc(){ return fInelasticXsc; }; 
00140   G4double GetProductionGlauberGribovXsc(){ return fProductionXsc; }; 
00141   G4double GetDiffractionGlauberGribovXsc(){ return fDiffractionXsc; }; 
00142   G4double GetRadiusConst()              { return fRadiusConst;  }; 
00143 
00144   G4double GetNucleusRadius(const G4DynamicParticle*, const G4Element*);
00145 
00146   G4double GetNucleusRadius(G4double Zt, G4double At);
00147   G4double GetNucleusRadiusGG(G4double At);
00148   G4double GetNucleusRadiusDE(G4double Zt, G4double At);
00149   G4double GetNucleusRadiusRMS(G4double Zt, G4double At);
00150 
00151   inline void SetEnergyLowerLimit(G4double E ){fLowerLimit=E;};
00152 
00153 private:
00154 
00155   const G4double fUpperLimit;
00156   G4double fLowerLimit; 
00157   const G4double fRadiusConst;
00158  
00159   G4double fTotalXsc, fElasticXsc, fInelasticXsc, fProductionXsc, fDiffractionXsc;
00160   G4double fHadronNucleonXsc;
00161  
00162   G4ParticleDefinition* theProton;
00163   G4ParticleDefinition* theNeutron;
00164   G4HadronNucleonXsc* hnXsc;
00165 };
00166 
00168 //
00169 // Inlines
00170 
00171 inline G4double
00172 G4ComponentGGNuclNuclXsc::GetElasticGlauberGribov(const G4DynamicParticle* dp,
00173                                                   G4int Z, G4int A)
00174 {
00175   GetZandACrossSection(dp, Z, A);
00176   return fElasticXsc;
00177 }
00178 
00180 
00181 inline G4double
00182 G4ComponentGGNuclNuclXsc::GetInelasticGlauberGribov(const G4DynamicParticle* dp,
00183                                                     G4int Z, G4int A)
00184 {
00185   GetZandACrossSection(dp, Z, A);
00186   return fInelasticXsc;
00187 }
00188 
00189 #endif

Generated on Mon May 27 17:47:56 2013 for Geant4 by  doxygen 1.4.7