G4RToEConvForGamma.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 // $Id$
00028 //
00029 //
00030 // --------------------------------------------------------------
00031 //      GEANT 4 class implementation file/  History:
00032 //    5 Oct. 2002, H.Kuirashige : Structure created based on object model
00033 // --------------------------------------------------------------
00034 
00035 #include "G4RToEConvForGamma.hh"
00036 #include "G4ParticleDefinition.hh"
00037 #include "G4ParticleTable.hh"
00038 #include "G4Material.hh"
00039 #include "G4PhysicsLogVector.hh"
00040 
00041 #include "G4ios.hh"
00042 #include "G4SystemOfUnits.hh"
00043 
00044 G4RToEConvForGamma::G4RToEConvForGamma() : G4VRangeToEnergyConverter()
00045 {    
00046   theParticle =  G4ParticleTable::GetParticleTable()->FindParticle("gamma");
00047   if (theParticle ==0) {
00048 #ifdef G4VERBOSE
00049     if (GetVerboseLevel()>0) {
00050       G4cout << " G4RToEConvForGamma::G4RToEConvForGamma() ";
00051       G4cout << " Gamma is not defined !!" << G4endl;
00052     }
00053 #endif
00054   } 
00055 }
00056 
00057 G4RToEConvForGamma::~G4RToEConvForGamma()
00058 { 
00059 }
00060 
00061 
00062 // ***********************************************************************
00063 // ******************* BuildAbsorptionLengthVector ***********************
00064 // ***********************************************************************
00065 void G4RToEConvForGamma::BuildAbsorptionLengthVector(
00066                             const G4Material* aMaterial,
00067                             G4RangeVector* absorptionLengthVector )
00068 {
00069   // fill the absorption length vector for this material
00070   // absorption length is defined here as
00071   //
00072   //    absorption length = 5./ macroscopic absorption cross section
00073   //
00074   const G4CrossSectionTable* aCrossSectionTable = (G4CrossSectionTable*)(theLossTable);
00075   const G4ElementVector* elementVector = aMaterial->GetElementVector();
00076   const G4double* atomicNumDensityVector = aMaterial->GetAtomicNumDensityVector();
00077 
00078   //  fill absorption length vector
00079   G4int NumEl = aMaterial->GetNumberOfElements();
00080   G4double absorptionLengthMax = 0.0;
00081   for (size_t ibin=0; ibin<size_t(TotBin); ibin++) {
00082     G4double SIGMA = 0. ;
00083     for (size_t iel=0; iel<size_t(NumEl); iel++) {
00084       G4int IndEl = (*elementVector)[iel]->GetIndex();
00085       SIGMA +=  atomicNumDensityVector[iel]*
00086                    (*((*aCrossSectionTable)[IndEl]))[ibin];
00087     }
00088     //  absorption length=5./SIGMA
00089     absorptionLengthVector->PutValue(ibin, 5./SIGMA);
00090     if (absorptionLengthMax < 5./SIGMA ) absorptionLengthMax = 5./SIGMA;
00091   }
00092 }
00093 
00094 
00095 
00096 // ***********************************************************************
00097 // ********************** ComputeCrossSection ****************************
00098 // ***********************************************************************
00099 G4double G4RToEConvForGamma::ComputeCrossSection(G4double AtomicNumber,
00100                                                  G4double KineticEnergy) const
00101 {
00102   //  Compute the "absorption" cross section of the photon "absorption"
00103   //  cross section means here the sum of the cross sections of the
00104   //  pair production, Compton scattering and photoelectric processes
00105   static G4double Z;  
00106   const  G4double t1keV = 1.*keV;
00107   const  G4double t200keV = 200.*keV;
00108   const  G4double t100MeV = 100.*MeV;
00109 
00110   static G4double s200keV, s1keV;
00111   static G4double tmin, tlow; 
00112   static G4double smin, slow;
00113   static G4double cmin, clow, chigh;
00114   //  compute Z dependent quantities in the case of a new AtomicNumber
00115   if(std::abs(AtomicNumber-Z)>0.1)  {
00116     Z = AtomicNumber;
00117     G4double Zsquare = Z*Z;
00118     G4double Zlog = std::log(Z);
00119     G4double Zlogsquare = Zlog*Zlog;
00120 
00121     s200keV = (0.2651-0.1501*Zlog+0.02283*Zlogsquare)*Zsquare;
00122     tmin = (0.552+218.5/Z+557.17/Zsquare)*MeV;
00123     smin = (0.01239+0.005585*Zlog-0.000923*Zlogsquare)*std::exp(1.5*Zlog);
00124     cmin=std::log(s200keV/smin)/(std::log(tmin/t200keV)*std::log(tmin/t200keV));
00125     tlow = 0.2*std::exp(-7.355/std::sqrt(Z))*MeV;
00126     slow = s200keV*std::exp(0.042*Z*std::log(t200keV/tlow)*std::log(t200keV/tlow));
00127     s1keV = 300.*Zsquare;
00128     clow =std::log(s1keV/slow)/std::log(tlow/t1keV);
00129 
00130     chigh=(7.55e-5-0.0542e-5*Z)*Zsquare*Z/std::log(t100MeV/tmin);
00131   }
00132 
00133   //  calculate the cross section (using an approximate empirical formula)
00134   G4double xs;
00135   if ( KineticEnergy<tlow ) {
00136     if(KineticEnergy<t1keV) xs = slow*std::exp(clow*std::log(tlow/t1keV));
00137     else                    xs = slow*std::exp(clow*std::log(tlow/KineticEnergy));
00138 
00139   } else if ( KineticEnergy<t200keV ) {
00140     xs = s200keV
00141          * std::exp(0.042*Z*std::log(t200keV/KineticEnergy)*std::log(t200keV/KineticEnergy));
00142 
00143   } else if( KineticEnergy<tmin ){
00144     xs = smin
00145          * std::exp(cmin*std::log(tmin/KineticEnergy)*std::log(tmin/KineticEnergy));
00146 
00147   } else {
00148     xs = smin + chigh*std::log(KineticEnergy/tmin);
00149 
00150   }
00151   return xs * barn;
00152 }
00153 

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