G4DNAScreenedRutherfordElasticModel Class Reference

#include <G4DNAScreenedRutherfordElasticModel.hh>

Inheritance diagram for G4DNAScreenedRutherfordElasticModel:

G4VEmModel

Public Member Functions

 G4DNAScreenedRutherfordElasticModel (const G4ParticleDefinition *p=0, const G4String &nam="DNAScreenedRutherfordElasticModel")
virtual ~G4DNAScreenedRutherfordElasticModel ()
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
virtual G4double CrossSectionPerVolume (const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax)
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
void SetKillBelowThreshold (G4double threshold)
G4double GetKillBelowThreshold ()

Protected Attributes

G4ParticleChangeForGammafParticleChangeForGamma

Detailed Description

Definition at line 39 of file G4DNAScreenedRutherfordElasticModel.hh.


Constructor & Destructor Documentation

G4DNAScreenedRutherfordElasticModel::G4DNAScreenedRutherfordElasticModel ( const G4ParticleDefinition p = 0,
const G4String nam = "DNAScreenedRutherfordElasticModel" 
)

Definition at line 41 of file G4DNAScreenedRutherfordElasticModel.cc.

References G4cout, and G4endl.

00042     :G4VEmModel(nam),isInitialised(false)
00043 {
00044     //  nistwater = G4NistManager::Instance()->FindOrBuildMaterial("G4_WATER");
00045     fpWaterDensity = 0;
00046 
00047     killBelowEnergy = 9*eV;
00048     lowEnergyLimit = 0 * eV;
00049     intermediateEnergyLimit = 200 * eV; // Switch between two final state models
00050     highEnergyLimit = 1. * MeV;
00051     SetLowEnergyLimit(lowEnergyLimit);
00052     SetHighEnergyLimit(highEnergyLimit);
00053 
00054     verboseLevel= 0;
00055     // Verbosity scale:
00056     // 0 = nothing
00057     // 1 = warning for energy non-conservation
00058     // 2 = details of energy budget
00059     // 3 = calculation of cross sections, file openings, sampling of atoms
00060     // 4 = entering in methods
00061 
00062     if( verboseLevel>0 )
00063     {
00064         G4cout << "Screened Rutherford Elastic model is constructed " << G4endl
00065                << "Energy range: "
00066                << lowEnergyLimit / eV << " eV - "
00067                << highEnergyLimit / MeV << " MeV"
00068                << G4endl;
00069     }
00070     fParticleChangeForGamma = 0;
00071 }

G4DNAScreenedRutherfordElasticModel::~G4DNAScreenedRutherfordElasticModel (  )  [virtual]

Definition at line 75 of file G4DNAScreenedRutherfordElasticModel.cc.

00076 {}


Member Function Documentation

