#include <G4LEAntiOmegaMinusInelastic.hh>
Inheritance diagram for G4LEAntiOmegaMinusInelastic:
Public Member Functions | |
G4LEAntiOmegaMinusInelastic () | |
~G4LEAntiOmegaMinusInelastic () | |
G4HadFinalState * | ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus) |
virtual void | ModelDescription (std::ostream &outFile) const |
Definition at line 45 of file G4LEAntiOmegaMinusInelastic.hh.
G4LEAntiOmegaMinusInelastic::G4LEAntiOmegaMinusInelastic | ( | ) | [inline] |
Definition at line 49 of file G4LEAntiOmegaMinusInelastic.hh.
References G4cout, G4endl, G4HadronicInteraction::SetMaxEnergy(), and G4HadronicInteraction::SetMinEnergy().
00049 : G4InelasticInteraction("G4LEAntiOmegaMinusInelastic") 00050 { 00051 SetMinEnergy(0.0); 00052 SetMaxEnergy(25.*CLHEP::GeV); 00053 G4cout << "WARNING: model G4LEAntiOmegaMinusInelastic is being deprecated and will\n" 00054 << "disappear in Geant4 version 10.0" << G4endl; 00055 }
G4LEAntiOmegaMinusInelastic::~G4LEAntiOmegaMinusInelastic | ( | ) | [inline] |
G4HadFinalState * G4LEAntiOmegaMinusInelastic::ApplyYourself | ( | const G4HadProjectile & | aTrack, | |
G4Nucleus & | targetNucleus | |||
) | [virtual] |
Implements G4HadronicInteraction.
Definition at line 55 of file G4LEAntiOmegaMinusInelastic.cc.
References G4InelasticInteraction::CalculateMomenta(), G4Nucleus::Cinema(), G4InelasticInteraction::DoIsotopeCounting(), G4Nucleus::EvaporationEffects(), G4cout, G4endl, G4UniformRand, G4HadProjectile::Get4Momentum(), G4HadProjectile::GetDefinition(), G4DynamicParticle::GetDefinition(), G4ReactionProduct::GetKineticEnergy(), G4HadProjectile::GetKineticEnergy(), G4HadProjectile::GetMaterial(), G4ReactionProduct::GetMomentum(), G4Material::GetName(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetPDGMass(), G4ReactionProduct::GetTotalMomentum(), 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.
00057 { 00058 const G4HadProjectile* originalIncident = &aTrack; 00059 if (originalIncident->GetKineticEnergy()<= 0.1*MeV) { 00060 theParticleChange.SetStatusChange(isAlive); 00061 theParticleChange.SetEnergyChange(aTrack.GetKineticEnergy()); 00062 theParticleChange.SetMomentumChange(aTrack.Get4Momentum().vect().unit()); 00063 return &theParticleChange; 00064 } 00065 00066 // create the target particle 00067 G4DynamicParticle* originalTarget = targetNucleus.ReturnTargetParticle(); 00068 00069 if (verboseLevel > 1) { 00070 const G4Material *targetMaterial = aTrack.GetMaterial(); 00071 G4cout << "kinetic energy = " << originalIncident->GetKineticEnergy()/MeV << "MeV, "; 00072 G4cout << "target material = " << targetMaterial->GetName() << ", "; 00073 G4cout << "target particle = " << originalTarget->GetDefinition()->GetParticleName() 00074 << G4endl; 00075 } 00076 00077 // Fermi motion and evaporation 00078 // As of Geant3, the Fermi energy calculation had not been Done 00079 G4double ek = originalIncident->GetKineticEnergy()/MeV; 00080 G4double amas = originalIncident->GetDefinition()->GetPDGMass()/MeV; 00081 G4ReactionProduct modifiedOriginal; 00082 modifiedOriginal = *originalIncident; 00083 00084 G4double tkin = targetNucleus.Cinema( ek ); 00085 ek += tkin; 00086 modifiedOriginal.SetKineticEnergy( ek*MeV ); 00087 G4double et = ek + amas; 00088 G4double p = std::sqrt( std::abs((et-amas)*(et+amas)) ); 00089 G4double pp = modifiedOriginal.GetMomentum().mag()/MeV; 00090 if (pp > 0.0) { 00091 G4ThreeVector momentum = modifiedOriginal.GetMomentum(); 00092 modifiedOriginal.SetMomentum( momentum * (p/pp) ); 00093 } 00094 00095 // calculate black track energies 00096 tkin = targetNucleus.EvaporationEffects( ek ); 00097 ek -= tkin; 00098 modifiedOriginal.SetKineticEnergy( ek*MeV ); 00099 et = ek + amas; 00100 p = std::sqrt( std::abs((et-amas)*(et+amas)) ); 00101 pp = modifiedOriginal.GetMomentum().mag()/MeV; 00102 if (pp > 0.0) { 00103 G4ThreeVector momentum = modifiedOriginal.GetMomentum(); 00104 modifiedOriginal.SetMomentum( momentum * (p/pp) ); 00105 } 00106 G4ReactionProduct currentParticle = modifiedOriginal; 00107 G4ReactionProduct targetParticle; 00108 targetParticle = *originalTarget; 00109 currentParticle.SetSide( 1 ); // incident always goes in forward hemisphere 00110 targetParticle.SetSide( -1 ); // target always goes in backward hemisphere 00111 G4bool incidentHasChanged = false; 00112 G4bool targetHasChanged = false; 00113 G4bool quasiElastic = false; 00114 G4FastVector<G4ReactionProduct,GHADLISTSIZE> vec; // vec will contain the secondary particles 00115 G4int vecLen = 0; 00116 vec.Initialize( 0 ); 00117 00118 const G4double cutOff = 0.1; 00119 const G4double anni = std::min( 1.3*currentParticle.GetTotalMomentum()/GeV, 0.4 ); 00120 00121 if ((currentParticle.GetKineticEnergy()/MeV > cutOff) || (G4UniformRand() > anni) ) 00122 Cascade(vec, vecLen, originalIncident, currentParticle, targetParticle, 00123 incidentHasChanged, targetHasChanged, quasiElastic); 00124 00125 CalculateMomenta(vec, vecLen, originalIncident, originalTarget, 00126 modifiedOriginal, targetNucleus, currentParticle, 00127 targetParticle, incidentHasChanged, targetHasChanged, 00128 quasiElastic); 00129 00130 SetUpChange(vec, vecLen, currentParticle, targetParticle, incidentHasChanged); 00131 00132 if (isotopeProduction) DoIsotopeCounting(originalIncident, targetNucleus); 00133 00134 delete originalTarget; 00135 return &theParticleChange; 00136 }
void G4LEAntiOmegaMinusInelastic::ModelDescription | ( | std::ostream & | outFile | ) | const [virtual] |
Reimplemented from G4HadronicInteraction.
Definition at line 41 of file G4LEAntiOmegaMinusInelastic.cc.
00042 { 00043 outFile << "G4LEAntiOmegaMinusInelastic is one of the Low Energy\n" 00044 << "Parameterized (LEP) models used to implement inelastic\n" 00045 << "antiOmega- scattering from nuclei. It is a re-engineered\n" 00046 << "version of the GHEISHA code of H. Fesefeldt. It divides the\n" 00047 << "initial collision products into backward- and forward-going\n" 00048 << "clusters which are then decayed into final state hadrons. The\n" 00049 << "model does not conserve energy on an event-by-event basis. It\n" 00050 << "may be applied to antiOmega- with initial energies between 0\n" 00051 << "and 25 GeV.\n"; 00052 }