Geant4_10
G4EnergyLossForExtrapolator.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: G4EnergyLossForExtrapolator.hh 74310 2013-10-03 06:44:36Z gcosmo $
27 //
28 //---------------------------------------------------------------------------
29 //
30 // ClassName: G4EnergyLossForExtrapolator
31 //
32 // Description: This class provide calculation of energy loss, fluctuation,
33 // and msc angle
34 //
35 // Author: 09.12.04 V.Ivanchenko
36 //
37 // Modification:
38 // 08-04-05 Rename Propogator -> Extrapolator
39 // 16-03-06 Add muon tables
40 // 21-03-06 Add verbosity defined in the constructor and Initialisation
41 // start only when first public method is called (V.Ivanchenko)
42 // 03-05-06 Remove unused pointer G4Material* from number of methods (VI)
43 // 28-07-07 Add maxEnergyTransfer for computation of energy loss (VI)
44 //
45 //----------------------------------------------------------------------------
46 //
47 
48 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
49 
50 #ifndef G4EnergyLossForExtrapolator_h
51 #define G4EnergyLossForExtrapolator_h 1
52 
53 #include <vector>
55 
56 #include "globals.hh"
57 #include "G4PhysicsTable.hh"
58 #include "G4DataVector.hh"
59 #include "G4Log.hh"
60 
62 class G4Material;
64 class G4ProductionCuts;
65 
66 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
67 
69 {
70 public:
71 
73 
75 
77 
79 
81 
82  G4double EnergyAfterStep(G4double kinEnergy, G4double step,
83  const G4Material*, const G4ParticleDefinition*);
84 
85  G4double EnergyBeforeStep(G4double kinEnergy, G4double step,
86  const G4Material*, const G4ParticleDefinition*);
87 
88  G4double TrueStepLength(G4double kinEnergy, G4double step,
89  const G4Material*, const G4ParticleDefinition* part);
90 
91  inline G4double EnergyAfterStep(G4double kinEnergy, G4double step,
92  const G4Material*, const G4String& particleName);
93 
94  inline G4double EnergyBeforeStep(G4double kinEnergy, G4double step,
95  const G4Material*, const G4String& particleName);
96 
97  inline G4double AverageScatteringAngle(G4double kinEnergy, G4double step,
98  const G4Material*,
99  const G4ParticleDefinition* part);
100 
101  inline G4double AverageScatteringAngle(G4double kinEnergy, G4double step,
102  const G4Material*,
103  const G4String& particleName);
104 
105  inline G4double ComputeTrueStep(const G4Material*,
106  const G4ParticleDefinition* part,
107  G4double kinEnergy, G4double stepLength);
108 
109  inline G4double EnergyDispersion(G4double kinEnergy, G4double step,
110  const G4Material*, const G4ParticleDefinition*);
111 
112  inline G4double EnergyDispersion(G4double kinEnergy, G4double step,
113  const G4Material*, const G4String& particleName);
114 
115  inline void SetVerbose(G4int val);
116 
117  inline void SetMinKinEnergy(G4double);
118 
119  inline void SetMaxKinEnergy(G4double);
120 
121  inline void SetMaxEnergyTransfer(G4double);
122 
123 private:
124 
125  void Initialisation();
126 
127  G4bool SetupKinematics(const G4ParticleDefinition*, const G4Material*,
128  G4double kinEnergy);
129 
130  G4PhysicsTable* PrepareTable();
131 
132  const G4ParticleDefinition* FindParticle(const G4String& name);
133 
134  void ComputeElectronDEDX(const G4ParticleDefinition* part,
135  G4PhysicsTable* table);
136 
137  void ComputeMuonDEDX(const G4ParticleDefinition* part,
138  G4PhysicsTable* table);
139 
140  void ComputeProtonDEDX(const G4ParticleDefinition* part,
141  G4PhysicsTable* table);
142 
143  void ComputeTrasportXS(const G4ParticleDefinition* part,
144  G4PhysicsTable* table);
145 
146  inline G4double ComputeValue(G4double x, const G4PhysicsTable* table);
147 
148  // hide assignment operator
151 
152  const G4ParticleDefinition* currentParticle;
153  const G4ParticleDefinition* electron;
154  const G4ParticleDefinition* positron;
155  const G4ParticleDefinition* muonPlus;
156  const G4ParticleDefinition* muonMinus;
157  const G4ParticleDefinition* proton;
158 
159  G4DataVector cuts;
160 
161  G4ProductionCuts* pcuts;
162  std::vector<const G4MaterialCutsCouple*> couples;
163 
164  G4String currentParticleName;
165 
166  G4PhysicsTable* dedxElectron;
167  G4PhysicsTable* dedxPositron;
168  G4PhysicsTable* dedxMuon;
169  G4PhysicsTable* dedxProton;
170  G4PhysicsTable* rangeElectron;
171  G4PhysicsTable* rangePositron;
172  G4PhysicsTable* rangeMuon;
173  G4PhysicsTable* rangeProton;
174  G4PhysicsTable* invRangeElectron;
175  G4PhysicsTable* invRangePositron;
176  G4PhysicsTable* invRangeMuon;
177  G4PhysicsTable* invRangeProton;
178  G4PhysicsTable* mscElectron;
179 
180  const G4Material* currentMaterial;
181  G4int index;
182  G4double electronDensity;
183  G4double radLength;
184  G4double mass;
185  G4double charge2;
186  G4double kineticEnergy;
187  G4double gam;
188  G4double bg2;
189  G4double beta2;
190  G4double tmax;
191 
192  G4double linLossLimit;
193  G4double emin;
194  G4double emax;
195  G4double maxEnergyTransfer;
196 
197  G4int nbins;
198  G4int nmat;
199  G4int verbose;
200  G4bool isInitialised;
201 };
202 
203 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
204 
205 inline G4double
207  G4double step,
208  const G4Material* mat,
209  const G4String& name)
210 {
211  return EnergyAfterStep(kinEnergy,step,mat,FindParticle(name));
212 }
213 
214 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
215 
216 inline G4double
218  G4double step,
219  const G4Material* mat,
220  const G4String& name)
221 {
222  return EnergyBeforeStep(kinEnergy,step,mat,FindParticle(name));
223 }
224 
225 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
226 
227 inline G4double
229  G4double step,
230  const G4Material* mat,
231  const G4String& name)
232 {
233  return AverageScatteringAngle(kinEnergy,step,mat,FindParticle(name));
234 }
235 
236 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
237 
238 inline G4double
240  G4double step,
241  const G4Material* mat,
242  const G4String& name)
243 {
244  return EnergyDispersion(kinEnergy,step,mat,FindParticle(name));
245 }
246 
247 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
248 
249 inline G4double
251  G4double stepLength,
252  const G4Material* mat,
253  const G4ParticleDefinition* part)
254 {
255  G4double theta = 0.0;
256  if(SetupKinematics(part, mat, kinEnergy)) {
257  G4double t = stepLength/radLength;
258  G4double y = std::max(0.001, t);
259  theta = 19.23*CLHEP::MeV*std::sqrt(charge2*t)*(1.0 + 0.038*G4Log(y))
260  /(beta2*gam*mass);
261  }
262  return theta;
263 }
264 
265 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
266 
267 inline G4double
269  const G4ParticleDefinition* part,
270  G4double kinEnergy,
271  G4double stepLength)
272 {
273  G4double theta = AverageScatteringAngle(kinEnergy,stepLength,mat,part);
274  return stepLength*std::sqrt(1.0 + 0.625*theta*theta);
275 }
276 
277 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
278 
279 inline G4double
281  G4double stepLength,
282  const G4Material* mat,
283  const G4ParticleDefinition* part)
284 {
285  G4double sig2 = 0.0;
286  if(SetupKinematics(part, mat, kinEnergy)) {
287  G4double step = ComputeTrueStep(mat,part,kinEnergy,stepLength);
288  sig2 = (1.0/beta2 - 0.5)*CLHEP::twopi_mc2_rcl2*tmax*step*electronDensity*charge2;
289  }
290  return sig2;
291 }
292 
293 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
294 
295 inline G4double
296 G4EnergyLossForExtrapolator::ComputeValue(G4double x,
297  const G4PhysicsTable* table)
298 {
299  G4double res = 0.0;
300  G4bool b;
301  if(table) res = ((*table)[index])->GetValue(x, b);
302  return res;
303 }
304 
305 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
306 
308 {
309  verbose = val;
310 }
311 
312 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
313 
315 {
316  emin = val;
317 }
318 
319 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
320 
322 {
323  emax = val;
324 }
325 
326 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
327 
329 {
330  maxEnergyTransfer = val;
331 }
332 
333 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
334 
335 #endif
336 
G4double EnergyDispersion(G4double kinEnergy, G4double step, const G4Material *, const G4ParticleDefinition *)
G4double AverageScatteringAngle(G4double kinEnergy, G4double step, const G4Material *, const G4ParticleDefinition *part)
const XML_Char * name
Definition: expat.h:151
tuple x
Definition: test.py:50
G4double TrueStepLength(G4double kinEnergy, G4double step, const G4Material *, const G4ParticleDefinition *part)
int G4int
Definition: G4Types.hh:78
Float_t mat
Definition: plot.C:40
Double_t y
Definition: plot.C:279
tuple b
Definition: test.py:12
G4double ComputeTrueStep(const G4Material *, const G4ParticleDefinition *part, G4double kinEnergy, G4double stepLength)
const G4ParticleDefinition const G4Material *G4double range
bool G4bool
Definition: G4Types.hh:79
TString part[npart]
Definition: Style.C:32
G4double ComputeDEDX(G4double kinEnergy, const G4ParticleDefinition *)
G4double ComputeRange(G4double kinEnergy, const G4ParticleDefinition *)
G4double EnergyBeforeStep(G4double kinEnergy, G4double step, const G4Material *, const G4ParticleDefinition *)
G4double EnergyAfterStep(G4double kinEnergy, G4double step, const G4Material *, const G4ParticleDefinition *)
G4double G4Log(G4double x)
Definition: G4Log.hh:227
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double ComputeEnergy(G4double range, const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76