G4PenelopeBremsstrahlungModel Class Reference

#include <G4PenelopeBremsstrahlungModel.hh>

Inheritance diagram for G4PenelopeBremsstrahlungModel:

G4VEmModel

Public Member Functions

 G4PenelopeBremsstrahlungModel (const G4ParticleDefinition *p=0, const G4String &processName="PenBrem")
virtual ~G4PenelopeBremsstrahlungModel ()
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *theParticle, G4double kinEnergy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX)
virtual G4double CrossSectionPerVolume (const G4Material *material, const G4ParticleDefinition *theParticle, G4double kineticEnergy, G4double cutEnergy, G4double maxEnergy=DBL_MAX)
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy)
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
void SetVerbosityLevel (G4int lev)
G4int GetVerbosityLevel ()

Protected Attributes

G4ParticleChangeForLoss * fParticleChange

Detailed Description

Definition at line 63 of file G4PenelopeBremsstrahlungModel.hh.

Constructor & Destructor Documentation

G4PenelopeBremsstrahlungModel::G4PenelopeBremsstrahlungModel ( const G4ParticleDefinition *  p = 0,
const G4String &  processName = "PenBrem" 
)

Definition at line 59 of file G4PenelopeBremsstrahlungModel.cc.

References G4PenelopeOscillatorManager::GetOscillatorManager(), G4VEmModel::SetDeexcitationFlag(), and G4VEmModel::SetHighEnergyLimit().

00061   :G4VEmModel(nam),fParticleChange(0),isInitialised(false),energyGrid(0),  
00062    XSTableElectron(0),XSTablePositron(0)
00063   
00064 {
00065   fIntrinsicLowEnergyLimit = 100.0*eV;
00066   fIntrinsicHighEnergyLimit = 100.0*GeV;
00067   nBins = 200;
00068 
00069   SetHighEnergyLimit(fIntrinsicHighEnergyLimit);
00070   //
00071   oscManager = G4PenelopeOscillatorManager::GetOscillatorManager();
00072   //
00073   verboseLevel= 0;
00074    
00075   // Verbosity scale:
00076   // 0 = nothing
00077   // 1 = warning for energy non-conservation
00078   // 2 = details of energy budget
00079   // 3 = calculation of cross sections, file openings, sampling of atoms
00080   // 4 = entering in methods
00081 
00082   fPenelopeFSHelper = new G4PenelopeBremsstrahlungFS();
00083   fPenelopeAngular = new G4PenelopeBremsstrahlungAngular();
00084 
00085   // Atomic deexcitation model activated by default
00086   SetDeexcitationFlag(true);
00087 
00088 }

G4PenelopeBremsstrahlungModel::~G4PenelopeBremsstrahlungModel (  )  [virtual]

Definition at line 92 of file G4PenelopeBremsstrahlungModel.cc.

00093 {
00094   ClearTables();
00095   if (fPenelopeFSHelper)
00096     delete fPenelopeFSHelper;
00097   if (fPenelopeAngular)
00098     delete fPenelopeAngular;
00099 }

Member Function Documentation

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

Reimplemented from G4VEmModel.

Definition at line 193 of file G4PenelopeBremsstrahlungModel.cc.

References G4cout, and G4endl.

00199 {
00200   G4cout << "*** G4PenelopeBremsstrahlungModel -- WARNING ***" << G4endl;
00201   G4cout << "Penelope Bremsstrahlung model v2008 does not calculate cross section _per atom_ " << G4endl;
00202   G4cout << "so the result is always zero. For physics values, please invoke " << G4endl;
00203   G4cout << "GetCrossSectionPerVolume() or GetMeanFreePath() via the G4EmCalculator" << G4endl;
00204   return 0;
00205 }

