source: trunk/examples/advanced/underground_physics/src/DMXPhysicsList.cc@ 1344

Last change on this file since 1344 was 1337, checked in by garnier, 15 years ago

tag geant4.9.4 beta 1 + modifs locales
File size: 29.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// GEANT 4 - Underground Dark Matter Detector Advanced Example
29//
30// For information related to this code contact: Alex Howard
31// e-mail: alexander.howard@cern.ch
32// --------------------------------------------------------------
33// Comments
34//
35// Underground Advanced
36// by A. Howard and H. Araujo
37// (27th November 2001)
38//
39// PhysicsList program
40//
41// Modified:
42//
43// 14-02-03 Fix bugs in msc and hIon instanciation + cut per region
44//
45// 05-02-05 AH - changes to G4Decay - added is not short lived protection
46// and redefined particles to allow non-static creation
47// i.e. changed construction to G4MesonConstructor, G4BaryonConstructor
48//
49// 23-10-09 LP - migrated EM physics from the LowEnergy processes (not supported) to
50// the new G4Livermore model implementation. Results unchanged.
51//
52// --------------------------------------------------------------
53
54#include "DMXPhysicsList.hh"
55
56#include "globals.hh"
57#include "G4ProcessManager.hh"
58#include "G4ProcessVector.hh"
59
60#include "G4ParticleDefinition.hh"
61#include "G4ParticleWithCuts.hh"
62#include "G4ParticleTypes.hh"
63#include "G4ParticleTable.hh"
64
65#include "G4ios.hh"
66#include <iomanip>
67
68#include "G4UserLimits.hh"
69
70
71// Constructor /////////////////////////////////////////////////////////////
72DMXPhysicsList::DMXPhysicsList() : G4VUserPhysicsList()
73{
74
75 defaultCutValue = 1.0*micrometer; //
76 cutForGamma = defaultCutValue;
77 cutForElectron = 1.0*nanometer;
78 cutForPositron = defaultCutValue;
79
80 VerboseLevel = 1;
81 OpVerbLevel = 0;
82
83 SetVerboseLevel(VerboseLevel);
84}
85
86
87// Destructor //////////////////////////////////////////////////////////////
88DMXPhysicsList::~DMXPhysicsList()
89{;}
90
91
92// Construct Particles /////////////////////////////////////////////////////
93void DMXPhysicsList::ConstructParticle()
94{
95
96 // In this method, static member functions should be called
97 // for all particles which you want to use.
98 // This ensures that objects of these particle types will be
99 // created in the program.
100
101 ConstructMyBosons();
102 ConstructMyLeptons();
103 ConstructMyHadrons();
104 ConstructMyShortLiveds();
105
106}
107
108
109// construct Bosons://///////////////////////////////////////////////////
110void DMXPhysicsList::ConstructMyBosons()
111{
112 // pseudo-particles
113 G4Geantino::GeantinoDefinition();
114 G4ChargedGeantino::ChargedGeantinoDefinition();
115
116 // gamma
117 G4Gamma::GammaDefinition();
118
119 //OpticalPhotons
120 G4OpticalPhoton::OpticalPhotonDefinition();
121
122}
123
124
125// construct Leptons://///////////////////////////////////////////////////
126void DMXPhysicsList::ConstructMyLeptons()
127{
128 // leptons
129 G4Electron::ElectronDefinition();
130 G4Positron::PositronDefinition();
131 G4MuonPlus::MuonPlusDefinition();
132 G4MuonMinus::MuonMinusDefinition();
133
134 G4NeutrinoE::NeutrinoEDefinition();
135 G4AntiNeutrinoE::AntiNeutrinoEDefinition();
136 G4NeutrinoMu::NeutrinoMuDefinition();
137 G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
138}
139
140
141#include "G4MesonConstructor.hh"
142#include "G4BaryonConstructor.hh"
143#include "G4IonConstructor.hh"
144
145// construct Hadrons://///////////////////////////////////////////////////
146void DMXPhysicsList::ConstructMyHadrons()
147{
148 // mesons
149 G4MesonConstructor mConstructor;
150 mConstructor.ConstructParticle();
151
152 // baryons
153 G4BaryonConstructor bConstructor;
154 bConstructor.ConstructParticle();
155
156 // ions
157 G4IonConstructor iConstructor;
158 iConstructor.ConstructParticle();
159
160}
161
162
163// construct Shortliveds://///////////////////////////////////////////////////
164void DMXPhysicsList::ConstructMyShortLiveds()
165{
166 // ShortLiveds
167 ;
168}
169
170
171
172
173// Construct Processes //////////////////////////////////////////////////////
174void DMXPhysicsList::ConstructProcess()
175{
176
177 AddTransportation();
178
179 ConstructEM();
180
181 ConstructOp();
182
183 ConstructHad();
184
185 ConstructGeneral();
186
187}
188
189
190// Transportation ///////////////////////////////////////////////////////////
191#include "DMXMaxTimeCuts.hh"
192#include "DMXMinEkineCuts.hh"
193#include "G4StepLimiter.hh"
194
195void DMXPhysicsList::AddTransportation() {
196
197 G4VUserPhysicsList::AddTransportation();
198
199 theParticleIterator->reset();
200 while( (*theParticleIterator)() ){
201 G4ParticleDefinition* particle = theParticleIterator->value();
202 G4ProcessManager* pmanager = particle->GetProcessManager();
203 G4String particleName = particle->GetParticleName();
204 // time cuts for ONLY neutrons:
205 if(particleName == "neutron")
206 pmanager->AddDiscreteProcess(new DMXMaxTimeCuts());
207 // Energy cuts to kill charged (embedded in method) particles:
208 pmanager->AddDiscreteProcess(new DMXMinEkineCuts());
209
210 // Step limit applied to all particles:
211 pmanager->AddProcess(new G4StepLimiter, -1,-1,1);
212
213 }
214}
215
216
217// Electromagnetic Processes ////////////////////////////////////////////////
218// all charged particles
219
220// gamma
221#include "G4PhotoElectricEffect.hh"
222#include "G4LivermorePhotoElectricModel.hh"
223
224#include "G4ComptonScattering.hh"
225#include "G4LivermoreComptonModel.hh"
226
227#include "G4GammaConversion.hh"
228#include "G4LivermoreGammaConversionModel.hh"
229
230#include "G4RayleighScattering.hh"
231#include "G4LivermoreRayleighModel.hh"
232
233
234// e-
235#include "G4eMultipleScattering.hh"
236
237#include "G4eIonisation.hh"
238#include "G4LivermoreIonisationModel.hh"
239
240#include "G4eBremsstrahlung.hh"
241#include "G4LivermoreBremsstrahlungModel.hh"
242
243
244// e+
245#include "G4eIonisation.hh"
246#include "G4eBremsstrahlung.hh"
247#include "G4eplusAnnihilation.hh"
248
249
250// alpha and GenericIon and deuterons, triton, He3:
251//hIonisation #include "G4hLowEnergyIonisation.hh" -> moved to G4hIonisation
252#include "G4EnergyLossTables.hh"
253// hLowEnergyIonisation uses Ziegler 1988 as the default
254
255
256//muon:
257#include "G4MuIonisation.hh"
258#include "G4MuBremsstrahlung.hh"
259#include "G4MuPairProduction.hh"
260#include "G4MuonMinusCaptureAtRest.hh"
261
262//OTHERS:
263#include "G4hIonisation.hh"
264#include "G4hMultipleScattering.hh"
265#include "G4hBremsstrahlung.hh"
266#include "G4ionIonisation.hh"
267#include "G4IonParametrisedLossModel.hh"
268
269//em process options to allow msc step-limitation to be switched off
270#include "G4EmProcessOptions.hh"
271
272void DMXPhysicsList::ConstructEM() {
273
274 //set a finer grid of the physic tables in order to improve precision
275 //former LowEnergy models have 200 bins up to 100 GeV
276 G4EmProcessOptions opt;
277 opt.SetMaxEnergy(100*GeV);
278 opt.SetDEDXBinning(200);
279 opt.SetLambdaBinning(200);
280
281 theParticleIterator->reset();
282 while( (*theParticleIterator)() ){
283 G4ParticleDefinition* particle = theParticleIterator->value();
284 G4ProcessManager* pmanager = particle->GetProcessManager();
285 G4String particleName = particle->GetParticleName();
286 G4String particleType = particle->GetParticleType();
287 G4double charge = particle->GetPDGCharge();
288
289 if (particleName == "gamma")
290 {
291 //gamma
292 G4RayleighScattering* theRayleigh = new G4RayleighScattering();
293 theRayleigh->SetModel(new G4LivermoreRayleighModel()); //not strictly necessary
294 pmanager->AddDiscreteProcess(theRayleigh);
295
296 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
297 thePhotoElectricEffect->SetModel(new G4LivermorePhotoElectricModel());
298 pmanager->AddDiscreteProcess(thePhotoElectricEffect);
299
300 G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
301 theComptonScattering->SetModel(new G4LivermoreComptonModel());
302 pmanager->AddDiscreteProcess(theComptonScattering);
303
304 G4GammaConversion* theGammaConversion = new G4GammaConversion();
305 theGammaConversion->SetModel(new G4LivermoreGammaConversionModel());
306 pmanager->AddDiscreteProcess(theGammaConversion);
307
308 }
309 else if (particleName == "e-")
310 {
311 //electron
312 // process ordering: AddProcess(name, at rest, along step, post step)
313 // Multiple scattering
314 G4eMultipleScattering* msc = new G4eMultipleScattering();
315 msc->SetStepLimitType(fUseDistanceToBoundary);
316 pmanager->AddProcess(msc,-1, 1, 1);
317
318 // Ionisation
319 G4eIonisation* eIonisation = new G4eIonisation();
320 eIonisation->SetEmModel(new G4LivermoreIonisationModel());
321 eIonisation->SetStepFunction(0.2, 100*um); //improved precision in tracking
322 pmanager->AddProcess(eIonisation,-1, 2, 2);
323
324 // Bremsstrahlung
325 G4eBremsstrahlung* eBremsstrahlung = new G4eBremsstrahlung();
326 eBremsstrahlung->SetEmModel(new G4LivermoreBremsstrahlungModel());
327 pmanager->AddProcess(eBremsstrahlung, -1,-3, 3);
328 }
329 else if (particleName == "e+")
330 {
331 //positron
332 G4eMultipleScattering* msc = new G4eMultipleScattering();
333 msc->SetStepLimitType(fUseDistanceToBoundary);
334 pmanager->AddProcess(msc,-1, 1, 1);
335
336 // Ionisation
337 G4eIonisation* eIonisation = new G4eIonisation();
338 eIonisation->SetStepFunction(0.2, 100*um); //
339 pmanager->AddProcess(eIonisation, -1, 2, 2);
340
341 //Bremsstrahlung (use default, no low-energy available)
342 pmanager->AddProcess(new G4eBremsstrahlung(), -1,-1, 3);
343
344 //Annihilation
345 pmanager->AddProcess(new G4eplusAnnihilation(),0,-1, 4);
346 }
347 else if( particleName == "mu+" ||
348 particleName == "mu-" )
349 {
350 //muon
351 pmanager->AddProcess(new G4eMultipleScattering, -1, 1, 1);
352 pmanager->AddProcess(new G4MuIonisation(), -1, 2, 2);
353 pmanager->AddProcess(new G4MuBremsstrahlung(), -1,-1, 3);
354 pmanager->AddProcess(new G4MuPairProduction(), -1,-1, 4);
355 if( particleName == "mu-" )
356 pmanager->AddProcess(new G4MuonMinusCaptureAtRest(), 0,-1,-1);
357 }
358 else if (particleName == "proton" ||
359 particleName == "pi+" ||
360 particleName == "pi-")
361 {
362 //multiple scattering
363 pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1);
364
365 //ionisation
366 G4hIonisation* hIonisation = new G4hIonisation();
367 hIonisation->SetStepFunction(0.2, 50*um);
368 pmanager->AddProcess(hIonisation, -1, 2, 2);
369
370 //bremmstrahlung
371 pmanager->AddProcess(new G4hBremsstrahlung, -1,-3, 3);
372 }
373 else if(particleName == "alpha" ||
374 particleName == "deuteron" ||
375 particleName == "triton" ||
376 particleName == "He3")
377 {
378 //multiple scattering
379 pmanager->AddProcess(new G4hMultipleScattering,-1,1,1);
380
381 //ionisation
382 G4ionIonisation* ionIoni = new G4ionIonisation();
383 ionIoni->SetStepFunction(0.1, 20*um);
384 pmanager->AddProcess(ionIoni, -1, 2, 2);
385 }
386 else if (particleName == "GenericIon")
387 {
388 // OBJECT may be dynamically created as either a GenericIon or nucleus
389 // G4Nucleus exists and therefore has particle type nucleus
390 // genericIon:
391
392 //multiple scattering
393 pmanager->AddProcess(new G4hMultipleScattering,-1,1,1);
394
395 //ionisation
396 G4ionIonisation* ionIoni = new G4ionIonisation();
397 ionIoni->SetEmModel(new G4IonParametrisedLossModel());
398 ionIoni->SetStepFunction(0.1, 20*um);
399 pmanager->AddProcess(ionIoni, -1, 2, 2);
400 }
401
402 else if ((!particle->IsShortLived()) &&
403 (charge != 0.0) &&
404 (particle->GetParticleName() != "chargedgeantino"))
405 {
406 //all others charged particles except geantino
407 G4hMultipleScattering* aMultipleScattering = new G4hMultipleScattering();
408 G4hIonisation* ahadronIon = new G4hIonisation();
409
410 //multiple scattering
411 pmanager->AddProcess(aMultipleScattering,-1,1,1);
412
413 //ionisation
414 pmanager->AddProcess(ahadronIon, -1,2,2);
415 }
416
417 }
418
419 // turn off msc step-limitation - especially as electron cut 1nm
420 opt.SetMscStepLimitation(fMinimal);
421
422}
423
424
425// Optical Processes ////////////////////////////////////////////////////////
426#include "G4Scintillation.hh"
427#include "G4OpAbsorption.hh"
428//#include "G4OpRayleigh.hh"
429#include "G4OpBoundaryProcess.hh"
430
431void DMXPhysicsList::ConstructOp()
432{
433 // default scintillation process
434 G4Scintillation* theScintProcessDef = new G4Scintillation("Scintillation");
435 // theScintProcessDef->DumpPhysicsTable();
436 theScintProcessDef->SetTrackSecondariesFirst(true);
437 theScintProcessDef->SetScintillationYieldFactor(1.0); //
438 theScintProcessDef->SetScintillationExcitationRatio(0.0); //
439 theScintProcessDef->SetVerboseLevel(OpVerbLevel);
440
441 // scintillation process for alpha:
442 G4Scintillation* theScintProcessAlpha = new G4Scintillation("Scintillation");
443 // theScintProcessNuc->DumpPhysicsTable();
444 theScintProcessAlpha->SetTrackSecondariesFirst(true);
445 theScintProcessAlpha->SetScintillationYieldFactor(1.1);
446 theScintProcessAlpha->SetScintillationExcitationRatio(1.0);
447 theScintProcessAlpha->SetVerboseLevel(OpVerbLevel);
448
449 // scintillation process for heavy nuclei
450 G4Scintillation* theScintProcessNuc = new G4Scintillation("Scintillation");
451 // theScintProcessNuc->DumpPhysicsTable();
452 theScintProcessNuc->SetTrackSecondariesFirst(true);
453 theScintProcessNuc->SetScintillationYieldFactor(0.2);
454 theScintProcessNuc->SetScintillationExcitationRatio(1.0);
455 theScintProcessNuc->SetVerboseLevel(OpVerbLevel);
456
457 // optical processes
458 G4OpAbsorption* theAbsorptionProcess = new G4OpAbsorption();
459 // G4OpRayleigh* theRayleighScatteringProcess = new G4OpRayleigh();
460 G4OpBoundaryProcess* theBoundaryProcess = new G4OpBoundaryProcess();
461 // theAbsorptionProcess->DumpPhysicsTable();
462 // theRayleighScatteringProcess->DumpPhysicsTable();
463 theAbsorptionProcess->SetVerboseLevel(OpVerbLevel);
464 // theRayleighScatteringProcess->SetVerboseLevel(OpVerbLevel);
465 theBoundaryProcess->SetVerboseLevel(OpVerbLevel);
466 G4OpticalSurfaceModel themodel = unified;
467 theBoundaryProcess->SetModel(themodel);
468
469 theParticleIterator->reset();
470 while( (*theParticleIterator)() )
471 {
472 G4ParticleDefinition* particle = theParticleIterator->value();
473 G4ProcessManager* pmanager = particle->GetProcessManager();
474 G4String particleName = particle->GetParticleName();
475 if (theScintProcessDef->IsApplicable(*particle)) {
476 // if(particle->GetPDGMass() > 5.0*GeV)
477 if(particle->GetParticleName() == "GenericIon") {
478 pmanager->AddProcess(theScintProcessNuc); // AtRestDiscrete
479 pmanager->SetProcessOrderingToLast(theScintProcessNuc,idxAtRest);
480 pmanager->SetProcessOrderingToLast(theScintProcessNuc,idxPostStep);
481 }
482 else if(particle->GetParticleName() == "alpha") {
483 pmanager->AddProcess(theScintProcessAlpha);
484 pmanager->SetProcessOrderingToLast(theScintProcessAlpha,idxAtRest);
485 pmanager->SetProcessOrderingToLast(theScintProcessAlpha,idxPostStep);
486 }
487 else {
488 pmanager->AddProcess(theScintProcessDef);
489 pmanager->SetProcessOrderingToLast(theScintProcessDef,idxAtRest);
490 pmanager->SetProcessOrderingToLast(theScintProcessDef,idxPostStep);
491 }
492 }
493
494 if (particleName == "opticalphoton") {
495 pmanager->AddDiscreteProcess(theAbsorptionProcess);
496 // pmanager->AddDiscreteProcess(theRayleighScatteringProcess);
497 pmanager->AddDiscreteProcess(theBoundaryProcess);
498 }
499 }
500}
501
502
503// Hadronic processes ////////////////////////////////////////////////////////
504
505// Elastic processes:
506#include "G4HadronElasticProcess.hh"
507
508// Inelastic processes:
509#include "G4PionPlusInelasticProcess.hh"
510#include "G4PionMinusInelasticProcess.hh"
511#include "G4KaonPlusInelasticProcess.hh"
512#include "G4KaonZeroSInelasticProcess.hh"
513#include "G4KaonZeroLInelasticProcess.hh"
514#include "G4KaonMinusInelasticProcess.hh"
515#include "G4ProtonInelasticProcess.hh"
516#include "G4AntiProtonInelasticProcess.hh"
517#include "G4NeutronInelasticProcess.hh"
518#include "G4AntiNeutronInelasticProcess.hh"
519#include "G4DeuteronInelasticProcess.hh"
520#include "G4TritonInelasticProcess.hh"
521#include "G4AlphaInelasticProcess.hh"
522
523// Low-energy Models: < 20GeV
524#include "G4LElastic.hh"
525#include "G4LEPionPlusInelastic.hh"
526#include "G4LEPionMinusInelastic.hh"
527#include "G4LEKaonPlusInelastic.hh"
528#include "G4LEKaonZeroSInelastic.hh"
529#include "G4LEKaonZeroLInelastic.hh"
530#include "G4LEKaonMinusInelastic.hh"
531#include "G4LEProtonInelastic.hh"
532#include "G4LEAntiProtonInelastic.hh"
533#include "G4LENeutronInelastic.hh"
534#include "G4LEAntiNeutronInelastic.hh"
535#include "G4LEDeuteronInelastic.hh"
536#include "G4LETritonInelastic.hh"
537#include "G4LEAlphaInelastic.hh"
538#include "G4HadronCaptureProcess.hh"
539// High-energy Models: >20 GeV
540#include "G4HEPionPlusInelastic.hh"
541#include "G4HEPionMinusInelastic.hh"
542#include "G4HEKaonPlusInelastic.hh"
543#include "G4HEKaonZeroInelastic.hh"
544#include "G4HEKaonZeroInelastic.hh"
545#include "G4HEKaonMinusInelastic.hh"
546#include "G4HEProtonInelastic.hh"
547#include "G4HEAntiProtonInelastic.hh"
548#include "G4HENeutronInelastic.hh"
549#include "G4HEAntiNeutronInelastic.hh"
550
551// Neutron high-precision models: <20 MeV
552#include "G4NeutronHPElastic.hh"
553#include "G4NeutronHPElasticData.hh"
554#include "G4NeutronHPCapture.hh"
555#include "G4NeutronHPCaptureData.hh"
556#include "G4NeutronHPInelastic.hh"
557#include "G4NeutronHPInelasticData.hh"
558#include "G4LCapture.hh"
559
560// Stopping processes
561#include "G4PiMinusAbsorptionAtRest.hh"
562#include "G4KaonMinusAbsorptionAtRest.hh"
563#include "G4AntiProtonAnnihilationAtRest.hh"
564#include "G4AntiNeutronAnnihilationAtRest.hh"
565
566
567// ConstructHad()
568// Makes discrete physics processes for the hadrons, at present limited
569// to those particles with GHEISHA interactions (INTRC > 0).
570// The processes are: Elastic scattering and Inelastic scattering.
571// F.W.Jones 09-JUL-1998
572void DMXPhysicsList::ConstructHad()
573{
574 G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess;
575 G4LElastic* theElasticModel = new G4LElastic;
576 theElasticProcess->RegisterMe(theElasticModel);
577
578 theParticleIterator->reset();
579 while ((*theParticleIterator)())
580 {
581 G4ParticleDefinition* particle = theParticleIterator->value();
582 G4ProcessManager* pmanager = particle->GetProcessManager();
583 G4String particleName = particle->GetParticleName();
584
585 if (particleName == "pi+")
586 {
587 pmanager->AddDiscreteProcess(theElasticProcess);
588 G4PionPlusInelasticProcess* theInelasticProcess =
589 new G4PionPlusInelasticProcess("inelastic");
590 G4LEPionPlusInelastic* theLEInelasticModel =
591 new G4LEPionPlusInelastic;
592 theInelasticProcess->RegisterMe(theLEInelasticModel);
593 G4HEPionPlusInelastic* theHEInelasticModel =
594 new G4HEPionPlusInelastic;
595 theInelasticProcess->RegisterMe(theHEInelasticModel);
596 pmanager->AddDiscreteProcess(theInelasticProcess);
597 }
598
599 else if (particleName == "pi-")
600 {
601 pmanager->AddDiscreteProcess(theElasticProcess);
602 G4PionMinusInelasticProcess* theInelasticProcess =
603 new G4PionMinusInelasticProcess("inelastic");
604 G4LEPionMinusInelastic* theLEInelasticModel =
605 new G4LEPionMinusInelastic;
606 theInelasticProcess->RegisterMe(theLEInelasticModel);
607 G4HEPionMinusInelastic* theHEInelasticModel =
608 new G4HEPionMinusInelastic;
609 theInelasticProcess->RegisterMe(theHEInelasticModel);
610 pmanager->AddDiscreteProcess(theInelasticProcess);
611 G4String prcNam;
612 pmanager->AddRestProcess(new G4PiMinusAbsorptionAtRest, ordDefault);
613 }
614
615 else if (particleName == "kaon+")
616 {
617 pmanager->AddDiscreteProcess(theElasticProcess);
618 G4KaonPlusInelasticProcess* theInelasticProcess =
619 new G4KaonPlusInelasticProcess("inelastic");
620 G4LEKaonPlusInelastic* theLEInelasticModel =
621 new G4LEKaonPlusInelastic;
622 theInelasticProcess->RegisterMe(theLEInelasticModel);
623 G4HEKaonPlusInelastic* theHEInelasticModel =
624 new G4HEKaonPlusInelastic;
625 theInelasticProcess->RegisterMe(theHEInelasticModel);
626 pmanager->AddDiscreteProcess(theInelasticProcess);
627 }
628
629 else if (particleName == "kaon0S")
630 {
631 pmanager->AddDiscreteProcess(theElasticProcess);
632 G4KaonZeroSInelasticProcess* theInelasticProcess =
633 new G4KaonZeroSInelasticProcess("inelastic");
634 G4LEKaonZeroSInelastic* theLEInelasticModel =
635 new G4LEKaonZeroSInelastic;
636 theInelasticProcess->RegisterMe(theLEInelasticModel);
637 G4HEKaonZeroInelastic* theHEInelasticModel =
638 new G4HEKaonZeroInelastic;
639 theInelasticProcess->RegisterMe(theHEInelasticModel);
640 pmanager->AddDiscreteProcess(theInelasticProcess);
641 }
642
643 else if (particleName == "kaon0L")
644 {
645 pmanager->AddDiscreteProcess(theElasticProcess);
646 G4KaonZeroLInelasticProcess* theInelasticProcess =
647 new G4KaonZeroLInelasticProcess("inelastic");
648 G4LEKaonZeroLInelastic* theLEInelasticModel =
649 new G4LEKaonZeroLInelastic;
650 theInelasticProcess->RegisterMe(theLEInelasticModel);
651 G4HEKaonZeroInelastic* theHEInelasticModel =
652 new G4HEKaonZeroInelastic;
653 theInelasticProcess->RegisterMe(theHEInelasticModel);
654 pmanager->AddDiscreteProcess(theInelasticProcess);
655 }
656
657 else if (particleName == "kaon-")
658 {
659 pmanager->AddDiscreteProcess(theElasticProcess);
660 G4KaonMinusInelasticProcess* theInelasticProcess =
661 new G4KaonMinusInelasticProcess("inelastic");
662 G4LEKaonMinusInelastic* theLEInelasticModel =
663 new G4LEKaonMinusInelastic;
664 theInelasticProcess->RegisterMe(theLEInelasticModel);
665 G4HEKaonMinusInelastic* theHEInelasticModel =
666 new G4HEKaonMinusInelastic;
667 theInelasticProcess->RegisterMe(theHEInelasticModel);
668 pmanager->AddDiscreteProcess(theInelasticProcess);
669 pmanager->AddRestProcess(new G4KaonMinusAbsorptionAtRest, ordDefault);
670 }
671
672 else if (particleName == "proton")
673 {
674 pmanager->AddDiscreteProcess(theElasticProcess);
675 G4ProtonInelasticProcess* theInelasticProcess =
676 new G4ProtonInelasticProcess("inelastic");
677 G4LEProtonInelastic* theLEInelasticModel = new G4LEProtonInelastic;
678 theInelasticProcess->RegisterMe(theLEInelasticModel);
679 G4HEProtonInelastic* theHEInelasticModel = new G4HEProtonInelastic;
680 theInelasticProcess->RegisterMe(theHEInelasticModel);
681 pmanager->AddDiscreteProcess(theInelasticProcess);
682 }
683
684 else if (particleName == "anti_proton")
685 {
686 pmanager->AddDiscreteProcess(theElasticProcess);
687 G4AntiProtonInelasticProcess* theInelasticProcess =
688 new G4AntiProtonInelasticProcess("inelastic");
689 G4LEAntiProtonInelastic* theLEInelasticModel =
690 new G4LEAntiProtonInelastic;
691 theInelasticProcess->RegisterMe(theLEInelasticModel);
692 G4HEAntiProtonInelastic* theHEInelasticModel =
693 new G4HEAntiProtonInelastic;
694 theInelasticProcess->RegisterMe(theHEInelasticModel);
695 pmanager->AddDiscreteProcess(theInelasticProcess);
696 }
697
698 else if (particleName == "neutron") {
699 // elastic scattering
700 G4HadronElasticProcess* theNeutronElasticProcess =
701 new G4HadronElasticProcess;
702 G4LElastic* theElasticModel1 = new G4LElastic;
703 G4NeutronHPElastic * theElasticNeutron = new G4NeutronHPElastic;
704 theNeutronElasticProcess->RegisterMe(theElasticModel1);
705 theElasticModel1->SetMinEnergy(19*MeV);
706 theNeutronElasticProcess->RegisterMe(theElasticNeutron);
707 G4NeutronHPElasticData * theNeutronData = new G4NeutronHPElasticData;
708 theNeutronElasticProcess->AddDataSet(theNeutronData);
709 pmanager->AddDiscreteProcess(theNeutronElasticProcess);
710 // inelastic scattering
711 G4NeutronInelasticProcess* theInelasticProcess =
712 new G4NeutronInelasticProcess("inelastic");
713 G4LENeutronInelastic* theInelasticModel = new G4LENeutronInelastic;
714 theInelasticModel->SetMinEnergy(19*MeV);
715 theInelasticProcess->RegisterMe(theInelasticModel);
716 G4NeutronHPInelastic * theLENeutronInelasticModel =
717 new G4NeutronHPInelastic;
718 theInelasticProcess->RegisterMe(theLENeutronInelasticModel);
719 G4NeutronHPInelasticData * theNeutronData1 =
720 new G4NeutronHPInelasticData;
721 theInelasticProcess->AddDataSet(theNeutronData1);
722 pmanager->AddDiscreteProcess(theInelasticProcess);
723 // capture
724 G4HadronCaptureProcess* theCaptureProcess =
725 new G4HadronCaptureProcess;
726 G4LCapture* theCaptureModel = new G4LCapture;
727 theCaptureModel->SetMinEnergy(19*MeV);
728 theCaptureProcess->RegisterMe(theCaptureModel);
729 G4NeutronHPCapture * theLENeutronCaptureModel = new G4NeutronHPCapture;
730 theCaptureProcess->RegisterMe(theLENeutronCaptureModel);
731 G4NeutronHPCaptureData * theNeutronData3 = new G4NeutronHPCaptureData;
732 theCaptureProcess->AddDataSet(theNeutronData3);
733 pmanager->AddDiscreteProcess(theCaptureProcess);
734 // G4ProcessManager* pmanager = G4Neutron::Neutron->GetProcessManager();
735 // pmanager->AddProcess(new G4UserSpecialCuts(),-1,-1,1);
736 }
737 else if (particleName == "anti_neutron")
738 {
739 pmanager->AddDiscreteProcess(theElasticProcess);
740 G4AntiNeutronInelasticProcess* theInelasticProcess =
741 new G4AntiNeutronInelasticProcess("inelastic");
742 G4LEAntiNeutronInelastic* theLEInelasticModel =
743 new G4LEAntiNeutronInelastic;
744 theInelasticProcess->RegisterMe(theLEInelasticModel);
745 G4HEAntiNeutronInelastic* theHEInelasticModel =
746 new G4HEAntiNeutronInelastic;
747 theInelasticProcess->RegisterMe(theHEInelasticModel);
748 pmanager->AddDiscreteProcess(theInelasticProcess);
749 }
750
751 else if (particleName == "deuteron")
752 {
753 pmanager->AddDiscreteProcess(theElasticProcess);
754 G4DeuteronInelasticProcess* theInelasticProcess =
755 new G4DeuteronInelasticProcess("inelastic");
756 G4LEDeuteronInelastic* theLEInelasticModel =
757 new G4LEDeuteronInelastic;
758 theInelasticProcess->RegisterMe(theLEInelasticModel);
759 pmanager->AddDiscreteProcess(theInelasticProcess);
760 }
761
762 else if (particleName == "triton")
763 {
764 pmanager->AddDiscreteProcess(theElasticProcess);
765 G4TritonInelasticProcess* theInelasticProcess =
766 new G4TritonInelasticProcess("inelastic");
767 G4LETritonInelastic* theLEInelasticModel =
768 new G4LETritonInelastic;
769 theInelasticProcess->RegisterMe(theLEInelasticModel);
770 pmanager->AddDiscreteProcess(theInelasticProcess);
771 }
772
773 else if (particleName == "alpha")
774 {
775 pmanager->AddDiscreteProcess(theElasticProcess);
776 G4AlphaInelasticProcess* theInelasticProcess =
777 new G4AlphaInelasticProcess("inelastic");
778 G4LEAlphaInelastic* theLEInelasticModel =
779 new G4LEAlphaInelastic;
780 theInelasticProcess->RegisterMe(theLEInelasticModel);
781 pmanager->AddDiscreteProcess(theInelasticProcess);
782 }
783
784 }
785}
786
787
788// Decays ///////////////////////////////////////////////////////////////////
789#include "G4Decay.hh"
790#include "G4RadioactiveDecay.hh"
791#include "G4IonTable.hh"
792#include "G4Ions.hh"
793
794void DMXPhysicsList::ConstructGeneral() {
795
796 // Add Decay Process
797 G4Decay* theDecayProcess = new G4Decay();
798 theParticleIterator->reset();
799 while( (*theParticleIterator)() )
800 {
801 G4ParticleDefinition* particle = theParticleIterator->value();
802 G4ProcessManager* pmanager = particle->GetProcessManager();
803
804 if (theDecayProcess->IsApplicable(*particle) && !particle->IsShortLived())
805 {
806 pmanager ->AddProcess(theDecayProcess);
807 // set ordering for PostStepDoIt and AtRestDoIt
808 pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
809 pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
810 }
811 }
812
813 // Declare radioactive decay to the GenericIon in the IonTable.
814 const G4IonTable *theIonTable =
815 G4ParticleTable::GetParticleTable()->GetIonTable();
816 G4RadioactiveDecay *theRadioactiveDecay = new G4RadioactiveDecay();
817
818 for (G4int i=0; i<theIonTable->Entries(); i++)
819 {
820 G4String particleName = theIonTable->GetParticle(i)->GetParticleName();
821 G4String particleType = theIonTable->GetParticle(i)->GetParticleType();
822
823 if (particleName == "GenericIon")
824 {
825 G4ProcessManager* pmanager =
826 theIonTable->GetParticle(i)->GetProcessManager();
827 pmanager->SetVerboseLevel(VerboseLevel);
828 pmanager ->AddProcess(theRadioactiveDecay);
829 pmanager ->SetProcessOrdering(theRadioactiveDecay, idxPostStep);
830 pmanager ->SetProcessOrdering(theRadioactiveDecay, idxAtRest);
831 }
832 }
833}
834
835// Cuts /////////////////////////////////////////////////////////////////////
836void DMXPhysicsList::SetCuts()
837{
838
839 if (verboseLevel >1)
840 G4cout << "DMXPhysicsList::SetCuts:";
841
842 if (verboseLevel>0){
843 G4cout << "DMXPhysicsList::SetCuts:";
844 G4cout << "CutLength : "
845 << G4BestUnit(defaultCutValue,"Length") << G4endl;
846 }
847
848 //special for low energy physics
849 G4double lowlimit=250*eV;
850 G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(lowlimit,100.*GeV);
851
852 // set cut values for gamma at first and for e- second and next for e+,
853 // because some processes for e+/e- need cut values for gamma
854 SetCutValue(cutForGamma, "gamma");
855 SetCutValue(cutForElectron, "e-");
856 SetCutValue(cutForPositron, "e+");
857
858 if (verboseLevel>0) DumpCutValuesTable();
859}
860
Note: See TracBrowser for help on using the repository browser.