G4LEPionPlusInelastic Class Reference

#include <G4LEPionPlusInelastic.hh>

Inheritance diagram for G4LEPionPlusInelastic:

G4InelasticInteraction G4HadronicInteraction

Public Member Functions

 G4LEPionPlusInelastic (const G4String &name="G4LEPionPlusInelastic")
 ~G4LEPionPlusInelastic ()
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
virtual void ModelDescription (std::ostream &outFile) const

Detailed Description

Definition at line 43 of file G4LEPionPlusInelastic.hh.


Constructor & Destructor Documentation

G4LEPionPlusInelastic::G4LEPionPlusInelastic ( const G4String name = "G4LEPionPlusInelastic"  ) 

Definition at line 41 of file G4LEPionPlusInelastic.cc.

References G4cout, G4endl, G4HadronicInteraction::SetMaxEnergy(), and G4HadronicInteraction::SetMinEnergy().

00042  :G4InelasticInteraction(name)
00043 {
00044   SetMinEnergy(0.0);
00045   SetMaxEnergy(55.*GeV);
00046   G4cout << "WARNING: model G4LEPionPlusInelastic is being deprecated and will\n"
00047          << "disappear in Geant4 version 10.0"  << G4endl;
00048 }

G4LEPionPlusInelastic::~G4LEPionPlusInelastic (  )  [inline]

Definition at line 49 of file G4LEPionPlusInelastic.hh.

00049 { }


Member Function Documentation

G4HadFinalState * G4LEPionPlusInelastic::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
) [virtual]

Implements G4HadronicInteraction.

Definition at line 66 of file G4LEPionPlusInelastic.cc.

