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

#include <G4CascadeInterface.hh>

Inheritance diagram for G4CascadeInterface:
Collaboration diagram for G4CascadeInterface:

Public Member Functions

 G4CascadeInterface (const G4String &name="BertiniCascade")
 
virtual ~G4CascadeInterface ()
 
G4ReactionProductVectorPropagate (G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus)
 
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
void SetVerboseLevel (G4int verbose)
 
G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
G4bool IsApplicable (const G4ParticleDefinition *aPD) const
 
void useCascadeDeexcitation ()
 
void usePreCompoundDeexcitation ()
 
virtual void ModelDescription (std::ostream &outFile) const
 
virtual void DumpConfiguration (std::ostream &outFile) 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
 
virtual void PropagateModelDescription (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)
 
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 ()
 

Static Public Member Functions

static void Initialize ()
 

Protected Member Functions

void clear ()
 
G4bool createBullet (const G4HadProjectile &aTrack)
 
G4bool createTarget (G4Nucleus &theNucleus)
 
G4bool createTarget (G4V3DNucleus *theNucleus)
 
G4bool createTarget (G4int A, G4int Z)
 
G4bool coulombBarrierViolation () const
 
G4bool retryInelasticProton () const
 
G4bool retryInelasticNucleus () const
 
void copyOutputToHadronicResult ()
 
G4ReactionProductVectorcopyOutputToReactionProducts ()
 
