MaCh3  2.4.2
Reference Guide
Functions | Variables
PlotLLH.cpp File Reference
#include <algorithm>
#include <iomanip>
#include "PlottingUtils/PlottingManager.h"
Include dependency graph for PlotLLH.cpp:

Go to the source code of this file.

Functions

void SetTPads (TPad *&LLHPad, TPad *&ratioPad)
 TPad is SUPER FRAGILE, it is safer to just make raw pointer, while ROOT behave weirdly with smart pointers. More...
 
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)
 
void drawRatioStack (THStack *ratioCompStack)
 
void makeLLHScanComparisons (const std::string &paramName, const std::string &LLHType, const std::string &outputFileName, const std::unique_ptr< TCanvas > &canv)
 
void makeSplitSampleLLHScanComparisons (const std::string &paramName, const std::string &outputFileName, const std::unique_ptr< TCanvas > &canv)
 
int PlotLLH ()
 
int main (int argc, char **argv)
 

Variables

double ratioPlotSplit
 
double yTitleOffset
 
double sampleLabelThreshold
 
int lineWidth
 
bool totalOnSplitPlots
 
bool sameAxis
 
double ratioLabelScaling
 
MaCh3Plotting::PlottingManager * PlotMan
 

Detailed Description

Author
Ewan Miller

Definition in file PlotLLH.cpp.

Function Documentation

◆ drawRatioStack()

void drawRatioStack ( THStack *  ratioCompStack)

Definition at line 117 of file PlotLLH.cpp.

117  {
118  // do this so 1.0 is in the middle of the plot vertically
119  double stackMax = ratioCompStack->GetMaximum("NOSTACK");
120  double stackMin = ratioCompStack->GetMinimum("NOSTACK");
121 
122  double stackLim = std::max(std::abs(1.0 - stackMax), std::abs(1.0 - stackMin));
123 
124  ratioCompStack->SetMinimum(1.0 - 1.05 * stackLim);
125  ratioCompStack->SetMaximum(1.0 + 1.05 * stackLim);
126 
127  // draw it
128  ratioCompStack->Draw(Form("NOSTACK%s", PlotMan->getDrawOptions().c_str()));
129 
130  // make it look a bit nicer
131  ratioCompStack->GetXaxis()->SetLabelSize(ratioLabelScaling *
132  ratioCompStack->GetXaxis()->GetLabelSize());
133  ratioCompStack->GetXaxis()->SetTitleSize(ratioLabelScaling *
134  ratioCompStack->GetXaxis()->GetLabelSize());
135  ratioCompStack->GetXaxis()->SetTitle("Parameter Variation");
136 
137  ratioCompStack->GetYaxis()->SetLabelSize(ratioLabelScaling *
138  ratioCompStack->GetYaxis()->GetLabelSize());
139  ratioCompStack->GetYaxis()->SetTitleSize(ratioLabelScaling *
140  ratioCompStack->GetYaxis()->GetLabelSize());
141  ratioCompStack->GetYaxis()->SetTitleOffset(yTitleOffset);
142  ratioCompStack->GetYaxis()->SetNdivisions(5, 2, 0);
143 }
double ratioLabelScaling
Definition: PlotLLH.cpp:24
MaCh3Plotting::PlottingManager * PlotMan
Definition: PlotLLH.cpp:27
double yTitleOffset
Definition: PlotLLH.cpp:18

◆ getSplitSampleStack()

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 at line 50 of file PlotLLH.cpp.

