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

#include <G4BinaryCascade.hh>

Inheritance diagram for G4BinaryCascade:
Collaboration diagram for G4BinaryCascade:

Public Member Functions

 G4BinaryCascade (G4VPreCompoundModel *ptr=0)
 
virtual ~G4BinaryCascade ()
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
virtual G4ReactionProductVectorPropagate (G4KineticTrackVector *, G4V3DNucleus *)
 
virtual void ModelDescription (std::ostream &) const
 
virtual void PropagateModelDescription (std::ostream &) const
 
- Public Member Functions inherited from G4VIntraNuclearTransportModel
 G4VIntraNuclearTransportModel (const G4String &mName="CascadeModel", G4VPreCompoundModel *ptr=nullptr)
 
virtual ~G4VIntraNuclearTransportModel ()
 
virtual G4ReactionProductVectorPropagateNuclNucl (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus, G4V3DNucleus *theProjectileNucleus)
 
void SetDeExcitation (G4VPreCompoundModel *ptr)
 
void Set3DNucleus (G4V3DNucleus *const value)
 
void SetPrimaryProjectile (const G4HadProjectile &aPrimary)
 
const G4StringGetModelName () 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 70 of file G4BinaryCascade.hh.

Constructor & Destructor Documentation

G4BinaryCascade::G4BinaryCascade ( G4VPreCompoundModel ptr = 0)

Definition at line 115 of file G4BinaryCascade.cc.

115  :
116 G4VIntraNuclearTransportModel("Binary Cascade", ptr)
117 {
118  // initialise the resonance sector
119  G4ShortLivedConstructor ShortLived;
120  ShortLived.ConstructParticle();
121 
122  theCollisionMgr = new G4CollisionManager;
123  theDecay=new G4BCDecay;
124  theImR.push_back(theDecay);
125  theLateParticle= new G4BCLateParticle;
127  theImR.push_back(aAb);
128  G4Scatterer * aSc=new G4Scatterer;
129  theH1Scatterer = new G4Scatterer;
130  theImR.push_back(aSc);
131 
132  thePropagator = new G4RKPropagation;
133  theCurrentTime = 0.;
134  theBCminP = 45*MeV;
135  theCutOnP = 90*MeV;
136  theCutOnPAbsorb= 0*MeV; // No Absorption of slow Mesons, other than above G4MesonAbsorption
137 
138  // reuse existing pre-compound model
139  if(!ptr) {
142  G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(p);
143  if(!pre) { pre = new G4PreCompoundModel(); }
144  SetDeExcitation(pre);
145  }
146  theExcitationHandler = GetDeExcitation()->GetExcitationHandler();
147  SetMinEnergy(0.0*GeV);
148  SetMaxEnergy(10.1*GeV);
149  //PrintWelcomeMessage();
150  thePrimaryEscape = true;
151  thePrimaryType = 0;
152 
154 
155  // init data members
156  currentA=currentZ=0;
157  lateA=lateZ=0;
158  initialA=initialZ=0;
159  projectileA=projectileZ=0;
160  currentInitialEnergy=initial_nuclear_mass=0.;
161  massInNucleus=0.;
162  theOuterRadius=0.;
163 }
G4VPreCompoundModel * GetDeExcitation() const
static constexpr double perCent
Definition: G4SIunits.hh:332
const char * p
Definition: xmltok.h:285
G4VIntraNuclearTransportModel(const G4String &mName="CascadeModel", G4VPreCompoundModel *ptr=nullptr)
void SetMinEnergy(G4double anEnergy)
G4ExcitationHandler * GetExcitationHandler() const
G4HadronicInteraction * FindModel(const G4String &name)
static G4HadronicInteractionRegistry * Instance()
static constexpr double GeV
Definition: G4SIunits.hh:217
void SetMaxEnergy(const G4double anEnergy)
void SetDeExcitation(G4VPreCompoundModel *ptr)
static constexpr double MeV
Definition: G4SIunits.hh:214
void SetEnergyMomentumCheckLevels(G4double relativeLevel, G4double absoluteLevel)

Here is the call graph for this function:

G4BinaryCascade::~G4BinaryCascade ( )
virtual

Definition at line 172 of file G4BinaryCascade.cc.

173 {
174  ClearAndDestroy(&theTargetList);
175  ClearAndDestroy(&theSecondaryList);
176  ClearAndDestroy(&theCapturedList);
177  delete thePropagator;
178  delete theCollisionMgr;
179  std::for_each(theImR.begin(), theImR.end(), Delete<G4BCAction>());
180  delete theLateParticle;
181  //delete theExcitationHandler;
182  delete theH1Scatterer;
183 }

