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

#include <G4GeneratorPrecompoundInterface.hh>

Inheritance diagram for G4GeneratorPrecompoundInterface:
Collaboration diagram for G4GeneratorPrecompoundInterface:

Public Member Functions

 G4GeneratorPrecompoundInterface (G4VPreCompoundModel *p=0)
 
virtual ~G4GeneratorPrecompoundInterface ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
 
virtual G4ReactionProductVectorPropagateNuclNucl (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus, G4V3DNucleus *theProjectileNucleus)
 
void SetCaptureThreshold (G4double)
 
virtual void PropagateModelDescription (std::ostream &) const
 
- Public Member Functions inherited from G4VIntraNuclearTransportModel
 G4VIntraNuclearTransportModel (const G4String &mName="CascadeModel", G4VPreCompoundModel *ptr=nullptr)
 
virtual ~G4VIntraNuclearTransportModel ()
 
void SetDeExcitation (G4VPreCompoundModel *ptr)
 
void Set3DNucleus (G4V3DNucleus *const value)
 
void SetPrimaryProjectile (const G4HadProjectile &aPrimary)
 
const G4StringGetModelName () const
 
virtual void ModelDescription (std::ostream &outFile) const
 
- 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 BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4VIntraNuclearTransportModel
G4V3DNucleusGet3DNucleus () const
 
G4VPreCompoundModelGetDeExcitation () const
 
const G4HadProjectileGetPrimaryProjectile () const
 
- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4VIntraNuclearTransportModel
G4String theTransportModelName
 
G4V3DNucleusthe3DNucleus
 
G4VPreCompoundModeltheDeExcitation
 
const G4HadProjectilethePrimaryProjectile
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 60 of file G4GeneratorPrecompoundInterface.hh.

Constructor & Destructor Documentation

G4GeneratorPrecompoundInterface::G4GeneratorPrecompoundInterface ( G4VPreCompoundModel p = 0)

Definition at line 79 of file G4GeneratorPrecompoundInterface.cc.

80 : CaptureThreshold(70*MeV) // Uzhi 1.05.2015 10 ->70
81 {
84 
87  He3 =G4He3::He3();
88  He4 =G4Alpha::Alpha();
89 
90  ANTIproton=G4AntiProton::AntiProton();
91  ANTIneutron=G4AntiNeutron::AntiNeutron();
92 
93  ANTIdeuteron=G4AntiDeuteron::AntiDeuteron();
94  ANTItriton =G4AntiTriton::AntiTriton();
95  ANTIHe3 =G4AntiHe3::AntiHe3();
96  ANTIHe4 =G4AntiAlpha::AntiAlpha();
97 
98  if(preModel) { SetDeExcitation(preModel); }
99  else {
100  G4HadronicInteraction* hadi =
102  G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(hadi);
103  if(!pre) { pre = new G4PreCompoundModel(); }
104  SetDeExcitation(pre);
105  }
106 }
static G4AntiHe3 * AntiHe3()
Definition: G4AntiHe3.cc:94
static G4AntiDeuteron * AntiDeuteron()
static G4AntiAlpha * AntiAlpha()
Definition: G4AntiAlpha.cc:89
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
static G4Triton * Triton()
Definition: G4Triton.cc:95
static G4Proton * Proton()
Definition: G4Proton.cc:93
static G4Neutron * Neutron()
Definition: G4Neutron.cc:104
static G4Deuteron * Deuteron()
Definition: G4Deuteron.cc:94
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
void SetDeExcitation(G4VPreCompoundModel *ptr)
static constexpr double MeV
Definition: G4SIunits.hh:214
static G4Alpha * Alpha()
Definition: G4Alpha.cc:89
static G4AntiTriton * AntiTriton()
Definition: G4AntiTriton.cc:94
static G4He3 * He3()
Definition: G4He3.cc:94
static G4AntiNeutron * AntiNeutron()

Here is the call graph for this function:

G4GeneratorPrecompoundInterface::~G4GeneratorPrecompoundInterface ( )
virtual

Definition at line 108 of file G4GeneratorPrecompoundInterface.cc.

109 {
110 }

Member Function Documentation

G4HadFinalState * G4GeneratorPrecompoundInterface::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtual

Implements G4HadronicInteraction.

Definition at line 326 of file G4GeneratorPrecompoundInterface.cc.

327 {
328  G4cout << "G4GeneratorPrecompoundInterface: ApplyYourself interface called stand-allone."
329  << G4endl;
330  G4cout << "This class is only a mediator between generator and precompound"<<G4endl;
331  G4cout << "Please remove from your physics list."<<G4endl;
332  throw G4HadronicException(__FILE__, __LINE__, "SEVERE: G4GeneratorPrecompoundInterface model interface called stand-allone.");
333  return new G4HadFinalState;
334 }
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61
G4ReactionProductVector * G4GeneratorPrecompoundInterface::Propagate ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus 
)
virtual

Implements G4VIntraNuclearTransportModel.

Definition at line 113 of file G4GeneratorPrecompoundInterface.cc.

114 {
115  #ifdef debugPrecoInt
116  G4cout<<G4endl<<"G4GeneratorPrecompoundInterface::Propagate"<<G4endl;
117  G4cout<<"Target A and Z "<<theNucleus->GetMassNumber()<<" "<<theNucleus->GetCharge()<<G4endl;
118  G4cout<<"Directly produced particles number "<<theSecondaries->size()<<G4endl;
119  #endif
120 
121  G4ReactionProductVector * theTotalResult = new G4ReactionProductVector;
122 
123  // decay the strong resonances
124  G4DecayKineticTracks decay(theSecondaries);
125  #ifdef debugPrecoInt
126  G4cout<<"Final stable particles number "<<theSecondaries->size()<<G4endl;
127  #endif
128 
129  // prepare the fragment
130  G4int anA=theNucleus->GetMassNumber();
131  G4int aZ=theNucleus->GetCharge();
132 // G4double TargetNucleusMass = G4NucleiProperties::GetNuclearMass(anA, aZ);
133 
134  G4int numberOfEx = 0;
135  G4int numberOfCh = 0;
136  G4int numberOfHoles = 0;
137 
138  G4double R = theNucleus->GetNuclearRadius();
139 
140  G4LorentzVector captured4Momentum(0.,0.,0.,0.);
141  G4LorentzVector Residual4Momentum(0.,0.,0.,0.); // TargetNucleusMass is not need at the moment
142  G4LorentzVector Secondary4Momentum(0.,0.,0.,0.);
143 
144  // loop over secondaries
145  G4KineticTrackVector::iterator iter;
146  for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
147  {
148  const G4ParticleDefinition* part = (*iter)->GetDefinition();
149  G4double e = (*iter)->Get4Momentum().e();
150  G4double mass = (*iter)->Get4Momentum().mag();
151  G4ThreeVector mom = (*iter)->Get4Momentum().vect();
152  if((part != proton && part != neutron) ||
153 // Uzhi 2.05.2015 (e > mass + CaptureThreshold) ||
154  ((*iter)->GetPosition().mag() > R)) {
155  G4ReactionProduct * theNew = new G4ReactionProduct(part);
156  theNew->SetMomentum(mom);
157  theNew->SetTotalEnergy(e);
158  theTotalResult->push_back(theNew);
159  Secondary4Momentum += (*iter)->Get4Momentum(); // Uzhi 29 April
160  #ifdef debugPrecoInt
161  G4cout<<"Secondary 4Mom "<<part->GetParticleName()<<" "<<(*iter)->Get4Momentum()<<" "
162  <<(*iter)->Get4Momentum().mag()<<G4endl;
163  #endif
164  } else {
165  if( e-mass > -CaptureThreshold*G4Log( G4UniformRand()) ) { // Added by Uzhi 2.05.2015
166  G4ReactionProduct * theNew = new G4ReactionProduct(part);
167  theNew->SetMomentum(mom);
168  theNew->SetTotalEnergy(e);
169  theTotalResult->push_back(theNew);
170  Secondary4Momentum += (*iter)->Get4Momentum(); // Uzhi 29 April
171  #ifdef debugPrecoInt
172  G4cout<<"Secondary 4Mom "<<part->GetParticleName()<<" "<<(*iter)->Get4Momentum()<<" "
173  <<(*iter)->Get4Momentum().mag()<<G4endl;
174  #endif
175  } else {
176  // within the nucleus, neutron or proton
177  // now calculate A, Z of the fragment, momentum, number of exciton states
178  ++anA;
179  ++numberOfEx;
180  G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
181  aZ += Z;
182  numberOfCh += Z;
183  captured4Momentum += (*iter)->Get4Momentum();
184  #ifdef debugPrecoInt
185  G4cout<<"Captured 4Mom "<<part->GetParticleName()<<(*iter)->Get4Momentum()<<G4endl;
186  #endif
187  }
188  }
189  delete (*iter);
190  }
191  delete theSecondaries;
192 
193  // loop over wounded nucleus
194  G4Nucleon * theCurrentNucleon =
195  theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
196  while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
197  {
198  if(theCurrentNucleon->AreYouHit()) {
199  ++numberOfHoles;
200  ++numberOfEx;
201  --anA;
202  aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/eplus + 0.1);
203 
204  Residual4Momentum -= theCurrentNucleon->Get4Momentum();
205  }
206  theCurrentNucleon = theNucleus->GetNextNucleon();
207  }
208 
209  #ifdef debugPrecoInt
210  G4cout<<G4endl;
211  G4cout<<"Secondary 4Mom "<<Secondary4Momentum<<G4endl;
212  G4cout<<"Captured 4Mom "<<captured4Momentum<<G4endl;
213  G4cout<<"Sec + Captured "<<Secondary4Momentum+captured4Momentum<<G4endl;
214  G4cout<<"Residual4Mom "<<Residual4Momentum<<G4endl;
215  G4cout<<"Sum 4 momenta "
216  <<Secondary4Momentum + captured4Momentum + Residual4Momentum <<G4endl;
217  #endif
218 
219  // Check that we use QGS model; loop over wounded nucleus
220  G4bool QGSM(false);
221  theCurrentNucleon = theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
222  while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
223  {
224  if(theCurrentNucleon->AreYouHit())
225  {
226  if(theCurrentNucleon->Get4Momentum().mag() <
227  theCurrentNucleon->GetDefinition()->GetPDGMass()) QGSM=true;
228  }
229  theCurrentNucleon = theNucleus->GetNextNucleon();
230  }
231 
232  #ifdef debugPrecoInt
233  if(!QGSM){
234  G4cout<<G4endl;
235  G4cout<<"Residual A and Z "<<anA<<" "<<aZ<<G4endl;
236  G4cout<<"Residual 4Mom "<<Residual4Momentum<<G4endl;
237  if(numberOfEx == 0)
238  {G4cout<<"Residual 4Mom = 0 means that there were not wounded and captured nucleons"<<G4endl;}
239  }
240  #endif
241 
242  if(anA == 0) return theTotalResult;
243 
244  G4LorentzVector exciton4Momentum(0.,0.,0.,0.); // Uzhi 29 April
245  if(anA >= aZ)
246  {
247  if(!QGSM)
248  { // FTF model was used
250 
251 // G4LorentzVector exciton4Momentum = Residual4Momentum + captured4Momentum;
252  exciton4Momentum = Residual4Momentum + captured4Momentum;
253 //exciton4Momentum.setE(std::sqrt(exciton4Momentum.vect().mag2()+sqr(fMass)));
254  G4double ActualMass = exciton4Momentum.mag();
255  if(ActualMass <= fMass ) { //E*<=0, Uzhi 5.05.2015
256  exciton4Momentum.setE(std::sqrt(exciton4Momentum.vect().mag2()+sqr(fMass))); // Uzhi 13.05.2015
257  }
258 
259  #ifdef debugPrecoInt
260  G4double exEnergy = 0.0;
261  if(ActualMass <= fMass ) {exEnergy = 0.;} // Uzhi 5.05.2015
262  else {exEnergy = ActualMass - fMass;}
263  G4cout<<"Ground state residual Mass "<<fMass<<" E* "<<exEnergy<<G4endl;
264  #endif
265  }
266  else
267  { // QGS model was used
268  G4double InitialTargetMass =
269  G4NucleiProperties::GetNuclearMass(theNucleus->GetMassNumber(), theNucleus->GetCharge());
270 
271  exciton4Momentum =
272  GetPrimaryProjectile()->Get4Momentum() + G4LorentzVector(0.,0.,0.,InitialTargetMass)
273  -Secondary4Momentum;
274 
276  G4double ActualMass = exciton4Momentum.mag();
277 
278  #ifdef debugPrecoInt
279  G4cout<<G4endl;
280  G4cout<<"Residual A and Z "<<anA<<" "<<aZ<<G4endl;
281  G4cout<<"Residual4Momentum "<<exciton4Momentum<<G4endl;
282  G4cout<<"ResidualMass, GroundStateMass and E* "<<ActualMass<<" "<<fMass<<" "
283  <<ActualMass - fMass<<G4endl;
284  #endif
285 
286  if(ActualMass - fMass < 0.)
287  {
288  G4double ResE = std::sqrt(exciton4Momentum.vect().mag2() + sqr(fMass+10*MeV));
289  exciton4Momentum.setE(ResE);
290  #ifdef debugPrecoInt
291  G4cout<<"ActualMass - fMass < 0. "<<ActualMass<<" "<<fMass<<" "<<ActualMass - fMass<<G4endl;
292  G4int Uzhi; G4cin>>Uzhi;
293  #endif
294  }
295  }
296  // Need to de-excite the remnant nucleus only if excitation energy > 0.
297  G4Fragment anInitialState(anA, aZ, exciton4Momentum);
298  anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
299  anInitialState.SetNumberOfCharged(numberOfCh);
300  anInitialState.SetNumberOfHoles(numberOfHoles);
301 
302  G4ReactionProductVector * aPrecoResult =
303  theDeExcitation->DeExcite(anInitialState);
304  // fill pre-compound part into the result, and return
305  #ifdef debugPrecoInt
306  G4cout<<"Target fragment number "<<aPrecoResult->size()<<G4endl;
307  #endif
308  for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
309  {
310  theTotalResult->push_back(aPrecoResult->operator[](ll));
311  #ifdef debugPrecoInt
312  G4cout<<"Fragment "<<ll<<" "
313  <<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()<<" "
314  <<aPrecoResult->operator[](ll)->GetMomentum()<<" "
315  <<aPrecoResult->operator[](ll)->GetTotalEnergy()<<" "
316  <<aPrecoResult->operator[](ll)->GetDefinition()->GetPDGMass()<<G4endl;
317  #endif
318  }
319  delete aPrecoResult;
320  }
321 
322  return theTotalResult;
323 }
static G4double GetNuclearMass(const G4double A, const G4double Z)
virtual G4int GetCharge()=0
const G4HadProjectile * GetPrimaryProjectile() const
virtual G4ReactionProductVector * DeExcite(G4Fragment &aFragment)=0
virtual G4double GetNuclearRadius()=0
virtual G4bool StartLoop()=0
void SetMomentum(const G4double x, const G4double y, const G4double z)
virtual G4int GetMassNumber()=0
virtual const G4LorentzVector & Get4Momentum() const
Definition: G4Nucleon.hh:72
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
virtual const G4ParticleDefinition * GetDefinition() const
Definition: G4Nucleon.hh:85
std::vector< G4ReactionProduct * > G4ReactionProductVector
#define G4cin
Definition: G4ios.hh:60
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
double mag() const
bool G4bool
Definition: G4Types.hh:79
G4bool AreYouHit() const
Definition: G4Nucleon.hh:97
void SetTotalEnergy(const G4double en)
static constexpr double eplus
Definition: G4SIunits.hh:199
const G4LorentzVector & Get4Momentum() const
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double GetPDGMass() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
virtual G4Nucleon * GetNextNucleon()=0
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
CLHEP::HepLorentzVector G4LorentzVector

