Message: Re: I125 Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None Re: I125 

Forum: Physics List
Re: Question I125 (Baljeet Seniwal)
Date: Jan 01, 14:05
From: Baljeet Seniwal <Baljeet Seniwal>

I modelled the electron decay spectrum of I125 with the same physics list that I am using to calculate the water sphere. I got this spectrum (attachment). When I compared it with literature (Howell 1992)(attached) I discovered that one peak is missing between 10keV and 100keV. Can you please guide me is it because of variation in spectras or my mistakes in PhysicsList?

This my PhysicsList that I am using

#include "PhysicsList.hh"

#include "G4EmDNAPhysics.hh"
#include "G4EmDNAPhysics_option1.hh"
#include "G4EmDNAPhysics_option2.hh"
#include "G4EmDNAPhysics_option3.hh"
#include "G4EmDNAPhysics_option4.hh"
#include "G4EmDNAPhysics_option5.hh"
#include "G4EmDNAPhysics_option6.hh"
#include "G4EmDNAPhysics_option7.hh"
#include "G4EmStandardPhysics_option3.hh"
#include<G4Ions.hh>
#include<G4DNAGenericIonsManager.hh>
#include "G4UnitsTable.hh"
#include "G4SystemOfUnits.hh"
#include "G4UserSpecialCuts.hh"

// particles

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

#include<G4IonTable.hh>
#include<G4EmParameters.hh>
#include<G4DecayPhysics.hh>
#include<G4RadioactiveDecayPhysics.hh>
#include<G4HadronElasticPhysicsHP.hh>
#include<G4HadronPhysicsFTFP_BERT_HP.hh>
#include<G4IonElasticPhysics.hh>
#include<G4IonPhysics.hh>
#include "G4EmExtraPhysics.hh"
#include<G4RadioactiveDecay.hh>
#include<G4LossTableManager.hh>
#include<G4UAtomicDeexcitation.hh>
#include<G4NuclideTable.hh>
#include<G4NuclearLevelData.hh>
#include "G4IonINCLXXPhysics.hh"
#include "G4UserSpecialCuts.hh"
#include "G4PhysicsListHelper.hh"
#include "G4Electron.hh"
#include<G4EmLivermorePhysics.hh>

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

PhysicsList::PhysicsList() : G4VModularPhysicsList()

{

    // EM physics
    RegisterPhysics( new G4EmDNAPhysics());
    //RegisterPhysics(new G4DecayPhysics());
    //RegisterPhysics(new G4RadioactiveDecayPhysics());

    // Ion Elastic scattering
   // RegisterPhysics( new G4IonElasticPhysics(verb));

    // Ion Inelastic physics
    //RegisterPhysics( new G4IonPhysics(verb));
    // RegisterPhysics(new G4EmLivermorePhysics());
    //  RegisterPhysics( new G4EmDNAPhysics_option1());
//    G4EmParameters* param = G4EmParameters::Instance();
//    param->SetAugerCascade(true);
//    param->SetStepFunction(0.1, 0.1*CLHEP::nm);
//    param->SetStepFunctionMuHad(0.1, 0.1*CLHEP::nm);
//    RegisterPhysics(new G4DecayPhysics());

}

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

PhysicsList::~PhysicsList() {

}

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

void PhysicsList::ConstructParticle()
{
    G4BosonConstructor  pBosonConstructor;
    pBosonConstructor.ConstructParticle();

    G4LeptonConstructor pLeptonConstructor;
    pLeptonConstructor.ConstructParticle();

    G4MesonConstructor pMesonConstructor;
    pMesonConstructor.ConstructParticle();

    G4BaryonConstructor pBaryonConstructor;
    pBaryonConstructor.ConstructParticle();

    G4IonConstructor pIonConstructor;
    pIonConstructor.ConstructParticle();

    G4ShortLivedConstructor pShortLivedConstructor;
    pShortLivedConstructor.ConstructParticle();

    G4DNAGenericIonsManager* genericIonsManager;
    genericIonsManager=G4DNAGenericIonsManager::Instance();
    genericIonsManager->GetIon("alpha++");
    genericIonsManager->GetIon("alpha+");
    genericIonsManager->GetIon("helium");
    genericIonsManager->GetIon("hydrogen");
    genericIonsManager->GetIon( "carbon");

}

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

#include "G4ProcessManager.hh"
#include "G4EmProcessOptions.hh"

#include<G4LossTableManager.hh>
#include<G4UAtomicDeexcitation.hh>
#include<G4RadioactiveDecay.hh>
#include<G4DecayProcessType.hh>
#include<G4Decay.hh>
#include "G4ComptonScattering.hh"
#include "G4GammaConversion.hh"
#include "G4PhotoElectricEffect.hh"

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

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

#include "G4hMultipleScattering.hh"
#include "G4hIonisation.hh"
#include "G4hBremsstrahlung.hh"
#include "G4hPairProduction.hh"

