Message: G4.9.4 vs G4.9.5 - Helium fragmentation tail disappeared Not Logged In (login)
 Next-in-Thread Next-in-Thread
 Next-in-Forum Next-in-Forum

None G4.9.4 vs G4.9.5 - Helium fragmentation tail disappeared 

Forum: Hadronic Processes
Date: 12 Jun, 2013
From: Hermann Fuchs <Hermann Fuchs>

Hello

I have recently upgraded from G4.9.4.p01 to G4.9.5.p02. I am using the same program with a self designed modular physics list.

Setup: Helium ion beam (100MeV/Nuc) in a water phantom. Longitudinal energy deposition is stored in root file.

Using G4.9.4.p01 the longitudinal energy deposition is correct, including a fragmentation tail after the Bragg Peak.

Using G4.9.5.p02 the position of the Bragg Peak is correct, however the longitudinal energy deposition misses the fragmentation tail after the Bragg Peak.

As I suppose it is related to hadronic processes I will post my physics list source code below. Is there a way to post the source code of the whole application?

Using one of the predefined recommended physics lists gives the correct behaviour. Unfortunately I have to use modular physics lists.

I would greatly appreciate your help.

Best regards, Hermann

#####################

    G4ProcessManager * pManager = 0;
   // Pi+ Physics
   pManager = G4PionPlus::PionPlus()->GetProcessManager();

   // add processes
   G4HadronElasticProcess* theppElasticProcess 
                         = new G4HadronElasticProcess();
   G4LElastic* theppElasticModel = new G4LElastic();
   theppElasticProcess->RegisterMe(theppElasticModel);
   pManager->AddDiscreteProcess(theppElasticProcess);

   G4PionPlusInelasticProcess* thePionPlusInelasticProcess 
                         = new G4PionPlusInelasticProcess(); 

   G4LEPionPlusInelastic* thePionPlusLEPModel = new G4LEPionPlusInelastic();
   G4HEPionPlusInelastic* thePionPlusHEPModel = new G4HEPionPlusInelastic();
   thePionPlusInelasticProcess->RegisterMe(thePionPlusLEPModel);
   thePionPlusInelasticProcess->RegisterMe(thePionPlusHEPModel);
   pManager->AddDiscreteProcess(thePionPlusInelasticProcess);

   G4VProcess* theppMultipleScattering = new G4hMultipleScattering();
   G4VProcess* theppIonisation        = new G4hIonisation();
   // 
   pManager->AddProcess(theppIonisation);
   pManager->AddProcess(theppMultipleScattering);
   // 
   // set ordering for AlongStepDoIt
   pManager->SetProcessOrdering(theppMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(theppIonisation,        idxAlongStep,2);
   // 
   // set ordering for PostStepDoIt
   pManager->SetProcessOrdering(theppMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(theppIonisation,        idxPostStep,2);

   // Pi- Physics
   pManager = G4PionMinus::PionMinus()->GetProcessManager();

   G4HadronElasticProcess* thepmElasticProcess 
                         = new G4HadronElasticProcess();
   G4LElastic* thepmElasticModel = new G4LElastic();
   thepmElasticProcess->RegisterMe(thepmElasticModel);
   pManager->AddDiscreteProcess(thepmElasticProcess);

   G4PionMinusInelasticProcess* thePionMinusInelasticProcess 
                         = new G4PionMinusInelasticProcess(); 

   G4LEPionMinusInelastic* thePionMinusLEPModel = new G4LEPionMinusInelastic();
   G4HEPionMinusInelastic* thePionMinusHEPModel = new G4HEPionMinusInelastic();
   thePionMinusInelasticProcess->RegisterMe(thePionMinusLEPModel);
   thePionMinusInelasticProcess->RegisterMe(thePionMinusHEPModel);
   pManager->AddDiscreteProcess(thePionMinusInelasticProcess);

   G4VProcess* thepmMultipleScattering = new G4hMultipleScattering();
   G4VProcess* thepmIonisation        = new G4hIonisation();
   // 
   // add processes
   pManager->AddProcess(thepmIonisation);
   pManager->AddProcess(thepmMultipleScattering);
   // 
   // set ordering for AlongStepDoIt
   pManager->SetProcessOrdering(thepmMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thepmIonisation,        idxAlongStep,2);
   // 
   // set ordering for PostStepDoIt
   pManager->SetProcessOrdering(thepmMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thepmIonisation,        idxPostStep,2);

   // Proton Physics
   pManager = G4Proton::Proton()->GetProcessManager();

   G4HadronElasticProcess* thepElasticProcess 
                         = new G4HadronElasticProcess();
   G4LElastic* thepElasticModel = new G4LElastic();
   thepElasticProcess->RegisterMe(thepElasticModel);
   pManager->AddDiscreteProcess(thepElasticProcess);

   G4ProtonInelasticProcess* thepInelasticProcess 
                         = new G4ProtonInelasticProcess(); 

   G4BinaryCascade* thepBCModel = new G4BinaryCascade();
   thepBCModel->SetMinEnergy(0. * MeV);
   thepBCModel->SetMaxEnergy(10. * GeV);

   thepInelasticProcess->RegisterMe(thepBCModel);

   pManager->AddDiscreteProcess(thepInelasticProcess);

   G4VProcess* thepMultipleScattering = new G4hMultipleScattering();
   G4VProcess* thepIonisation        = new G4hIonisation();

   pManager->AddProcess(thepIonisation);
   pManager->AddProcess(thepMultipleScattering);

   // set ordering for AlongStepDoIt
   pManager->SetProcessOrdering(thepMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thepIonisation,        idxAlongStep,2);
   // 
   // set ordering for PostStepDoIt
   pManager->SetProcessOrdering(thepMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thepIonisation,        idxPostStep,2);

   // anti-proton Physics
   pManager = G4AntiProton::AntiProton()->GetProcessManager();

   // add process
   G4HadronElasticProcess* theapElasticProcess 
                         = new G4HadronElasticProcess();
   G4LElastic* theapElasticModel = new G4LElastic();
   theapElasticProcess->RegisterMe(theapElasticModel);
   pManager->AddDiscreteProcess(theapElasticProcess);

   G4AntiProtonInelasticProcess* theapInelasticProcess 
                         = new G4AntiProtonInelasticProcess(); 

   G4LEAntiProtonInelastic* theapLEPModel = new G4LEAntiProtonInelastic();
   G4HEAntiProtonInelastic* theapHEPModel = new G4HEAntiProtonInelastic();
   theapInelasticProcess->RegisterMe(theapLEPModel);
   theapInelasticProcess->RegisterMe(theapHEPModel);
   pManager->AddDiscreteProcess(theapInelasticProcess);

   G4AntiProtonAnnihilationAtRest* theapAnnihilation
                            =  new G4AntiProtonAnnihilationAtRest();
   pManager->AddRestProcess(theapAnnihilation);

   G4VProcess* theapMultipleScattering = new G4hMultipleScattering();
   G4VProcess* theapIonisation        = new G4hIonisation();

   pManager->AddProcess(theapIonisation);
   pManager->AddProcess(theapMultipleScattering);
   //
   pManager->SetProcessOrdering(theapMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(theapIonisation,        idxAlongStep,2);
   // 
   pManager->SetProcessOrdering(theapMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(theapIonisation,        idxPostStep,2);

   // neutron Physics
   pManager = G4Neutron::Neutron()->GetProcessManager();

   G4HadronElasticProcess* thenElasticProcess 
                         = new G4HadronElasticProcess();
   G4LElastic* thenElasticModel = new G4LElastic();
   thenElasticProcess->RegisterMe(thenElasticModel);
   pManager->AddDiscreteProcess(thenElasticProcess);

   G4NeutronInelasticProcess* thenInelasticProcess 
                         = new G4NeutronInelasticProcess(); 

   G4NeutronHPInelasticData* neutronHPInelasticCrossSection 
                         = new G4NeutronHPInelasticData();

   G4NeutronHPInelastic* thenHPIModel 
                         = new G4NeutronHPInelastic();   
   G4BinaryCascade* thenBCModel = new G4BinaryCascade();
   G4double neutronHPInelasticMaxEnergy = 20. * MeV;
   G4double neutronBinaryMinEnergy = 19.0 * MeV;
   G4double neutronBinaryMaxEnergy = 10. * GeV;
   thenHPIModel->SetMaxEnergy(neutronHPInelasticMaxEnergy);
   thenBCModel->SetMinEnergy(neutronBinaryMinEnergy);
   thenBCModel->SetMaxEnergy(neutronBinaryMaxEnergy);

   thenInelasticProcess->RegisterMe(thenBCModel);
   thenInelasticProcess->RegisterMe(thenHPIModel);
   thenInelasticProcess->AddDataSet(neutronHPInelasticCrossSection);

   pManager->AddDiscreteProcess(thenInelasticProcess);

// fission

   G4HadronFissionProcess* thenFission
                         = new G4HadronFissionProcess();

   G4NeutronHPFissionData* thenHPFissionCrossSection 
                         = new G4NeutronHPFissionData();

   G4LFission* thenLFissionModel = new G4LFission();
    G4NeutronHPFission* thenHPFissionModel = new G4NeutronHPFission();

   G4double thenLFissionMinEnergy = 19.8 * MeV;
   G4double thenLFissionMaxEnergy = 100. * TeV;
   G4double thenHPFissionMinEnergy = 0. * MeV;
   G4double thenHPFissionMaxEnergy = 20. * MeV;

   thenLFissionModel->SetMinEnergy(thenLFissionMinEnergy);
   thenLFissionModel->SetMaxEnergy(thenLFissionMaxEnergy);
   thenHPFissionModel->SetMinEnergy(thenHPFissionMinEnergy);
   thenHPFissionModel->SetMaxEnergy(thenHPFissionMaxEnergy); 

   thenFission->RegisterMe(thenLFissionModel);
   thenFission->RegisterMe(thenHPFissionModel);
   thenFission->AddDataSet(thenHPFissionCrossSection);

   pManager->AddDiscreteProcess(thenFission);

// capture

   G4HadronCaptureProcess* thenCapture
                         = new G4HadronCaptureProcess();

   G4NeutronHPCaptureData* thenHPCaptureCrossSection 
                         = new G4NeutronHPCaptureData();

   G4LCapture* thenLCaptureModel = new G4LCapture();
   G4NeutronHPCapture* thenHPCaptureModel = new G4NeutronHPCapture();

   G4double thenLCaptureMinEnergy = 19.8 * MeV;
   G4double thenLCaptureMaxEnergy = 100. * TeV;
   G4double thenHPCaptureMinEnergy = 0. * MeV;
   G4double thenHPCaptureMaxEnergy = 20. * MeV;

   thenLCaptureModel->SetMinEnergy(thenLCaptureMinEnergy);
   thenLCaptureModel->SetMaxEnergy(thenLCaptureMaxEnergy);
   thenHPCaptureModel->SetMinEnergy(thenHPCaptureMinEnergy);
   thenHPCaptureModel->SetMaxEnergy(thenHPCaptureMaxEnergy);

   thenCapture->RegisterMe(thenLCaptureModel);
   thenCapture->RegisterMe(thenHPCaptureModel);
   thenCapture->AddDataSet(thenHPCaptureCrossSection );

   pManager->AddDiscreteProcess(thenCapture);

   // anti-neutron Physics
   pManager = G4AntiNeutron::AntiNeutron()->GetProcessManager();

   G4HadronElasticProcess* theanElasticProcess 
                         = new G4HadronElasticProcess();
   G4LElastic* theanElasticModel = new G4LElastic();
   theanElasticProcess->RegisterMe(theanElasticModel);
   pManager->AddDiscreteProcess(theanElasticProcess);

   G4AntiNeutronInelasticProcess* theanInelasticProcess 
                         = new G4AntiNeutronInelasticProcess(); 

   G4LEAntiNeutronInelastic* theanLEPModel = new G4LEAntiNeutronInelastic();
   G4HEAntiNeutronInelastic* theanHEPModel = new G4HEAntiNeutronInelastic();
   theanInelasticProcess->RegisterMe(theanLEPModel);
   theanInelasticProcess->RegisterMe(theanHEPModel);
   pManager->AddDiscreteProcess(theanInelasticProcess);

   G4AntiNeutronAnnihilationAtRest* theanAnnihilation
                            =  new G4AntiNeutronAnnihilationAtRest();
   pManager->AddRestProcess(theanAnnihilation);

// Generic ion

   pManager = G4GenericIon::GenericIon()->GetProcessManager();

   G4VProcess* thegionMultipleScattering = new G4hMultipleScattering();
   G4VProcess* thegionIonisation        = new G4ionIonisation();

   // Inelastic process
   G4HadronInelasticProcess* thegionInelastic =
       new G4HadronInelasticProcess("IonInelastic",G4GenericIon::GenericIon());
   thegionInelastic->AddDataSet(new G4TripathiCrossSection);
   thegionInelastic->AddDataSet(new G4IonsShenCrossSection);
   G4BinaryLightIonReaction* thegionBCModel = new G4BinaryLightIonReaction;

   thegionInelastic->AddDataSet(new G4TripathiLightCrossSection);
   thegionInelastic->AddDataSet(new G4TripathiCrossSection);
   thegionInelastic->AddDataSet(new G4IonsShenCrossSection);

   thegionInelastic->RegisterMe(thegionBCModel);

// Elastic process
   G4HadronElasticProcess* thegionElasticProcess = 
       new G4HadronElasticProcess();
   G4LElastic* thegionElasticModel = new G4LElastic();
   thegionElasticProcess->RegisterMe(thegionElasticModel);

   //
   pManager->AddProcess(thegionIonisation);
   pManager->AddProcess(thegionMultipleScattering);
   pManager->AddDiscreteProcess(thegionInelastic);
   pManager->AddDiscreteProcess(thegionElasticProcess);
   //
   pManager->SetProcessOrdering(thegionMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thegionIonisation,        idxAlongStep,2);
   //
   pManager->SetProcessOrdering(thegionMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thegionIonisation,        idxPostStep,2);

// radioactive decay

   G4RadioactiveDecay* radDecayProc = new G4RadioactiveDecay();
   pManager-> AddDiscreteProcess(radDecayProc);

// Deuteron

/*   G4double ionLEPMinEnergy = 0. * MeV;
   G4double ionLEPMaxEnergy = 100. * MeV;
   G4double ionBinaryMinEnergy = 80. * MeV;
   G4double ionBinaryMaxEnergy = 40. * GeV; */

   pManager = G4Deuteron::Deuteron()->GetProcessManager();

   G4HadronElasticProcess* thedeuElasticProcess
                         = new G4HadronElasticProcess();
   G4LElastic* thedeuElasticModel = new G4LElastic();
   thedeuElasticProcess->RegisterMe(thedeuElasticModel);
   pManager->AddDiscreteProcess(thedeuElasticProcess);

   G4DeuteronInelasticProcess* theDeuteronInelasticProcess
                         = new G4DeuteronInelasticProcess();

/*   G4LEDeuteronInelastic* theDeuteronLEPModel = new G4LEDeuteronInelastic();
   theDeuteronLEPModel->SetMinEnergy(ionLEPMinEnergy);
   theDeuteronLEPModel->SetMaxEnergy(ionLEPMaxEnergy); */

   G4BinaryLightIonReaction* theDeuteronBLIModel = new G4BinaryLightIonReaction;

   theDeuteronInelasticProcess->AddDataSet(new G4TripathiLightCrossSection);
   theDeuteronInelasticProcess->AddDataSet(new G4TripathiCrossSection);
   theDeuteronInelasticProcess->AddDataSet(new G4IonsShenCrossSection);

/*   theDeuteronBLIModel->SetMinEnergy(ionBinaryMinEnergy);
   theDeuteronBLIModel->SetMaxEnergy(ionBinaryMaxEnergy); */

//   theDeuteronInelasticProcess->RegisterMe(theDeuteronLEPModel);
   theDeuteronInelasticProcess->RegisterMe(theDeuteronBLIModel);

   pManager->AddDiscreteProcess(theDeuteronInelasticProcess);

   G4VProcess* thedeuMultipleScattering = new G4hMultipleScattering();
   G4VProcess* thedeuIonisation        = new G4hIonisation();
   //
   pManager->AddProcess(thedeuIonisation);
   pManager->AddProcess(thedeuMultipleScattering);
   //
   pManager->SetProcessOrdering(thedeuMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thedeuIonisation,        idxAlongStep,2);
   //
   pManager->SetProcessOrdering(thedeuMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thedeuIonisation,        idxPostStep,2);

   // Triton

   pManager = G4Triton::Triton()->GetProcessManager();

   G4HadronElasticProcess* thetriElasticProcess
                         = new G4HadronElasticProcess();
   G4LElastic* thetriElasticModel = new G4LElastic();
   thetriElasticProcess->RegisterMe(thetriElasticModel);
   pManager->AddDiscreteProcess(thetriElasticProcess);

   G4TritonInelasticProcess* thetriInelasticProcess
                         = new G4TritonInelasticProcess();

/*   G4LETritonInelastic* thetriLEPModel = new G4LETritonInelastic();
   thetriLEPModel->SetMinEnergy(ionLEPMinEnergy);
   thetriLEPModel->SetMaxEnergy(ionLEPMaxEnergy); */

   G4BinaryLightIonReaction* thetriBLIModel = new G4BinaryLightIonReaction;

   thetriInelasticProcess->AddDataSet(new G4TripathiLightCrossSection);
   thetriInelasticProcess->AddDataSet(new G4TripathiCrossSection);
   thetriInelasticProcess->AddDataSet(new G4IonsShenCrossSection);

/*   thetriBLIModel->SetMinEnergy(ionBinaryMinEnergy);
   thetriBLIModel->SetMaxEnergy(ionBinaryMaxEnergy); */

//   thetriInelasticProcess->RegisterMe(thetriLEPModel);
   thetriInelasticProcess->RegisterMe(thetriBLIModel);

   pManager->AddDiscreteProcess(thetriInelasticProcess);

   G4VProcess* thetriMultipleScattering = new G4hMultipleScattering();
   G4VProcess* thetriIonisation        = new G4hIonisation();
   //
   pManager->AddProcess(thetriIonisation);
   pManager->AddProcess(thetriMultipleScattering);
   //
   pManager->SetProcessOrdering(thetriMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thetriIonisation,        idxAlongStep,2);
   //
   pManager->SetProcessOrdering(thetriMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thetriIonisation,        idxPostStep,2);

   // Alpha
   pManager = G4Alpha::Alpha()->GetProcessManager();

   G4HadronElasticProcess* thealElasticProcess
                         = new G4HadronElasticProcess();
   G4LElastic* thealElasticModel = new G4LElastic();
   thealElasticProcess->RegisterMe(thealElasticModel);
   pManager->AddDiscreteProcess(thealElasticProcess);

   G4AlphaInelasticProcess* thealInelasticProcess
                         = new G4AlphaInelasticProcess();

/*   G4LEAlphaInelastic* thealLEPModel = new G4LEAlphaInelastic();
   thealLEPModel->SetMinEnergy(ionLEPMinEnergy);
   thealLEPModel->SetMaxEnergy(ionLEPMaxEnergy); */

   G4BinaryLightIonReaction* thealBLIModel = new G4BinaryLightIonReaction;

   thealInelasticProcess->AddDataSet(new G4TripathiLightCrossSection);
   thealInelasticProcess->AddDataSet(new G4TripathiCrossSection);
   thealInelasticProcess->AddDataSet(new G4IonsShenCrossSection);

/*   thealBLIModel->SetMinEnergy(ionBinaryMinEnergy);
   thealBLIModel->SetMaxEnergy(ionBinaryMaxEnergy); */

//   thealInelasticProcess->RegisterMe(thealLEPModel);
   thealInelasticProcess->RegisterMe(thealBLIModel);

   pManager->AddDiscreteProcess(thealInelasticProcess);

   G4VProcess* thealpMultipleScattering = new G4hMultipleScattering();
   G4VProcess* thealpIonisation        = new G4ionIonisation();

   //
   pManager->AddProcess(thealpIonisation);
   pManager->AddProcess(thealpMultipleScattering);
   //
   pManager->SetProcessOrdering(thealpMultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thealpIonisation,        idxAlongStep,2);
   //
   pManager->SetProcessOrdering(thealpMultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thealpIonisation,        idxPostStep,2);

   // He3
   pManager = G4He3::He3()->GetProcessManager();

   G4HadronElasticProcess* thehe3ElasticProcess
                         = new G4HadronElasticProcess();
   G4LElastic* thehe3ElasticModel = new G4LElastic();
   thehe3ElasticProcess->RegisterMe(thehe3ElasticModel);
   pManager->AddDiscreteProcess(thehe3ElasticProcess);

   G4VProcess* thehe3MultipleScattering = new G4hMultipleScattering();
   G4VProcess* thehe3Ionisation        = new G4ionIonisation();
   //
   pManager->AddProcess(thehe3Ionisation);
   pManager->AddProcess(thehe3MultipleScattering);
   //
   pManager->SetProcessOrdering(thehe3MultipleScattering, idxAlongStep,1);
   pManager->SetProcessOrdering(thehe3Ionisation,        idxAlongStep,2);
   //
   pManager->SetProcessOrdering(thehe3MultipleScattering, idxPostStep,1);
   pManager->SetProcessOrdering(thehe3Ionisation,        idxPostStep,2);

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

1 None: Re: G4.9.4 vs G4.9.5 - Helium fragmentation tail disappeared   (Dennis H. Wright - 14 Jun, 2013)
2 Idea: Re: G4.9.4 vs G4.9.5 - Helium fragmentation tail disappeared   (Vladimir Ivanchenko - 20 Jun, 2013)
 Add Message Add Message
to: "G4.9.4 vs G4.9.5 - Helium fragmentation tail disappeared"

 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 ]