Member Function Documentation

G4HadFinalState * G4BinaryCascade::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus theNucleus 
)
virtual

Implements G4HadronicInteraction.

Definition at line 253 of file G4BinaryCascade.cc.

256 {
257  if(getenv("BCDEBUG") ) G4cerr << " ######### Binary Cascade Reaction starts ######### "<< G4endl;
258 
259  G4LorentzVector initial4Momentum = aTrack.Get4Momentum();
260  const G4ParticleDefinition * definition = aTrack.GetDefinition();
261 
262  if(initial4Momentum.e()-initial4Momentum.m()<theBCminP &&
263  ( definition==G4Neutron::NeutronDefinition() || definition==G4Proton::ProtonDefinition() ) )
264  {
265  return theDeExcitation->ApplyYourself(aTrack, aNucleus);
266  }
267 
269  // initialize the G4V3DNucleus from G4Nucleus
271 
272  // Build a KineticTrackVector with the G4Track
273  G4KineticTrackVector * secondaries;// = new G4KineticTrackVector;
274  G4ThreeVector initialPosition(0., 0., 0.); // will be set later
275 
276  if(!getenv("I_Am_G4BinaryCascade_Developer") )
277  {
278  if(definition!=G4Neutron::NeutronDefinition() &&
279  definition!=G4Proton::ProtonDefinition() &&
280  definition!=G4PionPlus::PionPlusDefinition() &&
281  definition!=G4PionMinus::PionMinusDefinition() )
282  {
283  G4cerr << "You are trying to use G4BinaryCascade with " <<definition->GetParticleName()<<" as projectile."<<G4endl;
284  G4cerr << "G4BinaryCascade should not be used for projectiles other than nucleons or pions."<<G4endl;
285  G4cerr << "If you want to continue, please switch on the developer environment: "<<G4endl;
286  G4cerr << "setenv I_Am_G4BinaryCascade_Developer 1 "<<G4endl<<G4endl;
287  throw G4HadronicException(__FILE__, __LINE__, "G4BinaryCascade - used for unvalid particle type - Fatal");
288  }
289  }
290 
291  // keep primary
292  thePrimaryType = definition;
293  thePrimaryEscape = false;
294 
295  // try until an interaction will happen
296  G4ReactionProductVector * products=0;
297  G4int interactionCounter = 0,collisionLoopMaxCount;
298  do
299  {
300  // reset status that could be changed in previous loop event
301  theCollisionMgr->ClearAndDestroy();
302 
303  if(products != 0)
304  { // free memory from previous loop event
305  ClearAndDestroy(products);
306  delete products;
307  products=0;
308  }
309 
310  G4int massNumber=aNucleus.GetA_asInt();
311  the3DNucleus->Init(massNumber, aNucleus.GetZ_asInt());
312  thePropagator->Init(the3DNucleus);
313  G4KineticTrack * kt;
314  collisionLoopMaxCount = 200;
315  do // sample impact parameter until collisions are found
316  {
317  theCurrentTime=0;
319  initialPosition=GetSpherePoint(1.1*radius, initial4Momentum); // get random position
320  kt = new G4KineticTrack(definition, 0., initialPosition, initial4Momentum);
322  // secondaries has been cleared by Propagate() in the previous loop event
323  secondaries= new G4KineticTrackVector;
324  secondaries->push_back(kt);
325  if(massNumber > 1) // 1H1 is special case
326  {
327  products = Propagate(secondaries, the3DNucleus);
328  } else {
329  products = Propagate1H1(secondaries,the3DNucleus);
330  }
331  // until we FIND a collision ... or give up
332  } while(! products && --collisionLoopMaxCount>0); /* Loop checking, 31.08.2015, G.Folger */
333 
334  if(++interactionCounter>99) break;
335  // ...until we find an ALLOWED collision ... or give up
336  } while(products && products->size() == 0); /* Loop checking, 31.08.2015, G.Folger */
337 
338  if(products && products->size()>0)
339  {
340  // G4cout << "BIC Applyyourself: number of products " << products->size() << G4endl;
341 
342  // Fill the G4ParticleChange * with products
344  G4ReactionProductVector::iterator iter;
345 
346  for(iter = products->begin(); iter != products->end(); ++iter)
347  {
348  G4DynamicParticle * aNew =
349  new G4DynamicParticle((*iter)->GetDefinition(),
350  (*iter)->GetTotalEnergy(),
351  (*iter)->GetMomentum());
353  }
354 
355  //DebugFinalEpConservation(aTrack, products);
356 
357 
358  } else { // no interaction, return primary
359  if(getenv("BCDEBUG") ) G4cerr << " ######### Binary Cascade Reaction void, return intial state ######### "<< G4endl;
363  }
364 
365  if ( products )
366  {
367  ClearAndDestroy(products);
368  delete products;
369  }
370 
371  delete the3DNucleus;
372  the3DNucleus = NULL;
373 
374  if(getenv("BCDEBUG") ) G4cerr << " ######### Binary Cascade Reaction ends ######### "<< G4endl;
375 
376  return &theParticleChange;
377 }
static G4Proton * ProtonDefinition()
Definition: G4Proton.cc:88
int G4int
Definition: G4Types.hh:78
void SetStatusChange(G4HadFinalStateStatus aS)
std::vector< G4ReactionProduct * > G4ReactionProductVector
virtual G4double GetOuterRadius()=0
Hep3Vector vect() const
CascadeState SetState(const CascadeState new_state)
const G4ParticleDefinition * GetDefinition() const
virtual void Init(G4int theA, G4int theZ)=0
static G4PionMinus * PionMinusDefinition()
Definition: G4PionMinus.cc:93
G4double GetKineticEnergy() const
static G4PionPlus * PionPlusDefinition()
Definition: G4PionPlus.cc:93
virtual void Init(G4V3DNucleus *theNucleus)=0
virtual G4ReactionProductVector * Propagate(G4KineticTrackVector *, G4V3DNucleus *)
const G4LorentzVector & Get4Momentum() const
void SetEnergyChange(G4double anEnergy)
Hep3Vector unit() const
#define G4endl
Definition: G4ios.hh:61
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
static constexpr double fermi
Definition: G4SIunits.hh:103
static G4Neutron * NeutronDefinition()
Definition: G4Neutron.cc:99
void SetMomentumChange(const G4ThreeVector &aV)
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