G4double G4DNAScreenedRutherfordElasticModel::CrossSectionPerVolume ( const G4Material material,
const G4ParticleDefinition p,
G4double  ekin,
G4double  emin,
G4double  emax 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 156 of file G4DNAScreenedRutherfordElasticModel.cc.

References DBL_MAX, G4cout, G4endl, G4Material::GetAtomicNumDensityVector(), G4Material::GetIndex(), G4ParticleDefinition::GetParticleName(), CLHEP::detail::n, and G4INCL::Math::pi.

00161 {
00162     if (verboseLevel > 3)
00163         G4cout << "Calling CrossSectionPerVolume() of G4DNAScreenedRutherfordElasticModel" << G4endl;
00164 
00165     // Calculate total cross section for model
00166 
00167     G4double sigma=0;
00168 
00169     G4double waterDensity = (*fpWaterDensity)[material->GetIndex()];
00170 
00171     if(waterDensity!= 0.0)
00172   //  if (material == nistwater || material->GetBaseMaterial() == nistwater)
00173     {
00174 
00175         if (ekin < highEnergyLimit)
00176         {
00177 
00178             if (ekin < killBelowEnergy) return DBL_MAX;
00179 
00180             G4double z = 10.;
00181             G4double n = ScreeningFactor(ekin,z);
00182             G4double crossSection = RutherfordCrossSection(ekin, z);
00183             sigma = pi *  crossSection / (n * (n + 1.));
00184         }
00185 
00186         if (verboseLevel > 2)
00187         {
00188             G4cout << "__________________________________" << G4endl;
00189             G4cout << "°°° G4DNAScreenedRutherfordElasticModel - XS INFO START" << G4endl;
00190             G4cout << "°°° Kinetic energy(eV)=" << ekin/eV << " particle : " << particleDefinition->GetParticleName() << G4endl;
00191             G4cout << "°°° Cross section per water molecule (cm^2)=" << sigma/cm/cm << G4endl;
00192             G4cout << "°°° Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./cm) << G4endl;
00193             //      G4cout << " - Cross section per water molecule (cm^-1)=" << sigma*material->GetAtomicNumDensityVector()[1]/(1./cm) << G4endl;
00194             G4cout << "°°° G4DNAScreenedRutherfordElasticModel - XS INFO END" << G4endl;
00195         }
00196 
00197     }
00198 
00199     return sigma*material->GetAtomicNumDensityVector()[1];
00200 }

G4double G4DNAScreenedRutherfordElasticModel::GetKillBelowThreshold (  )  [inline]

Definition at line 64 of file G4DNAScreenedRutherfordElasticModel.hh.

00064 { return killBelowEnergy; }              

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

Implements G4VEmModel.

Definition at line 80 of file G4DNAScreenedRutherfordElasticModel.cc.

References fParticleChangeForGamma, G4cout, G4endl, G4Material::GetMaterial(), G4DNAMolecularMaterial::GetNumMolPerVolTableFor(), G4VEmModel::GetParticleChangeForGamma(), G4VEmModel::HighEnergyLimit(), G4DNAMolecularMaterial::Instance(), G4VEmModel::LowEnergyLimit(), G4VEmModel::SetHighEnergyLimit(), and G4VEmModel::SetLowEnergyLimit().

00082 {
00083 
00084     if (verboseLevel > 3)
00085         G4cout << "Calling G4DNAScreenedRutherfordElasticModel::Initialise()" << G4endl;
00086 
00087     // Energy limits
00088 
00089     if (LowEnergyLimit() < lowEnergyLimit)
00090     {
00091         G4cout << "G4DNAScreenedRutherfordElasticModel: low energy limit increased from " <<
00092                   LowEnergyLimit()/eV << " eV to " << lowEnergyLimit/eV << " eV" << G4endl;
00093         SetLowEnergyLimit(lowEnergyLimit);
00094     }
00095 
00096     if (HighEnergyLimit() > highEnergyLimit)
00097     {
00098         G4cout << "G4DNAScreenedRutherfordElasticModel: high energy limit decreased from " <<
00099                   HighEnergyLimit()/MeV << " MeV to " << highEnergyLimit/MeV << " MeV" << G4endl;
00100         SetHighEnergyLimit(highEnergyLimit);
00101     }
00102 
00103     // Constants for final stae by Brenner & Zaider
00104 
00105     betaCoeff.push_back(7.51525);
00106     betaCoeff.push_back(-0.41912);
00107     betaCoeff.push_back(7.2017E-3);
00108     betaCoeff.push_back(-4.646E-5);
00109     betaCoeff.push_back(1.02897E-7);
00110 
00111     deltaCoeff.push_back(2.9612);
00112     deltaCoeff.push_back(-0.26376);
00113     deltaCoeff.push_back(4.307E-3);
00114     deltaCoeff.push_back(-2.6895E-5);
00115     deltaCoeff.push_back(5.83505E-8);
00116 
00117     gamma035_10Coeff.push_back(-1.7013);
00118     gamma035_10Coeff.push_back(-1.48284);
00119     gamma035_10Coeff.push_back(0.6331);
00120     gamma035_10Coeff.push_back(-0.10911);
00121     gamma035_10Coeff.push_back(8.358E-3);
00122     gamma035_10Coeff.push_back(-2.388E-4);
00123 
00124     gamma10_100Coeff.push_back(-3.32517);
00125     gamma10_100Coeff.push_back(0.10996);
00126     gamma10_100Coeff.push_back(-4.5255E-3);
00127     gamma10_100Coeff.push_back(5.8372E-5);
00128     gamma10_100Coeff.push_back(-2.4659E-7);
00129 
00130     gamma100_200Coeff.push_back(2.4775E-2);
00131     gamma100_200Coeff.push_back(-2.96264E-5);
00132     gamma100_200Coeff.push_back(-1.20655E-7);
00133 
00134     //
00135 
00136     if( verboseLevel>0 )
00137     {
00138         G4cout << "Screened Rutherford elastic model is initialized " << G4endl
00139                << "Energy range: "
00140                << LowEnergyLimit() / eV << " eV - "
00141                << HighEnergyLimit() / MeV << " MeV"
00142                << G4endl;
00143     }
00144 
00145     // Initialize water density pointer
00146     fpWaterDensity = G4DNAMolecularMaterial::Instance()->GetNumMolPerVolTableFor(G4Material::GetMaterial("G4_WATER"));
00147 
00148     if (isInitialised) { return; }
00149     fParticleChangeForGamma = GetParticleChangeForGamma();
00150     isInitialised = true;
00151 
00152 }

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

Implements G4VEmModel.

Definition at line 255 of file G4DNAScreenedRutherfordElasticModel.cc.

References fParticleChangeForGamma, fStopAndKill, G4cout, G4endl, G4UniformRand, G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), G4INCL::Math::pi, G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), and G4ParticleChangeForGamma::SetProposedKineticEnergy().

