G4LivermorePolarizedPhotoElectricModel Class Reference

#include <G4LivermorePolarizedPhotoElectricModel.hh>

Inheritance diagram for G4LivermorePolarizedPhotoElectricModel:

G4VEmModel

Public Member Functions

 G4LivermorePolarizedPhotoElectricModel (const G4String &nam="LivermorePolarizedPhotoElectric")
virtual ~G4LivermorePolarizedPhotoElectricModel ()
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX)
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)

Protected Attributes

G4ParticleChangeForGammafParticleChange

Detailed Description

Definition at line 41 of file G4LivermorePolarizedPhotoElectricModel.hh.


Constructor & Destructor Documentation

G4LivermorePolarizedPhotoElectricModel::G4LivermorePolarizedPhotoElectricModel ( const G4String nam = "LivermorePolarizedPhotoElectric"  ) 

Definition at line 49 of file G4LivermorePolarizedPhotoElectricModel.cc.

References G4Electron::Electron(), G4cout, G4endl, G4Gamma::Gamma(), G4VEmModel::SetAngularDistribution(), and G4VEmModel::SetDeexcitationFlag().

00051   :G4VEmModel(nam),fParticleChange(0),
00052    crossSectionHandler(0), shellCrossSectionHandler(0)
00053 {
00054   lowEnergyLimit  = 10 * eV; // SI - Could be 10 eV ?
00055   highEnergyLimit = 100 * GeV;
00056   
00057   verboseLevel= 0;
00058   // Verbosity scale:
00059   // 0 = nothing 
00060   // 1 = warning for energy non-conservation 
00061   // 2 = details of energy budget
00062   // 3 = calculation of cross sections, file openings, sampling of atoms
00063   // 4 = entering in methods
00064 
00065   theGamma    = G4Gamma::Gamma();
00066   theElectron = G4Electron::Electron();
00067 
00068   // use default generator
00069   SetAngularDistribution(new G4SauterGavrilaAngularDistribution());
00070 
00071   // polarized generator needs fix
00072   //SetAngularDistribution(new G4PhotoElectricAngularGeneratorPolarized());
00073   SetDeexcitationFlag(true);
00074   fAtomDeexcitation = 0; 
00075   fDeexcitationActive = false; 
00076 
00077   if (verboseLevel > 1) {
00078     G4cout << "Livermore Polarized PhotoElectric is constructed " << G4endl
00079            << "Energy range: "
00080            << lowEnergyLimit / keV << " keV - "
00081            << highEnergyLimit / GeV << " GeV"
00082            << G4endl;
00083   }
00084 }

G4LivermorePolarizedPhotoElectricModel::~G4LivermorePolarizedPhotoElectricModel (  )  [virtual]

Definition at line 88 of file G4LivermorePolarizedPhotoElectricModel.cc.

00089 {  
00090   delete crossSectionHandler;
00091   delete shellCrossSectionHandler;
00092 }


Member Function Documentation

G4double G4LivermorePolarizedPhotoElectricModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition ,
G4double  kinEnergy,
G4double  Z,
G4double  A = 0,
G4double  cut = 0,
G4double  emax = DBL_MAX 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 146 of file G4LivermorePolarizedPhotoElectricModel.cc.

References G4VCrossSectionHandler::FindValue(), G4cout, G4endl, and G4lrint().

00151 {
00152   G4double cs = crossSectionHandler->FindValue(G4lrint(Z), GammaEnergy);
00153 
00154   if (verboseLevel > 3) {
00155     G4cout << "G4LivermorePolarizedPhotoElectricModel::ComputeCrossSectionPerAtom()" 
00156            << G4endl;
00157     G4cout << "  E(keV)= " << GammaEnergy/keV << "  Z= " << Z 
00158            << "  CrossSection(barn)= "
00159            << cs/barn << G4endl; 
00160   }
00161   return cs;
00162 }

