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

Last change on this file since 846 was 819, checked in by garnier, 16 years ago

import all except CVS
File size: 20.1 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.43 2007/10/29 08:38:58 vnivanch Exp $
27// GEANT4 tag $Name: geant4-09-01-patch-02 $
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//
58// Class Description:
59//
60// It is the unified Discrete process
61
62// -------------------------------------------------------------------
63//
64
65#ifndef G4VEmProcess_h
66#define G4VEmProcess_h 1
67
68#include "G4VDiscreteProcess.hh"
69#include "globals.hh"
70#include "G4Material.hh"
71#include "G4MaterialCutsCouple.hh"
72#include "G4Track.hh"
73#include "G4EmModelManager.hh"
74#include "G4UnitsTable.hh"
75#include "G4ParticleDefinition.hh"
76#include "G4ParticleChangeForGamma.hh"
77
78class G4Step;
79class G4VEmModel;
80class G4DataVector;
81class G4VParticleChange;
82class G4PhysicsTable;
83class G4PhysicsVector;
84
85//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
86
87class G4VEmProcess : public G4VDiscreteProcess
88{
89public:
90
91  G4VEmProcess(const G4String& name,
92                     G4ProcessType type = fElectromagnetic);
93
94  virtual ~G4VEmProcess();
95
96  //------------------------------------------------------------------------
97  // Virtual methods to be implemented in concrete processes
98  //------------------------------------------------------------------------
99
100  virtual G4bool IsApplicable(const G4ParticleDefinition& p) = 0;
101
102  virtual void PrintInfo() = 0;
103
104protected:
105
106  virtual void InitialiseProcess(const G4ParticleDefinition*) = 0;
107
108  //------------------------------------------------------------------------
109  // Methods with standard implementation; may be overwritten if needed
110  //------------------------------------------------------------------------
111
112  inline G4double RecalculateLambda(G4double kinEnergy,
113                                    const G4MaterialCutsCouple* couple);
114
115  //------------------------------------------------------------------------
116  // Generic methods common to all Discrete processes
117  //------------------------------------------------------------------------
118
119public:
120
121  void PrintInfoDefinition();
122
123  G4VParticleChange* PostStepDoIt(const G4Track&, const G4Step&);
124
125  void PreparePhysicsTable(const G4ParticleDefinition&);
126  // Initialise for build of tables
127
128  void BuildPhysicsTable(const G4ParticleDefinition&);
129  // Build physics table during initialisation
130
131  G4bool StorePhysicsTable(const G4ParticleDefinition*,
132                           const G4String& directory,
133                           G4bool ascii = false);
134  // Store PhysicsTable in a file.
135  // Return false in case of failure at I/O
136
137  G4bool RetrievePhysicsTable(const G4ParticleDefinition*,
138                              const G4String& directory,
139                              G4bool ascii);
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
146  //------------------------------------------------------------------------
147  // Specific methods for Discrete EM post step simulation
148  //------------------------------------------------------------------------
149
150  inline G4double MicroscopicCrossSection(G4double kineticEnergy,
151                                          const G4MaterialCutsCouple* couple);
152  // It returns the cross section of the process for energy/ material
153
154  inline G4double ComputeCrossSectionPerAtom(G4double kineticEnergy, 
155                                             G4double Z, G4double A=0., 
156                                             G4double cut=0.0);
157  // It returns the cross section of the process per atom
158
159  inline G4double MeanFreePath(const G4Track& track);
160
161  virtual G4double PostStepGetPhysicalInteractionLength(
162                             const G4Track& track,
163                             G4double   previousStepSize,
164                             G4ForceCondition* condition
165                            );
166
167  inline G4VEmModel* SelectModelForMaterial(G4double kinEnergy, 
168                                            size_t& idxRegion) const;
169
170  inline G4double GetLambda(G4double& kinEnergy, 
171                            const G4MaterialCutsCouple* couple);
172  // It returns the Lambda of the process
173
174  //------------------------------------------------------------------------
175  // Specific methods to build and access Physics Tables
176  //------------------------------------------------------------------------
177
178  inline void SetLambdaBinning(G4int nbins);
179  inline G4int LambdaBinning() const;
180  // Binning for lambda table
181
182  inline void SetMinKinEnergy(G4double e);
183  inline G4double MinKinEnergy() const;
184  // Min kinetic energy for tables
185
186  inline void SetMaxKinEnergy(G4double e);
187  inline G4double MaxKinEnergy() const;
188  // Max kinetic energy for tables
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
201  //------------------------------------------------------------------------
202
203  inline void AddEmModel(G4int, G4VEmModel*, const G4Region* region = 0);
204  // Add EM model coupled for the region
205   
206  inline void SetModel(G4VEmModel*);
207  // Assign a model to a process
208 
209  inline G4VEmModel* Model();
210  // return the assigned model
211   
212  inline void UpdateEmModel(const G4String&, G4double, G4double);
213  // Define new energy range for the model identified by the name
214
215  // Access to models
216  inline G4VEmModel* GetModelByIndex(G4int idx = 0);
217
218  //------------------------------------------------------------------------
219  // Get/set parameters used for simulation of energy loss
220  //------------------------------------------------------------------------
221
222  inline void ActivateDeexcitation(G4bool, const G4Region* r = 0);
223
224  inline void SetLambdaFactor(G4double val);
225
226  inline void SetIntegral(G4bool val);
227  inline G4bool IsIntegral() const;
228
229  inline void SetApplyCuts(G4bool val);
230 
231protected:
232
233  G4double GetMeanFreePath(const G4Track& track,
234                           G4double previousStepSize,
235                           G4ForceCondition* condition);
236
237  G4PhysicsVector* LambdaPhysicsVector(const G4MaterialCutsCouple*);
238
239  inline void SetParticle(const G4ParticleDefinition* p);
240 
241  inline void SetSecondaryParticle(const G4ParticleDefinition* p);
242
243  inline G4VEmModel* SelectModel(G4double& kinEnergy);
244
245  inline size_t CurrentMaterialCutsCoupleIndex() const;
246
247  inline G4double GetGammaEnergyCut();
248
249  inline G4double GetElectronEnergyCut();
250
251  inline void SetBuildTableFlag(G4bool val);
252
253  inline void SetStartFromNullFlag(G4bool val);
254
255private:
256
257  void Clear();
258
259  void BuildLambdaTable();
260
261  void FindLambdaMax();
262
263  inline void InitialiseStep(const G4Track&);
264
265  inline void DefineMaterial(const G4MaterialCutsCouple* couple);
266
267  inline void ComputeIntegralLambda(G4double kinEnergy);
268
269  inline G4double GetLambdaFromTable(G4double kinEnergy);
270
271  inline G4double GetCurrentLambda(G4double kinEnergy);
272
273  inline G4double ComputeCurrentLambda(G4double kinEnergy);
274
275  // hide  assignment operator
276
277  G4VEmProcess(G4VEmProcess &);
278  G4VEmProcess & operator=(const G4VEmProcess &right);
279
280// =====================================================================
281
282protected:
283
284  G4ParticleChangeForGamma     fParticleChange;
285
286private:
287
288  std::vector<G4DynamicParticle*> secParticles;
289
290  G4EmModelManager*            modelManager;
291  G4VEmModel*                  selectedModel; 
292
293  // tables and vectors
294  G4PhysicsTable*              theLambdaTable;
295  G4double*                    theEnergyOfCrossSectionMax;
296  G4double*                    theCrossSectionMax;
297
298  const G4ParticleDefinition*  particle;
299  const G4ParticleDefinition*  secondaryParticle;
300  const G4ParticleDefinition*  theGamma;
301  const G4ParticleDefinition*  theElectron;
302  const G4ParticleDefinition*  thePositron;
303
304  const std::vector<G4double>* theCuts;
305  const std::vector<G4double>* theCutsGamma;
306  const std::vector<G4double>* theCutsElectron;
307  const std::vector<G4double>* theCutsPositron;
308
309  G4int                        nLambdaBins;
310
311  G4double                     minKinEnergy;
312  G4double                     maxKinEnergy;
313  G4double                     lambdaFactor;
314
315  // cash
316  const G4Material*            currentMaterial;
317  const G4MaterialCutsCouple*  currentCouple;
318  size_t                       currentMaterialIndex;
319
320  G4double                     mfpKinEnergy;
321  G4double                     preStepKinEnergy;
322  G4double                     preStepLambda;
323
324  G4bool                       integral;
325  G4bool                       buildLambdaTable;
326  G4bool                       applyCuts;
327  G4bool                       startFromNull;
328
329  G4int                        nRegions;
330  std::vector<G4Region*>       regions;
331  std::vector<G4bool>          flagsDeexcitation;
332
333};
334
335//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
336//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
337
338inline void G4VEmProcess::DefineMaterial(const G4MaterialCutsCouple* couple)
339{
340  if(couple != currentCouple) {
341    currentCouple   = couple;
342    currentMaterial = couple->GetMaterial();
343    currentMaterialIndex = couple->GetIndex();
344    mfpKinEnergy = DBL_MAX;
345  }
346}
347
348//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
349
350inline void G4VEmProcess::InitialiseStep(const G4Track& track)
351{
352  preStepKinEnergy = track.GetKineticEnergy();
353  DefineMaterial(track.GetMaterialCutsCouple());
354  if (theNumberOfInteractionLengthLeft < 0.0) mfpKinEnergy = DBL_MAX;
355}
356
357//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
358
359inline G4double G4VEmProcess::GetLambda(G4double& kineticEnergy,
360                                  const G4MaterialCutsCouple* couple)
361{
362  DefineMaterial(couple);
363  return GetCurrentLambda(kineticEnergy);
364}
365
366//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
367
368inline G4double G4VEmProcess::GetCurrentLambda(G4double e)
369{
370  G4double x = 0.0;
371  if(theLambdaTable) x = GetLambdaFromTable(e);
372  else               x = ComputeCurrentLambda(e);
373  return x;
374}
375
376//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
377
378inline G4double G4VEmProcess::RecalculateLambda(G4double e, 
379                                                const G4MaterialCutsCouple* couple)
380{
381  DefineMaterial(couple);
382  return ComputeCurrentLambda(e);
383}
384
385//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
386
387inline G4double G4VEmProcess::ComputeCurrentLambda(G4double e)
388{
389  G4VEmModel* currentModel = SelectModel(e);
390  G4double x = 0.0;
391  if(currentModel) 
392    x = currentModel->CrossSectionPerVolume(currentMaterial,particle,
393                                            e,(*theCuts)[currentMaterialIndex]);
394  return x;
395}
396
397//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
398
399inline G4double G4VEmProcess::GetLambdaFromTable(G4double e)
400{
401  G4bool b;
402  return (((*theLambdaTable)[currentMaterialIndex])->GetValue(e, b));
403}
404
405//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
406
407inline void G4VEmProcess::ComputeIntegralLambda(G4double e)
408{
409  mfpKinEnergy  = theEnergyOfCrossSectionMax[currentMaterialIndex];
410  if (e <= mfpKinEnergy) {
411    preStepLambda = GetLambdaFromTable(e);
412
413  } else {
414    G4double e1 = e*lambdaFactor;
415    if(e1 > mfpKinEnergy) {
416      preStepLambda  = GetLambdaFromTable(e);
417      G4double preStepLambda1 = GetLambdaFromTable(e1);
418      if(preStepLambda1 > preStepLambda) {
419        mfpKinEnergy = e1;
420        preStepLambda = preStepLambda1;
421      }
422    } else {
423      preStepLambda = theCrossSectionMax[currentMaterialIndex];
424    }
425  }
426}
427
428//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
429
430inline G4double G4VEmProcess::MeanFreePath(const G4Track& track)
431{
432  DefineMaterial(track.GetMaterialCutsCouple());
433  preStepLambda = GetCurrentLambda(track.GetKineticEnergy());
434  G4double x = DBL_MAX;
435  if(DBL_MIN < preStepLambda) x = 1.0/preStepLambda;
436  return x;
437}
438
439//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
440
441inline G4VEmModel* G4VEmProcess::SelectModel(G4double& kinEnergy)
442{
443  return modelManager->SelectModel(kinEnergy, currentMaterialIndex);
444}
445
446//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
447
448inline G4VEmModel* G4VEmProcess::SelectModelForMaterial(
449                                   G4double kinEnergy, size_t& idxRegion) const
450{
451  return modelManager->SelectModel(kinEnergy, idxRegion);
452}
453
454//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
455
456inline const G4ParticleDefinition* G4VEmProcess::Particle() const
457{
458  return particle;
459}
460
461//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
462
463inline const G4ParticleDefinition* G4VEmProcess::SecondaryParticle() const
464{
465  return secondaryParticle;
466}
467
468//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
469
470inline G4double G4VEmProcess::GetGammaEnergyCut()
471{
472  return (*theCutsGamma)[currentMaterialIndex];
473}
474
475//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
476
477inline G4double G4VEmProcess::GetElectronEnergyCut()
478{
479  return (*theCutsElectron)[currentMaterialIndex];
480}
481
482//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
483
484inline void G4VEmProcess::SetLambdaFactor(G4double val)
485{
486  if(val > 0.0 && val <= 1.0) lambdaFactor = val;
487}
488
489//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
490
491inline G4VEmModel* G4VEmProcess::GetModelByIndex(G4int idx)
492{
493  return modelManager->GetModel(idx);
494}
495
496//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
497
498inline void G4VEmProcess::SetParticle(const G4ParticleDefinition* p)
499{
500  particle = p;
501}
502
503//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
504
505inline void G4VEmProcess::SetSecondaryParticle(const G4ParticleDefinition* p)
506{
507  secondaryParticle = p;
508}
509
510//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
511
512inline void G4VEmProcess::AddEmModel(G4int order, G4VEmModel* p, 
513                                     const G4Region* region)
514{
515  G4VEmFluctuationModel* fm = 0;
516  modelManager->AddEmModel(order, p, fm, region);
517  if(p) p->SetParticleChange(pParticleChange);
518}
519
520//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
521
522inline void G4VEmProcess::SetModel(G4VEmModel* model)
523{
524  selectedModel = model;
525}
526
527//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
528
529inline G4VEmModel* G4VEmProcess::Model()
530{
531  return selectedModel;
532}
533
534//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
535
536inline void G4VEmProcess::UpdateEmModel(const G4String& nam, 
537                                        G4double emin, G4double emax)
538{
539  modelManager->UpdateEmModel(nam, emin, emax);
540}
541
542//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
543
544inline G4double G4VEmProcess::ComputeCrossSectionPerAtom(
545                 G4double kineticEnergy, G4double Z, G4double A, G4double cut)
546{
547  G4VEmModel* model = SelectModel(kineticEnergy);
548  return model->ComputeCrossSectionPerAtom(particle,kineticEnergy,Z,A,cut);
549}
550
551//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
552
553inline void G4VEmProcess::SetLambdaBinning(G4int nbins)
554{
555  nLambdaBins = nbins;
556}
557
558//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
559
560inline G4int G4VEmProcess::LambdaBinning() const
561{
562  return nLambdaBins;
563}
564
565//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
566
567inline void G4VEmProcess::SetMinKinEnergy(G4double e)
568{
569  minKinEnergy = e;
570}
571
572//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
573
574inline G4double G4VEmProcess::MinKinEnergy() const
575{
576  return minKinEnergy;
577}
578
579//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
580
581inline void G4VEmProcess::SetMaxKinEnergy(G4double e)
582{
583  maxKinEnergy = e;
584}
585
586//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
587
588inline G4double G4VEmProcess::MaxKinEnergy() const
589{
590  return maxKinEnergy;
591}
592
593//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
594
595inline void G4VEmProcess::ActivateDeexcitation(G4bool, const G4Region*)
596{}
597
598//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
599
600inline const G4PhysicsTable* G4VEmProcess::LambdaTable() const
601{
602  return theLambdaTable;
603}
604
605//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
606
607inline void G4VEmProcess::SetIntegral(G4bool val)
608{
609  if(particle && particle != theGamma) integral = val;
610  if(integral) buildLambdaTable = true;
611}
612
613//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
614
615inline G4bool G4VEmProcess::IsIntegral() const
616{
617  return integral;
618}
619
620//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
621
622inline void G4VEmProcess::SetBuildTableFlag(G4bool val)
623{
624  buildLambdaTable = val;
625  if(!val) integral = false;
626}
627
628//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
629
630inline void G4VEmProcess::SetStartFromNullFlag(G4bool val)
631{
632  startFromNull = val;
633}
634
635//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
636
637inline void G4VEmProcess::SetApplyCuts(G4bool val)
638{
639  applyCuts = val;
640}
641
642//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
643
644inline size_t G4VEmProcess::CurrentMaterialCutsCoupleIndex() const 
645{
646  return currentMaterialIndex;
647}
648
649//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
650
651#endif
Note: See TracBrowser for help on using the repository browser.