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

#include <G4NeutronRadCapture.hh>

Inheritance diagram for G4NeutronRadCapture:
Collaboration diagram for G4NeutronRadCapture:

Public Member Functions

 G4NeutronRadCapture ()
 
virtual ~G4NeutronRadCapture ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus) final
 
virtual void InitialiseModel () final
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair
< G4double, G4double
GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double,
G4double
GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 52 of file G4NeutronRadCapture.hh.

Constructor & Destructor Documentation

G4NeutronRadCapture::G4NeutronRadCapture ( )
explicit

Definition at line 53 of file G4NeutronRadCapture.cc.

54  : G4HadronicInteraction("nRadCapture"),
55  photonEvaporation(nullptr),lab4mom(0.,0.,0.,0.)
56 {
57  lowestEnergyLimit = 10*CLHEP::eV;
58  minExcitation = 0.1*CLHEP::keV;
59  SetMinEnergy( 0.0*CLHEP::GeV );
60  SetMaxEnergy( 100.*CLHEP::TeV );
61 
62  theTableOfIons = G4ParticleTable::GetParticleTable()->GetIonTable();
63 }
static constexpr double keV
void SetMinEnergy(G4double anEnergy)
G4IonTable * GetIonTable() const
G4HadronicInteraction(const G4String &modelName="HadronicModel")
static constexpr double eV
static constexpr double GeV
static G4ParticleTable * GetParticleTable()
void SetMaxEnergy(const G4double anEnergy)
static constexpr double TeV

Here is the call graph for this function:

G4NeutronRadCapture::~G4NeutronRadCapture ( )
virtual

Definition at line 65 of file G4NeutronRadCapture.cc.

66 {
67  delete photonEvaporation;
68 }

Member Function Documentation

G4HadFinalState * G4NeutronRadCapture::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
finalvirtual

Implements G4HadronicInteraction.

Definition at line 82 of file G4NeutronRadCapture.cc.

