ATLAS Offline Software
Functions
JMR_MakeUncertaintyPlots.cxx File Reference
#include <filesystem>
#include <iostream>
#include <string>
#include <fstream>
#include <vector>
#include <cmath>
#include <memory>
#include <TROOT.h>
#include <TStyle.h>
#include <TFile.h>
#include <TTree.h>
#include <TCanvas.h>
#include <TLegend.h>
#include <TGraphErrors.h>
#include <TLorentzVector.h>
#include <TMath.h>
#include <TH1D.h>
#include <TH2D.h>
#include <TH2F.h>
#include <TH3D.h>
#include <TH3F.h>
#include <TRandom3.h>
#include <TString.h>
#include <TLatex.h>
#include <TLine.h>
#include <TPave.h>
#include <TPad.h>
#include <TMarker.h>
#include <TSystem.h>
#include <Math/Interpolator.h>
#include <JetUncertainties/JetHelpers.h>

Go to the source code of this file.

Functions

double Interpolate2D (const TH2 *histo, double x, double y)
 
double Read3DHistogram (const TH3 *histo, double x, double y, double z)
 
std::vector< double > makeLogBins (int nbins, double min, double max)
 
void MakeJMRPlot (int atlas_approved, const std::string &histofile, const std::string &uncmodel, const std::string &kind_of_mass, const std::string &kind_of_mc, bool dointerpol, int mass_bin, double eta_bin, bool pTbin=false)
 
int main (int argc, char *argv[])
 

Function Documentation

◆ Interpolate2D()

double Interpolate2D ( const TH2 *  histo,
double  x,
double  y 
)

Definition at line 51 of file JMR_MakeUncertaintyPlots.cxx.

51  {
52 
53  if (not histo){
54  std::cout << "Histogram pointer is null in FFJetSmearingTool::Interpolate2D" << std::endl;
55  return 0.;
56  }
57 
58  double aux_x = x;
59  double aux_y = y;
60 
62  // If the value is outside the histogram region, we take the closest value to that one
63 
64  double xMax = histo->GetXaxis()->GetBinLowEdge(histo->GetNbinsX()+1);
65  double xMin = histo->GetXaxis()->GetBinLowEdge(1);
66  double yMax = histo->GetYaxis()->GetBinLowEdge(histo->GetNbinsY()+1);
67  double yMin = histo->GetYaxis()->GetBinLowEdge(1);
68 
69  if (x >= xMax) aux_x = xMax-1e-6 ; //so it fits the up-most x-bin
70  if (x <= xMin) aux_x = xMin+1e-6 ; //so it fits the low-most x-bin
71  if (y >= yMax) aux_y = yMax-1e-6 ; //so it fits the up-most y-bin
72  if (y <= yMin) aux_y = yMin+1e-6 ; //so it fits the low-most y-bin
73 
74  Int_t bin_x = histo->GetXaxis()->FindFixBin(aux_x);
75  Int_t bin_y = histo->GetYaxis()->FindFixBin(aux_y);
76  if (bin_x<1 || bin_x>histo->GetNbinsX() || bin_y<1 || bin_y>histo->GetNbinsY()){
77  std::cout << "The point is outside the histogram domain." << std::endl;
78  return 0.;
79  }
80 
81  double interpolated_value = JetHelpers::Interpolate(histo, aux_x, aux_y);
82  return interpolated_value;
83 
84 }

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 568 of file JMR_MakeUncertaintyPlots.cxx.

