#include <G4LEDeuteronInelastic.hh>
Inheritance diagram for G4LEDeuteronInelastic:
Public Member Functions | |
G4LEDeuteronInelastic () | |
~G4LEDeuteronInelastic () | |
G4HadFinalState * | ApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus) |
virtual void | ModelDescription (std::ostream &outFile) const |
Definition at line 45 of file G4LEDeuteronInelastic.hh.
G4LEDeuteronInelastic::G4LEDeuteronInelastic | ( | ) | [inline] |
Definition at line 49 of file G4LEDeuteronInelastic.hh.
References G4cout, G4endl, G4HadronicInteraction::SetMaxEnergy(), and G4HadronicInteraction::SetMinEnergy().
00049 : G4InelasticInteraction("G4LEDeuteronInelastic") 00050 { 00051 SetMinEnergy( 0.0 ); 00052 // SetMaxEnergy( 100.*CLHEP::MeV ); // NUCREC only worked for energies < 100MeV 00053 // Work around to avoid exception in G4EnergyRangeManager 00054 SetMaxEnergy( 10.*CLHEP::TeV ); // NUCREC only worked for energies < 100MeV 00055 G4cout << "WARNING: model G4LEDeuteronInelastic is being deprecated and will\n" 00056 << "disappear in Geant4 version 10.0" << G4endl; 00057 }
G4LEDeuteronInelastic::~G4LEDeuteronInelastic | ( | ) | [inline] |
G4HadFinalState * G4LEDeuteronInelastic::ApplyYourself | ( | const G4HadProjectile & | aTrack, | |
G4Nucleus & | targetNucleus | |||
) | [virtual] |
Implements G4HadronicInteraction.
Definition at line 51 of file G4LEDeuteronInelastic.cc.
References G4HadFinalState::AddSecondary(), G4Nucleus::AtomicMass(), G4HadFinalState::Clear(), G4InelasticInteraction::DoIsotopeCounting(), G4cout, G4endl, G4HadProjectile::Get4Momentum(), G4Nucleus::GetA_asInt(), G4HadProjectile::GetKineticEnergy(), G4HadProjectile::GetMaterial(), G4Material::GetName(), G4Nucleus::GetZ_asInt(), G4FastVector< Type, N >::Initialize(), isAlive, G4InelasticInteraction::isotopeProduction, G4ReactionDynamics::NuclearReaction(), G4DynamicParticle::SetDefinition(), G4HadFinalState::SetEnergyChange(), G4DynamicParticle::SetMomentum(), G4HadFinalState::SetMomentumChange(), G4HadFinalState::SetStatusChange(), G4HadronicInteraction::theParticleChange, G4InelasticInteraction::theReactionDynamics, and G4HadronicInteraction::verboseLevel.
00053 { 00054 theParticleChange.Clear(); 00055 const G4HadProjectile* originalIncident = &aTrack; 00056 00057 if (verboseLevel > 1) { 00058 const G4Material *targetMaterial = aTrack.GetMaterial(); 00059 G4cout << "G4LEDeuteronInelastic::ApplyYourself called" << G4endl; 00060 G4cout << "kinetic energy = " << originalIncident->GetKineticEnergy()/MeV << "MeV, "; 00061 G4cout << "target material = " << targetMaterial->GetName() << ", "; 00062 } 00063 00064 // Work-around for lack of model above 100 MeV 00065 if (originalIncident->GetKineticEnergy()/MeV > 100. || 00066 originalIncident->GetKineticEnergy() <= 0.1*MeV) { 00067 theParticleChange.SetStatusChange(isAlive); 00068 theParticleChange.SetEnergyChange(aTrack.GetKineticEnergy()); 00069 theParticleChange.SetMomentumChange(aTrack.Get4Momentum().vect().unit()); 00070 return &theParticleChange; 00071 } 00072 00073 G4double A = targetNucleus.GetA_asInt(); 00074 G4double Z = targetNucleus.GetZ_asInt(); 00075 G4double theAtomicMass = targetNucleus.AtomicMass(A, Z); 00076 G4double massVec[9]; 00077 massVec[0] = targetNucleus.AtomicMass( A+2.0, Z+1.0 ); 00078 massVec[1] = targetNucleus.AtomicMass( A+1.0, Z+1.0 ); 00079 massVec[2] = targetNucleus.AtomicMass( A+1.0, Z ); 00080 massVec[3] = theAtomicMass; 00081 massVec[4] = 0.; 00082 if (A > 1.0 && A-1.0 > Z) 00083 massVec[4] = targetNucleus.AtomicMass(A-1.0, Z); 00084 massVec[5] = 0.; 00085 if (A > 2.0 && Z > 1.0 && A-2.0 > Z-1.0) 00086 massVec[5] = targetNucleus.AtomicMass(A-2.0, Z-1.0); 00087 massVec[6] = 0.; 00088 if (A > Z+1.0) 00089 massVec[6] = targetNucleus.AtomicMass(A, Z+1.0); 00090 massVec[7] = massVec[3]; 00091 massVec[8] = 0.; 00092 if (Z > 1.0) massVec[8] = targetNucleus.AtomicMass(A,Z-1.0); 00093 00094 G4FastVector<G4ReactionProduct,4> vec; // vec will contain the secondary particles 00095 G4int vecLen = 0; 00096 vec.Initialize( 0 ); 00097 00098 theReactionDynamics.NuclearReaction(vec, vecLen, originalIncident, 00099 targetNucleus, theAtomicMass, massVec); 00100 00101 G4double p = vec[0]->GetMomentum().mag(); 00102 theParticleChange.SetMomentumChange( vec[0]->GetMomentum() * (1.0/p) ); 00103 theParticleChange.SetEnergyChange( vec[0]->GetKineticEnergy() ); 00104 delete vec[0]; 00105 00106 if (vecLen <= 1) 00107 { 00108 theParticleChange.SetStatusChange(isAlive); 00109 theParticleChange.SetEnergyChange(aTrack.GetKineticEnergy()); 00110 theParticleChange.SetMomentumChange(aTrack.Get4Momentum().vect().unit()); 00111 if (isotopeProduction) DoIsotopeCounting(originalIncident, targetNucleus); 00112 return &theParticleChange; 00113 } 00114 00115 G4DynamicParticle* pd; 00116 for (G4int i=1; i<vecLen; ++i) { 00117 pd = new G4DynamicParticle(); 00118 pd->SetDefinition( vec[i]->GetDefinition() ); 00119 pd->SetMomentum( vec[i]->GetMomentum() ); 00120 theParticleChange.AddSecondary( pd ); 00121 delete vec[i]; 00122 } 00123 00124 if (isotopeProduction) DoIsotopeCounting(originalIncident, targetNucleus); 00125 return &theParticleChange; 00126 }
void G4LEDeuteronInelastic::ModelDescription | ( | std::ostream & | outFile | ) | const [virtual] |
Reimplemented from G4HadronicInteraction.
Definition at line 36 of file G4LEDeuteronInelastic.cc.
00037 { 00038 outFile << "G4LEDeuteronInelastic is one of the Low Energy Parameterized\n" 00039 << "(LEP) models used to implement inelastic deuteron scattering\n" 00040 << "from nuclei. It is a re-engineered version of the GHEISHA\n" 00041 << "code of H. Fesefeldt. It divides the initial collision\n" 00042 << "products into backward- and forward-going clusters which are\n" 00043 << "then decayed into final state hadrons. The model does not\n" 00044 << "conserve energy on an event-by-event basis. It may be\n" 00045 << "applied to deuterons with initial energies between 0 and 10\n" 00046 << "TeV.\n"; 00047 }