diff --git a/NPLib/Detectors/PISTA/CMakeLists.txt b/NPLib/Detectors/PISTA/CMakeLists.txt
index 41848764389a46d338eae0aa8d60b5e782d863c2..b0b2938ab4cecab4d566f682e9c56311ac2999f0 100644
--- a/NPLib/Detectors/PISTA/CMakeLists.txt
+++ b/NPLib/Detectors/PISTA/CMakeLists.txt
@@ -21,5 +21,6 @@ add_library(NPPISTA SHARED TPISTASpectra.cxx TPISTAData.cxx TPISTAPhysics.cxx
target_link_libraries(NPPISTA ${ROOT_LIBRARIES} NPCore)
-install(FILES TPISTAData.h TPISTAPhysics.h TPISTASpectra.h TFPMWData.h TFPMWPhysics.h TICData.h TICPhysics.h TVamosReconstruction.h TTimeData.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
+install(FILES TPISTAData.h TPISTAPhysics.h TPISTASpectra.h TFPMWData.h TFPMWPhysics.h TICData.h TICPhysics.h
+ TVamosReconstruction.h TTimeData.h TProfileEvaluator.h DESTINATION ${CMAKE_INCLUDE_OUTPUT_DIRECTORY})
diff --git a/NPLib/Detectors/PISTA/TICPhysics.cxx b/NPLib/Detectors/PISTA/TICPhysics.cxx
index c905cf16a25a1716af8a7b21396c65017eebfacb..45a2a22eccf667d5a76117b816c51f51a7518e03 100644
--- a/NPLib/Detectors/PISTA/TICPhysics.cxx
+++ b/NPLib/Detectors/PISTA/TICPhysics.cxx
@@ -9,10 +9,10 @@
* Original Author: P. Morfouace contact address: pierre.morfouace@cea.fr *
* *
* Creation Date : October 2023 *
- * Last update : *
+ * Last update : 22/01/24 *
*---------------------------------------------------------------------------*
* Decription: *
- * This class hold IC Treated data *
+ * This class hold IC Treated data *
* *
*---------------------------------------------------------------------------*
* Comment: *
@@ -38,6 +38,7 @@ using namespace std;
// ROOT
#include "TChain.h"
+#include "TKey.h"
ClassImp(TICPhysics);
@@ -47,7 +48,6 @@ TICPhysics::TICPhysics()
m_PreTreatedData(new TICData),
m_EventPhysics(this),
m_FPMW_Section(-1),
- m_number_of_splines(34),
m_Eres_Threshold(3000),
m_Z_SPLINE_CORRECTION(false),
m_NumberOfDetectors(0){
@@ -94,11 +94,26 @@ void TICPhysics::BuildPhysicalEvent() {
double GainInit = Cal->GetValue("IC/INIT_SEG"+NPL::itoa(segment)+"_ALIGN",0);
fIC_raw[i] = m_PreTreatedData->GetIC_Charge(i);
+
fIC[i] = gain*m_PreTreatedData->GetIC_Charge(i);
+
+ if(i < 4){
+ fIC_PID[i] = 0.5* m_PreTreatedData->GetIC_Charge(i);
+ }
+ else if(i >= 6){
+ fIC_PID[i] = 2* m_PreTreatedData->GetIC_Charge(i);
+ }
+ else {
+ fIC_PID[i] = m_PreTreatedData->GetIC_Charge(i);
+ }
+
fIC_Init[i] = GainInit * m_PreTreatedData->GetIC_Charge(i);
fIC_TS.push_back(m_PreTreatedData->GetIC_TS(i));
}
+ if(m_Y_SPLINE_CORRECTION && m_XY0_SPLINE_CORRECTION){
+ ApplyXYCorrections();
+ }
if (fIC_Init[1]>0 && fIC_Init[5]>0) {
EtotInit = 0 ;
@@ -121,10 +136,26 @@ void TICPhysics::BuildPhysicalEvent() {
EtotInit = -100 ;
}
- if(fIC[1]>0 && fIC[5]>0){
- DE = 0.5*(fIC_raw[0] + fIC_raw[1] + fIC_raw[2] + fIC_raw[3]) + fIC_raw[4];
- Eres = fIC_raw[5] + fIC_raw[6] + fIC_raw[7] + fIC_raw[8] + fIC_raw[9];
+ DE = (fIC_PID[0] + fIC_PID[1] + fIC_PID[2] + fIC_PID[3]) + fIC_PID[4];
+ Eres = fIC_PID[5] + (fIC_PID[6] + fIC_PID[7] + fIC_PID[8] + fIC_PID[9]);
+
+ if (m_DE_SPLINE_CORRECTION && m_Y_SPLINE_CORRECTION){
+ if (m_TimeData->GetMWPC13Mult() ==1 && fIC_TS.size()>=8 ){ //only mult 1 event
+ UShort_t FPMW_Section = m_FPMW_Section;
+ double TempY;
+
+ //Data year sensitive loading
+ if (m_Data_Year == 2024){
+ TempY = 10* (fIC_TS.at(0) - m_TimeData->GetTS_MWPC13(0)) - ((m_TimeData->GetTime_MWPC13(0)+m_TimeData->GetToff_DT13(FPMW_Section))) ;
+ }
+ else if (m_Data_Year == 2023){
+ TempY = m_Y ;
+ }
+ DE = DE * m_DEspline.at(0)->Eval(0) / m_DEspline.at(0)->Eval(TempY) ;
+ } // end if mult 1
+ } // end DE correction
+ if(fIC[1]>0 && fIC[5]>0){
double scalor = Cal->GetValue("IC/ETOT_SCALING_SEC"+NPL::itoa(m_FPMW_Section),0);
for(int i=0; i<10; i++){
@@ -232,17 +263,39 @@ void TICPhysics::ReadAnalysisConfig() {
m_Data_Year = atof(DataBuffer.c_str());
cout << whatToDo << " " << m_Data_Year << endl;
}
+ else if (whatToDo=="LOAD_Y_SPLINE") {
+ AnalysisConfigFile >> DataBuffer;
+ m_Y_SPLINE_PATH = DataBuffer;
+ cout << "*** Loading Y spline ***" << endl;
+ m_Y_SPLINE_CORRECTION = LoadSpline(m_Yspline,m_number_Y_spline,m_Y_SPLINE_PATH);
+ }
+ else if (whatToDo=="LOAD_XY0_PROFILE") {
+ AnalysisConfigFile >> DataBuffer;
+ m_XY0_PROFILE_PATH = DataBuffer;
+ cout << "*** Loading XY0 profile ***" << endl;
+ TString PathTemp = m_XY0_PROFILE_PATH;
+ TFile* ifile = new TFile(PathTemp,"read");
+ if(ifile->IsOpen() && !ifile->IsZombie()){
+ m_XY0_SPLINE_CORRECTION = true;
+ m_IC0_Profile.LoadProfile(PathTemp,"ICOneZeroProfile");
+ }
+ else {
+ m_XY0_SPLINE_CORRECTION = false ;
+ }
+ ifile->Close();
+ }
- else if(whatToDo=="NUMBER_OF_SPLINES"){
+ else if (whatToDo=="LOAD_DE_SPLINE") {
AnalysisConfigFile >> DataBuffer;
- m_number_of_splines = atoi(DataBuffer.c_str());
- cout << whatToDo << " " << m_number_of_splines << endl;
+ m_DE_SPLINE_PATH = DataBuffer;
+ cout << "*** Loading DE spline ***" << endl;
+ m_DE_SPLINE_CORRECTION = LoadSpline(m_DEspline,m_number_DE_spline,m_DE_SPLINE_PATH);
}
else if (whatToDo=="LOAD_Z_SPLINE"){
AnalysisConfigFile >> DataBuffer;
m_Z_SPLINE_PATH = DataBuffer;
cout << "*** Loading Z spline ***" << endl;
- LoadZSpline();
+ m_Z_SPLINE_CORRECTION = LoadSpline(m_Zspline,m_number_zspline,m_Z_SPLINE_PATH);
}
else {
ReadingStatus = false;
@@ -253,35 +306,51 @@ void TICPhysics::ReadAnalysisConfig() {
///////////////////////////////////////////////////////////////////////////
-void TICPhysics::LoadZSpline(){
- TString filename = m_Z_SPLINE_PATH;
+bool TICPhysics::LoadSpline(vector<TSpline3*> &iSpline, int &NSpline , string Path){
+ TString filename = Path;
TFile* ifile = new TFile(filename,"read");
+ if(ifile->IsOpen() && !ifile->IsZombie()){
+
+ // Get number of spline
+ TIter next(ifile->GetListOfKeys());
+ TKey* key;
+ NSpline = 0 ;
- if(ifile->IsOpen()){
- m_Z_SPLINE_CORRECTION = true;
- for(int i=0; i<m_number_of_splines; i++){
- m_Zspline[i] = (TSpline3*) ifile->FindObjectAny(Form("fspline_%d",i+1));
- cout << "*** " << m_Zspline[i]->GetName() << " is loaded!" << endl;
+ while ((key=(TKey*)next())){
+ if (std::string(key->GetClassName()) == "TSpline3"){
+ NSpline ++;
+ }
+ }
+ cout << "This file contains " << NSpline << " splines " << endl;
+ // Load Spline
+ for(int i=0; i<NSpline; i++){
+ iSpline.at(i) = (TSpline3*) ifile->FindObjectAny(Form("fspline_%d",i+1));
+ iSpline.at(i)->SetName(Form("fspline_%s_%d",Path.c_str(),i+1));
+ cout << iSpline.at(i)->GetName() << " is loaded!" << endl;
}
- }
- else
- cout << "File " << filename << " not found!" << endl;
ifile->Close();
+ return true;
+ }
+ else{
+ cout << "File " << filename << " not found!" << endl;
+ ifile->Close();
+ return false;
+ }
}
///////////////////////////////////////////////////////////////////////////
double TICPhysics::ApplyZSpline(){
double DEcorr;
- double FF_DEcorr0[m_number_of_splines];
- for(int i=0; i<m_number_of_splines; i++){
+ double FF_DEcorr0[m_number_zspline];
+ for(int i=0; i<m_number_zspline; i++){
FF_DEcorr0[i] = m_Zspline[i]->Eval(8500);
}
double DEspline0;
double Eval_DEspline;
int index=0;
- for(int i=0; i<m_number_of_splines; i++){
+ for(int i=0; i<m_number_zspline; i++){
Eval_DEspline = m_Zspline[i]->Eval(Eres);
if(DE<Eval_DEspline) break;
index = i;
@@ -290,7 +359,7 @@ double TICPhysics::ApplyZSpline(){
Eval_DEspline = m_Zspline[index]->Eval(Eres);
DEspline0 = FF_DEcorr0[index];
double dmin, dsup;
- if(index<(m_number_of_splines-1) && DE>m_Zspline[0]->Eval(Eres)){
+ if(index<(m_number_zspline-1) && DE>m_Zspline[0]->Eval(Eres)){
dmin = DE - m_Zspline[index]->Eval(Eres);
dsup = m_Zspline[index+1]->Eval(Eres) - DE;
@@ -299,7 +368,7 @@ double TICPhysics::ApplyZSpline(){
DEcorr = DEspline0 * DE / Eval_DEspline;
}
- else if(index==m_number_of_splines-1){
+ else if(index==m_number_zspline-1){
Eval_DEspline = m_Zspline[index]->Eval(Eres);
DEspline0 = FF_DEcorr0[index];
@@ -308,6 +377,66 @@ double TICPhysics::ApplyZSpline(){
return DEcorr;
}
+///////////////////////////////////////////////////////////////////////////
+void TICPhysics::ApplyXYCorrections(){
+ vector<double> ICcorr_Y(11), ICcorr_X(11);
+ double FF_DriftTime ;
+ if (m_TimeData->GetMWPC13Mult() ==1 && fIC_TS.size()>=8){ //only mult 1 event
+ UShort_t FPMW_Section = m_FPMW_Section;
+
+ // ***************************Different def of correction depending on year**********************************
+ if (m_Data_Year == 2024){
+ FF_DriftTime = 10* (fIC_TS.at(0) - m_TimeData->GetTS_MWPC13(0)) - ((m_TimeData->GetTime_MWPC13(0)+m_TimeData->GetToff_DT13(FPMW_Section))) ;
+ }
+ else if (m_Data_Year == 2023){
+ FF_DriftTime = m_Y;
+ }
+ else {
+ return ;
+ }
+
+ //************************** Correction of section 1 to 4 ***************************************************
+ for (int seg = 2; seg < fIC_TS.size() ; seg++) { // loop on multiplicity of event
+
+ if (m_Yspline.at(seg-2)==0){
+ ICcorr_Y.at(seg) = fIC_PID[seg];
+ }
+
+ else {
+ ICcorr_Y.at(seg) = fIC_PID[seg] * m_Yspline.at(seg-2)->Eval(0)/ m_Yspline.at(seg-2)->Eval(FF_DriftTime);
+ if (!(ICcorr_Y.at(seg)==ICcorr_Y.at(seg))) ICcorr_Y.at(seg) = 0;
+ } //end if non empty
+ if (seg == 0) break;
+ }//endloop seg
+
+
+ //************************** Correction of section 0 ***************************************************
+
+ Double_t PolX = 1.37622;
+
+ Double_t ICRatio = fIC_PID[1]/fIC_PID[0];
+ Double_t ICRatioCorr = ICRatio * PolX / m_IC0_Profile.Evaluate(m_X,FF_DriftTime,false);
+ Double_t ICcorr_Y0 = fIC_PID[0] / PolX * m_IC0_Profile.Evaluate(m_X,FF_DriftTime,false);
+
+ if ( ICRatioCorr<1.4 || ICRatioCorr >1.5){
+ ICRatioCorr = ICRatio * PolX / m_IC0_Profile.Evaluate(m_X,FF_DriftTime,false);
+ ICcorr_Y0 = fIC_PID[0] / PolX * m_IC0_Profile.Evaluate(m_X,FF_DriftTime,false);
+ if (ICRatioCorr >100) {
+ ICcorr_Y0 = fIC_PID[0];
+ }
+ }
+
+
+ //************************** Overwrite ICRAW ***************************************************
+ //Overwrite ic_raw
+ for (int i = 1 ; i<fIC_TS.size() ;i++){
+ if (ICcorr_Y.at(i) != 0 ){
+ fIC_PID[i] = ICcorr_Y.at(i);
+ }
+ }
+ fIC_PID[0] = ICcorr_Y0;
+ }
+}
///////////////////////////////////////////////////////////////////////////
void TICPhysics::Clear() {
@@ -320,6 +449,7 @@ void TICPhysics::Clear() {
fIC[i] = 0;
fIC_raw[i] = 0;
fIC_Init[i] = 0;
+ fIC_PID[i] = 0;
}
fIC_TS.clear();
diff --git a/NPLib/Detectors/PISTA/TICPhysics.h b/NPLib/Detectors/PISTA/TICPhysics.h
index c02c89ab1ccfbd4a6c7ef8fee2935fb7f8f7405b..2656c0361e46b4c25709f00699ebe351c99ea3ed 100644
--- a/NPLib/Detectors/PISTA/TICPhysics.h
+++ b/NPLib/Detectors/PISTA/TICPhysics.h
@@ -34,8 +34,11 @@ using namespace std;
#include "TH1.h"
#include "TVector3.h"
#include "TSpline.h"
+#include "TFile.h"
// NPTool headers
#include "TICData.h"
+#include "TTimeData.h"
+#include "TProfileEvaluator.h"
#include "NPCalibrationManager.h"
#include "NPVDetector.h"
#include "NPInputParser.h"
@@ -71,6 +74,7 @@ class TICPhysics : public TObject, public NPL::VDetector {
double fIC[11];
double fIC_raw[11];
+ double fIC_PID[11];
double fIC_Init[11];//!
private:
@@ -132,10 +136,22 @@ class TICPhysics : public TObject, public NPL::VDetector {
// needed for online analysis for example
void SetRawDataPointer(TICData* rawDataPointer) {m_EventData = rawDataPointer;}
+ // Setter and getter for section of the FPMW.
void SetFPMWSection(int section) {m_FPMW_Section = section;}
int GetFPMWSection() {return m_FPMW_Section;}
- void LoadZSpline();
+
+ //Setter and Getter for TTimeData used for drift time calculation and XY
+ void SetTTimeData(TTimeData *newTime) {m_TimeData = newTime;}
+ TTimeData* GetTTimeData() {return m_TimeData;}
+
+ void SetX(double PosX) {m_X = PosX;}
+ double GetX() {return m_X;}
+ void SetY(double PosY) {m_Y = PosY;}
+ double GetY() {return m_Y;}
+
+ bool LoadSpline(vector<TSpline3*> &iSpline, int &NSpline, string Path);
double ApplyZSpline();
+ void ApplyXYCorrections();
// objects are not written in the TTree
private:
@@ -152,14 +168,36 @@ class TICPhysics : public TObject, public NPL::VDetector {
private:
int m_NumberOfDetectors; //!
int m_FPMW_Section; //!
- int m_number_of_splines; //!
+ int m_number_zspline; //!
string m_Z_SPLINE_PATH; //!
bool m_Z_SPLINE_CORRECTION; //!
- TSpline3* m_Zspline[50]; //!
+ vector<TSpline3*> m_Zspline = vector<TSpline3*>(50); //!
double m_Eres_Threshold; //!
- //Year of acquisition
+ //Time temporary variable
+ TTimeData* m_TimeData; //!
+ //Correction in XY of IC
+ string m_Y_SPLINE_PATH; //!
+ string m_XY0_PROFILE_PATH; //!
+ string m_DE_SPLINE_PATH; //!
+
+ bool m_Y_SPLINE_CORRECTION = false; //!
+ bool m_DE_SPLINE_CORRECTION = false; //!
+ bool m_XY0_SPLINE_CORRECTION = false; //!
+
+ int m_number_Y_spline; //!
+ int m_number_XY0_spline; //!
+ int m_number_DE_spline; //!
+
+ double m_X; //!
+ double m_Y; //!
+ // Be careful this array contains only the correction from IC1 to 4
+ vector<TSpline3*> m_Yspline = vector<TSpline3*>(11); //!
+ vector<TSpline3*> m_DEspline = vector<TSpline3*>(11); //!
+ // The following contains correction for IC0
+ ProfileEvaluator m_IC0_Profile;//!
private:
+ //Year of acquisition
double m_Data_Year; //!
// Static constructor to be passed to the Detector Factory
public:
diff --git a/NPLib/Detectors/PISTA/TTimeData.cxx b/NPLib/Detectors/PISTA/TTimeData.cxx
index cb8713f660e51856e867ed046a518bcc9d047484..b3dab37f54621e72a2069e2ada357bd6f07c3aab 100644
--- a/NPLib/Detectors/PISTA/TTimeData.cxx
+++ b/NPLib/Detectors/PISTA/TTimeData.cxx
@@ -54,6 +54,9 @@ void TTimeData::Clear() {
fTime_MWPC23.clear();
fTime_MWPC24.clear();
+ fToff_DT13.clear();
+ fToff_DT14.clear();
+
fSection_MWPC3.clear();
fSection_MWPC4.clear();
diff --git a/NPLib/Detectors/PISTA/TTimeData.h b/NPLib/Detectors/PISTA/TTimeData.h
index d56982141903c3cad350f77b0d3193f5f893f3dd..5a6dad3c73408a5fca8bda46e541014a0b6b1f5e 100644
--- a/NPLib/Detectors/PISTA/TTimeData.h
+++ b/NPLib/Detectors/PISTA/TTimeData.h
@@ -44,6 +44,9 @@ class TTimeData : public TObject {
vector<float> fTime_MWPC23;
vector<float> fTime_MWPC24;
+ vector<double> fToff_DT13;
+ vector<double> fToff_DT14;
+
vector<short> fSection_MWPC3;
vector<short> fSection_MWPC4;
@@ -82,6 +85,9 @@ class TTimeData : public TObject {
inline void SetTime_MWPC23(float time ){fTime_MWPC23.push_back(time);};//!
inline void SetTime_MWPC24(float time ){fTime_MWPC24.push_back(time);};//!
+ inline void SetToff_DT13(double offset ){fToff_DT13.push_back(offset);};//!
+ inline void SetToff_DT14(double offset ){fToff_DT14.push_back(offset);};//!
+
inline void SetSection_MWPC3(short section ){fSection_MWPC3.push_back(section);};//!
inline void SetSection_MWPC4(short section ){fSection_MWPC4.push_back(section);};//!
@@ -107,6 +113,11 @@ class TTimeData : public TObject {
inline float GetTime_MWPC24(const unsigned int &i) const
{return fTime_MWPC24.at(i) ;}//!
+ inline double GetToff_DT13(const unsigned int &i) const
+ {return fToff_DT13.at(i) ;}//!
+ inline double GetToff_DT14(const unsigned int &i) const
+ {return fToff_DT14.at(i) ;}//!
+
inline short GetSection_MWPC3(const unsigned int &i) const
{return fSection_MWPC3.at(i) ;}//!
inline short GetSection_MWPC4(const unsigned int &i) const
@@ -123,7 +134,7 @@ class TTimeData : public TObject {
//////////////////////////////////////////////////////////////
// Required for ROOT dictionnary
- ClassDef(TTimeData,2) // TimeData structure
+ ClassDef(TTimeData,3) // TimeData structure
};
#endif