MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
PlotLLH.cpp
Go to the documentation of this file.
1// C++
2#include <algorithm>
3#include <iomanip>
4
5// MACH3 PLOTTING
7
8//this file has lots of usage of the ROOT plotting interface that only takes floats, turn this warning off for this CU for now
9#pragma GCC diagnostic ignored "-Wfloat-conversion"
10#pragma GCC diagnostic ignored "-Wconversion"
11
14
15// some options for the plots
22
24
26
27void getSplitSampleStack(int fileIdx, std::string parameterName, TH1D LLH_allSams,
28 std::vector<float> &cumSums, std::vector<bool> &drawLabel,
29 THStack *sampleStack, TLegend *splitSamplesLegend,
30 float baselineLLH_main = 0.00001)
31 {
32 std::vector<std::string> sampNames = man->input().getTaggedSamples(man->getOption<std::vector<std::string>>("sampleTags"));
33 size_t nSamples = sampNames.size();
34
35 cumSums.resize(nSamples);
36 drawLabel.resize(nSamples);
37
38 float LLH_main_integ = LLH_allSams.Integral();
39 float cumSum = 0.0;
40 int nBins = LLH_allSams.GetNbinsX();
41
42 MACH3LOG_DEBUG("getting split sample THStack for {} known samples", nSamples);
43 for (uint i = 0; i < nSamples; i++)
44 {
45 std::string sampName = sampNames[i];
46
47 MACH3LOG_DEBUG(" on sample {}/{}: {}", i, nSamples, sampName);
48
49 TH1D *LLH_indivSam =
50 new TH1D(man->input().getSampleSpecificLLHScan_TH1D(fileIdx, parameterName, sampName));
51 LLH_indivSam->SetName(Form("%i_%s_%s", fileIdx, parameterName.c_str(), sampName.c_str()));
52 LLH_indivSam->SetBit(kCanDelete);
53
54 // the hist for this sample didn't exist so move on
55 if (LLH_indivSam->GetNbinsX() == 1)
56 {
57 MACH3LOG_DEBUG(" sample hist had only 1 bin - assuming it doesn't exist");
58 delete LLH_indivSam;
59 drawLabel[i] = false;
60 cumSums[i] = -999.9;
61 continue;
62 }
63
64 LLH_indivSam->SetStats(0);
65 LLH_indivSam->SetLineColor(TColor::GetColorPalette(
66 floor(static_cast<float>(i) * TColor::GetNumberOfColors() / static_cast<float>(nSamples))));
67 LLH_indivSam->SetFillColor(TColor::GetColorPalette(
68 floor(static_cast<float>(i) * TColor::GetNumberOfColors() / static_cast<float>(nSamples))));
69 sampleStack->Add(LLH_indivSam);
70 splitSamplesLegend->AddEntry(LLH_indivSam, man->style().prettifySampleName(sampName).c_str(),
71 "lf");
72
73 float lastBinLLH = LLH_indivSam->GetBinContent(nBins);
74 cumSum += lastBinLLH;
75
76 MACH3LOG_DEBUG(" Last bin LLH = {} :: cumulative LLH = {}", lastBinLLH, cumSum);
77 MACH3LOG_DEBUG(" LLH fraction = {} / {} = {}", LLH_indivSam->Integral(), LLH_main_integ, LLH_indivSam->Integral() / LLH_main_integ);
78
79 cumSums[i] = cumSum;
80
81 // dont draw a label if the likelihood contribution is less than threshold%
82 if ((LLH_indivSam->Integral() / LLH_main_integ > sampleLabelThreshold) &&
83 (LLH_indivSam->Integral() / baselineLLH_main > sampleLabelThreshold))
84 {
85 drawLabel[i] = true;
86 } else
87 {
88 drawLabel[i] = false;
89 }
90 MACH3LOG_DEBUG(" drawLabel = {}", drawLabel.back());
92 }
93}
94
95// handy dandy helper function for drawing and nicely formatting a stack of ratio plots
96void drawRatioStack(THStack *ratioCompStack) {
97 // do this so 1.0 is in the middle of the plot vertically
98 double stackMax = ratioCompStack->GetMaximum("NOSTACK");
99 double stackMin = ratioCompStack->GetMinimum("NOSTACK");
100
101 double stackLim = std::max(std::abs(1.0 - stackMax), std::abs(1.0 - stackMin));
102
103 ratioCompStack->SetMinimum(1.0 - 1.05 * stackLim);
104 ratioCompStack->SetMaximum(1.0 + 1.05 * stackLim);
105
106 // draw it
107 ratioCompStack->Draw(Form("NOSTACK%s", man->getDrawOptions().c_str()));
108
109 // make it look a bit nicer
110 ratioCompStack->GetXaxis()->SetLabelSize(ratioLabelScaling *
111 ratioCompStack->GetXaxis()->GetLabelSize());
112 ratioCompStack->GetXaxis()->SetTitleSize(ratioLabelScaling *
113 ratioCompStack->GetXaxis()->GetLabelSize());
114 ratioCompStack->GetXaxis()->SetTitle("Parameter Variation");
115
116 ratioCompStack->GetYaxis()->SetLabelSize(ratioLabelScaling *
117 ratioCompStack->GetYaxis()->GetLabelSize());
118 ratioCompStack->GetYaxis()->SetTitleSize(ratioLabelScaling *
119 ratioCompStack->GetYaxis()->GetLabelSize());
120 ratioCompStack->GetYaxis()->SetTitleOffset(yTitleOffset);
121 ratioCompStack->GetYaxis()->SetNdivisions(5, 2, 0);
122}
123
124void makeLLHScanComparisons(std::string paramName, std::string LLHType, std::string outputFileName,
125 TCanvas *canv, TPad *LLHPad, TPad *ratioPad) {
126 // will use these to make comparisons of likelihoods
127 THStack *compStack = new THStack("LLH_Stack", "");
128 THStack *ratioCompStack = new THStack("LLH_Ratio_Stack", "");
129 auto legend = std::make_unique<TLegend>(0.3, 0.6, 0.7, 0.8);
130
131 // get the sample reweight hist from the main file
132 TH1D LLH_main = man->input().getLLHScan_TH1D(0, paramName, LLHType);
133 LLH_main.SetStats(0);
134
135 LLH_main.SetLineColor(kBlack);
136 compStack->Add(&LLH_main);
137 legend->AddEntry(&LLH_main, man->getFileLabel(0).c_str(), "l");
138
139 int nBins = LLH_main.GetNbinsX();
140
141 // go through the other files
142 for (unsigned int extraFileIdx = 1; extraFileIdx < man->input().getNInputFiles(); extraFileIdx++)
143 {
144 TH1D *compHist = new TH1D(man->input().getLLHScan_TH1D(extraFileIdx, paramName, LLHType));
145 compHist->SetBit(kCanDelete); // <- will allow this to be deleted by root once done plotting
146 if (compHist->GetNbinsX() == 0)
147 continue;
148
149 // make them look different to each other
150 compHist->SetLineColor(
151 TColor::GetColorPalette(floor(static_cast<float>(extraFileIdx) * TColor::GetNumberOfColors() /
152 static_cast<float>(man->input().getNInputFiles()))));
153 compHist->SetLineStyle(2 + extraFileIdx % 9);
154 compHist->SetLineWidth(lineWidth);
155
156 TH1D *divHist = static_cast<TH1D*>(compHist->Clone(Form("RatioHist_%i", extraFileIdx)));
157 divHist->SetBit(kCanDelete);
158
159 if (man->getPlotRatios())
160 {
161 // get the ratio hist
162 divHist->Divide(compHist, &LLH_main);
163 ratioCompStack->Add(divHist);
164 }
165
166 // add it to the comparisson hstack and legend
167 compStack->Add(compHist);
168 legend->AddEntry(compHist, man->getFileLabel(extraFileIdx).c_str(), "l");
169 }
170
171 // draw the log likelihoods
172 canv->cd();
173 canv->Draw();
174 LLHPad->Draw();
175 LLHPad->cd();
176 if (man->getDrawGrid())
177 LLHPad->SetGrid();
178 compStack->Draw(Form("NOSTACK%s", man->getDrawOptions().c_str()));
179 if (!man->getPlotRatios())
180 compStack->GetXaxis()->SetTitle("Parameter Variation");
181 compStack->GetYaxis()->SetTitle(Form("-2LLH_{%s}", LLHType.c_str()));
182 compStack->SetTitle(man->style().prettifyParamName(paramName).c_str());
183 compStack->GetYaxis()->SetTitleOffset(yTitleOffset);
184 legend->Draw();
185
186 // add the ratio plot if specified
187 if (man->getPlotRatios())
188 {
189 canv->cd();
190 ratioPad->Draw();
191 ratioPad->cd();
192 if (man->getDrawGrid())
193 ratioPad->SetGrid();
194
195 drawRatioStack(ratioCompStack);
196
197 // add horizontal line at 1 for reference
198 TLine line = TLine();
199 line.SetLineColor(kBlack);
200 line.SetLineWidth(lineWidth);
201 line.DrawLine(LLH_main.GetBinLowEdge(1), 1.0, LLH_main.GetBinLowEdge(nBins + 1), 1.0);
202 }
203
204 // save to the output file
205 canv->SaveAs(outputFileName.c_str());
206
207 // will use these to make comparisons of likelihoods
208 delete compStack;
209 delete ratioCompStack;
210}
211
212void makeSplitSampleLLHScanComparisons(std::string paramName, std::string outputFileName,
213 TCanvas *canv, TPad *LLHPad, TPad *ratioPad) {
214 MACH3LOG_DEBUG(" Making split sample LLH comparison");
215 canv->Clear();
216 canv->Draw();
217
218 canv->Divide(man->getNFiles());
219
220 // get the sample hist from the main file
221 TH1D LLH_main = man->input().getLLHScan_TH1D(0, paramName, "sample");
222 if (LLH_main.GetNbinsX() == 1)
223 {
224 MACH3LOG_DEBUG(" Main LLH had only 1 bin, assuming it doesn't exist");
225 return;
226 }
227
228 THStack *baseSplitSamplesStack = new THStack(
229 paramName.c_str(), Form("%s - %s", paramName.c_str(), man->getFileLabel(0).c_str()));
230 TLegend *baseSplitSamplesLegend = new TLegend(0.37, 0.475, 0.63, 0.9);
231
232 if (man->getDrawGrid())
233 canv->cd(1)->SetGrid();
234 else
235 canv->cd(1);
236
237 std::vector<float> cumSums;
238 std::vector<bool> drawLabel;
239
240 getSplitSampleStack(0, paramName, LLH_main, cumSums, drawLabel, baseSplitSamplesStack,
241 baseSplitSamplesLegend);
242
243 baseSplitSamplesStack->Draw(man->getDrawOptions().c_str());
244
246 {
247 LLH_main.SetLineWidth(1); // undo SetLineWidth that was done above
248 LLH_main.Draw(Form("same%s", man->getDrawOptions().c_str()));
249 baseSplitSamplesLegend->AddEntry(&LLH_main, "All Samples", "l");
250 }
251
252 baseSplitSamplesLegend->Draw();
253
254 TLatex *label = new TLatex;
255 // format the label
256 label->SetTextAlign(11);
257 label->SetTextAngle(-55);
258 label->SetTextSize(0.012);
259
260 // need to draw the labels after other stuff or they dont show up
261 for (uint i = 0; i < man->input().getTaggedSamples(man->getOption<std::vector<std::string>>("sampleTags")).size(); i++)
262 {
263 MACH3LOG_DEBUG(" Will I draw the label for sample {}??", i);
264 std::string sampName = man->input().getTaggedSamples(man->getOption<std::vector<std::string>>("sampleTags"))[i];
265 if (!drawLabel[i])
266 {
267 MACH3LOG_DEBUG(" - Not drawing label");
268 continue;
269 }
270 MACH3LOG_DEBUG(" - Drawing label");
271
272 label->DrawLatex(LLH_main.GetBinLowEdge(LLH_main.GetNbinsX() + 1), cumSums[i],
273 Form("#leftarrow%s", man->style().prettifySampleName(sampName).c_str()));
274 MACH3LOG_DEBUG(" I drew the label!");
275 }
276
277 // now we plot the comparisson file plots
278 for (unsigned int extraFileIdx = 1; extraFileIdx < man->getNFiles(); extraFileIdx++)
279 {
280 MACH3LOG_DEBUG(" - Adding plot for additional file {}", extraFileIdx);
281 canv->cd(1 + extraFileIdx);
282
283 std::vector<float> extraCumSums;
284 std::vector<bool> extraDrawLabel;
285
286 THStack *splitSamplesStack =
287 new THStack(paramName.c_str(),
288 Form("%s - %s", paramName.c_str(), man->getFileLabel(extraFileIdx).c_str()));
289 TLegend *splitSamplesLegend = new TLegend(0.37, 0.475, 0.63, 0.9);
290
291 splitSamplesStack->SetBit(kCanDelete);
292 splitSamplesLegend->SetBit(kCanDelete);
293
294 TH1D compLLH_main = man->input().getLLHScan_TH1D(extraFileIdx, paramName, "sample");
295 if (compLLH_main.GetNbinsX() == 1)
296 {
297 delete splitSamplesStack;
298 delete splitSamplesLegend;
299 continue;
300 }
301
302 // if on the same y axis, also check that the contribution of the sample compared to the
303 // baseline LLH integral is above the threshold otherwise the labels might get very crowded if
304 // the comparisson LLH is much smaller than the baseline one
305 if (sameAxis)
306 getSplitSampleStack(extraFileIdx, paramName, compLLH_main, extraCumSums, extraDrawLabel,
307 splitSamplesStack, splitSamplesLegend, LLH_main.Integral());
308 else
309 getSplitSampleStack(extraFileIdx, paramName, compLLH_main, extraCumSums, extraDrawLabel,
310 splitSamplesStack, splitSamplesLegend);
311
312 // go to the pad for the histograms
313 if (man->getDrawGrid())
314 LLHPad->SetGrid();
315 LLHPad->Draw();
316 LLHPad->cd();
317
318 splitSamplesStack->Draw(man->getDrawOptions().c_str());
319 splitSamplesStack->GetXaxis()->SetTitle("Parameter Variation");
320 splitSamplesStack->GetYaxis()->SetTitle("-2LLH_{sam}");
321 if (sameAxis)
322 splitSamplesStack->SetMaximum(baseSplitSamplesStack->GetMaximum());
323
325 {
326 compLLH_main.Draw(Form("same%s", man->getDrawOptions().c_str()));
327 compLLH_main.SetLineColor(kBlack);
328 splitSamplesLegend->AddEntry(&compLLH_main, "All Samples", "l");
329 }
330 splitSamplesLegend->Draw();
331
332 // need to draw the labels after other stuff or they dont show up
333 for (uint i = 0; i < man->input().getTaggedSamples(man->getOption<std::vector<std::string>>("sampleTags")).size(); i++)
334 {
335 std::string sampName = man->input().getTaggedSamples(man->getOption<std::vector<std::string>>("sampleTags"))[i];
336 if (!drawLabel[i])
337 continue;
338 label->DrawLatex(compLLH_main.GetBinLowEdge(compLLH_main.GetNbinsX() + 1), extraCumSums[i],
339 Form("#leftarrow%s", man->style().prettifySampleName(sampName).c_str()));
340 }
341
342 if (man->getPlotRatios())
343 {
344 THStack *splitSamplesStackRatios = new THStack(paramName.c_str(), "");
345
346 TList *baselineHistList = baseSplitSamplesStack->GetHists();
347 TList *compHistList = splitSamplesStack->GetHists();
348
349 for (uint sampleIdx = 0; sampleIdx < man->input().getTaggedSamples(man->getOption<std::vector<std::string>>("sampleTags")).size(); sampleIdx++)
350 {
351 TH1D *divHist = new TH1D(
352 Form("%s_%s_splitDiv_%i", paramName.c_str(), man->getFileLabel(extraFileIdx).c_str(),
353 sampleIdx),
354 Form("%s_%s_splitDiv", paramName.c_str(), man->getFileLabel(extraFileIdx).c_str()),
355 compLLH_main.GetNbinsX(), compLLH_main.GetBinLowEdge(1),
356 compLLH_main.GetBinLowEdge(compLLH_main.GetNbinsX() + 1));
357 divHist->Divide(static_cast<TH1D*>(compHistList->At(sampleIdx)),
358 static_cast<TH1D*>(baselineHistList->At(sampleIdx)));
359 splitSamplesStackRatios->Add(divHist);
360 divHist->SetLineColor((static_cast<TH1D*>(compHistList->At(sampleIdx))->GetLineColor()));
361 }
362
363 canv->cd(2 + extraFileIdx);
364 if (man->getDrawGrid())
365 ratioPad->SetGrid();
366 ratioPad->Draw();
367 ratioPad->cd();
368
369 drawRatioStack(splitSamplesStackRatios);
370 }
371 }
372
373 canv->SaveAs(outputFileName.c_str());
374
375 delete label;
376 delete baseSplitSamplesStack;
377 delete baseSplitSamplesLegend;
378}
379
380int PlotLLH() {
381
382 man->style().setPalette(man->getOption<std::string>("colorPalette"));
383
384 // make the canvas and other plotting stuff
385 TCanvas *canv = new TCanvas("canv", "", 1024, 1024);
386 gStyle->SetOptTitle(2);
387 // split up the canvas so can have side by side plots, one for each file
388 TCanvas *splitSamplesCanv = new TCanvas("splitSampCanv", "", 4096 * man->getNFiles(), 4096);
389
390 // split the canvas if plotting ratios
391 TPad *LLHPad, *ratioPad;
392 if (man->getPlotRatios())
393 {
394 LLHPad = new TPad("LLHPad", "LLHPad", 0.0, ratioPlotSplit, 1.0, 1.0);
395 LLHPad->SetBottomMargin(0.0);
396 ratioPad = new TPad("ratioPad", "ratioPad", 0.0, 0.0, 1.0, ratioPlotSplit);
397 ratioPad->SetTopMargin(0.0);
398 ratioPad->SetBottomMargin(0.3);
399 } else
400 {
401 LLHPad = new TPad("AllSampPad", "AllSampPad", 0.0, 0.0, 1.0, 1.0);
402 ratioPad = new TPad("AllSampRatioPad", "AllSampRatioPad", 0.0, 0.0, 0.0, 0.0);
403 }
404
405 canv->SaveAs((man->getOutputName("_Sample") + "[").c_str());
406 canv->SaveAs((man->getOutputName("_Penalty") + "[").c_str());
407 canv->SaveAs((man->getOutputName("_Total") + "[").c_str());
408
409 if (man->getSplitBySample())
410 canv->SaveAs((man->getOutputName("_bySample") + "[").c_str());
411
412 for( std::string par: man->getOption<std::vector<std::string>>("parameterTags")) std::cout << par << ", ";
413 // loop over the spline parameters
414 for (std::string paramName : man->input().getTaggedParameters(man->getOption<std::vector<std::string>>("parameterTags")))
415 {
416 MACH3LOG_DEBUG("working on parameter {}", paramName);
417 // ###############################################################
418 // First lets do just the straight up likelihoods from all samples
419 // ###############################################################
420
421 makeLLHScanComparisons(paramName, "sample", man->getOutputName("_Sample"), canv, LLHPad,
422 ratioPad);
423 makeLLHScanComparisons(paramName, "penalty", man->getOutputName("_Penalty"), canv, LLHPad,
424 ratioPad);
425 makeLLHScanComparisons(paramName, "total", man->getOutputName("_Total"), canv, LLHPad,
426 ratioPad);
427
428 // #########################################
429 // ## now lets make plots split by sample ##
430 // #########################################
431 if (man->getSplitBySample())
432 {
433 makeSplitSampleLLHScanComparisons(paramName, man->getOutputName("_bySample"),
434 splitSamplesCanv, LLHPad, ratioPad);
435 }
436 }
437
438 canv->SaveAs((man->getOutputName("_Sample") + "]").c_str());
439 canv->SaveAs((man->getOutputName("_Penalty") + "]").c_str());
440 canv->SaveAs((man->getOutputName("_Total") + "]").c_str());
441 if (man->getSplitBySample())
442 canv->SaveAs((man->getOutputName("_bySample") + "]").c_str());
443
444 delete LLHPad;
445 delete ratioPad;
446 delete canv;
447 if(man != nullptr) delete man;
448
449 return 0;
450}
451
452int main(int argc, char **argv) {
454
456 man->parseInputs(argc, argv);
457
458 man->setExec("PlotLLH");
459
460 ratioPlotSplit = man->getOption<double>("ratioPlotSplit");
461 yTitleOffset = man->getOption<double>("yTitleOffset");
462 sampleLabelThreshold = man->getOption<double>("sampleLabelThreshold");
463 lineWidth = man->getOption<int>("lineWidth");
464 totalOnSplitPlots = man->getOption<bool>("totalOnSplitPlots");
465 sameAxis = man->getOption<bool>("sameAxis");
466
467 // scale the ratio plot labels by this much to make them same size as the normal plot
468 ratioLabelScaling = (1.0 / ratioPlotSplit - 1.0);
469
470 return PlotLLH();
471}
TCanvas * canv
#define MACH3LOG_DEBUG
Definition: MaCh3Logger.h:22
void SetMaCh3LoggerFormat()
Set messaging format of the logger.
Definition: MaCh3Logger.h:30
bool sameAxis
Definition: PlotLLH.cpp:21
double ratioLabelScaling
Definition: PlotLLH.cpp:23
void drawRatioStack(THStack *ratioCompStack)
Definition: PlotLLH.cpp:96
int main(int argc, char **argv)
Definition: PlotLLH.cpp:452
double ratioPlotSplit
Definition: PlotLLH.cpp:16
double sampleLabelThreshold
Definition: PlotLLH.cpp:18
int PlotLLH()
Definition: PlotLLH.cpp:380
bool totalOnSplitPlots
Definition: PlotLLH.cpp:20
void makeLLHScanComparisons(std::string paramName, std::string LLHType, std::string outputFileName, TCanvas *canv, TPad *LLHPad, TPad *ratioPad)
Definition: PlotLLH.cpp:124
int lineWidth
Definition: PlotLLH.cpp:19
double yTitleOffset
Definition: PlotLLH.cpp:17
MaCh3Plotting::PlottingManager * man
Definition: PlotLLH.cpp:25
void makeSplitSampleLLHScanComparisons(std::string paramName, std::string outputFileName, TCanvas *canv, TPad *LLHPad, TPad *ratioPad)
Definition: PlotLLH.cpp:212
void getSplitSampleStack(int fileIdx, std::string parameterName, TH1D LLH_allSams, std::vector< float > &cumSums, std::vector< bool > &drawLabel, THStack *sampleStack, TLegend *splitSamplesLegend, float baselineLLH_main=0.00001)
Definition: PlotLLH.cpp:27
TH1D getLLHScan_TH1D(int fileNum, std::string paramName, std::string LLHType) const
Definition: inputManager.h:366
std::vector< std::string > getTaggedParameters(const std::vector< std::string > &tags, std::string checkType="all") const
Get all parameters which have some set of tags.
Definition: inputManager.h:480
TH1D getSampleSpecificLLHScan_TH1D(int fileNum, std::string paramName, std::string sample) const
Definition: inputManager.h:403
std::vector< std::string > getTaggedSamples(const std::vector< std::string > &tags, std::string checkType="all") const
Get all samples which have some set of tags.
Definition: inputManager.h:491
size_t getNInputFiles() const
Definition: inputManager.h:471
The main class to be used in plotting scripts.
const std::string getOutputName()
Get the straight up output file name with no bells or whistles, just the file extension.
void setExec(std::string execName)
Internally set the name of the executable that manager is being used in.
void parseInputs(int argc, char *const *argv)
Parse command line arguments.
const std::string getDrawOptions()
const std::string getFileLabel(int i)
const StyleManager & style()
Get the StyleManager contained within this PlottingManager, for doing style related things.
T getOption(std::string option)
Get a specific option from the config for this executable.
const InputManager & input()
Get the InputManager contained within this PlottingManager, for doing input related things.
std::string prettifyParamName(const std::string &origName) const
Convert hideous and vulgar internal parameter name into a beautiful presentable name.
Definition: styleManager.h:44
void setPalette(int rootPlotStyle) const
Set the root colour palette to one of the default root pallettes as defined in (root docs)[https://ro...
std::string prettifySampleName(const std::string &origName) const
Convert hideous and vulgar internal sample name into a beautiful presentable name.
Definition: styleManager.h:53