G4RPGKZeroInelastic Class Reference

#include <G4RPGKZeroInelastic.hh>

Inheritance diagram for G4RPGKZeroInelastic:

G4RPGInelastic G4HadronicInteraction

Public Member Functions

 G4RPGKZeroInelastic ()
 ~G4RPGKZeroInelastic ()
G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)

Detailed Description

Definition at line 43 of file G4RPGKZeroInelastic.hh.


Constructor & Destructor Documentation

G4RPGKZeroInelastic::G4RPGKZeroInelastic (  )  [inline]

Definition at line 47 of file G4RPGKZeroInelastic.hh.

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

00047                           : G4RPGInelastic("G4RPGKZeroInelastic")
00048     {
00049       SetMinEnergy( 0.0 );
00050       SetMaxEnergy( 25.*CLHEP::GeV );
00051     }

G4RPGKZeroInelastic::~G4RPGKZeroInelastic (  )  [inline]

Definition at line 53 of file G4RPGKZeroInelastic.hh.

00054     { }


Member Function Documentation

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

Implements G4HadronicInteraction.

Definition at line 35 of file G4RPGKZeroInelastic.cc.

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

Referenced by G4RPGKShortInelastic::ApplyYourself(), and G4RPGKLongInelastic::ApplyYourself().

00037 {
00038   const G4HadProjectile *originalIncident = &aTrack;
00039 
00040   // create the target particle
00041 
00042   G4DynamicParticle *originalTarget = targetNucleus.ReturnTargetParticle();
00043     
00044   if( verboseLevel > 1 )
00045   {
00046     const G4Material *targetMaterial = aTrack.GetMaterial();
00047     G4cout << "G4RPGKZeroInelastic::ApplyYourself called" << G4endl;
00048     G4cout << "kinetic energy = " << originalIncident->GetKineticEnergy()/MeV << "MeV, ";
00049     G4cout << "target material = " << targetMaterial->GetName() << ", ";
00050     G4cout << "target particle = " << originalTarget->GetDefinition()->GetParticleName()
00051            << G4endl;
00052   }
00053 
00054   // Fermi motion and evaporation
00055   // As of Geant3, the Fermi energy calculation had not been Done
00056 
00057     G4double ek = originalIncident->GetKineticEnergy()/MeV;
00058     G4double amas = originalIncident->GetDefinition()->GetPDGMass()/MeV;
00059     G4ReactionProduct modifiedOriginal;
00060     modifiedOriginal = *originalIncident;
00061     
00062     G4double tkin = targetNucleus.Cinema( ek );
00063     ek += tkin;
00064     modifiedOriginal.SetKineticEnergy( ek*MeV );
00065     G4double et = ek + amas;
00066     G4double p = std::sqrt( std::abs((et-amas)*(et+amas)) );
00067     G4double pp = modifiedOriginal.GetMomentum().mag()/MeV;
00068     if( pp > 0.0 )
00069     {
00070       G4ThreeVector momentum = modifiedOriginal.GetMomentum();
00071       modifiedOriginal.SetMomentum( momentum * (p/pp) );
00072     }
00073     //
00074     // calculate black track energies
00075     //
00076     tkin = targetNucleus.EvaporationEffects( ek );
00077     ek -= tkin;
00078     modifiedOriginal.SetKineticEnergy( ek*MeV );
00079     et = ek + amas;
00080     p = std::sqrt( std::abs((et-amas)*(et+amas)) );
00081     pp = modifiedOriginal.GetMomentum().mag()/MeV;
00082     if( pp > 0.0 )
00083     {
00084       G4ThreeVector momentum = modifiedOriginal.GetMomentum();
00085       modifiedOriginal.SetMomentum( momentum * (p/pp) );
00086     }
00087     G4ReactionProduct currentParticle = modifiedOriginal;
00088     G4ReactionProduct targetParticle;
00089     targetParticle = *originalTarget;
00090     currentParticle.SetSide( 1 ); // incident always goes in forward hemisphere
00091     targetParticle.SetSide( -1 );  // target always goes in backward hemisphere
00092     G4bool incidentHasChanged = false;
00093     G4bool targetHasChanged = false;
00094     G4bool quasiElastic = false;
00095     G4FastVector<G4ReactionProduct,GHADLISTSIZE> vec;  // vec will contain the secondary particles
00096     G4int vecLen = 0;
00097     vec.Initialize( 0 );
00098     
00099     const G4double cutOff = 0.1;
00100     if( currentParticle.GetKineticEnergy()/MeV > cutOff )
00101       Cascade( vec, vecLen,
00102                originalIncident, currentParticle, targetParticle,
00103                incidentHasChanged, targetHasChanged, quasiElastic );
00104     
00105     CalculateMomenta( vec, vecLen,
00106                       originalIncident, originalTarget, modifiedOriginal,
00107                       targetNucleus, currentParticle, targetParticle,
00108                       incidentHasChanged, targetHasChanged, quasiElastic );
00109     
00110     SetUpChange( vec, vecLen,
00111                  currentParticle, targetParticle,
00112                  incidentHasChanged );
00113     
00114     delete originalTarget;
00115     return &theParticleChange;
00116 }


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