G4double G4PenelopeBremsstrahlungModel::ComputeDEDXPerVolume ( const G4Material *  ,
const G4ParticleDefinition *  ,
G4double  kineticEnergy,
G4double  cutEnergy 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 209 of file G4PenelopeBremsstrahlungModel.cc.

References G4cout, G4endl, G4PenelopeOscillatorManager::GetAtomsPerMolecule(), G4PenelopeCrossSection::GetSoftStoppingPower(), and G4Material::GetTotNbOfAtomsPerVolume().

00213 {
00214   if (verboseLevel > 3)
00215     G4cout << "Calling ComputeDEDX() of G4PenelopeBremsstrahlungModel" << G4endl;
00216  
00217   G4PenelopeCrossSection* theXS = GetCrossSectionTableForCouple(theParticle,material,
00218                                                                 cutEnergy);
00219   G4double sPowerPerMolecule = 0.0;
00220   if (theXS)
00221     sPowerPerMolecule = theXS->GetSoftStoppingPower(kineticEnergy);
00222 
00223   
00224   G4double atomDensity = material->GetTotNbOfAtomsPerVolume();
00225   G4double atPerMol =  oscManager->GetAtomsPerMolecule(material);
00226                                                                                         
00227   G4double moleculeDensity = 0.; 
00228   if (atPerMol)
00229     moleculeDensity = atomDensity/atPerMol;
00230  
00231   G4double sPowerPerVolume = sPowerPerMolecule*moleculeDensity;
00232 
00233   if (verboseLevel > 2)
00234     {
00235       G4cout << "G4PenelopeBremsstrahlungModel " << G4endl;
00236       G4cout << "Stopping power < " << cutEnergy/keV << " keV at " <<
00237         kineticEnergy/keV << " keV = " << 
00238         sPowerPerVolume/(keV/mm) << " keV/mm" << G4endl;
00239     }
00240   return sPowerPerVolume;
00241 }

G4double G4PenelopeBremsstrahlungModel::CrossSectionPerVolume ( const G4Material *  material,
const G4ParticleDefinition *  theParticle,
G4double  kineticEnergy,
G4double  cutEnergy,
G4double  maxEnergy = DBL_MAX 
) [virtual]

Reimplemented from G4VEmModel.

Definition at line 145 of file G4PenelopeBremsstrahlungModel.cc.

References G4cout, G4endl, G4PenelopeOscillatorManager::GetAtomsPerMolecule(), G4PenelopeCrossSection::GetHardCrossSection(), G4Material::GetName(), G4Material::GetTotNbOfAtomsPerVolume(), and G4VEmModel::SetupForMaterial().

00150 {
00151   //
00152   if (verboseLevel > 3)
00153     G4cout << "Calling CrossSectionPerVolume() of G4PenelopeBremsstrahlungModel" << G4endl;
00154  
00155   SetupForMaterial(theParticle, material, energy);
00156 
00157   G4double crossPerMolecule = 0.;
00158 
00159   G4PenelopeCrossSection* theXS = GetCrossSectionTableForCouple(theParticle,material,
00160                                                                 cutEnergy);
00161 
00162   if (theXS)
00163     crossPerMolecule = theXS->GetHardCrossSection(energy);
00164 
00165   G4double atomDensity = material->GetTotNbOfAtomsPerVolume();
00166   G4double atPerMol =  oscManager->GetAtomsPerMolecule(material);
00167  
00168   if (verboseLevel > 3)
00169     G4cout << "Material " << material->GetName() << " has " << atPerMol <<
00170       "atoms per molecule" << G4endl;
00171 
00172   G4double moleculeDensity = 0.; 
00173   if (atPerMol)
00174     moleculeDensity = atomDensity/atPerMol;
00175  
00176   G4double crossPerVolume = crossPerMolecule*moleculeDensity;
00177 
00178   if (verboseLevel > 2)
00179   {
00180     G4cout << "G4PenelopeBremsstrahlungModel " << G4endl;
00181     G4cout << "Mean free path for gamma emission > " << cutEnergy/keV << " keV at " <<
00182       energy/keV << " keV = " << (1./crossPerVolume)/mm << " mm" << G4endl;
00183   }
00184   
00185   return crossPerVolume;
00186 }

G4int G4PenelopeBremsstrahlungModel::GetVerbosityLevel (  )  [inline]

Definition at line 107 of file G4PenelopeBremsstrahlungModel.hh.

00107 {return verboseLevel;};

void G4PenelopeBremsstrahlungModel::Initialise ( const G4ParticleDefinition *  ,
const G4DataVector &   
) [virtual]

Implements G4VEmModel.

Definition at line 103 of file G4PenelopeBremsstrahlungModel.cc.

References fParticleChange, G4cout, G4endl, G4VEmModel::GetParticleChangeForLoss(), G4VEmModel::HighEnergyLimit(), G4PenelopeBremsstrahlungAngular::Initialize(), and G4VEmModel::LowEnergyLimit().

00105 {
00106   if (verboseLevel > 3)
00107     G4cout << "Calling G4PenelopeBremsstrahlungModel::Initialise()" << G4endl;
00108  
00109   //Clear and re-build the tables 
00110   ClearTables();
00111   
00112   //forces the cleaning of tables, in this specific case
00113   if (fPenelopeAngular)
00114     fPenelopeAngular->Initialize();
00115     
00116 
00117   //Set the number of bins for the tables. 20 points per decade
00118   nBins = (size_t) (20*std::log10(HighEnergyLimit()/LowEnergyLimit()));
00119   nBins = std::max(nBins,(size_t)100);
00120   energyGrid = new G4PhysicsLogVector(LowEnergyLimit(),
00121                                       HighEnergyLimit(), 
00122                                       nBins-1); //one hidden bin is added
00123  
00124 
00125   XSTableElectron = new 
00126     std::map< std::pair<const G4Material*,G4double>, G4PenelopeCrossSection*>;
00127   XSTablePositron = new 
00128     std::map< std::pair<const G4Material*,G4double>, G4PenelopeCrossSection*>;
00129 
00130   if (verboseLevel > 2) {
00131     G4cout << "Penelope Bremsstrahlung model v2008 is initialized " << G4endl
00132            << "Energy range: "
00133            << LowEnergyLimit() / keV << " keV - "
00134            << HighEnergyLimit() / GeV << " GeV." 
00135            << G4endl;
00136   }
00137  
00138   if(isInitialised) return;
00139   fParticleChange = GetParticleChangeForLoss();
00140   isInitialised = true;
00141 }

G4double G4PenelopeBremsstrahlungModel::MinEnergyCut ( const G4ParticleDefinition *  ,
const G4MaterialCutsCouple *   
) [virtual]

Definition at line 386 of file G4PenelopeBremsstrahlungModel.cc.

00388 {
00389   return fIntrinsicLowEnergyLimit;
00390 }

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

Implements G4VEmModel.

Definition at line 245 of file G4PenelopeBremsstrahlungModel.cc.

References fParticleChange, G4cout, G4endl, G4Gamma::Gamma(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentum(), G4DynamicParticle::GetMomentumDirection(), G4Material::GetName(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForLoss::ProposeMomentumDirection(), G4PenelopeBremsstrahlungAngular::SampleDirection(), G4PenelopeBremsstrahlungFS::SampleGammaEnergy(), and G4ParticleChangeForLoss::SetProposedKineticEnergy().

00250 {
00251   if (verboseLevel > 3)
00252     G4cout << "Calling SampleSecondaries() of G4PenelopeBremsstrahlungModel" << G4endl;
00253 
00254   G4double kineticEnergy = aDynamicParticle->GetKineticEnergy();
00255   const G4Material* material = couple->GetMaterial();
00256 
00257   if (kineticEnergy <= fIntrinsicLowEnergyLimit)
00258    {
00259      fParticleChange->SetProposedKineticEnergy(0.);
00260      fParticleChange->ProposeLocalEnergyDeposit(kineticEnergy);
00261      return ;
00262    }
00263 
00264   G4ParticleMomentum particleDirection0 = aDynamicParticle->GetMomentumDirection();
00265   //This is the momentum
00266   G4ThreeVector initialMomentum =  aDynamicParticle->GetMomentum();
00267 
00268   //Not enough energy to produce a secondary! Return with nothing happened
00269   if (kineticEnergy < cutG)
00270     return;
00271 
00272   if (verboseLevel > 3)
00273     G4cout << "Going to sample gamma energy for: " <<material->GetName() << " " << 
00274       "energy = " << kineticEnergy/keV << ", cut = " << cutG/keV << G4endl;
00275 
00276    //Sample gamma's energy according to the spectrum
00277   G4double gammaEnergy = 
00278     fPenelopeFSHelper->SampleGammaEnergy(kineticEnergy,material,cutG);
00279 
00280   if (verboseLevel > 3)
00281     G4cout << "Sampled gamma energy: " << gammaEnergy/keV << " keV" << G4endl;
00282   
00283   //Now sample the direction for the Gamma. Notice that the rotation is already done
00284   //Z is unused here, I plug 0. The information is in the material pointer
00285    G4ThreeVector gammaDirection1 = 
00286      fPenelopeAngular->SampleDirection(aDynamicParticle,gammaEnergy,0,material);
00287   
00288   if (verboseLevel > 3)
00289     G4cout << "Sampled cosTheta for e-: " << gammaDirection1.cosTheta() << G4endl;
00290 
00291   G4double residualPrimaryEnergy = kineticEnergy-gammaEnergy;
00292   if (residualPrimaryEnergy < 0)
00293     {
00294       //Ok we have a problem, all energy goes with the gamma
00295       gammaEnergy += residualPrimaryEnergy;
00296       residualPrimaryEnergy = 0.0;
00297     }
00298   
00299   //Produce final state according to momentum conservation
00300   G4ThreeVector particleDirection1 = initialMomentum - gammaEnergy*gammaDirection1;
00301   particleDirection1 = particleDirection1.unit(); //normalize
00302 
00303   //Update the primary particle
00304   if (residualPrimaryEnergy > 0.)
00305     {
00306       fParticleChange->ProposeMomentumDirection(particleDirection1);
00307       fParticleChange->SetProposedKineticEnergy(residualPrimaryEnergy);
00308     }
00309   else
00310     {
00311       fParticleChange->SetProposedKineticEnergy(0.);
00312     }
00313   
00314   //Now produce the photon
00315   G4DynamicParticle* theGamma = new G4DynamicParticle(G4Gamma::Gamma(),
00316                                                       gammaDirection1,
00317                                                       gammaEnergy);
00318   fvect->push_back(theGamma);
00319   
00320   if (verboseLevel > 1)
00321     {
00322       G4cout << "-----------------------------------------------------------" << G4endl;
00323       G4cout << "Energy balance from G4PenelopeBremsstrahlung" << G4endl;
00324       G4cout << "Incoming primary energy: " << kineticEnergy/keV << " keV" << G4endl;
00325       G4cout << "-----------------------------------------------------------" << G4endl;
00326       G4cout << "Outgoing primary energy: " << residualPrimaryEnergy/keV << " keV" << G4endl;
00327       G4cout << "Bremsstrahlung photon " << gammaEnergy/keV << " keV" << G4endl;
00328       G4cout << "Total final state: " << (residualPrimaryEnergy+gammaEnergy)/keV 
00329              << " keV" << G4endl;
00330       G4cout << "-----------------------------------------------------------" << G4endl;
00331     }
00332 
00333   if (verboseLevel > 0)
00334     {
00335       G4double energyDiff = std::fabs(residualPrimaryEnergy+gammaEnergy-kineticEnergy);
00336       if (energyDiff > 0.05*keV)
00337         G4cout << "Warning from G4PenelopeBremsstrahlung: problem with energy conservation: " 
00338                <<
00339           (residualPrimaryEnergy+gammaEnergy)/keV <<
00340           " keV (final) vs. " <<
00341           kineticEnergy/keV << " keV (initial)" << G4endl;
00342     }  
00343   return;
00344 }

void G4PenelopeBremsstrahlungModel::SetVerbosityLevel ( G4int  lev  )  [inline]

Definition at line 106 of file G4PenelopeBremsstrahlungModel.hh.

00106 {verboseLevel = lev;};

Field Documentation

G4ParticleChangeForLoss* G4PenelopeBremsstrahlungModel::fParticleChange [protected]

Definition at line 107 of file G4PenelopeBremsstrahlungModel.hh.

Referenced by Initialise(), and SampleSecondaries().

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