Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4INCL::StandardPropagationModel Class Reference

#include <G4INCLStandardPropagationModel.hh>

Inheritance diagram for G4INCL::StandardPropagationModel:
Collaboration diagram for G4INCL::StandardPropagationModel:

Public Member Functions

 StandardPropagationModel (LocalEnergyType localEnergyType, LocalEnergyType localEnergyDeltaType, const G4double hTime=0.0)
 
virtual ~StandardPropagationModel ()
 
G4double getCurrentTime ()
 
void setNucleus (G4INCL::Nucleus *nucleus)
 
G4INCL::NucleusgetNucleus ()
 
G4double shoot (ParticleSpecies const &projectileSpecies, const G4double kineticEnergy, const G4double impactParameter, const G4double phi)
 
G4double shootParticle (ParticleType const t, const G4double kineticEnergy, const G4double impactParameter, const G4double phi)
 
G4double shootComposite (ParticleSpecies const &s, const G4double kineticEnergy, const G4double impactParameter, const G4double phi)
 
void setStoppingTime (G4double)
 
G4double getStoppingTime ()
 
void registerAvatar (G4INCL::IAvatar *anAvatar)
 
IAvatargenerateBinaryCollisionAvatar (Particle *const p1, Particle *const p2)
 Generate a two-particle avatar. More...
 
G4double getReflectionTime (G4INCL::Particle const *const aParticle)
 Get the reflection time. More...
 
G4double getTime (G4INCL::Particle const *const particleA, G4INCL::Particle const *const particleB, G4double *minDistOfApproach) const
 
void generateUpdatedCollisions (const ParticleList &updatedParticles, const ParticleList &particles)
 Generate and register collisions between a list of updated particles and all the other particles. More...
 
void generateCollisions (const ParticleList &particles)
 Generate and register collisions among particles in a list, except between those in another list. More...
 
void generateCollisions (const ParticleList &particles, const ParticleList &except)
 Generate and register collisions among particles in a list, except between those in another list. More...
 
void generateDecays (const ParticleList &particles)
 Generate decays for particles that can decay. More...
 
void updateAvatars (const ParticleList &particles)
 
void generateAllAvatars ()
 (Re)Generate all possible avatars. More...
 
G4INCL::IAvatarpropagate (FinalState const *const fs)
 
- Public Member Functions inherited from G4INCL::IPropagationModel
 IPropagationModel ()
 
virtual ~IPropagationModel ()
 

Additional Inherited Members

Detailed Description

Standard INCL4 particle propagation and avatar prediction

This class implements the standard INCL4 avatar prediction and particle propagation logic. The main idea is to predict all collisions between particles and their reflections from the potential wall. After this we select the avatar with the smallest time, propagate all particles to their positions at that time and return the avatar to the INCL kernel

See Also
G4INCL::Kernel.

The particle trajectories in this propagation model are straight lines and all particles are assumed to move with constant velocity.

Definition at line 69 of file G4INCLStandardPropagationModel.hh.

Constructor & Destructor Documentation

G4INCL::StandardPropagationModel::StandardPropagationModel ( LocalEnergyType  localEnergyType,
LocalEnergyType  localEnergyDeltaType,
const G4double  hTime = 0.0 
)

Definition at line 64 of file G4INCLStandardPropagationModel.cc.

65  :theNucleus(0), maximumTime(70.0), currentTime(0.0),
66  hadronizationTime(hTime),
67  firstAvatar(true),
68  theLocalEnergyType(localEnergyType),
69  theLocalEnergyDeltaType(localEnergyDeltaType)
70  {
71  }
G4INCL::StandardPropagationModel::~StandardPropagationModel ( )
virtual

Definition at line 73 of file G4INCLStandardPropagationModel.cc.

74  {
75  delete theNucleus;
76  }

Member Function Documentation

void G4INCL::StandardPropagationModel::generateAllAvatars ( )

(Re)Generate all possible avatars.

Definition at line 436 of file G4INCLStandardPropagationModel.cc.

