Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4eplusAnnihilation.cc
이 파일의 문서화 페이지로 가기
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: G4eplusAnnihilation.cc 109177 2018-04-03 06:55:14Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4eplusAnnihilation
34 //
35 // Author: Vladimir Ivanchenko on base of Michel Maire code
36 //
37 // Creation date: 02.08.2004
38 //
39 // Modifications:
40 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivanchenko)
41 // 08-04-05 Major optimisation of internal interfaces (V.Ivanchenko)
42 // 03-05-05 suppress Integral option (mma)
43 // 04-05-05 Make class to be default (V.Ivanchenko)
44 // 25-01-06 remove cut dependance in AtRestDoIt (mma)
45 // 09-08-06 add SetModel(G4VEmModel*) (mma)
46 // 12-09-06, move SetModel(G4VEmModel*) in G4VEmProcess (mma)
47 // 30-05-12 propagate parent weight to secondaries (D. Sawkey)
48 //
49 
50 //
51 // -------------------------------------------------------------------
52 //
53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
55 
56 #include "G4eplusAnnihilation.hh"
57 #include "G4PhysicalConstants.hh"
58 #include "G4MaterialCutsCouple.hh"
59 #include "G4Gamma.hh"
60 #include "G4Positron.hh"
61 #include "G4eeToTwoGammaModel.hh"
62 #include "G4EmBiasingManager.hh"
63 
64 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
65 
66 using namespace std;
67 
69  : G4VEmProcess(name), isInitialised(false)
70 {
72  SetIntegral(true);
73  SetBuildTableFlag(false);
74  SetStartFromNullFlag(false);
77  enableAtRestDoIt = true;
78  mainSecondaries = 2;
79 }
80 
81 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
82 
84 {}
85 
86 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
87 
89 {
90  return (&p == G4Positron::Positron());
91 }
92 
93 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
94 
97 {
98  *condition = NotForced;
99  return 0.0;
100 }
101 
102 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
103 
105 {
106  if(!isInitialised) {
107  isInitialised = true;
108  if(!EmModel(0)) { SetEmModel(new G4eeToTwoGammaModel()); }
111  AddEmModel(1, EmModel(0));
112  }
113 }
114 
115 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
116 
118  G4String) const
119 {}
120 
121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
122 
124  const G4Step& step)
125 // Performs the e+ e- annihilation when both particles are assumed at rest.
126 {
128  size_t idx = CurrentMaterialCutsCoupleIndex();
129  G4double ene(0.0);
130  G4VEmModel* model = SelectModel(ene, idx);
131 
132  // define new weight for primary and secondaries
134 
135  // sample secondaries
136  secParticles.clear();
137  G4double gammaCut = GetGammaEnergyCut();
139  track.GetDynamicParticle(), gammaCut);
140 
141  G4int num0 = secParticles.size();
142 
143  // splitting or Russian roulette
144  if(biasManager) {
146  G4double eloss = 0.0;
148  secParticles, track, model, &fParticleChange, eloss,
149  idx, gammaCut, step.GetPostStepPoint()->GetSafety());
150  if(eloss > 0.0) {
153  }
154  }
155  }
156  // save secondaries
157  G4int num = secParticles.size();
158  if(num > 0) {
159 
162  G4double time = track.GetGlobalTime();
163 
164  for (G4int i=0; i<num; ++i) {
165  if (secParticles[i]) {
168  G4double e = dp->GetKineticEnergy();
169  G4bool good = true;
170  if(ApplyCuts()) {
171  if (p == theGamma) {
172  if (e < gammaCut) { good = false; }
173  } else if (p == theElectron) {
174  if (e < GetElectronEnergyCut()) { good = false; }
175  }
176  // added secondary if it is good
177  }
178  if (good) {
179  G4Track* t = new G4Track(dp, time, track.GetPosition());
181  t->SetWeight(weight);
183 
184  // define type of secondary
186  else if(i < num0) {
187  if(p == theGamma) {
189  } else {
191  }
192  } else {
194  }
195  /*
196  G4cout << "Secondary(post step) has weight " << t->GetWeight()
197  << ", Ekin= " << t->GetKineticEnergy()/MeV << " MeV "
198  << GetProcessName() << " fluoID= " << fluoID
199  << " augerID= " << augerID <<G4endl;
200  */
201  } else {
202  delete dp;
203  edep += e;
204  }
205  }
206  }
208  }
209  return &fParticleChange;
210 }
211 
212 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
213 
214 void G4eplusAnnihilation::ProcessDescription(std::ostream& out) const
215 {
216  out << "<strong>Positron annihilation</strong>";
218 }
219 
220 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
virtual G4VParticleChange * AtRestDoIt(const G4Track &track, const G4Step &stepData) override
void AddSecondary(G4Track *aSecondary)
const XML_Char * name
Definition: expat.h:151
void SetWeight(G4double aValue)
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:707
void SetCreatorModelIndex(G4int idx)
const char * p
Definition: xmltok.h:285
void SetBuildTableFlag(G4bool val)
G4int mainSecondaries
G4VEmModel * EmModel(size_t index=0) const
const G4TouchableHandle & GetTouchableHandle() const
const G4MaterialCutsCouple * MaterialCutsCouple() const
virtual void ProcessDescription(std::ostream &) const override
G4double ApplySecondaryBiasing(std::vector< G4DynamicParticle * > &, const G4Track &track, G4VEmModel *currentModel, G4ParticleChangeForGamma *pParticleChange, G4double &eloss, G4int coupleIdx, G4double tcut, G4double safety=0.0)
G4double condition(const G4ErrorSymMatrix &m)
G4bool SecondaryBiasingRegion(G4int coupleIdx)
G4VEmModel * SelectModel(G4double &kinEnergy, size_t index)
const G4ParticleDefinition * theGamma
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
G4double GetLocalEnergyDeposit() const
size_t CurrentMaterialCutsCoupleIndex() const
G4bool enableAtRestDoIt
Definition: G4VProcess.hh:350
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
virtual void ProcessDescription(std::ostream &outFile) const override
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
void SetEmModel(G4VEmModel *, G4int index=0)
G4double GetElectronEnergyCut()
G4ParticleChangeForGamma fParticleChange
G4double MinKinEnergy() const
const G4ParticleDefinition * theElectron
G4double GetGlobalTime() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
std::vector< G4DynamicParticle * > secParticles
const G4ThreeVector & GetPosition() const
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &track, G4ForceCondition *condition) override
static G4Positron * Positron()
Definition: G4Positron.cc:94
G4EmBiasingManager * biasManager
void SetSecondaryParticle(const G4ParticleDefinition *p)
Definition: G4Step.hh:76
G4StepPoint * GetPostStepPoint() const
G4eplusAnnihilation(const G4String &name="annihil")
Double_t edep
void InitializeForPostStep(const G4Track &)
G4double MaxKinEnergy() const
double weight
Definition: plottest35.C:25
int G4int
Definition: G4Types.hh:78
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:714
G4double GetSafety() const
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:435
G4ForceCondition
G4bool ApplyCuts() const
virtual void StreamProcessInfo(std::ostream &outFile, G4String endOfLine=G4String("\n")) const override
G4double GetKineticEnergy() const
virtual G4bool IsApplicable(const G4ParticleDefinition &p) final
void AddEmModel(G4int, G4VEmModel *, const G4Region *region=nullptr)
const G4ParticleDefinition * GetParticleDefinition() const
virtual void InitialiseProcess(const G4ParticleDefinition *) override
void SetStartFromNullFlag(G4bool val)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
void SetIntegral(G4bool val)
G4double GetGammaEnergyCut()
const XML_Char XML_Content * model
Definition: expat.h:151
const G4DynamicParticle * GetDynamicParticle() const
G4double GetParentWeight() const
void SetNumberOfSecondaries(G4int totSecondaries)