Message: Re: Why an Ion of H3 does not generate any secondary? Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: Why an Ion of H3 does not generate any secondary? 

Forum: Hadronic Processes
Re: None Why an Ion of H3 does not generate any secondary? (Jorge Cabello)
Re: None Re: Why an Ion of H3 does not generate any secondary? (Vladimir IVANTCHENKO )
Re: None Re: Why an Ion of H3 does not generate any secondary? (Jorge Cabello)
Date: 01 Dec, 2006
From: Vladimir IVANTCHENKO <vnivanch@mail.cern.ch>

On Fri, 1 Dec 2006, Jorge wrote:

> *** Discussion title: Hadronic Processes
> Email replies to PublicHyperNews@slac.stanford.edu must include:
>   In-Reply-To: <"/hadronprocess/604/1/1"@geant4-hn.slac.stanford.edu>
>   Subject: ...change this to be about your reply.
> 
> Dear Vladimir,
> 
> I have introduced the changes the you suggested but the result is still
> the same:
> 
> ***************************************************************************
> * G4Track Information:   Particle = triton,   Track ID = 1,   Parent ID = 0
> ***************************************************************************
> 
> Step#    X(mm)    Y(mm)    Z(mm) KinE(MeV)  dE(MeV) StepLeng TrackLeng  NextVolume ProcName
>     0        0        0        0         0        0        0         0   Epitaxial initStep
>     1        0        0        0         0        0        0         0   Epitaxial User Limit
>   trackID: 1  energy deposit: 0 eV   position: 0 0 0 fm
> 
> Event no 1. Total energy deposit: 0 keV
> 
> As you can see the process in step 1 is "User Limit" while if I run the
> same code with C14 (higher energy) the process is "radioactive decay":
> 
> *****************************************************************************
> * G4Track Information:   Particle = C14[0.0],   Track ID = 1,   Parent ID = 0
> *****************************************************************************
> 
> Step#    X(mm)    Y(mm)    Z(mm) KinE(MeV)  dE(MeV) StepLeng TrackLeng  NextVolume ProcName
>     0        0        0        0         0        0        0         0   Epitaxial initStep
> G4ProcessManageer:: copy constructor
>     1        0        0        0         0        0        0         0   Epitaxial RadioactiveDecay
>     :----- List of 2ndaries - #SpawnInStep=  3(Rest= 3,Along= 0,Post= 0), #SpawnTotal=  3 ---------------
>     :         0         0         0  3.45e-06           N14[0.0]
>     :         0         0         0    0.0945          anti_nu_e
>     :         0         0         0     0.062                 e-
>     :----------------------------------------------------------------- EndOf2ndaries Info ---------------
>   trackID: 1  energy deposit: 0 eV   position: 0 0 0 fm
> 
> Event no 1. Total energy deposit: 45.60589138679122 keV
>

I would understand the cout as following: C14 has zero kinetic energy and 
radioactive decay. The triton has zero kinetic energy and is killed 
because of G4UserLimit, which is defined by you. Is it limit in time 
or in kinetic energy? Nothing to do with cuts.

Below in the Physics List everything is ok except one minor detail: the
radioactive decay is assigned to GenericIons. The triton (as alpha and
He3) are not Generic Ions - they are separate particles. So, if you need
to have radiactivity for the triton then add it additionally.

