ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4UrbanMscModel.hh
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4UrbanMscModel.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: G4UrbanMscModel
34 //
35 // Author: Laszlo Urban
36 //
37 // Creation date: 19.02.2013
38 //
39 // Created from G4UrbanMscModel96
40 //
41 // New parametrization for theta0
42 // Correction for very small step length
43 //
44 // Class Description:
45 //
46 // Implementation of the model of multiple scattering based on
47 // H.W.Lewis Phys Rev 78 (1950) 526 and L.Urban model
48 
49 // -------------------------------------------------------------------
50 //
51 
52 #ifndef G4UrbanMscModel_h
53 #define G4UrbanMscModel_h 1
54 
55 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56 
58 
59 #include "G4VMscModel.hh"
60 #include "G4MscStepLimitType.hh"
61 #include "G4Log.hh"
62 #include "G4Exp.hh"
63 
65 class G4SafetyHelper;
66 class G4LossTableManager;
67 
68 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
69 
71 {
72 
73 public:
74 
75  explicit G4UrbanMscModel(const G4String& nam = "UrbanMsc");
76 
77  ~G4UrbanMscModel() override;
78 
79  void Initialise(const G4ParticleDefinition*,
80  const G4DataVector&) override;
81 
82  void StartTracking(G4Track*) override;
83 
84  G4double
86  G4double KineticEnergy,
87  G4double AtomicNumber,
88  G4double AtomicWeight=0.,
89  G4double cut =0.,
90  G4double emax=DBL_MAX) override;
91 
93  G4double safety) override;
94 
96  G4double& currentMinimalStep) override;
97 
98  G4double ComputeGeomPathLength(G4double truePathLength) override;
99 
100  G4double ComputeTrueStepLength(G4double geomStepLength) override;
101 
102  G4double ComputeTheta0(G4double truePathLength, G4double KineticEnergy);
103 
104 private:
105 
106  G4double SampleCosineTheta(G4double trueStepLength, G4double KineticEnergy);
107 
108  void SampleDisplacement(G4double sinTheta, G4double phi);
109 
110  void SampleDisplacementNew(G4double sinTheta, G4double phi);
111 
112  inline void SetParticle(const G4ParticleDefinition*);
113 
114  inline void UpdateCache();
115 
116  inline G4double Randomizetlimit();
117 
118  inline G4double SimpleScattering(G4double xmeanth, G4double x2meanth);
119 
120  // hide assignment operator
121  G4UrbanMscModel & operator=(const G4UrbanMscModel &right) = delete;
122  G4UrbanMscModel(const G4UrbanMscModel&) = delete;
123 
125 
129 
132 
136 
145 
151 
153 
159 
162 
168 
170 
175 
178 
181 
184 
185 };
186 
187 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
188 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
189 
190 inline
192 {
193  if (p != particle) {
194  particle = p;
195  mass = p->GetPDGMass();
198  }
199 }
200 
201 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
202 
204 {
205  G4double res = tlimitmin;
206  if(tlimit > tlimitmin)
207  {
209  res = std::max(res, tlimitmin);
210  }
211  return res;
212 }
213 
214 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
215 
217 {
218  lnZ = G4Log(Zeff);
219  // correction in theta0 formula
220  G4double w = G4Exp(lnZ/6.);
221  G4double facz = 0.990395+w*(-0.168386+w*0.093286) ;
222  coeffth1 = facz*(1. - 8.7780e-2/Zeff);
223  coeffth2 = facz*(4.0780e-2 + 1.7315e-4*Zeff);
224 
225  // tail parameters
226  G4double Z13 = w*w;
227  coeffc1 = 2.3785 - Z13*(4.1981e-1 - Z13*6.3100e-2);
228  coeffc2 = 4.7526e-1 + Z13*(1.7694 - Z13*3.3885e-1);
229  coeffc3 = 2.3683e-1 - Z13*(1.8111 - Z13*3.2774e-1);
230  coeffc4 = 1.7888e-2 + Z13*(1.9659e-2 - Z13*2.6664e-3);
231 
232  Z2 = Zeff*Zeff;
233  Z23 = Z13*Z13;
234 
235  stepmina = 15.99/(1.+0.119*Zeff);
236  stepminb = 4.39/(1.+0.079*Zeff);
237 
238  Zold = Zeff;
239 }
240 
241 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
242 
243 inline
245 {
246  // 'large angle scattering'
247  // 2 model functions with correct xmean and x2mean
248  G4double a = (2.*xmeanth+9.*x2meanth-3.)/(2.*xmeanth-3.*x2meanth+1.);
249  G4double prob = (a+2.)*xmeanth/a;
250 
251  // sampling
252  G4double rdm = rndmEngineMod->flat();
253  G4double cth = (rndmEngineMod->flat() < prob)
254  ? -1.+2.*G4Exp(G4Log(rdm)/(a+1.)) : -1.+2.*rdm;
255  return cth;
256 }
257 
258 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
259 
260 
261 #endif
262