From 4a60422ad2d41fc13da037681c6d7cd7daf956f4 Mon Sep 17 00:00:00 2001
From: matta <matta@npt>
Date: Fri, 8 Feb 2013 16:01:08 +0000
Subject: [PATCH] * Fix bug in NPReaction (shoot option) * clean up the Physics
list file (remove trailling comment and reindent) * fix cout issue in event
action
---
Inputs/EventGenerator/10He.reaction | 24 +-
NPLib/Physics/NPReaction.cxx | 4 +-
NPSimulation/src/EventAction.cc | 3 +-
NPSimulation/src/PhysicsList.cc | 382 +++++++++++++---------------
4 files changed, 196 insertions(+), 217 deletions(-)
diff --git a/Inputs/EventGenerator/10He.reaction b/Inputs/EventGenerator/10He.reaction
index e2b06f82d..dee2a70a3 100644
--- a/Inputs/EventGenerator/10He.reaction
+++ b/Inputs/EventGenerator/10He.reaction
@@ -13,11 +13,11 @@
ExcitationEnergy= 1.4
Beam
- Particle= 11Li
- Energy= 550
- SigmaEnergy= 0
- SigmaThetaX= 0.01
- SigmaPhiY= 0.01
+ Particle= 7Li
+ Energy= 20
+ SigmaEnergy= 20
+ SigmaThetaX= 25
+ SigmaPhiY= 25
SigmaX= 5
SigmaY= 5
MeanThetaX= 0
@@ -30,20 +30,20 @@ Beam
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-TwoBodyReaction
- Beam= 11Li
+%TwoBodyReaction
+ Beam= 9Li
Target= 2H
Light= 3He
- Heavy= 10He
+ Heavy= 8He
ExcitationEnergyLight= 0.0
- ExcitationEnergyHeavy= 10.0
- CrossSectionPath= 11Li(d,3He)10He.txt CS10He
- ShootLight= 1
+ ExcitationEnergyHeavy= 0.0
+ CrossSectionPath= flat.txt CS10He
+ ShootLight= 0
ShootHeavy= 1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
-ParticleDecay 10He
+%ParticleDecay 10He
Daughter= 9He n
ExcitationEnergy= 0 0
DifferentialCrossSection= 11Li(d,3He)10He.txt particle9He
diff --git a/NPLib/Physics/NPReaction.cxx b/NPLib/Physics/NPReaction.cxx
index 727e6aa5d..1b06f63ff 100644
--- a/NPLib/Physics/NPReaction.cxx
+++ b/NPLib/Physics/NPReaction.cxx
@@ -320,7 +320,7 @@ void Reaction::ReadConfigurationFile(string Path){
ReactionFile >> DataBuffer;
if(atoi(DataBuffer.c_str()) == 0 )
- fshoot3 = true;
+ fshoot3 = false;
if(fVerboseLevel==1 && fshoot3) cout << "Shoot 3 : Yes " << endl;
else if (fVerboseLevel==1 ) cout << "Shoot 3 : No " << endl;
@@ -331,7 +331,7 @@ void Reaction::ReadConfigurationFile(string Path){
ReactionFile >> DataBuffer;
if(atoi(DataBuffer.c_str()) == 0 )
- fshoot4 = true;
+ fshoot4 = false;
if(fVerboseLevel==1 && fshoot4) cout << "Shoot 4 : Yes " << endl;
else if (fVerboseLevel==1 ) cout << "Shoot 4 : No " << endl;
diff --git a/NPSimulation/src/EventAction.cc b/NPSimulation/src/EventAction.cc
index 45322725f..b2a0dd506 100644
--- a/NPSimulation/src/EventAction.cc
+++ b/NPSimulation/src/EventAction.cc
@@ -32,6 +32,7 @@
#include "DetectorConstruction.hh"
#include "RootOutput.h"
+#include<iostream>
using namespace std;
@@ -52,7 +53,7 @@ void EventAction::BeginOfEventAction(const G4Event* event){
if ((event->GetEventID() + 1) % m_printModulo == 0)
// G4cout << "Event: " << event->GetEventID() + 1 << G4endl;
- G4cout << "\rEvent: " << event->GetEventID() + 1 << flush;
+ cout << "\rEvent: " << event->GetEventID() + 1 << flush;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/NPSimulation/src/PhysicsList.cc b/NPSimulation/src/PhysicsList.cc
index 78975bc4e..3e59f42e6 100644
--- a/NPSimulation/src/PhysicsList.cc
+++ b/NPSimulation/src/PhysicsList.cc
@@ -20,7 +20,7 @@
* A good improvement should be a modular physicis list in order to deal *
* accuratly with different physics cases. *
*****************************************************************************/
- #include "PhysicsList.hh"
+#include "PhysicsList.hh"
// I/O
#include "G4ios.hh"
@@ -45,7 +45,7 @@
#include "G4eBremsstrahlung.hh"
#include "G4eplusAnnihilation.hh"
-// Penelope
+// Penelope
//#include "G4PenelopeIonisation.hh"
//#include "G4PenelopeBremsstrahlung.hh"
//#include "G4PenelopeAnnihilation.hh"
@@ -98,74 +98,67 @@
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-PhysicsList::PhysicsList()
-{
- // ie: no secondaries
- defaultCutValue = 1000 * pc;
+PhysicsList::PhysicsList(){
+ // ie: no secondaries
+ defaultCutValue = 1000 * pc;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-PhysicsList::~PhysicsList()
-{
- ;
+PhysicsList::~PhysicsList(){
+ ;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructParticle()
-{
- // In this method, static member functions should be called
- // for all particles which you want to use.
- // This ensures that objects of these particle types will be
- // created in the program.
-
- //Usefull to test geometry
- G4Geantino::GeantinoDefinition();
-
- //Usefull for Gamma
- ConstructBosons();
-
- //Usefull for betaDecay
- ConstructLeptons();
-
- //Needed by G4 (4.9.2) to run on mac os X ;-)
- ConstructMesons();
-
- //usefull for p and n
- ConstructBaryons();
-
- //Usefull of course :p
- ConstructIons();
+void PhysicsList::ConstructParticle(){
+ // In this method, static member functions should be called
+ // for all particles which you want to use.
+ // This ensures that objects of these particle types will be
+ // created in the program.
+
+ //Usefull to test geometry
+ G4Geantino::GeantinoDefinition();
+
+ //Usefull for Gamma
+ ConstructBosons();
+
+ //Usefull for betaDecay
+ ConstructLeptons();
+
+ //Needed by G4 (4.9.2) to run on mac os X ;-)
+ ConstructMesons();
+
+ //usefull for p and n
+ ConstructBaryons();
+
+ //Usefull of course :p
+ ConstructIons();
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructBosons()
-{
- G4Geantino::GeantinoDefinition() ;
- G4ChargedGeantino::ChargedGeantinoDefinition() ;
- G4Gamma::GammaDefinition() ;
+void PhysicsList::ConstructBosons(){
+ G4Geantino::GeantinoDefinition() ;
+ G4ChargedGeantino::ChargedGeantinoDefinition() ;
+ G4Gamma::GammaDefinition() ;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructLeptons()
-{
- G4Electron::ElectronDefinition() ;
- G4Positron::PositronDefinition() ;
- G4NeutrinoE::NeutrinoEDefinition() ;
- G4AntiNeutrinoE::AntiNeutrinoEDefinition() ;
- //G4NeutrinoMu::NeutrinoMuDefinition() ;
- //G4AntiNeutrinoMu::AntiNeutrinoMuDefinition() ;
- //G4MuonPlus::MuonPlusDefinition() ;
- //G4MuonMinus::MuonMinusDefinition() ;
+void PhysicsList::ConstructLeptons(){
+ G4Electron::ElectronDefinition() ;
+ G4Positron::PositronDefinition() ;
+ G4NeutrinoE::NeutrinoEDefinition() ;
+ G4AntiNeutrinoE::AntiNeutrinoEDefinition() ;
+ //G4NeutrinoMu::NeutrinoMuDefinition() ;
+ //G4AntiNeutrinoMu::AntiNeutrinoMuDefinition() ;
+ //G4MuonPlus::MuonPlusDefinition() ;
+ //G4MuonMinus::MuonMinusDefinition() ;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructBaryons()
-{
- G4Proton::ProtonDefinition() ;
- G4AntiProton::AntiProtonDefinition() ;
- G4Neutron::NeutronDefinition();
+void PhysicsList::ConstructBaryons(){
+ G4Proton::ProtonDefinition() ;
+ G4AntiProton::AntiProtonDefinition() ;
+ G4Neutron::NeutronDefinition();
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructMesons()
-{
+void PhysicsList::ConstructMesons(){
// mesons
G4PionPlus ::PionPlusDefinition();
G4PionMinus ::PionMinusDefinition();
@@ -183,166 +176,151 @@ void PhysicsList::ConstructMesons()
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructIons()
-{
- G4He3::He3Definition() ;
- G4Deuteron::DeuteronDefinition() ;
- G4Triton::TritonDefinition() ;
- G4Alpha::AlphaDefinition() ;
-
- G4IonConstructor iConstructor ;
- iConstructor.ConstructParticle() ;
-
+void PhysicsList::ConstructIons(){
+ G4He3::He3Definition() ;
+ G4Deuteron::DeuteronDefinition() ;
+ G4Triton::TritonDefinition() ;
+ G4Alpha::AlphaDefinition() ;
+
+ G4IonConstructor iConstructor ;
+ iConstructor.ConstructParticle() ;
+
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::ConstructProcess()
-{
- AddTransportation() ;
- ConstructEM() ;
-
- SetCuts() ;
+void PhysicsList::ConstructProcess(){
+ AddTransportation() ;
+ ConstructEM() ;
+
+ SetCuts() ;
}
//....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
- //standard Geant4
- pmanager->AddDiscreteProcess(new G4PhotoElectricEffect) ;
- pmanager->AddDiscreteProcess(new G4ComptonScattering) ;
- pmanager->AddDiscreteProcess(new G4GammaConversion) ;
- //Low energy
- //pmanager->AddDiscreteProcess(new G4LowEnergyPhotoElectric) ;
- //pmanager->AddDiscreteProcess(new G4LowEnergyCompton) ;
- //pmanager->AddDiscreteProcess(new G4LowEnergyGammaConversion) ;
- //pmanager->AddDiscreteProcess(new G4LowEnergyRayleigh) ;
- // Penelope
- //pmanager->AddDiscreteProcess(new G4PenelopePhotoElectric) ;
- //pmanager->AddDiscreteProcess(new G4PenelopeCompton) ;
- //pmanager->AddDiscreteProcess(new G4PenelopeGammaConversion) ;
- //pmanager->AddDiscreteProcess(new G4PenelopeRayleigh) ;
-
- } else if (particleName == "e-") {
- //electron
- pmanager->AddProcess(new G4eMultipleScattering , -1, 1, 1) ;
- //standard geant4:
- pmanager->AddProcess(new G4eIonisation , -1, 2, 2) ;
- pmanager->AddProcess(new G4eBremsstrahlung , -1, -1, 3) ;
- // Low energy:
- //pmanager->AddProcess(new G4LowEnergyIonisation , -1, 2, 2) ;
- //pmanager->AddProcess(new G4LowEnergyBremsstrahlung , -1, -1, 3) ;
- // Penelope:
- // pmanager->AddProcess(new G4PenelopeIonisation , -1, 2, 2) ;
- // pmanager->AddProcess(new G4PenelopeBremsstrahlung , -1, -1, 3) ;
-
-
- } else if (particleName == "e+") {
- //positron
- pmanager->AddProcess(new G4eMultipleScattering , -1, 1, 1 );
- // standard Geant4 and Low energy
- pmanager->AddProcess(new G4eIonisation , -1, 2, 2 );
- pmanager->AddProcess(new G4eBremsstrahlung , -1, -1, 3 );
- pmanager->AddProcess(new G4eplusAnnihilation , 0, -1, 4 );
- //Penelope:
- //pmanager->AddProcess(new G4PenelopeIonisation , -1, 2, 2 );
- //pmanager->AddProcess(new G4PenelopeBremsstrahlung , -1, -1, 3 );
- //pmanager->AddProcess(new G4PenelopeAnnihilation , 0, -1, 4 );
-
-
-
- } else if (particleName == "mu+" ||
- particleName == "mu-") {
+void PhysicsList::ConstructEM(){
- //muon
- /* pmanager->AddProcess(new G4muMultipleScattering , -1, 1, 1 ) ;
- pmanager->AddProcess(new G4MuIonisation , -1, 2, 2 ) ;
- pmanager->AddProcess(new G4MuBremsstrahlung , -1, -1, 3 ) ;
- pmanager->AddProcess(new G4MuPairProduction , -1, -1, 4 ) ;*/
+ theParticleIterator->reset();
- } else if (particleName == "GenericIon") {
- pmanager->AddProcess(new G4hMultipleScattering(), -1, 1, 1) ;
- G4ionIonisation* iI = new G4ionIonisation ;
- // mod by Nicolas [07/05/09]
-// iI->ActivateNuclearStopping(true) ;
- iI->ActivateStoppingData(true) ;
- pmanager->AddProcess(iI , -1, 2, 2) ;
-
- //all others charged particles except geantino
- } else if ((!particle->IsShortLived()) &&
- (particle->GetPDGCharge() != 0.0) &&
- (particleName != "chargedgeantino")) {
-
- G4hIonisation* hI = new G4hIonisation ;
- // mod by Nicolas [07/05/09]
-// hI->ActivateNuclearStopping(true) ;
- pmanager->AddProcess(new G4hMultipleScattering , -1, 1, 1) ;
- pmanager->AddProcess(hI , -1, 2, 2) ;
-
-
- }//end else if
- }//end while particle
-
- G4EmProcessOptions opt ;
- opt.SetSubCutoff(true) ;
- opt.SetMinEnergy(0.001*eV) ;
- opt.SetMaxEnergy(1000.*MeV) ;
- opt.SetDEDXBinning(1000) ;
- opt.SetLambdaBinning(1000) ;
- // mod by Nicolas [07/05/09]
-// opt.SetLinearLossLimit(1.e-3) ;
-
+ while ((*theParticleIterator)()) {
+ G4ParticleDefinition* particle = theParticleIterator->value() ;
+ G4ProcessManager* pmanager = particle->GetProcessManager() ;
+ G4String particleName = particle->GetParticleName() ;
+
+ if (particleName == "gamma") {
+ // gamma
+ //standard Geant4
+ pmanager->AddDiscreteProcess(new G4PhotoElectricEffect) ;
+ pmanager->AddDiscreteProcess(new G4ComptonScattering) ;
+ pmanager->AddDiscreteProcess(new G4GammaConversion) ;
+ //Low energy
+ //pmanager->AddDiscreteProcess(new G4LowEnergyPhotoElectric) ;
+ //pmanager->AddDiscreteProcess(new G4LowEnergyCompton) ;
+ //pmanager->AddDiscreteProcess(new G4LowEnergyGammaConversion) ;
+ //pmanager->AddDiscreteProcess(new G4LowEnergyRayleigh) ;
+ // Penelope
+ //pmanager->AddDiscreteProcess(new G4PenelopePhotoElectric) ;
+ //pmanager->AddDiscreteProcess(new G4PenelopeCompton) ;
+ //pmanager->AddDiscreteProcess(new G4PenelopeGammaConversion) ;
+ //pmanager->AddDiscreteProcess(new G4PenelopeRayleigh) ;
+
+ }
+ else if (particleName == "e-") {
+ //electron
+ pmanager->AddProcess(new G4eMultipleScattering , -1, 1, 1) ;
+ //standard geant4:
+ pmanager->AddProcess(new G4eIonisation , -1, 2, 2) ;
+ pmanager->AddProcess(new G4eBremsstrahlung , -1, -1, 3) ;
+ // Low energy:
+ //pmanager->AddProcess(new G4LowEnergyIonisation , -1, 2, 2) ;
+ //pmanager->AddProcess(new G4LowEnergyBremsstrahlung , -1, -1, 3) ;
+ // Penelope:
+ // pmanager->AddProcess(new G4PenelopeIonisation , -1, 2, 2) ;
+ // pmanager->AddProcess(new G4PenelopeBremsstrahlung , -1, -1, 3) ;
+ }
+
+ else if (particleName == "e+") {
+ //positron
+ pmanager->AddProcess(new G4eMultipleScattering , -1, 1, 1 );
+ // standard Geant4 and Low energy
+ pmanager->AddProcess(new G4eIonisation , -1, 2, 2 );
+ pmanager->AddProcess(new G4eBremsstrahlung , -1, -1, 3 );
+ pmanager->AddProcess(new G4eplusAnnihilation , 0, -1, 4 );
+ //Penelope:
+ //pmanager->AddProcess(new G4PenelopeIonisation , -1, 2, 2 );
+ //pmanager->AddProcess(new G4PenelopeBremsstrahlung , -1, -1, 3 );
+ //pmanager->AddProcess(new G4PenelopeAnnihilation , 0, -1, 4 );
+
+ }
+
+ else if (particleName == "mu+" ||
+ particleName == "mu-") {
+ }
+
+ else if (particleName == "GenericIon") {
+ pmanager->AddProcess(new G4hMultipleScattering(), -1, 1, 1) ;
+ G4ionIonisation* iI = new G4ionIonisation ;
+ // mod by Nicolas [07/05/09]
+ iI->ActivateStoppingData(true) ;
+ pmanager->AddProcess(iI , -1, 2, 2) ;
+ }
+
+ else if ((!particle->IsShortLived()) &&
+ (particle->GetPDGCharge() != 0.0) &&
+ (particleName != "chargedgeantino")) {
+
+ G4hIonisation* hI = new G4hIonisation ;
+ // mod by Nicolas [07/05/09]
+ pmanager->AddProcess(new G4hMultipleScattering , -1, 1, 1) ;
+ pmanager->AddProcess(hI , -1, 2, 2) ;
+
+
+ }
+ }
+
+ G4EmProcessOptions opt ;
+ opt.SetSubCutoff(true) ;
+ opt.SetMinEnergy(0.001*eV) ;
+ opt.SetMaxEnergy(1000.*MeV) ;
+ opt.SetDEDXBinning(1000) ;
+ opt.SetLambdaBinning(1000) ;
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-void PhysicsList::SetCuts()
-{
- // uppress error messages even in case e/gamma/proton do not exist
- G4int temp = GetVerboseLevel();
- SetVerboseLevel(0);
- // " G4VUserPhysicsList::SetCutsWithDefault" method sets
- // the default cut value for all particle types
- SetCutsWithDefault();
-
- // for gamma-rays
-// SetCutValue(0.01*mm, "gamma");
-// SetCutValue(0.01*mm, "e-");
-// SetCutValue(0.01*mm, "e+");
-
- // Retrieve verbose level
- SetVerboseLevel(temp);
+void PhysicsList::SetCuts(){
+ // uppress error messages even in case e/gamma/proton do not exist
+ G4int temp = GetVerboseLevel();
+ SetVerboseLevel(0);
+ // " G4VUserPhysicsList::SetCutsWithDefault" method sets
+ // the default cut value for all particle types
+ SetCutsWithDefault();
+
+ // for gamma-rays
+ // SetCutValue(0.01*mm, "gamma");
+ // SetCutValue(0.01*mm, "e-");
+ // SetCutValue(0.01*mm, "e+");
+
+ // Retrieve verbose level
+ SetVerboseLevel(temp);
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
-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);
- }
- }
-//end Add Decay Process
-
+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);
+ }
+ }
+ //end Add Decay Process
+
}
-void PhysicsList::MyOwnConstruction()
-{
- ConstructDecay();
+void PhysicsList::MyOwnConstruction(){
+ ConstructDecay();
}
--
GitLab