54  {
55  std::vector<std::string> sampNames = PlotMan->input().getTaggedSamples(PlotMan->getOption<std::vector<std::string>>("sampleTags"));
56  size_t nSamples = sampNames.size();
57 
58  cumSums.resize(nSamples);
59  drawLabel.resize(nSamples);
60 
61  float LLH_main_integ = LLH_allSams.Integral();
62  float cumSum = 0.0;
63  int nBins = LLH_allSams.GetNbinsX();
64 
65  MACH3LOG_DEBUG("getting split sample THStack for {} known samples", nSamples);
66  for (uint i = 0; i < nSamples; i++)
67  {
68  std::string sampName = sampNames[i];
69 
70  MACH3LOG_DEBUG(" on sample {}/{}: {}", i, nSamples, sampName);
71 
72  TH1D *LLH_indivSam =
73  new TH1D(PlotMan->input().getSampleSpecificLLHScan_TH1D(fileIdx, parameterName, sampName));
74  LLH_indivSam->SetName(Form("%i_%s_%s", fileIdx, parameterName.c_str(), sampName.c_str()));
75  LLH_indivSam->SetBit(kCanDelete);
76 
77  // the hist for this sample didn't exist so move on
78  if (LLH_indivSam->GetNbinsX() == 1)
79  {
80  MACH3LOG_DEBUG(" sample hist had only 1 bin - assuming it doesn't exist");
81  delete LLH_indivSam;
82  drawLabel[i] = false;
83  cumSums[i] = -999.9;
84  continue;
85  }
86 
87  LLH_indivSam->SetStats(0);
88  LLH_indivSam->SetLineColor(TColor::GetColorPalette(
89  floor(static_cast<float>(i) * TColor::GetNumberOfColors() / static_cast<float>(nSamples))));
90  LLH_indivSam->SetFillColor(TColor::GetColorPalette(
91  floor(static_cast<float>(i) * TColor::GetNumberOfColors() / static_cast<float>(nSamples))));
92  sampleStack->Add(LLH_indivSam);
93  splitSamplesLegend->AddEntry(LLH_indivSam, PlotMan->style().prettifySampleName(sampName).c_str(), "lf");
94 
95  float lastBinLLH = LLH_indivSam->GetBinContent(nBins);
96  cumSum += lastBinLLH;
97 
98  MACH3LOG_DEBUG(" Last bin LLH = {} :: cumulative LLH = {}", lastBinLLH, cumSum);
99  MACH3LOG_DEBUG(" LLH fraction = {} / {} = {}", LLH_indivSam->Integral(), LLH_main_integ, LLH_indivSam->Integral() / LLH_main_integ);
100 
101  cumSums[i] = cumSum;
102 
103  // dont draw a label if the likelihood contribution is less than threshold%
104  if ((LLH_indivSam->Integral() / LLH_main_integ > sampleLabelThreshold) &&
105  (LLH_indivSam->Integral() / baselineLLH_main > sampleLabelThreshold))
106  {
107  drawLabel[i] = true;
108  } else {
109  drawLabel[i] = false;
110  }
111  MACH3LOG_DEBUG(" drawLabel = {}", drawLabel.back());
112  MACH3LOG_DEBUG("");
113  }
114 }
#define MACH3LOG_DEBUG
Definition: MaCh3Logger.h:34
double sampleLabelThreshold
Definition: PlotLLH.cpp:19

◆ main()

int main ( int  argc,
char **  argv 
)

Definition at line 469 of file PlotLLH.cpp.

469  {
472 
473  PlotMan = new MaCh3Plotting::PlottingManager();
474  PlotMan->parseInputs(argc, argv);
475 
476  PlotMan->setExec("PlotLLH");
477 
478  ratioPlotSplit = PlotMan->getOption<double>("ratioPlotSplit");
479  yTitleOffset = PlotMan->getOption<double>("yTitleOffset");
480  sampleLabelThreshold = PlotMan->getOption<double>("sampleLabelThreshold");
481  lineWidth = PlotMan->getOption<int>("lineWidth");
482  totalOnSplitPlots = PlotMan->getOption<bool>("totalOnSplitPlots");
483  sameAxis = PlotMan->getOption<bool>("sameAxis");
484 
485  // scale the ratio plot labels by this much to make them same size as the normal plot
486  ratioLabelScaling = (1.0 / ratioPlotSplit - 1.0);
487 
488  if(PlotMan->getPlotRatios() && PlotMan->input().getNInputFiles() == 1){
489  MACH3LOG_ERROR("Can't plot ratio with single file...");
490  throw MaCh3Exception(__FILE__, __LINE__);
491  }
492 
493  return PlotLLH();
494 }
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:37
void SetMaCh3LoggerFormat()
Set messaging format of the logger.
Definition: MaCh3Logger.h:61
bool sameAxis
Definition: PlotLLH.cpp:22
double ratioPlotSplit
Definition: PlotLLH.cpp:17
int PlotLLH()
Definition: PlotLLH.cpp:422
bool totalOnSplitPlots
Definition: PlotLLH.cpp:21
int lineWidth
Definition: PlotLLH.cpp:20
Custom exception class used throughout MaCh3.
void MaCh3Welcome()
KS: Prints welcome message with MaCh3 logo.
Definition: Monitor.cpp:12

