G4LivermoreRayleighModel Class Reference

#include <G4LivermoreRayleighModel.hh>

Inheritance diagram for G4LivermoreRayleighModel:

G4VEmModel

Public Member Functions

 G4LivermoreRayleighModel ()
virtual ~G4LivermoreRayleighModel ()
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)
void SetLowEnergyThreshold (G4double)

Detailed Description

Definition at line 39 of file G4LivermoreRayleighModel.hh.


Constructor & Destructor Documentation

G4LivermoreRayleighModel::G4LivermoreRayleighModel (  ) 

Definition at line 40 of file G4LivermoreRayleighModel.cc.

References G4cout, G4endl, and G4VEmModel::SetAngularDistribution().

00041   :G4VEmModel("LivermoreRayleigh"),isInitialised(false),maxZ(100)
00042 {
00043   fParticleChange = 0;
00044   lowEnergyLimit  = 10 * eV; 
00045   
00046   dataCS.resize(maxZ+1,0); 
00047 
00048   SetAngularDistribution(new G4RayleighAngularGenerator());
00049   
00050   verboseLevel= 0;
00051   // Verbosity scale for debugging purposes:
00052   // 0 = nothing 
00053   // 1 = calculation of cross sections, file openings...
00054   // 2 = entering in methods
00055 
00056   if(verboseLevel > 0) 
00057   {
00058     G4cout << "G4LivermoreRayleighModel is constructed " << G4endl;
00059   }
00060 
00061 }

G4LivermoreRayleighModel::~G4LivermoreRayleighModel (  )  [virtual]

Definition at line 65 of file G4LivermoreRayleighModel.cc.

00066 {  
00067   for(G4int i=0; i<=maxZ; ++i) { delete dataCS[i]; }
00068 }


Member Function Documentation

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

Reimplemented from G4VEmModel.

Definition at line 178 of file G4LivermoreRayleighModel.cc.

References G4PhysicsVector::Energy(), G4cout, G4endl, G4lrint(), G4PhysicsVector::GetVectorLength(), CLHEP::detail::n, and G4PhysicsVector::Value().

00183 {
00184   if (verboseLevel > 1) 
00185   {
00186     G4cout << "Calling ComputeCrossSectionPerAtom() of G4LivermoreRayleighModel" 
00187            << G4endl;
00188   }
00189 
00190   if(GammaEnergy < lowEnergyLimit) { return 0.0; }
00191   
00192   G4double xs = 0.0;
00193   
00194   G4int intZ = G4lrint(Z);
00195 
00196   if(intZ < 1 || intZ > maxZ) { return xs; }
00197 
00198   G4LPhysicsFreeVector* pv = dataCS[intZ];
00199 
00200   // element was not initialised
00201   if(!pv) 
00202   {
00203     char* path = getenv("G4LEDATA");
00204     ReadData(intZ, path);
00205     pv = dataCS[intZ];
00206     if(!pv) { return xs; }
00207   }
00208 
00209   G4int n = pv->GetVectorLength() - 1;
00210   G4double e = GammaEnergy/MeV;
00211   if(e >= pv->Energy(n)) {
00212     xs = (*pv)[n]/(e*e);  
00213   } else if(e >= pv->Energy(0)) {
00214     xs = pv->Value(e)/(e*e);  
00215   }
00216 
00217   if(verboseLevel > 0)
00218   {
00219     G4cout  <<  "****** DEBUG: tcs value for Z=" << Z << " at energy (MeV)=" << e << G4endl;
00220     G4cout  <<  "  cs (Geant4 internal unit)=" << xs << G4endl;
00221     G4cout  <<  "    -> first E*E*cs value in CS data file (iu) =" << (*pv)[0] << G4endl;
00222     G4cout  <<  "    -> last  E*E*cs value in CS data file (iu) =" << (*pv)[n] << G4endl;
00223     G4cout  <<  "*********************************************************" << G4endl;
00224   }
00225   return xs;
00226 }

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

Implements G4VEmModel.

Definition at line 72 of file G4LivermoreRayleighModel.cc.

References G4cout, G4endl, G4lrint(), G4Material::GetElementVector(), G4MaterialCutsCouple::GetMaterial(), G4ProductionCutsTable::GetMaterialCutsCouple(), G4Material::GetNumberOfElements(), G4VEmModel::GetParticleChangeForGamma(), G4ProductionCutsTable::GetProductionCutsTable(), G4ProductionCutsTable::GetTableSize(), G4VEmModel::HighEnergyLimit(), G4VEmModel::InitialiseElementSelectors(), and G4VEmModel::LowEnergyLimit().

00074 {
00075   if (verboseLevel > 1) 
00076   {
00077     G4cout << "Calling Initialise() of G4LivermoreRayleighModel." << G4endl
00078            << "Energy range: "
00079            << LowEnergyLimit() / eV << " eV - "
00080            << HighEnergyLimit() / GeV << " GeV"
00081            << G4endl;
00082   }
00083 
00084   // Initialise element selector
00085   InitialiseElementSelectors(particle, cuts);
00086 
00087   // Access to elements
00088   
00089   char* path = getenv("G4LEDATA");
00090   G4ProductionCutsTable* theCoupleTable =
00091     G4ProductionCutsTable::GetProductionCutsTable();
00092   G4int numOfCouples = theCoupleTable->GetTableSize();
00093   
00094   for(G4int i=0; i<numOfCouples; ++i) 
00095   {
00096     const G4MaterialCutsCouple* couple = 
00097       theCoupleTable->GetMaterialCutsCouple(i);
00098     const G4Material* material = couple->GetMaterial();
00099     const G4ElementVector* theElementVector = material->GetElementVector();
00100     G4int nelm = material->GetNumberOfElements();
00101     
00102     for (G4int j=0; j<nelm; ++j) 
00103     {
00104       G4int Z = G4lrint((*theElementVector)[j]->GetZ());
00105       if(Z < 1)          { Z = 1; }
00106       else if(Z > maxZ)  { Z = maxZ; }
00107       if( (!dataCS[Z]) ) { ReadData(Z, path); }
00108     }
00109   }
00110   //
00111   
00112   //fGenerator->PrintGeneratorInformation();
00113   
00114   if(isInitialised) { return; }
00115   fParticleChange = GetParticleChangeForGamma();
00116   isInitialised = true;
00117 
00118 }

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

Implements G4VEmModel.

Definition at line 230 of file G4LivermoreRayleighModel.cc.

References fStopAndKill, G4cout, G4endl, G4lrint(), G4VEmModel::GetAngularDistribution(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4Element::GetZ(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), G4VEmAngularDistribution::SampleDirection(), G4VEmModel::SelectRandomAtom(), and G4ParticleChangeForGamma::SetProposedKineticEnergy().

00235 {
00236   if (verboseLevel > 1) {
00237     G4cout << "Calling SampleSecondaries() of G4LivermoreRayleighModel" 
00238            << G4endl;
00239   }
00240   G4double photonEnergy0 = aDynamicGamma->GetKineticEnergy();
00241 
00242   // absorption of low-energy gamma  
00243   if (photonEnergy0 <= lowEnergyLimit)
00244     {
00245       fParticleChange->ProposeTrackStatus(fStopAndKill);
00246       fParticleChange->SetProposedKineticEnergy(0.);
00247       fParticleChange->ProposeLocalEnergyDeposit(photonEnergy0);
00248       return ;
00249     }
00250 
00251   // Select randomly one element in the current material
00252   const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
00253   const G4Element* elm = SelectRandomAtom(couple,particle,photonEnergy0);
00254   G4int Z = G4lrint(elm->GetZ());
00255 
00256   // Sample the angle of the scattered photon
00257   
00258   G4ThreeVector photonDirection = 
00259     GetAngularDistribution()->SampleDirection(aDynamicGamma, 
00260                                               photonEnergy0, Z, couple->GetMaterial());
00261   fParticleChange->ProposeMomentumDirection(photonDirection);
00262 }

void G4LivermoreRayleighModel::SetLowEnergyThreshold ( G4double   )  [inline]

Definition at line 85 of file G4LivermoreRayleighModel.hh.

00086 {
00087   lowEnergyLimit = val;
00088 }


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