MBDLaser/macros/plot_deltT.C

0001
0002 // I need to KNOW how to  do somthing like TH1D* h1_mbdtt = hist->ProjectionY();
0003 // nevents_tt += h1_mbdtt->Integral();
0004
0005
0006
0007 #include <fstream>
0008 #include <TCanvas.h>
0009 #include <TFile.h>
0010 #include <TGraph.h>
0011 #include <TH2F.h>
0012
0013 // varaibles needed
0014 const int MAXRUNS = 20000;
0015 TFile *savefile;
0016 TFile *tfile[MAXRUNS];
0017 TTree *t;
0018 Double_t runindex[MAXRUNS];
0019 TCanvas *ac[100];
0020 TString dir;
0021 TFile *refFile_all;
0022
0023 //TGraph
0024 TGraph *graph;
0025 TGraph *gt_chBych;
0026 TGraph *gq_chBych;
0027
0028 //Time in all runs
0029 TH2 *h2_mbdtt[MAXRUNS];
0030 Double_t T_current;
0031 TFile *refFile;
0032 TH2 *h2_ref;
0033 Double_t T_ref;
0034 Double_t deltaT;
0035
0036 //charge in each ch
0037 TH1 *h_mbdq[MAXRUNS][128];
0038 Double_t Q_current[MAXRUNS][128];
0039 TH1 *h_ref_Q[128];
0040 Double_t Q_ref[128];
0041 Double_t Q_ratio[MAXRUNS][128];
0042
0043
0044 //Time in each ch
0045 TH1 *h_mbdt[MAXRUNS][128];
0046 //Double_t btmean[128][MAXRUNS];        // bt mean in each ch
0047 //Double_t btmeanerr[128][MAXRUNS];
0048 Double_t T_current_all[MAXRUNS][128];
0049 TH1 *h_ref_all[128];
0050 Double_t T_ref_all[128];
0051 Double_t deltaT_all[MAXRUNS][128];
0052
0053
0054 //Canvas
0055 TCanvas *g_Canvas;
0056 TCanvas *gt_Canvas;
0057 TCanvas *gq_Canvas;
0058
0059
0060
0061
0062 //int Colors_maker[MAXRUNS]={kRed, kBlue, kGreen,kYellow, kMagenta, kCyan, kOrange, kSpring, kTeal, kAzure, kViolet, kPink, kGray, kBlack, kRed-1, kBlue-1, kGreen-1, kYellow-1, kMagenta-1, kCyan-1, kOrange-1, kSpring-1, kTeal-1, kAzure-1, kViolet-1};
0063
0064
0065 void plot_deltT(const char *flistname = "analyzed_Laser_data2023_To_current.list")
0066 {
0067     ifstream flist(flistname);
0068     TString dstfname;
0069     int nruns = 0;
0070     TString name;
0071     TString title;
0072     int icv = 0;
0073     savefile = new TFile(dir + "Tlaser.root", "RECREATE");
0074
0075     //initialize TGraph
0076     graph = new TGraph();
0077     gt_chBych = new TGraph();
0078     gq_chBych = new TGraph();
0079
0080           // get the time of the laser in the reference run
0081           refFile = new TFile("LASER_UNCALMBD-00040759-0000.root", "READ");
0082           h2_ref = (TH2 *)refFile->Get("h2_mbdtt");
0083           T_ref = h2_ref->GetMean(1);
0084
0085     while (flist >> dstfname)
0086     {
0087         tfile[nruns] = new TFile(dstfname, "READ");
0088         title = dstfname;
0089         cout << "Title = " << dstfname << endl;
0090
0091         h2_mbdtt[nruns] = (TH2 *)tfile[nruns]->Get("h2_mbdtt");
0092
0093         if (h2_mbdtt[nruns] != nullptr)
0094         {
0095           // get the time of the laser in the current run
0096           Double_t T_current = h2_mbdtt[nruns]->GetMean(1);
0097
0098
0099           // calculate delta-T
0100           deltaT = T_current - T_ref; // cout<< "deltaT = " << deltaT<< endl;
0101
0102           // fill
0103           graph->SetPoint(nruns, nruns, deltaT);
0104
0105
0106           //here you should write h_mbdt[128] and calculate the deltaT channels by channels
0107           // to get the meean and the mean error from  North MBD/MBD time && South MBD/MBD time
0108           for (int ipmt=0; ipmt<128; ipmt++)
0109           {
0110             // appending the value of ipmt to the string h_mbdq and then assigns it to the varaible name
0111             name = "h_mbdt"; name += ipmt; // h_mbdt0, h_mbdt1, ... so on
0112
0113             h_mbdt[nruns][ipmt] = (TH1*)tfile[nruns]->Get(name);
0114
0115             // btmean[ipmt][nruns] = h_mbdt[nruns][ipmt]->GetMean();
0116             // btmeanerr[ipmt][nruns] = h_mbdt[nruns][ipmt]->GetMeanError();
0117
0118             name = "h_mbdq"; name += ipmt;
0119             h_mbdq[nruns][ipmt] = (TH1*)tfile[nruns]->Get(name);
0120
0121
0122             //cout<< name<<" \t"<< btmean[ipmt][nruns]<< "\t"<< btmeanerr[ipmt][nruns]<<endl;
0123
0124             if (h_mbdt[nruns][ipmt] != nullptr)
0125             {
0126
0127               // get the time of the laser in the current run for each channls (we should have 128 values )
0128               T_current_all[nruns][ipmt] = h_mbdt[nruns][ipmt]->GetMean();  //cout<< T_current_all << "\t"<< name<< " \t" << ipmt<<endl;
0129
0130               Q_current[nruns][ipmt] = h_mbdq[nruns][ipmt]->GetMean();
0131
0132               // get the time of the laser in the reference run for each ch
0133               //if ( refFile_all==nullptr ) refFile_all = new TFile("LASER_UNCALMBD-00041145-0000.root", "READ");
0134               if ( refFile_all==nullptr ) refFile_all = new TFile("LASER_UNCALMBD-00044961-0000.root", "READ");
0135
0136               name = "h_mbdt"; name += ipmt;
0137               h_ref_all[ipmt] = (TH1*)refFile_all->Get(name);
0138               T_ref_all[ipmt] = h_ref_all[ipmt]->GetMean();
0139
0140
0141               name = "h_mbdq"; name += ipmt;
0142               h_ref_Q[ipmt] = (TH1*)refFile_all->Get(name);
0143               Q_ref[ipmt] = h_ref_Q[ipmt]->GetMean();
0144
0145
0146
0147               //calculation
0148               deltaT_all[nruns][ipmt] = T_current_all[nruns][ipmt] - T_ref_all[ipmt];
0149
0150               Q_ratio[nruns][ipmt] = Q_current[nruns][ipmt] / Q_ref[ipmt];
0151
0152
0153               //cout<< "deltaT = " << deltaT_all[nruns][ipmt]<< " \t" <<"T_current in each ch = " << T_current_all[nruns][ipmt] << "\t" << "T_ref in each ch = "<<T_ref_all[ipmt] << " \t" << name <<"\t"<< nruns<< endl;
0154
0155               if ( ipmt==0 || ipmt==127 ) cout<< "Q_ratio = " << Q_ratio[nruns][ipmt]<< " \t" <<"Q_current in each ch = " << Q_current[nruns][ipmt] << "\t" << "Q_ref in each ch = "<<Q_ref[ipmt] << " \t" << name <<"\t"<< nruns<< endl;
0156
0157               //Fill
0158               int n = gt_chBych->GetN();
0159               gt_chBych->SetPoint(n,nruns, deltaT_all[nruns][ipmt]);
0160
0161               n = gq_chBych->GetN();
0162               gq_chBych->SetPoint(n,nruns, Q_ratio[nruns][ipmt] );
0163
0164               // set the color maker
0165               //int makerColor = kRed + nruns % 24 ; // to use different color for each run ; 24 different colors
0166               //int colors_maker_index = nruns % MAXRUNS;
0167               //gt_chBych->SetMarkerColor(Colors_maker[(ipmt) % 25]); // the remainder is used as index to slect the colors from color maker[] so the rang is 0 to 24 (to find the avalible colors)
0168
0169               // maker color directly for each ch
0170               /*   for ( int ipmt_col = 0; ipmt_col<128; ipmt_col++)
0171                    {
0172                    int channelColor = kRed +(nruns *128 + ipmt_col) % 24;
0173               //gt_chBych->SetMarkerColor(channelColor);
0174
0175               }
0176               */
0177               //refFile_all->Close(); // I had to add it avoiding error like too many files open
0178               //delete refFile_all;
0179
0180
0181             }
0182
0183           }
0184
0185         }
0186         else
0187         {
0188           cout << "Histogram 'h2_mbdtt' not found in file: " << dstfname << endl;
0189         }
0190
0191         ac[icv] = new TCanvas(Form("h2_mbdtt_%d", nruns), Form("btt_vs_ch_%d", nruns), 1000, 500);
0192         h2_mbdtt[nruns]->SetTitle(title + nruns);
0193         h2_mbdtt[nruns]->GetXaxis()->SetTitle("ch");
0194         h2_mbdtt[nruns]->GetYaxis()->SetTitle("Time[ADC]");
0195         //h2_mbdtt[nruns]->Draw("colz");
0196         gPad->SetLogz(1);
0197         ac[icv]->SaveAs(Form("h2_mbdtt_%d.png", nruns));
0198
0199         nruns++;
0200         icv++;
0201
0202
0203
0204     }
0205     if(tfile[nruns] != nullptr) // to avoid such error like opening too many files
0206     {
0207       tfile[nruns]->Close();
0208       delete tfile[nruns];
0209     }
0210
0211
0212     g_Canvas = new TCanvas("deltaT_graph", "Delta-T vs Run Index ", 1000, 500);
0213     graph->SetTitle("Delta-T vs Run Index");
0214     graph->GetXaxis()->SetTitle("Run Index");
0215     graph->GetYaxis()->SetTitle("Delta-T[ns]");
0216     graph->SetMarkerSize(1);
0217     graph->SetMarkerStyle(20);
0218     graph->SetMarkerColor(kRed);
0219     graph->Draw("AP");
0220     g_Canvas->SaveAs("Delta-T_vs_Run_Index.png");
0221
0222
0223     gt_Canvas= new TCanvas("deltaT_graph_128ch", "Delta-T vs Run Index128 ch ", 1000, 500);
0224     gt_chBych->SetTitle("Delta-T vs Run Index");
0225     gt_chBych->GetXaxis()->SetTitle("Run Index");
0226     gt_chBych->GetYaxis()->SetTitle("Delta-T[ns]");
0227     gt_chBych->SetMarkerSize(0.2);
0228     gt_chBych->SetMarkerStyle(5);
0229     gt_chBych->SetMarkerColor(2.);
0230     gt_chBych->Draw("AP");
0231     gt_Canvas->SaveAs("Delta_T_vs_Run_Index128ch.png");
0232
0233
0234     gq_Canvas= new TCanvas("Q_ratio_graph_128ch", "Q_ratio vs Run Index128 ", 1000, 500);
0235     gq_chBych->SetTitle("Q_ratio vs Run Index");
0236     gq_chBych->GetXaxis()->SetTitle("Run Index");
0237     gq_chBych->GetYaxis()->SetTitle("Q_ratio[ADC]");
0238     gq_chBych->SetMarkerSize(0.2);
0239     gq_chBych->SetMarkerStyle(5);
0240     gq_chBych->SetMarkerColor(kBlue);
0241     gq_chBych->Draw("AP");
0242     gq_Canvas->SaveAs("Q_ratio_vs_Run_Index128ch.png");
0243
0244
0245
0246
0247
0248     savefile->Write();
0249     savefile->Close();
0250
0251     cout << "Processed " << nruns << " runs." << endl;
0252 }