436  {
437  ParticleList const &particles = theNucleus->getStore()->getParticles();
438 // assert(!particles.empty());
439  G4double time;
440  for(ParticleIter i=particles.begin(), e=particles.end(); i!=e; ++i) {
441  time = this->getReflectionTime(*i);
442  if(time <= maximumTime) registerAvatar(new SurfaceAvatar(*i, time, theNucleus));
443  }
444  generateCollisions(particles);
445  generateDecays(particles);
446  }
void registerAvatar(G4INCL::IAvatar *anAvatar)
ParticleList const & getParticles() const
Definition: G4INCLStore.hh:253
void generateCollisions(const ParticleList &particles)
Generate and register collisions among particles in a list, except between those in another list...
Store * getStore() const
G4double getReflectionTime(G4INCL::Particle const *const aParticle)
Get the reflection time.
void generateDecays(const ParticleList &particles)
Generate decays for particles that can decay.
double G4double
Definition: G4Types.hh:76
ParticleList::const_iterator ParticleIter

Here is the call graph for this function:

Here is the caller graph for this function:

IAvatar * G4INCL::StandardPropagationModel::generateBinaryCollisionAvatar ( Particle *const  p1,
Particle *const  p2 
)

Generate a two-particle avatar.

Generate a two-particle avatar, if all the appropriate conditions are met.

Definition at line 264 of file G4INCLStandardPropagationModel.cc.

264  {
265  // Is either particle a participant?
266  if(!p1->isParticipant() && !p2->isParticipant() && p1->getParticipantType()==p2->getParticipantType()) return NULL;
267 
268  // Is it a pi-resonance collision (we don't treat them)?
269  if((p1->isResonance() && p2->isPion()) || (p1->isPion() && p2->isResonance()))
270  return NULL;
271 
272  // Will the avatar take place between now and the end of the cascade?
273  G4double minDistOfApproachSquared = 0.0;
274  G4double t = getTime(p1, p2, &minDistOfApproachSquared);
275  if(t>maximumTime || t<currentTime+hadronizationTime) return NULL;
276 
277  // Local energy. Jump through some hoops to calculate the cross section
278  // at the collision point, and clean up after yourself afterwards.
279  G4bool hasLocalEnergy;
280  if(p1->isPion() || p2->isPion())
281  hasLocalEnergy = ((theLocalEnergyDeltaType == FirstCollisionLocalEnergy &&
282  theNucleus->getStore()->getBook().getAcceptedCollisions()==0) ||
283  theLocalEnergyDeltaType == AlwaysLocalEnergy);
284  else
285  hasLocalEnergy = ((theLocalEnergyType == FirstCollisionLocalEnergy &&
286  theNucleus->getStore()->getBook().getAcceptedCollisions()==0) ||
287  theLocalEnergyType == AlwaysLocalEnergy);
288  const G4bool p1HasLocalEnergy = (hasLocalEnergy && !p1->isPion());
289  const G4bool p2HasLocalEnergy = (hasLocalEnergy && !p2->isPion());
290 
291  if(p1HasLocalEnergy) {
292  backupParticle1 = *p1;
293  p1->propagate(t - currentTime);
294  if(p1->getPosition().mag() > theNucleus->getSurfaceRadius(p1)) {
295  *p1 = backupParticle1;
296  return NULL;
297  }
299  }
300  if(p2HasLocalEnergy) {
301  backupParticle2 = *p2;
302  p2->propagate(t - currentTime);
303  if(p2->getPosition().mag() > theNucleus->getSurfaceRadius(p2)) {
304  *p2 = backupParticle2;
305  if(p1HasLocalEnergy) {
306  *p1 = backupParticle1;
307  }
308  return NULL;
309  }
311  }
312 
313  // Compute the total cross section
314  const G4double totalCrossSection = CrossSections::total(p1, p2);
316 
317  // Restore particles to their state before the local-energy tweak
318  if(p1HasLocalEnergy) {
319  *p1 = backupParticle1;
320  }
321  if(p2HasLocalEnergy) {
322  *p2 = backupParticle2;
323  }
324 
325  // Is the CM energy > cutNN? (no cutNN on the first collision)
326  if(theNucleus->getStore()->getBook().getAcceptedCollisions()>0
327  && p1->isNucleon() && p2->isNucleon()
328  && squareTotalEnergyInCM < BinaryCollisionAvatar::getCutNNSquared()) return NULL;
329 
330  // Do the particles come close enough to each other?
331  if(Math::tenPi*minDistOfApproachSquared > totalCrossSection) return NULL;
332 
333  // Bomb out if the two collision partners are the same particle
334 // assert(p1->getID() != p2->getID());
335 
336  // Return a new avatar, then!
337  return new G4INCL::BinaryCollisionAvatar(t, totalCrossSection, theNucleus, p1, p2);
338  }
void transformToLocalEnergyFrame(Nucleus const *const n, Particle *const p)
G4int getAcceptedCollisions() const
Definition: G4INCLBook.hh:100
G4double squareTotalEnergyInCM(Particle const *const p1, Particle const *const p2)
const G4double tenPi
Store * getStore() const
void propagate(G4double step)
G4double getSurfaceRadius(Particle const *const particle) const
Get the maximum allowed radius for a given particle.
G4double getTime(G4INCL::Particle const *const particleA, G4INCL::Particle const *const particleB, G4double *minDistOfApproach) const
Book & getBook()
Definition: G4INCLStore.hh:259
bool G4bool
Definition: G4Types.hh:79
G4double total(Particle const *const p1, Particle const *const p2)
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

