Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4EmParameters.cc
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: G4EmParameters.cc 69320 2013-04-30 15:59:36Z vnivanch $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class file
31 //
32 //
33 // File name: G4EmParameters
34 //
35 // Author: Vladimir Ivanchenko
36 //
37 // Creation date: 18.05.2013
38 //
39 // Modifications:
40 //
41 //
42 //
43 // -------------------------------------------------------------------
44 //
45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47 
48 #include "G4EmParameters.hh"
49 #include "G4PhysicalConstants.hh"
50 #include "G4UnitsTable.hh"
51 #include "G4SystemOfUnits.hh"
52 #include "G4VEmProcess.hh"
53 #include "G4VEnergyLossProcess.hh"
54 #include "G4VAtomDeexcitation.hh"
56 #include "G4NistManager.hh"
57 #include "G4RegionStore.hh"
58 #include "G4Region.hh"
59 #include "G4ApplicationState.hh"
60 #include "G4StateManager.hh"
61 #include "G4Threading.hh"
62 
63 G4EmParameters* G4EmParameters::theInstance = nullptr;
64 
65 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
66 
68 {
69  if(nullptr == theInstance) {
70  static G4EmParameters manager;
71  theInstance = &manager;
72  }
73  return theInstance;
74 }
75 
76 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
77 
79 {
80  delete theMessenger;
81  delete emSaturation;
82 }
83 
84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
85 
86 G4EmParameters::G4EmParameters()
87 {
89  theMessenger = new G4EmParametersMessenger(this);
90 
91  fStateManager = G4StateManager::GetStateManager();
92  Initialise();
93  emSaturation = nullptr;
94 }
95 
97 {
98  if(!IsLocked()) { Initialise(); }
99 }
100 
101 void G4EmParameters::Initialise()
102 {
103  lossFluctuation = true;
104  buildCSDARange = false;
105  flagLPM = true;
106  spline = true;
107  cutAsFinalRange = false;
108  applyCuts = false;
109  fluo = false;
110  beardenFluoDir = false;
111  auger = false;
112  augerCascade = false;
113  pixe = false;
114  deexIgnoreCut = false;
115  lateralDisplacement = true;
116  muhadLateralDisplacement = false;
117  latDisplacementBeyondSafety = false;
118  useAngGeneratorForIonisation = false;
119  useMottCorrection = false;
120  integral = true;
121  birks = false;
122 
123  minSubRange = 1.0;
124  minKinEnergy = 0.1*CLHEP::keV;
125  maxKinEnergy = 100.0*CLHEP::TeV;
126  maxKinEnergyCSDA = 1.0*CLHEP::GeV;
127  lowestElectronEnergy = 1.0*CLHEP::keV;
128  lowestMuHadEnergy = 1.0*CLHEP::keV;
129  linLossLimit = 0.01;
130  bremsTh = maxKinEnergy;
131  lambdaFactor = 0.8;
132  factorForAngleLimit = 1.0;
133  thetaLimit = CLHEP::pi;
134  rangeFactor = 0.04;
135  rangeFactorMuHad = 0.2;
136  geomFactor = 2.5;
137  skin = 1.0;
138  dRoverRange = 0.2;
139  finalRange = CLHEP::mm;
140  dRoverRangeMuHad = 0.2;
141  finalRangeMuHad = 0.1*CLHEP::mm;
142 
143  nbins = 77;
144  nbinsPerDecade = 7;
145  verbose = 1;
146  workerVerbose = 0;
147 
148  mscStepLimit = fUseSafety;
149  mscStepLimitMuHad = fMinimal;
150  nucFormfactor = fExponentialNF;
151 
152  namePIXE = "Empirical";
153  nameElectronPIXE = "Livermore";
154 }
155 
157 {
158  if(IsLocked()) { return; }
159  lossFluctuation = val;
160 }
161 
163 {
164  return lossFluctuation;
165 }
166 
168 {
169  if(IsLocked()) { return; }
170  buildCSDARange = val;
171 }
172 
174 {
175  return buildCSDARange;
176 }
177 
179 {
180  if(IsLocked()) { return; }
181  flagLPM = val;
182 }
183 
185 {
186  return flagLPM;
187 }
188 
190 {
191  if(IsLocked()) { return; }
192  spline = val;
193 }
194 
196 {
197  return spline;
198 }
199 
201 {
202  if(IsLocked()) { return; }
203  cutAsFinalRange = val;
204 }
205 
207 {
208  return cutAsFinalRange;
209 }
210 
212 {
213  if(IsLocked()) { return; }
214  applyCuts = val;
215 }
216 
218 {
219  return applyCuts;
220 }
221 
223 {
224  if(IsLocked()) { return; }
225  fluo = val;
226 }
227 
229 {
230  return fluo;
231 }
232 
234 {
235  if(IsLocked()) { return; }
236  beardenFluoDir = val;
237 }
238 
240 {
241  return beardenFluoDir;
242 }
243 
245 {
246  if(IsLocked()) { return; }
247  auger = val;
248  if(val) { fluo = true; }
249 }
250 
252 {
253  return auger;
254 }
255 
257 {
258  if(IsLocked()) { return; }
259  augerCascade = val;
260  if(val) { fluo = true; auger = true; }
261 }
262 
264 {
265  return augerCascade;
266 }
267 
269 {
270  if(IsLocked()) { return; }
271  pixe = val;
272  if(val) { fluo = true; }
273 }
274 
276 {
277  return pixe;
278 }
279 
281 {
282  if(IsLocked()) { return; }
283  deexIgnoreCut = val;
284 }
285 
287 {
288  return deexIgnoreCut;
289 }
290 
292 {
293  if(IsLocked()) { return; }
294  lateralDisplacement = val;
295 }
296 
298 {
299  return lateralDisplacement;
300 }
301 
303 {
304  if(IsLocked()) { return; }
305  muhadLateralDisplacement = val;
306 }
307 
309 {
310  return muhadLateralDisplacement;
311 }
312 
314 {
315  if(IsLocked()) { return; }
316  latDisplacementBeyondSafety = val;
317 }
318 
320 {
321  return latDisplacementBeyondSafety;
322 }
323 
325 {
326  if(IsLocked()) { return; }
327  useAngGeneratorForIonisation = val;
328 }
329 
331 {
332  return useAngGeneratorForIonisation;
333 }
334 
336 {
337  if(IsLocked()) { return; }
338  useMottCorrection = val;
339 }
340 
342 {
343  return useMottCorrection;
344 }
345 
347 {
348  if(IsLocked()) { return; }
349  integral = val;
350 }
351 
353 {
354  return integral;
355 }
356 
358 {
359  if(IsLocked()) { return; }
360  birks = val;
361  if(birks) {
362  if(!emSaturation) { emSaturation = new G4EmSaturation(1); }
363  emSaturation->InitialiseG4Saturation();
364  }
365 }
366 
368 {
369  return birks;
370 }
371 
373 {
374  if(emSaturation != ptr) {
375  delete emSaturation;
376  emSaturation = ptr;
377  }
378 }
379 
381 {
382  if(!emSaturation) { SetBirksActive(true); }
383  return emSaturation;
384 }
385 
387 {
388  if(IsLocked()) { return; }
389  if(val > 0.0 && val < 1.0) {
390  minSubRange = val;
391  } else {
393  ed << "Value of MinSubRange is out of range (0 - 1): " << val
394  << " is ignored";
395  PrintWarning(ed);
396  }
397 }
398 
400 {
401  return minSubRange;
402 }
403 
405 {
406  if(IsLocked()) { return; }
407  if(val > 1.e-3*eV && val < maxKinEnergy) {
408  minKinEnergy = val;
409  nbins = nbinsPerDecade*G4lrint(std::log10(maxKinEnergy/minKinEnergy));
410  } else {
412  ed << "Value of MinKinEnergy is out of range: " << val/MeV
413  << " MeV is ignored";
414  PrintWarning(ed);
415  }
416 }
417 
419 {
420  return minKinEnergy;
421 }
422 
424 {
425  if(IsLocked()) { return; }
426  if(val > minKinEnergy && val < 1.e+7*TeV) {
427  maxKinEnergy = val;
428  nbins = nbinsPerDecade*G4lrint(std::log10(maxKinEnergy/minKinEnergy));
429  } else {
431  ed << "Value of MaxKinEnergy is out of range: "
432  << val/GeV << " GeV is ignored";
433  PrintWarning(ed);
434  }
435 }
436 
438 {
439  return maxKinEnergy;
440 }
441 
443 {
444  if(IsLocked()) { return; }
445  if(val > minKinEnergy && val <= 100*TeV) {
446  maxKinEnergyCSDA = val;
447  } else {
449  ed << "Value of MaxKinEnergyCSDA is out of range: "
450  << val/GeV << " GeV is ignored";
451  PrintWarning(ed);
452  }
453 }
454 
456 {
457  return maxKinEnergyCSDA;
458 }
459 
461 {
462  if(IsLocked()) { return; }
463  if(val >= 0.0) {
464  lowestElectronEnergy = val;
465  } else {
467  ed << "Value of lowestElectronEnergy is out of range: "
468  << val/MeV << " MeV is ignored";
469  PrintWarning(ed);
470  }
471 }
472 
474 {
475  return lowestElectronEnergy;
476 }
477 
479 {
480  if(IsLocked()) { return; }
481  if(val >= 0.0) {
482  lowestMuHadEnergy = val;
483  } else {
485  ed << "Value of lowestMuHadEnergy is out of range: "
486  << val/MeV << " MeV is ignored";
487  PrintWarning(ed);
488  }
489 }
490 
492 {
493  return lowestMuHadEnergy;
494 }
495 
497 {
498  if(IsLocked()) { return; }
499  if(val > 0.0 && val < 0.5) {
500  linLossLimit = val;
501  } else {
503  ed << "Value of linLossLimit is out of range: " << val
504  << " is ignored";
505  PrintWarning(ed);
506  }
507 }
508 
510 {
511  return linLossLimit;
512 }
513 
515 {
516  if(IsLocked()) { return; }
517  if(val > 0.0) {
518  bremsTh = val;
519  } else {
521  ed << "Value of bremsstrahlung threshold is out of range: "
522  << val/GeV << " GeV is ignored";
523  PrintWarning(ed);
524  }
525 }
526 
528 {
529  return bremsTh;
530 }
531 
533 {
534  if(IsLocked()) { return; }
535  if(val > 0.0 && val < 1.0) {
536  lambdaFactor = val;
537  } else {
539  ed << "Value of lambda factor is out of range: " << val
540  << " is ignored";
541  PrintWarning(ed);
542  }
543 }
544 
546 {
547  return lambdaFactor;
548 }
549 
551 {
552  if(IsLocked()) { return; }
553  if(val > 0.0) {
554  factorForAngleLimit = val;
555  } else {
557  ed << "Value of factor for enegry limit is out of range: "
558  << val << " is ignored";
559  PrintWarning(ed);
560  }
561 }
562 
564 {
565  return factorForAngleLimit;
566 }
567 
569 {
570  if(IsLocked()) { return; }
571  if(val >= 0.0 && val <= pi) {
572  thetaLimit = val;
573  } else {
575  ed << "Value of polar angle limit is out of range: "
576  << val << " is ignored";
577  PrintWarning(ed);
578  }
579 }
580 
582 {
583  return thetaLimit;
584 }
585 
587 {
588  if(IsLocked()) { return; }
589  if(val > 0.0 && val < 1.0) {
590  rangeFactor = val;
591  } else {
593  ed << "Value of rangeFactor is out of range: "
594  << val << " is ignored";
595  PrintWarning(ed);
596  }
597 }
598 
600 {
601  return rangeFactor;
602 }
603 
605 {
606  if(IsLocked()) { return; }
607  if(val > 0.0 && val < 1.0) {
608  rangeFactorMuHad = val;
609  } else {
611  ed << "Value of rangeFactorMuHad is out of range: "
612  << val << " is ignored";
613  PrintWarning(ed);
614  }
615 }
616 
618 {
619  return rangeFactorMuHad;
620 }
621 
623 {
624  if(IsLocked()) { return; }
625  if(val >= 1.0) {
626  geomFactor = val;
627  } else {
629  ed << "Value of geomFactor is out of range: "
630  << val << " is ignored";
631  PrintWarning(ed);
632  }
633 }
634 
636 {
637  return geomFactor;
638 }
639 
641 {
642  if(IsLocked()) { return; }
643  if(val >= 1.0) {
644  skin = val;
645  } else {
647  ed << "Value of skin is out of range: "
648  << val << " is ignored";
649  PrintWarning(ed);
650  }
651 }
652 
654 {
655  return skin;
656 }
657 
659 {
660  if(IsLocked()) { return; }
661  if(v1 > 0.0 && v1 <= 1.0 && v2 > 0.0) {
662  dRoverRange = v1;
663  finalRange = v2;
664  } else {
666  ed << "Values of step function are out of range: "
667  << v1 << ", " << v2/CLHEP::mm << " mm - are ignored";
668  PrintWarning(ed);
669  }
670 }
671 
673 {
674  if(IsLocked()) { return; }
675  if(v1 > 0.0 && v1 <= 1.0 && v2 > 0.0) {
676  dRoverRangeMuHad = v1;
677  finalRangeMuHad = v2;
678  } else {
680  ed << "Values of step function are out of range: "
681  << v1 << ", " << v2/CLHEP::mm << " mm - are ignored";
682  PrintWarning(ed);
683  }
684 }
685 
687 {
688  if(IsLocked()) { return; }
689  if(val >= 5 && val < 10000000) {
690  nbins = val;
691  nbinsPerDecade = G4lrint(nbins/std::log10(maxKinEnergy/minKinEnergy));
692  } else {
694  ed << "Value of number of bins is out of range: "
695  << val << " is ignored";
696  PrintWarning(ed);
697  }
698 }
699 
701 {
702  return nbins;
703 }
704 
706 {
707  if(IsLocked()) { return; }
708  if(val >= 5 && val < 1000000) {
709  nbinsPerDecade = val;
710  nbins = nbinsPerDecade*G4lrint(std::log10(maxKinEnergy/minKinEnergy));
711  } else {
713  ed << "Value of number of bins per decade is out of range: "
714  << val << " is ignored";
715  PrintWarning(ed);
716  }
717 }
718 
720 {
721  return nbinsPerDecade;
722 }
723 
725 {
726  if(IsLocked()) { return; }
727  verbose = val;
728  workerVerbose = std::min(workerVerbose, verbose);
729 }
730 
732 {
733  return verbose;
734 }
735 
737 {
738  if(IsLocked()) { return; }
739  workerVerbose = val;
740 }
741 
743 {
744  return workerVerbose;
745 }
746 
748 {
749  if(IsLocked()) { return; }
750  mscStepLimit = val;
751 }
752 
754 {
755  return mscStepLimit;
756 }
757 
759 {
760  if(IsLocked()) { return; }
761  mscStepLimitMuHad = val;
762 }
763 
765 {
766  return mscStepLimitMuHad;
767 }
768 
769 void
771 {
772  if(IsLocked()) { return; }
773  nucFormfactor = val;
774 }
775 
777 {
778  return nucFormfactor;
779 }
780 
782 {
783  G4cout << "G4EmParameters::SetPIXECrossSectionModel " << sss << G4endl;
784  if(IsLocked()) { return; }
785  namePIXE = sss;
786 }
787 
789 {
790  return namePIXE;
791 }
792 
794 {
795  if(IsLocked()) { return; }
796  nameElectronPIXE = sss;
797 }
798 
800 {
801  return nameElectronPIXE;
802 }
803 
804 void G4EmParameters::PrintWarning(G4ExceptionDescription& ed) const
805 {
806  G4Exception("G4EmParameters", "em0044", JustWarning, ed);
807 }
808 
809 G4String G4EmParameters::CheckRegion(const G4String& reg) const
810 {
811  G4String r = reg;
812  if(r == "" || r == "world" || r == "World") {
813  r = "DefaultRegionForTheWorld";
814  }
815  return r;
816 }
817 
819  const G4String& region,
820  const G4String& type)
821 {
822  G4String r = CheckRegion(region);
823  G4int nreg = m_regnamesPAI.size();
824  for(G4int i=0; i<nreg; ++i) {
825  if((m_particlesPAI[i] == particle ||
826  m_particlesPAI[i] == "all" ||
827  particle == "all") &&
828  (m_regnamesPAI[i] == r ||
829  m_regnamesPAI[i] == "DefaultRegionForTheWorld" ||
830  r == "DefaultRegionForTheWorld") ) {
831 
832  m_typesPAI[i] = type;
833  if(particle == "all") { m_particlesPAI[i] = particle; }
834  if(r == "DefaultRegionForTheWorld") { m_regnamesPAI[i] = r; }
835  return;
836  }
837  }
838  m_particlesPAI.push_back(particle);
839  m_regnamesPAI.push_back(r);
840  m_typesPAI.push_back(type);
841 }
842 
843 const std::vector<G4String>& G4EmParameters::ParticlesPAI() const
844 {
845  return m_particlesPAI;
846 }
847 
848 const std::vector<G4String>& G4EmParameters::RegionsPAI() const
849 {
850  return m_regnamesPAI;
851 }
852 
853 const std::vector<G4String>& G4EmParameters::TypesPAI() const
854 {
855  return m_typesPAI;
856 }
857 
859 {
860  G4String r = CheckRegion(region);
861  G4int nreg = m_regnamesME.size();
862  for(G4int i=0; i<nreg; ++i) {
863  if(r == m_regnamesME[i]) { return; }
864  }
865  m_regnamesME.push_back(r);
866 }
867 
868 const std::vector<G4String>& G4EmParameters::RegionsMicroElec() const
869 {
870  return m_regnamesME;
871 }
872 
873 void G4EmParameters::AddDNA(const G4String& region, const G4String& type)
874 {
875  G4String r = CheckRegion(region);
876  G4int nreg = m_regnamesDNA.size();
877  for(G4int i=0; i<nreg; ++i) {
878  if(r == m_regnamesDNA[i]) { return; }
879  }
880  m_regnamesDNA.push_back(r);
881  m_typesDNA.push_back(type);
882 }
883 
884 const std::vector<G4String>& G4EmParameters::RegionsDNA() const
885 {
886  return m_regnamesDNA;
887 }
888 
889 const std::vector<G4String>& G4EmParameters::TypesDNA() const
890 {
891  return m_typesDNA;
892 }
893 
894 void G4EmParameters::AddMsc(const G4String& region, const G4String& type)
895 {
896  G4String r = CheckRegion(region);
897  G4int nreg = m_regnamesMsc.size();
898  for(G4int i=0; i<nreg; ++i) {
899  if(r == m_regnamesMsc[i]) { return; }
900  }
901  m_regnamesMsc.push_back(r);
902  m_typesMsc.push_back(type);
903 }
904 
905 const std::vector<G4String>& G4EmParameters::RegionsMsc() const
906 {
907  return m_regnamesMsc;
908 }
909 
910 const std::vector<G4String>& G4EmParameters::TypesMsc() const
911 {
912  return m_typesMsc;
913 }
914 
916 {
917  if(IsLocked()) { return; }
918  G4String r = CheckRegion(region);
919  G4int nreg = m_regnamesSubCut.size();
920  for(G4int i=0; i<nreg; ++i) {
921  if(r == m_regnamesSubCut[i]) {
922  m_subCuts[i] = val;
923  return;
924  }
925  }
926  m_regnamesSubCut.push_back(r);
927  m_subCuts.push_back(val);
928 }
929 
930 void
932  G4bool fauger, G4bool fpixe)
933 {
934  if(IsLocked()) { return; }
935  if(fdeex) { fluo = true; }
936  G4String r = CheckRegion(region);
937  G4int nreg = m_regnamesDeex.size();
938  if(0 == nreg && r != "DefaultRegionForTheWorld") {
939  m_regnamesDeex.push_back("DefaultRegionForTheWorld");
940  m_fluo.push_back(false);
941  m_auger.push_back(false);
942  m_pixe.push_back(false);
943  nreg = 1;
944  }
945  for(G4int i=0; i<nreg; ++i) {
946  if(r == m_regnamesDeex[i]) {
947  m_fluo[i] = fdeex;
948  m_auger[i]= fauger;
949  m_pixe[i] = fpixe;
950  return;
951  }
952  }
953  m_regnamesDeex.push_back(r);
954  m_fluo.push_back(fdeex);
955  m_auger.push_back(fauger);
956  m_pixe.push_back(fpixe);
957 }
958 
959 void
961  G4double val, G4bool wflag)
962 {
963  if(IsLocked()) { return; }
964  if(val > 0.0) {
965  G4int n = m_procBiasedXS.size();
966  for(G4int i=0; i<n; ++i) {
967  if(procname == m_procBiasedXS[i]) {
968  m_factBiasedXS[i] = val;
969  m_weightBiasedXS[i]= wflag;
970  return;
971  }
972  }
973  m_procBiasedXS.push_back(procname);
974  m_factBiasedXS.push_back(val);
975  m_weightBiasedXS.push_back(wflag);
976  } else {
978  ed << "Process: " << procname << " XS biasing factor "
979  << val << " is negative - ignored";
980  PrintWarning(ed);
981  }
982 }
983 
984 void
986  const G4String& region,
987  G4double length,
988  G4bool wflag)
989 {
990  if(IsLocked()) { return; }
991  G4String r = CheckRegion(region);
992  if(length >= 0.0) {
993  G4int n = m_procForced.size();
994  for(G4int i=0; i<n; ++i) {
995  if(procname == m_procForced[i] && r == m_regnamesForced[i] ) {
996  m_lengthForced[i] = length;
997  m_weightForced[i]= wflag;
998  return;
999  }
1000  }
1001  m_regnamesForced.push_back(r);
1002  m_procForced.push_back(procname);
1003  m_lengthForced.push_back(length);
1004  m_weightForced.push_back(wflag);
1005  } else {
1007  ed << "Process: " << procname << " in region " << r
1008  << " : forced interacttion length= "
1009  << length << " is negative - ignored";
1010  PrintWarning(ed);
1011  }
1012 }
1013 
1014 void
1016  const G4String& region,
1017  G4double factor,
1018  G4double energyLimit)
1019 {
1020  if(IsLocked()) { return; }
1021  G4String r = CheckRegion(region);
1022  if(factor >= 0.0 && energyLimit >= 0.0) {
1023  G4int n = m_procBiasedSec.size();
1024  for(G4int i=0; i<n; ++i) {
1025  if(procname == m_procBiasedSec[i] && r == m_regnamesBiasedSec[i] ) {
1026  m_factBiasedSec[i] = factor;
1027  m_elimBiasedSec[i] = energyLimit;
1028  return;
1029  }
1030  }
1031  m_regnamesBiasedSec.push_back(r);
1032  m_procBiasedSec.push_back(procname);
1033  m_factBiasedSec.push_back(factor);
1034  m_elimBiasedSec.push_back(energyLimit);
1035  } else {
1037  ed << "Process: " << procname << " in region " << r
1038  << " : secondary bised factor= "
1039  << factor << ", Elim= " << energyLimit << " - ignored";
1040  PrintWarning(ed);
1041  }
1042 }
1043 
1045  G4bool isElectron) const
1046 {
1047  if(isElectron) { ptr->SetStepFunction(dRoverRange, finalRange, false); }
1048  else { ptr->SetStepFunction(dRoverRangeMuHad, finalRangeMuHad, false); }
1049 
1050  G4RegionStore* regionStore = G4RegionStore::GetInstance();
1051  G4int n = m_regnamesSubCut.size();
1052  for(G4int i=0; i<n; ++i) {
1053  const G4Region* reg = regionStore->GetRegion(m_regnamesSubCut[i], false);
1054  if(reg) { ptr->ActivateSubCutoff(m_subCuts[i], reg); }
1055  }
1056  n = m_procBiasedXS.size();
1057  for(G4int i=0; i<n; ++i) {
1058  if(ptr->GetProcessName() == m_procBiasedXS[i]) {
1059  ptr->SetCrossSectionBiasingFactor(m_factBiasedXS[i],
1060  m_weightBiasedXS[i]);
1061  break;
1062  }
1063  }
1064  n = m_procForced.size();
1065  for(G4int i=0; i<n; ++i) {
1066  if(ptr->GetProcessName() == m_procForced[i]) {
1067  ptr->ActivateForcedInteraction(m_lengthForced[i],
1068  m_regnamesForced[i],
1069  m_weightForced[i]);
1070  break;
1071  }
1072  }
1073  n = m_procBiasedSec.size();
1074  for(G4int i=0; i<n; ++i) {
1075  if(ptr->GetProcessName() == m_procBiasedSec[i]) {
1076  ptr->ActivateSecondaryBiasing(m_regnamesBiasedSec[i],
1077  m_factBiasedSec[i],
1078  m_elimBiasedSec[i]);
1079  break;
1080  }
1081  }
1082 }
1083 
1085 {
1086  G4int n = m_procBiasedXS.size();
1087  for(G4int i=0; i<n; ++i) {
1088  if(ptr->GetProcessName() == m_procBiasedXS[i]) {
1089  ptr->SetCrossSectionBiasingFactor(m_factBiasedXS[i],
1090  m_weightBiasedXS[i]);
1091  break;
1092  }
1093  }
1094  n = m_procForced.size();
1095  for(G4int i=0; i<n; ++i) {
1096  if(ptr->GetProcessName() == m_procForced[i]) {
1097  ptr->ActivateForcedInteraction(m_lengthForced[i],
1098  m_regnamesForced[i],
1099  m_weightForced[i]);
1100  break;
1101  }
1102  }
1103  n = m_procBiasedSec.size();
1104  for(G4int i=0; i<n; ++i) {
1105  if(ptr->GetProcessName() == m_procBiasedSec[i]) {
1106  ptr->ActivateSecondaryBiasing(m_regnamesBiasedSec[i],
1107  m_factBiasedSec[i],
1108  m_elimBiasedSec[i]);
1109  break;
1110  }
1111  }
1112 }
1113 
1115 {
1116  G4int n = m_regnamesDeex.size();
1117  for(G4int i=0; i<n; ++i) {
1118  ptr->SetDeexcitationActiveRegion(m_regnamesDeex[i],
1119  m_fluo[i], m_auger[i], m_pixe[i]);
1120  }
1121 }
1122 
1123 std::ostream& G4EmParameters::StreamInfo(std::ostream& os) const
1124 {
1125  G4int prec = os.precision(5);
1126  os << "=======================================================================" << "\n";
1127  os << "====== Electromagnetic Physics Parameters ========" << "\n";
1128  os << "=======================================================================" << "\n";
1129  os << "Fluctuations of dE/dx are enabled " <<lossFluctuation << "\n";
1130  os << "Build CSDA range enabled " <<buildCSDARange << "\n";
1131  os << "LPM effect enabled " <<flagLPM << "\n";
1132  os << "Spline of EM tables enabled " <<spline << "\n";
1133  os << "Use cut as a final range enabled " <<finalRange << "\n";
1134  os << "Apply cuts on all EM processes " <<applyCuts << "\n";
1135  os << "Fluorescence enabled " <<fluo << "\n";
1136  os << "Fluorescence Bearden data files enabled " <<beardenFluoDir << "\n";
1137  os << "Auger electron production enabled " <<auger << "\n";
1138  os << "Auger cascade enabled " <<augerCascade << "\n";
1139  os << "PIXE atomic de-excitation enabled " <<pixe << "\n";
1140  os << "De-excitation module ignores cuts " <<deexIgnoreCut << "\n";
1141  os << "Msc lateraral displacement for e+- enabled " <<lateralDisplacement << "\n";
1142  os << "Msc lateraral displacement for muons and hadrons " <<muhadLateralDisplacement << "\n";
1143  os << "Msc lateraral displacement beyond geometry safety " <<latDisplacementBeyondSafety << "\n";
1144  os << "Enable angular generator interface "
1145  <<useAngGeneratorForIonisation << "\n";
1146  os << "Use Mott correction for e- scattering "
1147  <<useMottCorrection << "\n";
1148  os << "Use integral approach for tracking "
1149  <<integral << "\n";
1150  os << "Use built-in Birks satuaration "
1151  << birks << "\n";
1152 
1153  os << "Factor of cut reduction for sub-cutoff method " <<minSubRange << "\n";
1154  os << "Min kinetic energy for tables "
1155  <<G4BestUnit(minKinEnergy,"Energy") << "\n";
1156  os << "Max kinetic energy for tables "
1157  <<G4BestUnit(maxKinEnergy,"Energy") << "\n";
1158  os << "Max kinetic energy for CSDA tables "
1159  <<G4BestUnit(maxKinEnergyCSDA,"Energy") << "\n";
1160  os << "Lowest e+e- kinetic energy "
1161  <<G4BestUnit(lowestElectronEnergy,"Energy") << "\n";
1162  os << "Lowest muon/hadron kinetic energy "
1163  <<G4BestUnit(lowestMuHadEnergy,"Energy") << "\n";
1164  os << "Linear loss limit " <<linLossLimit << "\n";
1165  os << "Bremsstrahlung energy threshold above which \n"
1166  << " primary is added to the list of secondary "
1167  <<G4BestUnit(bremsTh,"Energy") << "\n";
1168  os << "X-section factor for integral approach " <<lambdaFactor << "\n";
1169  os << "Factor used for dynamic computation of angular \n"
1170  << " limit between single and multiple scattering " << factorForAngleLimit << "\n";
1171  os << "Fixed angular limit between single \n"
1172  << " and multiple scattering "
1173  <<thetaLimit/rad << " rad" << "\n";
1174  os << "Range factor for msc step limit for e+- " <<rangeFactor << "\n";
1175  os << "Range factor for msc step limit for muons/hadrons " <<rangeFactorMuHad << "\n";
1176  os << "Geometry factor for msc step limitation of e+- " <<geomFactor << "\n";
1177  os << "Skin parameter for msc step limitation of e+- " <<skin << "\n";
1178  os << "Step function for e+- " <<"("<< dRoverRange
1179  << ", " << finalRange << " mm)\n";
1180  os << "Step function for muons/hadrons " <<"("<< dRoverRangeMuHad
1181  << ", " << finalRangeMuHad << " mm)\n";
1182 
1183  os << "Number of bins in tables " <<nbins << "\n";
1184  os << "Number of bins per decade of a table " <<nbinsPerDecade << "\n";
1185  os << "Verbose level " <<verbose << "\n";
1186  os << "Verbose level for worker thread " <<workerVerbose << "\n";
1187 
1188  os << "Type of msc step limit algorithm for e+- " <<mscStepLimit << "\n";
1189  os << "Type of msc step limit algorithm for muons/hadrons " <<mscStepLimitMuHad << "\n";
1190  os << "Type of nuclear form-factor " <<nucFormfactor << "\n";
1191 
1192  os << "Type of PIXE cross section for hadrons " <<namePIXE << "\n";
1193  os << "Type of PIXE cross section for e+- " <<nameElectronPIXE << "\n";
1194  os << "=======================================================================" << "\n";
1195  os.precision(prec);
1196  return os;
1197 }
1198 
1200 {
1201  StreamInfo(G4cout);
1202 }
1203 
1204 std::ostream& operator<< (std::ostream& os, const G4EmParameters& par)
1205 {
1206  return par.StreamInfo(os);
1207 }
1208 
1209 G4bool G4EmParameters::IsLocked() const
1210 {
1211  return (!G4Threading::IsMasterThread() ||
1212  (fStateManager->GetCurrentState() != G4State_PreInit &&
1213  fStateManager->GetCurrentState() != G4State_Init &&
1214  fStateManager->GetCurrentState() != G4State_Idle));
1215 }
1216 
1217 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.....
1218 
void ActivateForcedInteraction(G4double length=0.0, const G4String &r="", G4bool flag=true)
G4bool UseCutAsFinalRange() const
void SetLossFluctuations(G4bool val)
G4int NumberOfBinsPerDecade() const
void SetApplyCuts(G4bool val)
void SetProcessBiasingFactor(const G4String &procname, G4double val, G4bool wflag)
G4int NumberOfBins() const
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
void SetVerbose(G4int val)
G4bool Spline() const
G4int WorkerVerbose() const
G4double MaxKinEnergy() const
void SetEmSaturation(G4EmSaturation *)
G4bool isElectron(G4int ityp)
void ActivateSecondaryBiasing(const G4String &name, const G4String &region, G4double factor, G4double energyLimit)
void SetDeexcitationIgnoreCut(G4bool val)
void SetUseMottCorrection(G4bool val)
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
void SetLowestElectronEnergy(G4double val)
void SetLatDisplacementBeyondSafety(G4bool val)
G4MscStepLimitType MscMuHadStepLimitType() const
void SetMscStepLimitType(G4MscStepLimitType val)
static constexpr double keV
G4bool AugerCascade() const
G4double LowestElectronEnergy() const
G4bool BirksActive() const
void SetBeardenFluoDir(G4bool val)
void SetLinearLossLimit(G4double val)
G4double MscGeomFactor() const
G4double MscMuHadRangeFactor() const
G4double MscThetaLimit() const
void SetAuger(G4bool val)
const std::vector< G4String > & ParticlesPAI() const
const std::vector< G4String > & RegionsMicroElec() const
G4bool Integral() const
void SetDeexcitationActiveRegion(const G4String &rname, G4bool valDeexcitation, G4bool valAuger, G4bool valPIXE)
std::ostream & StreamInfo(std::ostream &os) const
void SetNumberOfBins(G4int val)
void SetMinSubRange(G4double val)
G4bool LPM() const
void SetMaxEnergyForCSDARange(G4double val)
const std::vector< G4String > & RegionsPAI() const
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
static constexpr double rad
Definition: G4SIunits.hh:149
G4bool ApplyCuts() const
void DefineRegParamForLoss(G4VEnergyLossProcess *, G4bool isElectron) const
void SetStepFunctionMuHad(G4double v1, G4double v2)
const std::vector< G4String > & TypesDNA() const
G4bool Fluo() const
void SetPIXEElectronCrossSectionModel(const G4String &)
int G4int
Definition: G4Types.hh:78
void SetDeexActiveRegion(const G4String &region, G4bool fdeex, G4bool fauger, G4bool fpixe)
const std::vector< G4String > & RegionsMsc() const
void SetMaxEnergy(G4double val)
void SetBremsstrahlungTh(G4double val)
void SetStepFunction(G4double v1, G4double v2, G4bool lock=true)
static G4NistManager * Instance()
G4double MinSubRange() const
G4bool DeexcitationIgnoreCut() const
void ActivateSecondaryBiasing(const G4String &region, G4double factor, G4double energyLimit)
static const G4double reg
G4bool LatDisplacementBeyondSafety() const
void SetBirksActive(G4bool val)
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
static constexpr double TeV
Definition: G4SIunits.hh:218
void AddPAIModel(const G4String &particle, const G4String &region, const G4String &type)
static G4RegionStore * GetInstance()
void SetLateralDisplacement(G4bool val)
G4bool BuildCSDARange() const
G4EmSaturation * GetEmSaturation()
void SetWorkerVerbose(G4int val)
static G4StateManager * GetStateManager()
G4double LinearLossLimit() const
void SetPIXECrossSectionModel(const G4String &)
G4double LambdaFactor() const
static const double prec
Definition: RanecuEngine.cc:58
void SetCrossSectionBiasingFactor(G4double f, G4bool flag=true)
G4GLOB_DLL std::ostream G4cout
G4int Verbose() const
G4bool MuHadLateralDisplacement() const
const G4String & PIXECrossSectionModel()
G4double LowestMuHadEnergy() const
bool G4bool
Definition: G4Types.hh:79
void SetMscRangeFactor(G4double val)
void InitialiseG4Saturation()
static constexpr double mm
Definition: SystemOfUnits.h:95
void SetAugerCascade(G4bool val)
void Dump() const
void SetNumberOfBinsPerDecade(G4int val)
G4bool LateralDisplacement() const
void ActivateForcedInteraction(G4double length, const G4String &region, G4bool flag=true)
static constexpr double eV
Definition: G4SIunits.hh:215
void SetNuclearFormfactorType(G4NuclearFormfactorType val)
void SetLowestMuHadEnergy(G4double val)
void SetMscGeomFactor(G4double val)
G4double MinKinEnergy() const
G4ApplicationState GetCurrentState() const
const G4int n
void SetMscMuHadStepLimitType(G4MscStepLimitType val)
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
void AddDNA(const G4String &region, const G4String &type)
void SetSubCutoff(G4bool val, const G4String &region="")
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static constexpr double GeV
G4NuclearFormfactorType NuclearFormfactorType() const
void SetBuildCSDARange(G4bool val)
void AddMicroElec(const G4String &region)
G4bool Auger() const
G4bool LossFluctuation() const
const std::vector< G4String > & TypesMsc() const
void SetMscMuHadRangeFactor(G4double val)
void SetPixe(G4bool val)
void SetSpline(G4bool val)
static const char sss[MAX_N_PAR+2]
Definition: Evaluator.cc:64
void SetMuHadLateralDisplacement(G4bool val)
int G4lrint(double ad)
Definition: templates.hh:163
void SetMinEnergy(G4double val)
void SetLPM(G4bool val)
G4bool UseAngularGeneratorForIonisation() const
G4NuclearFormfactorType
G4MscStepLimitType
G4double BremsstrahlungTh() const
static G4EmParameters * Instance()
void DefineRegParamForDeex(G4VAtomDeexcitation *) const
std::ostream & operator<<(std::ostream &, const BasicVector3D< float > &)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double MscRangeFactor() const
void SetIntegral(G4bool val)
static constexpr double GeV
Definition: G4SIunits.hh:217
void AddMsc(const G4String &region, const G4String &type)
void ActivateAngularGeneratorForIonisation(G4bool val)
void ActivateForcedInteraction(const G4String &procname, const G4String &region, G4double length, G4bool wflag)
void SetUseCutAsFinalRange(G4bool val)
G4bool UseMottCorrection() const
const std::vector< G4String > & RegionsDNA() const
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
G4bool Pixe() const
G4bool IsMasterThread()
Definition: G4Threading.cc:146
static constexpr double pi
Definition: G4SIunits.hh:75
static constexpr double TeV
G4MscStepLimitType MscStepLimitType() const
void ActivateSubCutoff(G4bool val, const G4Region *region=nullptr)
void SetMscThetaLimit(G4double val)
void SetLambdaFactor(G4double val)
void SetFactorForAngleLimit(G4double val)
double G4double
Definition: G4Types.hh:76
void DefineRegParamForEM(G4VEmProcess *) const
G4double MaxEnergyForCSDARange() const
G4double MscSkin() const
const G4String & PIXEElectronCrossSectionModel()
static constexpr double pi
Definition: SystemOfUnits.h:54
const std::vector< G4String > & TypesPAI() const
void SetStepFunction(G4double v1, G4double v2)
void SetFluo(G4bool val)
G4double FactorForAngleLimit() const
G4bool BeardenFluoDir() const
void SetMscSkin(G4double val)