568  {
569 
570  // Check if enough arguments are passed
571  if (argc < 5){
572  std::cerr << "Usage: " << argv[0] << " <input1: unc. model> <input2: kind_of_mc> <input3: dointerpol> <input4: eta_bin> <input5: histograms file>" << std::endl;
573  return 1; // Exit with an error code
574  }
575 
576  std::cout << "\n Inputs instructions:" << std::endl;
577  std::cout << " * <input1: unc. model> : can be 'SimpleJMR' or 'FullJMR' " << std::endl;
578  std::cout << " * <input2: kind_of_mc> : can be 'MC20', 'MC21', 'MC23', 'MC20AF3' " << std::endl;
579  std::cout << " * <input3: dointerpol> : can be 'true' or 'false' " << std::endl;
580  std::cout << " * <input4: eta_bin> : can be any positive decimal number below 2.0 " << std::endl;
581  std::cout << " * <input5: histograms file> : your '.root' conatining the uncertainties " << std::endl;
582 
583  // Example: "JMR_MakeUncertaintyPlots SimpleJMR MC20 true 0.1 ./source/JetUncertainties/share/DCM_240502_new/R10_JMR_AllComponents_Phase1.root"
584 
586  // Daniel's settings - JMR Phase-I R22 2024 recommendation //
588 
593 
594  // Type of label ---> 0="Internal", 1= "Preliminar", 2="ATLASLabel"(only "ATLAS")
595  int atlas_approved = 0;
596 
597  // Define uncertainty model and conditions
598  std::string kind_of_mass = "UFO";
599  std::string uncmodel = argv[1]; // "SimpleJMR" or "FullJMR"
600  std::string kind_of_mc = argv[2]; // "MC20", "MC21", "MC23", or "MC20AF3"
601  std::string stinterpol = argv[3]; // "true", "false"
602 
603  bool dointerpol = true;
604  if (stinterpol == "true"){
605  dointerpol = true;
606  }else{
607  dointerpol = false;
608  }
609 
610  // - The |eta| value between [the input has 16 bins from 0 and 3.0]
611  double eta_bin = std::stod(argv[4]); //0.1, 0.5, 1.05, 1.55, 1.99
612 
614  // - If false (true) the plot is done as a function of jet pT (mass)
615  bool pTbin = false;
616  // - The mass bin (FF has 2 measurements with 50 < mass < 120 GeV or 120 < mass < 300 GeV)
617  int mass_bin = 0;
618 
620  // Note: this has to be your custom input root file, user and calibration dependent.
621  std::string histofile = argv[5];
622 
624  // Plotting functions //
626 
627  // Plots as a function of jet pT
628  pTbin = false;
629  mass_bin = 40;
630  MakeJMRPlot(atlas_approved, histofile, uncmodel, kind_of_mass, kind_of_mc, dointerpol, mass_bin, eta_bin, pTbin);
631 
632  // Plots as a function of jet mass
633  pTbin = true;
634  mass_bin = 200;
635  MakeJMRPlot(atlas_approved, histofile, uncmodel, kind_of_mass, kind_of_mc, dointerpol, mass_bin, eta_bin, pTbin);
636 
637  return 0;
638 
639 }

◆ MakeJMRPlot()

void MakeJMRPlot ( int  atlas_approved,
const std::string &  histofile,
const std::string &  uncmodel,
const std::string &  kind_of_mass,
const std::string &  kind_of_mc,
bool  dointerpol,
int  mass_bin,
double  eta_bin,
bool  pTbin = false 
)

Definition at line 133 of file JMR_MakeUncertaintyPlots.cxx.

