37 std::pair<TH1*,TH1*>
takeRatios(TH1 *h_MCData_Num, TH1 *h_MCData_Denom, TH1 *h_MCMC_Num, TH1 *h_MCMC_Denom)
46 TH1 *h_MCData_Rat = (TH1*)h_MCData_Num->Clone();
47 h_MCData_Rat->Reset();
48 TH1 *h_MCMC_Rat = (TH1*)h_MCMC_Num->Clone();
51 h_MCData_Rat->Divide(h_MCData_Num, h_MCData_Denom, 1, 1,
"B");
52 h_MCMC_Rat ->Divide(h_MCMC_Num, h_MCMC_Denom, 1, 1,
"B");
54 std::pair<TH1*,TH1*> h_Rats = {h_MCData_Rat, h_MCMC_Rat};
61 TCanvas*
makeCanvas(TH1 *h_MCData, TH1 *h_MCMC, TString s_Obj)
65 h_MCData->SetLineColor(kBlack);
66 h_MCData->SetLineWidth(2);
67 h_MCData->GetXaxis()->SetTitle(
"");
68 h_MCData->GetYaxis()->SetTitle(
"Efficiency");
69 h_MCData->GetYaxis()->SetTitleSize(0.05);
70 h_MCData->GetYaxis()->SetTitleOffset(1.2);
71 h_MCData->GetYaxis()->CenterTitle();
72 h_MCMC->SetLineColor(
kRed);
73 h_MCMC->SetLineWidth(2);
74 double ymax =
std::max(h_MCMC->GetMaximum(), h_MCData->GetMaximum());
76 TCanvas *
c =
new TCanvas(
"c",
"c", 800, 800);
77 TPad *
p1 =
new TPad(
"p1",
"p1", 0, 0, 1, 1);
79 p1->SetBottomMargin(0.3);
80 p1->SetTopMargin(0.01);
81 p1->SetRightMargin(0.05);
82 p1->SetLeftMargin(0.13);
83 TPad *
p2 =
new TPad(
"p2",
"p2", 0, 0, 1, 1);
85 p2->SetTopMargin(0.7);
86 p2->SetRightMargin(0.05);
87 p2->SetBottomMargin(0.1);
88 p2->SetLeftMargin(0.13);
94 h_MCData->GetYaxis()->SetRangeUser(0, ymax + ymax*0.3);
98 TH1D *h_Rat = (TH1D*)h_MCData->Clone();
99 h_MCData->GetXaxis()->SetLabelSize(0);
100 h_Rat->Divide(h_MCMC);
101 h_Rat->SetAxisRange(0.89,1.11,
"Y");
102 h_Rat->SetLineColor(
kBlue);
104 h_Rat->GetYaxis()->SetTitle(
"MCData/MCMC");
105 h_Rat->GetYaxis()->SetTitleSize(0.05);
106 h_Rat->GetXaxis()->CenterTitle();
107 h_Rat->GetYaxis()->CenterTitle();
110 double x[2] = {h_Rat->GetXaxis()->GetXmin(), h_Rat->GetXaxis()->GetXmax()};
111 double y[2] = {1.0, 1.0};
112 TGraph* l_Xaxis =
new TGraph(2, x, y);
113 l_Xaxis->SetLineColor(kGray);
114 l_Xaxis->SetLineWidth(3);
115 l_Xaxis->SetLineStyle(2);
116 l_Xaxis->Draw(
"L SAME");
118 double x1_leg;
double x2_leg;
119 s_Obj.Contains(
"kPDG") ? x1_leg = 0.4 : x1_leg = 0.65;
120 s_Obj.Contains(
"kPDG") ? x2_leg = 0.6 : x2_leg = 0.89;
121 TLegend *
leg =
new TLegend(x1_leg, 0.84, x2_leg, 0.98);
122 leg->AddEntry(h_MCData,
"MC singles in Data",
"L");
123 leg->AddEntry(h_MCMC,
"MC singles in MC",
"L");
136 double effDiff = map_HistPair[
cut].first ->Integral()/map_HistPair[s_Denom].first->Integral() - map_HistPair[
cut].second->Integral()/map_HistPair[s_Denom].second->Integral();
137 effDiff /= map_HistPair[
cut].second->Integral()/map_HistPair[s_Denom].second->Integral();
147 double effDiffErrPc = 0.;
177 double nSelectedData = map_HistPair[
cut].first ->Integral();
178 double nSelectedMC = map_HistPair[
cut].second->Integral();
180 double pData = nSelectedData/map_HistPair[s_Denom].first ->Integral();
181 double pMC = nSelectedMC /map_HistPair[s_Denom].second->Integral();
183 double pVarData = pData*(1-pData)/map_HistPair[s_Denom].first ->
Integral();
184 double pVarMC = pMC *(1-pMC) /map_HistPair[s_Denom].
second->Integral();
189 double varDelta = pSdData*pSdData + pSdMC*pSdMC;
192 double sdDeltaFrac = sdDelta/pMC;
195 effDiffErrPc = sdDeltaFrac*100.;
204 gStyle->SetOptStat(0);
206 TFile *f_InFileMCData =
new TFile(s_InFileMCData.c_str(),
"READ");
207 TFile *f_InFileMCMC =
new TFile(s_InFileMCMC .c_str(),
"READ");
209 std::map<std::string,const Var> map_Vars;
210 std::map<std::string,const NuTruthVar> map_NuVars;
213 std::vector<std::string> vec_Cuts = help.
CutFlow(s_FileAppend, useAnaTrueECut, getDetailed);
215 std::string lastEntry = vec_Cuts.at(vec_Cuts.size()-1);
216 std::map<std::string,std::pair<TH1*,TH1*>> map_HistPair;
217 std::map<std::string,std::pair<TH1*,TH1*>> map_HistPair_AllVars;
218 TIter iter(f_InFileMCData->GetListOfKeys());
220 while((key=(TKey*)iter())){
221 TString s_Obj = ((TObject*)key->ReadObj())->
GetName();
223 for(
auto var : map_Vars){
225 if(!
makePlots) s_FileAppend.find(
"nue")!=std::string::npos ? s_EffVar =
"kTrueEben" : s_EffVar =
"kTrueE";
226 for(
auto cut : vec_Cuts){
228 std::string s_NuSpec =
"dir_-CUT-" +
cut +
"-VAR-"+ s_EffVar +
"_Nu";
235 TH1 *h_MCData = spec_MCData->ToTH1(spec_MCData->POT());
236 h_MCData->SetName(s_Obj);
237 TH1 *h_MCMC = spec_MCMC ->ToTH1(spec_MCMC ->
POT());
238 h_MCMC->SetName(s_Obj);
239 map_HistPair[
cut] = {h_MCData, h_MCMC};
240 if((
cut==firstEntry) || (
cut==lastEntry)){
241 map_HistPair_AllVars[(
std::string)s_Obj] = {h_MCData, h_MCMC};
253 std::map<std::string,std::pair<std::string,std::string>> map_SpecialPairs;
254 std::map<std::string,std::vector<TH1*>> map_SpecialHists;
255 map_SpecialPairs[
"kSpillPOT"] = {
"dir_-CUT-" + lastEntry +
"-VAR-kSpillPOT",
"dir_-CUT-" + firstEntry +
"-VAR-kSpillPOTInNuTruth_Nu"};
257 for(
auto pair : map_SpecialPairs){
258 std::unique_ptr<Spectrum> spec_MCData_Num =
Spectrum::LoadFrom(f_InFileMCData, (TString)pair.second.first );
259 std::unique_ptr<Spectrum> spec_MCData_Denom =
Spectrum::LoadFrom(f_InFileMCData, (TString)pair.second.second);
260 std::unique_ptr<Spectrum> spec_MCMC_Num =
Spectrum::LoadFrom(f_InFileMCMC , (TString)pair.second.first );
261 std::unique_ptr<Spectrum> spec_MCMC_Denom =
Spectrum::LoadFrom(f_InFileMCMC , (TString)pair.second.second);
263 std::vector<TH1*> vec_Temp;
264 vec_Temp.push_back(spec_MCData_Num ->ToTH1(spec_MCData_Num ->
POT()));
265 vec_Temp.at(vec_Temp.size()-1)->SetName(pair.second.first.c_str());
266 vec_Temp.push_back(spec_MCData_Denom->ToTH1(spec_MCData_Denom->POT()));
267 vec_Temp.at(vec_Temp.size()-1)->SetName(pair.second.second.c_str());
268 vec_Temp.push_back(spec_MCMC_Num ->ToTH1(spec_MCMC_Num ->
POT()));
269 vec_Temp.at(vec_Temp.size()-1)->SetName(pair.second.first.c_str());
270 vec_Temp.push_back(spec_MCMC_Denom ->ToTH1(spec_MCMC_Denom ->
POT()));
271 vec_Temp.at(vec_Temp.size()-1)->SetName(pair.second.second.c_str());
273 map_SpecialHists[pair.first] = vec_Temp;
277 if(useAnaTrueECut) s_FileAppend +=
"_ATE";
278 if(getDetailed) s_FileAppend +=
"_Detailed";
284 f_Out.open(s_OutDir +
"/CutFlowTable_" + s_FileAppend +
".tex");
286 f_Out <<
"\\documentclass{article}" <<
std::endl;
287 f_Out <<
"\\usepackage[a4paper,margin=1in,landscape]{geometry}"<<
std::endl;
288 f_Out <<
"\\usepackage{multirow}" <<
std::endl;
289 f_Out <<
"\\usepackage{array}" <<
std::endl;
290 f_Out <<
"\\usepackage{amsmath}" <<
std::endl;
291 f_Out <<
"\\newcolumntype{A}{>{\\centering\\arraybackslash} m{3cm}}" <<
std::endl;
292 f_Out <<
"\\newcolumntype{B}{>{\\centering\\arraybackslash} m{4.2cm}}" <<
std::endl;
294 f_Out <<
"\\begin{document}" <<
std::endl;
296 f_Out <<
"\\begin{tabular}{|A|A|A|A|A|B|}" <<
std::endl;
298 f_Out <<
"\\multirow{2}{*}{}&\\multicolumn{2}{c|}{\\textbf{MC in Data}} & \\multicolumn{2}{c|}{\\textbf{MC in MC}} &" <<
std::endl;
299 f_Out <<
"\\multirow{2}{*}{$\\frac{\\left(\\textbf{Data}_\\textbf{eff} - \\textbf{MC}_\\textbf{eff}\\right)}{\\textbf{MC}_\\textbf{eff}}$, \\%}\\\\" <<
std::endl;
301 f_Out <<
"& \\textbf{Events} & \\textbf{Eff, \\%} & \\textbf{Events} & \\textbf{Eff, \\%} & \\\\" <<
std::endl;
305 for(
unsigned int cut = 0;
cut < vec_Cuts.size();
cut++){
306 f_Out <<
"\\textbf{" << help.
getNiceCutLabel(vec_Cuts.at(
cut)) <<
"} & & & & & \\\\" << std::endl;
307 f_Out <<
" & " << map_HistPair[vec_Cuts.at(
cut)].first->Integral();
308 f_Out <<
" & " << map_HistPair[vec_Cuts.at(
cut)].first->Integral()*100./map_HistPair[s_Denom].first->Integral();
309 f_Out <<
" & " << map_HistPair[vec_Cuts.at(
cut)].second->Integral();
310 f_Out <<
" & " << map_HistPair[vec_Cuts.at(
cut)].second->Integral()*100./map_HistPair[s_Denom].second->Integral();
311 f_Out <<
" & \\large " <<
getEffDiffPc(map_HistPair, vec_Cuts.at(
cut), s_Denom);
312 if(
cut==vec_Cuts.size()-1){
325 for(
auto hist : map_HistPair_AllVars){
327 if(
hist.first.find(
"CUT-"+firstEntry+
"-VAR")!=std::string::npos){
330 s_Num.replace(s_Denom.find(firstEntry), firstEntry.length(), lastEntry);
331 s_Num = s_Num.substr(0, s_Num.length()-3);
333 TH1 *h_MCData_Num = map_HistPair_AllVars[s_Num].first;
334 TH1 *h_MCData_Denom =
hist.second.first;
335 TH1 *h_MCMC_Num = map_HistPair_AllVars[s_Num].second;
336 TH1 *h_MCMC_Denom =
hist.second.second;
338 std::pair<TH1*,TH1*> h_Rats =
takeRatios(h_MCData_Num, h_MCData_Denom, h_MCMC_Num, h_MCMC_Denom);
339 TCanvas *
c =
makeCanvas(h_Rats.first, h_Rats.second, s_Num);
342 c->SaveAs(s_OutName.c_str());
349 for(
auto hists : map_SpecialHists){
351 TCanvas *
c =
makeCanvas(h_Rats.first, h_Rats.second,
hists.second.at(1)->GetName());
353 if(((TString)
hists.second.at(1)->GetName()).Contains(
"kSpillPOT")){
354 TBox *box =
new TBox(0.955, 0.1, 1.0, 0.2);
355 box->SetFillColor(kWhite);
360 c->SaveAs(s_OutName.c_str());
void Analyse_GetEfficiency(std::string s_OutDir="", std::string s_FileAppend="", std::string s_InFileMCData="", std::string s_InFileMCMC="", bool useAnaTrueECut=true, bool makePlots=false, bool getDetailed=false)
T max(const caf::Proxy< T > &a, T b)
std::vector< std::string > CutFlow(std::string const &s_Selection, bool useAnaTrueECut, bool getDetailed)
Cuts and Vars for the 2020 FD DiF Study.
double Integral(const Spectrum &s, const double pot, int cvnregion)
double getEffDiffPc(std::map< std::string, std::pair< TH1 *, TH1 * >> map_HistPair, std::string const &cut, std::string const &s_Denom)
std::string getNiceCutLabel(std::string const &s_Cut)
void scheduleVars(std::map< std::string, const Var > &map_Vars, std::map< std::string, const NuTruthVar > &map_NuVars, std::string s_Selection)
std::string GetName(int i)
static std::unique_ptr< Spectrum > LoadFrom(TDirectory *dir, const std::string &name)
std::pair< TH1 *, TH1 * > takeRatios(TH1 *h_MCData_Num, TH1 *h_MCData_Denom, TH1 *h_MCMC_Num, TH1 *h_MCMC_Denom)
std::string VarLabel(std::string s_DirName, bool isNice)
double getEffDiffErrPc(std::map< std::string, std::pair< TH1 *, TH1 * >> map_HistPair, std::string const &cut, std::string const &s_Denom)
std::vector< double > POT
TCanvas * makeCanvas(TH1 *h_MCData, TH1 *h_MCMC, TString s_Obj)