Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4EmModelActivator.cc
이 파일의 문서화 페이지로 가기
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: G4EmModelActivator.cc 1651 2015-05-02 16:40:24Z vnivanch $
27 // GEANT4 tag $Name$
28 //
29 //---------------------------------------------------------------------------
30 //
31 // ClassName: G4EmModelActivator
32 //
33 // Author: V.Ivanchenko 01.06.2015
34 //
35 // Organisation: G4AI
36 // Contract: ESA contract 4000107387/12/NL/AT
37 // Customer: ESA/ESTEC
38 //
39 // Modified:
40 //
41 //----------------------------------------------------------------------------
42 //
43 
44 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
45 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
46 
47 #include "G4EmModelActivator.hh"
48 #include "G4EmParameters.hh"
49 #include "G4ParticleDefinition.hh"
50 #include "G4ParticleTable.hh"
51 #include "G4RegionStore.hh"
52 #include "G4Region.hh"
53 #include "G4VEnergyLossProcess.hh"
54 #include "G4VMscModel.hh"
55 #include "G4LossTableManager.hh"
56 #include "G4EmConfigurator.hh"
57 #include "G4PAIModel.hh"
58 #include "G4PAIPhotModel.hh"
59 #include "G4Gamma.hh"
60 #include "G4Electron.hh"
61 #include "G4Positron.hh"
62 #include "G4MuonPlus.hh"
63 #include "G4MuonMinus.hh"
64 #include "G4Proton.hh"
65 #include "G4GenericIon.hh"
66 #include "G4Alpha.hh"
67 #include "G4ProcessManager.hh"
68 #include "G4DummyModel.hh"
69 #include "G4EmProcessSubType.hh"
70 #include "G4PhysicsListHelper.hh"
71 #include "G4EmParticleList.hh"
72 
73 #include "G4MicroElecElastic.hh"
75 
76 #include "G4MicroElecInelastic.hh"
78 
79 #include "G4BraggModel.hh"
80 #include "G4BraggIonModel.hh"
81 #include "G4BetheBlochModel.hh"
82 #include "G4UrbanMscModel.hh"
83 #include "G4MollerBhabhaModel.hh"
84 #include "G4IonFluctuations.hh"
86 #include "G4LowECapture.hh"
87 #include "G4hMultipleScattering.hh"
90 #include "G4ionIonisation.hh"
91 #include "G4KleinNishinaModel.hh"
92 
93 #include "G4CoulombScattering.hh"
95 #include "G4WentzelVIModel.hh"
97 #include "G4RayleighScattering.hh"
98 #include "G4UrbanMscModel.hh"
100 #include "G4LowEPComptonModel.hh"
101 
108 
110 #include "G4PenelopeComptonModel.hh"
116 
117 #include "G4SystemOfUnits.hh"
118 #include "G4Threading.hh"
119 #include <vector>
120 
121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
122 
124  : baseName(emphys)
125 {
127 
128  const std::vector<G4String> regnamesPAI = theParameters->RegionsPAI();
129  if(regnamesPAI.size() > 0)
130  {
131  ActivatePAI();
132  }
133  const std::vector<G4String> regnamesME = theParameters->RegionsMicroElec();
134  if(regnamesME.size() > 0)
135  {
137  }
138  const std::vector<G4String> regnamesMSC = theParameters->RegionsMsc();
139  if(regnamesMSC.size() > 0)
140  {
142  }
143 }
144 
145 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
146 
148 {
149  const std::vector<G4String> regnamesPhys = theParameters->RegionsPhysics();
150  G4int nreg = regnamesPhys.size();
151  if(0 == nreg) { return; }
152  G4int verbose = theParameters->Verbose() - 1;
153  if(verbose > 0) {
154  G4cout << "### G4EmModelActivator::ActivateEmOptions for " << nreg << " regions"
155  << G4endl;
156  }
157  const std::vector<G4String> typesPhys = theParameters->TypesMsc();
158 
159  // start configuration of models
162  const G4ParticleDefinition* phot = G4Gamma::Gamma();
163  const G4ParticleDefinition* prot = G4Proton::Proton();
165  G4EmConfigurator* em_config = man->EmConfigurator();
166  G4VAtomDeexcitation* adeexc = man->AtomDeexcitation();
168  G4VEmModel* mod;
169 
170  // high energy limit for low-energy e+- model of msc
171  G4double mscEnergyLimit = G4EmParameters::Instance()->MscEnergyLimit();
172 
173  // high energy limit for Livermore and Penelope models
174  G4double highEnergyLimit = 1*GeV;
175 
176  // general high energy limit
177  G4double highEnergy = 100*TeV;
178 
179  for(G4int i=0; i<nreg; ++i) {
180  G4String reg = regnamesPhys[i];
181  if(verbose > 0) {
182  G4cout << i << ". region <" << reg << ">; type <" << typesPhys[i] << "> "
183  << G4endl;
184  }
185 
186  if(baseName == typesPhys[i]) { continue; }
187 
188  if("G4EmStandard" == typesPhys[i]) {
189  G4UrbanMscModel* msc = new G4UrbanMscModel();
190  msc->SetLocked(true);
191  AddStandardScattering(elec, em_config, msc, reg, mscEnergyLimit, highEnergy);
192 
193  msc = new G4UrbanMscModel();
194  msc->SetLocked(true);
195  AddStandardScattering(posi, em_config, msc, reg, mscEnergyLimit, highEnergy);
196 
197  } else if("G4EmStandard_opt1" == typesPhys[i] ||
198  "G4EmStandard_opt2" == typesPhys[i]) {
199  G4UrbanMscModel* msc = new G4UrbanMscModel();
201  msc->SetRangeFactor(0.2);
202  msc->SetLocked(true);
203  AddStandardScattering(elec, em_config, msc, reg, mscEnergyLimit, highEnergy);
204 
205  msc = new G4UrbanMscModel();
207  msc->SetRangeFactor(0.2);
208  msc->SetLocked(true);
209  AddStandardScattering(posi, em_config, msc, reg, mscEnergyLimit, highEnergy);
210 
211  } else if("G4EmStandard_opt3" == typesPhys[i]) {
212 
213  G4DummyModel* dummy = new G4DummyModel();
214  G4UrbanMscModel* msc = new G4UrbanMscModel();
216  msc->SetLocked(true);
217  em_config->SetExtraEmModel("e-", "msc", msc, reg);
218  FindOrAddProcess(elec, "CoulombScat");
219  em_config->SetExtraEmModel("e-", "CoulombScat", dummy, reg);
220 
221  msc = new G4UrbanMscModel();
223  msc->SetLocked(true);
224  em_config->SetExtraEmModel("e+", "msc", msc, reg);
225  FindOrAddProcess(posi, "CoulombScat");
226  em_config->SetExtraEmModel("e+", "CoulombScat", dummy, reg);
227 
228  msc = new G4UrbanMscModel();
230  msc->SetRangeFactor(0.2);
231  msc->SetLocked(true);
232  em_config->SetExtraEmModel("proton", "msc", msc, reg);
233  FindOrAddProcess(prot, "CoulombScat");
234  em_config->SetExtraEmModel("proton", "CoulombScat", dummy, reg);
235 
237  theParameters->SetDeexActiveRegion(reg, true, false, false);
238  }
240  FindOrAddProcess(phot, "Rayl");
241  mod = new G4LivermoreRayleighModel();
242  em_config->SetExtraEmModel("gamma", "Rayl", mod, reg);
243  FindOrAddProcess(phot, "phot");
244  mod = new G4LivermorePhotoElectricModel();
245  em_config->SetExtraEmModel("gamma", "phot", mod, reg);
246  FindOrAddProcess(phot, "compt");
247  mod = new G4KleinNishinaModel();
248  em_config->SetExtraEmModel("gamma", "compt", mod, reg);
249 
250  } else if("G4EmStandard_opt4" == typesPhys[i]) {
253  msc->SetRangeFactor(0.1);
254  msc->SetLocked(true);
255  AddStandardScattering(elec, em_config, msc, reg, mscEnergyLimit, highEnergy);
256 
257  msc = new G4GoudsmitSaundersonMscModel();
259  msc->SetRangeFactor(0.1);
260  msc->SetLocked(true);
261  AddStandardScattering(posi, em_config, msc, reg, mscEnergyLimit, highEnergy);
262 
264  theParameters->SetDeexActiveRegion(reg, true, false, false);
265  theParameters->SetStepFunction(0.2, 10*um);
267  }
269  FindOrAddProcess(phot, "Rayl");
270  mod = new G4LivermoreRayleighModel();
271  em_config->SetExtraEmModel("gamma", "Rayl", mod, reg);
272  FindOrAddProcess(phot, "compt");
273  mod = new G4KleinNishinaModel();
274  em_config->SetExtraEmModel("gamma", "compt", mod, reg);
275  mod = new G4LowEPComptonModel();
276  mod->SetHighEnergyLimit(20*MeV);
277  em_config->SetExtraEmModel("gamma", "compt", mod, reg);
278 
279  } else if("G4EmStandardGS" == typesPhys[i]) {
281  msc->SetRangeFactor(0.06);
282  msc->SetLocked(true);
283  AddStandardScattering(elec, em_config, msc, reg, mscEnergyLimit, highEnergy);
284 
285  msc = new G4GoudsmitSaundersonMscModel();
286  msc->SetRangeFactor(0.06);
287  msc->SetLocked(true);
288  AddStandardScattering(posi, em_config, msc, reg, mscEnergyLimit, highEnergy);
289 
290  } else if("G4EmStandardWVI" == typesPhys[i]) {
291  G4WentzelVIModel* msc = new G4WentzelVIModel();
292  msc->SetLocked(true);
293  AddStandardScattering(elec, em_config, msc, reg, mscEnergyLimit, highEnergy);
294 
295  msc = new G4WentzelVIModel();
296  msc->SetLocked(true);
297  AddStandardScattering(posi, em_config, msc, reg, mscEnergyLimit, highEnergy);
298 
299  FindOrAddProcess(elec, "CoulombScat");
300  FindOrAddProcess(posi, "CoulombScat");
301  mod = new G4eCoulombScatteringModel();
302  em_config->SetExtraEmModel("e-", "CoulombScat", mod, reg, 0.0, mscEnergyLimit);
303  mod = new G4eCoulombScatteringModel();
304  em_config->SetExtraEmModel("e+", "CoulombScat", mod, reg, 0.0, mscEnergyLimit);
305 
307  theParameters->SetDeexActiveRegion(regnamesPhys[i], true, false, false);
308  }
310 
311  } else if("G4EmStandardSS" == typesPhys[i] &&
312  baseName != "G4EmStandard_opt3") {
313  G4EmParticleList emList;
314  for(const auto& particleName : emList.PartNames()) {
315  G4ParticleDefinition* particle = table->FindParticle(particleName);
316  if(particle && 0.0 != particle->GetPDGCharge()) {
317  FindOrAddProcess(particle, "CoulombScat");
319  sc->SetPolarAngleLimit(0.0);
320  sc->SetLocked(true);
321  em_config->SetExtraEmModel(particleName, "CoulombScat", sc, reg);
322  if(particleName == "mu+" || particleName == "mu-") {
323  em_config->SetExtraEmModel(particleName, "muMsc",
324  new G4DummyModel(), reg);
325  } else {
326  em_config->SetExtraEmModel(particleName, "msc",
327  new G4DummyModel(), reg);
328  }
329  }
330  }
332  theParameters->SetDeexActiveRegion(regnamesPhys[i], true, true, true);
333  }
335 
336  } else if("G4EmLivermore" == typesPhys[i]) {
337 
338  mod = new G4LivermorePhotoElectricModel();
339  em_config->SetExtraEmModel("gamma", "phot", mod, reg);
340  mod = new G4LivermoreComptonModel();
341  em_config->SetExtraEmModel("gamma", "compt", mod, reg);
343  em_config->SetExtraEmModel("gamma", "conv", mod, reg);
344 
345  FindOrAddProcess(phot, "Rayl");
346  mod = new G4LivermoreRayleighModel();
347  em_config->SetExtraEmModel("gamma", "Rayl", mod, reg);
348 
349  mod = new G4LivermoreIonisationModel();
351  em_config->SetExtraEmModel("e-", "eIoni", mod, reg, 0.0, 0.1*MeV, uf);
352  mod = new G4LivermoreBremsstrahlungModel();
353  em_config->SetExtraEmModel("e-", "eBrem", mod, reg, 0.0, highEnergyLimit);
354 
356  theParameters->SetDeexActiveRegion(regnamesPhys[i], true, false, false);
357  }
359 
360  } else if("G4EmPenelope" == typesPhys[i]) {
361  mod = new G4PenelopePhotoElectricModel();
362  em_config->SetExtraEmModel("gamma", "phot", mod, reg);
363  mod = new G4PenelopeComptonModel();
364  em_config->SetExtraEmModel("gamma", "compt", mod, reg);
365  mod = new G4PenelopeGammaConversionModel();
366  em_config->SetExtraEmModel("gamma", "conv", mod, reg);
367 
368  FindOrAddProcess(phot, "Rayl");
369  mod = new G4PenelopeRayleighModel();
370  em_config->SetExtraEmModel("gamma", "Rayl", mod, reg);
371 
372  mod = new G4PenelopeIonisationModel();
374  em_config->SetExtraEmModel("e-", "eIoni", mod, reg, 0.0, highEnergyLimit, uf);
375  mod = new G4PenelopeBremsstrahlungModel();
376  em_config->SetExtraEmModel("e-", "eBrem", mod, reg, 0.0, highEnergyLimit);
377 
378  mod = new G4PenelopeIonisationModel();
379  uf = new G4UniversalFluctuation();
380  em_config->SetExtraEmModel("e+", "eIoni", mod, reg, 0.0, highEnergyLimit, uf);
381  mod = new G4PenelopeBremsstrahlungModel();
382  em_config->SetExtraEmModel("e+", "eBrem", mod, reg, 0.0, highEnergyLimit);
383  mod = new G4PenelopeAnnihilationModel();
384  em_config->SetExtraEmModel("e+", "annihil", mod, reg, 0.0, highEnergyLimit);
385 
387  theParameters->SetDeexActiveRegion(regnamesPhys[i], true, false, false);
388  }
390 
391  } else if("G4RadioactiveDecay" == typesPhys[i]) {
392 
396  }
398 
399  } else {
400  if(verbose > 0 && G4Threading::IsMasterThread()) {
401  G4cout << "### G4EmModelActivator::ActivateEmOptions WARNING: \n"
402  << " EM Physics configuration name <" << typesPhys[i]
403  << "> is not known - ignored" << G4endl;
404  }
405  }
406  }
407 }
408 
409 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
410 
412 {
413  const std::vector<G4String> regnamesPAI = theParameters->RegionsPAI();
414  G4int nreg = regnamesPAI.size();
415  if(0 == nreg) { return; }
416  G4int verbose = theParameters->Verbose() - 1;
417  if(verbose > 0) {
418  G4cout << "### G4EmModelActivator::ActivatePAI for " << nreg << " regions"
419  << G4endl;
420  }
421  const std::vector<G4String> particlesPAI = theParameters->ParticlesPAI();
422  const std::vector<G4String> typesPAI = theParameters->TypesPAI();
423 
424  const std::vector<G4VEnergyLossProcess*>& v = G4LossTableManager::Instance()
426  std::vector<G4VEnergyLossProcess*>::const_iterator itr;
427  G4RegionStore* regionStore = G4RegionStore::GetInstance();
428 
434 
435  for(G4int i = 0; i < nreg; ++i) {
436  const G4ParticleDefinition* p = nullptr;
437  if(particlesPAI[i] != "all") {
438  p = G4ParticleTable::GetParticleTable()->FindParticle(particlesPAI[i]);
439  if(!p) {
440  G4cout << "### WARNING: ActivatePAI::FindParticle fails to find "
441  << particlesPAI[i] << G4endl;
442  continue;
443  }
444  }
445  const G4Region* r = regionStore->GetRegion(regnamesPAI[i], false);
446  if(!r) {
447  G4cout << "### WARNING: ActivatePAI::GetRegion fails to find "
448  << regnamesPAI[i] << G4endl;
449  continue;
450  }
451 
452  G4String name = "hIoni";
453  if(p == elec || p == posi)
454  { name = "eIoni"; }
455  else if (p == mupl || p == mumi)
456  { name = "muIoni"; }
457  else if (p == gion)
458  { name = "ionIoni"; }
459 
460  //for(itr = v.begin(); itr != v.end(); ++itr) {
461  for(auto proc : v) {
462 
463  if(!proc->IsIonisationProcess()) { continue; }
464 
465  G4String namep = proc->GetProcessName();
466  if(p) {
467  if(name != namep) { continue; }
468  } else {
469  if(namep != "hIoni" && namep != "muIoni" &&
470  namep != "eIoni" && namep != "ionIoni")
471  { continue; }
472  }
473  G4VEmModel* em = nullptr;
474  G4VEmFluctuationModel* fm = nullptr;
475  if(typesPAI[i] == "PAIphoton" || typesPAI[i] == "pai_photon") {
476  G4PAIPhotModel* mod = new G4PAIPhotModel(p,"PAIPhotModel");
477  em = mod;
478  fm = mod;
479  } else {
480  G4PAIModel* mod = new G4PAIModel(p,"PAIModel");
481  em = mod;
482  fm = mod;
483  }
484  proc->AddEmModel(0, em, fm, r);
485  if(verbose > 0) {
486  G4cout << "### G4EmModelActivator: add <" << typesPAI[i]
487  << "> model for " << particlesPAI[i]
488  << " in the " << regnamesPAI[i] << G4endl;
489  }
490  }
491  }
492 }
493 
494 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
495 
497 {
498  const std::vector<G4String> regnamesME = theParameters->RegionsMicroElec();
499  G4int nreg = regnamesME.size();
500  if(0 == nreg)
501  {
502  return;
503  }
504  G4int verbose = theParameters->Verbose() - 1;
505  if(verbose > 0)
506  {
507  G4cout << "### G4EmModelActivator::ActivateMicroElec for " << nreg
508  << " regions" << G4endl;
509  }
511 
513  const G4ParticleDefinition* prot = G4Proton::Proton();
515  G4ProcessManager* eman = elec->GetProcessManager();
516  G4ProcessManager* pman = prot->GetProcessManager();
517  G4ProcessManager* iman = gion->GetProcessManager();
518 
519  G4bool emsc = HasMsc(eman);
520 
521  // MicroElec elastic is not active in the world
522  G4MicroElecElastic* eElasProc =
523  new G4MicroElecElastic("e-G4MicroElecElastic");
524  eman->AddDiscreteProcess(eElasProc);
525 
526  G4MicroElecInelastic* eInelProc =
527  new G4MicroElecInelastic("e-G4MicroElecInelastic");
528  eman->AddDiscreteProcess(eInelProc);
529 
530  G4MicroElecInelastic* pInelProc =
531  new G4MicroElecInelastic("p_G4MicroElecInelastic");
532  pman->AddDiscreteProcess(pInelProc);
533 
534  G4MicroElecInelastic* iInelProc =
535  new G4MicroElecInelastic("ion_G4MicroElecInelastic");
536  iman->AddDiscreteProcess(iInelProc);
537 
538  // start configuration of models
539  G4EmConfigurator* em_config = man->EmConfigurator();
540  G4VEmModel* mod;
541 
542  // limits for MicroElec applicability
543  G4double elowest = 16.7 * eV;
544  G4double elimel = 100 * MeV;
545  G4double elimin = 9 * MeV;
546  G4double pmin = 50 * keV;
547  G4double pmax = 99.9 * MeV;
548 
549  G4LowECapture* ecap = new G4LowECapture(elowest);
550  eman->AddDiscreteProcess(ecap);
551 
552  for(G4int i = 0; i < nreg; ++i)
553  {
554 
555  G4String reg = regnamesME[i];
556  G4cout << "### MicroElec models are activated for G4Region " << reg
557  << G4endl
558  << " Energy limits for e- elastic: " << elowest/eV << " eV - "
559  << elimel/MeV << " MeV"
560  << G4endl
561  << " Energy limits for e- inelastic: " << elowest/eV << " eV - "
562  << elimin/MeV << " MeV"
563  << G4endl
564  << " Energy limits for hadrons/ions: " << pmin/MeV << " MeV - "
565  << pmax/MeV << " MeV"
566  << G4endl;
567 
568  // e-
569  if(emsc)
570  {
571  G4UrbanMscModel* msc = new G4UrbanMscModel();
572  msc->SetActivationLowEnergyLimit(elimel);
573  em_config->SetExtraEmModel("e-", "msc", msc, reg);
574  }
575  else
576  {
577  mod = new G4DummyModel();
578  em_config->SetExtraEmModel("e-", "CoulombScat", mod, reg, 0.0, elimel);
579  }
580 
581  mod = new G4MicroElecElasticModel();
582  em_config->SetExtraEmModel("e-",
583  "e-G4MicroElecElastic",
584  mod,
585  reg,
586  elowest,
587  elimin);
588 
589  mod = new G4MollerBhabhaModel();
590  mod->SetActivationLowEnergyLimit(elimin);
591  em_config->SetExtraEmModel("e-",
592  "eIoni",
593  mod,
594  reg,
595  0.0,
596  10 * TeV,
597  new G4UniversalFluctuation());
598 
599  mod = new G4MicroElecInelasticModel();
600  em_config->SetExtraEmModel("e-",
601  "e-G4MicroElecInelastic",
602  mod,
603  reg,
604  elowest,
605  elimin);
606 
607  // proton
608  mod = new G4BraggModel();
609  mod->SetActivationHighEnergyLimit(pmin);
610  em_config->SetExtraEmModel("proton",
611  "hIoni",
612  mod,
613  reg,
614  0.0,
615  2 * MeV,
616  new G4UniversalFluctuation());
617 
618  mod = new G4BetheBlochModel();
619  mod->SetActivationLowEnergyLimit(pmax);
620  em_config->SetExtraEmModel("proton",
621  "hIoni",
622  mod,
623  reg,
624  2 * MeV,
625  10 * TeV,
626  new G4UniversalFluctuation());
627 
628  mod = new G4MicroElecInelasticModel();
629  em_config->SetExtraEmModel("proton",
630  "p_G4MicroElecInelastic",
631  mod,
632  reg,
633  pmin,
634  pmax);
635 
636  // ions
637  mod = new G4BraggIonModel();
638  mod->SetActivationHighEnergyLimit(pmin);
639  em_config->SetExtraEmModel("GenericIon",
640  "ionIoni",
641  mod,
642  reg,
643  0.0,
644  2 * MeV,
645  new G4IonFluctuations());
646 
647  mod = new G4BetheBlochModel();
648  mod->SetActivationLowEnergyLimit(pmax);
649  em_config->SetExtraEmModel("GenericIon",
650  "ionIoni",
651  mod,
652  reg,
653  2 * MeV,
654  10 * TeV,
655  new G4IonFluctuations());
656 
657  mod = new G4MicroElecInelasticModel();
658  em_config->SetExtraEmModel("GenericIon",
659  "ion_G4MicroElecInelastic",
660  mod,
661  reg,
662  pmin,
663  pmax);
664  }
665 }
666 
667 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
668 
670 {
671  G4bool res = false;
672  G4ProcessVector* pv = pm->GetProcessList();
673  G4int nproc = pm->GetProcessListLength();
674  for(G4int i = 0; i < nproc; ++i)
675  {
676  if(((*pv)[i])->GetProcessSubType() == fMultipleScattering)
677  {
678  res = true;
679  break;
680  }
681  }
682  return res;
683 }
684 
686  G4EmConfigurator* em_config,
687  G4VMscModel* mscmod,
688  const G4String& reg,
689  G4double e1, G4double e2)
690 {
691  G4String pname = part->GetParticleName();
692 
693  // low-energy msc model
694  em_config->SetExtraEmModel(pname, "msc", mscmod, reg, 0.0, e1);
695 
696  // high energy msc model
697  G4WentzelVIModel* msc = new G4WentzelVIModel();
698  msc->SetLocked(true);
699  em_config->SetExtraEmModel(pname, "msc", msc, reg, e1, e2);
700 
701  // high energy single scattering model
702  FindOrAddProcess(part, "CoulombScat");
705  mod->SetLocked(true);
706  em_config->SetExtraEmModel(pname, "CoulombScat", mod, reg, 0.0, e2);
707 }
708 
710  const G4String& name)
711 {
712  G4ProcessManager* pm = part->GetProcessManager();
713  G4ProcessVector* pv = pm->GetProcessList();
714  G4int nproc = pm->GetProcessListLength();
715  for(G4int i = 0; i<nproc; ++i) {
716  if(((*pv)[i])->GetProcessName() == name) { return; }
717  }
718  if(name == "CoulombScat") {
720  cs->SetEmModel(new G4DummyModel());
721  pm->AddDiscreteProcess(cs);
722  } else if(name == "Rayl") {
724  rs->SetEmModel(new G4DummyModel());
725  pm->AddDiscreteProcess(rs);
726  }
727 }
728 
729 
730 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
void SetDeexcitationIgnoreCut(G4bool val)
const std::vector< G4String > & PartNames() const
const XML_Char * name
Definition: expat.h:151
const std::vector< G4String > & RegionsPhysics() const
static G4ParticleTable * GetParticleTable()
G4VAtomDeexcitation * AtomDeexcitation()
static constexpr double MeV
Definition: G4SIunits.hh:214
void SetHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:707
G4int GetProcessListLength() const
G4int AddDiscreteProcess(G4VProcess *aProcess, G4int ord=ordDefault)
static constexpr double keV
Definition: G4SIunits.hh:216
void SetLocked(G4bool)
Definition: G4VEmModel.hh:822
#define G4endl
Definition: G4ios.hh:61
const char * p
Definition: xmltok.h:285
const G4String & GetParticleName() const
void SetActivationHighEnergyLimit(G4double)
Definition: G4VEmModel.hh:721
G4double GetPDGCharge() const
G4bool IsMasterThread()
Definition: G4Threading.cc:130
void SetStepFunctionMuHad(G4double v1, G4double v2)
static G4MuonMinus * MuonMinus()
Definition: G4MuonMinus.cc:100
static G4Proton * Proton()
Definition: G4Proton.cc:93
void SetExtraEmModel(const G4String &particleName, const G4String &processName, G4VEmModel *, const G4String &regionName="", G4double emin=0.0, G4double emax=DBL_MAX, G4VEmFluctuationModel *fm=nullptr)
G4ProcessVector * GetProcessList() const
const std::vector< G4VEnergyLossProcess * > & GetEnergyLossProcessVector()
const std::vector< G4String > & RegionsPAI() const
static constexpr double um
Definition: G4SIunits.hh:113
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
const std::vector< G4String > & RegionsMicroElec() const
static G4Gamma * Gamma()
Definition: G4Gamma.cc:86
void AddStandardScattering(const G4ParticleDefinition *, G4EmConfigurator *, G4VMscModel *, const G4String &, G4double, G4double)
static constexpr double TeV
Definition: G4SIunits.hh:218
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
TString part[npart]
Definition: Style.C:32
G4EmParameters * theParameters
void SetEmModel(G4VEmModel *, G4int index=0)
void SetActivationLowEnergyLimit(G4double)
Definition: G4VEmModel.hh:728
void SetStepFunction(G4double v1, G4double v2)
const std::vector< G4String > & TypesPAI() const
void SetStepLimitType(G4MscStepLimitType)
Definition: G4VMscModel.hh:223
void SetPolarAngleLimit(G4double)
Definition: G4VEmModel.hh:742
static G4Positron * Positron()
Definition: G4Positron.cc:94
static constexpr double eV
Definition: G4SIunits.hh:215
G4int Verbose() const
static G4Electron * Electron()
Definition: G4Electron.cc:94
void SetRangeFactor(G4double)
Definition: G4VMscModel.hh:209
G4bool HasMsc(G4ProcessManager *) const
void FindOrAddProcess(const G4ParticleDefinition *, const G4String &)
const std::vector< G4String > & TypesMsc() const
void SetAugerCascade(G4bool val)
static G4RegionStore * GetInstance()
int G4int
Definition: G4Types.hh:78
G4ProcessManager * GetProcessManager() const
static G4MuonPlus * MuonPlus()
Definition: G4MuonPlus.cc:99
G4EmModelActivator(const G4String &emphys="")
G4GLOB_DLL std::ostream G4cout
void DefineRegParamForDeex(G4VAtomDeexcitation *) const
static G4LossTableManager * Instance()
static const G4double reg
void SetDeexActiveRegion(const G4String &region, G4bool fdeex, G4bool fauger, G4bool fpixe)
G4EmConfigurator * EmConfigurator()
G4double MscEnergyLimit() const
const std::vector< G4String > & RegionsMsc() const
static constexpr double GeV
Definition: G4SIunits.hh:217
const std::vector< G4String > & ParticlesPAI() const
static G4EmParameters * Instance()