void G4BinaryCascade::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 185 of file G4BinaryCascade.cc.

186 {
187  outFile << "G4BinaryCascade is an intra-nuclear cascade model in which\n"
188  << "an incident hadron collides with a nucleon, forming two\n"
189  << "final-state particles, one or both of which may be resonances.\n"
190  << "The resonances then decay hadronically and the decay products\n"
191  << "are then propagated through the nuclear potential along curved\n"
192  << "trajectories until they re-interact or leave the nucleus.\n"
193  << "This model is valid for incident pions up to 1.5 GeV and\n"
194  << "nucleons up to 10 GeV.\n"
195  << "The remaining excited nucleus is handed on to ";
196  if (theDeExcitation) // pre-compound
197  {
198  outFile << theDeExcitation->GetModelName() << " : \n ";
200  }
201  else if (theExcitationHandler) // de-excitation
202  {
203  outFile << "G4ExcitationHandler"; //theExcitationHandler->GetModelName();
204  theExcitationHandler->ModelDescription(outFile);
205  }
206  else
207  {
208  outFile << "void.\n";
209  }
210  outFile<< " \n";
211 }
void ModelDescription(std::ostream &outFile) const
const G4String & GetModelName() const
virtual void DeExciteModelDescription(std::ostream &outFile) const =0

Here is the call graph for this function:

G4ReactionProductVector * G4BinaryCascade::Propagate ( G4KineticTrackVector secondaries,
G4V3DNucleus aNucleus 
)
virtual

Implements G4VIntraNuclearTransportModel.

Definition at line 379 of file G4BinaryCascade.cc.

