EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
QA_Draw_Utility.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file QA_Draw_Utility.C
1 // $Id: $
2 
11 #ifndef QA_Draw_Utility_C
12 #define QA_Draw_Utility_C
13 
14 #include <TCanvas.h>
15 #include <TClass.h>
16 #include <TDirectory.h>
17 #include <TFile.h>
18 #include <TSystem.h>
19 
20 #include <TGraphErrors.h>
21 #include <TH1F.h>
22 #include <TList.h>
23 #include <TObject.h>
24 #include <TString.h>
25 #include <TStyle.h>
26 #include <TTree.h>
27 
28 #include <TLegend.h>
29 #include <TLegendEntry.h>
30 #include <TLine.h>
31 #include <TPad.h>
32 
33 #include <cassert>
34 #include <cmath>
35 
36 #include <iostream>
37 
38 using namespace std;
39 
41 void SavePad(TPad *p)
42 {
43  if (!p)
44  return;
45 
46  TList *l = p->GetListOfPrimitives();
47  // l->Print();
48 
49  TIter next(l);
50  TObject *obj = NULL;
51  while ((obj = next()))
52  {
53  if (obj->IsA()->GetBaseClassOffset(TClass::GetClass("TPad")) >= 0)
54  {
55  if ((TPad *) obj != p)
56  SavePad((TPad *) obj);
57  }
58  else if (obj->IsA()->GetBaseClassOffset(TClass::GetClass("TH1")) >= 0)
59  {
60  cout << "Save TH1 " << obj->GetName() << endl;
61  obj->Clone()->Write(obj->GetName(), TObject::kOverwrite);
62  }
63  else if (obj->IsA()->GetBaseClassOffset(TClass::GetClass("TF1")) >= 0)
64  {
65  cout << "Save TF1 " << obj->GetName() << endl;
66  obj->Clone()->Write(obj->GetName(), TObject::kOverwrite);
67  }
68  else if (obj->IsA()->GetBaseClassOffset(TClass::GetClass("TGraph")) >= 0)
69  {
70  cout << "Save TGraph " << obj->GetName() << endl;
71  obj->Clone()->Write(obj->GetName(), TObject::kOverwrite);
72  }
73  }
74 }
75 
77 
82 void SaveCanvas(TCanvas *c, TString name = "", Bool_t bEPS = kTRUE)
83 {
84  if (name.Length() == 0)
85  name = c->GetName();
86 
87  c->Print(name + ".png");
88 
89  TDirectory *oldd = gDirectory;
90 
91  TString rootfilename;
92 
93  c->Print(rootfilename = name + ".root");
94 
95  TFile f(rootfilename, "update");
96 
97  SavePad(c);
98 
99  f.Close();
100 
101  oldd->cd();
102 
103  if (bEPS)
104  {
105  // c->Print(name + ".pdf");
106 
107  float x = 20;
108  float y = 20;
109  gStyle->GetPaperSize(x, y);
110 
111  gStyle->SetPaperSize(c->GetWindowWidth() / 72 * 2.54,
112  c->GetWindowHeight() / 72 * 2.54);
113  // c->Print(name + ".eps");
114  c->Print(name + ".svg");
115  gSystem->Exec("rsvg-convert -f pdf -o " + name + ".pdf " + name + ".svg");
116  gSystem->Exec("rm -fv " + name + ".svg");
117 
118  gStyle->SetPaperSize(x, y);
119  }
120  // c->Print(name+".C");
121 }
122 
125 double DrawReference(TH1 *hnew, TH1 *href, bool draw_href_error = false)
126 {
127  hnew->SetLineColor(kBlue + 3);
128  hnew->SetMarkerColor(kBlue + 3);
129  // hnew->SetLineWidth(2);
130  hnew->SetMarkerStyle(kFullCircle);
131  // hnew->SetMarkerSize(1);
132 
133  if (href)
134  {
135  if (draw_href_error)
136  {
137  href->SetLineColor(kGreen + 1);
138  href->SetFillColor(kGreen + 1);
139  href->SetLineStyle(kSolid);
140  href->SetMarkerColor(kGreen + 1);
141  // href->SetLineWidth(2);
142  href->SetMarkerStyle(kDot);
143  href->SetMarkerSize(0);
144  }
145  else
146  {
147  href->SetLineColor(kGreen + 1);
148  href->SetFillColor(kGreen + 1);
149  href->SetLineStyle(0);
150  href->SetMarkerColor(kGreen + 1);
151  href->SetLineWidth(0);
152  href->SetMarkerStyle(kDot);
153  href->SetMarkerSize(0);
154  }
155  }
156 
157  hnew->Draw(); // set scale
158 
159  double ks_test = numeric_limits<double>::signaling_NaN();
160 
161  if (href)
162  {
163  if (draw_href_error)
164  {
165  href->DrawClone("E2 same");
166  href->SetFillStyle(0);
167  href->SetLineWidth(8);
168  href->DrawClone("HIST same ][");
169  }
170  else
171  href->Draw("HIST same");
172  hnew->Draw("same"); // over lay data points
173 
174  ks_test = hnew->KolmogorovTest(href, "D");
175  }
176 
177  // ---------------------------------
178  // now, make summary header
179  // ---------------------------------
180  if (href)
181  {
182  gPad->SetTopMargin(.14);
183  TLegend *legend = new TLegend(0, .93, 0, 1, hnew->GetTitle(), "NB NDC");
184  legend->Draw();
185  legend = new TLegend(0, .86, .3, .93, NULL, "NB NDC");
186  legend->AddEntry(href, Form("Reference"), "f");
187  legend->Draw();
188  legend = new TLegend(0.3, .86, 1, .93, NULL, "NB NDC");
189  TLegendEntry *le = legend->AddEntry(hnew, Form("New: KS-Test P=%.3f", ks_test), "lpe");
190  if (ks_test>=1)
191  le->SetTextColor(kBlue+1);
192  else if (ks_test>=.2)
193  le->SetTextColor(kGreen+2);
194  else if (ks_test>=.05)
195  le->SetTextColor(kYellow+1);
196  else
197  le->SetTextColor(kRed+1);
198  legend->Draw();
199  }
200  else
201  {
202  gPad->SetTopMargin(.07);
203  TLegend *legend = new TLegend(0, .93, 0, 1, hnew->GetTitle(), "NB NDC");
204  legend->Draw();
205  }
206 
207  return ks_test;
208 }
209 
212 void DrawReference(TGraph *hnew, TGraph *href, bool draw_href_error = true)
213 {
214  hnew->SetLineColor(kBlue + 3);
215  hnew->SetMarkerColor(kBlue + 3);
216  hnew->SetLineWidth(2);
217  hnew->SetMarkerStyle(kFullCircle);
218  hnew->SetMarkerSize(1);
219 
220  if (href)
221  {
222  if (draw_href_error)
223  {
224  href->SetLineColor(kGreen + 1);
225  href->SetFillColor(kGreen + 1);
226  href->SetFillStyle(0);
227  href->SetLineStyle(kSolid);
228  href->SetMarkerColor(kGreen + 1);
229  href->SetLineWidth(4);
230  href->SetMarkerStyle(kDot);
231  href->SetMarkerSize(0);
232  }
233  else
234  {
235  href->SetLineColor(kGreen + 1);
236  href->SetFillColor(kGreen + 1);
237  href->SetLineStyle(0);
238  href->SetMarkerColor(kGreen + 1);
239  href->SetLineWidth(0);
240  href->SetMarkerStyle(kDot);
241  href->SetMarkerSize(0);
242  }
243  }
244 
245  if (href)
246  {
247  if (draw_href_error)
248  {
249  href->DrawClone("E2");
250  }
251  else
252  href->Draw("HIST same");
253  hnew->Draw("p e"); // over lay data points
254  }
255 }
256 
258 TGraphErrors *
259 FitResolution(const TH2F *h2, const bool normalize_mean = true)
260 {
261  TProfile *p2 = h2->ProfileX();
262 
263  int n = 0;
264  double x[1000];
265  double ex[1000];
266  double y[1000];
267  double ey[1000];
268 
269  for (int i = 1; i <= h2->GetNbinsX(); i++)
270  {
271  TH1D *h1 = h2->ProjectionY(Form("htmp_%d", rand()), i, i);
272 
273  if (h1->GetSum() < 10)
274  continue;
275 
276  TF1 fgaus("fgaus", "gaus", -p2->GetBinError(i) * 4,
277  p2->GetBinError(i) * 4);
278 
279  TF1 f2(Form("dgaus"), "gaus + [3]*exp(-0.5*((x-[1])/[4])**2) + [5]",
280  -p2->GetBinError(i) * 4, p2->GetBinError(i) * 4);
281 
282  fgaus.SetParameter(1, p2->GetBinContent(i));
283  fgaus.SetParameter(2, p2->GetBinError(i));
284 
285  h1->Fit(&fgaus, "MQ0");
286 
287  x[n] = p2->GetBinCenter(i);
288  ex[n] = (p2->GetBinCenter(2) - p2->GetBinCenter(1)) / 2;
289 
290  const double norm = normalize_mean ? fgaus.GetParameter(1) : 1;
291 
292  y[n] = fgaus.GetParameter(2) / norm;
293  ey[n] = fgaus.GetParError(2) / norm;
294 
295  n++;
296  delete h1;
297  }
298 
299  TGraphErrors *ge = new TGraphErrors(n, x, y, 0, ey);
300  ge->SetName(TString(h2->GetName()) + "_FitResolution");
301 
302  ge->SetLineColor(kBlue + 3);
303  ge->SetMarkerColor(kBlue + 3);
304  ge->SetLineWidth(2);
305  ge->SetMarkerStyle(kFullCircle);
306  ge->SetMarkerSize(1);
307  return ge;
308 }
309 
311 TGraphErrors *
312 FitProfile(const TH2F *h2)
313 {
314  TProfile *p2 = h2->ProfileX();
315 
316  int n = 0;
317  double x[1000];
318  double ex[1000];
319  double y[1000];
320  double ey[1000];
321 
322  for (int i = 1; i <= h2->GetNbinsX(); i++)
323  {
324  TH1D *h1 = h2->ProjectionY(Form("htmp_%d", rand()), i, i);
325 
326  if (h1->GetSum() < 10)
327  continue;
328 
329  TF1 fgaus("fgaus", "gaus", -p2->GetBinError(i) * 4,
330  p2->GetBinError(i) * 4);
331 
332  TF1 f2(Form("dgaus"), "gaus + [3]*exp(-0.5*((x-[1])/[4])**2) + [5]",
333  -p2->GetBinError(i) * 4, p2->GetBinError(i) * 4);
334 
335  fgaus.SetParameter(1, p2->GetBinContent(i));
336  fgaus.SetParameter(2, p2->GetBinError(i));
337 
338  h1->Fit(&fgaus, "MQ0");
339 
340  // f2.SetParameters(fgaus.GetParameter(0) / 2, fgaus.GetParameter(1),
341  // fgaus.GetParameter(2), fgaus.GetParameter(0) / 2,
342  // fgaus.GetParameter(2) / 4, 0);
343  //
344  // h1->Fit(&f2, "MQ0");
345 
346  // new TCanvas;
347  // h1->Draw();
348  // fgaus.Draw("same");
349  // break;
350 
351  x[n] = p2->GetBinCenter(i);
352  ex[n] = (p2->GetBinCenter(2) - p2->GetBinCenter(1)) / 2;
353  y[n] = fgaus.GetParameter(1);
354  ey[n] = fgaus.GetParameter(2);
355 
356  // p2->SetBinContent(i, fgaus.GetParameter(1));
357  // p2->SetBinError(i, fgaus.GetParameter(2));
358 
359  n++;
360  delete h1;
361  }
362 
363  TGraphErrors *ge = new TGraphErrors(n, x, y, ex, ey);
364  ge->SetName(TString(h2->GetName()) + "_FitProfile");
365  ge->SetLineColor(kBlue + 3);
366  ge->SetMarkerColor(kBlue + 3);
367  ge->SetLineWidth(2);
368  ge->SetMarkerStyle(kFullCircle);
369  ge->SetMarkerSize(1);
370  return ge;
371 }
372 
374 TH1 *GetBinominalRatio(TH1 *h_pass, TH1 *h_n_trial, bool process_zero_bins = false)
375 {
376  assert(h_pass);
377  assert(h_n_trial);
378 
379  assert(h_pass->GetNbinsX() == h_n_trial->GetNbinsX());
380  assert(h_pass->GetNbinsY() == h_n_trial->GetNbinsY());
381  assert(h_pass->GetNbinsZ() == h_n_trial->GetNbinsZ());
382 
383  TH1 *h_ratio = (TH1 *) h_pass->Clone(TString(h_pass->GetName()) + "_Ratio");
384  assert(h_ratio);
385  h_ratio->Divide(h_n_trial); // a rough estimation first, also taking care of the overflow bins and zero bins
386 
387  for (int x = 1; x <= h_n_trial->GetNbinsX(); ++x)
388  for (int y = 1; y <= h_n_trial->GetNbinsY(); ++y)
389  for (int z = 1; z <= h_n_trial->GetNbinsZ(); ++z)
390  {
391  const double n_trial = h_n_trial->GetBinContent(x, y, z);
392 
393  if (n_trial > 0)
394  {
395  const double p = h_pass->GetBinContent(x, y, z) / n_trial;
396 
397  // Wilson score interval
398  h_ratio->SetBinContent(x, y, z, //
399  (p + 1 / (2 * n_trial)) / (1 + 1 / n_trial));
400  h_ratio->SetBinError(x, y,
401  z, //
402  TMath::Sqrt(
403  1. / n_trial * p * (1 - p) + 1. / (4 * n_trial * n_trial)) /
404  (1 + 1 / n_trial));
405  }
406  else if (process_zero_bins)
407  {
408  h_ratio->SetBinContent(x, y, z, 0.5);
409  h_ratio->SetBinError(x, y, z, 0.5);
410  }
411  }
412 
413  return h_ratio;
414 }
415 
421 void PutInputFileName(TCanvas *c1, const double height, const char *new_file_name, const char *ref_file_name)
422 {
423  c1->cd();
424  TPad *pad = new TPad("PutInputFileName", "PutInputFileName", 0, 0, 1, height, 18);
425  pad->Draw();
426  pad->cd();
427 
428  if (new_file_name)
429  {
430  TText *t = new TText(0.05, 0.75, TString("New QA File: ") + TString(new_file_name));
431  t->SetTextAlign(12);
432  t->SetTextColor(kBlue + 3);
433  t->SetTextSize(.4);
434  t->Draw();
435  }
436  if (ref_file_name)
437  {
438  TText *t = new TText(0.05, 0.25, TString("Reference QA File: ") + TString(ref_file_name));
439  t->SetTextAlign(12);
440  t->SetTextColor(kGreen + 1);
441  t->SetTextSize(.4);
442  t->Draw();
443  }
444 }
445 
446 #endif