AliPhysics  e34b7ac (e34b7ac)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
ComparepPb.C
Go to the documentation of this file.
1 #include <TMultiGraph.h>
2 #include <TString.h>
3 #include <TCanvas.h>
4 #include <TROOT.h>
5 #include <TGraphAsymmErrors.h>
6 #include <TGraphErrors.h>
7 #include <TError.h>
8 #include <TLegend.h>
9 #include <TLegendEntry.h>
10 #include <TH1.h>
11 #include <TAxis.h>
12 #include <THStack.h>
13 
14 //____________________________________________________________________
23 Int_t Marker(Int_t sys, Int_t style)
24 {
25  return (sys == 4 ? (style == 23 ? 9 : 4) : 0) + style;
26 }
27 
28 //____________________________________________________________________
38 TMultiGraph*
39 GetOneFwd(UShort_t sys, const char* meth, Int_t style)
40 {
41 
42 
43  TString what(meth);
44  what.ReplaceAll("X", (sys == 3 ? "A" : "C"));
45  TString sysN(sys == 3 ? "pPb" : "Pbp");
46 
47  TString cmd;
48  cmd = Form("Drawer::GetStack(0,\"%s\",5023,\"CENT%s\")",
49  sysN.Data(), meth);
50 
51  Long_t ret = gROOT->ProcessLine(cmd);
52  if (!ret) {
53  Warning("GetOne", "didn't get anything back for %s", what.Data());
54  return 0;
55  }
56  THStack* stack = reinterpret_cast<THStack*>(ret);
57  if (!stack) {
58  Warning("GetOne", "didn't get anything back for %s", what.Data());
59  return 0;
60  }
61  Int_t sty = Marker(sys, style);
62  TMultiGraph* mg = new TMultiGraph;
63  TIter nextH(stack->GetHists());
64  TH1* hist = 0;
65  while ((hist = static_cast<TH1*>(nextH()))) {
66  ret = gROOT->ProcessLine(Form("Drawer::H2G((TH1*)%p,0)",hist));
67  if (!ret) {
68  Warning("GetOne", "Couldn't convert %s", hist->GetName());
69  continue;
70  }
71  TGraphErrors* g = reinterpret_cast<TGraphErrors*>(ret);
72  g->SetMarkerStyle(sty);
73  mg->Add(g);
74  }
75 
76  Info("", "Got %d graphs", mg->GetListOfGraphs()->GetEntries());
77  mg->GetListOfGraphs()->Sort();
78  mg->SetName(Form("fwd_%s_%s", sysN.Data(), meth));
79  return mg;
80 }
81 //____________________________________________________________________
91 TMultiGraph*
92 GetOneCen(UShort_t sys, const char* meth, Int_t style)
93 {
94  TString what(meth);
95  what.ReplaceAll("X", (sys == 3 ? "A" : "C"));
96  UShort_t cnt = 0;
97  if (what.Contains("V0M")) cnt = (0x01 << 4);
98  else if (what.Contains("V0A")) cnt = (0x02 << 4);
99  else if (what.Contains("ZNA")) cnt = (0x04 << 4);
100  else if (what.Contains("ZNC")) cnt = (0x08 << 4);
101  else if (what.Contains("V0C")) cnt = (0x10 << 4);
102  else if (what.Contains("CL1")) cnt = (0x20 << 4);
103  Info("GetOne", "%d/%s -> 0x%04x", sys, meth, cnt);
104 
105  Long_t ret =
106  gROOT->ProcessLine(Form("RefData::GetData(%d,5023,0x%x,0,100,0xf);",
107  sys, cnt));
108  if (!ret) {
109  Warning("GetOne", "didn't get anything back for %s", what.Data());
110  return 0;
111  }
112  TMultiGraph* mg = reinterpret_cast<TMultiGraph*>(ret);
113  if (!mg) {
114  Warning("GetOne", "didn't get anything back for %s", what.Data());
115  return 0;
116  }
117 
118  Int_t sty = Marker(sys, style);
119  TIter next(mg->GetListOfGraphs());
120  TGraph* g = 0;
121  while ((g = static_cast<TGraph*>(next())))
122  g->SetMarkerStyle(sty);
123  Info("", "Got %d graphs", mg->GetListOfGraphs()->GetEntries());
124  mg->GetListOfGraphs()->Sort();
125  mg->SetName(Form("fwd_%s_%s", (sys == 3 ? "pPb" : "Pbp"), meth));
126  return mg;
127 }
128 
129 //____________________________________________________________________
135 void
137 {
138  Int_t n = g->GetN();
139  for (Int_t i = 0; i < n/2; i++) {
140  Int_t i1 = i;
141  Int_t i2 = n-i-1;
142  if (i1 == i2) continue;
143  Double_t x1 = g->GetX()[i1];
144  Double_t x2 = g->GetX()[i2];
145  Double_t y1 = g->GetY()[i1];
146  Double_t y2 = g->GetY()[i2];
147  Double_t ex1 = g->GetEX()[i1];
148  Double_t ex2 = g->GetEX()[i2];
149  Double_t ey1 = g->GetEY()[i1];
150  Double_t ey2 = g->GetEY()[i2];
151  g->SetPoint(i1, -x2, y2);
152  g->SetPointError(i1, ex2, ey2);
153  g->SetPoint(i2, -x1, y1);
154  g->SetPointError(i2, ex1, ey1);
155  }
156 }
157 //____________________________________________________________________
163 void
165 {
166  Int_t n = g->GetN();
167  for (Int_t i = 0; i < n/2; i++) {
168  Int_t i1 = i;
169  Int_t i2 = n-i-1;
170  if (i1 == i2) continue;
171  Double_t x1 = g->GetX()[i1];
172  Double_t x2 = g->GetX()[i2];
173  Double_t y1 = g->GetY()[i1];
174  Double_t y2 = g->GetY()[i2];
175  Double_t exl1 = g->GetEXlow()[i1];
176  Double_t exl2 = g->GetEXlow()[i2];
177  Double_t exh1 = g->GetEXhigh()[i1];
178  Double_t exh2 = g->GetEXhigh()[i2];
179  Double_t eyl1 = g->GetEYlow()[i1];
180  Double_t eyl2 = g->GetEYlow()[i2];
181  Double_t eyh1 = g->GetEYhigh()[i1];
182  Double_t eyh2 = g->GetEYhigh()[i2];
183  g->SetPoint(i1, -x2, y2);
184  g->SetPointError(i1, exl2, exh2, eyl2, eyh2);
185  g->SetPoint(i2, -x1, y1);
186  g->SetPointError(i2, exl1, exh1, eyl1, eyh1);
187  }
188 }
189 //____________________________________________________________________
195 void
197 {
198  if (g->IsA()->InheritsFrom(TGraphAsymmErrors::Class()))
199  FlipG(static_cast<TGraphAsymmErrors*>(g));
200  else if (g->IsA()->InheritsFrom(TGraphErrors::Class()))
201  FlipG(static_cast<TGraphErrors*>(g));
202  else
203  Warning("FlipG", "Don't know how to flip a %s", g->ClassName());
204 }
205 
206 //____________________________________________________________________
214 TMultiGraph*
215 FlipMG(TMultiGraph* mg)
216 {
217  TMultiGraph* ret = static_cast<TMultiGraph*>(mg->Clone());
218  TIter next(ret->GetListOfGraphs());
219  TGraph* g = 0;
220  while ((g = static_cast<TGraph*>(next())))
221  FlipG(g);
222  return ret;
223 }
224 
225 //____________________________________________________________________
232 void
234 {
235  Info("RatioGG(GE)", "%s (%s:%d)/ %s (%s:%d)",
236  num->GetName(), num->ClassName(), num->GetN(),
237  denom->GetName(), denom->ClassName(), denom->GetN());
238  Int_t n = num->GetN();
239  Double_t dh = denom->GetX()[n-1]+denom->GetEX()[n-1];
240  Double_t dl = denom->GetX()[0] +denom->GetEX()[0];
241  TBits clean(n);
242  for (Int_t i = 0; i < n; i++) {
243  Double_t x = num->GetX()[i];
244  if (x > dh || x < dl) {
245  clean.SetBitNumber(i);
246  continue;
247  }
248  Double_t d = denom->Eval(x);
249  Double_t y = num->GetY()[i];
250  Double_t ex = num->GetEX()[i];
251  Double_t ey = num->GetEY()[i];
252  num->SetPoint(i, x, y/d);
253  num->SetPointError(i, ex, ey/d);
254  }
255  for (Int_t i = clean.GetNbits()-1; i >= 0; i--) {
256  if (clean.TestBitNumber(i)) num->RemovePoint(i);
257  }
258 }
259 //____________________________________________________________________
266 void
268 {
269  Info("RatioGG(GA)", "%s (%s:%d)/ %s (%s:%d)",
270  num->GetName(), num->ClassName(), num->GetN(),
271  denom->GetName(), denom->ClassName(), denom->GetN());
272  Int_t n = num->GetN();
273  for (Int_t i = 0; i < n; i++) {
274  Double_t x = num->GetX()[i];
275  Double_t d = denom->Eval(x);
276  Double_t y = num->GetY()[i];
277  Double_t exl = num->GetEXlow()[i];
278  Double_t exh = num->GetEXhigh()[i];
279  Double_t eyl = num->GetEYlow()[i];
280  Double_t eyh = num->GetEYhigh()[i];
281  num->SetPoint(i, x, y/d);
282  num->SetPointError(i, exl, exh, eyl/d, eyh/d);
283  }
284 }
285 //____________________________________________________________________
292 void
293 RatioGG(TGraph* num, TGraph* denom)
294 {
295  if (num->IsA()->InheritsFrom(TGraphAsymmErrors::Class()))
296  RatioGG(static_cast<TGraphAsymmErrors*>(num),denom);
297  else if (num->IsA()->InheritsFrom(TGraphErrors::Class()))
298  RatioGG(static_cast<TGraphErrors*>(num),denom);
299  else
300  Warning("RatioGG", "Don't know how to ratio a %s to a %s",
301  num->ClassName(), denom->ClassName());
302 }
303 
304 
305 //____________________________________________________________________
314 TMultiGraph*
315 RatioMG(TMultiGraph* num, TMultiGraph* denom)
316 {
317  TMultiGraph* ret = static_cast<TMultiGraph*>(num->Clone());
318  TIter nextN(ret->GetListOfGraphs());
319  TIter nextD(denom->GetListOfGraphs());
320  TGraph* n = 0;
321  TGraph* d = 0;
322  while ((n = static_cast<TGraph*>(nextN())) &&
323  (d = static_cast<TGraph*>(nextD())))
324  RatioGG(n, d);
325  return ret;
326 }
327 
328 //____________________________________________________________________
337 void
339  Bool_t hasFwd,
340  Double_t fac=1.2,
341  const char* ytitle="1/#it{N} d#it{N}_{ch}/d#it{#eta}")
342 {
343  static Int_t cnt = 0;
344  TMultiGraph* mg = static_cast<TMultiGraph*>(o);
345  TH1* h = mg->GetHistogram();
346  mg->SetTitle("");
347  if (!h) return;
348 
349 
350  Double_t etaMin = (hasFwd ? -6 : -2.5);
351  Double_t etaMax = (hasFwd ? +6 : +2.5);
352  h->GetXaxis()->Set((etaMax-etaMin)/0.5, etaMin, etaMax);
353  h->Rebuild();
354  h->SetName(Form("%s_%d", o->GetName(), cnt++));
355  h->SetTitle("");
356  h->SetXTitle("#it{#eta}");
357  h->SetYTitle(ytitle);
358  h->GetXaxis()->SetNdivisions(210);
359  h->GetYaxis()->SetNdivisions(210);
360  h->GetXaxis()->SetLabelFont(42);
361  h->GetXaxis()->SetLabelSize(0.05);
362  h->GetXaxis()->SetTitleFont(42);
363  h->GetXaxis()->SetTitleSize(0.05);
364  h->GetXaxis()->SetTitleOffset(0.7);
365  h->GetYaxis()->SetLabelFont(42);
366  h->GetYaxis()->SetLabelSize(0.05);
367  h->GetYaxis()->SetTitleFont(42);
368  h->GetYaxis()->SetTitleSize(0.05);
369  h->GetYaxis()->SetTitleOffset(0.7);
370  h->SetMaximum(fac*h->GetMaximum());
371  h->SetMinimum(0.01);
372 }
373 
374 //____________________________________________________________________
382 const char*
383 ExtractCent(const char* name)
384 {
385  TString n(name);
386  TObjArray* a = n.Tokenize("_");
387  TObjString* oc1 = static_cast<TObjString*>(a->At(1));
388  TObjString* oc2 = static_cast<TObjString*>(a->At(2));
389  TString sc1 = oc1->String().Strip(TString::kLeading, '0');
390  TString sc2 = oc2->String().Strip(TString::kLeading, '0');
391  delete a;
392  return Form("%3d-%3d%%", sc1.Atoi(), sc2.Atoi());
393 }
394 
395 //____________________________________________________________________
402 void
403 CompareOne(const char* meth, Int_t style)
404 {
405  TMultiGraph* fwdpPb = GetOneFwd(3, meth, style);
406  TMultiGraph* fwdPbp = GetOneFwd(4, meth, style);
407  Bool_t hasFwd = true;
408  if (!fwdpPb || !fwdPbp) {
409  Warning("CompareOne", "FWD pPb=%p Pbp=%p", fwdpPb, fwdPbp);
410  hasFwd = false;
411  }
412  TMultiGraph* cenpPb = GetOneCen(3, meth, style);
413  TMultiGraph* cenPbp = GetOneCen(4, meth, style);
414  Bool_t hasCen = true;
415  if (!cenpPb || !cenPbp) {
416  Warning("CompareOne", "CEN pPb=%p Pbp=%p", cenpPb, cenPbp);
417  hasCen = false;
418  }
419  if (!hasCen && !hasFwd) return;
420 
421  TCanvas* c = new TCanvas(meth, meth, 900, 1000);
422  c->SetTopMargin(0.01);
423  c->SetRightMargin(0.01);
424  c->cd();
425  c->Divide(1,3,0,0);
426  Double_t lx = 0.8;
427  c->GetPad(1)->SetRightMargin(1-lx+.01);
428  c->GetPad(2)->SetRightMargin(1-lx+.01);
429  c->GetPad(3)->SetRightMargin(1-lx+.01);
430  c->GetPad(1)->SetGridx();
431  c->GetPad(2)->SetGridx();
432  c->GetPad(3)->SetGridx();
433  // c->GetPad(2)->SetBottomMargin(0.1);
434 
435  c->cd(1);
436  TMultiGraph* frame = 0;
437  if (hasFwd) {
438  fwdpPb->Draw("pa");
439  fwdPbp->Draw("p");
440  frame = fwdpPb;
441  }
442  if (hasCen) {
443  cenpPb->Draw(frame ? "p" : "pa");
444  cenPbp->Draw("p");
445  if (!frame) frame = cenpPb;
446  }
447  FixFrame(frame, hasFwd, 1.2);
448 
449  TLegend* l = new TLegend(lx, 0.01, 0.99, 0.99,
450  Form("%s estimator", meth));
451  l->SetFillStyle(0);
452  l->SetBorderSize(0);
453  l->SetTextFont(42);
454  l->SetNColumns(1);
455  TLegendEntry* e = l->AddEntry("", "p-Pb", "p");
456  e->SetMarkerStyle(Marker(3, style));
457  e = l->AddEntry("", "Pb-p", "p");
458  e->SetMarkerStyle(Marker(4, style));
459  l->Draw();
460 
461  c->cd(2);
462  TMultiGraph* ffwdPbp = 0;
463  if (hasFwd) {
464  ffwdPbp = FlipMG(fwdPbp);
465  ffwdPbp->SetName(Form("fwd_%s_fPbp", meth));
466  fwdpPb->Draw("ap");
467  ffwdPbp->Draw("p");
468  }
469  TMultiGraph* fcenPbp = 0;
470  if (hasCen) {
471  fcenPbp = FlipMG(cenPbp);
472  fcenPbp->SetName(Form("cen_%s_fPbp", meth));
473  cenpPb->Draw(hasFwd ? "p" : "ap");
474  fcenPbp->Draw("p");
475  }
476 
477  TLegend* ll = new TLegend(lx, 0.01, 0.99, 0.9,
478  "Pb-p flipped around #it{#eta}=0");
479  ll->SetFillStyle(0);
480  ll->SetBorderSize(0);
481  ll->SetTextFont(42);
482  ll->SetNColumns(1);
483  TIter next(frame->GetListOfGraphs());
484  TGraph* g = 0;
485  while ((g = static_cast<TGraph*>(next()))) {
486  e = ll->AddEntry("", ExtractCent(g->GetName()), "f");
487  e->SetFillColor(g->GetMarkerColor());
488  e->SetFillStyle(1001);
489  e->SetLineColor(g->GetMarkerColor());
490  }
491  ll->Draw();
492 
493  c->cd(3);
494  frame = 0;
495 
496  Double_t minEta = (hasFwd ? -3.5 : -2.);
497  Double_t maxEta = (hasFwd ? +5.0 : +2.);
498  TGraphErrors* band = new TGraphErrors(2);
499  band->SetPoint(0, minEta, 1);
500  band->SetPoint(1, maxEta, 1);
501  band->SetPointError(0, 0, .05);
502  band->SetPointError(1, 0, .05);
503  band->SetFillColor(kYellow-8);
504  band->SetFillStyle(3001);
505  band->SetLineColor(kYellow-8);
506 
507  if (hasFwd) {
508  TMultiGraph* fwdRat = RatioMG(ffwdPbp, fwdpPb);
509  fwdRat->SetName(Form("fwd_%s_rat", meth));
510  fwdRat->GetListOfGraphs()->AddFirst(band, "3");
511  frame = fwdRat;
512  fwdRat->Draw("ap");
513  }
514  if (hasCen) {
515  TMultiGraph* cenRat = RatioMG(fcenPbp, cenpPb);
516  cenRat->SetName(Form("cen_%s_rat", meth));
517  if (!frame) {
518  cenRat->GetListOfGraphs()->AddFirst(band, "3");
519  frame = cenRat;
520  }
521  cenRat->Draw(hasFwd ? "p" : "ap");
522  }
523  FixFrame(frame, hasFwd, 1, "Flipped Pb-p over p-Pb");
524  if (meth[0] != 'C') {
525  frame->GetHistogram()->SetMinimum(0.74);
526  frame->GetHistogram()->SetMaximum(1.26);
527  }
528  TLegend* lll = new TLegend(lx, 0.4, 0.99, 0.5);
529  lll->SetFillStyle(0);
530  lll->SetBorderSize(0);
531  lll->SetTextFont(42);
532  lll->SetNColumns(1);
533  lll->AddEntry(band, "5% error band", "F");
534  lll->Draw();
535 
536  c->Print(Form("%s.pdf", meth));
537 }
538 
539 //____________________________________________________________________
540 void
542 {
543  gROOT->LoadMacro("$ALICE_PHYSICS/PWGLF/FORWARD/analysis2/dndeta/Drawer.C+g");
544  gROOT->LoadMacro("$ALICE_PHYSICS/PWGLF/FORWARD/analysis2/OtherData.C+g");
545 
546  const char* meths[] = { "V0M", "CL1", "V0X", "ZNX", 0 };
547  Int_t styles[] = { 20, 21, 22, 23, 0 };
548 
549  const char** pmeth = meths;
550  Int_t* pstyle = styles;
551 
552  while (*pmeth && *pstyle) {
553  CompareOne(*pmeth, *pstyle);
554  pmeth++;
555  pstyle++;
556  }
557 }
558 //____________________________________________________________________
559 //
560 // EOF
561 //
double Double_t
Definition: External.C:58
TMultiGraph * FlipMG(TMultiGraph *mg)
Definition: ComparepPb.C:215
TCanvas * c
Definition: TestFitELoss.C:172
TMultiGraph * RatioMG(TMultiGraph *num, TMultiGraph *denom)
Definition: ComparepPb.C:315
void FixFrame(TObject *o, Bool_t hasFwd, Double_t fac=1.2, const char *ytitle="1/#it{N} d#it{N}_{ch}/d#it{#eta}")
Definition: ComparepPb.C:338
int Int_t
Definition: External.C:63
void FlipG(TGraphErrors *g)
Definition: ComparepPb.C:136
TMultiGraph * GetOneCen(UShort_t sys, const char *meth, Int_t style)
Definition: ComparepPb.C:92
Int_t Marker(Int_t sys, Int_t style)
Definition: ComparepPb.C:23
TMultiGraph * GetOneFwd(UShort_t sys, const char *meth, Int_t style)
Definition: ComparepPb.C:39
unsigned short UShort_t
Definition: External.C:28
bool Bool_t
Definition: External.C:53
void CompareOne(const char *meth, Int_t style)
Definition: ComparepPb.C:403
const char * ExtractCent(const char *name)
Definition: ComparepPb.C:383
Definition: External.C:196
void RatioGG(TGraphErrors *num, TGraph *denom)
Definition: ComparepPb.C:233
void ComparepPb()
Definition: ComparepPb.C:541