Geant4  9.6.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VEnergyLossProcess.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$
27 // GEANT4 tag $Name:
28 //
29 // -------------------------------------------------------------------
30 //
31 // GEANT4 Class header file
32 //
33 //
34 // File name: G4VEnergyLossProcess
35 //
36 // Author: Vladimir Ivanchenko on base of Laszlo Urban code
37 //
38 // Creation date: 03.01.2002
39 //
40 // Modifications:
41 //
42 // 26-12-02 Secondary production moved to derived classes (V.Ivanchenko)
43 // 20-01-03 Migrade to cut per region (V.Ivanchenko)
44 // 24-01-03 Make models region aware (V.Ivanchenko)
45 // 05-02-03 Fix compilation warnings (V.Ivanchenko)
46 // 13-02-03 SubCutoffProcessors defined for regions (V.Ivanchenko)
47 // 17-02-03 Fix problem of store/restore tables (V.Ivanchenko)
48 // 26-02-03 Region dependent step limit (V.Ivanchenko)
49 // 26-03-03 Add GetDEDXDispersion (V.Ivanchenko)
50 // 09-04-03 Fix problem of negative range limit for non integral (V.Ivanchenko)
51 // 13-05-03 Add calculation of precise range (V.Ivanchenko)
52 // 21-07-03 Add UpdateEmModel method (V.Ivanchenko)
53 // 12-11-03 G4EnergyLossSTD -> G4EnergyLossProcess (V.Ivanchenko)
54 // 14-01-04 Activate precise range calculation (V.Ivanchenko)
55 // 10-03-04 Fix problem of step limit calculation (V.Ivanchenko)
56 // 30-06-04 make destructor virtual (V.Ivanchenko)
57 // 05-07-04 fix problem of GenericIons seen at small cuts (V.Ivanchenko)
58 // 03-08-04 Add DEDX table to all processes for control on integral range(VI)
59 // 06-08-04 Clear up names of member functions (V.Ivanchenko)
60 // 27-08-04 Add NeedBuildTables method (V.Ivanchneko)
61 // 09-09-04 Bug fix for the integral mode with 2 peaks (V.Ivanchneko)
62 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivanchenko)
63 // 08-04-05 Major optimisation of internal interfaces (V.Ivanchenko)
64 // 11-04-05 Use MaxSecondaryEnergy from a model (V.Ivanchenko)
65 // 10-01-05 Remove SetStepLimits (V.Ivanchenko)
66 // 10-01-06 PreciseRange -> CSDARange (V.Ivantchenko)
67 // 13-01-06 Remove AddSubCutSecondaries and cleanup (V.Ivantchenko)
68 // 20-01-06 Introduce G4EmTableType and reducing number of methods (VI)
69 // 26-01-06 Add public method GetCSDARange (V.Ivanchenko)
70 // 22-03-06 Add SetDynamicMassCharge (V.Ivanchenko)
71 // 23-03-06 Use isIonisation flag (V.Ivanchenko)
72 // 13-05-06 Add method to access model by index (V.Ivanchenko)
73 // 14-01-07 add SetEmModel(index) and SetFluctModel() (mma)
74 // 15-01-07 Add separate ionisation tables and reorganise get/set methods for
75 // dedx tables (V.Ivanchenko)
76 // 13-03-07 use SafetyHelper instead of navigator (V.Ivanchenko)
77 // 27-07-07 use stl vector for emModels instead of C-array (V.Ivanchenko)
78 // 25-09-07 More accurate handling zero xsect in
79 // PostStepGetPhysicalInteractionLength (V.Ivanchenko)
80 // 27-10-07 Virtual functions moved to source (V.Ivanchenko)
81 // 15-07-08 Reorder class members for further multi-thread development (VI)
82 //
83 // Class Description:
84 //
85 // It is the unified energy loss process it calculates the continuous
86 // energy loss for charged particles using a set of Energy Loss
87 // models valid for different energy regions. There are a possibility
88 // to create and access to dE/dx and range tables, or to calculate
89 // that information on fly.
90 
91 // -------------------------------------------------------------------
92 //
93 
94 #ifndef G4VEnergyLossProcess_h
95 #define G4VEnergyLossProcess_h 1
96 
97 #include "G4VContinuousDiscreteProcess.hh"
98 #include "globals.hh"
99 #include "G4Material.hh"
100 #include "G4MaterialCutsCouple.hh"
101 #include "G4Track.hh"
102 #include "G4EmModelManager.hh"
103 #include "G4UnitsTable.hh"
104 #include "G4ParticleChangeForLoss.hh"
105 #include "G4EmTableType.hh"
106 #include "G4PhysicsTable.hh"
107 #include "G4PhysicsVector.hh"
108 
109 class G4Step;
110 class G4ParticleDefinition;
111 class G4VEmModel;
112 class G4VEmFluctuationModel;
113 class G4DataVector;
114 class G4Region;
115 class G4SafetyHelper;
116 class G4VAtomDeexcitation;
117 class G4EmBiasingManager;
118 
119 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
120 
121 class G4VEnergyLossProcess : public G4VContinuousDiscreteProcess
122 {
123 public:
124 
125  G4VEnergyLossProcess(const G4String& name = "EnergyLoss",
126  G4ProcessType type = fElectromagnetic);
127 
128  virtual ~G4VEnergyLossProcess();
129 
130 private:
131  // clean vectors and arrays
132  void Clean();
133 
134  //------------------------------------------------------------------------
135  // Virtual methods to be implemented in concrete processes
136  //------------------------------------------------------------------------
137 
138 public:
139  virtual G4bool IsApplicable(const G4ParticleDefinition& p) = 0;
140 
141  virtual void PrintInfo() = 0;
142 
143 protected:
144 
145  virtual void InitialiseEnergyLossProcess(const G4ParticleDefinition*,
146  const G4ParticleDefinition*) = 0;
147 
148  //------------------------------------------------------------------------
149  // Methods with standard implementation; may be overwritten if needed
150  //------------------------------------------------------------------------
151 
152  virtual G4double MinPrimaryEnergy(const G4ParticleDefinition*,
153  const G4Material*, G4double cut);
154 
155  //------------------------------------------------------------------------
156  // Virtual methods implementation common to all EM ContinuousDiscrete
157  // processes. Further inheritance is not assumed
158  //------------------------------------------------------------------------
159 
160 public:
161 
162  // prepare all tables
163  void PreparePhysicsTable(const G4ParticleDefinition&);
164 
165  // build all tables
166  void BuildPhysicsTable(const G4ParticleDefinition&);
167 
168  // build a table
169  G4PhysicsTable* BuildDEDXTable(G4EmTableType tType = fRestricted);
170 
171  // build a table
172  G4PhysicsTable* BuildLambdaTable(G4EmTableType tType = fRestricted);
173 
174  // summary printout after initialisation
175  void PrintInfoDefinition();
176 
177  // Called before tracking of each new G4Track
178  void StartTracking(G4Track*);
179 
180  // Step limit from AlongStep
181  G4double AlongStepGetPhysicalInteractionLength(const G4Track&,
182  G4double previousStepSize,
183  G4double currentMinimumStep,
184  G4double& currentSafety,
185  G4GPILSelection* selection);
186 
187  // Step limit from cross section
188  G4double PostStepGetPhysicalInteractionLength(const G4Track& track,
189  G4double previousStepSize,
190  G4ForceCondition* condition);
191 
192  // AlongStep computations
193  G4VParticleChange* AlongStepDoIt(const G4Track&, const G4Step&);
194 
195  // Sampling of secondaries in vicinity of geometrical boundary
196  // Return sum of secodaries energy
197  G4double SampleSubCutSecondaries(std::vector<G4Track*>&, const G4Step&,
198  G4VEmModel* model, G4int matIdx);
199 
200  // PostStep sampling of secondaries
201  G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&);
202 
203  // Store all PhysicsTable in files.
204  // Return false in case of any fatal failure at I/O
205  G4bool StorePhysicsTable(const G4ParticleDefinition*,
206  const G4String& directory,
207  G4bool ascii = false);
208 
209  // Retrieve all Physics from a files.
210  // Return true if all the Physics Table are built.
211  // Return false if any fatal failure.
212  G4bool RetrievePhysicsTable(const G4ParticleDefinition*,
213  const G4String& directory,
214  G4bool ascii);
215 
216 private:
217  // store a table
218  G4bool StoreTable(const G4ParticleDefinition* p,
219  G4PhysicsTable*, G4bool ascii,
220  const G4String& directory,
221  const G4String& tname);
222 
223  // retrieve a table
224  G4bool RetrieveTable(const G4ParticleDefinition* p,
225  G4PhysicsTable*, G4bool ascii,
226  const G4String& directory,
227  const G4String& tname,
228  G4bool mandatory);
229 
230  //------------------------------------------------------------------------
231  // Public interface to cross section, mfp and sampling of fluctuations
232  // These methods are not used in run time
233  //------------------------------------------------------------------------
234 
235 public:
236  // access to dispersion of restricted energy loss
237  G4double GetDEDXDispersion(const G4MaterialCutsCouple *couple,
238  const G4DynamicParticle* dp,
239  G4double length);
240 
241  // Access to cross section table
242  G4double CrossSectionPerVolume(G4double kineticEnergy,
243  const G4MaterialCutsCouple* couple);
244 
245  // access to cross section
246  G4double MeanFreePath(const G4Track& track);
247 
248  // access to step limit
249  G4double ContinuousStepLimit(const G4Track& track,
250  G4double previousStepSize,
251  G4double currentMinimumStep,
252  G4double& currentSafety);
253 
254 protected:
255 
256  // implementation of the pure virtual method
257  G4double GetMeanFreePath(const G4Track& track,
258  G4double previousStepSize,
259  G4ForceCondition* condition);
260 
261  // implementation of the pure virtual method
262  G4double GetContinuousStepLimit(const G4Track& track,
263  G4double previousStepSize,
264  G4double currentMinimumStep,
265  G4double& currentSafety);
266 
267  //------------------------------------------------------------------------
268  // Run time method which may be also used by derived processes
269  //------------------------------------------------------------------------
270 
271  // creeation of an empty vector for cross section
272  G4PhysicsVector* LambdaPhysicsVector(const G4MaterialCutsCouple*,
273  G4double cut);
274 
275  inline size_t CurrentMaterialCutsCoupleIndex() const;
276 
277  inline G4double GetCurrentRange() const;
278 
279  //------------------------------------------------------------------------
280  // Specific methods to set, access, modify models
281  //------------------------------------------------------------------------
282 
283  // Select model in run time
284  inline void SelectModel(G4double kinEnergy);
285 
286 public:
287  // Select model by energy and region index
288  inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy,
289  size_t& idx) const;
290 
291  // Add EM model coupled with fluctuation model for region, smaller value
292  // of order defines which pair of models will be selected for a given
293  // energy interval
294  void AddEmModel(G4int, G4VEmModel*,
295  G4VEmFluctuationModel* fluc = 0,
296  const G4Region* region = 0);
297 
298  // Define new energy range for the model identified by the name
299  void UpdateEmModel(const G4String&, G4double, G4double);
300 
301  // Assign a model to a process
302  void SetEmModel(G4VEmModel*, G4int index=1);
303 
304  // return the assigned model
305  G4VEmModel* EmModel(G4int index=1);
306 
307  // Access to models
308  G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false);
309 
310  G4int NumberOfModels();
311 
312  // Assign a fluctuation model to a process
313  void SetFluctModel(G4VEmFluctuationModel*);
314 
315  // return the assigned fluctuation model
316  inline G4VEmFluctuationModel* FluctModel();
317 
318  //------------------------------------------------------------------------
319  // Define and access particle type
320  //------------------------------------------------------------------------
321 
322 protected:
323  inline void SetParticle(const G4ParticleDefinition* p);
324  inline void SetSecondaryParticle(const G4ParticleDefinition* p);
325 
326 public:
327  inline void SetBaseParticle(const G4ParticleDefinition* p);
328  inline const G4ParticleDefinition* Particle() const;
329  inline const G4ParticleDefinition* BaseParticle() const;
330  inline const G4ParticleDefinition* SecondaryParticle() const;
331 
332  //------------------------------------------------------------------------
333  // Get/set parameters to configure the process at initialisation time
334  //------------------------------------------------------------------------
335 
336  // Add subcutoff option for the region
337  void ActivateSubCutoff(G4bool val, const G4Region* region = 0);
338 
339  // Activate biasing
340  void SetCrossSectionBiasingFactor(G4double f, G4bool flag = true);
341 
342  void ActivateForcedInteraction(G4double length = 0.0,
343  const G4String& region = "",
344  G4bool flag = true);
345 
346  void ActivateSecondaryBiasing(const G4String& region, G4double factor,
347  G4double energyLimit);
348 
349  // Add subcutoff process (bremsstrahlung) to sample secondary
350  // particle production in vicinity of the geometry boundary
351  void AddCollaborativeProcess(G4VEnergyLossProcess*);
352 
353  inline void SetLossFluctuations(G4bool val);
354  inline void SetRandomStep(G4bool val);
355 
356  inline void SetIntegral(G4bool val);
357  inline G4bool IsIntegral() const;
358 
359  // Set/Get flag "isIonisation"
360  inline void SetIonisation(G4bool val);
361  inline G4bool IsIonisationProcess() const;
362 
363  // Redefine parameteters for stepping control
364  inline void SetLinearLossLimit(G4double val);
365  inline void SetMinSubRange(G4double val);
366  inline void SetLambdaFactor(G4double val);
367  inline void SetStepFunction(G4double v1, G4double v2);
368  inline void SetLowestEnergyLimit(G4double);
369 
370  inline G4int NumberOfSubCutoffRegions() const;
371 
372  //------------------------------------------------------------------------
373  // Specific methods to path Physics Tables to the process
374  //------------------------------------------------------------------------
375 
376  void SetDEDXTable(G4PhysicsTable* p, G4EmTableType tType);
377  void SetCSDARangeTable(G4PhysicsTable* pRange);
378  void SetRangeTableForLoss(G4PhysicsTable* p);
379  void SetSecondaryRangeTable(G4PhysicsTable* p);
380  void SetInverseRangeTable(G4PhysicsTable* p);
381 
382  void SetLambdaTable(G4PhysicsTable* p);
383  void SetSubLambdaTable(G4PhysicsTable* p);
384 
385  // Binning for dEdx, range, inverse range and labda tables
386  inline void SetDEDXBinning(G4int nbins);
387  inline void SetLambdaBinning(G4int nbins);
388 
389  // Binning for dEdx, range, and inverse range tables
390  inline void SetDEDXBinningForCSDARange(G4int nbins);
391 
392  // Min kinetic energy for tables
393  inline void SetMinKinEnergy(G4double e);
394  inline G4double MinKinEnergy() const;
395 
396  // Max kinetic energy for tables
397  inline void SetMaxKinEnergy(G4double e);
398  inline G4double MaxKinEnergy() const;
399 
400  // Max kinetic energy for tables
401  inline void SetMaxKinEnergyForCSDARange(G4double e);
402 
403  // Biasing parameters
404  inline G4double CrossSectionBiasingFactor() const;
405 
406  // Return values for given G4MaterialCutsCouple
407  inline G4double GetDEDX(G4double& kineticEnergy, const G4MaterialCutsCouple*);
408  inline G4double GetDEDXForSubsec(G4double& kineticEnergy,
409  const G4MaterialCutsCouple*);
410  inline G4double GetRange(G4double& kineticEnergy, const G4MaterialCutsCouple*);
411  inline G4double GetCSDARange(G4double& kineticEnergy, const G4MaterialCutsCouple*);
412  inline G4double GetRangeForLoss(G4double& kineticEnergy, const G4MaterialCutsCouple*);
413  inline G4double GetKineticEnergy(G4double& range, const G4MaterialCutsCouple*);
414  inline G4double GetLambda(G4double& kineticEnergy, const G4MaterialCutsCouple*);
415 
416  inline G4bool TablesAreBuilt() const;
417 
418  // Access to specific tables
419  inline G4PhysicsTable* DEDXTable() const;
420  inline G4PhysicsTable* DEDXTableForSubsec() const;
421  inline G4PhysicsTable* DEDXunRestrictedTable() const;
422  inline G4PhysicsTable* IonisationTable() const;
423  inline G4PhysicsTable* IonisationTableForSubsec() const;
424  inline G4PhysicsTable* CSDARangeTable() const;
425  inline G4PhysicsTable* RangeTableForLoss() const;
426  inline G4PhysicsTable* InverseRangeTable() const;
427  inline G4PhysicsTable* LambdaTable();
428  inline G4PhysicsTable* SubLambdaTable();
429 
430  //------------------------------------------------------------------------
431  // Run time method for simulation of ionisation
432  //------------------------------------------------------------------------
433 
434  // access atom on which interaction happens
435  const G4Element* GetCurrentElement() const;
436 
437  // sample range at the end of a step
438  // inline G4double SampleRange();
439 
440  // Set scaling parameters for ions is needed to G4EmCalculator
441  inline void SetDynamicMassCharge(G4double massratio, G4double charge2ratio);
442 
443 private:
444 
445  void FillSecondariesAlongStep(G4double& eloss, G4double& weight);
446 
447  // define material and indexes
448  inline void DefineMaterial(const G4MaterialCutsCouple* couple);
449 
450  //------------------------------------------------------------------------
451  // Compute values using scaling relation, mass and charge of based particle
452  //------------------------------------------------------------------------
453 
454  inline G4double GetDEDXForScaledEnergy(G4double scaledKinEnergy);
455  inline G4double GetSubDEDXForScaledEnergy(G4double scaledKinEnergy);
456  inline G4double GetIonisationForScaledEnergy(G4double scaledKinEnergy);
457  inline G4double GetSubIonisationForScaledEnergy(G4double scaledKinEnergy);
458  inline G4double GetScaledRangeForScaledEnergy(G4double scaledKinEnergy);
459  inline G4double GetLimitScaledRangeForScaledEnergy(G4double scaledKinEnergy);
460  inline G4double ScaledKinEnergyForLoss(G4double range);
461  inline G4double GetLambdaForScaledEnergy(G4double scaledKinEnergy);
462  inline void ComputeLambdaForScaledEnergy(G4double scaledKinEnergy);
463 
464  // hide assignment operator
465  G4VEnergyLossProcess(G4VEnergyLossProcess &);
466  G4VEnergyLossProcess & operator=(const G4VEnergyLossProcess &right);
467 
468  // ======== Parameters of the class fixed at construction =========
469 
470  G4EmModelManager* modelManager;
471  G4EmBiasingManager* biasManager;
472  G4SafetyHelper* safetyHelper;
473 
474  const G4ParticleDefinition* secondaryParticle;
475  const G4ParticleDefinition* theElectron;
476  const G4ParticleDefinition* thePositron;
477  const G4ParticleDefinition* theGamma;
478  const G4ParticleDefinition* theGenericIon;
479 
480  // G4PhysicsVector* vstrag;
481 
482  // ======== Parameters of the class fixed at initialisation =======
483 
484  std::vector<G4VEmModel*> emModels;
485  G4VEmFluctuationModel* fluctModel;
486  G4VAtomDeexcitation* atomDeexcitation;
487  std::vector<const G4Region*> scoffRegions;
488  G4int nSCoffRegions;
489  G4bool* idxSCoffRegions;
490 
491  std::vector<G4VEnergyLossProcess*> scProcesses;
492  G4int nProcesses;
493 
494  // tables and vectors
495  G4PhysicsTable* theDEDXTable;
496  G4PhysicsTable* theDEDXSubTable;
497  G4PhysicsTable* theDEDXunRestrictedTable;
498  G4PhysicsTable* theIonisationTable;
499  G4PhysicsTable* theIonisationSubTable;
500  G4PhysicsTable* theRangeTableForLoss;
501  G4PhysicsTable* theCSDARangeTable;
502  G4PhysicsTable* theSecondaryRangeTable;
503  G4PhysicsTable* theInverseRangeTable;
504  G4PhysicsTable* theLambdaTable;
505  G4PhysicsTable* theSubLambdaTable;
506 
507  std::vector<G4double> theDEDXAtMaxEnergy;
508  std::vector<G4double> theRangeAtMaxEnergy;
509  std::vector<G4double> theEnergyOfCrossSectionMax;
510  std::vector<G4double> theCrossSectionMax;
511 
512  const std::vector<G4double>* theDensityFactor;
513  const std::vector<G4int>* theDensityIdx;
514 
515  const G4DataVector* theCuts;
516  const G4DataVector* theSubCuts;
517 
518  const G4ParticleDefinition* baseParticle;
519 
520  G4int nBins;
521  G4int nBinsCSDA;
522 
523  G4double lowestKinEnergy;
524  G4double minKinEnergy;
525  G4double maxKinEnergy;
526  G4double maxKinEnergyCSDA;
527 
528  G4double linLossLimit;
529  G4double minSubRange;
530  G4double dRoverRange;
531  G4double finalRange;
532  G4double lambdaFactor;
533  G4double biasFactor;
534 
535  G4bool lossFluctuationFlag;
536  G4bool rndmStepFlag;
537  G4bool tablesAreBuilt;
538  G4bool integral;
539  G4bool isIon;
540  G4bool isIonisation;
541  G4bool useSubCutoff;
542  G4bool useDeexcitation;
543  G4bool biasFlag;
544  G4bool weightFlag;
545 
546 protected:
547 
548  G4ParticleChangeForLoss fParticleChange;
549 
550  // ======== Cashed values - may be state dependent ================
551 
552 private:
553 
554  std::vector<G4DynamicParticle*> secParticles;
555  std::vector<G4Track*> scTracks;
556 
557  const G4ParticleDefinition* particle;
558 
559  G4VEmModel* currentModel;
560  const G4Material* currentMaterial;
561  const G4MaterialCutsCouple* currentCouple;
562  size_t currentCoupleIndex;
563  size_t basedCoupleIndex;
564 
565  G4int nWarnings;
566 
567  G4double massRatio;
568  G4double fFactor;
569  G4double reduceFactor;
570  G4double chargeSqRatio;
571 
572  G4double preStepLambda;
573  G4double fRange;
574  G4double preStepKinEnergy;
575  G4double preStepScaledEnergy;
576  G4double mfpKinEnergy;
577 
578  G4GPILSelection aGPILSelection;
579 
580 };
581 
582 // ======== Run time inline methods ================
583 
584 inline size_t G4VEnergyLossProcess::CurrentMaterialCutsCoupleIndex() const
585 {
586  return currentCoupleIndex;
587 }
588 
589 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
590 
591 inline G4double G4VEnergyLossProcess::GetCurrentRange() const
592 {
593  return fRange;
594 }
595 
596 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
597 
598 inline void G4VEnergyLossProcess::SelectModel(G4double kinEnergy)
599 {
600  currentModel = modelManager->SelectModel(kinEnergy, currentCoupleIndex);
601  currentModel->SetCurrentCouple(currentCouple);
602 }
603 
604 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
605 
606 inline G4VEmModel* G4VEnergyLossProcess::SelectModelForMaterial(
607  G4double kinEnergy, size_t& idx) const
608 {
609  return modelManager->SelectModel(kinEnergy, idx);
610 }
611 
612 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
613 
614 inline void
615 G4VEnergyLossProcess::DefineMaterial(const G4MaterialCutsCouple* couple)
616 {
617  if(couple != currentCouple) {
618  currentCouple = couple;
619  currentMaterial = couple->GetMaterial();
620  currentCoupleIndex = couple->GetIndex();
621  basedCoupleIndex = (*theDensityIdx)[currentCoupleIndex];
622  fFactor = chargeSqRatio*biasFactor*(*theDensityFactor)[currentCoupleIndex];
623  reduceFactor = 1.0/(fFactor*massRatio);
624  mfpKinEnergy = DBL_MAX;
625  }
626 }
627 
628 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
629 
630 inline void G4VEnergyLossProcess::SetDynamicMassCharge(G4double massratio,
631  G4double charge2ratio)
632 {
633  massRatio = massratio;
634  fFactor *= charge2ratio/chargeSqRatio;
635  chargeSqRatio = charge2ratio;
636  reduceFactor = 1.0/(fFactor*massRatio);
637 }
638 
639 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
640 
641 inline G4double G4VEnergyLossProcess::GetDEDXForScaledEnergy(G4double e)
642 {
643  //G4cout << "G4VEnergyLossProcess::GetDEDX: Idx= "
644  // << basedCoupleIndex << " E(MeV)= " << e << G4endl;
645  G4double x = fFactor*(*theDEDXTable)[basedCoupleIndex]->Value(e);
646  if(e < minKinEnergy) { x *= std::sqrt(e/minKinEnergy); }
647  return x;
648 }
649 
650 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
651 
652 inline G4double G4VEnergyLossProcess::GetSubDEDXForScaledEnergy(G4double e)
653 {
654  G4double x = fFactor*(*theDEDXSubTable)[basedCoupleIndex]->Value(e);
655  if(e < minKinEnergy) { x *= std::sqrt(e/minKinEnergy); }
656  return x;
657 }
658 
659 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
660 
661 inline G4double G4VEnergyLossProcess::GetIonisationForScaledEnergy(G4double e)
662 {
663  G4double x = fFactor*(*theIonisationTable)[basedCoupleIndex]->Value(e);
664  if(e < minKinEnergy) { x *= std::sqrt(e/minKinEnergy); }
665  return x;
666 }
667 
668 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
669 
670 inline
671 G4double G4VEnergyLossProcess::GetSubIonisationForScaledEnergy(G4double e)
672 {
673  G4double x = fFactor*(*theIonisationSubTable)[basedCoupleIndex]->Value(e);
674  if(e < minKinEnergy) { x *= std::sqrt(e/minKinEnergy); }
675  return x;
676 }
677 
678 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
679 
680 inline G4double G4VEnergyLossProcess::GetScaledRangeForScaledEnergy(G4double e)
681 {
682  //G4cout << "G4VEnergyLossProcess::GetRange: Idx= "
683  // << basedCoupleIndex << " E(MeV)= " << e << G4endl;
684  G4double x = ((*theRangeTableForLoss)[basedCoupleIndex])->Value(e);
685  if(e < minKinEnergy) { x *= std::sqrt(e/minKinEnergy); }
686  return x;
687 }
688 
689 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
690 
691 inline G4double
692 G4VEnergyLossProcess::GetLimitScaledRangeForScaledEnergy(G4double e)
693 {
694  G4double x;
695  if (e < maxKinEnergyCSDA) {
696  x = ((*theCSDARangeTable)[basedCoupleIndex])->Value(e);
697  if(e < minKinEnergy) { x *= std::sqrt(e/minKinEnergy); }
698  } else {
699  x = theRangeAtMaxEnergy[basedCoupleIndex] +
700  (e - maxKinEnergyCSDA)/theDEDXAtMaxEnergy[basedCoupleIndex];
701  }
702  return x;
703 }
704 
705 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
706 
707 inline G4double G4VEnergyLossProcess::ScaledKinEnergyForLoss(G4double r)
708 {
709  //G4cout << "G4VEnergyLossProcess::GetEnergy: Idx= "
710  // << basedCoupleIndex << " R(mm)= " << r << G4endl;
711  G4PhysicsVector* v = (*theInverseRangeTable)[basedCoupleIndex];
712  G4double rmin = v->Energy(0);
713  G4double e = 0.0;
714  if(r >= rmin) { e = v->Value(r); }
715  else if(r > 0.0) {
716  G4double x = r/rmin;
717  e = minKinEnergy*x*x;
718  }
719  return e;
720 }
721 
722 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
723 
724 inline G4double G4VEnergyLossProcess::GetLambdaForScaledEnergy(G4double e)
725 {
726  return fFactor*((*theLambdaTable)[basedCoupleIndex])->Value(e);
727 }
728 
729 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
730 
731 inline G4double
732 G4VEnergyLossProcess::GetDEDX(G4double& kineticEnergy,
733  const G4MaterialCutsCouple* couple)
734 {
735  DefineMaterial(couple);
736  return GetDEDXForScaledEnergy(kineticEnergy*massRatio);
737 }
738 
739 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
740 
741 inline G4double
742 G4VEnergyLossProcess::GetDEDXForSubsec(G4double& kineticEnergy,
743  const G4MaterialCutsCouple* couple)
744 {
745  DefineMaterial(couple);
746  return GetSubDEDXForScaledEnergy(kineticEnergy*massRatio);
747 }
748 
749 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
750 
751 inline G4double
752 G4VEnergyLossProcess::GetRange(G4double& kineticEnergy,
753  const G4MaterialCutsCouple* couple)
754 {
755  G4double x = fRange;
756  if(couple != currentCouple || kineticEnergy != preStepKinEnergy) {
757  DefineMaterial(couple);
758  if(theCSDARangeTable) {
759  x = GetLimitScaledRangeForScaledEnergy(kineticEnergy*massRatio)
760  * reduceFactor;
761  } else if(theRangeTableForLoss) {
762  x = GetScaledRangeForScaledEnergy(kineticEnergy*massRatio)*reduceFactor;
763  }
764  }
765  return x;
766 }
767 
768 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
769 
770 inline G4double
771 G4VEnergyLossProcess::GetCSDARange(G4double& kineticEnergy,
772  const G4MaterialCutsCouple* couple)
773 {
774  DefineMaterial(couple);
775  G4double x = DBL_MAX;
776  if(theCSDARangeTable) {
777  x = GetLimitScaledRangeForScaledEnergy(kineticEnergy*massRatio)*reduceFactor;
778  }
779  return x;
780 }
781 
782 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
783 
784 inline G4double
785 G4VEnergyLossProcess::GetRangeForLoss(G4double& kineticEnergy,
786  const G4MaterialCutsCouple* couple)
787 {
788  G4double x = fRange;
789  if(couple != currentCouple || kineticEnergy != preStepKinEnergy) {
790  DefineMaterial(couple);
791  x = GetScaledRangeForScaledEnergy(kineticEnergy*massRatio)*reduceFactor;
792  }
793  // G4cout << "Range from " << GetProcessName()
794  // << " e= " << kineticEnergy << " r= " << x << G4endl;
795  return x;
796 }
797 
798 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
799 
800 inline G4double
801 G4VEnergyLossProcess::GetKineticEnergy(G4double& range,
802  const G4MaterialCutsCouple* couple)
803 {
804  DefineMaterial(couple);
805  return ScaledKinEnergyForLoss(range/reduceFactor)/massRatio;
806 }
807 
808 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
809 
810 inline G4double
811 G4VEnergyLossProcess::GetLambda(G4double& kineticEnergy,
812  const G4MaterialCutsCouple* couple)
813 {
814  DefineMaterial(couple);
815  G4double x = 0.0;
816  if(theLambdaTable) { x = GetLambdaForScaledEnergy(kineticEnergy*massRatio); }
817  return x;
818 }
819 
820 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
821 
822 inline void G4VEnergyLossProcess::ComputeLambdaForScaledEnergy(G4double e)
823 {
824  mfpKinEnergy = theEnergyOfCrossSectionMax[currentCoupleIndex];
825  if (e <= mfpKinEnergy) {
826  preStepLambda = GetLambdaForScaledEnergy(e);
827 
828  } else {
829  G4double e1 = e*lambdaFactor;
830  if(e1 > mfpKinEnergy) {
831  preStepLambda = GetLambdaForScaledEnergy(e);
832  G4double preStepLambda1 = GetLambdaForScaledEnergy(e1);
833  if(preStepLambda1 > preStepLambda) {
834  mfpKinEnergy = e1;
835  preStepLambda = preStepLambda1;
836  }
837  } else {
838  preStepLambda = fFactor*theCrossSectionMax[currentCoupleIndex];
839  }
840  }
841 }
842 
843 // ======== Get/Set inline methods used at initialisation ================
844 
845 inline void G4VEnergyLossProcess::SetFluctModel(G4VEmFluctuationModel* p)
846 {
847  fluctModel = p;
848 }
849 
850 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
851 
852 inline G4VEmFluctuationModel* G4VEnergyLossProcess::FluctModel()
853 {
854  return fluctModel;
855 }
856 
857 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
858 
859 inline void G4VEnergyLossProcess::SetParticle(const G4ParticleDefinition* p)
860 {
861  particle = p;
862 }
863 
864 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
865 
866 inline void G4VEnergyLossProcess::SetSecondaryParticle(const G4ParticleDefinition* p)
867 {
868  secondaryParticle = p;
869 }
870 
871 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
872 
873 inline void G4VEnergyLossProcess::SetBaseParticle(const G4ParticleDefinition* p)
874 {
875  baseParticle = p;
876 }
877 
878 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
879 
880 inline const G4ParticleDefinition* G4VEnergyLossProcess::Particle() const
881 {
882  return particle;
883 }
884 
885 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
886 
887 inline const G4ParticleDefinition* G4VEnergyLossProcess::BaseParticle() const
888 {
889  return baseParticle;
890 }
891 
892 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
893 
894 inline const G4ParticleDefinition* G4VEnergyLossProcess::SecondaryParticle() const
895 {
896  return secondaryParticle;
897 }
898 
899 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
900 
901 inline void G4VEnergyLossProcess::SetLossFluctuations(G4bool val)
902 {
903  lossFluctuationFlag = val;
904 }
905 
906 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
907 
908 inline void G4VEnergyLossProcess::SetRandomStep(G4bool val)
909 {
910  rndmStepFlag = val;
911 }
912 
913 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
914 
915 inline void G4VEnergyLossProcess::SetIntegral(G4bool val)
916 {
917  integral = val;
918 }
919 
920 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
921 
922 inline G4bool G4VEnergyLossProcess::IsIntegral() const
923 {
924  return integral;
925 }
926 
927 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
928 
929 inline void G4VEnergyLossProcess::SetIonisation(G4bool val)
930 {
931  isIonisation = val;
932  if(val) { aGPILSelection = CandidateForSelection; }
933  else { aGPILSelection = NotCandidateForSelection; }
934 }
935 
936 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
937 
938 inline G4bool G4VEnergyLossProcess::IsIonisationProcess() const
939 {
940  return isIonisation;
941 }
942 
943 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
944 
945 inline void G4VEnergyLossProcess::SetLinearLossLimit(G4double val)
946 {
947  linLossLimit = val;
948 }
949 
950 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
951 
952 inline void G4VEnergyLossProcess::SetMinSubRange(G4double val)
953 {
954  minSubRange = val;
955 }
956 
957 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
958 
959 inline void G4VEnergyLossProcess::SetLambdaFactor(G4double val)
960 {
961  if(val > 0.0 && val <= 1.0) { lambdaFactor = val; }
962 }
963 
964 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
965 
966 void G4VEnergyLossProcess::SetStepFunction(G4double v1, G4double v2)
967 {
968  dRoverRange = v1;
969  finalRange = v2;
970  if (dRoverRange > 0.999) { dRoverRange = 1.0; }
971  currentCouple = 0;
972  mfpKinEnergy = DBL_MAX;
973 }
974 
975 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
976 
977 inline void G4VEnergyLossProcess::SetLowestEnergyLimit(G4double val)
978 {
979  lowestKinEnergy = val;
980 }
981 
982 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
983 
984 inline G4int G4VEnergyLossProcess::NumberOfSubCutoffRegions() const
985 {
986  return nSCoffRegions;
987 }
988 
989 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
990 
991 inline void G4VEnergyLossProcess::SetDEDXBinning(G4int nbins)
992 {
993  nBins = nbins;
994 }
995 
996 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
997 
998 inline void G4VEnergyLossProcess::SetLambdaBinning(G4int nbins)
999 {
1000  nBins = nbins;
1001 }
1002 
1003 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1004 
1005 inline void G4VEnergyLossProcess::SetDEDXBinningForCSDARange(G4int nbins)
1006 {
1007  nBinsCSDA = nbins;
1008 }
1009 
1010 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1011 
1012 inline void G4VEnergyLossProcess::SetMinKinEnergy(G4double e)
1013 {
1014  minKinEnergy = e;
1015 }
1016 
1017 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1018 
1019 inline G4double G4VEnergyLossProcess::MinKinEnergy() const
1020 {
1021  return minKinEnergy;
1022 }
1023 
1024 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1025 
1026 inline void G4VEnergyLossProcess::SetMaxKinEnergy(G4double e)
1027 {
1028  maxKinEnergy = e;
1029  if(e < maxKinEnergyCSDA) { maxKinEnergyCSDA = e; }
1030 }
1031 
1032 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1033 
1034 inline G4double G4VEnergyLossProcess::MaxKinEnergy() const
1035 {
1036  return maxKinEnergy;
1037 }
1038 
1039 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1040 
1041 inline void G4VEnergyLossProcess::SetMaxKinEnergyForCSDARange(G4double e)
1042 {
1043  maxKinEnergyCSDA = e;
1044 }
1045 
1046 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1047 
1048 inline G4double G4VEnergyLossProcess::CrossSectionBiasingFactor() const
1049 {
1050  return biasFactor;
1051 }
1052 
1053 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1054 
1055 inline G4bool G4VEnergyLossProcess::TablesAreBuilt() const
1056 {
1057  return tablesAreBuilt;
1058 }
1059 
1060 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1061 
1062 inline G4PhysicsTable* G4VEnergyLossProcess::DEDXTable() const
1063 {
1064  return theDEDXTable;
1065 }
1066 
1067 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1068 
1069 inline G4PhysicsTable* G4VEnergyLossProcess::DEDXTableForSubsec() const
1070 {
1071  return theDEDXSubTable;
1072 }
1073 
1074 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1075 
1076 inline G4PhysicsTable* G4VEnergyLossProcess::DEDXunRestrictedTable() const
1077 {
1078  return theDEDXunRestrictedTable;
1079 }
1080 
1081 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1082 
1083 inline G4PhysicsTable* G4VEnergyLossProcess::IonisationTable() const
1084 {
1085  G4PhysicsTable* t = theDEDXTable;
1086  if(theIonisationTable) { t = theIonisationTable; }
1087  return t;
1088 }
1089 
1090 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1091 
1092 inline G4PhysicsTable* G4VEnergyLossProcess::IonisationTableForSubsec() const
1093 {
1094  G4PhysicsTable* t = theDEDXSubTable;
1095  if(theIonisationSubTable) { t = theIonisationSubTable; }
1096  return t;
1097 }
1098 
1099 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1100 
1101 inline G4PhysicsTable* G4VEnergyLossProcess::CSDARangeTable() const
1102 {
1103  return theCSDARangeTable;
1104 }
1105 
1106 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1107 
1108 inline G4PhysicsTable* G4VEnergyLossProcess::RangeTableForLoss() const
1109 {
1110  return theRangeTableForLoss;
1111 }
1112 
1113 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1114 
1115 inline G4PhysicsTable* G4VEnergyLossProcess::InverseRangeTable() const
1116 {
1117  return theInverseRangeTable;
1118 }
1119 
1120 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1121 
1122 inline G4PhysicsTable* G4VEnergyLossProcess::LambdaTable()
1123 {
1124  return theLambdaTable;
1125 }
1126 
1127 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1128 
1129 inline G4PhysicsTable* G4VEnergyLossProcess::SubLambdaTable()
1130 {
1131  return theSubLambdaTable;
1132 }
1133 
1134 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
1135 
1136 #endif
Generated on Sat May 25 2013 14:33:38 for Geant4 by   doxygen 1.8.4