Geant4  10.03
G4MuonMinusAtomicCapture.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4MuonMinusAtomicCapture.cc $
27 //
28 //---------------------------------------------------------------------
29 //
30 // GEANT4 Class
31 //
32 // GEANT4 Class header file
33 //
34 // File name: G4MuonMinusAtomicCapture
35 //
36 // 20160701 K.L. Genser - New process using G4MuonicAtom somewhat based on G4HadronStoppingProcess
37 //
38 // Class Description:
39 //
40 // Stopping of mu-
41 //
42 // G4VParticleChange will contain gammas from G4EmCaptureCascade and
43 // resulting G4MuonicAtom
44 //
45 //
46 //------------------------------------------------------------------------
47 
49 #include "G4ParticleDefinition.hh"
50 #include "G4HadronicProcessType.hh"
51 #include "G4MuonMinusBoundDecay.hh"
52 #include "G4HadronicInteraction.hh"
54 #include "G4EmCaptureCascade.hh"
55 #include "G4MuonMinus.hh"
56 #include "G4IonTable.hh"
57 #include "G4RandomDirection.hh"
58 #include "G4HadSecondary.hh"
59 
60 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
61 
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 }
72 
73 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
74 
76 {}
77 
78 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
79 
81 {
82  return (&p == G4MuonMinus::MuonMinus());
83 }
84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
85 
86 void
88 {
90 }
91 
92 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
93 
95 {
97 }
98 
99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
100 
103 {
104  *condition = NotForced;
105  return 0.0;
106 }
107 
108 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
109 
112 {
113  *condition = NotForced;
114  return DBL_MAX;
115 }
116 
117 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
118 
120  const G4Step&)
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 }
207 
208 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
209 
210 void G4MuonMinusAtomicCapture::ProcessDescription(std::ostream& outFile) const
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 }
217 
218 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
G4double condition(const G4ErrorSymMatrix &m)
G4int GetA_asInt() const
Definition: G4Nucleus.hh:109
virtual G4Element * SelectZandA(const G4Track &track, G4Nucleus *)
G4HadSecondary * GetSecondary(size_t i)
static G4HadronicProcessStore * Instance()
G4int GetPDGcode() const
const G4ThreeVector & GetPosition() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
void SetGlobalTime(G4double t)
G4ThreeVector G4RandomDirection()
G4HadProjectile thePro
G4ParticleDefinition * GetDefinition() const
const char * name(G4int ptype)
virtual G4VParticleChange * AtRestDoIt(const G4Track &, const G4Step &)
G4double GetTime() const
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
void SetWeight(G4double aValue)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void ProposeWeight(G4double finalWeight)
void RegisterParticleForExtraProcess(G4VProcess *, const G4ParticleDefinition *)
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &track, G4ForceCondition *condition)
G4GLOB_DLL std::ostream G4cout
G4ParticleChange * theTotalResult
bool G4bool
Definition: G4Types.hh:79
Definition: G4Step.hh:76
G4bool enablePostStepDoIt
Definition: G4VProcess.hh:352
void ProcessDescription(std::ostream &outFile) const
G4double GetGlobalTime() const
const G4TouchableHandle & GetTouchableHandle() const
void Initialise(const G4Track &aT)
static G4IonTable * GetIonTable()
Definition: G4IonTable.hh:78
void RegisterExtraProcess(G4VProcess *)
virtual void Initialize(const G4Track &)
void SetTime(G4double aT)
virtual void PreparePhysicsTable(const G4ParticleDefinition &)
G4Nucleus * GetTargetNucleusPointer()
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
void SetNumberOfSecondaries(G4int totSecondaries)
G4DynamicParticle * GetParticle()
G4ParticleDefinition * GetMuonicAtom(G4Ions const *)
Definition: G4IonTable.cc:1871
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
G4bool IsApplicable(const G4ParticleDefinition &)
void SetBoundEnergy(G4double e)
void AddSecondary(G4Track *aSecondary)
G4double GetWeight() const
static G4MuonMinus * MuonMinus()
Definition: G4MuonMinus.cc:100
#define G4endl
Definition: G4ios.hh:61
void AddSecondary(G4DynamicParticle *aP, G4int mod=-1)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
G4ForceCondition
G4HadronicInteraction * fEmCascade
#define DBL_MAX
Definition: templates.hh:83
G4int GetNumberOfSecondaries() const
virtual G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &targetNucleus)=0
G4double GetLocalEnergyDeposit() const
void PrintInfo(const G4ParticleDefinition *)
G4MuonMinusAtomicCapture(const G4String &name="muMinusAtomicCaptureAtRest")
G4double GetWeight() const