Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4VMultipleScattering.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 // $Id: G4VMultipleScattering.hh 106714 2017-10-20 09:38:06Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class header file
31 //
32 //
33 // File name: G4VMultipleScattering
34 //
35 // Author: Vladimir Ivanchenko on base of Laszlo Urban code
36 //
37 // Creation date: 12.03.2002
38 //
39 // Modifications:
40 //
41 // 16-07-03 Update GetRange interface (V.Ivanchenko)
42 //
43 //
44 // Class Description:
45 //
46 // It is the generic process of multiple scattering it includes common
47 // part of calculations for all charged particles
48 //
49 // 26-11-03 bugfix in AlongStepDoIt (L.Urban)
50 // 25-05-04 add protection against case when range is less than steplimit (VI)
51 // 30-06-04 make destructor virtual (V.Ivanchenko)
52 // 27-08-04 Add InitialiseForRun method (V.Ivanchneko)
53 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivanchenko)
54 // 15-04-05 optimize internal interfaces (V.Ivanchenko)
55 // 15-04-05 remove boundary flag (V.Ivanchenko)
56 // 07-10-05 error in a protection in GetContinuousStepLimit corrected (L.Urban)
57 // 27-10-05 introduce virtual function MscStepLimitation() (V.Ivanchenko)
58 // 26-01-06 Rename GetRange -> GetRangeFromRestricteDEDX (V.Ivanchenko)
59 // 17-02-06 Save table of transport cross sections not mfp (V.Ivanchenko)
60 // 07-03-06 Move step limit calculation to model (V.Ivanchenko)
61 // 13-05-06 Add method to access model by index (V.Ivanchenko)
62 // 12-02-07 Add get/set skin (V.Ivanchenko)
63 // 27-10-07 Virtual functions moved to source (V.Ivanchenko)
64 // 15-07-08 Reorder class members for further multi-thread development (VI)
65 // 07-04-09 Moved msc methods from G4VEmModel to G4VMscModel (VI)
66 //
67 
68 // -------------------------------------------------------------------
69 //
70 
71 #ifndef G4VMultipleScattering_h
72 #define G4VMultipleScattering_h 1
73 
75 #include "globals.hh"
76 #include "G4Material.hh"
78 #include "G4Track.hh"
79 #include "G4Step.hh"
80 #include "G4EmModelManager.hh"
81 #include "G4VMscModel.hh"
82 #include "G4EmParameters.hh"
83 #include "G4MscStepLimitType.hh"
84 
87 class G4LossTableManager;
88 class G4SafetyHelper;
89 
90 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
91 
93 {
94 public:
95 
96  G4VMultipleScattering(const G4String& name = "msc",
98 
99  virtual ~G4VMultipleScattering();
100 
101  //------------------------------------------------------------------------
102  // Virtual methods to be implemented for the concrete model
103  //------------------------------------------------------------------------
104 
105  virtual G4bool IsApplicable(const G4ParticleDefinition& p) override = 0;
106 
107  // obsolete
108  virtual void PrintInfo() {};
109 
110  virtual void ProcessDescription(std::ostream& outFile) const override;
111 
112 protected:
113 
114  virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
115 
116  virtual void StreamProcessInfo(std::ostream&, G4String) const {};
117 
118 public:
119 
120  //------------------------------------------------------------------------
121  // Generic methods common to all ContinuousDiscrete processes
122  //------------------------------------------------------------------------
123 
124  // Initialise for build of tables
125  void PreparePhysicsTable(const G4ParticleDefinition&) override;
126 
127  // Build physics table during initialisation
128  void BuildPhysicsTable(const G4ParticleDefinition&) override;
129 
130  // Store PhysicsTable in a file.
131  // Return false in case of failure at I/O
133  const G4String& directory,
134  G4bool ascii = false) override;
135 
136  // Retrieve Physics from a file.
137  // (return true if the Physics Table can be build by using file)
138  // (return false if the process has no functionality or in case of failure)
139  // File name should is constructed as processName+particleName and the
140  // should be placed under the directory specifed by the argument.
142  const G4String& directory,
143  G4bool ascii) override;
144 
145  // This is called in the beginning of tracking for a new track
146  void StartTracking(G4Track*) override;
147 
148  // The function overloads the corresponding function of the base
149  // class.It limits the step near to boundaries only
150  // and invokes the method GetMscContinuousStepLimit at every step.
152  const G4Track&,
153  G4double previousStepSize,
154  G4double currentMinimalStep,
155  G4double& currentSafety,
156  G4GPILSelection* selection) override;
157 
158  // The function overloads the corresponding function of the base
159  // class.
161  const G4Track&,
162  G4double previousStepSize,
163  G4ForceCondition* condition) override;
164 
165  // Along step actions
166  G4VParticleChange* AlongStepDoIt(const G4Track&, const G4Step&) override;
167 
168  // Post step actions
169  G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&) override;
170 
171  // This method does not used for tracking, it is intended only for tests
173  G4double previousStepSize,
174  G4double currentMinimalStep,
175  G4double& currentSafety);
176 
177  //------------------------------------------------------------------------
178  // Specific methods to set, access, modify models
179  //------------------------------------------------------------------------
180 
181  // Select model in run time
182  inline G4VEmModel* SelectModel(G4double kinEnergy, size_t idx);
183 
184 public:
185 
186  // Add model for region, smaller value of order defines which
187  // model will be selected for a given energy interval
188  void AddEmModel(G4int order, G4VEmModel*, const G4Region* region = nullptr);
189 
190  // Assign a model to a process local list, to enable the list in run time
191  // the derived process should execute AddEmModel(..) for all such models
192  void SetEmModel(G4VMscModel*, size_t index = 0);
193 
194  // return a model from the local list
195  G4VMscModel* EmModel(size_t index = 0) const;
196 
197  // Access to run time models by index
198  inline G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false) const;
199 
200  //------------------------------------------------------------------------
201  // Get/Set parameters for simulation of multiple scattering
202  //------------------------------------------------------------------------
203 
205 
206  inline G4bool LateralDisplasmentFlag() const;
207  inline void SetLateralDisplasmentFlag(G4bool val);
208 
209  inline G4double Skin() const;
210  inline void SetSkin(G4double val);
211 
212  inline G4double RangeFactor() const;
213  inline void SetRangeFactor(G4double val);
214 
215  inline G4double GeomFactor() const;
216 
217  inline G4double PolarAngleLimit() const;
218 
219  inline G4MscStepLimitType StepLimitType() const;
220  inline void SetStepLimitType(G4MscStepLimitType val);
221 
222  inline G4double LowestKinEnergy() const;
223  inline void SetLowestKinEnergy(G4double val);
224 
225  inline const G4ParticleDefinition* FirstParticle() const;
226 
227  //------------------------------------------------------------------------
228  // Run time methods
229  //------------------------------------------------------------------------
230 
231 protected:
232 
233  // This method is not used for tracking, it returns mean free path value
235  G4double,
236  G4ForceCondition* condition) override;
237 
238  // This method is not used for tracking, it returns step limit
240  G4double previousStepSize,
241  G4double currentMinimalStep,
242  G4double& currentSafety) override ;
243 
244  // return number of already added models
245  inline G4int NumberOfModels() const;
246 
247 private:
248 
249  // hide assignment operator
252  operator=(const G4VMultipleScattering &right) = delete;
253 
254  // Print out of generic class parameters
255  void StreamInfo(std::ostream& outFile, const G4ParticleDefinition&,
256  G4String endOfLine=G4String("\n")) const;
257 
258  // ======== Parameters of the class fixed at construction =========
259 
263 
264  // ======== Parameters of the class fixed at initialisation =======
265 
267 
268  std::vector<G4VMscModel*> mscModels;
270 
273 
275 
278 
282 
283  // ======== Cached values - may be state dependent ================
284 
285 protected:
286 
288 
289 private:
290 
291  // cache
294 
300 
305 };
306 
307 // ======== Run time inline methods ================
308 
309 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
310 
311 inline G4VEmModel*
312 G4VMultipleScattering::SelectModel(G4double kinEnergy, size_t coupleIndex)
313 {
314  return modelManager->SelectModel(kinEnergy, coupleIndex);
315 }
316 
317 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
318 
320 {
321  return latDisplacement;
322 }
323 
324 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
325 
327 {
328  latDisplacement = val;
329 }
330 
331 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
332 
334 {
335  return theParameters->MscSkin();
336 }
337 
338 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
339 
341 {
343 }
344 
345 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
346 
348 {
349  return facrange;
350 }
351 
352 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
353 
355 {
356  if(val > 0.0 && val < 1.0) {
357  facrange = val;
359  }
360 }
361 
362 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
363 
365 {
366  return theParameters->MscGeomFactor();
367 }
368 
369 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
370 
372 {
373  return theParameters->MscThetaLimit();
374 }
375 
376 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
377 
379 {
380  return stepLimit;
381 }
382 
383 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
384 
386 {
387  stepLimit = val;
388  if(val == fMinimal) { SetRangeFactor(0.2); }
390 }
391 
392 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
393 
395 {
396  return lowestKinEnergy;
397 }
398 
399 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
400 
402 {
403  lowestKinEnergy = val;
404 }
405 
406 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
407 
409 {
410  return firstParticle;
411 }
412 
413 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
414 
416 {
417  return modelManager->NumberOfModels();
418 }
419 
420 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
421 
422 inline G4VEmModel*
424 {
425  return modelManager->GetModel(idx, ver);
426 }
427 
428 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
429 
430 #endif
void SetRangeFactor(G4double val)
const XML_Char * name
Definition: expat.h:151
G4VMultipleScattering(const G4String &name="msc", G4ProcessType type=fElectromagnetic)
G4ParticleChangeForMSC fParticleChange
G4MscStepLimitType
G4bool StorePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii=false) override
G4MscStepLimitType StepLimitType() const
G4double LowestKinEnergy() const
G4bool RetrievePhysicsTable(const G4ParticleDefinition *, const G4String &directory, G4bool ascii) override
const char * p
Definition: xmltok.h:285
G4MscStepLimitType stepLimit
G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &) override
G4double condition(const G4ErrorSymMatrix &m)
void BuildPhysicsTable(const G4ParticleDefinition &) override
void SetLowestKinEnergy(G4double val)
G4VEmModel * GetModel(G4int idx, G4bool ver=false)
G4VMultipleScattering & operator=(const G4VMultipleScattering &right)=delete
const G4ParticleDefinition * currParticle
G4double ContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety)
const G4ParticleDefinition * firstParticle
G4bool LateralDisplasmentFlag() const
G4int NumberOfModels() const
void SetMscSkin(G4double val)
G4double PolarAngleLimit() const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
void SetMscRangeFactor(G4double val)
G4ProcessType
G4double MscSkin() const
G4GPILSelection
G4VEmModel * SelectModel(G4double kinEnergy, size_t idx)
G4VEmModel * SelectModel(G4double &energy, size_t &index)
G4double GetContinuousStepLimit(const G4Track &track, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety) override
G4VMscModel * EmModel(size_t index=0) const
G4EmModelManager * modelManager
G4VEmModel * GetModelByIndex(G4int idx=0, G4bool ver=false) const
void StartTracking(G4Track *) override
G4LossTableManager * emManager
void StreamInfo(std::ostream &outFile, const G4ParticleDefinition &, G4String endOfLine=G4String("\n")) const
Definition: G4Step.hh:76
void SetIonisation(G4VEnergyLossProcess *)
G4VParticleChange * AlongStepDoIt(const G4Track &, const G4Step &) override
G4double GetMeanFreePath(const G4Track &track, G4double, G4ForceCondition *condition) override
G4double PostStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4ForceCondition *condition) override
void PreparePhysicsTable(const G4ParticleDefinition &) override
void SetMscStepLimitType(G4MscStepLimitType val)
virtual G4bool IsApplicable(const G4ParticleDefinition &p) override=0
int G4int
Definition: G4Types.hh:78
G4double MscThetaLimit() const
void SetEmModel(G4VMscModel *, size_t index=0)
G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double previousStepSize, G4double currentMinimalStep, G4double &currentSafety, G4GPILSelection *selection) override
std::vector< G4VMscModel * > mscModels
virtual void ProcessDescription(std::ostream &outFile) const override
virtual void StreamProcessInfo(std::ostream &, G4String) const
void SetStepLimitType(G4MscStepLimitType val)
void SetLateralDisplasmentFlag(G4bool val)
G4ForceCondition
G4double MscGeomFactor() const
virtual void InitialiseProcess(const G4ParticleDefinition *)=0
void AddEmModel(G4int order, G4VEmModel *, const G4Region *region=nullptr)
const G4ParticleDefinition * FirstParticle() const
G4VEnergyLossProcess * fIonisation