◆ makeLLHScanComparisons()

void makeLLHScanComparisons ( const std::string &  paramName,
const std::string &  LLHType,
const std::string &  outputFileName,
const std::unique_ptr< TCanvas > &  canv 
)

Definition at line 145 of file PlotLLH.cpp.

148  {
149  // will use these to make comparisons of likelihoods
150  auto compStack = std::make_unique<THStack>("LLH_Stack", "");
151  auto ratioCompStack = std::make_unique<THStack>("LLH_Ratio_Stack", "");
152  auto legend = std::make_unique<TLegend>(0.3, 0.6, 0.7, 0.8);
153 
154  TPad* LLHPad = nullptr;
155  TPad* ratioPad = nullptr;
156  SetTPads(LLHPad, ratioPad);
157 
158  // get the sample reweight hist from the main file
159  TH1D LLH_main = PlotMan->input().getLLHScan_TH1D(0, paramName, LLHType);
160  LLH_main.SetStats(0);
161 
162  LLH_main.SetLineColor(kBlack);
163  compStack->Add(&LLH_main);
164  legend->AddEntry(&LLH_main, PlotMan->getFileLabel(0).c_str(), "l");
165 
166  int nBins = LLH_main.GetNbinsX();
167 
168  // go through the other files
169  for (unsigned int extraFileIdx = 1; extraFileIdx < PlotMan->input().getNInputFiles(); extraFileIdx++)
170  {
171  int originalErrorLevel = gErrorIgnoreLevel;
172  gErrorIgnoreLevel = kFatal;
173 
174  TH1D *compHist = new TH1D(PlotMan->input().getLLHScan_TH1D(extraFileIdx, paramName, LLHType));
175  compHist->SetName(Form("LLHScan_%s_%s_%d", paramName.c_str(), LLHType.c_str(), extraFileIdx));
176  compHist->SetBit(kCanDelete); // <- will allow this to be deleted by root once done plotting
177  gErrorIgnoreLevel = originalErrorLevel;
178  if (compHist->GetNbinsX() == 0)
179  continue;
180 
181  // make them look different to each other
182  compHist->SetLineColor(
183  TColor::GetColorPalette(floor(static_cast<float>(extraFileIdx) * TColor::GetNumberOfColors() /
184  static_cast<float>(PlotMan->input().getNInputFiles()))));
185  compHist->SetLineStyle(2 + extraFileIdx % 9);
186  compHist->SetLineWidth(lineWidth);
187 
188  TH1D *divHist = static_cast<TH1D*>(compHist->Clone(Form("RatioHist_%i", extraFileIdx)));
189  divHist->SetBit(kCanDelete);
190 
191  if (PlotMan->getPlotRatios())
192  {
193  // get the ratio hist
194  divHist->Divide(compHist, &LLH_main);
195  ratioCompStack->Add(divHist);
196  }
197 
198  // add it to the comparison stack and legend
199  compStack->Add(compHist);
200  legend->AddEntry(compHist, PlotMan->getFileLabel(extraFileIdx).c_str(), "l");
201  }
202 
203  // draw the log likelihoods
204  canv->cd();
205  canv->Draw();
206  LLHPad->Draw();
207  LLHPad->cd();
208  if (PlotMan->getDrawGrid())
209  LLHPad->SetGrid();
210  compStack->Draw(Form("NOSTACK%s", PlotMan->getDrawOptions().c_str()));
211  if (!PlotMan->getPlotRatios())
212  compStack->GetXaxis()->SetTitle("Parameter Variation");
213  compStack->GetYaxis()->SetTitle(Form("-2LLH_{%s}", LLHType.c_str()));
214  compStack->SetTitle(PlotMan->style().prettifyParamName(paramName).c_str());
215  compStack->GetYaxis()->SetTitleOffset(yTitleOffset);
216  legend->Draw();
217 
218  // add the ratio plot if specified
219  if (PlotMan->getPlotRatios())
220  {
221  canv->cd();
222  ratioPad->Draw();
223  ratioPad->cd();
224  if (PlotMan->getDrawGrid())
225  ratioPad->SetGrid();
226 
227  drawRatioStack(ratioCompStack.get());
228 
229  // add horizontal line at 1 for reference
230  TLine line = TLine();
231  line.SetLineColor(kBlack);
232  line.SetLineWidth(lineWidth);
233  line.DrawLine(LLH_main.GetBinLowEdge(1), 1.0, LLH_main.GetBinLowEdge(nBins + 1), 1.0);
234  }
235 
236  // save to the output file
237  canv->SaveAs(outputFileName.c_str());
238  delete LLHPad;
239  delete ratioPad;
240 }
TCanvas * canv
void drawRatioStack(THStack *ratioCompStack)
Definition: PlotLLH.cpp:117
void SetTPads(TPad *&LLHPad, TPad *&ratioPad)
TPad is SUPER FRAGILE, it is safer to just make raw pointer, while ROOT behave weirdly with smart poi...
Definition: PlotLLH.cpp:30

