G4LivermorePolarizedGammaConversionModel Class Reference

#include <G4LivermorePolarizedGammaConversionModel.hh>

Inheritance diagram for G4LivermorePolarizedGammaConversionModel:

G4VEmModel

Public Member Functions

 G4LivermorePolarizedGammaConversionModel (const G4ParticleDefinition *p=0, const G4String &nam="LivermorePolarizedGammaConversion")
virtual ~G4LivermorePolarizedGammaConversionModel ()
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 Member Functions

G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)

Protected Attributes

G4ParticleChangeForGammafParticleChange

Detailed Description

Definition at line 45 of file G4LivermorePolarizedGammaConversionModel.hh.


Constructor & Destructor Documentation

G4LivermorePolarizedGammaConversionModel::G4LivermorePolarizedGammaConversionModel ( const G4ParticleDefinition p = 0,
const G4String nam = "LivermorePolarizedGammaConversion" 
)

Definition at line 43 of file G4LivermorePolarizedGammaConversionModel.cc.

References G4cout, G4endl, G4VCrossSectionHandler::Initialise(), G4VEmModel::SetHighEnergyLimit(), and G4VEmModel::SetLowEnergyLimit().

00045   :G4VEmModel(nam),fParticleChange(0),
00046    isInitialised(false),meanFreePathTable(0),crossSectionHandler(0)
00047 {
00048   lowEnergyLimit = 2*electron_mass_c2;
00049   highEnergyLimit = 100 * GeV;
00050   SetLowEnergyLimit(lowEnergyLimit);
00051   SetHighEnergyLimit(highEnergyLimit);
00052   smallEnergy = 2.*MeV;
00053 
00054   Phi=0.;
00055   Psi=0.;
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, samping of atoms
00063   // 4 = entering in methods
00064 
00065   if(verboseLevel > 0) {
00066     G4cout << "Livermore Polarized GammaConversion is constructed " << G4endl
00067            << "Energy range: "
00068            << lowEnergyLimit / keV << " keV - "
00069            << highEnergyLimit / GeV << " GeV"
00070            << G4endl;
00071   }
00072 
00073   crossSectionHandler = new G4CrossSectionHandler();
00074   crossSectionHandler->Initialise(0,lowEnergyLimit,highEnergyLimit,400);
00075 }

G4LivermorePolarizedGammaConversionModel::~G4LivermorePolarizedGammaConversionModel (  )  [virtual]

Definition at line 79 of file G4LivermorePolarizedGammaConversionModel.cc.

00080 {  
00081   delete crossSectionHandler;
00082 }


Member Function Documentation

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

Reimplemented from G4VEmModel.

Definition at line 149 of file G4LivermorePolarizedGammaConversionModel.cc.

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

00154 {
00155   if (verboseLevel > 3) {
00156     G4cout << "G4LivermorePolarizedGammaConversionModel::ComputeCrossSectionPerAtom()" 
00157            << G4endl;
00158   }
00159   if(Z < 0.9 || GammaEnergy <= lowEnergyLimit) { return 0.0; }
00160   G4double cs = crossSectionHandler->FindValue(G4int(Z), GammaEnergy);
00161   return cs;
00162 }

G4double G4LivermorePolarizedGammaConversionModel::GetMeanFreePath ( const G4Track aTrack,
G4double  previousStepSize,
G4ForceCondition condition 
) [protected]

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

Implements G4VEmModel.

Definition at line 88 of file G4LivermorePolarizedGammaConversionModel.cc.

References G4VCrossSectionHandler::Clear(), fParticleChange, G4cout, G4endl, G4VEmModel::GetParticleChangeForGamma(), G4VEmModel::HighEnergyLimit(), G4VEmModel::InitialiseElementSelectors(), G4VCrossSectionHandler::LoadData(), and G4VEmModel::LowEnergyLimit().

00090 {
00091   if (verboseLevel > 3)
00092     G4cout << "Calling G4LivermorePolarizedGammaConversionModel::Initialise()" << G4endl;
00093 
00094   if (crossSectionHandler)
00095   {
00096     crossSectionHandler->Clear();
00097     delete crossSectionHandler;
00098   }
00099 
00100   // Energy limits
00101   /*
00102   // V.Ivanchenko: this was meanless check  
00103   if (LowEnergyLimit() < lowEnergyLimit)
00104     {
00105       G4cout << "G4LivermorePolarizedGammaConversionModel: low energy limit increased from " << 
00106         LowEnergyLimit()/eV << " eV to " << lowEnergyLimit << " eV" << G4endl;
00107       //      SetLowEnergyLimit(lowEnergyLimit);
00108     }
00109   */
00110   if (HighEnergyLimit() > highEnergyLimit)
00111     {
00112       G4cout << "G4LivermorePolarizedGammaConversionModel: high energy limit decreased from " << 
00113         HighEnergyLimit()/GeV << " GeV to " << highEnergyLimit << " GeV" << G4endl;
00114       // V.Ivanchenko: this is forbidden 
00115      // SetHighEnergyLimit(highEnergyLimit);
00116     }
00117 
00118   // Reading of data files - all materials are read
00119   
00120   crossSectionHandler = new G4CrossSectionHandler;
00121   crossSectionHandler->Clear();
00122   G4String crossSectionFile = "pair/pp-cs-";
00123   crossSectionHandler->LoadData(crossSectionFile);
00124 
00125   //
00126   if (verboseLevel > 2) {
00127     G4cout << "Loaded cross section files for Livermore Polarized GammaConversion model" 
00128            << G4endl;
00129   }
00130   InitialiseElementSelectors(particle,cuts);
00131 
00132   if(verboseLevel > 0) {
00133     G4cout << "Livermore Polarized GammaConversion model is initialized " << G4endl
00134            << "Energy range: "
00135            << LowEnergyLimit() / keV << " keV - "
00136            << HighEnergyLimit() / GeV << " GeV"
00137            << G4endl;
00138   }
00139 
00140   //    
00141   if(!isInitialised) { 
00142     isInitialised = true;
00143     fParticleChange = GetParticleChangeForGamma();
00144   }
00145 }

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

Implements G4VEmModel.

Definition at line 167 of file G4LivermorePolarizedGammaConversionModel.cc.

References G4Electron::Electron(), fParticleChange, fStopAndKill, G4cout, G4endl, G4UniformRand, G4DynamicParticle::GetDefinition(), G4Element::GetfCoulomb(), G4Element::GetIonisation(), G4DynamicParticle::GetKineticEnergy(), G4IonisParamElm::GetlogZ3(), G4DynamicParticle::GetMomentumDirection(), G4DynamicParticle::GetPolarization(), G4IonisParamElm::GetZ3(), G4Positron::Positron(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), G4VEmModel::SelectRandomAtom(), and G4ParticleChangeForGamma::SetProposedKineticEnergy().

