Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-09 08:12:15 

0001 #include "vtxZDist_comp.h"
0002 
0003 vtxZDist_comp::vtxZDist_comp(
0004         std::string sPH_labeling_in,
0005         std::string MC_labeling_in,
0006         std::vector<std::pair<std::string, std::pair<std::string, std::string>>> data_input_directory_pair_vec_in, // note: directory, and description
0007         std::vector<std::pair<std::string, std::pair<std::string, std::string>>> MC_input_directory_pair_vec_in, // note: directory, and description
0008         std::map<std::string, std::vector<std::tuple<double,double,std::string>>> labelling_vec_map_in, // note: starting XY, and description 
0009         std::string output_directory_in,
0010 
0011         bool WithVtxZReWeighting_in // note: MC with reweighting
0012 ):
0013     sPH_labeling(sPH_labeling_in),
0014     MC_labeling(MC_labeling_in),
0015     data_input_directory_pair_vec(data_input_directory_pair_vec_in),
0016     MC_input_directory_pair_vec(MC_input_directory_pair_vec_in),
0017     labelling_vec_map(labelling_vec_map_in),
0018     output_directory(output_directory_in),
0019 
0020     WithVtxZReWeighting(WithVtxZReWeighting_in),
0021 
0022     logy_max_factor(1000),
0023     lineary_max_factor(1.5)
0024 {
0025     system(Form("if [ ! -d %s ]; then mkdir -p %s; fi", output_directory.c_str(), output_directory.c_str()));    
0026 
0027     for (auto &pair : data_input_directory_pair_vec) // note : {direcotry -> <description, folder_name>}
0028     {
0029         std::cout<<"data_input_directory_pair_vec: "<<pair.first<<" "<<pair.second.first<<" "<<pair.second.second<<std::endl;
0030         system(Form("if [ ! -d %s/Data_%s ]; then mkdir -p %s/Data_%s; fi", output_directory.c_str(), pair.second.second.c_str(), output_directory.c_str(), pair.second.second.c_str()));    
0031     }
0032 
0033     for (auto &pair : MC_input_directory_pair_vec)
0034     {
0035         std::cout<<"MC_input_directory_pair_vec: "<<pair.first<<" "<<pair.second.first<<" "<<pair.second.second<<std::endl;
0036         system(Form("if [ ! -d %s/MC_%s ]; then mkdir -p %s/MC_%s; fi", output_directory.c_str(), pair.second.second.c_str(), output_directory.c_str(), pair.second.second.c_str()));    
0037     }
0038 
0039     data_input_directory_map.clear();
0040     for (auto &pair : data_input_directory_pair_vec)
0041     {
0042         data_input_directory_map.insert(
0043             std::make_pair(
0044                 pair.first,
0045                 pair.second
0046             )
0047         );
0048     }
0049 
0050     MC_input_directory_map.clear();
0051     for (auto &pair : MC_input_directory_pair_vec)
0052     {
0053         MC_input_directory_map.insert(
0054             std::make_pair(
0055                 pair.first,
0056                 pair.second
0057             )
0058         );
0059     }
0060 
0061     SetsPhenixStyle();
0062     TGaxis::SetMaxDigits(4);
0063     c1 = new TCanvas("c1","c1",950,800);
0064 
0065     ltx = new TLatex();
0066     ltx->SetNDC();
0067     ltx->SetTextSize(0.045);
0068     ltx->SetTextAlign(31);
0069 
0070     All_leg = new TLegend(0.4, 0.82, 0.8, 0.88);
0071     All_leg -> SetBorderSize(0);
0072     All_leg -> SetTextSize(0.025);
0073 
0074     All_leg_long = new TLegend(0.7, 0.5, 0.8, 0.88);
0075     All_leg_long -> SetBorderSize(0);
0076 
0077     sub_ltx = new TLatex();
0078     sub_ltx->SetNDC();
0079     sub_ltx->SetTextSize(0.03);
0080     // ltx->SetTextAlign(31);
0081 
0082     PrepareInputFiles();
0083 
0084     nCentrality_bin = Constants::centrality_edges.size() - 1;
0085 }
0086 
0087 std::string vtxZDist_comp::string_convertor(std::string input_string)
0088 {
0089     std::string output_string;
0090 
0091     for (int i = 0; i < input_string.size(); i++)
0092     {
0093         if (input_string[i] == ' ')
0094         {
0095             output_string += "_";
0096         }
0097         else if (input_string[i] == '.')
0098         {
0099             output_string += "p";
0100         }
0101         else if (input_string[i] == ',')
0102         {
0103             output_string += "_";
0104         }
0105         else
0106         {
0107             output_string += input_string[i];
0108         }
0109     }
0110 
0111     // std::cout<<"before: "<<input_string<<" after: "<<output_string<<std::endl;
0112 
0113     return output_string;
0114 }
0115 
0116 void vtxZDist_comp::PrepareInputFiles()
0117 {
0118     data_h1D_map_map.clear();
0119     data_h2D_map_map.clear();
0120     MC_h1D_map_map.clear();
0121     MC_h2D_map_map.clear();
0122 
0123     // note : for data ----------------------------------------------------------------------------------------------------
0124     for (auto &pair : data_input_directory_pair_vec) // note : <file_name, description>
0125     {
0126         file_data_vec.push_back(
0127             TFile::Open(Form("%s",pair.first.c_str()))
0128         );
0129 
0130         if (file_data_vec.back() == nullptr)
0131         {
0132             std::cout << "Error: file_in can not be opened" << std::endl;
0133             exit(1);
0134         }
0135 
0136         // note : initialize the vector
0137         data_h2D_map_map.insert( // note : file_name, hist_name, hist
0138             std::make_pair(
0139                 pair.first,
0140                 std::map<std::string, TH2D*>{} // note : also map<string,TH2D*>
0141             )
0142         ); // note : the size is number of files
0143         data_h1D_map_map.insert(
0144             std::make_pair(
0145                 pair.first,
0146                 std::map<std::string, TH1D*>{} // note : also map<string,TH1D*>
0147             )
0148         ); // note : the size is number of files
0149 
0150         for (TObject* keyAsObj : *file_data_vec.back()->GetListOfKeys())
0151         {
0152             auto key = dynamic_cast<TKey*>(keyAsObj);
0153             std::string hist_name  = key->GetName();
0154             std::string class_name = key->GetClassName();
0155 
0156             if (class_name == "TH2D")
0157             {
0158                 data_h2D_map_map[pair.first].insert(
0159                     std::make_pair(
0160                         hist_name.c_str(),
0161                         (TH2D*) file_data_vec.back() -> Get( hist_name.c_str() )
0162                     )
0163                 );
0164             }
0165             else if (class_name == "TH1D")
0166             {
0167                 data_h1D_map_map[pair.first].insert(
0168                     std::make_pair(
0169                         hist_name.c_str(),
0170                         (TH1D*) file_data_vec.back() -> Get( hist_name.c_str() )
0171                     )
0172                 );
0173 
0174                 data_h1D_map_map[pair.first][hist_name] -> SetLineColor(TColor::GetColor(color_code[data_h1D_map_map.size()-1].c_str()));
0175                 data_h1D_map_map[pair.first][hist_name] -> SetMarkerColor(TColor::GetColor(color_code[data_h1D_map_map.size()-1].c_str()));
0176                 data_h1D_map_map[pair.first][hist_name] -> SetMarkerStyle(20);
0177                 data_h1D_map_map[pair.first][hist_name] -> SetMarkerSize(0.8);
0178                 data_h1D_map_map[pair.first][hist_name] -> Sumw2(true);
0179                 data_h1D_map_map[pair.first][hist_name] -> Scale(1.0/data_h1D_map_map[pair.first][hist_name]->Integral());
0180 
0181                 std::string Yaxis_title = data_h1D_map_map[pair.first][hist_name] -> GetYaxis() -> GetTitle();
0182                 Yaxis_title += "(A.U.)";
0183                 data_h1D_map_map[pair.first][hist_name] -> GetYaxis() -> SetTitle(Yaxis_title.c_str());
0184             }
0185         }
0186     }
0187 
0188     // note : for MC ----------------------------------------------------------------------------------------------------
0189     for (auto &pair : MC_input_directory_pair_vec)
0190     {
0191         file_MC_vec.push_back(
0192             TFile::Open(Form("%s",pair.first.c_str()))
0193         );
0194 
0195         if (file_MC_vec.back() == nullptr)
0196         {
0197             std::cout << "Error: file_in can not be opened" << std::endl;
0198             exit(1);
0199         }
0200 
0201         // note : initialize the vector
0202         MC_h2D_map_map.insert(
0203             std::make_pair(
0204                 pair.first,
0205                 std::map<std::string, TH2D*>{}
0206             )
0207         );
0208         MC_h1D_map_map.insert(
0209             std::make_pair(
0210                 pair.first,
0211                 std::map<std::string, TH1D*>{}
0212             )
0213         );
0214 
0215         for (TObject* keyAsObj : *file_MC_vec.back()->GetListOfKeys())
0216         {
0217             auto key = dynamic_cast<TKey*>(keyAsObj);
0218             std::string hist_name  = key->GetName();
0219             std::string class_name = key->GetClassName();
0220 
0221             if (class_name == "TH2D")
0222             {
0223                 MC_h2D_map_map[pair.first].insert(
0224                     std::make_pair(
0225                         hist_name.c_str(),
0226                         (TH2D*) file_MC_vec.back() -> Get( hist_name.c_str() )
0227                     )
0228                 );
0229             }
0230             else if (class_name == "TH1D")
0231             {
0232                 MC_h1D_map_map[pair.first].insert(
0233                     std::make_pair(
0234                         hist_name.c_str(),
0235                         (TH1D*) file_MC_vec.back() -> Get( hist_name.c_str() )
0236                     )
0237                 );
0238 
0239                 MC_h1D_map_map[pair.first][hist_name] -> SetLineColor(TColor::GetColor( color_code[MC_h1D_map_map.size() - 1 + data_h1D_map_map.size()].c_str() ));
0240                 MC_h1D_map_map[pair.first][hist_name] -> SetMarkerColor(TColor::GetColor( color_code[MC_h1D_map_map.size() - 1 + data_h1D_map_map.size()].c_str() ));
0241                 MC_h1D_map_map[pair.first][hist_name] -> SetMarkerStyle(20);
0242                 MC_h1D_map_map[pair.first][hist_name] -> Sumw2(true);
0243                 MC_h1D_map_map[pair.first][hist_name] -> Scale(1.0/MC_h1D_map_map[pair.first][hist_name]->Integral());
0244 
0245                 std::string Yaxis_title = MC_h1D_map_map[pair.first][hist_name] -> GetYaxis() -> GetTitle();
0246                 Yaxis_title += "(A.U.)";
0247                 MC_h1D_map_map[pair.first][hist_name] -> GetYaxis() -> SetTitle(Yaxis_title.c_str());
0248             }
0249         }
0250     }
0251 
0252 }
0253 
0254 void vtxZDist_comp::MakePlots(std::string draw_method, bool isData)
0255 {
0256     std::string sub_label_text = (isData) ? "Data" : "MC";
0257     
0258     std::map<std::string, std::pair<std::string, std::string>> input_directory_map;
0259     std::vector<std::string> h1_only_logy_plot_vec;
0260     std::map<std::string, std::map<std::string, TH1D*>> *h1D_map_map ;
0261     std::map<std::string, std::map<std::string, TH2D*>> *h2D_map_map ;
0262 
0263     if (isData)
0264     {
0265         input_directory_map = data_input_directory_map;
0266         h1_only_logy_plot_vec = h1data_only_logy_plot_vec;
0267         h1D_map_map = &data_h1D_map_map;
0268         h2D_map_map = &data_h2D_map_map;
0269     }
0270     else
0271     {
0272         input_directory_map = MC_input_directory_map;
0273         h1_only_logy_plot_vec = h1MC_only_logy_plot_vec;
0274         h1D_map_map = &MC_h1D_map_map;
0275         h2D_map_map = &MC_h2D_map_map;
0276     }
0277 
0278     for (auto &pair : *(h1D_map_map)) // note : pair.first : directory name
0279     {
0280         for (auto &pair2 : pair.second) // note : pair2.first : hist name
0281         {   
0282             c1 -> cd();
0283 
0284             bool set_log_y = false;
0285             if (std::find(h1comp_logy_plot_vec.begin(), h1comp_logy_plot_vec.end(), pair2.first) != h1comp_logy_plot_vec.end())
0286             {
0287                 set_log_y = true;
0288             }
0289             if (std::find(h1_only_logy_plot_vec.begin(), h1_only_logy_plot_vec.end(), pair2.first) != h1_only_logy_plot_vec.end())
0290             {
0291                 set_log_y = true;
0292             }
0293 
0294             c1 -> SetLogy(set_log_y);
0295 
0296             std::string log_text;
0297 
0298             if (set_log_y){
0299                 pair2.second -> SetMaximum( pair2.second->GetBinContent(pair2.second->GetMaximumBin()) * logy_max_factor );
0300                 log_text = "_logy";
0301             }
0302             else {
0303                 pair2.second -> SetMaximum( pair2.second->GetBinContent(pair2.second->GetMaximumBin()) * lineary_max_factor );
0304                 log_text = "";
0305             }
0306 
0307             
0308             pair2.second -> GetXaxis() -> SetNdivisions(505);
0309             pair2.second -> Draw("hist");
0310 
0311             std::string final_labelling = (isData) ? sPH_labeling : MC_labeling;
0312 
0313             ltx->DrawLatex(1 - gPad->GetRightMargin(), 1 - gPad->GetTopMargin() + 0.01, Form("#it{#bf{sPHENIX}} %s", final_labelling.c_str()));        
0314 
0315             std::string leg_label = (draw_method == "hist") ? "f" : "ep";
0316 
0317             All_leg -> Clear();
0318             All_leg -> AddEntry(pair2.second, input_directory_map[pair.first].first.c_str(), leg_label.c_str());
0319             All_leg -> Draw("same");
0320 
0321             // note : draw labelling
0322             if (labelling_vec_map.find(sub_label_text) != labelling_vec_map.end())
0323             {
0324                 for (auto &tuple : labelling_vec_map[sub_label_text])
0325                 {
0326                     double x1 = std::get<0>(tuple);
0327                     double y1 = std::get<1>(tuple);
0328                     std::string text = std::get<2>(tuple);
0329 
0330                     sub_ltx->DrawLatex(x1, y1, text.c_str());
0331                 }
0332 
0333             }
0334 
0335             // todo: no sensor plot
0336             if (pair2.first.find("_l") != std::string::npos && pair2.first.find("_phi") != std::string::npos && pair2.first.find("_z") != std::string::npos)
0337             {
0338                 std::cout<<"what? : "<<pair2.first<<std::endl;
0339                 // c1 -> Print(Form("%s/data_%s_Draw%s.pdf",(output_directory + "/" + sensor_folder).c_str(),pair.first.c_str(),draw_method.c_str()));
0340             }
0341             else
0342             {
0343                 c1 -> Print(
0344                     Form(
0345                         "%s/%s_%s_Draw%s%s.pdf",
0346                         (output_directory + "/" + sub_label_text + "_" + input_directory_map[pair.first].second).c_str(),
0347                         sub_label_text.c_str(),
0348                         pair2.first.c_str(),
0349                         draw_method.c_str(),
0350                         log_text.c_str()
0351                     )
0352                 );
0353             }
0354 
0355             c1 -> Clear();
0356 
0357         } // note : end of all the histograms {h1D}
0358     } // note : end of all the data files {h1D}
0359 
0360     for (auto &pair : *(h2D_map_map)) // note : pair.first : directory name
0361     {
0362         for (auto &pair2 : pair.second) // note : pair2.first : hist name
0363         {   
0364             c1 -> cd();
0365             c1 -> SetLogy(false);
0366 
0367             std::string log_text = "";
0368             if (std::find(h2_logz_plot_vec.begin(), h2_logz_plot_vec.end(), pair2.first) != h2_logz_plot_vec.end())
0369             {
0370                 c1 -> SetLogz();
0371                 log_text = "_logz";
0372             }
0373 
0374             pair2.second -> Draw("colz0");
0375 
0376             std::string final_labelling = (isData) ? sPH_labeling : MC_labeling;
0377 
0378             ltx->DrawLatex(1 - gPad->GetRightMargin(), 1 - gPad->GetTopMargin() + 0.01, Form("#it{#bf{sPHENIX}} %s", final_labelling.c_str()));        
0379 
0380             // note : draw labelling
0381             if (labelling_vec_map.find(sub_label_text) != labelling_vec_map.end())
0382             {
0383                 for (auto &tuple : labelling_vec_map[sub_label_text])
0384                 {
0385                     double x1 = std::get<0>(tuple);
0386                     double y1 = std::get<1>(tuple);
0387                     std::string text = std::get<2>(tuple);
0388 
0389                     sub_ltx->DrawLatex(x1, y1, text.c_str());
0390                 }
0391 
0392             }
0393 
0394             c1 -> Print(
0395                 Form(
0396                     "%s/%s_%s%s.pdf",
0397                     (output_directory + "/" + sub_label_text + "_" + input_directory_map[pair.first].second).c_str(),
0398                     sub_label_text.c_str(),
0399                     pair2.first.c_str(), 
0400                     log_text.c_str()
0401                 )
0402             );
0403             c1 -> Clear();
0404             c1 -> SetLogz(false);
0405 
0406         } // note : end of all the histograms {h2D}
0407     } // note : end of all the data files {h2D}
0408 }
0409 
0410 void vtxZDist_comp::MakeDataPlot(std::string draw_method)
0411 {
0412     MakePlots(draw_method, true);
0413 }
0414 
0415 void vtxZDist_comp::MakeMCPlot(std::string draw_method)
0416 {
0417     MakePlots(draw_method, false);
0418 }
0419 
0420 void vtxZDist_comp::MakeVtxZCheckPlot()
0421 {
0422     // note : use the Inclusive70 as the reference
0423     for (auto &pair : data_h1D_map_map) // note : files
0424     {
0425         std::vector<std::pair<std::string, TH1D*>> all_Mbin_vtxZ_map; all_Mbin_vtxZ_map.clear();
0426 
0427         for (int i = 0; i < nCentrality_bin; i++)
0428         {
0429             make_comparison(
0430                 {
0431                     {Form("Centrality [%s]%%", Constants::centrality_text[i].c_str()), (TH1D*)(pair.second[Form("h1D_INTTz_Mbin%d", i)])->Clone(Form("h1D_INTTz_Mbin%d", i))}
0432                 },
0433                 {
0434                     {"Inclusive70", (TH1D*)(pair.second["h1D_INTTz_Inclusive70"])->Clone("h1D_INTTz_Inclusive70")}
0435                 },
0436                 output_directory + "/" + "Data_"+data_input_directory_map[pair.first].second,
0437                 Form("Comp_Inclusive70_INTTz_Mbin%d", i),
0438                 "Data",
0439                 false,
0440                 false, // note : isData_more
0441                 0,
0442                 -99
0443             );
0444 
0445             all_Mbin_vtxZ_map.push_back(
0446                 std::make_pair(
0447                     Form("Centrality [%s]%%", Constants::centrality_text[i].c_str()),
0448                     (TH1D*) pair.second[Form("h1D_INTTz_Mbin%d", i)] -> Clone(Form("h1D_INTTz_Mbin%d", i))
0449                 )
0450             );
0451         }
0452 
0453         make_comparison(
0454             all_Mbin_vtxZ_map,
0455             {
0456                 {"Inclusive70", (TH1D*)(pair.second["h1D_INTTz_Inclusive70"])->Clone("h1D_INTTz_Inclusive70")}
0457             },
0458             output_directory + "/" + "Data_"+data_input_directory_map[pair.first].second,
0459             Form("Comp_Inclusive70_INTTz_AllMbin"),
0460             "Data",
0461             false,
0462             true, // note : isData_more
0463             0,
0464             -99
0465         );
0466     }
0467 
0468 
0469     for (auto &pair : MC_h1D_map_map) // note : files
0470     {
0471         std::vector<std::pair<std::string, TH1D*>> all_Mbin_vtxZ_map; all_Mbin_vtxZ_map.clear();
0472 
0473         for (int i = 0; i < nCentrality_bin; i++)
0474         {
0475             make_comparison(
0476                 {
0477                     {Form("Centrality [%s]%%", Constants::centrality_text[i].c_str()), (TH1D*)(pair.second[Form("h1D_INTTz_Mbin%d", i)])->Clone(Form("h1D_INTTz_Mbin%d", i))}
0478                 },
0479                 {
0480                     {"Inclusive70", (TH1D*)(pair.second["h1D_INTTz_Inclusive70"])->Clone("h1D_INTTz_Inclusive70")}
0481                 },
0482                 output_directory + "/" + "MC_"+MC_input_directory_map[pair.first].second,
0483                 Form("Comp_Inclusive70_INTTz_Mbin%d", i),
0484                 "MC",
0485                 false,
0486                 false, // note : isData_more
0487                 0,
0488                 -99
0489             );
0490 
0491             all_Mbin_vtxZ_map.push_back(
0492                 std::make_pair(
0493                     Form("Centrality [%s]%%", Constants::centrality_text[i].c_str()),
0494                     (TH1D*) pair.second[Form("h1D_INTTz_Mbin%d", i)] -> Clone(Form("h1D_INTTz_Mbin%d", i))
0495                 )
0496             );
0497         }
0498 
0499         make_comparison(
0500             all_Mbin_vtxZ_map,
0501             {
0502                 {"Inclusive70", (TH1D*)(pair.second["h1D_INTTz_Inclusive70"])->Clone("h1D_INTTz_Inclusive70")}
0503             },
0504             output_directory + "/" + "MC_"+MC_input_directory_map[pair.first].second,
0505             Form("Comp_Inclusive70_INTTz_AllMbin"),
0506             "MC",
0507             false,
0508             true, // note : isData_more
0509             0,
0510             -99
0511         );
0512     }
0513 }
0514 
0515 void vtxZDist_comp::MakeComparisonPlot()
0516 {
0517     std::vector<std::pair<std::string, TH1D*>> data_hist_vec; data_hist_vec.clear();
0518     std::vector<std::pair<std::string, TH1D*>> MC_hist_vec; MC_hist_vec.clear();
0519 
0520     std::vector<std::string> plot_name; plot_name.clear();
0521     for (auto &pair : data_h1D_map_map){ // note : std::map<std::string, std::map<std::string, TH1D*>>
0522         for (auto &pair2 : pair.second){
0523             plot_name.push_back(pair2.first);
0524         }
0525         break;
0526     }
0527 
0528     for (int i = 0; i < plot_name.size(); i++)
0529     {
0530         for (auto &pair : data_h1D_map_map){ // note : std::map<std::string, std::map<std::string, TH1D*>>
0531             data_hist_vec.push_back(
0532                 {
0533                     data_input_directory_map[pair.first].first, // note : the description of that root file
0534                     pair.second[plot_name[i]]
0535                 }
0536             );
0537 
0538             // std::cout<<"Index: "<<i<<", for plot: "<<plot_name[i]<<", data: "<<pair.first<<", hist: "<<pair.second[plot_name[i]]->GetName()<<", NbinsX"<<pair.second[plot_name[i]]->GetNbinsX()<<std::endl;
0539         }
0540 
0541         for (auto &pair : MC_h1D_map_map){ // note : std::map<std::string, std::map<std::string, TH1D*>>
0542             MC_hist_vec.push_back(
0543                 {
0544                     MC_input_directory_map[pair.first].first, // note : the description of that root file
0545                     pair.second[plot_name[i]]
0546                 }
0547             );
0548 
0549             // std::cout<<"Index: "<<i<<", for plot: "<<plot_name[i]<<", MC: "<<pair.first<<", hist: "<<pair.second[plot_name[i]]->GetName()<<", NbinsX"<<pair.second[plot_name[i]]->GetNbinsX()<<std::endl;
0550         }
0551 
0552         bool m_set_log_y = false;
0553 
0554         if (std::find(h1comp_logy_plot_vec.begin(), h1comp_logy_plot_vec.end(), plot_name[i]) != h1comp_logy_plot_vec.end())
0555         {
0556             m_set_log_y = true;
0557         }
0558 
0559         make_comparison(
0560             data_hist_vec,
0561             MC_hist_vec,
0562             output_directory,
0563             Form("Comp_%s", plot_name[i].c_str()),
0564             "Comp",
0565             m_set_log_y,
0566             false, // note : isData_more
0567             0,
0568             -99
0569         );
0570 
0571         data_hist_vec.clear();
0572         MC_hist_vec.clear();
0573 
0574     }
0575 
0576 
0577 }
0578 
0579 void vtxZDist_comp::PrepareINTTvtxZReWeight()
0580 {
0581     TFile * file_out = new TFile(Form("%s/INTTvtxZReWeight.root", output_directory.c_str()), "RECREATE");
0582     
0583     // todo : only the first data file is used as the reference 
0584     // todo : only based on the Inclusive70, and Inclusive100
0585 
0586     TH1D * data_reference_Inclusive70 = (TH1D*) data_h1D_map_map[data_input_directory_pair_vec[0].first]["h1D_INTTz_Inclusive70"] -> Clone("h1D_INTTz_Inclusive70");
0587     data_reference_Inclusive70 -> Sumw2(true);
0588 
0589     TH1D * data_reference_Inclusive100 = (TH1D*) data_h1D_map_map[data_input_directory_pair_vec[0].first]["h1D_INTTz_Inclusive100"] -> Clone("h1D_INTTz_Inclusive100");
0590     data_reference_Inclusive100 -> Sumw2(true);
0591 
0592     TH1D * data_reference_Inclusive70_tight = (TH1D*) data_h1D_map_map[data_input_directory_pair_vec[0].first]["h1D_INTTz_Inclusive70_tight"] -> Clone("h1D_INTTz_Inclusive70_tight");
0593     data_reference_Inclusive70_tight -> Sumw2(true);
0594 
0595     TH1D * data_reference_Inclusive100_tight = (TH1D*) data_h1D_map_map[data_input_directory_pair_vec[0].first]["h1D_INTTz_Inclusive100_tight"] -> Clone("h1D_INTTz_Inclusive100_tight");
0596     data_reference_Inclusive100_tight -> Sumw2(true);
0597 
0598     std::map<std::string, TH1D*> reweighting_ratio; reweighting_ratio.clear();
0599 
0600     int count = 0; // note : should correspond to file index
0601 
0602     for (auto &pair : MC_h1D_map_map){
0603         reweighting_ratio.insert(
0604             std::make_pair(
0605                 Form("%d_Inclusive70", count),
0606                 (TH1D*) pair.second["h1D_INTTz_Inclusive70"] -> Clone((MC_input_directory_map[pair.first].second + "_Inclusive70").c_str())
0607             )
0608         );
0609 
0610         reweighting_ratio[Form("%d_Inclusive70", count)] -> Sumw2(true);
0611         reweighting_ratio[Form("%d_Inclusive70", count)] -> Divide(data_reference_Inclusive70, pair.second["h1D_INTTz_Inclusive70"], 1, 1);
0612         reweighting_ratio[Form("%d_Inclusive70", count)] -> GetYaxis() -> SetTitle("Data/MC");
0613         reweighting_ratio[Form("%d_Inclusive70", count)] -> Write((MC_input_directory_map[pair.first].second + "_Inclusive70").c_str());
0614 
0615 
0616         reweighting_ratio.insert(
0617             std::make_pair(
0618                 Form("%d_Inclusive100", count),
0619                 (TH1D*) pair.second["h1D_INTTz_Inclusive100"] -> Clone((MC_input_directory_map[pair.first].second + "_Inclusive100").c_str())
0620             )
0621         );
0622 
0623         reweighting_ratio[Form("%d_Inclusive100", count)] -> Sumw2(true);
0624         reweighting_ratio[Form("%d_Inclusive100", count)] -> Divide(data_reference_Inclusive100, pair.second["h1D_INTTz_Inclusive100"], 1, 1);
0625         reweighting_ratio[Form("%d_Inclusive100", count)] -> GetYaxis() -> SetTitle("Data/MC");
0626         reweighting_ratio[Form("%d_Inclusive100", count)] -> Write((MC_input_directory_map[pair.first].second + "_Inclusive100").c_str());
0627 
0628 
0629         // Division ---for tight-----------------------------------------------------------------------------------------------------------------------------
0630         reweighting_ratio.insert(
0631             std::make_pair(
0632                 Form("%d_Inclusive70_tight", count),
0633                 (TH1D*) pair.second["h1D_INTTz_Inclusive70_tight"] -> Clone((MC_input_directory_map[pair.first].second + "_Inclusive70_tight").c_str())
0634             )
0635         );
0636 
0637         reweighting_ratio[Form("%d_Inclusive70_tight", count)] -> Sumw2(true);
0638         reweighting_ratio[Form("%d_Inclusive70_tight", count)] -> Divide(data_reference_Inclusive70_tight, pair.second["h1D_INTTz_Inclusive70_tight"], 1, 1);
0639         reweighting_ratio[Form("%d_Inclusive70_tight", count)] -> GetYaxis() -> SetTitle("Data/MC");
0640         reweighting_ratio[Form("%d_Inclusive70_tight", count)] -> Write((MC_input_directory_map[pair.first].second + "_Inclusive70_tight").c_str());
0641 
0642 
0643         reweighting_ratio.insert(
0644             std::make_pair(
0645                 Form("%d_Inclusive100_tight", count),
0646                 (TH1D*) pair.second["h1D_INTTz_Inclusive100_tight"] -> Clone((MC_input_directory_map[pair.first].second + "_Inclusive100_tight").c_str())
0647             )
0648         );
0649 
0650         reweighting_ratio[Form("%d_Inclusive100_tight", count)] -> Sumw2(true);
0651         reweighting_ratio[Form("%d_Inclusive100_tight", count)] -> Divide(data_reference_Inclusive100_tight, pair.second["h1D_INTTz_Inclusive100_tight"], 1, 1);
0652         reweighting_ratio[Form("%d_Inclusive100_tight", count)] -> GetYaxis() -> SetTitle("Data/MC");
0653         reweighting_ratio[Form("%d_Inclusive100_tight", count)] -> Write((MC_input_directory_map[pair.first].second + "_Inclusive100_tight").c_str());
0654 
0655         count++;
0656     }
0657 
0658     file_out -> Close();
0659 
0660 }
0661 
0662 double vtxZDist_comp::GetNonZeroMinimalBin(TH1D * hist_in)
0663 {
0664     double min = 10000000;
0665     for (int i = 1; i <= hist_in->GetNbinsX(); i++)
0666     {
0667         if (hist_in->GetBinContent(i) > 0 && hist_in->GetBinContent(i) < min)
0668         {
0669             min = hist_in->GetBinContent(i);
0670         }
0671     }
0672 
0673     return min;
0674 }
0675 
0676 std::pair<double,double> vtxZDist_comp::make_comparison(
0677     std::vector<std::pair<std::string, TH1D*>> data_hist_vec_in,  // note : description, hist
0678     std::vector<std::pair<std::string, TH1D*>> MC_hist_vec_in,    // note : description, hist
0679     std::string output_directory_in, 
0680     std::string output_filename_in, 
0681     std::string sub_label_text_in, // note : "Data", "MC", "Comp"
0682     bool set_log_y, 
0683     bool isData_more,
0684     double MainPlot_y_low, 
0685     double MainPlot_y_top, 
0686     double ratio_y_low,
0687     double ratio_y_top
0688 )
0689 {
0690 
0691     std::cout<<"In make_comparison, working on "<<output_filename_in<<std::endl;
0692 
0693     if (MainPlot_y_top == -99)
0694     {
0695         for (auto &pair : data_hist_vec_in)
0696         {
0697             MainPlot_y_top = (pair.second->GetBinContent(pair.second->GetMaximumBin()) > MainPlot_y_top) ? pair.second->GetBinContent(pair.second->GetMaximumBin()) : MainPlot_y_top;
0698         }
0699 
0700         for (auto &pair : MC_hist_vec_in)
0701         {
0702             MainPlot_y_top = (pair.second->GetBinContent(pair.second->GetMaximumBin()) > MainPlot_y_top) ? pair.second->GetBinContent(pair.second->GetMaximumBin()) : MainPlot_y_top;
0703         }
0704 
0705 
0706         if (set_log_y){            
0707             MainPlot_y_top = MainPlot_y_top * logy_max_factor;
0708         }
0709         else{
0710             MainPlot_y_top = MainPlot_y_top * lineary_max_factor;
0711         }
0712     }
0713 
0714     if (MainPlot_y_low == 0 && set_log_y)
0715     {
0716         MainPlot_y_low = 1000000;
0717 
0718         for (auto &pair : data_hist_vec_in)
0719         {
0720             double min = GetNonZeroMinimalBin(pair.second);
0721 
0722             MainPlot_y_low = (min < MainPlot_y_low) ? min : MainPlot_y_low;
0723         }
0724 
0725         for (auto &pair : MC_hist_vec_in)
0726         {
0727             double min = GetNonZeroMinimalBin(pair.second);
0728 
0729             MainPlot_y_low = (min < MainPlot_y_low) ? min : MainPlot_y_low;
0730         }
0731     }
0732 
0733     std::vector<int> data_hist_nbinsX; data_hist_nbinsX.clear();
0734     std::vector<int> MC_hist_nbinsX; MC_hist_nbinsX.clear();
0735 
0736     int color_index = 0;
0737 
0738     for (auto &pair : data_hist_vec_in)
0739     {
0740 
0741         std::cout<<"Data: "<<pair.first<<", "<<pair.second->GetName()<<", "<<pair.second->GetNbinsX()<<std::endl;
0742 
0743         data_hist_nbinsX.push_back(pair.second->GetNbinsX());
0744         pair.second->SetLineColor(TColor::GetColor(color_code[color_index].c_str()));
0745         pair.second->SetMarkerColor(TColor::GetColor(color_code[color_index].c_str()));
0746         color_index += 1;
0747     }
0748 
0749     for (auto &pair : MC_hist_vec_in)
0750     {
0751         std::cout<<"MC: "<<pair.first<<", "<<pair.second->GetName()<<", "<<pair.second->GetNbinsX()<<std::endl;
0752 
0753         MC_hist_nbinsX.push_back(pair.second->GetNbinsX());
0754         pair.second->SetLineColor(TColor::GetColor(color_code[color_index].c_str()));
0755         pair.second->SetMarkerColor(TColor::GetColor(color_code[color_index].c_str()));
0756         color_index += 1;
0757     }
0758 
0759     for (int &data_nbin : data_hist_nbinsX){
0760         for (int &MC_nbin : MC_hist_nbinsX){
0761             if (data_nbin != MC_nbin){
0762                 std::cout<<"Error: Different binning for data and MC, data: "<<data_nbin<<", MC: "<<MC_nbin<<std::endl;
0763                 exit(1);
0764             }
0765         }
0766     }
0767 
0768     bool bins_are_same = true;
0769     for (int i = 1; i <= data_hist_vec_in[0].second->GetNbinsX(); i++)
0770     {
0771         if (data_hist_vec_in[0].second->GetBinCenter(i) != MC_hist_vec_in[0].second->GetBinCenter(i))
0772         {
0773             bins_are_same = false;
0774             break;
0775         }
0776     }
0777 
0778     TH1D * h1D_pad1 = (TH1D*)data_hist_vec_in.front().second->Clone("h1D_pad1");
0779     TH1D * h1D_pad2 = (TH1D*)data_hist_vec_in.front().second->Clone("h1D_pad2");
0780     h1D_pad1 -> Reset("ICESM");
0781     h1D_pad2 -> Reset("ICESM");
0782     h1D_pad1 -> SetLineColorAlpha(1,0.);
0783     h1D_pad2 -> SetLineColorAlpha(1,0.);
0784     h1D_pad1 -> SetFillColorAlpha(1,0.);
0785     h1D_pad2 -> SetFillColorAlpha(1,0.);
0786     h1D_pad1 -> SetMarkerColorAlpha(1,0.);
0787     h1D_pad2 -> SetMarkerColorAlpha(1,0.);
0788 
0789 
0790     SetsPhenixStyle();
0791 
0792     TCanvas * c1 = new TCanvas("", "c1", 950, 800);
0793 
0794     TPad * pad2 = new TPad("", "pad2", 0.0, 0.0, 1.0, 0.29);
0795     pad2 -> SetBottomMargin(0.3);
0796     pad2->Draw();
0797     
0798     TPad * pad1 = new TPad("", "pad1", 0.0, 0.29, 1.0, 1.0);
0799     pad1 -> SetTopMargin(0.065);
0800     pad1 -> SetBottomMargin(0.02);
0801     pad1->Draw();
0802 
0803     int total_hist = data_hist_vec_in.size() + MC_hist_vec_in.size();
0804 
0805     // todo: set the size fo the legend
0806     double leg_x1 = (total_hist <= 4) ? 0.45 : 0.65;
0807     double leg_y1 = (total_hist <= 4) ? 0.76 : 0.4;
0808     double leg_x2 = (total_hist <= 4) ? 0.8 : 0.8;
0809     double leg_y2 = (total_hist <= 4) ? 0.88 : 0.88;
0810 
0811     TLegend * leg = new TLegend(leg_x1, leg_y1, leg_x2, leg_y2);
0812     leg -> SetBorderSize(0);
0813     leg -> SetTextSize(0.025);
0814 
0815     h1D_pad1->GetXaxis()->SetLabelOffset(999);
0816     h1D_pad1->GetXaxis()->SetLabelSize(0);
0817     h1D_pad1->GetXaxis()->SetTitleOffset(999);
0818     h1D_pad1->GetYaxis()->SetLabelSize(0.045);
0819     h1D_pad1->GetYaxis()->SetTitleOffset(0.8);
0820     h1D_pad1->GetYaxis()->SetTitleSize(0.07);    
0821     h1D_pad1->GetXaxis()->SetNdivisions();
0822 
0823     for (auto &pair : MC_hist_vec_in){
0824         leg -> AddEntry(pair.second, pair.first.c_str(), "f");
0825     }
0826     
0827     for (auto &pair : data_hist_vec_in){
0828         // pair.second->GetXaxis()->SetNdivisions();
0829         leg -> AddEntry(pair.second, pair.first.c_str(), "lep");
0830     }
0831     
0832     pad1->cd();
0833 
0834     std::cout<<"For "<<output_filename_in<<", MainPlot_y_low: "<<MainPlot_y_low<<", MainPlot_y_top: "<<MainPlot_y_top<<std::endl;
0835 
0836     h1D_pad1 -> SetMinimum(MainPlot_y_low);
0837     h1D_pad1 -> SetMaximum(MainPlot_y_top);
0838     h1D_pad1 -> Draw("hist");
0839     
0840     for (int i = 0; i < MC_hist_vec_in.size(); i++){
0841         MC_hist_vec_in[i].second -> Draw("hist same");
0842     }
0843 
0844     for (int i = 0; i < data_hist_vec_in.size(); i++){
0845         data_hist_vec_in[i].second -> Draw("ep same");
0846     }
0847 
0848     leg -> Draw("same");
0849 
0850     pad1->SetLogy(set_log_y);
0851 
0852     ltx->DrawLatex(1 - gPad->GetRightMargin(), 1 - gPad->GetTopMargin() + 0.01, Form("#it{#bf{sPHENIX}} %s", sPH_labeling.c_str()));
0853 
0854     if (labelling_vec_map.find(sub_label_text_in) != labelling_vec_map.end())
0855     {
0856         for (auto &tuple : labelling_vec_map[sub_label_text_in])
0857         {
0858             double x1 = std::get<0>(tuple);
0859             double y1 = std::get<1>(tuple);
0860             std::string text = std::get<2>(tuple);
0861 
0862             sub_ltx->DrawLatex(x1, y1, text.c_str());
0863         }
0864     }
0865 
0866     TLine * line = new TLine();
0867     line -> SetLineStyle(2);
0868     line -> SetLineWidth(2);
0869     line -> SetLineColor(28);
0870 
0871     pad2->cd();
0872 
0873     std::vector<TH1D*> ratio_hist_vec_in; ratio_hist_vec_in.clear();
0874 
0875     std::vector<std::pair<std::string, TH1D*>> hist_vec_in = (isData_more) ? data_hist_vec_in : MC_hist_vec_in;
0876 
0877     for (int i = 0; i < hist_vec_in.size(); i++){
0878         
0879         std::pair<std::string, TH1D*> pair = hist_vec_in[i];
0880 
0881         // todo: can only do the ratio w.r.t the first data hist
0882         if (isData_more){
0883             ratio_hist_vec_in.push_back( // note : hist_vec_in: data
0884                 (TH1D*) MC_hist_vec_in[0].second->Clone(Form("ratio_hist_%s", pair.first.c_str()))
0885             );
0886 
0887             ratio_hist_vec_in.back() -> Sumw2(true);
0888             if (bins_are_same) {ratio_hist_vec_in.back() -> Divide(pair.second, ratio_hist_vec_in.back());}
0889 
0890         }
0891         else { // note : hist_vec_in: MC
0892             ratio_hist_vec_in.push_back(
0893                 (TH1D*) data_hist_vec_in[0].second->Clone(Form("ratio_hist_%s", pair.first.c_str()))
0894             );
0895 
0896             ratio_hist_vec_in.back() -> Sumw2(true);
0897             if (bins_are_same) {ratio_hist_vec_in.back() -> Divide(pair.second);}
0898         }
0899 
0900 
0901         // if (ratio_hist_vec_in.back()->GetXaxis()->GetLabelOffset() > 10)
0902         // {
0903         //     ratio_hist_vec_in.back()->GetXaxis()->SetLabelOffset(0.01);
0904         // }
0905 
0906         ratio_hist_vec_in.back()->SetMarkerStyle(20);
0907         ratio_hist_vec_in.back()->SetMarkerSize(0.8);
0908         ratio_hist_vec_in.back()->SetMarkerColor(pair.second->GetLineColor());
0909         ratio_hist_vec_in.back()->SetLineColor(pair.second->GetLineColor());
0910 
0911         if (i == 0){
0912             h1D_pad2->GetXaxis()->SetLabelSize(0.11);
0913             h1D_pad2->GetXaxis()->SetTitleOffset(0.8);
0914             h1D_pad2->GetXaxis()->SetTitleSize(0.14);
0915             h1D_pad2->GetXaxis()->SetNdivisions();
0916 
0917             h1D_pad2->GetYaxis()->SetLabelSize(0.11);
0918             h1D_pad2->GetYaxis()->SetTitleOffset(0.4);
0919             h1D_pad2->GetYaxis()->SetTitleSize(0.14);
0920             h1D_pad2->GetYaxis()->SetNdivisions(505);
0921 
0922             h1D_pad2 -> SetMinimum(ratio_y_low);
0923             h1D_pad2 -> SetMaximum(ratio_y_top);
0924 
0925             h1D_pad2->GetYaxis()->SetTitle("Data/MC");
0926 
0927 
0928             h1D_pad2->Draw("ep");
0929             line -> DrawLine(h1D_pad2->GetXaxis()->GetXmin(), 1, h1D_pad2->GetXaxis()->GetXmax(), 1);
0930 
0931             ratio_hist_vec_in.back()->Draw("ep same");
0932         }
0933         else {
0934             ratio_hist_vec_in.back()->Draw("ep same");
0935         }
0936     }
0937 
0938 
0939     
0940     // std::string data_hist_name = data_hist->GetName(); 
0941 
0942 
0943     // if (data_hist_name.find("_l") != std::string::npos && data_hist_name.find("_phi") != std::string::npos && data_hist_name.find("_z") != std::string::npos)
0944     // {
0945     //     std::cout<<"what? : "<<data_hist_name<<std::endl;
0946 
0947     //     // c1 -> Print(Form("%s/Comp_%s.pdf", (output_directory + "/" + sensor_folder).c_str(), output_filename.c_str()));
0948     // }
0949     if (true)
0950     {
0951         c1 -> Print(Form("%s/Comp_%s.pdf", output_directory_in.c_str(), output_filename_in.c_str()));
0952     }
0953 
0954     return {};
0955     
0956     // if (data_hist_name.find("_l") == std::string::npos || data_hist_name.find("_phi") == std::string::npos || data_hist_name.find("_z") == std::string::npos) {return {};}
0957 
0958     // return GetHistYInfo(ratio_hist);
0959 
0960     
0961     // pad1 -> Clear();
0962     // pad2 -> Clear();
0963     // c1 -> Clear();
0964 
0965     // delete pad1;
0966     // delete pad2;
0967     // delete c1;
0968     // delete ratio_hist;
0969 
0970     // return;
0971 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!