G4HadFinalStateNoInteraction (const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
 
void checkFinalResult ()
 
void throwNonConservationFailure ()
 
G4DynamicParticlemakeDynamicParticle (const G4InuclElementaryParticle &iep) const
 
G4DynamicParticlemakeDynamicParticle (const G4InuclNuclei &inuc) const
 
- 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 ()
 

Additional Inherited Members

- 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 88 of file G4CascadeInterface.hh.

Constructor & Destructor Documentation

G4CascadeInterface::G4CascadeInterface ( const G4String name = "BertiniCascade")

Definition at line 146 of file G4CascadeInterface.cc.

148  randomFile(G4CascadeParameters::randomFile()),
149  maximumTries(20), numberOfTries(0),
150  collider(new G4InuclCollider), balance(new G4CascadeCheckBalance(name)),
151  bullet(0), target(0), output(new G4CollisionOutput) {
152  // Set up global objects for master thread or sequential build
154 
156  balance->setLimits(5*perCent, 10*MeV/GeV); // Bertini internal units
158 
161  else
163 }
const XML_Char * target
Definition: expat.h:268
static constexpr double perCent
Definition: G4SIunits.hh:332
G4VIntraNuclearTransportModel(const G4String &mName="CascadeModel", G4VPreCompoundModel *ptr=nullptr)
void setLimits(G4double relative, G4double absolute)
static const G4String & randomFile()
static constexpr double GeV
Definition: G4SIunits.hh:217
static constexpr double MeV
Definition: G4SIunits.hh:214
G4bool IsMasterThread()
Definition: G4Threading.cc:146
void SetVerboseLevel(G4int verbose)
void SetEnergyMomentumCheckLevels(G4double relativeLevel, G4double absoluteLevel)
static G4bool usePreCompound()

Here is the call graph for this function:

G4CascadeInterface::~G4CascadeInterface ( )
virtual

Definition at line 165 of file G4CascadeInterface.cc.

165  {
166  clear();
167  delete collider; collider=0;
168  delete balance; balance=0;
169  delete output; output=0;
170 }

Here is the call graph for this function:

Member Function Documentation

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

Implements G4HadronicInteraction.

Definition at line 248 of file G4CascadeInterface.cc.

249  {
250  if (verboseLevel)
251  G4cout << " >>> G4CascadeInterface::ApplyYourself" << G4endl;
252 
253  if (aTrack.GetKineticEnergy() < 0.) {
254  G4cerr << " >>> G4CascadeInterface got negative-energy track: "
255  << aTrack.GetDefinition()->GetParticleName() << " Ekin = "
256  << aTrack.GetKineticEnergy() << G4endl;
257  }
258 
259 #ifdef G4CASCADE_DEBUG_INTERFACE
260  static G4int counter(0);
261  counter++;
262  G4cerr << "Reaction number "<< counter << " "
263  << aTrack.GetDefinition()->GetParticleName() << " "
264  << aTrack.GetKineticEnergy() << " MeV" << G4endl;
265 #endif
266 
267  if (!randomFile.empty()) { // User requested random-seed capture
268  if (verboseLevel>1)
269  G4cout << " Saving random engine state to " << randomFile << G4endl;
271  }
272 
274  clear();
275 
276  // Abort processing if no interaction is possible
277  if (!IsApplicable(aTrack, theNucleus)) {
278  if (verboseLevel) G4cerr << " No interaction possible " << G4endl;
279  return NoInteraction(aTrack, theNucleus);
280  }
281 
282  // Make conversion between native Geant4 and Bertini cascade classes.
283  if (!createBullet(aTrack)) {
284  if (verboseLevel) G4cerr << " Unable to create usable bullet" << G4endl;
285  return NoInteraction(aTrack, theNucleus);
286  }
287 
288  if (!createTarget(theNucleus)) {
289  if (verboseLevel) G4cerr << " Unable to create usable target" << G4endl;
290  return NoInteraction(aTrack, theNucleus);
291  }
292 
293  // Different retry conditions for proton target vs. nucleus
294  const G4bool isHydrogen = (theNucleus.GetA_asInt() == 1);
295 
296  numberOfTries = 0;
297  do { // we try to create inelastic interaction
298  if (verboseLevel > 1)
299  G4cout << " Generating cascade attempt " << numberOfTries << G4endl;
300 
301  output->reset();
302  collider->collide(bullet, target, *output);
303  balance->collide(bullet, target, *output);
304 
305  numberOfTries++;
306  /* Loop checking 08.06.2015 MHK */
307  } while ( isHydrogen ? retryInelasticProton() : retryInelasticNucleus() );
308 
309  // Null event if unsuccessful
310  if (numberOfTries >= maximumTries) {
311  if (verboseLevel)
312  G4cout << " Cascade aborted after trials " << numberOfTries << G4endl;
313  return NoInteraction(aTrack, theNucleus);
314  }
315 
316  // Abort job if energy or momentum are not conserved
317  if (!balance->okay()) {
319  return NoInteraction(aTrack, theNucleus);
320  }
321 
322  // Successful cascade -- clean up and return
323  if (verboseLevel) {
324  G4cout << " Cascade output after trials " << numberOfTries << G4endl;
325  if (verboseLevel > 1) output->printCollisionOutput();
326  }
327 
328  // Rotate event to put Z axis along original projectile direction
329  output->rotateEvent(bulletInLabFrame);
330 
332 
333  // Report violations of conservation laws in original frame
335 
336  // Clean up and return final result;
337  clear();
338 /*
339  G4int nSec = theParticleChange.GetNumberOfSecondaries();
340  for (G4int i = 0; i < nSec; i++) {
341  G4HadSecondary* sec = theParticleChange.GetSecondary(i);
342  G4DynamicParticle* dp = sec->GetParticle();
343  if (dp->GetDefinition()->GetParticleName() == "neutron")
344  G4cout << dp->GetDefinition()->GetParticleName() << " has "
345  << dp->GetKineticEnergy()/MeV << " MeV " << G4endl;
346  }
347 */
348  return &theParticleChange;
349 }
G4bool retryInelasticNucleus() const
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
const XML_Char * target
Definition: expat.h:268
void collide(G4InuclParticle *bullet, G4InuclParticle *target, G4CollisionOutput &output)
void printCollisionOutput(std::ostream &os=G4cout) const
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * GetDefinition() const
bool G4bool
Definition: G4Types.hh:79
G4double GetKineticEnergy() const
G4bool createTarget(G4Nucleus &theNucleus)
void collide(G4InuclParticle *bullet, G4InuclParticle *target, G4CollisionOutput &globalOutput)
void rotateEvent(const G4LorentzRotation &rotate)
static void saveEngineStatus(const char filename[]="Config.conf")
Definition: Random.cc:275
#define G4endl
Definition: G4ios.hh:61
G4bool createBullet(const G4HadProjectile &aTrack)
G4bool retryInelasticProton() const
G4GLOB_DLL std::ostream G4cerr
G4bool IsApplicable(const G4HadProjectile &aTrack, G4Nucleus &theNucleus)
G4HadFinalState * NoInteraction(const G4HadProjectile &aTrack, G4Nucleus &theNucleus)

Here is the call graph for this function:

void G4CascadeInterface::checkFinalResult ( )
protected

Definition at line 626 of file G4CascadeInterface.cc.

626  {
627  balance->collide(bullet, target, *output);
628 
629  if (verboseLevel > 2) {
630  if (!balance->baryonOkay()) {
631  G4cerr << "ERROR: no baryon number conservation, sum of baryons = "
632  << balance->deltaB() << G4endl;
633  }
634 
635  if (!balance->chargeOkay()) {
636  G4cerr << "ERROR: no charge conservation, sum of charges = "
637  << balance->deltaQ() << G4endl;
638  }
639 
640  if (std::abs(balance->deltaKE()) > 0.01 ) { // GeV
641  G4cerr << "Kinetic energy conservation violated by "
642  << balance->deltaKE() << " GeV" << G4endl;
643  }
644 
645  G4double eInit = bullet->getEnergy() + target->getEnergy();
646  G4double eFinal = eInit + balance->deltaE();
647 
648  G4cout << "Initial energy " << eInit << " final energy " << eFinal
649  << "\nTotal energy conservation at level "
650  << balance->deltaE() * GeV << " MeV" << G4endl;
651 
652  if (balance->deltaKE() > 5.0e-5 ) { // 0.05 MeV
653  G4cerr << "FATAL ERROR: kinetic energy created "
654  << balance->deltaKE() * GeV << " MeV" << G4endl;
655  }
656  }
657 }
const XML_Char * target
Definition: expat.h:268
void collide(G4InuclParticle *bullet, G4InuclParticle *target, G4CollisionOutput &output)
G4double getEnergy() const
G4GLOB_DLL std::ostream G4cout
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

Here is the caller graph for this function:

void G4CascadeInterface::clear ( void  )
protected

Definition at line 190 of file G4CascadeInterface.cc.

190  {
191  bullet=0;
192  target=0;
193 }
const XML_Char * target
Definition: expat.h:268

Here is the caller graph for this function:

void G4CascadeInterface::copyOutputToHadronicResult ( )
protected

Definition at line 557 of file G4CascadeInterface.cc.

557  {
558  if (verboseLevel > 1)
559  G4cout << " >>> G4CascadeInterface::copyOutputToHadronicResult" << G4endl;
560 
561  const std::vector<G4InuclNuclei>& outgoingNuclei = output->getOutgoingNuclei();
562  const std::vector<G4InuclElementaryParticle>& particles = output->getOutgoingParticles();
563 
566 
567  // Get outcoming particles
568  if (!particles.empty()) {
569  particleIterator ipart = particles.begin();
570  for (; ipart != particles.end(); ipart++) {
572  }
573  }
574 
575  // get nuclei fragments
576  if (!outgoingNuclei.empty()) {
577  nucleiIterator ifrag = outgoingNuclei.begin();
578  for (; ifrag != outgoingNuclei.end(); ifrag++) {
580  }
581  }
582 }
void SetStatusChange(G4HadFinalStateStatus aS)
G4GLOB_DLL std::ostream G4cout
G4DynamicParticle * makeDynamicParticle(const G4InuclElementaryParticle &iep) const
void SetEnergyChange(G4double anEnergy)
const std::vector< G4InuclNuclei > & getOutgoingNuclei() const
const std::vector< G4InuclElementaryParticle > & getOutgoingParticles() const
std::vector< G4InuclNuclei >::const_iterator nucleiIterator
std::vector< G4InuclElementaryParticle >::iterator particleIterator
Definition: G4BigBanger.cc:65
#define G4endl
Definition: G4ios.hh:61
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)

Here is the call graph for this function:

Here is the caller graph for this function:

G4ReactionProductVector * G4CascadeInterface::copyOutputToReactionProducts ( )
protected

Definition at line 584 of file G4CascadeInterface.cc.

584  {
585  if (verboseLevel > 1)
586  G4cout << " >>> G4CascadeInterface::copyOutputToReactionProducts" << G4endl;
587 
588  const std::vector<G4InuclElementaryParticle>& particles = output->getOutgoingParticles();
589  const std::vector<G4InuclNuclei>& fragments = output->getOutgoingNuclei();
590 
592 
593  G4ReactionProduct* rp = 0; // Buffers to create outgoing tracks
594  G4DynamicParticle* dp = 0;
595 
596  // Get outcoming particles
597  if (!particles.empty()) {
598  particleIterator ipart = particles.begin();
599  for (; ipart != particles.end(); ipart++) {
600  rp = new G4ReactionProduct;
601  dp = makeDynamicParticle(*ipart);
602  (*rp) = (*dp); // This does all the necessary copying
603  propResult->push_back(rp);
604  delete dp;
605  }
606  }
607 
608  // get nuclei fragments
609  if (!fragments.empty()) {
610  nucleiIterator ifrag = fragments.begin();
611  for (; ifrag != fragments.end(); ifrag++) {
612  rp = new G4ReactionProduct;
613  dp = makeDynamicParticle(*ifrag);
614  (*rp) = (*dp); // This does all the necessary copying
615  propResult->push_back(rp);
616  delete dp;
617  }
618  }
619 
620  return propResult;
621 }
std::vector< G4ReactionProduct * > G4ReactionProductVector
G4GLOB_DLL std::ostream G4cout
G4DynamicParticle * makeDynamicParticle(const G4InuclElementaryParticle &iep) const
const std::vector< G4InuclNuclei > & getOutgoingNuclei() const
const std::vector< G4InuclElementaryParticle > & getOutgoingParticles() const
std::vector< G4InuclNuclei >::const_iterator nucleiIterator
std::vector< G4InuclElementaryParticle >::iterator particleIterator
Definition: G4BigBanger.cc:65
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::coulombBarrierViolation ( ) const
protected