#include "G4ionIonisation.hh"
#include "G4IonParametrisedLossModel.hh"
#include "G4NuclearStopping.hh"
#include<G4NucleusLimits.hh>
void PhysicsList::ConstructProcess() {

    // transportation
    //
    AddTransportation();
    G4EmDNAPhysics().ConstructProcess();
    G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
    G4RadioactiveDecay* radioactiveDecay = new G4RadioactiveDecay();
     ph->RegisterProcess(radioactiveDecay, G4GenericIon::GenericIon());
    AddTrackingCut();
}

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

void PhysicsList::AddTrackingCut() {

    defaultCutValue = 1. * nm;
    SetCutValue(defaultCutValue, "gamma");
    SetCutValue(defaultCutValue, "e-");
    SetCutValue(defaultCutValue, "e+");
    G4double lowLimit = 250. * eV;
    G4double highLimit = 100. * GeV;

    G4ProductionCutsTable::GetProductionCutsTable()->SetEnergyRange(lowLimit,
                                                                    highLimit);

    // Print the cuts
    if (verboseLevel>0) DumpCutValuesTable();

}

PrimaryGeneratorAction.cc is

#include "PrimaryGeneratorAction.hh"
#include<G4ParticleTable.hh>
#include<G4Event.hh>
#include<G4SystemOfUnits.hh>
#include<Randomize.hh>
#include<G4ParticleGun.hh>
#include<G4ParticleDefinition.hh>
#include<G4PhysicalConstants.hh>
#include<G4RandomDirection.hh>
#include<G4Isotope.hh>
#include<G4ParticleTable.hh>
#include<G4IonTable.hh>
#include<G4Geantino.hh>

#include"DetectorConstruction.hh"
#include<G4RunManager.hh>
#include<G4StateManager.hh>
//#include<G4VStateDependent.hh>

PrimaryGeneratorAction::PrimaryGeneratorAction() :
    G4VUserPrimaryGeneratorAction(), fParticleGun(0),fDetector(0)
{
    fDetector = dynamic_cast<const DetectorConstruction*>(G4RunManager::GetRunManager()->GetUserDetectorConstruction());

    G4int n_particle = 1;
    fParticleGun  = new G4ParticleGun(n_particle);

    fParticleGun->SetParticleEnergy(0*eV);
    // vertex volume
    //

            fRmin= 0. *micrometer;
            fRmax = 4. *micrometer;

            //opening angle
            //
            G4double alphaMin =  0.*deg;
            G4double alphaMax = 360.*deg;
            fCosAlphaMin = std::cos(alphaMin);
            fCosAlphaMax = std::cos(alphaMax);

}

PrimaryGeneratorAction::~PrimaryGeneratorAction() {

    if(fParticleGun) delete fParticleGun;

}

void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) {

    //vertex position uniform in spherical shell

            G4double cosTheta = 2*G4UniformRand() - 1;  //cosTheta uniform in [0, pi]
            G4double sinTheta = std::sqrt(1. - cosTheta*cosTheta);
            G4double phi      = twopi*G4UniformRand();  //phi uniform in [0, 2*pi]
            G4ThreeVector ur(sinTheta*std::cos(phi),sinTheta*std::sin(phi),cosTheta);

            G4double R3 = fRmin + G4UniformRand()*(fRmax - fRmin);

            fParticleGun->SetParticlePosition(R3*ur);

            //particle direction uniform around ur
            //
            //1- in World frame
            //cosAlpha uniform in [cos(alphaMin), cos(alphaMax)]
            G4double cosAlpha = fCosAlphaMin-G4UniformRand()*(fCosAlphaMin-fCosAlphaMax);
            G4double sinAlpha = std::sqrt(1. - cosAlpha*cosAlpha);
            G4double psi      = twopi*G4UniformRand();  //psi uniform in (0,2*pi)
            G4ThreeVector dir(sinAlpha*std::cos(psi),sinAlpha*std::sin(psi),cosAlpha);

            //2- rotate dir   (rotateUz transforms uz to ur)
            dir.rotateUz(ur);

            fParticleGun->SetParticleMomentumDirection(dir);
    if (fParticleGun->GetParticleDefinition() == G4Geantino::Geantino()) {
        G4int Z = 10, A = 24;
        G4double ionCharge   = 0.*eplus;
        G4double excitEnergy = 0.*keV;

        G4ParticleDefinition* ion
                = G4IonTable::GetIonTable()->GetIon(Z,A,excitEnergy);
        fParticleGun->SetParticleDefinition(ion);
        fParticleGun->SetParticleCharge(ionCharge);

    }

    //fParticleGun->SetParticlePosition(G4ThreeVector());

    //create vertex
    //
    fParticleGun->GeneratePrimaryVertex(anEvent);
}

   Attachment:
      http://hypernews.slac.stanford.edu/HyperNews/geant4/get/AUX/2019/01/01/14.01-18902-AAPM_SPECTRUM.png
      http://hypernews.slac.stanford.edu/HyperNews/geant4/get/AUX/2019/01/01/14.01-64191-lectron_spectrum_I123.png

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

1 Idea: Re: I125   (Vladimir Ivanchenko - Jan 03, 02:16)
2 Idea: Re: I125   (Vladimir Ivanchenko - Jan 05, 03:10)
1 Idea: Re: I125   (Vladimir Ivanchenko - Jan 05, 03:46)
 Add Message Add Message
to: "Re: I125"

 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 ]