G4PEEffectFluoModel Class Reference

#include <G4PEEffectFluoModel.hh>

Inheritance diagram for G4PEEffectFluoModel:

G4VEmModel G4PolarizedPEEffectModel

Public Member Functions

 G4PEEffectFluoModel (const G4String &nam="PhotoElectric")
virtual ~G4PEEffectFluoModel ()
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A, G4double, G4double)
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy)
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)

Detailed Description

Definition at line 57 of file G4PEEffectFluoModel.hh.


Constructor & Destructor Documentation

G4PEEffectFluoModel::G4PEEffectFluoModel ( const G4String nam = "PhotoElectric"  ) 

Definition at line 65 of file G4PEEffectFluoModel.cc.

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

00066   : G4VEmModel(nam)
00067 {
00068   theGamma    = G4Gamma::Gamma();
00069   theElectron = G4Electron::Electron();
00070   fminimalEnergy = 1.0*eV;
00071   SetDeexcitationFlag(true);
00072   fParticleChange = 0;
00073   fAtomDeexcitation = 0;
00074 
00075   // default generator
00076   SetAngularDistribution(new G4SauterGavrilaAngularDistribution());
00077 }

G4PEEffectFluoModel::~G4PEEffectFluoModel (  )  [virtual]

Definition at line 81 of file G4PEEffectFluoModel.cc.

00082 {}


Member Function Documentation

G4double G4PEEffectFluoModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition ,
G4double  kinEnergy,
G4double  Z,
G4double  A,
G4double  ,
G4double   
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 96 of file G4PEEffectFluoModel.cc.

References G4SandiaTable::GetSandiaCofPerAtom().

00100 {
00101   G4double* SandiaCof = G4SandiaTable::GetSandiaCofPerAtom((G4int)Z, energy);
00102 
00103   G4double energy2 = energy*energy;
00104   G4double energy3 = energy*energy2;
00105   G4double energy4 = energy2*energy2;
00106 
00107   return SandiaCof[0]/energy  + SandiaCof[1]/energy2 +
00108     SandiaCof[2]/energy3 + SandiaCof[3]/energy4;
00109 }

G4double G4PEEffectFluoModel::CrossSectionPerVolume ( const G4Material ,
const G4ParticleDefinition ,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 114 of file G4PEEffectFluoModel.cc.

References G4SandiaTable::GetSandiaCofForMaterial(), and G4Material::GetSandiaTable().

00118 {
00119   G4double* SandiaCof = 
00120     material->GetSandiaTable()->GetSandiaCofForMaterial(energy);
00121                                 
00122   G4double energy2 = energy*energy;
00123   G4double energy3 = energy*energy2;
00124   G4double energy4 = energy2*energy2;
00125           
00126   return SandiaCof[0]/energy  + SandiaCof[1]/energy2 +
00127     SandiaCof[2]/energy3 + SandiaCof[3]/energy4; 
00128 }

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

Implements G4VEmModel.

Reimplemented in G4PolarizedPEEffectModel.

Definition at line 86 of file G4PEEffectFluoModel.cc.

References G4LossTableManager::AtomDeexcitation(), G4VEmModel::GetParticleChangeForGamma(), and G4LossTableManager::Instance().

Referenced by G4PolarizedPEEffectModel::Initialise().

00088 {
00089   fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation();
00090   if(!fParticleChange) { fParticleChange = GetParticleChangeForGamma(); }
00091 }

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

Implements G4VEmModel.

Reimplemented in G4PolarizedPEEffectModel.

Definition at line 133 of file G4PEEffectFluoModel.cc.

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

Referenced by G4PolarizedPEEffectModel::SampleSecondaries().

00138 {
00139   const G4Material* aMaterial = couple->GetMaterial();
00140 
00141   G4double energy = aDynamicPhoton->GetKineticEnergy();
00142 
00143   // select randomly one element constituing the material.
00144   const G4Element* anElement = SelectRandomAtom(aMaterial,theGamma,energy);
00145   
00146   //
00147   // Photo electron
00148   //
00149 
00150   // Select atomic shell
00151   G4int nShells = anElement->GetNbOfAtomicShells();
00152   G4int i = 0;  
00153   for(; i<nShells; ++i) {
00154     /*
00155     G4cout << "i= " << i << " E(eV)= " << energy/eV 
00156            << " Eb(eV)= " << anElement->GetAtomicShell(i)/eV
00157            << "  " << anElement->GetName() 
00158            << G4endl;
00159     */
00160     if(energy >= anElement->GetAtomicShell(i)) { break; }
00161   }
00162 
00163   G4double edep = energy;
00164 
00165   // Normally one shell is available 
00166   if (i < nShells) { 
00167   
00168     G4double bindingEnergy  = anElement->GetAtomicShell(i);
00169     G4double elecKineEnergy = energy - bindingEnergy;
00170 
00171     // create photo electron
00172     //
00173     if (elecKineEnergy > fminimalEnergy) {
00174       edep = bindingEnergy;
00175       G4ThreeVector elecDirection =
00176         GetAngularDistribution()->SampleDirection(aDynamicPhoton, 
00177                                                   elecKineEnergy,
00178                                                   i, 
00179                                                   couple->GetMaterial());
00180    
00181       G4DynamicParticle* aParticle = 
00182         new G4DynamicParticle(theElectron, elecDirection, elecKineEnergy);
00183       fvect->push_back(aParticle);
00184     } 
00185 
00186     // sample deexcitation
00187     //
00188     if(fAtomDeexcitation) {
00189       G4int index = couple->GetIndex();
00190       if(fAtomDeexcitation->CheckDeexcitationActiveRegion(index)) {
00191         G4int Z = G4lrint(anElement->GetZ());
00192         G4AtomicShellEnumerator as = G4AtomicShellEnumerator(i);
00193         const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as);
00194         size_t nbefore = fvect->size();
00195         fAtomDeexcitation->GenerateParticles(fvect, shell, Z, index);
00196         size_t nafter = fvect->size();
00197         if(nafter > nbefore) {
00198           for (size_t j=nbefore; j<nafter; ++j) {
00199             edep -= ((*fvect)[j])->GetKineticEnergy();
00200           } 
00201         }
00202       }
00203     }
00204   }
00205 
00206   // kill primary photon
00207   fParticleChange->SetProposedKineticEnergy(0.);
00208   fParticleChange->ProposeTrackStatus(fStopAndKill);
00209   if(edep > 0.0) {
00210     fParticleChange->ProposeLocalEnergyDeposit(edep);
00211   }
00212 }


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