Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4PreCompoundFragment Class Referenceabstract

#include <G4PreCompoundFragment.hh>

Inheritance diagram for G4PreCompoundFragment:
Collaboration diagram for G4PreCompoundFragment:

Public Member Functions

 G4PreCompoundFragment (const G4ParticleDefinition *, G4VCoulombBarrier *aCoulombBarrier)
 
virtual ~G4PreCompoundFragment ()
 
G4double CalcEmissionProbability (const G4Fragment &aFragment)
 
G4double SampleKineticEnergy (const G4Fragment &aFragment)
 
- Public Member Functions inherited from G4VPreCompoundFragment
 G4VPreCompoundFragment (const G4ParticleDefinition *, G4VCoulombBarrier *aCoulombBarrier)
 
virtual ~G4VPreCompoundFragment ()
 
void Initialize (const G4Fragment &aFragment)
 
G4bool IsItPossible (const G4Fragment &aFragment) const
 
G4ReactionProductGetReactionProduct () const
 
G4int GetA () const
 
G4int GetZ () const
 
G4int GetRestA () const
 
G4int GetRestZ () const
 
G4double GetBindingEnergy () const
 
G4double GetEnergyThreshold () const
 
G4double GetEmissionProbability () const
 
G4double GetNuclearMass () const
 
G4double GetRestNuclearMass () const
 
const G4LorentzVectorGetMomentum () const
 
void SetMomentum (const G4LorentzVector &value)
 
void SetOPTxs (G4int)
 
void UseSICB (G4bool)
 

Protected Member Functions

virtual G4double GetAlpha () const =0
 
virtual G4double GetBeta () const =0
 
G4double CrossSection (G4double ekin) const
 
virtual G4double ProbabilityDistributionFunction (G4double K, const G4Fragment &aFragment)=0
 

Additional Inherited Members

- Protected Attributes inherited from G4VPreCompoundFragment
G4DeexPrecoParameterstheParameters
 
G4Powg4calc
 
G4int theA
 
G4int theZ
 
G4int theResA
 
G4int theResZ
 
G4int theFragA
 
G4int theFragZ
 
G4double theResA13
 
G4double theBindingEnergy
 
G4double theMinKinEnergy
 
G4double theMaxKinEnergy
 
G4double theResMass
 
G4double theReducedMass
 
G4double theMass
 
G4double theEmissionProbability
 
G4double theCoulombBarrier
 
G4int OPTxs
 
G4bool useSICB
 

Detailed Description

Definition at line 44 of file G4PreCompoundFragment.hh.

Constructor & Destructor Documentation

G4PreCompoundFragment::G4PreCompoundFragment ( const G4ParticleDefinition p,
G4VCoulombBarrier aCoulombBarrier 
)

Definition at line 42 of file G4PreCompoundFragment.cc.

44  : G4VPreCompoundFragment(p, aCoulBarrier)
45 {
46  muu = probmax = 0.0;
47  if(0 == theZ) { index = 0; }
48  else if(1 == theZ) { index = theA; }
49  else { index = theA + 1; }
50 }
G4VPreCompoundFragment(const G4ParticleDefinition *, G4VCoulombBarrier *aCoulombBarrier)
G4PreCompoundFragment::~G4PreCompoundFragment ( )
virtual

Definition at line 52 of file G4PreCompoundFragment.cc.

53 {}

Member Function Documentation

G4double G4PreCompoundFragment::CalcEmissionProbability ( const G4Fragment aFragment)
virtual

Implements G4VPreCompoundFragment.

Definition at line 56 of file G4PreCompoundFragment.cc.