382 {
383  G4ping debug("debug_G4BinaryCascade");
384 #ifdef debug_BIC_Propagate
385  G4cout << "G4BinaryCascade Propagate starting -------------------------------------------------------" <<G4endl;
386 #endif
387 
388  the3DNucleus=aNucleus;
390  theOuterRadius = the3DNucleus->GetOuterRadius();
391  theCurrentTime=0;
392  theProjectile4Momentum=G4LorentzVector(0,0,0,0);
393  theMomentumTransfer=G4ThreeVector(0,0,0);
394  // build theSecondaryList, theProjectileList and theCapturedList
395  ClearAndDestroy(&theCapturedList);
396  ClearAndDestroy(&theSecondaryList);
397  theSecondaryList.clear();
398  ClearAndDestroy(&theFinalState);
399  std::vector<G4KineticTrack *>::iterator iter;
400  theCollisionMgr->ClearAndDestroy();
401 
402  theCutOnP=90*MeV;
403  if(the3DNucleus->GetMass()>30) theCutOnP = 70*MeV;
404  if(the3DNucleus->GetMass()>60) theCutOnP = 50*MeV;
405  if(the3DNucleus->GetMass()>120) theCutOnP = 45*MeV;
406 
407 
408  BuildTargetList();
409 
410 #ifdef debug_BIC_GetExcitationEnergy
411  G4cout << "ExcitationEnergy0 " << GetExcitationEnergy() << G4endl;
412 #endif
413 
414  thePropagator->Init(the3DNucleus);
415 
416  G4bool success = BuildLateParticleCollisions(secondaries);
417  if (! success ) // fails if no excitation energy left....
418  {
419  products=HighEnergyModelFSProducts(products, secondaries);
420  ClearAndDestroy(secondaries);
421  delete secondaries;
422 
423 #ifdef debug_G4BinaryCascade
424  G4cout << "G4BinaryCascade::Propagate: warning - high energy model failed energy conservation, returning unchanged high energy final state" << G4endl;
425 #endif
426 
427  return products;
428  }
429  // check baryon and charge ...
430 
431  _CheckChargeAndBaryonNumber_("lateparticles");
432  _DebugEpConservation(" be4 findcollisions");
433 
434  // if called stand alone find first collisions
435  FindCollisions(&theSecondaryList);
436 
437 
438  if(theCollisionMgr->Entries() == 0 ) //late particles ALWAYS create Entries
439  {
440  //G4cout << " no collsions -> return 0" << G4endl;
441  delete products;
442 #ifdef debug_BIC_return
443  G4cout << "return @ begin2, no collisions "<< G4endl;
444 #endif
445  return 0;
446  }
447 
448  // end of initialization: do the job now
449  // loop until there are no more collisions
450 
451 
452  G4bool haveProducts = false;
453  G4int collisionCount=0;
454  G4int collisionLoopMaxCount=1000000;
455  while(theCollisionMgr->Entries() > 0 && currentZ && --collisionLoopMaxCount>0) /* Loop checking, 31.08.2015, G.Folger */
456  {
457  if(Absorb()) { // absorb secondaries, pions only
458  haveProducts = true;
459  }
460  if(Capture()) { // capture secondaries, nucleons only
461  haveProducts = true;
462  }
463 
464  // propagate to the next collision if any (collisions could have been deleted
465  // by previous absorption or capture)
466  if(theCollisionMgr->Entries() > 0)
467  {
469  nextCollision = theCollisionMgr->GetNextCollision();
470 #ifdef debug_BIC_Propagate_Collisions
471  G4cout << " NextCollision * , Time, curtime = " << nextCollision << " "
472  <<nextCollision->GetCollisionTime()<< " " <<
473  theCurrentTime<< G4endl;
474 #endif
475  if (!DoTimeStep(nextCollision->GetCollisionTime()-theCurrentTime) )
476  {
477  // Check if nextCollision is still valid, ie. particle did not leave nucleus
478  if (theCollisionMgr->GetNextCollision() != nextCollision )
479  {
480  nextCollision = 0;
481  }
482  }
483  //_DebugEpConservation("Stepped");
484 
485  if( nextCollision )
486  {
487  if (ApplyCollision(nextCollision))
488  {
489  //G4cerr << "ApplyCollision success " << G4endl;
490  haveProducts = true;
491  collisionCount++;
492  //_CheckChargeAndBaryonNumber_("ApplyCollision");
493  //_DebugEpConservation("ApplyCollision");
494 
495  } else {
496  //G4cerr << "ApplyCollision failure " << G4endl;
497  theCollisionMgr->RemoveCollision(nextCollision);
498  }
499  }
500  }
501  }
502 
503  //--------- end of on Collisions
504  //G4cout << "currentZ @ end loop " << currentZ << G4endl;
505  G4int nProtons(0);
506  for(iter = theTargetList.begin(); iter != theTargetList.end(); ++iter)
507  {
508  if ( (*iter)->GetDefinition() == G4Proton::Proton() ) ++nProtons;
509  }
510  if ( ! theTargetList.size() || ! nProtons ){
511  // nucleus completely destroyed, fill in ReactionProductVector
512  products = FillVoidNucleusProducts(products);
513 #ifdef debug_BIC_return
514  G4cout << "return @ Z=0 after collision loop "<< G4endl;
515  PrintKTVector(&theSecondaryList,std::string(" theSecondaryList"));
516  G4cout << "theTargetList size: " << theTargetList.size() << G4endl;
517  PrintKTVector(&theTargetList,std::string(" theTargetList"));
518  PrintKTVector(&theCapturedList,std::string(" theCapturedList"));
519 
520  G4cout << " ExcitE be4 Correct : " <<GetExcitationEnergy() << G4endl;
521  G4cout << " Mom Transfered to nucleus : " << theMomentumTransfer << " " << theMomentumTransfer.mag() << G4endl;
522  PrintKTVector(&theFinalState,std::string(" FinalState uncorrected"));
523  G4cout << "returned products: " << products->size() << G4endl;
524  _CheckChargeAndBaryonNumber_("destroyed Nucleus");
525  _DebugEpConservation("destroyed Nucleus");
526 #endif
527 
528  return products;
529  }
530 
531  // No more collisions: absorb, capture and propagate the secondaries out of the nucleus
532  if(Absorb()) {
533  haveProducts = true;
534  // G4cout << "Absorb sucess " << G4endl;
535  }
536 
537  if(Capture()) {
538  haveProducts = true;
539  // G4cout << "Capture sucess " << G4endl;
540  }
541 
542  if(!haveProducts) // no collisions happened. Return an empty vector.
543  {
544 #ifdef debug_BIC_return
545  G4cout << "return 3, no products "<< G4endl;
546 #endif
547  return products;
548  }
549 
550 
551 #ifdef debug_BIC_Propagate
552  G4cout << " Momentum transfer to Nucleus " << theMomentumTransfer << " " << theMomentumTransfer.mag() << G4endl;
553  G4cout << " Stepping particles out...... " << G4endl;
554 #endif
555 
556  StepParticlesOut();
557  _DebugEpConservation("stepped out");
558 
559 
560  if ( theSecondaryList.size() > 0 )
561  {
562 #ifdef debug_G4BinaryCascade
563  G4cerr << "G4BinaryCascade: Warning, have active particles at end" << G4endl;
564  PrintKTVector(&theSecondaryList, "active particles @ end added to theFinalState");
565 #endif
566  // add left secondaries to FinalSate
567  for ( iter =theSecondaryList.begin(); iter != theSecondaryList.end(); ++iter)
568  {
569  theFinalState.push_back(*iter);
570  }
571  theSecondaryList.clear();
572 
573  }
574  while ( theCollisionMgr->Entries() > 0 ) /* Loop checking, 31.08.2015, G.Folger */
575  {
576 #ifdef debug_G4BinaryCascade
577  G4cerr << " Warning: remove left over collision(s) " << G4endl;
578 #endif
579  theCollisionMgr->RemoveCollision(theCollisionMgr->GetNextCollision());
580  }
581 
582 #ifdef debug_BIC_Propagate_Excitation
583 
584  PrintKTVector(&theSecondaryList,std::string(" theSecondaryList"));
585  G4cout << "theTargetList size: " << theTargetList.size() << G4endl;
586  // PrintKTVector(&theTargetList,std::string(" theTargetList"));
587  PrintKTVector(&theCapturedList,std::string(" theCapturedList"));
588 
589  G4cout << " ExcitE be4 Correct : " <<GetExcitationEnergy() << G4endl;
590  G4cout << " Mom Transfered to nucleus : " << theMomentumTransfer << " " << theMomentumTransfer.mag() << G4endl;
591  PrintKTVector(&theFinalState,std::string(" FinalState uncorrected"));
592 #endif
593 
594  //
595 
596 
597  G4double ExcitationEnergy=GetExcitationEnergy();
598 
599 #ifdef debug_BIC_Propagate_finals
600  PrintKTVector(&theFinalState,std::string(" FinalState be4 corr"));
601  G4cout << " Excitation Energy prefinal, #collisions:, out, captured "
602  << ExcitationEnergy << " "
603  << collisionCount << " "
604  << theFinalState.size() << " "
605  << theCapturedList.size()<<G4endl;
606 #endif
607 
608  if (ExcitationEnergy < 0 )
609  {
610  G4int maxtry=5, ntry=0;
611  do {
612  CorrectFinalPandE();
613  ExcitationEnergy=GetExcitationEnergy();
614  } while ( ++ntry < maxtry && ExcitationEnergy < 0 ); /* Loop checking, 31.08.2015, G.Folger */
615  }
616  _DebugEpConservation("corrected");
617 
618 #ifdef debug_BIC_Propagate_finals
619  PrintKTVector(&theFinalState,std::string(" FinalState corrected"));
620  G4cout << " Excitation Energy final, #collisions:, out, captured "
621  << ExcitationEnergy << " "
622  << collisionCount << " "
623  << theFinalState.size() << " "
624  << theCapturedList.size()<<G4endl;
625 #endif
626 
627 
628  if ( ExcitationEnergy < 0. )
629  {
630  #ifdef debug_G4BinaryCascade
631  G4cerr << "G4BinaryCascade-Warning: negative excitation energy ";
632  G4cerr <<ExcitationEnergy<<G4endl;
633  PrintKTVector(&theFinalState,std::string("FinalState"));
634  PrintKTVector(&theCapturedList,std::string("captured"));
635  G4cout << "negative ExE:Final 4Momentum .mag: " << GetFinal4Momentum()
636  << " "<< GetFinal4Momentum().mag()<< G4endl
637  << "negative ExE:FinalNucleusMom .mag: " << GetFinalNucleusMomentum()
638  << " "<< GetFinalNucleusMomentum().mag()<< G4endl;
639  #endif
640  #ifdef debug_BIC_return
641  G4cout << " negative Excitation E return empty products " << products << G4endl;
642  G4cout << "return 4, excit < 0 "<< G4endl;
643  #endif
644 
645  ClearAndDestroy(products);
646  return products; // return empty products- FIXME
647  }
648 
649  G4ReactionProductVector * precompoundProducts=DeExcite();
650 
651 
652  G4DecayKineticTracks decay(&theFinalState);
653  _DebugEpConservation("decayed");
654 
655  products= ProductsAddFinalState(products, theFinalState);
656 
657  products= ProductsAddPrecompound(products, precompoundProducts);
658 
659 // products=ProductsAddFakeGamma(products);
660 
661 
662  thePrimaryEscape = true;
663 
664  #ifdef debug_BIC_return
665  G4cout << "BIC: return @end, all ok "<< G4endl;
666  //G4cout << " return products " << products << G4endl;
667  #endif
668 
669  return products;
670 }
Definition: G4ping.hh:33
CLHEP::Hep3Vector G4ThreeVector
void RemoveCollision(G4CollisionInitialState *collision)
G4CollisionInitialState * GetNextCollision()
#define _DebugEpConservation(val)
#define _CheckChargeAndBaryonNumber_(val)
int G4int
Definition: G4Types.hh:78
std::vector< G4ReactionProduct * > G4ReactionProductVector
virtual G4double GetOuterRadius()=0
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
double mag() const
bool G4bool
Definition: G4Types.hh:79
static G4Proton * Proton()
Definition: G4Proton.cc:93
virtual void Init(G4V3DNucleus *theNucleus)=0
#define debug
virtual G4double GetMass()=0
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76
double mag() const
G4GLOB_DLL std::ostream G4cerr
CLHEP::HepLorentzVector G4LorentzVector

