Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4MuonVDNuclearModel Class Reference

#include <G4MuonVDNuclearModel.hh>

Inheritance diagram for G4MuonVDNuclearModel:
Collaboration diagram for G4MuonVDNuclearModel:

Public Member Functions

 G4MuonVDNuclearModel ()
 
virtual ~G4MuonVDNuclearModel ()
 
virtual G4HadFinalStateApplyYourself (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
virtual void ModelDescription (std::ostream &outFile) const
 
- Public Member Functions inherited from G4HadronicInteraction
 G4HadronicInteraction (const G4String &modelName="HadronicModel")
 
virtual ~G4HadronicInteraction ()
 
virtual G4double SampleInvariantT (const G4ParticleDefinition *p, G4double plab, G4int Z, G4int A)
 
virtual G4bool IsApplicable (const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)
 
G4double GetMinEnergy () const
 
G4double GetMinEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMinEnergy (G4double anEnergy)
 
void SetMinEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMinEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4double GetMaxEnergy () const
 
G4double GetMaxEnergy (const G4Material *aMaterial, const G4Element *anElement) const
 
void SetMaxEnergy (const G4double anEnergy)
 
void SetMaxEnergy (G4double anEnergy, const G4Element *anElement)
 
void SetMaxEnergy (G4double anEnergy, const G4Material *aMaterial)
 
G4int GetVerboseLevel () const
 
void SetVerboseLevel (G4int value)
 
const G4StringGetModelName () const
 
void DeActivateFor (const G4Material *aMaterial)
 
void ActivateFor (const G4Material *aMaterial)
 
void DeActivateFor (const G4Element *anElement)
 
void ActivateFor (const G4Element *anElement)
 
G4bool IsBlocked (const G4Material *aMaterial) const
 
G4bool IsBlocked (const G4Element *anElement) const
 
void SetRecoilEnergyThreshold (G4double val)
 
G4double GetRecoilEnergyThreshold () const
 
virtual const std::pair
< G4double, G4double
GetFatalEnergyCheckLevels () const
 
virtual std::pair< G4double,
G4double
GetEnergyMomentumCheckLevels () const
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void InitialiseModel ()
 

Additional Inherited Members

- Protected Member Functions inherited from G4HadronicInteraction
void SetModelName (const G4String &nam)
 
G4bool IsBlocked () const
 
void Block ()
 
- Protected Attributes inherited from G4HadronicInteraction
G4HadFinalState theParticleChange
 
G4int verboseLevel
 
G4double theMinEnergy
 
G4double theMaxEnergy
 
G4bool isBlocked
 

Detailed Description

Definition at line 50 of file G4MuonVDNuclearModel.hh.

Constructor & Destructor Documentation

G4MuonVDNuclearModel::G4MuonVDNuclearModel ( )
explicit

Definition at line 72 of file G4MuonVDNuclearModel.cc.

73  : G4HadronicInteraction("G4MuonVDNuclearModel"),isMaster(false)
74 {
76  GetCrossSectionDataSet(G4KokoulinMuonNuclearXS::Default_Name());
77 
78  SetMinEnergy(0.0);
80  CutFixed = 0.2*CLHEP::GeV;
81 
82  if(!fElementData && G4Threading::IsMasterThread()) {
83  fElementData = new G4ElementData();
84  MakeSamplingTable();
85  isMaster = true;
86  }
87 
88  // reuse existing pre-compound model
89  G4GeneratorPrecompoundInterface* precoInterface
93  G4VPreCompoundModel* pre = static_cast<G4VPreCompoundModel*>(p);
94  if(!pre) { pre = new G4PreCompoundModel(); }
95  precoInterface->SetDeExcitation(pre);
96 
97  // Build FTFP model
98  ftfp = new G4TheoFSGenerator();
99  ftfp->SetTransport(precoInterface);
100  theFragmentation = new G4LundStringFragmentation();
101  theStringDecay = new G4ExcitedStringDecay(theFragmentation);
102  G4FTFModel* theStringModel = new G4FTFModel;
103  theStringModel->SetFragmentationModel(theStringDecay);
104  ftfp->SetHighEnergyGenerator(theStringModel);
105 
106  // Build Bertini cascade
107  bert = new G4CascadeInterface();
108 }
void SetFragmentationModel(G4VStringFragmentation *aModel)
const char * p
Definition: xmltok.h:285
void SetHighEnergyGenerator(G4VHighEnergyGenerator *const value)
void SetMinEnergy(G4double anEnergy)
static const char * Default_Name()
static G4CrossSectionDataSetRegistry * Instance()
G4HadronicInteraction(const G4String &modelName="HadronicModel")
G4HadronicInteraction * FindModel(const G4String &name)
static constexpr double GeV
static G4HadronicInteractionRegistry * Instance()
void SetMaxEnergy(const G4double anEnergy)
void SetDeExcitation(G4VPreCompoundModel *ptr)
void SetTransport(G4VIntraNuclearTransportModel *const value)
static constexpr double PeV
G4bool IsMasterThread()
Definition: G4Threading.cc:146

Here is the call graph for this function:

G4MuonVDNuclearModel::~G4MuonVDNuclearModel ( )
virtual

Definition at line 110 of file G4MuonVDNuclearModel.cc.

111 {
112  delete theFragmentation;
113  delete theStringDecay;
114 
115  if(isMaster) {
116  delete fElementData;
117  fElementData = nullptr;
118  }
119 }

Member Function Documentation

G4HadFinalState * G4MuonVDNuclearModel::ApplyYourself ( const G4HadProjectile aTrack,
G4Nucleus targetNucleus 
)
virtual

Implements G4HadronicInteraction.

Definition at line 122 of file G4MuonVDNuclearModel.cc.

124 {
126 
127  // For very low energy, return initial track
128  G4double epmax = aTrack.GetTotalEnergy() - 0.5*proton_mass_c2;
129  if (epmax <= CutFixed) {
133  return &theParticleChange;
134  }
135 
136  // Produce recoil muon and transferred photon
137  G4DynamicParticle* transferredPhoton = CalculateEMVertex(aTrack, targetNucleus);
138 
139  // Interact the gamma with the nucleus
140  CalculateHadronicVertex(transferredPhoton, targetNucleus);
141  return &theParticleChange;
142 }
void SetStatusChange(G4HadFinalStateStatus aS)
Hep3Vector vect() const
G4double GetKineticEnergy() const
float proton_mass_c2
Definition: hepunit.py:275
const G4LorentzVector & Get4Momentum() const
void SetEnergyChange(G4double anEnergy)
Hep3Vector unit() const
double G4double
Definition: G4Types.hh:76
void SetMomentumChange(const G4ThreeVector &aV)
G4double GetTotalEnergy() const

Here is the call graph for this function:

void G4MuonVDNuclearModel::ModelDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronicInteraction.

Definition at line 389 of file G4MuonVDNuclearModel.cc.

390 {
391  outFile << "G4MuonVDNuclearModel handles the inelastic scattering\n"
392  << "of mu- and mu+ from nuclei using the equivalent photon\n"
393  << "approximation in which the incoming lepton generates a\n"
394  << "virtual photon at the electromagnetic vertex, and the\n"
395  << "virtual photon is converted to a real photon. At low\n"
396  << "energies, the photon interacts directly with the nucleus\n"
397  << "using the Bertini cascade. At high energies the photon\n"
398  << "is converted to a pi0 which interacts using the FTFP\n"
399  << "model. The muon-nuclear cross sections of R. Kokoulin \n"
400  << "are used to generate the virtual photon spectrum\n";
401 }

The documentation for this class was generated from the following files: