Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4WentzelVIRelXSection.hh
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: G4WentzelVIRelXSection.hh 96934 2016-05-18 09:10:41Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 //
31 // GEANT4 Class header file
32 //
33 //
34 // File name: G4WentzelVIRelXSection
35 //
36 // Authors: V.Ivanchenko
37 //
38 // Creation date: 08.06.2012 from G4WentzelOKandVIxSection
39 //
40 // Modifications:
41 //
42 //
43 // Class Description:
44 //
45 // Implementation of the computation of total and transport cross sections,
46 // sample scattering angle for the single scattering case.
47 // to be used by single and multiple scattering models. References:
48 // 1) G.Wentzel, Z. Phys. 40 (1927) 590.
49 // 2) J.M. Fernandez-Varea et al., NIM B73 (1993) 447.
50 //
51 // -------------------------------------------------------------------
52 //
53 
54 #ifndef G4WentzelVIRelXSection_h
55 #define G4WentzelVIRelXSection_h 1
56 
57 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
58 
59 #include "globals.hh"
60 #include "G4Material.hh"
61 #include "G4Element.hh"
62 #include "G4ElementVector.hh"
63 #include "G4NistManager.hh"
64 #include "G4ThreeVector.hh"
65 #include "G4Pow.hh"
66 
68 
69 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
70 
72 {
73 
74 public:
75 
76  explicit G4WentzelVIRelXSection(G4bool combined = true);
77 
78  virtual ~G4WentzelVIRelXSection();
79 
80  void Initialise(const G4ParticleDefinition*, G4double CosThetaLim);
81 
83 
84  // return cos(ThetaMax) for msc and cos(thetaMin) for single scattering
85  // cut = DBL_MAX means no scattering off electrons
87 
89 
91  G4double CosThetaMax,
92  G4double elecRatio);
93 
94  inline G4double ComputeNuclearCrossSection(G4double CosThetaMin,
95  G4double CosThetaMax);
96 
98  G4double CosThetaMax);
99 
100  inline G4double SetupKinematic(G4double kinEnergy,
101  const G4Material* mat,
102  G4double cut,
103  G4double tmass);
104 
105  inline G4double GetMomentumSquare() const;
106 
107  inline G4double GetCosThetaNuc() const;
108 
109  inline G4double GetCosThetaElec() const;
110 
111 private:
112 
113  void ComputeMaxElectronScattering(G4double cut);
114 
115  // hide assignment operator
116  G4WentzelVIRelXSection & operator=
117  (const G4WentzelVIRelXSection &right) = delete;
119 
120  const G4ParticleDefinition* theProton;
121  const G4ParticleDefinition* theElectron;
122  const G4ParticleDefinition* thePositron;
123  const G4Material* currentMaterial;
124 
125  G4NistManager* fNistManager;
126  G4Pow* fG4pow;
127 
128  G4ThreeVector temp;
129 
130  G4double numlimit;
131 
132  // integer parameters
133  G4int nwarnings;
134  G4int nwarnlimit;
135 
136  G4bool isCombined;
137 
138  // single scattering parameters
139  G4double coeff;
140  G4double cosTetMaxElec;
141  G4double cosTetMaxNuc;
142  G4double cosThetaMax;
143  G4double alpha2;
144 
145  // projectile
146  const G4ParticleDefinition* particle;
147 
148  G4double chargeSquare;
149  G4double charge3;
150  G4double spin;
151  G4double mass;
152  G4double tkin;
153  G4double mom2;
154  G4double momCM2;
155  G4double invbeta2;
156  G4double kinFactor;
157  G4double etag;
158  G4double ecut;
159  G4double lowEnergyLimit;
160 
161  // target
162  G4int targetZ;
163  G4double targetMass;
164  G4double screenZ;
165  G4double formfactA;
166  G4double factorA2;
167  G4double factB;
168  G4double factB1;
169  G4double factD;
170  G4double gam0pcmp;
171  G4double pcmp2;
172 
173  static G4double ScreenRSquare[100];
174  static G4double FormFactor[100];
175 };
176 
177 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
178 
179 inline G4double
181  const G4Material* mat,
182  G4double cut,
183  G4double tmass)
184 {
185  if(kinEnergy != tkin || mat != currentMaterial ||
186  ecut != cut || tmass != targetMass) {
187 
188  currentMaterial = mat;
189  ecut = cut;
190  tkin = kinEnergy;
191  G4double momLab2 = tkin*(tkin + 2.0*mass);
192 
193  G4double etot = tkin + mass;
194  G4double ptot = std::sqrt(momLab2);
195  G4double m12 = mass*mass;
196 
197  targetMass = tmass;
198 
199  // relativistic reduced mass from publucation
200  // A.P. Martynenko, R.N. Faustov, Teoret. mat. Fiz. 64 (1985) 179
201 
202  //incident particle & target nucleus
203  G4double Ecm = std::sqrt(m12 + targetMass*targetMass + 2.0*etot*targetMass);
204  G4double mu_rel = mass*targetMass/Ecm;
205  G4double momCM = ptot*targetMass/Ecm;
206  // relative system
207  mom2 = momCM*momCM;
208  invbeta2 = 1.0 + mu_rel*mu_rel/mom2;
209 
210  factB = spin/invbeta2;
211  factD = std::sqrt(mom2)/tmass;
212  if(isCombined) {
213  G4double cost = 1.-factorA2*mat->GetIonisation()->GetInvA23()/mom2;
214  if(cost > cosTetMaxNuc) { cosTetMaxNuc = cost; }
215  }
216  }
217  return cosTetMaxNuc;
218 
219 }
220 
221 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
222 
224 {
225  return mom2;
226 }
227 
228 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
229 
231 {
232  return cosTetMaxNuc;
233 }
234 
235 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
236 
238 {
239  return cosTetMaxElec;
240 }
241 
242 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
243 
244 inline G4double
246  G4double cosTMax)
247 {
248  G4double xsec = 0.0;
249  if(cosTMax < cosTMin) {
250  xsec = targetZ*kinFactor*(cosTMin - cosTMax)/
251  ((1.0 - cosTMin + screenZ)*(1.0 - cosTMax + screenZ));
252  }
253  return xsec;
254 }
255 
256 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
257 
258 inline G4double
260  G4double cosTMax)
261 {
262  G4double xsec = 0.0;
263  G4double cost1 = std::max(cosTMin,cosTetMaxElec);
264  G4double cost2 = std::max(cosTMax,cosTetMaxElec);
265  if(cost1 > cost2) {
266  xsec = kinFactor*(cost1 - cost2)/
267  ((1.0 - cost1 + screenZ)*(1.0 - cost2 + screenZ));
268  }
269  return xsec;
270 }
271 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
272 
273 #endif
274 
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:226
G4double ComputeNuclearCrossSection(G4double CosThetaMin, G4double CosThetaMax)
G4WentzelVIRelXSection(G4bool combined=true)
Definition: G4Pow.hh:56
G4double ComputeTransportCrossSectionPerAtom(G4double CosThetaMax)
int G4int
Definition: G4Types.hh:78
void SetupParticle(const G4ParticleDefinition *)
G4double GetMomentumSquare() const
G4double SetupKinematic(G4double kinEnergy, const G4Material *mat, G4double cut, G4double tmass)
G4double SetupTarget(G4int Z, G4double cut)
bool G4bool
Definition: G4Types.hh:79
G4double GetInvA23() const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double ComputeElectronCrossSection(G4double CosThetaMin, G4double CosThetaMax)
void Initialise(const G4ParticleDefinition *, G4double CosThetaLim)
double G4double
Definition: G4Types.hh:76
G4ThreeVector & SampleSingleScattering(G4double CosThetaMin, G4double CosThetaMax, G4double elecRatio)