void G4INCL::StandardPropagationModel::generateCollisions ( const ParticleList particles)

Generate and register collisions among particles in a list, except between those in another list.

This method generates all possible collisions among the particles. Each collision is generated only once.

Parameters
particleslist of particles

Definition at line 396 of file G4INCLStandardPropagationModel.cc.

396  {
397  // Loop over all the particles
398  for(ParticleIter p1=particles.begin(), e=particles.end(); p1!=e; ++p1) {
399  // Loop over the rest of the particles
400  for(ParticleIter p2 = p1 + 1; p2 != particles.end(); ++p2) {
402  }
403  }
404  }
void registerAvatar(G4INCL::IAvatar *anAvatar)
IAvatar * generateBinaryCollisionAvatar(Particle *const p1, Particle *const p2)
Generate a two-particle avatar.
ParticleList::const_iterator ParticleIter

Here is the call graph for this function:

Here is the caller graph for this function:

void G4INCL::StandardPropagationModel::generateCollisions ( const ParticleList particles,
const ParticleList except 
)

Generate and register collisions among particles in a list, except between those in another list.

This method generates all possible collisions among the particles. Each collision is generated only once. The collision is NOT generated if BOTH collision partners belong to the except list.

You should pass an empty list as the except parameter if you want to generate all possible collisions among particles.

Parameters
particleslist of particles
exceptlist of excluded particles

Definition at line 406 of file G4INCLStandardPropagationModel.cc.

406  {
407 
408  const G4bool haveExcept = !except.empty();
409 
410  // Loop over all the particles
411  for(ParticleIter p1=particles.begin(), e=particles.end(); p1!=e; ++p1)
412  {
413  // Loop over the rest of the particles
414  ParticleIter p2 = p1;
415  for(++p2; p2 != particles.end(); ++p2)
416  {
417  // Skip the collision if both particles must be excluded
418  if(haveExcept && except.contains(*p1) && except.contains(*p2)) continue;
419 
421  }
422  }
423 
424  }
void registerAvatar(G4INCL::IAvatar *anAvatar)
bool G4bool
Definition: G4Types.hh:79
IAvatar * generateBinaryCollisionAvatar(Particle *const p1, Particle *const p2)
Generate a two-particle avatar.
ParticleList::const_iterator ParticleIter

Here is the call graph for this function:

void G4INCL::StandardPropagationModel::generateDecays ( const ParticleList particles)

Generate decays for particles that can decay.

The list of particles given as an argument is allowed to contain also stable particles.

Parameters
particleslist of particles to (possibly) generate decays for

Definition at line 464 of file G4INCLStandardPropagationModel.cc.

464  {
465  for(ParticleIter i=particles.begin(), e=particles.end(); i!=e; ++i) {
466  if((*i)->isDelta()) {
468  G4double time = currentTime + decayTime;
469  if(time <= maximumTime) {
470  registerAvatar(new DecayAvatar((*i), time, theNucleus));
471  }
472  }
473 /* if((*i)->isOmega()) {
474  G4double decayTimeOmega = PionResonanceDecayChannel::computeDecayTime((*i));
475  G4double timeOmega = currentTime + decayTimeOmega;
476  if(timeOmega <= maximumTime) {
477  registerAvatar(new DecayAvatar((*i), timeOmega, theNucleus));
478  }
479  }*/
480  }
481  }
void registerAvatar(G4INCL::IAvatar *anAvatar)
static G4double computeDecayTime(Particle *p)
double G4double
Definition: G4Types.hh:76
ParticleList::const_iterator ParticleIter

