G4BraggModel.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:     G4BraggModel
00034 //
00035 // Author:        Vladimir Ivanchenko
00036 //
00037 // Creation date: 03.01.2002
00038 //
00039 // Modifications:
00040 // 23-12-02 V.Ivanchenko change interface in order to moveto cut per region
00041 // 24-01-03 Make models region aware (V.Ivanchenko)
00042 // 13-02-03 Add name (V.Ivanchenko)
00043 // 12-11-03 Fix for GenericIons (V.Ivanchenko)
00044 // 11-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
00045 // 15-02-06 ComputeCrossSectionPerElectron, ComputeCrossSectionPerAtom (mma)
00046 // 25-04-06 Added stopping data from PSTAR (V.Ivanchenko)
00047 // 12-08-08 Added methods GetParticleCharge, GetChargeSquareRatio, 
00048 //          CorrectionsAlongStep needed for ions(V.Ivanchenko)
00049 
00050 //
00051 // Class Description:
00052 //
00053 // Implementation of energy loss and delta-electron production
00054 // by heavy slow charged particles using ICRU'49 and NIST evaluated data 
00055 // for protons
00056 
00057 // -------------------------------------------------------------------
00058 //
00059 
00060 #ifndef G4BraggModel_h
00061 #define G4BraggModel_h 1
00062 
00063 #include <CLHEP/Units/PhysicalConstants.h>
00064 
00065 #include "G4VEmModel.hh"
00066 #include "G4PSTARStopping.hh"
00067 
00068 class G4ParticleChangeForLoss;
00069 class G4EmCorrections;
00070 
00071 class G4BraggModel : public G4VEmModel
00072 {
00073 
00074 public:
00075 
00076   G4BraggModel(const G4ParticleDefinition* p = 0,
00077                const G4String& nam = "Bragg");
00078 
00079   virtual ~G4BraggModel();
00080 
00081   virtual void Initialise(const G4ParticleDefinition*, const G4DataVector&);
00082 
00083   virtual G4double ComputeCrossSectionPerElectron(
00084                                  const G4ParticleDefinition*,
00085                                  G4double kineticEnergy,
00086                                  G4double cutEnergy,
00087                                  G4double maxEnergy);
00088                                  
00089   virtual G4double ComputeCrossSectionPerAtom(
00090                                  const G4ParticleDefinition*,
00091                                  G4double kineticEnergy,
00092                                  G4double Z, G4double A,
00093                                  G4double cutEnergy,
00094                                  G4double maxEnergy);
00095                                                                  
00096   virtual G4double CrossSectionPerVolume(const G4Material*,
00097                                  const G4ParticleDefinition*,
00098                                  G4double kineticEnergy,
00099                                  G4double cutEnergy,
00100                                  G4double maxEnergy);
00101                                  
00102   virtual G4double ComputeDEDXPerVolume(const G4Material*,
00103                                 const G4ParticleDefinition*,
00104                                 G4double kineticEnergy,
00105                                 G4double cutEnergy);
00106 
00107   virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*,
00108                                  const G4MaterialCutsCouple*,
00109                                  const G4DynamicParticle*,
00110                                  G4double tmin,
00111                                  G4double maxEnergy);
00112 
00113   // Compute ion charge 
00114   virtual G4double GetChargeSquareRatio(const G4ParticleDefinition*,
00115                                         const G4Material*,
00116                                         G4double kineticEnergy);
00117 
00118   virtual G4double GetParticleCharge(const G4ParticleDefinition* p,
00119                                      const G4Material* mat,
00120                                      G4double kineticEnergy);
00121 
00122 protected:
00123 
00124   virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition*,
00125                                       G4double kinEnergy);
00126 
00127   inline G4double GetChargeSquareRatio() const;
00128 
00129   inline void SetChargeSquareRatio(G4double val);
00130 
00131 private:
00132 
00133   inline void SetParticle(const G4ParticleDefinition* p);
00134 
00135   G4bool HasMaterial(const G4Material* material);
00136 
00137   G4double StoppingPower(const G4Material* material,
00138                                G4double kineticEnergy);
00139 
00140   G4double ElectronicStoppingPower(G4double z,
00141                                    G4double kineticEnergy) const;
00142 
00143   G4double DEDX(const G4Material* material, G4double kineticEnergy);
00144 
00145   G4bool MolecIsInZiegler1988(const G4Material* material);
00146 
00147   G4double ChemicalFactor(G4double kineticEnergy, G4double eloss125) const;
00148 
00149   void SetExpStopPower125(G4double value) {expStopPower125 = value;};
00150 
00151   // hide assignment operator
00152   G4BraggModel & operator=(const  G4BraggModel &right);
00153   G4BraggModel(const  G4BraggModel&);
00154 
00155 
00156   G4EmCorrections*            corr;
00157 
00158   const G4ParticleDefinition* particle;
00159   G4ParticleDefinition*       theElectron;
00160   G4ParticleChangeForLoss*    fParticleChange;
00161   G4PSTARStopping             pstar;
00162 
00163   const G4Material*           currentMaterial;
00164 
00165   G4double mass;
00166   G4double spin;
00167   G4double chargeSquare;
00168   G4double massRate;
00169   G4double ratio;
00170   G4double lowestKinEnergy;
00171   G4double protonMassAMU;
00172   G4double theZieglerFactor;
00173   G4double expStopPower125;    // Experimental Stopping power at 125keV
00174 
00175   G4int    iMolecula;          // index in the molecula's table
00176   G4int    iPSTAR;             // index in PSTAR
00177   G4bool   isIon;
00178   G4bool   isInitialised;
00179 };
00180 
00181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
00182 
00183 inline void G4BraggModel::SetParticle(const G4ParticleDefinition* p)
00184 {
00185   particle = p;
00186   mass = particle->GetPDGMass();
00187   spin = particle->GetPDGSpin();
00188   G4double q = particle->GetPDGCharge()/CLHEP::eplus;
00189   chargeSquare = q*q;
00190   massRate     = mass/CLHEP::proton_mass_c2;
00191   ratio = CLHEP::electron_mass_c2/mass;
00192 }
00193 
00194 inline G4double G4BraggModel::GetChargeSquareRatio() const
00195 {
00196   return chargeSquare;
00197 }
00198 
00199 inline void G4BraggModel::SetChargeSquareRatio(G4double val)
00200 {
00201   chargeSquare = val;
00202 }
00203 
00204 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
00205 
00206 #endif

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