Definition at line 662 of file G4CascadeInterface.cc.

662  {
663  G4bool violated = false; // by default coulomb analysis is OK
664 
665  const G4double coulumbBarrier = 8.7 * MeV/GeV; // Bertini uses GeV
666 
667  const std::vector<G4InuclElementaryParticle>& p =
668  output->getOutgoingParticles();
669 
670  for (particleIterator ipart=p.begin(); ipart != p.end(); ipart++) {
671  if (ipart->type() == proton) {
672  violated |= (ipart->getKineticEnergy() < coulumbBarrier);
673  }
674  }
675 
676  return violated;
677 }
const char * p
Definition: xmltok.h:285
bool G4bool
Definition: G4Types.hh:79
const std::vector< G4InuclElementaryParticle > & getOutgoingParticles() const
static constexpr double GeV
Definition: G4SIunits.hh:217
std::vector< G4InuclElementaryParticle >::iterator particleIterator
Definition: G4BigBanger.cc:65
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::createBullet ( const G4HadProjectile aTrack)
protected

Definition at line 444 of file G4CascadeInterface.cc.

444  {
445  const G4ParticleDefinition* trkDef = aTrack.GetDefinition();
446  G4int bulletType = 0; // For elementary particles
447  G4int bulletA = 0, bulletZ = 0; // For nucleus projectile
448 
449  if (trkDef->GetAtomicMass() <= 1) {
450  bulletType = G4InuclElementaryParticle::type(trkDef);
451  } else {
452  bulletA = trkDef->GetAtomicMass();
453  bulletZ = trkDef->GetAtomicNumber();
454  }
455 
456  if (0 == bulletType && 0 == bulletA*bulletZ) {
457  if (verboseLevel) {
458  G4cerr << " G4CascadeInterface: " << trkDef->GetParticleName()
459  << " not usable as bullet." << G4endl;
460  }
461  bullet = 0;
462  return false;
463  }
464 
465  // Code momentum and energy -- Bertini wants z-axis and GeV units
466  G4LorentzVector projectileMomentum = aTrack.Get4Momentum()/GeV;
467 
468  // Rotation/boost to get from z-axis back to original frame
469  bulletInLabFrame = G4LorentzRotation::IDENTITY; // Initialize
470  bulletInLabFrame.rotateZ(-projectileMomentum.phi());
471  bulletInLabFrame.rotateY(-projectileMomentum.theta());
472  bulletInLabFrame.invert();
473 
474  G4LorentzVector momentumBullet(0., 0., projectileMomentum.rho(),
475  projectileMomentum.e());
476 
477  if (G4InuclElementaryParticle::valid(bulletType)) {
478  hadronBullet.fill(momentumBullet, bulletType);
479  bullet = &hadronBullet;
480  } else {
481  nucleusBullet.fill(momentumBullet, bulletA, bulletZ);
482  bullet = &nucleusBullet;
483  }
484 
485  if (verboseLevel > 2) G4cout << "Bullet: \n" << *bullet << G4endl;
486 
487  return true;
488 }
void fill(G4int a, G4int z, G4double exc=0., Model model=DefaultModel)
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
HepLorentzRotation & rotateY(double delta)
G4int GetAtomicNumber() const
double phi() const
G4GLOB_DLL std::ostream G4cout
const G4ParticleDefinition * GetDefinition() const
double theta() const
double rho() const
const G4LorentzVector & Get4Momentum() const
G4int GetAtomicMass() const
HepLorentzRotation & rotateZ(double delta)
HepLorentzRotation & invert()
static constexpr double GeV
Definition: G4SIunits.hh:217
void fill(G4int ityp, Model model=DefaultModel)
#define G4endl
Definition: G4ios.hh:61
static DLL_API const HepLorentzRotation IDENTITY
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::createTarget ( G4Nucleus theNucleus)
protected

Definition at line 493 of file G4CascadeInterface.cc.

493  {
494  return createTarget(theNucleus.GetA_asInt(), theNucleus.GetZ_asInt());
495 }
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
G4bool createTarget(G4Nucleus &theNucleus)
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::createTarget ( G4V3DNucleus theNucleus)
protected

Definition at line 497 of file G4CascadeInterface.cc.

497  {
498  return createTarget(theNucleus->GetMassNumber(), theNucleus->GetCharge());
499 }
virtual G4int GetCharge()=0
virtual G4int GetMassNumber()=0
G4bool createTarget(G4Nucleus &theNucleus)

Here is the call graph for this function:

G4bool G4CascadeInterface::createTarget ( G4int  A,
G4int  Z 
)
protected

Definition at line 501 of file G4CascadeInterface.cc.

501  {
502  if (A > 1) {
503  nucleusTarget.fill(A, Z);
504  target = &nucleusTarget;
505  } else {
506  hadronTarget.fill(0., (Z==1?proton:neutron));
507  target = &hadronTarget;
508  }
509 
510  if (verboseLevel > 2) G4cout << "Target: \n" << *target << G4endl;
511 
512  return true; // Right now, target never fails
513 }
void fill(G4int a, G4int z, G4double exc=0., Model model=DefaultModel)
const XML_Char * target
Definition: expat.h:268
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
void fill(G4int ityp, Model model=DefaultModel)
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

void G4CascadeInterface::DumpConfiguration ( std::ostream &  outFile) const
virtual

Definition at line 186 of file G4CascadeInterface.cc.

186  {
188 }
static void DumpConfiguration(std::ostream &os)

Here is the call graph for this function:

void G4CascadeInterface::Initialize ( )
static

Definition at line 198 of file G4CascadeInterface.cc.

198  {
203  if (!ch || !pn || !nn || !pp) return; // Avoid "unused variables"
204 }
static const G4CascadeChannel * GetTable(G4int initialState)
static G4Diproton * Definition()
Definition: G4Diproton.cc:68
static G4Dineutron * Definition()
Definition: G4Dineutron.cc:68
static G4UnboundPN * Definition()
Definition: G4UnboundPN.cc:67

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::IsApplicable ( const G4HadProjectile aTrack,
G4Nucleus theNucleus 
)
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 231 of file G4CascadeInterface.cc.

232  {
233  return IsApplicable(aTrack.GetDefinition());
234 }
const G4ParticleDefinition * GetDefinition() const
G4bool IsApplicable(const G4HadProjectile &aTrack, G4Nucleus &theNucleus)

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::IsApplicable ( const G4ParticleDefinition aPD) const

Definition at line 236 of file G4CascadeInterface.cc.

236  {
237  if (aPD->GetAtomicMass() > 1) return true; // Nuclei are okay
238 
239  // Valid particle and have interactions available
241  return (G4CascadeChannelTables::GetTable(type));
242 }
static const G4CascadeChannel * GetTable(G4int initialState)
int G4int
Definition: G4Types.hh:78
G4int GetAtomicMass() const

Here is the call graph for this function:

G4DynamicParticle * G4CascadeInterface::makeDynamicParticle ( const G4InuclElementaryParticle iep) const
protected

Definition at line 519 of file G4CascadeInterface.cc.

519  {
520  G4int outgoingType = iep.type();
521 
522  if (iep.quasi_deutron()) {
523  G4cerr << " ERROR: G4CascadeInterface incompatible particle type "
524  << outgoingType << G4endl;
525  return 0;
526  }
527 
528  // Copy local G4DynPart to public output (handle kaon mixing specially)
529  if (outgoingType == kaonZero || outgoingType == kaonZeroBar) {
530  G4ThreeVector momDir = iep.getMomentum().vect().unit();
531  G4double ekin = iep.getKineticEnergy()*GeV; // Bertini -> G4 units
532 
534  if (G4UniformRand() > 0.5) pd = G4KaonZeroLong::Definition();
535 
536  return new G4DynamicParticle(pd, momDir, ekin);
537  } else {
538  return new G4DynamicParticle(iep.getDynamicParticle());
539  }
540 
541  return 0; // Should never get here!
542 }
const G4DynamicParticle & getDynamicParticle() const
G4LorentzVector getMomentum() const
int G4int
Definition: G4Types.hh:78
G4double getKineticEnergy() const
Hep3Vector vect() const
#define G4UniformRand()
Definition: Randomize.hh:97
static G4KaonZeroLong * Definition()
Hep3Vector unit() const
static G4KaonZeroShort * Definition()
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

Here is the caller graph for this function:

G4DynamicParticle * G4CascadeInterface::makeDynamicParticle ( const G4InuclNuclei inuc) const
protected

Definition at line 545 of file G4CascadeInterface.cc.

545  {
546  if (verboseLevel > 2) {
547  G4cout << " Nuclei fragment: \n" << inuc << G4endl;
548  }
549 
550  // Copy local G4DynPart to public output
551  return new G4DynamicParticle(inuc.getDynamicParticle());
552 }
const G4DynamicParticle & getDynamicParticle() const
G4GLOB_DLL std::ostream G4cout
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

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

Reimplemented from G4VIntraNuclearTransportModel.

Definition at line 172 of file G4CascadeInterface.cc.

173 {
174  outFile << "The Bertini-style cascade implements the inelastic scattering\n"
175  << "of hadrons by nuclei. Nucleons, pions, kaons and hyperons\n"
176  << "from 0 to 15 GeV may be used as projectiles in this model.\n"
177  << "Final state hadrons are produced by a classical cascade\n"
178  << "consisting of individual hadron-nucleon scatterings which use\n"
179  << "free-space partial cross sections, corrected for various\n"
180  << "nuclear medium effects. The target nucleus is modeled as a\n"
181  << "set of 1, 3 or 6 spherical shells, in which scattered hadrons\n"
182  << "travel in straight lines until they are reflected from or\n"
183  << "transmitted through shell boundaries.\n";
184 }
G4HadFinalState * G4CascadeInterface::NoInteraction ( const G4HadProjectile aTrack,
G4Nucleus theNucleus 
)
protected

Definition at line 427 of file G4CascadeInterface.cc.

428  {
429  if (verboseLevel)
430  G4cout << " >>> G4CascadeInterface::NoInteraction" << G4endl;
431 
434 
435  G4double ekin = aTrack.GetKineticEnergy()>0. ? aTrack.GetKineticEnergy() : 0.;
436  theParticleChange.SetEnergyChange(ekin); // Protect against rounding
437 
438  return &theParticleChange;
439 }
void SetStatusChange(G4HadFinalStateStatus aS)
G4GLOB_DLL std::ostream G4cout
G4double GetKineticEnergy() const
void SetEnergyChange(G4double anEnergy)
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4ReactionProductVector * G4CascadeInterface::Propagate ( G4KineticTrackVector theSecondaries,
G4V3DNucleus theNucleus 
)
virtual

Implements G4VIntraNuclearTransportModel.

Definition at line 352 of file G4CascadeInterface.cc.

353  {
354  if (verboseLevel) G4cout << " >>> G4CascadeInterface::Propagate" << G4endl;
355 
356 #ifdef G4CASCADE_DEBUG_INTERFACE
357  if (verboseLevel>1) {
358  G4cout << " G4V3DNucleus A " << theNucleus->GetMassNumber()
359  << " Z " << theNucleus->GetCharge()
360  << "\n " << theSecondaries->size() << " secondaries:" << G4endl;
361  for (size_t i=0; i<theSecondaries->size(); i++) {
362  G4KineticTrack* kt = (*theSecondaries)[i];
363  G4cout << " " << i << ": " << kt->GetDefinition()->GetParticleName()
364  << " p " << kt->Get4Momentum() << " @ " << kt->GetPosition()
365  << " t " << kt->GetFormationTime() << G4endl;
366  }
367  }
368 #endif
369 
370  if (!randomFile.empty()) { // User requested random-seed capture
371  if (verboseLevel>1)
372  G4cout << " Saving random engine state to " << randomFile << G4endl;
374  }
375 
377  clear();
378 
379  // Process input secondaries list to eliminate resonances
380  G4DecayKineticTracks decay(theSecondaries);
381 
382  // NOTE: Requires 9.4-ref-03 mods to base class and G4TheoFSGenerator
383  const G4HadProjectile* projectile = GetPrimaryProjectile();
384  if (projectile) createBullet(*projectile);
385 
386  if (!createTarget(theNucleus)) {
387  if (verboseLevel) G4cerr << " Unable to create usable target" << G4endl;
388  return 0; // FIXME: This will cause a segfault later
389  }
390 
391  numberOfTries = 0;
392  do {
393  if (verboseLevel > 1)
394  G4cout << " Generating rescatter attempt " << numberOfTries << G4endl;
395 
396  output->reset();
397  collider->rescatter(bullet, theSecondaries, theNucleus, *output);
398  balance->collide(bullet, target, *output);
399 
400  numberOfTries++;
401  // FIXME: retry checks will SEGFAULT until we can define the bullet!
402  } while (retryInelasticNucleus()); /* Loop checking 08.06.2015 MHK */
403 
404  // Check whether repeated attempts have all failed; report and exit
405  if (numberOfTries >= maximumTries && !balance->okay()) {
406  throwNonConservationFailure(); // This terminates the job
407  }
408 
409  // Successful cascade -- clean up and return
410  if (verboseLevel) {
411  G4cout << " Cascade rescatter after trials " << numberOfTries << G4endl;
412  if (verboseLevel > 1) output->printCollisionOutput();
413  }
414 
415  // Does calling code take ownership? I hope so!
417 
418  // Clean up and and return final result
419  clear();
420  return propResult;
421 }
G4bool retryInelasticNucleus() const
void rescatter(G4InuclParticle *bullet, G4KineticTrackVector *theSecondaries, G4V3DNucleus *theNucleus, G4CollisionOutput &globalOutput)
const XML_Char * target
Definition: expat.h:268
virtual G4int GetCharge()=0
const G4HadProjectile * GetPrimaryProjectile() const
const G4ThreeVector & GetPosition() const
void collide(G4InuclParticle *bullet, G4InuclParticle *target, G4CollisionOutput &output)
void printCollisionOutput(std::ostream &os=G4cout) const
virtual G4int GetMassNumber()=0
const G4String & GetParticleName() const
std::vector< G4ReactionProduct * > G4ReactionProductVector
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
G4double GetFormationTime() const
G4bool createTarget(G4Nucleus &theNucleus)
G4ReactionProductVector * copyOutputToReactionProducts()
static void saveEngineStatus(const char filename[]="Config.conf")
Definition: Random.cc:275
#define G4endl
Definition: G4ios.hh:61
G4bool createBullet(const G4HadProjectile &aTrack)
const G4LorentzVector & Get4Momentum() const
const G4ParticleDefinition * GetDefinition() const
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

