AliPhysics  aaf9c62 (aaf9c62)
Compare.C
Go to the documentation of this file.
1 #ifndef __CINT__
2 # include <TFile.h>
3 # include <THStack.h>
4 # include <TH1.h>
5 # include <TError.h>
6 # include <TMath.h>
7 # include <TClass.h>
8 # include <TCanvas.h>
9 # include <TSystem.h>
10 #else
11 class TFile;
12 class THStack;
13 class TH1;
14 class TCanvas;
15 #endif
16 
18 Int_t cW = 1200;
19 Int_t cH = 800;
20 const char* obs = "\\mathrm{d}N_{\\mathrm{ch}}/\\mathrm{d}\\eta";
21 
22 //____________________________________________________________________
29 void PrintH(TH1* h, Int_t prec=2)
30 {
31  Printf("Content of %s - %s", h->GetName(), h->GetTitle());
32  for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
33  if (h->GetBinContent(i) <= 1e-6) continue;
34  printf("%3d (%+6.3f): %.*f+/-%.*f\n", i,
35  h->GetXaxis()->GetBinCenter(i),
36  prec, h->GetBinContent(i),
37  prec, h->GetBinError(i));
38  }
39 }
40 
41 TObject* GetO(TDirectory* dir, const char* name="result", TClass* cls=0)
42 {
43  if (!dir) {
44  Warning("GetHS", "No directory");
45  return 0;
46  }
47  TString par; par = gSystem->DirName(name);
48  TString bse; bse = gSystem->BaseName(name);
49 
50  if (!par.EqualTo(".")) {
51  TDirectory* save = dir;
52  dir = save->GetDirectory(par);
53  if (!dir) {
54  Warning("", "Failed to get directory %s", par.Data());
55  // save->ls();
56  return 0;
57  }
58  }
59  TObject* o = dir->Get(bse);
60  if (!o) {
61  Warning("GetHS", "%s not found in %s", name, dir->GetName());
62  return 0;
63  }
64  if (!cls) return o;
65  if (!o->IsA()->InheritsFrom(cls)) {
66  Warning("GetHS", "%s is not a %s!", name, cls->GetName());
67  return 0;
68  }
69  return o;
70 }
71 
72 THStack* GetHS(TDirectory* dir, const char* name="result")
73 {
74  return static_cast<THStack*>(GetO(dir,name,THStack::Class()));
75 }
76 TH1* GetH1(TDirectory* dir, const char* name="result")
77 {
78  return static_cast<TH1*>(GetO(dir,name,TH1::Class()));
79 }
80 
82  const char* var,
83  const char* stackName="result",
84  const char* sub="")
85 {
86  const char* filename = Form("results/combine_%s_0x%x.root", var, flags);
87  TFile* file = TFile::Open(filename, "READ");
88  if (!file) return 0;
89 
90  THStack* stack = GetHS(file, Form("%s%s",sub,stackName));
91  if (!stack) return 0;
92 
93  return stack->GetHists();
94 }
95 
97  const char* var)
98 {
99  const char* filename = Form("results/combine_%s_0x%x.root", var, flags);
100  TFile* file = TFile::Open(filename, "READ");
101  if (!file) return 0;
102 
103  return GetH1(file, "mid");
104 }
106  const char* var)
107 {
108  const char* filename = Form("results/combine_%s_0x%x.root", var, flags);
109  TFile* file = TFile::Open(filename, "READ");
110  if (!file) return 0;
111 
112  return GetH1(file, "cent");
113 }
114 
115 void AddLine(TH1* h)
116 {
117  TLine* l = new TLine(h->GetXaxis()->GetXmin(), 1,
118  h->GetXaxis()->GetXmax(), 1);;
119  l->SetLineStyle(7);
120  l->SetLineWidth(2);
121  h->GetListOfFunctions()->Add(l);
122 }
123 
124 TH1* Compare(TH1* def, TH1* oth)
125 {
126  // Printf("dividing %s by %s", oth->GetName(), def->GetName());
127  TH1* r = static_cast<TH1*>(oth->Clone());
128  r->SetDirectory(0);
129  r->Divide(def);
130  for (Int_t i = 1; i <= def->GetNbinsX(); i++) {
131  Double_t cd = def->GetBinContent(i);
132  Double_t co = oth->GetBinContent(i);
133  Double_t ed = def->GetBinError(i);
134  Double_t eo = oth->GetBinError(i);
135  if (cd < 1e-6 || co < 1e-6) {
136  r->SetBinContent(i, 0);
137  r->SetBinError (i, 0);
138  }
139 
140  Double_t ch = r->GetBinContent(i);
141  r->SetBinError(i, TMath::Max(ed/cd*ch,eo/co*ch));
142  }
143  // Printf("Returning %p (%s)", r, r->GetName());
144  return r;
145 }
146 
147 
149  const char* var,
150  Bool_t noTruth=true)
151 {
152  TList* defs = GetHists(flags, "none");
153  TList* others = GetHists(flags, var);
154  if (!defs || !others) return;
155 
156  TIter nextDef(defs);
157  TIter nextOth(others);
158  THStack* result = new THStack(Form("compare_%s_0x%x",
159  var, flags),
160  Form("\\hbox{%s vs. default }%s\\hbox{ %s}",
161  var, obs,
162  (flags & 0x3) == 0x3 ?
163  "combinatorics" :
164  "injection"));
165  TH1* def = 0;
166  TH1* oth = 0;
167  while ((def = static_cast<TH1*>(nextDef())) &&
168  (oth = static_cast<TH1*>(nextOth()))) {
169  if (noTruth && TString(def->GetName()).EqualTo("truth")) continue;
170  def->SetName("default");
171  oth->SetName(var);
172  TH1* h = Compare(def,oth);
173  h->SetYTitle(Form("%s / default", var));
174  result->Add(h);
175  }
176  TCanvas* c = new TCanvas(result->GetName(),result->GetTitle(), cW, cH);
177  c->SetTopMargin(0.01);
178  c->SetRightMargin(0.01);
179  c->SetTicks();
180  result->Draw("nostack");
181  result->GetHistogram()->SetMinimum(0.9);
182  result->GetHistogram()->SetMaximum(1.2);
183  result->GetHistogram()->GetYaxis()->SetTitleOffset(1.5);
184  result->GetHistogram()->SetXTitle("\\eta");
185  result->GetHistogram()->SetYTitle(Form("%s / default", var));
186  result->SetMinimum(0.9);
187  result->SetMaximum(1.1);
188  AddLine(result->GetHistogram());
189  c->Modified();
190  c->Print(Form("plots/%s.png", result->GetName()));
191 }
192 
194 {
195  Int_t nbin = 9;
196  Double_t bins[] = { 0., 5., 10., 20., 30., 40., 50., 60., 70., 80., };
197  Double_t vals[] = { 1948, 1590, 1180, 786, 512, 318, 183, 96.3, 44.9, };
198  Double_t errs[] = { 38, 32, 31, 20, 15, 12, 8, 5.8, 3.4, };
199  TH1* ret = new TH1D("published", "dNch/deta in |eta|<0.5",
200  nbin, bins);
201  for (Int_t i = 0; i < nbin; i++) {
202  ret->SetBinContent(i+1, vals[i]);
203  ret->SetBinError (i+1, errs[i]);
204  }
205  return ret;
206 }
207 void PrintAxis(TAxis* axis)
208 {
209  printf("%3d bins: ");
210  for (Int_t i = 1; i <= axis->GetNbins(); i++)
211  printf("%6.3f-", axis->GetBinLowEdge(i));
212  Printf("%6.3f", axis->GetBinUpEdge(i));
213 }
215  const char* var)
216 {
217  TH1* def = GetPubl(); // GetMid(flags, "none");
218  TH1* other = GetMid(flags, var);
219  if (!def || !other) return;
220  // PrintAxis(def->GetXaxis());
221  // PrintAxis(other->GetXaxis());
222 
223  def->SetName("default");
224  other->SetName(var);
225  TH1* hs = static_cast<TH1*>(def->Clone("ratioSys"));
226  TH1* ht = static_cast<TH1*>(def->Clone("ratioStat"));
227  hs->SetTitle("Syst.uncer.");
228  ht->SetTitle("Stat.error");
229  hs->SetStats(0);
230  ht->SetStats(0);
231  ht->SetMarkerStyle(20);
232  ht->SetMarkerSize(2);
233  ht->SetMarkerColor(kBlack);
234  ht->SetLineColor(kBlack);
235  ht->SetLineWidth(2);
236  hs->SetFillColor(kBlue-10);
237  hs->SetFillStyle(1001);
238  hs->SetYTitle(Form("%s / published", var));
239  ht->SetYTitle(Form("%s / published", var));
240  hs->Reset();
241  ht->Reset();
242  for (Int_t i = 1; i <= def->GetNbinsX(); i++) {
243  Double_t oc = other->GetBinContent(i);
244  Double_t oe = other->GetBinError (i);
245  Double_t dc = def ->GetBinContent(i);
246  Double_t de = def ->GetBinError (i);
247  Double_t y = oc / dc;
248  Double_t st = oe/oc*y;
249  Double_t sy = de/dc*y;
250  hs->SetBinContent(i, y);
251  ht->SetBinContent(i, y);
252  hs->SetBinError (i, sy);
253  ht->SetBinError (i, st);
254  }
255  THStack* stack = new THStack(Form("mid_%s_0x%x", var, flags), "");
256  stack->Add(hs, "e2");
257  stack->Add(ht);
258  AddLine(ht);
259  // PrintH(def);
260  // PrintH(other);
261  // PrintH(h);
262 
263  TCanvas* c = new TCanvas(ht->GetName(),hs->GetTitle(), cW, cH);
264  c->SetTopMargin(0.01);
265  c->SetRightMargin(0.01);
266  c->SetTicks();
267  stack->SetMinimum(0.9);
268  stack->SetMaximum(1.1);
269  stack->Draw("nostack");
270  TH1* h = stack->GetHistogram();
271  // h->SetMinimum(0.9);
272  // h->SetMaximum(1.1);
273  h->GetYaxis()->SetTitleOffset(1.5);
274  h->SetXTitle("Centrality [%]");
275  h->SetYTitle(Form("%s / published", var));
276 
277  TLegend* l = c->BuildLegend(c->GetLeftMargin(),
278  c->GetBottomMargin(),
279  .7,
280  .3);
281  l->SetFillStyle(0);
282  l->SetBorderSize(0);
283  c->Modified();
284  c->Print(Form("plots/%s.png", h->GetName()));
285 }
286 
287 
289  UShort_t f2,
290  const char* var,
291  Bool_t noTruth=true)
292 {
293  TList* defs = GetHists(f1, var);
294  TList* oths = GetHists(f2, var);
295 
296  TIter nextDef(defs);
297  TIter nextOth(oths);
298  THStack* result = new THStack(Form("compare_%s", var),
299  Form("\\hbox{combinatorics vs. injection }"
300  "%s\\hbox{ %s}", obs, var));
301  TH1* def = 0;
302  TH1* oth = 0;
303  while ((def = static_cast<TH1*>(nextDef())) &&
304  (oth = static_cast<TH1*>(nextOth()))) {
305  if (noTruth && TString(def->GetName()).EqualTo("truth")) continue;
306  def->SetName(f1 == 0x3 ? "combinatorics" : "injection");
307  oth->SetName(f2 == 0x3 ? "combinatorics" : "injection");
308  result->Add(Compare(def,oth));
309  }
310  TCanvas* c = new TCanvas(result->GetName(),result->GetTitle(), cW, cH);
311  c->SetTopMargin(0.01);
312  c->SetRightMargin(0.01);
313  c->SetTicks();
314  result->Draw("nostack");
315  result->GetHistogram()->SetMinimum(.9);
316  result->GetHistogram()->SetMaximum(1.1);
317  result->GetHistogram()->GetYaxis()->SetTitleOffset(1.5);
318  result->GetHistogram()->SetXTitle("\\eta");
319  result->GetHistogram()->SetYTitle(Form("%s / %s",
320  f2==0?"injection":"combinatorial",
321  f1==0?"injection":"combinatorial"));
322  AddLine(result->GetHistogram());
323  result->SetMinimum(.9);
324  result->SetMaximum(1.1);
325  c->Modified();
326  c->Print(Form("plots/%s.png", result->GetName()));
327 }
328 
329 void CompareSubs(UShort_t flgs, const char* var, Bool_t noTruth=false)
330 {
331  TFile* file = TFile::Open(Form("results/combine_none_0x%x.root",flgs),"READ");
332  if (!file) return;
333  TH1* cent = GetH1(file, "cent");
334  if (!cent) return;
335 
336  for (Int_t i = 1; i <= cent->GetNbinsX(); i++) {
337  TString name;
338  Double_t c1 = cent->GetXaxis()->GetBinLowEdge(i);
339  Double_t c2 = cent->GetXaxis()->GetBinUpEdge (i);
340  name.Form("cent%03d%02d_%03d%02d",
341  Int_t(c1), Int_t(c1*100)%100, Int_t(c2), Int_t(c2*100)%100);
342  TString sname(name); sname.Append("/summary");
343  TString pname(sname);
344  pname.ReplaceAll("_", "");
345  pname.ReplaceAll("cent","");
346  pname.ReplaceAll("00", "");
347  TList* defs = GetHists(flgs, "none", sname);
348  TList* oths = GetHists(flgs, var, sname);
349  if (!defs || !oths) break;
350 
351  TIter nextDef(defs);
352  TIter nextOth(oths);
353  THStack* result = new THStack(Form("compare_%s_0x%x_%s",
354  var, flgs, name.Data()),
355  Form("\\hbox{%s vs. default }%s"
356  "\\hbox{ %s %s}",
357  var, obs,
358  (flgs & 0x3) == 0x3 ?
359  "combinatorics" :
360  "injection",
361  pname.Data()));
362  TH1* def = 0;
363  TH1* oth = 0;
364  while ((def = static_cast<TH1*>(nextDef())) &&
365  (oth = static_cast<TH1*>(nextOth()))) {
366  if (noTruth && TString(def->GetName()).EqualTo("truth")) continue;
367  result->Add(Compare(def,oth));
368  }
369  TCanvas* c = new TCanvas(result->GetName(),result->GetTitle(), cW, cH);
370  c->SetTopMargin(0.01);
371  c->SetRightMargin(0.01);
372  result->Draw("nostack");
373  result->GetHistogram()->SetMinimum(.9);
374  result->GetHistogram()->SetMaximum(1.1);
375  result->GetHistogram()->GetYaxis()->SetTitleOffset(1.5);
376  result->GetHistogram()->SetXTitle("\\eta");
377  result->GetHistogram()->SetYTitle(Form("%s / default", var));
378  result->SetMinimum(.9);
379  result->SetMaximum(1.1);
380  AddLine(result->GetHistogram());
381  c->Modified();
382  c->Print(Form("plots/%s.png", result->GetName()));
383  }
384 }
385 
386 void CentPlot(UShort_t flags)
387 {
388  TList* defs = GetHists(flags, "none");
389  TH1* cent = GetCent(flags, "none");
390  // defs->ls();
391  TCanvas* c = new TCanvas("legend", "Legend", cW, cH);
392  c->SetTopMargin(0.01);
393  c->SetRightMargin(0.01);
394  TLegend* l = new TLegend(0.01, 0.01, .99, .99);
395  l->SetBorderSize(0);
396  l->SetFillStyle(0);
397  // l->SetNColumns(2);
398  TIter next(defs);
399  TH1* h = 0;
400  Int_t bin = 1;
401  while ((h = static_cast<TH1*>(next()))) {
402  if (TString(h->GetName()).Contains("truth")) continue;
403  TString nme; nme.Form("%2d - %3d%%",
404  Int_t(cent->GetXaxis()->GetBinLowEdge(bin)),
405  Int_t(cent->GetXaxis()->GetBinUpEdge(bin)));
406  bin++;
407  TLegendEntry* e = l->AddEntry("", nme, "F");
408  e->SetFillColor(h->GetMarkerColor());
409  e->SetFillStyle(1001);
410  }
411  c->cd();
412  l->Draw();
413  c->Modified();
414  c->Print(Form("plots/%s.png", c->GetName()));
415 }
416 
417 
418 
419 
420 void
421 Compare(const char* var)
422 {
423  if (!gROOT->GetClass("GraphSysErr"))
424  gROOT->LoadMacro("$HOME/GraphSysErr/GraphSysErr.C+g");
425  gSystem->mkdir("plots",1);
426  TString w(var);
427  // CompareVars(0x0,var);
428  if (!w.EqualTo("none")) CompareVars(0x3,var);
429  // CompareFlgs(0x0,0x3,var);
430  // CompareMids(0x0,var);
431  CompareMids(0x3,var);
432  if (w.EqualTo("none")) CentPlot(0x3);
433  // CompareSubs(0x3,var);
434 }
TH1 * GetPubl()
Definition: Compare.C:193
const char * filename
Definition: TestFCM.C:1
double Double_t
Definition: External.C:58
TSystem * gSystem
void CompareFlgs(UShort_t f1, UShort_t f2, const char *var, Bool_t noTruth=true)
Definition: Compare.C:288
TCanvas * c
Definition: TestFitELoss.C:172
const char * obs
Definition: Compare.C:20
AliStack * stack
void CompareMids(UShort_t flags, const char *var)
Definition: Compare.C:214
TH1 * GetH1(TDirectory *dir, const char *name="result")
Definition: Compare.C:76
void PrintAxis(TAxis *axis)
Definition: Compare.C:207
int Int_t
Definition: External.C:63
Int_t cW
Definition: Compare.C:18
const Bool_t kCompareVarLoaded
Definition: Compare.C:17
Definition: External.C:212
TH1 * Compare(TH1 *def, TH1 *oth)
Definition: Compare.C:124
TH1 * GetCent(UShort_t flags, const char *var)
Definition: Compare.C:105
void CentPlot(UShort_t flags)
Definition: Compare.C:386
Int_t cH
Definition: Compare.C:19
void AddLine(TH1 *h)
Definition: Compare.C:115
TFile * file
TList with histograms for a given trigger.
void PrintH(TH1 *h, Int_t prec=2)
Definition: Compare.C:29
unsigned short UShort_t
Definition: External.C:28
THStack * GetHS(TDirectory *dir, const char *name="result")
Definition: Compare.C:72
TList * GetHists(UShort_t flags, const char *var, const char *stackName="result", const char *sub="")
Definition: Compare.C:81
bool Bool_t
Definition: External.C:53
TH1 * GetMid(UShort_t flags, const char *var)
Definition: Compare.C:96
TObject * GetO(TDirectory *dir, const char *name="result", TClass *cls=0)
Definition: Compare.C:41
Definition: External.C:196
void CompareVars(UShort_t flags, const char *var, Bool_t noTruth=true)
Definition: Compare.C:148
TDirectoryFile * dir
void CompareSubs(UShort_t flgs, const char *var, Bool_t noTruth=false)
Definition: Compare.C:329