◆ makeSplitSampleLLHScanComparisons()

void makeSplitSampleLLHScanComparisons ( const std::string &  paramName,
const std::string &  outputFileName,
const std::unique_ptr< TCanvas > &  canv 
)

Definition at line 242 of file PlotLLH.cpp.

244  {
245  MACH3LOG_DEBUG(" Making split sample LLH comparison");
246  // split the canvas if plotting ratios
247  TPad* LLHPad = nullptr;
248  TPad* ratioPad = nullptr;
249  SetTPads(LLHPad, ratioPad);
250  canv->Clear();
251  canv->Draw();
252 
253  canv->Divide(PlotMan->getNFiles());
254 
255  // get the sample hist from the main file
256  TH1D LLH_main = PlotMan->input().getLLHScan_TH1D(0, paramName, "sample");
257  if (LLH_main.GetNbinsX() == 1)
258  {
259  MACH3LOG_DEBUG(" Main LLH had only 1 bin, assuming it doesn't exist");
260  return;
261  }
262 
263  auto baseSplitSamplesStack = std::make_unique<THStack>(
264  paramName.c_str(), Form("%s - %s", paramName.c_str(), PlotMan->getFileLabel(0).c_str()));
265 
266  auto baseSplitSamplesLegend = std::make_unique<TLegend>(0.37, 0.475, 0.63, 0.9);
267 
268  if (PlotMan->getDrawGrid())
269  canv->cd(1)->SetGrid();
270  else
271  canv->cd(1);
272 
273  canv->cd(1)->SetLeftMargin(0.15);
274  std::vector<float> cumSums;
275  std::vector<bool> drawLabel;
276 
277  getSplitSampleStack(0, paramName, LLH_main, cumSums, drawLabel, baseSplitSamplesStack.get(),
278  baseSplitSamplesLegend.get());
279  baseSplitSamplesStack->Draw(PlotMan->getDrawOptions().c_str());
280  // KS: Not sure why but need to plot after it's drawn otherwise ROOT throws segfault...
281  baseSplitSamplesStack->GetYaxis()->SetTitle("-2LLH_{sam}");
282  if (PlotMan->getPlotRatios() == false) baseSplitSamplesStack->GetXaxis()->SetTitle("Parameter Variation");
283  gPad->Modified();
284  gPad->Update();
285 
286  if (totalOnSplitPlots)
287  {
288  LLH_main.SetLineWidth(1); // undo SetLineWidth that was done above
289  LLH_main.Draw(Form("same%s", PlotMan->getDrawOptions().c_str()));
290  baseSplitSamplesLegend->AddEntry(&LLH_main, "All Samples", "l");
291  }
292 
293  baseSplitSamplesLegend->Draw();
294 
295  auto label = std::make_unique<TLatex>();
296  // format the label
297  label->SetTextAlign(11);
298  label->SetTextAngle(-55);
299  label->SetTextSize(0.012);
300 
301  // need to draw the labels after other stuff or they don't show up
302  for (uint i = 0; i < PlotMan->input().getTaggedSamples(PlotMan->getOption<std::vector<std::string>>("sampleTags")).size(); i++)
303  {
304  MACH3LOG_DEBUG(" Will I draw the label for sample {}??", i);
305  std::string sampName = PlotMan->input().getTaggedSamples(PlotMan->getOption<std::vector<std::string>>("sampleTags"))[i];
306  if (!drawLabel[i])
307  {
308  MACH3LOG_DEBUG(" - Not drawing label");
309  continue;
310  }
311  MACH3LOG_DEBUG(" - Drawing label");
312 
313  label->DrawLatex(LLH_main.GetBinLowEdge(LLH_main.GetNbinsX() + 1), cumSums[i],
314  Form("#leftarrow%s", PlotMan->style().prettifySampleName(sampName).c_str()));
315  MACH3LOG_DEBUG(" I drew the label!");
316  }
317 
318  // now we plot the comparison file plots
319  for (unsigned int extraFileIdx = 1; extraFileIdx < PlotMan->getNFiles(); extraFileIdx++)
320  {
321  MACH3LOG_DEBUG(" - Adding plot for additional file {}", extraFileIdx);
322  canv->cd(1 + extraFileIdx);
323 
324  std::vector<float> extraCumSums;
325  std::vector<bool> extraDrawLabel;
326 
327  THStack *splitSamplesStack =
328  new THStack(paramName.c_str(),
329  Form("%s - %s", paramName.c_str(), PlotMan->getFileLabel(extraFileIdx).c_str()));
330 
331  auto splitSamplesLegend = std::make_unique<TLegend>(0.37, 0.475, 0.63, 0.9);
332  splitSamplesStack->SetBit(kCanDelete);
333  splitSamplesLegend->SetBit(kCanDelete);
334 
335  int originalErrorLevel = gErrorIgnoreLevel;
336  gErrorIgnoreLevel = kFatal;
337  TH1D compLLH_main = PlotMan->input().getLLHScan_TH1D(extraFileIdx, paramName, "sample");
338  compLLH_main.SetName((compLLH_main.GetName() + std::string("_file_") + extraFileIdx).Data());
339  gErrorIgnoreLevel = originalErrorLevel;
340  if (compLLH_main.GetNbinsX() == 1)
341  {
342  delete splitSamplesStack;
343  continue;
344  }
345 
346  // if on the same y axis, also check that the contribution of the sample compared to the
347  // baseline LLH integral is above the threshold otherwise the labels might get very crowded if
348  // the comparisson LLH is much smaller than the baseline one
349  if (sameAxis)
350  getSplitSampleStack(extraFileIdx, paramName, compLLH_main, extraCumSums, extraDrawLabel,
351  splitSamplesStack, splitSamplesLegend.get(), LLH_main.Integral());
352  else
353  getSplitSampleStack(extraFileIdx, paramName, compLLH_main, extraCumSums, extraDrawLabel,
354  splitSamplesStack, splitSamplesLegend.get());
355 
356  // go to the pad for the histograms
357  if (PlotMan->getDrawGrid())
358  LLHPad->SetGrid();
359  LLHPad->Draw();
360  LLHPad->cd();
361 
362  splitSamplesStack->Draw(PlotMan->getDrawOptions().c_str());
363  if (sameAxis)
364  splitSamplesStack->SetMaximum(baseSplitSamplesStack->GetMaximum());
365 
366  if (totalOnSplitPlots)
367  {
368  compLLH_main.Draw(Form("same%s", PlotMan->getDrawOptions().c_str()));
369  compLLH_main.SetLineColor(kBlack);
370  splitSamplesLegend->AddEntry(&compLLH_main, "All Samples", "l");
371  }
372  splitSamplesLegend->Draw();
373 
374  // KS: Not sure why but need to plot after it's drawn otherwise ROOT throws segfault...
375  baseSplitSamplesStack->GetYaxis()->SetTitle("-2LLH_{sam}");
376  if (PlotMan->getPlotRatios() == false) baseSplitSamplesStack->GetXaxis()->SetTitle("Parameter Variation");
377  gPad->Modified();
378  gPad->Update();
379 
380  // need to draw the labels after other stuff or they dont show up
381  for (uint i = 0; i < PlotMan->input().getTaggedSamples(PlotMan->getOption<std::vector<std::string>>("sampleTags")).size(); i++)
382  {
383  std::string sampName = PlotMan->input().getTaggedSamples(PlotMan->getOption<std::vector<std::string>>("sampleTags"))[i];
384  if (!drawLabel[i])
385  continue;
386  label->DrawLatex(compLLH_main.GetBinLowEdge(compLLH_main.GetNbinsX() + 1), extraCumSums[i],
387  Form("#leftarrow%s", PlotMan->style().prettifySampleName(sampName).c_str()));
388  }
389 
390  if (PlotMan->getPlotRatios())
391  {
392  THStack *splitSamplesStackRatios = new THStack(paramName.c_str(), "");
393 
394  TList *baselineHistList = baseSplitSamplesStack->GetHists();
395  TList *compHistList = splitSamplesStack->GetHists();
396  for (uint sampleIdx = 0; sampleIdx < PlotMan->input().getTaggedSamples(PlotMan->getOption<std::vector<std::string>>("sampleTags")).size(); sampleIdx++)
397  {
398  TH1D *divHist = new TH1D(
399  Form("%s_%s_splitDiv_%i", paramName.c_str(), PlotMan->getFileLabel(extraFileIdx).c_str(),
400  sampleIdx),
401  Form("%s_%s_splitDiv", paramName.c_str(), PlotMan->getFileLabel(extraFileIdx).c_str()),
402  compLLH_main.GetNbinsX(), compLLH_main.GetBinLowEdge(1),
403  compLLH_main.GetBinLowEdge(compLLH_main.GetNbinsX() + 1));
404  divHist->Divide(static_cast<TH1D*>(compHistList->At(sampleIdx)),
405  static_cast<TH1D*>(baselineHistList->At(sampleIdx)));
406  splitSamplesStackRatios->Add(divHist);
407  divHist->SetLineColor((static_cast<TH1D*>(compHistList->At(sampleIdx))->GetLineColor()));
408  }
409 
410  canv->cd(2 + extraFileIdx);
411  if (PlotMan->getDrawGrid())
412  ratioPad->SetGrid();
413  ratioPad->Draw();
414  ratioPad->cd();
415  drawRatioStack(splitSamplesStackRatios);
416  }
417  }
418  canv->SaveAs(outputFileName.c_str());
419  delete LLHPad;
420 }
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:50

