G4hNuclearStoppingModel.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 //
00027 // -------------------------------------------------------------------
00028 //
00029 // GEANT4 Class file
00030 //
00031 //
00032 // File name:     G4hNuclearStoppingModel
00033 //
00034 // Author:        V.Ivanchenko (Vladimir.Ivanchenko@cern.ch)
00035 // 
00036 // Creation date: 20 July 2000
00037 //
00038 // Modifications:
00039 // 20/07/2000  V.Ivanchenko First implementation
00040 // 22/08/2000  V.Ivanchenko Bug fixed in call of a model
00041 // 03/10/2000  V.Ivanchenko CodeWizard clean up
00042 //
00043 // Class Description: 
00044 //
00045 // Low energy protons/ions nuclear stopping parametrisation
00046 //
00047 // Class Description: End 
00048 //
00049 // -------------------------------------------------------------------
00050 //
00051 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00052 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00053 
00054 #include "G4hNuclearStoppingModel.hh" 
00055 
00056 #include "globals.hh"
00057 #include "G4PhysicalConstants.hh"
00058 #include "G4SystemOfUnits.hh"
00059 #include "G4UnitsTable.hh"
00060 #include "G4hZiegler1977Nuclear.hh"
00061 #include "G4hZiegler1985Nuclear.hh"
00062 #include "G4hICRU49Nuclear.hh"
00063 #include "G4DynamicParticle.hh"
00064 #include "G4ParticleDefinition.hh"
00065 #include "G4ElementVector.hh"
00066 #include "G4Material.hh"
00067 
00068 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00069 
00070 G4hNuclearStoppingModel::G4hNuclearStoppingModel(const G4String& name)
00071   :G4VLowEnergyModel(name), modelName(name)
00072 {
00073   InitializeMe() ;
00074 }
00075 
00076 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00077 
00078 void G4hNuclearStoppingModel::InitializeMe()
00079 {
00080   // Constants
00081   highEnergyLimit = 100.*MeV ;
00082   lowEnergyLimit  = 1.*eV ;
00083   factorPDG2AMU   = 1.007276/proton_mass_c2 ;
00084   theZieglerFactor= eV*cm2*1.0e-15 ; 
00085 
00086   // Registration of parametrisation models of nuclear energy losses
00087   G4String blank  = G4String(" ") ;
00088   G4String zi77   = G4String("Ziegler1977") ;
00089   G4String ir49   = G4String("ICRU_R49") ;
00090   G4String zi85   = G4String("Ziegler1985") ;
00091   if(zi77 == modelName) { 
00092       nStopingPowerTable = new G4hZiegler1977Nuclear();
00093     
00094   } else if(ir49 == modelName || blank == modelName) {
00095       nStopingPowerTable = new G4hICRU49Nuclear();
00096 
00097   } else if(zi85 == modelName) {
00098       nStopingPowerTable = new G4hZiegler1985Nuclear();
00099         
00100   } else {
00101     G4cout << 
00102     "G4hLowEnergyIonisation warning: There is no table with the modelName <" 
00103  << modelName << ">" 
00104  << " for nuclear stopping, <ICRU_R49> is applied " 
00105  << G4endl; 
00106     nStopingPowerTable = new G4hICRU49Nuclear();
00107   }
00108 
00109   // Default is nuclear stopping fluctuations On
00110   //  nStopingPowerTable->SetNuclearStoppingFluctuationsOn();  
00111   nStopingPowerTable->SetNuclearStoppingFluctuationsOff();  
00112 }
00113 
00114 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00115 
00116 G4hNuclearStoppingModel::~G4hNuclearStoppingModel() 
00117 {
00118   delete nStopingPowerTable;
00119 }
00120 
00121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00122 
00123 G4double G4hNuclearStoppingModel::TheValue(
00124                           const G4DynamicParticle* particle,
00125                           const G4Material* material) 
00126 {
00127   // Projectile nucleus
00128   G4double energy = particle->GetKineticEnergy() ;
00129   G4double z1 = std::abs((particle->GetCharge())/eplus) ;
00130   G4double m1 = (particle->GetMass())*factorPDG2AMU ;
00131 
00132   G4double nloss = StoppingPower(material, energy, z1, m1) * theZieglerFactor; 
00133 
00134   return nloss;
00135 }
00136 
00137 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00138 
00139 G4double G4hNuclearStoppingModel::TheValue(
00140                           const G4ParticleDefinition* aParticle,
00141                           const G4Material* material,
00142                                 G4double kineticEnergy)  
00143 {
00144   // Projectile nucleus
00145   G4double z1 = std::abs((aParticle->GetPDGCharge())/eplus) ;
00146   G4double m1 = (aParticle->GetPDGMass())*factorPDG2AMU ;
00147 
00148   G4double nloss = StoppingPower(material, kineticEnergy, z1, m1)
00149                  * theZieglerFactor; 
00150 
00151   return nloss;
00152 }
00153 
00154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00155 
00156 G4double G4hNuclearStoppingModel::StoppingPower(
00157                           const G4Material* material,
00158                                 G4double kineticEnergy,
00159                                 G4double z1, G4double m1) const
00160 {
00161   // Target nucleus
00162   G4int NumberOfElements = material->GetNumberOfElements() ;
00163   if(0 == NumberOfElements) return 0.0 ;
00164 
00165   const G4ElementVector* theElementVector = 
00166                                  material->GetElementVector() ;
00167   const G4double* theAtomicNumDensityVector = 
00168                                  material->GetAtomicNumDensityVector() ;
00169 
00170   //  loop for the elements in the material
00171 
00172   G4double nloss = 0.0;
00173 
00174   for (G4int iel=0; iel<NumberOfElements; iel++) {
00175     const G4Element* element = (*theElementVector)[iel] ;
00176     G4double z2 = element->GetZ();
00177     G4double m2Local = element->GetA()*mole/g ;
00178     nloss += (nStopingPowerTable->
00179               NuclearStoppingPower(kineticEnergy, z1, z2, m1, m2Local))
00180            * theAtomicNumDensityVector[iel] ;    
00181   }
00182 
00183   return nloss;
00184 }
00185 
00186 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00187 
00188 
00189 
00190 
00191 
00192 
00193 
00194 

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