G4bool G4CascadeInterface::retryInelasticNucleus ( ) const
protected

Definition at line 713 of file G4CascadeInterface.cc.

713  {
714  // Quantities necessary for conditional block below
715  G4int npart = output->numberOfOutgoingParticles();
716  G4int nfrag = output->numberOfOutgoingNuclei();
717 
718  const G4ParticleDefinition* firstOut = (npart == 0) ? 0 :
719  output->getOutgoingParticles().begin()->getDefinition();
720 
721 #ifdef G4CASCADE_DEBUG_INTERFACE
722  // Report on all retry conditions, in order of return logic
723  G4cout << " retryInelasticNucleus: numberOfTries "
724  << ((numberOfTries < maximumTries) ? "RETRY (t)" : "EXIT (f)")
725  << "\n retryInelasticNucleus: AND outputParticles "
726  << ((npart != 0) ? "NON-ZERO (t)" : "EMPTY (f)")
727 #ifdef G4CASCADE_COULOMB_DEV
728  << "\n retryInelasticNucleus: AND coulombBarrier (COULOMB_DEV) "
729  << (coulombBarrierViolation() ? "VIOLATED (t)" : "PASSED (f)")
730  << "\n retryInelasticNucleus: AND collision type (COULOMB_DEV) "
731  << ((npart+nfrag > 2) ? "INELASTIC (t)" : "ELASTIC (f)")
732 #else
733  << "\n retryInelasticNucleus: AND collsion type "
734  << ((npart+nfrag < 3) ? "ELASTIC (t)" : "INELASTIC (f)")
735  << "\n retryInelasticNucleus: AND Leading particle bullet "
736  << ((firstOut == bullet->getDefinition()) ? "YES (t)" : "NO (f)")
737 #endif
738  << "\n retryInelasticNucleus: OR conservation "
739  << (!balance->okay() ? "FAILED (t)" : "PASSED (f)")
740  << G4endl;
741 #endif
742 
743  return ( (numberOfTries < maximumTries) &&
744  ( ((npart != 0) &&
745 #ifdef G4CASCADE_COULOMB_DEV
746  (coulombBarrierViolation() && npart+nfrag > 2)
747 #else
748  (npart+nfrag < 3 && firstOut == bullet->getDefinition())
749 #endif
750  )
751 #ifndef G4CASCADE_SKIP_ECONS
752  || (!balance->okay())
753 #endif
754  )
755  );
756 }
G4bool coulombBarrierViolation() const
const G4ParticleDefinition * getDefinition() const
int G4int
Definition: G4Types.hh:78
G4GLOB_DLL std::ostream G4cout
G4int numberOfOutgoingParticles() const
G4int numberOfOutgoingNuclei() const
const std::vector< G4InuclElementaryParticle > & getOutgoingParticles() const
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

G4bool G4CascadeInterface::retryInelasticProton ( ) const
protected

Definition at line 681 of file G4CascadeInterface.cc.

