G4eBremsstrahlungModel.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 // $Id$
00027 //
00028 // -------------------------------------------------------------------
00029 //
00030 // GEANT4 Class header file
00031 //
00032 //
00033 // File name:     G4eBremsstrahlungModel
00034 //
00035 // Author:        Vladimir Ivanchenko on base of Laszlo Urban code
00036 //
00037 // Creation date: 07.01.2002
00038 //
00039 // Modifications:
00040 //
00041 // 23-12-02 Change interface in order to move to cut per region (V.Ivanchenko)
00042 // 24-01-03 Make models region aware (V.Ivanchenko)
00043 // 13-02-03 Add name (V.Ivanchenko)
00044 // 08-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
00045 // 07-02-06  public function ComputeCrossSectionPerAtom() (mma)
00046 //
00047 //
00048 // Class Description:
00049 //
00050 // Implementation of energy loss for gamma emission by electrons and
00051 // positrons
00052 
00053 // -------------------------------------------------------------------
00054 //
00055 
00056 #ifndef G4eBremsstrahlungModel_h
00057 #define G4eBremsstrahlungModel_h 1
00058 
00059 #include "G4VEmModel.hh"
00060 
00061 class G4Element;
00062 class G4ParticleChangeForLoss;
00063 
00064 class G4eBremsstrahlungModel : public G4VEmModel
00065 {
00066 
00067 public:
00068 
00069   G4eBremsstrahlungModel(const G4ParticleDefinition* p = 0, 
00070                          const G4String& nam = "eBrem");
00071 
00072   virtual ~G4eBremsstrahlungModel();
00073 
00074   virtual void Initialise(const G4ParticleDefinition*, const G4DataVector&);
00075 
00076   virtual G4double ComputeDEDXPerVolume(const G4Material*,
00077                                         const G4ParticleDefinition*,
00078                                         G4double kineticEnergy,
00079                                         G4double cutEnergy);
00080                                         
00081   virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition*,
00082                                                   G4double tkin, 
00083                                                   G4double Z,   G4double,
00084                                                   G4double cut,
00085                                                   G4double maxE = DBL_MAX);
00086   
00087   virtual G4double CrossSectionPerVolume(const G4Material*,
00088                                          const G4ParticleDefinition*,
00089                                          G4double kineticEnergy,
00090                                          G4double cutEnergy,
00091                                          G4double maxEnergy);
00092 
00093   virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*,
00094                                  const G4MaterialCutsCouple*,
00095                                  const G4DynamicParticle*,
00096                                  G4double tmin,
00097                                  G4double maxEnergy);
00098 
00099 protected:
00100 
00101   const G4Element* SelectRandomAtom(const G4MaterialCutsCouple* couple);
00102 
00103 private:
00104 
00105   void SetParticle(const G4ParticleDefinition* p);
00106 
00107   G4double ComputeBremLoss(G4double Z, G4double tkin, G4double cut);
00108 
00109   G4double PositronCorrFactorLoss(G4double Z, G4double tkin, G4double cut);
00110 
00111   G4double PositronCorrFactorSigma(G4double Z, G4double tkin, G4double cut);
00112 
00113   G4DataVector* ComputePartialSumSigma(const G4Material* material,
00114                                              G4double tkin, G4double cut);
00115 
00116   G4double SupressionFunction(const G4Material* material, G4double tkin,
00117                                     G4double gammaEnergy);
00118 
00119   inline G4double ScreenFunction1(G4double ScreenVariable);
00120 
00121   inline G4double ScreenFunction2(G4double ScreenVariable);
00122 
00123   // hide assignment operator
00124   G4eBremsstrahlungModel & operator=(const  G4eBremsstrahlungModel &right);
00125   G4eBremsstrahlungModel(const  G4eBremsstrahlungModel&);
00126 
00127 protected:
00128 
00129   const G4ParticleDefinition* particle;
00130   G4ParticleDefinition*       theGamma;
00131   G4ParticleChangeForLoss*    fParticleChange;
00132 
00133   G4bool   isElectron;
00134 
00135 private:
00136 
00137   G4double highKinEnergy;
00138   G4double lowKinEnergy;
00139   G4double probsup;
00140   G4double MigdalConstant;
00141   G4double LPMconstant;
00142   G4bool   isInitialised;
00143 
00144   std::vector<G4DataVector*> partialSumSigma;
00145 
00146 };
00147 
00148 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00149 
00150 inline G4double G4eBremsstrahlungModel::ScreenFunction1(G4double ScreenVariable)
00151 
00152 // compute the value of the screening function 3*PHI1 - PHI2
00153 
00154 {
00155   G4double screenVal;
00156 
00157   if (ScreenVariable > 1.)
00158     screenVal = 42.24 - 8.368*std::log(ScreenVariable+0.952);
00159   else
00160     screenVal = 42.392 - ScreenVariable* (7.796 - 1.961*ScreenVariable);
00161 
00162   return screenVal;
00163 } 
00164 
00165 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00166 
00167 inline
00168 G4double G4eBremsstrahlungModel::ScreenFunction2(G4double ScreenVariable)
00169 
00170 // compute the value of the screening function 1.5*PHI1 - 0.5*PHI2
00171 
00172 {
00173   G4double screenVal;
00174 
00175   if (ScreenVariable > 1.)
00176     screenVal = 42.24 - 8.368*std::log(ScreenVariable+0.952);
00177   else
00178     screenVal = 41.734 - ScreenVariable* (6.484 - 1.250*ScreenVariable);
00179 
00180   return screenVal;
00181 } 
00182 
00183 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00184 
00185 #endif

Generated on Mon May 27 17:48:05 2013 for Geant4 by  doxygen 1.4.7