00172 {
00173 
00174   // Fluorescence generated according to:
00175   // J. Stepanek ,"A program to determine the radiation spectra due to a single atomic
00176   // subshell ionisation by a particle or due to deexcitation or decay of radionuclides",
00177   // Comp. Phys. Comm. 1206 pp 1-1-9 (1997)
00178 
00179   if (verboseLevel > 3)
00180     G4cout << "Calling SampleSecondaries() of G4LivermorePolarizedGammaConversionModel" << G4endl;
00181 
00182   G4double photonEnergy = aDynamicGamma->GetKineticEnergy();
00183   // Within energy limit?
00184 
00185   if(photonEnergy <= lowEnergyLimit)
00186     {
00187       fParticleChange->ProposeTrackStatus(fStopAndKill);
00188       fParticleChange->SetProposedKineticEnergy(0.);
00189       return;
00190     }
00191 
00192 
00193   G4ThreeVector gammaPolarization0 = aDynamicGamma->GetPolarization();
00194   G4ThreeVector gammaDirection0 = aDynamicGamma->GetMomentumDirection();
00195 
00196   // Make sure that the polarization vector is perpendicular to the
00197   // gamma direction. If not
00198 
00199   if(!(gammaPolarization0.isOrthogonal(gammaDirection0, 1e-6))||(gammaPolarization0.mag()==0))
00200     { // only for testing now
00201       gammaPolarization0 = GetRandomPolarization(gammaDirection0);
00202     }
00203   else
00204     {
00205       if ( gammaPolarization0.howOrthogonal(gammaDirection0) != 0)
00206         {
00207           gammaPolarization0 = GetPerpendicularPolarization(gammaDirection0, gammaPolarization0);
00208         }
00209     }
00210 
00211   // End of Protection
00212 
00213 
00214   G4double epsilon ;
00215   G4double epsilon0Local = electron_mass_c2 / photonEnergy ;
00216 
00217   // Do it fast if photon energy < 2. MeV
00218 
00219   if (photonEnergy < smallEnergy )
00220     {
00221       epsilon = epsilon0Local + (0.5 - epsilon0Local) * G4UniformRand();
00222     }
00223   else
00224     {
00225 
00226  // Select randomly one element in the current material
00227 
00228       //     G4int Z = crossSectionHandler->SelectRandomAtom(couple,photonEnergy);
00229       //const G4Element* element = crossSectionHandler->SelectRandomElement(couple,photonEnergy);
00230 
00231       const G4ParticleDefinition* particle =  aDynamicGamma->GetDefinition();
00232       const G4Element* element = SelectRandomAtom(couple,particle,photonEnergy);
00233 
00234       /*
00235       if (element == 0)
00236         {
00237           G4cout << "G4LivermorePolarizedGammaConversionModel::PostStepDoIt - element = 0" << G4endl;
00238         }
00239       */
00240       
00241       G4IonisParamElm* ionisation = element->GetIonisation();
00242       
00243       /*
00244       if (ionisation == 0)
00245         {
00246           G4cout << "G4LivermorePolarizedGammaConversionModel::PostStepDoIt - ionisation = 0" << G4endl;
00247         }
00248       */
00249 
00250 
00251       // Extract Coulomb factor for this Element
00252 
00253       G4double fZ = 8. * (ionisation->GetlogZ3());
00254       if (photonEnergy > 50. * MeV) fZ += 8. * (element->GetfCoulomb());
00255 
00256       // Limits of the screening variable
00257       G4double screenFactor = 136. * epsilon0Local / (element->GetIonisation()->GetZ3()) ;
00258       G4double screenMax = exp ((42.24 - fZ)/8.368) - 0.952 ;
00259       G4double screenMin = std::min(4.*screenFactor,screenMax) ;
00260 
00261       // Limits of the energy sampling
00262       G4double epsilon1 = 0.5 - 0.5 * sqrt(1. - screenMin / screenMax) ;
00263       G4double epsilonMin = std::max(epsilon0Local,epsilon1);
00264       G4double epsilonRange = 0.5 - epsilonMin ;
00265 
00266       // Sample the energy rate of the created electron (or positron)
00267       G4double screen;
00268       G4double gReject ;
00269 
00270       G4double f10 = ScreenFunction1(screenMin) - fZ;
00271       G4double f20 = ScreenFunction2(screenMin) - fZ;
00272       G4double normF1 = std::max(f10 * epsilonRange * epsilonRange,0.);
00273       G4double normF2 = std::max(1.5 * f20,0.);
00274 
00275       do {
00276         if (normF1 / (normF1 + normF2) > G4UniformRand() )
00277           {
00278             epsilon = 0.5 - epsilonRange * pow(G4UniformRand(), 0.3333) ;
00279             screen = screenFactor / (epsilon * (1. - epsilon));
00280             gReject = (ScreenFunction1(screen) - fZ) / f10 ;
00281           }
00282         else
00283           {
00284             epsilon = epsilonMin + epsilonRange * G4UniformRand();
00285             screen = screenFactor / (epsilon * (1 - epsilon));
00286             gReject = (ScreenFunction2(screen) - fZ) / f20 ;
00287 
00288 
00289           }
00290       } while ( gReject < G4UniformRand() );
00291 
00292     }   //  End of epsilon sampling
00293   
00294   // Fix charges randomly
00295   
00296   G4double electronTotEnergy;
00297   G4double positronTotEnergy;
00298 
00299 
00300   if (G4int(2*G4UniformRand()))  
00301     {
00302       electronTotEnergy = (1. - epsilon) * photonEnergy;
00303       positronTotEnergy = epsilon * photonEnergy;
00304     }
00305   else
00306     {
00307       positronTotEnergy = (1. - epsilon) * photonEnergy;
00308       electronTotEnergy = epsilon * photonEnergy;
00309     }
00310 
00311   // Scattered electron (positron) angles. ( Z - axis along the parent photon)
00312   // Universal distribution suggested by L. Urban (Geant3 manual (1993) Phys211),
00313   // derived from Tsai distribution (Rev. Mod. Phys. 49, 421 (1977)
00314 
00315 /*
00316   G4double u;
00317   const G4double a1 = 0.625;
00318   G4double a2 = 3. * a1;
00319 
00320   if (0.25 > G4UniformRand())
00321     {
00322       u = - log(G4UniformRand() * G4UniformRand()) / a1 ;
00323     }
00324   else
00325     {
00326       u = - log(G4UniformRand() * G4UniformRand()) / a2 ;
00327     }
00328 */
00329 
00330   G4double Ene = electronTotEnergy/electron_mass_c2; // Normalized energy
00331 
00332   G4double cosTheta = 0.;
00333   G4double sinTheta = 0.;
00334 
00335   SetTheta(&cosTheta,&sinTheta,Ene);
00336 
00337   //  G4double theta = u * electron_mass_c2 / photonEnergy ;
00338   //  G4double phi  = twopi * G4UniformRand() ;
00339 
00340   G4double phi,psi=0.;
00341 
00342   //corrected e+ e- angular angular distribution //preliminary!
00343 
00344   //  if(photonEnergy>50*MeV)
00345   // {
00346   phi = SetPhi(photonEnergy);
00347   psi = SetPsi(photonEnergy,phi);
00348   //  }
00349   //else
00350   // {
00351   //psi = G4UniformRand()*2.*pi;
00352   //phi = pi; // coplanar
00353   // }
00354 
00355   Psi = psi;
00356   Phi = phi;
00357   //G4cout << "PHI " << phi << G4endl;
00358   //G4cout << "PSI " << psi << G4endl;
00359 
00360   G4double phie = psi; //azimuthal angle for the electron
00361 
00362   G4double dirX = sinTheta*cos(phie);
00363   G4double dirY = sinTheta*sin(phie);
00364   G4double dirZ = cosTheta;
00365   G4ThreeVector electronDirection(dirX,dirY,dirZ);
00366   // Kinematics of the created pair:
00367   // the electron and positron are assumed to have a symetric angular
00368   // distribution with respect to the Z axis along the parent photon
00369 
00370   //G4double localEnergyDeposit = 0. ;
00371 
00372   G4double electronKineEnergy = std::max(0.,electronTotEnergy - electron_mass_c2) ;
00373 
00374   SystemOfRefChange(gammaDirection0,electronDirection,
00375                     gammaPolarization0);
00376 
00377   G4DynamicParticle* particle1 = new G4DynamicParticle (G4Electron::Electron(),
00378                                                         electronDirection,
00379                                                         electronKineEnergy);
00380 
00381   // The e+ is always created (even with kinetic energy = 0) for further annihilation
00382 
00383   Ene = positronTotEnergy/electron_mass_c2; // Normalized energy
00384 
00385   cosTheta = 0.;
00386   sinTheta = 0.;
00387 
00388   SetTheta(&cosTheta,&sinTheta,Ene);
00389   G4double phip = phie+phi; //azimuthal angle for the positron
00390 
00391   dirX = sinTheta*cos(phip);
00392   dirY = sinTheta*sin(phip);
00393   dirZ = cosTheta;
00394   G4ThreeVector positronDirection(dirX,dirY,dirZ);
00395 
00396   G4double positronKineEnergy = std::max(0.,positronTotEnergy - electron_mass_c2) ;
00397   SystemOfRefChange(gammaDirection0,positronDirection,
00398                     gammaPolarization0);
00399 
00400   // Create G4DynamicParticle object for the particle2
00401   G4DynamicParticle* particle2 = new G4DynamicParticle(G4Positron::Positron(),
00402                                                        positronDirection, positronKineEnergy);
00403 
00404 
00405   fvect->push_back(particle1);
00406   fvect->push_back(particle2);
00407 
00408 
00409 
00410   // Kill the incident photon
00411 
00412 
00413 
00414   // Create lists of pointers to DynamicParticles (photons and electrons)
00415   // (Is the electron vector necessary? To be checked)
00416   //  std::vector<G4DynamicParticle*>* photonVector = 0;
00417   //std::vector<G4DynamicParticle*> electronVector;
00418 
00419   fParticleChange->ProposeMomentumDirection( 0., 0., 0. );
00420   fParticleChange->SetProposedKineticEnergy(0.);
00421   fParticleChange->ProposeTrackStatus(fStopAndKill);
00422 
00423 }


Field Documentation

G4ParticleChangeForGamma* G4LivermorePolarizedGammaConversionModel::fParticleChange [protected]

Definition at line 73 of file G4LivermorePolarizedGammaConversionModel.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