LightFlavorRatios/swimming_correction/getCutEfficiency_systematics.C

0001 #include "../util/HistogramTools.h"
0002 #include "../config/binning.h"
0003
0004 void getCutEfficiency_systematics(std::string central_filename = "/sphenix/tg/tg01/hf/mjpeters/LightFlavorProduction/cutEfficiencyCorrection/cut_efficiency_correction.root",
0005                               std::string min_filename = "/sphenix/tg/tg01/hf/mjpeters/LightFlavorProduction/cutEfficiencyCorrection/systematics_min/cut_efficiency_correction.root",
0006                               std::string max_filename = "/sphenix/tg/tg01/hf/mjpeters/LightFlavorProduction/cutEfficiencyCorrection/systematics_max/cut_efficiency_correction.root",
0007                               std::string outfilename = "./syserr_test.root", std::string numerator_name = "Lambda0", std::string denominator_name = "K_S0",
0008                               std::string numerator_title = "(#Lambda+#bar{#Lambda})", std::string denominator_title = "K_{S}^{0}")
0009 {
0010   TFile* f_central = TFile::Open(central_filename.c_str());
0011   TFile* f_minthreshold = TFile::Open(min_filename.c_str());
0012   TFile* f_maxthreshold = TFile::Open(max_filename.c_str());
0013   TFile* f_out = new TFile(outfilename.c_str(),"RECREATE");
0014
0015   std::vector<HistogramInfo> variables =
0016   {
0017     BinInfo::final_pt_bins,
0018     BinInfo::final_eta_bins,
0019     BinInfo::final_phi_bins,
0020     BinInfo::final_rapidity_bins
0021   };
0022
0023   for(HistogramInfo& h : variables)
0024   {
0025     TH1F* central_numerator_cuteff = (TH1F*)f_central->Get((numerator_name+"_cuteff_vs"+h.name).c_str());
0026     TH1F* central_denominator_cuteff = (TH1F*)f_central->Get((denominator_name+"_cuteff_vs"+h.name).c_str());
0027     TH1F* central_correction = (TH1F*)f_central->Get((numerator_name+"_over_"+denominator_name+"_cuteff_correction_vs"+h.name).c_str());
0028     TH1F* central_inverse_correction = (TH1F*)f_central->Get((denominator_name+"_over_"+numerator_name+"_cuteff_correction_vs"+h.name).c_str());
0029
0030     TH1F* min_numerator_cuteff = (TH1F*)f_minthreshold->Get((numerator_name+"_cuteff_vs"+h.name).c_str());
0031     TH1F* min_denominator_cuteff = (TH1F*)f_minthreshold->Get((denominator_name+"_cuteff_vs"+h.name).c_str());
0032     TH1F* min_correction = (TH1F*)f_minthreshold->Get((numerator_name+"_over_"+denominator_name+"_cuteff_correction_vs"+h.name).c_str());
0033     TH1F* min_inverse_correction = (TH1F*)f_minthreshold->Get((denominator_name+"_over_"+numerator_name+"_cuteff_correction_vs"+h.name).c_str());
0034
0035     TH1F* max_numerator_cuteff = (TH1F*)f_maxthreshold->Get((numerator_name+"_cuteff_vs"+h.name).c_str());
0036     TH1F* max_denominator_cuteff = (TH1F*)f_maxthreshold->Get((denominator_name+"_cuteff_vs"+h.name).c_str());
0037     TH1F* max_correction = (TH1F*)f_maxthreshold->Get((numerator_name+"_over_"+denominator_name+"_cuteff_correction_vs"+h.name).c_str());
0038     TH1F* max_inverse_correction = (TH1F*)f_maxthreshold->Get((denominator_name+"_over_"+numerator_name+"_cuteff_correction_vs"+h.name).c_str());
0039
0040     f_out->cd();
0041
0042     std::vector<std::string> syserr_names =
0043     {
0044       numerator_name+"_cuteff_syserr",
0045       denominator_name+"_cuteff_syserr",
0046       numerator_name+"_over_"+denominator_name+"_cuteff_correction_syserr",
0047       denominator_name+"_over_"+numerator_name+"_cuteff_correction_syserr"
0048     };
0049
0050     std::vector<std::string> syserr_titles =
0051     {
0052       numerator_title+" cut efficiency systematic uncertainty",
0053       denominator_title+" cut efficiency systematic uncertainty",
0054       numerator_title+" / "+denominator_title+" cut efficiency correction systematic uncertainty",
0055       denominator_title+" / "+numerator_title+" cut efficiency correction systematic uncertainty"
0056     };
0057
0058     std::vector<std::array<TH1F*,3>> syserr_triplets =
0059     {
0060       {central_numerator_cuteff,min_numerator_cuteff,max_numerator_cuteff},
0061       {central_denominator_cuteff,min_denominator_cuteff,max_denominator_cuteff},
0062       {central_correction,min_correction,max_correction},
0063       {central_inverse_correction,min_inverse_correction,max_inverse_correction}
0064     };
0065
0066     for(int i=0;i<syserr_triplets.size();i++)
0067     {
0068       TH1F* syserr = makeHistogram(syserr_names[i].c_str(),syserr_titles[i].c_str(),h);
0069       for(int bin=1;bin<=syserr_triplets[i][0]->GetNbinsX();bin++)
0070       {
0071         std::vector<double> vals;
0072         for(TH1F* th : syserr_triplets[i])
0073         {
0074           vals.push_back(th->GetBinContent(bin));
0075         }
0076         std::sort(vals.begin(),vals.end());
0077         syserr->SetBinContent(bin,fabs(vals.back()-vals.front())/2.);
0078       }
0079       syserr->Write();
0080     }
0081   }
0082 }