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

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