Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
draw_hists_3D.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file draw_hists_3D.C
1 #include "sPhenixStyle.C"
2 #include "sPhenixStyle.h"
3 
4 TPad* getPad(int j, int k, int l, float cw, float ch, const int nEBins, const int nCutBins);
5 
6 //void drawCanvas_invMass(TCanvas *c, TH1F *hInvMass, int pad_x, int pad_y, TPad *pad, int isEta);
7 void drawCanvas_invMass(TCanvas *c, TH1F *hCorr, int pad_x, int pad_y, TPad *pad, int isEta, float peakPos, float peakW);
8 
9 void SetHistoStyle(TH1F *histo, int cutBin, int eBin, float low, float hi, int isEta);
10 
11 void DrawCutLabel(int x, int y, int isEta);
12 
13 void GetSubtractedDist(TH1F *histOrig, TH1F *histSub, TF1 *invMassFit, TF1 *invMassBG);
14 
15 //float eBins[] = {1,3,5,7,9,11,13,15,17,20};
16 float eBins[] = {1,2,3,4,5,6,7,8,9,10,11,12,13};
17 float eBinsEta[] = {12,14,16,18,20};
18 float eCuts[] = {0.5,0.6,0.7,0.8,0.9,1.,1.1/*,1.2,1.3,1.4,1.5*/};
19 float eCutsEta[] = {3,3.5,4,4.5,5,5.5,6};
20 const int nEtaBins = 5;
21 void draw_hists_3D(const char *input)
22 {
23 
24  SetsPhenixStyle();
25  TFile *fin = new TFile(input);
26 
27  const int nEBins = (int)sizeof(eBins)/sizeof(eBins[0]);
28  const int nCutBins = (int)sizeof(eCuts)/sizeof(eCuts[0]);
29  const int nCutBinsEta = (int)sizeof(eCutsEta)/sizeof(eCutsEta[0]);
30  const int nEBinsEta = (int)sizeof(eBinsEta)/sizeof(eBinsEta[0]);
31 
32  TH3F *hist[nEtaBins];
33  TH3F *histMatched[nEtaBins];
34  float cw = 2*966, ch = 1.3*637;
35 
36  TCanvas *cPi0Cut[nEtaBins];
37  TPad *pad_invMass[nEBins-1][nCutBins][nEtaBins];
38 
39  TCanvas *cSubMass[nEtaBins];
40  TPad *pad_invMassSub[nEBins-1][nCutBins][nEtaBins];
41 
42  TCanvas *cMatchMass[nEtaBins];
43  TPad *pad_matchMass[nEBins-1][nCutBins][nEtaBins];
44 
45 
46  TGraphErrors *mass[nCutBins][nEtaBins];
47  TGraphErrors *width[nCutBins][nEtaBins];
48  TGraphErrors *yields[nCutBins][nEtaBins];
49  TGraphErrors *yieldsMatched[nCutBins][nEtaBins];
50  TGraphErrors *matchZeroFrac[nCutBins][nEtaBins];
51  TF1 *invMassFit, *invMassFitBG;
52  int colors[] = {1,2,4,kGreen+2, kViolet,kCyan,kOrange+2,kMagenta+2,kAzure-2};
53  TLegend *leggy = new TLegend(.161, .729, .722, .939 );
54  leggy -> SetNColumns(3);
55  leggy -> SetFillStyle(0);
56  leggy -> SetBorderSize(0);
57 
58  float xEndFit = 0.25;
59  float xStartFit = 0.07;
60  float xEndHist = 0.25;
61  float xStartHist = 0.07;
62  float invMassWindow[2] = {0.112, 0.162};
63  float peakPos[nEBins - 1][nCutBins][nEtaBins];
64  float peakW[nEBins - 1][nCutBins][nEtaBins];
65  float sigma = 3.;
66  TLatex marks;
67 
68  for(int k = 0; k < nEtaBins; k++)
69  {
70 
71 
72  cSubMass[k]= new TCanvas(Form("cMassSubPi0_Eta%d",k),"Subtracted Mass Dist",cw,ch);
73  cSubMass[k] -> Range(0,0,1,1);
74 
75  cPi0Cut[k] = new TCanvas(Form("cMassPi0_Eta%d",k),"Invariant Mass",cw,ch);
76  cPi0Cut[k] -> Range(0,0,1,1);
77 
78  hist[k] = (TH3F*)fin -> Get(Form("ePi0InvMassEcut_Eta%d",k));
79 
80  for(int i = 0; i < nEBins - 1; i++)
81  {
82  for(int j = 0; j < nCutBins; j++)
83  {
84  if(i == 0)
85  {
86  mass[j][k] = new TGraphErrors();
87  mass[j][k] -> SetMarkerColor(colors[j]);
88 
89  width[j][k] = new TGraphErrors();
90  width[j][k] -> SetMarkerColor(colors[j]);
91 
92  yields[j][k] = new TGraphErrors();
93  yields[j][k] -> SetMarkerColor(colors[j]);
94 
95 
96  }
97 
98  pad_invMass[i][j][k] = getPad(i,j,k,cw,ch,nEBins , nCutBins);//for unsubtracted invariant mass distributions
99  pad_invMass[i][j][k] -> cd();
100 
101  pad_invMassSub[i][j][k] = getPad(i,j,k,cw,ch,nEBins , nCutBins);//for subtracted mass distributions
102 
103 
104  TH1F *hist1D = (TH1F*)hist[k] -> ProjectionY(Form("invMass_E%d_Cut%d_Eta%d",i,j,k), hist[k] -> GetXaxis() -> FindBin(eBins[i] + 0.001), hist[k] -> GetXaxis() -> FindBin(eBins[i+1] - 0.0001), hist[k] -> GetZaxis() -> FindBin(eCuts[j] + 0.0001), hist[k] -> GetNbinsZ());
105 
106 
107 
108  TH1F *invMassSub = (TH1F*)hist1D -> Clone();
109  if(eBins[i+1] <= 13)
110  {
111  invMassFit = new TF1(Form("invMassFit_E%d_Cut%d_Eta%d",i,j,k),"[5]*exp(-(([7]-x)^2)/(2*[6]^2))");
112 
113 
114  invMassFit -> SetLineWidth(1);
115  invMassFit -> SetParameter(6,0.015);
116  invMassFit -> SetParLimits(6,0.005,0.1);
117  invMassFit -> SetLineColor(2);
118  invMassFit -> SetParameter(7,0.135);
119  invMassFit -> SetParLimits(7, 0.135 - 0.07, 0.135 + 0.07);
120  if((eCuts[j] >= 0.9) && (eBins[i] == 3 || eBins[i] == 4))
121  {
122  invMassFit -> SetParameter(5, hist1D -> GetMaximum()*2);
123 
124  }
125  else invMassFit -> SetParameter(5,hist1D -> GetMaximum());
126  hist1D -> Fit(invMassFit,"Q","",xStartFit, xEndFit);
127  // GetSubtractedDist(hist1D, invMassSub, invMassFit, invMassFitBG);
128 
129  float peakPosition = invMassFit -> GetParameter(7);
130  peakPos[i][j][k] = peakPosition;
131 
132  float peakWidth = invMassFit -> GetParameter(6);
133  peakW[i][j][k] = peakWidth;
134  Double_t yieldErr;
135 
136  float yield;
137  if(hist1D -> Integral() == 0) yield = 0;
138  else yield = invMassSub -> IntegralAndError(invMassSub -> GetXaxis() -> FindBin((peakPosition - sigma*peakWidth)+0.0001), invMassSub -> GetXaxis() -> FindBin((peakPosition + sigma*peakWidth)+0.0001), yieldErr, "");
139 
140  if(yield != 0)
141  {
142  mass[j][k] -> SetPoint(i, (eBins[i]+eBins[i+1])/2, invMassFit -> GetParameter(7));
143  mass[j][k] -> SetPointError(i, 0, invMassFit -> GetParError(7));
144 
145  width[j][k] -> SetPoint(i, (eBins[i]+eBins[i+1])/2, invMassFit -> GetParameter(6));
146  width[j][k] -> SetPointError(i, 0, invMassFit -> GetParError(6));
147 
148  yields[j][k] -> SetPoint(i, (eBins[i]+eBins[i+1])/2, yield);
149  yields[j][k] -> SetPointError(i, 0, yieldErr);
150 
151 
152  }
153 
154  }
155  hist1D -> GetXaxis() -> SetRangeUser(xStartHist,xEndHist);
156 
157  TCanvas *cMassFitCheck = new TCanvas();
158  hist1D -> SetTitle(";M_{#gamma#gamma} [GeV/c^{2}];");
159  hist1D -> Draw("epx0");
160 
161 
162  marks.DrawLatexNDC(0.7,0.8,Form("#splitline{E_{Min} > %g GeV}{%g > E_{#pi^{0}} > %g}",eCuts[j], eBins[i], eBins[i+1]));
163  cMassFitCheck -> SaveAs(Form("plots/invMassFitCheck_%dE_%dCut_Sigma%g_Eta%d.pdf",i,j,sigma,k));
164 
165  drawCanvas_invMass(cPi0Cut[k], hist1D, j, i, pad_invMass[i][j][k], 0, peakPos[i][j][k], sigma*peakW[i][j][k]);
166  pad_invMassSub[i][j][k] -> cd();
167 
168  invMassSub -> GetXaxis() -> SetRangeUser(xStartHist, xEndHist);
169  drawCanvas_invMass(cSubMass[k], invMassSub, j, i, pad_invMassSub[i][j][k], 0, peakPos[i][j][k], sigma*peakW[i][j][k]);
170 
171 
172  }
173  }
174 
175  cPi0Cut[k] -> SaveAs(Form("plots/invMassCuts_Sigma%g_Eta%d.pdf",sigma, k));
176  cSubMass[k] -> SaveAs(Form("plots/invMassCutsSubtracted_Sigma%g_Eta%d.pdf",sigma, k));
177  }
178 
179 
180  TCanvas *cMass[nEtaBins];
181  for(int k = 0; k < nEtaBins; k++)
182  {
183  cMass[k] = new TCanvas();
184  for(int i = 0; i < nCutBins; i++)
185  {
186  if(i == 0)
187  {
188  mass[i][k] -> Draw("ap");
189  mass[i][k] -> GetYaxis() -> SetRangeUser(0.130,0.180);
190  mass[i][k] -> GetXaxis() -> SetLimits(0,13);
191  mass[i][k] -> GetXaxis() -> SetTitle("E_{#pi^{0}}^{Reco} [GeV]");
192  mass[i][k] -> GetYaxis() -> SetTitle("Mass Peak Location [GeV/c^{2}]");
193  }
194  else mass[i][k] -> Draw("samep");
195  if(k==0)leggy -> AddEntry(mass[i][k], Form("E_{min} > %g",eCuts[i]),"p");
196  }
197 
198  leggy -> Draw("same");
199  cMass[k] -> SaveAs(Form("plots/invMassGraphPeakPos_Sigma%g_Eta%d.pdf",sigma,k));
200  }
201 
202  TCanvas *cWidth[nEtaBins];
203  for(int k = 0; k < nEtaBins; k++)
204  {
205  cWidth[k] = new TCanvas();
206 
207  for(int i = 0; i < nCutBins; i++)
208  {
209  if(i == 0)
210  {
211  width[i][k] -> Draw("ap");
212  width[i][k] -> GetYaxis() -> SetRangeUser(0.013, 0.028);
213  width[i][k] -> GetXaxis() -> SetLimits(0,13);
214  width[i][k] -> GetXaxis() -> SetTitle("E_{#pi^{0}}^{Reco} [GeV]");
215  width[i][k] -> GetYaxis() -> SetTitle("Mass Peak Width [GeV/c^{2}]");
216  }
217  else width[i][k] -> Draw("samep");
218  }
219 
220  leggy -> Draw("same");
221  cWidth[k] -> SaveAs(Form("plots/invMassGraphPeakWidth_Sigma%g_Eta%d.pdf",sigma,k));
222  }
223 
224  TCanvas *cYields[nEtaBins];
225  for(int k = 0; k < nEtaBins; k++)
226  {
227  cYields[k] = new TCanvas();
228  for(int i = 0; i < nCutBins; i++)
229  {
230  if(i == 0)
231  {
232  yields[i][k] -> Draw("ap");
233  yields[i][k] -> GetYaxis() -> SetRangeUser(1e3,1e7);
234  yields[i][k] -> GetXaxis() -> SetTitle("E_{#pi^{0}}^{Reco} [GeV]");
235  yields[i][k] -> GetYaxis() -> SetTitle("N_{#pi^{0}}");
236  }
237  else yields[i][k] -> Draw("samep");
238  gPad -> SetLogy();
239  }
240 
241  leggy -> Draw("same");
242  cYields[k] -> SaveAs(Form("plots/invMassYield_Sigma%g_Eta%d.pdf",sigma,k));
243  }
244 
245 
246  TH1F *pi0TruSpec[nEtaBins];
247  for(int i = 0; i < nEtaBins; i++) pi0TruSpec[i] = (TH1F*)fin -> Get(Form("primPi0E_Eta%d",i));
248  TGraphErrors *pi0Eff[nCutBins][nEtaBins];
249  TCanvas *cEff[nEtaBins];
250  float scaleMax;
251  float scaleMin;
252 
253  for(int k = 0; k < nEtaBins; k++)
254  {
255  cEff[k] = new TCanvas();
256  scaleMax = -9999;
257  scaleMin = 9999;
258  cEff[k] -> cd();
259  for(int i = 0; i < nCutBins ; i++)
260  {
261  pi0Eff[i][k] = new TGraphErrors();
262  for(int j = 0; j < nEBins -1; j++)
263  {
264  Double_t pi0TruthYieldErr;
265  float pi0TruthYield = pi0TruSpec[k] -> IntegralAndError(pi0TruSpec[k] -> GetXaxis() -> FindBin(eBins[j]+0.001), pi0TruSpec[k] -> GetXaxis() -> FindBin(eBins[j+1]-0.001), pi0TruthYieldErr, "");
266 
267 
268  Double_t x, pi0YieldMeas, pi0YieldErr;
269  yields[i][k] -> GetPoint(j, x, pi0YieldMeas);
270  pi0YieldErr = yields[i][k] -> GetErrorY(j);
271  float eff = pi0YieldMeas/pi0TruthYield;
272  float err = eff*sqrt(pow(pi0YieldErr/pi0YieldMeas,2) + pow(pi0TruthYieldErr/pi0TruthYield,2));
273 
274  if(eff < scaleMin) scaleMin = eff;
275  if(eff > scaleMax) scaleMax = eff;
276  pi0Eff[i][k] -> SetPoint(j, (eBins[j]+eBins[j+1])/2, eff);
277  pi0Eff[i][k] -> SetPointError(j, 0, err);
278  }
279  pi0Eff[i][k] -> SetMarkerColor(colors[i]);
280  pi0Eff[i][k] -> SetName(Form("pi0Eff_Eta%d_ECut%d",k,i));
281 
282  if(i == 0)
283  {
284  pi0Eff[i][k] -> Draw("ap");
285  pi0Eff[i][k] -> GetYaxis() -> SetRangeUser(-0.05,scaleMax*1.4);
286  pi0Eff[i][k] -> GetXaxis() -> SetTitle("E_{#pi^{0}}^{Reco} [GeV]");
287  pi0Eff[i][k] -> GetYaxis() -> SetTitle("#epsilon");
288  }
289  else pi0Eff[i][k] -> Draw("samep");
290 
291  }
292  leggy -> Draw("same");
293  cEff[k] -> SaveAs(Form("plots/pi0Eff_Sigma%g_Eta%d.pdf",sigma,k));
294  }
295 
296 
297  TCanvas *cAngle2D = new TCanvas();
298  TH2F *deltaR2D = (TH2F*)fin -> Get("truthPi0EDeltaR");
299 
300  deltaR2D -> SetTitle(";E_{#pi^{0}}^{Truth} [GeV];#Delta R");
301  deltaR2D -> GetYaxis() -> SetRangeUser(0,0.5);
302  deltaR2D -> Draw("colz");
303  cAngle2D -> SaveAs("plots/deltaR.pdf");
304 
305  TH1F *deltaRproj[nEBins-1];
306  TCanvas *cAngleProj = new TCanvas();
307 
308  for(int i = 0; i < nEBins - 1; i++)
309  {
310  deltaRproj[i] = (TH1F*)deltaR2D -> ProjectionY(Form("deltaRProj_%d",i),deltaR2D -> GetXaxis() -> FindBin(eBins[i]+0.001), deltaR2D -> GetXaxis() -> FindBin(eBins[i+1]-0.001));
311  deltaRproj[i] -> Draw();
312  cAngleProj -> SaveAs(Form("plots/deltaR_proj%d.pdf",i));
313  }
314 
315  TCanvas *cAsym = new TCanvas();
316  TH2F *asym = (TH2F*)fin -> Get("truthPi0EAsym");
317 
318  asym -> SetTitle(";E_{#pi^{0}}^{Truth [GeV];#frac{|E^{#gamma_{1}} - E^{#gamma_{2}}|}{#lower[0.2]{E^{#gamma_{1}} + E^{#gamma_{2}}}}");
319  asym -> GetYaxis() -> SetTitleOffset(1.3);
320  asym -> Draw("colz");
321  cAsym -> SaveAs("plots/asymE.pdf");
322 
323 }
324 
325 TPad* getPad(int j, int k, int l, float cw, float ch, const int nEBins, const int nCutBins)
326 {
327  float xpos[nCutBins+1], ypos[nEBins];
328  int nCutOffset = 1;
329  float marginLeft = 0.07;
330  float marginRight = 0.035;
331 
332  float marginTop = 0.03;
333  float marginBottom = 0.06;
334 
335  float marginColumn = 0.04;
336 
337  float w = (1 - marginLeft - marginRight - marginColumn*(nCutBins-2 + nCutOffset))/(nCutBins - 1+ nCutOffset);
338  float h = (1 - marginTop - marginBottom)/(nEBins - 1);
339 
340  for(int m = 0; m < nCutBins + nCutOffset; m++)
341  {
342  if(m == 0) xpos[m] = 0;
343  else if(m == nCutBins - 1+ nCutOffset) xpos[m] = 1.;
344  else xpos[m] = marginLeft + m*w + (m-1)*marginColumn;
345  }
346 
347  for(int m = 0; m < nEBins; m++)
348  {
349  if(m==0) ypos[nEBins -1 -m] = 0;
350  else if (m == nEBins - 1) ypos[nEBins-1-m] = 1.;
351  else ypos[nEBins - 1 -m] = marginBottom +m*h;
352  }
353 
354  TPad *myPad = new TPad(Form("pad%d%d%d",j,k,l),Form("pad%d%d%d",j,k,l),xpos[k],ypos[j+1],xpos[k+1],ypos[j]);
355 
356  if(k==0) myPad->SetLeftMargin(marginLeft/(xpos[k+1] - xpos[k]));
357  else myPad -> SetLeftMargin(marginColumn/(xpos[k+1] - xpos[k]));
358 
359  if(k == nCutBins - 2+ nCutOffset) myPad -> SetRightMargin(marginRight/(xpos[k+1]-xpos[k]));
360  else myPad -> SetRightMargin(0.012);
361 
362  if(j == 0) myPad -> SetTopMargin(marginTop/(ypos[j] - ypos[j+1]));
363  else myPad -> SetTopMargin(0);
364 
365  if(j == nEBins - 2) myPad -> SetBottomMargin(marginBottom/(ypos[j] - ypos[j+1]));
366  else myPad -> SetBottomMargin(0);
367 
368  return myPad;
369 }
370 
371 void drawCanvas_invMass(TCanvas *c, TH1F *hCorr, int pad_x, int pad_y, TPad *pad, int mode, float peakPos, float peakW)
372 {
373 
374  float range = hCorr -> GetMaximum() - hCorr -> GetMinimum();
375  SetHistoStyle(hCorr, pad_x, pad_y, hCorr -> GetMinimum(), hCorr->GetMaximum() + 0.2*range, mode);
376  c -> cd(0);
377  /*TLine *pi0Mass = new TLine(0.135,0,0.135,hCorr->GetMaximum());
378  pi0Mass -> SetLineStyle(7);
379  pi0Mass -> SetLineColor(2);
380  TLine *etaMass = new TLine(0.548,0,0.548,hCorr->GetMaximum());
381  etaMass -> SetLineStyle(7);
382  etaMass -> SetLineColor(2);*/
383  TLine *massWindow[2];
384  //massWindow[0] = new TLine(0.112,0,0.112,hCorr->GetMaximum());
385  pad -> Draw();
386  pad -> cd();
387 
388  hCorr -> Draw("epX0");
389  DrawCutLabel(pad_x,pad_y,mode);
390  if(mode == 1) gPad -> SetLogy();
391  if(peakPos)
392  {
393  massWindow[0] = new TLine(peakPos - peakW,0,peakPos - peakW,hCorr->GetMaximum());
394  massWindow[0] -> SetLineStyle(7);
395  massWindow[0] -> SetLineColor(2);
396 
397  //massWindow[1] = new TLine(0.162,0,0.162,hCorr->GetMaximum());
398  massWindow[1] = new TLine(peakPos + peakW,0,peakPos + peakW,hCorr->GetMaximum());
399  massWindow[1] -> SetLineStyle(7);
400  massWindow[1] -> SetLineColor(2);
401  massWindow[0] -> Draw("same");
402  massWindow[1] -> Draw("same");
403  }
404 
405 
406  //etaMass -> Draw("same");
407  //pi0Mass -> Draw("same");
408 
409 }
410 
411 void SetHistoStyle(TH1F *histo, int cutBin, int eBin, float low, float hi, int mode)
412 {
413  if(hi > low) histo -> GetYaxis() -> SetRangeUser(low,hi);
414  if(mode == 1)
415  {
416  histo -> GetYaxis() -> SetRangeUser(0.1,10*hi);
417  }
418 
419  histo -> SetTitle(";;");
420  int xTitleBin[3][2] = {{11,3},{11,3},{3,3}};
421  int yTitleBin[3][2] = {{5,0},{5,0},{2,0}};
422 
423  float yTitleOffset[3] = {0.3, 0.4, 1};
424  float xTitleOffset[3] = {0.7, 0.7, 0.7};
425 
426  float yTitleSize[3] = {0.5, 0.5, 0.15};
427  float xTitleSize[3] = {0.28, 0.28, 0.13};
428 
429  if(eBin == xTitleBin[mode][0] && cutBin == xTitleBin[mode][1]) histo -> GetXaxis() -> SetTitle("M_{#gamma#gamma}");
430  if(eBin == yTitleBin[mode][0] && cutBin == yTitleBin[mode][1]) histo -> GetYaxis() -> SetTitle("N_{#gamma#gamma}");
431 
432  histo -> GetXaxis() -> SetTitleSize(xTitleSize[mode]);
433  histo -> GetXaxis() -> SetTitleOffset(xTitleOffset[mode]);
434  histo -> GetXaxis() -> CenterTitle();
435 
436  histo -> GetXaxis() -> SetLabelFont(43);
437  histo -> GetXaxis() -> SetLabelSize(18);
438  histo -> GetXaxis() -> SetLabelOffset(0.0);
439  histo -> SetTickLength(0.05,"X");
440  histo -> SetNdivisions(505,"x");
441 
442 
443  histo -> GetYaxis() -> SetNdivisions(505);
444  histo -> GetYaxis() -> SetTitleSize(yTitleSize[mode]);
445  histo -> GetYaxis() -> SetTitleOffset(yTitleOffset[mode]);
446  histo -> GetYaxis() -> CenterTitle();
447 
448  histo -> GetYaxis() -> SetLabelFont(43);
449  histo -> GetYaxis() -> SetLabelSize(19);
450  histo -> GetYaxis() -> SetLabelOffset(0.003);
451  histo -> GetYaxis() -> SetNoExponent();
452  histo -> SetTickLength(0.04,"Y");
453 
454  histo -> SetStats(0);
455  //histo -> SetLineWidth(1.4);
456  histo -> SetLineColor(1);
457  histo -> SetMarkerColor(1);
458  histo -> SetMarkerStyle(20);
459  histo -> SetMarkerSize(0.4);
460 }
461 
462 void DrawCutLabel(int x, int y, int mode)
463 {
464  gPad -> cd();
465  gPad -> Update();
466  float Y2 = gPad -> GetUymax();
467  float Y1 = gPad -> GetUymin();
468  float rangeY = Y2 - Y1;
469  float xtop[3] = {0.8, -0.1, 0};
470  float xbottom[3] = {0.273152, 1.01, 1.1};
471  double ytop[3] = {Y2+0.01*rangeY, Y2+1*rangeY, Y2+0.01*rangeY};
472  double ybottom[3] = { 0.98*Y2, 0.98*Y2, 0.9*Y2};
473 
474  TText *textCut;
475  if(y == 0)
476  {
477  // if(mode < 2)
478  // {
479  // if(x <= 4)textCut = new TText(0.08, Y2+0.01*rangeY, Form("Min E > %g GeV", eCuts[x]));
480  // else textCut = new TText(0.08, Y2+0.01*rangeY, Form("Min E > %g GeV", eCuts[x]));
481  // }
482  //else textCut = new TText(0, Y2+0.01*rangeY, Form("Min E > %g GeV", eCutsEta[x]));
483  textCut = new TText(xtop[mode], ytop[mode], Form("Min E > %g GeV", eCuts[x]));
484  textCut -> SetTextFont(43);
485  textCut -> SetTextSize(18);
486  textCut -> Draw();
487  }
488  if(x == sizeof(eCuts)/sizeof(eCuts[0]) - 1)
489  {
490  //if(mode < 2)textCut = new TText(0.273152, 0.98*Y2, Form("%g-%gGeV",eBins[y], eBins[y+1]));
491  //else textCut = new TText(1.1, 0.9*Y2, Form("%g-%gGeV",eBinsEta[y], eBinsEta[y+1]));
492  if(mode < 2)textCut = new TText(xbottom[mode], ybottom[mode], Form("%g-%gGeV",eBins[y], eBins[y+1]));
493  else textCut = new TText(xbottom[mode], ybottom[mode], Form("%g-%gGeV",eBinsEta[y], eBinsEta[y+1]));
494  textCut -> SetTextAngle(270);
495  textCut -> SetTextFont(43);
496  textCut -> SetTextSize(18);
497  textCut -> Draw();
498  }
499 
500 }
501 
502 void GetSubtractedDist(TH1F *histOrig, TH1F *histSub, TF1 *invMassFit, TF1 *invMassBG)
503 {
504 
505 
506  for(int i = 0; i < invMassFit -> GetNpar(); i++)
507  {
508  invMassBG -> SetParameter(i, invMassFit -> GetParameter(i));
509  }
510 
511  for(int i = 1; i < histSub -> GetNbinsX() + 1; i++)
512  {
513  histSub -> SetBinContent(i, histOrig -> GetBinContent(i) - invMassBG -> Eval(histOrig -> GetBinCenter(i)));
514  }
515 
516 
517 }