analysis_INTT/plot/Hitmap_evtsel.C

0001 #include <TCanvas.h>
0002 #include <TFile.h>
0003 #include <TH2F.h>
0004 #include <TStyle.h>
0005 #include <TSystem.h>
0006 #include <TTree.h>
0007 #include <algorithm>
0008 #include <cstdlib>
0009 #include <iostream>
0010 #include <string>
0011 #include <vector>
0012
0013 float TickSize = 0.03;
0014 float AxisTitleSize = 0.05;
0015 float AxisLabelSize = 0.04;
0016 float LeftMargin = 0.15;
0017 float RightMargin = 0.14;
0018 float TopMargin = 0.08;
0019 float BottomMargin = 0.13;
0020
0021 int m_inttFeeOffset = 23; // 23 is the offset for INTT in streaming mode; https://github.com/sPHENIX-Collaboration/coresoftware/blob/master/offline/packages/intt/InttCombinedRawDataDecoder.h#L67C3-L67C76
0022
0023 std::vector<int> colorlist = {kTVibBlue, kTVibTeal, kTVibOrange, kTBriPurple};
0024
0025 // void Draw_2Dhist(TH2F *hist, bool IsData, bool logz, bool norm1, float rmargin, const char *XaxisName, const char *YaxisName, const char *ZaxisName, const char *drawopt, const char *outname)
0026 // {
0027 //     // get the INTT BCO number from the histogram name: split the name by "_", the bco number is the last element
0028 //     std::string histname = hist->GetName();
0029 //     std::string bco = histname.substr(histname.find_last_of("_") + 1);
0030
0031 //     TCanvas *c = new TCanvas("c", "c", 1000, 700);
0032 //     if (logz)
0033 //     {
0034 //         c->SetLogz();
0035 //     }
0036 //     c->cd();
0037 //     if (strcmp(ZaxisName, "") == 0)
0038 //     {
0039 //         gPad->SetRightMargin(rmargin);
0040 //     }
0041 //     else
0042 //     {
0043 //         gPad->SetRightMargin(rmargin + 0.03);
0044 //     }
0045
0046 //     gPad->SetTopMargin(TopMargin);
0047 //     gPad->SetLeftMargin(LeftMargin);
0048 //     gPad->SetBottomMargin(BottomMargin);
0049 //     if (norm1)
0050 //     {
0051 //         hist->Scale(1. / hist->Integral(-1, -1, -1, -1));
0052 //     }
0053 //     hist->GetXaxis()->SetTitle(XaxisName);
0054 //     hist->GetYaxis()->SetTitle(YaxisName);
0055 //     hist->GetXaxis()->SetTickSize(TickSize);
0056 //     hist->GetYaxis()->SetTickSize(TickSize);
0057 //     hist->GetXaxis()->SetTitleSize(AxisTitleSize);
0058 //     hist->GetYaxis()->SetTitleSize(AxisTitleSize);
0059 //     hist->GetXaxis()->SetLabelSize(AxisLabelSize);
0060 //     hist->GetYaxis()->SetLabelSize(AxisLabelSize);
0061 //     if (strcmp(ZaxisName, "") != 0)
0062 //     {
0063 //         hist->GetZaxis()->SetTitle(ZaxisName);
0064 //         hist->GetZaxis()->SetTitleSize(AxisTitleSize);
0065 //         hist->GetZaxis()->SetTitleOffset(1.1);
0066 //     }
0067
0068 //     hist->GetXaxis()->SetTitleOffset(1.1);
0069 //     hist->GetYaxis()->SetTitleOffset(1.3);
0070 //     hist->GetZaxis()->SetLabelSize(AxisLabelSize);
0071 //     // hist->GetZaxis()->SetRangeUser(hist->GetMinimum(0), hist->GetMaximum());
0072 //     hist->SetContour(1000);
0073 //     hist->Draw(drawopt);
0074
0075 //     // Draw the INTT BCO number
0076 //     TLatex *tex = new TLatex();
0077 //     tex->SetNDC();
0078 //     // tex->SetTextFont(43);
0079 //     tex->SetTextSize(0.035);
0080 //     tex->DrawLatex(0.2, 0.8, Form("INTT BCO: %s", bco.c_str()));
0081
0082 //     c->RedrawAxis();
0083 //     c->Draw();
0084 //     c->SaveAs((std::string(outname) + ".pdf").c_str());
0085 //     c->SaveAs((std::string(outname) + ".png").c_str());
0086 //     if (c)
0087 //     {
0088 //         c->Close();
0089 //         gSystem->ProcessEvents();
0090 //         delete c;
0091 //         c = nullptr;
0092 //         gSystem->ProcessEvents();
0093 //     }
0094 // }
0095
0096 TCanvas *set2Dhist(TH2F *hist, bool IsBadCorr, bool IsData, bool logz, bool norm1, float rmargin, const char *XaxisName, const char *YaxisName, const char *ZaxisName, const char *drawopt, std::vector<std::string> plotinfo)
0097 {
0098     // get the INTT BCO number from the histogram name: split the name by "_", the bco number is the last element
0099     std::string histname = hist->GetName();
0100     std::string bco = histname.substr(histname.find_last_of("_") + 1);
0101
0102     TCanvas *c = new TCanvas(Form("c_%s", histname.c_str()), Form("c_%s", histname.c_str()), 1000, 700);
0103     if (logz)
0104     {
0105         c->SetLogz();
0106     }
0107     c->cd();
0108     if (strcmp(ZaxisName, "") == 0)
0109     {
0110         gPad->SetRightMargin(rmargin);
0111     }
0112     else
0113     {
0114         gPad->SetRightMargin(rmargin + 0.03);
0115     }
0116
0117     gPad->SetTopMargin(TopMargin);
0118     gPad->SetLeftMargin(LeftMargin);
0119     gPad->SetBottomMargin(BottomMargin);
0120     if (norm1)
0121     {
0122         hist->Scale(1. / hist->Integral(-1, -1, -1, -1));
0123     }
0124     hist->GetXaxis()->SetTitle(XaxisName);
0125     hist->GetYaxis()->SetTitle(YaxisName);
0126     hist->GetXaxis()->SetTickSize(TickSize);
0127     hist->GetYaxis()->SetTickSize(TickSize);
0128     hist->GetXaxis()->SetTitleSize(AxisTitleSize);
0129     hist->GetYaxis()->SetTitleSize(AxisTitleSize);
0130     hist->GetXaxis()->SetLabelSize(AxisLabelSize);
0131     hist->GetYaxis()->SetLabelSize(AxisLabelSize);
0132     if (strcmp(ZaxisName, "") != 0)
0133     {
0134         hist->GetZaxis()->SetTitle(ZaxisName);
0135         hist->GetZaxis()->SetTitleSize(AxisTitleSize);
0136         hist->GetZaxis()->SetTitleOffset(1.1);
0137     }
0138
0139     hist->GetXaxis()->SetTitleOffset(1.1);
0140     hist->GetYaxis()->SetTitleOffset(1.3);
0141     hist->GetZaxis()->SetLabelSize(AxisLabelSize);
0142     // hist->GetZaxis()->SetRangeUser(hist->GetMinimum(0), hist->GetMaximum());
0143     hist->SetContour(1000);
0144     hist->Draw(drawopt);
0145
0146     // Draw the plotinfo strings
0147     TLatex *tex = new TLatex();
0148     tex->SetNDC();
0149     // tex->SetTextFont(43);
0150     tex->SetTextSize(0.03);
0151     if (IsBadCorr)
0152     {
0153         tex->SetTextColor(kRed);
0154         for (int i = 0; i < plotinfo.size(); i++)
0155         {
0156             tex->DrawLatex(0.18, 0.87 - 0.03 * i, plotinfo[i].c_str());
0157         }
0158     }
0159     else
0160     {
0161         for (int i = 0; i < plotinfo.size(); i++)
0162         {
0163             tex->DrawLatex(0.18, 0.87 - 0.03 * i, plotinfo[i].c_str());
0164         }
0165     }
0166
0167     c->RedrawAxis();
0168
0169     return c;
0170 }
0171
0172 void Hitmap_evtsel()
0173 {
0174     bool plothitmap = true;
0175
0176     gStyle->SetPalette(kRainBow);
0177
0178     // TString filename = "/sphenix/tg/tg01/commissioning/INTT/work/cwshih/seflgendata/run_54280/From_official_INTTRAWHIT_DST_ColdHotChannelMask/completed/Trkr_DST_physics_intt-54280_HotChannel_BCOCut_CaloCombined_00000.root";
0179     // TString filename = "/sphenix/tg/tg01/commissioning/INTT/work/cwshih/seflgendata/run_54280/From_official_INTTRAWHIT_DST_ColdHotChannelMask/completed/Trkr_DST_physics_intt-54280_HotChannel_BCOCut_CaloCombined_merged.root";
0180     std::string plotdir = "./TrkrHit/Hitmap_InnerOuterCorrBranch_perBCO_Run54280/";
0181     system(("mkdir -p " + plotdir).c_str());
0182
0183     int z_offset[4] = {5, 0, 13, 21};
0184     int layer_offset[4] = {0, 12, 24, 40};
0185     int layerhitset_offset[4] = {0, 48, 96, 160};
0186
0187     // integrated hit map for non-corr branch
0188     TH2F *hM_hitmap_RowColumnFlatten_noncorr = new TH2F("hM_hitmap_RowColumnFlatten_noncorr", "hM_hitmap_RowColumnFlatten_noncorr", 4 * 32, 0, 4 * 32, 360 * 16, 0, 360 * 16);
0189     TH2F *hM_hitmap_RowColumnFlatten_noncorr_weightAdc = new TH2F("hM_hitmap_RowColumnFlatten_noncorr_weightAdc", "hM_hitmap_RowColumnFlatten_noncorr_weightAdc", 4 * 32, 0, 4 * 32, 360 * 16, 0, 360 * 16);
0190     TH2F *hM_hitmap_RowColumnFlatten_pm4noncorr = new TH2F("hM_hitmap_RowColumnFlatten_pm4noncorr", "hM_hitmap_RowColumnFlatten_pm4noncorr", 4 * 32, 0, 4 * 32, 360 * 16, 0, 360 * 16);
0191     TH2F *hM_hitmap_RowColumnFlatten_pm4noncorr_weightAdc = new TH2F("hM_hitmap_RowColumnFlatten_pm4noncorr_weightAdc", "hM_hitmap_RowColumnFlatten_pm4noncorr_weightAdc", 4 * 32, 0, 4 * 32, 360 * 16, 0, 360 * 16);
0192
0193     // hit adc distribution
0194     TH1F *hM_TrkrhitADC_noncorr = new TH1F("hM_TrkrhitADC_noncorr", "hM_TrkrhitADC_noncorr", 251, -0.5, 250.5);
0195     TH1F *hM_TrkrhitADC_pm4noncorr = new TH1F("hM_TrkrhitADC_pm4noncorr", "hM_TrkrhitADC_pm4noncorr", 251, -0.5, 250.5);
0196
0197     // event bco difference = BCO of (i+1)th event - BCO of ith event
0198     TH1F *hM_BCOdiff_noncorrnext = new TH1F("hM_BCOdiff_noncorrnext", "hM_BCOdiff_noncorrnext", 131, -0.5, 130.5);
0199
0200     // TH2F *hM_HitTimeBucket_Layer_noncorrnext = new TH2F("hM_HitTimeBucket_Layer_noncorrnext", "hM_HitTimeBucket_Layer_noncorrnext", 131, -0.5, 130.5, 4, 2.5, 6.5);
0201     std::vector<TH1F *> list_hM_HitTimeBucket_Layer_noncorrnext;
0202     for (int l = 0; l < 4; l++)
0203     {
0204         list_hM_HitTimeBucket_Layer_noncorrnext.push_back(new TH1F(Form("hM_HitTimeBucket_Layer_noncorrnext_Layer%d", l + 3), Form("hM_HitTimeBucket_Layer_noncorrnext_Layer%d", l + 3), 131, -0.5, 130.5));
0205     }
0206
0207     // output text file
0208     std::ofstream outtxtfile(Form("%s/EventList.txt", plotdir.c_str()), std::ios::out);
0209
0210     // TFile *fin = new TFile(filename, "READ");
0211     // TTree *t = (TTree *)fin->Get("EventTree");
0212     TChain *t = new TChain("EventTree");
0213     t->Add("/sphenix/tg/tg01/hf/hjheng/ppg02/dst/Data_Run54280_20241127_ProdA2024/ntuple_00*.root");
0214     int event;
0215     float MBD_z_vtx;
0216     int NTrkrhits, NClus_Layer1, NClus;
0217     uint64_t INTT_BCO;
0218     bool is_min_bias;
0219     vector<int> *firedTriggers = 0;
0220     std::vector<uint64_t> *InttRawHit_bco = 0;
0221     std::vector<uint16_t> *InttRawHit_FPHX_BCO = 0;
0222     std::vector<uint16_t> *TrkrHitRow = 0, *TrkrHitColumn = 0, *TrkrHitADC = 0;
0223     std::vector<uint8_t> *TrkrHitLadderZId = 0, *TrkrHitLadderPhiId = 0, *TrkrHitLayer = 0;
0224     std::vector<float> *ClusPhiSize = 0;
0225     t->SetBranchAddress("event", &event);
0226     t->SetBranchAddress("is_min_bias", &is_min_bias);
0227     t->SetBranchAddress("MBD_z_vtx", &MBD_z_vtx);
0228     t->SetBranchAddress("INTT_BCO", &INTT_BCO);
0229     t->SetBranchAddress("NTrkrhits", &NTrkrhits);
0230     t->SetBranchAddress("NClus_Layer1", &NClus_Layer1);
0231     t->SetBranchAddress("NClus", &NClus);
0232     t->SetBranchAddress("firedTriggers", &firedTriggers);
0233     if (plothitmap)
0234     {
0235         t->SetBranchAddress("InttRawHit_bco", &InttRawHit_bco);
0236         t->SetBranchAddress("InttRawHit_FPHX_BCO", &InttRawHit_FPHX_BCO);
0237         t->SetBranchAddress("TrkrHitRow", &TrkrHitRow);
0238         t->SetBranchAddress("TrkrHitColumn", &TrkrHitColumn);
0239         t->SetBranchAddress("TrkrHitADC", &TrkrHitADC);
0240         t->SetBranchAddress("TrkrHitLadderZId", &TrkrHitLadderZId);
0241         t->SetBranchAddress("TrkrHitLadderPhiId", &TrkrHitLadderPhiId);
0242         t->SetBranchAddress("TrkrHitLayer", &TrkrHitLayer);
0243         t->SetBranchAddress("ClusPhiSize", &ClusPhiSize);
0244     }
0245
0246     // Pre-allocate vectors to avoid reallocations
0247     int Nbco2plot = 1500;
0248     std::map<int, uint64_t> map_evtbco;
0249
0250     bool fillbcodiffnext = false;
0251     uint64_t bco_ofinterest = 0;
0252
0253     // Get total entries once instead of calling each time
0254     Long64_t nEntries = t->GetEntries();
0255     std::cout << "Total entries: " << nEntries << std::endl;
0256
0257     // Cache divisions and comparisons
0258     const float outer_inner_ratio_threshold = 1.2f;
0259     const float mbd_z_vtx_threshold = 10.0f;
0260     const int nclus_layer1_threshold = 500;
0261
0262     for (Long64_t ev = 0; ev < nEntries; ++ev)
0263     {
0264         if (map_evtbco.size() >= Nbco2plot)
0265             break;
0266
0267         t->GetEntry(ev);
0268
0269         if (fillbcodiffnext)
0270         {
0271             hM_BCOdiff_noncorrnext->Fill(INTT_BCO - bco_ofinterest);
0272             fillbcodiffnext = false;
0273         }
0274
0275         // Combine conditions to reduce branches
0276         // bool passes_cuts = (is_min_bias) && (std::abs(MBD_z_vtx) <= mbd_z_vtx_threshold) && (NClus_Layer1 >= nclus_layer1_threshold) && ((static_cast<float>(NClus - NClus_Layer1) / NClus_Layer1) >= outer_inner_ratio_threshold);
0277         // (NClus_Layer1<(NClus-NClus_Layer1)-5*TMath::Sqrt(NClus-NClus_Layer1))
0278         bool passes_cuts = (is_min_bias && NClus_Layer1 < (NClus - NClus_Layer1) - 5 * TMath::Sqrt(NClus - NClus_Layer1));
0279
0280         if (!passes_cuts)
0281             continue;
0282
0283         cout << "Event: " << ev << " , map_evtbco.size(): " << map_evtbco.size() << ", INTT_BCO: " << INTT_BCO << " NClus (inner): " << NClus_Layer1 << " NClus (outer): " << (NClus - NClus_Layer1) << " Outer/Inner: " << static_cast<float>(NClus - NClus_Layer1) / static_cast<float>(NClus_Layer1) << endl;
0284         cout << "Event: " << ev << " has trigger bit: ";
0285         for (int j = 0; j < firedTriggers->size(); j++)
0286         {
0287             cout << firedTriggers->at(j) << " ";
0288         }
0289         cout << endl;
0290
0291         if (plothitmap)
0292             map_evtbco[ev] = INTT_BCO;
0293
0294         bco_ofinterest = INTT_BCO;
0295         fillbcodiffnext = true;
0296     }
0297
0298     if (plothitmap)
0299     {
0300         for (auto [ev, inttbco] : map_evtbco)
0301         {
0302             uint64_t bco_ofinterest = inttbco;
0303
0304             TH1F *hM_hitTimebucket_noncorr = new TH1F(Form("hM_hitTimebucket_noncorr_%lu", bco_ofinterest), Form("hM_hitTimebucket_noncorr_%lu", bco_ofinterest), 131, -0.5, 130.5);
0305             TH1F *hM_hitTimebucket_noncorrnext = new TH1F(Form("hM_hitTimebucket_noncorrnext_%lu", bco_ofinterest), Form("hM_hitTimebucket_noncorrnext_%lu", bco_ofinterest), 131, -0.5, 130.5);
0306             TH1F *hM_hitTimebucket_noncorrnext_wbcointerest = new TH1F(Form("hM_hitTimebucket_noncorrnext_wbcointerest_%lu", bco_ofinterest), Form("hM_hitTimebucket_noncorrnext_wbcointerest_%lu", bco_ofinterest), 131, -0.5, 130.5);
0307
0308             // look at +-5 events around the selected event
0309             for (int i = -5; i <= 5; i++)
0310             {
0311                 if (ev + i < 0 || ev + i >= t->GetEntriesFast())
0312                     continue;
0313
0314                 std::cout << "Processing event, BCO " << inttbco << std::endl;
0315                 t->GetEntry(ev + i);
0316
0317                 TH2F *hM_hitmap_RowColumnFlatten = new TH2F(Form("hM_hitmap_RowColumnFlatten_%lu", INTT_BCO), Form("hM_hitmap_RowColumnFlatten_%lu", INTT_BCO), 4 * 32, 0, 4 * 32, 360 * 16, 0, 360 * 16);
0318                 TH2F *hM_hitmap_RowColumnFlatten_weighttimebucket = new TH2F(Form("hM_hitmap_RowColumnFlatten_weighttimebucket_%lu", INTT_BCO), Form("hM_hitmap_RowColumnFlatten_weighttimebucket_%lu", INTT_BCO), 4 * 32, 0, 4 * 32, 360 * 16, 0, 360 * 16);
0319
0320                 float avg_trkrhitadc = 0;
0321                 for (int h = 0; h < NTrkrhits; h++)
0322                 {
0323                     int flat_x = TrkrHitColumn->at(h) + z_offset[TrkrHitLadderZId->at(h)] + 32 * (TrkrHitLayer->at(h) - 3);
0324                     int flat_y = TrkrHitRow->at(h) + 360 * TrkrHitLadderPhiId->at(h);
0325                     int time_bucket = (InttRawHit_FPHX_BCO->at(h) - (InttRawHit_bco->at(h) & 0x7fU) - m_inttFeeOffset + 128) % 128;
0326
0327                     hM_hitmap_RowColumnFlatten->Fill(flat_x, flat_y, TrkrHitADC->at(h) + 1);
0328                     hM_hitmap_RowColumnFlatten_weighttimebucket->Fill(flat_x, flat_y, time_bucket);
0329                     avg_trkrhitadc += (TrkrHitADC->at(h) + 1);
0330
0331                     if (i == 0)
0332                     {
0333                         if (h == 0)
0334                         {
0335                             // print to output text file
0336                             outtxtfile << "EOI, Event: " << event << ", INTT_BCO: " << INTT_BCO << " has trigger bit: ";
0337                             for (int j = 0; j < firedTriggers->size(); j++)
0338                             {
0339                                 outtxtfile << firedTriggers->at(j) << " ";
0340                             }
0341                             outtxtfile << endl;
0342                         }
0343
0344                         if (is_min_bias)
0345                         {
0346                             hM_hitmap_RowColumnFlatten_noncorr->Fill(flat_x, flat_y);
0347                             hM_hitmap_RowColumnFlatten_noncorr_weightAdc->Fill(flat_x, flat_y, TrkrHitADC->at(h) + 1);
0348                             hM_TrkrhitADC_noncorr->Fill(TrkrHitADC->at(h) + 1);
0349                         }
0350
0351                         int time_bucket = (InttRawHit_FPHX_BCO->at(h) - (InttRawHit_bco->at(h) & 0x7fU) - m_inttFeeOffset + 128) % 128;
0352                         hM_hitTimebucket_noncorr->Fill(time_bucket);
0353                     }
0354                     else
0355                     {
0356                         if (i == 1)
0357                         {
0358                             if (h == 0)
0359                             {
0360                                 outtxtfile << "Event after EOI, Event: " << event << ", INTT_BCO: " << INTT_BCO << " has trigger bit: ";
0361                                 for (int j = 0; j < firedTriggers->size(); j++)
0362                                 {
0363                                     outtxtfile << firedTriggers->at(j) << " ";
0364                                 }
0365                                 outtxtfile << endl;
0366                             }
0367
0368                             int time_bucket = (InttRawHit_FPHX_BCO->at(h) - (InttRawHit_bco->at(h) & 0x7fU) - m_inttFeeOffset + 128) % 128;
0369                             hM_hitTimebucket_noncorrnext->Fill(time_bucket);
0370                             int time_bucket_wbcointerest = (InttRawHit_FPHX_BCO->at(h) - (bco_ofinterest & 0x7fU) - m_inttFeeOffset + 128) % 128;
0371                             hM_hitTimebucket_noncorrnext_wbcointerest->Fill(time_bucket_wbcointerest);
0372
0373                             // hM_HitTimeBucket_Layer_noncorrnext->Fill(time_bucket, TrkrHitLayer->at(h));
0374                             list_hM_HitTimeBucket_Layer_noncorrnext[TrkrHitLayer->at(h) - 3]->Fill(time_bucket);
0375                         }
0376
0377                         if (is_min_bias)
0378                         {
0379                             hM_hitmap_RowColumnFlatten_pm4noncorr->Fill(flat_x, flat_y);
0380                             hM_hitmap_RowColumnFlatten_pm4noncorr_weightAdc->Fill(flat_x, flat_y, TrkrHitADC->at(h) + 1);
0381                             hM_TrkrhitADC_pm4noncorr->Fill(TrkrHitADC->at(h) + 1);
0382                         }
0383                     }
0384                 }
0385
0386                 avg_trkrhitadc /= NTrkrhits;
0387
0388                 // plotinfo strings
0389                 std::string evtinfo = Form("Event Id: %d", event);
0390                 std::string bcoinfo = (INTT_BCO == bco_ofinterest) ? Form("BCO*: %lu", INTT_BCO) : Form("BCO: %lu", INTT_BCO);
0391                 std::string nclusinfo1 = Form("Number of clusters (inner): %d", NClus_Layer1);
0392                 std::string nclusinfo2 = Form("Number of clusters (outer): %d", NClus - NClus_Layer1);
0393                 std::string nclusinfo3 = Form("Outer/Inner: %.4f", static_cast<float>(NClus - NClus_Layer1) / static_cast<float>(NClus_Layer1));
0394                 std::string avgadcinfo = Form("Average TrkrHit ADC: %.2f", avg_trkrhitadc);
0395                 std::string minbiasinfo = (is_min_bias) ? "Is Min Bias" : "Is not Min Bias";
0396                 std::vector<std::string> plotinfo = {evtinfo, bcoinfo, nclusinfo1, nclusinfo2, nclusinfo3, minbiasinfo};
0397
0398                 TCanvas *c = set2Dhist(hM_hitmap_RowColumnFlatten, (INTT_BCO == bco_ofinterest), true, false, false, RightMargin, "(z and layer) index", "phi index", "TrkrHit ADC", "colz", plotinfo);
0399                 hM_hitmap_RowColumnFlatten_weighttimebucket->GetZaxis()->SetRangeUser(0, 128);
0400                 // hM_hitmap_RowColumnFlatten_weighttimebucket->SetContour(6);
0401                 TCanvas *c_weighttimebucket = set2Dhist(hM_hitmap_RowColumnFlatten_weighttimebucket, (INTT_BCO == bco_ofinterest), true, false, false, RightMargin, "(z and layer) index", "phi index", "Time bucket", "colz", plotinfo);
0402
0403                 if (i == -5)
0404                 {
0405                     c->Print((plotdir + Form("Hitmap_evtsel_%lu.pdf[", bco_ofinterest)).c_str());
0406                     c_weighttimebucket->Print((plotdir + Form("Hitmap_evtsel_weighttimebucket_%lu.pdf[", bco_ofinterest)).c_str());
0407                 }
0408                 else if (i == 5)
0409                 {
0410                     c->Print((plotdir + Form("Hitmap_evtsel_%lu.pdf]", bco_ofinterest)).c_str());
0411                     c_weighttimebucket->Print((plotdir + Form("Hitmap_evtsel_weighttimebucket_%lu.pdf]", bco_ofinterest)).c_str());
0412                 }
0413                 else
0414                 {
0415                     c->Print((plotdir + Form("Hitmap_evtsel_%lu.pdf", bco_ofinterest)).c_str());
0416                     c_weighttimebucket->Print((plotdir + Form("Hitmap_evtsel_weighttimebucket_%lu.pdf", bco_ofinterest)).c_str());
0417                 }
0418
0419                 delete hM_hitmap_RowColumnFlatten;
0420                 delete hM_hitmap_RowColumnFlatten_weighttimebucket;
0421
0422                 if (c)
0423                 {
0424                     c->Close();
0425                     gSystem->ProcessEvents();
0426                     delete c;
0427                     c = nullptr;
0428                     gSystem->ProcessEvents();
0429                 }
0430
0431                 if (c_weighttimebucket)
0432                 {
0433                     c_weighttimebucket->Close();
0434                     gSystem->ProcessEvents();
0435                     delete c_weighttimebucket;
0436                     c_weighttimebucket = nullptr;
0437                     gSystem->ProcessEvents();
0438                 }
0439             }
0440
0441             TCanvas *c_timebucket = new TCanvas("c_timebucket", "c_timebucket", 800, 700);
0442             c_timebucket->cd();
0443             // c_timebucket->SetLogy();
0444             hM_hitTimebucket_noncorr->GetXaxis()->SetTitle("Time bucket");
0445             hM_hitTimebucket_noncorr->GetYaxis()->SetTitle("Entries");
0446             hM_hitTimebucket_noncorr->GetXaxis()->SetTitleOffset(1.3);
0447             hM_hitTimebucket_noncorr->GetYaxis()->SetTitleOffset(1.5);
0448             hM_hitTimebucket_noncorr->GetYaxis()->SetRangeUser(0, hM_hitTimebucket_noncorr->GetMaximum() * 1.4);
0449             hM_hitTimebucket_noncorr->SetLineColor(kTBriBlue);
0450             hM_hitTimebucket_noncorr->SetLineWidth(2);
0451             hM_hitTimebucket_noncorr->SetMarkerColor(kTBriBlue);
0452             hM_hitTimebucket_noncorr->SetMarkerStyle(20);
0453             hM_hitTimebucket_noncorr->SetMarkerSize(1);
0454             hM_hitTimebucket_noncorrnext->SetLineColor(kTBriRed);
0455             hM_hitTimebucket_noncorrnext->SetLineWidth(2);
0456             hM_hitTimebucket_noncorrnext->SetMarkerColor(kTBriRed);
0457             hM_hitTimebucket_noncorrnext->SetMarkerStyle(20);
0458             hM_hitTimebucket_noncorrnext->SetMarkerSize(1);
0459             hM_hitTimebucket_noncorrnext_wbcointerest->SetLineColor(kTBriGreen);
0460             hM_hitTimebucket_noncorrnext_wbcointerest->SetLineWidth(2);
0461             hM_hitTimebucket_noncorrnext_wbcointerest->SetMarkerColor(kTBriGreen);
0462             hM_hitTimebucket_noncorrnext_wbcointerest->SetMarkerStyle(20);
0463             hM_hitTimebucket_noncorrnext_wbcointerest->SetMarkerSize(1);
0464             hM_hitTimebucket_noncorr->Draw("PE");
0465             hM_hitTimebucket_noncorrnext->Draw("PE same");
0466             hM_hitTimebucket_noncorrnext_wbcointerest->Draw("PE same");
0467             c_timebucket->RedrawAxis();
0468
0469             static TLegend *leg = nullptr;
0470             if (!leg)
0471             {
0472                 leg = new TLegend(0.35, 0.8, 0.65, 0.9);
0473                 leg->AddEntry(hM_hitTimebucket_noncorr, "Event of interest (EOI)", "pel");
0474                 leg->AddEntry(hM_hitTimebucket_noncorrnext, "Next event", "pel");
0475                 leg->AddEntry(hM_hitTimebucket_noncorrnext_wbcointerest, "Next event, w.r.t the BCO of EOI", "pel");
0476                 leg->SetBorderSize(0);
0477                 leg->SetFillStyle(0);
0478                 leg->SetTextSize(0.04);
0479             }
0480             leg->Draw();
0481
0482             std::string pdfPath = plotdir + Form("HitTimebucket_evtsel_%lu.pdf", bco_ofinterest);
0483             std::string pngPath = plotdir + Form("HitTimebucket_evtsel_%lu.png", bco_ofinterest);
0484             c_timebucket->SaveAs(pdfPath.c_str());
0485             c_timebucket->SaveAs(pngPath.c_str());
0486
0487             if (c_timebucket)
0488             {
0489                 c_timebucket->Close();
0490                 gSystem->ProcessEvents();
0491                 delete c_timebucket;
0492                 c_timebucket = nullptr;
0493                 gSystem->ProcessEvents();
0494             }
0495         }
0496
0497         // Create and save hit map plots
0498         auto saveHitMapPlot = [&](TH2F *hist, const std::string &name, const std::vector<std::string> &labels)
0499         {
0500             TCanvas *c = set2Dhist(hist, false, true, false, false, RightMargin, "(z and layer) index", "phi index", "Average Trkrhit ADC", "colz", labels);
0501             std::string pdfPath = plotdir + name + ".pdf";
0502             std::string pngPath = plotdir + name + ".png";
0503             c->SaveAs(pdfPath.c_str());
0504             c->SaveAs(pngPath.c_str());
0505             delete c;
0506         };
0507
0508         TH2F *hM_hitmap_RowColumnFlatten_noncorr_avgAdc = (TH2F *)hM_hitmap_RowColumnFlatten_noncorr_weightAdc->Clone("hM_hitmap_RowColumnFlatten_noncorr_avgAdc");
0509         hM_hitmap_RowColumnFlatten_noncorr_avgAdc->Divide(hM_hitmap_RowColumnFlatten_noncorr);
0510         hM_hitmap_RowColumnFlatten_noncorr_avgAdc->GetZaxis()->SetRangeUser(0, 210);
0511         saveHitMapPlot(hM_hitmap_RowColumnFlatten_noncorr_avgAdc, "Hitmap_evtsel_noncorr", {"Integrated hit map", "|MBD Zvtx|#leq10cm", "NClus_{inner}#geq1000", "Outer/Inner#geq1.2", "Is Min Bias"});
0512
0513         TH2F *hM_hitmap_RowColumnFlatten_pm4noncorr_avgAdc = (TH2F *)hM_hitmap_RowColumnFlatten_pm4noncorr_weightAdc->Clone("hM_hitmap_RowColumnFlatten_pm4noncorr_avgAdc");
0514         hM_hitmap_RowColumnFlatten_pm4noncorr_avgAdc->Divide(hM_hitmap_RowColumnFlatten_pm4noncorr);
0515         hM_hitmap_RowColumnFlatten_pm4noncorr_avgAdc->GetZaxis()->SetRangeUser(0, 210);
0516         saveHitMapPlot(hM_hitmap_RowColumnFlatten_pm4noncorr_avgAdc, "Hitmap_evtsel_pm4noncorr", {"Integrated hit map", "[-4, -1] and [+1, +4] events", "Is Min Bias"});
0517
0518         // Normalize and save Trkrhit ADC plots
0519         hM_TrkrhitADC_noncorr->Scale(1. / hM_TrkrhitADC_noncorr->Integral(-1, -1));
0520         hM_TrkrhitADC_pm4noncorr->Scale(1. / hM_TrkrhitADC_pm4noncorr->Integral(-1, -1));
0521         float binmax_adc = std::max(hM_TrkrhitADC_noncorr->GetMaximum(), hM_TrkrhitADC_pm4noncorr->GetMaximum());
0522
0523         TCanvas *c_adc = new TCanvas("c_adc", "c_adc", 800, 700);
0524         c_adc->cd();
0525         hM_TrkrhitADC_noncorr->GetXaxis()->SetTitle("TrkrHit ADC");
0526         hM_TrkrhitADC_noncorr->GetYaxis()->SetTitle("Normalized entries");
0527         hM_TrkrhitADC_noncorr->GetXaxis()->SetTitleOffset(1.1);
0528         hM_TrkrhitADC_noncorr->GetYaxis()->SetTitleOffset(1.4);
0529         hM_TrkrhitADC_noncorr->GetYaxis()->SetRangeUser(0, binmax_adc * 1.5);
0530         hM_TrkrhitADC_noncorr->SetLineColor(kTVibBlue);
0531         hM_TrkrhitADC_noncorr->SetLineWidth(2);
0532         hM_TrkrhitADC_noncorr->SetMarkerColor(kTVibBlue);
0533         hM_TrkrhitADC_noncorr->SetMarkerStyle(20);
0534         hM_TrkrhitADC_noncorr->SetMarkerSize(1);
0535         hM_TrkrhitADC_pm4noncorr->SetLineColor(kTVibRed);
0536         hM_TrkrhitADC_pm4noncorr->SetLineWidth(2);
0537         hM_TrkrhitADC_pm4noncorr->SetMarkerColor(kTVibRed);
0538         hM_TrkrhitADC_pm4noncorr->SetMarkerStyle(20);
0539         hM_TrkrhitADC_pm4noncorr->SetMarkerSize(1);
0540         hM_TrkrhitADC_noncorr->Draw("PE");
0541         hM_TrkrhitADC_pm4noncorr->Draw("PE same");
0542         TLegend *leg = new TLegend(0.3, 0.8, 0.7, 0.9);
0543         leg->AddEntry(hM_TrkrhitADC_noncorr, "Events of interest (Is MB)", "pel");
0544         leg->AddEntry(hM_TrkrhitADC_pm4noncorr, "[-4, -1] and [+1, +4] events (Is MB)", "pel");
0545         leg->SetBorderSize(0);
0546         leg->SetFillStyle(0);
0547         leg->SetTextSize(0.04);
0548         leg->Draw();
0549         c_adc->SaveAs((plotdir + "TrkrHitADC_noncorr.pdf").c_str());
0550         c_adc->SaveAs((plotdir + "TrkrHitADC_noncorr.png").c_str());
0551
0552         // hM_HitTimeBucket_Layer_noncorrnext->Scale(1. / hM_HitTimeBucket_Layer_noncorrnext->Integral(-1, -1, -1, -1));
0553         float binmax_timebucket = 0;
0554         for (int l = 0; l < 4; l++)
0555         {
0556             // list_hM_HitTimeBucket_Layer_noncorrnext[l]->Scale(1. / list_hM_HitTimeBucket_Layer_noncorrnext[l]->Integral(-1, -1));
0557             if (list_hM_HitTimeBucket_Layer_noncorrnext[l]->GetMaximum() > binmax_timebucket)
0558             {
0559                 binmax_timebucket = list_hM_HitTimeBucket_Layer_noncorrnext[l]->GetMaximum();
0560             }
0561         }
0562         TCanvas *c_timebucketlayer = new TCanvas("c_timebucketlayer", "c_timebucketlayer", 800, 700);
0563         TLegend *leg_timebucketlayer = new TLegend(0.6, 0.7, 0.85, 0.9);
0564         leg_timebucketlayer->SetHeader("Next event of EOI");
0565         c_timebucketlayer->cd();
0566         c_timebucketlayer->SetTopMargin(0.07);
0567         // hM_HitTimeBucket_Layer_noncorrnext->GetXaxis()->SetTitle("Time bucket");
0568         // hM_HitTimeBucket_Layer_noncorrnext->GetYaxis()->SetTitle("Layer");
0569         // hM_HitTimeBucket_Layer_noncorrnext->GetXaxis()->SetTitleOffset(1.3);
0570         // hM_HitTimeBucket_Layer_noncorrnext->GetYaxis()->SetTitleOffset(1.4);
0571         // hM_HitTimeBucket_Layer_noncorrnext->GetYaxis()->SetNdivisions(4);
0572         // hM_HitTimeBucket_Layer_noncorrnext->Draw("colz");
0573         list_hM_HitTimeBucket_Layer_noncorrnext[0]->GetXaxis()->SetTitle("Time bucket");
0574         list_hM_HitTimeBucket_Layer_noncorrnext[0]->GetYaxis()->SetTitle("Entries");
0575         list_hM_HitTimeBucket_Layer_noncorrnext[0]->GetXaxis()->SetTitleOffset(1.3);
0576         list_hM_HitTimeBucket_Layer_noncorrnext[0]->GetYaxis()->SetTitleOffset(1.4);
0577         list_hM_HitTimeBucket_Layer_noncorrnext[0]->GetYaxis()->SetRangeUser(0, binmax_timebucket * 1.3);
0578         list_hM_HitTimeBucket_Layer_noncorrnext[0]->SetLineColor(colorlist[0]);
0579         list_hM_HitTimeBucket_Layer_noncorrnext[0]->SetLineWidth(1);
0580         list_hM_HitTimeBucket_Layer_noncorrnext[0]->SetMarkerColor(colorlist[0]);
0581         list_hM_HitTimeBucket_Layer_noncorrnext[0]->SetMarkerStyle(20);
0582         list_hM_HitTimeBucket_Layer_noncorrnext[0]->SetMarkerSize(0.1);
0583         list_hM_HitTimeBucket_Layer_noncorrnext[0]->Draw("E hist");
0584         leg_timebucketlayer->AddEntry(list_hM_HitTimeBucket_Layer_noncorrnext[0], "Layer 3", "pel");
0585         for (int l = 1; l < 4; l++)
0586         {
0587             // TLine *line = new TLine(-0.5, 2.5+i, 130.5, 2.5+i);
0588             // line->SetLineColor(kBlack);
0589             // line->SetLineWidth(1);
0590             // line->Draw();
0591             list_hM_HitTimeBucket_Layer_noncorrnext[l]->SetLineWidth(1);
0592             list_hM_HitTimeBucket_Layer_noncorrnext[l]->SetLineColor(colorlist[l]);
0593             list_hM_HitTimeBucket_Layer_noncorrnext[l]->SetMarkerSize(0.1);
0594             list_hM_HitTimeBucket_Layer_noncorrnext[l]->SetMarkerStyle(20);
0595             list_hM_HitTimeBucket_Layer_noncorrnext[l]->SetMarkerColor(colorlist[l]);
0596             list_hM_HitTimeBucket_Layer_noncorrnext[l]->Draw("E hist same");
0597             leg_timebucketlayer->AddEntry(list_hM_HitTimeBucket_Layer_noncorrnext[l], Form("Layer %d", l + 3), "pel");
0598         }
0599         leg_timebucketlayer->SetBorderSize(0);
0600         leg_timebucketlayer->SetFillStyle(0);
0601         leg_timebucketlayer->SetTextSize(0.04);
0602         leg_timebucketlayer->Draw();
0603         c_timebucketlayer->RedrawAxis();
0604         c_timebucketlayer->SaveAs((plotdir + "HitTimebucketLayer_evtsel_noncorrnext.pdf").c_str());
0605         c_timebucketlayer->SaveAs((plotdir + "HitTimebucketLayer_evtsel_noncorrnext.png").c_str());
0606     }
0607
0608     outtxtfile.close();
0609
0610     TCanvas *c_bcodiff = new TCanvas("c_bcodiff", "c_bcodiff", 800, 700);
0611     c_bcodiff->cd();
0612     hM_BCOdiff_noncorrnext->GetXaxis()->SetTitle("BCO_{next} - BCO_{of interest}");
0613     hM_BCOdiff_noncorrnext->GetYaxis()->SetTitle("Entries");
0614     // Draw overflow bin
0615     hM_BCOdiff_noncorrnext->GetXaxis()->SetRange(1, hM_BCOdiff_noncorrnext->GetNbinsX());
0616     hM_BCOdiff_noncorrnext->GetXaxis()->SetTitleOffset(1.3);
0617     hM_BCOdiff_noncorrnext->GetYaxis()->SetTitleOffset(1.4);
0618     hM_BCOdiff_noncorrnext->SetLineColor(kBlack);
0619     hM_BCOdiff_noncorrnext->SetLineWidth(2);
0620     hM_BCOdiff_noncorrnext->Draw("hist");
0621     c_bcodiff->SaveAs((plotdir + "BCOdiff_noncorrnext.pdf").c_str());
0622     c_bcodiff->SaveAs((plotdir + "BCOdiff_noncorrnext.png").c_str());
0623 }