ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4EmExtraPhysics.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4EmExtraPhysics.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 //
27 //---------------------------------------------------------------------------
28 //
29 // ClassName: G4EmExtraPhysics
30 //
31 // Author: 2002 J.P. Wellisch
32 //
33 // Modified:
34 //
35 // 10.11.2005 V.Ivanchenko edit to provide a standard
36 // 19.06.2006 V.Ivanchenko add mu-nuclear process
37 // 16.10.2012 A.Ribon: renamed G4EmExtraBertiniPhysics as G4EmExtraPhysics
38 // 10.04.2014 A.Dotti: Add MT functionality for messenger
39 // 24.04.2014 A.Ribon: switched on muon-nuclear by default
40 // 29.01.2018 V.Grichine, adding neutrinos
41 // 07.05.2019 V.Grichine, adding muon neutrino nucleus interactions
42 //
43 //
45 
46 #include "G4EmExtraPhysics.hh"
47 
48 #include "G4SystemOfUnits.hh"
49 
50 #include "G4ParticleDefinition.hh"
51 #include "G4ParticleTable.hh"
52 #include "G4Gamma.hh"
53 #include "G4Electron.hh"
54 #include "G4Positron.hh"
55 #include "G4MuonPlus.hh"
56 #include "G4MuonMinus.hh"
57 #include "G4AntiNeutrinoE.hh"
58 #include "G4NeutrinoE.hh"
59 #include "G4AntiNeutrinoMu.hh"
60 #include "G4NeutrinoMu.hh"
61 #include "G4AntiNeutrinoTau.hh"
62 #include "G4NeutrinoTau.hh"
63 
64 #include "G4SynchrotronRadiation.hh"
65 #include "G4MuonNuclearProcess.hh"
66 #include "G4MuonVDNuclearModel.hh"
67 #include "G4ElectroVDNuclearModel.hh"
68 #include "G4TheoFSGenerator.hh"
69 #include "G4GeneratorPrecompoundInterface.hh"
70 #include "G4QGSModel.hh"
71 #include "G4GammaParticipants.hh"
72 #include "G4QGSMFragmentation.hh"
73 #include "G4ExcitedStringDecay.hh"
74 #include "G4CascadeInterface.hh"
75 
76 #include "G4LENDorBERTModel.hh"
77 #include "G4LENDCombinedCrossSection.hh"
78 
79 #include "G4GammaConversionToMuons.hh"
80 #include "G4AnnihiToMuPair.hh"
81 #include "G4eeToHadrons.hh"
82 
83 #include "G4PhotoNuclearProcess.hh"
84 #include "G4ElectronNuclearProcess.hh"
85 #include "G4PositronNuclearProcess.hh"
86 
87 #include "G4NeutrinoElectronProcess.hh"
88 #include "G4NeutrinoElectronTotXsc.hh"
89 #include "G4NeutrinoElectronCcModel.hh"
90 #include "G4NeutrinoElectronNcModel.hh"
91 #include "G4MuNeutrinoNucleusProcess.hh"
92 #include "G4MuNeutrinoNucleusTotXsc.hh"
93 #include "G4NuMuNucleusCcModel.hh"
94 #include "G4NuMuNucleusNcModel.hh"
95 #include "G4GammaGeneralProcess.hh"
96 #include "G4LossTableManager.hh"
97 
98 #include "G4HadronicParameters.hh"
99 #include "G4PhysicsListHelper.hh"
100 #include "G4BuilderType.hh"
101 
102 // factory
103 #include "G4PhysicsConstructorFactory.hh"
104 //
105 G4_DECLARE_PHYSCONSTR_FACTORY(G4EmExtraPhysics);
106 
108 
109 G4EmExtraPhysics::G4EmExtraPhysics(G4int ver):
110  G4VPhysicsConstructor("G4GammaLeptoNuclearPhys"),
111  gnActivated (true),
112  eActivated (true),
113  gLENDActivated(false),
114  munActivated(true),
115  synActivated(false),
116  synActivatedForAll(false),
117  gmumuActivated(false),
118  pmumuActivated(false),
119  phadActivated (false),
120  fNuActivated (false),
121  fNuETotXscActivated (false),
122  gmumuFactor (1.0),
123  pmumuFactor (1.0),
124  phadFactor (1.0),
125  fNuEleCcBias(1.0),
126  fNuEleNcBias(1.0),
127  fNuNucleusBias(1.0),
128  fNuDetectorName("0"),
129  verbose(ver)
130 {
131  theMessenger = new G4EmMessenger(this);
132  SetPhysicsType(bEmExtra);
133  if(verbose > 1) G4cout << "### G4EmExtraPhysics" << G4endl;
134 }
135 
136 G4EmExtraPhysics::G4EmExtraPhysics(const G4String&)
137  : G4EmExtraPhysics(1)
138 {}
139 
140 G4EmExtraPhysics::~G4EmExtraPhysics()
141 {
142  delete theMessenger;
143  theMessenger = nullptr;
144 }
145 
146 void G4EmExtraPhysics::Synch(G4bool val)
147 {
148  synActivated = val;
149 }
150 
151 void G4EmExtraPhysics::SynchAll(G4bool val)
152 {
153  synActivatedForAll = val;
154  if(synActivatedForAll) { synActivated = true; }
155 }
156 
157 void G4EmExtraPhysics::GammaNuclear(G4bool val)
158 {
159  gnActivated = val;
160 }
161 
162 void G4EmExtraPhysics::LENDGammaNuclear(G4bool val)
163 {
164  gLENDActivated = val;
165 }
166 
167 void G4EmExtraPhysics::ElectroNuclear(G4bool val)
168 {
169  eActivated = val;
170 }
171 
172 void G4EmExtraPhysics::MuonNuclear(G4bool val)
173 {
174  munActivated = val;
175 }
176 
177 void G4EmExtraPhysics::GammaToMuMu(G4bool val)
178 {
179  gmumuActivated = val;
180 }
181 
182 void G4EmExtraPhysics::PositronToMuMu(G4bool val)
183 {
184  pmumuActivated = val;
185 }
186 
187 void G4EmExtraPhysics::PositronToHadrons(G4bool val)
188 {
189  phadActivated = val;
190 }
191 
192 void G4EmExtraPhysics::GammaToMuMuFactor(G4double val)
193 {
194  if(val > 0.0) gmumuFactor = val;
195 }
196 
197 void G4EmExtraPhysics::PositronToMuMuFactor(G4double val)
198 {
199  if(val > 0.0) pmumuFactor = val;
200 }
201 
202 void G4EmExtraPhysics::PositronToHadronsFactor(G4double val)
203 {
204  if(val > 0.0) phadFactor = val;
205 }
206 
208 
209 void G4EmExtraPhysics::NeutrinoActivated(G4bool val)
210 {
211  fNuActivated = val;
212 }
213 
214 void G4EmExtraPhysics::NuETotXscActivated(G4bool val)
215 {
216  fNuETotXscActivated = val;
217 }
218 
219 void G4EmExtraPhysics::SetNuEleCcBias(G4double bf)
220 {
221  if(bf > 0.0) fNuEleCcBias = bf;
222 }
223 
224 void G4EmExtraPhysics::SetNuEleNcBias(G4double bf)
225 {
226  if(bf > 0.0) fNuEleNcBias = bf;
227 }
228 
229 void G4EmExtraPhysics::SetNuNucleusBias(G4double bf)
230 {
231  if(bf > 0.0) fNuNucleusBias = bf;
232 }
233 
234 void G4EmExtraPhysics::SetNuDetectorName(const G4String& dn)
235 {
236  fNuDetectorName = dn;
237 }
238 
240 
241 void G4EmExtraPhysics::ConstructParticle()
242 {
243  G4Gamma::Gamma();
244  G4Electron::Electron();
245  G4Positron::Positron();
246  G4MuonPlus::MuonPlus();
247  G4MuonMinus::MuonMinus();
248 
249  G4AntiNeutrinoE::AntiNeutrinoE();
250  G4NeutrinoE::NeutrinoE();
251  G4AntiNeutrinoMu::AntiNeutrinoMu();
252  G4NeutrinoMu::NeutrinoMu();
253  G4AntiNeutrinoTau::AntiNeutrinoTau();
254  G4NeutrinoTau::NeutrinoTau();
255 }
256 
257 void G4EmExtraPhysics::ConstructProcess()
258 {
259  G4ParticleDefinition* gamma = G4Gamma::Gamma();
260  G4ParticleDefinition* electron = G4Electron::Electron();
261  G4ParticleDefinition* positron = G4Positron::Positron();
262  G4ParticleDefinition* muonplus = G4MuonPlus::MuonPlus();
263  G4ParticleDefinition* muonminus = G4MuonMinus::MuonMinus();
264 
265  G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
266  G4LossTableManager* emManager = G4LossTableManager::Instance();
267 
268  if(gnActivated) { ConstructGammaElectroNuclear(); }
269 
270  if(munActivated) {
271  G4MuonNuclearProcess* muNucProcess = new G4MuonNuclearProcess();
272  G4MuonVDNuclearModel* muNucModel = new G4MuonVDNuclearModel();
273  muNucProcess->RegisterMe(muNucModel);
274  ph->RegisterProcess( muNucProcess, muonplus);
275  ph->RegisterProcess( muNucProcess, muonminus);
276  }
277  if(gmumuActivated) {
278  G4GammaConversionToMuons* theGammaToMuMu = new G4GammaConversionToMuons();
279  theGammaToMuMu->SetCrossSecFactor(gmumuFactor);
280  G4GammaGeneralProcess* sp =
281  (G4GammaGeneralProcess*)emManager->GetGammaGeneralProcess();
282  if(sp) {
283  sp->AddMMProcess(theGammaToMuMu);
284  } else {
285  ph->RegisterProcess(theGammaToMuMu, gamma);
286  }
287  }
288  if(pmumuActivated) {
289  G4AnnihiToMuPair* thePosiToMuMu = new G4AnnihiToMuPair();
290  thePosiToMuMu->SetCrossSecFactor(pmumuFactor);
291  ph->RegisterProcess(thePosiToMuMu, positron);
292  }
293  if(phadActivated) {
294  G4eeToHadrons* thePosiToHadrons = new G4eeToHadrons();
295  thePosiToHadrons->SetCrossSecFactor(phadFactor);
296  ph->RegisterProcess(thePosiToHadrons, positron);
297  }
298  if(synActivated) {
299  G4SynchrotronRadiation* theSynchRad = new G4SynchrotronRadiation();
300  ph->RegisterProcess( theSynchRad, electron);
301  ph->RegisterProcess( theSynchRad, positron);
302  if(synActivatedForAll) {
303  auto myParticleIterator=GetParticleIterator();
304  myParticleIterator->reset();
305  G4ParticleDefinition* particle = nullptr;
306 
307  while( (*myParticleIterator)() ) {
308  particle = myParticleIterator->value();
309  if( particle->GetPDGStable() && particle->GetPDGCharge() != 0.0) {
310  if(verbose > 1) {
311  G4cout << "### G4SynchrotronRadiation for "
312  << particle->GetParticleName() << G4endl;
313  }
314  ph->RegisterProcess( theSynchRad, particle);
315  }
316  }
317  }
318  }
319  if( fNuActivated )
320  {
321  G4ParticleDefinition* anuelectron = G4AntiNeutrinoE::AntiNeutrinoE();
322  G4ParticleDefinition* nuelectron = G4NeutrinoE::NeutrinoE();
323  G4ParticleDefinition* anumuon = G4AntiNeutrinoMu::AntiNeutrinoMu();
324  G4ParticleDefinition* numuon = G4NeutrinoMu::NeutrinoMu();
325  G4ParticleDefinition* anutau = G4AntiNeutrinoTau::AntiNeutrinoTau();
326  G4ParticleDefinition* nutau = G4NeutrinoTau::NeutrinoTau();
327 
328  G4NeutrinoElectronProcess* theNuEleProcess =
329  new G4NeutrinoElectronProcess(fNuDetectorName);
330  G4NeutrinoElectronTotXsc* theNuEleTotXsc = new G4NeutrinoElectronTotXsc();
331 
332  if(fNuETotXscActivated)
333  {
334  G4double bftot = std::max(fNuEleCcBias,fNuEleNcBias);
335  theNuEleProcess->SetBiasingFactor(bftot);
336  }
337  else
338  {
339  theNuEleProcess->SetBiasingFactors(fNuEleCcBias,fNuEleNcBias);
340  theNuEleTotXsc->SetBiasingFactors(fNuEleCcBias,fNuEleNcBias);
341  }
342  theNuEleProcess->AddDataSet(theNuEleTotXsc);
343 
344  G4NeutrinoElectronCcModel* ccModel = new G4NeutrinoElectronCcModel();
345  G4NeutrinoElectronNcModel* ncModel = new G4NeutrinoElectronNcModel();
346  theNuEleProcess->RegisterMe(ccModel);
347  theNuEleProcess->RegisterMe(ncModel);
348 
349  ph->RegisterProcess(theNuEleProcess, anuelectron);
350  ph->RegisterProcess(theNuEleProcess, nuelectron);
351  ph->RegisterProcess(theNuEleProcess, anumuon);
352  ph->RegisterProcess(theNuEleProcess, numuon);
353  ph->RegisterProcess(theNuEleProcess, anutau);
354  ph->RegisterProcess(theNuEleProcess, nutau);
355 
356  // nu_mu nucleus interactions
357  G4MuNeutrinoNucleusProcess* theNuMuNucleusProcess = new G4MuNeutrinoNucleusProcess(fNuDetectorName);
358  G4MuNeutrinoNucleusTotXsc* theNuMuNucleusTotXsc = new G4MuNeutrinoNucleusTotXsc();
359 
360  if(fNuETotXscActivated)
361  {
362  theNuMuNucleusProcess->SetBiasingFactor(fNuNucleusBias);
363  }
364  theNuMuNucleusProcess->AddDataSet(theNuMuNucleusTotXsc);
365 
366  G4NuMuNucleusCcModel* numunuclcc = new G4NuMuNucleusCcModel();
367  G4NuMuNucleusNcModel* numunuclnc = new G4NuMuNucleusNcModel();
368 
369  theNuMuNucleusProcess->RegisterMe(numunuclcc);
370  theNuMuNucleusProcess->RegisterMe(numunuclnc);
371 
372  ph->RegisterProcess(theNuMuNucleusProcess, anumuon);
373  ph->RegisterProcess(theNuMuNucleusProcess, numuon);
374  }
375 }
376 
377 void G4EmExtraPhysics::ConstructGammaElectroNuclear()
378 {
379  G4LossTableManager* emManager = G4LossTableManager::Instance();
380  G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
381 
382  G4PhotoNuclearProcess* gnuc = new G4PhotoNuclearProcess();
383 
384  G4QGSModel< G4GammaParticipants >* theStringModel =
385  new G4QGSModel< G4GammaParticipants >;
386  G4QGSMFragmentation* theFrag = new G4QGSMFragmentation();
387  G4ExcitedStringDecay* theStringDecay = new G4ExcitedStringDecay(theFrag);
388  theStringModel->SetFragmentationModel(theStringDecay);
389 
390  G4GeneratorPrecompoundInterface* theCascade =
391  new G4GeneratorPrecompoundInterface;
392 
393  G4TheoFSGenerator* theModel = new G4TheoFSGenerator();
394  theModel->SetTransport(theCascade);
395  theModel->SetHighEnergyGenerator(theStringModel);
396 
397  G4HadronicParameters* param = G4HadronicParameters::Instance();
398 
399  G4CascadeInterface* cascade = new G4CascadeInterface;
400  cascade->SetMaxEnergy(param->GetMaxEnergyTransitionFTF_Cascade());
401  gnuc->RegisterMe(cascade);
402  theModel->SetMinEnergy(param->GetMinEnergyTransitionFTF_Cascade());
403  theModel->SetMaxEnergy(param->GetMaxEnergy());
404  gnuc->RegisterMe(theModel);
405 
406  G4GammaGeneralProcess* sp =
407  (G4GammaGeneralProcess*)emManager->GetGammaGeneralProcess();
408  if(sp) {
409  sp->AddHadProcess(gnuc);
410  } else {
411  // LEND may be activated if the general process is not activated
412  ph->RegisterProcess(gnuc, G4Gamma::Gamma());
413  if(gLENDActivated) { ConstructLENDGammaNuclear(cascade, gnuc); }
414  }
415 
416  if(eActivated) {
417  G4ElectronNuclearProcess* enuc = new G4ElectronNuclearProcess;
418  G4PositronNuclearProcess* pnuc = new G4PositronNuclearProcess;
419  G4ElectroVDNuclearModel* eModel = new G4ElectroVDNuclearModel;
420 
421  enuc->RegisterMe(eModel);
422  ph->RegisterProcess(enuc, G4Electron::Electron());
423 
424  pnuc->RegisterMe(eModel);
425  ph->RegisterProcess(enuc, G4Positron::Positron());
426  }
427 }
428 
429 void G4EmExtraPhysics::ConstructLENDGammaNuclear(
430  G4CascadeInterface* cascade, G4PhotoNuclearProcess* gnuc)
431 {
432  if (std::getenv("G4LENDDATA") == nullptr ) {
433  G4String message = "\n Skipping activation of Low Energy Nuclear Data (LEND) model for gamma nuclear interactions.\n The LEND model needs data files and they are available from ftp://gdo-nuclear.ucllnl.org/GND_after2013/GND_v1.3.tar.gz.\n Please set the environment variable G4LENDDATA to point to the directory named v1.3 extracted from the archive file.\n";
434  G4Exception( "G4EmExtraPhysics::ConstructLENDGammaNuclear()"
435  , "G4LENDBertiniGammaElectroNuclearBuilder001"
436  , JustWarning , message);
437  return;
438  }
439 
440  cascade->SetMinEnergy(19.9*MeV);
441  G4LENDorBERTModel* theGammaReactionLowE =
442  new G4LENDorBERTModel( G4Gamma::Gamma() );
443  theGammaReactionLowE->DumpLENDTargetInfo(true);
444  G4LENDCombinedCrossSection* theGammaCrossSectionLowE =
445  new G4LENDCombinedCrossSection( G4Gamma::Gamma() );
446  theGammaReactionLowE->SetMaxEnergy(20*MeV);
447  gnuc->RegisterMe(theGammaReactionLowE);
448  gnuc->AddDataSet(theGammaCrossSectionLowE);
449 }