diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/CorrectXY.C b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/CorrectXY.C
new file mode 100644
index 0000000000000000000000000000000000000000..68d3a24a0d5893da0f3e42a45fbcdae50171ca2f
--- /dev/null
+++ b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/CorrectXY.C
@@ -0,0 +1,163 @@
+#include <TCanvas.h>
+#include <TChain.h>
+#include "ProfileEvaluator.h"
+#include "TICPhysics.h"
+#include "TTimeData.h"
+#include <TStyle.h>
+#include <chrono>
+#include <TH3.h>
+#include <TH2.h>
+
+
+void EvaluateProfileFromFile(const char* filename = "input.root", const char* profname = "prof", Double_t x = 0, Double_t y = 0) ;
+
+vector<double> TxtToVector(const char *Path);
+
+
+void CorrectXY(){
+
+ ProfileEvaluator evaluator;
+ evaluator.LoadProfile("Output/RatioProfile.root","ICOneZeroProfile");
+
+ TChain* chain = new TChain("PhysicsTree");
+ chain->Add("../../../root/analysis/VamosCalib247.root");
+
+ TTimeData *Time = new TTimeData();
+ TICPhysics* IC = new TICPhysics() ;
+ double FF_IC_X, FF_IC_Y, FF_V13, FF_DriftTime;
+
+
+
+ vector<double> Toff13 , Toff14, Toff23, Toff24;
+ const char* Path13 = "../../mwpc/Toff/output/Toff13.txt";
+ const char* Path14 = "../../mwpc/Toff/output/Toff14.txt";
+ const char* Path23 = "../../mwpc/Toff/output/Toff23.txt";
+ const char* Path24 = "../../mwpc/Toff/output/Toff24.txt";
+
+ Toff13 = TxtToVector(Path13);
+ Toff14 = TxtToVector(Path14);
+ Toff23 = TxtToVector(Path23);
+ Toff24 = TxtToVector(Path24);
+
+
+ chain->SetBranchStatus("FF_IC_X", true);
+ chain->SetBranchAddress("FF_IC_X", &FF_IC_X);
+ chain->SetBranchStatus("FF_IC_Y", true);
+ chain->SetBranchAddress("FF_IC_Y", &FF_IC_Y);
+
+ chain->SetBranchStatus("IC", true);
+ chain->SetBranchAddress("IC", &IC);
+
+ chain->SetBranchStatus("Time", true);
+ chain->SetBranchAddress("Time", &Time);
+
+ TProfile2D *Pxyz = new TProfile2D("ICOneZeroProfile","ICOneZeroProfile",100,-1000.0,1000.0,100,-1200.0,1000.0);
+ TProfile2D *PxyzCorr = new TProfile2D("PxyzCorr","PxyzCorr",100,-1000.0,1000.0,100,-1200.0,1000.0);
+
+
+ //int Nentries = chain->GetEntries();
+ int Nentries = 1000000;
+ auto start = std::chrono::high_resolution_clock::now();
+
+ for (int e = 0; e < Nentries; e++) {
+ if (e % 100000 == 0 && e > 0 ) {
+ auto now = std::chrono::high_resolution_clock::now();
+ std::chrono::duration<double> elapsed = now - start;
+ double avgTimePerIteration = elapsed.count() / e;
+ double timeLeft = avgTimePerIteration * (Nentries - e);
+
+ std::cout << "********** Estimated time left: " << int(timeLeft) << " seconds **********" << "\r" << flush;
+ }
+
+ chain->GetEntry(e);
+
+ if(Time->GetMWPC13Mult()==1 && IC->fIC_raw[8] >0){
+ UShort_t FPMW_Section = Time->GetSection_MWPC3(0);
+ FF_DriftTime = 10* (IC->fIC_TS.at(0) - Time->GetTS_MWPC13(0)) - ((Time->GetTime_MWPC13(0)+Toff13.at(FPMW_Section))) ;
+
+
+ Double_t ICRatio = IC->fIC_raw[1]/IC->fIC_raw[0];
+
+ if(IC->fIC_raw[8] >0 && ICRatio <= 8 && ICRatio >0.1 ) {
+
+ Double_t ICRatioCorr = ICRatio * 1.3 / evaluator.Evaluate(FF_IC_X,FF_DriftTime,true);
+ if (ICRatioCorr >100 ) ICRatioCorr = 1.2;
+ Pxyz->Fill(FF_IC_X,FF_DriftTime,ICRatio);
+ PxyzCorr->Fill(FF_IC_X,FF_DriftTime,ICRatioCorr);
+ }
+ }
+
+ }
+
+ gStyle->SetCanvasPreferGL(true);
+ TCanvas *c1 = new TCanvas("c1","c1");
+ gStyle->SetPalette(kRainBow);
+ Pxyz->SetXTitle("X");
+ Pxyz->SetYTitle("DriftTime");
+ Pxyz->SetZTitle("RatioIC");
+ Pxyz->SetMinimum(1.0);
+ Pxyz->SetMaximum(1.5);
+ Pxyz->Draw("GL SURF2"); // Configure the 3D view
+ c1->Update();
+
+ TCanvas *c2 = new TCanvas("c2","c2");
+ gStyle->SetPalette(kRainBow);
+ PxyzCorr->SetXTitle("X");
+ PxyzCorr->SetYTitle("DriftTime");
+ PxyzCorr->SetZTitle("RatioIC");
+ PxyzCorr->SetMinimum(1.0);
+ PxyzCorr->SetMaximum(1.5);
+ PxyzCorr->Draw("GL SURF2"); // Configure the 3D view
+
+ TFile *ofile = new TFile("Output/RatioProfile.root" ,"recreate" );
+ Pxyz ->Write();
+
+}
+
+
+
+void EvaluateProfileFromFile(const char* filename ,
+ const char* profname ,
+ Double_t x ,
+ Double_t y ) {
+ ProfileEvaluator evaluator;
+ if (!evaluator.LoadProfile(filename, profname)) {
+ printf("Failed to load profile\n");
+ return;
+ }
+
+ evaluator.PrintInfo();
+ Double_t z = evaluator.Evaluate(x, y, true);
+ printf("z(%.2f, %.2f) = %.3f\n", x, y, z);
+
+ // Quick visualization
+ TCanvas* c = new TCanvas("c", "Profile View", 800, 400);
+ c->Divide(2,1);
+
+ c->cd(1);
+ evaluator.GetProfile()->Draw("colz");
+
+ c->cd(2);
+ evaluator.GetProfile()->Draw("surf1");
+}
+
+vector<double> TxtToVector(const char *Path){
+ string line;
+ vector<double> values;
+ ifstream file(Path);
+
+ if (file.is_open()) {
+ while (std::getline(file, line)) {
+ try {
+ values.push_back(std::stod(line));
+ } catch (const std::invalid_argument& e) {
+ std::cerr << "Invalid number in line: " << line << '\n';
+ }
+ }
+ file.close();
+ } else {
+ std::cerr << "Error opening file.\n";
+ }
+
+ return values;
+}
diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/DECorr.C b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/DECorr.C
new file mode 100644
index 0000000000000000000000000000000000000000..21101091a8d15f3b790ce6d385e3b33e95ab8b50
--- /dev/null
+++ b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/DECorr.C
@@ -0,0 +1,456 @@
+#include "TTimeData.h"
+#include <TICPhysics.h>
+#include <TCanvas.h>
+#include <TChain.h>
+#include <TF1.h>
+#include <TFile.h>
+#include <TH2.h>
+#include <TSpline.h>
+#include <TCutG.h>
+#include <TStyle.h>
+#include <fstream>
+#include <vector>
+#include "ProfileEvaluator.h"
+
+using namespace std;
+vector<double> TxtToVector(const char *Path);
+int GetNumberKey(TFile *infile ,const char* ClassName);
+
+TSpline3* MakeSpline(TH2F* hInput, Int_t NCallee, double XMin , double
+ XMax, double YMin, double YMax);
+//////////////////////////////////////////////////////////////////////////////////////////
+void DECorr(bool cut = false, bool spline = false) {
+ //===========================================================================================================
+ // Loading var
+ //===========================================================================================================
+ TChain* chain = new TChain("PhysicsTree");
+ chain->Add("../../../root/analysis/VamosCalib247.root");
+
+ TICPhysics* IC = new TICPhysics() ;
+ TTimeData *Time = new TTimeData();
+ double FF_IC_X, FF_IC_Y, FF_V13, Y;
+
+ chain->SetBranchStatus("FF_IC_X", true);
+ chain->SetBranchAddress("FF_IC_X", &FF_IC_X);
+
+ chain->SetBranchStatus("FF_IC_Y", true);
+ chain->SetBranchAddress("FF_IC_Y", &FF_IC_Y);
+
+ chain->SetBranchStatus("IC", true);
+ chain->SetBranchAddress("IC", &IC);
+
+ chain->SetBranchStatus("FF_V13", true);
+ chain->SetBranchAddress("FF_V13", &FF_V13);
+
+ chain->SetBranchStatus("Time", true);
+ chain->SetBranchAddress("Time", &Time);
+
+ vector<double> Toff13 , Toff14, Toff23, Toff24;
+ const char* Path13 = "../../mwpc/Toff/output/Toff13.txt";
+ const char* Path14 = "../../mwpc/Toff/output/Toff14.txt";
+ const char* Path23 = "../../mwpc/Toff/output/Toff23.txt";
+ const char* Path24 = "../../mwpc/Toff/output/Toff24.txt";
+
+ Toff13 = TxtToVector(Path13);
+ Toff14 = TxtToVector(Path14);
+ Toff23 = TxtToVector(Path23);
+ Toff24 = TxtToVector(Path24);
+
+ ProfileEvaluator Profile;
+ Profile.LoadProfile("Output/RatioProfile.root","ICOneZeroProfile");
+
+ TCutG *CutZ;
+ TCutG *CutZbis;
+ TFile *fCut;
+ if (cut) {
+ fCut = TFile::Open("Output/CutDeCorr.root");
+ CutZbis = (TCutG*)fCut->Get("CutZbis") ;
+ CutZ = (TCutG*)fCut->Get("CutZ") ;
+ }
+
+ TFile *fSpline;
+ TSpline3 *splineDEcorr, *splineDEbis;
+ fSpline= new TFile("Output/splineDE.root","open");
+ splineDEcorr = (TSpline3*) fSpline->Get("SplineDe_0");
+ splineDEbis = (TSpline3*) fSpline->Get("SplineDebis");
+
+ vector<TSpline3*> vsplineDE;
+ // Get number of spline
+ int SplineCountDe = 0 ;
+ TIter next1(fSpline->GetListOfKeys());
+ TKey* key1;
+
+ while ((key1=(TKey*)next1())){
+ if (std::string(key1->GetClassName()) == "TSpline3"){
+ SplineCountDe ++;
+ }
+ }
+
+ //Fill spline vector
+ for (int i =0 ; i<SplineCountDe ; i++){
+ vsplineDE.push_back((TSpline3*) fSpline->Get(Form("SplineDe_%d",i)));
+ }
+
+ //empty backup one
+ for (int i =0 ; i<SplineCountDe ; i++){
+ if (vsplineDE[i]==nullptr) vsplineDE.erase(vsplineDE.begin()+i);
+ }
+
+
+ //Defining output histo
+ TH2F* hDE_E = new TH2F("DE_E","DE_E",1000,0,22000,1000,0,28000);
+ TH2F* hDE_E_splined = new TH2F("DE_E_splined","DE_E_splined",1000,0,22000,1000,0,28000);
+ TH2F* hDE_Ebis = new TH2F("DE_Ebis","DE_Ebis",1000,0,22000,1000,0,24000);
+
+
+ TH2F* hDE_Y = new TH2F("DE_Y","DE_Y",1000,-1000,1000,1000,0,22000);
+ TH2F* hDE_Y_splined = new TH2F("DE_Ysplined","DE_Ysplined",1000,-1000,1000,1000,0,22000);
+ TH2F* hDE_Y_splined_nocut = new TH2F("DE_Ysplined_nocut","DE_Ysplined_nocut",1000,-1000,1000,1000,0,22000);
+ TH2F* hIC_Y = new TH2F("IC_Y","IC_Y",1000,-1000,1000,1000,0,22000);
+ TH2F* hIC_Y_splined = new TH2F("IC_Ysplined","IC_Ysplined",1000,-1000,1000,1000,0,22000);
+
+
+ TH2F* hDE_Y_corr = new TH2F("DE_Ycorr","DE_Ycorr",1000,-1000,1000,1000,0,22000);
+ TH2F* hDE_Y_corr_nocut = new TH2F("DE_Ycorr_nocut","DE_Ycorr_nocut",1000,-1000,1000,1000,0,22000);
+
+ TH2F* hDE_Y_bis = new TH2F("DE_Ybis","DE_Ybis",1000,-1000,1000,1000,0,22000);
+ TH2F* hDEbis_Y_corr = new TH2F("DE_Ybiscorr","DE_Ybiscorr",1000,-1000,1000,1000,0,22000);
+
+ TH2F* hDE_E_corr = new TH2F("DE_Ecorr","DE_Ecorr",1000,0,22000,1000,0,28000);
+ TH2F* hDE_Ebis_corr = new TH2F("DE_Ebiscorr","DE_Ebiscorr",1000,0,22000,1000,0,24000);
+
+ TH2F* hDE_V = new TH2F("DE_V","DE_V",1000,0,5,1000,0,28000);
+
+ // Def Lim spline
+ double Ymin = -1000 , Ymax =1000;
+ //===========================================================================================================
+ // Event Loop
+ //===========================================================================================================
+ int Nentries = chain->GetEntries();
+ //int Nentries = 1000000;
+ auto start = std::chrono::high_resolution_clock::now();
+
+ for (int e = 0; e < Nentries; e++) {
+ if (e % 100000 == 0 && e > 0 ) {
+ auto now = std::chrono::high_resolution_clock::now();
+ std::chrono::duration<double> elapsed = now - start;
+ double avgTimePerIteration = elapsed.count() / e;
+ double timeLeft = avgTimePerIteration * (Nentries - e);
+
+ std::cout << "********** Estimated time left: " << int(timeLeft) << " seconds **********" << "\r" << flush;
+ }
+
+ chain->GetEntry(e);
+
+
+ if (Time->GetMWPC13Mult() ==1 && IC->fIC_TS.size()>=8){ //only mult 1 event
+ UShort_t FPMW_Section = Time->GetSection_MWPC3(0);
+ double FF_DriftTime = 10* (IC->fIC_TS.at(0) - Time->GetTS_MWPC13(0)) - ((Time->GetTime_MWPC13(0)+Toff13.at(FPMW_Section))) ;
+
+ Double_t ICRatio = IC->fIC_raw[1]/IC->fIC_raw[0];
+ Double_t ICRatioCorr = ICRatio *1.3 / Profile.Evaluate(FF_IC_X,FF_DriftTime,true);
+ Double_t ICcorr_Y = IC->fIC_raw[0] / 1.3 * Profile.Evaluate(FF_IC_X,FF_DriftTime,true);
+
+
+ double DE = 0 , DE_splined=0 ,DE_corr=0, DE_Ybis=0;
+ double E =0;
+
+ for (int seg = 0 ; seg < sizeof(IC->fIC_raw)/sizeof(IC->fIC_raw[0]) ; seg++ ){
+ if (seg >4){
+ E += double(IC->fIC_raw[seg]) ;
+ }
+ }
+
+
+
+ // ********************* DE setter *********************
+ double IC1corr = (IC->fIC_raw[1]*(1-0.000686068*FF_IC_Y))*(1-4.88238e-05*FF_IC_Y+7.40395e-06*FF_IC_Y*FF_IC_Y);
+ double DE_Bis = 0.5*(IC1corr+IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4];
+
+ DE = 0.5*(double(IC->fIC_raw[0])+ double(IC->fIC_raw[1]) +double(IC->fIC_raw[2])+ double(IC->fIC_raw[3])) + double(IC->fIC_raw[4]);
+ DE_splined = 0.5*(ICcorr_Y+ IC->fIC_raw[1] +IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4];
+
+ //DE_splined = 0.5*( ICcorr_Y[1] +IC->fIC_raw[2]+IC->fIC_raw[3])+IC->fIC_raw[4];
+
+
+ //************* Cut ******************
+ bool CutV = abs( FF_V13-3.94)<0.1;
+ bool CutVbis = abs( FF_V13-3.94)<0.1;
+
+ //bool CutV = true;
+ //bool CutVbis = true;
+
+ //***********Fill histo*************
+
+ if (FF_DriftTime >Ymin && FF_DriftTime <Ymax ){
+
+ //*************Init DE_E*****************
+ hDE_E->Fill(E,DE);
+ hDE_E_splined->Fill(E,DE_splined);
+ hDE_Ebis->Fill(E,DE_Bis);
+
+ hDE_V->Fill(FF_V13,DE_splined);
+
+ //******* Coor DE ************
+ if (spline== true ){
+ DE_corr = DE_splined * splineDEcorr->Eval(0) / splineDEcorr->Eval(FF_DriftTime);
+ hDE_E_corr->Fill(E,DE_corr);
+
+ DE_Ybis = DE_Bis * splineDEbis->Eval(0) / splineDEbis->Eval(FF_DriftTime);
+ hDE_Ebis_corr->Fill(E,DE_Ybis);
+
+ }
+
+ else {
+ DE_corr = DE_splined;
+ hDE_E_corr->Fill(E,DE_corr);
+
+ DE_Ybis = DE_Bis ;
+ hDE_Ebis_corr->Fill(E,DE_Ybis);
+ }
+
+ bool CutCharge ;
+
+
+ if (cut == true) CutCharge = CutZ->IsInside(FF_V13,DE_splined);
+ //bool CutChargeBis = CutZ->IsInside(E,DE_Y_Recur);
+ //bool CutVbis = CutZbis->IsInside(E,DE_Bis);
+ // DE_Y in a cut
+
+
+ if (cut == true && CutV) {
+ if ( CutV) {
+ hDE_Y_corr_nocut->Fill(FF_DriftTime,DE_corr);
+ hDE_Y_splined_nocut->Fill(FF_DriftTime,DE_splined);
+ hDE_Y->Fill(FF_DriftTime,DE);
+ hDE_Y_splined->Fill(FF_DriftTime,DE_splined);
+ hDE_Y_corr->Fill(FF_DriftTime,DE_corr);
+
+ hIC_Y->Fill(FF_DriftTime,IC->fIC_raw[0]);
+ hIC_Y_splined->Fill(FF_DriftTime,ICcorr_Y);
+ }
+ }//end cut
+ else {
+ hDE_Y->Fill(FF_DriftTime,DE);
+ hDE_Y_splined->Fill(FF_DriftTime,DE_splined);
+ hDE_Y_corr->Fill(FF_DriftTime,DE_corr);
+
+ hIC_Y->Fill(FF_DriftTime,IC->fIC_raw[0]);
+ hIC_Y_splined->Fill(FF_DriftTime,ICcorr_Y);
+ }
+
+ if (cut == true && CutVbis) {
+ hDE_Y_bis->Fill(Y,DE_Bis);
+ hDEbis_Y_corr->Fill(Y,DE_Ybis);
+ }
+ else{;}
+ //end cut
+ }//end cond Y
+ } // End mult 1
+ } //endl loop event
+
+ // making spline
+
+ TSpline3 *splineDE, *osplineDEbis;
+ if (spline == true){
+ splineDE = MakeSpline(hDE_Y_splined,0,Ymin,Ymax,0,23000);
+ osplineDEbis = MakeSpline(hDE_Y_bis,1,Ymin,Ymax,0,23000);
+ TFile *fSpline = new TFile("Output/splineDE.root","recreate");
+ splineDE->SetName("SplineDe_0");
+ splineDE->Write();
+ osplineDEbis->SetName("SplineDebis");
+ osplineDEbis->Write();
+ fSpline->Close();
+ }
+ //===========================================================================================================
+ // Drawing histo
+ //===========================================================================================================
+
+ gStyle->SetPalette(kRainBow);
+
+ TCanvas* c = new TCanvas("c","c",1500,1000);
+ c->Divide(2,2);
+ c->cd(1);
+ gPad->SetLogz();
+ hDE_Y->ProfileX()->Draw();
+ c->cd(2);
+ gPad->SetLogz();
+ hDE_Y_splined->ProfileX()->Draw();
+ c->cd(3);
+ gPad->SetLogz();
+ hIC_Y->ProfileX()->Draw();
+ c->cd(4);
+ gPad->SetLogz();
+ hIC_Y_splined->ProfileX()->Draw();
+
+ TCanvas *c4 = new TCanvas("c4","c4",1500,1000);
+ c4->Divide(2);
+ c4->cd(1);
+ hDE_Y_splined_nocut->Draw();
+ c4->cd(2);
+ hDE_Y_corr_nocut->Draw();
+
+ TCanvas *c5 = new TCanvas("c5","c5",1500,1000);
+ hDE_E->Draw("colz");
+
+
+
+ TCanvas *c1 = new TCanvas("c1","c1",1500,1000);
+ hDE_V->Draw();
+
+ TCanvas* c2 = new TCanvas("c2","c2",1500,1000);
+ c2->Divide(3 + 1*spline);
+
+ c2->cd(1);
+ hDE_E->SetTitle("Unmodified");
+ gPad->Modified();
+ gPad->Update();
+ hDE_E->Draw("colz");
+ hDE_E->GetXaxis()->SetTitle("E");
+ hDE_E->GetYaxis()->SetTitle("DE");
+
+ c2->cd(2);
+ hDE_Ebis->SetTitle("Online");
+ gPad->Modified();
+ gPad->Update();
+ hDE_Ebis->Draw("colz");
+ hDE_Ebis->GetXaxis()->SetTitle("E");
+ hDE_Ebis->GetYaxis()->SetTitle("DE");
+
+ c2->cd(3);
+ hDE_E_splined->SetTitle("IC Corrected in DT");
+ gPad->Modified();
+ gPad->Update();
+ hDE_E_splined->Draw("colz");
+ hDE_E_splined->GetXaxis()->SetTitle("E");
+ hDE_E_splined->GetYaxis()->SetTitle("DE");
+
+ if ( spline == true){
+ c2->cd(4);
+ hDE_E_corr->SetTitle("DE Corrected");
+ gPad->Modified();
+ gPad->Update();
+ hDE_E_corr->Draw("colz");
+ hDE_E_corr->GetXaxis()->SetTitle("E");
+ hDE_E_corr->GetYaxis()->SetTitle("DE");
+ }
+} // End spline chio XY
+
+
+int GetNumberKey(TFile *infile ,const char* ClassName){
+ // Get number of spline
+ int KeyCount = 0 ;
+ TIter next(infile->GetListOfKeys());
+ TKey* keyFit;
+
+ while ((keyFit=(TKey*)next())){
+ if (std::string(keyFit->GetClassName()) == ClassName){
+ KeyCount ++;
+ }
+ }
+
+ return KeyCount;
+
+}
+
+TSpline3* MakeSpline(TH2F* hInput, Int_t NCallee, double XMin = 0 , double
+ XMax = 1000 , double YMin=0 , double YMax = 100) {
+
+ TSpline3* gspline;
+
+ TH1F* hProfile;
+ TH1F* pfx; // extended hProfile
+
+ TCanvas* canfit = new TCanvas(Form("CanSpline_%02d",NCallee),
+ Form("canfit_%02d",NCallee), 0, 0, 1000, 500);
+
+
+ //===========================================================================================================
+ // Init profile
+ //===========================================================================================================
+
+ // Create the TProfile
+ hProfile = (TH1F*)hInput->ProfileX();
+ hProfile->SetLineWidth(2);
+ hProfile->SetDirectory(0);
+ canfit->cd();
+
+ hProfile->GetYaxis()->SetRangeUser(YMin, YMax);
+ hProfile->Draw();
+ //===========================================================================================================
+ // First and last bin to get to 6k
+ // event get promoted
+ //===========================================================================================================
+ int FirstBin, LastBin;
+ double Treshold = 0;
+ for (int bin =1 ; bin<hProfile->GetNbinsX(); bin++) {
+ FirstBin = bin;
+ if (hProfile->GetBinLowEdge(bin)> XMin) break;
+ }
+ for (int bin = hProfile->GetNbinsX(); bin>1 ; bin--) {
+ LastBin = bin;
+ if (hProfile->GetBinLowEdge(bin)< XMax) break;
+ }
+ //===========================================================================================================
+ // Init Extended profile function
+ //===========================================================================================================
+ // Create the extended TProfile
+ pfx = new TH1F(
+ Form("pfx_%02d", Int_t(NCallee/2)), Form("pfx_%02d", Int_t(NCallee/2)), hProfile->GetNbinsX(), hProfile->GetBinLowEdge(1),
+ hProfile->GetBinLowEdge(hProfile->GetNbinsX()) + hProfile->GetBinWidth(hProfile->GetNbinsX()));
+ pfx->SetLineColor(8);
+ pfx->SetDirectory(0);
+ //===========================================================================================================
+ // Fill extended profile
+ //===========================================================================================================
+
+ // fill the extended TProfile
+ float newval, lastval, lasterr;
+ for (int bin = 1; bin <= FirstBin; bin++) {
+ newval = 0;
+ pfx->SetBinContent(bin, newval);
+ }
+
+
+ for (int bin = FirstBin; bin <= LastBin; bin++) {
+ newval = hProfile->GetBinContent(bin);
+ if (newval != 0) {
+ pfx->SetBinContent(bin, newval);
+ pfx->SetBinError(bin, hProfile->GetBinError(bin));
+ lastval = newval;
+ lasterr = hProfile->GetBinError(bin);
+ }
+ else {
+ pfx->SetBinContent(bin, lastval);
+ pfx->SetBinError(bin, lasterr);
+ }
+ }
+ pfx->Draw("same");
+
+ gspline = new TSpline3(pfx);
+ gspline->SetName(Form("fspline_%d", NCallee + 1));
+
+ return gspline;
+
+} // end makespline
+
+vector<double> TxtToVector(const char *Path){
+ string line;
+ vector<double> values;
+ ifstream file(Path);
+
+ if (file.is_open()) {
+ while (std::getline(file, line)) {
+ try {
+ values.push_back(std::stod(line));
+ } catch (const std::invalid_argument& e) {
+ std::cerr << "Invalid number in line: " << line << '\n';
+ }
+ }
+ file.close();
+ } else {
+ std::cerr << "Error opening file.\n";
+ }
+
+ return values;
+}
diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/RatioProfile.root b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/RatioProfile.root
new file mode 100644
index 0000000000000000000000000000000000000000..e181a9879429be4a95e61fd1d388837306f078a2
Binary files /dev/null and b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/RatioProfile.root differ
diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/ReadMe.md b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/ReadMe.md
new file mode 100644
index 0000000000000000000000000000000000000000..ae505d9092a2a100ef1a1e6f5b3368d68e9fb6f8
--- /dev/null
+++ b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/ReadMe.md
@@ -0,0 +1 @@
+# Put here the output profile
diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/splineDE.root b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/splineDE.root
new file mode 100644
index 0000000000000000000000000000000000000000..330cfde1c7a3895b11fcdc11f57c399d884723ca
Binary files /dev/null and b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/Output/splineDE.root differ
diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/ProfileEvaluator.h b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/ProfileEvaluator.h
new file mode 100644
index 0000000000000000000000000000000000000000..0994ba971a8eb0dc2bb445420355b322baf2af80
--- /dev/null
+++ b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/ProfileEvaluator.h
@@ -0,0 +1,127 @@
+#ifndef PROFILE_EVAL
+#define PROFILE_EVAL
+
+#include <TFile.h>
+#include <TProfile2D.h>
+#include <iostream>
+using namespace std;
+
+class ProfileEvaluator {
+ private:
+ TFile* fFile;
+ TProfile2D* fProfile;
+
+ public:
+
+ //***********************************************************//
+
+ ProfileEvaluator() : fFile(nullptr), fProfile(nullptr) {}
+
+ //***********************************************************//
+
+ bool LoadProfile(const char* filename, const char* profname) {
+ fFile = TFile::Open(filename, "READ");
+ if (!fFile || fFile->IsZombie()) return false;
+
+ fProfile = (TProfile2D*)fFile->Get(profname)->Clone("Profile");
+ fProfile->SetDirectory(0);
+ return fProfile != nullptr;
+ }
+
+ //***********************************************************//
+
+ Double_t Evaluate(Double_t x, Double_t y, Bool_t useInterpolation = true) {
+ if (!fProfile) return 0.0;
+
+ Int_t binx = fProfile->GetXaxis()->FindBin(x);
+ Int_t biny = fProfile->GetYaxis()->FindBin(y);
+
+ if (!useInterpolation) return fProfile->GetBinContent(binx, biny);
+
+ //Def Bin dist
+ Double_t BinDistX = fProfile->GetXaxis()->GetBinCenter(2) - fProfile->GetXaxis()->GetBinCenter(1);
+ Double_t BinDistY = fProfile->GetYaxis()->GetBinCenter(2) - fProfile->GetYaxis()->GetBinCenter(1);
+
+ // Normalisation var
+ Double_t norma = 0;
+
+ // Fetch Bin Center for surrounding 3 by 3 bin and value
+ vector<Double_t> Cx(3) , Cy(3) ;
+ vector<vector<Double_t>> Q(3, vector<Double_t>(3)), W(3,vector<Double_t>(3));
+
+ for (int i = -1 ; i<=1 ; i++){
+ // Loading nearby bins
+ int TempX = binx +i ;
+ int PosX = i + 1 ;
+
+ // Edge cases handling
+ if (TempX < 1) TempX = 1 ;
+ if (TempX > fProfile->GetNbinsX()) TempX = fProfile->GetNbinsX() ;
+
+ // Loading center
+ Cx.at(PosX) = fProfile->GetXaxis()->GetBinCenter(TempX);
+ //Normat that we use POSX here
+
+ // Loading value and weight
+ for (int j = -1 ; j<=1 ; j++){
+ // Need to translate to acess the vector
+ int PosY = j + 1 ;
+ int TempY = biny + j ;
+
+ if (TempY < 1) TempY = 1 ;
+ if (TempY > fProfile->GetNbinsY()) TempY = fProfile->GetNbinsY() ;
+
+ Cy.at(PosY) = fProfile->GetYaxis()->GetBinCenter(TempY);
+
+ // Loading value
+ Q.at(PosX).at(PosY) = fProfile->GetBinContent(TempX,TempY);
+
+ // Calculating weight
+ Double_t Wx = 0 , Wy = 0 ;
+ Wx = 1 - ( abs( x - Cx.at(PosX)) / (BinDistX) );
+ Wy = 1 - ( abs( y - Cy.at(PosY)) / (BinDistY) );
+
+ //handling max range
+ if(Wx < 0 || Wy < 0) W.at(PosX).at(PosY) = 0;
+ else W.at(PosX).at(PosY) = Wx * Wy;
+ // Verifiy that sum of the weight is 1
+ norma += W.at(PosX).at(PosY) ;
+ }
+ }
+
+ // now just calculate the return value
+
+ Double_t InterPol = 0 ;
+
+ for (int PosX = 0 ; PosX <= 2 ; PosX ++ ){
+ for (int PosY = 0 ; PosY <= 2 ; PosY ++ ){
+ InterPol += Q.at(PosX).at(PosY) * W.at(PosX).at(PosY);
+ }
+ }
+ return InterPol;
+ }
+
+ //***********************************************************//
+
+ void PrintInfo() {
+ if (!fProfile) return;
+ printf("Profile: %s\n", fProfile->GetName());
+ printf("X bins: %d, range: [%.2f, %.2f]\n",
+ fProfile->GetNbinsX(),
+ fProfile->GetXaxis()->GetXmin(),
+ fProfile->GetXaxis()->GetXmax());
+ printf("Y bins: %d, range: [%.2f, %.2f]\n",
+ fProfile->GetNbinsY(),
+ fProfile->GetYaxis()->GetXmin(),
+ fProfile->GetYaxis()->GetXmax());
+ }
+
+ //***********************************************************//
+
+ TProfile2D* GetProfile() { return fProfile; }
+
+ //***********************************************************//
+ ~ProfileEvaluator() {
+ }
+};
+#endif
diff --git a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/VerifCorrelation.C b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/VerifCorrelation.C
index c386bbdcb379fd8c03b08ebfcb3c739389d12f48..b8fda7cd0ad7bd7aef63feefbaef0ff10fb1058e 100644
--- a/Projects/AlPhaPha/2024/macro/chio/XYCalibration/VerifCorrelation.C
+++ b/Projects/AlPhaPha/2024/macro/chio/XYCalibration/VerifCorrelation.C
@@ -1,9 +1,12 @@
#include "TICPhysics.h"
#include "TTimeData.h"
#include <TChain.h>
+#include <TFile.h>
#include <TH2.h>
#include <TH3.h>
-
+#include <TStyle.h>
+#include <TProfile2D.h>
+#include <chrono>
vector<double> TxtToVector(const char *Path);
void VerifCorrelation(){
@@ -20,10 +23,10 @@ void VerifCorrelation(){
vector<double> Toff13 , Toff14, Toff23, Toff24;
- const char* Path13 = "../../../mwpc/Toff/output/Toff13.txt";
- const char* Path14 = "../../../mwpc/Toff/output/Toff14.txt";
- const char* Path23 = "../../../mwpc/Toff/output/Toff23.txt";
- const char* Path24 = "../../../mwpc/Toff/output/Toff24.txt";
+ const char* Path13 = "../../mwpc/Toff/output/Toff13.txt";
+ const char* Path14 = "../../mwpc/Toff/output/Toff14.txt";
+ const char* Path23 = "../../mwpc/Toff/output/Toff23.txt";
+ const char* Path24 = "../../mwpc/Toff/output/Toff24.txt";
Toff13 = TxtToVector(Path13);
Toff14 = TxtToVector(Path14);
@@ -43,7 +46,8 @@ void VerifCorrelation(){
chain->SetBranchAddress("Time", &Time);
- TH2F *ICXY = new TH2F("ICXY","ICXY",1000,-1000,1000,1000,-1000,1000);
+ TH3F *ICXY = new TH3F("ICXY","ICXY",100,-550,400,100,-1000,1000,100,0,2);
+ TProfile2D *Pxyz = new TProfile2D("ICOneZeroProfile","ICOneZeroProfile",100,-550.0,400.0,100,-900.0,500.0);
int Nentries = chain->GetEntries();
//int Nentries = 1000000;
@@ -67,35 +71,58 @@ void VerifCorrelation(){
FF_DriftTime = 10* (IC->fIC_TS.at(0) - Time->GetTS_MWPC13(0)) - ((Time->GetTime_MWPC13(0)+Toff13.at(FPMW_Section))) ;
- double ICRatio = IC->fIC_raw[0]/IC->fIC_raw[1];
+ Double_t ICRatio = IC->fIC_raw[1]/IC->fIC_raw[0];
- if(IC->fIC_raw[8] >0 ) {
- ICXY->Fill(FF_IC_X,FF_IC_Y,ICRatio);
+ if(IC->fIC_raw[8] >0 && ICRatio <= 3 && ICRatio >0.1 ) {
+ ICXY->Fill(FF_IC_X,FF_DriftTime,ICRatio);
+ Pxyz->Fill(FF_IC_X,FF_DriftTime,ICRatio);
}
}
+
}
+ gStyle->SetCanvasPreferGL(true);
TCanvas *c = new TCanvas("c","c");
- ICXY->SetMinimum(0);
- ICXY->Draw("surf");
+ gStyle->SetPalette(kRainBow);
+ ICXY->SetXTitle("X");
+ ICXY->SetYTitle("DriftTime");
+ ICXY->SetZTitle("RatioIC");
+ ICXY->SetMinimum(0);
+ ICXY->Draw("GL"); // Configure the 3D view
+
+ // Update the canvas
+ c->Update();
+
+ TCanvas *c1 = new TCanvas("c1","c1");
+ gStyle->SetPalette(kRainBow);
+ Pxyz->SetXTitle("X");
+ Pxyz->SetYTitle("DriftTime");
+ Pxyz->SetZTitle("RatioIC");
+ Pxyz->SetMinimum(1.0);
+ Pxyz->SetMaximum(1.5);
+ Pxyz->Draw("GL SURF2"); // Configure the 3D view
+
+
+ TFile *ofile = new TFile("Output/RatioProfile.root","recreate");
+ Pxyz->Write();
}
vector<double> TxtToVector(const char *Path){
- string line;
- vector<double> values;
- ifstream file(Path);
-
- if (file.is_open()) {
- while (std::getline(file, line)) {
- try {
- values.push_back(std::stod(line));
- } catch (const std::invalid_argument& e) {
- std::cerr << "Invalid number in line: " << line << '\n';
- }
- }
- file.close();
- } else {
- std::cerr << "Error opening file.\n";
- }
-
- return values;
+ string line;
+ vector<double> values;
+ ifstream file(Path);
+
+ if (file.is_open()) {
+ while (std::getline(file, line)) {
+ try {
+ values.push_back(std::stod(line));
+ } catch (const std::invalid_argument& e) {
+ std::cerr << "Invalid number in line: " << line << '\n';
+ }
+ }
+ file.close();
+ } else {
+ std::cerr << "Error opening file.\n";
+ }
+
+ return values;
}