source: trunk/source/processes/electromagnetic/utils/include/G4VEmProcess.hh @ 1315

Last change on this file since 1315 was 1315, checked in by garnier, 14 years ago

update geant4-09-04-beta-cand-01 interfaces-V09-03-09 vis-V09-03-08
File size: 19.6 KB
Line 
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: G4VEmProcess.hh,v 1.60 2010/04/28 14:43:13 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-04-beta-cand-01 $
28//
29// -------------------------------------------------------------------
30//
31// GEANT4 Class header file
32//
33//
34// File name:     G4VEmProcess
35//
36// Author:        Vladimir Ivanchenko
37//
38// Creation date: 01.10.2003
39//
40// Modifications:
41// 30-06-04 make destructor virtual (V.Ivanchenko)
42// 09-08-04 optimise integral option (V.Ivanchenko)
43// 11-08-04 add protected methods to access cuts (V.Ivanchenko)
44// 09-09-04 Bug fix for the integral mode with 2 peaks (V.Ivanchneko)
45// 16-09-04 Add flag for LambdaTable and method RecalculateLambda (VI)
46// 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivantchenko)
47// 08-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
48// 18-04-05 Use G4ParticleChangeForGamma (V.Ivantchenko)
49// 09-05-05 Fix problem in logic when path boundary between materials (VI)
50// 11-01-06 add A to parameters of ComputeCrossSectionPerAtom (VI)
51// 01-02-06 put default value A=0. to keep compatibility with v5.2 (mma)
52// 13-05-06 Add method to access model by index (V.Ivanchenko)
53// 12-09-06 add SetModel() (mma)
54// 25-09-07 More accurate handling zero xsect in
55//          PostStepGetPhysicalInteractionLength (V.Ivanchenko)
56// 27-10-07 Virtual functions moved to source (V.Ivanchenko)
57// 15-07-08 Reorder class members for further multi-thread development (VI)
58// 17-02-10 Added pointer currentParticle (VI)
59//
60// Class Description:
61//
62// It is the unified Discrete process
63
64// -------------------------------------------------------------------
65//
66
67#ifndef G4VEmProcess_h
68#define G4VEmProcess_h 1
69
70#include "G4VDiscreteProcess.hh"
71#include "globals.hh"
72#include "G4Material.hh"
73#include "G4MaterialCutsCouple.hh"
74#include "G4Track.hh"
75#include "G4EmModelManager.hh"
76#include "G4UnitsTable.hh"
77#include "G4ParticleDefinition.hh"
78#include "G4ParticleChangeForGamma.hh"
79
80class G4Step;
81class G4VEmModel;
82class G4DataVector;
83class G4VParticleChange;
84class G4PhysicsTable;
85class G4PhysicsVector;
86
87//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88
89class G4VEmProcess : public G4VDiscreteProcess
90{
91public:
92
93  G4VEmProcess(const G4String& name,
94               G4ProcessType type = fElectromagnetic);
95
96  virtual ~G4VEmProcess();
97
98  //------------------------------------------------------------------------
99  // Virtual methods to be implemented in concrete processes
100  //------------------------------------------------------------------------
101
102  virtual G4bool IsApplicable(const G4ParticleDefinition& p) = 0;
103
104  virtual void PrintInfo() = 0;
105
106protected:
107
108  virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
109
110  //------------------------------------------------------------------------
111  // Implementation of virtual methods common to all Discrete processes
112  //------------------------------------------------------------------------
113
114public:
115
116  // Initialise for build of tables
117  void PreparePhysicsTable(const G4ParticleDefinition&);
118
119  // Build physics table during initialisation
120  void BuildPhysicsTable(const G4ParticleDefinition&);
121
122  void PrintInfoDefinition();
123
124  // implementation of virtual method, specific for G4VEmProcess
125  G4double PostStepGetPhysicalInteractionLength(
126                             const G4Track& track,
127                             G4double   previousStepSize,
128                             G4ForceCondition* condition
129                            );
130
131  // implementation of virtual method, specific for G4VEmProcess
132  G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&);
133
134  // Store PhysicsTable in a file.
135  // Return false in case of failure at I/O
136  G4bool StorePhysicsTable(const G4ParticleDefinition*,
137                           const G4String& directory,
138                           G4bool ascii = false);
139
140  // Retrieve Physics from a file.
141  // (return true if the Physics Table can be build by using file)
142  // (return false if the process has no functionality or in case of failure)
143  // File name should is constructed as processName+particleName and the
144  // should be placed under the directory specifed by the argument.
145  G4bool RetrievePhysicsTable(const G4ParticleDefinition*,
146                              const G4String& directory,
147                              G4bool ascii);
148
149  // deexcitation activated per G4Region
150  void ActivateDeexcitation(G4bool, const G4Region* r = 0);
151
152  //------------------------------------------------------------------------
153  // Specific methods for Discrete EM post step simulation
154  //------------------------------------------------------------------------
155
156  // It returns the cross section per volume for energy/ material
157  G4double CrossSectionPerVolume(G4double kineticEnergy,
158                                 const G4MaterialCutsCouple* couple);
159
160  // It returns the cross section of the process per atom
161  G4double ComputeCrossSectionPerAtom(G4double kineticEnergy, 
162                                      G4double Z, G4double A=0., 
163                                      G4double cut=0.0);
164
165  G4double MeanFreePath(const G4Track& track);
166
167  // It returns cross section per volume
168  inline G4double GetLambda(G4double& kinEnergy, 
169                            const G4MaterialCutsCouple* couple);
170
171  //------------------------------------------------------------------------
172  // Specific methods to build and access Physics Tables
173  //------------------------------------------------------------------------
174
175  // Binning for lambda table
176  inline void SetLambdaBinning(G4int nbins);
177  inline G4int LambdaBinning() const;
178
179  // Min kinetic energy for tables
180  inline void SetMinKinEnergy(G4double e);
181  inline G4double MinKinEnergy() const;
182
183  // Max kinetic energy for tables
184  inline void SetMaxKinEnergy(G4double e);
185  inline G4double MaxKinEnergy() const;
186
187  inline void SetPolarAngleLimit(G4double a);
188  inline G4double PolarAngleLimit() const;
189
190  inline const G4PhysicsTable* LambdaTable() const;
191
192  //------------------------------------------------------------------------
193  // Define and access particle type
194  //------------------------------------------------------------------------
195
196  inline const G4ParticleDefinition* Particle() const;
197  inline const G4ParticleDefinition* SecondaryParticle() const;
198
199  //------------------------------------------------------------------------
200  // Specific methods to set, access, modify models and basic parameters
201  //------------------------------------------------------------------------
202
203protected:
204  // Select model in run time
205  inline G4VEmModel* SelectModel(G4double& kinEnergy, size_t index); 
206
207public:
208  // Select model by energy and region index
209  inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy, 
210                                            size_t& idxRegion) const;
211   
212  // Add model for region, smaller value of order defines which
213  // model will be selected for a given energy interval 
214  void AddEmModel(G4int, G4VEmModel*, const G4Region* region = 0);
215
216  // Assign a model to a process
217  void SetModel(G4VEmModel*, G4int index = 1);
218 
219  // return the assigned model
220  G4VEmModel* Model(G4int index = 1);
221   
222  // Define new energy range for the model identified by the name
223  void UpdateEmModel(const G4String&, G4double, G4double);
224
225  // Access to models
226  G4VEmModel* GetModelByIndex(G4int idx = 0, G4bool ver = false);
227
228  // access atom on which interaction happens
229  const G4Element* GetCurrentElement() const;
230
231  inline void SetLambdaFactor(G4double val);
232
233  inline void SetIntegral(G4bool val);
234  inline G4bool IsIntegral() const;
235
236  inline void SetApplyCuts(G4bool val);
237
238  inline void SetBuildTableFlag(G4bool val);
239
240  //------------------------------------------------------------------------
241  // Other generic methods
242  //------------------------------------------------------------------------
243 
244protected:
245
246  G4double GetMeanFreePath(const G4Track& track,
247                           G4double previousStepSize,
248                           G4ForceCondition* condition);
249
250  G4PhysicsVector* LambdaPhysicsVector(const G4MaterialCutsCouple*);
251
252  inline G4double RecalculateLambda(G4double kinEnergy,
253                                    const G4MaterialCutsCouple* couple);
254
255  inline G4ParticleChangeForGamma* GetParticleChange();
256
257  inline void SetParticle(const G4ParticleDefinition* p);
258 
259  inline void SetSecondaryParticle(const G4ParticleDefinition* p);
260
261  inline size_t CurrentMaterialCutsCoupleIndex() const;
262
263  inline G4double GetGammaEnergyCut();
264
265  inline G4double GetElectronEnergyCut();
266
267  inline void SetStartFromNullFlag(G4bool val);
268
269private:
270
271  void Clear();
272
273  void BuildLambdaTable();
274
275  void FindLambdaMax();
276
277  inline void InitialiseStep(const G4Track&);
278
279  inline void DefineMaterial(const G4MaterialCutsCouple* couple);
280
281  inline void ComputeIntegralLambda(G4double kinEnergy);
282
283  inline G4double GetLambdaFromTable(G4double kinEnergy);
284
285  inline G4double GetCurrentLambda(G4double kinEnergy);
286
287  inline G4double ComputeCurrentLambda(G4double kinEnergy);
288
289  // copy constructor and hide assignment operator
290  G4VEmProcess(G4VEmProcess &);
291  G4VEmProcess & operator=(const G4VEmProcess &right);
292
293  // ======== Parameters of the class fixed at construction =========
294
295  G4EmModelManager*            modelManager;
296  const G4ParticleDefinition*  theGamma;
297  const G4ParticleDefinition*  theElectron;
298  const G4ParticleDefinition*  thePositron;
299  const G4ParticleDefinition*  secondaryParticle;
300
301  G4bool                       buildLambdaTable;
302
303  // ======== Parameters of the class fixed at initialisation =======
304
305  std::vector<G4VEmModel*>     emModels;
306
307  // tables and vectors
308  G4PhysicsTable*              theLambdaTable;
309  G4double*                    theEnergyOfCrossSectionMax;
310  G4double*                    theCrossSectionMax;
311
312  const std::vector<G4double>* theCuts;
313  const std::vector<G4double>* theCutsGamma;
314  const std::vector<G4double>* theCutsElectron;
315  const std::vector<G4double>* theCutsPositron;
316
317  G4int                        nLambdaBins;
318
319  G4double                     minKinEnergy;
320  G4double                     maxKinEnergy;
321  G4double                     lambdaFactor;
322  G4double                     polarAngleLimit;
323
324  G4bool                       integral;
325  G4bool                       applyCuts;
326  G4bool                       startFromNull;
327  G4bool                       useDeexcitation;
328
329  G4int                        nDERegions;
330  std::vector<const G4Region*> deRegions;
331  G4bool*                      idxDERegions;
332
333  // ======== Cashed values - may be state dependent ================
334
335protected:
336
337  G4ParticleChangeForGamma     fParticleChange;
338
339private:
340
341  std::vector<G4DynamicParticle*> secParticles;
342
343  G4VEmModel*                  currentModel; 
344
345  const G4ParticleDefinition*  particle;
346  const G4ParticleDefinition*  currentParticle;
347
348  // cash
349  const G4Material*            currentMaterial;
350  const G4MaterialCutsCouple*  currentCouple;
351  size_t                       currentCoupleIndex;
352
353  G4double                     mfpKinEnergy;
354  G4double                     preStepKinEnergy;
355  G4double                     preStepLambda;
356
357};
358
359// ======== Run time inline methods ================
360
361inline size_t G4VEmProcess::CurrentMaterialCutsCoupleIndex() const 
362{
363  return currentCoupleIndex;
364}
365
366//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
367
368inline G4double G4VEmProcess::GetGammaEnergyCut()
369{
370  return (*theCutsGamma)[currentCoupleIndex];
371}
372
373//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
374
375inline G4double G4VEmProcess::GetElectronEnergyCut()
376{
377  return (*theCutsElectron)[currentCoupleIndex];
378}
379
380//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
381
382inline void G4VEmProcess::DefineMaterial(const G4MaterialCutsCouple* couple)
383{
384  if(couple != currentCouple) {
385    currentCouple   = couple;
386    currentMaterial = couple->GetMaterial();
387    currentCoupleIndex = couple->GetIndex();
388    mfpKinEnergy = DBL_MAX;
389  }
390}
391
392//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
393
394inline 
395G4VEmModel* G4VEmProcess::SelectModel(G4double& kinEnergy, size_t index)
396{
397  currentModel = modelManager->SelectModel(kinEnergy, index);
398  currentModel->SetCurrentCouple(currentCouple);
399  return currentModel;
400}
401
402//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
403
404inline 
405G4VEmModel* G4VEmProcess::SelectModelForMaterial(G4double kinEnergy, 
406                                                 size_t& idxRegion) const
407{
408  return modelManager->SelectModel(kinEnergy, idxRegion);
409}
410
411//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
412
413inline void G4VEmProcess::InitialiseStep(const G4Track& track)
414{
415  currentParticle = track.GetDefinition();
416  preStepKinEnergy = track.GetKineticEnergy();
417  DefineMaterial(track.GetMaterialCutsCouple());
418  SelectModel(preStepKinEnergy, currentCoupleIndex);
419  if (theNumberOfInteractionLengthLeft < 0.0) mfpKinEnergy = DBL_MAX;
420}
421
422//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
423
424inline G4double G4VEmProcess::GetLambdaFromTable(G4double e)
425{
426  return (((*theLambdaTable)[currentCoupleIndex])->Value(e));
427}
428
429//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
430
431inline G4double G4VEmProcess::ComputeCurrentLambda(G4double e)
432{
433  SelectModel(e, currentCoupleIndex);
434  return currentModel->CrossSectionPerVolume(currentMaterial,currentParticle,
435                                             e,(*theCuts)[currentCoupleIndex]);
436}
437
438//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
439
440inline G4double G4VEmProcess::GetCurrentLambda(G4double e)
441{
442  G4double x = 0.0;
443  if(theLambdaTable) { x = GetLambdaFromTable(e); }
444  else               { x = ComputeCurrentLambda(e); }
445  return x;
446}
447
448//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
449
450inline G4double G4VEmProcess::GetLambda(G4double& kineticEnergy,
451                                        const G4MaterialCutsCouple* couple)
452{
453  DefineMaterial(couple);
454  return GetCurrentLambda(kineticEnergy);
455}
456
457//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
458
459inline G4double G4VEmProcess::RecalculateLambda(G4double e, 
460                                                const G4MaterialCutsCouple* couple)
461{
462  DefineMaterial(couple);
463  return ComputeCurrentLambda(e);
464}
465
466//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
467
468inline void G4VEmProcess::ComputeIntegralLambda(G4double e)
469{
470  mfpKinEnergy  = theEnergyOfCrossSectionMax[currentCoupleIndex];
471  if (e <= mfpKinEnergy) {
472    preStepLambda = GetLambdaFromTable(e);
473
474  } else {
475    G4double e1 = e*lambdaFactor;
476    if(e1 > mfpKinEnergy) {
477      preStepLambda  = GetLambdaFromTable(e);
478      G4double preStepLambda1 = GetLambdaFromTable(e1);
479      if(preStepLambda1 > preStepLambda) {
480        mfpKinEnergy = e1;
481        preStepLambda = preStepLambda1;
482      }
483    } else {
484      preStepLambda = theCrossSectionMax[currentCoupleIndex];
485    }
486  }
487}
488
489// ======== Get/Set inline methods used at initialisation ================
490
491inline void G4VEmProcess::SetLambdaBinning(G4int nbins)
492{
493  nLambdaBins = nbins;
494}
495
496//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
497
498inline G4int G4VEmProcess::LambdaBinning() const
499{
500  return nLambdaBins;
501}
502
503//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
504
505inline void G4VEmProcess::SetMinKinEnergy(G4double e)
506{
507  minKinEnergy = e;
508}
509
510//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
511
512inline G4double G4VEmProcess::MinKinEnergy() const
513{
514  return minKinEnergy;
515}
516
517//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
518
519inline void G4VEmProcess::SetMaxKinEnergy(G4double e)
520{
521  maxKinEnergy = e;
522}
523
524//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
525
526inline G4double G4VEmProcess::MaxKinEnergy() const
527{
528  return maxKinEnergy;
529}
530
531//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
532
533inline void G4VEmProcess::SetPolarAngleLimit(G4double val)
534{
535  if(val < 0.0)     polarAngleLimit = 0.0;
536  else if(val > pi) polarAngleLimit = pi;
537  else              polarAngleLimit = val;
538}
539
540//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
541
542inline G4double G4VEmProcess::PolarAngleLimit() const
543{
544  return polarAngleLimit;
545}
546
547//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
548
549inline const G4PhysicsTable* G4VEmProcess::LambdaTable() const
550{
551  return theLambdaTable;
552}
553
554//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
555
556inline const G4ParticleDefinition* G4VEmProcess::Particle() const
557{
558  return particle;
559}
560
561//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
562
563inline const G4ParticleDefinition* G4VEmProcess::SecondaryParticle() const
564{
565  return secondaryParticle;
566}
567
568//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
569
570inline void G4VEmProcess::SetLambdaFactor(G4double val)
571{
572  if(val > 0.0 && val <= 1.0) { lambdaFactor = val; }
573}
574
575//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
576
577inline void G4VEmProcess::SetIntegral(G4bool val)
578{
579  if(particle && particle != theGamma) { integral = val; }
580  if(integral) { buildLambdaTable = true; }
581}
582
583//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
584
585inline G4bool G4VEmProcess::IsIntegral() const
586{
587  return integral;
588}
589
590//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
591
592inline void G4VEmProcess::SetApplyCuts(G4bool val)
593{
594  applyCuts = val;
595}
596
597//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
598
599inline void G4VEmProcess::SetBuildTableFlag(G4bool val)
600{
601  buildLambdaTable = val;
602  if(!val) { integral = false; }
603}
604
605//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
606
607inline G4ParticleChangeForGamma* G4VEmProcess::GetParticleChange()
608{
609  return &fParticleChange;
610}
611
612//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
613
614inline void G4VEmProcess::SetParticle(const G4ParticleDefinition* p)
615{
616  particle = p;
617  currentParticle = p;
618}
619
620//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
621
622inline void G4VEmProcess::SetSecondaryParticle(const G4ParticleDefinition* p)
623{
624  secondaryParticle = p;
625}
626
627//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
628
629inline void G4VEmProcess::SetStartFromNullFlag(G4bool val)
630{
631  startFromNull = val;
632}
633
634//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
635
636#endif
Note: See TracBrowser for help on using the repository browser.