AliPhysics  fceccc5 (fceccc5)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
ComputeDmesonYield.C
Go to the documentation of this file.
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include "TFile.h"
3 #include "TNtuple.h"
4 #include "TString.h"
5 #include "TGraphAsymmErrors.h"
6 #include "TH1.h"
7 #include "TH2.h"
8 #include "TLine.h"
9 #include "TMarker.h"
10 #include "TStyle.h"
11 #include "TLatex.h"
12 #include "TCanvas.h"
13 #include "TMath.h"
14 #include "TLegend.h"
15 #include "AliHFSystErr.h"
16 #endif
17 
18 // input files
19 TString filnamPPref="~/alice/Charm/PbYield/2011/Final/ppref/D0Kpi_276TeV_FONLLExtrapolationAndExtrapolation_from7TeVAlicedata_combinedFD_141211_010312_160712_rebin_1_2_3_4_5_6_8_12_16_24.root";
20 
21 TString filnamSpectrumNb="~/alice/Charm/PbYield/2011/Final/010/Dzero/HFPtSpectrum_D0toKpi_010_Nb_160615.root";
22 TString filnamSpectrumFc="~/alice/Charm/PbYield/2011/Final/010/Dzero/HFPtSpectrum_D0toKpi_010_fc_160615.root";
23 TString filnamRaaNb="~/alice/Charm/PbYield/2011/Final/010/Dzero/HFPtSpectrumRaa_D0toKpi_010_Nb_160615.root";
24 TString filnamRaaFc="~/alice/Charm/PbYield/2011/Final/010/Dzero/HFPtSpectrumRaa_D0toKpi_010_fc_160615.root";
25 
26 //configuration
27 enum dspec{kDzero,kDplus,kDstar,kDs};
28 Int_t mesonSpecie=kDzero;
29 TString centrality="0-10";
30 const Int_t nPtBins=9;
31 Double_t binlim[nPtBins+1]={1.,2.,3.,4.,5.,6.,8.,12.,16.,24.};
32 Int_t method=2; // 1=fc; 2=Nb --> default=2
33 Int_t optErrFD=2; // 0=from histos, not combined;
34  // 1= from ntuple with Rbc hypo, not combining Nb and fc;
35  // 2= from ntuple with Rbc hypo, combining Nb and fc
36  // --> default=2; change to 1 only to produce files for D-meson ratios
37 Float_t centHypoFdOverPr=2.;
38 Float_t lowHypoFdOverPr=1.;
39 Float_t highHypoFdOverPr=3.;
40 
41 
42 // Graphical styles
43 Bool_t draw[nPtBins]={1,1,0,1,0,1,0,0,1};
44 Int_t colors[nPtBins]={kGray+2,kMagenta+1,kMagenta,kBlue,kCyan,kGreen+2,kYellow+2,kOrange+1,kRed+1};
45 Int_t lstyle[nPtBins]={9,10,3,5,7,1,3,6,2};
46 
47 
48 Bool_t PbPbDataSyst(AliHFSystErr *syst, TH1D* heff, Double_t pt, Double_t &dataSystUp, Double_t &dataSystDown);
49 
50 
51 void ComputeDmesonYield(){
52 
53  TString mesName="Dzero";
54  Int_t mesCode=1;
55  Double_t brat=0.0388;
56  TString mesSymb="D^{0}";
57  Float_t ptForExtrap=16.;
58  if(mesonSpecie==kDplus){
59  mesName="Dplus";
60  mesCode=2;
61  brat=0.0913;
62  mesSymb="D^{+}";
63  ptForExtrap=24.;
64  }else if(mesonSpecie==kDstar){
65  mesName="Dstar";
66  mesCode=3;
67  brat=0.0388*0.677;
68  mesSymb="D^{*+}";
69  ptForExtrap=24.;
70  }else if(mesonSpecie==kDs){
71  mesName="Ds";
72  mesCode=4;
73  brat=0.0224;
74  mesSymb="D_{s}^{+}";
75  Float_t centHypoFdOverPr=1.;
76  Float_t lowHypoFdOverPr=1./3.;
77  Float_t highHypoFdOverPr=3.;
78  ptForExtrap=16.;
79  }
80 
81  TString centralityNoDash=centrality.Data();
82  centralityNoDash.ReplaceAll("-","");
83  TString filnamChi;
84  TString filnamCnt;
85  TString filnamCntS;
86  Int_t optCombineFDEloss=1; // 0=enevlope, 1= sum in quadrature
87  Bool_t correctForBR=kTRUE;
88  Bool_t showRcbSystNorm=kTRUE;
89 
90  if(method==1){
91  filnamChi=filnamSpectrumFc.Data();
92  filnamCnt=filnamRaaFc.Data();
93  filnamCntS=filnamRaaNb.Data();
94  }
95  else if(method==2){
96  filnamChi=filnamSpectrumNb.Data();
97  filnamCnt=filnamRaaNb.Data();
98  filnamCntS=filnamRaaFc.Data();
99  }
100 
101  Int_t colorSystFD=kGray+1;
102  Int_t colorSystRb=kOrange;
103  Int_t colorSystDatappC=kBlue+1;//kMagenta+1;
104  Int_t colorSystDataaaC=kRed+1;//kMagenta+1;
105  Int_t colorppC=kBlue+1;
106  Int_t coloraaC=kRed+1;
107  Int_t linestppC=1;
108  Int_t linestaaC=1;
109  Int_t linewippC=2;
110  Int_t linewiaaC=2;
111  Int_t markerppC=21;
112  Int_t markeraaC=23;
113  Int_t msizppC=1.2;
114  Int_t msizaaC=1.2;
115  Float_t sizesystdata=0.4;
116  Float_t sizesystfd=0.3;
117  Float_t sizesystrb=0.2;
118  Float_t sizesysttot=0.15;
119 
120 
121  // pp reference
122 
123  TFile *filPP=new TFile(filnamPPref.Data());
124  TH1D *hppRef = (TH1D*)filPP->Get("hReference");//fhScaledData");
125  TH1D *hppRefSystData = (TH1D*)filPP->Get("fhScaledSystData");
126 
127  Float_t relsysterrLowDatapp[nPtBins],relsysterrHiDatapp[nPtBins];
128  TGraphAsymmErrors* gppRefSyst=(TGraphAsymmErrors*)filPP->Get("gReferenceSyst");
129  for(Int_t i=0; i<gppRefSyst->GetN(); i++){
130  Double_t x,y;
131  gppRefSyst->GetPoint(i,x,y);
132  Int_t hBin=TMath::BinarySearch(nPtBins,binlim,x);
133  if(x>ptForExtrap){
134  relsysterrLowDatapp[hBin]=gppRefSyst->GetErrorYlow(i)/y;
135  relsysterrHiDatapp[hBin]=gppRefSyst->GetErrorYhigh(i)/y;
136  }
137  }
138  Float_t relstaterrpp[nPtBins];
139  for(Int_t ib=0; ib<nPtBins; ib++){
140  Int_t hBin=hppRef->FindBin(0.5*(binlim[ib]+binlim[ib+1]));
141  relstaterrpp[ib]=hppRef->GetBinError(hBin)/hppRef->GetBinContent(hBin);
142  printf("-- Relative stat errors PP: Bin %d(%f)--- histobin=%d Err %f\n",ib,0.5*(binlim[ib]+binlim[ib+1]),hBin,relstaterrpp[ib]);
143  Int_t hBin2=hppRefSystData->FindBin(0.5*(binlim[ib]+binlim[ib+1]));
144  if(hppRefSystData->GetBinCenter(hBin2)<ptForExtrap){
145  Float_t relsysterrDatapp=hppRefSystData->GetBinError(hBin2)/hppRefSystData->GetBinContent(hBin2);
146  relsysterrLowDatapp[ib]=relsysterrDatapp;
147  relsysterrHiDatapp[ib]=relsysterrDatapp;
148  }
149  printf(" ---- Check SYST err DATA PP Bin %d CS=%f Err+%f -%f\n",ib, hppRef->GetBinContent(hBin),relsysterrHiDatapp[ib],relsysterrLowDatapp[ib]);
150  }
151 
152  Float_t relsysterrLowEnScalpp[nPtBins];
153  Float_t relsysterrHiEnScalpp[nPtBins];
154  for(Int_t ib=0; ib<nPtBins; ib++){
155  relsysterrLowEnScalpp[ib]=0.;
156  relsysterrHiEnScalpp[ib]=0.;
157  }
158  TGraphAsymmErrors* gsystppEnSc=(TGraphAsymmErrors*)filPP->Get("gScaledDataSystExtrap");
159  for(Int_t i=0; i<gsystppEnSc->GetN(); i++){
160  Double_t x,y;
161  gsystppEnSc->GetPoint(i,x,y);
162  Int_t hBin=TMath::BinarySearch(nPtBins,binlim,x);
163  if(hBin>=0 && hBin<nPtBins && x<binlim[nPtBins]){
164  if(x<ptForExtrap){
165  relsysterrLowEnScalpp[hBin]=gsystppEnSc->GetErrorYlow(i)/y;
166  relsysterrHiEnScalpp[hBin]=gsystppEnSc->GetErrorYhigh(i)/y;
167  }
168  }
169  }
170 
171  TGraphAsymmErrors* gsystppFD=(TGraphAsymmErrors*)filPP->Get("gScaledDataSystFeedDown");
172  Float_t relsysterrLowFDpp[nPtBins];
173  Float_t relsysterrHiFDpp[nPtBins];
174  for(Int_t i=0; i<nPtBins; i++){
175  relsysterrHiFDpp[i]=0.;
176  relsysterrLowFDpp[i]=0.;
177  }
178  for(Int_t i=0; i<gsystppFD->GetN(); i++){
179  Double_t x,y;
180  gsystppFD->GetPoint(i,x,y);
181  Int_t hBin=TMath::BinarySearch(nPtBins,binlim,x);
182  if(hBin>=0 && hBin<nPtBins && x<binlim[nPtBins]){
183  relsysterrLowFDpp[hBin]=gsystppFD->GetErrorYlow(i)/y;
184  relsysterrHiFDpp[hBin]=gsystppFD->GetErrorYhigh(i)/y;
185  printf(" ---- Check syst err FD pp Bin %d Err+%f -%f\n",hBin,relsysterrHiFDpp[hBin],relsysterrLowFDpp[hBin]);
186  // printf("%d %f %f\n",hBin,relsysterrLowFDpp[hBin],relsysterrHiFDpp[hBin]);
187  }
188  }
189 
190  // A-A events
191 
192  printf("--- Pb-Pb events ---\n");
193 
194  TFile *filChi= new TFile(filnamChi.Data());
195  TH1D *hSpC = (TH1D*)filChi->Get("hRECpt");
196  Float_t relstaterrPbPb[nPtBins];
197  for(Int_t ib=0; ib<nPtBins; ib++){
198  Int_t hBin=hSpC->FindBin(0.5*(binlim[ib]+binlim[ib+1]));
199  relstaterrPbPb[ib]=hSpC->GetBinError(hBin)/hSpC->GetBinContent(hBin);
200  printf("-- Relative stat errors AA from yield: Bin %d(%d) pt=%f Err %f\n",ib,hBin,hSpC->GetBinCenter(hBin),relstaterrPbPb[ib]);
201  }
202 
203  AliHFSystErr *systematicsABcent = new AliHFSystErr();
204  systematicsABcent->SetRunNumber(11);
205  systematicsABcent->SetCollisionType(1);
206  systematicsABcent->SetCentrality(centralityNoDash.Data());
207  systematicsABcent->Init(mesCode);
208  TH1D* heffC=(TH1D*)filChi->Get("hDirectEffpt");
209 
210 
211  TGraphAsymmErrors* gsystaaFDc=(TGraphAsymmErrors*)filChi->Get("gSigmaCorrConservative"); // FD due to scales
212  Float_t relsysterrLowFDaa[nPtBins];
213  Float_t relsysterrHiFDaa[nPtBins];
214  for(Int_t i=0; i<nPtBins; i++){
215  relsysterrLowFDaa[i]=0.;
216  relsysterrHiFDaa[i]=0.;
217  }
218  for(Int_t i=0; i<gsystaaFDc->GetN(); i++){
219  Double_t x,y;
220  gsystaaFDc->GetPoint(i,x,y);
221  Int_t hBin=TMath::BinarySearch(nPtBins,binlim,x);
222  if(hBin>=0 && hBin<nPtBins && x<binlim[nPtBins]){
223  relsysterrLowFDaa[hBin]=gsystaaFDc->GetErrorYlow(i)/y;
224  relsysterrHiFDaa[hBin]=gsystaaFDc->GetErrorYhigh(i)/y;
225  printf(" ---- Check syst err FD AA Bin %d Err+%f -%f\n",hBin,relsysterrHiFDaa[hBin],relsysterrLowFDaa[hBin]);
226  }
227  }
228 
229  TFile *filCnt=new TFile(filnamCnt.Data());
230  Float_t relstaterrPbPb2[nPtBins];
231  TH1D* hraaCcheck2=(TH1D*)filCnt->Get("hRABvsPt");
232  for(Int_t ib=0; ib<nPtBins; ib++){
233  Int_t hBin=hraaCcheck2->FindBin(0.5*(binlim[ib]+binlim[ib+1]));
234  Double_t aux=hraaCcheck2->GetBinError(hBin)/hraaCcheck2->GetBinContent(hBin);
235  relstaterrPbPb2[ib]=TMath::Sqrt(aux*aux-relstaterrpp[ib]*relstaterrpp[ib]);
236 
237  printf("-- Relative stat errors AA from RAA-PP: Bin %d(%f)---%d Err %f\n",ib,0.5*(binlim[ib]+binlim[ib+1]),hBin,relstaterrPbPb2[ib]);
238  }
239 
240  TNtuple* ntC=(TNtuple*)filCnt->Get("ntupleRAB");
241  Float_t pt,TAB,sigmaPP,invyieldAB,RABCharm,RABBeauty;
242  Float_t invyieldABFDHigh,invyieldABFDLow,fprompt;
243  ntC->SetBranchAddress("pt",&pt);
244  ntC->SetBranchAddress("TAB",&TAB);
245  ntC->SetBranchAddress("sigmaPP",&sigmaPP);
246  ntC->SetBranchAddress("invyieldAB",&invyieldAB);
247  ntC->SetBranchAddress("RABCharm",&RABCharm);
248  ntC->SetBranchAddress("RABBeauty",&RABBeauty);
249  ntC->SetBranchAddress("invyieldABFDHigh",&invyieldABFDHigh);
250  ntC->SetBranchAddress("invyieldABFDLow",&invyieldABFDLow);
251  ntC->SetBranchAddress("fc",&fprompt);
252 
253  Float_t raac[nPtBins],invypp[nPtBins],invyPbPb[nPtBins],invyPbPbLo[nPtBins],invyPbPbHi[nPtBins],minval[nPtBins];
254  Float_t raacLowHyp[nPtBins],invyPbPbLowHyp[nPtBins],minvalLowHyp[nPtBins];
255  Float_t raacHigHyp[nPtBins],invyPbPbHigHyp[nPtBins],minvalHigHyp[nPtBins];
256  Float_t fPromptCent[nPtBins],fPromptHigHyp[nPtBins],fPromptLowHyp[nPtBins];
257  Float_t invyPbPbLoSingleSyst[nPtBins],invyPbPbHiSingleSyst[nPtBins];
258 
259  for(Int_t ib=0; ib<nPtBins; ib++){
260  minval[ib]=9999.;
261  invypp[ib]=-999.;
262  invyPbPb[ib]=-999.;
263  invyPbPbLo[ib]=-999.;
264  invyPbPbHi[ib]=-999.;
265  invyPbPbLoSingleSyst[ib]=-999.;
266  invyPbPbHiSingleSyst[ib]=-999.;
267  raac[ib]=-999.;
268  raacLowHyp[ib]=-999.;
269  invyPbPbLowHyp[ib]=-999.;
270  minvalLowHyp[ib]=9999.;
271  raacHigHyp[ib]=-999.;
272  invyPbPbHigHyp[ib]=-999.;
273  minvalHigHyp[ib]=9999.;
274  fPromptCent[ib]=-9999.;
275  fPromptHigHyp[ib]=-9999.;
276  fPromptLowHyp[ib]=-9999.;
277  }
278 
279  TGraph** gcb=new TGraph*[nPtBins];
280  TGraph** gcrbc=new TGraph*[nPtBins];
281  for(Int_t ib=0; ib<nPtBins; ib++){
282  gcb[ib]=new TGraph(0);
283  gcrbc[ib]=new TGraph(0);
284  }
285 
286  for(Int_t ie=0; ie<ntC->GetEntries(); ie++){
287  ntC->GetEvent(ie);
288  if(correctForBR){
289  invyieldAB/=brat;
290  invyieldABFDLow/=brat;
291  invyieldABFDHigh/=brat;
292  sigmaPP/=brat;
293  }
294  if(pt>binlim[nPtBins]) continue;
295  Int_t theBin=TMath::BinarySearch(nPtBins,binlim,(Double_t)pt);
296  if(theBin<0 || theBin>=nPtBins) continue;
297  Float_t rFdPr=RABBeauty/RABCharm;
298  Float_t dist=TMath::Abs(rFdPr-centHypoFdOverPr);
299  if(dist<minval[theBin]){
300  raac[theBin]=RABCharm;
301  minval[theBin]=dist;
302  invyPbPb[theBin]=invyieldAB;
303  invyPbPbLo[theBin]=invyieldABFDLow;
304  invyPbPbHi[theBin]=invyieldABFDHigh;
305  invypp[theBin]=TAB*sigmaPP;
306  fPromptCent[theBin]=fprompt;
307  }
308  Float_t distLowHyp=TMath::Abs(rFdPr-lowHypoFdOverPr);
309  if(distLowHyp<minvalLowHyp[theBin]){
310  // LowHyp -> lower Raa(feeddown) -> less Fd to be subtracted -> higher fprompt -> higher prompt yield -> higher Raa(prompt)
311  raacLowHyp[theBin]=RABCharm;
312  invyPbPbLowHyp[theBin]=invyieldAB;
313  invyPbPbHiSingleSyst[theBin]=invyieldABFDHigh;
314  minvalLowHyp[theBin]=distLowHyp;
315  fPromptLowHyp[theBin]=fprompt;
316  }
317  Float_t distHigHyp=TMath::Abs(rFdPr-highHypoFdOverPr);
318  if(distHigHyp<minvalHigHyp[theBin]){
319  // HigHyp -> higher Raa(feeddown) -> more Fd to be subtracted -> lower fprompt -> lower prompt yield -> lower Raa(prompt)
320  raacHigHyp[theBin]=RABCharm;
321  invyPbPbHigHyp[theBin]=invyieldAB;
322  invyPbPbLoSingleSyst[theBin]=invyieldABFDLow;
323  minvalHigHyp[theBin]=distHigHyp;
324  fPromptHigHyp[theBin]=fprompt;
325  }
326  if(theBin<nPtBins && theBin>=0 && pt<binlim[nPtBins]){
327  gcb[theBin]->SetPoint(gcb[theBin]->GetN(),RABBeauty,RABCharm);
328  gcrbc[theBin]->SetPoint(gcrbc[theBin]->GetN(),rFdPr,RABCharm);
329  }
330  }
331 
332  TH1D* hcheck[nPtBins];
333  for(Int_t i=0; i<nPtBins; i++){
334  hcheck[i]=(TH1D*)filCnt->Get(Form("hRCharmVsElossHypo_%d",i+1));
335  }
336 
337  TMarker** mC=new TMarker*[nPtBins];
338  for(Int_t ib=0; ib<nPtBins; ib++){
339  mC[ib]=new TMarker(1.,raac[ib],20);
340  mC[ib]->SetMarkerSize(1.2);
341  if(showRcbSystNorm){
342  Double_t auxx,auxy;
343  for(Int_t ip=0; ip<gcrbc[ib]->GetN(); ip++){
344  gcrbc[ib]->GetPoint(ip,auxx,auxy);
345  auxy/=raac[ib];
346  gcrbc[ib]->SetPoint(ip,auxx,(auxy-1)*100.);
347  }
348  }
349  }
350 
351  TFile *filCntS=new TFile(filnamCntS.Data());
352  TNtuple* ntCS=(TNtuple*)filCntS->Get("ntupleRAB");
353  ntCS->SetBranchAddress("pt",&pt);
354  ntCS->SetBranchAddress("TAB",&TAB);
355  ntCS->SetBranchAddress("sigmaPP",&sigmaPP);
356  ntCS->SetBranchAddress("invyieldAB",&invyieldAB);
357  ntCS->SetBranchAddress("RABCharm",&RABCharm);
358  ntCS->SetBranchAddress("RABBeauty",&RABBeauty);
359  ntCS->SetBranchAddress("invyieldABFDHigh",&invyieldABFDHigh);
360  ntCS->SetBranchAddress("invyieldABFDLow",&invyieldABFDLow);
361  ntCS->SetBranchAddress("fc",&fprompt);
362 
363  Float_t invyPbPbS[nPtBins],invyPbPbLoS[nPtBins],invyPbPbHiS[nPtBins],minvalS[nPtBins];
364  Float_t invyPbPbLoSingleSystS[nPtBins],invyPbPbHiSingleSystS[nPtBins],minvalLowHypS[nPtBins],minvalHigHypS[nPtBins];
365  Float_t fPromptLowHypS[nPtBins],fPromptHigHypS[nPtBins];
366  for(Int_t ib=0; ib<nPtBins; ib++){
367  minvalS[ib]=9999.;
368  minvalLowHypS[ib]=9999.;
369  minvalHigHypS[ib]=9999.;
370  invyPbPbS[ib]=-999.;
371  invyPbPbLoS[ib]=-999.;
372  invyPbPbHiS[ib]=-999.;
373  invyPbPbLoSingleSystS[ib]=-999.;
374  invyPbPbHiSingleSystS[ib]=-999.;
375  fPromptLowHypS[ib]=-999.;
376  fPromptHigHypS[ib]=-999.;
377  }
378  for(Int_t ie=0; ie<ntCS->GetEntries(); ie++){
379  ntCS->GetEvent(ie);
380  if(correctForBR){
381  invyieldAB/=brat;
382  invyieldABFDLow/=brat;
383  invyieldABFDHigh/=brat;
384  sigmaPP/=brat;
385  }
386  if(pt>binlim[nPtBins]) continue;
387  Int_t theBin=TMath::BinarySearch(nPtBins,binlim,(Double_t)pt);
388  if(theBin<0 || theBin>=nPtBins) continue;
389  Float_t rFdPr=RABBeauty/RABCharm;
390  Float_t dist=TMath::Abs(rFdPr-centHypoFdOverPr);
391  if(dist<minvalS[theBin]){
392  minvalS[theBin]=dist;
393  invyPbPbS[theBin]=invyieldAB;
394  invyPbPbLoS[theBin]=invyieldABFDLow;
395  invyPbPbHiS[theBin]=invyieldABFDHigh;
396  }
397  Float_t distLowHyp=TMath::Abs(rFdPr-lowHypoFdOverPr);
398  if(distLowHyp<minvalLowHypS[theBin]){
399  // LowHyp -> lower Raa(feeddown) -> less Fd to be subtracted -> higher fprompt -> higher prompt yield -> higher Raa(prompt)
400  invyPbPbHiSingleSystS[theBin]=invyieldABFDHigh;
401  minvalLowHypS[theBin]=distLowHyp;
402  fPromptLowHypS[theBin]=fprompt;
403  }
404  Float_t distHigHyp=TMath::Abs(rFdPr-highHypoFdOverPr);
405  if(distHigHyp<minvalHigHypS[theBin]){
406  // HigHyp -> higher Raa(feeddown) -> more Fd to be subtracted -> lower fprompt -> lower prompt yield -> lower Raa(prompt)
407  invyPbPbLoSingleSystS[theBin]=invyieldABFDLow;
408  minvalHigHypS[theBin]=distHigHyp;
409  fPromptHigHypS[theBin]=fprompt;
410  }
411  }
412 
413  TH1F* hfPromptCent=new TH1F("hfPromptCent"," ; p_{T} (Gev/c) ; f_{prompt}",nPtBins,binlim);
414  TH1F* hfPromptMinNb=new TH1F("hfPromptMinNb"," ; p_{T} (Gev/c) ; f_{prompt}",nPtBins,binlim);
415  TH1F* hfPromptMinfc=new TH1F("hfPromptMinfc"," ; p_{T} (Gev/c) ; f_{prompt}",nPtBins,binlim);
416  TH1F* hfPromptMaxNb=new TH1F("hfPromptMaxNb"," ; p_{T} (Gev/c) ; f_{prompt}",nPtBins,binlim);
417  TH1F* hfPromptMaxfc=new TH1F("hfPromptMaxfc"," ; p_{T} (Gev/c) ; f_{prompt}",nPtBins,binlim);
418  hfPromptCent->SetStats(0);
419  hfPromptMinNb->SetStats(0);
420  hfPromptMaxNb->SetStats(0);
421  hfPromptMinfc->SetStats(0);
422  hfPromptMaxfc->SetStats(0);
423 
424  for(Int_t ib=0; ib<nPtBins; ib++){
425  printf("Bin %d\n",ib);
426  printf(" fprompt central=%f --- Nb fpromptmin=%f fpromptmax=%f --- fc fpromptmin=%f fpromptmax=%f\n",fPromptCent[ib],fPromptHigHyp[ib],fPromptLowHyp[ib],fPromptHigHypS[ib],fPromptLowHypS[ib]);
427  hfPromptCent->SetBinContent(ib+1,fPromptCent[ib]);
428  hfPromptMinNb->SetBinContent(ib+1,fPromptHigHyp[ib]);
429  hfPromptMaxNb->SetBinContent(ib+1,fPromptLowHyp[ib]);
430  hfPromptMinfc->SetBinContent(ib+1,fPromptHigHypS[ib]);
431  hfPromptMaxfc->SetBinContent(ib+1,fPromptLowHypS[ib]);
432  }
433 
434  TH1F* hppC=new TH1F("hppC",Form("pp reference for %s%% CC",centrality.Data()),nPtBins,binlim);
435  TGraphAsymmErrors *gppCsystdata=new TGraphAsymmErrors(0);
436  gppCsystdata->SetName("gppCsystdata");
437  gppCsystdata->SetTitle(Form("Data Syst. Err. pp, scaled to %s%% CC",centrality.Data()));
438  TGraphAsymmErrors *gppCsystFD=new TGraphAsymmErrors(0);
439  gppCsystFD->SetName("gppCsystFD");
440  gppCsystFD->SetTitle(Form("B feed-down Syst. Err. pp, scaled to %s%% CC",centrality.Data()));
441  TH1F* hAAC=new TH1F("hAAC",Form("PbPb %s%% CC",centrality.Data()),nPtBins,binlim);
442  TGraphAsymmErrors *gaaCsystdata=new TGraphAsymmErrors(0);
443  gaaCsystdata->SetName("gaaCsystdata");
444  gaaCsystdata->SetTitle(Form("Data Syst. Err. PbPb %s%% CC",centrality.Data()));
445  TGraphAsymmErrors *gaaCsystFD=new TGraphAsymmErrors(0);
446  gaaCsystFD->SetName("gaaCsystFD");
447  gaaCsystFD->SetTitle(Form("B feed-down Syst. Err. PbPb %s%% CC",centrality.Data()));
448  TGraphAsymmErrors *gaaCsystRb=new TGraphAsymmErrors(0);
449  gaaCsystRb->SetName("gaaCsystRb");
450  gaaCsystRb->SetTitle(Form("Raa(B) Syst. Err. PbPb %s%% CC",centrality.Data()));
451  TGraphAsymmErrors *gaaCsystB=new TGraphAsymmErrors(0);
452  gaaCsystB->SetName("gaaCsystB");
453  gaaCsystB->SetTitle(Form("B Syst. Err. PbPb %s%% CC",centrality.Data()));
454  TGraphAsymmErrors *gaaCsystTot=new TGraphAsymmErrors(0);
455  gaaCsystTot->SetName("gaaCsystTot");
456  gaaCsystTot->SetTitle(Form("Tot Syst. Err. PbPb %s%% CC",centrality.Data()));
457  TH1F* hRAAC=new TH1F("hRAAC","",nPtBins,binlim);
458  TGraphAsymmErrors *graaC=new TGraphAsymmErrors(0);
459 
460  for(Int_t ib=0; ib<nPtBins; ib++){
461  printf("Bin %d\n",ib);
462  printf("raa(centHyp)=%f raa(lowHyp)=%f raa(highHyp)=%f RelDiff=%f\n",
463  raac[ib],raacLowHyp[ib],raacHigHyp[ib],(raacHigHyp[ib]-raacLowHyp[ib])/raac[ib]);
464  if(raac[ib]>0.){
465  Float_t relstaterrRaa=TMath::Sqrt(relstaterrpp[ib]*relstaterrpp[ib]+relstaterrPbPb2[ib]*relstaterrPbPb2[ib]);
466  hRAAC->SetBinContent(ib+1,raac[ib]);
467  hRAAC->SetBinError(ib+1,relstaterrRaa*raac[ib]);
468  Int_t nP=graaC->GetN();
469  graaC->SetPoint(nP,hRAAC->GetBinCenter(ib+1),raac[ib]);
470  graaC->SetPointEXlow(nP,hRAAC->GetBinCenter(ib+1)-hRAAC->GetBinLowEdge(ib+1));
471  graaC->SetPointEXhigh(nP,hRAAC->GetBinLowEdge(ib+2)-hRAAC->GetBinCenter(ib+1));
472  graaC->SetPointEYlow(nP,raac[ib]-raacHigHyp[ib]);
473  graaC->SetPointEYhigh(nP,raacLowHyp[ib]-raac[ib]);
474  }
475  if(invypp[ib]>0.){
476  hppC->SetBinContent(ib+1,invypp[ib]);
477  hppC->SetBinError(ib+1,relstaterrpp[ib]*invypp[ib]);
478  Int_t nP=gppCsystdata->GetN();
479  gppCsystdata->SetPoint(nP,hppC->GetBinCenter(ib+1),invypp[ib]);
480  gppCsystdata->SetPointEXlow(nP,sizesystdata);
481  gppCsystdata->SetPointEXhigh(nP,sizesystdata);
482  Double_t edatl=relsysterrLowDatapp[ib]*invypp[ib];
483  Double_t edath=relsysterrHiDatapp[ib]*invypp[ib];
484  Double_t eenscl=relsysterrLowEnScalpp[ib]*invypp[ib];
485  Double_t eensch=relsysterrHiEnScalpp[ib]*invypp[ib];
486  Double_t elow=TMath::Sqrt(edatl*edatl+eenscl*eenscl);
487  Double_t ehig=TMath::Sqrt(edath*edath+eensch*eensch);
488  gppCsystdata->SetPointEYlow(nP,elow);
489  gppCsystdata->SetPointEYhigh(nP,ehig);
490  nP=gppCsystFD->GetN();
491  gppCsystFD->SetPoint(nP,hppC->GetBinCenter(ib+1),invypp[ib]);
492  gppCsystFD->SetPointEXlow(nP,sizesystfd);
493  gppCsystFD->SetPointEXhigh(nP,sizesystfd);
494  gppCsystFD->SetPointEYlow(nP,relsysterrLowFDpp[ib]*invypp[ib]);
495  gppCsystFD->SetPointEYhigh(nP,relsysterrHiFDpp[ib]*invypp[ib]);
496  }
497  if(invyPbPb[ib]>0.){
498  printf("yielddAA(centHyp)=%f yielddAA(lowHyp)=%f yielddAA(highHyp)=%f RelDiff=%f\n",
499  invyPbPb[ib],invyPbPbLowHyp[ib],invyPbPbHigHyp[ib],(invyPbPbHigHyp[ib]-invyPbPbLowHyp[ib])/invyPbPb[ib]);
500  hAAC->SetBinContent(ib+1,invyPbPb[ib]);
501  hAAC->SetBinError(ib+1,relstaterrPbPb2[ib]*invyPbPb[ib]);
502  Int_t nP=gaaCsystdata->GetN();
503  gaaCsystdata->SetPoint(nP,hAAC->GetBinCenter(ib+1),invyPbPb[ib]);
504  gaaCsystdata->SetPointEXlow(nP,sizesystdata);
505  gaaCsystdata->SetPointEXhigh(nP,sizesystdata);
506  Double_t systundatalow,systundatahigh;
507  Bool_t isOk=PbPbDataSyst(systematicsABcent,heffC,hAAC->GetBinCenter(ib+1),systundatahigh,systundatalow);
508  printf("Check PID syst: %f %f %f\n",systematicsABcent->GetTotalSystErr(hAAC->GetBinCenter(ib+1)),systundatahigh,systundatalow);
509  systundatalow*=invyPbPb[ib];
510  systundatahigh*=invyPbPb[ib];
511  gaaCsystdata->SetPointEYlow(nP,systundatalow);
512  gaaCsystdata->SetPointEYhigh(nP,systundatahigh);
513 
514  nP=gaaCsystRb->GetN();
515  gaaCsystRb->SetPoint(nP,hAAC->GetBinCenter(ib+1),invyPbPb[ib]);
516  gaaCsystRb->SetPointEXlow(nP,sizesystrb);
517  gaaCsystRb->SetPointEXhigh(nP,sizesystrb);
518  gaaCsystRb->SetPointEYlow(nP,invyPbPb[ib]-invyPbPbHigHyp[ib]);
519  gaaCsystRb->SetPointEYhigh(nP,invyPbPbLowHyp[ib]-invyPbPb[ib]);
520  nP=gaaCsystFD->GetN();
521  gaaCsystFD->SetPoint(nP,hAAC->GetBinCenter(ib+1),invyPbPb[ib]);
522  gaaCsystFD->SetPointEXlow(nP,sizesystfd);
523  gaaCsystFD->SetPointEXhigh(nP,sizesystfd);
524  gaaCsystB->SetPoint(nP,hAAC->GetBinCenter(ib+1),invyPbPb[ib]);
525  gaaCsystB->SetPointEXlow(nP,sizesystfd);
526  gaaCsystB->SetPointEXhigh(nP,sizesystfd);
527  gaaCsystTot->SetPoint(nP,hAAC->GetBinCenter(ib+1),invyPbPb[ib]);
528  gaaCsystTot->SetPointEXlow(nP,sizesysttot);
529  gaaCsystTot->SetPointEXhigh(nP,sizesysttot);
530  if(optErrFD==0){
531  gaaCsystFD->SetPointEYlow(nP,relsysterrLowFDaa[ib]*invyPbPb[ib]);
532  gaaCsystFD->SetPointEYhigh(nP,relsysterrHiFDaa[ib]*invyPbPb[ib]);
533  }else{
534  Double_t minyield,maxyield;
535  Double_t minyieldsing,maxyieldsing;
536  if(optErrFD==1){
537  minyield=invyPbPbLo[ib];
538  maxyield=invyPbPbHi[ib];
539  minyieldsing=invyPbPbLoSingleSyst[ib];
540  maxyieldsing=invyPbPbHiSingleSyst[ib];
541  }else{
542  minyield=TMath::Min(invyPbPbLo[ib],invyPbPbLoS[ib]);
543  maxyield=TMath::Max(invyPbPbHi[ib],invyPbPbHiS[ib]);
544  minyieldsing=TMath::Min(invyPbPbLoSingleSyst[ib],invyPbPbLoSingleSystS[ib]);
545  maxyieldsing=TMath::Max(invyPbPbHiSingleSyst[ib],invyPbPbHiSingleSystS[ib]);
546  }
547  gaaCsystFD->SetPointEYlow(nP,invyPbPb[ib]-minyield);
548  gaaCsystFD->SetPointEYhigh(nP,maxyield-invyPbPb[ib]);
549  printf("Relative syst from FD (FONLL scales+masses)=-%f +%f\n",gaaCsystFD->GetErrorYlow(nP)/invyPbPb[ib],gaaCsystFD->GetErrorYhigh(nP)/invyPbPb[ib]);
550  Double_t systunBlow=0.;
551  Double_t systunBhigh=0.;
552  if(optCombineFDEloss==0){ // envelope
553  systunBlow=invyPbPb[ib]-minyieldsing;
554  systunBhigh=maxyieldsing-invyPbPb[ib];
555  }else{ // sum in quadrature
556  Double_t e1l=gaaCsystFD->GetErrorYlow(nP);
557  Double_t e1h=gaaCsystFD->GetErrorYhigh(nP);
558  Double_t e2l=gaaCsystRb->GetErrorYlow(nP);
559  Double_t e2h=gaaCsystRb->GetErrorYhigh(nP);
560  systunBlow=TMath::Sqrt(e1l*e1l+e2l*e2l);
561  systunBhigh=TMath::Sqrt(e1h*e1h+e2h*e2h);
562  }
563  gaaCsystB->SetPointEYlow(nP,systunBlow);
564  gaaCsystB->SetPointEYhigh(nP,systunBhigh);
565 
566  Double_t totSystlow=TMath::Sqrt(systunBlow*systunBlow+systundatalow*systundatalow);
567  Double_t totSysthigh=TMath::Sqrt(systunBhigh*systunBhigh+systundatahigh*systundatahigh);
568  gaaCsystTot->SetPointEYlow(nP,totSystlow);
569  gaaCsystTot->SetPointEYhigh(nP,totSysthigh);
570  }
571  }
572  }
573 
574  TH1F* hcheckRAAC=(TH1F*)hAAC->Clone("hcheckRAAC");
575  hcheckRAAC->Divide(hppC);
576 
577 
578  // Plots
579 
580  gStyle->SetOptTitle(0);
581 
582  TCanvas* c1=new TCanvas("c1","RcVsRb");
583  TLegend* leg= new TLegend(0.19,0.18,0.55,0.41);
584  leg->SetFillColor(0);
585  leg->SetFillStyle(0);
586  leg->SetBorderSize(0);
587  leg->SetTextFont(43);
588  leg->SetTextSize(28);
589  leg->SetMargin(0.32);
590  TLegendEntry* ent;
591  Bool_t first=kTRUE;
592 
593  for(Int_t ib=0; ib<nPtBins; ib++){
594  if(draw[ib] && gcb[ib]->GetN()>0){
595  gcb[ib]->SetLineColor(colors[ib]);
596  gcb[ib]->SetLineWidth(3);
597  gcb[ib]->SetLineStyle(lstyle[ib]);
598  if(first){
599  gcb[ib]->Draw("AL");
600  gcb[ib]->GetXaxis()->SetLimits(0.,1.);
601  gcb[ib]->GetXaxis()->SetTitle("#it{R}_{AA} feed-down");
602  gcb[ib]->GetYaxis()->SetTitle("#it{R}_{AA} prompt");
603  gcb[ib]->SetMinimum(0.);
604  gcb[ib]->SetMaximum(1.);
605  first=kFALSE;
606  }else{
607  gcb[ib]->Draw("lsame");
608  }
609  ent=leg->AddEntry(gcb[ib],Form("%d<#it{p}_{T}<%d GeV/#it{c}",(Int_t)binlim[ib],(Int_t)binlim[ib+1]),"L");
610  }
611  }
612  leg->Draw();
613 
614  first=kTRUE;
615 
616  TCanvas* c2=new TCanvas("c2x","RcVsRcb",700,700);
617  c2->SetTickx();
618  c2->SetTicky();
619  // c2->SetGridx();
620  // c2->SetGridy();
621  c2->SetLeftMargin(0.14);
622  c2->SetBottomMargin(0.13);
623  c2->SetTopMargin(0.035);
624  c2->SetRightMargin(0.045);
625  for(Int_t ib=0; ib<nPtBins; ib++){
626  if(draw[ib] && gcrbc[ib]->GetN()>0){
627  gcrbc[ib]->SetLineColor(colors[ib]);
628  gcrbc[ib]->SetLineWidth(3);
629  gcrbc[ib]->SetLineStyle(lstyle[ib]);
630  if(first){
631  gcrbc[ib]->Draw("AL");
632  gcrbc[ib]->GetXaxis()->SetLimits(lowHypoFdOverPr,highHypoFdOverPr);
633  gcrbc[ib]->GetXaxis()->SetTitle("Hypothesis on (#it{R}_{AA} feed-down)/(#it{R}_{AA} prompt)");
634  gcrbc[ib]->GetYaxis()->SetTitleOffset(1.2);
635  gcrbc[ib]->GetXaxis()->SetTitleOffset(1.2);
636  gcrbc[ib]->GetYaxis()->SetTitleFont(43);
637  gcrbc[ib]->GetXaxis()->SetTitleFont(43);
638  gcrbc[ib]->GetYaxis()->SetTitleSize(30);
639  gcrbc[ib]->GetXaxis()->SetTitleSize(30);
640  gcrbc[ib]->GetYaxis()->SetLabelFont(43);
641  gcrbc[ib]->GetXaxis()->SetLabelFont(43);
642  gcrbc[ib]->GetYaxis()->SetLabelSize(28);
643  gcrbc[ib]->GetXaxis()->SetLabelSize(28);
644  if(!showRcbSystNorm){
645  gcrbc[ib]->GetYaxis()->SetTitle("#it{R}_{AA} prompt");
646  gcrbc[ib]->SetMinimum(0.);
647  gcrbc[ib]->SetMaximum(0.8);
648  }else{
649  gcrbc[ib]->GetYaxis()->SetTitle("Relative variation of #it{R}_{AA} prompt (%)");
650  gcrbc[ib]->SetMinimum(-20);
651  gcrbc[ib]->SetMaximum(20);
652  }
653  first=kFALSE;
654  }else{
655  gcrbc[ib]->Draw("lsame");
656  }
657  if(!showRcbSystNorm){
658  mC[ib]->SetMarkerColor(colors[ib]);
659  mC[ib]->Draw("same");
660  }
661  }
662  }
663 // TBox* b=new TBox(0.5,0.,2.,1.);
664 // b->SetFillColor(kYellow);
665 // b->SetFillStyle(3003);
666  // b->Draw("same");
667  leg->Draw();
668  TLatex* tali0=new TLatex(0.18,0.89,"ALICE");
669  tali0->SetNDC();
670  tali0->SetTextFont(43);
671  tali0->SetTextSize(28);
672  tali0->Draw();
673  TLatex* tpbpb=new TLatex(0.41,0.89,Form("Pb-Pb, %s%%, #sqrt{#it{s}_{NN}} = 2.76 TeV",centrality.Data()));
674  tpbpb->SetNDC();
675  tpbpb->SetTextFont(43);
676  tpbpb->SetTextSize(28);
677  tpbpb->Draw();
678  // TLatex* tdmes=new TLatex(0.63,0.85,Form("%s meson",mesSymb.Data()));
679  TLatex* tdmes=new TLatex(0.3,0.43,Form("%s meson",mesSymb.Data()));
680  tdmes->SetNDC();
681  tdmes->SetTextFont(43);
682  tdmes->SetTextSize(28);
683  tdmes->Draw();
684 
685  TLine* lin0=new TLine(lowHypoFdOverPr,0.,highHypoFdOverPr,0.);
686  lin0->SetLineStyle(2);
687  lin0->SetLineColor(kGray+1);
688  lin0->Draw();
689  if(method==2 && optErrFD==2){
690  c2->SaveAs(Form("%s-RcVsRcb_method%d_optErrFD%d_br%d.eps",mesName.Data(),method,optErrFD,correctForBR));
691  }
692 
693 
694  TCanvas* cfp=new TCanvas("cfp","fprompt",700,700);
695  gPad->SetTickx();
696  gPad->SetTicky();
697  hfPromptCent->SetMinimum(0);
698  hfPromptCent->SetMaximum(1.05);
699  hfPromptCent->SetLineWidth(3);
700  hfPromptCent->Draw();
701  hfPromptMinNb->SetLineWidth(2);
702  hfPromptMaxNb->SetLineWidth(2);
703  hfPromptMinNb->SetLineStyle(2);
704  hfPromptMaxNb->SetLineStyle(2);
705  hfPromptMinNb->Draw("same");
706  hfPromptMaxNb->Draw("same");
707  hfPromptMinfc->SetLineWidth(2);
708  hfPromptMaxfc->SetLineWidth(2);
709  hfPromptMinfc->SetLineStyle(2);
710  hfPromptMaxfc->SetLineStyle(2);
711  hfPromptMinfc->SetLineColor(2);
712  hfPromptMaxfc->SetLineColor(2);
713  hfPromptMinfc->Draw("same");
714  hfPromptMaxfc->Draw("same");
715  TLatex* t2=new TLatex(0.17,0.2,"Nb");
716  t2->SetNDC();
717  t2->Draw();
718  TLatex* t1=new TLatex(0.17,0.15,"fc");
719  t1->SetNDC();
720  t1->SetTextColor(2);
721  t1->Draw();
722 
723  hppC->SetMarkerStyle(markerppC);
724  hppC->SetMarkerSize(msizppC);
725  hppC->SetMarkerColor(colorppC);
726  hppC->SetLineColor(colorppC);
727  hppC->SetLineWidth(linewippC);
728  hppC->SetLineStyle(linestppC);
729  gppCsystdata->SetLineColor(colorSystDatappC);
730  gppCsystdata->SetLineWidth(2);
731  gppCsystdata->SetFillStyle(0);
732  gppCsystdata->SetFillColor(colorSystDatappC);
733  gppCsystFD->SetLineColor(colorSystFD);
734  gppCsystFD->SetFillColor(colorSystFD);
735  hAAC->SetMarkerStyle(markeraaC);
736  hAAC->SetMarkerSize(msizaaC);
737  hAAC->SetMarkerColor(coloraaC);
738  hAAC->SetLineColor(coloraaC);
739  hAAC->SetLineWidth(linewiaaC);
740  hAAC->SetLineStyle(linestaaC);
741  gaaCsystdata->SetLineColor(colorSystDataaaC);
742  gaaCsystdata->SetLineWidth(2);
743  gaaCsystdata->SetFillStyle(0);
744 
745  gaaCsystB->SetFillColor(colorSystRb);
746  gaaCsystRb->SetFillColor(colorSystRb);
747  gaaCsystFD->SetFillColor(colorSystFD);
748 
749  gaaCsystTot->SetLineColor(1);
750  gaaCsystTot->SetLineWidth(3);
751  gaaCsystTot->SetFillStyle(0);
752 
753 
754  Float_t ymax=4.7*hppC->GetMaximum();
755  Float_t ymin=999.;
756  Float_t y16C=hAAC->GetBinContent(hAAC->FindBin(14.));
757  if(y16C>0. && y16C<ymin) ymin=y16C;
758  Float_t y12C=hAAC->GetBinContent(hAAC->FindBin(10.));
759  if(y12C>0. && y12C<ymin) ymin=y12C;
760  ymin*=0.255;
761 
762  ymax=19.6;
763  ymin=7E-7;
764 
765  TH2F *hempty=new TH2F("hempty","",100,0.,38.998,100,ymin,ymax);
766  hempty->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
767  hempty->GetYaxis()->SetTitle("d#it{N}/d#it{p}_{T}_{ }|_{ |#it{y}|<0.5} (1/GeV/#it{c})");
768  hempty->GetYaxis()->SetTitleSize(0.05);
769  hempty->GetXaxis()->SetTitleSize(0.05);
770  hempty->GetYaxis()->SetTitleOffset(1.3);
771  hempty->GetYaxis()->SetTitleFont(42);
772  hempty->GetXaxis()->SetTitleFont(42);
773  hempty->GetYaxis()->SetLabelFont(42);
774  hempty->GetXaxis()->SetLabelFont(42);
775  hempty->SetStats(0);
776  TLatex* tdec =new TLatex(0.18,0.89,Form("ALICE %s #sqrt{#it{s}_{NN}} = 2.76 TeV",mesSymb.Data()));
777  tdec->SetNDC();
778  tdec->SetTextSize(0.038);
779  tdec->SetTextFont(42);
780 
781  gStyle->SetMarkerSize(1.8);
782 
783  TCanvas* c3=new TCanvas("c3","Yield1Pad",750,800);
784  c3->SetLeftMargin(0.15);
785  c3->SetRightMargin(0.07);
786  c3->SetBottomMargin(0.12);
787  c3->SetTopMargin(0.07);
788  c3->SetLogy();
789  c3->SetTicky();
790  hempty->Draw();
791  gppCsystFD->Draw("e2same");
792  gppCsystdata->DrawClone("e2same");
793  gaaCsystFD->Draw("E2same");
794  gaaCsystRb->Draw("E2same");
795  gaaCsystdata->DrawClone("e2same");
796  hppC->DrawClone("Psame");
797  hAAC->DrawClone("Psame");
798  // TLegend* legC= new TLegend(0.59,0.77,0.89,0.91);
799  TLegend* legC= new TLegend(0.53,0.55,0.89,0.69);
800  //TLegend* legC= new TLegend(0.59,0.37,0.89,0.51);
801  legC->SetHeader(Form("Centrality %s",centrality.Data()));
802  legC->SetFillColor(0);
803  legC->SetTextFont(42);
804  ent=legC->AddEntry(hppC,"p-p rescaled reference","PL");
805  ent=legC->AddEntry("","(#pm 6% norm. unc. not shown)","");
806  //ent->SetTextColor(hppC->GetLineColor());
807  ent=legC->AddEntry(hAAC,"Pb-Pb","PL");
808  //ent->SetTextColor(hAAC->GetLineColor());
809  legC->Draw();
810  //TLegend* legSy= new TLegend(0.18,0.16,0.45,0.32);
811  //TLegend* legSy= new TLegend(0.62,0.57,0.89,0.75);
812  TLegend* legSy= new TLegend(0.566,0.715,0.89,0.86);
813  legSy->SetFillColor(0);
814  legSy->SetLineColor(kGray+2);
815  legSy->SetTextFont(42);
816  ent=legSy->AddEntry(gppCsystdata,"Syst. unc. from Data","F");
817  ent=legSy->AddEntry(gppCsystFD,"Syst. unc. from FONLL feed-down corr.","F");
818  ent=legSy->AddEntry(gaaCsystRb,"Syst. unc. from #it{R}_{AA} feed-down","F");
819  legSy->Draw();
820  tdec->Draw();
821  c3->SaveAs(Form("%s-Yields_1pad_method%d_optErrFD%d_br%d.eps",mesName.Data(),method,optErrFD,correctForBR));
822 
823 
824 
825  TH2F *hemptyr=new TH2F("hemptyr","",100,0.,38.998,100,0.,1.2);
826  hemptyr->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
827  hemptyr->GetYaxis()->SetTitle(Form("#it{R}_{AA} (prompt %s)",mesSymb.Data()));
828  hemptyr->GetYaxis()->SetTitleOffset(1.2);
829  hemptyr->SetStats(0);
830 
831  TCanvas* c5=new TCanvas("c5","RAAcheck");
832  c5->SetGridy();
833  hRAAC->SetLineWidth(2);
834  hRAAC->SetMarkerStyle(20);
835  hRAAC->SetMarkerColor(1);
836  hRAAC->SetLineColor(1);
837  hRAAC->SetMinimum(0.);
838  hRAAC->SetMaximum(1.2);
839  hcheckRAAC->SetMarkerStyle(24);
840  hcheckRAAC->SetMarkerSize(1.2);
841  hcheckRAAC->SetMarkerColor(2);
842  hraaCcheck2->SetMarkerStyle(25);
843  hraaCcheck2->SetMarkerColor(6);
844  hraaCcheck2->SetMarkerSize(1.5);
845 
846  graaC->SetFillColor(kGray);
847  hemptyr->Draw();
848  graaC->Draw("E2same");
849  hRAAC->Draw("PEsame");
850  hraaCcheck2->Draw("PSAME");
851  hcheckRAAC->Draw("PSAME");
852  TLegend* legr=new TLegend(0.7,0.7,0.89,0.89);
853  ent=legr->AddEntry(hRAAC,Form("%s%%",centrality.Data()),"PL");
854  legr->SetFillColor(0);
855  legr->SetTextFont(42);
856  legr->Draw();
857  c5->Update();
858  // c5->SaveAs(Form("%s-RAA_check_method%d.eps",mesName.Data(),method));
859 
860 
861  // TCanvas* c6=new TCanvas("c6","Rccheck");
862  // c6->Divide(3,2);
863  // for(Int_t ib=0; ib<nPtBins; ib++){
864  // c6->cd(ib+1);
865  // if(gcrbc[ib]->GetN()>0){
866  // hcheck[ib]->Draw("col");
867  // gcrbc[ib]->SetLineColor(colors[ib]);
868  // gcrbc[ib]->Draw("lsame");
869  // }
870  // }
871 
872 
873  TFile* outfil=new TFile(Form("%sYield_method%d_fd%d_br%d.root",mesName.Data(),method,optErrFD,correctForBR),"recreate");
874  hAAC->Write();
875  gppCsystFD->Write();
876  gppCsystdata->Write();
877  gaaCsystFD->Write();
878  gaaCsystRb->Write();
879  gaaCsystB->Write();
880  gaaCsystdata->Write();
881  gaaCsystTot->Write();
882  hppC->Write();
883  hfPromptCent->Write();
884  hfPromptMinNb->Write();
885  hfPromptMaxNb->Write();
886  hfPromptMinfc->Write();
887  hfPromptMaxfc->Write();
888  outfil->Close();
889 }
890 
891 
892 //____________________________________________________________
893 Bool_t PbPbDataSyst(AliHFSystErr *syst, TH1D* heff, Double_t pt, Double_t &dataSystUp, Double_t &dataSystDown)
894 {
895 
896  Double_t err = syst->GetTotalSystErr(pt)*syst->GetTotalSystErr(pt);
897  Int_t theBin=heff->FindBin(pt);
898  Double_t errRel=heff->GetBinError(theBin)/heff->GetBinContent(theBin);
899 
900  err += (errRel*errRel);
901 
902  Double_t errDown = err ;
903  Double_t errUp = err ;
904 
905 
906  dataSystUp = TMath::Sqrt(errUp);
907  dataSystDown = TMath::Sqrt(errDown);
908 
909  printf("Pt %f Bin %d Eff %f RelErrEff %f TotSyst +%f -%f\n",pt,theBin,heff->GetBinContent(theBin),errRel,dataSystUp,dataSystDown);
910 
911  return kTRUE;
912 }
913 
914 
kDstar
Definition: ComputeDmesonYield.C:27
ymax
const Double_t ymax
Definition: AddTaskCFDStar.C:7
Double_t
double Double_t
Definition: External.C:58
TH2F
Definition: External.C:236
filnamSpectrumFc
TString filnamSpectrumFc
Definition: ComputeDmesonYield.C:22
mesonSpecie
Int_t mesonSpecie
Definition: ComputeDmesonYield.C:28
filnamSpectrumNb
TString filnamSpectrumNb
Definition: ComputeDmesonYield.C:21
filnamRaaFc
TString filnamRaaFc
Definition: ComputeDmesonYield.C:24
centrality
centrality
Definition: AverageDmesonRaa.C:71
AliHFSystErr::SetCentrality
void SetCentrality(TString centrality)
Definition: AliHFSystErr.h:60
kDs
Definition: ComputeDmesonYield.C:27
lowHypoFdOverPr
Float_t lowHypoFdOverPr
Definition: ComputeDmesonYield.C:38
kDzero
Definition: ComputeDmesonYield.C:27
nPtBins
const Int_t nPtBins
Definition: ComputeDmesonYield.C:30
Int_t
int Int_t
Definition: External.C:63
kDplus
Definition: ComputeDmesonYield.C:27
Float_t
float Float_t
Definition: External.C:68
method
Int_t method
Definition: ComputeDmesonYield.C:32
TH1D
Definition: External.C:212
ComputeDmesonYield
void ComputeDmesonYield()
Definition: ComputeDmesonYield.C:51
AliHFSystErr::GetTotalSystErr
Double_t GetTotalSystErr(Double_t pt, Double_t feeddownErr=0) const
Definition: AliHFSystErr.cxx:5669
colors
Int_t colors[nPtBins]
Definition: ComputeDmesonYield.C:44
binlim
Double_t binlim[nPtBins+1]
Definition: ComputeDmesonYield.C:31
AliHFSystErr::Init
void Init(Int_t decay)
Function to initialize the variables/histograms.
Definition: AliHFSystErr.cxx:97
draw
Bool_t draw[nPtBins]
Definition: ComputeDmesonYield.C:43
filnamPPref
TString filnamPPref
Definition: ComputeDmesonYield.C:19
AliHFSystErr.h
ymin
const Double_t ymin
Definition: AddTaskCFDStar.C:6
AliHFSystErr::SetCollisionType
void SetCollisionType(Int_t type)
Definition: AliHFSystErr.h:51
centHypoFdOverPr
Float_t centHypoFdOverPr
Definition: ComputeDmesonYield.C:37
highHypoFdOverPr
Float_t highHypoFdOverPr
Definition: ComputeDmesonYield.C:39
Bool_t
bool Bool_t
Definition: External.C:53
optErrFD
Int_t optErrFD
Definition: ComputeDmesonYield.C:33
dspec
dspec
Definition: ComputeDmesonYield.C:27
AliHFSystErr::SetRunNumber
void SetRunNumber(Int_t number)
Definition: AliHFSystErr.h:44
lstyle
Int_t lstyle[nPtBins]
Definition: ComputeDmesonYield.C:45
PbPbDataSyst
Bool_t PbPbDataSyst(AliHFSystErr *syst, TH1D *heff, Double_t pt, Double_t &dataSystUp, Double_t &dataSystDown)
Definition: ComputeDmesonYield.C:893
filnamRaaNb
TString filnamRaaNb
Definition: ComputeDmesonYield.C:23