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

Last change on this file since 1321 was 1230, checked in by garnier, 14 years ago

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