133  {
134 
135  // Make a log-binned histogram: used to read the uncertainties originally
136  std::vector<double> npTbinning;
137  double original_ptmin = 0, original_ptmax = 0;
138  if (pTbin){
139  original_ptmin = 40, original_ptmax = 600;
140  npTbinning = makeLogBins(560,original_ptmin,original_ptmax);
141  }else{
142  original_ptmin = 200, original_ptmax = 3000;
143  npTbinning = makeLogBins(2800,original_ptmin,original_ptmax);
144  }
145  int npTbins = -1+npTbinning.size();
146 
147  std::unique_ptr<TFile> f_JMR_uncertainties_histograms = std::make_unique<TFile>(histofile.c_str(), "READ");
148  if (f_JMR_uncertainties_histograms == 0){
149  std::cout << "Can't find input file.\n";
150  return;
151  }
152 
153  std::cout << "\n > Input root file opened \n" << std::endl;
154 
155  std::vector<std::string> UFO_JMR_uncertainties_histo_names = {};
156  std::vector<std::string> uncertainty_names = {};
157  std::vector<std::string> uncertainty_treatment = {};
158 
159  if (uncmodel == "SimpleJMR"){
160 
161  UFO_JMR_uncertainties_histo_names = {"JMRUnc_mc20vsmc21_"+kind_of_mc+"_PreRec", "JMRUnc_mc20vsmc23_"+kind_of_mc+"_PreRec", "JMRUnc_mc20FullSimvsmc20AF3_"+kind_of_mc+"_PreRec", "JMRUnc_Noise_PreRec", "JMRUnc_mc20UFOvsmc16UFO_PreRec", "JMR_FlatSmearing20_AntiKt10UFOCSSKSoftDropBeta100Zcut10"};
162 
163  uncertainty_names = {"MC20 vs MC21 ("+kind_of_mc+")", "MC20 vs MC23 ("+kind_of_mc+")", "FullSim vs AF3 ("+kind_of_mc+")", "Noise", "MC20 UFO vs MC16 UFO", "Flat smearing 20%"};
164 
165  uncertainty_treatment = {"eLOGmOeAbsEta", "eLOGmOeAbsEta", "PtAbsMass", "eLOGmOeAbsEta", "eLOGmOeAbsEta", "PtAbsMass"};
166 
167  }else if (uncmodel == "FullJMR"){
168 
169  // Not including the Topology as it is only applied for QCD-like jets
170 
171  UFO_JMR_uncertainties_histo_names = {"JMRUnc_mc20vsmc21_"+kind_of_mc+"_PreRec", "JMRUnc_mc20vsmc23_"+kind_of_mc+"_PreRec", "JMRUnc_mc20FullSimvsmc20AF3_"+kind_of_mc+"_PreRec", "JMRUnc_Noise_PreRec", "JMRUnc_mc20UFOvsmc16UFO_PreRec", "JMR_UncertaintiesOutsideCalibratedRegion_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_Smoothing_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_Nom_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_Stat_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_Det_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_MCisr_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_MCfsr_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_MCps_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_MCcr_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_MCrecoil_AntiKt10UFOCSSKSoftDropBeta100Zcut10", "JMR_FF_MCpthard_AntiKt10UFOCSSKSoftDropBeta100Zcut10"};
172 
173  uncertainty_names = {"MC20 vs MC21 ("+kind_of_mc+")", "MC20 vs MC23 ("+kind_of_mc+")", "FullSim vs AF3 ("+kind_of_mc+")", "Noise", "MC20 UFO vs MC16 UFO", "Flat smearing 20% (outside)", "FF smoothing", "FF nom.", "FF stat.", "FF Detector", "FF ISR", "FF FSR", "FF PS", "FF Color recon.", "FF recoil", "FF pthard"};
174 
175  uncertainty_treatment = {"eLOGmOeAbsEta", "eLOGmOeAbsEta", "PtAbsMass", "eLOGmOeAbsEta", "eLOGmOeAbsEta", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass", "PtAbsMass"};
176 
177  }
178 
179  std::cout << " > Uncertainty model: " << uncmodel << "\n" << std::endl;
180 
181  if (uncertainty_names.size() != UFO_JMR_uncertainties_histo_names.size()){
182  std::cout << "List of syst histograms and syst names of different sizes.\n";
183  return;
184  }
185 
186  std::cout << " > Input uncertainty names and histograms provided \n" << std::endl;
187 
188  // A vector with all the histograms for the large-R jet JMR uncertainties (from the tool)
189  std::vector<TH3F*> JMR_uncertainties_hists_3d;
190  int histocount3d = 0;
191  std::vector<TH2D*> JMR_uncertainties_hists_2d;
192  int histocount2d = 0;
193  int iprime[50] = {};
194 
195  for (int i=0; i < std::ssize(uncertainty_names); i++){
196 
197  TString histoname = UFO_JMR_uncertainties_histo_names[i].c_str();
198  TString unctreatment = uncertainty_treatment[i].c_str();
199  if (f_JMR_uncertainties_histograms->Get(histoname) != 0){
200 
201  if (unctreatment == "eLOGmOeAbsEta"){
202  JMR_uncertainties_hists_3d.push_back((TH3F*)f_JMR_uncertainties_histograms->Get(histoname));
203  iprime[i] = histocount3d;
204  histocount3d += 1;
205  }else if (unctreatment == "PtAbsMass"){
206  JMR_uncertainties_hists_2d.push_back((TH2D*)f_JMR_uncertainties_histograms->Get(histoname));
207  iprime[i] = histocount2d;
208  histocount2d += 1;
209  }
210  std::cout << " >> NP read '" << UFO_JMR_uncertainties_histo_names[i] << "' with iprime["<<i<<"] = " <<iprime[i]<< std::endl;
211  std::cout << " >>> Uncertainty treatment = " << unctreatment << std::endl;
212 
213  }else{
214  std::cout << "Can't read input histo." << UFO_JMR_uncertainties_histo_names[i] << "\n";
215  return;
216  }
217 
218  }
219 
220  std::cout << "\n > Input 2D and 3D histograms are read now \n" << std::endl;
221 
222  // Make histogram in 1Dim for the plots
223  std::vector<std::unique_ptr<TH1D>> JMR_uncertainties_1Dhists;
224 
225  for (int i=0; i < std::ssize(uncertainty_names); i++){
226 
227  JMR_uncertainties_1Dhists.push_back(std::make_unique<TH1D>(("UFO_"+std::to_string(i)).c_str(),("UFO_"+std::to_string(i)).c_str(), npTbins, &npTbinning[0]));
228  auto* hist0 = JMR_uncertainties_1Dhists.at(i).get(); // Cache the pointer
229 
230  TString unctreatment = uncertainty_treatment[i].c_str();
231  float uncertainty_value = 0;
232 
233  if (unctreatment == "eLOGmOeAbsEta"){
234 
235  auto* hist1 = JMR_uncertainties_hists_3d.at(iprime[i]); // Cache the pointer
236 
237  for (int bin=1; bin<=npTbins; ++bin){
238  double j = hist0->GetBinCenter(bin);
239 
240  // Calculate the hyperbolic cosine of eta and the energy
241  double cosh_eta = std::cosh(eta_bin);
242  double energy_bin = sqrt((mass_bin*mass_bin) + (j*j)*pow(cosh_eta,2));
243 
244  // If we are out of the range of the input histograms we just take the last value
245  double x = energy_bin, y = TMath::Log(mass_bin/energy_bin), z = eta_bin;
246 
247  Int_t bin_energy = hist1->GetXaxis()->FindBin(x);
248  Int_t bin_mass = hist1->GetYaxis()->FindBin(y);
249  Int_t bin_eta = hist1->GetZaxis()->FindBin(z);
250 
251  // If true, the plot is done as a function of mass instead, and we interprete "mass_bin" as "pT_bin"
252  if (pTbin){
253  energy_bin = sqrt((j*j) + (mass_bin*mass_bin)*pow(cosh_eta,2));
254 
255  // If we are out of the range of the input histograms we just take the last value
256  x = energy_bin, y = TMath::Log(j/energy_bin), z = eta_bin;
257 
258  bin_energy = hist1->GetXaxis()->FindBin(x);
259  bin_mass = hist1->GetYaxis()->FindBin(y);
260  bin_eta = hist1->GetZaxis()->FindBin(z);
261  }
262 
263  // In principle we should not take the fabs(), since there's UP and DOWN, but this is just for illustration
264  if (dointerpol){
265  uncertainty_value = fabs( Read3DHistogram(hist1, x, y, z) ); // Interpolating as in the tool
266  }else{
267  uncertainty_value = fabs( hist1->GetBinContent(bin_energy, bin_mass, bin_eta) ); // No interpolation
268  }
269 
270  // Save the systematic uncertainty
271  hist0->SetBinContent(bin, uncertainty_value);
272 
273  }
274 
275  }else if (unctreatment == "PtAbsMass"){
276 
277  auto* hist2 = JMR_uncertainties_hists_2d.at(iprime[i]); // Cache the pointer
278 
279  for (int bin=1; bin<=npTbins; ++bin){
280 
281  double j = hist0->GetBinCenter(bin);
282 
283  // If we are out of the range of the input histograms we just take the last value
284  double x = j, y = mass_bin;
285 
286  Int_t bin_pt = hist2->GetXaxis()->FindBin(x);
287  Int_t bin_mass = hist2->GetYaxis()->FindBin(y);
288 
289  // If true, the plot is done as a function of mass instead and the value provided as "mass_bin" is pT
290  if (pTbin){
291  x = mass_bin, y = j;
292 
293  bin_pt = hist2->GetXaxis()->FindBin(x);
294  bin_mass = hist2->GetYaxis()->FindBin(y);
295  }
296 
297  // In principle we should not take the fabs(), since there's UP and DOWN, but this is just for illustration
298  if (dointerpol){
299  uncertainty_value = fabs( Interpolate2D(hist2, x, y) ); // Interpolating as in the tool
300  }else{
301  uncertainty_value = fabs( hist2->GetBinContent(bin_pt, bin_mass) ); // No interpolation
302  }
303 
304  // Save the systematic uncertainty
305  hist0->SetBinContent(bin, uncertainty_value);
306 
307  }
308 
309  }
310 
311  }
312 
313  std::cout << " > Individual uncertainties evaluated \n" << std::endl;
314 
316 
317  std::unique_ptr<TH1D> JMR_All_uncertainties_1Dhist = std::make_unique<TH1D>("total_uncertainty", "total_uncertainty", npTbins, &npTbinning[0]);
318 
319  for (int j=1; j<=JMR_All_uncertainties_1Dhist->GetXaxis()->GetNbins(); j++){
320 
321  double total_uncertainty = 0;
322 
323  for (int i=0; i < std::ssize(uncertainty_names); i++){
324 
325  double uncval = JMR_uncertainties_1Dhists.at(i)->GetBinContent(j);
326 
327  total_uncertainty += (uncval) * (uncval);
328  //std::cout << " total error bin (comp: i = "<<i<<", j = "<<j<<") = " << uncval << std::endl;
329 
330  }
331 
332  total_uncertainty = sqrt(total_uncertainty);
333  JMR_All_uncertainties_1Dhist->SetBinContent(j, total_uncertainty);
334  //std::cout << " total error bin (j = "<<j<<") = " << total_uncertainty << std::endl;
335 
336  }
337 
338  // The constant plot with the old uncertainty (0.2 for all pt bins)
339 
340  std::unique_ptr<TH1D> JMR_old_uncertainties_1Dhist = std::make_unique<TH1D>("Old_uncertainty","Old_uncertainty", npTbins, &npTbinning[0]);
341 
342  for (int j=1; j<=JMR_old_uncertainties_1Dhist->GetXaxis()->GetNbins(); j++){
343  JMR_old_uncertainties_1Dhist->SetBinContent(j, 0.2);
344  }
345 
346  std::cout << " > Total uncertainty calculated \n" << std::endl;
347 
348  // We will also need a total graph with the sum
349 
350  gStyle->SetOptStat(0);
351  gROOT->Reset();
352  gROOT->SetStyle("ATLAS");
353  gROOT->ForceStyle();
354  gStyle->SetPadLeftMargin(0.15);
355  gStyle->SetPadRightMargin(0.08);
356  gStyle->SetEndErrorSize(4);
357  gStyle->SetPaintTextFormat("4.3f");
358  gStyle->SetPalette(8,0, 0.3);
359 
360  std::unique_ptr<TCanvas> canvas = std::make_unique<TCanvas>("c1", "c1", 1200, 800);
361 
362  // An invisible graph to set the limits of the canvas
363  int nbins = (original_ptmax - original_ptmin)/(1.);
364  std::unique_ptr<TH1D> DrawHistogram = std::make_unique<TH1D>("draw_hist","draw_hist_title",nbins, original_ptmin, original_ptmax);
365 
366  DrawHistogram->SetMinimum(0);
367  DrawHistogram->SetMaximum(1.0);
368  DrawHistogram->SetLineWidth(0);//invisible histogram to set the limits of the plot
369  if (!pTbin){
370  canvas->SetLogx();//Logaritmic scale in the X axis
371  }
372 
373  DrawHistogram->GetXaxis()->SetTitle("p_{T}^{jet} [GeV]");
374  DrawHistogram->GetXaxis()->SetLabelOffset(0.0125);
375  DrawHistogram->GetXaxis()->SetNdivisions(510);
376  DrawHistogram->GetXaxis()->SetMoreLogLabels(true);
377  if (pTbin){
378  DrawHistogram->GetXaxis()->SetTitle("m^{jet} [GeV]");
379  DrawHistogram->GetXaxis()->SetNdivisions(505);
380  DrawHistogram->GetXaxis()->SetMoreLogLabels(false);
381  }
382  DrawHistogram->GetXaxis()->SetTitleOffset(1.50);
383 
384  DrawHistogram->GetYaxis()->SetTitle("Fractional JMR uncertainty");
385  DrawHistogram->GetYaxis()->SetLabelOffset(0.0145);
386  DrawHistogram->GetYaxis()->SetTitleOffset(1.40);
387 
388  // The numbers in the axis are still wrong
389  DrawHistogram->Draw();
390 
392 
393  JMR_old_uncertainties_1Dhist->SetLineColor(kTeal+5);
394  JMR_old_uncertainties_1Dhist->SetLineStyle(1);
395  JMR_old_uncertainties_1Dhist->SetLineWidth(2);
396  JMR_old_uncertainties_1Dhist->SetFillStyle(4000);
397  JMR_old_uncertainties_1Dhist->SetFillColor(kTeal-9);
398  JMR_old_uncertainties_1Dhist->Draw("hist same");
399 
400  JMR_All_uncertainties_1Dhist->SetLineColor(kBlue+2);
401  JMR_All_uncertainties_1Dhist->SetLineStyle(1);
402  JMR_All_uncertainties_1Dhist->SetLineWidth(2);
403  JMR_All_uncertainties_1Dhist->SetFillStyle(4000);
404  JMR_All_uncertainties_1Dhist->SetFillColor(590);
405  JMR_All_uncertainties_1Dhist->Draw("hist same");
406 
407  for (int index=0; index < std::ssize(uncertainty_names); index++){
408  int indexp = index + 2;
409  if (index == 8) indexp = 14;
410  if (index == 9) indexp = 30;
411  if (index == 10) indexp = 43;
412  if (index == 11) indexp = 46;
413 
414  auto* hist0 = JMR_uncertainties_1Dhists.at(index).get(); // Cache the pointer
415  hist0->SetLineColor(indexp);
416  hist0->SetLineStyle(index+2);
417  hist0->SetLineWidth(2);
418  hist0->Draw("hist same");
419  }
420 
422 
423  // ATLAS label
424 
425  TString ATLAS_LABELLING="";
426  if (atlas_approved==0){
427  ATLAS_LABELLING = "Internal";
428  }else if (atlas_approved==1)
429  ATLAS_LABELLING = "Preliminary";
430  else if (atlas_approved==2){
431  ATLAS_LABELLING = "";
432  }
433 
434  TLatex p;
435  TLatex l;
436  l.SetNDC();
437  l.SetTextFont(72);
438  l.SetTextSize(0.05);
439  //l.SetTextColor(color);
440  p.SetNDC();
441  p.SetTextFont(42);
442  p.SetTextSize(0.05);
443  //ATLAS_LABEL(0.2,0.6,kBlack);
444 
445  l.DrawLatex(0.67,0.89,"ATLAS");
446  p.DrawLatex(0.67+0.115,0.89,ATLAS_LABELLING);
447 
448  // Other information
449 
450  TLatex data_text;
451  data_text.SetNDC();
452  data_text.SetTextFont(42);
453  data_text.SetTextSize(0.04);
454  data_text.DrawLatex(0.18,0.89,"Data 2015-2017, #sqrt{s} = 13 TeV");
455 
456  TLatex jet_info_1_text;
457  jet_info_1_text.SetNDC();
458  jet_info_1_text.SetTextFont(42);
459  jet_info_1_text.SetTextSize(0.04);
460  jet_info_1_text.DrawLatex(0.18,0.89-0.06,"Phase-I unc. for UFO R = 1.0 jets");
461 
462  TLatex jet_info_2_text;
463  jet_info_2_text.SetNDC();
464  jet_info_2_text.SetTextFont(42);
465  jet_info_2_text.SetTextSize(0.04);
466 
467  std::string jet_info_2_string = "";
468  if (kind_of_mass == "UFO"){
469 
470  jet_info_2_string = "|#eta^{jet}| = "+TString(Form("%f",eta_bin))+", m_{UFO}^{jet} = "+TString(Form("%i",mass_bin))+" GeV";
471 
472  if (pTbin){
473  jet_info_2_string = "|#eta^{jet}| = "+TString(Form("%f",eta_bin))+", p_{T}^{jet} = "+TString(Form("%i",mass_bin))+" GeV";
474  }
475 
476  }
477  jet_info_2_text.DrawLatex(0.18,0.89-2*0.06,jet_info_2_string.c_str());
478 
480 
481  std::unique_ptr<TLegend> legend0_1 = std::make_unique<TLegend>(0.17,0.89-3.33*0.06,0.45,0.89-2.33*0.06);
482  legend0_1->SetTextSize(0.035);
483  legend0_1->SetTextFont(42);
484  legend0_1->SetBorderSize(0);
485  legend0_1->SetFillStyle(0);
486  legend0_1->AddEntry(JMR_All_uncertainties_1Dhist.get(), "Tot. uncert.","f");
487  legend0_1->Draw("same");
488 
489  std::unique_ptr<TLegend> legend0_2 = std::make_unique<TLegend>(0.17,0.89-4.33*0.06,0.45,0.89-3.33*0.06);
490  legend0_2->SetTextSize(0.035);
491  legend0_2->SetTextFont(42);
492  legend0_2->SetBorderSize(0);
493  legend0_2->SetFillStyle(0);
494  legend0_2->AddEntry(JMR_old_uncertainties_1Dhist.get(), "Tot. uncert. (previous)","f");
495  legend0_2->Draw("same");
496 
497  // List all the uncertainty components
498  std::unique_ptr<TLegend> legend1;
499  if (uncmodel == "SimpleJMR"){
500  legend1 = std::make_unique<TLegend>(0.55,0.70,0.85,0.88);
501  }else if (uncmodel == "FullJMR"){
502  legend1 = std::make_unique<TLegend>(0.55,0.55,0.85,0.88);
503  }
504  legend1->SetTextSize(0.035);
505  legend1->SetTextFont(42);
506  legend1->SetBorderSize(0);
507  legend1->SetFillStyle(0);
508  for (int index=0; index < std::ssize(uncertainty_names); index++){
509  legend1->AddEntry(JMR_uncertainties_1Dhists.at(index).get(), uncertainty_names[index].c_str(),"l");
510  }
511  legend1->Draw("same");
512 
513  DrawHistogram->Draw("axis same");
514 
515  std::cout << " > Plots are done now \n" << std::endl;
516 
517  std::string outputdir = "./source/JetUncertainties/util/output/";
518 
519  try{
520  // Check if the directory exists; if not, create it
521  if (!std::filesystem::exists(outputdir)) {
522  if (std::filesystem::create_directory(outputdir)){
523  std::cout << "Output directory created successfully: " << outputdir << std::endl;
524  }else{
525  std::cerr << "Failed to create directory: " << outputdir << std::endl;
526  return; // Exit with error code
527  }
528  } else {
529  std::cout << "Output directory already exists: " << outputdir << std::endl;
530  }
531  }catch (const std::filesystem::filesystem_error& e){
532  std::cerr << "Filesystem error: " << e.what() << std::endl;
533  return; // Exit with error code
534  }
535 
536  std::string string_write;
537  std::string string_write_pdf;
538  std::string string_write_png;
539 
540  TString mBinStr = "mass"+std::to_string(mass_bin);
541  if (pTbin){
542  mBinStr = "pT"+std::to_string(mass_bin);
543  }
544 
545  TString EtaBinStr = "Eta"+std::to_string(eta_bin);
546 
547  if (atlas_approved == 0){
548  string_write = outputdir+"FractionalJMRuncertainty_"+kind_of_mc+"_"+kind_of_mass+"_"+uncmodel+"_"+EtaBinStr+"_"+mBinStr+ "GeV_ATLASint";
549  string_write_pdf = string_write+".pdf";
550  string_write_png = string_write+".png";
551  }else if (atlas_approved == 1){
552  string_write = outputdir+"FractionalJMRuncertainty_"+kind_of_mc+"_"+kind_of_mass+"_"+uncmodel+"_"+EtaBinStr+"_"+mBinStr+ "GeV_ATLASprel";
553  string_write_pdf = string_write+".pdf";
554  string_write_png = string_write+".png";
555  }else if (atlas_approved == 2){
556  string_write = outputdir+"FractionalJMRuncertainty_"+kind_of_mc+"_"+kind_of_mass+"_"+uncmodel+"_"+EtaBinStr+"_"+mBinStr+ "GeV_ATLAS";
557  string_write_pdf = string_write+".pdf";
558  string_write_png = string_write+".png";
559  }
560 
561  canvas->Print(string_write_pdf.c_str());
562  //canvas->Print(string_write_png.c_str());
563 
564  return;
565 
566 }

◆ makeLogBins()

std::vector<double> makeLogBins ( int  nbins,
double  min,
double  max 
)

Definition at line 123 of file JMR_MakeUncertaintyPlots.cxx.

123  {
124  std::vector<double> bins;
125  double lmin = log(min), k=(log(max)-log(min))/nbins;
126  for (int i=0;i<=nbins;++i){
127  bins.push_back(exp(lmin+k*i));
128  //std::cout << " bins.at("<<i<<") = " << exp(lmin+k*i) << std::endl;
129  }
130  return bins;
131 }

◆ Read3DHistogram()

double Read3DHistogram ( const TH3 *  histo,
double  x,
double  y,
double  z 
)

Definition at line 87 of file JMR_MakeUncertaintyPlots.cxx.

87  {
88 
89  if (not histo){
90  std::cout << "Histogram pointer is null in FFJetSmearingTool::Read3DHistogram" << std::endl;
91  return 0.;
92  }
93 
94  double aux_x = x;
95  double aux_y = y;
96  double aux_z = z;
97 
99  // If the value is outside the histogram region, we take the closest value to that one
100 
101  double xMax = histo->GetXaxis()->GetBinLowEdge(histo->GetNbinsX()+1);
102  double xMin = histo->GetXaxis()->GetBinLowEdge(1);
103  double yMax = histo->GetYaxis()->GetBinLowEdge(histo->GetNbinsY()+1);
104  double yMin = histo->GetYaxis()->GetBinLowEdge(1);
105  double zMax = histo->GetZaxis()->GetBinLowEdge(histo->GetNbinsZ()+1);
106  double zMin = histo->GetZaxis()->GetBinLowEdge(1);
107 
108  if (x >= xMax) aux_x = xMax-1e-6 ; //so it fits the up-most x-bin
109  if (x <= xMin) aux_x = xMin+1e-6 ; //so it fits the low-most x-bin
110  if ( std::isnan(y)) return 0; // no weight if the input is NaN, can happen for log(X)
111  if (y >= yMax) aux_y = yMax-1e-6 ; //so it fits the up-most y-bin
112  if (y <= yMin) aux_y = yMin+1e-6 ; //so it fits the low-most y-bin
113  if (z >= zMax) aux_z = zMax-1e-6 ; //so it fits the up-most z-bin
114  if (z <= zMin) aux_z = zMin+1e-6 ; //so it fits the low-most z-bin
115 
116  //Use the interpolate function from JetHelpers.cxx
117  double interpolated_value = JetHelpers::Interpolate(histo, aux_x, aux_y, aux_z);
118  return interpolated_value;
119 
120 }
AllowedVariables::e
e
Definition: AsgElectronSelectorTool.cxx:37
index
Definition: index.py:1
max
constexpr double max()
Definition: ap_fixedTest.cxx:33
min
constexpr double min()
Definition: ap_fixedTest.cxx:26
bin
Definition: BinsDiffFromStripMedian.h:43
UploadAMITag.l
list l
Definition: UploadAMITag.larcaf.py:158
Interpolate2D
double Interpolate2D(const TH2 *histo, double x, double y)
Definition: JMR_MakeUncertaintyPlots.cxx:51
drawFromPickle.exp
exp
Definition: drawFromPickle.py:36
x
#define x
JetHelpers::Interpolate
double Interpolate(const TH1 *histo, const double x)
Definition: JetHelpers.cxx:16
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:210
lumiFormat.i
int i
Definition: lumiFormat.py:85
z
#define z
LArCellNtuple.argv
argv
Definition: LArCellNtuple.py:152
DeMoAtlasDataLoss.canvas
dictionary canvas
Definition: DeMoAtlasDataLoss.py:187
plotting.yearwise_luminosity_vs_mu.bins
bins
Definition: yearwise_luminosity_vs_mu.py:30
DQHistogramMergeRegExp.argc
argc
Definition: DQHistogramMergeRegExp.py:20
doL1CaloHVCorrections.eta_bin
eta_bin
Definition: doL1CaloHVCorrections.py:368
makeLogBins
std::vector< double > makeLogBins(int nbins, double min, double max)
Definition: JMR_MakeUncertaintyPlots.cxx:123
PixelAthHitMonAlgCfg.histoname
string histoname
Definition: PixelAthHitMonAlgCfg.py:139
CompareRootFiles.hist1
hist1
Definition: CompareRootFiles.py:36
ActsTrk::to_string
std::string to_string(const DetectorType &type)
Definition: GeometryDefs.h:34
plotBeamSpotVxVal.bin
int bin
Definition: plotBeamSpotVxVal.py:83
SCT_CalibAlgs::nbins
@ nbins
Definition: SCT_CalibNumbers.h:10
DeMoScan.index
string index
Definition: DeMoScan.py:364
MakeJMRPlot
void MakeJMRPlot(int atlas_approved, const std::string &histofile, const std::string &uncmodel, const std::string &kind_of_mass, const std::string &kind_of_mc, bool dointerpol, int mass_bin, double eta_bin, bool pTbin=false)
Definition: JMR_MakeUncertaintyPlots.cxx:133
y
#define y
Read3DHistogram
double Read3DHistogram(const TH3 *histo, double x, double y, double z)
Definition: JMR_MakeUncertaintyPlots.cxx:87
python.CaloCondTools.log
log
Definition: CaloCondTools.py:20
python.dummyaccess.exists
def exists(filename)
Definition: dummyaccess.py:9
plotBeamSpotCompare.histo
histo
Definition: plotBeamSpotCompare.py:415
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
fitman.k
k
Definition: fitman.py:528
CompareRootFiles.hist2
hist2
Definition: CompareRootFiles.py:37