Here is the call graph for this function:

Here is the caller graph for this function:

void G4INCL::StandardPropagationModel::generateUpdatedCollisions ( const ParticleList updatedParticles,
const ParticleList particles 
)

Generate and register collisions between a list of updated particles and all the other particles.

This method does not generate collisions among the particles in updatedParticles; in other words, it generates a collision between one of the updatedParticles and one of the particles ONLY IF the latter does not belong to updatedParticles.

If you intend to generate all possible collisions among particles in a list, use generateCollisions().

Parameters
updatedParticleslist of updated particles
particleslist of particles

Definition at line 377 of file G4INCLStandardPropagationModel.cc.

377  {
378 
379  // Loop over all the updated particles
380  for(ParticleIter updated=updatedParticles.begin(), e=updatedParticles.end(); updated!=e; ++updated)
381  {
382  // Loop over all the particles
383  for(ParticleIter particle=particles.begin(), end=particles.end(); particle!=end; ++particle)
384  {
385  /* Consider the generation of a collision avatar only if (*particle)
386  * is not one of the updated particles.
387  * The criterion makes sure that you don't generate avatars between
388  * updated particles. */
389  if(updatedParticles.contains(*particle)) continue;
390 
391  registerAvatar(generateBinaryCollisionAvatar(*particle,*updated));
392  }
393  }
394  }
void registerAvatar(G4INCL::IAvatar *anAvatar)
IAvatar * generateBinaryCollisionAvatar(Particle *const p1, Particle *const p2)
Generate a two-particle avatar.
ParticleList::const_iterator ParticleIter

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4INCL::StandardPropagationModel::getCurrentTime ( )
virtual

Returns the current global time of the system.

Implements G4INCL::IPropagationModel.

Definition at line 250 of file G4INCLStandardPropagationModel.cc.

250  {
251  return currentTime;
252  }
G4INCL::Nucleus * G4INCL::StandardPropagationModel::getNucleus ( )
virtual

Get the nucleus.

Implements G4INCL::IPropagationModel.

Definition at line 78 of file G4INCLStandardPropagationModel.cc.

79  {
80  return theNucleus;
81  }
G4double G4INCL::StandardPropagationModel::getReflectionTime ( G4INCL::Particle const *const  aParticle)

Get the reflection time.

Returns the reflection time of a particle on the potential wall.

Parameters
aParticlepointer to the particle

Definition at line 340 of file G4INCLStandardPropagationModel.cc.

340  {
341  Intersection theIntersection(
343  aParticle->getPosition(),
344  aParticle->getPropagationVelocity(),
345  theNucleus->getSurfaceRadius(aParticle)));
346  G4double time;
347  if(theIntersection.exists) {
348  time = currentTime + theIntersection.time;
349  } else {
350  INCL_ERROR("Imaginary reflection time for particle: " << '\n'
351  << aParticle->print());
352  time = 10000.0;
353  }
354  return time;
355  }
#define INCL_ERROR(x)
G4double getSurfaceRadius(Particle const *const particle) const
Get the maximum allowed radius for a given particle.
Intersection getLaterTrajectoryIntersection(const ThreeVector &x0, const ThreeVector &p, const G4double r)
Compute the second intersection of a straight particle trajectory with a sphere.
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4INCL::StandardPropagationModel::getStoppingTime ( )
virtual

Get the current stopping time.

Implements G4INCL::IPropagationModel.

Definition at line 241 of file G4INCLStandardPropagationModel.cc.

241  {
242  return maximumTime;
243  }
G4double G4INCL::StandardPropagationModel::getTime ( G4INCL::Particle const *const  particleA,
G4INCL::Particle const *const  particleB,
G4double minDistOfApproach 
) const

Get the predicted time of the collision between two particles.

Definition at line 357 of file G4INCLStandardPropagationModel.cc.

