Geant4  10.01
G4VEmModel.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: G4VEmModel.hh 84661 2014-10-17 14:30:16Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class header file
31 //
32 //
33 // File name: G4VEmModel
34 //
35 // Author: Vladimir Ivanchenko
36 //
37 // Creation date: 03.01.2002
38 //
39 // Modifications:
40 //
41 // 23-12-02 V.Ivanchenko change interface before move to cut per region
42 // 24-01-03 Cut per region (V.Ivanchenko)
43 // 13-02-03 Add name (V.Ivanchenko)
44 // 25-02-03 Add sample theta and displacement (V.Ivanchenko)
45 // 23-07-03 Replace G4Material by G4MaterialCutCouple in dE/dx and CrossSection
46 // calculation (V.Ivanchenko)
47 // 01-03-04 L.Urban signature changed in SampleCosineTheta
48 // 23-04-04 L.urban signature of SampleCosineTheta changed back
49 // 17-11-04 Add method CrossSectionPerAtom (V.Ivanchenko)
50 // 14-03-05 Reduce number of pure virtual methods and make inline part
51 // separate (V.Ivanchenko)
52 // 24-03-05 Remove IsInCharge and add G4VParticleChange in the constructor (VI)
53 // 08-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
54 // 15-04-05 optimize internal interface for msc (V.Ivanchenko)
55 // 08-05-05 A -> N (V.Ivanchenko)
56 // 25-07-05 Move constructor and destructor to the body (V.Ivanchenko)
57 // 02-02-06 ComputeCrossSectionPerAtom: default value A=0. (mma)
58 // 06-02-06 add method ComputeMeanFreePath() (mma)
59 // 07-03-06 Optimize msc methods (V.Ivanchenko)
60 // 29-06-06 Add member currentElement and Get/Set methods (V.Ivanchenko)
61 // 29-10-07 Added SampleScattering (V.Ivanchenko)
62 // 15-07-08 Reorder class members and improve comments (VI)
63 // 21-07-08 Added vector of G4ElementSelector and methods to use it (VI)
64 // 12-09-08 Added methods GetParticleCharge, GetChargeSquareRatio,
65 // CorrectionsAlongStep, ActivateNuclearStopping (VI)
66 // 16-02-09 Moved implementations of virtual methods to source (VI)
67 // 07-04-09 Moved msc methods from G4VEmModel to G4VMscModel (VI)
68 // 13-10-10 Added G4VEmAngularDistribution (VI)
69 //
70 // Class Description:
71 //
72 // Abstract interface to energy loss models
73 
74 // -------------------------------------------------------------------
75 //
76 
77 #ifndef G4VEmModel_h
78 #define G4VEmModel_h 1
79 
80 #include "globals.hh"
81 #include "G4DynamicParticle.hh"
82 #include "G4ParticleDefinition.hh"
83 #include "G4MaterialCutsCouple.hh"
84 #include "G4Material.hh"
85 #include "G4Element.hh"
86 #include "G4ElementVector.hh"
87 #include "G4Isotope.hh"
88 #include "G4DataVector.hh"
89 #include "G4VEmFluctuationModel.hh"
91 #include "G4EmElementSelector.hh"
92 #include <CLHEP/Random/RandomEngine.h>
93 #include <vector>
94 
95 class G4ElementData;
96 class G4PhysicsTable;
97 class G4Region;
98 class G4VParticleChange;
101 class G4Track;
102 class G4LossTableManager;
103 
105 {
106 
107 public:
108 
109  G4VEmModel(const G4String& nam);
110 
111  virtual ~G4VEmModel();
112 
113  //------------------------------------------------------------------------
114  // Virtual methods to be implemented for any concrete model
115  //------------------------------------------------------------------------
116 
117  virtual void Initialise(const G4ParticleDefinition*,
118  const G4DataVector&) = 0;
119 
120  virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*,
121  const G4MaterialCutsCouple*,
122  const G4DynamicParticle*,
123  G4double tmin = 0.0,
124  G4double tmax = DBL_MAX) = 0;
125 
126  //------------------------------------------------------------------------
127  // Methods for initialisation of MT; may be overwritten if needed
128  //------------------------------------------------------------------------
129 
130  // initilisation in local thread
131  virtual void InitialiseLocal(const G4ParticleDefinition*,
132  G4VEmModel* masterModel);
133 
134  // initilisation of a new material at run time
135  virtual void InitialiseForMaterial(const G4ParticleDefinition*,
136  const G4Material*);
137 
138  // initilisation of a new element at run time
139  virtual void InitialiseForElement(const G4ParticleDefinition*,
140  G4int Z);
141 
142  //------------------------------------------------------------------------
143  // Methods with standard implementation; may be overwritten if needed
144  //------------------------------------------------------------------------
145 
146  // main method to compute dEdx
148  const G4ParticleDefinition*,
149  G4double kineticEnergy,
150  G4double cutEnergy = DBL_MAX);
151 
152  // main method to compute cross section per Volume
154  const G4ParticleDefinition*,
155  G4double kineticEnergy,
156  G4double cutEnergy = 0.0,
157  G4double maxEnergy = DBL_MAX);
158 
159  // main method to compute cross section per atom
161  G4double kinEnergy,
162  G4double Z,
163  G4double A = 0., /* amu */
164  G4double cutEnergy = 0.0,
165  G4double maxEnergy = DBL_MAX);
166 
167  // Compute effective ion charge square
168  virtual G4double ChargeSquareRatio(const G4Track&);
169 
170  // Compute effective ion charge square
172  const G4Material*,
173  G4double kineticEnergy);
174 
175  // Compute ion charge
177  const G4Material*,
178  G4double kineticEnergy);
179 
180  // Initialisation for a new track
181  virtual void StartTracking(G4Track*);
182 
183  // add correction to energy loss and compute non-ionizing energy loss
184  virtual void CorrectionsAlongStep(const G4MaterialCutsCouple*,
185  const G4DynamicParticle*,
186  G4double& eloss,
187  G4double& niel,
188  G4double length);
189 
190  // value which may be tabulated (by default cross section)
191  virtual G4double Value(const G4MaterialCutsCouple*,
192  const G4ParticleDefinition*,
193  G4double kineticEnergy);
194 
195  // threshold for zero value
196  virtual G4double MinPrimaryEnergy(const G4Material*,
197  const G4ParticleDefinition*,
198  G4double cut = 0.0);
199 
200  // model can define low-energy limit for the cut
202  const G4MaterialCutsCouple*);
203 
204  // initilisation at run time for a given material
205  virtual void SetupForMaterial(const G4ParticleDefinition*,
206  const G4Material*,
207  G4double kineticEnergy);
208 
209  // add a region for the model
210  virtual void DefineForRegion(const G4Region*);
211 
212 protected:
213 
214  // initialisation of the ParticleChange for the model
216 
217  // initialisation of the ParticleChange for the model
219 
220  // kinematically allowed max kinetic energy of a secondary
222  G4double kineticEnergy);
223 
224 public:
225 
226  //------------------------------------------------------------------------
227  // Generic methods common to all models
228  //------------------------------------------------------------------------
229 
230  // should be called at initialisation to build element selectors
232  const G4DataVector&);
233 
234  // should be called at initialisation to access element selectors
235  inline std::vector<G4EmElementSelector*>* GetElementSelectors();
236 
237  // should be called at initialisation to set element selectors
238  inline void SetElementSelectors(std::vector<G4EmElementSelector*>*);
239 
240  // dEdx per unit length
242  const G4ParticleDefinition*,
243  G4double kineticEnergy,
244  G4double cutEnergy = DBL_MAX);
245 
246  // cross section per volume
248  const G4ParticleDefinition*,
249  G4double kineticEnergy,
250  G4double cutEnergy = 0.0,
251  G4double maxEnergy = DBL_MAX);
252 
253  // compute mean free path via cross section per volume
255  G4double kineticEnergy,
256  const G4Material*,
257  G4double cutEnergy = 0.0,
258  G4double maxEnergy = DBL_MAX);
259 
260  // generic cross section per element
262  const G4Element*,
263  G4double kinEnergy,
264  G4double cutEnergy = 0.0,
265  G4double maxEnergy = DBL_MAX);
266 
267  // select isotope in order to have precise mass of the nucleus
268  inline G4int SelectIsotopeNumber(const G4Element*);
269 
270  // atom can be selected effitiantly if element selectors are initialised
271  inline const G4Element* SelectRandomAtom(const G4MaterialCutsCouple*,
272  const G4ParticleDefinition*,
273  G4double kineticEnergy,
274  G4double cutEnergy = 0.0,
275  G4double maxEnergy = DBL_MAX);
276 
277  // to select atom cross section per volume is recomputed for each element
278  const G4Element* SelectRandomAtom(const G4Material*,
279  const G4ParticleDefinition*,
280  G4double kineticEnergy,
281  G4double cutEnergy = 0.0,
282  G4double maxEnergy = DBL_MAX);
283 
284  // to select atom if cross section is proportional number of electrons
285  inline G4int SelectRandomAtomNumber(const G4Material*);
286 
287  //------------------------------------------------------------------------
288  // Get/Set methods
289  //------------------------------------------------------------------------
290 
292 
294 
295  inline G4ElementData* GetElementData();
296 
298 
300 
302 
304 
305  inline G4double HighEnergyLimit() const;
306 
307  inline G4double LowEnergyLimit() const;
308 
309  inline G4double HighEnergyActivationLimit() const;
310 
311  inline G4double LowEnergyActivationLimit() const;
312 
313  inline G4double PolarAngleLimit() const;
314 
315  inline G4double SecondaryThreshold() const;
316 
317  inline G4bool LPMFlag() const;
318 
319  inline G4bool DeexcitationFlag() const;
320 
321  inline G4bool ForceBuildTableFlag() const;
322 
323  inline G4bool UseAngularGeneratorFlag() const;
324 
325  inline void SetAngularGeneratorFlag(G4bool);
326 
327  inline void SetHighEnergyLimit(G4double);
328 
329  inline void SetLowEnergyLimit(G4double);
330 
332 
334 
335  inline G4bool IsActive(G4double kinEnergy);
336 
337  inline void SetPolarAngleLimit(G4double);
338 
339  inline void SetSecondaryThreshold(G4double);
340 
341  inline void SetLPMFlag(G4bool val);
342 
343  inline void SetDeexcitationFlag(G4bool val);
344 
345  inline void SetForceBuildTable(G4bool val);
346 
347  inline void SetMasterThread(G4bool val);
348 
349  inline G4bool IsMaster() const;
350 
351  inline G4double MaxSecondaryKinEnergy(const G4DynamicParticle* dynParticle);
352 
353  inline const G4String& GetName() const;
354 
355  inline void SetCurrentCouple(const G4MaterialCutsCouple*);
356 
357  inline const G4Element* GetCurrentElement() const;
358 
359  inline const G4Isotope* GetCurrentIsotope() const;
360 
361 protected:
362 
363  inline const G4MaterialCutsCouple* CurrentCouple() const;
364 
365  inline void SetCurrentElement(const G4Element*);
366 
367 private:
368 
369  // hide assignment operator
371  G4VEmModel(const G4VEmModel&);
372 
373  // ======== Parameters of the class fixed at construction =========
374 
377  const G4String name;
378 
379  // ======== Parameters of the class fixed at initialisation =======
380 
391 
396  std::vector<G4EmElementSelector*>* elmSelectors;
397 
398 protected:
399 
400  CLHEP::HepRandomEngine* rndmEngineMod;
401 
405  const std::vector<G4double>* theDensityFactor;
406  const std::vector<G4int>* theDensityIdx;
407  size_t idxTable;
408  const static G4double inveplus;
409 
410  // ======== Cashed values - may be state dependent ================
411 
412 private:
413 
418 
420  std::vector<G4double> xsec;
421 
422 };
423 
424 // ======== Run time inline methods ================
425 
427 {
428  fCurrentCouple = p;
429 }
430 
431 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
432 
434 {
435  return fCurrentCouple;
436 }
437 
438 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
439 
441 {
442  fCurrentElement = elm;
443  fCurrentIsotope = 0;
444 }
445 
446 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
447 
449 {
450  return fCurrentElement;
451 }
452 
453 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
454 
456 {
457  return fCurrentIsotope;
458 }
459 
460 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
461 
462 inline
464 {
465  return MaxSecondaryEnergy(dynPart->GetParticleDefinition(),
466  dynPart->GetKineticEnergy());
467 }
468 
469 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
470 
472  const G4ParticleDefinition* part,
473  G4double kinEnergy,
474  G4double cutEnergy)
475 {
476  SetCurrentCouple(couple);
477  return ComputeDEDXPerVolume(couple->GetMaterial(),part,kinEnergy,cutEnergy);
478 }
479 
480 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
481 
483  const G4ParticleDefinition* part,
484  G4double kinEnergy,
485  G4double cutEnergy,
486  G4double maxEnergy)
487 {
488  SetCurrentCouple(couple);
489  return CrossSectionPerVolume(couple->GetMaterial(),part,kinEnergy,
490  cutEnergy,maxEnergy);
491 }
492 
493 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
494 
495 inline
497  G4double ekin,
498  const G4Material* material,
499  G4double emin,
500  G4double emax)
501 {
502  G4double mfp = DBL_MAX;
503  G4double cross = CrossSectionPerVolume(material,part,ekin,emin,emax);
504  if (cross > DBL_MIN) { mfp = 1./cross; }
505  return mfp;
506 }
507 
508 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
509 
511  const G4ParticleDefinition* part,
512  const G4Element* elm,
513  G4double kinEnergy,
514  G4double cutEnergy,
515  G4double maxEnergy)
516 {
517  SetCurrentElement(elm);
518  return ComputeCrossSectionPerAtom(part,kinEnergy,elm->GetZ(),elm->GetN(),
519  cutEnergy,maxEnergy);
520 }
521 
522 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
523 
524 inline const G4Element*
526  const G4ParticleDefinition* part,
527  G4double kinEnergy,
528  G4double cutEnergy,
529  G4double maxEnergy)
530 {
531  fCurrentCouple = couple;
532  if(nSelectors > 0) {
533  fCurrentElement =
534  ((*elmSelectors)[couple->GetIndex()])->SelectRandomAtom(kinEnergy);
535  } else {
536  fCurrentElement = SelectRandomAtom(couple->GetMaterial(),part,kinEnergy,
537  cutEnergy,maxEnergy);
538  }
539  fCurrentIsotope = 0;
540  return fCurrentElement;
541 }
542 
543 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
544 
546 {
547  // this algorith assumes that cross section is proportional to
548  // number electrons multiplied by number of atoms
549  size_t nn = mat->GetNumberOfElements();
550  const G4ElementVector* elmv = mat->GetElementVector();
551  G4int Z = G4lrint((*elmv)[0]->GetZ());
552  if(1 < nn) {
553  const G4double* at = mat->GetVecNbOfAtomsPerVolume();
554  G4double tot = mat->GetTotNbOfAtomsPerVolume()*rndmEngineMod->flat();
555  for( size_t i=0; i<nn; ++i) {
556  Z = G4lrint((*elmv)[0]->GetZ());
557  tot -= Z*at[i];
558  if(tot <= 0.0) { break; }
559  }
560  }
561  return Z;
562 }
563 
564 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
565 
567 {
568  SetCurrentElement(elm);
569  G4int N = G4lrint(elm->GetN());
570  G4int ni = elm->GetNumberOfIsotopes();
571  fCurrentIsotope = 0;
572  if(ni > 0) {
573  G4int idx = 0;
574  if(ni > 1) {
576  G4double x = rndmEngineMod->flat();
577  for(; idx<ni; ++idx) {
578  x -= ab[idx];
579  if (x <= 0.0) { break; }
580  }
581  if(idx >= ni) { idx = ni - 1; }
582  }
583  fCurrentIsotope = elm->GetIsotope(idx);
584  N = fCurrentIsotope->GetN();
585  }
586  return N;
587 }
588 
589 // ======== Get/Set inline methods used at initialisation ================
590 
592 {
593  return flucModel;
594 }
595 
596 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
597 
599 {
600  return anglModel;
601 }
602 
603 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
604 
606 {
607  if(p != anglModel) {
608  delete anglModel;
609  anglModel = p;
610  }
611 }
612 
613 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
614 
616 {
617  return highLimit;
618 }
619 
620 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
621 
623 {
624  return lowLimit;
625 }
626 
627 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
628 
630 {
631  return eMaxActive;
632 }
633 
634 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
635 
637 {
638  return eMinActive;
639 }
640 
641 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
642 
644 {
645  return polarAngleLimit;
646 }
647 
648 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
649 
651 {
652  return secondaryThreshold;
653 }
654 
655 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
656 
658 {
659  return theLPMflag;
660 }
661 
662 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
663 
665 {
666  return flagDeexcitation;
667 }
668 
669 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
670 
672 {
673  return flagForceBuildTable;
674 }
675 
676 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
677 
679 {
680  return useAngularGenerator;
681 }
682 
683 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
684 
686 {
687  useAngularGenerator = val;
688 }
689 
690 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
691 
693 {
694  isMaster = val;
695 }
696 
697 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
698 
700 {
701  return isMaster;
702 }
703 
704 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
705 
707 {
708  highLimit = val;
709 }
710 
711 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
712 
714 {
715  lowLimit = val;
716 }
717 
718 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
719 
721 {
722  eMaxActive = val;
723 }
724 
725 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
726 
728 {
729  eMinActive = val;
730 }
731 
732 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
733 
735 {
736  return (kinEnergy >= eMinActive && kinEnergy <= eMaxActive);
737 }
738 
739 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
740 
742 {
743  polarAngleLimit = val;
744 }
745 
746 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
747 
749 {
750  secondaryThreshold = val;
751 }
752 
753 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
754 
756 {
757  theLPMflag = val;
758 }
759 
760 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
761 
763 {
764  flagDeexcitation = val;
765 }
766 
767 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
768 
770 {
771  flagForceBuildTable = val;
772 }
773 
774 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
775 
776 inline const G4String& G4VEmModel::GetName() const
777 {
778  return name;
779 }
780 
781 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
782 
783 inline std::vector<G4EmElementSelector*>* G4VEmModel::GetElementSelectors()
784 {
785  return elmSelectors;
786 }
787 
788 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
789 
790 inline void
791 G4VEmModel::SetElementSelectors(std::vector<G4EmElementSelector*>* p)
792 {
793  elmSelectors = p;
794  if(elmSelectors) { nSelectors = elmSelectors->size(); }
795  localElmSelectors = false;
796 }
797 
798 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
799 
801 {
802  return fElementData;
803 }
804 
805 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
806 
808 {
809  return xSectionTable;
810 }
811 
812 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
813 
814 #endif
815 
size_t GetNumberOfIsotopes() const
Definition: G4Element.hh:158
G4double LowEnergyActivationLimit() const
Definition: G4VEmModel.hh:636
G4double HighEnergyActivationLimit() const
Definition: G4VEmModel.hh:629
G4bool localTable
Definition: G4VEmModel.hh:392
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:463
size_t idxTable
Definition: G4VEmModel.hh:407
void SetActivationHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:720
G4bool flagForceBuildTable
Definition: G4VEmModel.hh:389
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:622
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:271
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
Definition: G4VEmModel.cc:257
virtual void StartTracking(G4Track *)
Definition: G4VEmModel.cc:297
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:133
G4double SecondaryThreshold() const
Definition: G4VEmModel.hh:650
std::vector< G4Element * > G4ElementVector
G4bool ForceBuildTableFlag() const
Definition: G4VEmModel.hh:671
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
Definition: G4VEmModel.cc:365
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
G4double GetKineticEnergy() const
void InitialiseElementSelectors(const G4ParticleDefinition *, const G4DataVector &)
Definition: G4VEmModel.cc:161
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:615
G4bool flagDeexcitation
Definition: G4VEmModel.hh:388
G4double GetN() const
Definition: G4Element.hh:134
G4bool theLPMflag
Definition: G4VEmModel.hh:387
virtual void DefineForRegion(const G4Region *)
Definition: G4VEmModel.cc:335
G4bool isMaster
Definition: G4VEmModel.hh:390
virtual G4double MinPrimaryEnergy(const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
Definition: G4VEmModel.cc:381
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:232
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:406
G4double GetZ() const
Definition: G4Element.hh:131
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:598
void SetSecondaryThreshold(G4double)
Definition: G4VEmModel.hh:748
G4bool LPMFlag() const
Definition: G4VEmModel.hh:657
G4bool IsMaster() const
Definition: G4VEmModel.hh:699
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:807
const G4String name
Definition: G4VEmModel.hh:377
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:591
virtual G4double MinEnergyCut(const G4ParticleDefinition *, const G4MaterialCutsCouple *)
Definition: G4VEmModel.cc:390
G4ElementData * GetElementData()
Definition: G4VEmModel.hh:800
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:188
int G4int
Definition: G4Types.hh:78
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:416
G4VEmFluctuationModel * flucModel
Definition: G4VEmModel.hh:375
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:76
G4int nSelectors
Definition: G4VEmModel.hh:395
G4LossTableManager * fManager
Definition: G4VEmModel.hh:414
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:706
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)=0
void SetAngularGeneratorFlag(G4bool)
Definition: G4VEmModel.hh:685
const std::vector< G4int > * theDensityIdx
Definition: G4VEmModel.hh:406
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:204
G4double eMaxActive
Definition: G4VEmModel.hh:384
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:433
virtual void InitialiseForMaterial(const G4ParticleDefinition *, const G4Material *)
Definition: G4VEmModel.cc:238
G4int GetN() const
Definition: G4Isotope.hh:94
G4double secondaryThreshold
Definition: G4VEmModel.hh:386
CLHEP::HepRandomEngine * rndmEngineMod
Definition: G4VEmModel.hh:400
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=0)
Definition: G4VEmModel.cc:413
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:415
G4double polarAngleLimit
Definition: G4VEmModel.hh:385
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:678
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:396
bool G4bool
Definition: G4Types.hh:79
G4VEmModel & operator=(const G4VEmModel &right)
G4double CrossSection(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:482
static const G4double inveplus
Definition: G4VEmModel.hh:408
G4double * GetRelativeAbundanceVector() const
Definition: G4Element.hh:166
virtual ~G4VEmModel()
Definition: G4VEmModel.cc:106
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:376
void SetCrossSectionTable(G4PhysicsTable *, G4bool isLocal)
Definition: G4VEmModel.cc:421
G4int SelectIsotopeNumber(const G4Element *)
Definition: G4VEmModel.hh:566
const G4ParticleDefinition * GetParticleDefinition() const
std::vector< G4EmElementSelector * > * GetElementSelectors()
Definition: G4VEmModel.hh:783
static const G4double A[nN]
G4int nsec
Definition: G4VEmModel.hh:419
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:326
virtual G4double ChargeSquareRatio(const G4Track &)
Definition: G4VEmModel.cc:340
const std::vector< G4double > * theDensityFactor
Definition: G4VEmModel.hh:405
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:734
void SetActivationLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:727
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:692
G4double eMinActive
Definition: G4VEmModel.hh:383
std::vector< G4double > xsec
Definition: G4VEmModel.hh:420
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:207
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:348
virtual G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:357
G4bool DeexcitationFlag() const
Definition: G4VEmModel.hh:664
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.cc:262
void SetElementSelectors(std::vector< G4EmElementSelector * > *)
Definition: G4VEmModel.hh:791
static const G4double ab
G4double ComputeDEDX(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.hh:471
int G4lrint(double ad)
Definition: templates.hh:163
G4bool localElmSelectors
Definition: G4VEmModel.hh:393
void SetLPMFlag(G4bool val)
Definition: G4VEmModel.hh:755
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:426
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:605
G4double highLimit
Definition: G4VEmModel.hh:382
G4bool useAngularGenerator
Definition: G4VEmModel.hh:394
G4double PolarAngleLimit() const
Definition: G4VEmModel.hh:643
#define DBL_MIN
Definition: templates.hh:75
G4double lowLimit
Definition: G4VEmModel.hh:381
void SetForceBuildTable(G4bool val)
Definition: G4VEmModel.hh:769
const G4Isotope * GetIsotope(G4int iso) const
Definition: G4Element.hh:169
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:403
const G4String & GetName() const
Definition: G4VEmModel.hh:776
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:398
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
G4PhysicsTable * xSectionTable
Definition: G4VEmModel.hh:404
double G4double
Definition: G4Types.hh:76
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:713
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:417
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:762
G4double ComputeMeanFreePath(const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:496
void SetCurrentElement(const G4Element *)
Definition: G4VEmModel.hh:440
#define DBL_MAX
Definition: templates.hh:83
G4ElementData * fElementData
Definition: G4VEmModel.hh:402
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:525
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:741
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:147
virtual G4double Value(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:372
const G4Material * GetMaterial() const
G4int SelectRandomAtomNumber(const G4Material *)
Definition: G4VEmModel.hh:545
const G4Isotope * GetCurrentIsotope() const
Definition: G4VEmModel.hh:455
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:448