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

Last change on this file since 819 was 819, checked in by garnier, 17 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.