84 {
87 
88  G4int A = theNucleus.GetA_asInt();
89  G4int Z = theNucleus.GetZ_asInt();
90 
91  G4double time = aTrack.GetGlobalTime();
92 
93  // Create initial state
94  lab4mom.set(0.,0.,0.,G4NucleiProperties::GetNuclearMass(A, Z));
95  lab4mom += aTrack.Get4Momentum();
96 
97  G4double M = lab4mom.mag();
98  ++A;
100  //G4cout << "Capture start: Z= " << Z << " A= " << A
101  // << " LabM= " << M << " Mcompound= " << mass << G4endl;
102 
103  // simplified method of 1 gamma emission
104  if(A <= 4) {
105 
106  G4ThreeVector bst = lab4mom.boostVector();
107 
108  if(M - mass <= lowestEnergyLimit) {
109  return &theParticleChange;
110  }
111 
112  if (verboseLevel > 1) {
113  G4cout << "G4NeutronRadCapture::DoIt: Eini(MeV)="
114  << aTrack.GetKineticEnergy()/MeV << " Eexc(MeV)= "
115  << (M - mass)/MeV
116  << " Z= " << Z << " A= " << A << G4endl;
117  }
118  G4double e1 = (M - mass)*(M + mass)/(2*M);
119 
120  G4double cost = 2.0*G4UniformRand() - 1.0;
121  if(cost > 1.0) {cost = 1.0;}
122  else if(cost < -1.0) {cost = -1.0;}
123  G4double sint = std::sqrt((1. - cost)*(1.0 + cost));
125 
126  G4LorentzVector lv2(e1*sint*std::cos(phi),e1*sint*std::sin(phi),
127  e1*cost,e1);
128  lv2.boost(bst);
129  G4HadSecondary* news =
131  news->SetTime(time);
133  delete news;
134 
135  const G4ParticleDefinition* theDef = 0;
136 
137  lab4mom -= lv2;
138  if (Z == 1 && A == 2) {theDef = G4Deuteron::Deuteron();}
139  else if (Z == 1 && A == 3) {theDef = G4Triton::Triton();}
140  else if (Z == 2 && A == 3) {theDef = G4He3::He3();}
141  else if (Z == 2 && A == 4) {theDef = G4Alpha::Alpha();}
142  else { theDef = theTableOfIons->GetIon(Z,A,0.0,noFloat,0); }
143 
144  if (verboseLevel > 1) {
145  G4cout << "Gamma 4-mom: " << lv2 << " "
146  << theDef->GetParticleName() << " " << lab4mom << G4endl;
147  }
148  if(theDef) {
149  news = new G4HadSecondary(new G4DynamicParticle(theDef, lab4mom));
150  news->SetTime(time);
152  delete news;
153  }
154 
155  // Use photon evaporation
156  } else {
157 
158  // protection against wrong kinematic
159  if(M < mass) {
160  G4double etot = std::max(mass, lab4mom.e());
161  G4double ptot = std::sqrt((etot - mass)*(etot + mass));
162  G4ThreeVector v = lab4mom.vect().unit();
163  lab4mom.set(v.x()*ptot,v.y()*ptot,v.z()*ptot,etot);
164  }
165 
166  G4Fragment* aFragment = new G4Fragment(A, Z, lab4mom);
167 
168  if (verboseLevel > 1) {
169  G4cout << "G4NeutronRadCapture::ApplyYourself initial G4Fragmet:"
170  << G4endl;
171  G4cout << aFragment << G4endl;
172  }
173 
174  //
175  // Sample final state
176  //
177  G4FragmentVector* fv = photonEvaporation->BreakUpFragment(aFragment);
178  if(!fv) { fv = new G4FragmentVector(); }
179  fv->push_back(aFragment);
180  size_t n = fv->size();
181 
182  if (verboseLevel > 1) {
183  G4cout << "G4NeutronRadCapture: " << n << " final particle" << G4endl;
184  }
185  for(size_t i=0; i<n; ++i) {
186 
187  G4Fragment* f = (*fv)[i];
188  G4double etot = f->GetMomentum().e();
189 
190  Z = f->GetZ_asInt();
191  A = f->GetA_asInt();
192 
193  const G4ParticleDefinition* theDef;
194  if(0 == Z && 0 == A) {theDef = f->GetParticleDefinition();}
195  else if (Z == 1 && A == 2) {theDef = G4Deuteron::Deuteron();}
196  else if (Z == 1 && A == 3) {theDef = G4Triton::Triton();}
197  else if (Z == 2 && A == 3) {theDef = G4He3::He3();}
198  else if (Z == 2 && A == 4) {theDef = G4Alpha::Alpha();}
199  else {
200  G4double eexc = f->GetExcitationEnergy();
201  if(eexc <= minExcitation) { eexc = 0.0; }
202  theDef = theTableOfIons->GetIon(Z, A, eexc, noFloat, 0);
203  /*
204  G4cout << "### NC Find ion Z= " << Z << " A= " << A
205  << " Eexc(MeV)= " << eexc/MeV << " "
206  << theDef << G4endl;
207  */
208  }
209  G4double ekin = std::max(0.0,etot - theDef->GetPDGMass());
210  if (verboseLevel > 1) {
211  G4cout << i << ". " << theDef->GetParticleName()
212  << " Ekin(MeV)= " << etot/MeV
213  << " p: " << f->GetMomentum().vect()
214  << G4endl;
215  }
216  G4HadSecondary* news = new G4HadSecondary(
217  new G4DynamicParticle(theDef,
218  f->GetMomentum().vect().unit(),
219  ekin));
220  G4double timeF = f->GetCreationTime();
221  if(timeF < 0.0) { timeF = 0.0; }
222  news->SetTime(time + timeF);
224  delete news;
225  delete f;
226  }
227  delete fv;
228  }
229  //G4cout << "Capture done" << G4endl;
230  return &theParticleChange;
231 }
Hep3Vector boostVector() const
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4FragmentVector * BreakUpFragment(G4Fragment *theNucleus)
double x() const
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:503
const G4ParticleDefinition * GetParticleDefinition() const
Definition: G4Fragment.hh:438
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
double z() const
void SetStatusChange(G4HadFinalStateStatus aS)
Hep3Vector vect() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
G4int GetA_asInt() const
Definition: G4Fragment.hh:266
G4double GetCreationTime() const
Definition: G4Fragment.hh:448
double mag() const
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:307
G4double GetKineticEnergy() const
G4double GetGlobalTime() const
std::vector< G4Fragment * > G4FragmentVector
Definition: G4Fragment.hh:63
static G4Triton * Triton()
Definition: G4Triton.cc:95
const G4int n
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
const G4LorentzVector & Get4Momentum() const
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
tuple v
Definition: test.py:18
void SetTime(G4double aT)
void set(double x, double y, double z, double t)
G4double GetPDGMass() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
Hep3Vector unit() const
double y() const
G4int GetZ_asInt() const
Definition: G4Fragment.hh:271
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
static G4He3 * He3()
Definition: G4He3.cc:94
static constexpr double twopi
Definition: SystemOfUnits.h:55
#define noFloat
Definition: G4Ions.hh:118
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:283

Here is the call graph for this function:

void G4NeutronRadCapture::InitialiseModel ( )
finalvirtual

Reimplemented from G4HadronicInteraction.

Definition at line 70 of file G4NeutronRadCapture.cc.

71 {
72  if(photonEvaporation != nullptr) { return; }
73  G4DeexPrecoParameters* param =
75  minExcitation = param->GetMinExcitation();
76 
77  photonEvaporation = new G4PhotonEvaporation();
78  photonEvaporation->Initialise();
79  photonEvaporation->SetICM(true);
80 }
virtual void SetICM(G4bool)
G4DeexPrecoParameters * GetParameters()
static G4NuclearLevelData * GetInstance()
G4double GetMinExcitation() const

Here is the call graph for this function:


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