359  {
360  G4double time;
361  G4INCL::ThreeVector t13 = particleA->getPropagationVelocity();
362  t13 -= particleB->getPropagationVelocity();
363  G4INCL::ThreeVector distance = particleA->getPosition();
364  distance -= particleB->getPosition();
365  const G4double t7 = t13.dot(distance);
366  const G4double dt = t13.mag2();
367  if(dt <= 1.0e-10) {
368  (*minDistOfApproach) = 100000.0;
369  return currentTime + 100000.0;
370  } else {
371  time = -t7/dt;
372  }
373  (*minDistOfApproach) = distance.mag2() + time * t7;
374  return currentTime + time;
375  }
G4double dot(const ThreeVector &v) const
G4double mag2() const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4INCL::IAvatar * G4INCL::StandardPropagationModel::propagate ( FinalState const *const  fs)
virtual

Propagate all particles and return the first avatar.

Implements G4INCL::IPropagationModel.

Definition at line 483 of file G4INCLStandardPropagationModel.cc.

484  {
485  if(fs) {
486  // We update only the information related to particles that were updated
487  // by the previous avatar.
488 #ifdef INCL_REGENERATE_AVATARS
489 #warning "The INCL_REGENERATE_AVATARS code has not been tested in a while. Use it at your peril."
490  if(!fs->getModifiedParticles().empty() || !fs->getEnteringParticles().empty() || !fs->getCreatedParticles().empty()) {
491  // Regenerates the entire avatar list, skipping collisions between
492  // updated particles
493  theNucleus->getStore()->clearAvatars();
494  theNucleus->getStore()->initialiseParticleAvatarConnections();
495  generateAllAvatarsExceptUpdated(fs);
496  }
497 #else
498  ParticleList const &updatedParticles = fs->getModifiedParticles();
499  if(fs->getValidity()==PauliBlockedFS) {
500  // This final state might represents the outcome of a Pauli-blocked delta
501  // decay
502 // assert(updatedParticles.empty() || (updatedParticles.size()==1 && updatedParticles.front()->isResonance()));
503 // assert(fs->getEnteringParticles().empty() && fs->getCreatedParticles().empty() && fs->getOutgoingParticles().empty() && fs->getDestroyedParticles().empty());
504  generateDecays(updatedParticles);
505  } else {
506  ParticleList const &entering = fs->getEnteringParticles();
507  generateDecays(updatedParticles);
508  generateDecays(entering);
509 
510  ParticleList const &created = fs->getCreatedParticles();
511  if(created.empty() && entering.empty())
512  updateAvatars(updatedParticles);
513  else {
514  ParticleList updatedParticlesCopy = updatedParticles;
515  updatedParticlesCopy.insert(updatedParticlesCopy.end(), entering.begin(), entering.end());
516  updatedParticlesCopy.insert(updatedParticlesCopy.end(), created.begin(), created.end());
517  updateAvatars(updatedParticlesCopy);
518  }
519  }
520 #endif
521  }
522 
523  G4INCL::IAvatar *theAvatar = theNucleus->getStore()->findSmallestTime();
524  if(theAvatar == 0) return 0; // Avatar list is empty
525  // theAvatar->dispose();
526 
527  if(theAvatar->getTime() < currentTime) {
528  INCL_ERROR("Avatar time = " << theAvatar->getTime() << ", currentTime = " << currentTime << '\n');
529  return 0;
530  } else if(theAvatar->getTime() > currentTime) {
531  theNucleus->getStore()->timeStep(theAvatar->getTime() - currentTime);
532 
533  currentTime = theAvatar->getTime();
534  theNucleus->getStore()->getBook().setCurrentTime(currentTime);
535  }
536 
537  return theAvatar;
538  }
void clearAvatars()
Definition: G4INCLStore.cc:193
#define INCL_ERROR(x)
IAvatar * findSmallestTime()
Definition: G4INCLStore.cc:142
Store * getStore() const
void updateAvatars(const ParticleList &particles)
G4double getTime() const
Book & getBook()
Definition: G4INCLStore.hh:259
void generateDecays(const ParticleList &particles)
Generate decays for particles that can decay.
void setCurrentTime(G4double t)
Definition: G4INCLBook.hh:97
void timeStep(G4double step)
Definition: G4INCLStore.cc:168

Here is the call graph for this function:

void G4INCL::StandardPropagationModel::registerAvatar ( G4INCL::IAvatar anAvatar)

Add an avatar to the storage.

Definition at line 259 of file G4INCLStandardPropagationModel.cc.

260  {
261  if(anAvatar) theNucleus->getStore()->add(anAvatar);
262  }
Store * getStore() const
void add(Particle *p)
Definition: G4INCLStore.cc:58

Here is the call graph for this function:

Here is the caller graph for this function:

void G4INCL::StandardPropagationModel::setNucleus ( G4INCL::Nucleus nucleus)
virtual

Set the nucleus for this propagation model.

Implements G4INCL::IPropagationModel.

Definition at line 254 of file G4INCLStandardPropagationModel.cc.

255  {
256  theNucleus = nucleus;
257  }
void G4INCL::StandardPropagationModel::setStoppingTime ( G4double  time)
virtual

Set the stopping time of the simulation.

Implements G4INCL::IPropagationModel.

Definition at line 245 of file G4INCLStandardPropagationModel.cc.

245  {
246 // assert(time>0.0);
247  maximumTime = time;
248  }
G4double G4INCL::StandardPropagationModel::shoot ( ParticleSpecies const &  projectileSpecies,
const G4double  kineticEnergy,
const G4double  impactParameter,
const G4double  phi 
)
virtual

Implements G4INCL::IPropagationModel.

Definition at line 83 of file G4INCLStandardPropagationModel.cc.

83  {
84  if(projectileSpecies.theType==Composite)
85  return shootComposite(projectileSpecies, kineticEnergy, impactParameter, phi);
86  else
87  return shootParticle(projectileSpecies.theType, kineticEnergy, impactParameter, phi);
88  }
G4double shootComposite(ParticleSpecies const &s, const G4double kineticEnergy, const G4double impactParameter, const G4double phi)
G4double shootParticle(ParticleType const t, const G4double kineticEnergy, const G4double impactParameter, const G4double phi)

Here is the call graph for this function:

G4double G4INCL::StandardPropagationModel::shootComposite ( ParticleSpecies const &  s,
const G4double  kineticEnergy,
const G4double  impactParameter,
const G4double  phi 
)
virtual

Implements G4INCL::IPropagationModel.

Definition at line 167 of file G4INCLStandardPropagationModel.cc.

167  {
168  theNucleus->setNucleusNucleusCollision();
169  currentTime = 0.0;
170 
171  // Create the ProjectileRemnant object
172  ProjectileRemnant *pr = new ProjectileRemnant(species, kineticEnergy);
173 
174  // Same stopping time as for nucleon-nucleus
175  maximumTime = 29.8 * std::pow(theNucleus->getA(), 0.16);
176 
177  // If the incoming cluster is slow, use a larger stopping time
178  const G4double rms = ParticleTable::getLargestNuclearRadius(pr->getA(), pr->getZ());
179  const G4double rMax = theNucleus->getUniverseRadius();
180  const G4double distance = 2.*rMax + 2.725*rms;
181  const G4double projectileVelocity = pr->boostVector().mag();
182  const G4double traversalTime = distance / projectileVelocity;
183  if(maximumTime < traversalTime)
184  maximumTime = traversalTime;
185  INCL_DEBUG("Cascade stopping time is " << maximumTime << '\n');
186 
187  // If Coulomb is activated, do not process events with impact
188  // parameter larger than the maximum impact parameter, taking into
189  // account Coulomb distortion.
190  if(impactParameter>CoulombDistortion::maxImpactParameter(pr,theNucleus)) {
191  INCL_DEBUG("impactParameter>CoulombDistortion::maxImpactParameter" << '\n');
192  delete pr;
193  return -1.;
194  }
195 
196  // Position the cluster at the right impact parameter
197  ThreeVector position(impactParameter * std::cos(phi),
198  impactParameter * std::sin(phi),
199  0.);
200  pr->setPosition(position);
201 
202  // Fill in the relevant kinematic variables
203  theNucleus->setIncomingAngularMomentum(pr->getAngularMomentum());
204  theNucleus->setIncomingMomentum(pr->getMomentum());
205  theNucleus->setInitialEnergy(pr->getEnergy()
206  + ParticleTable::getTableMass(theNucleus->getA(),theNucleus->getZ()));
207 
209  firstAvatar = false;
210 
211  // Get the entry avatars from Coulomb
212  IAvatarList theAvatarList
213  = CoulombDistortion::bringToSurface(pr, theNucleus);
214 
215  if(theAvatarList.empty()) {
216  INCL_DEBUG("No ParticleEntryAvatar found, transparent event" << '\n');
217  delete pr;
218  return -1.;
219  }
220 
221  /* Store the internal kinematics of the projectile remnant.
222  *
223  * Note that this is at variance with the Fortran version, which stores
224  * the initial kinematics of the particles *after* putting the spectators
225  * on mass shell, but *before* removing the same energy from the
226  * participants. Due to the different code flow, doing so in the C++
227  * version leads to wrong excitation energies for the forced compound
228  * nucleus.
229  */
230  pr->storeComponents();
231 
232  // Tell the Nucleus about the ProjectileRemnant
233  theNucleus->setProjectileRemnant(pr);
234 
235  // Register the ParticleEntryAvatars
236  theNucleus->getStore()->addParticleEntryAvatars(theAvatarList);
237 
238  return pr->getTransversePosition().mag();
239  }
G4int getA() const
Returns the baryon number.
ParticleEntryAvatar * bringToSurface(Particle *p, Nucleus *const n)
Modify the momentum of an incoming particle and position it on the surface of the nucleus...
void generateAllAvatars()
(Re)Generate all possible avatars.
void setIncomingAngularMomentum(const ThreeVector &j)
Set the incoming angular-momentum vector.
void setIncomingMomentum(const ThreeVector &p)
Set the incoming momentum vector.
void setInitialEnergy(const G4double e)
Set the initial energy.
Store * getStore() const
G4double maxImpactParameter(ParticleSpecies const &p, const G4double kinE, Nucleus const *const n)
Return the maximum impact parameter for Coulomb-distorted trajectories.
void addParticleEntryAvatars(IAvatarList const &al)
Add one ParticleEntry avatar.
Definition: G4INCLStore.cc:78
#define position
Definition: xmlparse.cc:622
G4double getLargestNuclearRadius(const G4int A, const G4int Z)
G4int getZ() const
Returns the charge number.
UnorderedVector< IAvatar * > IAvatarList
G4ThreadLocal NuclearMassFn getTableMass
Static pointer to the mass function for nuclei.
void setProjectileRemnant(ProjectileRemnant *const c)
Set the projectile remnant.
G4double getUniverseRadius() const
Getter for theUniverseRadius.
double G4double
Definition: G4Types.hh:76
#define INCL_DEBUG(x)
void setNucleusNucleusCollision()
Set a nucleus-nucleus collision.

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4INCL::StandardPropagationModel::shootParticle ( ParticleType const  t,
const G4double  kineticEnergy,
const G4double  impactParameter,
const G4double  phi 
)
virtual

Implements G4INCL::IPropagationModel.

Definition at line 90 of file G4INCLStandardPropagationModel.cc.

