PWG3/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 {
  64   //
  65   // Default Constructor
  66   //
  67 }
  68
  69 //--------------------------------------------------------------------------
  70 AliHFSystErr::~AliHFSystErr() {
  71   //
  72   // Default Destructor
  73   //
  74   /*
  75
  76   if(fNorm)         { delete fNorm; fNorm=0; }
  77   if(fRawYield)     { delete fRawYield; fRawYield=0; }
  78   if(fTrackingEff)  { delete fTrackingEff; fTrackingEff=0; }
  79   if(fBR)           { delete fBR; fBR=0; }
  80   if(fCutsEff)      { delete fCutsEff; fCutsEff=0; }
  81   if(fPIDEff)       { delete fPIDEff; fPIDEff=0; }
  82   if(fMCPtShape)    { delete fMCPtShape; fMCPtShape=0; }
  83   if(fPartAntipart) { delete fPartAntipart; fPartAntipart=0; }
  84   */
  85 }
  86
  87 //--------------------------------------------------------------------------
  88 void AliHFSystErr::Init(Int_t decay){
  89   //
  90   // Variables/histos initialization
  91   //
  92
  93   if (fRunNumber!=10 && fIsLowEnergy==false) {
  94     AliError("Only settings for 2010 and the low energy runs are implemented so far");
  95   }
  96   if (fCentralityClass!="020" && fCentralityClass!="4080" && fCentralityClass!="0100"){
  97     AliError("Only settings for MB2010 are implemented so far");
  98   }
  99
 100   switch(decay) {
 101   case 1: // D0->Kpi
 102     if (fCollisionType==0) {
 103       if (fIsLowEnergy) InitD0toKpi2010ppLowEn();
 104       else InitD0toKpi2010pp();
 105     } else if (fCollisionType==1) {
 106       if (fCentralityClass=="020") InitD0toKpi2010PbPb020();
 107       else if (fCentralityClass=="4080") InitD0toKpi2010PbPb4080();
 108       else AliError("Not yet implemented");
 109     }
 110     //    else if (fCollisionType==2) InitD0toKpi2010ppLowEn();
 111     break;
 112   case 2: // D+->Kpipi
 113     if (fCollisionType==0) {
 114       if (fIsLowEnergy) InitDplustoKpipi2010ppLowEn();
 115       else InitDplustoKpipi2010pp();
 116     } else if (fCollisionType==1) {
 117       if (fCentralityClass=="020") InitDplustoKpipi2010PbPb020();
 118       else if (fCentralityClass=="4080") InitDplustoKpipi2010PbPb4080();
 119       else AliError("Not yet implemented");
 120     }
 121     break;
 122   case 3: // D*->D0pi
 123     if (fCollisionType==0) {
 124       if(fIsLowEnergy)  InitDstartoD0pi2010ppLowEn();
 125       else InitDstartoD0pi2010pp();
 126     }else if (fCollisionType==1) {
 127       if (fCentralityClass=="020")  InitDstartoD0pi2010PbPb020();
 128       if (fCentralityClass=="2040") InitDstartoD0pi2010PbPb2040();
 129       if (fCentralityClass=="4080") InitDstartoD0pi2010PbPb4080();
 130       if (fCentralityClass!="4080" && fCentralityClass!="2040" && fCentralityClass!="020")  AliError("Not yet implemented");
 131     }
 132     break;
 133   default:
 134     printf("Invalid decay type: %d\n",decay);
 135     break;
 136   }
 137
 138 }
 139
 140 //--------------------------------------------------------------------------
 141 void AliHFSystErr::InitD0toKpi2010pp() {
 142   //
 143   // D0->Kpi syst errors. Responsible: A. Rossi
 144   //   2010 pp sample
 145   //
 146
 147   // Normalization
 148   fNorm = new TH1F("fNorm","fNorm",24,0,24);
 149   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 150
 151   // Branching ratio
 152   fBR = new TH1F("fBR","fBR",24,0,24);
 153   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
 154
 155   // Tracking efficiency
 156   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);
 157   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.02); // 2% (1% per track)
 158
 159   // Raw yield extraction
 160   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
 161   fRawYield->SetBinContent(1,1);
 162   fRawYield->SetBinContent(2,0.2);
 163   for(Int_t i=3;i<=16;i++) fRawYield->SetBinContent(i,0.05);
 164   for(Int_t i=13;i<=16;i++) fRawYield->SetBinContent(i,0.10);
 165   for(Int_t i=17;i<=24;i++) fRawYield->SetBinContent(i,1);
 166
 167   // Cuts efficiency (from cuts variation)
 168   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
 169   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 170
 171   // PID efficiency (from PID/noPID)
 172   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
 173   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.03); // 3%
 174   fPIDEff->SetBinContent(2,0.10); // 10%
 175
 176   // MC dN/dpt
 177   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
 178   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0);
 179   fMCPtShape->SetBinContent(1,0.03);
 180   fMCPtShape->SetBinContent(2,0.03);
 181
 182   // particle-antiparticle
 183   //  fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",24,0,24);
 184   //  fPartAntipart->SetBinContent(1,1);
 185   //  for(Int_t i=2;i<=24;i++) fPartAntipart->SetBinContent(i,0.05);
 186
 187   return;
 188 }
 189 //--------------------------------------------------------------------------
 190 void AliHFSystErr::InitD0toKpi2010PbPb020() {
 191   //
 192   // D0->Kpi syst errors. Responsible: A. Rossi
 193   //   2010 PbPb sample, 0-20 CC
 194   //
 195
 196   // Normalization
 197   fNorm = new TH1F("fNorm","fNorm",20,0,20);
 198   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 199
 200   // Branching ratio
 201   fBR = new TH1F("fBR","fBR",20,0,20);
 202   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
 203
 204   // Tracking efficiency
 205   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
 206   for(Int_t i=2;i<=12;i++) fTrackingEff->SetBinContent(i,0.05+0.005*(Float_t)i);
 207
 208
 209   // Raw yield extraction
 210   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 211   fRawYield->SetBinContent(1,0);
 212   fRawYield->SetBinContent(2,0);
 213   fRawYield->SetBinContent(3,0.05);
 214   fRawYield->SetBinContent(4,0.05);
 215   fRawYield->SetBinContent(5,0.10);
 216   fRawYield->SetBinContent(6,0.10);
 217   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05);
 218  // for(Int_t i=13;i<=20;i++) fRawYield->SetBinContent(i,0);
 219
 220   // Cuts efficiency (from cuts variation)
 221   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 222   fCutsEff->SetBinContent(1,0.);
 223   fCutsEff->SetBinContent(2,0.);
 224   fCutsEff->SetBinContent(3,0.18);
 225   fCutsEff->SetBinContent(4,0.18);
 226   fCutsEff->SetBinContent(5,0.18);
 227   fCutsEff->SetBinContent(6,0.15);
 228   fCutsEff->SetBinContent(7,0.15);
 229   fCutsEff->SetBinContent(8,0.15);
 230 fCutsEff->SetBinContent(9,0.15);
 231  fCutsEff->SetBinContent(10,0.15);
 232   fCutsEff->SetBinContent(11,0.15);
 233   fCutsEff->SetBinContent(12,0.15);
 234  // for(Int_t i=13;i<=20;i++) fCutsEff->SetBinContent(i,0.);
 235
 236   // PID efficiency (from PID/noPID)
 237   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
 238   for(Int_t i=3;i<=12;i++) fPIDEff->SetBinContent(i,0.05);
 239   fPIDEff->SetBinContent(5,0.10);
 240   fPIDEff->SetBinContent(6,0.10);
 241
 242   // MC dN/dpt
 243   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 244   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0);
 245   fMCPtShape->SetBinContent(1,0.03);
 246   fMCPtShape->SetBinContent(2,0.03);
 247
 248
 249   // particle-antiparticle
 250   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
 251   for(Int_t i=3;i<=12;i++) fPartAntipart->SetBinContent(i,0.05);
 252   fPartAntipart->SetBinContent(3,0.10);
 253   fPartAntipart->SetBinContent(4,0.10);
 254   fPartAntipart->SetBinContent(7,0.10);
 255   fPartAntipart->SetBinContent(8,0.10);
 256
 257   return;
 258 }
 259 //--------------------------------------------------------------------------
 260 void AliHFSystErr::InitD0toKpi2010PbPb4080() {
 261   //
 262   // D0->Kpi syst errors. Responsible: A. Rossi
 263   //   2010 PbPb sample, 40-80 CC
 264   //
 265
 266   // Normalization
 267   fNorm = new TH1F("fNorm","fNorm",20,0,20);
 268   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 269
 270   // Branching ratio
 271   fBR = new TH1F("fBR","fBR",20,0,20);
 272   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
 273
 274   // Tracking efficiency
 275   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
 276   for(Int_t i=2;i<=12;i++) fTrackingEff->SetBinContent(i,0.5*(0.05+0.005*(Float_t)i));
 277
 278
 279   // Raw yield extraction
 280   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 281   fRawYield->SetBinContent(1,0);
 282   fRawYield->SetBinContent(2,0);
 283   for(Int_t i=3;i<=12;i++) fRawYield->SetBinContent(i,0.05);
 284   //for(Int_t i=13;i<=20;i++) fRawYield->SetBinContent(i,0);
 285
 286   // Cuts efficiency (from cuts variation)
 287   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 288   fCutsEff->SetBinContent(1,0.);
 289   fCutsEff->SetBinContent(2,0.);
 290   fCutsEff->SetBinContent(3,0.15);
 291   fCutsEff->SetBinContent(4,0.15);
 292   fCutsEff->SetBinContent(5,0.15);
 293   fCutsEff->SetBinContent(6,0.15);
 294   fCutsEff->SetBinContent(7,0.15);
 295   fCutsEff->SetBinContent(8,0.15);
 296   fCutsEff->SetBinContent(9,0.15);
 297   fCutsEff->SetBinContent(10,0.15);
 298   fCutsEff->SetBinContent(11,0.15);
 299   fCutsEff->SetBinContent(12,0.15);
 300  // for(Int_t i=13;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<=12;i++) fPIDEff->SetBinContent(i,0.05);
 305
 306   // MC dN/dpt
 307   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 308   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0);
 309   fMCPtShape->SetBinContent(1,0.03);
 310   fMCPtShape->SetBinContent(2,0.03);
 311
 312   // particle-antiparticle
 313   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
 314   for(Int_t i=3;i<=12;i++) fPartAntipart->SetBinContent(i,0.05);
 315
 316   return;
 317 }
 318
 319 //--------------------------------------------------------------------------
 320 void AliHFSystErr::InitD0toKpi2010ppLowEn() {
 321   //
 322   // D0->Kpi syst errors. Low energy run
 323   //   2011 2.76 TeV pp sample
 324   //
 325   AliInfo(" Settings for D0 --> K pi, p-p collisions at 2.76 TeV");
 326
 327   // Normalization
 328   fNorm = new TH1F("fNorm","fNorm",20,0,20);
 329   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 330
 331   // Branching ratio
 332   fBR = new TH1F("fBR","fBR",20,0,20);
 333   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
 334
 335   // Tracking efficiency
 336   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
 337   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.10); //10% (5% per track)
 338
 339   // Raw yield extraction
 340   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 341   fRawYield->SetBinContent(1,1);
 342   fRawYield->SetBinContent(2,1);
 343   fRawYield->SetBinContent(3,0.2);
 344   fRawYield->SetBinContent(4,0.2);
 345   fRawYield->SetBinContent(5,0.1);
 346   fRawYield->SetBinContent(6,0.1);
 347   fRawYield->SetBinContent(7,0.2);
 348   fRawYield->SetBinContent(8,0.2);
 349   for(Int_t i=9;i<=20;i++) fRawYield->SetBinContent(i,0.065);
 350
 351   // Cuts efficiency (from cuts variation)
 352   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 353   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.18); // 10%
 354
 355   // PID efficiency (from PID/noPID)
 356   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
 357   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.20); // 10%
 358   fPIDEff->SetBinContent(3,0.25); // 10%
 359   fPIDEff->SetBinContent(4,0.25); // 10%
 360   fPIDEff->SetBinContent(7,0.25); // 10%
 361   fPIDEff->SetBinContent(8,0.25); // 10%
 362
 363   // MC dN/dpt
 364   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 365   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0);
 366   fMCPtShape->SetBinContent(1,0.03);
 367   fMCPtShape->SetBinContent(2,0.03);
 368
 369   // particle-antiparticle
 370   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
 371   fPartAntipart->SetBinContent(1,1);
 372   fPartAntipart->SetBinContent(2,1);
 373   for(Int_t i=3;i<=6;i++) fPartAntipart->SetBinContent(i,0.08);
 374
 375   return;
 376 }
 377
 378 //--------------------------------------------------------------------------
 379 void AliHFSystErr::InitDplustoKpipi2010pp() {
 380   //
 381   // D+->Kpipi syst errors. Responsible: R. Bala
 382   //  2010 pp sample
 383   //
 384
 385
 386 // Normalization
 387   fNorm = new TH1F("fNorm","fNorm",24,0,24);
 388   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.07); // 10% error on sigmaV0and
 389
 390   // Branching ratio
 391   fBR = new TH1F("fBR","fBR",24,0,24);
 392   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.04); // 4% PDG2010
 393
 394   // Tracking efficiency
 395   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);
 396   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.03); // 3% (1% per track)
 397
 398
 399   // Raw yield extraction
 400   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
 401   fRawYield->SetBinContent(1,0.25);
 402   fRawYield->SetBinContent(2,0.25);
 403   fRawYield->SetBinContent(3,0.25);
 404   fRawYield->SetBinContent(4,0.25);
 405   fRawYield->SetBinContent(5,0.09);
 406   fRawYield->SetBinContent(6,0.09);
 407   for(Int_t i=7;i<=12;i++) fRawYield->SetBinContent(i,0.05);  //5 to 10%
 408   for(Int_t i=12;i<=24;i++) fRawYield->SetBinContent(i,0.10);  // 15%
 409
 410   // Cuts efficiency (from cuts variation)
 411   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
 412   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 413
 414   // PID efficiency (from PID/noPID)
 415   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
 416   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.05); // 5%
 417   fPIDEff->SetBinContent(1,0.15); // 15%
 418   fPIDEff->SetBinContent(2,0.15); // 15%
 419   fPIDEff->SetBinContent(3,0.15); // 15%
 420   fPIDEff->SetBinContent(4,0.15); // 15%
 421   for(Int_t i=12;i<=16;i++) fPIDEff->SetBinContent(i,0.10); // 5%
 422
 423   // MC dN/dpt
 424   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
 425   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0);
 426   fMCPtShape->SetBinContent(1,0.03);
 427   fMCPtShape->SetBinContent(2,0.03);
 428
 429
 430   // particle-antiparticle
 431   /*
 432   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
 433   fPartAntipart->SetBinContent(1,1);
 434   fPartAntipart->SetBinContent(2,1);
 435   fPartAntipart->SetBinContent(3,0.12);
 436   for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
 437   */
 438   return;
 439 }
 440
 441 //--------------------------------------------------------------------------
 442 void AliHFSystErr::InitDplustoKpipi2010PbPb020() {
 443   //
 444   // D+->Kpipi syst errors. Responsible: ??
 445   //  2010 PbPb sample, 0-20 CC
 446   //
 447
 448  // Normalization
 449   fNorm = new TH1F("fNorm","fNorm",20,0,20);
 450   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 451
 452   // Branching ratio
 453   fBR = new TH1F("fBR","fBR",20,0,20);
 454   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.04); // 4% PDG2010
 455
 456   // Tracking efficiency
 457   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
 458   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,i,0.075+0.005*(Float_t)i);
 459
 460
 461   // Raw yield extraction
 462   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 463   for(Int_t i=1;i<=20;i++) fRawYield->SetBinContent(i,1);  //5 to 10%
 464 //   for(Int_t i=6;i<=7;i++) fRawYield->SetBinContent(i,0.32);
 465 //   for(Int_t i=7;i<=9;i++) fRawYield->SetBinContent(i,0.32);
 466 //   for(Int_t i=9;i<=12;i++) fRawYield->SetBinContent(i,0.05);//31./156.);
 467   fRawYield->SetBinContent(6,0.12);
 468   fRawYield->SetBinContent(7,0.17);
 469   fRawYield->SetBinContent(8,0.17);
 470   fRawYield->SetBinContent(9,0.08);
 471   fRawYield->SetBinContent(10,0.08);
 472   fRawYield->SetBinContent(11,0.08);
 473   fRawYield->SetBinContent(12,0.08);
 474
 475   // Cuts efficiency (from cuts variation)
 476   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 477   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.15); // 10%
 478
 479   // PID efficiency (from PID/noPID)
 480   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
 481   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
 482   fPIDEff->SetBinContent(3,0.13); // 13%
 483
 484
 485   // MC dN/dpt  (copied from D0 : will update later)
 486   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 487   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.);
 488   fMCPtShape->SetBinContent(1,0.03);
 489   fMCPtShape->SetBinContent(2,0.03);
 490
 491
 492   // particle-antiparticle
 493   /*
 494   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
 495   fPartAntipart->SetBinContent(1,1);
 496   fPartAntipart->SetBinContent(2,1);
 497   fPartAntipart->SetBinContent(3,0.12);
 498   for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
 499   */
 500   return;
 501 }
 502
 503 //--------------------------------------------------------------------------
 504 void AliHFSystErr::InitDplustoKpipi2010PbPb4080() {
 505   //
 506   // D+->Kpipi syst errors. Responsible: ??
 507   //  2010 PbPb sample, 40-80 CC
 508   //
 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.10); // 10% error on sigmaV0and
 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.04); // 4% 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.03); // 3% (1% per track)
 522
 523
 524   // Raw yield extraction
 525   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 526   fRawYield->SetBinContent(1,1);
 527   fRawYield->SetBinContent(2,1);
 528   fRawYield->SetBinContent(3,1);
 529   fRawYield->SetBinContent(4,1.);
 530   fRawYield->SetBinContent(5,4./23.);
 531   fRawYield->SetBinContent(6,4./23.);
 532   fRawYield->SetBinContent(7,5./20.);
 533   fRawYield->SetBinContent(8,5./20.);
 534   for(Int_t i=9;i<=12;i++) fRawYield->SetBinContent(i,1./16.);
 535   for(Int_t i=13;i<=20;i++) fRawYield->SetBinContent(i,1);  //5 to 10%
 536
 537   // Cuts efficiency (from cuts variation)
 538   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 539   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 540
 541   // PID efficiency (from PID/noPID)
 542   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
 543   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
 544   fPIDEff->SetBinContent(3,0.13); // 13%
 545
 546
 547   // MC dN/dpt  (copied from D0 : will update later)
 548   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 549   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0);
 550   fMCPtShape->SetBinContent(1,0.03);
 551   fMCPtShape->SetBinContent(2,0.03);
 552
 553
 554   // particle-antiparticle
 555   /*
 556   fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
 557   fPartAntipart->SetBinContent(1,1);
 558   fPartAntipart->SetBinContent(2,1);
 559   fPartAntipart->SetBinContent(3,0.12);
 560   for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);   //5 to 12%
 561   */
 562   return;
 563 }
 564
 565 //--------------------------------------------------------------------------
 566 void AliHFSystErr::InitDplustoKpipi2010ppLowEn() {
 567
 568   //
 569   // D+->Kpipi syst errors. Responsible: R. Bala
 570   //  2011 2.76 TeV pp sample
 571   //
 572   AliInfo(" Settings for D+ --> K pi pi p-p collisions at 2.76 TeV");
 573
 574   // Normalization
 575   fNorm = new TH1F("fNorm","fNorm",20,0,20);
 576   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 577
 578   // Branching ratio
 579   fBR = new TH1F("fBR","fBR",20,0,20);
 580   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.04); // 4% PDG2010
 581
 582   // Tracking efficiency
 583   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
 584   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.15); // 3% (1% per track)
 585
 586   // Raw yield extraction
 587   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 588   fRawYield->SetBinContent(1,1);
 589   fRawYield->SetBinContent(2,1);
 590   fRawYield->SetBinContent(3,0.20);
 591   fRawYield->SetBinContent(4,0.20);
 592   for(Int_t i=5;i<=20;i++) fRawYield->SetBinContent(i,0.055);  //5 to 10%
 593
 594   // Cuts efficiency (from cuts variation)
 595   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 596   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 597
 598   // PID efficiency (from PID/noPID)
 599   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
 600   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
 601   fPIDEff->SetBinContent(3,0.15); // 13%
 602   fPIDEff->SetBinContent(4,0.15); // 13%
 603
 604   // MC dN/dpt  (copied from D0 : will update later)
 605   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 606   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.);
 607   fMCPtShape->SetBinContent(1,0.03);
 608   fMCPtShape->SetBinContent(2,0.03);
 609
 610   return;
 611 }
 612
 613 //--------------------------------------------------------------------------
 614 void AliHFSystErr::InitDstartoD0pi2010pp() {
 615   //
 616   // D*+->D0pi syst errors. Responsible: A. Grelli, Y. Wang
 617   //  2010 pp sample
 618   //
 619
 620  // Normalization
 621   fNorm = new TH1F("fNorm","fNorm",24,0,24);
 622   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 623
 624   // Branching ratio
 625   fBR = new TH1F("fBR","fBR",24,0,24);
 626   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
 627
 628   // Tracking efficiency
 629   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);
 630   fTrackingEff->SetBinContent(1,1.0);
 631   fTrackingEff->SetBinContent(2,0.12);
 632   fTrackingEff->SetBinContent(3,0.08);
 633   fTrackingEff->SetBinContent(3,0.05);
 634   for(Int_t i=4;i<=24;i++) fTrackingEff->SetBinContent(i,0.03); // 3% (1% per track)
 635
 636
 637   // Raw yield extraction
 638   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
 639   fRawYield->SetBinContent(1,1.0);
 640   fRawYield->SetBinContent(2,0.09);
 641   fRawYield->SetBinContent(3,0.04);
 642   fRawYield->SetBinContent(4,0.02);
 643   fRawYield->SetBinContent(5,0.03);
 644   fRawYield->SetBinContent(6,0.06);
 645   fRawYield->SetBinContent(7,0.04);
 646   fRawYield->SetBinContent(8,0.11);
 647   for(Int_t i=5;i<=24;i++) fRawYield->SetBinContent(i,0.03);  //4%
 648   fRawYield->SetBinContent(13,0.09);
 649   fRawYield->SetBinContent(14,0.09);
 650   fRawYield->SetBinContent(15,0.09);
 651   fRawYield->SetBinContent(16,0.09);
 652   fRawYield->SetBinContent(17,0.24);
 653   fRawYield->SetBinContent(18,0.24);
 654   fRawYield->SetBinContent(19,0.24);
 655   fRawYield->SetBinContent(20,0.24);
 656   fRawYield->SetBinContent(21,0.24);
 657   fRawYield->SetBinContent(22,0.24);
 658   fRawYield->SetBinContent(23,0.24);
 659   fRawYield->SetBinContent(24,0.24);
 660
 661   // Cuts efficiency (from cuts variation)
 662   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
 663   fCutsEff->SetBinContent(2,0.22);
 664   for(Int_t i=3;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 665
 666   // PID efficiency (from PID/noPID)
 667   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
 668   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.04); // 3%
 669
 670
 671   // MC dN/dpt  (copied from D0 : will update later)
 672   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
 673   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0);
 674   fMCPtShape->SetBinContent(1,0.03);
 675   fMCPtShape->SetBinContent(2,0.03);
 676
 677   return;
 678
 679
 680 }
 681 //--------------------------------------------------------------------------
 682 void AliHFSystErr::InitDstartoD0pi2010ppLowEn() {
 683
 684   //
 685   // D+->Kpipi syst errors. Responsible: A. Grelli
 686   //  2011 2.76 TeV pp sample
 687   //
 688   AliInfo(" Settings for D*+ --> D0 pi p-p collisions at 2.76 TeV");
 689
 690 // Normalization
 691   fNorm = new TH1F("fNorm","fNorm",20,0,20);
 692   for(Int_t i=1;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 693
 694   // Branching ratio
 695   fBR = new TH1F("fBR","fBR",20,0,20);
 696   for(Int_t i=1;i<=20;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
 697
 698   // Tracking efficiency
 699   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
 700   for(Int_t i=1;i<=20;i++) fTrackingEff->SetBinContent(i,0.15); //10% (to be checked!!)
 701
 702   // Raw yield extraction
 703   fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
 704   fRawYield->SetBinContent(1,1);
 705   fRawYield->SetBinContent(2,1);
 706   fRawYield->SetBinContent(3,0.2);
 707   fRawYield->SetBinContent(4,0.2);
 708   fRawYield->SetBinContent(5,0.08);
 709   fRawYield->SetBinContent(6,0.08);
 710   fRawYield->SetBinContent(7,0.1);
 711   fRawYield->SetBinContent(8,0.1);
 712   fRawYield->SetBinContent(9,0.2);
 713   fRawYield->SetBinContent(10,0.2);
 714   for(Int_t i=9;i<=20;i++) fRawYield->SetBinContent(i,0.065);
 715
 716   // Cuts efficiency (from cuts variation)
 717   fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
 718   for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.18);
 719   fCutsEff->SetBinContent(3,0.2);
 720   fCutsEff->SetBinContent(4,0.2);
 721   fCutsEff->SetBinContent(5,0.1);
 722   fCutsEff->SetBinContent(6,0.1);
 723   fCutsEff->SetBinContent(7,0.12);
 724   fCutsEff->SetBinContent(8,0.12);
 725   fCutsEff->SetBinContent(9,0.2);
 726   fCutsEff->SetBinContent(10,0.2);
 727   fCutsEff->SetBinContent(11,0.2);
 728   fCutsEff->SetBinContent(12,0.2);
 729
 730   // PID efficiency (from PID/noPID)
 731   fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
 732   for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 10%
 733
 734   // MC dN/dpt
 735   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
 736   for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,0.);
 737   fMCPtShape->SetBinContent(1,0.03);
 738   fMCPtShape->SetBinContent(2,0.03);
 739
 740
 741   return;
 742 }
 743 //------------------------------------------------------------------------
 744 void AliHFSystErr::InitDstartoD0pi2010PbPb020() {
 745   //
 746   // D*+->D0pi syst errors. Responsible: A. Grelli
 747   //  2010 PbPb sample, 0-20 CC
 748   //
 749
 750   AliInfo(" Settings for D*+ --> D0pi Pb-Pb collisions at 2.76 TeV - 0-20 centrality - DUMMY");
 751
 752  // Normalization
 753   fNorm = new TH1F("fNorm","fNorm",24,0,24);
 754   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 755
 756   // Branching ratio
 757   fBR = new TH1F("fBR","fBR",24,0,24);
 758   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
 759
 760   // Tracking efficiency
 761   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);;
 762   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.12); // 3% (1% per track)
 763
 764
 765   // Raw yield extraction
 766   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
 767   for(Int_t i=1;i<=24;i++) fRawYield->SetBinContent(i,0.15);  //4%
 768
 769   // Cuts efficiency (from cuts variation)
 770   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
 771   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 772
 773   // PID efficiency (from PID/noPID)
 774   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
 775   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.04); // 3%
 776
 777
 778   // MC dN/dpt  (copied from D0 : will update later)
 779   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
 780   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0.);
 781   fMCPtShape->SetBinContent(1,0.03);
 782   fMCPtShape->SetBinContent(2,0.03);
 783
 784   return;
 785
 786 }
 787 //-------------------------------------------------------------------------
 788 void AliHFSystErr::InitDstartoD0pi2010PbPb2040() {
 789   //
 790   // D*+->D0pi syst errors. Responsible: A. Grelli
 791   //  2010 PbPb sample, 20-40 CC
 792   //
 793
 794   AliInfo(" Settings for D*+ --> D0pi Pb-Pb collisions at 2.76 TeV - 20-40 centrality - DUMMY");
 795
 796  // Normalization
 797   fNorm = new TH1F("fNorm","fNorm",24,0,24);
 798   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 799
 800   // Branching ratio
 801   fBR = new TH1F("fBR","fBR",24,0,24);
 802   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
 803
 804   // Tracking efficiency
 805   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);;
 806   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.12); // 3% (1% per track)
 807
 808
 809   // Raw yield extraction
 810   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
 811   for(Int_t i=1;i<=24;i++) fRawYield->SetBinContent(i,0.15);  //4%
 812
 813   // Cuts efficiency (from cuts variation)
 814   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
 815   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 816
 817   // PID efficiency (from PID/noPID)
 818   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
 819   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.04); // 3%
 820
 821
 822   // MC dN/dpt  (copied from D0 : will update later)
 823   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
 824   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0.);
 825   fMCPtShape->SetBinContent(1,0.03);
 826   fMCPtShape->SetBinContent(2,0.03);
 827
 828   return;
 829
 830 }
 831 //--------------------------------------------------------------------------
 832 void AliHFSystErr::InitDstartoD0pi2010PbPb4080() {
 833   //
 834   // D*+->D0pi syst errors. Responsible: A. Grelli
 835   //  2010 PbPb sample, 40-80 CC
 836   //
 837
 838   AliInfo(" Settings for D*+ --> D0pi Pb-Pb collisions at 2.76 TeV - 40-80 centrality - DUMMY");
 839
 840  // Normalization
 841   fNorm = new TH1F("fNorm","fNorm",24,0,24);
 842   for(Int_t i=1;i<=24;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
 843
 844   // Branching ratio
 845   fBR = new TH1F("fBR","fBR",24,0,24);
 846   for(Int_t i=1;i<=24;i++) fBR->SetBinContent(i,0.015); // 1.5% PDG2010
 847
 848   // Tracking efficiency
 849   fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",24,0,24);;
 850   for(Int_t i=1;i<=24;i++) fTrackingEff->SetBinContent(i,0.12); // 3% (1% per track)
 851
 852
 853   // Raw yield extraction
 854   fRawYield = new TH1F("fRawYield","fRawYield",24,0,24);
 855   for(Int_t i=1;i<=24;i++) fRawYield->SetBinContent(i,0.15);  //4%
 856
 857   // Cuts efficiency (from cuts variation)
 858   fCutsEff = new TH1F("fCutsEff","fCutsEff",24,0,24);
 859   for(Int_t i=1;i<=24;i++) fCutsEff->SetBinContent(i,0.10); // 10%
 860
 861   // PID efficiency (from PID/noPID)
 862   fPIDEff = new TH1F("fPIDEff","fPIDEff",24,0,24);
 863   for(Int_t i=1;i<=24;i++) fPIDEff->SetBinContent(i,0.04); // 3%
 864
 865
 866   // MC dN/dpt  (copied from D0 : will update later)
 867   fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",24,0,24);
 868   for(Int_t i=1;i<=24;i++) fMCPtShape->SetBinContent(i,0.);
 869   fMCPtShape->SetBinContent(1,0.03);
 870   fMCPtShape->SetBinContent(2,0.03);
 871
 872   return;
 873
 874 }
 875 //--------------------------------------------------------------------------
 876 Double_t AliHFSystErr::GetCutsEffErr(Double_t pt) const {
 877   //
 878   // Get error
 879   //
 880
 881   Int_t bin=fCutsEff->FindBin(pt);
 882
 883   return fCutsEff->GetBinContent(bin);
 884 }
 885 //--------------------------------------------------------------------------
 886 Double_t AliHFSystErr::GetMCPtShapeErr(Double_t pt) const {
 887   //
 888   // Get error
 889   //
 890
 891   Int_t bin=fMCPtShape->FindBin(pt);
 892
 893   return fMCPtShape->GetBinContent(bin);
 894 }
 895 //--------------------------------------------------------------------------
 896 Double_t AliHFSystErr::GetSeleEffErr(Double_t pt) const {
 897   //
 898   // Get error
 899   //
 900
 901   Double_t err=GetCutsEffErr(pt)*GetCutsEffErr(pt)+GetMCPtShapeErr(pt)*GetMCPtShapeErr(pt);
 902
 903   return TMath::Sqrt(err);
 904 }
 905 //--------------------------------------------------------------------------
 906 Double_t AliHFSystErr::GetPIDEffErr(Double_t pt) const {
 907   //
 908   // Get error
 909   //
 910
 911   Int_t bin=fPIDEff->FindBin(pt);
 912
 913   return fPIDEff->GetBinContent(bin);
 914 }
 915 //--------------------------------------------------------------------------
 916 Double_t AliHFSystErr::GetTrackingEffErr(Double_t pt) const {
 917   //
 918   // Get error
 919   //
 920
 921   Int_t bin=fTrackingEff->FindBin(pt);
 922
 923   return fTrackingEff->GetBinContent(bin);
 924 }
 925 //--------------------------------------------------------------------------
 926 Double_t AliHFSystErr::GetRawYieldErr(Double_t pt) const {
 927   //
 928   // Get error
 929   //
 930
 931   Int_t bin=fRawYield->FindBin(pt);
 932
 933   return fRawYield->GetBinContent(bin);
 934 }
 935 //--------------------------------------------------------------------------
 936 Double_t AliHFSystErr::GetPartAntipartErr(Double_t pt) const {
 937   //
 938   // Get error
 939   //
 940
 941   Int_t bin=fPartAntipart->FindBin(pt);
 942
 943   return fPartAntipart->GetBinContent(bin);
 944 }
 945 //--------------------------------------------------------------------------
 946 Double_t AliHFSystErr::GetTotalSystErr(Double_t pt,Double_t feeddownErr) const {
 947   //
 948   // Get total syst error (except norm. error)
 949   //
 950
 951   Double_t err=0.;
 952
 953   if(fRawYield) err += GetRawYieldErr(pt)*GetRawYieldErr(pt);
 954   if(fTrackingEff) err += GetTrackingEffErr(pt)*GetTrackingEffErr(pt);
 955   if(fBR) err += GetBRErr()*GetBRErr();
 956   if(fCutsEff) err += GetCutsEffErr(pt)*GetCutsEffErr(pt);
 957   if(fPIDEff) err += GetPIDEffErr(pt)*GetPIDEffErr(pt);
 958   if(fMCPtShape) err += GetMCPtShapeErr(pt)*GetMCPtShapeErr(pt);
 959   if(fPartAntipart) err += GetPartAntipartErr(pt)*GetPartAntipartErr(pt);
 960
 961   err += feeddownErr*feeddownErr;
 962
 963   return TMath::Sqrt(err);
 964 }
 965 //---------------------------------------------------------------------------
 966 void AliHFSystErr::DrawErrors(TGraphAsymmErrors *grErrFeeddown) const {
 967   //
 968   // Draw errors
 969   //
 970   gStyle->SetOptStat(0);
 971
 972   TCanvas *cSystErr = new TCanvas("cSystErr","Systematic Errors",300,80,640,500);
 973   cSystErr->Range(0.20,-0.5,18.4,0.34);
 974   cSystErr->SetRightMargin(0.318);
 975   cSystErr->SetFillColor(0);
 976
 977   TH2F *hFrame = new TH2F("hFrame","Systematic errors; p_{t} [GeV/c]; Relative Error",20,0,20,100,-1,+1);
 978   hFrame->SetAxisRange(1.,15.9,"X");
 979   hFrame->SetAxisRange(-0.5,0.5,"Y");
 980   hFrame->Draw();
 981
 982   TLegend *leg = new TLegend(0.69,0.44,0.98,0.86,NULL,"brNDC");
 983   leg->SetTextSize(0.03601695);
 984   leg->SetFillStyle(0);
 985   leg->SetBorderSize(0);
 986
 987   TH1F *hTotErr=new TH1F("hTotErr","",24,0,24);
 988   Int_t nbins = fNorm->GetNbinsX();
 989   TGraphAsymmErrors *gTotErr = new TGraphAsymmErrors(nbins);
 990   for(Int_t i=1;i<=24;i++) {
 991     Double_t pt = hTotErr->GetBinCenter(i);
 992     Double_t ptwidth = hTotErr->GetBinWidth(i);
 993
 994     if(grErrFeeddown) {
 995       Double_t x=0., y=0., errxl=0., errxh=0., erryl=0., erryh=0.;
 996       Double_t toterryl=0., toterryh=0.;
 997       for(Int_t j=0; j<grErrFeeddown->GetN(); j++) {
 998         grErrFeeddown->GetPoint(j,x,y);
 999         errxh = grErrFeeddown->GetErrorXhigh(j);
1000         errxl = grErrFeeddown->GetErrorXlow(j);
1001         if ( ( (x-errxl) <= pt) && ( (x+errxl) >= pt) ) {
1002           erryh = grErrFeeddown->GetErrorYhigh(j);
1003           erryl = grErrFeeddown->GetErrorYlow(j);
1004         }
1005       }
1006       if (erryl>=1e-3) toterryl = GetTotalSystErr(pt,erryl);
1007       else toterryl = GetTotalSystErr(pt);
1008       if (erryh>=1e-3) toterryh = GetTotalSystErr(pt,erryh);
1009       else toterryh = GetTotalSystErr(pt);
1010
1011       hTotErr->SetBinContent(i,toterryh);
1012       gTotErr->SetPoint(i,pt,0.);
1013       gTotErr->SetPointError(i,ptwidth/2.,ptwidth/2.,toterryl,toterryh); // i, exl, exh, eyl, eyh
1014     }
1015     else {
1016       hTotErr->SetBinContent(i,GetTotalSystErr(pt));
1017       gTotErr->SetPoint(i,pt,0.);
1018       gTotErr->SetPointError(i,ptwidth/2.,ptwidth/2.,GetTotalSystErr(pt),GetTotalSystErr(pt)); // i, exl, exh, eyl, eyh
1019     }
1020
1021   }
1022   gTotErr->SetLineColor(kBlack);
1023   gTotErr->SetFillColor(kRed);
1024   gTotErr->SetFillStyle(3002);
1025   gTotErr->Draw("2");
1026   leg->AddEntry(gTotErr,"Total (excl. norm.)","f");
1027 //   hTotErr->SetLineColor(1);
1028 //   hTotErr->SetLineWidth(3);
1029 //   hTotErr->Draw("same");
1030 //   leg->AddEntry(hTotErr,"Total (excl. norm.)","l");
1031
1032
1033   fNorm->SetFillColor(1);
1034   fNorm->SetFillStyle(3002);
1035   //fNorm->Draw("same");
1036   //TH1F *hNormRefl = ReflectHisto(fNorm);
1037   //hNormRefl->Draw("same");
1038   leg->AddEntry(fNorm,"Normalization (10%)","");
1039
1040   if(grErrFeeddown) {
1041     grErrFeeddown->SetFillColor(kTeal-8);
1042     grErrFeeddown->SetFillStyle(3001);
1043     grErrFeeddown->Draw("2");
1044     leg->AddEntry(grErrFeeddown,"Feed-down from B","f");
1045   }
1046   if(fTrackingEff) {
1047     fTrackingEff->SetFillColor(4);
1048     fTrackingEff->SetFillStyle(3006);
1049     fTrackingEff->Draw("same");
1050     TH1F *hTrackingEffRefl = ReflectHisto(fTrackingEff);
1051     hTrackingEffRefl->Draw("same");
1052     leg->AddEntry(fTrackingEff,"Tracking efficiency","f");
1053   }
1054   if(fBR) {
1055     fBR->SetFillColor(6);
1056     fBR->SetFillStyle(3005);
1057     //fBR->SetFillStyle(3020);
1058     fBR->Draw("same");
1059     TH1F *hBRRefl = ReflectHisto(fBR);
1060     hBRRefl->Draw("same");
1061     leg->AddEntry(fBR,"Branching ratio","f");
1062   }
1063   if(fRawYield) {
1064     Int_t ci;   // for color index setting
1065     ci = TColor::GetColor("#00cc00");
1066     fRawYield->SetLineColor(ci);
1067     //    fRawYield->SetLineColor(3);
1068     fRawYield->SetLineWidth(3);
1069     fRawYield->Draw("same");
1070     TH1F *hRawYieldRefl = ReflectHisto(fRawYield);
1071     hRawYieldRefl->Draw("same");
1072     leg->AddEntry(fRawYield,"Yield extraction","l");
1073   }
1074   if(fCutsEff) {
1075     fCutsEff->SetLineColor(4);
1076     fCutsEff->SetLineWidth(3);
1077     fCutsEff->Draw("same");
1078     TH1F *hCutsEffRefl = ReflectHisto(fCutsEff);
1079     hCutsEffRefl->Draw("same");
1080     leg->AddEntry(fCutsEff,"Cuts efficiency","l");
1081   }
1082   if(fPIDEff) {
1083     fPIDEff->SetLineColor(7);
1084     fPIDEff->SetLineWidth(3);
1085     fPIDEff->Draw("same");
1086     TH1F *hPIDEffRefl = ReflectHisto(fPIDEff);
1087     hPIDEffRefl->Draw("same");
1088     leg->AddEntry(fPIDEff,"PID efficiency","l");
1089   }
1090   if(fMCPtShape) {
1091     Int_t ci = TColor::GetColor("#9933ff");
1092     fMCPtShape->SetLineColor(ci);
1093     //    fMCPtShape->SetLineColor(8);
1094     fMCPtShape->SetLineWidth(3);
1095     fMCPtShape->Draw("same");
1096     TH1F *hMCPtShapeRefl = ReflectHisto(fMCPtShape);
1097     hMCPtShapeRefl->Draw("same");
1098     leg->AddEntry(fMCPtShape,"MC p_{t} shape","l");
1099   }
1100   if(fPartAntipart) {
1101     Int_t ci = TColor::GetColor("#ff6600");
1102     fPartAntipart->SetLineColor(ci);
1103     //    fPartAntipart->SetLineColor(9);
1104     fPartAntipart->SetLineWidth(3);
1105     fPartAntipart->Draw("same");
1106     TH1F *hPartAntipartRefl = ReflectHisto(fPartAntipart);
1107     hPartAntipartRefl->Draw("same");
1108     leg->AddEntry(fPartAntipart,"D = #bar{D}","l");
1109   }
1110
1111
1112   leg->Draw();
1113
1114   cSystErr->SaveAs("RelativeSystematics.eps");
1115
1116   return;
1117 }
1118 //-------------------------------------------------------------------------
1119 TH1F* AliHFSystErr::ReflectHisto(TH1F *hin) const {
1120   //
1121   // Clones and reflects histogram
1122   //
1123   TH1F *hout=(TH1F*)hin->Clone("hout");
1124   hout->Scale(-1.);
1125
1126   return hout;
1127 }