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

#include <G4MuonMinusAtomicCapture.hh>

Inheritance diagram for G4MuonMinusAtomicCapture:
Collaboration diagram for G4MuonMinusAtomicCapture:

Public Member Functions

 G4MuonMinusAtomicCapture (const G4String &name="muMinusAtomicCaptureAtRest")
 
 ~G4MuonMinusAtomicCapture ()
 
G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
void ProcessDescription (std::ostream &outFile) const
 
void SetElementSelector (G4ElementSelector *ptr)
 
void SetEmCascade (G4HadronicInteraction *ptr)
 
- Public Member Functions inherited from G4HadronicProcess
 G4HadronicProcess (const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)
 
 G4HadronicProcess (const G4String &processName, G4HadronicProcessType subType)
 
virtual ~G4HadronicProcess ()
 
void RegisterMe (G4HadronicInteraction *a)
 
G4double GetElementCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
 
G4double GetMicroscopicCrossSection (const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &aTrack, const G4Step &aStep)
 
void DumpPhysicsTable (const G4ParticleDefinition &p)
 
void AddDataSet (G4VCrossSectionDataSet *aDataSet)
 
std::vector
< G4HadronicInteraction * > & 
GetHadronicInteractionList ()
 
G4double GetMeanFreePath (const G4Track &aTrack, G4double, G4ForceCondition *)
 
const G4NucleusGetTargetNucleus () const
 
const G4IsotopeGetTargetIsotope ()
 
void BiasCrossSectionByFactor (G4double aScale)
 
void SetIntegral (G4bool val)
 
void SetEpReportLevel (G4int level)
 
void SetEnergyMomentumCheckLevels (G4double relativeLevel, G4double absoluteLevel)
 
std::pair< G4double, G4doubleGetEnergyMomentumCheckLevels () const
 
G4CrossSectionDataStoreGetCrossSectionDataStore ()
 
void MultiplyCrossSectionBy (G4double factor)
 
- Public Member Functions inherited from G4VDiscreteProcess
 G4VDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VDiscreteProcess (G4VDiscreteProcess &)
 
virtual ~G4VDiscreteProcess ()
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

G4double GetMeanLifeTime (const G4Track &, G4ForceCondition *)
 
- Protected Member Functions inherited from G4HadronicProcess
G4HadronicInteractionChooseHadronicInteraction (const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, G4Material *aMaterial, G4Element *anElement)
 
G4NucleusGetTargetNucleusPointer ()
 
void DumpState (const G4Track &, const G4String &, G4ExceptionDescription &)
 
G4HadronicInteractionGetHadronicInteraction () const
 
G4double GetLastCrossSection ()
 
void FillResult (G4HadFinalState *aR, const G4Track &aT)
 