void G4LivermorePolarizedPhotoElectricModel::Initialise ( const G4ParticleDefinition ,
const G4DataVector  
) [virtual]

Implements G4VEmModel.

Definition at line 97 of file G4LivermorePolarizedPhotoElectricModel.cc.

References G4LossTableManager::AtomDeexcitation(), G4VCrossSectionHandler::Clear(), fParticleChange, G4cout, G4endl, G4VEmModel::GetParticleChangeForGamma(), G4VEmModel::HighEnergyLimit(), G4LossTableManager::Instance(), G4VAtomDeexcitation::IsFluoActive(), G4VCrossSectionHandler::LoadData(), G4VCrossSectionHandler::LoadShellData(), and G4VEmModel::LowEnergyLimit().

00100 {
00101   if (verboseLevel > 3) {
00102     G4cout << "Calling G4LivermorePolarizedPhotoElectricModel::Initialise()" << G4endl;
00103   }
00104   if (crossSectionHandler)
00105   {
00106     crossSectionHandler->Clear();
00107     delete crossSectionHandler;
00108   }
00109 
00110   if (shellCrossSectionHandler)
00111     {
00112       shellCrossSectionHandler->Clear();
00113       delete shellCrossSectionHandler;
00114     }
00115 
00116   
00117   // Reading of data files - all materials are read  
00118   crossSectionHandler = new G4CrossSectionHandler;
00119   crossSectionHandler->Clear();
00120   G4String crossSectionFile = "phot/pe-cs-";
00121   crossSectionHandler->LoadData(crossSectionFile);
00122 
00123   shellCrossSectionHandler = new G4CrossSectionHandler();
00124   shellCrossSectionHandler->Clear();
00125   G4String shellCrossSectionFile = "phot/pe-ss-cs-";
00126   shellCrossSectionHandler->LoadShellData(shellCrossSectionFile);
00127 
00128   fParticleChange = GetParticleChangeForGamma();
00129 
00130   fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
00131   if(fAtomDeexcitation) { 
00132     fDeexcitationActive = fAtomDeexcitation->IsFluoActive(); 
00133   }
00134   if (verboseLevel > 1) {
00135     G4cout << "Livermore Polarized PhotoElectric model is initialized " 
00136            << G4endl
00137            << "Energy range: "
00138            << LowEnergyLimit() / keV << " keV - "
00139            << HighEnergyLimit() / GeV << " GeV"
00140            << G4endl;
00141   }
00142 }

void G4LivermorePolarizedPhotoElectricModel::SampleSecondaries ( std::vector< G4DynamicParticle * > *  ,
const G4MaterialCutsCouple ,
const G4DynamicParticle ,
G4double  tmin,
G4double  maxEnergy 
) [virtual]

Implements G4VEmModel.

Definition at line 166 of file G4LivermorePolarizedPhotoElectricModel.cc.

References G4AtomicShell::BindingEnergy(), G4InuclSpecialFunctions::bindingEnergy(), G4VAtomDeexcitation::CheckDeexcitationActiveRegion(), fParticleChange, fStopAndKill, G4cout, G4endl, G4lrint(), G4VAtomDeexcitation::GenerateParticles(), G4VEmModel::GetAngularDistribution(), G4Element::GetAtomicShell(), G4VAtomDeexcitation::GetAtomicShell(), G4MaterialCutsCouple::GetIndex(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4Element::GetNbOfAtomicShells(), G4Element::GetZ(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4VParticleChange::ProposeTrackStatus(), G4VEmAngularDistribution::SampleDirection(), G4VEmModel::SelectRandomAtom(), G4VCrossSectionHandler::SelectRandomShell(), and G4ParticleChangeForGamma::SetProposedKineticEnergy().

00172 {
00173   if (verboseLevel > 3) {
00174     G4cout << "Calling SampleSecondaries() of G4LivermorePolarizedPhotoElectricModel" 
00175            << G4endl;
00176   }
00177 
00178   G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
00179   
00180   // kill incident photon  
00181   fParticleChange->SetProposedKineticEnergy(0.);
00182   fParticleChange->ProposeTrackStatus(fStopAndKill);
00183   
00184   // low-energy gamma is absorpted by this process
00185 
00186   if (photonEnergy <= lowEnergyLimit)
00187     {
00188       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy);
00189       return;
00190     }
00191 
00192   // Select randomly one element in the current material
00193   //G4cout << "Select random atom Egamma(keV)= " << photonEnergy/keV << G4endl;
00194   const G4Element* elm = 
00195     SelectRandomAtom(couple->GetMaterial(),theGamma,photonEnergy);
00196   G4int Z = G4lrint(elm->GetZ());
00197 
00198   // Select the ionised shell in the current atom according to shell cross sections
00199   //G4cout << "Select random shell Z= " << Z << G4endl;
00200   size_t shellIndex = shellCrossSectionHandler->SelectRandomShell(Z,photonEnergy);
00201   //G4cout << "Shell index= " << shellIndex 
00202   //     << " nShells= " << elm->GetNbOfAtomicShells() << G4endl;
00203   G4double bindingEnergy;
00204   const G4AtomicShell* shell = 0;
00205   if(fDeexcitationActive) {
00206     G4AtomicShellEnumerator as = G4AtomicShellEnumerator(shellIndex);
00207     shell = fAtomDeexcitation->GetAtomicShell(Z, as);
00208     bindingEnergy = shell->BindingEnergy();
00209   } else {
00210     G4int nshells = elm->GetNbOfAtomicShells() - 1;
00211     if(G4int(shellIndex) > nshells) { shellIndex = std::max(0, nshells); }
00212     bindingEnergy = elm->GetAtomicShell(shellIndex);
00213   }
00214 
00215   // There may be cases where the binding energy of the selected 
00216   // shell is > photon energy
00217   // In such cases do not generate secondaries
00218   if(photonEnergy < bindingEnergy) {
00219     fParticleChange->ProposeLocalEnergyDeposit(photonEnergy);
00220     return;
00221   }
00222   //G4cout << "Z= " << Z <<  "  shellIndex= " << shellIndex 
00223   //     << " Ebind(keV)= " <<  bindingEnergy/keV << G4endl; 
00224 
00225 
00226   // Primary outcoming electron
00227   G4double eKineticEnergy = photonEnergy - bindingEnergy;
00228   G4double edep = bindingEnergy;
00229 
00230   // Calculate direction of the photoelectron
00231   G4ThreeVector electronDirection = 
00232     GetAngularDistribution()->SampleDirection(aDynamicGamma, 
00233                                               eKineticEnergy,
00234                                               shellIndex, 
00235                                               couple->GetMaterial());
00236 
00237   // The electron is created ...
00238   G4DynamicParticle* electron = new G4DynamicParticle (theElectron,
00239                                                        electronDirection,
00240                                                        eKineticEnergy);
00241   fvect->push_back(electron);
00242 
00243   // Sample deexcitation
00244   if(fDeexcitationActive) {
00245     G4int index = couple->GetIndex();
00246     if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
00247       size_t nbefore = fvect->size();
00248 
00249       fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
00250       size_t nafter = fvect->size();
00251       if(nafter > nbefore) {
00252         for (size_t j=nbefore; j<nafter; ++j) {
00253           edep -= ((*fvect)[j])->GetKineticEnergy();
00254         } 
00255       }
00256     }
00257   }
00258 
00259   // energy balance - excitation energy left
00260   if(edep > 0.0) {
00261     fParticleChange->ProposeLocalEnergyDeposit(edep);
00262   }
00263 }


Field Documentation

G4ParticleChangeForGamma* G4LivermorePolarizedPhotoElectricModel::fParticleChange [protected]

Definition at line 69 of file G4LivermorePolarizedPhotoElectricModel.hh.

Referenced by Initialise(), and SampleSecondaries().


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:26 2013 for Geant4 by  doxygen 1.4.7