AliPhysics  5eaf189 (5eaf189)
ExtractGSE2.C
Go to the documentation of this file.
1 
11 #ifndef __CINT__
12 # include "GraphSysErr.C"
13 # include <TFile.h>
14 # include <TList.h>
15 # include <TMath.h>
16 # include <TCanvas.h>
17 # include <TLegend.h>
18 # include <TSeqCollection.h>
19 # include <TLegendEntry.h>
20 # include <THStack.h>
21 #else
22 class TSeqCollection;
23 class TCanvas;
24 class TGraphErrors;
25 class THStack;
26 class TLegend;
27 class TH1;
28 #endif
29 
30 //====================================================================
40 {
41  Double_t c = (c1+c2) / 2;
42  if (c < 5) return 0;
43  else if (c < 10) return 1;
44  else if (c < 20) return 2;
45  else if (c < 30) return 3;
46  else if (c < 40) return 4;
47  else if (c < 50) return 5;
48  else if (c < 60) return 6;
49  else if (c < 70) return 7;
50  else if (c < 80) return 8;
51  else if (c < 90) return 9;
52  return 10;
53 }
54 
55 //____________________________________________________________________
64 static Color_t PbPbColor(Double_t c1, Double_t c2)
65 {
66  const Color_t cc[] = { kMagenta+2,
67  kBlue+2,
68  kAzure-1, // 10,
69  kCyan+2,
70  kGreen+1,
71  kSpring+5,//+10,
72  kYellow+1,
73  kOrange+5,//+10,
74  kRed+1,
75  kPink+5,//+10,
76  kBlack };
77  Int_t bin = PbPbBin(c1,c2);
78  return cc[bin];
79 }
80 //____________________________________________________________________
89 void MakeCommon(TObject* o, const char* name, Double_t val, Color_t c)
90 {
91  GraphSysErr* gse = static_cast<GraphSysErr*>(o);
92  Int_t id = gse->DefineCommon(name,true,val,GraphSysErr::kBox);
93  gse->SetSysFillColor(id, c);
94  gse->SetSysLineColor(id, c);
95 }
96 //____________________________________________________________________
106 Int_t MakeP2P(TObject* o, const char* name, Color_t c)
107 {
108  GraphSysErr* gse = static_cast<GraphSysErr*>(o);
109  Int_t id = gse->DeclarePoint2Point(name,true,GraphSysErr::kBox);
110  gse->SetSysFillColor(id, c);
111  gse->SetSysLineColor(id, c);
112  return id;
113 }
114 
115 //====================================================================
117 {
118  return sMin + TMath::Power(x/xMax, 2)*(sMax-sMin);
119 }
120 
122 {
123  return SysEval(x, sMin, sMax, 80);
124 }
126 {
127  return SysEval(x, sMin, sMax, 2);
128 }
129 
130 //____________________________________________________________________
142 TObject* MakeGSE(TDirectory* d,
143  Int_t dimen,
144  UShort_t sNN,
145  Double_t c1,
146  Double_t c2)
147 {
148  if (!gROOT->GetClass("GraphSysErr")) return 0;
149  TString bin; bin.Form("%03dd%02d_%03dd%02d",
150  Int_t(c1), Int_t(c1*100)%100,
151  Int_t(c2), Int_t(c2*100)%100);
152  Bool_t isAll = (c1+1.e-9 >= c2);
153  if (isAll) bin = "all";
154  TString sub(bin);
155  if (!isAll) sub.Prepend("cent");
156  sub.Append(Form("/results%dd/result", dimen));
157  TString nme(bin); nme.Prepend("CENT_");
158  TH1* g = GetH1(d, sub);
159  if (!g) return 0;
160 
161  Double_t eff = g->GetBinContent(0);
162  Printf("Trigger efficiency: %6.4f", eff);
163  if (eff < 1e-6) eff = 1;
164  g->Scale(eff);
165  g->SetBinContent(0,0);
166 
167  Color_t col = PbPbColor(c1,c2); // g->GetMarkerColor();
168  // Double_t bg = (1-c1/100)*(2-0.1)+0.1;
169  // Double_t c = TMath::Power(c1/100,2)*(6.2-0.4)+0.4;
170  // Preliminary centrality sys error for Xe-Xe
171  Double_t aMin = sNN == 5440 ? -0.09 : 0.02;
172  Double_t pComp = sNN == 5440 ? 0.002 : 0.01;
173  Double_t wDcy = sNN == 5440 ? 0.004 : 0.01;
174  Double_t ptInt = sNN == 5440 ? 0.004 : 0.02;
175  Double_t ptMax = sNN == 5440 ? 0.010 : 0.00;
176  Double_t egDep = sNN == 5440 ? 0.000 : 0.02;
177  Double_t cMin = sNN == 5440 ? 0.001 : sNN == 5023 ? 0.005 : 0.004;
178  Double_t cMax = sNN == 5440 ? 0.048 : sNN == 5023 ? 0.075 : 0.062;
179  Double_t tMax = sNN == 5440 ? 0.040 : 0;
180  Double_t tMin = sNN == 5440 ? 0.003 : 0;
181  Double_t bgMax = sNN == 5440 ? 0.011 : 0.02;
182  Double_t bgMin = sNN == 5440 ? 0.005 : 0.001;
183  Double_t bg = CSysEval(c2, bgMax, bgMin);
184  Double_t c = CSysEval(c2, cMin, cMax);
185  Double_t em1 = sNN == 5440 ? 0.012 : sNN == 5023 ? 0.04 : 0.02;
186  Double_t em2 = sNN == 5440 ? 0.040 : 5023 ? 0.04 : 0.02;
187  Double_t em = c1 >= 80 ? em2 : c1 >= 70 ? em1 : 0;
188  Double_t tkl = CSysEval(c2, tMax, tMin);
189  GraphSysErr* gse = new GraphSysErr(g->GetNbinsX());
190  gse->SetName(nme);
191  gse->SetTitle(Form("%5.1f - %5.1f%%", c1, c2));
192  //gse->SetKey("author",(sNN == 5023?"PREGHENELLA : 2015":"SHAHOYAN : 2013"));
193  gse->SetKey("title", Form("dNch/deta in PbPb at %d GeV", sNN));
194  gse->SetKey("obskey", "DN/DETARAP");
195  gse->SetKey("reackey", "PB PB --> CHARGED X");
196  gse->SetKey("laboratory", "CERN");
197  gse->SetKey("accelerator", "LHC");
198  gse->SetKey("detector", "TRACKLETS");
199  gse->SetKey("reference", sNN==5023 ? "ALICE-AN-2830" : "ALICE-AN-2180");
200  if (!isAll)gse->AddQualifier("CENTRALITY IN PCT", Form("%.1f TO %.1f",c1,c2));
201  gse->AddQualifier("SQRT(S)/NUCLEON IN GEV", Form("%d", sNN));
202  gse->SetXTitle("ETARAP");
203  gse->SetYTitle("DN/DETARAP");
204  gse->SetMarkerStyle(g->GetMarkerStyle());
205  gse->SetMarkerSize(g->GetMarkerSize());
207  gse->SetMarkerColor(col);
208  gse->SetLineColor(col);
209  gse->SetFillColor(col);
210  gse->SetSumFillColor(col);
211  gse->SetSumLineColor(col);
213  gse->SetCommonSumFillColor(col);
214  gse->SetCommonSumLineColor(col);
216  MakeCommon(gse, "Particle composition", pComp, col);
217  MakeCommon(gse, "Weak decay", wDcy, col);
218  Int_t pte = -1;
219  if (ptMax < 1e-3)
220  MakeCommon(gse, "pT extrapolation", ptInt, col);
221  else
222  pte = MakeP2P(gse, "pT extrapolation", col);
223  MakeCommon(gse, "Background subtraction", bg, col);
224  MakeCommon(gse, "Centrality", c, col);
225  MakeCommon(gse, "EM contamination", em, col);
226  if (egDep > 1e-6)
227  MakeCommon(gse, "EG dependence", egDep, col);
228  if (sNN==5440)
229  MakeCommon(gse, "Material budget", 0.001, col);
230  if (tkl > 1e-6)
231  MakeCommon(gse, "Tracklet selection", tkl, col);
232 
233  if (eff) {
234  MakeCommon(gse, "TRIGGER", 0.02, col);
235  }
236  Int_t acc = -1;
237  if (aMin > 0)
238  acc = MakeP2P(gse, "Acceptance", col);
239  else
240  MakeCommon(gse, "Acceptance", -aMin, col);
241 
242 
243  Int_t j = 0;
244  for (Int_t i = 1; i <= g->GetNbinsX(); i++) {
245  Double_t eta = g->GetXaxis()->GetBinCenter(i);
246  Double_t eEta = g->GetXaxis()->GetBinWidth(i)/2;
247  Double_t xo = TMath::Abs(eta)+eEta;
248  if (xo > 2) continue;
249  Double_t ea = aMin*TMath::Power(xo/2,2);
250  Double_t ep = EtaSysEval(xo,ptInt, ptMax);
251  gse->SetPoint(j, eta, g->GetBinContent(i));
252  gse->SetPointError(j, eEta, eEta);
253  gse->SetStatError(j, g->GetBinError(i),g->GetBinError(i));
254  if (acc > 0) gse->SetSysError(acc, j, eEta, eEta, ea/100, ea/100);
255  if (pte > 0) gse->SetSysError(pte, j, eEta, eEta, ep/100, ep/100);
256  j++;
257  }
258  return gse;
259 
260 
261 }
262 
263 //____________________________________________________________________
275 TObject* MakeTGSE(TDirectory* d,
276  Int_t dimen,
277  UShort_t sNN,
278  Double_t c1,
279  Double_t c2)
280 {
281  if (!gROOT->GetClass("GraphSysErr")) return 0;
282  TString bin; bin.Form("%03dd%02d_%03dd%02d",
283  Int_t(c1), Int_t(c1*100)%100,
284  Int_t(c2), Int_t(c2*100)%100);
285  Bool_t isAll = (c1+1.e-9 >= c2);
286  if (isAll) bin = "all";
287  TString sub(bin);
288  if (!isAll) sub.Prepend("cent");
289  sub.Append(Form("/results%dd/simG",dimen));
290  TString nme(bin); nme.Prepend("CENTT_");
291  TH1* g = GetH1(d, sub);
292  if (!g) return 0;
293 
294  Color_t col = PbPbColor(c1,c2); // g->GetMarkerColor();
295  // Double_t bg = (1-c1/100)*(2-0.1)+0.1;
296  // Double_t c = TMath::Power(c1/100,2)*(6.2-0.4)+0.4;
297  GraphSysErr* gse = new GraphSysErr(g->GetNbinsX());
298  gse->SetName(nme);
299  gse->SetTitle(Form("%5.1f - %5.1f%%", c1, c2));
300  gse->SetKey("author", (sNN == 5023 ? "PREGHENELLA : 2015":"SHAHOYAN : 2013"));
301  gse->SetKey("title", Form("dNch/deta in PbPb at %d GeV", sNN));
302  gse->SetKey("obskey", "DN/DETARAP");
303  gse->SetKey("reackey", "PB PB --> CHARGED X");
304  gse->SetKey("laboratory", "CERN");
305  gse->SetKey("accelerator", "LHC");
306  gse->SetKey("detector", "TRACKLETS");
307  gse->SetKey("reference", sNN==5023 ? "ALICE-AN-2830" : "ALICE-AN-2180");
308  if (!isAll)gse->AddQualifier("CENTRALITY IN PCT", Form("%.1f TO %.1f",c1,c2));
309  gse->AddQualifier("SQRT(S)/NUCLEON IN GEV", Form("%d", sNN));
310  gse->SetXTitle("ETARAP");
311  gse->SetYTitle("DN/DETARAP");
312  gse->SetMarkerStyle(g->GetMarkerStyle());
313  gse->SetMarkerSize(g->GetMarkerSize());
315  gse->SetMarkerColor(col);
316  gse->SetLineColor(col);
317  gse->SetFillColor(col);
318  gse->SetSumFillColor(col);
319  gse->SetSumLineColor(col);
321  gse->SetCommonSumFillColor(col);
322  gse->SetCommonSumLineColor(col);
324 
325  Int_t j = 0;
326  for (Int_t i = 1; i <= g->GetNbinsX(); i++) {
327  Double_t eta = g->GetXaxis()->GetBinCenter(i);
328  Double_t eEta = g->GetXaxis()->GetBinWidth(i)/2;
329  Double_t xo = TMath::Abs(eta)+eEta;
330  if (xo > 2) continue;
331  Double_t ea = 0.02*TMath::Power(xo/2,2);
332  gse->SetPoint(j, eta, g->GetBinContent(i));
333  gse->SetPointError(j, eEta, eEta);
334  gse->SetStatError(j, g->GetBinError(i),g->GetBinError(i));
335  j++;
336  }
337  return gse;
338 
339 
340 }
341 
342 //____________________________________________________________________
343 TObject*
344 GetO(TDirectory* dir, const char* name, TClass* cls=0)
345 {
346  if (!dir) {
347  Warning("GetO", "No directory passed");
348  return 0;
349  }
350 
351  TObject* o = dir->Get(name);
352  if (!o) {
353  Warning("GetO", "object %s not found in %s",
354  name, dir->GetPath());
355  return 0;
356  }
357  if (!cls) return o;
358  if (!o->IsA()->InheritsFrom(cls)) {
359  Warning("GetO", "Object %s in %s is not a %s, but a %s",
360  name, dir->GetPath(), cls->GetName(), o->ClassName());
361  return 0;
362  }
363  return o;
364 }
365 //____________________________________________________________________
366 TDirectory* GetD(TDirectory* dir, const char* name)
367 {
368  return static_cast<TDirectory*>(GetO(dir,name,TDirectory::Class()));
369 }
370 
371 //____________________________________________________________________
372 TH1* GetH1(TDirectory* dir, const char* name)
373 {
374  return static_cast<TH1*>(GetO(dir,name,TH1::Class()));
375 }
376 //____________________________________________________________________
381 void
382 ExtractGSE2(const char* input, UShort_t sNN=5023)
383 {
384  if (!gROOT->GetClass("GraphSysErr"))
385  gROOT->LoadMacro("$HOME/GraphSysErr/GraphSysErr.C+g");
386 
387  TString base = input; //gSystem->DirName(input); base.ReplaceAll(".root","");
388  TFile* file = TFile::Open(Form("%s/result.root",input), "READ");
389  if (!file) {
390  Warning("ExtractGSE2", "Failed to open %s/result.root", input);
391  return;
392  }
393 
394  Int_t dimen = -1;
395  if (base.EndsWith("unit")) dimen = 0;
396  else if (base.EndsWith("const")) dimen = 1;
397  else if (base.EndsWith("eta")) dimen = 2;
398  else if (base.EndsWith("etaipz")) dimen = 3;
399  if (dimen < 0) {
400  Error("ExtractGSE2", "Don't know how to extract dimension from %s",base);
401  return;
402  }
403 
404  TH1* frame = 0;
405  TH1* cent = GetH1(file, "realCent");
406 
407  Bool_t first = true;
408  TList* stack = new TList;
409  TList* truths = new TList;
410  for (Int_t i = 1; i <= cent->GetNbinsX(); i++) {
411  Double_t c1 = cent->GetXaxis()->GetBinLowEdge(i);
412  Double_t c2 = cent->GetXaxis()->GetBinUpEdge(i);
413  TObject* g = MakeGSE(file, dimen, sNN, c1, c2);
414  TObject* t = MakeTGSE(file, dimen, sNN, c1, c2);
415  Double_t gmin, gmax, tmin, tmax;
416  if (!g) continue;
417  stack->Add(g);
418  truths->Add(t);
419  if (first) g->Draw("quad stat combine axis");
420  else g->Draw("quad stat combine");
421  if (t) t->Draw("quad");
422  first = false;
423  if (!frame) {
424  GraphSysErr* gse = static_cast<GraphSysErr*>(g);
425  if (gse->GetMulti())
426  frame = gse->GetMulti()->GetHistogram();
427  }
428  }
429  if (frame) frame->SetMinimum(1);
430 
431  // TString obase(base); obase.Prepend("GSE_");
432 
433  std::ofstream out(Form("%s/gse.input", base.Data()));
434  GraphSysErr::Export(stack, out, "HFC", 2);
435  out << "*E" << std::endl;
436  out.close();
437 
438  TFile* rout = TFile::Open(Form("%s/gse.root", base.Data()), "RECREATE");
439  stack->AddAll(truths);
440  stack->Write("container", TObject::kSingleKey);
441  rout->Write();
442 
443 }
444 
445 void
447 {
448  ExtractGSE(Form("MiddNdeta_0x%x.root", flags));
449 }
450 //____________________________________________________________________
void Export(std::ostream &out=std::cout, Option_t *option="", Int_t nsign=2)
Definition: GraphSysErr.C:2612
TDirectory * GetD(TDirectory *dir, const char *name)
Definition: ExtractGSE2.C:366
void ExtractGSE2(const char *input, UShort_t sNN=5023)
Definition: ExtractGSE2.C:382
void SetSysLineColor(Int_t id, Color_t color)
Definition: GraphSysErr.C:4359
void SetYTitle(const char *title)
Definition: GraphSysErr.C:3598
const Color_t cc[]
Definition: DrawKs.C:1
double Double_t
Definition: External.C:58
TObject * GetO(TDirectory *dir, const char *name, TClass *cls=0)
Definition: ExtractGSE2.C:344
void SetKey(const char *key, const char *value, Bool_t replace=false)
Definition: GraphSysErr.C:4587
TObject * MakeGSE(TDirectory *d, Int_t dimen, UShort_t sNN, Double_t c1, Double_t c2)
Definition: ExtractGSE2.C:142
void SetPointError(Int_t i, Double_t ex)
Definition: GraphSysErr.C:3468
void SetSysError(Int_t id, Double_t eyl, Double_t eyh)
Definition: GraphSysErr.C:3528
void SetCommonSumLineColor(Color_t color)
Definition: GraphSysErr.C:4503
Double_t EtaSysEval(Double_t x, Double_t sMin, Double_t sMax)
Definition: ExtractGSE2.C:125
TCanvas * c
Definition: TestFitELoss.C:172
void AddQualifier(const TString &key, const TString &value, Bool_t replace=false)
Definition: GraphSysErr.C:4759
void SetCommonSumOption(EDrawOption_t opt)
Definition: GraphSysErr.C:4491
AliStack * stack
void SetStatError(Int_t i, Double_t ey)
Definition: GraphSysErr.C:3505
TH1 * GetH1(TDirectory *dir, const char *name)
Definition: ExtractGSE2.C:372
void SetPoint(Int_t i, Double_t x, Double_t y)
Definition: GraphSysErr.C:3457
void SetTitle(const char *name)
Definition: GraphSysErr.C:3576
UInt_t DefineCommon(const char *title, Bool_t relative, Double_t ey, EDrawOption_t option=kFill)
Definition: GraphSysErr.C:3350
Int_t MakeP2P(TObject *o, const char *name, Color_t c)
Definition: ExtractGSE2.C:106
int Int_t
Definition: External.C:63
UInt_t DeclarePoint2Point(const char *title, Bool_t relative, EDrawOption_t option=kBar)
Definition: GraphSysErr.C:3393
Double_t SysEval(Double_t x, Double_t sMin, Double_t sMax, Double_t xMax)
Definition: ExtractGSE2.C:116
Double_t CSysEval(Double_t x, Double_t sMin, Double_t sMax)
Definition: ExtractGSE2.C:121
void SetDataOption(EDrawOption_t opt)
Definition: GraphSysErr.C:3586
void SetSysFillColor(Int_t id, Color_t color)
Definition: GraphSysErr.C:4395
static Color_t PbPbColor(Double_t c1, Double_t c2)
Definition: ExtractGSE2.C:64
An (X,Y) graph with configurable errors.
TMultiGraph * GetMulti(Option_t *option="")
Definition: GraphSysErr.C:924
void SetSumOption(EDrawOption_t opt)
Definition: GraphSysErr.C:4449
void MakeCommon(TObject *o, const char *name, Double_t val, Color_t c)
Definition: ExtractGSE2.C:89
TObject * MakeTGSE(TDirectory *d, Int_t dimen, UShort_t sNN, Double_t c1, Double_t c2)
Definition: ExtractGSE2.C:275
TFile * file
TList with histograms for a given trigger.
unsigned short UShort_t
Definition: External.C:28
void SetSumFillColor(Color_t color)
Definition: GraphSysErr.C:4479
bool Bool_t
Definition: External.C:53
Double_t ptMax
void ExtractGSE(Int_t flags)
Definition: ExtractGSE2.C:446
void SetSumLineColor(Color_t color)
Definition: GraphSysErr.C:4461
static Int_t PbPbBin(Double_t c1, Double_t c2)
Definition: ExtractGSE2.C:39
void SetXTitle(const char *title)
Definition: GraphSysErr.C:3592
void SetCommonSumFillColor(Color_t color)
Definition: GraphSysErr.C:4521
Definition: External.C:196
TDirectoryFile * dir