G4HadFinalStateCheckResult (const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
 
void CheckEnergyMomentumConservation (const G4Track &, const G4Nucleus &)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4HadronicProcess
G4HadProjectile thePro
 
G4ParticleChangetheTotalResult
 
G4int epReportLevel
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 

Detailed Description

Definition at line 61 of file G4MuonMinusAtomicCapture.hh.

Constructor & Destructor Documentation

G4MuonMinusAtomicCapture::G4MuonMinusAtomicCapture ( const G4String name = "muMinusAtomicCaptureAtRest")
explicit

Definition at line 62 of file G4MuonMinusAtomicCapture.cc.

63  : G4HadronicProcess(name, fHadronAtRest),// name, process type
64  fElementSelector(new G4ElementSelector()),
65  fEmCascade(new G4EmCaptureCascade()) // Owned by InteractionRegistry
66 {
67  // Modify G4VProcess flags to emulate G4VRest instead of G4VDiscrete
68  enableAtRestDoIt = true;
69  enablePostStepDoIt = false;
71 }
static G4HadronicProcessStore * Instance()
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
G4bool enablePostStepDoIt
Definition: G4VProcess.hh:352
void RegisterExtraProcess(G4VProcess *)
G4HadronicProcess(const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)

Here is the call graph for this function:

G4MuonMinusAtomicCapture::~G4MuonMinusAtomicCapture ( )

Definition at line 75 of file G4MuonMinusAtomicCapture.cc.

76 {}

Member Function Documentation

G4VParticleChange * G4MuonMinusAtomicCapture::AtRestDoIt ( const G4Track track,
const G4Step  
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 119 of file G4MuonMinusAtomicCapture.cc.

121 {
122  // if primary is not Alive then do nothing (how?)
123  theTotalResult->Initialize(track);
124 
125  G4Nucleus* nucleus = GetTargetNucleusPointer();
126  // the call below actually sets the nucleus params;
127  // G4Nucleus targetNucleus; is a member of G4HadronicProcess
128  // G4Element* elm =
129  fElementSelector->SelectZandA(track, nucleus);
130 
131  G4HadFinalState* result = 0;
132  thePro.Initialise(track); // thePro is G4HadProjectile from G4HadronicProcess
133 
134  // save track time an dstart capture from zero time
135  thePro.SetGlobalTime(0.0);
136  G4double time0 = track.GetGlobalTime();
137 
138  // Do the electromagnetic cascade in the nuclear field.
139  // EM cascade should keep G4HadFinalState object,
140  // because it will not be deleted at the end of this method
141  //
142  result = fEmCascade->ApplyYourself(thePro, *nucleus);
143  G4double ebound = result->GetLocalEnergyDeposit(); // may need to carry this over; review
144  G4double edep = 0.0;
145  G4int nSecondaries = result->GetNumberOfSecondaries();
146  thePro.SetBoundEnergy(ebound);
147 
148  // creating the muonic atom
149  ++nSecondaries;
150 
152  G4ParticleDefinition* muonicAtom = itp->GetMuonicAtom(nucleus->GetZ_asInt(),
153  nucleus->GetA_asInt());
154 
155  G4DynamicParticle* dp = new G4DynamicParticle(muonicAtom,G4RandomDirection(),0.);
156  G4HadSecondary hadSec(dp);
157  hadSec.SetTime(time0);
158  result->AddSecondary(hadSec);
159 
160  // Fill results
161  //
164  theTotalResult->SetNumberOfSecondaries(nSecondaries);
165  G4double w = track.GetWeight();
167 
168  G4cout << __func__
169  << " nSecondaries "
170  << nSecondaries
171  << G4endl;
172 
173  for(G4int i=0; i<nSecondaries; ++i) {
174  G4HadSecondary* sec = result->GetSecondary(i);
175 
176  // add track global time to the reaction time
177  G4double time = sec->GetTime();
178  if(time < 0.0) { time = 0.0; }
179  time += time0;
180 
181  G4cout << __func__
182  << " "
183  << i
184  << " Resulting secondary "
185  << sec->GetParticle()->GetPDGcode()
186  << " "
188  << G4endl;
189 
190  // create secondary track
191  G4Track* t = new G4Track(sec->GetParticle(),
192  time,
193  track.GetPosition());
194  t->SetWeight(w*sec->GetWeight());
195 
198  }
199  result->Clear();
200 
201  // fixme: needs to be done at the MuonicAtom level
202  // if (epReportLevel != 0) { // G4HadronicProcess::
203  // CheckEnergyMomentumConservation(track, *nucleus);
204  // }
205  return theTotalResult;
206 }
G4double G4ParticleHPJENDLHEData::G4double result
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
virtual G4Element * SelectZandA(const G4Track &track, G4Nucleus *)
G4HadSecondary * GetSecondary(size_t i)
G4int GetPDGcode() const
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
void SetGlobalTime(G4double t)
G4ThreeVector G4RandomDirection()
G4HadProjectile thePro
G4ParticleDefinition * GetDefinition() const
G4double GetTime() const
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void ProposeWeight(G4double finalWeight)
G4GLOB_DLL std::ostream G4cout
G4ParticleChange * theTotalResult
G4double GetGlobalTime() const
const G4TouchableHandle & GetTouchableHandle() const
void Initialise(const G4Track &aT)
static G4IonTable * GetIonTable()
Definition: G4IonTable.hh:78
virtual void Initialize(const G4Track &)
G4Nucleus * GetTargetNucleusPointer()
void SetNumberOfSecondaries(G4int totSecondaries)
G4DynamicParticle * GetParticle()
G4ParticleDefinition * GetMuonicAtom(G4Ions const *)
Definition: G4IonTable.cc:1871
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void SetBoundEnergy(G4double e)
void AddSecondary(G4Track *aSecondary)
G4double GetWeight() const
#define G4endl
Definition: G4ios.hh:61
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
G4int GetNumberOfSecondaries() const
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
G4double GetLocalEnergyDeposit() const
G4double GetWeight() const

Here is the call graph for this function:

G4double G4MuonMinusAtomicCapture::AtRestGetPhysicalInteractionLength ( const G4Track track,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 101 of file G4MuonMinusAtomicCapture.cc.

103 {
104  *condition = NotForced;
105  return 0.0;
106 }
G4double condition(const G4ErrorSymMatrix &m)
void G4MuonMinusAtomicCapture::BuildPhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 94 of file G4MuonMinusAtomicCapture.cc.

95 {
97 }
static G4HadronicProcessStore * Instance()
void PrintInfo(const G4ParticleDefinition *)

Here is the call graph for this function:

G4double G4MuonMinusAtomicCapture::GetMeanLifeTime ( const G4Track ,
G4ForceCondition  
)
inlineprotected

Definition at line 95 of file G4MuonMinusAtomicCapture.hh.

96  { return -1.0; }
G4bool G4MuonMinusAtomicCapture::IsApplicable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4VProcess.

Definition at line 80 of file G4MuonMinusAtomicCapture.cc.

81 {
82  return (&p == G4MuonMinus::MuonMinus());
83 }
static G4MuonMinus * MuonMinus()
Definition: G4MuonMinus.cc:100

Here is the call graph for this function:

G4double G4MuonMinusAtomicCapture::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VDiscreteProcess.

Definition at line 110 of file G4MuonMinusAtomicCapture.cc.

112 {
113  *condition = NotForced;
114  return DBL_MAX;
115 }
G4double condition(const G4ErrorSymMatrix &m)
#define DBL_MAX
Definition: templates.hh:83
void G4MuonMinusAtomicCapture::PreparePhysicsTable ( const G4ParticleDefinition p)
virtual

Reimplemented from G4HadronicProcess.

Definition at line 87 of file G4MuonMinusAtomicCapture.cc.

88 {
90 }
static G4HadronicProcessStore * Instance()
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)

Here is the call graph for this function:

void G4MuonMinusAtomicCapture::ProcessDescription ( std::ostream &  outFile) const
virtual

Reimplemented from G4HadronicProcess.

Definition at line 210 of file G4MuonMinusAtomicCapture.cc.

211 {
212  outFile << "Stopping of mu- using default element selector, EM cascade"
213  << " sampling and bound decay sampling.\n"
214  << "Bertini model is used for nuclear capture\n"
215  << "G4MuonicAtom is created\n";
216 }
void G4MuonMinusAtomicCapture::SetElementSelector ( G4ElementSelector ptr)
inline
void G4MuonMinusAtomicCapture::SetEmCascade ( G4HadronicInteraction ptr)
inline

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