From eef51f28d0fd3383fe79d789322841a43b93c57b Mon Sep 17 00:00:00 2001
From: "audrey.chatillon" <audrey.chatillon@gmail.com>
Date: Fri, 24 Jan 2025 17:32:28 +0100
Subject: [PATCH] [Epic] ON-GOING need to have one input per TrackID for 1s
history data
---
NPLib/Detectors/Epic/TEpicData.cxx | 2 -
NPLib/Detectors/Epic/TEpicData.h | 28 +----
NPLib/Detectors/Epic/TEpicPhysics.cxx | 17 +--
NPSimulation/Detectors/Epic/Epic.cc | 105 ++++++++++--------
.../Scorers/GaseousDetectorScorers.cc | 24 ++--
.../Scorers/GaseousDetectorScorers.hh | 10 +-
6 files changed, 79 insertions(+), 107 deletions(-)
diff --git a/NPLib/Detectors/Epic/TEpicData.cxx b/NPLib/Detectors/Epic/TEpicData.cxx
index a97f178a8..955f35ec2 100644
--- a/NPLib/Detectors/Epic/TEpicData.cxx
+++ b/NPLib/Detectors/Epic/TEpicData.cxx
@@ -52,8 +52,6 @@ void TEpicData::Dump() const {
for (size_t i = 0 ; i < mysize ; i++){
cout << "AnodeNbr: " << GetAnodeNbr(i)
<< " Q1: " << GetQ1(i)
- << " Q2: " << GetQ2(i)
- << " Qmax: " << GetQmax(i)
<< " Time: " << GetTime(i) << endl;
}
}
diff --git a/NPLib/Detectors/Epic/TEpicData.h b/NPLib/Detectors/Epic/TEpicData.h
index ab3bc580d..4d612196b 100644
--- a/NPLib/Detectors/Epic/TEpicData.h
+++ b/NPLib/Detectors/Epic/TEpicData.h
@@ -36,15 +36,9 @@ class TEpicData : public TObject {
struct EpicAnodeData {
UShort_t AnodeNbr;
Double_t Q1;
- Double_t Q2;
- Double_t Qmax;
Double_t Time;
- Double_t Time_HF;
- Double_t ToF;
- Bool_t isFakeFission;
// vector of {Zstep, DEstep, DTstep}
std::vector<string> ParticleName;
- std::vector<int> ParentID;
std::vector<int> TrackID;
std::vector<Double_t> PosZ;
std::vector<Double_t> DE;
@@ -53,14 +47,8 @@ class TEpicData : public TObject {
// Getters
UShort_t GetAnodeNbr() const { return AnodeNbr; }
Double_t GetQ1() const { return Q1; }
- Double_t GetQ2() const { return Q2; }
- Double_t GetQmax() const { return Qmax; }
Double_t GetTime() const { return Time; }
- Double_t GetTimeHF() const { return Time_HF; }
- Double_t GetToF() const { return ToF; }
- Bool_t IsFakeFission() const { return isFakeFission; }
const std::vector<string>& GetParticleName() const { return ParticleName; }
- const std::vector<int>& GetParentID() const { return ParentID; }
const std::vector<int>& GetTrackID() const { return TrackID; }
const std::vector<Double_t>& GetPosZ() const { return PosZ; }
const std::vector<Double_t>& GetDE() const { return DE; }
@@ -95,11 +83,11 @@ class TEpicData : public TObject {
// add //! to avoid ROOT creating dictionnary for the methods
public:
////////////////////// SETTERS ////////////////////////
- void Set(const UShort_t& n, const Double_t& q1, const Double_t& q2, const Double_t& qmax,
- const Double_t& t, const Double_t& thf, const Double_t& tof, const Bool_t& fake,
- const std::vector<string>& name, const std::vector<int>& pid, const std::vector<int>& tid,
- const std::vector<double>& z, const std::vector<double>& de, const std::vector<double>& dt) {
- fEpic_Data.push_back({n, q1, q2, qmax, t, thf, tof, fake, name, pid, tid, z, de, dt});
+ void Set(const UShort_t& n, const Double_t& q1, const Double_t& t,
+ const std::vector<string>& name, const std::vector<int>& tid,
+ const std::vector<double>& z, const std::vector<double>& de,
+ const std::vector<double>& dt) {
+ fEpic_Data.push_back({n, q1, t, name, tid, z, de, dt});
}
const EpicAnodeData& operator[](const unsigned int& i) const {return fEpic_Data[i];}
@@ -107,14 +95,8 @@ class TEpicData : public TObject {
inline UShort_t GetMultiplicity() const {return fEpic_Data.size();}
UShort_t GetAnodeNbr(const unsigned short& i) const { return fEpic_Data[i].AnodeNbr; }
Double_t GetQ1(const unsigned int& i) const { return fEpic_Data[i].Q1; }
- Double_t GetQ2(const unsigned int& i) const { return fEpic_Data[i].Q2; }
- Double_t GetQmax(const unsigned int& i) const { return fEpic_Data[i].Qmax; }
Double_t GetTime(const unsigned int& i) const { return fEpic_Data[i].Time; }
- Double_t GetTimeHF(const unsigned int& i) const { return fEpic_Data[i].Time_HF; }
- Double_t GetToF(const unsigned int& i) const { return fEpic_Data[i].ToF; }
- Bool_t GetFakeFissionStatus(const unsigned int& i) const { return fEpic_Data[i].isFakeFission; }
const std::vector<string>& GetParticleName(const unsigned int& i) const { return fEpic_Data[i].ParticleName; }
- const std::vector<int>& GetParentID(const unsigned int& i) const { return fEpic_Data[i].ParentID; }
const std::vector<int>& GetTrackID(const unsigned int& i) const { return fEpic_Data[i].TrackID; }
const std::vector<Double_t>& GetPosZ(const unsigned int& i) const { return fEpic_Data[i].PosZ; }
const std::vector<Double_t>& GetDE(const unsigned int& i) const { return fEpic_Data[i].DE; }
diff --git a/NPLib/Detectors/Epic/TEpicPhysics.cxx b/NPLib/Detectors/Epic/TEpicPhysics.cxx
index 1d09de705..1a3015acd 100644
--- a/NPLib/Detectors/Epic/TEpicPhysics.cxx
+++ b/NPLib/Detectors/Epic/TEpicPhysics.cxx
@@ -88,13 +88,8 @@ void TEpicPhysics::BuildPhysicalEvent() {
int A = m_EventData->GetAnodeNbr(e);
AnodeNumber.push_back(m_PreTreatedData->GetAnodeNbr(e));
Q1.push_back(m_PreTreatedData->GetQ1(e));
- Q2.push_back(m_PreTreatedData->GetQ2(e));
- Qmax.push_back(m_PreTreatedData->GetQmax(e));
Time.push_back(m_PreTreatedData->GetTime(e));
- isFakeFission.push_back(m_PreTreatedData->GetFakeFissionStatus(e));
- Time_HF.push_back(m_PreTreatedData->GetTimeHF(e));
- ToF.push_back(m_PreTreatedData->GetToF(e));
}
}
@@ -118,18 +113,8 @@ void TEpicPhysics::PreTreat() {
double TimeOffset = 0;
TimeOffset = Cal->GetValue("Epic/ANODE"+NPL::itoa(AnodeNumber)+"_TIMEOFFSET",0);
double Time = m_EventData->GetTime(i);
- double TimeHF = m_EventData->GetTimeHF(i);
- double tof = Time - TimeHF - TimeOffset;
- if(tof < 0){
- tof += Cal->GetValue("Epic/PULSE_TIMEOFFSET",0) ;
- }
- m_PreTreatedData->Set(AnodeNumber, Q1,
- m_EventData->GetQ2(i),
- m_EventData->GetQmax(i),
- Time, TimeHF, tof,
- m_EventData->GetFakeFissionStatus(i),
+ m_PreTreatedData->Set(AnodeNumber, Q1, Time,
m_EventData->GetParticleName(i),
- m_EventData->GetParentID(i),
m_EventData->GetTrackID(i),
m_EventData->GetPosZ(i),
m_EventData->GetDE(i),
diff --git a/NPSimulation/Detectors/Epic/Epic.cc b/NPSimulation/Detectors/Epic/Epic.cc
index a87b151f3..539aa377a 100644
--- a/NPSimulation/Detectors/Epic/Epic.cc
+++ b/NPSimulation/Detectors/Epic/Epic.cc
@@ -360,6 +360,8 @@ G4AssemblyVolume* Epic::BuildEpic(){
double posZ_first_cathode = posZ_anode - 0.5*thickness_A*mm - m_Distance_AK*mm - 0.5*m_Thickness_K*mm;
// build the stack of cathodes / actinide samples / anodes to do build sensitive gas
+ //double step_length_limit = m_Distance_AK*mm / 100.;
+ //G4UserLimits * user_limits = new G4UserLimits(step_length_limit*mm, 33.*mm, 1.*us, 0, 0);
G4UserLimits * user_limits = new G4UserLimits(25.*um, 33.*mm, 1.*us, 0, 0);
for(int i=0; i<m_nA; i++){
@@ -535,60 +537,75 @@ void Epic::ReadSensitive(const G4Event* ){
// GaseousDetector scorer
GaseousDetectorScorers::PS_GaseousDetector* Scorer= (GaseousDetectorScorers::PS_GaseousDetector*) m_EpicScorer->GetPrimitive(0);
+ int previous_trackID = -1;
unsigned int size = Scorer->GetMult();
+ cout << endl << " ================= Epic::ReadSensitive() Scorer size = " << size << endl;
for(unsigned int i = 0 ; i < size ; i++){
vector<unsigned int> level = Scorer->GetLevel(i);
+ //double Time = RandGauss::shoot(Scorer->GetTime(i),Epic_NS::ResoTime);
//double Energy = RandGauss::shoot(Scorer->GetEnergy(i),Epic_NS::ResoEnergy);
double Energy = Scorer->GetEnergy(i); // Attention Energy of all particles (FF, alpha, e-, ...)
vector<string> step_name = Scorer->GetParticleName(i);
- vector<int> step_parentid = Scorer->GetParentID(i);
- vector<int> step_trackid = Scorer->GetTrackID(i);
+ vector<int> step_trackID = Scorer->GetTrackID(i);
vector<double> step_posZ = Scorer->GetStepPosZ(i);
vector<double> step_DE = Scorer->GetEnergyLossPerStep(i);
vector<double> step_time = Scorer->GetStepTime(i);
- double influence = 0 ;
- if(Energy>Epic_NS::EnergyThreshold){
- double Time = RandGauss::shoot(Scorer->GetTime(i),Epic_NS::ResoTime);
- int Anode = level[0];
-
- double thickness_A = 2*Epic_NS::Cu_Thickness*mm + Epic_NS::Kapton_Thickness*mm;
- double posZ_anode = (1-m_nA)*(m_Distance_AK*mm + m_Thickness_K*mm) - std::trunc(0.5*m_nA)*(m_InterDistance_KK*mm + thickness_A);
- posZ_anode += m_mapping_A[Anode] * (2. * (m_Distance_AK*mm + m_Thickness_K*mm) + m_InterDistance_KK*mm + thickness_A) ;
- vector<string> name;
- vector<int> parentid;
- vector<int> trackid;
- vector<double> z;
- vector<double> de;
- vector<double> dt;
-
- for(int j=0; j<step_name.size(); j++){
- //cout << "step #" << j << " , particle name = " << step_name.at(j) << " , z = " << step_posZ.at(j) << " , dE = " << step_DE.at(j) << " , t = " << step_time.at(j) << endl;
- if(step_name.at(j)=="e-" || step_name.at(j)=="anti_nu_e") {
- continue;
- }
- else{ // FF or alpha
- //cout << "[FF or alpha]step #" << j << " , particle name = " << step_name.at(j) << ", z = " << step_posZ.at(j) << " , dE = " << step_DE.at(j) << " , t = " << step_time.at(j) << endl;
- influence += step_DE.at(j) * TMath::Abs(step_posZ.at(j) - posZ_anode);
- name.push_back(step_name.at(j));
- parentid.push_back(step_parentid.at(j));
- trackid.push_back(step_trackid.at(j));
- z.push_back(step_posZ.at(j));
- de.push_back(step_DE.at(j));
- dt.push_back(step_time.at(j));
- }
- }
- m_Event->Set(m_mapping_A[Anode], // set anode number
- influence / m_Distance_AK*mm, // set Q1 FF or alpha only
- Energy, // set Q2
- Energy, // set Qmax
- Time, // set Time
- 0., // set TimeHF
- 0., // set ToF
- 0, // status FakeFission
- name,parentid,trackid,z, de, dt);
- }
- }
+ int Anode = level[0];
+ double thickness_A = 2*Epic_NS::Cu_Thickness*mm + Epic_NS::Kapton_Thickness*mm;
+ double posZ_anode = (1-m_nA)*(m_Distance_AK*mm + m_Thickness_K*mm) - std::trunc(0.5*m_nA)*(m_InterDistance_KK*mm + thickness_A);
+ posZ_anode += m_mapping_A[Anode] * (2. * (m_Distance_AK*mm + m_Thickness_K*mm) + m_InterDistance_KK*mm + thickness_A) ;
+
+ //// FIXME : need a fixed binning in dz of 25*um
+ //// FIXME : one input per TrackID -> change the TEpicData TrackID and name should not be a vector anymore
+ //vector<string> name;
+ //vector<int> trackID;
+ //vector<double> dz;
+ //vector<double> de;
+ //vector<double> dt;
+
+ //double influence_pertrackID = 0 ;
+ //bool end_of_new_trackID = false ;
+ //cout << "step_name.size() = " << step_name.size() << endl;
+ //for(int j=0; j<step_name.size(); j++){
+ // cout << "j = " << j << " : TrackID = " << step_trackID.at(j) << " [previous is " << previous_trackID << "]. Particle = " << step_name.at(j) << " Actual m_Event mult : " << m_Event->GetMultiplicity() << endl;
+ // if(step_name.at(j)!="e-" || step_name.at(j)=="anti_nu_e"){
+ // continue;
+ // }
+ // else if(step_name.at(j)=="alpha"){ // FF or alpha // to be fixed !
+ //
+ // cout << "j = " << j << " : TrackID = " << step_trackID.at(j) << " [previous is " << previous_trackID << "]. Particle = " << step_name.at(j) << " Actual m_Event mult : " << m_Event->GetMultiplicity() << endl;
+ //
+ // if(step_trackID.at(j) != previous_trackID){
+ // //if(m_Event->GetMultiplicity()>0) end_of_new_trackID = true;
+ // if(previous_trackID != -1) end_of_new_trackID = true;
+ // previous_trackID = step_trackID.at(j) ;
+ // }
+ //
+ // if(end_of_new_trackID == true) {
+ // m_Event->Set(m_mapping_A[Anode], // set anode number
+ // influence_pertrackID, // set Q1 FF or alpha only
+ // step_time.at(j), // set Time
+ // name,trackID, dz, de, dt);
+ // name.clear();
+ // trackID.clear();
+ // dz.clear();
+ // de.clear();
+ // dt.clear();
+ // end_of_new_trackID = false;
+ // }
+
+ // double dz_anode = TMath::Abs(step_posZ.at(j) - posZ_anode) - 0.5*thickness_A ; // need to remove half thickness of Anode
+ // influence_pertrackID += (step_DE.at(j) * dz_anode) / m_Distance_AK*mm;
+ //
+ // name.push_back(step_name.at(j));
+ // trackID.push_back(step_trackID.at(j));
+ // dz.push_back(dz_anode);
+ // de.push_back(step_DE.at(j));
+ // dt.push_back(step_time.at(j));
+ // }// end of else
+ //}// end of loop over j
+ }// end of loop over i
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
diff --git a/NPSimulation/Scorers/GaseousDetectorScorers.cc b/NPSimulation/Scorers/GaseousDetectorScorers.cc
index 7cf09f69d..cc1bc635a 100644
--- a/NPSimulation/Scorers/GaseousDetectorScorers.cc
+++ b/NPSimulation/Scorers/GaseousDetectorScorers.cc
@@ -48,10 +48,10 @@ vector<GaseousDetectorData>::iterator
GaseousDetectorDataVector::find(const unsigned int& index) {
for (vector<GaseousDetectorData>::iterator it = m_Data.begin();
it != m_Data.end(); it++) {
- if ((*it).GetIndex() == index)
- return it;
+ if ((*it).GetIndex() == index) // comparison of the current index with the argument index
+ return it; // return immediatly the iterator
}
- return m_Data.end();
+ return m_Data.end(); // return m_Data.end() if index is not found
}
//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
@@ -73,46 +73,45 @@ G4bool PS_GaseousDetector::ProcessHits(G4Step* aStep, G4TouchableHistory*) {
// Contain Particle Id, Energy, Time and as many copy number as nested volume
unsigned int mysize = m_NestingLevel.size();
string particlename = aStep->GetTrack()->GetParticleDefinition()->GetParticleName();
- int parentID = aStep->GetTrack()->GetParentID();
+ //if (particlename != "alpha") return TRUE;
+
int trackID = aStep->GetTrack()->GetTrackID();
double step_posZ = aStep->GetPreStepPoint()->GetPosition().z();
-
t_Energy = aStep->GetTotalEnergyDeposit();
t_Time = aStep->GetPreStepPoint()->GetGlobalTime(); // DeltaT [ns] since the begining of the simulated event up to the begining of the step
-
t_ParticleName.push_back(particlename);
- t_ParentID.push_back(parentID);
t_TrackID.push_back(trackID);
t_StepPosZ.push_back(step_posZ);
t_EnergyLossPerStep.push_back(t_Energy);
t_StepTime.push_back(t_Time);
t_Level.clear();
-
for (unsigned int i = 0; i < mysize; i++) {
t_Level.push_back(
aStep->GetPreStepPoint()->GetTouchableHandle()->GetCopyNumber(
m_NestingLevel[i]));
}
+
+
// Check if the particle has interact before, if yes, add up the energies and update the track
vector<GaseousDetectorData>::iterator it;
it = m_Data.find(GaseousDetectorData::CalculateIndex(t_Level));
if (it != m_Data.end()) {
- // update the m_Data vector at each step
+ // update the m_Data vector at each step corresponding at the same sensitive volume
it->Add(t_Energy);
it->SetParticleName(particlename);
- it->SetParentID(parentID);
it->SetTrackID(trackID);
it->SetStepPosZ(step_posZ);
it->SetEnergyLossPerStep(t_Energy);
it->SetStepTime(t_Time);
}
else {
- // first step: initialize a new m_Data vector
+ // first step: initialize a new m_Data vector<GaseousDetectorData>
// t_Time [ns] is the time at first step
- m_Data.Set(t_Energy, t_Time, t_Level, t_ParticleName, t_ParentID, t_TrackID, t_StepPosZ, t_EnergyLossPerStep, t_StepTime);
+ m_Data.Set(t_Energy, t_Time, t_Level, t_ParticleName, t_TrackID, t_StepPosZ, t_EnergyLossPerStep, t_StepTime);
}
+
return TRUE;
}
@@ -127,7 +126,6 @@ void PS_GaseousDetector::clear() {
m_Data.clear();
t_Level.clear();
t_ParticleName.clear();
- t_ParentID.clear();
t_TrackID.clear();
t_StepPosZ.clear();
t_EnergyLossPerStep.clear();
diff --git a/NPSimulation/Scorers/GaseousDetectorScorers.hh b/NPSimulation/Scorers/GaseousDetectorScorers.hh
index 21e813707..f403a4d92 100644
--- a/NPSimulation/Scorers/GaseousDetectorScorers.hh
+++ b/NPSimulation/Scorers/GaseousDetectorScorers.hh
@@ -39,7 +39,6 @@ namespace GaseousDetectorScorers {
const double& Time,
const vector<unsigned int>& Nesting,
const vector<string>& name,
- const vector<int>& pID,
const vector<int>& tID,
const vector<double>& posZ,
const vector<double>& DE,
@@ -49,7 +48,6 @@ namespace GaseousDetectorScorers {
m_Energy=Energy;
m_Time=Time;
m_ParticleName=name;
- m_ParentID=pID;
m_TrackID=tID;
m_StepPosZ=posZ;
m_EnergyLossPerStep=DE;
@@ -63,7 +61,6 @@ namespace GaseousDetectorScorers {
double m_Energy;
double m_Time;
vector<string> m_ParticleName;
- vector<int> m_ParentID;
vector<int> m_TrackID;
vector<double> m_StepPosZ;
vector<double> m_EnergyLossPerStep;
@@ -78,7 +75,6 @@ namespace GaseousDetectorScorers {
inline double GetEnergy() const {return m_Energy;}
inline double GetTime() const {return m_Time;}
inline vector<string> GetParticleName() const {return m_ParticleName;}
- inline vector<int> GetParentID() const {return m_ParentID;}
inline vector<int> GetTrackID() const {return m_TrackID;}
inline vector<double> GetStepPosZ() const {return m_StepPosZ;}
inline vector<double> GetEnergyLossPerStep() const {return m_EnergyLossPerStep;}
@@ -91,7 +87,6 @@ namespace GaseousDetectorScorers {
void Add(const double& DE) {m_Energy+=DE;};
// add specific information at each step
inline void SetParticleName(const string& name){m_ParticleName.push_back(name);}
- inline void SetParentID(const int& id){m_ParentID.push_back(id);}
inline void SetTrackID(const int& id){m_TrackID.push_back(id);}
inline void SetStepPosZ(const double& z){m_StepPosZ.push_back(z);}
inline void SetEnergyLossPerStep(const double& DE){m_EnergyLossPerStep.push_back(DE);}
@@ -120,12 +115,11 @@ namespace GaseousDetectorScorers {
const double& Time,
const vector<unsigned int>& Nesting,
const vector<string>& name,
- const vector<int>& pID,
const vector<int>& tID,
const vector<double>& z,
const vector<double>& DE,
const vector<double>& t) {
- m_Data.push_back(GaseousDetectorData(Energy,Time,Nesting,name,pID,tID,z,DE,t));
+ m_Data.push_back(GaseousDetectorData(Energy,Time,Nesting,name,tID,z,DE,t));
};
const GaseousDetectorData* operator[](const unsigned int& i) const {return &m_Data[i];};
};
@@ -161,7 +155,6 @@ namespace GaseousDetectorScorers {
double t_Time; // time at first step
// all these vectors are used to accumulate info at each step
vector<string> t_ParticleName;
- vector<int> t_ParentID;
vector<int> t_TrackID;
vector<double> t_StepPosZ;
vector<double> t_EnergyLossPerStep;
@@ -173,7 +166,6 @@ namespace GaseousDetectorScorers {
inline double GetEnergy(const unsigned int& i) {return m_Data[i]->GetEnergy();};
inline double GetTime(const unsigned int& i) {return m_Data[i]->GetTime();};
inline vector<string> GetParticleName(const unsigned int& i) const {return m_Data[i]->GetParticleName();}
- inline vector<int> GetParentID(const unsigned int& i) const {return m_Data[i]->GetParentID();}
inline vector<int> GetTrackID(const unsigned int& i) const {return m_Data[i]->GetTrackID();}
inline vector<double> GetStepPosZ(const unsigned int& i) const {return m_Data[i]->GetStepPosZ();}
inline vector<double> GetEnergyLossPerStep(const unsigned int& i) const {return m_Data[i]->GetEnergyLossPerStep();}
--
GitLab