G4E1SingleProbability1 Class Reference

#include <G4E1SingleProbability1.hh>

Inheritance diagram for G4E1SingleProbability1:

G4VEmissionProbability

Public Member Functions

 G4E1SingleProbability1 ()
virtual ~G4E1SingleProbability1 ()
G4double EmissionProbability (const G4Fragment &frag, G4double excite)
G4double EmissionProbDensity (const G4Fragment &frag, G4double ePhoton)

Detailed Description

Definition at line 38 of file G4E1SingleProbability1.hh.


Constructor & Destructor Documentation

G4E1SingleProbability1::G4E1SingleProbability1 (  ) 

Definition at line 41 of file G4E1SingleProbability1.cc.

00042 {}

G4E1SingleProbability1::~G4E1SingleProbability1 (  )  [virtual]

Definition at line 44 of file G4E1SingleProbability1.cc.

00045 {}


Member Function Documentation

G4double G4E1SingleProbability1::EmissionProbability ( const G4Fragment frag,
G4double  excite 
) [virtual]

Implements G4VEmissionProbability.

Definition at line 121 of file G4E1SingleProbability1.cc.

References G4Fragment::GetExcitationEnergy().

00123 {
00124 
00125   // From nuclear fragment properties and the excitation energy, calculate
00126   // the probability for photon evaporation down to the level
00127   // Uexcite-exciteE.
00128   // fragment = nuclear fragment BEFORE de-excitation
00129 
00130   G4double theProb = 0.0;
00131 
00132   G4double ScaleFactor = 1.0;     // playing with scale factors
00133 
00134   const G4double Uexcite = frag.GetExcitationEnergy();
00135   G4double Uafter = Uexcite - exciteE;
00136 
00137   G4double normC = 3.0;
00138 
00139   const G4double upperLim = Uexcite;
00140   const G4double lowerLim = Uafter;
00141   const G4int numIters = 25;
00142 
00143   // Need to integrate EmissionProbDensity from lowerLim to upperLim 
00144   // and multiply by normC
00145 
00146   G4double integ = normC *
00147            EmissionIntegration(frag,exciteE,lowerLim,upperLim,numIters);
00148 
00149   if(integ > 0.0) theProb = integ;
00150 
00151   return theProb * ScaleFactor;
00152 
00153 }

G4double G4E1SingleProbability1::EmissionProbDensity ( const G4Fragment frag,
G4double  ePhoton 
)

Definition at line 50 of file G4E1SingleProbability1.cc.

References G4Fragment::GetA_asInt(), G4Fragment::GetExcitationEnergy(), G4Pow::GetInstance(), G4Fragment::GetZ_asInt(), G4ConstantLevelDensityParameter::LevelDensityParameter(), G4INCL::Math::pi, and G4Pow::powZ().

00052 {
00053 
00054   // Calculate the probability density here
00055 
00056   // From nuclear fragment properties and the excitation energy, calculate
00057   // the probability density for photon evaporation from U to U - exciteE
00058   // (U = nucleus excitation energy, exciteE = total evaporated photon
00059   // energy).
00060   // fragment = nuclear fragment BEFORE de-excitation
00061 
00062   G4double theProb = 0.0;
00063 
00064   G4int Afrag = frag.GetA_asInt();
00065   G4int Zfrag = frag.GetZ_asInt();
00066   G4double Uexcite = frag.GetExcitationEnergy();
00067 
00068   if( (Uexcite-exciteE) < 0.0 || exciteE < 0 || Uexcite <= 0) return theProb;
00069 
00070   // Need a level density parameter.
00071   // For now, just use the constant approximation (not reliable near magic
00072   // nuclei).
00073 
00074   G4ConstantLevelDensityParameter a;
00075   G4double aLevelDensityParam = a.LevelDensityParameter(Afrag,Zfrag,Uexcite);
00076 
00077   G4double levelDensBef = std::exp(2.0*std::sqrt(aLevelDensityParam*Uexcite));
00078   G4double levelDensAft = std::exp(2.0*std::sqrt(aLevelDensityParam*(Uexcite-exciteE)));
00079 
00080   // Now form the probability density
00081 
00082   // Define constants for the photoabsorption cross-section (the reverse
00083   // process of our de-excitation)
00084 
00085   G4double sigma0 = 2.5 * Afrag * millibarn;  // millibarns
00086 
00087   G4double Egdp = (40.3 / G4Pow::GetInstance()->powZ(Afrag,0.2) )*MeV;
00088   G4double GammaR = 0.30 * Egdp;
00089  
00090   const G4double normC = 1.0 / ((pi * hbarc)*(pi * hbarc));
00091 
00092   // CD
00093   //cout<<"  PROB TESTS "<<G4endl;
00094   //cout<<" hbarc = "<<hbarc<<G4endl;
00095   //cout<<" pi = "<<pi<<G4endl;
00096   //cout<<" Uexcite, exciteE = "<<Uexcite<<"  "<<exciteE<<G4endl;
00097   //cout<<" Uexcite, exciteE = "<<Uexcite*MeV<<"  "<<exciteE*MeV<<G4endl;
00098   //cout<<" lev density param = "<<aLevelDensityParam<<G4endl;
00099   //cout<<" level densities = "<<levelDensBef<<"  "<<levelDensAft<<G4endl;
00100   //cout<<" sigma0 = "<<sigma0<<G4endl;
00101   //cout<<" Egdp, GammaR = "<<Egdp<<"  "<<GammaR<<G4endl;
00102   //cout<<" normC = "<<normC<<G4endl;
00103 
00104   G4double numerator = sigma0 * exciteE*exciteE * GammaR*GammaR;
00105   G4double denominator = (exciteE*exciteE - Egdp*Egdp)*
00106            (exciteE*exciteE - Egdp*Egdp) + GammaR*GammaR*exciteE*exciteE;
00107 
00108   G4double sigmaAbs = numerator/denominator; 
00109 
00110   theProb = normC * sigmaAbs * exciteE*exciteE *
00111             levelDensAft/levelDensBef;
00112 
00113   // CD
00114   //cout<<" sigmaAbs = "<<sigmaAbs<<G4endl;
00115   //cout<<" Probability = "<<theProb<<G4endl;
00116 
00117   return theProb;
00118 
00119 }


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