analysis_INTT/plot/quickdraw_mbdvtx.C

0001 #include <filesystem>
0002
0003 // std::string ntpdir = "/sphenix/tg/tg01/hf/hjheng/ppg02/dst/Run54280_combinedNtuple/";
0004 std::string ntpdir = "/sphenix/tg/tg01/commissioning/INTT/work/cwshih/seflgendata/run_54280/From_official_INTTRAWHIT_DST_ColdHotChannelMask/completed/";
0005 std::string vtxmapdir = "/sphenix/user/hjheng/sPHENIXRepo/analysis/dNdEta_Run2023/analysis_INTT/minitree/VtxEvtMap_Data_CombinedNtuple_Run54280_HotChannel_BCO/";
0006 std::string plotdir = "./RecoPV_ana/Run54280_trigcomp";
0007
0008 std::map<uint64_t, vector<float>> EvtVtx_map_inttbco(const char *vtxfname)
0009 {
0010     std::map<uint64_t, vector<float>> EvtVtx_map;
0011
0012     TFile *f = new TFile(vtxfname, "READ");
0013     TTree *t = (TTree *)f->Get("minitree");
0014     uint64_t intt_bco;
0015     float PV_x, PV_y, PV_z;
0016     t->SetBranchAddress("INTT_BCO", &intt_bco);
0017     t->SetBranchAddress("PV_x", &PV_x);
0018     t->SetBranchAddress("PV_y", &PV_y);
0019     t->SetBranchAddress("PV_z", &PV_z);
0020     for (int ev = 0; ev < t->GetEntriesFast(); ev++)
0021     {
0022         t->GetEntry(ev);
0023         EvtVtx_map[intt_bco] = {PV_x, PV_y, PV_z};
0024     }
0025
0026     return EvtVtx_map;
0027 }
0028
0029 void quickdraw_mbdvtx()
0030 {
0031     // mkdir if not exist
0032     system(Form("mkdir -p %s", plotdir.c_str()));
0033
0034     TH1F *hM_MbdZvtx_trig12 = new TH1F("hM_MbdZvtx_trig12", "M_MbdZvtx_trig12", 160, -40, 40);
0035     TH1F *hM_MbdZvtx_trig12_InttZvtxle10 = new TH1F("hM_MbdZvtx_trig12_InttZvtxle10", "M_MbdZvtx_trig12_InttZvtxle10", 160, -40, 40);
0036     TH1F *hM_MbdZvtx_trig13 = new TH1F("hM_MbdZvtx_trig13", "hM_MbdZvtx_trig13", 160, -40, 40);
0037     TH1F *hM_MbdZvtx_trig13_InttZvtxle10 = new TH1F("hM_MbdZvtx+trig13_InttZvtxle10", "hM_MbdZvtx+trig13_InttZvtxle10", 160, -40, 40);
0038
0039     // loop over all files in the directory, sort by filename
0040     int nfilestorun = 10;
0041     for (int i = 0; i < nfilestorun; i++)
0042     {
0043         std::string formatidx = std::string(TString::Format("%05d", i).Data());
0044         std::string filename = Form("%s/Trkr_DST_physics_intt-54280_HotChannel_BCOCut_CaloCombined_%s.root", ntpdir.c_str(), formatidx.c_str());
0045         // check if the file exists
0046         std::filesystem::path entry(filename);
0047         if (!std::filesystem::exists(entry))
0048             continue;
0049
0050         std::cout << "Processing " << filename << std::endl;
0051
0052         std::map<uint64_t, vector<float>> evtvtxmap = EvtVtx_map_inttbco(Form("%s/minitree_%s.root", vtxmapdir.c_str(), formatidx.c_str()));
0053
0054         TFile *f = new TFile(filename.c_str(), "READ");
0055         TTree *t = (TTree *)f->Get("EventTree");
0056         uint64_t INTT_BCO;
0057         std::vector<int> *firedTriggers = 0;
0058         float MBD_z_vtx, MBD_south_charge_sum, MBD_north_charge_sum;
0059         t->SetBranchAddress("LiveTrigger_Vec", &firedTriggers);
0060         t->SetBranchAddress("INTT_BCO", &INTT_BCO);
0061         t->SetBranchAddress("MBD_z_vtx", &MBD_z_vtx);
0062         t->SetBranchAddress("MBD_south_charge_sum", &MBD_south_charge_sum);
0063         t->SetBranchAddress("MBD_north_charge_sum", &MBD_north_charge_sum);
0064         for (int ev = 0; ev < t->GetEntriesFast(); ev++)
0065         {
0066             t->GetEntry(ev);
0067
0068             // try to get the vertex from the map. If not found, set to nan
0069             float INTT_z_vtx = std::numeric_limits<float>::quiet_NaN();
0070             if (evtvtxmap.find(INTT_BCO) != evtvtxmap.end())
0071                 INTT_z_vtx = evtvtxmap[INTT_BCO][2];
0072
0073             bool validInttVtx = (INTT_z_vtx == INTT_z_vtx);
0074             bool validMbdVtx = (MBD_z_vtx == MBD_z_vtx);
0075             bool MbdNScharge = (MBD_north_charge_sum > 0 || MBD_south_charge_sum > 0);
0076             bool firedTrig12 = (std::find(firedTriggers->begin(), firedTriggers->end(), 12) != firedTriggers->end());
0077             bool firedTrig13 = (std::find(firedTriggers->begin(), firedTriggers->end(), 13) != firedTriggers->end());
0078
0079             if (!validMbdVtx || !MbdNScharge || !validInttVtx)
0080                 continue;
0081
0082             if (firedTrig12)
0083             {
0084                 hM_MbdZvtx_trig12->Fill(MBD_z_vtx);
0085                 if (fabs(INTT_z_vtx) < 10)
0086                     hM_MbdZvtx_trig12_InttZvtxle10->Fill(MBD_z_vtx);
0087             }
0088
0089             if (firedTrig13)
0090             {
0091                 hM_MbdZvtx_trig13->Fill(MBD_z_vtx);
0092                 if (fabs(INTT_z_vtx) < 10)
0093                     hM_MbdZvtx_trig13_InttZvtxle10->Fill(MBD_z_vtx);
0094             }
0095         }
0096
0097         f->Close();
0098     }
0099
0100     TCanvas *c = new TCanvas("c", "c", 800, 700);
0101     c->cd();
0102     hM_MbdZvtx_trig13->GetXaxis()->SetTitle("MBD Z vertex [cm]");
0103     hM_MbdZvtx_trig13->GetYaxis()->SetTitle("Counts");
0104     hM_MbdZvtx_trig13->GetYaxis()->SetTitleOffset(1.5);
0105     hM_MbdZvtx_trig13->GetYaxis()->SetRangeUser(0, hM_MbdZvtx_trig13->GetMaximum() * 1.5);
0106     hM_MbdZvtx_trig13->SetLineWidth(2);
0107     hM_MbdZvtx_trig13->SetLineColor(kTBriRed);
0108     hM_MbdZvtx_trig13->SetMarkerColor(kTBriRed);
0109     hM_MbdZvtx_trig13->SetMarkerSize(0.7);
0110     hM_MbdZvtx_trig13->Draw("PE");
0111     hM_MbdZvtx_trig12->SetLineWidth(2);
0112     hM_MbdZvtx_trig12->SetLineColor(kTBriBlue);
0113     hM_MbdZvtx_trig12->SetMarkerColor(kTBriBlue);
0114     hM_MbdZvtx_trig12->SetMarkerSize(0.7);
0115     hM_MbdZvtx_trig12->Draw("PE same");
0116     hM_MbdZvtx_trig12_InttZvtxle10->SetLineWidth(2);
0117     hM_MbdZvtx_trig12_InttZvtxle10->SetLineColor(kTBriCyan);
0118     hM_MbdZvtx_trig12_InttZvtxle10->SetMarkerColor(kTBriCyan);
0119     hM_MbdZvtx_trig12_InttZvtxle10->SetMarkerSize(0.7);
0120     hM_MbdZvtx_trig12_InttZvtxle10->Draw("PE same");
0121     hM_MbdZvtx_trig13_InttZvtxle10->SetLineWidth(2);
0122     hM_MbdZvtx_trig13_InttZvtxle10->SetLineColor(kTBriPurple);
0123     hM_MbdZvtx_trig13_InttZvtxle10->SetMarkerColor(kTBriPurple);
0124     hM_MbdZvtx_trig13_InttZvtxle10->SetMarkerSize(0.7);
0125     hM_MbdZvtx_trig13_InttZvtxle10->Draw("PE same");
0126     TLegend *leg = new TLegend(0.2, 0.72, 0.5, 0.9);
0127     leg->AddEntry(hM_MbdZvtx_trig13, "Fired trigger 13 (MBD S&N#geq2, |vtx|<30cm)", "pel");
0128     leg->AddEntry(hM_MbdZvtx_trig12, "Fired trigger 12 (MBD S&N#geq2, |vtx|<10cm)", "pel");
0129     leg->AddEntry(hM_MbdZvtx_trig13_InttZvtxle10, "Fired trigger 13, |INTT Z vtx|<10cm", "pel");
0130     leg->AddEntry(hM_MbdZvtx_trig12_InttZvtxle10, "Fired trigger 12, |INTT Z vtx|<10cm", "pel");
0131     leg->SetTextSize(0.04);
0132     leg->SetBorderSize(0);
0133     leg->SetFillStyle(0);
0134     leg->Draw();
0135     c->RedrawAxis();
0136     c->SaveAs(Form("%s/quickdraw_mbdvtx.png", plotdir.c_str()));
0137     c->SaveAs(Form("%s/quickdraw_mbdvtx.pdf", plotdir.c_str()));
0138 }