ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4WentzelVIModel.hh
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4WentzelVIModel.hh
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 //
27 // -------------------------------------------------------------------
28 //
29 //
30 // GEANT4 Class header file
31 //
32 //
33 // File name: G4WentzelVIModel
34 //
35 // Author: V.Ivanchenko
36 //
37 // Creation date: 09.04.2008 from G4MuMscModel
38 //
39 // Modifications:
40 // 27-05-2010 V.Ivanchenko added G4WentzelOKandVIxSection class to
41 // compute cross sections and sample scattering angle
42 //
43 // Class Description:
44 //
45 // Implementation of the model of multiple scattering based on
46 // G.Wentzel, Z. Phys. 40 (1927) 590.
47 // H.W.Lewis, Phys Rev 78 (1950) 526.
48 // J.M. Fernandez-Varea et al., NIM B73 (1993) 447.
49 // L.Urban, CERN-OPEN-2006-077.
50 
51 // -------------------------------------------------------------------
52 //
53 
54 #ifndef G4WentzelVIModel_h
55 #define G4WentzelVIModel_h 1
56 
57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
58 
59 #include "G4VMscModel.hh"
60 #include "G4MaterialCutsCouple.hh"
62 
63 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
64 
66 {
67 
68 public:
69 
70  explicit G4WentzelVIModel(G4bool comb=true, const G4String& nam = "WentzelVIUni");
71 
72  virtual ~G4WentzelVIModel();
73 
74  virtual void Initialise(const G4ParticleDefinition*,
75  const G4DataVector&) override;
76 
77  virtual void InitialiseLocal(const G4ParticleDefinition*,
78  G4VEmModel* masterModel) override;
79 
80  virtual void StartTracking(G4Track*) override;
81 
83  G4double KineticEnergy,
84  G4double AtomicNumber,
85  G4double AtomicWeight=0.,
86  G4double cut = DBL_MAX,
87  G4double emax= DBL_MAX) override;
88 
90  G4double safety) override;
91 
92  virtual G4double
94  G4double& currentMinimalStep) override;
95 
96  virtual G4double ComputeGeomPathLength(G4double truePathLength) override;
97 
98  virtual G4double ComputeTrueStepLength(G4double geomStepLength) override;
99 
100  // defines low energy limit on energy transfer to atomic electron
101  inline void SetFixedCut(G4double);
102 
103  // low energy limit on energy transfer to atomic electron
104  inline G4double GetFixedCut() const;
105 
106  // access to cross section class
108 
110 
111  inline void SetUseSecondMoment(G4bool);
112 
113  inline G4bool UseSecondMoment() const;
114 
116 
118  const G4MaterialCutsCouple*,
119  G4double kineticEnergy);
120 
122 
124 
125 protected:
126 
128 
129  inline void SetupParticle(const G4ParticleDefinition*);
130 
131 private:
132 
134  G4double kineticEnergy);
135 
136  // hide assignment operator
138  G4WentzelVIModel(const G4WentzelVIModel&) = delete;
139 
140 protected:
141 
143 
147 
148  // cache kinematics
155 
156  // cache material
160 
164 
167 
168  // cache kinematics
170 
171  // single scattering parameters
174 
176  size_t idx2;
177 
178  // data for single scattering mode
181  std::vector<G4double> xsecn;
182  std::vector<G4double> prob;
183 
186 
187  // projectile
189 
190  // flags
194 };
195 
196 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
197 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
198 
200 {
201  // Initialise mass and charge
202  if(p != particle) {
203  particle = p;
204  wokvi->SetupParticle(p);
205  }
206 }
207 
208 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
209 
211 {
212  fixedCut = val;
213 }
214 
215 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
216 
218 {
219  return fixedCut;
220 }
221 
222 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
223 
225 {
226  if(ptr != wokvi) {
227  delete wokvi;
228  wokvi = ptr;
229  }
230 }
231 
232 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
233 
235 {
236  return wokvi;
237 }
238 
239 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
240 
242 {
243  useSecondMoment = val;
244 }
245 
246 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
247 
249 {
250  return useSecondMoment;
251 }
252 
253 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
254 
256 {
257  return fSecondMoments;
258 }
259 
260 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
261 
262 inline G4double
264  const G4MaterialCutsCouple* couple,
265  G4double ekin)
266 {
267  G4double x = 0.0;
268  if(useSecondMoment) {
269  DefineMaterial(couple);
270  x = (fSecondMoments) ?
271  (*fSecondMoments)[(*theDensityIdx)[currentMaterialIndex]]->Value(ekin, idx2)
272  *(*theDensityFactor)[currentMaterialIndex]/(ekin*ekin)
273  : ComputeSecondMoment(part, ekin);
274  }
275  return x;
276 }
277 
278 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
279 
280 #endif
281