G4BremsstrahlungCrossSectionHandler.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 // -------------------------------------------------------------------
00029 //
00030 // GEANT4 Class file
00031 //
00032 //
00033 // File name:     G4BremsstrahlungCrossSectionHandler
00034 //
00035 // Author:        V.Ivanchenko (Vladimir.Ivanchenko@cern.ch)
00036 // 
00037 // Creation date: 25 September 2001
00038 //
00039 // Modifications:
00040 // 
00041 // 10.10.2001 MGP Revision to improve code quality and consistency with design
00042 // 21.01.2003 VI  cut per region
00043 // 03.03.2009 LP  Added public method to make a easier migration of
00044 //                G4LowEnergyBremsstrahlung to G4LivermoreBremsstrahlungModel
00045 //
00046 // 15 Jul 2009   Nicolas A. Karakatsanis
00047 //
00048 //                           - BuildCrossSectionForMaterials method was revised in order to calculate the 
00049 //                             logarithmic values of the loaded data. 
00050 //                             It retrieves the data values from the G4EMLOW data files but, then, calculates the
00051 //                             respective log values and loads them to seperate data structures.
00052 //                             The EM data sets, initialized this way, contain both non-log and log values.
00053 //                             These initialized data sets can enhance the computing performance of data interpolation
00054 //                             operations
00055 //
00056 //
00057 //
00058 //
00059 // -------------------------------------------------------------------
00060 
00061 #include "G4BremsstrahlungCrossSectionHandler.hh"
00062 #include "G4eBremsstrahlungSpectrum.hh"
00063 #include "G4DataVector.hh"
00064 #include "G4CompositeEMDataSet.hh"
00065 #include "G4VDataSetAlgorithm.hh"
00066 #include "G4SemiLogInterpolation.hh"
00067 #include "G4VEMDataSet.hh"
00068 #include "G4EMDataSet.hh"
00069 #include "G4Material.hh"
00070 #include "G4ProductionCutsTable.hh"
00071 
00072 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00073 
00074 G4BremsstrahlungCrossSectionHandler::G4BremsstrahlungCrossSectionHandler(const G4VEnergySpectrum* spec,
00075                                                                          G4VDataSetAlgorithm* )
00076   : theBR(spec)
00077 {
00078   interp = new G4SemiLogInterpolation();
00079 }
00080 
00081 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00082 
00083 G4BremsstrahlungCrossSectionHandler::~G4BremsstrahlungCrossSectionHandler()
00084 {
00085   delete interp;
00086 }
00087 
00088 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00089 
00090 std::vector<G4VEMDataSet*>*
00091 G4BremsstrahlungCrossSectionHandler::BuildCrossSectionsForMaterials(const G4DataVector& energyVector,
00092                                                                     const G4DataVector* energyCuts)
00093 {
00094   std::vector<G4VEMDataSet*>* set = new std::vector<G4VEMDataSet*>;
00095 
00096   G4DataVector* energies;
00097   G4DataVector* cs;
00098 
00099   G4DataVector* log_energies;
00100   G4DataVector* log_cs;
00101 
00102   G4int nOfBins = energyVector.size();
00103 
00104   const G4ProductionCutsTable* theCoupleTable=
00105         G4ProductionCutsTable::GetProductionCutsTable();
00106   size_t numOfCouples = theCoupleTable->GetTableSize();
00107 
00108   for (size_t mLocal=0; mLocal<numOfCouples; mLocal++) {
00109 
00110     const G4MaterialCutsCouple* couple = theCoupleTable->GetMaterialCutsCouple(mLocal);
00111     const G4Material* material= couple->GetMaterial();
00112     const G4ElementVector* elementVector = material->GetElementVector();
00113     const G4double* nAtomsPerVolume = material->GetVecNbOfAtomsPerVolume();
00114     G4int nElements = material->GetNumberOfElements();
00115 
00116     G4double tcut  = (*energyCuts)[mLocal];
00117 
00118     G4VDataSetAlgorithm* algo = interp->Clone();
00119     G4VEMDataSet* setForMat = new G4CompositeEMDataSet(algo,1.,1.);
00120 
00121     for (G4int i=0; i<nElements; i++) {
00122 
00123       G4int Z = (G4int) ((*elementVector)[i]->GetZ());
00124 
00125       energies = new G4DataVector;
00126       cs       = new G4DataVector;
00127 
00128       log_energies = new G4DataVector;
00129       log_cs       = new G4DataVector;
00130 
00131       G4double density = nAtomsPerVolume[i];
00132 
00133       for (G4int bin=0; bin<nOfBins; bin++) {
00134 
00135         G4double e = energyVector[bin];
00136         energies->push_back(e);
00137         if (e==0.) e=1e-300;
00138         log_energies->push_back(std::log10(e));
00139         G4double value = 0.0;
00140 
00141         if(e > tcut) {
00142           G4double elemCs = FindValue(Z, e);
00143 
00144           value  = theBR->Probability(Z, tcut, e, e);
00145 
00146           value *= elemCs*density;
00147         }
00148         cs->push_back(value);
00149 
00150         if (value==0.) value=1e-300;
00151         log_cs->push_back(std::log10(value));
00152       }
00153       G4VDataSetAlgorithm* algol = interp->Clone();
00154 
00155       //G4VEMDataSet* elSet = new G4EMDataSet(i,energies,cs,algol,1.,1.);
00156 
00157       G4VEMDataSet* elSet = new G4EMDataSet(i,energies,cs,log_energies,log_cs,algol,1.,1.);
00158 
00159       setForMat->AddComponent(elSet);
00160     }
00161     set->push_back(setForMat);
00162   }
00163 
00164   return set;
00165 }
00166 
00167 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
00168 
00169 G4double G4BremsstrahlungCrossSectionHandler::GetCrossSectionAboveThresholdForElement(G4double energy,
00170                                                                                       G4double cutEnergy,
00171                                                                                       G4int Z)
00172 {
00173   G4double value = 0.;
00174   if(energy > cutEnergy)
00175     {
00176       G4double elemCs = FindValue(Z, energy);
00177       value  = theBR->Probability(Z,cutEnergy, energy, energy);   
00178       value *= elemCs;
00179     }
00180   return value;
00181 }

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