]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGHF/vertexingHF/AliHFSystErr.cxx
updates for train
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliHFSystErr.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 /////////////////////////////////////////////////////////////
19 //
20 // Class to handle systematic errors for charm hadrons
21 //
22 // Usage:
23 // AliHFSystEff syst;           // DECAY = 1 for D0, 2, for D+, 3 for D* 
24 // syst.SetRunNumber(YEAR);     // YEAR = two last numbers of the year (is 10 for 2010)
25 // syst.SetCollisionType(TYPE);  // TYPE =  0 is pp, 1 is PbPb
26 // syst.SetCentrality(CENT);     // CENT is centrality, 0100 for MB, 020 (4080) for 0-20 (40-80) CC...
27 // syst.Init(DECAY);             // DECAY = 1 for D0, 2, for D+, 3 for D* 
28 // syst.DrawErrors(); // to see a plot of the error contributions
29 // syst.GetTotalSystErr(pt); // to get the total err at pt 
30 //
31 // Author: A.Dainese, andrea.dainese@pd.infn.it
32 /////////////////////////////////////////////////////////////
33
34 #include <TStyle.h>
35 #include <TGraphAsymmErrors.h>
36 #include <TMath.h>
37 #include <TCanvas.h>
38 #include <TH2F.h>
39 #include <TLegend.h>
40 #include <TColor.h>
41
42 #include "AliLog.h"
43 #include "AliHFSystErr.h"
44
45
46 ClassImp(AliHFSystErr)
47
48 //--------------------------------------------------------------------------
49 AliHFSystErr::AliHFSystErr(const Char_t* name, const Char_t* title) : 
50 TNamed(name,title),
51 fNorm(0),
52 fRawYield(0),
53 fTrackingEff(0),
54 fBR(0),
55 fCutsEff(0),
56 fPIDEff(0),
57 fMCPtShape(0),
58 fPartAntipart(0),
59 fRunNumber(10),
60 fCollisionType(0),
61 fCentralityClass("0100"),
62 fIsLowEnergy(false),
63 fIsCentScan(false)
64 {
65   //
66   // Default Constructor
67   //
68 }
69
70 //--------------------------------------------------------------------------
71 AliHFSystErr::~AliHFSystErr() {
72   //  
73   // Default Destructor
74   //
75   /*
76
77   if(fNorm)         { delete fNorm; fNorm=0; }
78   if(fRawYield)     { delete fRawYield; fRawYield=0; }
79   if(fTrackingEff)  { delete fTrackingEff; fTrackingEff=0; }
80   if(fBR)           { delete fBR; fBR=0; }
81   if(fCutsEff)      { delete fCutsEff; fCutsEff=0; }
82   if(fPIDEff)       { delete fPIDEff; fPIDEff=0; }
83   if(fMCPtShape)    { delete fMCPtShape; fMCPtShape=0; }
84   if(fPartAntipart) { delete fPartAntipart; fPartAntipart=0; }
85   */
86 }
87
88 //--------------------------------------------------------------------------
89 void AliHFSystErr::Init(Int_t decay){
90   //
91   // Variables/histos initialization
92   //
93
94   if (fRunNumber!=10 && fIsLowEnergy==false) { 
95     AliFatal("Only settings for 2010 and the low energy runs are implemented so far");
96   }
97
98   switch(decay) {
99   case 1: // D0->Kpi
100     if (fCollisionType==0) {
101       if (fIsLowEnergy) InitD0toKpi2010ppLowEn();
102       else InitD0toKpi2010pp();
103     } else if (fCollisionType==1) {
104       if (fCentralityClass=="010") InitD0toKpi2010PbPb010CentScan();
105       else if (fCentralityClass=="1020") InitD0toKpi2010PbPb1020CentScan();
106       else if (fCentralityClass=="020") InitD0toKpi2010PbPb020();
107       else if (fCentralityClass=="2040") InitD0toKpi2010PbPb2040CentScan();
108       else if (fCentralityClass=="4060") InitD0toKpi2010PbPb4060CentScan();
109       else if (fCentralityClass=="6080") InitD0toKpi2010PbPb6080CentScan();
110       else if (fCentralityClass=="4080") InitD0toKpi2010PbPb4080();
111       else AliFatal("Not yet implemented");
112     }
113     //    else if (fCollisionType==2) InitD0toKpi2010ppLowEn();
114     break;
115   case 2: // D+->Kpipi
116     if (fCollisionType==0) {
117       if (fIsLowEnergy) InitDplustoKpipi2010ppLowEn();
118       else InitDplustoKpipi2010pp();
119     } else if (fCollisionType==1) {
120       if (fCentralityClass=="010") InitDplustoKpipi2010PbPb010CentScan();
121       else if (fCentralityClass=="1020") InitDplustoKpipi2010PbPb1020CentScan();
122       else if (fCentralityClass=="020") InitDplustoKpipi2010PbPb020();
123       else if (fCentralityClass=="2040") InitDplustoKpipi2010PbPb2040CentScan();
124       else if (fCentralityClass=="4060") InitDplustoKpipi2010PbPb4060CentScan();
125       else if (fCentralityClass=="6080") InitDplustoKpipi2010PbPb6080CentScan();
126       else if (fCentralityClass=="4080") InitDplustoKpipi2010PbPb4080();
127       else AliFatal("Not yet implemented");
128     }
129     break;
130   case 3: // D*->D0pi
131     if (fCollisionType==0) {
132       if(fIsLowEnergy)  InitDstartoD0pi2010ppLowEn();
133       else InitDstartoD0pi2010pp();
134     }else if (fCollisionType==1) {
135       if (fCentralityClass=="010") InitDstartoD0pi2010PbPb010CentScan();
136       else if (fCentralityClass=="1020") InitDstartoD0pi2010PbPb1020CentScan();
137       else if (fCentralityClass=="020") InitDstartoD0pi2010PbPb020();
138       else if (fCentralityClass=="2040" && fIsCentScan) InitDstartoD0pi2010PbPb2040CentScan();
139       else if (fCentralityClass=="2040") InitDstartoD0pi2010PbPb2040();
140       else if (fCentralityClass=="4060") InitDstartoD0pi2010PbPb4060CentScan();
141       else if (fCentralityClass=="6080") InitDstartoD0pi2010PbPb6080CentScan();
142       else if (fCentralityClass=="4080") InitDstartoD0pi2010PbPb4080();
143       else AliFatal("Not yet implemented");
144     }
145     break;
146   case 4: // D+s->KKpi
147     if (fCollisionType==0) InitDstoKKpi2010pp();
148     else AliFatal("Not yet implemented");
149     break;
150     
151   default:
152     printf("Invalid decay type: %d\n",decay);
153     break;
154   }
155
156 }
157   
158 //--------------------------------------------------------------------------
159 void AliHFSystErr::InitD0toKpi2010pp() {
160   // 
161   // D0->Kpi syst errors. Responsible: A. Rossi
162   //   2010 pp sample
163   //
164
165   // Normalization
166   fNorm = new TH1F("fNorm","fNorm",24,0,24);
167   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.035); // 4% error on sigmaV0and
168
169   // Branching ratio 
170   fBR = new TH1F("fBR","fBR",24,0,24);
171   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
172
173   // Tracking efficiency
174   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);
175   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.08); // 8% (4% per track)
176
177   // Raw yield extraction
178   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
179   fRawYield->SetBinContent(1,1);
180   fRawYield->SetBinContent(2,0.22);
181   fRawYield->SetBinContent(3,0.1);
182   for(Int_t i=4;i<=7;i++) fRawYield->SetBinContent(i,0.04);
183   for(Int_t i=8;i<=12;i++) fRawYield->SetBinContent(i,0.07);
184   for(Int_t i=13;i<=16;i++) fRawYield->SetBinContent(i,0.10);
185   for(Int_t i=17;i<=24;i++) fRawYield->SetBinContent(i,1);
186
187   // Cuts efficiency (from cuts variation)
188   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
189   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
190
191   // PID efficiency (from PID/noPID)
192   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
193   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.03); // 3%
194   fPIDEff->SetBinContent(2,0.05); // 5%
195
196   // MC dN/dpt
197   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
198   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0);
199   fMCPtShape->SetBinContent(1,0.03);
200   fMCPtShape->SetBinContent(2,0.03);
201
202   // particle-antiparticle
203   //  fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",24,0,24);
204   //  fPartAntipart->SetBinContent(1,1); 
205   //  for(Int_t i=2;i<=24;i++) fPartAntipart->SetBinContent(i,0.05);
206   
207   return;
208 }
209 //--------------------------------------------------------------------------
210 void AliHFSystErr::InitD0toKpi2010PbPb020() {
211   // 
212   // D0->Kpi syst errors. Responsible: A. Rossi
213   //   2010 PbPb sample, 0-20 CC
214   //
215
216   // Normalization
217   fNorm = new TH1F("fNorm","fNorm",20,0,20);
218   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.05); // TAA and pp norm
219
220   // Branching ratio 
221   fBR = new TH1F("fBR","fBR",20,0,20);
222   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
223
224   // Tracking efficiency
225   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
226   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.10);// Jacek, 5% per track
227
228   // Raw yield extraction
229   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
230   fRawYield->SetBinContent(1,0);
231   fRawYield->SetBinContent(2,0);
232   fRawYield->SetBinContent(3,0.08);
233   for(Int_t i=4;i<=12;i++) fRawYield->SetBinContent(i,0.06);
234   for(Int_t i=13;i<=16;i++) fRawYield->SetBinContent(i,0.10);
235
236   // Cuts efficiency (from cuts variation)
237   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
238   fCutsEff->SetBinContent(1,0.);
239   fCutsEff->SetBinContent(2,0.);
240   fCutsEff->SetBinContent(3,0.13);
241   fCutsEff->SetBinContent(4,0.11);  
242   for(Int_t i=5;i<=16;i++) fCutsEff->SetBinContent(i,0.10);
243   for(Int_t i=17;i<=20;i++) fCutsEff->SetBinContent(i,0.);
244
245   // PID efficiency (from PID/noPID)
246   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
247   for(Int_t i=3;i<=16;i++) fPIDEff->SetBinContent(i,0.05);
248
249   // MC dN/dpt
250   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
251   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.01);
252   fMCPtShape->SetBinContent(3,0.04);
253   fMCPtShape->SetBinContent(4,0.02);
254   for(Int_t i=13;i<=16;i++) fMCPtShape->SetBinContent(i,0.03); 
255
256 //   // particle-antiparticle
257 //   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
258 //   for(Int_t i=3;i<=12;i++) fPartAntipart->SetBinContent(i,0.05);
259 //   fPartAntipart->SetBinContent(3,0.10);
260 //   fPartAntipart->SetBinContent(4,0.10);
261 //   fPartAntipart->SetBinContent(7,0.10);
262 //   fPartAntipart->SetBinContent(8,0.10);
263   
264   return;
265 }
266 //--------------------------------------------------------------------------
267 void AliHFSystErr::InitD0toKpi2010PbPb4080() {
268   // 
269   // D0->Kpi syst errors. Responsible: A. Rossi
270   //   2010 PbPb sample, 40-80 CC
271   //
272
273   // Normalization
274   fNorm = new TH1F("fNorm","fNorm",20,0,20);
275   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.07); // TAA and pp norm
276
277   // Branching ratio 
278   fBR = new TH1F("fBR","fBR",20,0,20);
279   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
280
281   // Tracking efficiency
282   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
283   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.10); // Jacek, 5% per track
284
285
286   // Raw yield extraction
287   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
288   fRawYield->SetBinContent(1,0);
289   fRawYield->SetBinContent(2,0);
290   for(Int_t i=3;i<=16;i++) fRawYield->SetBinContent(i,0.05);
291   //for(Int_t i=13;i<=20;i++) fRawYield->SetBinContent(i,0);
292
293   // Cuts efficiency (from cuts variation)
294   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
295   fCutsEff->SetBinContent(1,0.);
296   fCutsEff->SetBinContent(2,0.);
297   fCutsEff->SetBinContent(3,0.13);
298   fCutsEff->SetBinContent(4,0.11);  
299   for(Int_t i=5;i<=16;i++) fCutsEff->SetBinContent(i,0.10);
300   for(Int_t i=17;i<=20;i++) fCutsEff->SetBinContent(i,0.);
301
302   // PID efficiency (from PID/noPID)
303   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
304 //   for(Int_t i=3;i<=6;i++) fPIDEff->SetBinContent(i,0.10);
305 //   for(Int_t i=7;i<=16;i++) fPIDEff->SetBinContent(i,0.05);
306   for(Int_t i=3;i<=16;i++) fPIDEff->SetBinContent(i,0.05);
307
308   // MC dN/dpt
309   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
310   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.01);
311   //  fMCPtShape->SetBinContent(3,0.04); Not set for peripherals (Raa Vs pt is flat)
312   //  fMCPtShape->SetBinContent(4,0.02);
313   for(Int_t i=13;i<=16;i++) fMCPtShape->SetBinContent(i,0.03); 
314   
315   //   // particle-antiparticle
316   //   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
317   //   for(Int_t i=3;i<=12;i++) fPartAntipart->SetBinContent(i,0.05);
318   
319   return;
320 }
321
322 //--------------------------------------------------------------------------
323 void AliHFSystErr::InitD0toKpi2010ppLowEn() {
324   // 
325   // D0->Kpi syst errors. Low energy run
326   //   2011 2.76 TeV pp sample
327   //
328   AliInfo(" Settings for D0 --> K pi, p-p collisions at 2.76 TeV"); 
329
330   // Normalization
331   fNorm = new TH1F("fNorm","fNorm",20,0,20);
332   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.031); // 4% error on sigmaV0and
333
334   // Branching ratio 
335   fBR = new TH1F("fBR","fBR",20,0,20);
336   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
337
338   // Tracking efficiency
339   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
340   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.10); //10% (5% per track)
341
342   // Raw yield extraction
343   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
344   fRawYield->SetBinContent(1,1);
345   for(Int_t i=1;i<=20;i++) fRawYield->SetBinContent(i,0.15);
346
347   // Cuts efficiency (from cuts variation)
348   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
349   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10); // 10% 
350   fCutsEff->SetBinContent(2,0.20);
351   for(Int_t i=7;i<=20;i++) fCutsEff->SetBinContent(i,0.15); // 10% 
352
353
354   // PID efficiency (from PID/noPID)
355   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
356   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.15); // 10%
357   //  fPIDEff->SetBinContent(2,0.20);
358   for(Int_t i=7;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 10%
359
360   // MC dN/dpt
361   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
362   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.01);
363   fMCPtShape->SetBinContent(1,0.03);
364   fMCPtShape->SetBinContent(2,0.03);
365
366 //   // particle-antiparticle
367 //   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
368 //   fPartAntipart->SetBinContent(1,1);
369 //   fPartAntipart->SetBinContent(2,1);
370 //   for(Int_t i=3;i<=6;i++) fPartAntipart->SetBinContent(i,0.08);
371 //   for(Int_t i=1;i<=20;i++) fPartAntipart->SetBinContent(i,0.);
372   
373   return;
374 }
375
376 //--------------------------------------------------------------------------
377 void AliHFSystErr::InitDplustoKpipi2010pp() {
378   // 
379   // D+->Kpipi syst errors. Responsible: R. Bala
380   //  2010 pp sample
381   //
382
383
384 // Normalization
385   fNorm = new TH1F("fNorm","fNorm",24,0,24);
386   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.035); // 4% error on sigmaV0and
387
388   // Branching ratio 
389   fBR = new TH1F("fBR","fBR",24,0,24);
390   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.021); // 2.1% PDG2010
391
392   // Tracking efficiency
393   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);
394   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.12); // 12% (4% per track)
395
396
397   // Raw yield extraction
398   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
399   fRawYield->SetBinContent(1,1);
400   fRawYield->SetBinContent(2,0.25);
401   fRawYield->SetBinContent(3,0.25);
402   fRawYield->SetBinContent(4,0.20);
403   fRawYield->SetBinContent(5,0.09);
404   fRawYield->SetBinContent(6,0.09);
405   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05);  
406   for(Int_t i=12;i<=24;i++) fRawYield->SetBinContent(i,0.10);  
407   
408   // Cuts efficiency (from cuts variation)
409   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
410   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
411
412   // PID efficiency (from PID/noPID)
413   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
414   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.05); // 5%
415   fPIDEff->SetBinContent(1,0.15); // 15%
416   fPIDEff->SetBinContent(2,0.15); // 15%
417   fPIDEff->SetBinContent(3,0.15); // 15%
418   fPIDEff->SetBinContent(4,0.15); // 15%
419   for(Int_t i=12;i<=16;i++) fPIDEff->SetBinContent(i,0.10); // 5%
420
421   // MC dN/dpt  
422   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
423   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0);
424   fMCPtShape->SetBinContent(1,0.03);
425   fMCPtShape->SetBinContent(2,0.03);
426
427
428   // particle-antiparticle
429   /*
430   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
431   fPartAntipart->SetBinContent(1,1);
432   fPartAntipart->SetBinContent(2,1);
433   fPartAntipart->SetBinContent(3,0.12);
434   for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
435   */
436   return;
437 }
438  
439 //--------------------------------------------------------------------------
440 void AliHFSystErr::InitDstoKKpi2010pp() {
441   // 
442   // D+s->KKpi syst errors. Responsible: G.M. Innocenti
443   //  2010 pp sample
444   //
445
446
447 // Normalization
448   fNorm = new TH1F("fNorm","fNorm",12,0,12);
449   for(Int_t i=1;i<=12;i++) fNorm->SetBinContent(i,0.035); // 3.5% error on sigmaV0and
450
451   // Branching ratio 
452   fBR = new TH1F("fBR","fBR",12,0,12);
453   for(Int_t i=1;i<=12;i++) fBR->SetBinContent(i,0.06); // 0.14/2.32 PDG2010
454
455   // Tracking efficiency
456   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",12,0,12);
457   for(Int_t i=1;i<=12;i++) fTrackingEff->SetBinContent(i,0.12); // 12% (4% per track)
458
459
460   // Raw yield extraction
461   fRawYield = new TH1F("fRawYield","fRawYield",12,0,12);
462   fRawYield->SetBinContent(1,1);
463   fRawYield->SetBinContent(2,1);
464   fRawYield->SetBinContent(3,0.20);
465   fRawYield->SetBinContent(4,0.20);
466   fRawYield->SetBinContent(5,0.15);
467   fRawYield->SetBinContent(6,0.15);
468   fRawYield->SetBinContent(7,0.15);
469   fRawYield->SetBinContent(8,0.15);
470   fRawYield->SetBinContent(9,0.20);
471   fRawYield->SetBinContent(10,0.20);
472   fRawYield->SetBinContent(11,0.20);
473   fRawYield->SetBinContent(12,0.20);
474   
475   // Cuts efficiency (from cuts variation)
476   fCutsEff = new TH1F("fCutsEff","fCutsEff",12,0,12);
477   for(Int_t i=1;i<=12;i++) fCutsEff->SetBinContent(i,0.15); // 15%
478
479   // PID efficiency (from PID/noPID)
480   fPIDEff = new TH1F("fPIDEff","fPIDEff",12,0,12);
481   for(Int_t i=1;i<=12;i++) fPIDEff->SetBinContent(i,0.07); // 7%
482
483   // MC dN/dpt 
484   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",12,0,12);
485   for(Int_t i=1; i<=2; i++) fMCPtShape->SetBinContent(i,1.);
486   for(Int_t i=3; i<=4; i++) fMCPtShape->SetBinContent(i,0.03);
487   for(Int_t i=5; i<=6; i++) fMCPtShape->SetBinContent(i,0.03);
488   for(Int_t i=7; i<=8; i++) fMCPtShape->SetBinContent(i,0.02);
489   for(Int_t i=9; i<=12; i++) fMCPtShape->SetBinContent(i,0.02);
490
491
492   // particle-antiparticle
493   /*
494   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",12,0,12);
495   fPartAntipart->SetBinContent(1,1);
496   fPartAntipart->SetBinContent(2,1);
497   fPartAntipart->SetBinContent(3,0.12);
498   for(Int_t i=4;i<=12;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
499   */
500   return;
501 }
502    
503  
504 //--------------------------------------------------------------------------
505 void AliHFSystErr::InitDplustoKpipi2010PbPb020() {
506   // 
507   // D+->Kpipi syst errors. Responsible: ??
508   //  2010 PbPb sample, 0-20 CC
509   //
510
511  // Normalization
512   fNorm = new TH1F("fNorm","fNorm",20,0,20);
513   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.05); // TAA and pp norm
514
515   // Branching ratio 
516   fBR = new TH1F("fBR","fBR",20,0,20);
517   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.021); // 2.1% PDG2010
518
519   // Tracking efficiency
520   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
521   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.15); // Jacek, 5% per track
522
523   // Raw yield extraction
524   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
525   for(Int_t i=1;i<=20;i++) fRawYield->SetBinContent(i,.10);  //5 to 10%
526   // fRawYield->SetBinContent(5,0.23);
527   //fRawYield->SetBinContent(6,0.23);
528   fRawYield->SetBinContent(7,0.20);
529   fRawYield->SetBinContent(8,0.20);
530   fRawYield->SetBinContent(9,0.15);
531   fRawYield->SetBinContent(10,0.15);
532   fRawYield->SetBinContent(11,0.15);
533   fRawYield->SetBinContent(12,0.15);
534
535   // Cuts efficiency (from cuts variation)
536   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
537   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.15); // 10%
538
539   // PID efficiency (from PID/noPID)
540   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
541   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
542
543   // MC dN/dpt  (2/2/2012)
544   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
545   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.);
546   for(Int_t iBin=7; iBin<=8; iBin++) fMCPtShape->SetBinContent(iBin,0.01);
547   for(Int_t iBin=9; iBin<=12; iBin++) fMCPtShape->SetBinContent(iBin,0.05);
548   for(Int_t iBin=13; iBin<=16; iBin++) fMCPtShape->SetBinContent(iBin,0.05);
549
550
551   // particle-antiparticle
552   /*
553   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
554   fPartAntipart->SetBinContent(1,1);
555   fPartAntipart->SetBinContent(2,1);
556   fPartAntipart->SetBinContent(3,0.12);
557   for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
558   */
559
560   return;
561 }
562
563 //--------------------------------------------------------------------------
564 void AliHFSystErr::InitDplustoKpipi2010PbPb4080() {
565   // 
566   // D+->Kpipi syst errors. Responsible: ??
567   //  2010 PbPb sample, 40-80 CC
568   //
569   
570
571  // Normalization
572   fNorm = new TH1F("fNorm","fNorm",20,0,20);
573   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.07); // TAA and pp norm
574
575   // Branching ratio 
576   fBR = new TH1F("fBR","fBR",20,0,20);
577   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.021); // 2.1% 
578
579   // Tracking efficiency
580   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
581   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.15); // Jacek, 5% per track
582
583
584   // Raw yield extraction
585   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
586   fRawYield->SetBinContent(1,1);
587   fRawYield->SetBinContent(2,1);
588   fRawYield->SetBinContent(3,1);
589   fRawYield->SetBinContent(4,0.15);
590   fRawYield->SetBinContent(5,0.05);
591   fRawYield->SetBinContent(6,0.05);
592   fRawYield->SetBinContent(7,0.15);
593   fRawYield->SetBinContent(8,0.15);
594   for(Int_t i=9;i<=12;i++) fRawYield->SetBinContent(i,0.15);
595   for(Int_t i=13;i<=20;i++) fRawYield->SetBinContent(i,1);  //5 to 10%
596
597   // Cuts efficiency (from cuts variation)
598   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
599   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10); // 10%
600
601   // PID efficiency (from PID/noPID)
602   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
603   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
604   fPIDEff->SetBinContent(3,0.13); // 13%
605  
606
607   // MC dN/dpt  (2/2/2012)
608   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
609   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0);
610   for(Int_t iBin=4; iBin<=8; iBin++) fMCPtShape->SetBinContent(iBin,0.01);
611   for(Int_t iBin=9; iBin<=12; iBin++) fMCPtShape->SetBinContent(iBin,0.03);
612   for(Int_t iBin=13; iBin<=16; iBin++) fMCPtShape->SetBinContent(iBin,0.03);
613
614
615   // particle-antiparticle
616   /*
617   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
618   fPartAntipart->SetBinContent(1,1);
619   fPartAntipart->SetBinContent(2,1);
620   fPartAntipart->SetBinContent(3,0.12);
621   for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
622   */
623   return;
624 }
625
626 //--------------------------------------------------------------------------
627 void AliHFSystErr::InitDplustoKpipi2010ppLowEn() {
628
629   // 
630   // D+->Kpipi syst errors. Responsible: R. Bala
631   //  2011 2.76 TeV pp sample
632   //
633   AliInfo(" Settings for D+ --> K pi pi p-p collisions at 2.76 TeV"); 
634
635   // Normalization
636   fNorm = new TH1F("fNorm","fNorm",20,0,20);
637   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.031); // 10% error on sigmaV0and
638
639   // Branching ratio 
640   fBR = new TH1F("fBR","fBR",20,0,20);
641   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.021); // 2.1% PDG2010
642
643   // Tracking efficiency
644   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
645   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.15); // 3% (1% per track)
646
647   // Raw yield extraction
648   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
649   fRawYield->SetBinContent(1,1);
650   fRawYield->SetBinContent(2,1);
651   for(Int_t i=3;i<=6;i++) fRawYield->SetBinContent(i,0.10);  //5 to 10%
652   fRawYield->SetBinContent(7,0.15);
653   fRawYield->SetBinContent(8,0.15); 
654   for(Int_t i=9;i<=20;i++) fRawYield->SetBinContent(i,0.055);  //5 to 10%
655
656   // Cuts efficiency (from cuts variation)
657   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
658   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.15); // 10%
659
660   // PID efficiency (from PID/noPID)
661   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
662   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
663   fPIDEff->SetBinContent(3,0.10); // 13%
664   fPIDEff->SetBinContent(4,0.10); // 13%
665  
666   // MC dN/dpt  (copied from D0 : will update later)
667   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
668   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.01);
669   fMCPtShape->SetBinContent(1,0.03);
670   fMCPtShape->SetBinContent(2,0.03);
671
672   return;
673 }
674
675 //--------------------------------------------------------------------------
676 void AliHFSystErr::InitDstartoD0pi2010pp() {
677   // 
678   // D*+->D0pi syst errors. Responsible: A. Grelli, Y. Wang
679   //  2010 pp sample
680   //
681
682  // Normalization
683   fNorm = new TH1F("fNorm","fNorm",24,0,24);
684   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.035); // 4% error on sigmaV0and
685
686   // Branching ratio 
687   fBR = new TH1F("fBR","fBR",24,0,24);
688   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
689
690   // Tracking efficiency
691   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);
692   fTrackingEff->SetBinContent(1,1.0);
693   fTrackingEff->SetBinContent(2,0.13); // 10% (ITSsa) \oplus 8% (4% per ITSTPC track)
694   fTrackingEff->SetBinContent(3,0.12);
695   fTrackingEff->SetBinContent(3,0.12);
696   for(Int_t i=4;i<=24;i++) fTrackingEff->SetBinContent(i,0.12); // 12% (4% per track)
697
698
699   // Raw yield extraction
700   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
701   fRawYield->SetBinContent(1,1.0);
702   fRawYield->SetBinContent(2,0.10);
703   fRawYield->SetBinContent(3,0.04);
704   fRawYield->SetBinContent(4,0.03);
705   fRawYield->SetBinContent(5,0.03);
706   fRawYield->SetBinContent(6,0.05);
707   fRawYield->SetBinContent(7,0.05);
708   fRawYield->SetBinContent(8,0.05);
709   for(Int_t i=9;i<=12;i++) fRawYield->SetBinContent(i,0.04);  //4%
710   for(Int_t i=13;i<=16;i++) fRawYield->SetBinContent(i,0.09);  //4%
711   for(Int_t i=17;i<=24;i++) fRawYield->SetBinContent(i,0.2);  //4%
712
713   // Cuts efficiency (from cuts variation)
714   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
715   fCutsEff->SetBinContent(2,0.22);
716   for(Int_t i=3;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
717
718   // PID efficiency (from PID/noPID)
719   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
720   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.04); // 3%
721  
722
723   // MC dN/dpt  (copied from D0 : will update later)
724   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
725   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0);
726   fMCPtShape->SetBinContent(1,0.03);
727   fMCPtShape->SetBinContent(2,0.03);
728
729   return;
730
731
732 }
733 //--------------------------------------------------------------------------
734 void AliHFSystErr::InitDstartoD0pi2010ppLowEn() {
735
736   // 
737   // D+->Kpipi syst errors. Responsible: A. Grelli
738   //  2011 2.76 TeV pp sample
739   //
740   AliInfo(" Settings for D*+ --> D0 pi p-p collisions at 2.76 TeV"); 
741
742 // Normalization
743   fNorm = new TH1F("fNorm","fNorm",20,0,20);
744   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.031); // 10% error on sigmaV0and
745
746   // Branching ratio 
747   fBR = new TH1F("fBR","fBR",20,0,20);
748   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
749
750   // Tracking efficiency
751   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
752   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.15); //10% (to be checked!!)
753
754   // Raw yield extraction
755   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
756   fRawYield->SetBinContent(1,1);
757   fRawYield->SetBinContent(2,1);
758   fRawYield->SetBinContent(3,0.14);
759   fRawYield->SetBinContent(4,0.14);
760   fRawYield->SetBinContent(5,0.12);
761   fRawYield->SetBinContent(6,0.12);
762   fRawYield->SetBinContent(7,0.06);
763   fRawYield->SetBinContent(8,0.06);
764   fRawYield->SetBinContent(9,0.08);
765   fRawYield->SetBinContent(10,0.08);
766   fRawYield->SetBinContent(11,0.08);
767   fRawYield->SetBinContent(12,0.08);
768   for(Int_t i=9;i<=20;i++) fRawYield->SetBinContent(i,0.065);
769
770   // Cuts efficiency (from cuts variation)
771   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
772   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10);  
773   fCutsEff->SetBinContent(3,0.15);
774   fCutsEff->SetBinContent(4,0.15);
775   fCutsEff->SetBinContent(5,0.15);
776   fCutsEff->SetBinContent(6,0.15);
777   fCutsEff->SetBinContent(7,0.10);
778   fCutsEff->SetBinContent(8,0.10);
779   fCutsEff->SetBinContent(9,0.10);
780   fCutsEff->SetBinContent(10,0.10);
781   fCutsEff->SetBinContent(11,0.10);
782   fCutsEff->SetBinContent(12,0.10);
783
784   // PID efficiency (from PID/noPID)
785   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
786   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 10%
787
788   // MC dN/dpt
789   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
790   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.01);
791   fMCPtShape->SetBinContent(1,0.03);
792   fMCPtShape->SetBinContent(2,0.03);
793
794   return;
795 }
796
797 //------------------------------------------------------------------------
798 void AliHFSystErr::InitDstartoD0pi2010PbPb020() {
799   // 
800   // D*+->D0pi syst errors. Responsible: A. Grelli
801   //  2010 PbPb sample, 0-20 CC
802   //
803
804   AliInfo(" Settings for D*+ --> D0pi Pb-Pb collisions at 2.76 TeV - 0-20 centrality - DUMMY"); 
805
806  // Normalization
807   fNorm = new TH1F("fNorm","fNorm",24,0,24);
808   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.05); // TAA and pp norm
809
810   // Branching ratio 
811   fBR = new TH1F("fBR","fBR",24,0,24);
812   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
813
814   // Tracking efficiency
815   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);;
816   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.15); // Jacek, 5% per track
817
818
819   // Raw yield extraction
820   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
821   for(Int_t i=1;i<=24;i++) fRawYield->SetBinContent(i,0.1);  //4%
822   fRawYield->SetBinContent(3,0.2);
823   fRawYield->SetBinContent(4,0.2);
824   fRawYield->SetBinContent(5,0.2);
825   fRawYield->SetBinContent(6,0.2);
826  
827   // Cuts efficiency (from cuts variation)
828   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
829   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
830   fCutsEff->SetBinContent(4,0.15);
831   fCutsEff->SetBinContent(5,0.15);
832   fCutsEff->SetBinContent(6,0.15);
833
834   // PID efficiency (from PID/noPID)
835   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
836   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.05); // 3%
837  
838
839   // MC dN/dpt  (from study on D* pt shape)
840   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
841   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0.045);
842   fMCPtShape->SetBinContent(4,0.025);
843   fMCPtShape->SetBinContent(5,0.025);
844   fMCPtShape->SetBinContent(6,0.025);
845   fMCPtShape->SetBinContent(7,0.04);
846   fMCPtShape->SetBinContent(8,0.04);
847   fMCPtShape->SetBinContent(9,0.03);
848   fMCPtShape->SetBinContent(10,0.03);
849   fMCPtShape->SetBinContent(11,0.03);
850   fMCPtShape->SetBinContent(12,0.03);
851   
852   
853
854   return;
855
856 }
857
858 //-------------------------------------------------------------------------
859 void AliHFSystErr::InitDstartoD0pi2010PbPb2040() {
860   // 
861   // D*+->D0pi syst errors. Responsible: A. Grelli
862   //  2010 PbPb sample, 20-40 CC
863   //
864
865   AliInfo(" Settings for D*+ --> D0pi Pb-Pb collisions at 2.76 TeV - 20-40 centrality - DUMMY"); 
866
867  // Normalization
868   fNorm = new TH1F("fNorm","fNorm",24,0,24);
869   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
870
871   // Branching ratio 
872   fBR = new TH1F("fBR","fBR",24,0,24);
873   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
874
875   // Tracking efficiency
876   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);;
877   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.15); // Jacek, 5% per track
878
879
880   // Raw yield extraction
881   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
882   for(Int_t i=1;i<=24;i++) fRawYield->SetBinContent(i,0.15);  //4%
883  
884   // Cuts efficiency (from cuts variation)
885   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
886   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
887
888   // PID efficiency (from PID/noPID)
889   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
890   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.04); // 3%
891  
892
893   // MC dN/dpt  (copied from D0 : will update later)
894   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
895   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0.);
896   fMCPtShape->SetBinContent(1,0.03);
897   fMCPtShape->SetBinContent(2,0.03);
898
899   return;
900
901 }
902
903 //--------------------------------------------------------------------------
904 void AliHFSystErr::InitDstartoD0pi2010PbPb4080() {
905   // 
906   // D*+->D0pi syst errors. Responsible: A. Grelli
907   //  2010 PbPb sample, 40-80 CC
908   //
909
910   AliInfo(" Settings for D*+ --> D0pi Pb-Pb collisions at 2.76 TeV - 40-80 centrality - DUMMY"); 
911
912  // Normalization
913   fNorm = new TH1F("fNorm","fNorm",24,0,24);
914   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.07); // TAA and pp norm
915
916   // Branching ratio 
917   fBR = new TH1F("fBR","fBR",24,0,24);
918   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
919
920   // Tracking efficiency
921   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);;
922   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.15); // Jacek, 5% per track
923
924
925   // Raw yield extraction
926   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
927   for(Int_t i=1;i<=24;i++) fRawYield->SetBinContent(i,0.2);  //4%
928   fRawYield->SetBinContent(1,1);
929   fRawYield->SetBinContent(2,0.15);
930   fRawYield->SetBinContent(3,0.15);
931   fRawYield->SetBinContent(4,0.15);
932   fRawYield->SetBinContent(5,0.15);
933   fRawYield->SetBinContent(6,0.10);
934   fRawYield->SetBinContent(7,0.10);
935   fRawYield->SetBinContent(8,0.10);
936   fRawYield->SetBinContent(9,0.11);
937   fRawYield->SetBinContent(10,0.11);
938   fRawYield->SetBinContent(11,0.11);
939   fRawYield->SetBinContent(12,0.11);
940   fRawYield->SetBinContent(13,0.08);
941   fRawYield->SetBinContent(14,0.08);
942   fRawYield->SetBinContent(15,0.08);
943   fRawYield->SetBinContent(16,0.08);
944
945
946   // Cuts efficiency (from cuts variation)
947   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
948   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
949
950   // PID efficiency (from PID/noPID)
951   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
952   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.05); // 3%
953  
954
955   // MC dN/dpt  (copied from D0 : will update later)
956   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
957   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0.01);
958   fMCPtShape->SetBinContent(2,0.05);
959   fMCPtShape->SetBinContent(3,0.05);
960   fMCPtShape->SetBinContent(4,0.05);
961   fMCPtShape->SetBinContent(5,0.04);
962   fMCPtShape->SetBinContent(6,0.02);
963   fMCPtShape->SetBinContent(7,0.04);
964   fMCPtShape->SetBinContent(8,0.04);
965  
966   return;
967
968 }
969
970 //--------------------------------------------------------------------------
971 void AliHFSystErr::InitD0toKpi2010PbPb010CentScan(){
972   // define errors for RAA vs. centrality
973   InitD0toKpi2010PbPb020();
974   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05);
975   for(Int_t i=3;i<=5;i++) fMCPtShape->SetBinContent(i,0.17);      
976   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08); 
977 }
978 //--------------------------------------------------------------------------
979 void AliHFSystErr::InitD0toKpi2010PbPb1020CentScan(){
980   // define errors for RAA vs. centrality
981   InitD0toKpi2010PbPb020();
982   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05);
983   for(Int_t i=3;i<=5;i++)  fMCPtShape->SetBinContent(i,0.17);    
984   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08); 
985 }
986 //--------------------------------------------------------------------------
987 void AliHFSystErr::InitD0toKpi2010PbPb2040CentScan(){
988   // define errors for RAA vs. centrality
989   InitD0toKpi2010PbPb4080();
990   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05);
991   for(Int_t i=3;i<=5;i++)  fMCPtShape->SetBinContent(i,0.14);
992   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08); 
993 }
994 //--------------------------------------------------------------------------
995 void AliHFSystErr::InitD0toKpi2010PbPb4060CentScan(){
996    // define errors for RAA vs. centrality
997   InitD0toKpi2010PbPb4080();
998   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.06);
999   for(Int_t i=3;i<=5;i++)  fMCPtShape->SetBinContent(i,0.11);  
1000   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08); 
1001 }
1002 //--------------------------------------------------------------------------
1003 void AliHFSystErr::InitD0toKpi2010PbPb6080CentScan(){
1004    // define errors for RAA vs. centrality
1005   InitD0toKpi2010PbPb4080();
1006   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.08);
1007   for(Int_t i=3;i<=5;i++) fMCPtShape->SetBinContent(i,0.08);
1008   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08); 
1009 }
1010 //--------------------------------------------------------------------------
1011 void AliHFSystErr::InitDplustoKpipi2010PbPb010CentScan(){
1012   // define errors for RAA vs. centrality
1013   InitDplustoKpipi2010PbPb020();
1014   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.18);
1015   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.09);
1016
1017 }
1018 //--------------------------------------------------------------------------
1019 void AliHFSystErr::InitDplustoKpipi2010PbPb1020CentScan(){
1020   // define errors for RAA vs. centrality
1021   InitDplustoKpipi2010PbPb020();
1022   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.23);
1023   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08);
1024 }
1025 //--------------------------------------------------------------------------
1026 void AliHFSystErr::InitDplustoKpipi2010PbPb2040CentScan(){
1027   // define errors for RAA vs. centrality
1028   InitDplustoKpipi2010PbPb020();
1029   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.08);
1030   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.095);
1031 }
1032 //--------------------------------------------------------------------------
1033 void AliHFSystErr::InitDplustoKpipi2010PbPb4060CentScan(){
1034   // define errors for RAA vs. centrality
1035   InitDplustoKpipi2010PbPb4080();
1036   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.08);
1037   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08);
1038 }
1039 //--------------------------------------------------------------------------
1040 void AliHFSystErr::InitDplustoKpipi2010PbPb6080CentScan(){
1041   // define errors for RAA vs. centrality
1042   InitDplustoKpipi2010PbPb4080();
1043   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.15);
1044   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.07);
1045 }
1046
1047 //--------------------------------------------------------------------------
1048 void AliHFSystErr::InitDstartoD0pi2010PbPb010CentScan(){
1049   // define errors for RAA vs. centrality
1050   InitDstartoD0pi2010PbPb020();
1051   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.16); 
1052   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.15);
1053 }
1054 //--------------------------------------------------------------------------
1055 void AliHFSystErr::InitDstartoD0pi2010PbPb1020CentScan(){
1056   // define errors for RAA vs. centrality
1057   InitDstartoD0pi2010PbPb020();
1058   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05); 
1059   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.15);
1060 }
1061 //--------------------------------------------------------------------------
1062 void AliHFSystErr::InitDstartoD0pi2010PbPb2040CentScan(){
1063   // define errors for RAA vs. centrality
1064   InitDstartoD0pi2010PbPb2040();
1065   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.10); 
1066   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.08);
1067 }
1068 //--------------------------------------------------------------------------
1069 void AliHFSystErr::InitDstartoD0pi2010PbPb4060CentScan(){
1070   // define errors for RAA vs. centrality
1071   InitDstartoD0pi2010PbPb4080();
1072   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.10); 
1073   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.045);
1074 }
1075 //--------------------------------------------------------------------------
1076 void AliHFSystErr::InitDstartoD0pi2010PbPb6080CentScan(){
1077   // define errors for RAA vs. centrality
1078   InitDstartoD0pi2010PbPb4080();
1079   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.10); 
1080   for(Int_t i=7;i<=12;i++) fMCPtShape->SetBinContent(i,0.045);
1081 }
1082
1083
1084 //--------------------------------------------------------------------------
1085 Double_t AliHFSystErr::GetCutsEffErr(Double_t pt) const {
1086   // 
1087   // Get error
1088   //
1089
1090   Int_t bin=fCutsEff->FindBin(pt);
1091
1092   return fCutsEff->GetBinContent(bin);
1093 }
1094 //--------------------------------------------------------------------------
1095 Double_t AliHFSystErr::GetMCPtShapeErr(Double_t pt) const {
1096   // 
1097   // Get error
1098   //
1099
1100   Int_t bin=fMCPtShape->FindBin(pt);
1101
1102   return fMCPtShape->GetBinContent(bin);
1103 }
1104 //--------------------------------------------------------------------------
1105 Double_t AliHFSystErr::GetSeleEffErr(Double_t pt) const {
1106   // 
1107   // Get error
1108   //
1109
1110   Double_t err=GetCutsEffErr(pt)*GetCutsEffErr(pt)+GetMCPtShapeErr(pt)*GetMCPtShapeErr(pt);
1111
1112   return TMath::Sqrt(err);
1113 }
1114 //--------------------------------------------------------------------------
1115 Double_t AliHFSystErr::GetPIDEffErr(Double_t pt) const {
1116   // 
1117   // Get error
1118   //
1119
1120   Int_t bin=fPIDEff->FindBin(pt);
1121
1122   return fPIDEff->GetBinContent(bin);
1123 }
1124 //--------------------------------------------------------------------------
1125 Double_t AliHFSystErr::GetTrackingEffErr(Double_t pt) const {
1126   // 
1127   // Get error
1128   //
1129
1130   Int_t bin=fTrackingEff->FindBin(pt);
1131
1132   return fTrackingEff->GetBinContent(bin);
1133 }
1134 //--------------------------------------------------------------------------
1135 Double_t AliHFSystErr::GetRawYieldErr(Double_t pt) const {
1136   // 
1137   // Get error
1138   //
1139
1140   Int_t bin=fRawYield->FindBin(pt);
1141
1142   return fRawYield->GetBinContent(bin);
1143 }
1144 //--------------------------------------------------------------------------
1145 Double_t AliHFSystErr::GetPartAntipartErr(Double_t pt) const {
1146   // 
1147   // Get error
1148   //
1149
1150   Int_t bin=fPartAntipart->FindBin(pt);
1151
1152   return fPartAntipart->GetBinContent(bin);
1153 }
1154 //--------------------------------------------------------------------------
1155 Double_t AliHFSystErr::GetTotalSystErr(Double_t pt,Double_t feeddownErr) const {
1156   // 
1157   // Get total syst error (except norm. error)
1158   //
1159
1160   Double_t err=0.;
1161
1162   if(fRawYield) err += GetRawYieldErr(pt)*GetRawYieldErr(pt);
1163   if(fTrackingEff) err += GetTrackingEffErr(pt)*GetTrackingEffErr(pt);
1164   //  if(fBR) err += GetBRErr()*GetBRErr();
1165   if(fCutsEff) err += GetCutsEffErr(pt)*GetCutsEffErr(pt);
1166   if(fPIDEff) err += GetPIDEffErr(pt)*GetPIDEffErr(pt);
1167   if(fMCPtShape) err += GetMCPtShapeErr(pt)*GetMCPtShapeErr(pt);
1168   if(fPartAntipart) err += GetPartAntipartErr(pt)*GetPartAntipartErr(pt);
1169
1170   err += feeddownErr*feeddownErr;
1171
1172   return TMath::Sqrt(err);
1173 }
1174 //---------------------------------------------------------------------------
1175 void AliHFSystErr::DrawErrors(TGraphAsymmErrors *grErrFeeddown) const {
1176   //
1177   // Draw errors
1178   //
1179   gStyle->SetOptStat(0);
1180
1181   TCanvas *cSystErr = new TCanvas("cSystErr","Systematic Errors",300,80,640,500);
1182   cSystErr->Range(0.20,-0.5,18.4,0.34);
1183   cSystErr->SetRightMargin(0.318);
1184   cSystErr->SetFillColor(0);
1185
1186   TH2F *hFrame = new TH2F("hFrame","Systematic errors; p_{t} (GeV/c); Relative Error",30,0,30,100,-1,+1);
1187   hFrame->SetAxisRange(1.,24.9,"X");
1188   hFrame->SetAxisRange(-0.5,0.5,"Y");
1189   hFrame->Draw();
1190
1191   TLegend *leg = new TLegend(0.69,0.44,0.98,0.86,NULL,"brNDC");
1192   leg->SetTextSize(0.03601695);
1193   leg->SetFillStyle(0);
1194   leg->SetBorderSize(0);
1195   
1196   TH1F *hTotErr=new TH1F("hTotErr","",24,0,24);
1197   Int_t nbins = fNorm->GetNbinsX();
1198   TGraphAsymmErrors *gTotErr = new TGraphAsymmErrors(nbins);
1199   for(Int_t i=1;i<=24;i++) {
1200     Double_t pt = hTotErr->GetBinCenter(i);
1201     Double_t ptwidth = hTotErr->GetBinWidth(i);
1202
1203     if(grErrFeeddown) {
1204       Double_t x=0., y=0., errxl=0., errxh=0., erryl=0., erryh=0.;
1205       Double_t toterryl=0., toterryh=0.;
1206       for(Int_t j=0; j<grErrFeeddown->GetN(); j++) {
1207         grErrFeeddown->GetPoint(j,x,y);
1208         errxh = grErrFeeddown->GetErrorXhigh(j);
1209         errxl = grErrFeeddown->GetErrorXlow(j);
1210         if ( ( (x-errxl) <= pt) && ( (x+errxl) >= pt) ) {
1211           erryh = grErrFeeddown->GetErrorYhigh(j);
1212           erryl = grErrFeeddown->GetErrorYlow(j);
1213         }
1214       }
1215       if (erryl>=1e-3) toterryl = GetTotalSystErr(pt,erryl);
1216       else toterryl = GetTotalSystErr(pt);
1217       if (erryh>=1e-3) toterryh = GetTotalSystErr(pt,erryh);
1218       else toterryh = GetTotalSystErr(pt);
1219
1220       hTotErr->SetBinContent(i,toterryh);
1221       gTotErr->SetPoint(i,pt,0.);
1222       gTotErr->SetPointError(i,ptwidth/2.,ptwidth/2.,toterryl,toterryh); // i, exl, exh, eyl, eyh
1223     }
1224     else {
1225       hTotErr->SetBinContent(i,GetTotalSystErr(pt));
1226       gTotErr->SetPoint(i,pt,0.);
1227       gTotErr->SetPointError(i,ptwidth/2.,ptwidth/2.,GetTotalSystErr(pt),GetTotalSystErr(pt)); // i, exl, exh, eyl, eyh
1228     }
1229
1230   }
1231   gTotErr->SetLineColor(kBlack);
1232   gTotErr->SetFillColor(kRed);
1233   gTotErr->SetFillStyle(3002);
1234   gTotErr->Draw("2");
1235   leg->AddEntry(gTotErr,"Total (excl. norm.)","f");
1236 //   hTotErr->SetLineColor(1);
1237 //   hTotErr->SetLineWidth(3);
1238 //   hTotErr->Draw("same");
1239 //   leg->AddEntry(hTotErr,"Total (excl. norm.)","l");
1240   
1241
1242   fNorm->SetFillColor(1);
1243   fNorm->SetFillStyle(3002);
1244   //fNorm->Draw("same");
1245   //TH1F *hNormRefl = ReflectHisto(fNorm);
1246   //hNormRefl->Draw("same");
1247   Double_t norm = fNorm->GetBinContent(1)*100;
1248   leg->AddEntry(fNorm,Form("Normalization (%.1f%s)",norm,"%"),"");
1249
1250   if(grErrFeeddown) {
1251     grErrFeeddown->SetFillColor(kTeal-8);
1252     grErrFeeddown->SetFillStyle(3001);
1253     grErrFeeddown->Draw("2");
1254     leg->AddEntry(grErrFeeddown,"Feed-down from B","f");
1255   }
1256   if(fTrackingEff) {
1257     fTrackingEff->SetFillColor(4);
1258     fTrackingEff->SetFillStyle(3006);
1259     fTrackingEff->Draw("same");
1260     TH1F *hTrackingEffRefl = ReflectHisto(fTrackingEff);
1261     hTrackingEffRefl->Draw("same");
1262     leg->AddEntry(fTrackingEff,"Tracking efficiency","f");
1263   }
1264   if(fBR) {
1265     fBR->SetFillColor(6);
1266     fBR->SetFillStyle(3005);
1267     //fBR->SetFillStyle(3020);
1268     fBR->Draw("same");
1269     TH1F *hBRRefl = ReflectHisto(fBR);
1270     hBRRefl->Draw("same");
1271     leg->AddEntry(fBR,"Branching ratio","f");
1272   }
1273   if(fRawYield) {
1274     Int_t ci;   // for color index setting
1275     ci = TColor::GetColor("#00cc00");
1276     fRawYield->SetLineColor(ci);
1277     //    fRawYield->SetLineColor(3);
1278     fRawYield->SetLineWidth(3);
1279     fRawYield->Draw("same");
1280     TH1F *hRawYieldRefl = ReflectHisto(fRawYield);
1281     hRawYieldRefl->Draw("same");
1282     leg->AddEntry(fRawYield,"Yield extraction","l");
1283   }
1284   if(fCutsEff) {
1285     fCutsEff->SetLineColor(4);
1286     fCutsEff->SetLineWidth(3);
1287     fCutsEff->Draw("same");
1288     TH1F *hCutsEffRefl = ReflectHisto(fCutsEff);
1289     hCutsEffRefl->Draw("same");
1290     leg->AddEntry(fCutsEff,"Cuts efficiency","l");
1291   }
1292   if(fPIDEff) {
1293     fPIDEff->SetLineColor(7);
1294     fPIDEff->SetLineWidth(3);
1295     fPIDEff->Draw("same");
1296     TH1F *hPIDEffRefl = ReflectHisto(fPIDEff);
1297     hPIDEffRefl->Draw("same");
1298     leg->AddEntry(fPIDEff,"PID efficiency","l");
1299   }
1300   if(fMCPtShape) {
1301     Int_t ci = TColor::GetColor("#9933ff");
1302     fMCPtShape->SetLineColor(ci);
1303     //    fMCPtShape->SetLineColor(8);
1304     fMCPtShape->SetLineWidth(3);
1305     fMCPtShape->Draw("same");
1306     TH1F *hMCPtShapeRefl = ReflectHisto(fMCPtShape);
1307     hMCPtShapeRefl->Draw("same");
1308     leg->AddEntry(fMCPtShape,"MC p_{t} shape","l");
1309   }
1310   if(fPartAntipart) {
1311     Int_t ci = TColor::GetColor("#ff6600");
1312     fPartAntipart->SetLineColor(ci);
1313     //    fPartAntipart->SetLineColor(9);
1314     fPartAntipart->SetLineWidth(3);
1315     fPartAntipart->Draw("same");
1316     TH1F *hPartAntipartRefl = ReflectHisto(fPartAntipart);
1317     hPartAntipartRefl->Draw("same");
1318     leg->AddEntry(fPartAntipart,"D = #bar{D}","l");
1319   }
1320
1321
1322   leg->Draw();
1323
1324   cSystErr->SaveAs("RelativeSystematics.eps");
1325
1326   return;
1327 }
1328 //-------------------------------------------------------------------------
1329 TH1F* AliHFSystErr::ReflectHisto(TH1F *hin) const {
1330   //
1331   // Clones and reflects histogram 
1332   // 
1333   TH1F *hout=(TH1F*)hin->Clone("hout");
1334   hout->Scale(-1.);
1335
1336   return hout;
1337 }