Here is the call graph for this function:

Here is the caller graph for this function:

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

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 212 of file G4BinaryCascade.cc.

213 {
214  outFile << "G4BinaryCascade propagtes secondaries produced by a high\n"
215  << "energy model through the wounded nucleus.\n"
216  << "Secondaries are followed after the formation time and if\n"
217  << "within the nucleus are propagated through the nuclear\n"
218  << "potential along curved trajectories until they interact\n"
219  << "with a nucleon, decay, or leave the nucleus.\n"
220  << "An interaction of a secondary with a nucleon produces two\n"
221  << "final-state particles, one or both of which may be resonances.\n"
222  << "Resonances decay hadronically and the decay products\n"
223  << "are in turn propagated through the nuclear potential along curved\n"
224  << "trajectories until they re-interact or leave the nucleus.\n"
225  << "This model is valid for pions up to 1.5 GeV and\n"
226  << "nucleons up to about 3.5 GeV.\n"
227  << "The remaining excited nucleus is handed on to ";
228  if (theDeExcitation) // pre-compound
229  {
230  outFile << theDeExcitation->GetModelName() << " : \n ";
232  }
233  else if (theExcitationHandler) // de-excitation
234  {
235  outFile << "G4ExcitationHandler"; //theExcitationHandler->GetModelName();
236  theExcitationHandler->ModelDescription(outFile);
237  }
238  else
239  {
240  outFile << "void.\n";
241  }
242 outFile<< " \n";
243 }
void ModelDescription(std::ostream &outFile) const
const G4String & GetModelName() const
virtual void DeExciteModelDescription(std::ostream &outFile) const =0

Here is the call graph for this function:


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