MaCh3 2.2.1
Reference Guide
Loading...
Searching...
No Matches
GetPostfitParamPlots.cpp
Go to the documentation of this file.
1#include <iostream>
2#include <sstream>
3#include <iomanip>
4#include <algorithm>
5
8
10#include "TROOT.h"
11#include "TGaxis.h"
12#include "TString.h"
13#include "TStyle.h"
14#include "TH1.h"
15#include "TH2.h"
16#include "TF1.h"
17#include "TLegend.h"
18#include "TPad.h"
19#include "TCanvas.h"
20#include "TTree.h"
21#include "TFile.h"
22#include "TVectorD.h"
23#include "TCandle.h"
24#include "TFrame.h"
25#include "TGraphAsymmErrors.h"
27
44
45//this file has lots of usage of the ROOT plotting interface that only takes floats, turn this warning off for this CU for now
46#pragma GCC diagnostic ignored "-Wfloat-conversion"
47#pragma GCC diagnostic ignored "-Wconversion"
48
50TH1D *Prefit;
51
52TH2D *Violin;
53TH2D *Violin2;
55
58
61
64
66
67std::vector<int> NDSamplesBins;
68std::vector<std::string> NDSamplesNames;
70TCanvas *canv;
71
72std::string SaveName;
73
74TPad *p3;
75TPad *p4;
76
77// KS: Color for 0 - prefit, 1 postfit, 2 another postfit, 3 you know the drill
78Color_t PlotColor[] = {kRed, kBlack, kBlue, kGreen};
79std::string plotType;
80
81std::vector <TH1D *> PostfitHistVec;
82
83void copyParToBlockHist(int localBin, std::string paramName, TH1D*blockHist, std::string type, int fileId, bool setLabels = true){
84 // Set the values in the sub-histograms
85 MACH3LOG_DEBUG("copying data from at local bin {}: for parameter {}", localBin, paramName);
86 MACH3LOG_DEBUG(" Fitter specific name: {}", man->input().translateName(fileId, MaCh3Plotting::kPostFit, paramName));
87 MACH3LOG_DEBUG(" value: {}", man->input().getPostFitValue(fileId, paramName, type));
88 MACH3LOG_DEBUG(" error: {}", man->input().getPostFitError(fileId, paramName, type));
89
90 blockHist->SetBinContent(localBin +1, man->input().getPostFitValue(fileId, paramName, type));
91 blockHist->SetBinError(localBin +1, man->input().getPostFitError(fileId, paramName, type));
92
93 if(setLabels){
94 blockHist->GetXaxis()->SetBinLabel(localBin +1, paramName.c_str());
95 blockHist->GetXaxis()->LabelsOption("v");
96 }
97}
98
99inline void InitializePads(TCanvas* canvas, TPad*& pad3, TPad*& pad4) {
100 // Initialize TPad p3
101 pad3 = new TPad("Top", "Top", 0.0, 0.4, 1.0, 1.0);
102 pad3->SetLeftMargin(canvas->GetLeftMargin());
103 pad3->SetRightMargin(canvas->GetRightMargin());
104 pad3->SetTopMargin(canvas->GetTopMargin());
105 pad3->SetBottomMargin(0);
106 pad3->SetGrid();
107
108 // Initialize TPad p4
109 pad4 = new TPad("Bottom", "Bottom", 0.0, 0.0, 1.0, 0.4);
110 pad4->SetLeftMargin(canvas->GetLeftMargin());
111 pad4->SetRightMargin(canvas->GetRightMargin());
112 pad4->SetTopMargin(0);
113 pad4->SetBottomMargin(0.75);
114 pad4->SetGrid();
115}
116
117void CopyViolinToBlock(TH2D* FullViolin, TH2D* ReducedViolin, const std::vector<std::string>& ParamNames) {
118 for(unsigned int i = 0; i < ParamNames.size(); i++)
119 {
120 int ParamBinId = -999;
121 for (int ix = 0; ix < FullViolin->GetXaxis()->GetNbins(); ++ix) {
122 if(FullViolin->GetXaxis()->GetBinLabel(ix+1) == ParamNames[i])
123 {
124 ParamBinId = ix+1;
125 break;
126 }
127 }
128 if(ParamBinId == -999) {
129 MACH3LOG_WARN("Didn't find param {}", ParamNames[i]);
130 continue;
131 }
132 //KS Fill content of reduced violin
133 for (int iy = 0; iy < FullViolin->GetYaxis()->GetNbins(); ++iy) {
134 ReducedViolin->SetBinContent(i+1, iy+1, FullViolin->GetBinContent(ParamBinId, iy+1));
135 ReducedViolin->GetXaxis()->SetBinLabel(i+1, ParamNames[i].c_str());
136 }
137 }
138 ReducedViolin->SetFillColor(FullViolin->GetFillColor());
139 ReducedViolin->SetFillColorAlpha(FullViolin->GetMarkerColor(), 0.35);
140 ReducedViolin->SetLineColor(FullViolin->GetMarkerColor());
141
142 ReducedViolin->SetMarkerColor(FullViolin->GetMarkerColor());
143 ReducedViolin->SetMarkerStyle(FullViolin->GetMarkerStyle());
144 ReducedViolin->SetMarkerSize(FullViolin->GetMarkerSize());
145
146 ReducedViolin->GetYaxis()->SetTitleOffset(FullViolin->GetTitleOffset());
147 ReducedViolin->GetYaxis()->SetTitle(FullViolin->GetYaxis()->GetTitle());
148 ReducedViolin->GetXaxis()->LabelsOption("v");
149}
150
151void PrettifyTitles(TH1D *Hist) {
152 for (int i = 0; i < Hist->GetXaxis()->GetNbins(); ++i)
153 {
154 std::string title = Hist->GetXaxis()->GetBinLabel(i+1);
155 title = man->style().prettifyParamName(title);
156 Hist->GetXaxis()->SetBinLabel(i+1, title.c_str());
157 }
158}
159
160void PrettifyTitles(TH2D *Hist) {
161 for (int i = 0; i < Hist->GetXaxis()->GetNbins(); ++i)
162 {
163 std::string title = Hist->GetXaxis()->GetBinLabel(i+1);
164
165 title = man->style().prettifyParamName(title);
166 Hist->GetXaxis()->SetBinLabel(i+1, title.c_str());
167 }
168}
169
170bool ReadSettings(std::shared_ptr<TFile> File1)
171{
172 MACH3LOG_DEBUG("Reading settings for file {}", File1->GetName());
173 #ifdef DEBUG
174 File1->ls();
175 #endif
176 TTree *Settings = (File1->Get<TTree>("Settings"));
177
178 // can't find settings tree :(
179 if (!Settings) return false;
180
181 MACH3LOG_DEBUG("Got settings tree");
182 Settings->Print();
183
184 Settings->SetBranchAddress("CrossSectionParameters", &CrossSectionParameters);
185 MACH3LOG_DEBUG("XSec params: {}", CrossSectionParameters);
186 Settings->SetBranchAddress("FluxParameters", &FluxParameters);
187 Settings->SetBranchAddress("CrossSectionParametersStartingPos", &XsecStartingPos);
188 Settings->SetBranchAddress("NDParameters", &NDParameters);
189 Settings->SetBranchAddress("NDParametersStartingPos", &NDParametersStartingPos);
190 Settings->SetBranchAddress("FDParameters", &FDParameters);
191 Settings->SetBranchAddress("FDParametersStartingPos", &FDParametersStartingPos);
192
193 std::vector<int> *NDSamples_Bins = 0;
194 std::vector<std::string> *NDSamples_Names = 0;
195 Settings->SetBranchAddress("NDSamplesNames", &NDSamples_Names);
196 Settings->SetBranchAddress("NDSamplesBins", &NDSamples_Bins);
197
198 Settings->GetEntry(0);
199
200 NDSamplesNames = *NDSamples_Names;
201 NDSamplesBins = *NDSamples_Bins;
202
203 MACH3LOG_DEBUG("Read settings tree successfully");
204 return true;
205}
206
207inline TH1D* makeRatio(TH1D *PrefitCopy, TH1D *PostfitCopy, bool setAxes){
208 // set up the ratio hist
209 TH1D* Ratio = static_cast<TH1D*>(PrefitCopy->Clone());
210 Ratio->GetYaxis()->SetTitle("(x_{Post}-#mu_{Prior})/#sigma_{Prior}");
211 Ratio->SetMinimum(-3.7);
212 Ratio->SetMaximum(3.7);
213
214 for (int j = 0; j < Ratio->GetXaxis()->GetNbins(); ++j)
215 {
216 if ( PrefitCopy->GetBinError(j+1) > 1.e-5 )
217 {
218 Ratio->SetBinContent(j+1, (PostfitCopy->GetBinContent(j+1)-PrefitCopy->GetBinContent(j+1))/PrefitCopy->GetBinError(j+1));
219
220 double up = (PostfitCopy->GetBinContent(j+1)+PostfitCopy->GetBinError(j+1)-PrefitCopy->GetBinContent(j+1))/PrefitCopy->GetBinError(j+1);
221 double down = (PostfitCopy->GetBinContent(j+1)-PostfitCopy->GetBinError(j+1)-PrefitCopy->GetBinContent(j+1))/PrefitCopy->GetBinError(j+1);
222
223 double maximum = up-Ratio->GetBinContent(j+1);
224 double minimum = Ratio->GetBinContent(j+1)-down;
225
226 Ratio->SetBinError(j+1, std::max(maximum, minimum));
227 }
228 //KS: Most likely flat prior
229 else {
230 Ratio->SetBinContent(j+1, (PostfitCopy->GetBinContent(j+1)-PrefitCopy->GetBinContent(j+1)));
231
232 double up = (PostfitCopy->GetBinContent(j+1)+PostfitCopy->GetBinError(j+1)-PrefitCopy->GetBinContent(j+1));
233 double down = (PostfitCopy->GetBinContent(j+1)-PostfitCopy->GetBinError(j+1)-PrefitCopy->GetBinContent(j+1));
234
235 double maximum = up-Ratio->GetBinContent(j+1);
236 double minimum = Ratio->GetBinContent(j+1)-down;
237
238 Ratio->SetBinError(j+1, std::max(maximum, minimum));
239 }
240 } //end loop over parameters
241
242 if(setAxes){
243 Ratio->SetFillStyle(0);
244 Ratio->SetFillColor(0);
245
246 Ratio->SetLineColor(PostfitCopy->GetLineColor());
247 if (Ratio->GetLineColor() == 0) Ratio->SetLineColor(kBlack);
248 Ratio->SetMarkerColor(PostfitCopy->GetMarkerColor());
249
250 Ratio->SetLineWidth(man->getOption<int>("plotLineWidth"));
251 Ratio->SetTitle("");
252
253 Ratio->SetMarkerSize(2);
254 Ratio->SetMarkerStyle(20);
255
256 Ratio->GetYaxis()->SetTitleSize(25);
257 Ratio->GetYaxis()->SetTitleFont(43);
258 Ratio->GetYaxis()->SetTitleOffset(2.0);
259 Ratio->GetYaxis()->SetLabelFont(43);
260 Ratio->GetYaxis()->SetLabelSize(25);
261 Ratio->GetYaxis()->CenterTitle();
262 Ratio->GetYaxis()->SetNdivisions(5,2,0);
263
264 Ratio->GetXaxis()->SetTitleSize(25);
265 Ratio->GetXaxis()->SetTitleFont(43);
266 Ratio->GetXaxis()->SetTitleOffset(4.0);
267 Ratio->GetXaxis()->SetLabelOffset(0.025);
268 Ratio->GetXaxis()->SetLabelFont(43);
269 Ratio->GetXaxis()->SetLabelSize(25);
270 }
271
272 return Ratio;
273}
274
275inline void DrawPlots(TCanvas *plotCanv, TH1D* PrefitCopy, const std::vector<TH1D*>& PostfitVec, TPad *mainPad, TPad *ratioPad) {
276 // Draw!
277 plotCanv->cd();
278 mainPad->Draw();
279 mainPad->cd();
280 PrefitCopy->GetYaxis()->SetTitle("Parameter Value");
281
282 PrefitCopy->GetYaxis()->SetLabelSize(0.);
283 PrefitCopy->GetYaxis()->SetTitleSize(0.05);
284 PrefitCopy->GetYaxis()->SetTitleOffset(1.3);
285 PrefitCopy->Draw("e2");
286
287 for (int fileId = 0; fileId < static_cast<int>(PostfitVec.size()); fileId++) {
288 TH1D *postFitHist = PostfitVec[fileId];
289
290 postFitHist->SetMarkerColor(TColor::GetColorPalette(fileId));
291 postFitHist->SetLineColor(TColor::GetColorPalette(fileId));
292 postFitHist->SetMarkerStyle(7);
293 postFitHist->SetLineStyle(1+fileId);
294 postFitHist->SetLineWidth(man->getOption<int>("plotLineWidth"));
295
296 postFitHist->Draw("e1, same");
297 }
298
299 plotCanv->Update();
300 TGaxis *axis = new TGaxis(PrefitCopy->GetXaxis()->GetBinLowEdge(PrefitCopy->GetXaxis()->GetFirst()), gPad->GetUymin()+0.01,
301 PrefitCopy->GetXaxis()->GetBinLowEdge(PrefitCopy->GetXaxis()->GetFirst()), gPad->GetUymax(),
302 gPad->GetUymin()+0.01, gPad->GetUymax(), 510, "");
303 axis->SetLabelFont(43);
304 axis->SetLabelSize(25);
305 axis->Draw();
306
307 plotCanv->cd();
308 ratioPad->Draw();
309 ratioPad->cd();
310
311 std::vector<TH1D*> ratioHists;
312
313 // save pointers to these so we can delete them once we are done
314 ratioHists.push_back(makeRatio(PrefitCopy, PostfitVec[0], true));
315
316 ratioHists[0]->Draw("p");
317 for(int postFitIdx = 1; postFitIdx < static_cast<int>(PostfitVec.size()); postFitIdx++){
318 ratioHists.push_back(makeRatio(PrefitCopy, PostfitVec[postFitIdx], true));
319
320 ratioHists[postFitIdx]->SetMarkerColor(TColor::GetColorPalette(postFitIdx));
321 ratioHists[postFitIdx]->SetLineColor(TColor::GetColorPalette(postFitIdx));
322 ratioHists[postFitIdx]->SetMarkerStyle(7);
323 ratioHists[postFitIdx]->SetLineStyle(1+postFitIdx);
324 ratioHists[postFitIdx]->SetLineWidth(man->getOption<int>("plotLineWidth"));
325
326 ratioHists[postFitIdx]->Draw("p same");
327 }
328
329 // draw lines across the plot at +-1 and 0
330 TLine line(ratioHists[0]->GetXaxis()->GetBinLowEdge(ratioHists[0]->GetXaxis()->GetFirst()), 0.0, ratioHists[0]->GetXaxis()->GetBinLowEdge(ratioHists[0]->GetXaxis()->GetLast()+1), 0.0);
331 TLine line2(ratioHists[0]->GetXaxis()->GetBinLowEdge(ratioHists[0]->GetXaxis()->GetFirst()), 1.0, ratioHists[0]->GetXaxis()->GetBinLowEdge(ratioHists[0]->GetXaxis()->GetLast()+1), 1.0);
332 TLine line3(ratioHists[0]->GetXaxis()->GetBinLowEdge(ratioHists[0]->GetXaxis()->GetFirst()), -1.0, ratioHists[0]->GetXaxis()->GetBinLowEdge(ratioHists[0]->GetXaxis()->GetLast()+1), -1.0);
333
334 line.SetLineColor(kRed);
335 line.SetLineStyle(kDashed);
336 line.SetLineWidth(man->getOption<int>("refLineWidth"));
337 line2.SetLineColor(kRed);
338 line2.SetLineStyle(kDashed);
339 line2.SetLineWidth(man->getOption<int>("refLineWidth"));
340 line3.SetLineColor(kRed);
341 line3.SetLineStyle(kDashed);
342 line3.SetLineWidth(man->getOption<int>("refLineWidth"));
343
344 line.Draw("same");
345 line2.Draw("same");
346 line3.Draw("same");
347
348 plotCanv->Print((SaveName).c_str());
349
350 ratioHists.clear();
351 delete axis;
352}
353
355{
356 // get the names of the blocks of parameters to group together
357 std::vector<std::string> const blockNames = man->getOption<std::vector<std::string>>("paramGroups");
358 const int XsecPlots = static_cast<int>(blockNames.size());
359
360 for (int i = 0; i < XsecPlots; i++)
361 {
362 // get the configuration for this parameter
363 std::string blockName = blockNames[i];
364 YAML::Node paramBlock = man->getOption(blockName);
365 std::string blockTitle = paramBlock[0].as<std::string>();
366 std::vector<double> blockLimits = paramBlock[1].as<std::vector<double>>();
367 std::vector<std::string> blockContents = paramBlock[2].as<std::vector<std::string>>();
368
369 // get num of params in the block
370 const int nParams = static_cast<int>(blockContents.size());
371
372 // set some plot things
373 TH1D *blockHist_prefit = new TH1D(blockName.c_str(),
374 blockTitle.c_str(), nParams, 0.0, static_cast<double>(nParams));
375
376 man->style().setTH1Style(blockHist_prefit, man->getOption<std::string>("prefitHistStyle"));
377
378 // set the errors for the prefit block hist
379 for(int localBin=0; localBin < nParams; localBin ++){
380 // the "local" bin is the params index within the group of parameters
381 std::string paramName = blockContents[localBin];
382 copyParToBlockHist(localBin, paramName, blockHist_prefit, "Prior", 0);
383 }
384
385 // now set for the postfit blocks for all files
386 std::vector <TH1D *> blockHist_postfit_Vec;
387 for(unsigned int fileId = 0; fileId < man->getNFiles(); fileId++){
388
389 TH1D *blockHist_postfit = new TH1D((blockName + man->getFileName(fileId)).c_str(),
390 blockTitle.c_str(), nParams, 0.0, static_cast<double>(nParams));
391
392 // loop through all the parameters in this block and set the contents in the blocks TH1
393 for(int localBin=0; localBin < nParams; localBin ++){
394 // the "local" bin is the params index within the group of parameters
395 std::string paramName = blockContents[localBin];
396 copyParToBlockHist(localBin, paramName, blockHist_postfit, "", fileId);
397 }
398
399 blockHist_postfit_Vec.push_back(blockHist_postfit);
400 }
401
402 // set the y axis limits we got from config
403 blockHist_prefit->GetYaxis()->SetRangeUser(blockLimits[0], blockLimits[1]);
404
405 // Do some fancy replacements
406 PrettifyTitles(blockHist_prefit);
407
408 DrawPlots(canv, blockHist_prefit, blockHist_postfit_Vec, p3, p4);
409 delete blockHist_prefit;
410 for (auto hist : blockHist_postfit_Vec) {
411 delete hist;
412 }
413 blockHist_postfit_Vec.clear();
414 }
415}
416
418{
419 // these for non named params where we dont need as much space
420 TPad* p1 = new TPad("p1", "p1", 0.0, 0.3, 1.0, 1.0);
421 TPad* p2 = new TPad("p2", "p2", 0.0, 0.0, 1.0, 0.3);
422 p1->SetLeftMargin(canv->GetLeftMargin());
423 p1->SetRightMargin(canv->GetRightMargin());
424 p1->SetTopMargin(canv->GetTopMargin());
425 p1->SetBottomMargin(0);
426 p1->SetGrid();
427 p2->SetLeftMargin(canv->GetLeftMargin());
428 p2->SetRightMargin(canv->GetRightMargin());
429 p2->SetTopMargin(0);
430 p2->SetBottomMargin(0.25);
431 p2->SetGrid();
432
433 p1->SetLogx(true);
434 p2->SetLogx(true);
435
436 // get the names of the blocks of parameters to group together
437 std::vector<std::string> const fluxBlockNames = man->getOption<std::vector<std::string>>("fluxGroups");
438 auto const fluxBinningTable = man->getOption("FluxBinning");
439
440 const int FluxPlots = static_cast<int>(fluxBlockNames.size());
441
442 for (int i = 0; i < FluxPlots; i++)
443 {
444 // get the configuration for this block
445 std::string fluxBlockName = fluxBlockNames[i];
446 YAML::Node paramBlock = man->getOption(fluxBlockName);
447 std::string blockTitle = paramBlock[0].as<std::string>();
448 std::vector<double> blockLimits = paramBlock[1].as<std::vector<double>>();
449 std::string blockBinningName = paramBlock[2].as<std::string>();
450 std::vector<int> blockContents = paramBlock[3].as<std::vector<int>>();
451
452 // get the binning for this block of flux params
453 std::vector<double> binning = fluxBinningTable[blockBinningName].as<std::vector<double>>();
454
455 // get num of params in the block
456 int nParams = blockContents[1] - blockContents[0] +1;
457 // check for sanity
458 if(nParams <= 0 || blockContents.size() > 2){
459 MACH3LOG_CRITICAL("Invalid flux parameter block endpoints specified for {}", fluxBlockName);
460 MACH3LOG_CRITICAL(" Should have the form [<low index>, <up index>]");
461 throw MaCh3Exception(__FILE__ , __LINE__ );
462 }
463 if (nParams != static_cast<int>(binning.size()) - 1) {
464 MACH3LOG_CRITICAL("Binning provided for flux param block {} does not match the number of parameters specified for the block", fluxBlockName);
465 MACH3LOG_CRITICAL(" Provided {} parameters but {} bins", nParams, binning.size() -1);
466 throw MaCh3Exception(__FILE__ , __LINE__ );
467 }
468
469 TH1D *blockHist_prefit = new TH1D(fluxBlockName.c_str(), blockTitle.c_str(), nParams, binning.data());
470 blockHist_prefit->GetYaxis()->SetTitle("Parameter Variation");
471 blockHist_prefit->GetXaxis()->SetTitle("E_{#nu} (GeV)");
472 blockHist_prefit->GetXaxis()->SetTitleOffset(blockHist_prefit->GetXaxis()->GetTitleOffset()*1.2);
473 man->style().setTH1Style(blockHist_prefit, man->getOption<std::string>("prefitHistStyle"));
474 // set the errors for the prefit block hist
475 for(int fluxParId = blockContents[0]; fluxParId <= blockContents[1]; fluxParId++){
476 int localBin = fluxParId - blockContents[0];
477 std::string paramName = "b_" + std::to_string(fluxParId);
478 copyParToBlockHist(localBin, paramName, blockHist_prefit, "Prior", 0, false);
479 }
480
481 // now set for the postfit blocks for all files
482 std::vector <TH1D *> blockHist_postfit_Vec;
483 for(unsigned int fileId = 0; fileId < man->getNFiles(); fileId++){
484 TH1D *blockHist_postfit = new TH1D(fluxBlockName.c_str(), blockTitle.c_str(), nParams, binning.data());
485
486 for(int fluxParId = blockContents[0]; fluxParId <= blockContents[1]; fluxParId++){
487 int localBin = fluxParId - blockContents[0];
488 std::string paramName = "b_" + std::to_string(fluxParId);
489
490 copyParToBlockHist(localBin, paramName, blockHist_postfit, "", fileId, false);
491 }
492
493 blockHist_postfit_Vec.push_back(blockHist_postfit);
494 }
495
496 // set the y axis limits we got from config
497 blockHist_prefit->GetYaxis()->SetRangeUser(blockLimits[0], blockLimits[1]);
498
499 DrawPlots(canv, blockHist_prefit, blockHist_postfit_Vec, p1, p2);
500 delete blockHist_prefit;
501 for (auto hist : blockHist_postfit_Vec) {
502 delete hist;
503 }
504 blockHist_postfit_Vec.clear();
505 }
506
507 canv->cd();
508 canv->SetLogx(false);
509 canv->SetBottomMargin(canv->GetBottomMargin()*1.7);
510 delete p1;
511 delete p2;
512}
513
515{
516 MACH3LOG_INFO("ND detector parameters: {}", NDParameters);
517 Prefit->GetYaxis()->SetTitleOffset(Prefit->GetYaxis()->GetTitleOffset()*1.2);
518
519 TPad* pTop = nullptr;
520 TPad* pDown = nullptr;
521 InitializePads(canv, pTop, pDown);
522
523 int NDbinCounter = NDParametersStartingPos;
524 int Start = NDbinCounter;
525
526 MACH3LOG_INFO("Running on {} samples", NDSamplesNames.size());
527
528 for (unsigned int i = 0; i < NDSamplesNames.size(); ++i)
529 {
530 MACH3LOG_DEBUG("--- On sample {}", NDSamplesNames[i]);
531 NDbinCounter += NDSamplesBins[i];
532
533 std::vector<TH1D*> PostfitNDDetHistVec(man->getNFiles());
534 TH1D *PreFitNDDetHist = man->input().getFile(0).file->Get<TH1D>(Form("param_%s_prefit", NDSamplesNames[i].c_str()));
535 man->style().setTH1Style(PreFitNDDetHist, man->getOption<std::string>("prefitHistStyle"));
536
537 std::string temp = NDSamplesNames[i].c_str();
538 while (temp.find("_") != std::string::npos) {
539 temp.replace(temp.find("_"), 1, std::string(" "));
540 }
541 PreFitNDDetHist->SetTitle(temp.c_str());
542 PreFitNDDetHist->GetXaxis()->SetRangeUser(Start, NDbinCounter);
543
544 MACH3LOG_DEBUG(" Start bin: {} :: End bin: {}", Start, NDbinCounter);
545 // set the x range for the postfits
546 for(unsigned int fileId = 0; fileId < man->getNFiles(); fileId++){
547 PostfitNDDetHistVec[fileId] = man->input().getFile(fileId).file->Get<TH1D>(Form("param_%s_%s", NDSamplesNames[i].c_str(), plotType.c_str()));
548 }
549
550 //KS: We dont' need name for every nd param
551 for(int j = 0; j < NDSamplesBins[i]; ++j)
552 {
553 bool ProductOfTen = false;
554 if(j % 10) ProductOfTen = true;
555 if(j != 0 && ProductOfTen) PreFitNDDetHist->GetXaxis()->SetBinLabel(Start+j+1, " ");
556 else{
557 PreFitNDDetHist->GetXaxis()->SetBinLabel(Start+j+1, Form("Det Variation Bin %i", Start+j));
558 }
559 }
560
561 PreFitNDDetHist->GetYaxis()->SetRangeUser(man->getOption<double>("detParYRange_low"), man->getOption<double>("detParYRange_high"));
562
563 Start += NDSamplesBins[i];
564
565 DrawPlots(canv, PreFitNDDetHist, PostfitNDDetHistVec, pTop, pDown);
566 canv->Update();
567 }
568 delete pTop;
569 delete pDown;
570}
571
573{
574 Prefit->GetYaxis()->SetTitleOffset(Prefit->GetYaxis()->GetTitleOffset()*1.2);
575
576 // these for named parameters where we need a nice big gap at the bottom to fit the names
577 TPad* pTop = nullptr;
578 TPad* pDown = nullptr;
579 InitializePads(canv, pTop, pDown);
580
581 int FDbinCounter = FDParametersStartingPos;
582 int Start = FDbinCounter;
583 //KS: WARNING This is hardcoded
584 double FDSamplesBins[7] = {12,18,30,36,44,57,58};
585 std::string FDSamplesNames[7] = {"FHC 1Re", "FHC 1R#mu", "RHC 1Re","RHC 1R#mu","FHC 1Re 1 d.e.","FHC MR#mu 1 or 2 d.e.", "Momentum Scale"};
586 for (unsigned int i = 0; i < 7; ++i)
587 {
588 FDbinCounter = FDParametersStartingPos + FDSamplesBins[i];
589
590 canv->cd();
591 Prefit->SetTitle(FDSamplesNames[i].c_str());
592 Prefit->GetXaxis()->SetRangeUser(Start, FDbinCounter);
593 //KS: We don't' need name for every FD param
594 for(int j = 0; j < FDSamplesBins[i]; ++j)
595 {
596 //bool ProductOfTen = false;
597 //if(j % 10) ProductOfTen = true;
598 //if(j != 0 && ProductOfTen) Prefit->GetXaxis()->SetBinLabel(Start+j+1, " ");
599 }
600 Prefit->GetYaxis()->SetRangeUser(man->getOption<double>("detParYRange_low"), man->getOption<double>("detParYRange_high"));
601 Prefit->Draw("e2");
602
603 Start = FDParametersStartingPos + FDSamplesBins[i];
604
605 DrawPlots(canv, Prefit, PostfitHistVec, pTop, pDown);
606 canv->Update();
607 }
608 delete pTop;
609 delete pDown;
610}
611
613{
614 gStyle->SetPalette(51);
615
616 TCanvas *blankCanv = new TCanvas ("blankCanv", "blankCanv", 2048, 2048);
617 blankCanv->SaveAs("RidgePlots.pdf[");
618
619 // get the names of the blocks of parameters to group together
620 std::vector<std::string> const blockNames = man->getOption<std::vector<std::string>>("paramGroups");
621 const int XsecPlots = static_cast<int>(blockNames.size());
622
623 double padTopMargin = 0.9;
624 double padBottomMargin = 0.1;
625 double padOverlap = 0.9;
626 double ridgeLineWidth = 1.0;
627
628 for (int i = 0; i < XsecPlots; i++)
629 {
630 // get the configuration for this parameter
631 std::string blockName = blockNames[i];
632 auto const &paramBlock = man->getOption(blockName);
633 auto blockTitle = paramBlock[0].as<std::string>();
634 auto blockLimits = paramBlock[1].as<std::vector<double>>();
635 auto blockContents = paramBlock[2].as<std::vector<std::string>>();
636
637 // the directory of histograms
638 TDirectoryFile *posteriorDir = man->input().getFile(0).file->Get<TDirectoryFile>("Post_1d_hists");
639
640 // get num of params in the block
641 int nParams = static_cast<int>(blockContents.size());
642 TCanvas *ridgeCanv = new TCanvas ("RidgePlotCanv", "RidgePlotCanv", 2048, 2048);
643 ridgeCanv->Divide(1,1+nParams, 0.01, 0.0);
644
645 auto title = std::make_unique<TLatex>();
646 title->SetTextAlign(21);
647 title->SetTextSize(0.03);
648 title->DrawLatex(0.5, 0.95, blockTitle.c_str());
649
650 auto label = std::make_unique<TLatex>();
651 label->SetTextAlign(31);
652 label->SetTextSize(0.02);
653
654 auto line = std::make_unique<TLine>();
655 line->SetLineColor(kBlack);
656 line->SetLineWidth(ridgeLineWidth);
657
658 // use this to set the limits and also to plot the x axis and grid
659 TH1D *axisPlot = new TH1D("axis plot", "", 1, blockLimits[0], blockLimits[1]);
660
661 for(int parId = 0; parId < nParams; parId++) {
662 std::string paramName = blockContents[parId];
663
664 TH1D *posteriorDist = nullptr;
665 // get the list of objects in the directory
666 TIter next(posteriorDir->GetListOfKeys());
667 while (TKey* key = static_cast<TKey*>(next())) {
668 // check if the end of the param name matches with the MaCh3 name, do this so we exclude things like nds_ at the start of the name
669 std::string str(key->GetTitle());
670 std::string name = man->input().translateName(0, MaCh3Plotting::kPostFit, paramName);
671 uint pos = str.find(name);
672 bool foundPar = (pos == str.length() - name.length());
673
674 if(foundPar){
675 posteriorDist = posteriorDir->Get<TH1D>(key->GetName());
676 }
677 }
678
679 if(posteriorDist == nullptr){
680 MACH3LOG_WARN("Couldn't find parameter {} when making ridgeline plots", paramName);
681 continue;
682 }
683
684 // EM: do some funky scaling so that we always get evenly spaced pads in the range [bottomMargin, TopMargin] with the specified overlap
685 double padAnchor = padBottomMargin + (static_cast<double>(nParams - parId - 1) /
686 static_cast<double>(nParams - 1)) * (padTopMargin - padBottomMargin);
687 double padWidth = (padTopMargin - padBottomMargin) / static_cast<double>(nParams);
688 double norm = (padTopMargin - padBottomMargin);
689
690 double padTop = padWidth * (1.0 + padOverlap) * (padTopMargin - padAnchor) / norm + padAnchor;
691 double padBottom = padAnchor - padWidth * (1.0 + padOverlap) * (padAnchor - padBottomMargin) / norm;
692
693 TPad *pad = new TPad(paramName.c_str(), "", 0.3, padBottom, 0.9, padTop, -1, 0, -1);
694 ridgeCanv->cd();
695
696 pad->SetBottomMargin(0.0);
697 pad->SetTopMargin(0.0);
698 pad->SetLeftMargin(0.0);
699 pad->SetRightMargin(0.0);
700
701 pad->Draw();
702 pad->cd();
703 pad->SetFillStyle(4000);
704
705 gPad->SetFrameFillStyle(4000);
706 posteriorDist->GetFunction("Gauss")->SetBit(TF1::kNotDraw);
707 posteriorDist->SetTitle("");
708 posteriorDist->SetLineWidth(ridgeLineWidth);
709
710 TH1D* axisPlot_tmp = static_cast<TH1D*>(axisPlot->Clone(Form("AxisPlot_%s", paramName.c_str())));
711 axisPlot_tmp->Draw("A");
712 posteriorDist->Draw("H SAME");
713
714 axisPlot_tmp->GetYaxis()->SetRangeUser(0.0, 0.7 *posteriorDist->GetMaximum());
715 posteriorDist->SetLineColor(kWhite);
716 posteriorDist->SetFillColorAlpha(TColor::GetColorPalette(floor(static_cast<float>(parId) *
717 TColor::GetNumberOfColors() / static_cast<float>(nParams))), 0.85);
718
719 posteriorDist->GetXaxis()->SetRangeUser(blockLimits[0], blockLimits[1]);
720 posteriorDist->GetYaxis()->SetTitle(paramName.c_str());
721
722 //EM: Have to get the frame drawn by the histogram and then set it to be transparent... it took me an hour to get rid of this one line on a plot
723 gPad->Modified(); gPad->Update();
724 TFrame *frame = gPad->GetFrame();
725 frame->SetLineColorAlpha(0, 0.0);
726
727 ridgeCanv->cd();
728 label->DrawLatexNDC(0.29, padBottom + 0.005, man->style().prettifyParamName(paramName).c_str());
729 line->DrawLine(0.1, padBottom, 0.9, padBottom);
730 }
731
732 ridgeCanv->cd();
733 ridgeCanv->SetGrid(1,1);
734 TPad *axisPad = new TPad("AxisPad", "", 0.3, 0.0, 0.9, 1.0, -1, 0, -1);
735 axisPad->SetLeftMargin(0.0);
736 axisPad->SetRightMargin(0.0);
737 axisPad->Draw();
738 axisPad->cd();
739 axisPad->SetGrid(1,1);
740 axisPad->SetFrameFillStyle(4000);
741
742 axisPlot->GetXaxis()->SetTickSize(0.01);
743 axisPlot->GetXaxis()->SetTitle("");
744 axisPlot->GetYaxis()->SetLabelOffset(9999);
745 axisPlot->GetYaxis()->SetLabelSize(0);
746 axisPlot->GetYaxis()->SetTickSize(0);
747 axisPlot->GetYaxis()->SetAxisColor(0,0.0);
748 axisPlot->Draw("AXIS");
749 axisPlot->Draw("AXIG SAME");
750
751 axisPlot->SetFillStyle(4000);
752 axisPad->SetFillStyle(4000);
753
754 axisPad->SetGrid(1,1);
755 gPad->Modified(); gPad->Update();
756 gPad->SetFrameFillStyle(4000);
757
758 gPad->Modified(); gPad->Update();
759 TFrame *frame = gPad->GetFrame();
760 frame->SetLineColorAlpha(0, 0.0);
761
762 ridgeCanv->SaveAs("RidgePlots.pdf");
763 delete ridgeCanv;
764 delete axisPlot;
765 }
766 blankCanv->SaveAs("RidgePlots.pdf]");
767 delete blankCanv;
768}
769
771{
773
774 //KS: By default we take HPD values, but by changing "plotType" you can use for example Gauss
775 plotType = "HPD";
776 //plotType = "gaus";
777
778 MACH3LOG_INFO("Plotting {} errors", plotType);
779
780 // if we have one MaCh3 nd file then we can get settings from it
781 bool plotNDDet = false;
782 for (size_t fileId = 0; fileId < man->input().getNInputFiles(); fileId++) {
783 if(!ReadSettings(man->input().getFile(0).file)) {
784 MACH3LOG_INFO("at least one file provided does not have 'settings' tree indicating it is not MaCh3 ND file");
785 MACH3LOG_INFO(" sadly this means I cannot plot ND Det parameters as this is only supported for MaCh3 ND files for now... sorry :(");
786 plotNDDet = false;
787 }
788 }
789
790 canv = new TCanvas("canv", "canv", 1024, 1024);
791 //gStyle->SetPalette(51);
792 gStyle->SetOptStat(0); //Set 0 to disable statystic box
793 canv->SetLeftMargin(0.12);
794 canv->SetBottomMargin(0.12);
795 canv->SetTopMargin(0.08);
796 canv->SetRightMargin(0.04);
797
798 canv->Print((SaveName+"[").c_str());
799
800 // these for named parameters where we need a nice big gap at the botto to fit the names
802
803 // Make a Legend page
804 auto leg = std::make_unique<TLegend>(0.0, 0.0, 1.0, 1.0);
805 // make a dummy TH1 to set out legend
806 Prefit = new TH1D();
807 man->style().setTH1Style(Prefit, man->getOption<std::string>("prefitHistStyle"));
808 leg->AddEntry(Prefit, "Prior", "lpf");
809
810 for(unsigned int fileId = 0; fileId < man->getNFiles(); fileId++){
811 TH1D *postFitHist_tmp = new TH1D();
812 postFitHist_tmp->SetBit(kCanDelete);
813
814 postFitHist_tmp->SetMarkerColor(TColor::GetColorPalette(fileId));
815 postFitHist_tmp->SetLineColor(TColor::GetColorPalette(fileId));
816 postFitHist_tmp->SetMarkerStyle(7);
817 postFitHist_tmp->SetLineStyle(1+fileId);
818 postFitHist_tmp->SetLineWidth(man->getOption<int>("plotLineWidth"));
819 leg->AddEntry(postFitHist_tmp, man->getFileLabel(fileId).c_str(), "lpf");
820 }
821
822 canv->cd();
823 canv->Clear();
824 leg->Draw();
825 canv->Print((SaveName).c_str());
826
828
830
831 //KS: By default we don't run ProcessMCMC with PlotDet as this take some time, in case we did let's make fancy plots
832 if(plotNDDet & (NDParameters > 0)) MakeNDDetPlots();
833
834 //KS: Same as above but for FD parameters,
836
837 canv->Print((SaveName+"]").c_str());
838
840
841 delete canv;
842 delete Prefit;
843}
844
845inline TGraphAsymmErrors* MakeTGraphAsymmErrors(std::shared_ptr<TFile> File, std::vector<int> Index = {})
846{
847 int GraphBins = Index.size() == 0 ? nBins : Index.size();
848 std::vector<double> x(GraphBins);
849 std::vector<double> y(GraphBins);
850 std::vector<double> exl(GraphBins);
851 std::vector<double> eyl(GraphBins);
852 std::vector<double> exh(GraphBins);
853 std::vector<double> eyh(GraphBins);
854
855 TH1D* PostHist = static_cast<TH1D*>(File->Get( ("param_xsec_"+plotType).c_str() ));
856
857 TVectorD* Errors_HPD_Positive = static_cast<TVectorD*>(File->Get( "Errors_HPD_Positive" ));
858 TVectorD* Errors_HPD_Negative = static_cast<TVectorD*>(File->Get( "Errors_HPD_Negative" ));
859 //KS: I am tempted to multithread this...
860 for(int i = 0; i < GraphBins; ++i)
861 {
862 int Counter = Index.size() == 0 ? i : Index[i];
863 //KS: We are extracting value from three object each having different numbering scheme, I have checked carefully so this is correct please don't change all these +1 +0.5 etc. it just work...
864 x[i] = i + 0.5;
865 y[i] = PostHist->GetBinContent(Counter+1);
866
867 //KS: We don't want x axis errors as they are confusing in Violin plot
868 exh[i] = 0.00001;
869 exl[i] = 0.00001;
870 eyh[i] = (*Errors_HPD_Positive)(Counter);
871 eyl[i] = (*Errors_HPD_Negative)(Counter);
872 }
873 TGraphAsymmErrors* PostGraph = new TGraphAsymmErrors(GraphBins, x.data(), y.data(), exl.data(), exh.data(), eyl.data(), eyh.data());
874 PostGraph->SetTitle("");
875
876 return PostGraph;
877}
878
879//KS: Make fancy violin plots
880void GetViolinPlots(std::string FileName1 = "", std::string FileName2 = "")
881{
882 //KS: Should be in some config... either way it control whether you plot symmetric or asymmetric error bars
883 bool PlotAssym = true;
884
885 //KS: No idea why but ROOT changed treatment of violin in R6. If you have non uniform binning this will results in very hard to see violin plots.
886 TCandle::SetScaledViolin(false);
887
888 MACH3LOG_INFO("Making Violin Plot");
889 if (!FileName1.empty()) MACH3LOG_INFO("File 1: {} ", FileName1);
890 if (!FileName2.empty()) MACH3LOG_INFO("File 2: {}", FileName2);
891
892 SaveName = FileName1;
893 SaveName = SaveName.substr(0, SaveName.find(".root"));
894 if(FileName2 != "") SaveName += FileName2;
895 if(FileName2 != "") SaveName = SaveName.substr(0, SaveName.find(".root"));
896 SaveName += "_Violin";
897 if(PlotAssym) SaveName += "_Assym";
898
899 std::shared_ptr<TFile> File1 = std::make_shared<TFile>(FileName1.c_str());
900 std::shared_ptr<TFile> File2 = nullptr;
901 if(FileName2 != "") File2 = std::make_shared<TFile>(FileName2.c_str());
902
903 canv = new TCanvas("canv", "canv", 1024, 1024);
904 canv->SetGrid();
905 gStyle->SetOptStat(0);
906 //KS: Remove errors on X axis as they are confusing in violin type of plot
907 if(!PlotAssym) gStyle->SetErrorX(0.0001);
908 canv->SetTickx();
909 canv->SetTicky();
910 canv->SetBottomMargin(0.25);
911 canv->SetTopMargin(0.08);
912 canv->SetRightMargin(0.03);
913 canv->SetLeftMargin(0.10);
914 canv->Print((SaveName+".pdf[").c_str());
915 canv->SetGrid();
916
917 Violin = nullptr;
918 ViolinPre = File1->Get<TH2D>( "param_violin_prior" );
919 Violin = File1->Get<TH2D>( "param_violin" );
920 if(Violin == nullptr)
921 {
922 MACH3LOG_ERROR("Couldn't find violin plot, make sure method from MCMCProcessor is being called");
923 return;
924 }
925 // Do some fancy replacements
926 ViolinPre->SetFillColor(kRed);
927 ViolinPre->SetFillColorAlpha(kRed, 0.35);
928 ViolinPre->SetMarkerColor(kRed);
929 ViolinPre->SetMarkerStyle(20);
930 ViolinPre->SetMarkerSize(0.5);
931
932 ViolinPre->GetYaxis()->SetTitleOffset(1.3);
933 ViolinPre->GetYaxis()->SetTitle("Parameter Value");
934 ViolinPre->GetXaxis()->LabelsOption("v");
935
936 Violin->SetFillColor(kBlue);
937 Violin->SetFillColorAlpha(kBlue, 0.35);
938 Violin->SetMarkerColor(kBlue);
939 Violin->SetMarkerStyle(20);
940 Violin->SetMarkerSize(0.5);
941
942 TH1D* Postfit = File1->Get<TH1D>( ("param_xsec_"+plotType).c_str() );
943 Postfit->SetMarkerColor(kRed);
944 Postfit->SetLineColor(kRed);
945 Postfit->SetMarkerStyle(7);
946
947 if(File2 != nullptr)
948 {
949 Violin2 = File2->Get<TH2D>( "param_violin" );
950 Violin2->SetMarkerColor(kGreen);
951 Violin2->SetLineColor(kGreen);
952 Violin2->SetFillColor(kGreen);
953 Violin2->SetFillColorAlpha(kGreen, 0.35);
954 }
955
956 TGraphAsymmErrors* PostGraphAll = MakeTGraphAsymmErrors(File1);
957
958 // Make a Legend page
959 auto leg = std::make_unique<TLegend>(0.0, 0.0, 1.0, 1.0);
960 if (ViolinPre != nullptr) leg->AddEntry(ViolinPre, "Prior", "lpf");
961 if (Violin != nullptr) leg->AddEntry(Violin, "Posterior", "lpf");
962 if (Violin2 != nullptr) leg->AddEntry(Violin2, "Second Violin", "lpf");
963 if(PlotAssym) leg->AddEntry(PostGraphAll, "HPD Assym", "lp");
964 else leg->AddEntry(Postfit, "HPD", "lpf");
965
966 canv->cd();
967 canv->Clear();
968 leg->Draw();
969 canv->Print((SaveName+".pdf").c_str());
970
971 // get the names of the blocks of parameters to group together
972 std::vector<std::string> const blockNames = man->getOption<std::vector<std::string>>("paramGroups");
973 const int XsecPlots = static_cast<int>(blockNames.size());
974
975 for (int i = 0; i < XsecPlots; i++)
976 {
977 // get the configuration for this parameter
978 std::string blockName = blockNames[i];
979 YAML::Node paramBlock = man->getOption(blockName);
980 std::string blockTitle = paramBlock[0].as<std::string>();
981 std::vector<double> blockLimits = paramBlock[1].as<std::vector<double>>();
982 std::vector<std::string> blockContents = paramBlock[2].as<std::vector<std::string>>();
983
984 // get num of params in the block
985 const int nParams = static_cast<int>(blockContents.size());
986
987 // set some plot things
988 auto blockHist_prefit = std::make_unique<TH2D>(blockName.c_str(), blockTitle.c_str(), nParams, 0.0, static_cast<double>(nParams),
989 ViolinPre->GetYaxis()->GetNbins(), ViolinPre->GetYaxis()->GetXmin(), ViolinPre->GetYaxis()->GetXmax());
990 CopyViolinToBlock(ViolinPre, blockHist_prefit.get(), blockContents);
991 // set the y axis limits we got from config
992 blockHist_prefit->GetYaxis()->SetRangeUser(blockLimits[0], blockLimits[1]);
993 auto blockHist_Violin1 = std::make_unique<TH2D>((blockTitle + "Violin1").c_str(), (blockTitle + "Violin1").c_str(), nParams, 0.0, static_cast<double>(nParams),
994 Violin->GetYaxis()->GetNbins(), Violin->GetYaxis()->GetXmin(), Violin->GetYaxis()->GetXmax());
995 CopyViolinToBlock(Violin, blockHist_Violin1.get(), blockContents);
996 blockHist_Violin1->GetYaxis()->SetRangeUser(blockLimits[0], blockLimits[1]);
997 std::unique_ptr<TH2D> blockHist_Violin2 = nullptr;
998 if(Violin2 != nullptr) {
999 blockHist_Violin2 = std::make_unique<TH2D>((blockTitle + "Violin2").c_str(), (blockTitle + "Violin2").c_str(), nParams, 0.0, static_cast<double>(nParams),
1000 Violin2->GetYaxis()->GetNbins(), Violin2->GetYaxis()->GetXmin(), Violin2->GetYaxis()->GetXmax());
1001 CopyViolinToBlock(Violin2, blockHist_Violin2.get(), blockContents);
1002 blockHist_Violin2->GetYaxis()->SetRangeUser(blockLimits[0], blockLimits[1]);
1003 }
1004 // Do some fancy replacements
1005 PrettifyTitles(blockHist_prefit.get());
1006
1007 // set some plot things
1008 auto blockHist_Best = std::make_unique<TH1D>(blockName.c_str(), blockTitle.c_str(),
1009 nParams, 0.0, static_cast<double>(nParams));
1010 // set the errors for the prefit block hist
1011 for(int localBin=0; localBin < nParams; localBin ++){
1012 // the "local" bin is the params index within the group of parameters
1013 std::string paramName = blockContents[localBin];
1014 copyParToBlockHist(localBin, paramName, blockHist_Best.get(), "", 0);
1015 }
1016
1017 std::vector<int> Index;
1018 for(unsigned int is = 0; is < blockContents.size(); is++) {
1019 int ParamBinId = -999;
1020 for (int ix = 0; ix < ViolinPre->GetXaxis()->GetNbins(); ++ix) {
1021 if(ViolinPre->GetXaxis()->GetBinLabel(ix+1) == blockContents[is]) {
1022 ParamBinId = ix;
1023 break;
1024 }
1025 }
1026 Index.push_back(ParamBinId);
1027 }
1028 TGraphAsymmErrors* PostGraph = MakeTGraphAsymmErrors(File1, Index);
1029 PostGraph->SetMarkerColor(kBlack);
1030 PostGraph->SetLineColor(kBlack);
1031 PostGraph->SetMarkerStyle(7);
1032 PostGraph->SetLineWidth(2);
1033 PostGraph->SetLineStyle(kSolid);
1034
1035 blockHist_prefit->Draw("violinX(03100300)");
1036 blockHist_Violin1->Draw("violinX(03100300) SAME");
1037 if(blockHist_Violin2 != nullptr) {
1038 blockHist_Violin2->Draw("violinX(03100300) SAME");
1039 }
1040 if(PlotAssym) PostGraph->Draw("P SAME");
1041 else Postfit->Draw("SAME");
1042 canv->Print((SaveName+".pdf").c_str());
1043 delete PostGraph;
1044 }
1045
1046 canv->Print((SaveName+".pdf]").c_str());
1047 delete canv;
1048 delete ViolinPre;
1049 delete Violin;
1050 delete PostGraphAll;
1051 if(Violin2 != nullptr) delete Violin2;
1052 delete Postfit;
1053 File1->Close();
1054 if(File2 != nullptr) {
1055 File2->Close();
1056 }
1057}
1058
1059int main(int argc, char *argv[])
1060{
1062 // Avoid Info in <TCanvas::Print>
1063 gErrorIgnoreLevel = kWarning;
1064
1066 man->parseInputs(argc, argv);
1067 #ifdef DEBUG
1068 std::cout << std::endl << std::endl << "====================" << std::endl;
1069 man->input().getFile(0).file->ls();
1070 #endif
1071 man->setExec("GetPostfitParamPlots");
1072
1073 man->style().setPalette(man->getOption<std::string>("colorPalette"));
1074
1076
1077 if (man->input().getNInputFiles() == 1)
1078 {
1079 std::string filename = man->getFileName(0);
1080 GetViolinPlots(filename);
1081 }
1082 else if (man->input().getNInputFiles() == 2)
1083 {
1084 std::string filename1 = man->getFileName(0);
1085 std::string filename2 = man->getFileName(1);
1086 GetViolinPlots(filename1, filename2);
1087 }
1088 else if (man->input().getNInputFiles() == 3)
1089 {
1090 std::string filename1 = man->getFileName(0);
1091 std::string filename2 = man->getFileName(1);
1092 std::string filename3 = man->getFileName(3);
1093 //KS: Violin plot currently not supported by three file version although it should be super easy to adapt
1094 }
1095
1096 delete man;
1097 return 0;
1098}
#define _MaCh3_Safe_Include_Start_
KS: Avoiding warning checking for headers.
Definition: Core.h:106
#define _MaCh3_Safe_Include_End_
KS: Restore warning checking after including external headers.
Definition: Core.h:117
std::vector< std::string > ParamNames
int main(int argc, char *argv[])
int FDParameters
void DrawPlots(TCanvas *plotCanv, TH1D *PrefitCopy, const std::vector< TH1D * > &PostfitVec, TPad *mainPad, TPad *ratioPad)
bool ReadSettings(std::shared_ptr< TFile > File1)
std::vector< int > NDSamplesBins
TH1D * Prefit
void InitializePads(TCanvas *canvas, TPad *&pad3, TPad *&pad4)
int XsecStartingPos
TH2D * ViolinPre
int CrossSectionParameters
void MakeNDDetPlots()
int FDParametersStartingPos
std::string SaveName
TPad * p3
TH1D * makeRatio(TH1D *PrefitCopy, TH1D *PostfitCopy, bool setAxes)
void GetViolinPlots(std::string FileName1="", std::string FileName2="")
void MakeFluxPlots()
TH2D * Violin
MaCh3Plotting::PlottingManager * man
int NDParametersStartingPos
Color_t PlotColor[]
std::vector< std::string > NDSamplesNames
TH2D * Violin2
TPad * p4
void PrettifyTitles(TH1D *Hist)
void GetPostfitParamPlots()
int NDParameters
std::vector< TH1D * > PostfitHistVec
int FluxParameters
void copyParToBlockHist(int localBin, std::string paramName, TH1D *blockHist, std::string type, int fileId, bool setLabels=true)
TGraphAsymmErrors * MakeTGraphAsymmErrors(std::shared_ptr< TFile > File, std::vector< int > Index={})
void CopyViolinToBlock(TH2D *FullViolin, TH2D *ReducedViolin, const std::vector< std::string > &ParamNames)
TCanvas * canv
std::string plotType
void MakeXsecPlots()
void MakeXsecRidgePlots()
void MakeFDDetPlots()
#define MACH3LOG_CRITICAL
Definition: MaCh3Logger.h:26
#define MACH3LOG_DEBUG
Definition: MaCh3Logger.h:22
#define MACH3LOG_ERROR
Definition: MaCh3Logger.h:25
#define MACH3LOG_INFO
Definition: MaCh3Logger.h:23
void SetMaCh3LoggerFormat()
Set messaging format of the logger.
Definition: MaCh3Logger.h:30
#define MACH3LOG_WARN
Definition: MaCh3Logger.h:24
Custom exception class for MaCh3 errors.
double getPostFitError(int fileNum, const std::string &paramName, std::string errorType="") const
Get the post fit error for a particular parameter from a particular input file.
InputFile const & getFile(int fileId) const
Definition: inputManager.h:498
std::string translateName(int fileId, fileTypeEnum fileType, std::string paramName) const
Definition: inputManager.h:499
double getPostFitValue(int fileNum, const std::string &paramName, std::string errorType="") const
Get the post fit value for a particular parameter from a particular input file.
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 getFileName(int i)
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 setTH1Style(TH1 *hist, std::string styleName) const
Set the style of a TH1 to one of the styles defined in the style config.
void setPalette(int rootPlotStyle) const
Set the root colour palette to one of the default root pallettes as defined in (root docs)[https://ro...
@ kPostFit
Processed post fit errors.
Definition: inputManager.h:18
std::shared_ptr< TFile > file
Pointer to the underlying file for this InputFile instance.
Definition: inputManager.h:111