Here is the call graph for this function:

void G4GeneratorPrecompoundInterface::PropagateModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 335 of file G4GeneratorPrecompoundInterface.cc.

336 {
337  outFile << "G4GeneratorPrecompoundInterface interfaces a high\n"
338  << "energy model through the wounded nucleus to precompound de-excition.\n"
339  << "Low energy protons and neutron present among secondaries produced by \n"
340  << "the high energy generator and within the nucleus are captured. The wounded\n"
341  << "nucleus and the captured particles form an excited nuclear fragment. This\n"
342  << "fragment is passed to the Geant4 pre-compound model for de-excitation.\n"
343  << "Nuclear de-excitation:\n";
344  // preco
345 
346 }
G4ReactionProductVector * G4GeneratorPrecompoundInterface::PropagateNuclNucl ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus,
G4V3DNucleus theProjectileNucleus 
)
virtual

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 350 of file G4GeneratorPrecompoundInterface.cc.

352 {
353 #ifdef debugPrecoInt
354  G4cout<<G4endl<<"G4GeneratorPrecompoundInterface::PropagateNuclNucl "<<G4endl;
355  G4cout<<"Projectile A and Z "<<theProjectileNucleus->GetMassNumber()<<" "
356  <<theProjectileNucleus->GetCharge()<<G4endl;
357  G4cout<<"Target A and Z "<<theNucleus->GetMassNumber()<<" "
358  <<theNucleus->GetCharge()<<G4endl;
359  G4cout<<"Directly produced particles number "<<theSecondaries->size()<<G4endl;
360  G4cout<<"Projectile 4Mom and mass "<<GetPrimaryProjectile()->Get4Momentum()<<" "
361  <<GetPrimaryProjectile()->Get4Momentum().mag()<<G4endl<<G4endl;
362 #endif
363 
364  // prepare the target residual
365  G4int anA=theNucleus->GetMassNumber();
366  G4int aZ=theNucleus->GetCharge();
367  G4int numberOfEx = 0;
368  G4int numberOfCh = 0;
369  G4int numberOfHoles = 0;
370  G4double exEnergy = 0.0;
371  G4double R = theNucleus->GetNuclearRadius();
372  G4LorentzVector Target4Momentum(0.,0.,0.,0.);
373 
374  // loop over wounded target nucleus
375  G4Nucleon * theCurrentNucleon =
376  theNucleus->StartLoop() ? theNucleus->GetNextNucleon() : 0;
377  while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
378  {
379  if(theCurrentNucleon->AreYouHit()) {
380  ++numberOfHoles;
381  ++numberOfEx;
382  --anA;
383  aZ -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
384  eplus + 0.1);
385  exEnergy += theCurrentNucleon->GetBindingEnergy();
386  Target4Momentum -=theCurrentNucleon->Get4Momentum();
387  }
388  theCurrentNucleon = theNucleus->GetNextNucleon();
389  }
390 
391 #ifdef debugPrecoInt
392  G4cout<<"Residual Target A Z E* 4mom "<<anA<<" "<<aZ<<" "<<exEnergy<<" "
393  <<Target4Momentum<<G4endl;
394 #endif
395 
396  // prepare the projectile residual
397 
398  G4bool ProjectileIsAntiNucleus=
400 
402 
403  G4int anAb=theProjectileNucleus->GetMassNumber();
404  G4int aZb=theProjectileNucleus->GetCharge();
405  G4int numberOfExB = 0;
406  G4int numberOfChB = 0;
407  G4int numberOfHolesB = 0;
408  G4double exEnergyB = 0.0;
409  G4double Rb = theProjectileNucleus->GetNuclearRadius();
410  G4LorentzVector Projectile4Momentum(0.,0.,0.,0.);
411 
412  // loop over wounded projectile nucleus
413  theCurrentNucleon =
414  theProjectileNucleus->StartLoop() ? theProjectileNucleus->GetNextNucleon() : 0;
415  while(theCurrentNucleon) /* Loop checking, 31.08.2015, G.Folger */
416  {
417  if(theCurrentNucleon->AreYouHit()) {
418  ++numberOfHolesB;
419  ++numberOfExB;
420  --anAb;
421  if(!ProjectileIsAntiNucleus) {
422  aZb -= G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
423  eplus + 0.1);
424  } else {
425  aZb += G4int(theCurrentNucleon->GetDefinition()->GetPDGCharge()/
426  eplus - 0.1);
427  }
428  exEnergyB += theCurrentNucleon->GetBindingEnergy();
429  Projectile4Momentum -=theCurrentNucleon->Get4Momentum();
430  }
431  theCurrentNucleon = theProjectileNucleus->GetNextNucleon();
432  }
433 
434  G4bool ExistTargetRemnant = G4double (numberOfHoles) <
435  0.3* G4double (numberOfHoles + anA);
436  G4bool ExistProjectileRemnant= G4double (numberOfHolesB) <
437  0.3*G4double (numberOfHolesB + anAb);
438 
439 #ifdef debugPrecoInt
440  G4cout<<"Projectile residual A Z E* 4mom "<<anAb<<" "<<aZb<<" "<<exEnergyB<<" "
441  <<Projectile4Momentum<<G4endl;
442  G4cout<<" ExistTargetRemnant ExistProjectileRemnant "
443  <<ExistTargetRemnant<<" "<< ExistProjectileRemnant<<G4endl;
444 #endif
445  //-----------------------------------------------------------------------------
446  // decay the strong resonances
447  G4ReactionProductVector * theTotalResult = new G4ReactionProductVector;
448  G4DecayKineticTracks decay(theSecondaries);
449  #ifdef debugPrecoInt
450  G4cout<<"Secondary stable particles number "<<theSecondaries->size()<<G4endl;
451  #endif
452 
453 #ifdef debugPrecoInt
454  G4LorentzVector secondary4Momemtum(0,0,0,0);
455  G4int SecondrNum(0);
456 #endif
457 
458  // loop over secondaries
459  G4KineticTrackVector::iterator iter;
460  for(iter=theSecondaries->begin(); iter !=theSecondaries->end(); ++iter)
461  {
462  const G4ParticleDefinition* part = (*iter)->GetDefinition();
463  G4LorentzVector aTrack4Momentum=(*iter)->Get4Momentum();
464 
465  if( part != proton && part != neutron &&
466  (part != ANTIproton && ProjectileIsAntiNucleus) &&
467  (part != ANTIneutron && ProjectileIsAntiNucleus) )
468  {
469  G4ReactionProduct * theNew = new G4ReactionProduct(part);
470  theNew->SetMomentum(aTrack4Momentum.vect());
471  theNew->SetTotalEnergy(aTrack4Momentum.e());
472  theTotalResult->push_back(theNew);
473 #ifdef debugPrecoInt
474  SecondrNum++;
475  secondary4Momemtum += (*iter)->Get4Momentum();
476  G4cout<<"Secondary "<<SecondrNum<<" "
477  <<theNew->GetDefinition()->GetParticleName()<<" "
478  <<theNew->GetMomentum()<<" "<<theNew->GetTotalEnergy()<<G4endl;
479 
480 #endif
481  delete (*iter);
482  continue;
483  }
484 
485  G4bool CanBeCapturedByTarget = false;
486  if( part == proton || part == neutron)
487  {
488  CanBeCapturedByTarget = ExistTargetRemnant &&
489  (-CaptureThreshold*G4Log( G4UniformRand()) >
490  (aTrack4Momentum + Target4Momentum).mag() -
491  aTrack4Momentum.mag() - Target4Momentum.mag()) &&
492  ((*iter)->GetPosition().mag() < R);
493  }
494  // ---------------------------
495  G4LorentzVector Position((*iter)->GetPosition(), (*iter)->GetFormationTime());
496  Position.boost(bst);
497 
498  G4bool CanBeCapturedByProjectile = false;
499 
500  if( !ProjectileIsAntiNucleus &&
501  ( part == proton || part == neutron))
502  {
503  CanBeCapturedByProjectile = ExistProjectileRemnant &&
504  (-CaptureThreshold*G4Log( G4UniformRand()) >
505  (aTrack4Momentum + Projectile4Momentum).mag() -
506  aTrack4Momentum.mag() - Projectile4Momentum.mag()) &&
507  (Position.vect().mag() < Rb);
508  }
509 
510  if( ProjectileIsAntiNucleus &&
511  ( part == ANTIproton || part == ANTIneutron))
512  {
513  CanBeCapturedByProjectile = ExistProjectileRemnant &&
514  (-CaptureThreshold*G4Log( G4UniformRand()) >
515  (aTrack4Momentum + Projectile4Momentum).mag() -
516  aTrack4Momentum.mag() - Projectile4Momentum.mag()) &&
517  (Position.vect().mag() < Rb);
518  }
519 
520  if(CanBeCapturedByTarget && CanBeCapturedByProjectile)
521  {
522  if(G4UniformRand() < 0.5)
523  { CanBeCapturedByTarget = true; CanBeCapturedByProjectile = false;}
524  else
525  { CanBeCapturedByTarget = false; CanBeCapturedByProjectile = true;}
526  }
527 
528  if(CanBeCapturedByTarget)
529  {
530  // within the target nucleus, neutron or proton
531  // now calculate A, Z of the fragment, momentum,
532  // number of exciton states
533 #ifdef debugPrecoInt
534  G4cout<<"Track is CapturedByTarget "<<" "<<part->GetParticleName()<<" "
535  <<aTrack4Momentum<<" "<<aTrack4Momentum.mag()<<G4endl;
536 #endif
537  ++anA;
538  ++numberOfEx;
539  G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
540  aZ += Z;
541  numberOfCh += Z;
542  Target4Momentum +=aTrack4Momentum;
543  delete (*iter);
544  } else if(CanBeCapturedByProjectile)
545  {
546  // within the projectile nucleus, neutron or proton
547  // now calculate A, Z of the fragment, momentum,
548  // number of exciton states
549 #ifdef debugPrecoInt
550  G4cout<<"Track is CapturedByProjectile"<<" "<<part->GetParticleName()<<" "
551  <<aTrack4Momentum<<" "<<aTrack4Momentum.mag()<<G4endl;
552 #endif
553  ++anAb;
554  ++numberOfExB;
555  G4int Z = G4int(part->GetPDGCharge()/eplus + 0.1);
556  if( ProjectileIsAntiNucleus ) Z=-Z;
557  aZb += Z;
558  numberOfChB += Z;
559  Projectile4Momentum +=aTrack4Momentum;
560  delete (*iter);
561  } else
562  { // the track is not captured
563  G4ReactionProduct * theNew = new G4ReactionProduct(part);
564  theNew->SetMomentum(aTrack4Momentum.vect());
565  theNew->SetTotalEnergy(aTrack4Momentum.e());
566  theTotalResult->push_back(theNew);
567 
568 #ifdef debugPrecoInt
569  SecondrNum++;
570  secondary4Momemtum += (*iter)->Get4Momentum();
571 /*
572  G4cout<<"Secondary "<<SecondrNum<<" "
573  <<theNew->GetDefinition()->GetParticleName()<<" "
574  <<secondary4Momemtum<<G4endl;
575 */
576 #endif
577  delete (*iter);
578  continue;
579  }
580  }
581  delete theSecondaries;
582  //-----------------------------------------------------
583 
584  #ifdef debugPrecoInt
585  G4cout<<"Final target residual A Z E* 4mom "<<anA<<" "<<aZ<<" "
586  <<exEnergy<<" "<<Target4Momentum<<G4endl;
587  #endif
588 
589  if(0!=anA )
590  {
592 
593  if((anA == theNucleus->GetMassNumber()) && (exEnergy <= 0.))
594  {Target4Momentum.setE(fMass);}
595 
596  G4double RemnMass=Target4Momentum.mag();
597 
598  if(RemnMass < fMass)
599  {
600  RemnMass=fMass + exEnergy;
601  Target4Momentum.setE(std::sqrt(Target4Momentum.vect().mag2() +
602  RemnMass*RemnMass));
603  } else
604  { exEnergy=RemnMass-fMass;}
605 
606  if( exEnergy < 0.) exEnergy=0.;
607 
608  // Need to de-excite the remnant nucleus
609  G4Fragment anInitialState(anA, aZ, Target4Momentum);
610  anInitialState.SetNumberOfParticles(numberOfEx-numberOfHoles);
611  anInitialState.SetNumberOfCharged(numberOfCh);
612  anInitialState.SetNumberOfHoles(numberOfHoles);
613 
614  G4ReactionProductVector * aPrecoResult =
615  theDeExcitation->DeExcite(anInitialState);
616 
617  #ifdef debugPrecoInt
618  G4cout<<"Target fragment number "<<aPrecoResult->size()<<G4endl;
619  #endif
620 
621  // fill pre-compound part into the result, and return
622  for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
623  {
624  theTotalResult->push_back(aPrecoResult->operator[](ll));
625  #ifdef debugPrecoInt
626  G4cout<<"Target fragment "<<ll<<" "
627  <<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()<<" "
628  <<aPrecoResult->operator[](ll)->GetMomentum()<<" "
629  <<aPrecoResult->operator[](ll)->GetTotalEnergy()<<" "
630  <<aPrecoResult->operator[](ll)->GetMass()<<G4endl;
631  #endif
632  }
633  delete aPrecoResult;
634  }
635 
636  //-----------------------------------------------------
637  if((anAb == theProjectileNucleus->GetMassNumber())&& (exEnergyB <= 0.))
638  {Projectile4Momentum = GetPrimaryProjectile()->Get4Momentum();}
639 
640  #ifdef debugPrecoInt
641  G4cout<<"Final projectile residual A Z E* Pmom Pmag2 "<<anAb<<" "<<aZb<<" "
642  <<exEnergyB<<" "<<Projectile4Momentum<<" "
643  <<Projectile4Momentum.mag2()<<G4endl;
644  #endif
645 
646  if(0!=anAb)
647  {
648  // G4ThreeVector bstToCM =Projectile4Momentum.findBoostToCM(); // Uzhi Apr. 2015
649  // Projectile4Momentum.boost(bstToCM); // Uzhi Apr. 2015
650 
651  G4double fMass = G4NucleiProperties::GetNuclearMass(anAb, aZb);
652  G4double RemnMass=Projectile4Momentum.mag();
653 
654  if(RemnMass < fMass)
655  {
656  RemnMass=fMass + exEnergyB;
657  Projectile4Momentum.setE(std::sqrt(Projectile4Momentum.vect().mag2() + // Uzhi 8.05.2015
658  RemnMass*RemnMass)); // Uzhi 8.05.2015
659  } else
660  { exEnergyB=RemnMass-fMass;}
661 
662  if( exEnergyB < 0.) exEnergyB=0.;
663 
664  G4ThreeVector bstToCM =Projectile4Momentum.findBoostToCM(); // Uzhi Apr. 2015
665  Projectile4Momentum.boost(bstToCM); // Uzhi Apr. 2015
666 
667  // Need to de-excite the remnant nucleus
668  G4Fragment anInitialState(anAb, aZb, Projectile4Momentum);
669  anInitialState.SetNumberOfParticles(numberOfExB-numberOfHolesB);
670  anInitialState.SetNumberOfCharged(numberOfChB);
671  anInitialState.SetNumberOfHoles(numberOfHolesB);
672 
673  G4ReactionProductVector * aPrecoResult =
674  theDeExcitation->DeExcite(anInitialState);
675 
676  #ifdef debugPrecoInt
677  G4cout<<"Projectile fragment number "<<aPrecoResult->size()<<G4endl;
678  #endif
679 
680  // fill pre-compound part into the result, and return
681  for(unsigned int ll=0; ll<aPrecoResult->size(); ++ll)
682  {
683  G4LorentzVector tmp=G4LorentzVector(aPrecoResult->operator[](ll)->GetMomentum(), // Uzhi 2015
684  aPrecoResult->operator[](ll)->GetTotalEnergy());// Uzhi 2015
685  tmp.boost(-bstToCM); // Transformation to the system of original remnant // Uzhi 2015
686  aPrecoResult->operator[](ll)->SetMomentum(tmp.vect()); // Uzhi 2015
687  aPrecoResult->operator[](ll)->SetTotalEnergy(tmp.e()); // Uzhi 2015
688 
689  if(ProjectileIsAntiNucleus)
690  {
691  const G4ParticleDefinition * aFragment=aPrecoResult->operator[](ll)->GetDefinition();
692  const G4ParticleDefinition * LastFragment=aFragment;
693  if (aFragment == proton) {LastFragment=G4AntiProton::AntiProtonDefinition();}
694  else if(aFragment == neutron) {LastFragment=G4AntiNeutron::AntiNeutronDefinition();}
695  else if(aFragment == deuteron){LastFragment=G4AntiDeuteron::AntiDeuteronDefinition();}
696  else if(aFragment == triton) {LastFragment=G4AntiTriton::AntiTritonDefinition();}
697  else if(aFragment == He3) {LastFragment=G4AntiHe3::AntiHe3Definition();}
698  else if(aFragment == He4) {LastFragment=G4AntiAlpha::AntiAlphaDefinition();}
699  else {}
700 
701  aPrecoResult->operator[](ll)->SetDefinitionAndUpdateE(LastFragment);
702  }
703 
704  #ifdef debugPrecoInt
705  G4cout<<"Projectile fragment "<<ll<<" "
706  <<aPrecoResult->operator[](ll)->GetDefinition()->GetParticleName()<<" "
707  <<aPrecoResult->operator[](ll)->GetMomentum()<<" "
708  <<aPrecoResult->operator[](ll)->GetTotalEnergy()<<" "
709  <<aPrecoResult->operator[](ll)->GetMass()<<G4endl;
710  #endif
711 //Uzhi
712  theTotalResult->push_back(aPrecoResult->operator[](ll));
713  }
714 
715  delete aPrecoResult;
716  }
717 
718  return theTotalResult;
719 }
static G4AntiTriton * AntiTritonDefinition()
Definition: G4AntiTriton.cc:89
Hep3Vector boostVector() const
static G4double GetNuclearMass(const G4double A, const G4double Z)
virtual G4int GetCharge()=0
const G4HadProjectile * GetPrimaryProjectile() const
virtual G4ReactionProductVector * DeExcite(G4Fragment &aFragment)=0
virtual G4double GetNuclearRadius()=0
virtual G4bool StartLoop()=0
void SetMomentum(const G4double x, const G4double y, const G4double z)
virtual G4int GetMassNumber()=0
static G4AntiDeuteron * AntiDeuteronDefinition()
virtual const G4LorentzVector & Get4Momentum() const
Definition: G4Nucleon.hh:72
static G4AntiProton * AntiProtonDefinition()
Definition: G4AntiProton.cc:88
int G4int
Definition: G4Types.hh:78
static G4AntiNeutron * AntiNeutronDefinition()
const G4String & GetParticleName() const
virtual const G4ParticleDefinition * GetDefinition() const
Definition: G4Nucleon.hh:85
std::vector< G4ReactionProduct * > G4ReactionProductVector
const G4ParticleDefinition * GetDefinition() const
Hep3Vector vect() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
const G4ParticleDefinition * GetDefinition() const
double mag() const
bool G4bool
Definition: G4Types.hh:79
HepLorentzVector & boost(double, double, double)
G4bool AreYouHit() const
Definition: G4Nucleon.hh:97
void SetTotalEnergy(const G4double en)
static constexpr double eplus
Definition: G4SIunits.hh:199
const G4LorentzVector & Get4Momentum() const
static G4AntiHe3 * AntiHe3Definition()
Definition: G4AntiHe3.cc:89
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double GetTotalEnergy() const
G4ThreeVector GetMomentum() const
static G4AntiAlpha * AntiAlphaDefinition()
Definition: G4AntiAlpha.cc:84
#define G4endl
Definition: G4ios.hh:61
virtual G4Nucleon * GetNextNucleon()=0
double G4double
Definition: G4Types.hh:76
G4double GetPDGCharge() const
G4double GetBindingEnergy() const
Definition: G4Nucleon.hh:75
CLHEP::HepLorentzVector G4LorentzVector

Here is the call graph for this function:

void G4GeneratorPrecompoundInterface::SetCaptureThreshold ( G4double  value)
inline

Definition at line 107 of file G4GeneratorPrecompoundInterface.hh.

108 {
109  CaptureThreshold=value;
110 }
const XML_Char int const XML_Char * value
Definition: expat.h:331

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