Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DoPercentDifference.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file DoPercentDifference.C
1 // ----------------------------------------------------------------------------
2 // 'DoPercentDifference.C'
3 // Derek Anderson
4 // 04.27.2023
5 //
6 // Calculates the percent difference between multiple
7 // histograms and plots them.
8 //
9 // NOTE: 'FSubMode' controls how the difference is
10 // calculated.
11 // FSubMode = 0: (ref - comp) / ref
12 // FSubMode = 1: (ln(ref) - ln(comp)) / ln(ref)
13 // FSubMode = 2: ref - comp
14 // FSubMode = 3: ln(ref) - ln(comp)
15 // ----------------------------------------------------------------------------
16 
17 // c includes
18 #include <iostream>
19 // root includes
20 #include "TH1.h"
21 #include "TPad.h"
22 #include "TFile.h"
23 #include "TMath.h"
24 #include "TLine.h"
25 #include "TError.h"
26 #include "TString.h"
27 #include "TCanvas.h"
28 #include "TLegend.h"
29 #include "TPaveText.h"
30 
31 using namespace std;
32 
33 // global constants
34 static const UInt_t NModes(4);
35 static const UInt_t NRange(2);
36 static const UInt_t NComp(7);
37 static const UInt_t NPlot(2);
38 static const UInt_t NPad(2);
39 static const UInt_t NVtx(4);
40 static const UInt_t NTxt(3);
41 
42 // calculation parameters
43 static const UInt_t FSubMode(0);
44 static const Bool_t DoRebin(false);
45 static const Bool_t DoIntNorm(true);
46 static const Bool_t DoFluctuationRemoval(true);
47 
48 
49 
50 void DoPercentDifference() {
51 
52  // lower verbosity
53  gErrorIgnoreLevel = kWarning;
54  cout << "\n Beginning percent difference script..." << endl;
55 
56  // output and reference parameters
57  const TString sOut("tuningDetEffects.testingPerDiffVsDiff_perDiffOverComp_pt15.pp200py8jet40.d27m4y2023.root");
58  const TString sRef("tuning_detector_effects/SphenixRecoUnscaledWithConstitHistos_output.root");
59  const TString sHistRef("EEC0");
60  const TString sNameRef("hRecoEEC_pt15");
61  const TString sLabelRef("Reco.");
62 
63  // comparison parameters
64  const TString sComp[NComp] = {"tuning_detector_effects/SphenixTruthShift0Check_output.root",
65  "tuning_detector_effects/SphenixTruthShift.001Check_output.root",
66  "tuning_detector_effects/SphenixTruthShift.005Check_output.root",
67  "tuning_detector_effects/SphenixTruthShift.01Check_output.root",
68  "tuning_detector_effects/SphenixTruthShift.015Check_output.root",
69  "tuning_detector_effects/SphenixTruthShift.02Check_output.root",
70  "tuning_detector_effects/SphenixTruthShift.025Check_output.root"};
71  const TString sHistComp[NComp] = {"EEC0",
72  "EEC0",
73  "EEC0",
74  "EEC0",
75  "EEC0",
76  "EEC0",
77  "EEC0"};
78  const TString sNameComp[NComp] = {"hTrueEEC_pt15_noSmear",
79  "hTrueEEC_pt15_smear001",
80  "hTrueEEC_pt15_smear005",
81  "hTrueEEC_pt15_smear01",
82  "hTrueEEC_pt15_smear015",
83  "hTrueEEC_pt15_smear02",
84  "hTrueEEC_pt15_smear025"};
85  const TString sNameSub[NComp] = {"hSub_pt15_noSmear",
86  "hSub_pt15_smear001",
87  "hSub_pt15_smear005",
88  "hSub_pt15_smear01",
89  "hSub_pt15_smear015",
90  "hSub_pt15_smear02",
91  "hSub_pt15_smear025"};
92  const TString sLabelComp[NComp] = {"Truth (no smearing)",
93  "Truth (smear width = 0.001)",
94  "Truth (smear width = 0.005)",
95  "Truth (smear width = 0.01)",
96  "Truth (smear width = 0.015)",
97  "Truth (smear width = 0.02)",
98  "Truth (smear width = 0.025)"};
99 
100  // calculation parameters
101  const UInt_t nRebin(2);
102  const Float_t rMainPeak[NRange] = {0.001, 1.};
103  const TString sTitleSubs[NModes] = {"\%-difference", "\%-diff. of logs", "difference", "diff. of logs"};
104 
105  // plot parameters
106  const TString sOptRef("");
107  const TString sOptComp[NComp] = {"SAME", "SAME", "SAME", "SAME", "SAME", "SAME", "SAME"};
108  const TString sOptSub[NComp] = {"HIST P", "SAME HIST P", "SAME HIST P", "SAME HIST P", "SAME HIST P", "SAME HIST P", "SAME HIST p"};
109  const Float_t xPlotRange[NPlot] = {0.0001, 1.};
110 
111  // style parameters
112  const TString sTitle("");
113  const TString sTitleX("#DeltaR");
114  const TString sTitleY("Normalized EEC");
115  const UInt_t fColDen(923);
116  const UInt_t fMarDen(20);
117  const UInt_t fColNum[NComp] = {921, 799, 819, 839, 859, 879, 899};
118  const UInt_t fMarNum[NComp] = {24, 26, 32, 25, 27, 28, 30};
119 
120  // text parameters
121  const TString sHeader("");
122  const TString sTxt[NTxt] = {"#bf{#it{sPHENIX}} Simulation [Run6]", "PYTHIA-8, JS 40 GeV jet sample", "#bf{p_{T}^{jet} #in (15, 20) GeV/c}"};
123 
124  // open output and reference files
125  TFile *fOut = new TFile(sOut.Data(), "recreate");
126  TFile *fRef = new TFile(sRef.Data(), "read");
127  if (!fOut || !fRef) {
128  cerr << "PANIC: couldn't open output or reference file!\n"
129  << " fOut = " << fOut << ", fRef = " << fRef
130  << endl;
131  return;
132  }
133  cout << " Opened output and reference files." << endl;
134 
135  // open comparison files
136  TFile *fComp[NComp];
137  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
138  fComp[iComp] = new TFile(sComp[iComp].Data(), "read");
139  if (!fComp[iComp]) {
140  cerr << "PANIC: couldn't open comparison file #" << iComp << "!" << endl;
141  return;
142  }
143  }
144  cout << " Opened comparison files." << endl;
145 
146  // grab reference histogram
147  TH1D *hRef = (TH1D*) fRef -> Get(sHistRef.Data());
148  if (!hRef) {
149  cerr << "PANIC: couldn't grab reference histogram!" << endl;
150  return;
151  }
152  hRef -> SetName(sNameRef.Data());
153  cout << " Grabbed reference histogram." << endl;
154 
155  // grab comparison histograms
156  TH1D *hComp[NComp];
157  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
158  hComp[iComp] = (TH1D*) fComp[iComp] -> Get(sHistComp[iComp]);
159  if (!hComp[iComp]) {
160  cerr << "PANIC: couldn't grab comparison histogram #" << iComp << "!" << endl;
161  return;
162  }
163  hComp[iComp] -> SetName(sNameComp[iComp].Data());
164  }
165  cout << " Grabbed comparison histograms." << endl;
166 
167  // rebin histograms (if needed)
168  if (DoRebin) {
169  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
170  hComp[iComp] -> Rebin(nRebin);
171  }
172  hRef -> Rebin(nRebin);
173  cout << " Rebinned histograms." << endl;
174  }
175 
176  // do normalization (if needed)
177  if (DoIntNorm) {
178 
179  // remove fluctuations (if needed)
180  if (DoFluctuationRemoval) {
181 
182  // remove fluctions for reference histogram
183  {
184 
185  // no. of bins
186  const UInt_t nBins = hRef -> GetNbinsX();
187 
188  // loop over bins and determine height of main peak
189  Double_t yMainPeak = 0.;
190  for (UInt_t iBin = 1; iBin < (nBins + 1); iBin++) {
191 
192  // check if bin is in range
193  const Double_t xBinCenter = hRef -> GetBinCenter(iBin);
194  const Bool_t isInMainPeakRange = ((xBinCenter >= rMainPeak[0]) && (xBinCenter < rMainPeak[1]));
195  if (!isInMainPeakRange) continue;
196 
197  // check if bigger than current max
198  const Double_t yBin = hRef -> GetBinContent(iBin);
199  if (yBin > yMainPeak) {
200  yMainPeak = yBin;
201  }
202  } // end bin loop
203 
204  // loop over bins below main peak and 0 out any which are unusually high
205  for (UInt_t iBin = 1; iBin < (nBins + 1); iBin++) {
206 
207  // check if bin is less than main peak range
208  const Double_t xBinCenter = hRef -> GetBinCenter(iBin);
209  const Bool_t isBelowMainPeak = (xBinCenter < rMainPeak[0]);
210  if (!isBelowMainPeak) continue;
211 
212  // 0 out bin if above main peak
213  const Double_t yBin = hRef -> GetBinContent(iBin);
214  if (yBin >= yMainPeak) {
215  hRef -> SetBinContent(iBin, 0.);
216  hRef -> SetBinError(iBin, 0.);
217  }
218  } // end bin loop
219  } // end reference fluctuation removal
220  cout << " Removed fluctuations for reference histogram." << endl;
221 
222  // loop over comparisons histograms
223  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
224 
225  // no. of bins
226  const UInt_t nBins = hRef -> GetNbinsX();
227 
228  // loop over bins and determine height of main peak
229  Double_t yMainPeak = 0.;
230  for (UInt_t iBin = 1; iBin < (nBins + 1); iBin++) {
231 
232  // check if bin is in range
233  const Double_t xBinCenter = hComp[iComp] -> GetBinCenter(iBin);
234  const Bool_t isInMainPeakRange = ((xBinCenter >= rMainPeak[0]) && (xBinCenter < rMainPeak[1]));
235  if (!isInMainPeakRange) continue;
236 
237  // check if bigger than current max
238  const Double_t yBin = hComp[iComp] -> GetBinContent(iBin);
239  if (yBin > yMainPeak) {
240  yMainPeak = yBin;
241  }
242  } // end bin loop
243 
244  // loop over bins below main peak and 0 out any which are unusually high
245  for (UInt_t iBin = 1; iBin < (nBins + 1); iBin++) {
246 
247  // check if bin is less than main peak range
248  const Double_t xBinCenter = hComp[iComp] -> GetBinCenter(iBin);
249  const Bool_t isBelowMainPeak = (xBinCenter < rMainPeak[0]);
250  if (!isBelowMainPeak) continue;
251 
252  // 0 out bin if above main peak
253  const Double_t yBin = hComp[iComp] -> GetBinContent(iBin);
254  if (yBin >= yMainPeak) {
255  hComp[iComp] -> SetBinContent(iBin, 0.);
256  hComp[iComp] -> SetBinError(iBin, 0.);
257  }
258  } // end bin loop
259  } // end comparison loop
260  cout << " Removed fluctuations for comparison histograms." << endl;
261  } // end if (DoFluctuationRemoval)
262 
263  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
264  const Double_t intComp = hComp[iComp] -> Integral();
265  if (intComp > 0.) hComp[iComp] -> Scale(1. / intComp);
266  }
267 
268  const Double_t intRef = hRef -> Integral();
269  if (intRef > 0.) hRef -> Scale(1. / intRef);
270  cout << " Normalized histograms." << endl;
271  }
272 
273  // do subtractions
274  TH1D *hSub[NComp];
275  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
276 
277  // initialize subtraction histograms
278  hSub[iComp] = (TH1D*) hRef -> Clone();
279  hSub[iComp] -> Reset("ICE");
280  hSub[iComp] -> SetName(sNameSub[iComp]);
281 
282  // loop over bins
283  const UInt_t nBins = hSub[iComp] -> GetNbinsX();
284  for (UInt_t iBin = 1; iBin < (nBins + 1); iBin++) {
285 
286  // get content and errors
287  const Double_t yRef = hRef -> GetBinContent(iBin);
288  const Double_t yComp = hComp[iComp] -> GetBinContent(iBin);
289  const Double_t eRef = hRef -> GetBinError(iBin);
290  const Double_t eComp = hComp[iComp] -> GetBinError(iBin);
291 
292  // do subtraction
293  Double_t ySub(-1.);
294  Double_t eSub(0.);
295  switch (FSubMode) {
296  case 0:
297  ySub = (yRef - yComp) / yComp;
298  eSub = 0; // FIXME do propagation of errors
299  break;
300  case 1:
301  ySub = (TMath::Log(yRef) - TMath::Log(yComp)) / TMath::Log(yRef);
302  eSub = 0; // FIXME do propagation of errors
303  break;
304  case 2:
305  ySub = yRef - yComp;
306  eSub = 0; // FIXME do propagation of errors
307  break;
308  case 3:
309  ySub = TMath::Log(yRef) - TMath::Log(yComp);
310  eSub = 0; // FIXME do propagation of errors
311  break;
312  default:
313  ySub = (yRef - yComp) / yRef;
314  eSub = 0; // FIXME do propagation of errors
315  break;
316  }
317 
318  // set bin content/errors
319  if (yRef > 0.) {
320  hSub[iComp] -> SetBinContent(iBin, ySub);
321  hSub[iComp] -> SetBinError(iBin, eSub);
322  } else {
323  hSub[iComp] -> SetBinContent(iBin, 0.);
324  hSub[iComp] -> SetBinError(iBin, 0.);
325  }
326  } // end bin loop
327  } // end comparison loop
328  cout << " Calculated ratios." << endl;
329 
330  // pick out subtraction title
331  TString sTitleS("");
332  switch (FSubMode) {
333  case 0:
334  sTitleS = sTitleSubs[0];
335  break;
336  case 1:
337  sTitleS = sTitleSubs[1];
338  break;
339  case 2:
340  sTitleS = sTitleSubs[2];
341  break;
342  case 3:
343  sTitleS = sTitleSubs[3];
344  break;
345  default:
346  sTitleS = sTitleSubs[0];
347  break;
348  }
349 
350  // set styles
351  const UInt_t fFil(0);
352  const UInt_t fLin(1);
353  const UInt_t fWid(1);
354  const UInt_t fTxt(42);
355  const UInt_t fAln(12);
356  const UInt_t fCnt(1);
357  const Float_t fLab[NPad] = {0.074, 0.04};
358  const Float_t fTit[NPad] = {0.074, 0.04};
359  const Float_t fOffX[NPad] = {1.1, 1.};
360  const Float_t fOffY[NPad] = {0.7, 1.3};
361  hRef -> SetMarkerColor(fColDen);
362  hRef -> SetMarkerStyle(fMarDen);
363  hRef -> SetFillColor(fColDen);
364  hRef -> SetFillStyle(fFil);
365  hRef -> SetLineColor(fColDen);
366  hRef -> SetLineStyle(fLin);
367  hRef -> SetLineWidth(fWid);
368  hRef -> SetTitle(sTitle.Data());
369  hRef -> SetTitleFont(fTxt);
370  hRef -> GetXaxis() -> SetRangeUser(xPlotRange[0], xPlotRange[1]);
371  hRef -> GetXaxis() -> SetTitle(sTitleX.Data());
372  hRef -> GetXaxis() -> SetTitleFont(fTxt);
373  hRef -> GetXaxis() -> SetTitleSize(fTit[1]);
374  hRef -> GetXaxis() -> SetTitleOffset(fOffX[1]);
375  hRef -> GetXaxis() -> SetLabelFont(fTxt);
376  hRef -> GetXaxis() -> SetLabelSize(fLab[1]);
377  hRef -> GetXaxis() -> CenterTitle(fCnt);
378  hRef -> GetYaxis() -> SetTitle(sTitleY.Data());
379  hRef -> GetYaxis() -> SetTitleFont(fTxt);
380  hRef -> GetYaxis() -> SetTitleSize(fTit[1]);
381  hRef -> GetYaxis() -> SetTitleOffset(fOffY[1]);
382  hRef -> GetYaxis() -> SetLabelFont(fTxt);
383  hRef -> GetYaxis() -> SetLabelSize(fLab[1]);
384  hRef -> GetYaxis() -> CenterTitle(fCnt);
385  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
386  hComp[iComp] -> SetMarkerColor(fColNum[iComp]);
387  hComp[iComp] -> SetMarkerStyle(fMarNum[iComp]);
388  hComp[iComp] -> SetFillColor(fColNum[iComp]);
389  hComp[iComp] -> SetFillStyle(fFil);
390  hComp[iComp] -> SetLineColor(fColNum[iComp]);
391  hComp[iComp] -> SetLineStyle(fLin);
392  hComp[iComp] -> SetLineWidth(fWid);
393  hComp[iComp] -> SetTitle(sTitle.Data());
394  hComp[iComp] -> SetTitleFont(fTxt);
395  hComp[iComp] -> GetXaxis() -> SetRangeUser(xPlotRange[0], xPlotRange[1]);
396  hComp[iComp] -> GetXaxis() -> SetTitle(sTitleX.Data());
397  hComp[iComp] -> GetXaxis() -> SetTitleFont(fTxt);
398  hComp[iComp] -> GetXaxis() -> SetTitleSize(fTit[1]);
399  hComp[iComp] -> GetXaxis() -> SetTitleOffset(fOffX[1]);
400  hComp[iComp] -> GetXaxis() -> SetLabelFont(fTxt);
401  hComp[iComp] -> GetXaxis() -> SetLabelSize(fLab[1]);
402  hComp[iComp] -> GetXaxis() -> CenterTitle(fCnt);
403  hComp[iComp] -> GetYaxis() -> SetTitle(sTitleY.Data());
404  hComp[iComp] -> GetYaxis() -> SetTitleFont(fTxt);
405  hComp[iComp] -> GetYaxis() -> SetTitleSize(fTit[1]);
406  hComp[iComp] -> GetYaxis() -> SetTitleOffset(fOffY[1]);
407  hComp[iComp] -> GetYaxis() -> SetLabelFont(fTxt);
408  hComp[iComp] -> GetYaxis() -> SetLabelSize(fLab[1]);
409  hComp[iComp] -> GetYaxis() -> CenterTitle(fCnt);
410  hSub[iComp] -> SetMarkerColor(fColNum[iComp]);
411  hSub[iComp] -> SetMarkerStyle(fMarNum[iComp]);
412  hSub[iComp] -> SetFillColor(fColNum[iComp]);
413  hSub[iComp] -> SetFillStyle(fFil);
414  hSub[iComp] -> SetLineColor(fColNum[iComp]);
415  hSub[iComp] -> SetLineStyle(fLin);
416  hSub[iComp] -> SetLineWidth(fWid);
417  hSub[iComp] -> SetTitle(sTitle.Data());
418  hSub[iComp] -> SetTitleFont(fTxt);
419  hSub[iComp] -> GetXaxis() -> SetRangeUser(xPlotRange[0], xPlotRange[1]);
420  hSub[iComp] -> GetXaxis() -> SetTitle(sTitleX.Data());
421  hSub[iComp] -> GetXaxis() -> SetTitleFont(fTxt);
422  hSub[iComp] -> GetXaxis() -> SetTitleSize(fTit[0]);
423  hSub[iComp] -> GetXaxis() -> SetTitleOffset(fOffX[0]);
424  hSub[iComp] -> GetXaxis() -> SetLabelFont(fTxt);
425  hSub[iComp] -> GetXaxis() -> SetLabelSize(fLab[0]);
426  hSub[iComp] -> GetXaxis() -> CenterTitle(fCnt);
427  hSub[iComp] -> GetYaxis() -> SetTitle(sTitleS.Data());
428  hSub[iComp] -> GetYaxis() -> SetTitleFont(fTxt);
429  hSub[iComp] -> GetYaxis() -> SetTitleSize(fTit[0]);
430  hSub[iComp] -> GetYaxis() -> SetTitleOffset(fOffY[0]);
431  hSub[iComp] -> GetYaxis() -> SetLabelFont(fTxt);
432  hSub[iComp] -> GetYaxis() -> SetLabelSize(fLab[0]);
433  hSub[iComp] -> GetYaxis() -> CenterTitle(fCnt);
434  }
435  cout << " Set styles." << endl;
436 
437  // make legend
438  const UInt_t fColLe = 0;
439  const UInt_t fFilLe = 0;
440  const UInt_t fLinLe = 0;
441  const Float_t yObjLe = 0.1 + ((NComp + 1) * 0.05);
442  const Float_t fLegXY[NVtx] = {0.1, 0.1, 0.3, yObjLe};
443 
444  TLegend *leg = new TLegend(fLegXY[0], fLegXY[1], fLegXY[2], fLegXY[3], sHeader.Data());
445  leg -> SetFillColor(fColLe);
446  leg -> SetFillStyle(fFilLe);
447  leg -> SetLineColor(fColLe);
448  leg -> SetLineStyle(fLinLe);
449  leg -> SetTextFont(fTxt);
450  leg -> SetTextAlign(fAln);
451  leg -> AddEntry(hRef, sLabelRef.Data(), "pf");
452  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
453  leg -> AddEntry(hComp[iComp], sLabelComp[iComp], "pf");
454  }
455  cout << " Made legend." << endl;
456 
457  // make text
458  const UInt_t fColTx = 0;
459  const UInt_t fFilTx = 0;
460  const UInt_t fLinTx = 0;
461  const Float_t yObjTx = 0.1 + (NTxt * 0.05);
462  const Float_t fTxtXY[NVtx] = {0.3, 0.1, 0.5, yObjTx};
463 
464  TPaveText *txt = new TPaveText(fTxtXY[0], fTxtXY[1], fTxtXY[2], fTxtXY[3], "NDC NB");
465  txt -> SetFillColor(fColTx);
466  txt -> SetFillStyle(fFilTx);
467  txt -> SetLineColor(fColTx);
468  txt -> SetLineStyle(fLinTx);
469  txt -> SetTextFont(fTxt);
470  txt -> SetTextAlign(fAln);
471  for (UInt_t iTxt = 0; iTxt < NTxt; iTxt++) {
472  txt -> AddText(sTxt[iTxt].Data());
473  }
474  cout << " Made text." << endl;
475 
476  // make line
477  const UInt_t fColLi(1);
478  const UInt_t fLinLi(9);
479  const UInt_t fWidLi(1);
480  const Float_t fLinXY[NVtx] = {xPlotRange[0], 0., xPlotRange[1], 0.};
481 
482  TLine *line = new TLine(fLinXY[0], fLinXY[1], fLinXY[2], fLinXY[3]);
483  line -> SetLineColor(fColLi);
484  line -> SetLineStyle(fLinLi);
485  line -> SetLineWidth(fWidLi);
486  cout << " Made line." << endl;
487 
488  // make plot
489  const UInt_t width(750);
490  const UInt_t height(950);
491  const UInt_t fMode(0);
492  const UInt_t fBord(2);
493  const UInt_t fGrid(0);
494  const UInt_t fTick(1);
495  const UInt_t fLogX(1);
496  const UInt_t fLogY1(0);
497  const UInt_t fLogY2(1);
498  const UInt_t fFrame(0);
499  const Float_t fMarginL(0.15);
500  const Float_t fMarginR(0.02);
501  const Float_t fMarginT1(0.005);
502  const Float_t fMarginT2(0.02);
503  const Float_t fMarginTNR(0.02);
504  const Float_t fMarginB1(0.25);
505  const Float_t fMarginB2(0.005);
506  const Float_t fMarginBNR(0.15);
507  const Float_t fPadXY1[NVtx] = {0., 0., 1., 0.35};
508  const Float_t fPadXY2[NVtx] = {0., 0.35, 1., 1.};
509 
510  TCanvas *cPlot = new TCanvas("cPlot", "", width, height);
511  TPad *pPad1 = new TPad("pPad1", "", fPadXY1[0], fPadXY1[1], fPadXY1[2], fPadXY1[3]);
512  TPad *pPad2 = new TPad("pPad2", "", fPadXY2[0], fPadXY2[1], fPadXY2[2], fPadXY2[3]);
513  cPlot -> SetGrid(fGrid, fGrid);
514  cPlot -> SetTicks(fTick, fTick);
515  cPlot -> SetBorderMode(fMode);
516  cPlot -> SetBorderSize(fBord);
517  pPad1 -> SetGrid(fGrid, fGrid);
518  pPad1 -> SetTicks(fTick, fTick);
519  pPad1 -> SetLogx(fLogX);
520  pPad1 -> SetLogy(fLogY1);
521  pPad1 -> SetBorderMode(fMode);
522  pPad1 -> SetBorderSize(fBord);
523  pPad1 -> SetFrameBorderMode(fFrame);
524  pPad1 -> SetLeftMargin(fMarginL);
525  pPad1 -> SetRightMargin(fMarginR);
526  pPad1 -> SetTopMargin(fMarginT1);
527  pPad1 -> SetBottomMargin(fMarginB1);
528  pPad2 -> SetGrid(fGrid, fGrid);
529  pPad2 -> SetTicks(fTick, fTick);
530  pPad2 -> SetLogx(fLogX);
531  pPad2 -> SetLogy(fLogY2);
532  pPad2 -> SetBorderMode(fMode);
533  pPad2 -> SetBorderSize(fBord);
534  pPad2 -> SetFrameBorderMode(fFrame);
535  pPad2 -> SetLeftMargin(fMarginL);
536  pPad2 -> SetRightMargin(fMarginR);
537  pPad2 -> SetTopMargin(fMarginT2);
538  pPad2 -> SetBottomMargin(fMarginB2);
539  cPlot -> cd();
540  pPad1 -> Draw();
541  pPad2 -> Draw();
542  pPad1 -> cd();
543  hSub[0] -> Draw(sOptSub[0].Data());
544  for (UInt_t iComp = 1; iComp < NComp; iComp++) {
545  hSub[iComp] -> Draw(sOptSub[iComp].Data());
546  }
547  line -> Draw();
548  pPad2 -> cd();
549  hRef -> Draw(sOptRef.Data());
550  for(UInt_t iComp = 0; iComp < NComp; iComp++) {
551  hComp[iComp] -> Draw(sOptComp[iComp].Data());
552  }
553  leg -> Draw();
554  txt -> Draw();
555  fOut -> cd();
556  cPlot -> Write();
557  cPlot -> Close();
558  cout << " Made plot." << endl;
559 
560  // save histograms
561  fOut -> cd();
562  hRef -> Write();
563  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
564  hComp[iComp] -> Write();
565  hSub[iComp] -> Write();
566  }
567  cout << " Saved histograms." << endl;
568 
569  // close files
570  fOut -> cd();
571  fOut -> Close();
572  fRef -> cd();
573  fRef -> Close();
574  for (UInt_t iComp = 0; iComp < NComp; iComp++) {
575  fComp[iComp] -> cd();
576  fComp[iComp] -> Close();
577  }
578  cout << " Finished plot!\n" << endl;
579 
580 }
581 
582 // end ------------------------------------------------------------------------