source: trunk/examples/advanced/cosmicray_charging/src/LISAPhysicsList.cc @ 1313

Last change on this file since 1313 was 807, checked in by garnier, 16 years ago

update
File size: 31.1 KB
Line 
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// * cosmicray_charging advanced example for Geant4                   *
29// * (adapted simulation of test-mass charging in the LISA mission)   *
30// *                                                                  *
31// * Henrique Araujo (h.araujo@imperial.ac.uk) & Peter Wass           *
32// * Imperial College London                                          *
33// *                                                                  *
34// * LISAPhysicsList class                                            *
35// *                                                                  *
36// ********************************************************************
37//
38// HISTORY
39// 22/02/2004: migrated from LISA-V04
40// 09/08/2004: Removed call by pointer of hadronics classes
41// 09/08/2004: Added MuNuclear interaction
42// 08/12/2005: changed particle construction
43//
44// ********************************************************************
45
46// Hadronics :
47//           : PreCo for p, n at      0 < E < 70 MeV
48//           : BiC   for p, n at 65 MeV < E < 6.1 GeV
49//           : QGSP  for p, n at  6 GeV < E < 100 TeV
50//           : BiC   for pi at        0 < E < 1.5 GeV
51//           : LEP   for pi at  1.4 GeV < E < 6.1 GeV
52//           : QGSP  for pi at    6 GeV < E < 100 TeV
53//           : LEP for H2, H3, He4 at      0 < E < 100 MeV
54//           : BiC for H2, H3, He4 at 80 MeV < E < 40 GeV
55//           : BiC for He3, GenIon at      0 < E < 30 GeV
56
57
58
59#include "LISAPhysicsList.hh"
60
61#include "G4ProcessManager.hh"
62#include "G4ProcessVector.hh"
63
64#include "G4ParticleDefinition.hh"
65#include "G4ParticleWithCuts.hh"
66#include "G4ParticleTypes.hh"
67#include "G4ParticleTable.hh"
68
69#include "globals.hh"
70#include "G4ios.hh"
71#include <iomanip>
72
73
74// Constructor /////////////////////////////////////////////////////////////
75LISAPhysicsList::LISAPhysicsList() : G4VUserPhysicsList() {
76
77  VerboseLevel = 1;
78  SetVerboseLevel(VerboseLevel);
79}
80
81
82// Destructor //////////////////////////////////////////////////////////////
83LISAPhysicsList::~LISAPhysicsList() 
84{;}
85
86
87
88////////////////////////////////////////////////////////////////////////////
89// Construct Particles /////////////////////////////////////////////////////
90
91#include "G4LeptonConstructor.hh"
92#include "G4BaryonConstructor.hh"
93#include "G4MesonConstructor.hh"
94#include "G4BosonConstructor.hh"
95#include "G4ShortLivedConstructor.hh"
96#include "G4IonConstructor.hh"
97
98void LISAPhysicsList::ConstructParticle() {
99
100  G4LeptonConstructor lepton;
101  lepton.ConstructParticle();
102 
103  G4BosonConstructor boson;
104  boson.ConstructParticle();
105 
106  G4MesonConstructor meson;
107  meson.ConstructParticle();
108 
109  G4BaryonConstructor baryon;
110  baryon.ConstructParticle();
111 
112  G4ShortLivedConstructor shortLived;
113  shortLived.ConstructParticle();
114 
115  G4IonConstructor ion;
116  ion.ConstructParticle();
117
118}
119
120/////////////////////////////////////////////////////////////////////////////
121// Construct Processes //////////////////////////////////////////////////////
122
123void LISAPhysicsList::ConstructProcess() {
124
125  AddTransportation();
126
127  ElectromagneticPhysics();
128
129  HadronicPhysics();
130
131  ElectroNuclearPhysics();
132
133  GeneralPhysics();
134
135}
136
137
138/////////////////////////////////////////////////////////////////////////////
139// Transportation ///////////////////////////////////////////////////////////
140
141void LISAPhysicsList::AddTransportation() {
142
143  G4VUserPhysicsList::AddTransportation();
144 
145}
146
147
148/////////////////////////////////////////////////////////////////////////////
149// Electromagnetic Processes ////////////////////////////////////////////////
150
151// all charged particles
152#include "G4MultipleScattering.hh"
153
154// gamma
155#include "G4LowEnergyRayleigh.hh"
156#include "G4LowEnergyPhotoElectric.hh"
157#include "G4LowEnergyCompton.hh" 
158#include "G4LowEnergyGammaConversion.hh"
159
160// e-
161#include "G4LowEnergyIonisation.hh"
162#include "G4LowEnergyBremsstrahlung.hh"
163
164// e+
165#include "G4eIonisation.hh"
166#include "G4eBremsstrahlung.hh"
167#include "G4eplusAnnihilation.hh"
168
169// muons
170#include "G4MuIonisation.hh"
171#include "G4MuBremsstrahlung.hh"
172#include "G4MuPairProduction.hh"
173#include "G4MuonMinusCaptureAtRest.hh"
174
175// hadrons and ions
176#include "G4hIonisation.hh"
177#include "G4ionIonisation.hh"
178
179
180void LISAPhysicsList::ElectromagneticPhysics() {
181
182
183  G4cout << "Electromagnetic Physics" << G4endl;
184
185
186   // processes
187
188  G4LowEnergyPhotoElectric*  lowePhot = new G4LowEnergyPhotoElectric();
189  G4LowEnergyIonisation*     loweIon  = new G4LowEnergyIonisation();
190  G4LowEnergyBremsstrahlung* loweBrem = new G4LowEnergyBremsstrahlung();
191  lowePhot->SetCutForLowEnSecPhotons(100*eV);
192  loweIon ->SetCutForLowEnSecPhotons(100*eV);
193  loweBrem->SetCutForLowEnSecPhotons(100*eV);
194
195  theParticleIterator->reset();
196  while( (*theParticleIterator)() ){
197    G4ParticleDefinition* particle = theParticleIterator->value();
198    G4ProcessManager* pmanager = particle->GetProcessManager();
199    G4String particleName      = particle->GetParticleName();
200    G4String particleType      = particle->GetParticleType();
201    G4double particleCharge    = particle->GetPDGCharge();
202   
203    if (particleName == "gamma") {
204      //gamma
205      pmanager->AddDiscreteProcess(new G4LowEnergyRayleigh());
206      pmanager->AddDiscreteProcess(lowePhot);
207      pmanager->AddDiscreteProcess(new G4LowEnergyCompton());
208      pmanager->AddDiscreteProcess(new G4LowEnergyGammaConversion());
209
210    } else if (particleName == "e-") {
211      //electron
212      G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
213      // Modifying Facrange from default value (0.199)
214      // to improve backscattering fraction for electrons
215      //      aMultipleScattering->SetRangeFactor(0.01);
216      pmanager->AddProcess(aMultipleScattering,      -1, 1, 1);
217      pmanager->AddProcess(loweIon,                  -1, 2, 2);
218      pmanager->AddProcess(loweBrem,                 -1,-1, 3);
219
220    } else if (particleName == "e+") {
221      //positron
222      G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
223      pmanager->AddProcess(aMultipleScattering,      -1, 1, 1);
224      pmanager->AddProcess(new G4eIonisation(),      -1, 2, 2);
225      pmanager->AddProcess(new G4eBremsstrahlung(),  -1,-1, 3);
226      pmanager->AddProcess(new G4eplusAnnihilation(), 0,-1, 4);     
227
228    } else if( particleName == "mu+" || 
229               particleName == "mu-"    ) {
230      //muon
231      G4MultipleScattering* aMultipleScattering = new G4MultipleScattering();
232      pmanager->AddProcess(aMultipleScattering,      -1, 1, 1);
233      pmanager->AddProcess(new G4MuIonisation(),     -1, 2, 2);
234      pmanager->AddProcess(new G4MuBremsstrahlung(), -1,-1, 3);
235      pmanager->AddProcess(new G4MuPairProduction(), -1,-1, 4);
236      if( particleName == "mu-" )
237        pmanager->AddProcess(new G4MuonMinusCaptureAtRest(),0,-1,-1);
238     
239    } else if( particleName == "GenericIon" || 
240               particleName == "alpha" ||
241               particleName == "He3") { 
242      // ions
243      pmanager->AddProcess(new G4MultipleScattering, -1, 1, 1);
244      pmanager->AddProcess(new G4ionIonisation,      -1, 2, 2);
245     
246    } else if (!particle->IsShortLived() &&
247               particleCharge != 0.0 && 
248               particleName   != "chargedgeantino") {
249      // all other stable charged particles
250      pmanager->AddProcess(new G4MultipleScattering, -1, 1, 1);
251      pmanager->AddProcess(new G4hIonisation,      -1, 2, 2);
252    }
253   
254  }
255
256}
257
258
259
260
261///////////////////////////////////////////////////////////////////////////
262// ElectroNuclear Physics /////////////////////////////////////////////////
263void LISAPhysicsList::ElectroNuclearPhysics() {
264
265  G4cout << "ElectroNuclear Physics" << G4endl;
266
267  // gamma
268  G4ProcessManager* pmanager = G4Gamma::Gamma()->GetProcessManager();
269  // low energy
270  theGammaReaction = new G4GammaNuclearReaction;
271  theGammaReaction->SetMaxEnergy(3.5*GeV);
272  thePhotoNuclearProcess.RegisterMe(theGammaReaction);
273  // high energy
274  theHEModel_PN = new G4TheoFSGenerator;
275  theCascade_PN = new G4StringChipsParticleLevelInterface;
276  theHEModel_PN->SetTransport(theCascade_PN);
277  theHEModel_PN->SetHighEnergyGenerator(theStringModel_PN);
278  theStringDecay_PN = new G4ExcitedStringDecay(&theFragmentation_PN);
279
280  theStringModel_PN = new G4QGSModel<G4GammaParticipants>;
281  theStringModel_PN -> SetFragmentationModel(theStringDecay_PN);
282  theHEModel_PN->SetMinEnergy(3.*GeV);
283  theHEModel_PN->SetMaxEnergy(100*TeV);
284  thePhotoNuclearProcess.RegisterMe(theHEModel_PN);
285  pmanager->AddDiscreteProcess(&thePhotoNuclearProcess);
286 
287  // e-
288  pmanager = G4Electron::Electron()->GetProcessManager();
289  // see G4ElectroNuclearReaction.hh for defaults
290  // 0 < E < 10 TeV
291  theElectroReaction = new G4ElectroNuclearReaction;
292  theElectroReaction->SetMaxEnergy(10*TeV);
293  theElectronNuclearProcess.RegisterMe(theElectroReaction);
294  pmanager->AddDiscreteProcess(&theElectronNuclearProcess);
295
296  // e+
297  pmanager = G4Positron::Positron()->GetProcessManager();
298  // see G4ElectroNuclearReaction.hh for defaults
299  // 0 < E < 10 TeV
300  thePositronNuclearProcess.RegisterMe(theElectroReaction);
301  pmanager->AddDiscreteProcess(&thePositronNuclearProcess);
302
303  // Mu-nuclear reaction
304  // mu-
305  pmanager = G4MuonMinus::MuonMinus()->GetProcessManager();
306  pmanager->AddDiscreteProcess(&theMuMinusNuclearInteraction);
307  // mu+
308  pmanager = G4MuonPlus::MuonPlus()->GetProcessManager();
309  pmanager->AddDiscreteProcess(&theMuPlusNuclearInteraction);
310
311
312}
313
314
315
316/////////////////////////////////////////////////////////////////////////////
317// Hadronic Physics /////////////////////////////////////////////////////////
318void LISAPhysicsList::HadronicPhysics() {
319
320
321  G4cout << "Hadronic Physics" << G4endl;
322
323
324  // **************************************************//
325  // *** preparing inelastic reactions for hadrons *** //
326  // **************************************************//
327  //
328  // high energy model for proton, neutron, pions and kaons
329  theHEModel = new G4TheoFSGenerator;
330  // all models for treatment of thermal nucleus
331  theEvaporation = new G4Evaporation;
332  theFermiBreakUp = new G4FermiBreakUp;
333  theMF = new G4StatMF;
334  // evaporation logic
335  theHandler = new G4ExcitationHandler;
336  theHandler->SetEvaporation(theEvaporation);
337  theHandler->SetFermiModel(theFermiBreakUp);
338  theHandler->SetMultiFragmentation(theMF);
339  theHandler->SetMaxAandZForFermiBreakUp(12, 6);
340  theHandler->SetMinEForMultiFrag(3.*MeV);
341
342  // pre-equilibrium stage
343  thePreEquilib = new G4PreCompoundModel(theHandler);
344  thePreEquilib->SetMaxEnergy(70*MeV);
345
346  // a no-cascade generator-precompound interaface
347  theCascade = new G4GeneratorPrecompoundInterface;
348  theCascade->SetDeExcitation(thePreEquilib);
349
350  // QGSP model
351  theStringModel = new G4QGSModel<G4QGSParticipants>;
352  theHEModel->SetTransport(theCascade);
353  theHEModel->SetHighEnergyGenerator(theStringModel);
354  theHEModel->SetMinEnergy(6*GeV);
355  theHEModel->SetMaxEnergy(100*TeV);
356  // Binary cascade for p, n
357  theCasc = new G4BinaryCascade;
358  theCasc->SetMinEnergy(65*MeV);
359  theCasc->SetMaxEnergy(6.1*GeV);
360  // fragmentation
361  theFragmentation = new G4QGSMFragmentation;
362  theStringDecay = new G4ExcitedStringDecay(theFragmentation);
363  theStringModel->SetFragmentationModel(theStringDecay);
364  //
365  // Binary Cascade for Pi
366  theCascForPi = new G4BinaryCascade;
367  theCascForPi->SetMinEnergy(0*MeV);
368  theCascForPi->SetMaxEnergy(1.5*GeV);
369  // LEP to fill the gap
370  theLEPionPlusInelasticModel = new G4LEPionPlusInelastic();
371  theLEPionPlusInelasticModel->SetMinEnergy(1.4*GeV);
372  theLEPionPlusInelasticModel->SetMaxEnergy(6.1*GeV);
373  theLEPionMinusInelasticModel = new G4LEPionMinusInelastic();
374  theLEPionMinusInelasticModel->SetMinEnergy(1.4*GeV);
375  theLEPionMinusInelasticModel->SetMaxEnergy(6.1*GeV);
376
377
378  // *******************************************************//
379  // *** preparing inelastic reactions for light nuclei *** //
380  // *******************************************************//
381  //
382  // binary cascade for light nuclei
383  // NOTE: Shen XS only up to 10 GeV/n;
384  theIonCascade= new G4BinaryLightIonReaction;
385  theIonCascade->SetMinEnergy(80*MeV);
386  theIonCascade->SetMaxEnergy(40*GeV);
387  theTripathiCrossSection = new G4TripathiCrossSection;
388  theShenCrossSection = new G4IonsShenCrossSection;
389  //
390  // deuteron
391  theLEDeuteronInelasticModel = new G4LEDeuteronInelastic();
392  theLEDeuteronInelasticModel->SetMaxEnergy(100*MeV);
393  //
394  // triton
395  theLETritonInelasticModel = new G4LETritonInelastic();
396  theLETritonInelasticModel->SetMaxEnergy(100*MeV);
397  //
398  // alpha
399  theLEAlphaInelasticModel = new G4LEAlphaInelastic();
400  theLEAlphaInelasticModel->SetMaxEnergy(100*MeV);
401  //
402  // Generic Ion and He3
403  // NOTE: Shen XS only up to 10 GeV/n;
404  theGenIonCascade = new G4BinaryLightIonReaction;
405  theGenIonCascade->SetMinEnergy(0*MeV);
406  theGenIonCascade->SetMaxEnergy(30*GeV);
407
408 
409  // ***************************//
410  // *** elastic scattering *** //
411  // ***************************//
412  //
413  theElasticModel = new G4LElastic();
414  theElasticProcess.RegisterMe(theElasticModel);
415
416
417  // *****************************************//
418  // *** attaching processes to particles *** //
419  // *****************************************//
420  //
421  theParticleIterator->reset();
422  while ((*theParticleIterator)()) {
423    G4ParticleDefinition* particle = theParticleIterator->value();
424    G4ProcessManager* pmanager = particle->GetProcessManager();
425    G4String particleName = particle->GetParticleName();
426   
427    if (particleName == "pi+") {
428      pmanager->AddDiscreteProcess(&theElasticProcess);
429      // NOTE: PreCo crahes for Pi+
430      // thePionPlusInelasticProcess.RegisterMe(thePreEquilib);
431      thePionPlusInelasticProcess.RegisterMe(theCascForPi);
432      thePionPlusInelasticProcess.RegisterMe(theLEPionPlusInelasticModel);
433      thePionPlusInelasticProcess.RegisterMe(theHEModel);
434      pmanager->AddDiscreteProcess(&thePionPlusInelasticProcess);
435
436    } else if (particleName == "pi-") {
437      pmanager->AddDiscreteProcess(&theElasticProcess);
438      // thePionMinusInelasticProcess.RegisterMe(thePreEquilib);
439      thePionMinusInelasticProcess.RegisterMe(theCascForPi);
440      thePionMinusInelasticProcess.RegisterMe(theLEPionMinusInelasticModel);
441      thePionMinusInelasticProcess.RegisterMe(theHEModel);
442      pmanager->AddDiscreteProcess(&thePionMinusInelasticProcess);
443      pmanager->AddRestProcess(&thePiMinusAbsorptionAtRest, ordDefault);
444
445    } else if (particleName == "kaon+") {
446      pmanager->AddDiscreteProcess(&theElasticProcess);
447      theLEKaonPlusInelasticModel = new G4LEKaonPlusInelastic();
448      theLEKaonPlusInelasticModel->SetMaxEnergy(25*GeV);
449      theKaonPlusInelasticProcess.RegisterMe(theLEKaonPlusInelasticModel);
450      theKaonPlusInelasticProcess.RegisterMe(theHEModel);
451      pmanager->AddDiscreteProcess(&theKaonPlusInelasticProcess);
452     
453    } else if (particleName == "kaon0S") {
454      pmanager->AddDiscreteProcess(&theElasticProcess);
455      theLEKaonZeroSInelasticModel = new G4LEKaonZeroSInelastic();
456      theLEKaonZeroSInelasticModel->SetMaxEnergy(25*GeV);
457      theKaonZeroSInelasticProcess.RegisterMe(theLEKaonZeroSInelasticModel);
458      theKaonZeroSInelasticProcess.RegisterMe(theHEModel);
459      pmanager->AddDiscreteProcess(&theKaonZeroSInelasticProcess);
460     
461    } else if (particleName == "kaon0L") {
462      pmanager->AddDiscreteProcess(&theElasticProcess);
463      theLEKaonZeroLInelasticModel = new G4LEKaonZeroLInelastic();
464      theLEKaonZeroLInelasticModel->SetMaxEnergy(25*GeV);
465      theKaonZeroLInelasticProcess.RegisterMe(theLEKaonZeroLInelasticModel);
466      theKaonZeroLInelasticProcess.RegisterMe(theHEModel);
467      pmanager->AddDiscreteProcess(&theKaonZeroLInelasticProcess);
468     
469    } else if (particleName == "kaon-") {
470      pmanager->AddDiscreteProcess(&theElasticProcess);
471      theLEKaonMinusInelasticModel = new G4LEKaonMinusInelastic();   
472      theLEKaonMinusInelasticModel->SetMaxEnergy(25*GeV);
473      theKaonMinusInelasticProcess.RegisterMe(theLEKaonMinusInelasticModel);
474      theKaonMinusInelasticProcess.RegisterMe(theHEModel);
475      pmanager->AddDiscreteProcess(&theKaonMinusInelasticProcess);
476      pmanager->AddRestProcess(&theKaonMinusAbsorptionAtRest, ordDefault);
477
478    } else if (particleName == "proton") {
479      pmanager->AddDiscreteProcess(&theElasticProcess);
480      theProtonInelasticProcess.RegisterMe(thePreEquilib);
481      theProtonInelasticProcess.RegisterMe(theCasc);
482      theProtonInelasticProcess.RegisterMe(theHEModel);
483      pmanager->AddDiscreteProcess(&theProtonInelasticProcess);
484
485    } else if (particleName == "anti_proton") {
486      pmanager->AddDiscreteProcess(&theElasticProcess);
487      theLEAntiProtonInelasticModel = new G4LEAntiProtonInelastic();
488      theHEAntiProtonInelasticModel = new G4HEAntiProtonInelastic();
489      theLEAntiProtonInelasticModel->SetMaxEnergy(25*GeV);
490      theAntiProtonInelasticProcess.RegisterMe(theLEAntiProtonInelasticModel);
491      theAntiProtonInelasticProcess.RegisterMe(theHEAntiProtonInelasticModel);
492      pmanager->AddDiscreteProcess(&theAntiProtonInelasticProcess);
493      pmanager->AddRestProcess(&theAntiProtonAnnihilationAtRest, ordDefault);
494
495    } else if (particleName == "neutron") {
496      // elastic scattering
497      // LEP
498      theNeutronElasticModel1 = new G4LElastic();
499      //   theNeutronElasticModel1->SetMinEnergy(19*MeV);
500      theNeutronElasticProcess.RegisterMe(theNeutronElasticModel1);
501      //   // HP
502      //   theNeutronElasticModel2 = new G4NeutronHPElastic();
503      //   theNeutronElasticModel2->SetMaxEnergy(19.1*MeV);
504      //   theNeutronElasticData = new G4NeutronHPElasticData;
505      //   theNeutronElasticProcess.AddDataSet(theNeutronElasticData);
506      //   theNeutronElasticProcess.RegisterMe(theNeutronElasticModel2);
507      pmanager->AddDiscreteProcess(&theNeutronElasticProcess);
508      // inelastic scattering
509      //   // HP
510      //   theNeutronInelasticModel1 = new G4NeutronHPInelastic();
511      //   theNeutronInelasticProcess.RegisterMe(theNeutronInelasticModel1);
512      //   theNeutronInelasticData1 = new G4NeutronHPInelasticData;
513      //   theNeutronInelasticProcess.AddDataSet(theNeutronInelasticData1);
514      // Preco_n + BiC + QGSP
515      theNeutronInelasticProcess.RegisterMe(thePreEquilib);
516      theNeutronInelasticProcess.RegisterMe(theCasc);
517      theNeutronInelasticProcess.RegisterMe(theHEModel);
518      pmanager->AddDiscreteProcess(&theNeutronInelasticProcess);
519      // capture
520      theNeutronCaptureModel1 = new G4LCapture();
521      //   theNeutronCaptureModel1->SetMinEnergy(19*MeV);
522      theNeutronCaptureProcess.RegisterMe(theNeutronCaptureModel1);
523      //   theNeutronCaptureModel2 = new G4NeutronHPCapture;
524      //   theNeutronCaptureProcess.RegisterMe(theNeutronCaptureModel2);
525      //   theNeutronCaptureData = new G4NeutronHPCaptureData;
526      //   theNeutronCaptureProcess.AddDataSet(theNeutronCaptureData);
527      pmanager->AddDiscreteProcess(&theNeutronCaptureProcess);
528      // fission
529      theNeutronFissionModel = new G4LFission();
530      theNeutronFissionProcess.RegisterMe(theNeutronFissionModel);
531      pmanager->AddDiscreteProcess(&theNeutronFissionProcess);
532
533    } else if (particleName == "anti_neutron") {
534      pmanager->AddDiscreteProcess(&theElasticProcess);
535      theLEAntiNeutronInelasticModel = new G4LEAntiNeutronInelastic();
536      theHEAntiNeutronInelasticModel = new G4HEAntiNeutronInelastic();
537      theLEAntiNeutronInelasticModel->SetMaxEnergy(25*GeV);
538      theAntiNeutronInelasticProcess.RegisterMe
539        (theLEAntiNeutronInelasticModel);
540      theAntiNeutronInelasticProcess.RegisterMe
541        (theHEAntiNeutronInelasticModel);
542      pmanager->AddDiscreteProcess(&theAntiNeutronInelasticProcess);
543      pmanager->AddRestProcess(&theAntiNeutronAnnihilationAtRest,ordDefault);
544
545  } else if (particleName == "deuteron") {
546      pmanager->AddDiscreteProcess(&theElasticProcess);
547      theDeuteronInelasticProcess = new G4DeuteronInelasticProcess;
548      theDeuteronInelasticProcess->AddDataSet(theTripathiCrossSection);
549      theDeuteronInelasticProcess->AddDataSet(theShenCrossSection);
550      theDeuteronInelasticProcess->RegisterMe(theLEDeuteronInelasticModel);
551      theDeuteronInelasticProcess->RegisterMe(theIonCascade);
552      pmanager->AddDiscreteProcess(theDeuteronInelasticProcess);
553
554    } else if (particleName == "triton") {
555      pmanager->AddDiscreteProcess(&theElasticProcess);
556      theTritonInelasticProcess = new G4TritonInelasticProcess;
557      theTritonInelasticProcess->AddDataSet(theTripathiCrossSection);
558      theTritonInelasticProcess->AddDataSet(theShenCrossSection);
559      theTritonInelasticProcess->RegisterMe(theLETritonInelasticModel);
560      theTritonInelasticProcess->RegisterMe(theIonCascade);
561      pmanager->AddDiscreteProcess(theTritonInelasticProcess);
562
563    } else if (particleName == "alpha") {
564      pmanager->AddDiscreteProcess(&theElasticProcess);
565      theAlphaInelasticProcess = new G4AlphaInelasticProcess;
566      theAlphaInelasticProcess->AddDataSet(theTripathiCrossSection);
567      theAlphaInelasticProcess->AddDataSet(theShenCrossSection);
568      theAlphaInelasticProcess->RegisterMe(theLEAlphaInelasticModel);
569      theAlphaInelasticProcess->RegisterMe(theIonCascade);
570      pmanager->AddDiscreteProcess(theAlphaInelasticProcess);
571     
572    } else if (particleName == "He3") {
573      // NOTE elastic scattering does not stick to He3!
574      pmanager->AddDiscreteProcess(&theElasticProcess);
575      theHe3InelasticProcess = new G4HadronInelasticProcess
576        ("He3Inelastic", G4He3::He3());
577      theHe3InelasticProcess->AddDataSet(theTripathiCrossSection);
578      theHe3InelasticProcess->AddDataSet(theShenCrossSection);
579      theHe3InelasticProcess->RegisterMe(theGenIonCascade);
580      pmanager->AddDiscreteProcess(theHe3InelasticProcess);
581
582    } else if (particleName == "GenericIon") {
583      pmanager->AddDiscreteProcess(&theElasticProcess);
584      theGenericIonInelasticProcess = new G4HadronInelasticProcess
585        ("IonInelastic", G4GenericIon::GenericIon());
586      theGenericIonInelasticProcess->AddDataSet(theTripathiCrossSection);
587      theGenericIonInelasticProcess->AddDataSet(theShenCrossSection);
588      theGenericIonInelasticProcess->RegisterMe(theGenIonCascade);
589      pmanager->AddDiscreteProcess(theGenericIonInelasticProcess);
590     
591    } else if (particleName == "lambda") {
592      pmanager->AddDiscreteProcess(&theElasticProcess);
593      theLELambdaInelasticModel = new G4LELambdaInelastic();
594      theHELambdaInelasticModel = new G4HELambdaInelastic();
595      theLELambdaInelasticModel->SetMaxEnergy(25*GeV);
596      theLambdaInelasticProcess.RegisterMe(theLELambdaInelasticModel);
597      theLambdaInelasticProcess.RegisterMe(theHELambdaInelasticModel);
598      pmanager->AddDiscreteProcess(&theLambdaInelasticProcess);
599
600    } else if (particleName == "anti_lambda") {
601      pmanager->AddDiscreteProcess(&theElasticProcess);
602      theLEAntiLambdaInelasticModel = new G4LEAntiLambdaInelastic();
603      theHEAntiLambdaInelasticModel = new G4HEAntiLambdaInelastic();
604      theLEAntiLambdaInelasticModel->SetMaxEnergy(25*GeV);
605      theAntiLambdaInelasticProcess.RegisterMe(theLEAntiLambdaInelasticModel);
606      theAntiLambdaInelasticProcess.RegisterMe(theHEAntiLambdaInelasticModel);
607      pmanager->AddDiscreteProcess(&theAntiLambdaInelasticProcess);
608
609    } else if (particleName == "omega-") {
610      pmanager->AddDiscreteProcess(&theElasticProcess);
611      theLEOmegaMinusInelasticModel = new G4LEOmegaMinusInelastic();
612      theHEOmegaMinusInelasticModel = new G4HEOmegaMinusInelastic();
613      theLEOmegaMinusInelasticModel->SetMaxEnergy(25*GeV);
614      theOmegaMinusInelasticProcess.RegisterMe(theLEOmegaMinusInelasticModel);
615      theOmegaMinusInelasticProcess.RegisterMe(theHEOmegaMinusInelasticModel);
616      pmanager->AddDiscreteProcess(&theOmegaMinusInelasticProcess);
617
618    } else if (particleName == "anti_omega-") {
619      pmanager->AddDiscreteProcess(&theElasticProcess);
620      theLEAntiOmegaMinusInelasticModel = new G4LEAntiOmegaMinusInelastic(); 
621      theHEAntiOmegaMinusInelasticModel = new G4HEAntiOmegaMinusInelastic(); 
622      theLEAntiOmegaMinusInelasticModel->SetMaxEnergy(25*GeV);
623      theAntiOmegaMinusInelasticProcess.RegisterMe
624        (theLEAntiOmegaMinusInelasticModel);
625      theAntiOmegaMinusInelasticProcess.RegisterMe
626        (theHEAntiOmegaMinusInelasticModel);
627      pmanager->AddDiscreteProcess(&theAntiOmegaMinusInelasticProcess);
628
629    } else if (particleName == "sigma-") {
630      pmanager->AddDiscreteProcess(&theElasticProcess);
631      theLESigmaMinusInelasticModel = new G4LESigmaMinusInelastic();
632      theHESigmaMinusInelasticModel = new G4HESigmaMinusInelastic();
633      theLESigmaMinusInelasticModel->SetMaxEnergy(25*GeV);
634      theSigmaMinusInelasticProcess.RegisterMe(theLESigmaMinusInelasticModel);
635      theSigmaMinusInelasticProcess.RegisterMe(theHESigmaMinusInelasticModel);
636      pmanager->AddDiscreteProcess(&theSigmaMinusInelasticProcess);
637
638    } else if (particleName == "anti_sigma-") {
639      pmanager->AddDiscreteProcess(&theElasticProcess);
640      theLEAntiSigmaMinusInelasticModel = new G4LEAntiSigmaMinusInelastic();
641      theHEAntiSigmaMinusInelasticModel = new G4HEAntiSigmaMinusInelastic();
642      theLEAntiSigmaMinusInelasticModel->SetMaxEnergy(25*GeV);
643      theAntiSigmaMinusInelasticProcess.RegisterMe
644        (theLEAntiSigmaMinusInelasticModel);
645      theAntiSigmaMinusInelasticProcess.RegisterMe
646        (theHEAntiSigmaMinusInelasticModel);
647      pmanager->AddDiscreteProcess(&theAntiSigmaMinusInelasticProcess);
648
649    } else if (particleName == "sigma+") {
650      pmanager->AddDiscreteProcess(&theElasticProcess);
651      theLESigmaPlusInelasticModel = new G4LESigmaPlusInelastic();
652      theHESigmaPlusInelasticModel = new G4HESigmaPlusInelastic();     
653      theLESigmaPlusInelasticModel->SetMaxEnergy(25*GeV);
654      theSigmaPlusInelasticProcess.RegisterMe(theLESigmaPlusInelasticModel);
655      theSigmaPlusInelasticProcess.RegisterMe(theHESigmaPlusInelasticModel);
656      pmanager->AddDiscreteProcess(&theSigmaPlusInelasticProcess);
657
658    } else if (particleName == "anti_sigma+") {
659      pmanager->AddDiscreteProcess(&theElasticProcess);
660      theLEAntiSigmaPlusInelasticModel = new G4LEAntiSigmaPlusInelastic();
661      theHEAntiSigmaPlusInelasticModel = new G4HEAntiSigmaPlusInelastic();
662      theLEAntiSigmaPlusInelasticModel->SetMaxEnergy(25*GeV);
663      theAntiSigmaPlusInelasticProcess.RegisterMe
664        (theLEAntiSigmaPlusInelasticModel);
665      theAntiSigmaPlusInelasticProcess.RegisterMe
666        (theHEAntiSigmaPlusInelasticModel);
667      pmanager->AddDiscreteProcess(&theAntiSigmaPlusInelasticProcess);
668
669    } else if (particleName == "xi0") {
670      pmanager->AddDiscreteProcess(&theElasticProcess);
671      theLEXiZeroInelasticModel = new G4LEXiZeroInelastic();
672      theHEXiZeroInelasticModel = new G4HEXiZeroInelastic();
673      theLEXiZeroInelasticModel->SetMaxEnergy(25*GeV);
674      theXiZeroInelasticProcess.RegisterMe(theLEXiZeroInelasticModel);
675      theXiZeroInelasticProcess.RegisterMe(theHEXiZeroInelasticModel);
676      pmanager->AddDiscreteProcess(&theXiZeroInelasticProcess);
677
678    } else if (particleName == "anti_xi0") {
679      pmanager->AddDiscreteProcess(&theElasticProcess);
680      theLEAntiXiZeroInelasticModel = new G4LEAntiXiZeroInelastic();
681      theHEAntiXiZeroInelasticModel = new G4HEAntiXiZeroInelastic();
682      theLEAntiXiZeroInelasticModel->SetMaxEnergy(25*GeV);
683      theAntiXiZeroInelasticProcess.RegisterMe(theLEAntiXiZeroInelasticModel);
684      theAntiXiZeroInelasticProcess.RegisterMe(theHEAntiXiZeroInelasticModel);
685      pmanager->AddDiscreteProcess(&theAntiXiZeroInelasticProcess);
686
687    } else if (particleName == "xi-") {
688      pmanager->AddDiscreteProcess(&theElasticProcess);
689      theLEXiMinusInelasticModel = new G4LEXiMinusInelastic();
690      theHEXiMinusInelasticModel = new G4HEXiMinusInelastic();
691      theLEXiMinusInelasticModel->SetMaxEnergy(25*GeV);
692      theXiMinusInelasticProcess.RegisterMe(theLEXiMinusInelasticModel);
693      theXiMinusInelasticProcess.RegisterMe(theHEXiMinusInelasticModel);
694      pmanager->AddDiscreteProcess(&theXiMinusInelasticProcess);
695
696    } else if (particleName == "anti_xi-") {
697      pmanager->AddDiscreteProcess(&theElasticProcess);
698      theLEAntiXiMinusInelasticModel = new G4LEAntiXiMinusInelastic();
699      theHEAntiXiMinusInelasticModel = new G4HEAntiXiMinusInelastic();
700      theLEAntiXiMinusInelasticModel->SetMaxEnergy(25*GeV);
701      theAntiXiMinusInelasticProcess.RegisterMe
702        (theLEAntiXiMinusInelasticModel);
703      theAntiXiMinusInelasticProcess.RegisterMe
704        (theHEAntiXiMinusInelasticModel);
705      pmanager->AddDiscreteProcess(&theAntiXiMinusInelasticProcess);
706    }
707  }
708}
709
710
711// Decays ///////////////////////////////////////////////////////////////////
712#include "G4Decay.hh"
713// #include "G4RadioactiveDecay.hh"
714#include "G4IonTable.hh"
715#include "G4Ions.hh"
716
717void LISAPhysicsList::GeneralPhysics() {
718
719  // Add Decay Process
720  G4Decay* theDecayProcess = new G4Decay("decay");
721  theParticleIterator->reset();
722  while( (*theParticleIterator)() ){
723    G4ParticleDefinition* particle = theParticleIterator->value();
724    G4ProcessManager* pmanager = particle->GetProcessManager();
725
726    if (theDecayProcess->IsApplicable(*particle)) { 
727      pmanager ->AddProcess(theDecayProcess);
728      pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
729      pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
730    }
731  }
732
733}
734
735
736
737// Cuts /////////////////////////////////////////////////////////////////////
738#include "G4Region.hh"
739#include "G4RegionStore.hh"
740#include "G4ProductionCuts.hh"
741
742void LISAPhysicsList::SetCuts() {
743
744  // low energy limit
745  G4double lowlimit=250*eV;
746  G4ProductionCutsTable::GetProductionCutsTable()
747    ->SetEnergyRange(lowlimit, 100.*GeV);
748
749  // default cuts for world volume
750  G4double cutValue = 2.0*mm;
751  SetCutValue(cutValue,"gamma");
752  SetCutValue(cutValue,"e-");
753  SetCutValue(cutValue,"e+");
754
755  // cuts per region: inertial sensor (250 eV for gamma, e-, e+)
756  G4Region* region = G4RegionStore::GetInstance()->GetRegion("sensor");
757  G4ProductionCuts* cuts = new G4ProductionCuts();
758  cuts->SetProductionCut(0.050*micrometer);
759  region->SetProductionCuts(cuts);
760 
761  if (verboseLevel>0) DumpCutValuesTable();
762
763}
764
Note: See TracBrowser for help on using the repository browser.