681  {
682  const std::vector<G4InuclElementaryParticle>& out =
683  output->getOutgoingParticles();
684 
685 #ifdef G4CASCADE_DEBUG_INTERFACE
686  // Report on all retry conditions, in order of return logic
687  G4cout << " retryInelasticProton: number of Tries "
688  << ((numberOfTries < maximumTries) ? "RETRY (t)" : "EXIT (f)")
689  << "\n retryInelasticProton: AND collision type ";
690  if (out.empty()) G4cout << "FAILED" << G4endl;
691  else {
692  G4cout << (out.size() == 2 ? "ELASTIC (t)" : "INELASTIC (f)")
693  << "\n retryInelasticProton: AND Leading particles bullet "
694  << (out.size() >= 2 &&
695  (out[0].getDefinition() == bullet->getDefinition() ||
696  out[1].getDefinition() == bullet->getDefinition())
697  ? "YES (t)" : "NO (f)")
698  << G4endl;
699  }
700 #endif
701 
702  return ( (numberOfTries < maximumTries) &&
703  (out.empty() ||
704  (out.size() == 2 &&
705  (out[0].getDefinition() == bullet->getDefinition() ||
706  out[1].getDefinition() == bullet->getDefinition())))
707  );
708 }
const G4ParticleDefinition * getDefinition() const
G4GLOB_DLL std::ostream G4cout
const std::vector< G4InuclElementaryParticle > & getOutgoingParticles() const
#define G4endl
Definition: G4ios.hh:61

Here is the call graph for this function:

Here is the caller graph for this function:

void G4CascadeInterface::SetVerboseLevel ( G4int  verbose)

Definition at line 221 of file G4CascadeInterface.cc.

221  {
223  collider->setVerboseLevel(verbose);
224  balance->setVerboseLevel(verbose);
225  output->setVerboseLevel(verbose);
226 }
void setVerboseLevel(G4int verbose)
virtual void setVerboseLevel(G4int verbose=0)
void setVerboseLevel(G4int verbose=0)
void SetVerboseLevel(G4int value)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4CascadeInterface::throwNonConservationFailure ( )
protected

Definition at line 762 of file G4CascadeInterface.cc.

762  {
763  // NOTE: Once G4HadronicException is changed, use the following line!
764  // G4ExceptionDescription errInfo;
765  std::ostream& errInfo = G4cerr;
766 
767  errInfo << " >>> G4CascadeInterface has non-conserving"
768  << " cascade after " << numberOfTries << " attempts." << G4endl;
769 
770  G4String throwMsg = "G4CascadeInterface - ";
771  if (!balance->energyOkay()) {
772  throwMsg += "Energy";
773  errInfo << " Energy conservation violated by " << balance->deltaE()
774  << " GeV (" << balance->relativeE() << ")" << G4endl;
775  }
776 
777  if (!balance->momentumOkay()) {
778  throwMsg += "Momentum";
779  errInfo << " Momentum conservation violated by " << balance->deltaP()
780  << " GeV/c (" << balance->relativeP() << ")" << G4endl;
781  }
782 
783  if (!balance->baryonOkay()) {
784  throwMsg += "Baryon number";
785  errInfo << " Baryon number violated by " << balance->deltaB() << G4endl;
786  }
787 
788  if (!balance->chargeOkay()) {
789  throwMsg += "Charge";
790  errInfo << " Charge conservation violated by " << balance->deltaQ()
791  << G4endl;
792  }
793 
794  errInfo << " Final event output, for debugging:\n"
795  << " Bullet: \n" << *bullet << G4endl
796  << " Target: \n" << *target << G4endl;
797  output->printCollisionOutput(errInfo);
798 
799  throwMsg += " non-conservation. More info in output.";
800  throw G4HadronicException(__FILE__, __LINE__, throwMsg); // Job ends here!
801 }
const XML_Char * target
Definition: expat.h:268
void printCollisionOutput(std::ostream &os=G4cout) const
#define G4endl
Definition: G4ios.hh:61
G4GLOB_DLL std::ostream G4cerr

Here is the call graph for this function:

Here is the caller graph for this function:

void G4CascadeInterface::useCascadeDeexcitation ( )

Definition at line 210 of file G4CascadeInterface.cc.

210  {
211  collider->useCascadeDeexcitation();
212 }
void useCascadeDeexcitation()

Here is the call graph for this function:

Here is the caller graph for this function:

void G4CascadeInterface::usePreCompoundDeexcitation ( )

Definition at line 214 of file G4CascadeInterface.cc.

214  {
215  collider->usePreCompoundDeexcitation();
216 }
void usePreCompoundDeexcitation()

Here is the call graph for this function:

Here is the caller graph for this function:


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