Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4VEmModel.hh
이 파일의 문서화 페이지로 가기
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 108386 2018-02-09 15:38:32Z 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"
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  explicit 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*, const G4DataVector&) = 0;
118 
119  virtual void SampleSecondaries(std::vector<G4DynamicParticle*>*,
120  const G4MaterialCutsCouple*,
121  const G4DynamicParticle*,
122  G4double tmin = 0.0,
123  G4double tmax = DBL_MAX) = 0;
124 
125  //------------------------------------------------------------------------
126  // Methods for initialisation of MT; may be overwritten if needed
127  //------------------------------------------------------------------------
128 
129  // initilisation in local thread
130  virtual void InitialiseLocal(const G4ParticleDefinition*,
131  G4VEmModel* masterModel);
132 
133  // initilisation of a new material at run time
134  virtual void InitialiseForMaterial(const G4ParticleDefinition*,
135  const G4Material*);
136 
137  // initilisation of a new element at run time
138  virtual void InitialiseForElement(const G4ParticleDefinition*,
139  G4int Z);
140 
141  //------------------------------------------------------------------------
142  // Methods with standard implementation; may be overwritten if needed
143  //------------------------------------------------------------------------
144 
145  // main method to compute dEdx
147  const G4ParticleDefinition*,
149  G4double cutEnergy = DBL_MAX);
150 
151  // main method to compute cross section per Volume
153  const G4ParticleDefinition*,
154  G4double kineticEnergy,
155  G4double cutEnergy = 0.0,
156  G4double maxEnergy = DBL_MAX);
157 
158  // method to get partial cross section
160  G4int level,
161  const G4ParticleDefinition*,
162  G4double kineticEnergy);
163 
164  // main method to compute cross section per atom
166  G4double kinEnergy,
167  G4double Z,
168  G4double A = 0., /* amu */
169  G4double cutEnergy = 0.0,
170  G4double maxEnergy = DBL_MAX);
171 
172  // main method to compute cross section per atomic shell
174  G4int Z, G4int shellIdx,
175  G4double kinEnergy,
176  G4double cutEnergy = 0.0,
177  G4double maxEnergy = DBL_MAX);
178 
179  // Compute effective ion charge square
180  virtual G4double ChargeSquareRatio(const G4Track&);
181 
182  // Compute effective ion charge square
184  const G4Material*,
185  G4double kineticEnergy);
186 
187  // Compute ion charge
189  const G4Material*,
190  G4double kineticEnergy);
191 
192  // Initialisation for a new track
193  virtual void StartTracking(G4Track*);
194 
195  // add correction to energy loss and compute non-ionizing energy loss
196  virtual void CorrectionsAlongStep(const G4MaterialCutsCouple*,
197  const G4DynamicParticle*,
198  G4double& eloss,
199  G4double& niel,
200  G4double length);
201 
202  // value which may be tabulated (by default cross section)
203  virtual G4double Value(const G4MaterialCutsCouple*,
204  const G4ParticleDefinition*,
205  G4double kineticEnergy);
206 
207  // threshold for zero value
208  virtual G4double MinPrimaryEnergy(const G4Material*,
209  const G4ParticleDefinition*,
210  G4double cut = 0.0);
211 
212  // model can define low-energy limit for the cut
214  const G4MaterialCutsCouple*);
215 
216  // initilisation at run time for a given material
217  virtual void SetupForMaterial(const G4ParticleDefinition*,
218  const G4Material*,
219  G4double kineticEnergy);
220 
221  // add a region for the model
222  virtual void DefineForRegion(const G4Region*);
223 
224  // for automatic documentation
225  virtual void ModelDescription(std::ostream& outFile) const;
226 
227  virtual void ModelDescription(std::ostream& outFile,
228  G4String endOfLine) const;
229 
230 protected:
231 
232  // initialisation of the ParticleChange for the model
234 
235  // initialisation of the ParticleChange for the model
237 
238  // kinematically allowed max kinetic energy of a secondary
240  G4double kineticEnergy);
241 
242 public:
243 
244  //------------------------------------------------------------------------
245  // Generic methods common to all models
246  //------------------------------------------------------------------------
247 
248  // should be called at initialisation to build element selectors
250  const G4DataVector&);
251 
252  // should be called at initialisation to access element selectors
253  inline std::vector<G4EmElementSelector*>* GetElementSelectors();
254 
255  // should be called at initialisation to set element selectors
256  inline void SetElementSelectors(std::vector<G4EmElementSelector*>*);
257 
258  // dEdx per unit length
259  virtual inline G4double ComputeDEDX(const G4MaterialCutsCouple*,
260  const G4ParticleDefinition*,
261  G4double kineticEnergy,
262  G4double cutEnergy = DBL_MAX);
263 
264  // cross section per volume
266  const G4ParticleDefinition*,
267  G4double kineticEnergy,
268  G4double cutEnergy = 0.0,
269  G4double maxEnergy = DBL_MAX);
270 
271  // compute mean free path via cross section per volume
273  G4double kineticEnergy,
274  const G4Material*,
275  G4double cutEnergy = 0.0,
276  G4double maxEnergy = DBL_MAX);
277 
278  // generic cross section per element
280  const G4Element*,
281  G4double kinEnergy,
282  G4double cutEnergy = 0.0,
283  G4double maxEnergy = DBL_MAX);
284 
285  // atom can be selected effitiantly if element selectors are initialised
286  inline const G4Element* SelectRandomAtom(const G4MaterialCutsCouple*,
287  const G4ParticleDefinition*,
288  G4double kineticEnergy,
289  G4double cutEnergy = 0.0,
290  G4double maxEnergy = DBL_MAX);
291 
292  // to select atom cross section per volume is recomputed for each element
293  const G4Element* SelectRandomAtom(const G4Material*,
294  const G4ParticleDefinition*,
295  G4double kineticEnergy,
296  G4double cutEnergy = 0.0,
297  G4double maxEnergy = DBL_MAX);
298 
299  // to select atom if cross section is proportional number of electrons
301 
302  // select isotope in order to have precise mass of the nucleus
304 
305  //------------------------------------------------------------------------
306  // Get/Set methods
307  //------------------------------------------------------------------------
308 
310 
312 
313  inline G4ElementData* GetElementData();
314 
316 
318 
320 
321  inline G4VEmModel* GetTripletModel();
322 
323  inline void SetTripletModel(G4VEmModel*);
324 
326 
327  inline G4double HighEnergyLimit() const;
328 
329  inline G4double LowEnergyLimit() const;
330 
331  inline G4double HighEnergyActivationLimit() const;
332 
333  inline G4double LowEnergyActivationLimit() const;
334 
335  inline G4double PolarAngleLimit() const;
336 
337  inline G4double SecondaryThreshold() const;
338 
339  inline G4bool LPMFlag() const;
340 
341  inline G4bool DeexcitationFlag() const;
342 
343  inline G4bool ForceBuildTableFlag() const;
344 
345  inline G4bool UseAngularGeneratorFlag() const;
346 
347  inline void SetAngularGeneratorFlag(G4bool);
348 
349  inline void SetHighEnergyLimit(G4double);
350 
351  inline void SetLowEnergyLimit(G4double);
352 
354 
356 
357  inline G4bool IsActive(G4double kinEnergy);
358 
359  inline void SetPolarAngleLimit(G4double);
360 
361  inline void SetSecondaryThreshold(G4double);
362 
363  inline void SetLPMFlag(G4bool val);
364 
365  inline void SetDeexcitationFlag(G4bool val);
366 
367  inline void SetForceBuildTable(G4bool val);
368 
369  inline void SetFluctuationFlag(G4bool val);
370 
371  inline void SetMasterThread(G4bool val);
372 
373  inline G4bool IsMaster() const;
374 
375  inline G4double MaxSecondaryKinEnergy(const G4DynamicParticle* dynParticle);
376 
377  inline const G4String& GetName() const;
378 
379  inline void SetCurrentCouple(const G4MaterialCutsCouple*);
380 
381  inline const G4Element* GetCurrentElement() const;
382 
383  inline const G4Isotope* GetCurrentIsotope() const;
384 
385  inline G4bool IsLocked() const;
386 
387  inline void SetLocked(G4bool);
388 
389 protected:
390 
391  inline const G4MaterialCutsCouple* CurrentCouple() const;
392 
393  inline void SetCurrentElement(const G4Element*);
394 
395 private:
396 
397  // hide assignment operator
398  G4VEmModel & operator=(const G4VEmModel &right) = delete;
399  G4VEmModel(const G4VEmModel&) = delete;
400 
401  // ======== Parameters of the class fixed at construction =========
402 
405  const G4String name;
406 
407  // ======== Parameters of the class fixed at initialisation =======
408 
419 
425  std::vector<G4EmElementSelector*>* elmSelectors;
427 
428 protected:
429 
433  const std::vector<G4double>* theDensityFactor;
434  const std::vector<G4int>* theDensityIdx;
435  size_t idxTable;
438 
439  // ======== Cached values - may be state dependent ================
440 
441 private:
442 
447 
449  std::vector<G4double> xsec;
450 
451 };
452 
453 // ======== Run time inline methods ================
454 
456 {
457  fCurrentCouple = p;
458 }
459 
460 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
461 
463 {
464  return fCurrentCouple;
465 }
466 
467 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
468 
470 {
471  fCurrentElement = elm;
472  fCurrentIsotope = nullptr;
473 }
474 
475 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
476 
478 {
479  return fCurrentElement;
480 }
481 
482 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
483 
485 {
486  return fCurrentIsotope;
487 }
488 
489 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
490 
491 inline
493 {
494  return MaxSecondaryEnergy(dynPart->GetParticleDefinition(),
495  dynPart->GetKineticEnergy());
496 }
497 
498 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
499 
501  const G4ParticleDefinition* part,
502  G4double kinEnergy,
503  G4double cutEnergy)
504 {
505  SetCurrentCouple(couple);
506  return ComputeDEDXPerVolume(couple->GetMaterial(),part,kinEnergy,cutEnergy);
507 }
508 
509 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
510 
512  const G4ParticleDefinition* part,
513  G4double kinEnergy,
514  G4double cutEnergy,
515  G4double maxEnergy)
516 {
517  SetCurrentCouple(couple);
518  return CrossSectionPerVolume(couple->GetMaterial(),part,kinEnergy,
519  cutEnergy,maxEnergy);
520 }
521 
522 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
523 
524 inline
526  G4double ekin,
527  const G4Material* material,
528  G4double emin,
529  G4double emax)
530 {
531  G4double cross = CrossSectionPerVolume(material,part,ekin,emin,emax);
532  return cross > 0.0 ? 1./cross : DBL_MAX;
533 }
534 
535 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
536 
537 inline G4double
539  const G4Element* elm,
540  G4double kinEnergy,
541  G4double cutEnergy,
542  G4double maxEnergy)
543 {
544  SetCurrentElement(elm);
545  return ComputeCrossSectionPerAtom(part,kinEnergy,elm->GetZ(),elm->GetN(),
546  cutEnergy,maxEnergy);
547 }
548 
549 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
550 
551 inline const G4Element*
553  const G4ParticleDefinition* part,
554  G4double kinEnergy,
555  G4double cutEnergy,
556  G4double maxEnergy)
557 {
558  fCurrentCouple = couple;
559  fCurrentElement = (nSelectors > 0) ?
560  ((*elmSelectors)[couple->GetIndex()])->SelectRandomAtom(kinEnergy) :
561  SelectRandomAtom(couple->GetMaterial(),part,kinEnergy,cutEnergy,maxEnergy);
562  fCurrentIsotope = nullptr;
563  return fCurrentElement;
564 }
565 
566 // ======== Get/Set inline methods used at initialisation ================
567 
569 {
570  return flucModel;
571 }
572 
573 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
574 
576 {
577  return anglModel;
578 }
579 
580 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
581 
583 {
584  if(p != anglModel) {
585  delete anglModel;
586  anglModel = p;
587  }
588 }
589 
590 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
591 
593 {
594  return fTripletModel;
595 }
596 
597 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
598 
600 {
601  if(p != fTripletModel) {
602  delete fTripletModel;
603  fTripletModel = p;
604  }
605 }
606 
607 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
608 
610 {
611  return highLimit;
612 }
613 
614 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
615 
617 {
618  return lowLimit;
619 }
620 
621 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
622 
624 {
625  return eMaxActive;
626 }
627 
628 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
629 
631 {
632  return eMinActive;
633 }
634 
635 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
636 
638 {
639  return polarAngleLimit;
640 }
641 
642 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
643 
645 {
646  return secondaryThreshold;
647 }
648 
649 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
650 
652 {
653  return theLPMflag;
654 }
655 
656 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
657 
659 {
660  return flagDeexcitation;
661 }
662 
663 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
664 
666 {
667  return flagForceBuildTable;
668 }
669 
670 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
671 
673 {
674  return useAngularGenerator;
675 }
676 
677 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
678 
680 {
681  useAngularGenerator = val;
682 }
683 
684 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
685 
687 {
688  lossFlucFlag = val;
689 }
690 
691 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
692 
694 {
695  isMaster = val;
696 }
697 
698 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
699 
701 {
702  return isMaster;
703 }
704 
705 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
706 
708 {
709  highLimit = val;
710 }
711 
712 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
713 
715 {
716  lowLimit = val;
717 }
718 
719 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
720 
722 {
723  eMaxActive = val;
724 }
725 
726 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
727 
729 {
730  eMinActive = val;
731 }
732 
733 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
734 
736 {
737  return (kinEnergy >= eMinActive && kinEnergy <= eMaxActive);
738 }
739 
740 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
741 
743 {
744  if(!isLocked) { polarAngleLimit = val; }
745 }
746 
747 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
748 
750 {
751  secondaryThreshold = val;
752 }
753 
754 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
755 
757 {
758  theLPMflag = val;
759 }
760 
761 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
762 
764 {
765  flagDeexcitation = val;
766 }
767 
768 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
769 
771 {
772  flagForceBuildTable = val;
773 }
774 
775 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
776 
777 inline const G4String& G4VEmModel::GetName() const
778 {
779  return name;
780 }
781 
782 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
783 
784 inline std::vector<G4EmElementSelector*>* G4VEmModel::GetElementSelectors()
785 {
786  return elmSelectors;
787 }
788 
789 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
790 
791 inline void
792 G4VEmModel::SetElementSelectors(std::vector<G4EmElementSelector*>* p)
793 {
794  elmSelectors = p;
795  nSelectors = (elmSelectors) ? elmSelectors->size() : 0;
796  localElmSelectors = false;
797 }
798 
799 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
800 
802 {
803  return fElementData;
804 }
805 
806 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
807 
809 {
810  return xSectionTable;
811 }
812 
813 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
814 
816 {
817  return isLocked;
818 }
819 
820 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
821 
823 {
824  isLocked = val;
825 }
826 
827 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
828 
829 #endif
830 
void SetCrossSectionTable(G4PhysicsTable *, G4bool isLocal)
Definition: G4VEmModel.cc:462
void SetAngularGeneratorFlag(G4bool)
Definition: G4VEmModel.hh:679
void SetElementSelectors(std::vector< G4EmElementSelector * > *)
Definition: G4VEmModel.hh:792
virtual G4double GetPartialCrossSection(const G4Material *, G4int level, const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:277
G4double lowLimit
Definition: G4VEmModel.hh:409
G4bool IsActive(G4double kinEnergy)
Definition: G4VEmModel.hh:735
virtual void SetupForMaterial(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:447
G4int SelectRandomAtomNumber(const G4Material *)
Definition: G4VEmModel.cc:313
G4double highLimit
Definition: G4VEmModel.hh:410
G4double eMinActive
Definition: G4VEmModel.hh:411
G4VEmModel * fTripletModel
Definition: G4VEmModel.hh:446
G4bool IsLocked() const
Definition: G4VEmModel.hh:815
const G4MaterialCutsCouple * fCurrentCouple
Definition: G4VEmModel.hh:443
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:707
G4bool ForceBuildTableFlag() const
Definition: G4VEmModel.hh:665
const G4Isotope * fCurrentIsotope
Definition: G4VEmModel.hh:445
G4double MaxSecondaryKinEnergy(const G4DynamicParticle *dynParticle)
Definition: G4VEmModel.hh:492
G4VEmAngularDistribution * anglModel
Definition: G4VEmModel.hh:404
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:389
void SetLocked(G4bool)
Definition: G4VEmModel.hh:822
G4bool isMaster
Definition: G4VEmModel.hh:418
G4int SelectIsotopeNumber(const G4Element *)
Definition: G4VEmModel.cc:335
const char * p
Definition: xmltok.h:285
virtual ~G4VEmModel()
Definition: G4VEmModel.cc:92
G4double eMaxActive
Definition: G4VEmModel.hh:412
void SetActivationHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:721
G4ElementData * fElementData
Definition: G4VEmModel.hh:430
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:552
static const G4double emax
std::vector< G4EmElementSelector * > * GetElementSelectors()
Definition: G4VEmModel.hh:784
void SetMasterThread(G4bool val)
Definition: G4VEmModel.hh:693
void SetCurrentCouple(const G4MaterialCutsCouple *)
Definition: G4VEmModel.hh:455
virtual void DefineForRegion(const G4Region *)
Definition: G4VEmModel.cc:376
G4double LowEnergyActivationLimit() const
Definition: G4VEmModel.hh:630
void SetCurrentElement(const G4Element *)
Definition: G4VEmModel.hh:469
G4double LowEnergyLimit() const
Definition: G4VEmModel.hh:616
G4double CrossSection(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:511
G4double GetN() const
Definition: G4Element.hh:135
Float_t Z
G4double secondaryThreshold
Definition: G4VEmModel.hh:414
void SetLPMFlag(G4bool val)
Definition: G4VEmModel.hh:756
G4ParticleChangeForLoss * GetParticleChangeForLoss()
Definition: G4VEmModel.cc:116
void SetDeexcitationFlag(G4bool val)
Definition: G4VEmModel.hh:763
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
TString part[npart]
Definition: Style.C:32
std::vector< G4double > xsec
Definition: G4VEmModel.hh:449
G4PhysicsTable * GetCrossSectionTable()
Definition: G4VEmModel.hh:808
G4double inveplus
Definition: G4VEmModel.hh:437
G4double PolarAngleLimit() const
Definition: G4VEmModel.hh:637
virtual void InitialiseForElement(const G4ParticleDefinition *, G4int Z)
Definition: G4VEmModel.cc:238
void InitialiseElementSelectors(const G4ParticleDefinition *, const G4DataVector &)
Definition: G4VEmModel.cc:146
void SetActivationLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:728
G4VEmModel * GetTripletModel()
Definition: G4VEmModel.hh:592
G4VEmFluctuationModel * GetModelOfFluctuations()
Definition: G4VEmModel.hh:568
G4bool lossFlucFlag
Definition: G4VEmModel.hh:436
G4bool flagDeexcitation
Definition: G4VEmModel.hh:416
G4VEmFluctuationModel * flucModel
Definition: G4VEmModel.hh:403
G4bool DeexcitationFlag() const
Definition: G4VEmModel.hh:658
virtual void CorrectionsAlongStep(const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
Definition: G4VEmModel.cc:406
const G4Isotope * GetCurrentIsotope() const
Definition: G4VEmModel.hh:484
G4bool flagForceBuildTable
Definition: G4VEmModel.hh:417
double A(double temperature)
virtual G4double ComputeDEDX(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.hh:500
G4bool localElmSelectors
Definition: G4VEmModel.hh:421
virtual void StartTracking(G4Track *)
Definition: G4VEmModel.cc:286
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:742
virtual G4double MaxSecondaryEnergy(const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:439
void SetAngularDistribution(G4VEmAngularDistribution *)
Definition: G4VEmModel.hh:582
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:357
virtual G4double ComputeCrossSectionPerShell(const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:367
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &)=0
virtual void InitialiseForMaterial(const G4ParticleDefinition *, const G4Material *)
Definition: G4VEmModel.cc:224
void SetSecondaryThreshold(G4double)
Definition: G4VEmModel.hh:749
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:462
virtual G4double GetParticleCharge(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:398
virtual G4double CrossSectionPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:252
int G4int
Definition: G4Types.hh:78
G4bool LPMFlag() const
Definition: G4VEmModel.hh:651
G4VParticleChange * pParticleChange
Definition: G4VEmModel.hh:431
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
void SetParticleChange(G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
Definition: G4VEmModel.cc:454
G4int nSelectors
Definition: G4VEmModel.hh:424
const G4String & GetName() const
Definition: G4VEmModel.hh:777
void SetLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:714
const G4String name
Definition: G4VEmModel.hh:405
void SetFluctuationFlag(G4bool val)
Definition: G4VEmModel.hh:686
G4bool UseAngularGeneratorFlag() const
Definition: G4VEmModel.hh:672
G4double SecondaryThreshold() const
Definition: G4VEmModel.hh:644
const std::vector< G4int > * theDensityIdx
Definition: G4VEmModel.hh:434
G4bool useAngularGenerator
Definition: G4VEmModel.hh:422
G4LossTableManager * fEmManager
Definition: G4VEmModel.hh:426
G4double GetKineticEnergy() const
G4double polarAngleLimit
Definition: G4VEmModel.hh:413
void SetForceBuildTable(G4bool val)
Definition: G4VEmModel.hh:770
G4bool IsMaster() const
Definition: G4VEmModel.hh:700
G4double HighEnergyActivationLimit() const
Definition: G4VEmModel.hh:623
const G4Element * fCurrentElement
Definition: G4VEmModel.hh:444
G4double ComputeMeanFreePath(const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:525
G4int nsec
Definition: G4VEmModel.hh:448
G4ElementData * GetElementData()
Definition: G4VEmModel.hh:801
G4VEmAngularDistribution * GetAngularDistribution()
Definition: G4VEmModel.hh:575
G4double GetZ() const
Definition: G4Element.hh:131
const G4Element * GetCurrentElement() const
Definition: G4VEmModel.hh:477
const G4Material * GetMaterial() const
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:131
virtual G4double MinPrimaryEnergy(const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
Definition: G4VEmModel.cc:422
const G4ParticleDefinition * GetParticleDefinition() const
virtual void ModelDescription(std::ostream &outFile) const
Definition: G4VEmModel.cc:476
virtual G4double ChargeSquareRatio(const G4Track &)
Definition: G4VEmModel.cc:381
G4bool localTable
Definition: G4VEmModel.hh:420
size_t idxTable
Definition: G4VEmModel.hh:435
G4bool isLocked
Definition: G4VEmModel.hh:423
G4PhysicsTable * xSectionTable
Definition: G4VEmModel.hh:432
G4double HighEnergyLimit() const
Definition: G4VEmModel.hh:609
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:66
virtual G4double ComputeDEDXPerVolume(const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
Definition: G4VEmModel.cc:243
#define DBL_MAX
Definition: templates.hh:83
virtual G4double MinEnergyCut(const G4ParticleDefinition *, const G4MaterialCutsCouple *)
Definition: G4VEmModel.cc:431
void SetTripletModel(G4VEmModel *)
Definition: G4VEmModel.hh:599
virtual G4double Value(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
Definition: G4VEmModel.cc:413
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel)
Definition: G4VEmModel.cc:218
G4VEmModel & operator=(const G4VEmModel &right)=delete
std::vector< G4EmElementSelector * > * elmSelectors
Definition: G4VEmModel.hh:425
G4bool theLPMflag
Definition: G4VEmModel.hh:415
const std::vector< G4double > * theDensityFactor
Definition: G4VEmModel.hh:433