8 #include "../testingMacros/atlasstyle/AtlasStyle.C"
9 #include "../testingMacros/atlasstyle/AtlasLabels.C"
35 if (fabs(numerator) < 1.
e-3)
38 std::cout <<
"Setting error code in steering" << std::endl;
48 static TLatex* tex = []()
52 tex->SetTextSize(0.045);
53 tex->SetTextColor(kBlack);
62 const TString config2 = uncToolDiff ? uncToolDiff->getConfigFile() :
"";
63 const TString
year = config1.Contains(
"_2011/") ?
"2011" : config1.Contains(
"_2012/") ?
"2012": config1.Contains(
"_2015/") ?
"2015" :
"UNKNOWN";
65 const TString scenario1 = config1.Contains(
"StrongerCorrelations") ?
"stronger" : (config1.Contains(
"WeakerCorrelations") ?
"weaker" :
"nominal");
66 const TString scenario2 = config2.Contains(
"StrongerCorrelations") ?
"stronger" : (config2.Contains(
"WeakerCorrelations") ?
"weaker" :
"nominal");
67 const TString scenarioDiff = Form(
"%s - %s correlation scenarios",scenario1.Data(),scenario2.Data());
68 const TString reduction1 = config1.Contains(
"_ByCategory") ?
"category reduction" : (config1.Contains(
"NP") ?
"global reduction" :
"original");
69 const TString reduction2 = config2.Contains(
"_ByCategory") ?
"category reduction" : (config2.Contains(
"NP") ?
"global reduction" :
"original");
70 const TString reductionDiff = Form(
"%s - %s",reduction1.Data(),reduction2.Data());
72 TString configDescString;
75 if (scenario1 == scenario2 && reduction1 == reduction2)
76 configDescString =
"UNKNOWN COMPARISON TYPE";
77 else if (scenario1 == scenario2 && reduction1 != reduction2)
78 configDescString = Form(
"%s correlation scenario, %s",scenario1.Data(),reductionDiff.Data());
79 else if (scenario1 != scenario2 && reduction1 == reduction2)
80 configDescString = Form(
"%s%s",scenarioDiff.Data(),reduction1 ==
"original" ?
"" : Form(
", %s",reduction1.Data()));
82 configDescString = Form(
"%s, %s",scenarioDiff.Data(),reductionDiff.Data());
86 configDescString = Form(
"%s correlation scenario%s",scenario1.Data(),reduction1 ==
"original" ?
"" : Form(
", %s",reduction1.Data()));
91 const TString fixedValString = TString(
histo->GetXaxis()->GetTitle()).Contains(
"#it{p}") ? Form(
"(#eta^{jet1},#eta^{jet2}) = (%.1f,%.1f)",fixedValue1,fixedValue2) : Form(
"(#it{p}_{T}^{jet1},#it{p}_{T}^{jet2}) = (%ld,%ld) GeV",lround(fixedValue1),lround(fixedValue2));
95 const TString reductionString =
"str.red";
97 const TString scenarioString = TString(
histo->GetName()).BeginsWith(
"diff_")?Form(
"Correlation difference, Rep_{full}^{%sJES} - Rep_{%s}^{%sJES}, %s",
selector.Data(),reductionString.Data(),
selector.Data(),fixedValString.Data()):Form(
"Correlation matrix, Rep_{%s}^{%sJES}, %s",TString(
histo->GetName()).Contains(
"_0_")?
"full":reductionString.Data(),
selector.Data(),fixedValString.Data());
98 tex->DrawLatex(0.13,0.8751,Form(
"anti-#it{k}_{t} #it{R} = %s, %s+JES %s",
jetDef.Contains(
"AntiKt4") ?
"0.4" :
jetDef.Contains(
"AntiKt6") ?
"0.6" :
"UNKNOWN",
jetDef.Contains(
"EM") ?
"EM" :
jetDef.Contains(
"LC") ?
"LCW" :
"UNKNOWN",
year.Data()));
99 tex->DrawLatex(0.13,0.9305,scenarioString.Data());
103 tex->DrawLatex(0.48,0.960,Form(
"anti-#it{k}_{t} #it{R} = %s, %s+JES + #it{in situ}",
jetDef.Contains(
"AntiKt4") ?
"0.4" :
jetDef.Contains(
"AntiKt6") ?
"0.6" :
"UNKNOWN",
jetDef.Contains(
"EM") ?
"EM" :
jetDef.Contains(
"LC") ?
"LCW" :
"UNKNOWN"));
104 tex->DrawLatex(0.48,0.905,fixedValString.Data());
105 tex->DrawLatex(0.13,0.905,Form(
"Data %s, #sqrt{s} = %d TeV",
year.Data(),
year==
"2012"?8:
year==
"2011"?7:
year==
"2015"||
year==
"2016"?13:0));
111 tex->DrawLatex(0.15,0.84,configDescString.Data());
112 tex->DrawLatex(0.15,0.79,TString(uncTool->
getAnalysisFile().c_str()).Contains(
"Unknown") ?
"unknown flavour composition" : TString(uncTool->
getAnalysisFile().c_str()).Contains(
"InclusiveJets") ?
"inclusive jets composition" :
"Unexpected flavour composition file");
118 if (TString(
histo->GetXaxis()->GetTitle()).Contains(
"#it{p}"))
119 tex->DrawLatex(0.01,0.015,Form(
"Mean value %.2f, max %.2f at (%ld,%ld) GeV",
round(100*
mean)/100.,
round(100*extremum)/100.,lround(
histo->GetXaxis()->GetBinCenter(extremeX)),lround(
histo->GetYaxis()->GetBinCenter(extremeY))));
121 tex->DrawLatex(0.01,0.015,Form(
"Mean value %.2f, max %.2f at (%.2f,%.2f)",
round(100*
mean)/100.,
round(100*extremum)/100.,
histo->GetXaxis()->GetBinCenter(extremeX),
histo->GetYaxis()->GetBinCenter(extremeY)));
127 histo->SetTitleSize(0.06,
"z");
128 histo->SetLabelSize(0.04,
"z");
129 histo->SetTitleOffset(1.225,
"x");
130 histo->SetTitleOffset(1.05,
"y");
131 histo->SetTitleOffset(0.76,
"z");
132 histo->SetLabelOffset(0.002,
"x");
133 histo->GetXaxis()->SetMoreLogLabels();
136 void PlotCorrelationHistos(
const TString&
outFile, TCanvas*
canvas, TFile* outHistFile,
const std::vector<JetUncertaintiesTool*>& providers,
const std::vector<TH2D*>& corrMats,
const double fixedValue1,
const double fixedValue2)
140 std::vector<TH2D*> corrMatDiffs;
141 std::vector< std::pair<JetUncertaintiesTool*,JetUncertaintiesTool*> > diffProviders;
142 for (
size_t iHisto = 0; iHisto + 1 < corrMats.size(); ++iHisto)
145 TH2D*
diff =
new TH2D(*corrMats.at(iHisto));
146 diff->SetName(Form(
"diff_%s",
diff->GetName()));
151 if (corrMats.at(iHisto)->GetBinContent(
binX,
binY) < -100)
165 corrMatDiffs.push_back(
diff);
166 diffProviders.emplace_back(providers.at(iHisto),providers.at(iHisto+1));
170 for (
size_t iHisto = 0; iHisto < corrMats.size(); ++iHisto)
172 TH2D*
histo = corrMats.at(iHisto);
173 histo->GetZaxis()->SetTitle(
"correlation");
174 histo->GetZaxis()->SetRangeUser(0.0,1.0);
193 for (
size_t iHisto = 0; iHisto < corrMatDiffs.size(); ++iHisto)
195 TH2D*
histo = corrMatDiffs.at(iHisto);
204 unsigned long long numValidBins = 0;
223 meanDiff /= numValidBins;
228 const double maxDiffRounded =
static_cast<int>(fabs(maxDiff)*20+1)/20.;
229 histo->GetZaxis()->SetRangeUser(-maxDiffRounded,maxDiffRounded);
235 histo->GetZaxis()->SetTitle(
"correlation difference");
238 DrawLabels(
histo,fixedValue1,fixedValue2,diffProviders.at(iHisto).first,diffProviders.at(iHisto).second,meanDiff,maxDiff,maxX,maxY);
245 if (outHistFile && corrMatDiffs.size() == 1)
254 void MakeCorrelationPlots(
const TString&
outFile, TCanvas*
canvas, TFile* outHistFile,
const std::vector<JetUncertaintiesTool*>& providers,
const std::vector< std::pair<double,double> >& fixedEta,
const std::vector< std::pair<double,double> >& fixedPt)
256 const TString
release = providers.at(0)->getRelease();
257 const double minPtVal =
release.EndsWith(
"TLA") ? 75 : 17;
258 const double maxPtVal =
testJER ? 1500 :
release.EndsWith(
"TLA") ? 1000 : (
release.BeginsWith(
"2011") ||
release.BeginsWith(
"2012")) ? 2500 : 3000;
259 const double minEtaVal = 0;
260 const double maxEtaVal =
release.EndsWith(
"TLA") ? 2.8 : 4.5;
265 for (
size_t iFixed = 0; iFixed < fixedEta.size(); ++iFixed)
267 const double etaVal1 = fixedEta.at(iFixed).first;
268 const double etaVal2 = fixedEta.at(iFixed).second;
270 printf(
"Processing (eta1,eta2) = (%.2f,%.2f)\n",etaVal1,etaVal2);
273 std::vector<TH2D*> corrPt;
274 for (
size_t iProv = 0; iProv < providers.size(); ++iProv)
276 corrPt.push_back(providers.at(iProv)->getPtCorrelationMatrix(100,minPtVal,maxPtVal,etaVal1,etaVal2));
277 corrPt.back()->GetXaxis()->SetTitle(
"#it{p}_{T}^{jet} [GeV]");
278 corrPt.back()->GetYaxis()->SetTitle(corrPt.back()->GetXaxis()->GetTitle());
282 for (
int binX = 1;
binX <= corrPt.back()->GetNbinsX(); ++
binX)
283 for (
int binY = 1;
binY <= corrPt.back()->GetNbinsY(); ++
binY)
284 if (corrPt.back()->GetBinContent(
binX,
binY) > -100) {
285 corrPt.back()->SetBinContent(
binX,
binY,1-corrPt.back()->GetBinContent(
binX,
binY));
297 for (
size_t iFixed = 0; iFixed < fixedPt.size(); ++iFixed)
299 const double ptVal1 = fixedPt.at(iFixed).first;
300 const double ptVal2 = fixedPt.at(iFixed).second;
302 printf(
"Processing (pT1,pT2) = (%.0f,%.0f)\n",ptVal1,ptVal2);
305 std::vector<TH2D*> corrEta;
306 for (
size_t iProv = 0; iProv < providers.size(); ++iProv)
308 corrEta.push_back(providers.at(iProv)->getEtaCorrelationMatrix(
round(10*(maxEtaVal-minEtaVal)),minEtaVal,maxEtaVal,ptVal1,ptVal2));
309 corrEta.back()->GetXaxis()->SetTitle(
"#eta^{jet}");
310 corrEta.back()->GetYaxis()->SetTitle(corrEta.back()->GetXaxis()->GetTitle());
314 for (
int binX = 1;
binX <= corrEta.back()->GetNbinsX(); ++
binX)
315 for (
int binY = 1;
binY <= corrEta.back()->GetNbinsY(); ++
binY)
316 if (corrEta.back()->GetBinContent(
binX,
binY) > -100)
317 corrEta.back()->SetBinContent(
binX,
binY,1-corrEta.back()->GetBinContent(
binX,
binY));
331 printf(
"Expected arguments:\n");
332 printf(
"\t1. Output file (ps, pdf, eps)\n");
333 printf(
"\t2. Jet definition(s), semicolon delimited\n");
334 printf(
"\t3. MC type\n");
335 printf(
"\t4. Config file(s), semicolon delimited\n");
336 printf(
"\t\tDifferences will be performed by left minus right\n");
337 printf(
"\t\tExample: A;B;C --> A,B,C,A-B,B-C\n");
338 printf(
"\t5. Fixed eta values to consider for pT correlations, semicolon delimited\n");
339 printf(
"\t Note: Comma-delimited pairs if specifying both eta1 and eta2\n");
340 printf(
"\t6. Fixed pT values to consider for eta correlations, semicolon delimited\n");
341 printf(
"\t Note: Comma-delimited pairs if specifying both pt1 and pt2\n");
342 printf(
"\t Note: Expected units are GeV, not MeV\n");
343 printf(
"\t7. Output histogram root file (OPTIONAL)\n");
344 printf(
"\t8. Whether to switch to relative differences (OPTIONAL, default \"false\")\n");
347 TFile* outHistFile =
nullptr;
349 std::vector<TString>
jetDefs = jet::utils::vectorize<TString>(
argv[2],
";");
350 TString mcType =
argv[3];
351 std::vector<TString> configs = jet::utils::vectorize<TString>(
argv[4],
";");
352 std::vector<TString> fixedEtaS = jet::utils::vectorize<TString>(
argv[5],
";");
353 std::vector<TString> fixedPtS = jet::utils::vectorize<TString>(
argv[6],
";");
355 outHistFile = TFile::Open(
argv[7],
"RECREATE");
359 #ifdef XAOD_STANDALONE
360 StatusCode::enableFailure();
361 #endif // XAOD_STANDALONE
365 printf(
"Failed to find any configs: %s\n",
argv[4]);
370 if (fixedEtaS.empty() && fixedPtS.empty())
372 printf(
"No fixed pT or eta values were specified, nothing to do.\n");
377 std::vector< std::pair<double,double> > fixedEta;
378 std::vector< std::pair<double,double> > fixedPt;
381 for (
size_t iEta1 = 0; iEta1 < fixedEtaS.size(); ++iEta1)
382 for (
size_t iEta2 = 0; iEta2 < fixedEtaS.size(); ++iEta2)
383 fixedEta.emplace_back(jet::utils::getTypeObjFromString<double>(fixedEtaS.at(iEta1)),jet::utils::getTypeObjFromString<double>(fixedEtaS.at(iEta2)));
389 std::vector<double> temp = jet::utils::vectorize<double>(fixedEtaS.at(
iEta),
",");
391 fixedEta.emplace_back(temp.at(0),temp.at(0));
392 else if (temp.size() == 2)
393 fixedEta.emplace_back(temp.at(0),temp.at(1));
396 printf(
"Specified a fixed eta term which is not 1 or 2 values: %s\n",fixedEtaS.at(
iEta).Data());
397 printf(
"Unsure of how to interpret this term, exiting.\n");
405 for (
size_t iPt1 = 0; iPt1 < fixedPtS.size(); ++iPt1)
406 for (
size_t iPt2 = 0; iPt2 < fixedPtS.size(); ++iPt2)
407 fixedPt.emplace_back(jet::utils::getTypeObjFromString<double>(fixedPtS.at(iPt1)),jet::utils::getTypeObjFromString<double>(fixedPtS.at(iPt2)));
411 for (
size_t iPt = 0; iPt < fixedPtS.size(); ++iPt)
413 std::vector<double> temp = jet::utils::vectorize<double>(fixedPtS.at(iPt),
",");
414 if (temp.size() == 1)
415 fixedPt.emplace_back(temp.at(0),temp.at(0));
416 else if (temp.size() == 2)
417 fixedPt.emplace_back(temp.at(0),temp.at(1));
420 printf(
"Specified a fixed pt term which is not 1 or 2 values: %s\n",fixedPtS.at(iPt).Data());
421 printf(
"Unsure of how to interpret this term, exiting.\n");
429 gStyle->SetPalette(1);
430 TCanvas*
canvas =
new TCanvas(
"canvas");
432 canvas->SetFillStyle(4000);
434 canvas->SetFrameBorderMode(0);
443 for (
size_t iJetDef = 0; iJetDef <
jetDefs.size(); ++iJetDef)
448 std::vector<JetUncertaintiesTool*> providers;
449 printf(
"Processing %zu config(s)...\n",configs.size());
450 for (
size_t iConfig = 0; iConfig < configs.size(); ++iConfig)
456 if (providers.back()->setProperty(
"JetDefinition",
jetDef.Data()).isFailure())
458 printf(
"Failed to set JetDefinition to %s\n",
jetDef.Data());
461 if (providers.back()->setProperty(
"MCType",mcType.Data()).isFailure())
463 printf(
"Failed to set MCType to %s\n",mcType.Data());
466 if (providers.back()->setProperty(
"ConfigFile",configs.at(iConfig).Data()).isFailure())
468 printf(
"Failed to set ConfigFile to %s\n",configs.at(iConfig).Data());
484 if (providers.back()->setScaleToGeV().isFailure())
486 printf(
"Failed to set tool scale to GeV for config: %s\n",configs.at(iConfig).Data());
491 if (providers.back()->initialize().isFailure())
493 printf(
"Failed to initialize tool for config: %s\n",configs.at(iConfig).Data());
502 for (
size_t iProv = 0; iProv < providers.size(); ++iProv)
503 delete providers.at(iProv);
514 outHistFile->Close();
515 outHistFile =
nullptr;