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

Last change on this file was 807, checked in by garnier, 17 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.