90  {
91  theNucleus->setParticleNucleusCollision();
92  currentTime = 0.0;
93 
94  // Create the projectile particle
95  const G4double projectileMass = ParticleTable::getTableParticleMass(type);
96  G4double energy = kineticEnergy + projectileMass;
97  G4double momentumZ = std::sqrt(energy*energy - projectileMass*projectileMass);
98  ThreeVector momentum(0.0, 0.0, momentumZ);
99  Particle *p= new G4INCL::Particle(type, energy, momentum, ThreeVector());
100 
101  G4double temfin;
102  G4double TLab;
103  if( p->isPion() ) {
104  temfin = 30.18 * std::pow(theNucleus->getA(), 0.17);
105  TLab = p->getKineticEnergy();
106  } else {
107  temfin = 29.8 * std::pow(theNucleus->getA(), 0.16);
108  TLab = p->getKineticEnergy()/p->getA();
109  }
110 
111  // energy-dependent stopping time above 2 AGeV
112  if(TLab>2000.)
113  temfin *= (5.8E4-TLab)/5.6E4;
114 
115  maximumTime = temfin;
116 
117  // If the incoming particle is slow, use a larger stopping time
118  const G4double rMax = theNucleus->getUniverseRadius();
119  const G4double distance = 2.*rMax;
120  const G4double projectileVelocity = p->boostVector().mag();
121  const G4double traversalTime = distance / projectileVelocity;
122  if(maximumTime < traversalTime)
123  maximumTime = traversalTime;
124  INCL_DEBUG("Cascade stopping time is " << maximumTime << '\n');
125 
126  // If Coulomb is activated, do not process events with impact
127  // parameter larger than the maximum impact parameter, taking into
128  // account Coulomb distortion.
129  if(impactParameter>CoulombDistortion::maxImpactParameter(p->getSpecies(), kineticEnergy, theNucleus)) {
130  INCL_DEBUG("impactParameter>CoulombDistortion::maxImpactParameter" << '\n');
131  delete p;
132  return -1.;
133  }
134 
135  ThreeVector position(impactParameter * std::cos(phi),
136  impactParameter * std::sin(phi),
137  0.);
138  p->setPosition(position);
139 
140  // Fill in the relevant kinematic variables
141  theNucleus->setIncomingAngularMomentum(p->getAngularMomentum());
142  theNucleus->setIncomingMomentum(p->getMomentum());
143  theNucleus->setInitialEnergy(p->getEnergy()
144  + ParticleTable::getTableMass(theNucleus->getA(),theNucleus->getZ()));
145 
146  // Reset the particle kinematics to the INCL values
147  p->setINCLMass();
148  p->setEnergy(p->getMass() + kineticEnergy);
149  p->adjustMomentumFromEnergy();
150 
151  p->makeProjectileSpectator();
153  firstAvatar = false;
154 
155  // Get the entry avatars from Coulomb and put them in the Store
156  ParticleEntryAvatar *theEntryAvatar = CoulombDistortion::bringToSurface(p, theNucleus);
157  if(theEntryAvatar) {
158  theNucleus->getStore()->addParticleEntryAvatar(theEntryAvatar);
159 
160  return p->getTransversePosition().mag();
161  } else {
162  delete p;
163  return -1.;
164  }
165  }
G4int getA() const
Returns the baryon number.
ParticleEntryAvatar * bringToSurface(Particle *p, Nucleus *const n)
Modify the momentum of an incoming particle and position it on the surface of the nucleus...
void generateAllAvatars()
(Re)Generate all possible avatars.
void setIncomingAngularMomentum(const ThreeVector &j)
Set the incoming angular-momentum vector.
void setIncomingMomentum(const ThreeVector &p)
Set the incoming momentum vector.
G4ThreadLocal ParticleMassFn getTableParticleMass
Static pointer to the mass function for particles.
void setInitialEnergy(const G4double e)
Set the initial energy.
const char * p
Definition: xmltok.h:285
Store * getStore() const
void setParticleNucleusCollision()
Set a particle-nucleus collision.
G4double maxImpactParameter(ParticleSpecies const &p, const G4double kinE, Nucleus const *const n)
Return the maximum impact parameter for Coulomb-distorted trajectories.
#define position
Definition: xmlparse.cc:622
G4int getZ() const
Returns the charge number.
void addParticleEntryAvatar(IAvatar *a)
Add one ParticleEntry avatar.
Definition: G4INCLStore.cc:66
G4double energy(const ThreeVector &p, const G4double m)
G4ThreadLocal NuclearMassFn getTableMass
Static pointer to the mass function for nuclei.
G4double getUniverseRadius() const
Getter for theUniverseRadius.
double G4double
Definition: G4Types.hh:76
#define INCL_DEBUG(x)

Here is the call graph for this function:

Here is the caller graph for this function:

void G4INCL::StandardPropagationModel::updateAvatars ( const ParticleList particles)

Update all avatars related to a particle.

Definition at line 426 of file G4INCLStandardPropagationModel.cc.

426  {
427 
428  for(ParticleIter iter=particles.begin(), e=particles.end(); iter!=e; ++iter) {
429  G4double time = this->getReflectionTime(*iter);
430  if(time <= maximumTime) registerAvatar(new SurfaceAvatar(*iter, time, theNucleus));
431  }
432  ParticleList const &p = theNucleus->getStore()->getParticles();
433  generateUpdatedCollisions(particles, p); // Predict collisions with spectators and participants
434  }
void registerAvatar(G4INCL::IAvatar *anAvatar)
ParticleList const & getParticles() const
Definition: G4INCLStore.hh:253
const char * p
Definition: xmltok.h:285
Store * getStore() const
G4double getReflectionTime(G4INCL::Particle const *const aParticle)
Get the reflection time.
double G4double
Definition: G4Types.hh:76
void generateUpdatedCollisions(const ParticleList &updatedParticles, const ParticleList &particles)
Generate and register collisions between a list of updated particles and all the other particles...
ParticleList::const_iterator ParticleIter

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: