src/PhysicsList.cc

Go to the documentation of this file.

//
// $Id: PhysicsList.cc 189 2008-10-23 16:32:49Z gunter $
//
// PhysicsList.cc
//
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#include "globals.hh"
#include "PhysicsList.hh"

#include "G4ProcessManager.hh"
#include "G4ParticleTypes.hh"

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

PhysicsList::PhysicsList():  G4VUserPhysicsList()
{
  defaultCutValue = 1.0*mm;
  SetVerboseLevel(1);
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

PhysicsList::~PhysicsList()
{}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#include "G4MesonConstructor.hh"
#include "G4BaryonConstructor.hh"
#include "G4ShortLivedConstructor.hh"
#include "G4BosonConstructor.hh"
#include "G4LeptonConstructor.hh"
#include "G4IonConstructor.hh"

void PhysicsList::ConstructParticle()
{
  // pseudo-particles
  G4Geantino::GeantinoDefinition();

  // gamma
  G4Gamma::GammaDefinition();
  
  // e+/-
  G4Electron::ElectronDefinition();
  G4Positron::PositronDefinition();
  
  // muon
  G4MuonPlus::MuonPlusDefinition();
  G4MuonMinus::MuonMinusDefinition();
  // neutrino
  G4NeutrinoE::NeutrinoEDefinition();
  G4AntiNeutrinoE::AntiNeutrinoEDefinition();
  G4NeutrinoMu::NeutrinoMuDefinition();
  G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
  
  G4Proton::ProtonDefinition();
  G4Neutron::NeutronDefinition();

  G4LeptonConstructor pLeptonConstructor;
  pLeptonConstructor.ConstructParticle();
  
  G4MesonConstructor pMesonConstructor;
  pMesonConstructor.ConstructParticle();

  G4BaryonConstructor pBaryonConstructor;
  pBaryonConstructor.ConstructParticle();

  G4ShortLivedConstructor pShortLivedConstructor;
  pShortLivedConstructor.ConstructParticle();  

  G4BosonConstructor  pBosonConstructor;
  pBosonConstructor.ConstructParticle();

  G4IonConstructor pConstructor;
  pConstructor.ConstructParticle();  

}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void PhysicsList::ConstructProcess()
{
  AddTransportation();
  ConstructEM();
  ConstructDecay();
  ConstructHad();

  // step limitation (as a full process)
  //  
  AddStepMax();      
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#include "G4ComptonScattering.hh"
#include "G4GammaConversion.hh"
#include "G4PhotoElectricEffect.hh"

#include "G4MultipleScattering.hh"

#include "G4eIonisation.hh"
#include "G4eBremsstrahlung.hh"
#include "G4eplusAnnihilation.hh"

#include "G4MuIonisation.hh"
#include "G4MuBremsstrahlung.hh"
#include "G4MuPairProduction.hh"

#include "G4hIonisation.hh"

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void PhysicsList::ConstructEM()
{
  theParticleIterator->reset();
  while( (*theParticleIterator)() ){
    G4ParticleDefinition* particle = theParticleIterator->value();
    G4ProcessManager* pmanager = particle->GetProcessManager();
    G4String particleName = particle->GetParticleName();
     
    if (particleName == "gamma") {
      // gamma         
      pmanager->AddDiscreteProcess(new G4PhotoElectricEffect);
      pmanager->AddDiscreteProcess(new G4ComptonScattering);
      pmanager->AddDiscreteProcess(new G4GammaConversion);
      
    } else if (particleName == "e-") {
      //electron
      pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
      pmanager->AddProcess(new G4eIonisation,        -1, 2,2);
      pmanager->AddProcess(new G4eBremsstrahlung,    -1, 3,3);      

    } else if (particleName == "e+") {
      //positron
      pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
      pmanager->AddProcess(new G4eIonisation,        -1, 2,2);
      pmanager->AddProcess(new G4eBremsstrahlung,    -1, 3,3);
      pmanager->AddProcess(new G4eplusAnnihilation,   0,-1,4);

    } else if( particleName == "mu+" || 
               particleName == "mu-"    ) {
      //muon  
      pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
      pmanager->AddProcess(new G4MuIonisation,       -1, 2,2);
      pmanager->AddProcess(new G4MuBremsstrahlung,   -1, 3,3);
      pmanager->AddProcess(new G4MuPairProduction,   -1, 4,4);
                  
    } else if (particle->GetPDGCharge() != 0.0 && 
               ! particle->IsShortLived()      &&
               particle->GetParticleName() != "chargedgeantino" ) {
      pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
      pmanager->AddProcess(new G4hIonisation,        -1, 2,2);
    }
  }
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#include "G4Decay.hh"

void PhysicsList::ConstructDecay()
{
  // Add Decay Process
  G4Decay* theDecayProcess = new G4Decay();
  theParticleIterator->reset();
  while( (*theParticleIterator)() ){
    G4ParticleDefinition* particle = theParticleIterator->value();
    G4ProcessManager* pmanager = particle->GetProcessManager();
    if (theDecayProcess->IsApplicable(*particle)) {
      pmanager ->AddProcess(theDecayProcess);
      // set ordering for PostStepDoIt and AtRestDoIt
      pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
      pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
    }
  }  
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

#include "G4HadronElasticProcess.hh"
#include "G4LElastic.hh"
#include "G4ProtonInelasticProcess.hh"
#include "G4NeutronInelasticProcess.hh"
#include "G4BinaryCascade.hh"
#include "G4HadronCaptureProcess.hh"
#include "G4LCapture.hh"
#include "G4HadronFissionProcess.hh"
#include "G4LFission.hh"
#include "G4NucleonNuclearCrossSection.hh"

// needed for higher energies 
#include "G4TheoFSGenerator.hh"
#include "G4FTFModel.hh"
#include "G4ExcitedStringDecay.hh"
#include "G4LundStringFragmentation.hh"

void PhysicsList::ConstructHad()
{
// add hadronic processes for Proton
//
  G4ParticleDefinition* proton = G4Proton::Proton();
  G4ProcessManager* proton_manager = proton->GetProcessManager();

  // process: elastic scattering
  //
  G4HadronElasticProcess* theElasticProcess = new G4HadronElasticProcess();
  
  // model:  use LElastic, i.e. parameterized model
  //
  G4LElastic* theElasticModel = new G4LElastic();
  theElasticProcess->RegisterMe(theElasticModel);
  proton_manager->AddDiscreteProcess(theElasticProcess);

  // process: inelastic interaction
  //   
  G4ProtonInelasticProcess* protonInelasticprocess 
                                    =  new G4ProtonInelasticProcess();
                                    
  // model:  use Binary cascade as model (good up to several GeV)
  //
  G4BinaryCascade* theModel = new G4BinaryCascade();  
  protonInelasticprocess->RegisterMe(theModel);
  
  // cross Section
  //
  protonInelasticprocess->AddDataSet(new G4NucleonNuclearCrossSection());
  
  // add process
  //
  proton_manager->AddDiscreteProcess(protonInelasticprocess);

// add hadronic processes for Neutron 
//  
  G4ParticleDefinition* neutron = G4Neutron::Neutron();
  G4ProcessManager* neutron_manager = neutron->GetProcessManager();

  // process:  elastic scattering
  //
  neutron_manager->AddDiscreteProcess(theElasticProcess);

  // process: inelastic interaction
  //     
  G4NeutronInelasticProcess* neutronInelasticprocess 
                                  =  new G4NeutronInelasticProcess();
  
  // model
  //                              
  neutronInelasticprocess->RegisterMe(theModel);
  
  // cross Section
  //
  neutronInelasticprocess->AddDataSet(new G4NucleonNuclearCrossSection());
  
  // add process
  //  
  neutron_manager->AddDiscreteProcess(neutronInelasticprocess);

// process : Capture
//  
  G4HadronCaptureProcess * theNeutronCapture = new G4HadronCaptureProcess();
  G4LCapture * theNeutronCaptureModel = new G4LCapture();
  theNeutronCapture->RegisterMe(theNeutronCaptureModel);
  neutron_manager->AddDiscreteProcess(theNeutronCapture);

// process : Fission
//    
  G4HadronFissionProcess * theNeutronFission = new G4HadronFissionProcess();
  G4LFission * theNeutronFissionModel = new G4LFission();
  theNeutronFission->RegisterMe(theNeutronFissionModel);
  neutron_manager->AddDiscreteProcess(theNeutronFission);

//  code for higher energies
/* 

   G4TheoFSGenerator * theHEModel = new G4TheoFSGenerator();
   G4FTFModel * theStringModel = new G4FTFModel();
   G4ExcitedStringDecay * theStringFrag = 
                    new G4ExcitedStringDecay(new G4LundStringFragmentation());
   theStringModel->SetFragmentationModel(theStringFrag);
        
   theHEModel->SetTransport(new G4BinaryCascade());
   theHEModel->SetMinEnergy(4.*GeV);
   theHEModel->SetMaxEnergy(100*TeV);
   theHEModel->SetHighEnergyGenerator(theStringModel);
   
   theModel->SetMaxEnergy(5*GeV);  // reduce use of casacde to below 5 GeV

   protonInelasticprocess->RegisterMe(theHEModel);
  
   neutronInelasticprocess->RegisterMe(theHEModel);
   
*/    
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......


#include "G4StepLimiter.hh"

void PhysicsList::AddStepMax()
{
  // Step limitation seen as a process
  G4StepLimiter* stepMaxProcess = new G4StepLimiter();

  theParticleIterator->reset();
  while ((*theParticleIterator)()){
      G4ParticleDefinition* particle = theParticleIterator->value();
      G4ProcessManager* pmanager = particle->GetProcessManager();
      if (particle->GetPDGCharge() != 0.0)
        {
          pmanager ->AddDiscreteProcess(stepMaxProcess);
        }
  }
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

void PhysicsList::SetCuts()
{
  //G4VUserPhysicsList::SetCutsWithDefault method sets 
  //the default cut value for all particle types 
  //
  SetCutsWithDefault();
     
  if (verboseLevel>0) DumpCutValuesTable();
}

//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......

---