References G4InelasticInteraction::CalculateMomenta(), G4Nucleus::Cinema(), G4InelasticInteraction::DoIsotopeCounting(), G4Nucleus::EvaporationEffects(), G4cout, G4endl, G4HadProjectile::Get4Momentum(), G4HadProjectile::GetDefinition(), G4DynamicParticle::GetDefinition(), G4ReactionProduct::GetKineticEnergy(), G4HadProjectile::GetKineticEnergy(), G4HadProjectile::GetMaterial(), G4ReactionProduct::GetMomentum(), G4Material::GetName(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGMass(), G4FastVector< Type, N >::Initialize(), isAlive, G4InelasticInteraction::isotopeProduction, G4InuclParticleNames::pp, G4Nucleus::ReturnTargetParticle(), G4HadFinalState::SetEnergyChange(), G4ReactionProduct::SetKineticEnergy(), G4ReactionProduct::SetMomentum(), G4HadFinalState::SetMomentumChange(), G4ReactionProduct::SetSide(), G4HadFinalState::SetStatusChange(), G4InelasticInteraction::SetUpChange(), G4HadronicInteraction::theParticleChange, and G4HadronicInteraction::verboseLevel.

00068 {
00069   const G4HadProjectile *originalIncident = &aTrack;
00070   if (originalIncident->GetKineticEnergy()<= 0.1*MeV) {
00071     theParticleChange.SetStatusChange(isAlive);
00072     theParticleChange.SetEnergyChange(aTrack.GetKineticEnergy());
00073     theParticleChange.SetMomentumChange(aTrack.Get4Momentum().vect().unit()); 
00074     return &theParticleChange;      
00075   }
00076 
00077   // create the target particle
00078     
00079   G4DynamicParticle *originalTarget = targetNucleus.ReturnTargetParticle();
00080   G4ReactionProduct targetParticle( originalTarget->GetDefinition() );
00081     
00082   if (verboseLevel > 1) {
00083     const G4Material* targetMaterial = aTrack.GetMaterial();
00084     G4cout << "G4LEPionPlusInelastic::ApplyYourself called" << G4endl;
00085     G4cout << "kinetic energy = " << originalIncident->GetKineticEnergy() << "MeV, ";
00086     G4cout << "target material = " << targetMaterial->GetName() << ", ";
00087     G4cout << "target particle = " << originalTarget->GetDefinition()->GetParticleName()
00088            << G4endl;
00089   }
00090   G4ReactionProduct currentParticle( 
00091   const_cast<G4ParticleDefinition *>(originalIncident->GetDefinition() ) );
00092   currentParticle.SetMomentum( originalIncident->Get4Momentum().vect() );
00093   currentParticle.SetKineticEnergy( originalIncident->GetKineticEnergy() );
00094     
00095   // Fermi motion and evaporation
00096   // As of Geant3, the Fermi energy calculation had not been Done  
00097   G4double ek = originalIncident->GetKineticEnergy();
00098   G4double amas = originalIncident->GetDefinition()->GetPDGMass();
00099     
00100   G4double tkin = targetNucleus.Cinema(ek);
00101   ek += tkin;
00102   currentParticle.SetKineticEnergy(ek);
00103   G4double et = ek + amas;
00104   G4double p = std::sqrt( std::abs((et-amas)*(et+amas)) );
00105   G4double pp = currentParticle.GetMomentum().mag();
00106   if (pp > 0.0) {
00107     G4ThreeVector momentum = currentParticle.GetMomentum();
00108     currentParticle.SetMomentum( momentum * (p/pp) );
00109   }
00110     
00111   // calculate black track energies  
00112   tkin = targetNucleus.EvaporationEffects(ek);
00113   ek -= tkin;
00114   currentParticle.SetKineticEnergy(ek);
00115   et = ek + amas;
00116   p = std::sqrt( std::abs((et-amas)*(et+amas)) );
00117   pp = currentParticle.GetMomentum().mag();
00118   if (pp > 0.0) {
00119     G4ThreeVector momentum = currentParticle.GetMomentum();
00120     currentParticle.SetMomentum( momentum * (p/pp) );
00121   }
00122 
00123   G4ReactionProduct modifiedOriginal = currentParticle;
00124 
00125   currentParticle.SetSide(1); // incident always goes in forward hemisphere
00126   targetParticle.SetSide(-1);  // target always goes in backward hemisphere
00127   G4bool incidentHasChanged = false;
00128   G4bool targetHasChanged = false;
00129   G4bool quasiElastic = false;
00130   G4FastVector<G4ReactionProduct,GHADLISTSIZE> vec;  // vec will contain the secondary particles
00131   G4int vecLen = 0;
00132   vec.Initialize(0);
00133     
00134   const G4double cutOff = 0.1*MeV;
00135   if (currentParticle.GetKineticEnergy() > cutOff)
00136     Cascade(vec, vecLen, originalIncident, currentParticle, targetParticle,
00137             incidentHasChanged, targetHasChanged, quasiElastic);
00138     
00139   CalculateMomenta(vec, vecLen,
00140                    originalIncident, originalTarget, modifiedOriginal,
00141                    targetNucleus, currentParticle, targetParticle,
00142                    incidentHasChanged, targetHasChanged, quasiElastic);
00143     
00144   SetUpChange(vec, vecLen, currentParticle, targetParticle, incidentHasChanged);
00145 
00146   if (isotopeProduction) DoIsotopeCounting(originalIncident, targetNucleus);
00147 
00148   delete originalTarget;
00149   return &theParticleChange;
00150 }

void G4LEPionPlusInelastic::ModelDescription ( std::ostream &  outFile  )  const [virtual]

Reimplemented from G4HadronicInteraction.

Definition at line 51 of file G4LEPionPlusInelastic.cc.

00052 {
00053   outFile << "G4LEPionPlusInelastic is one of the Low Energy Parameterized\n"
00054           << "(LEP) models used to implement inelastic pi+ scattering\n"
00055           << "from nuclei.  It is a re-engineered version of the GHEISHA\n"
00056           << "code of H. Fesefeldt.  It divides the initial collision\n"
00057           << "products into backward- and forward-going clusters which are\n"
00058           << "then decayed into final state hadrons.  The model does not\n"
00059           << "conserve energy on an event-by-event basis.  It may be\n"
00060           << "applied to pions with initial energies between 0 and 25\n"
00061           << "GeV.\n";
00062 }


The documentation for this class was generated from the following files:
Generated on Mon May 27 17:52:23 2013 for Geant4 by  doxygen 1.4.7