ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4EmDNAPhysicsActivator.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4EmDNAPhysicsActivator.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 // add elastic scattering processes of proton, hydrogen, helium, alpha+, alpha++
27 
29 
30 #include "G4EmParameters.hh"
31 #include "G4ParticleDefinition.hh"
32 #include "G4ParticleTable.hh"
33 #include "G4RegionStore.hh"
34 #include "G4Region.hh"
35 #include "G4VEnergyLossProcess.hh"
36 #include "G4LossTableManager.hh"
37 #include "G4EmConfigurator.hh"
38 
39 #include "G4Gamma.hh"
40 #include "G4Electron.hh"
41 #include "G4Positron.hh"
42 #include "G4MuonPlus.hh"
43 #include "G4MuonMinus.hh"
44 #include "G4Proton.hh"
45 #include "G4GenericIon.hh"
46 #include "G4Alpha.hh"
47 
48 #include "G4ProcessManager.hh"
49 #include "G4DummyModel.hh"
50 #include "G4EmProcessSubType.hh"
51 #include "G4PhysicsListHelper.hh"
52 
53 #include "G4BraggModel.hh"
54 #include "G4BraggIonModel.hh"
55 #include "G4BetheBlochModel.hh"
56 #include "G4UrbanMscModel.hh"
57 #include "G4WentzelVIModel.hh"
58 #include "G4MollerBhabhaModel.hh"
59 #include "G4IonFluctuations.hh"
61 #include "G4LowECapture.hh"
62 #include "G4hMultipleScattering.hh"
65 #include "G4hIonisation.hh"
67 
68 // Processes and models for Geant4-DNA
70 
71 #include "G4DNAElastic.hh"
73 //#include "G4DNAScreenedRutherfordElasticModel.hh"
75 #include "G4DNAIonElasticModel.hh"
76 
77 #include "G4DNAExcitation.hh"
78 #include "G4DNAAttachment.hh"
79 #include "G4DNAVibExcitation.hh"
80 #include "G4DNAIonisation.hh"
81 #include "G4DNAChargeDecrease.hh"
82 #include "G4DNAChargeIncrease.hh"
83 
84 #include "G4SystemOfUnits.hh"
85 #include <vector>
86 
87 #include "G4DNAChemistryManager.hh"
91 
92 #include "G4Threading.hh"
93 
94 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
95 
97  : G4VPhysicsConstructor("G4EmDNAPhysicsActivator"), verbose(ver)
98 {
101 }
102 
103 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
104 
106 {}
107 
109 {
110  return (0 < verbose && G4Threading::IsMasterThread());
111 }
112 
113 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
114 
116 {
117 // bosons
118  G4Gamma::Gamma();
119 
120 // leptons
123 
124 // baryons
126 
128  G4Alpha::Alpha();
129 
130  G4DNAGenericIonsManager * genericIonsManager;
131  genericIonsManager=G4DNAGenericIonsManager::Instance();
132  genericIonsManager->GetIon("alpha+");
133  genericIonsManager->GetIon("helium");
134  genericIonsManager->GetIon("hydrogen");
135  //genericIonsManager->GetIon("carbon");
136  //genericIonsManager->GetIon("nitrogen");
137  //genericIonsManager->GetIon("oxygen");
138  //genericIonsManager->GetIon("iron");
139 
140 }
141 
142 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
143 
145 {
146  const std::vector<G4String>& regnamesDNA = theParameters->RegionsDNA();
147  G4int nreg = regnamesDNA.size();
148  if(0 == nreg)
149  {
150  return;
151  }
152  const std::vector<G4String>& typesDNA = theParameters->TypesDNA();
153 
154  if(IsVerbose()) {
155  G4cout << "### G4EmDNAPhysicsActivator::ConstructProcess for " << nreg
156  << " regions; DNA physics type " << typesDNA[0] << G4endl;
157  }
158 
159  // list of particles
161  const G4ParticleDefinition* prot = G4Proton::Proton();
163 
164  G4DNAGenericIonsManager * genericIonsManager =
166  const G4ParticleDefinition* alpha2 = G4Alpha::Alpha();
167  const G4ParticleDefinition* alpha1 = genericIonsManager->GetIon("alpha+");
168  const G4ParticleDefinition* alpha0 = genericIonsManager->GetIon("helium");
169  const G4ParticleDefinition* h0 = genericIonsManager->GetIon("hydrogen");
170 
171  G4ProcessManager* eman = elec->GetProcessManager();
172  G4ProcessManager* pman = prot->GetProcessManager();
173  G4ProcessManager* iman = gion->GetProcessManager();
174  G4ProcessManager* a2man = alpha2->GetProcessManager();
175  G4ProcessManager* a1man = alpha1->GetProcessManager();
176  G4ProcessManager* a0man = alpha0->GetProcessManager();
177  G4ProcessManager* h0man = h0->GetProcessManager();
178 
179  // alpha+ standard processes
181  G4ParticleDefinition* alpha11 = const_cast<G4ParticleDefinition*>(alpha1);
182  ph->RegisterProcess(new G4hMultipleScattering(), alpha11);
183  ph->RegisterProcess(new G4hIonisation(), alpha11);
184 
185  G4bool emsc = HasMsc(eman);
186  G4bool pmsc = HasMsc(pman);
187  G4bool a2msc = HasMsc(a2man);
188  G4bool a1msc = HasMsc(a1man);
189  // G4bool imsc = HasMsc(iman);
190 
191  // processes are defined with dummy models for the world
192  // elastic scatetring
193  G4DNAElastic* theDNAeElasticProcess = new G4DNAElastic("e-_G4DNAElastic");
194  theDNAeElasticProcess->SetEmModel(new G4DummyModel());
195  eman->AddDiscreteProcess(theDNAeElasticProcess);
196 
197  G4DNAElastic* theDNApElasticProcess =
198  new G4DNAElastic("proton_G4DNAElastic");
199  theDNApElasticProcess->SetEmModel(new G4DummyModel());
200  pman->AddDiscreteProcess(theDNApElasticProcess);
201 
202  G4DNAElastic* theDNAa2ElasticProcess =
203  new G4DNAElastic("alpha_G4DNAElastic");
204  theDNAa2ElasticProcess->SetEmModel(new G4DummyModel());
205  a2man->AddDiscreteProcess(theDNAa2ElasticProcess);
206 
207  G4DNAElastic* theDNAa1ElasticProcess =
208  new G4DNAElastic("alpha+_G4DNAElastic");
209  theDNAa1ElasticProcess->SetEmModel(new G4DummyModel());
210  a1man->AddDiscreteProcess(theDNAa1ElasticProcess);
211 
212  G4DNAElastic* theDNAa0ElasticProcess =
213  new G4DNAElastic("helium_G4DNAElastic");
214  theDNAa0ElasticProcess->SetEmModel(new G4DummyModel());
215  a0man->AddDiscreteProcess(theDNAa0ElasticProcess);
216 
217  G4DNAElastic* theDNAh0ElasticProcess =
218  new G4DNAElastic("hydrogen_G4DNAElastic");
219  theDNAh0ElasticProcess->SetEmModel(new G4DummyModel());
220  h0man->AddDiscreteProcess(theDNAh0ElasticProcess);
221 
222  // excitation
223  G4DNAExcitation* theDNAeExcProcess =
224  new G4DNAExcitation("e-_G4DNAExcitation");
225  theDNAeExcProcess->SetEmModel(new G4DummyModel());
226  eman->AddDiscreteProcess(theDNAeExcProcess);
227 
228  G4DNAExcitation* theDNApExcProcess =
229  new G4DNAExcitation("proton_G4DNAExcitation");
230  theDNApExcProcess->SetEmModel(new G4DummyModel());
231  pman->AddDiscreteProcess(theDNApExcProcess);
232 
233  G4DNAExcitation* theDNAa2ExcProcess =
234  new G4DNAExcitation("alpha_G4DNAExcitation");
235  theDNAa2ExcProcess->SetEmModel(new G4DummyModel());
236  a2man->AddDiscreteProcess(theDNAa2ExcProcess);
237 
238  G4DNAExcitation* theDNAa1ExcProcess =
239  new G4DNAExcitation("alpha+_G4DNAExcitation");
240  theDNAa1ExcProcess->SetEmModel(new G4DummyModel());
241  a1man->AddDiscreteProcess(theDNAa1ExcProcess);
242 
243  G4DNAExcitation* theDNAa0ExcProcess =
244  new G4DNAExcitation("helium_G4DNAExcitation");
245  theDNAa0ExcProcess->SetEmModel(new G4DummyModel());
246  a0man->AddDiscreteProcess(theDNAa0ExcProcess);
247 
248  G4DNAExcitation* theDNAh0ExcProcess =
249  new G4DNAExcitation("hydrogen_G4DNAExcitation");
250  theDNAh0ExcProcess->SetEmModel(new G4DummyModel());
251  h0man->AddDiscreteProcess(theDNAh0ExcProcess);
252 
253  // vibration excitation
254  G4DNAVibExcitation* theDNAeVibExcProcess =
255  new G4DNAVibExcitation("e-_G4DNAVibExcitation");
256  theDNAeVibExcProcess->SetEmModel(new G4DummyModel());
257  eman->AddDiscreteProcess(theDNAeVibExcProcess);
258 
259  // ionisation
260  G4DNAIonisation* theDNAeIoniProcess =
261  new G4DNAIonisation("e-_G4DNAIonisation");
262  theDNAeIoniProcess->SetEmModel(new G4DummyModel());
263  eman->AddDiscreteProcess(theDNAeIoniProcess);
264 
265  G4DNAIonisation* theDNApIoniProcess =
266  new G4DNAIonisation("proton_G4DNAIonisation");
267  theDNApIoniProcess->SetEmModel(new G4DummyModel());
268  pman->AddDiscreteProcess(theDNApIoniProcess);
269 
270  G4DNAIonisation* theDNAa2IoniProcess =
271  new G4DNAIonisation("alpha_G4DNAIonisation");
272  theDNAa2IoniProcess->SetEmModel(new G4DummyModel());
273  a2man->AddDiscreteProcess(theDNAa2IoniProcess);
274 
275  G4DNAIonisation* theDNAa1IoniProcess =
276  new G4DNAIonisation("alpha+_G4DNAIonisation");
277  theDNAa1IoniProcess->SetEmModel(new G4DummyModel());
278  a1man->AddDiscreteProcess(theDNAa1IoniProcess);
279 
280  G4DNAIonisation* theDNAa0IoniProcess =
281  new G4DNAIonisation("helium_G4DNAIonisation");
282  theDNAa0IoniProcess->SetEmModel(new G4DummyModel());
283  a0man->AddDiscreteProcess(theDNAa0IoniProcess);
284 
285  G4DNAIonisation* theDNAh0IoniProcess =
286  new G4DNAIonisation("hydrogen_G4DNAIonisation");
287  theDNAh0IoniProcess->SetEmModel(new G4DummyModel());
288  h0man->AddDiscreteProcess(theDNAh0IoniProcess);
289 
290  G4DNAIonisation* theDNAiIoniProcess =
291  new G4DNAIonisation("GenericIon_G4DNAIonisation");
292  theDNAiIoniProcess->SetEmModel(new G4DummyModel());
293  iman->AddDiscreteProcess(theDNAiIoniProcess);
294 
295  // attachment
296  G4DNAAttachment* theDNAAttachProcess =
297  new G4DNAAttachment("e-_G4DNAAttachment");
298  theDNAAttachProcess->SetEmModel(new G4DummyModel());
299  eman->AddDiscreteProcess(theDNAAttachProcess);
300 
301  // charge exchange
302  G4DNAChargeDecrease* theDNApChargeDecreaseProcess =
303  new G4DNAChargeDecrease("proton_G4DNAChargeDecrease");
304  theDNApChargeDecreaseProcess->SetEmModel(new G4DummyModel());
305  pman->AddDiscreteProcess(theDNApChargeDecreaseProcess);
306 
307  G4DNAChargeDecrease* theDNAa2ChargeDecreaseProcess =
308  new G4DNAChargeDecrease("alpha_G4DNAChargeDecrease");
309  theDNAa2ChargeDecreaseProcess->SetEmModel(new G4DummyModel());
310  a2man->AddDiscreteProcess(theDNAa2ChargeDecreaseProcess);
311 
312  G4DNAChargeDecrease* theDNAa1ChargeDecreaseProcess =
313  new G4DNAChargeDecrease("alpha+_G4DNAChargeDecrease");
314  theDNAa1ChargeDecreaseProcess->SetEmModel(new G4DummyModel());
315  a1man->AddDiscreteProcess(theDNAa1ChargeDecreaseProcess);
316 
317  G4DNAChargeIncrease* theDNAa1ChargeIncreaseProcess =
318  new G4DNAChargeIncrease("alpha+_G4DNAChargeIncrease");
319  theDNAa1ChargeIncreaseProcess->SetEmModel(new G4DummyModel());
320  a1man->AddDiscreteProcess(theDNAa1ChargeIncreaseProcess);
321 
322  G4DNAChargeIncrease* theDNAa0ChargeIncreaseProcess =
323  new G4DNAChargeIncrease("helium_G4DNAChargeIncrease");
324  theDNAa0ChargeIncreaseProcess->SetEmModel(new G4DummyModel());
325  a0man->AddDiscreteProcess(theDNAa0ChargeIncreaseProcess);
326 
327  G4DNAChargeIncrease* theDNAh0ChargeIncreaseProcess =
328  new G4DNAChargeIncrease("hydrogen_G4DNAChargeIncrease");
329  theDNAh0ChargeIncreaseProcess->SetEmModel(new G4DummyModel());
330  h0man->AddDiscreteProcess(theDNAh0ChargeIncreaseProcess);
331 
332  // limits for DNA model applicability
333  static const G4double elowest= 7.4 * eV;
334  static const G4double elimel = 1 * MeV;
335  static const G4double pminbb = 2 * MeV;
336  static const G4double pmin = 0.1 * keV;
337  static const G4double pmax = 100 * MeV;
338  static const G4double hemin = 1 * keV;
339  static const G4double ionmin = 0.5 * MeV;
340 
341  // low-energy capture
342  G4LowECapture* ecap = nullptr;
344  {
345  // Note: G4DNAElectronSolvation could also be used instead of G4LowECapture
346  ecap = new G4LowECapture(elowest);
347  // ecap->SetVerboseLevel(1);
348  eman->AddDiscreteProcess(ecap);
349  }
350  else
351  {
352  // When chemistry is activated: G4DNAElectronSolvation turns the electron
353  // to a solvated electron, otherwise it kills the electron at the
354  // corresponding high energy limit of the model
355  G4DNAElectronSolvation* pSolvatation =
356  new G4DNAElectronSolvation("e-_G4DNAElectronSolvation");
358  eman->AddDiscreteProcess(pSolvatation);
359  }
360 
361  G4LowECapture* pcap = new G4LowECapture(pmin);
362  pman->AddDiscreteProcess(pcap);
363  G4LowECapture* icap = new G4LowECapture(ionmin);
364  iman->AddDiscreteProcess(icap);
365  G4LowECapture* a2cap = new G4LowECapture(hemin);
366  a2man->AddDiscreteProcess(a2cap);
367  G4LowECapture* a1cap = new G4LowECapture(hemin);
368  a1man->AddDiscreteProcess(a1cap);
369  G4LowECapture* a0cap = new G4LowECapture(hemin);
370  a0man->AddDiscreteProcess(a0cap);
371  G4LowECapture* h0cap = new G4LowECapture(ionmin);
372  h0man->AddDiscreteProcess(h0cap);
373 
374  // loop over regions
375  for(G4int i = 0; i < nreg; ++i)
376  {
377  G4String reg = regnamesDNA[i];
378  if(IsVerbose())
379  {
380  G4cout << "### DNA models type " << typesDNA[i]
381  << " are activated for G4Region " << reg << G4endl;
382  }
383 
384  // type of DNA physics
385  G4int itype = 0;
386 
387  if(0 == itype) {
388  AddElectronModels0(reg, ecap, emsc, elowest, elimel);
389  AddProtonModels0(reg, pmsc, elimel, pminbb, pmax);
390  AddHeliumModels0(reg, a1msc, a2msc, elimel, pminbb, pmax);
391  AddGenericIonModels0(reg, pminbb);
392  DeactivateNuclearStopping(pman, elimel);
393  DeactivateNuclearStopping(a1man, elimel);
394  DeactivateNuclearStopping(a2man, elimel);
395  }
396  }
398 }
399 
400 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
401 
403  G4LowECapture* ecap,
404  G4bool emsc,
405  G4double elowest,
406  G4double elimel)
407 {
408  G4EmConfigurator* em_config =
410  G4VEmModel* mod;
411 
412  static const G4double elimin = 1 * MeV;
413  static const G4double elimvb = 100 * eV;
414  static const G4double elimat = 13 * eV;
415  static const G4double elim1 = 10 * keV;
416 
417  // for e- 100 MeV is a limit between different msc models
419 
420  if(emsc) {
421  G4UrbanMscModel* msc = new G4UrbanMscModel();
422  msc->SetActivationLowEnergyLimit(elimel);
423  G4double emaxmsc = std::min(100*MeV, emax);
424  em_config->SetExtraEmModel("e-", "msc", msc, reg, 0.0, emaxmsc);
425  } else {
426  mod = new G4eCoulombScatteringModel();
427  mod->SetActivationLowEnergyLimit(elimel);
428  em_config->SetExtraEmModel("e-", "CoulombScat", mod, reg, 0.0, emax);
429  }
430 
431  // cuts and solvation
432  if(ecap) {
433  ecap->AddRegion(reg);
434  } else {
435  // For this condition to work, the chemistry list has to be called
436  // before any standard EM physics list
438  em_config->SetExtraEmModel("e-", "e-_G4DNAElectronSolvation",
439  mod, reg, 0., elowest);
440  }
441 
442  // elastic
443  mod = new G4DNAChampionElasticModel();
444  em_config->SetExtraEmModel("e-", "e-_G4DNAElastic",
445  mod, reg, 0.0, elimel);
446  // ionisation
447  mod = new G4MollerBhabhaModel();
448  mod->SetActivationLowEnergyLimit(elimin);
449  em_config->SetExtraEmModel("e-", "eIoni",
450  mod, reg, 0.0, emax,
451  new G4UniversalFluctuation());
452 
453  mod = new G4DNABornIonisationModel();
454  em_config->SetExtraEmModel("e-", "e-_G4DNAIonisation",
455  mod, reg, elim1, elimin);
456 
458  em_config->SetExtraEmModel("e-", "e-_G4DNAIonisation",
459  mod, reg, 0.0, elim1);
460 
461  // exc
463  em_config->SetExtraEmModel("e-", "e-_G4DNAExcitation",
464  mod, reg, 0.0, elim1);
465 
466  mod = new G4DNABornExcitationModel();
467  em_config->SetExtraEmModel("e-", "e-_G4DNAExcitation",
468  mod, reg, elim1, elimin);
469 
470  mod = new G4DNASancheExcitationModel();
471  em_config->SetExtraEmModel("e-", "e-_G4DNAVibExcitation",
472  mod, reg, 0.0, elimvb);
473 
474  mod = new G4DNAMeltonAttachmentModel();
475  em_config->SetExtraEmModel("e-", "e-_G4DNAAttachment",
476  mod, reg, 0.0, elimat);
477 
478 }
479 
480 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
481 
483  G4bool pmsc, G4double elimel,
484  G4double pminbb, G4double pmax)
485 {
486  G4EmConfigurator* em_config =
488  G4VEmModel* mod;
489 
490  static const G4double gmmax = 500 * keV;
491 
493 
494  // proton
495 
496  // if SS physics list msc process does not exist
497  if(pmsc) {
498  G4WentzelVIModel* msc = new G4WentzelVIModel();
499  msc->SetActivationLowEnergyLimit(elimel);
500  em_config->SetExtraEmModel("proton", "msc", msc, reg, 0.0, emax);
501  }
502  // single scattering always applied
503  mod = new G4eCoulombScatteringModel();
504  mod->SetActivationLowEnergyLimit(elimel);
505  em_config->SetExtraEmModel("proton", "CoulombScat", mod, reg, 0.0, emax);
506 
507  mod = new G4BraggModel();
508  mod->SetActivationLowEnergyLimit(std::min(pminbb, pmax));
509  em_config->SetExtraEmModel("proton", "hIoni",
510  mod, reg, 0.0, pminbb,
511  new G4UniversalFluctuation());
512 
513  mod = new G4BetheBlochModel();
514  mod->SetActivationLowEnergyLimit(pmax);
515  em_config->SetExtraEmModel("proton", "hIoni",
516  mod, reg, pminbb, emax,
517  new G4UniversalFluctuation());
518 
519  mod = new G4DNARuddIonisationModel();
520  em_config->SetExtraEmModel("proton", "proton_G4DNAIonisation",
521  mod, reg, 0.0, gmmax);
522 
523  mod = new G4DNABornIonisationModel();
524  em_config->SetExtraEmModel("proton", "proton_G4DNAIonisation",
525  mod, reg, gmmax, pmax);
526 
528  em_config->SetExtraEmModel("proton", "proton_G4DNAExcitation",
529  mod, reg, 0.0, gmmax);
530 
531  mod = new G4DNABornExcitationModel();
532  em_config->SetExtraEmModel("proton", "proton_G4DNAExcitation",
533  mod, reg, gmmax, pmax);
534 
536  em_config->SetExtraEmModel("proton", "proton_G4DNAChargeDecrease",
537  mod, reg, 0.0, pmax);
538 
539  mod = new G4DNAIonElasticModel();
540  em_config->SetExtraEmModel("proton", "proton_G4DNAElastic",
541  mod, reg, 0.0, elimel);
542 
543  // hydrogen
544  mod = new G4DNARuddIonisationModel();
545  em_config->SetExtraEmModel("hydrogen", "hydrogen_G4DNAIonisation",
546  mod, reg, 0.0, pmax);
547 
549  em_config->SetExtraEmModel("hydrogen", "hydrogen_G4DNAExcitation",
550  mod, reg, 0.0, gmmax);
551 
553  em_config->SetExtraEmModel("hydrogen", "hydrogen_G4DNAChargeIncrease",
554  mod, reg, 0.0, pmax);
555 
556  mod = new G4DNAIonElasticModel();
557  em_config->SetExtraEmModel("hydrogen", "hydrogen_G4DNAElastic",
558  mod, reg, 0.0, elimel);
559 
560 }
561 
562 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
563 
565  G4double pminbb)
566 {
567  G4EmConfigurator* em_config =
569  G4VEmModel* mod;
570 
572  G4double iemax = std::min(10*MeV, emax);
573  //G4double iemin = 100*eV;
574 
575  mod = new G4BraggIonModel();
576  mod->SetActivationLowEnergyLimit(iemax);
577  em_config->SetExtraEmModel("GenericIon", "ionIoni",
578  mod, reg, 0.0, pminbb,
579  new G4IonFluctuations());
580 
581  mod = new G4BetheBlochModel();
582  mod->SetActivationLowEnergyLimit(iemax);
583  em_config->SetExtraEmModel("GenericIon", "ionIoni",
584  mod, reg, pminbb, emax,
585  new G4IonFluctuations());
586 
588  em_config->SetExtraEmModel("GenericIon", "GenericIon_G4DNAIonisation",
589  mod, reg, 0.0, iemax);
590 
591 }
592 
593 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
594 
596  G4bool a1msc,
597  G4bool a2msc, G4double elimel,
598  G4double pminbb, G4double)
599 {
600  G4EmConfigurator* em_config =
602  G4VEmModel* mod;
603 
604  static const G4double hemax = 400 * MeV;
605  static const G4double massRatio = G4Alpha::Alpha()->GetPDGMass()/CLHEP::proton_mass_c2;
606 
608  G4double pminbba = massRatio*pminbb;
609  if(IsVerbose()) {
610  G4cout << "AddHeliumModels0 for <" << reg << "> a1msc: " << a1msc <<" a2msc: " << a2msc
611  << " elimel= " << elimel << " pminbba= " << pminbba << G4endl;
612  }
613  // alpha++
614  if(elimel < emax) {
615  if(a2msc) {
616  G4UrbanMscModel* msc = new G4UrbanMscModel();
617  msc->SetActivationLowEnergyLimit(elimel);
618  em_config->SetExtraEmModel("alpha", "msc", msc, reg, 0.0, emax);
619  } else {
620  mod = new G4IonCoulombScatteringModel();
621  mod->SetActivationLowEnergyLimit(elimel);
622  em_config->SetExtraEmModel("alpha", "CoulombScat", mod, reg, 0.0, emax);
623  }
624  }
625 
626  mod = new G4BraggIonModel();
627  mod->SetActivationLowEnergyLimit(hemax/massRatio);
628  em_config->SetExtraEmModel("alpha", "ionIoni",
629  mod, reg, 0.0, pminbba,
630  new G4IonFluctuations());
631 
632  mod = new G4BetheBlochModel();
633  mod->SetActivationLowEnergyLimit(hemax/massRatio);
634  em_config->SetExtraEmModel("alpha", "ionIoni",
635  mod, reg, pminbba, emax,
636  new G4IonFluctuations());
637 
638  mod = new G4DNARuddIonisationModel();
639  em_config->SetExtraEmModel("alpha", "alpha_G4DNAIonisation",
640  mod, reg, 0.0, hemax);
641 
643  em_config->SetExtraEmModel("alpha", "alpha_G4DNAExcitation",
644  mod, reg, 0.0, hemax);
645 
647  em_config->SetExtraEmModel("alpha", "alpha_G4DNAChargeDecrease",
648  mod, reg, 0.0, hemax);
649 
650  mod = new G4DNAIonElasticModel();
651  em_config->SetExtraEmModel("alpha", "alpha_G4DNAElastic",
652  mod, reg, 0.0, elimel);
653 
654  // ---
655  // alpha+
656  if(elimel < emax) {
657  if(a1msc) {
658  G4UrbanMscModel* msc = new G4UrbanMscModel();
659  msc->SetActivationLowEnergyLimit(elimel);
660  em_config->SetExtraEmModel("alpha+", "msc", msc, reg, 0.0, emax);
661  } else {
662  mod = new G4IonCoulombScatteringModel();
663  mod->SetActivationLowEnergyLimit(elimel);
664  em_config->SetExtraEmModel("alpha+", "CoulombScat", mod, reg, 0.0, emax);
665  }
666  }
667 
668  mod = new G4BraggIonModel();
669  mod->SetActivationLowEnergyLimit(hemax/massRatio);
670  em_config->SetExtraEmModel("alpha+", "hIoni",
671  mod, reg, 0.0, pminbba,
672  new G4IonFluctuations());
673 
674  mod = new G4BetheBlochModel();
675  mod->SetActivationLowEnergyLimit(hemax/massRatio);
676  em_config->SetExtraEmModel("alpha+", "hIoni",
677  mod, reg, pminbba, emax,
678  new G4IonFluctuations());
679 
680  mod = new G4DNARuddIonisationModel();
681  em_config->SetExtraEmModel("alpha+", "alpha+_G4DNAIonisation",
682  mod, reg, 0.0, hemax);
683 
685  em_config->SetExtraEmModel("alpha+", "alpha+_G4DNAExcitation",
686  mod, reg, 0.0, hemax);
687 
689  em_config->SetExtraEmModel("alpha+", "alpha+_G4DNAChargeDecrease",
690  mod, reg, 0.0, hemax);
691 
693  em_config->SetExtraEmModel("alpha+", "alpha+_G4DNAChargeIncrease",
694  mod, reg, 0.0, hemax);
695 
696  mod = new G4DNAIonElasticModel();
697  em_config->SetExtraEmModel("alpha+", "alpha+_G4DNAElastic",
698  mod, reg, 0.0, elimel);
699 
700  // ---
701  // helium
702  mod = new G4DNARuddIonisationModel();
703  em_config->SetExtraEmModel("helium", "helium_G4DNAIonisation",
704  mod, reg, 0.0, hemax);
705 
707  em_config->SetExtraEmModel("helium", "helium_G4DNAExcitation",
708  mod, reg, 0.0, hemax);
709 
711  em_config->SetExtraEmModel("helium", "helium_G4DNAChargeIncrease",
712  mod, reg, 0.0, hemax);
713 
714  mod = new G4DNAIonElasticModel();
715  em_config->SetExtraEmModel("helium", "helium_G4DNAElastic",
716  mod, reg, 0.0, elimel);
717 }
718 
719 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
720 
722  G4double elimel)
723 {
724  G4ProcessVector* pv = pman->GetProcessList();
725  G4int nproc = pman->GetProcessListLength();
726  for(G4int i = 0; i < nproc; ++i) {
727  if(((*pv)[i])->GetProcessSubType() == fNuclearStopping) {
728  G4VEmProcess* proc = static_cast<G4VEmProcess*>((*pv)[i]);
729  if(proc) {
731  mod->SetActivationLowEnergyLimit(elimel);
732  proc->SetEmModel(mod);
733  }
734  break;
735  }
736  }
737 }
738 
739 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
740 
742 {
743  G4bool res = false;
744  G4ProcessVector* pv = pman->GetProcessList();
745  G4int nproc = pman->GetProcessListLength();
746  for(G4int i = 0; i < nproc; ++i)
747  {
748  if(((*pv)[i])->GetProcessSubType() == fMultipleScattering)
749  {
750  res = true;
751  break;
752  }
753  }
754  return res;
755 }
756 
757 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......