From 8a7c458769f28037a4a1c58025e699f85d5e8a95 Mon Sep 17 00:00:00 2001
From: Greg <gchristian@tamu.edu>
Date: Tue, 11 Jul 2017 13:19:33 -0500
Subject: [PATCH] Added coincidence source ability to EventGeneratorIsotropic
---
.../EventGenerator/EventGeneratorIsotropic.cc | 115 ++++++++++--------
.../EventGenerator/EventGeneratorIsotropic.hh | 34 +++---
Projects/T40/configs/ConfigGeTamu.dat | 2 +-
3 files changed, 83 insertions(+), 68 deletions(-)
diff --git a/NPSimulation/EventGenerator/EventGeneratorIsotropic.cc b/NPSimulation/EventGenerator/EventGeneratorIsotropic.cc
index bd3e604b6..96397b44b 100644
--- a/NPSimulation/EventGenerator/EventGeneratorIsotropic.cc
+++ b/NPSimulation/EventGenerator/EventGeneratorIsotropic.cc
@@ -40,15 +40,7 @@ using namespace CLHEP;
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
EventGeneratorIsotropic::EventGeneratorIsotropic(){
- m_EnergyLow = 0 ;
- m_EnergyHigh = 0 ;
- m_x0 = 0 ;
- m_y0 = 0 ;
- m_z0 = 0 ;
- m_SigmaX = 0;
- m_SigmaY = 0;
m_ParticleStack = ParticleStack::getInstance();
- m_particle = NULL;
}
@@ -58,99 +50,118 @@ EventGeneratorIsotropic::~EventGeneratorIsotropic(){
}
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+EventGeneratorIsotropic::SourceParameters::SourceParameters(){
+ m_EnergyLow = 0 ;
+ m_EnergyHigh = 0 ;
+ m_x0 = 0 ;
+ m_y0 = 0 ;
+ m_z0 = 0 ;
+ m_SigmaX = 0 ;
+ m_SigmaY = 0 ;
+ m_particle = NULL;
+}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void EventGeneratorIsotropic::ReadConfiguration(NPL::InputParser parser){
vector<NPL::InputBlock*> blocks = parser.GetAllBlocksWithToken("Isotropic");
+ m_Parameters.reserve(blocks.size());
if(NPOptionManager::getInstance()->GetVerboseLevel())
cout << endl << "\033[1;35m//// Isotropic reaction found " << endl;
-
+
vector<string> token = {"EnergyLow","EnergyHigh","HalfOpenAngleMin","HalfOpenAngleMax","x0","y0","z0","Particle"};
for(unsigned int i = 0 ; i < blocks.size() ; i++){
if(blocks[i]->HasTokenList(token)){
- m_EnergyLow =blocks[i]->GetDouble("EnergyLow","MeV");
- m_EnergyHigh =blocks[i]->GetDouble("EnergyHigh","MeV");
- m_HalfOpenAngleMin =blocks[i]->GetDouble("HalfOpenAngleMin","deg");
- m_HalfOpenAngleMax =blocks[i]->GetDouble("HalfOpenAngleMax","deg");
- m_x0 =blocks[i]->GetDouble("x0","mm");
- m_y0 =blocks[i]->GetDouble("y0","mm");
- m_z0 =blocks[i]->GetDouble("z0","mm");
+ m_Parameters.push_back(SourceParameters());
+ const vector<SourceParameters>::reverse_iterator it = m_Parameters.rbegin();
+
+ it->m_EnergyLow =blocks[i]->GetDouble("EnergyLow","MeV");
+ it->m_EnergyHigh =blocks[i]->GetDouble("EnergyHigh","MeV");
+ it->m_HalfOpenAngleMin =blocks[i]->GetDouble("HalfOpenAngleMin","deg");
+ it->m_HalfOpenAngleMax =blocks[i]->GetDouble("HalfOpenAngleMax","deg");
+ it->m_x0 =blocks[i]->GetDouble("x0","mm");
+ it->m_y0 =blocks[i]->GetDouble("y0","mm");
+ it->m_z0 =blocks[i]->GetDouble("z0","mm");
vector<string> particleName =blocks[i]->GetVectorString("Particle");
for(unsigned int j = 0 ; j < particleName.size() ; j++){
- if(particleName[j]=="proton"){ m_particleName.push_back("1H") ;}
- else if(particleName[j]=="deuton"){ m_particleName.push_back("2H") ; }
- else if(particleName[j]=="triton"){ m_particleName.push_back("3H") ; }
- else if(particleName[j]=="3He" || particleName[j]=="He3") { m_particleName.push_back("3He") ; }
- else if(particleName[j]=="alpha") { m_particleName.push_back("4He") ; }
- else if(particleName[j]=="gamma") { m_particleName.push_back("gamma") ;}
- else if(particleName[j]=="neutron") {m_particleName.push_back("neutron") ;}
- else m_particleName.push_back(particleName[j]);
+ if(particleName[j]=="proton"){ it->m_particleName.push_back("1H") ;}
+ else if(particleName[j]=="deuton"){ it->m_particleName.push_back("2H") ; }
+ else if(particleName[j]=="triton"){ it->m_particleName.push_back("3H") ; }
+ else if(particleName[j]=="3He" || particleName[j]=="He3") { it->m_particleName.push_back("3He") ; }
+ else if(particleName[j]=="alpha") { it->m_particleName.push_back("4He") ; }
+ else if(particleName[j]=="gamma") { it->m_particleName.push_back("gamma") ;}
+ else if(particleName[j]=="neutron") {it->m_particleName.push_back("neutron") ;}
+ else it->m_particleName.push_back(particleName[j]);
}
if(blocks[i]->HasToken("ExcitationEnergy"))
- m_ExcitationEnergy =blocks[i]->GetVectorDouble("ExcitationEnergy","MeV");
+ it->m_ExcitationEnergy =blocks[i]->GetVectorDouble("ExcitationEnergy","MeV");
if(blocks[i]->HasToken("SigmaX"))
- m_SigmaX=blocks[i]->GetDouble("SigmaX","mm");
+ it->m_SigmaX=blocks[i]->GetDouble("SigmaX","mm");
if(blocks[i]->HasToken("SigmaY"))
- m_SigmaX=blocks[i]->GetDouble("SigmaY","mm");
+ it->m_SigmaY=blocks[i]->GetDouble("SigmaY","mm");
if(blocks[i]->HasToken("Multiplicity"))
- m_Multiplicty=blocks[i]->GetVectorInt("Multiplicity");
+ it->m_Multiplicty=blocks[i]->GetVectorInt("Multiplicity");
}
else{
cout << "ERROR: check your input file formatting \033[0m" << endl;
exit(1);
}
- }
-
- if(m_ExcitationEnergy.size()==0)
- m_ExcitationEnergy.resize(m_particleName.size(),0);
- if(m_Multiplicty.size()==0)
- m_Multiplicty.resize(m_particleName.size(),1);
-
+ }
+
+ for(auto& par : m_Parameters) {
+ if(par.m_ExcitationEnergy.size()==0)
+ par.m_ExcitationEnergy.resize(par.m_particleName.size(),0);
+ if(par.m_Multiplicty.size()==0)
+ par.m_Multiplicty.resize(par.m_particleName.size(),1);
+ }
}
- //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
- void EventGeneratorIsotropic::GenerateEvent(G4Event*){
- for(unsigned int p=0; p<m_particleName.size(); p++){
- for(int i=0; i<m_Multiplicty[p]; i++){
- m_particle=NULL;
- if(m_particle==NULL){
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+void EventGeneratorIsotropic::GenerateEvent(G4Event*){
+
+ for(auto& par : m_Parameters) {
+ for(unsigned int p=0; p<par.m_particleName.size(); p++){
+ for(int i=0; i<par.m_Multiplicty[p]; i++){
+ par.m_particle=NULL;
+ if(par.m_particle==NULL){
- if(m_particleName[p]=="gamma" || m_particleName[p]=="neutron" || m_particleName[p]=="opticalphoton"){
- m_particle = G4ParticleTable::GetParticleTable()->FindParticle(m_particleName[p].c_str());
+ if(par.m_particleName[p]=="gamma" || par.m_particleName[p]=="neutron" || par.m_particleName[p]=="opticalphoton"){
+ par.m_particle = G4ParticleTable::GetParticleTable()->FindParticle(par.m_particleName[p].c_str());
}
else{
- NPL::Nucleus* N = new NPL::Nucleus(m_particleName[p]);
- m_particle = G4ParticleTable::GetParticleTable()->GetIonTable()->GetIon(N->GetZ(), N->GetA(),m_ExcitationEnergy[p]);
+ NPL::Nucleus* N = new NPL::Nucleus(par.m_particleName[p]);
+ par.m_particle = G4ParticleTable::GetParticleTable()->GetIonTable()->GetIon(N->GetZ(), N->GetA(),par.m_ExcitationEnergy[p]);
delete N;
}
}
- G4double cos_theta_min = cos(m_HalfOpenAngleMin);
- G4double cos_theta_max = cos(m_HalfOpenAngleMax);
+ G4double cos_theta_min = cos(par.m_HalfOpenAngleMin);
+ G4double cos_theta_max = cos(par.m_HalfOpenAngleMax);
G4double cos_theta = cos_theta_min + (cos_theta_max - cos_theta_min) * RandFlat::shoot();
G4double theta = acos(cos_theta) ;
G4double phi = RandFlat::shoot() * 2 * pi ;
- G4double particle_energy = m_EnergyLow + RandFlat::shoot() * (m_EnergyHigh - m_EnergyLow) ;
+ G4double particle_energy = par.m_EnergyLow + RandFlat::shoot() * (par.m_EnergyHigh - par.m_EnergyLow) ;
// Direction of particle, energy and laboratory angle
G4double momentum_x = sin(theta) * cos(phi) ;
G4double momentum_y = sin(theta) * sin(phi) ;
G4double momentum_z = cos(theta) ;
- G4double x0 = RandGauss::shoot(m_x0,m_SigmaX);
- G4double y0 = RandGauss::shoot(m_y0,m_SigmaY);
+ G4double x0 = RandGauss::shoot(par.m_x0,par.m_SigmaX);
+ G4double y0 = RandGauss::shoot(par.m_y0,par.m_SigmaY);
- Particle particle(m_particle, theta,particle_energy,G4ThreeVector(momentum_x, momentum_y, momentum_z),G4ThreeVector(x0, y0, m_z0));
+ Particle particle(par.m_particle, theta,particle_energy,G4ThreeVector(momentum_x, momentum_y, momentum_z),G4ThreeVector(x0, y0, par.m_z0));
m_ParticleStack->AddParticleToStack(particle);
}
}
}
-
+
+}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/NPSimulation/EventGenerator/EventGeneratorIsotropic.hh b/NPSimulation/EventGenerator/EventGeneratorIsotropic.hh
index 4388bac99..d509f8409 100644
--- a/NPSimulation/EventGenerator/EventGeneratorIsotropic.hh
+++ b/NPSimulation/EventGenerator/EventGeneratorIsotropic.hh
@@ -44,22 +44,26 @@ public: // Inherit from VEventGenerator Class
void ReadConfiguration(NPL::InputParser) ;
void GenerateEvent(G4Event*) ;
void InitializeRootOutput() ;
-
+
private: // Source parameter from input file
- G4double m_EnergyLow ; // Lower limit of energy range
- G4double m_EnergyHigh ; // Upper limit of energy range
- G4double m_HalfOpenAngleMin ; // Min Half open angle of the source
- G4double m_HalfOpenAngleMax ; // Max Half open angle of the source
- G4double m_x0 ; // Vertex Position X
- G4double m_y0 ; // Vertex Position Y
- G4double m_z0 ; // Vertex Position Z
- G4double m_SigmaX ;
- G4double m_SigmaY ;
- G4ParticleDefinition* m_particle ; // Kind of particle to shoot isotropically
- vector<G4double> m_ExcitationEnergy ; // Excitation energy of the emitted particle
- vector<string> m_particleName ;
- ParticleStack* m_ParticleStack ;
- vector<G4int> m_Multiplicty;
+ struct SourceParameters {
+ SourceParameters() ;
+ G4double m_EnergyLow ; // Lower limit of energy range
+ G4double m_EnergyHigh ; // Upper limit of energy range
+ G4double m_HalfOpenAngleMin ; // Min Half open angle of the source
+ G4double m_HalfOpenAngleMax ; // Max Half open angle of the source
+ G4double m_x0 ; // Vertex Position X
+ G4double m_y0 ; // Vertex Position Y
+ G4double m_z0 ; // Vertex Position Z
+ G4double m_SigmaX ;
+ G4double m_SigmaY ;
+ G4ParticleDefinition* m_particle ; // Kind of particle to shoot isotropically
+ vector<G4double> m_ExcitationEnergy ; // Excitation energy of the emitted particle
+ vector<string> m_particleName ;
+ vector<G4int> m_Multiplicty ;
+ };
+ vector<SourceParameters> m_Parameters ;
+ ParticleStack* m_ParticleStack ;
};
#endif
diff --git a/Projects/T40/configs/ConfigGeTamu.dat b/Projects/T40/configs/ConfigGeTamu.dat
index 42066dea4..2cdacb530 100644
--- a/Projects/T40/configs/ConfigGeTamu.dat
+++ b/Projects/T40/configs/ConfigGeTamu.dat
@@ -15,4 +15,4 @@ ConfigGeTamu
%ADD_BACK_ARRAY
%ADD_BACK_FACING
%%%%%%%% Randomise raw energy within ADC bin
- ADC_RANDOM_BIN
+% ADC_RANDOM_BIN
--
GitLab