G4LivermoreIonisationCrossSection.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 // $Id$
00027 //
00028 // Author: Vladimir Ivanchenko
00029 //
00030 // History:
00031 // --------
00032 // 31 May 2011 V.Ivanchenko  Created  
00033 // 09 Mar 2012 L.Pandola     update methods
00034 // 
00035 //
00036 
00037 #include "G4LivermoreIonisationCrossSection.hh"
00038 #include "G4SystemOfUnits.hh"
00039 #include "G4AtomicTransitionManager.hh"
00040 #include "G4VCrossSectionHandler.hh"
00041 #include "G4eCrossSectionHandler.hh"
00042 #include "G4SemiLogInterpolation.hh"
00043 
00044 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00045 
00046 
00047 G4LivermoreIonisationCrossSection::G4LivermoreIonisationCrossSection(
00048   const G4String& nam) : G4VhShellCrossSection(nam), crossSectionHandler(0)
00049 {
00050   fLowEnergyLimit  = 10.0*eV;
00051   fHighEnergyLimit = 100.0*GeV;
00052 
00053   transitionManager = G4AtomicTransitionManager::Instance();
00054 
00055   verboseLevel = 0;
00056   
00057   Initialise();
00058 }
00059 
00060 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00061 
00062 G4LivermoreIonisationCrossSection::~G4LivermoreIonisationCrossSection()
00063 {
00064   delete crossSectionHandler;
00065 }
00066 
00067 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00068 
00069 void G4LivermoreIonisationCrossSection::Initialise()
00070 {
00071   const G4int binForFluo = 20;
00072   G4int nbin = G4int(std::log10(fHighEnergyLimit/fLowEnergyLimit) + 0.5);
00073   if(nbin <= 0) { nbin = 1; }
00074   nbin *= binForFluo;
00075 
00076   // Data on shell ionisation x-sections
00077   if (crossSectionHandler) { 
00078     crossSectionHandler->Clear();
00079     delete crossSectionHandler; 
00080   }
00081 
00082   G4VDataSetAlgorithm* inter = new G4SemiLogInterpolation();
00083   crossSectionHandler = 
00084     new G4eCrossSectionHandler(inter,fLowEnergyLimit,fHighEnergyLimit,nbin);
00085   crossSectionHandler->LoadShellData("ioni/ion-ss-cs-");
00086   //G4cout << "!!!  G4LivermoreIonisationCrossSection::Initialise()" << G4endl;
00087 }
00088 
00089 G4double 
00090 G4LivermoreIonisationCrossSection::CrossSection(G4int Z, G4AtomicShellEnumerator shell,
00091                                                 G4double kinEnergy, G4double,
00092                                                 const G4Material*) 
00093 {
00094   G4double cross = 0.0;
00095   G4int n = G4int(shell);
00096   G4int nmax = std::min(9,transitionManager->NumberOfShells(Z));
00097   if(Z > 6 && Z < 93 && n < nmax && 
00098      kinEnergy >= fLowEnergyLimit && kinEnergy <= fHighEnergyLimit) {
00099     //G4cout << "Z= " << Z << "  n= " << n << " E(MeV)= " << kinEnergy/MeV << G4endl;
00100     cross = crossSectionHandler->FindValue(Z, kinEnergy, n);
00101   }
00102   return cross;
00103 }
00104 
00105 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00106 
00107 std::vector<G4double> 
00108 G4LivermoreIonisationCrossSection::GetCrossSection(G4int Z,
00109                                                    G4double kinEnergy,
00110                                                    G4double, G4double, 
00111                                                    const G4Material*)
00112 {
00113   G4int nmax = std::min(9,transitionManager->NumberOfShells(Z));
00114   std::vector<G4double> vec(nmax,0.0); 
00115   for(G4int i=0; i<nmax; ++i) {
00116     vec[i] = CrossSection(Z, G4AtomicShellEnumerator(i), kinEnergy); 
00117   }
00118   return vec;
00119 }
00120 
00121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00122 
00123 std::vector<G4double> 
00124 G4LivermoreIonisationCrossSection::Probabilities(G4int Z,
00125                                                  G4double kinEnergy,
00126                                                  G4double,
00127                                                  G4double,
00128                                                  const G4Material*)
00129 {
00130   std::vector<G4double> vec = GetCrossSection(Z, kinEnergy);
00131   size_t n = vec.size();
00132   size_t i;
00133   G4double sum = 0.0;
00134   for(i=0; i<n; ++i) { sum += vec[i]; }
00135   if(sum > 0.0) { 
00136     sum = 1.0/sum; 
00137     for(i=0; i<n; ++i) { vec[i] = vec[i]*sum; }
00138   }
00139   return vec;
00140 }
00141 
00142 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....

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