| Geant4 Cross Reference |
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 // GEANT4 Class header file
29 //
30 //
31 // File name: G4VMscModel
32 //
33 // Author: Vladimir Ivanchenko
34 //
35 // Creation date: 07.03.2008
36 //
37 // Modifications:
38 // 07.04.2009 V.Ivanchenko moved msc methods from G4VEmModel to G4VMscModel
39 // 26.03.2012 V.Ivanchenko added transport x-section pointer and Get?Set methods
40 //
41 // Class Description:
42 //
43 // General interface to msc models
44
45 // -------------------------------------------------------------------
46 //
47 #ifndef G4VMscModel_h
48 #define G4VMscModel_h 1
49
50 #include <CLHEP/Units/SystemOfUnits.h>
51
52 #include "G4VEmModel.hh"
53 #include "G4MscStepLimitType.hh"
54 #include "globals.hh"
55 #include "G4ThreeVector.hh"
56 #include "G4Track.hh"
57 #include "G4SafetyHelper.hh"
58 #include "G4VEnergyLossProcess.hh"
59 #include "G4PhysicsTable.hh"
60 #include "G4ThreeVector.hh"
61 #include <vector>
62
63 class G4ParticleChangeForMSC;
64 class G4ParticleDefinition;
65
66 class G4VMscModel : public G4VEmModel
67 {
68
69 public:
70
71 explicit G4VMscModel(const G4String& nam);
72
73 ~G4VMscModel() override;
74
75 virtual G4double ComputeTruePathLengthLimit(const G4Track& track,
76 G4double& stepLimit) = 0;
77
78 virtual G4double ComputeGeomPathLength(G4double truePathLength) = 0;
79
80 virtual G4double ComputeTrueStepLength(G4double geomPathLength) = 0;
81
82 virtual G4ThreeVector& SampleScattering(const G4ThreeVector&,
83 G4double safety) = 0;
84
85 void InitialiseParameters(const G4ParticleDefinition*);
86
87 void DumpParameters(std::ostream& out) const;
88
89 // empty method
90 void SampleSecondaries(std::vector<G4DynamicParticle*>*,
91 const G4MaterialCutsCouple*,
92 const G4DynamicParticle*,
93 G4double tmin, G4double tmax) override;
94
95 //================================================================
96 // Set parameters of multiple scattering models
97 //================================================================
98
99 inline void SetStepLimitType(G4MscStepLimitType);
100
101 inline void SetLateralDisplasmentFlag(G4bool val);
102
103 inline void SetRangeFactor(G4double);
104
105 inline void SetGeomFactor(G4double);
106
107 inline void SetSkin(G4double);
108
109 inline void SetLambdaLimit(G4double);
110
111 inline void SetSafetyFactor(G4double);
112
113 inline void SetSampleZ(G4bool);
114
115 //================================================================
116 // Get/Set access to Physics Tables
117 //================================================================
118
119 inline G4VEnergyLossProcess* GetIonisation() const;
120
121 inline void SetIonisation(G4VEnergyLossProcess*,
122 const G4ParticleDefinition* part);
123
124 //================================================================
125 // Run time methods
126 //================================================================
127
128 protected:
129
130 // initialisation of the ParticleChange for the model
131 // initialisation of interface with geometry and ionisation
132 G4ParticleChangeForMSC*
133 GetParticleChangeForMSC(const G4ParticleDefinition* p = nullptr);
134
135 // convert true length to geometry length
136 inline G4double ConvertTrueToGeom(G4double& tLength, G4double& gLength);
137
138 public:
139
140 // compute safety
141 inline G4double ComputeSafety(const G4ThreeVector& position,
142 G4double limit= DBL_MAX);
143
144 // compute linear distance to a geometry boundary
145 inline G4double ComputeGeomLimit(const G4Track&, G4double& presafety,
146 G4double limit);
147
148 inline G4double GetDEDX(const G4ParticleDefinition* part,
149 G4double kineticEnergy,
150 const G4MaterialCutsCouple* couple);
151 inline G4double GetDEDX(const G4ParticleDefinition* part,
152 G4double kineticEnergy,
153 const G4MaterialCutsCouple* couple,
154 G4double logKineticEnergy);
155
156 inline G4double GetRange(const G4ParticleDefinition* part,
157 G4double kineticEnergy,
158 const G4MaterialCutsCouple* couple);
159 inline G4double GetRange(const G4ParticleDefinition* part,
160 G4double kineticEnergy,
161 const G4MaterialCutsCouple* couple,
162 G4double logKineticEnergy);
163
164 inline G4double GetEnergy(const G4ParticleDefinition* part,
165 G4double range,
166 const G4MaterialCutsCouple* couple);
167
168 // G4MaterialCutsCouple should be defined before call to this method
169 inline
170 G4double GetTransportMeanFreePath(const G4ParticleDefinition* part,
171 G4double kinEnergy);
172 inline
173 G4double GetTransportMeanFreePath(const G4ParticleDefinition* part,
174 G4double kinEnergy,
175 G4double logKinEnergy);
176
177 // hide assignment operator
178 G4VMscModel & operator=(const G4VMscModel &right) = delete;
179 G4VMscModel(const G4VMscModel&) = delete;
180
181 private:
182
183 G4SafetyHelper* safetyHelper;
184 G4VEnergyLossProcess* ionisation;
185 const G4ParticleDefinition* currentPart;
186
187 G4double dedx;
188 G4double localtkin;
189 G4double localrange;
190
191 protected:
192
193 G4double facrange;
194 G4double facgeom;
195 G4double facsafety;
196 G4double skin;
197 G4double dtrl;
198 G4double lambdalimit;
199 G4double geomMin;
200 G4double geomMax;
201
202 G4ThreeVector fDisplacement;
203 G4MscStepLimitType steppingAlgorithm;
204
205 G4bool samplez;
206 G4bool latDisplasment;
207
208 };
209
210 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
211 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
212
213 inline void G4VMscModel::SetLateralDisplasmentFlag(G4bool val)
214 {
215 if(!IsLocked()) { latDisplasment = val; }
216 }
217
218 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
219
220 inline void G4VMscModel::SetSkin(G4double val)
221 {
222 if(!IsLocked()) { skin = val; }
223 }
224
225 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
226
227 inline void G4VMscModel::SetRangeFactor(G4double val)
228 {
229 if(!IsLocked()) { facrange = val; }
230 }
231
232 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
233
234 inline void G4VMscModel::SetGeomFactor(G4double val)
235 {
236 if(!IsLocked()) { facgeom = val; }
237 }
238
239 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
240
241 inline void G4VMscModel::SetLambdaLimit(G4double val)
242 {
243 if(!IsLocked()) { lambdalimit = val; }
244 }
245
246 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
247
248 inline void G4VMscModel::SetSafetyFactor(G4double val)
249 {
250 if(!IsLocked()) { facsafety = val; }
251 }
252
253 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
254
255 inline void G4VMscModel::SetStepLimitType(G4MscStepLimitType val)
256 {
257 if(!IsLocked()) { steppingAlgorithm = val; }
258 }
259
260 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
261
262 inline void G4VMscModel::SetSampleZ(G4bool val)
263 {
264 if(!IsLocked()) { samplez = val; }
265 }
266
267 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
268
269 inline G4double G4VMscModel::ComputeSafety(const G4ThreeVector& position,
270 G4double limit)
271 {
272 return safetyHelper->ComputeSafety(position, limit);
273 }
274
275 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
276
277 inline G4double G4VMscModel::ConvertTrueToGeom(G4double& tlength,
278 G4double& glength)
279 {
280 glength = ComputeGeomPathLength(tlength);
281 // should return true length
282 return tlength;
283 }
284
285 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
286
287 inline G4double G4VMscModel::ComputeGeomLimit(const G4Track& track,
288 G4double& presafety,
289 G4double limit)
290 {
291 return safetyHelper->CheckNextStep(
292 track.GetStep()->GetPreStepPoint()->GetPosition(),
293 track.GetMomentumDirection(),
294 limit, presafety);
295 }
296
297 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
298
299 inline G4double
300 G4VMscModel::GetDEDX(const G4ParticleDefinition* part, G4double kinEnergy,
301 const G4MaterialCutsCouple* couple)
302 {
303 G4double x;
304 if (ionisation) {
305 x = ionisation->GetDEDX(kinEnergy, couple);
306 } else {
307 const G4double q = part->GetPDGCharge()*inveplus;
308 x = dedx*q*q;
309 }
310 return x;
311 }
312
313 inline G4double
314 G4VMscModel::GetDEDX(const G4ParticleDefinition* part, G4double kinEnergy,
315 const G4MaterialCutsCouple* couple, G4double logKinEnergy)
316 {
317 G4double x;
318 if (ionisation) {
319 x = ionisation->GetDEDX(kinEnergy, couple, logKinEnergy);
320 } else {
321 const G4double q = part->GetPDGCharge()*inveplus;
322 x = dedx*q*q;
323 }
324 return x;
325 }
326
327 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
328
329 inline G4double
330 G4VMscModel::GetRange(const G4ParticleDefinition* part,G4double kinEnergy,
331 const G4MaterialCutsCouple* couple)
332 {
333 //G4cout << "G4VMscModel::GetRange E(MeV)= " << kinEnergy << " "
334 // << ionisation << " " << part->GetParticleName()
335 // << G4endl;
336 localtkin = kinEnergy;
337 if (ionisation) {
338 localrange = ionisation->GetRangeForLoss(kinEnergy, couple);
339 } else {
340 const G4double q = part->GetPDGCharge()*inveplus;
341 localrange = kinEnergy/(dedx*q*q*couple->GetMaterial()->GetDensity());
342 }
343 //G4cout << "R(mm)= " << localrange << " " << ionisation << G4endl;
344 return localrange;
345 }
346
347 inline G4double
348 G4VMscModel::GetRange(const G4ParticleDefinition* part,G4double kinEnergy,
349 const G4MaterialCutsCouple* couple, G4double logKinEnergy)
350 {
351 //G4cout << "G4VMscModel::GetRange E(MeV)= " << kinEnergy << " "
352 // << ionisation << " " << part->GetParticleName()
353 // << G4endl;
354 localtkin = kinEnergy;
355 if (ionisation) {
356 localrange = ionisation->GetRangeForLoss(kinEnergy, couple, logKinEnergy);
357 } else {
358 const G4double q = part->GetPDGCharge()*inveplus;
359 localrange = kinEnergy/(dedx*q*q*couple->GetMaterial()->GetDensity());
360 }
361 //G4cout << "R(mm)= " << localrange << " " << ionisation << G4endl;
362 return localrange;
363 }
364
365 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
366
367 inline G4double
368 G4VMscModel::GetEnergy(const G4ParticleDefinition* part,
369 G4double range, const G4MaterialCutsCouple* couple)
370 {
371 G4double e;
372 //G4cout << "G4VMscModel::GetEnergy R(mm)= " << range << " " << ionisation
373 // << " Rlocal(mm)= " << localrange << " Elocal(MeV)= " << localtkin
374 // << G4endl;
375 if(ionisation) { e = ionisation->GetKineticEnergy(range, couple); }
376 else {
377 e = localtkin;
378 if(localrange > range) {
379 G4double q = part->GetPDGCharge()*inveplus;
380 e -= (localrange - range)*dedx*q*q*couple->GetMaterial()->GetDensity();
381 }
382 }
383 return e;
384 }
385
386 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
387
388 inline G4VEnergyLossProcess* G4VMscModel::GetIonisation() const
389 {
390 return ionisation;
391 }
392
393 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
394
395 inline void G4VMscModel::SetIonisation(G4VEnergyLossProcess* p,
396 const G4ParticleDefinition* part)
397 {
398 ionisation = p;
399 currentPart = part;
400 }
401
402 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
403
404 inline G4double
405 G4VMscModel::GetTransportMeanFreePath(const G4ParticleDefinition* part,
406 G4double ekin)
407 {
408 G4double x;
409 if (xSectionTable) {
410 const G4int idx = CurrentCouple()->GetIndex();
411 x = (*xSectionTable)[idx]->Value(ekin, idxTable)/(ekin*ekin);
412 } else {
413 x = CrossSectionPerVolume(CurrentCouple()->GetMaterial(), part, ekin,
414 0.0, DBL_MAX);
415 }
416 return (x > 0.0) ? 1.0/x : DBL_MAX;
417 }
418
419 inline G4double
420 G4VMscModel::GetTransportMeanFreePath(const G4ParticleDefinition* part,
421 G4double ekin, G4double logekin)
422 {
423 G4double x;
424 if (xSectionTable) {
425 const G4int idx = CurrentCouple()->GetIndex();
426 x = (*xSectionTable)[idx]->LogVectorValue(ekin, logekin)/(ekin*ekin);
427 } else {
428 x = CrossSectionPerVolume(CurrentCouple()->GetMaterial(), part, ekin,
429 0.0, DBL_MAX);
430 }
431 return (x > 0.0) ? 1.0/x : DBL_MAX;
432 }
433
434 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
435
436 #endif
437
438