00260 {
00261 
00262     if (verboseLevel > 3)
00263         G4cout << "Calling SampleSecondaries() of G4DNAScreenedRutherfordElasticModel" << G4endl;
00264 
00265     G4double electronEnergy0 = aDynamicElectron->GetKineticEnergy();
00266 
00267     if (electronEnergy0 < killBelowEnergy)
00268     {
00269         fParticleChangeForGamma->SetProposedKineticEnergy(0.);
00270         fParticleChangeForGamma->ProposeTrackStatus(fStopAndKill);
00271         fParticleChangeForGamma->ProposeLocalEnergyDeposit(electronEnergy0);
00272         return ;
00273     }
00274 
00275     G4double cosTheta = 0.;
00276 
00277     if (electronEnergy0>= killBelowEnergy && electronEnergy0 < highEnergyLimit)
00278     {
00279         if (electronEnergy0<intermediateEnergyLimit)
00280         {
00281             if (verboseLevel > 3) G4cout << "---> Using Brenner & Zaider model" << G4endl;
00282             cosTheta = BrennerZaiderRandomizeCosTheta(electronEnergy0);
00283         }
00284 
00285         if (electronEnergy0>=intermediateEnergyLimit)
00286         {
00287             if (verboseLevel > 3) G4cout << "---> Using Screened Rutherford model" << G4endl;
00288             G4double z = 10.;
00289             cosTheta = ScreenedRutherfordRandomizeCosTheta(electronEnergy0,z);
00290         }
00291 
00292         G4double phi = 2. * pi * G4UniformRand();
00293 
00294         G4ThreeVector zVers = aDynamicElectron->GetMomentumDirection();
00295         G4ThreeVector xVers = zVers.orthogonal();
00296         G4ThreeVector yVers = zVers.cross(xVers);
00297 
00298         G4double xDir = std::sqrt(1. - cosTheta*cosTheta);
00299         G4double yDir = xDir;
00300         xDir *= std::cos(phi);
00301         yDir *= std::sin(phi);
00302 
00303         G4ThreeVector zPrimeVers((xDir*xVers + yDir*yVers + cosTheta*zVers));
00304 
00305         fParticleChangeForGamma->ProposeMomentumDirection(zPrimeVers.unit()) ;
00306 
00307         fParticleChangeForGamma->SetProposedKineticEnergy(electronEnergy0);
00308     }
00309 
00310 }

void G4DNAScreenedRutherfordElasticModel::SetKillBelowThreshold ( G4double  threshold  )  [inline]

Definition at line 108 of file G4DNAScreenedRutherfordElasticModel.hh.

References G4Exception(), and JustWarning.

00109 { 
00110     killBelowEnergy = threshold; 
00111     if (threshold < 9*CLHEP::eV)
00112      G4Exception ("*** WARNING : the G4DNAScreenedRutherfordElasticModel class is not validated below 9 eV !","",JustWarning,"") ;   
00113 }                


Field Documentation

G4ParticleChangeForGamma* G4DNAScreenedRutherfordElasticModel::fParticleChangeForGamma [protected]

Definition at line 68 of file G4DNAScreenedRutherfordElasticModel.hh.

Referenced by Initialise(), and SampleSecondaries().


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