◆ PlotLLH()

int PlotLLH ( )

Definition at line 422 of file PlotLLH.cpp.

422  {
423  PlotMan->style().setPalette(PlotMan->getOption<std::string>("colorPalette"));
424 
425  // make the canvas and other plotting stuff
426  auto canv = std::make_unique<TCanvas>("canv", "", 1024, 1024);
427  gStyle->SetOptTitle(2);
428  // split up the canvas so can have side by side plots, one for each file
429  auto splitSamplesCanv = std::make_unique<TCanvas>("splitSampCanv", "", 4096 * PlotMan->getNFiles(), 4096);
430 
431  canv->SaveAs((PlotMan->getOutputName("_Sample") + "[").c_str());
432  canv->SaveAs((PlotMan->getOutputName("_Penalty") + "[").c_str());
433  canv->SaveAs((PlotMan->getOutputName("_Total") + "[").c_str());
434 
435  if (PlotMan->getSplitBySample())
436  canv->SaveAs((PlotMan->getOutputName("_bySample") + "[").c_str());
437 
438  for( std::string par: PlotMan->getOption<std::vector<std::string>>("parameterTags")) std::cout << par << ", ";
439  // loop over the spline parameters
440  for (std::string paramName : PlotMan->input().getTaggedParameters(PlotMan->getOption<std::vector<std::string>>("parameterTags")))
441  {
442  MACH3LOG_DEBUG("working on parameter {}", paramName);
443  // ###############################################################
444  // First lets do just the straight up likelihoods from all samples
445  // ###############################################################
446  makeLLHScanComparisons(paramName, "sample", PlotMan->getOutputName("_Sample"), canv);
447  makeLLHScanComparisons(paramName, "penalty", PlotMan->getOutputName("_Penalty"), canv);
448  makeLLHScanComparisons(paramName, "total", PlotMan->getOutputName("_Total"), canv);
449  // #########################################
450  // ## now lets make plots split by sample ##
451  // #########################################
452  if (PlotMan->getSplitBySample())
453  {
454  makeSplitSampleLLHScanComparisons(paramName, PlotMan->getOutputName("_bySample"), splitSamplesCanv);
455  }
456  }
457 
458  canv->SaveAs((PlotMan->getOutputName("_Sample") + "]").c_str());
459  canv->SaveAs((PlotMan->getOutputName("_Penalty") + "]").c_str());
460  canv->SaveAs((PlotMan->getOutputName("_Total") + "]").c_str());
461  if (PlotMan->getSplitBySample())
462  canv->SaveAs((PlotMan->getOutputName("_bySample") + "]").c_str());
463 
464  if(PlotMan != nullptr) delete PlotMan;
465 
466  return 0;
467 }
void makeSplitSampleLLHScanComparisons(const std::string &paramName, const std::string &outputFileName, const std::unique_ptr< TCanvas > &canv)
Definition: PlotLLH.cpp:242
void makeLLHScanComparisons(const std::string &paramName, const std::string &LLHType, const std::string &outputFileName, const std::unique_ptr< TCanvas > &canv)
Definition: PlotLLH.cpp:145

