Message: Re: Problems with low energy physics processes Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Problems with low energy physics processes 

Forum: Physics List
Re: Question Problems with low energy physics processes (Hans Boie)
Date: 24 Dec, 2004
From: Vladimir IVANTCHENKO <vnivanch@mail.cern.ch>

On Mon, 20 Dec 2004, Hans Boie wrote:

> *** Discussion title: Physics List
> Email replies to PublicHyperNews@slac.stanford.edu must include:
>   In-Reply-To: <"/phys-list/113"@geant4-hn.slac.stanford.edu>
>   Subject: ...change this to be about your reply.
> 
> Hi!
> 
> I'm simulating my experiment in Geant4 (mainly three germanium and two
> silicon detectors) - which works great so far - thanks for your good
> work!!
> 
> Now i included the muonic background (cosmic muons) but, as described in
> the manual, i do not get the X-rays from the lead shield of my
> experiment. So i tried to include low energy physics procresses by
> looking at one of the Geant4 examples and copying the physics list.
> 
> In fact it seems that i get the X-rays but the simulation slowed down
> very much and i loose all other lines in my germanium detector. I tried
> to vary the cuts but also this did not help very much.
> 
> I would be glad if someone who is more familar with the low energy
> processes could give me hint, what i did wrong.
> 
> This is the physics list class i use:
> 
> BraZPhysicsList::BraZPhysicsList(BraZData* Data)
>   : G4VUserPhysicsList(), data(Data)
> {
>    defaultCutValue = data->default_cut_length;
>    SetVerboseLevel(1);
> }
> 
> BraZPhysicsList::~BraZPhysicsList() {}
> 
> void BraZPhysicsList::ConstructParticle()
> {
>   ConstructBosons();
>   ConstructLeptons();
>   ConstructBaryons();
>   ConstructIons();
> }
> 
> void BraZPhysicsList::ConstructBosons() {
> 
>   // gamma
>   G4Gamma::GammaDefinition();
>   // optical photon
>   G4OpticalPhoton::OpticalPhotonDefinition();
> }
> 
> void BraZPhysicsList::ConstructLeptons()
> {
>   // leptons
>   //  e+/-
>   G4Electron::ElectronDefinition();
>   G4Positron::PositronDefinition();
>   // mu+/-
>   G4MuonPlus::MuonPlusDefinition();
>   G4MuonMinus::MuonMinusDefinition();
>   // nu_e
>   G4NeutrinoE::NeutrinoEDefinition();
>   G4AntiNeutrinoE::AntiNeutrinoEDefinition();
>   // nu_mu
>   G4NeutrinoMu::NeutrinoMuDefinition();
>   G4AntiNeutrinoMu::AntiNeutrinoMuDefinition();
> }
> 
> void BraZPhysicsList::ConstructBaryons()
> {
>   //  barions
>   G4Proton::ProtonDefinition();
>   G4AntiProton::AntiProtonDefinition();
> 
>   G4Neutron::NeutronDefinition();
>   G4AntiNeutron::AntiNeutronDefinition();
> }
> 
> void BraZPhysicsList::ConstructProcess()
> {
>   AddTransportation();
>   ConstructEM();
>   ConstructGeneral();
> }
> 
> void BraZPhysicsList::ConstructIons()
> {
>   G4Alpha::AlphaDefinition();
> }
> 
> void BraZPhysicsList::ConstructEM()
> {
>   theParticleIterator->reset();
> 
>   while( (*theParticleIterator)() )
>     {
>       G4ParticleDefinition* particle = theParticleIterator->value();
>       G4ProcessManager* pmanager = particle->GetProcessManager();
>       G4String particleName = particle->GetParticleName();
> 
>       if(particleName == "gamma")
>         {
>           if(!data->uselowenphysics)
>             {
>               // gamma
>               pmanager->AddDiscreteProcess(new G4GammaConversion());
>               pmanager->AddDiscreteProcess(new G4ComptonScattering());      
>               pmanager->AddDiscreteProcess(new G4PhotoElectricEffect());
>             }
>           else
>             {
>               G4LowEnergyPhotoElectric* theLowEnPhoto
>                 = new G4LowEnergyPhotoElectric();
> 
>               theLowEnPhoto->ActivateAuger(true);
>               theLowEnPhoto->SetCutForLowEnSecPhotons(data->fluorescence_cut);
>               theLowEnPhoto->SetCutForLowEnSecElectrons(data->fluorescence_cut);
>               G4VProcess* theLowEnCompton = new G4LowEnergyCompton();
>               G4VProcess* theLowEnPair = new G4LowEnergyGammaConversion();
>               G4VProcess* theLowEnRayleigh = new G4LowEnergyRayleigh();
> 
>               pmanager->AddDiscreteProcess(theLowEnPhoto);
>               pmanager->AddDiscreteProcess(theLowEnCompton);
>               pmanager->AddDiscreteProcess(theLowEnPair);
>               pmanager->AddDiscreteProcess(theLowEnRayleigh);
>             }     
>         }
>       else if(particleName == "e-")
>         {
>           if(!data->uselowenphysics)
>             {
>               //electron
>               G4VProcess* theeminusMultipleScattering
>                 = new G4MultipleScattering();
>               G4VProcess* theeminusBremsstrahlung = new G4eBremsstrahlung();
>               G4VProcess* theeminusIonisation = new G4eIonisation();
> 
>               // add processes
>               pmanager->AddProcess(theeminusMultipleScattering);
>               pmanager->AddProcess(theeminusBremsstrahlung);
>               pmanager->AddProcess(theeminusIonisation);
> 
>               // set ordering for AlongStepDoIt
>               pmanager->SetProcessOrdering(theeminusMultipleScattering,
>                                            idxAlongStep,1);
>               pmanager->SetProcessOrdering(theeminusIonisation,
>                                            idxAlongStep,2);
> 
>               // set ordering for PostStepDoIt
>               pmanager->SetProcessOrdering(theeminusMultipleScattering,
>                                            idxPostStep,1);
>               pmanager->SetProcessOrdering(theeminusBremsstrahlung,
>                                            idxPostStep,3);
>               pmanager->SetProcessOrdering(theeminusIonisation,
>                                            idxPostStep,2);
>             }
>           else
>             {
>               G4LowEnergyBremsstrahlung* theLowEnBremss
>                 = new G4LowEnergyBremsstrahlung();
>               theLowEnBremss->SetCutForLowEnSecPhotons(data->fluorescence_cut);
>               G4MultipleScattering* aMultipleScattering
>                 = new G4MultipleScattering();
>               G4LowEnergyIonisation* theLowEnIon
>                 = new G4LowEnergyIonisation();
>               theLowEnIon->ActivateAuger(true);
>               theLowEnIon->SetCutForLowEnSecPhotons(data->fluorescence_cut);
>               theLowEnIon->SetCutForLowEnSecElectrons(data->fluorescence_cut);
> 
>               pmanager->AddProcess(aMultipleScattering, -1, 1, 1);
>               pmanager->AddProcess(theLowEnBremss, -1, -1, 3);
>               pmanager->AddProcess(theLowEnIon, -1, 2, 2);
>             }
>         }
>       else if(particleName == "e+")
>         {
>           //positron
>           G4VProcess* theeplusMultipleScattering = new G4MultipleScattering();
>           G4VProcess* theeplusIonisation = new G4eIonisation();
>           G4VProcess* theeplusBremsstrahlung = new G4eBremsstrahlung();
>           G4VProcess* theeplusAnnihilation = new G4eplusAnnihilation();
> 
>           // add processes
>           pmanager->AddProcess(theeplusMultipleScattering);
>           pmanager->AddProcess(theeplusIonisation);
>           pmanager->AddProcess(theeplusBremsstrahlung);
>           pmanager->AddProcess(theeplusAnnihilation);
> 
>           // set ordering for AtRestDoIt
>           pmanager->SetProcessOrderingToFirst(theeplusAnnihilation, idxAtRest);
> 
>           // set ordering for AlongStepDoIt
>           pmanager->SetProcessOrdering(theeplusMultipleScattering,
>                                        idxAlongStep,1);
>           pmanager->SetProcessOrdering(theeplusIonisation,
>                                        idxAlongStep,2);
> 
>           // set ordering for PostStepDoIt
>           pmanager->SetProcessOrdering(theeplusMultipleScattering,
>                                        idxPostStep,1);
>           pmanager->SetProcessOrdering(theeplusIonisation,
>                                        idxPostStep,2);
>           pmanager->SetProcessOrdering(theeplusBremsstrahlung,
>                                        idxPostStep,3);
>           pmanager->SetProcessOrdering(theeplusAnnihilation,
>                                        idxPostStep,4);
>         }
>       else if( particleName == "mu+" || particleName == "mu-"    )
>         {
>           //muon  
>           G4VProcess* aMultipleScattering = new G4MultipleScattering();
>           G4VProcess* aBremsstrahlung = new G4MuBremsstrahlung();
>           G4VProcess* aPairProduction = new G4MuPairProduction();
>           G4VProcess* anIonisation = new G4MuIonisation();
> 
>           // add processes
>           pmanager->AddProcess(anIonisation);
>           pmanager->AddProcess(aMultipleScattering);
>           pmanager->AddProcess(aBremsstrahlung);
>           pmanager->AddProcess(aPairProduction);
> 
>           // set ordering for AlongStepDoIt
>           pmanager->SetProcessOrdering(aMultipleScattering, idxAlongStep,1);
>           pmanager->SetProcessOrdering(anIonisation, idxAlongStep,2);
> 
>           // set ordering for PostStepDoIt
>           pmanager->SetProcessOrdering(aMultipleScattering, idxPostStep,1);
>           pmanager->SetProcessOrdering(anIonisation, idxPostStep,2);
>           pmanager->SetProcessOrdering(aBremsstrahlung, idxPostStep,3);
>           pmanager->SetProcessOrdering(aPairProduction, idxPostStep,4);
>         }
>       else if((particleName == "proton"     ||
>                particleName == "alpha"      ||
>                particleName == "deuteron"   ||
>                particleName == "triton"     ||
>                particleName == "He3"        ||
>                particleName == "GenericIon" )
>               && data->uselowenphysics )
>         {
>           G4MultipleScattering* aMultipleScattering
>             = new G4MultipleScattering();
>           G4hLowEnergyIonisation* ahadronLowEIon
>             = new G4hLowEnergyIonisation();
>           pmanager->AddProcess(aMultipleScattering,-1,1,1);
>           pmanager->AddProcess(ahadronLowEIon,-1,2,2);
>         }
>       else if ((!particle->IsShortLived()) &&
>                (particle->GetPDGCharge() != 0.0) && 
>                (particle->GetParticleName() != "chargedgeantino"))
>         {
>           // all others charged particles except geantino
>           G4VProcess* aMultipleScattering = new G4MultipleScattering();
>           G4VProcess* anIonisation = new G4hIonisation();
> 
>           // add processes
>           pmanager->AddProcess(anIonisation);
>           pmanager->AddProcess(aMultipleScattering);
> 
>           // set ordering for AlongStepDoIt
>           pmanager->SetProcessOrdering(aMultipleScattering, idxAlongStep,1);
>           pmanager->SetProcessOrdering(anIonisation, idxAlongStep,2);
> 
>           // set ordering for PostStepDoIt
>           pmanager->SetProcessOrdering(aMultipleScattering, idxPostStep,1);
>           pmanager->SetProcessOrdering(anIonisation, idxPostStep,2);
>         }
>     }
> }
> 
> #include "G4Decay.hh"
> void BraZPhysicsList::ConstructGeneral()
> {
>   // 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);
>         }
>     }
> }
> 
> void BraZPhysicsList::SetCuts()
> {
>   SetCutsWithDefault();
> 
>   if(data->uselowenphysics)
>     {
>       SetCutValue(data->cut_for_gamma,"gamma");
>       SetCutValue(data->cut_for_electron,"e-");
>       SetCutValue(data->cut_for_electron,"e+");
>       SetCutValue(data->cut_for_proton, "proton");
>       SetCutValue(data->cut_for_proton, "alpha");
>     }
> 
>   if (verboseLevel>0)
>     DumpCutValuesTable();
> }
> 
> 

Hello,

In general PhysicsList is correct. Some comments:

- cut values for proton and alha have no effect

- muon ionisation does not provide PIXE 

- if you really need it, then the only way is to define the ionisation for
muons as for hadrons and forget about high energy muons processes

In any case, simulation of deexcitations currently is slow because it
requires simulation of the deexcitation cascade in an atom. The only
advice for the time being is to set cuts on Auger and fluorescence to be
as big as possible. If Auger is not interesting in your application -
switch it off, that will improve performance.

VI

Inline Depth:
 1 1
 All All
Outline Depth:
 1 1
 2 2
 All All
Add message: (add)

 Add Message Add Message
to: "Re: Problems with low energy physics processes"

 Subscribe Subscribe

This site runs SLAC HyperNews version 1.11-slac-98, derived from the original HyperNews


[ Geant 4 Home | Geant 4 HyperNews | Search | Request New Forum | Feedback ]