diff --git a/NPLib/Detectors/ChiNu/TChiNuData.cxx b/NPLib/Detectors/ChiNu/TChiNuData.cxx
index 93290619cafb5eb0dbd74ee5f8716ea514a8b213..d7487596938067a57ed1ecdf698243e667ec8265 100755
--- a/NPLib/Detectors/ChiNu/TChiNuData.cxx
+++ b/NPLib/Detectors/ChiNu/TChiNuData.cxx
@@ -50,6 +50,9 @@ void TChiNuData::Clear() {
// Time
fChiNu_T_DetectorNbr.clear();
fChiNu_Time.clear();
+ // PArticle type
+ fChiNu_IG_DetectorNbr.clear();
+ fChiNu_IsGamma.clear();
}
diff --git a/NPLib/Detectors/ChiNu/TChiNuData.h b/NPLib/Detectors/ChiNu/TChiNuData.h
index c4b01d371fe718e5dff20a9bb40a1f2c2f4fc030..24cb64fbdce1d3186503949a3cbbdea86d4e3c66 100755
--- a/NPLib/Detectors/ChiNu/TChiNuData.h
+++ b/NPLib/Detectors/ChiNu/TChiNuData.h
@@ -42,6 +42,10 @@ class TChiNuData : public TObject {
vector<UShort_t> fChiNu_T_DetectorNbr;
vector<Double_t> fChiNu_Time;
+ // Particle Type
+ vector<UShort_t> fChiNu_IG_DetectorNbr;
+ vector<Bool_t> fChiNu_IsGamma;
+
//////////////////////////////////////////////////////////////
// Constructor and destructor
@@ -77,6 +81,12 @@ class TChiNuData : public TObject {
fChiNu_Time.push_back(Time);
};//!
+ // Particle Type
+ inline void SetIsGamma(const UShort_t& DetNbr,const Bool_t& Type) {
+ fChiNu_IG_DetectorNbr.push_back(DetNbr);
+ fChiNu_IsGamma.push_back(Type);
+ };//!
+
////////////////////// GETTERS ////////////////////////
// Energy
@@ -95,6 +105,14 @@ class TChiNuData : public TObject {
inline Double_t Get_Time(const unsigned int &i) const
{return fChiNu_Time[i];}//!
+ // Particle Type
+ inline UShort_t GetMultIsGamma() const
+ {return fChiNu_IG_DetectorNbr.size();}
+ inline UShort_t GetIG_DetectorNbr(const unsigned int &i) const
+ {return fChiNu_IG_DetectorNbr[i];}//!
+ inline Double_t Get_IsGamma(const unsigned int &i) const
+ {return fChiNu_IsGamma[i];}//!
+
//////////////////////////////////////////////////////////////
// Required for ROOT dictionnary
diff --git a/NPLib/Detectors/ChiNu/TChiNuPhysics.cxx b/NPLib/Detectors/ChiNu/TChiNuPhysics.cxx
index bb1d4c3e8023a0a19155eb0bc9cac0bca5d0493b..49c445858ed4e0a39be89ab38221d2147819997d 100755
--- a/NPLib/Detectors/ChiNu/TChiNuPhysics.cxx
+++ b/NPLib/Detectors/ChiNu/TChiNuPhysics.cxx
@@ -92,6 +92,7 @@ void TChiNuPhysics::BuildPhysicalEvent() {
DetectorNumber.push_back(m_PreTreatedData->GetE_DetectorNbr(e));
Energy.push_back(m_PreTreatedData->Get_Energy(e));
Time.push_back(m_PreTreatedData->Get_Time(t));
+ IsGamma.push_back(m_PreTreatedData->Get_IsGamma(e));
}
}
}
@@ -125,6 +126,14 @@ void TChiNuPhysics::PreTreat() {
Double_t Time= Cal->ApplyCalibration("ChiNu/TIME"+NPL::itoa(m_EventData->GetT_DetectorNbr(i)),m_EventData->Get_Time(i));
m_PreTreatedData->SetTime(m_EventData->GetT_DetectorNbr(i), Time);
}
+
+ // IsGamma
+ mysize = m_EventData->GetMultIsGamma();
+ for (UShort_t i = 0; i < mysize; ++i) {
+ Bool_t IsGamma= m_EventData->Get_IsGamma(i);
+ m_PreTreatedData->SetIsGamma(m_EventData->GetIG_DetectorNbr(i), IsGamma);
+ }
+
}
@@ -198,6 +207,7 @@ void TChiNuPhysics::Clear() {
DetectorNumber.clear();
Energy.clear();
Time.clear();
+ IsGamma.clear();
}
diff --git a/NPLib/Detectors/ChiNu/TChiNuPhysics.h b/NPLib/Detectors/ChiNu/TChiNuPhysics.h
index 6afd1af4956f7936823b704fdf265f066c3cefa4..0c837db9fd8638d18ac44f4d01b5bab2eaa1973f 100755
--- a/NPLib/Detectors/ChiNu/TChiNuPhysics.h
+++ b/NPLib/Detectors/ChiNu/TChiNuPhysics.h
@@ -65,6 +65,7 @@ class TChiNuPhysics : public TObject, public NPL::VDetector {
vector<int> DetectorNumber;
vector<double> Energy;
vector<double> Time;
+ vector<bool> IsGamma;
/// A usefull method to bundle all operation to add a detector
void AddDetector(TVector3 POS);
diff --git a/NPSimulation/Detectors/ChiNu/ChiNu.cc b/NPSimulation/Detectors/ChiNu/ChiNu.cc
index 6328e77a107dad6248f65eff4d818ef39841c021..082f7ec88366feb4225705e31e00b7dcaf75f7ef 100755
--- a/NPSimulation/Detectors/ChiNu/ChiNu.cc
+++ b/NPSimulation/Detectors/ChiNu/ChiNu.cc
@@ -23,6 +23,7 @@
#include <sstream>
#include <cmath>
#include <limits>
+#include <regex>
//G4 Geometry object
#include "G4Tubs.hh"
#include "G4Box.hh"
@@ -45,6 +46,7 @@
// NPTool header
#include "ChiNu.hh"
#include "CalorimeterScorers.hh"
+#include "ProcessScorers.hh"
#include "InteractionScorers.hh"
#include "RootOutput.h"
#include "MaterialManager.hh"
@@ -318,31 +320,89 @@ void ChiNu::InitializeRootOutput(){
void ChiNu::ReadSensitive(const G4Event* ){
m_Event->Clear();
+
+ ///////////
+ // Process scorer
+ ProcessScorers::PS_Process* Process_scorer = (ProcessScorers::PS_Process*) m_ChiNuScorer->GetPrimitive(2);
+ unsigned int ProcessMult = Process_scorer->GetMult();
+ vector<int> IsGamma;
+
+ int m_VID=-1;
+ int g_hit=0;
+
+ for(unsigned int i = 0 ; i < ProcessMult ; i++){
+ string particle_name = Process_scorer->GetParticleName(i);
+ string volume_name = Process_scorer->GetVolumeName(i);
+ string process_name = Process_scorer->GetProcessName(i);
+
+
+ regex pattern("impr_(\\d+)_logic");
+ smatch match;
+ int VolumeNbr;
+
+ if(regex_search( volume_name,match,pattern))
+ VolumeNbr = stoi(match[1]);
+ else VolumeNbr = -1;
+
+ if(VolumeNbr!=m_VID && particle_name=="gamma")
+ {// if the first particle in the volume is a gamma, it will give a wrong time
+ // cout << "!!!!!!!! First particle is gamma !!!!!!!!!!" << endl;
+ // cout << "i=" << i << " particle_name=" << particle_name << endl;
+ // cout << "volume_name=" << volume_name << endl;
+ // cout << "process_name=" << process_name << endl;
+ g_hit=1;
+ }
+ else if(VolumeNbr!=m_VID)
+ g_hit=0;
+ else if(g_hit==1){
+ // Only if an electron ionizes the liquid scintillator will there be a measured energy deposit and associated time.
+ // cout << "//////////" << endl;
+ // cout << "i=" << i << " particle_name=" << particle_name << endl;
+ // cout << "volume_name=" << volume_name << endl;
+ // cout << "process_name=" << process_name << endl;
+ if(particle_name=="e-")
+ {
+ g_hit++;
+ IsGamma.push_back(m_VID);
+ }
+ }
+
+ m_VID=VolumeNbr;
+ }
+
///////////
// Calorimeter scorer
CalorimeterScorers::PS_Calorimeter* Scorer= (CalorimeterScorers::PS_Calorimeter*) m_ChiNuScorer->GetPrimitive(0);
- unsigned int size = Scorer->GetMult();
+ unsigned int size = Scorer->GetMult();
+ vector<int> DetID;
for(unsigned int i = 0 ; i < size ; i++){
vector<unsigned int> level = Scorer->GetLevel(i);
int DetectorNbr = level[0];
-
+ DetID.push_back(DetectorNbr);
// Light output equation:
// L(Ep)=A*Ep - B*(1-exp(-C*pow(Ep,D))) from F. Pino et al. Applied Radiation and Isotopes 89 (2014) 79-84.
double A=0.62, B=1.3, C=0.39, D=0.97;
double Ep=RandGauss::shoot(Scorer->GetEnergy(i),ChiNu_NS::ResoEnergyOffset+ChiNu_NS::ResoEnergySlope*Scorer->GetEnergy(i));
double Energy = A*Ep - B*(1-exp(-C*pow(Ep,D)));
////////////////////////
-
- // Energy = RandGauss::shoot(Energy,ChiNu_NS::ResoEnergyOffset+ChiNu_NS::ResoEnergySlope*Energy);
+
if(Energy>ChiNu_NS::EnergyThreshold){
// double Time = RandGauss::shoot(Scorer->GetTime(i),ChiNu_NS::ResoTime); // One time resolution for all detectors
double Time = RandGauss::shoot(Scorer->GetTime(i),m_ResoTime.at(DetectorNbr-1)); // Time resolution that can be set different for each detector
m_Event->SetEnergy(DetectorNbr,Energy);
- m_Event->SetTime(DetectorNbr,Time);
+ m_Event->SetTime(DetectorNbr,Time);
+
+ // Check if the detector was first lighted by a gamma ray, which would have given a wrong time.
+ auto it = find(IsGamma.begin(), IsGamma.end(), DetectorNbr);
+ if(it != IsGamma.end())
+ m_Event->SetIsGamma(DetectorNbr,true);
+ else
+ m_Event->SetIsGamma(DetectorNbr,false);
}
}
+ IsGamma.clear();
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
@@ -359,9 +419,11 @@ void ChiNu::InitializeScorers() {
vector<int> level; level.push_back(0);
G4VPrimitiveScorer* Calorimeter= new CalorimeterScorers::PS_Calorimeter("Calorimeter",level, 0) ;
G4VPrimitiveScorer* Interaction= new InteractionScorers::PS_Interactions("Interaction",ms_InterCoord, 0) ;
+ G4VPrimitiveScorer* Process= new ProcessScorers::PS_Process("Process", 0) ;
//and register it to the multifunctionnal detector
m_ChiNuScorer->RegisterPrimitive(Calorimeter);
m_ChiNuScorer->RegisterPrimitive(Interaction);
+ m_ChiNuScorer->RegisterPrimitive(Process);
G4SDManager::GetSDMpointer()->AddNewDetector(m_ChiNuScorer) ;
}
diff --git a/NPSimulation/Scorers/ProcessScorers.cc b/NPSimulation/Scorers/ProcessScorers.cc
index 8d8c4fd7327ec106d196b5642e9d322ee009a9ef..5ff6d44dd9589a4f0ee7d72330b37db94fec50fd 100644
--- a/NPSimulation/Scorers/ProcessScorers.cc
+++ b/NPSimulation/Scorers/ProcessScorers.cc
@@ -29,6 +29,7 @@ PS_Process::PS_Process(G4String name, G4int depth) : G4VPrimitiveScorer(name, de
// m_NestingLevel = NestingLevel;
auto tree = RootOutput::getInstance()->GetTree();
tree->Branch("ProcScor_Particle", &t_particlename);
+ tree->Branch("ProcScor_Volume", &t_volumename);
tree->Branch("ProcScor_Process", &t_processname);
tree->Branch("ProcScor_Time", &t_processtime);
tree->Branch("ProcScor_Efficiency", &t_efficiency);
@@ -71,6 +72,7 @@ G4bool PS_Process::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
//step_length = aStep->GetTrack()->GetStepLength();
track_kineE = aStep->GetTrack()->GetKineticEnergy();
t_particlename.push_back(particlename);
+ t_volumename.push_back(volumename);
t_kinetic_energy.push_back(track_kineE);
//Looking after the neutrons
@@ -143,6 +145,7 @@ void PS_Process::EndOfEvent(G4HCofThisEvent*) {}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void PS_Process::clear() {
t_particlename.clear();
+ t_volumename.clear();
t_processname.clear();
t_processtime.clear();
t_gamma_energy.clear();
diff --git a/NPSimulation/Scorers/ProcessScorers.hh b/NPSimulation/Scorers/ProcessScorers.hh
index 52a54275570cfa3796ba05ff89e97869bd46d28b..e8768df6ee14f9f8862442a5d1a3f37b7fccdb55 100644
--- a/NPSimulation/Scorers/ProcessScorers.hh
+++ b/NPSimulation/Scorers/ProcessScorers.hh
@@ -59,6 +59,7 @@ namespace ProcessScorers {
//double t_Time;
//vector<unsigned int> t_Level;
vector<string> t_particlename;
+ vector<string> t_volumename;
vector<string> t_processname;
vector<double> t_processtime;
vector<double> t_gamma_energy;
@@ -74,6 +75,7 @@ namespace ProcessScorers {
inline unsigned int GetMult() {return t_processname.size();};
inline string GetProcessName(const unsigned int& i) {return t_processname[i];};
inline string GetParticleName(const unsigned int& i) {return t_particlename[i];};
+ inline string GetVolumeName(const unsigned int& i) {return t_volumename[i];};
inline double GetProcessTime(const unsigned int& i) {return t_processtime[i];};
inline vector<double> GetGammaEnergy() {return t_gamma_energy;};
inline vector<double> GetKinEnergy() {return t_kinetic_energy;};