◆ SetTPads()

void SetTPads ( TPad *&  LLHPad,
TPad *&  ratioPad 
)

TPad is SUPER FRAGILE, it is safer to just make raw pointer, while ROOT behave weirdly with smart pointers.

Definition at line 30 of file PlotLLH.cpp.

31 {
32  int originalErrorLevel = gErrorIgnoreLevel;
33  gErrorIgnoreLevel = kFatal;
34  if (PlotMan->getPlotRatios())
35  {
36  LLHPad = new TPad("LLHPad", "LLHPad", 0.0, ratioPlotSplit, 1.0, 1.0);
37  LLHPad->SetBottomMargin(0.0);
38  ratioPad = new TPad("ratioPad", "ratioPad", 0.0, 0.0, 1.0, ratioPlotSplit);
39  ratioPad->SetTopMargin(0.0);
40  ratioPad->SetBottomMargin(0.3);
41  } else {
42  LLHPad = new TPad("AllSampPad", "AllSampPad", 0.0, 0.0, 1.0, 1.0);
43  ratioPad = new TPad("AllSampRatioPad", "AllSampRatioPad", 0.0, 0.0, 0.0, 0.0);
44  }
45  LLHPad->AppendPad();
46  ratioPad->AppendPad();
47  gErrorIgnoreLevel = originalErrorLevel;
48 }

Variable Documentation

◆ lineWidth

int lineWidth

Definition at line 20 of file PlotLLH.cpp.

◆ PlotMan

MaCh3Plotting::PlottingManager* PlotMan
Warning
KS: keep raw pointer or ensure manual delete of PlotMan. If spdlog in automatically deleted before PlotMan then destructor has some spdlog and this could cause segfault

Definition at line 27 of file PlotLLH.cpp.

◆ ratioLabelScaling

double ratioLabelScaling

Definition at line 24 of file PlotLLH.cpp.

◆ ratioPlotSplit

double ratioPlotSplit

Definition at line 17 of file PlotLLH.cpp.

◆ sameAxis

bool sameAxis

Definition at line 22 of file PlotLLH.cpp.

◆ sampleLabelThreshold

double sampleLabelThreshold

Definition at line 19 of file PlotLLH.cpp.

◆ totalOnSplitPlots

bool totalOnSplitPlots

Definition at line 21 of file PlotLLH.cpp.

◆ yTitleOffset

double yTitleOffset

Definition at line 18 of file PlotLLH.cpp.