57 {
58  //G4cout << theCoulombBarrier << " " << GetMaximalKineticEnergy() << G4endl;
59  // If theCoulombBarrier effect is included in the emission probabilities
60  // Coulomb barrier is the lower limit of integration over kinetic energy
61 
63 
64  if (theMaxKinEnergy <= theMinKinEnergy) { return 0.0; }
65 
66  // compute power once
67  if(OPTxs <= 2) {
69  } else {
71  }
72 
74  IntegrateEmissionProbability(theMinKinEnergy,theMaxKinEnergy,aFragment);
75  /*
76  G4cout << "## G4PreCompoundFragment::CalcEmisProb "
77  << "Z= " << aFragment.GetZ_asInt()
78  << " A= " << aFragment.GetA_asInt()
79  << " Elow= " << LowerLimit/MeV
80  << " Eup= " << UpperLimit/MeV
81  << " prob= " << theEmissionProbability
82  << G4endl;
83  */
85 }
static G4double ComputePowerParameter(G4int resA, G4int idx)
static G4double ComputePowerParameter(G4int resA, G4int idx)

Here is the call graph for this function:

G4double G4PreCompoundFragment::CrossSection ( G4double  ekin) const
protected

Definition at line 116 of file G4PreCompoundFragment.cc.

117 {
118  G4double res;
119  if(OPTxs == 0) {
120  res = GetOpt0(ekin);
121 
122  } else if(OPTxs <= 2) {
124  theResA13, muu,
125  index, theZ, theResA);
126 
127  } else {
129  theResA13, muu,
130  index, theZ, theA, theResA);
131  }
132  return res;
133 }
static G4double ComputeCrossSection(G4double K, G4double cb, G4double resA13, G4double amu1, G4int idx, G4int Z, G4int resA)
static G4double ComputeCrossSection(G4double K, G4double cb, G4double resA13, G4double amu1, G4int idx, G4int Z, G4int A, G4int resA)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

virtual G4double G4PreCompoundFragment::GetAlpha ( ) const
protectedpure virtual
virtual G4double G4PreCompoundFragment::GetBeta ( ) const
protectedpure virtual
virtual G4double G4PreCompoundFragment::ProbabilityDistributionFunction ( G4double  K,
const G4Fragment aFragment 
)
protectedpure virtual

Implemented in G4PreCompoundIon, and G4PreCompoundNucleon.

Here is the caller graph for this function:

G4double G4PreCompoundFragment::SampleKineticEnergy ( const G4Fragment aFragment)
virtual

Implements G4VPreCompoundFragment.

Definition at line 144 of file G4PreCompoundFragment.cc.

145 {
147  static const G4double toler = 1.25;
148  probmax *= toler;
149  G4double prob, T(0.0);
150  CLHEP::HepRandomEngine* rndm = G4Random::getTheEngine();
151  G4int i;
152  for(i=0; i<100; ++i) {
153  T = theMinKinEnergy + delta*rndm->flat();
154  prob = ProbabilityDistributionFunction(T, fragment);
155  /*
156  if(prob > probmax) {
157  G4cout << "G4PreCompoundFragment WARNING: prob= " << prob
158  << " probmax= " << probmax << G4endl;
159  G4cout << "i= " << i << " Z= " << theZ << " A= " << theA
160  << " resZ= " << theResZ << " resA= " << theResA << "\n"
161  << " T= " << T << " Tmax= " << theMaxKinEnergy
162  << " Tmin= " << limit
163  << G4endl;
164  for(G4int i=0; i<N; ++i) { G4cout << " " << probability[i]; }
165  G4cout << G4endl;
166  }
167  */
168  // Loop checking, 05-Aug-2015, Vladimir Ivanchenko
169  if(probmax*rndm->flat() <= prob) { break; }
170  }
171  /*
172  G4cout << "G4PreCompoundFragment: i= " << i << " Z= " << theZ << " A= " << theA
173  <<" T(MeV)= " << T << " Emin(MeV)= " << theMinKinEnergy << " Emax= "
174  << theMaxKinEnergy << G4endl;
175  */
176  return T;
177 }
virtual double flat()=0
int G4int
Definition: G4Types.hh:78
virtual G4double ProbabilityDistributionFunction(G4double K, const G4Fragment &aFragment)=0
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:


The documentation for this class was generated from the following files: