analysis_INTT/plot/makePaperPlot.C

0001 #include "./plotutil.h"
0002
0003 float TickSize = 0.03;
0004 float AxisTitleSize = 0.06;
0005 float AxisLabelSize = 0.05;
0006 float LeftMargin = 0.15;
0007 float RightMargin = 0.05;
0008 float TopMargin = 0.08;
0009 float BottomMargin = 0.13;
0010 float textscale_pad1 = 1.2;
0011 float textscale_pad2 = 2.6;
0012
0013 std::vector<std::string> histfnames = {"./hists/Data_Run54280_20250210_ProdA2024/dRcut0p5_NominalVtxZ_RandomClusSet0_clusAdcCutSet0_clusPhiSizeCutSet0/hists_merged.root", //
0014                                        "./hists/Sim_HIJING_MDC2_ana472_20250307/dRcut0p5_NominalVtxZ_RandomClusSet0_clusAdcCutSet0_clusPhiSizeCutSet0/hists_merged.root",  //
0015                                        "./hists/Sim_EPOS_MDC2_ana472_20250310/dRcut0p5_NominalVtxZ_RandomClusSet0_clusAdcCutSet0_clusPhiSizeCutSet0/hists_merged.root",    //
0016                                        "./hists/Sim_AMPT_MDC2_ana472_20250310/dRcut0p5_NominalVtxZ_RandomClusSet0_clusAdcCutSet0_clusPhiSizeCutSet0/hists_merged.root"};
0017
0018 std::vector<TH1 *> hists;
0019
0020 std::vector<std::string> histlegend = {"Data", "HIJING", "EPOS", "AMPT"};
0021
0022 std::vector<int> histlinestyle = {1, 1, 6, 7};
0023
0024 std::vector<std::string> histcolor = {"#09122C", "#99cc99", "#6699cc", "#f2777a"};
0025
0026 // std::string prelimtext = "Internal";
0027 std::string prelimtext = "";
0028
0029 void drawcompare_ratio(std::vector<TH1 *> hists, //
0030                        bool logy,                //
0031                        float ymaxscale,          //
0032                        std::string xtitle,       //
0033                        std::string ytitle,       //
0034                        std::string legtitle,     //
0035                        std::string plotName      //
0036 )
0037 {
0038     float ymin = 1E10, ymax = -1E10;
0039     float ratio_min = 999., ratio_max = -1.;
0040     // normalize the histogram and take ratio
0041     std::vector<TH1 *> v_ratio;
0042     v_ratio.clear();
0043
0044     for (size_t i = 0; i < hists.size(); ++i)
0045     {
0046         if (hists[i])
0047         {
0048             hists[i]->Scale(1. / hists[i]->Integral(-1, -1), "width");
0049
0050             if (hists[i]->GetMinimum(0) < ymin)
0051                 ymin = hists[i]->GetMinimum(0);
0052
0053             if (hists[i]->GetMaximum() > ymax)
0054                 ymax = hists[i]->GetMaximum();
0055
0056             hists[i]->SetLineColor(TColor::GetColor(histcolor[i].c_str()));
0057             hists[i]->SetLineWidth(2);
0058             hists[i]->SetLineStyle(histlinestyle[i]);
0059
0060             if (i == 0)
0061             {
0062                 hists[i]->SetMarkerStyle(20);
0063                 hists[i]->SetMarkerSize(0.8);
0064             }
0065             else
0066             {
0067                 hists[i]->SetMarkerSize(0);
0068                 // calculate ratio to the first histogram
0069                 TH1 *h_ratio = (TH1 *)hists[i]->Clone(Form("h_ratio_%s", histlegend[i].c_str()));
0070                 h_ratio->Divide(hists[0]);
0071                 v_ratio.push_back(h_ratio);
0072                 h_ratio->SetLineColor(TColor::GetColor(histcolor[i].c_str()));
0073                 h_ratio->SetLineWidth(2);
0074                 h_ratio->SetLineStyle(histlinestyle[i]);
0075                 h_ratio->SetMarkerSize(0);
0076
0077                 if (h_ratio->GetMinimum(0) < ratio_min)
0078                     ratio_min = h_ratio->GetMinimum(0);
0079
0080                 if (h_ratio->GetMaximum() > ratio_max)
0081                     ratio_max = h_ratio->GetMaximum();
0082             }
0083         }
0084     }
0085
0086     TCanvas *c = new TCanvas("c", "c", 800, 700);
0087     TPad *pad1 = new TPad("pad1", "pad1", 0, 0.3, 1, 1);
0088     pad1->SetBottomMargin(0);
0089     pad1->SetTopMargin(0.08);
0090     pad1->SetRightMargin(RightMargin);
0091     pad1->Draw();
0092     TPad *pad2 = new TPad("pad2", "pad2", 0, 0, 1, 0.3);
0093     pad2->SetRightMargin(RightMargin);
0094     pad2->SetTopMargin(0);
0095     pad2->SetBottomMargin(0.4);
0096     pad2->Draw();
0097     pad1->cd();
0098     pad1->SetLogy(logy);
0099     hists[0]->GetYaxis()->SetTitle(ytitle.c_str());
0100     hists[0]->GetYaxis()->SetTitleOffset(1.);
0101     hists[0]->GetYaxis()->SetTitleSize(AxisTitleSize * textscale_pad1);
0102     hists[0]->GetYaxis()->SetLabelSize(AxisLabelSize * textscale_pad1);
0103     hists[0]->GetYaxis()->SetRangeUser(ymin * 0.7, ymax * ymaxscale);
0104     // hists[0]->GetYaxis()->SetMoreLogLabels();
0105     hists[0]->GetXaxis()->SetTitleSize(0);
0106     hists[0]->GetXaxis()->SetLabelSize(0);
0107     hists[0]->Draw("PE");
0108     for (size_t i = 1; i < hists.size(); ++i)
0109     {
0110         if (hists[i])
0111         {
0112             hists[i]->Draw("hist same");
0113         }
0114     }
0115
0116     TLegend *sphnxleg = new TLegend(1 - pad1->GetRightMargin() - 0.47, //
0117                                     1 - pad1->GetTopMargin() - 0.23,    //
0118                                     1 - pad1->GetRightMargin() - 0.1,  //
0119                                     1 - pad1->GetTopMargin() - 0.07    //
0120     );
0121     sphnxleg->SetTextAlign(kHAlignLeft + kVAlignCenter);
0122     sphnxleg->SetTextSize(0.06);
0123     sphnxleg->SetFillStyle(0);
0124     sphnxleg->AddEntry("", Form("#it{#bf{sPHENIX}} %s", prelimtext.c_str()), "");
0125     sphnxleg->AddEntry("", "Au+Au #sqrt{s_{NN}}=200 GeV", "");
0126     sphnxleg->Draw();
0127     TLegend *histleg = new TLegend(1 - pad1->GetRightMargin() - 0.382, //
0128                                    1 - pad1->GetTopMargin() - 0.47,   //
0129                                    1 - pad1->GetRightMargin(),        //
0130                                    1 - pad1->GetTopMargin() - 0.25    //
0131     );
0132     histleg->SetHeader(Form("%s", legtitle.c_str()));
0133     histleg->SetNColumns(2);
0134     histleg->AddEntry(hists[0], histlegend[0].c_str(), "pe");
0135     for (size_t i = 1; i < hists.size(); ++i)
0136     {
0137         histleg->AddEntry(hists[i], histlegend[i].c_str(), "l");
0138     }
0139     histleg->SetBorderSize(0);
0140     histleg->SetFillStyle(0);
0141     histleg->SetTextSize(0.055);
0142     histleg->Draw();
0143
0144     pad2->cd();
0145     // pad2->SetLogy(logy);
0146     for (size_t i = 0; i < v_ratio.size(); ++i)
0147     {
0148         if (i == 0)
0149         {
0150             v_ratio[i]->GetYaxis()->SetTitle("Model/data");
0151             v_ratio[i]->GetYaxis()->SetTitleOffset(0.5);
0152             v_ratio[i]->GetYaxis()->SetTitleSize(AxisTitleSize * textscale_pad2);
0153             v_ratio[i]->GetYaxis()->SetLabelSize(AxisLabelSize * textscale_pad2);
0154             v_ratio[i]->GetYaxis()->SetRangeUser(0.501, 1.499);
0155             v_ratio[i]->GetYaxis()->SetNdivisions(5, 3, 0);
0156             v_ratio[i]->GetXaxis()->SetTitle(xtitle.c_str());
0157             v_ratio[i]->GetXaxis()->SetTitleSize(AxisTitleSize * textscale_pad2);
0158             v_ratio[i]->GetXaxis()->SetLabelSize(AxisLabelSize * textscale_pad2);
0159             v_ratio[i]->GetXaxis()->SetTitleOffset(1.25);
0160             v_ratio[i]->Draw("hist l");
0161         }
0162         else
0163         {
0164             v_ratio[i]->Draw("hist l same");
0165         }
0166     }
0167
0168     TLine *line = new TLine(v_ratio[0]->GetXaxis()->GetXmin(), 1, v_ratio[0]->GetXaxis()->GetXmax(), 1);
0169     line->SetLineColor(kBlack);
0170     line->SetLineWidth(2);
0171     line->SetLineStyle(2);
0172     line->Draw("same");
0173     c->RedrawAxis();
0174     c->SaveAs(Form("%s.pdf", plotName.c_str()));
0175     c->SaveAs(Form("%s.png", plotName.c_str()));
0176 }
0177
0178 void makePaperPlot()
0179 {
0180     std::string plotdir = "./PaperPlot";
0181
0182     for (auto &fname : histfnames)
0183     {
0184         TH1 *h = getHist<TH1>(fname, "hM_dR_reco_Centrality0to70_MBDAsymLe0p75_VtxZm20to20");
0185         if (h)
0186         {
0187             hists.push_back(h);
0188         }
0189     }
0190
0191     drawcompare_ratio(hists, true, 2, "Tracklet #DeltaR", "#frac{1}{N} #frac{dN}{d(#DeltaR)}", "Centrality 0-70%", Form("%s/Paperplot_RecotklDeltaR", plotdir.c_str()));
0192 }