VI

 
> I don't know if I am not applying the cuts in the right way. Here is my
> whole code in the Physics list, I think these is going to be easier to
> detect the problem:
> 
> #include "globals.hh"
> #include "Exp2PhysicsList.hh"
> 
> #include "G4ProcessManager.hh"
> #include "G4ParticleTypes.hh"
> 
> #include "G4IonConstructor.hh"
> #include "G4ParticleDefinition.hh"
> 
> Exp2PhysicsList::Exp2PhysicsList():  G4VUserPhysicsList()
> {
>   defaultCutValue = .005*um;  
>    SetVerboseLevel(1);
> }
> 
> Exp2PhysicsList::~Exp2PhysicsList() {}
> 
> void Exp2PhysicsList::ConstructParticle()
> {
>   ConstructBosons();
>   ConstructLeptons();
>   ConstructMesons();
>   ConstructBaryons();
>   ConstructIons();
> }
> 
> void Exp2PhysicsList::ConstructBosons()
> {
>   // pseudo-particles
>   G4Geantino::GeantinoDefinition();
>   G4ChargedGeantino::ChargedGeantinoDefinition();
> 
>   // gamma
>   G4Gamma::GammaDefinition();
> }
> 
> void Exp2PhysicsList::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 Exp2PhysicsList::ConstructMesons()
> {
>   //  mesons
>   //    light mesons
>   G4PionPlus::PionPlusDefinition();
>   G4PionMinus::PionMinusDefinition();
>   G4PionZero::PionZeroDefinition();
>   G4Eta::EtaDefinition();
>   G4EtaPrime::EtaPrimeDefinition();
>   G4KaonPlus::KaonPlusDefinition();
>   G4KaonMinus::KaonMinusDefinition();
>   G4KaonZero::KaonZeroDefinition();
>   G4AntiKaonZero::AntiKaonZeroDefinition();
>   G4KaonZeroLong::KaonZeroLongDefinition();
>   G4KaonZeroShort::KaonZeroShortDefinition();
> }
> 
> void Exp2PhysicsList::ConstructBaryons()
> {
>   //  barions
>   G4Proton::ProtonDefinition();
>   G4AntiProton::AntiProtonDefinition();
> 
>   G4Neutron::NeutronDefinition();
>   G4AntiNeutron::AntiNeutronDefinition();
> }
> 
> void Exp2PhysicsList::ConstructIons()
> {
>   //  ions
>   G4IonConstructor iConstructor;
>   iConstructor.ConstructParticle();
> }
> 
> void Exp2PhysicsList::ConstructProcess()
> {
>   AddTransportation();
>   ConstructEM();  
>   ConstructGeneral();
> }
> 
> #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"
> #include "G4ionIonisation.hh"
> #include "G4RadioactiveDecay.hh"
> #include "G4StepLimiter.hh"
> #include "G4UserSpecialCuts.hh"
> 
> // low energy e-
> 
> #include "G4LowEnergyIonisation.hh" 
> #include "G4LowEnergyBremsstrahlung.hh" 
> 
> void Exp2PhysicsList::ConstructEM()
> {
>   theParticleIterator->reset();
>   while( (*theParticleIterator)() ){
>     G4ParticleDefinition* particle = theParticleIterator->value();
>     G4ProcessManager* pmanager = particle->GetProcessManager();
>     G4String particleName = particle->GetParticleName();
> 
>     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);
>       //pmanager->AddProcess(new G4UserSpecialCuts,   -1,-1,4);
> 
>       /*loweIon  = new G4LowEnergyIonisation("LowEnergyIoni");
>       loweBrem = new G4LowEnergyBremsstrahlung("LowEnBrem");
>       loweBrem->SetAngularGenerator("tsai");
> 
>       pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
>       pmanager->AddProcess(loweIon,     -1, 2,2);
>       pmanager->AddProcess(loweBrem,    -1,-3,3);*/
> 
>     } else if (particleName == "GenericIon") {  
>       // ion         
>       pmanager->AddProcess(new G4MultipleScattering,-1,1,1);
>       pmanager->AddProcess(new G4hIonisation,-1,2,2);
>       pmanager->AddProcess(new G4RadioactiveDecay,-1,3,3);
> 
>     } else if (particleName == "gamma") {
>       // gamma         
>       pmanager->AddDiscreteProcess(new G4PhotoElectricEffect);
>       pmanager->AddDiscreteProcess(new G4ComptonScattering);
>       pmanager->AddDiscreteProcess(new G4GammaConversion);
> 
>     }  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->IsShortLived()) &&
>                (particle->GetPDGCharge() != 0.0) && 
>                (particle->GetParticleName() != "chargedgeantino")) {
>       //all others charged particles except geantino
>       pmanager->AddProcess(new G4MultipleScattering,-1, 1,1);
>       pmanager->AddProcess(new G4hIonisation,       -1, 2,2); 
>     }
>   }
> }
> 
> #include "G4Decay.hh"
> 
> void Exp2PhysicsList::ConstructGeneral() {
> 
>   // Declare radioactive decay to the GenericIon in the IonTable.
>   //
>   const G4IonTable *theIonTable = G4ParticleTable::GetParticleTable()->GetIonTable();  
>   G4GenericIon* ion = G4GenericIon::GenericIon();
>   //Initialize these members to pass the "IsApplicable" method to run "G4RadioactiveDecay"
>   ion->SetAtomicNumber(1);
>   ion->SetAtomicMass(1);
>   G4RadioactiveDecay*  theRadioactiveDecay = new G4RadioactiveDecay("RadioactiveDecay");
>   ion->GetProcessManager()->AddProcess(theRadioactiveDecay, 0, -1, 3);
> 
>   // 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);
>     }
>   }
> 
>   for (G4int i=0; i<theIonTable->Entries(); i++)
>     {
>       G4String particleName = theIonTable->GetParticle(i)->GetParticleName();      
>       if (particleName == "GenericIon")
>         {
>           G4ProcessManager* pmanager = theIonTable->GetParticle(i)->GetProcessManager();
>           pmanager->SetVerboseLevel(3);
>           //pmanager ->AddProcess(theRadioactiveDecay);   
>           pmanager ->SetProcessOrdering(theRadioactiveDecay, idxPostStep);
>           pmanager ->SetProcessOrdering(theRadioactiveDecay, idxAtRest);
>         }
>     }
> 
> }
> 
> void Exp2PhysicsList::SetCuts()
> {      
>   defaultCutValue = 0.001*um;
>   SetCutsWithDefault();
> 
>   G4double lowlimit=1.*eV; 
>   G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(lowlimit,100*GeV);
> }
> 
> Another thing about the low energy package is that I don't know if I
> HAVE to use it or I can use it if I want to obtain a better accuracy for
> low energy electrons (as you can see I have commented that piece of
> code). At the moment I am more concerned with making the code works than
> with the accuracy so maybe it's better not to use it, at least till the
> program is working.
> 
> Thanks a lot,
> 

-- 
    Vladimir Ivanchenko  
    ### CERN PH SFT Tel: +41-22-76-78-871  Vladimir.Ivantchenko@cern.ch ###

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

1 None: Re: Why an Ion of H3 does not generate any secondary?   (Vladimir IVANTCHENKO - 01 Dec, 2006)
1 None: Re: Why an Ion of H3 does not generate any secondary?   (Jorge - 01 Dec, 2006)
3 None: Re: Why an Ion of H3 does not generate any secondary?   (Jorge - 05 Dec, 2006)
 Add Message Add Message
to: "Re: Why an Ion of H3 does not generate any secondary?"

 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 ]