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