1 /**************************************************************************
2 * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 /////////////////////////////////////////////////////////////
18 // Class to handle systematic errors for charm hadrons
21 // AliHFSystEff syst(DECAY); // DECAY = 1 for D0, 2, for D+, 3 for D*
22 // syst.DrawErrors(); // to see a plot of the error contributions
23 // syst.GetTotalSystErr(pt); // to get the total err at pt
25 // Author: A.Dainese, andrea.dainese@pd.infn.it
26 /////////////////////////////////////////////////////////////
29 #include <TGraphAsymmErrors.h>
36 #include "AliHFSystErr.h"
39 ClassImp(AliHFSystErr)
41 //--------------------------------------------------------------------------
42 AliHFSystErr::AliHFSystErr(const Char_t* name, const Char_t* title) :
54 // Default Constructor
57 //--------------------------------------------------------------------------
58 AliHFSystErr::AliHFSystErr(Int_t decay,const Char_t* name, const Char_t* title) :
70 // Standard Constructor
84 printf("Invalid decay type: %d\n",decay);
90 //--------------------------------------------------------------------------
91 AliHFSystErr::~AliHFSystErr() {
96 if(fNorm) { delete fNorm; fNorm=0; }
97 if(fRawYield) { delete fRawYield; fRawYield=0; }
98 if(fTrackingEff) { delete fTrackingEff; fTrackingEff=0; }
99 if(fBR) { delete fBR; fBR=0; }
100 if(fCutsEff) { delete fCutsEff; fCutsEff=0; }
101 if(fPIDEff) { delete fPIDEff; fPIDEff=0; }
102 if(fMCPtShape) { delete fMCPtShape; fMCPtShape=0; }
103 if(fPartAntipart) { delete fPartAntipart; fPartAntipart=0; }
106 //--------------------------------------------------------------------------
107 void AliHFSystErr::InitD0toKpi() {
109 // D0->Kpi syst errors. Responsible: A. Rossi
113 fNorm = new TH1F("fNorm","fNorm",20,0,20);
114 for(Int_t i=4;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
117 fBR = new TH1F("fBR","fBR",20,0,20);
118 for(Int_t i=4;i<=20;i++) fBR->SetBinContent(i,0.012); // 1.2% PDG2010
120 // Tracking efficiency
121 fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
122 for(Int_t i=4;i<=20;i++) fTrackingEff->SetBinContent(i,0.02); // 2% (1% per track)
124 // Raw yield extraction
125 fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
126 fRawYield->SetBinContent(1,1);
127 fRawYield->SetBinContent(2,1);
128 fRawYield->SetBinContent(3,0.15);
129 for(Int_t i=4;i<=20;i++) fRawYield->SetBinContent(i,0.065);
131 // Cuts efficiency (from cuts variation)
132 fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
133 for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10); // 10%
135 // PID efficiency (from PID/noPID)
136 fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
137 for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.03); // 3%
138 fPIDEff->SetBinContent(4,0.15); // 15%
141 fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
142 for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,(Float_t)i*0.006);
144 // particle-antiparticle
145 fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
146 fPartAntipart->SetBinContent(1,1);
147 fPartAntipart->SetBinContent(2,1);
148 fPartAntipart->SetBinContent(3,0.20);
149 for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05);
153 //--------------------------------------------------------------------------
154 void AliHFSystErr::InitDplustoKpipi() {
156 // D+->Kpipi syst errors. Responsible: R. Bala
160 fNorm = new TH1F("fNorm","fNorm",20,0,20);
161 for(Int_t i=4;i<=20;i++) fNorm->SetBinContent(i,0.10); // 10% error on sigmaV0and
164 fBR = new TH1F("fBR","fBR",20,0,20);
165 for(Int_t i=4;i<=20;i++) fBR->SetBinContent(i,0.04); // 4% PDG2010
167 // Tracking efficiency
168 fTrackingEff = new TH1F("fTrackingEff","fTrackingEff",20,0,20);
169 for(Int_t i=4;i<=20;i++) fTrackingEff->SetBinContent(i,0.03); // 3% (1% per track)
172 // Raw yield extraction
173 fRawYield = new TH1F("fRawYield","fRawYield",20,0,20);
174 fRawYield->SetBinContent(1,1);
175 fRawYield->SetBinContent(2,1);
176 fRawYield->SetBinContent(3,0.10);
177 for(Int_t i=4;i<=20;i++) fRawYield->SetBinContent(i,0.055); //5 to 10%
179 // Cuts efficiency (from cuts variation)
180 fCutsEff = new TH1F("fCutsEff","fCutsEff",20,0,20);
181 for(Int_t i=1;i<=20;i++) fCutsEff->SetBinContent(i,0.10); // 10%
183 // PID efficiency (from PID/noPID)
184 fPIDEff = new TH1F("fPIDEff","fPIDEff",20,0,20);
185 for(Int_t i=1;i<=20;i++) fPIDEff->SetBinContent(i,0.05); // 5%
186 fPIDEff->SetBinContent(3,0.13); // 13%
189 // MC dN/dpt (copied from D0 : will update later)
190 fMCPtShape = new TH1F("fMCPtShape","fMCPtShape",20,0,20);
191 for(Int_t i=1;i<=20;i++) fMCPtShape->SetBinContent(i,(Float_t)i*0.006);
194 // particle-antiparticle
196 fPartAntipart = new TH1F("fPartAntipart","fPartAntipart",20,0,20);
197 fPartAntipart->SetBinContent(1,1);
198 fPartAntipart->SetBinContent(2,1);
199 fPartAntipart->SetBinContent(3,0.12);
200 for(Int_t i=4;i<=20;i++) fPartAntipart->SetBinContent(i,0.05); //5 to 12%
204 //--------------------------------------------------------------------------
205 void AliHFSystErr::InitDstartoD0pi() {
207 // D*+->D0pi syst errors. Responsible: tbd
212 //--------------------------------------------------------------------------
213 Double_t AliHFSystErr::GetCutsEffErr(Double_t pt) const {
218 Int_t bin=fCutsEff->FindBin(pt);
220 return fCutsEff->GetBinContent(bin);
222 //--------------------------------------------------------------------------
223 Double_t AliHFSystErr::GetMCPtShapeErr(Double_t pt) const {
228 Int_t bin=fMCPtShape->FindBin(pt);
230 return fMCPtShape->GetBinContent(bin);
232 //--------------------------------------------------------------------------
233 Double_t AliHFSystErr::GetSeleEffErr(Double_t pt) const {
238 Double_t err=GetCutsEffErr(pt)*GetCutsEffErr(pt)+GetMCPtShapeErr(pt)*GetMCPtShapeErr(pt);
240 return TMath::Sqrt(err);
242 //--------------------------------------------------------------------------
243 Double_t AliHFSystErr::GetPIDEffErr(Double_t pt) const {
248 Int_t bin=fPIDEff->FindBin(pt);
250 return fPIDEff->GetBinContent(bin);
252 //--------------------------------------------------------------------------
253 Double_t AliHFSystErr::GetTrackingEffErr(Double_t pt) const {
258 Int_t bin=fTrackingEff->FindBin(pt);
260 return fTrackingEff->GetBinContent(bin);
262 //--------------------------------------------------------------------------
263 Double_t AliHFSystErr::GetRawYieldErr(Double_t pt) const {
268 Int_t bin=fRawYield->FindBin(pt);
270 return fRawYield->GetBinContent(bin);
272 //--------------------------------------------------------------------------
273 Double_t AliHFSystErr::GetPartAntipartErr(Double_t pt) const {
278 Int_t bin=fPartAntipart->FindBin(pt);
280 return fPartAntipart->GetBinContent(bin);
282 //--------------------------------------------------------------------------
283 Double_t AliHFSystErr::GetTotalSystErr(Double_t pt,Double_t feeddownErr) const {
285 // Get total syst error (except norm. error)
290 if(fRawYield) err += GetRawYieldErr(pt)*GetRawYieldErr(pt);
291 if(fTrackingEff) err += GetTrackingEffErr(pt)*GetTrackingEffErr(pt);
292 if(fBR) err += GetBRErr()*GetBRErr();
293 if(fCutsEff) err += GetCutsEffErr(pt)*GetCutsEffErr(pt);
294 if(fPIDEff) err += GetPIDEffErr(pt)*GetPIDEffErr(pt);
295 if(fMCPtShape) err += GetMCPtShapeErr(pt)*GetMCPtShapeErr(pt);
296 if(fPartAntipart) err += GetPartAntipartErr(pt)*GetPartAntipartErr(pt);
298 err += feeddownErr*feeddownErr;
300 return TMath::Sqrt(err);
302 //---------------------------------------------------------------------------
303 void AliHFSystErr::DrawErrors(TGraphAsymmErrors *grErrFeeddown) const {
307 gStyle->SetOptStat(0);
309 TCanvas *cSystErr = new TCanvas("cSystErr","Systematic Errors",0,0,500,500);
310 cSystErr->SetFillColor(0);
312 TH2F *hFrame = new TH2F("hFrame","Systematic errors; p_{t} [GeV/c]; Relative Error",20,0,20,100,0,+1);
313 hFrame->SetAxisRange(2.,11.9,"X");
314 hFrame->SetAxisRange(0.,0.5,"Y");
317 TLegend *leg=new TLegend(0.5,0.5,0.9,0.9);
318 leg->SetFillStyle(0);
319 leg->SetBorderSize(0);
322 fNorm->SetFillColor(1);
323 fNorm->SetFillStyle(3002);
325 leg->AddEntry(fNorm,"Normalization","f");
328 grErrFeeddown->SetFillColor(kTeal-8);
329 grErrFeeddown->SetFillStyle(3001);
330 grErrFeeddown->Draw("2");
331 leg->AddEntry(grErrFeeddown,"Feed-down from B","f");
334 fTrackingEff->SetFillColor(4);
335 fTrackingEff->SetFillStyle(3006);
336 fTrackingEff->Draw("same");
337 leg->AddEntry(fTrackingEff,"Tracking efficiency","f");
340 fBR->SetFillColor(6);
341 fBR->SetFillStyle(3005);
342 //fBR->SetFillStyle(3020);
344 leg->AddEntry(fBR,"Branching ratio","f");
347 Int_t ci; // for color index setting
348 ci = TColor::GetColor("#00cc00");
349 fRawYield->SetLineColor(ci);
350 // fRawYield->SetLineColor(3);
351 fRawYield->SetLineWidth(3);
352 fRawYield->Draw("same");
353 leg->AddEntry(fRawYield,"Inv. mass analysis","l");
356 fCutsEff->SetLineColor(4);
357 fCutsEff->SetLineWidth(3);
358 fCutsEff->Draw("same");
359 leg->AddEntry(fCutsEff,"Cuts efficiency","l");
362 fPIDEff->SetLineColor(7);
363 fPIDEff->SetLineWidth(3);
364 fPIDEff->Draw("same");
365 leg->AddEntry(fPIDEff,"PID efficiency","l");
368 Int_t ci = TColor::GetColor("#9933ff");
369 fMCPtShape->SetLineColor(ci);
370 // fMCPtShape->SetLineColor(8);
371 fMCPtShape->SetLineWidth(3);
372 fMCPtShape->Draw("same");
373 leg->AddEntry(fMCPtShape,"MC p_{t} shape","l");
376 Int_t ci = TColor::GetColor("#ff6600");
377 fPartAntipart->SetLineColor(ci);
378 // fPartAntipart->SetLineColor(9);
379 fPartAntipart->SetLineWidth(3);
380 fPartAntipart->Draw("same");
381 leg->AddEntry(fPartAntipart,"D = #bar{D}","l");
385 TH1F *hTotErr=new TH1F("hTotErr","",20,0,20);
386 for(Int_t i=1;i<=20;i++) {
387 Double_t pt=hTotErr->GetBinCenter(i);
389 Double_t x=0., y=0., erryh = 0.;
390 for(Int_t j=0; j<grErrFeeddown->GetN(); j++) {
391 grErrFeeddown->GetPoint(j,x,y);
392 if ( TMath::Abs(x - pt) < 0.5) erryh = grErrFeeddown->GetErrorYhigh(j);
394 hTotErr->SetBinContent(i,GetTotalSystErr(pt,erryh));
396 hTotErr->SetBinContent(i,GetTotalSystErr(pt));
399 hTotErr->SetLineColor(1);
400 hTotErr->SetLineWidth(3);
401 hTotErr->Draw("same");
402 leg->AddEntry(hTotErr,"Total (excl. norm.)","l");
407 cSystErr->SaveAs("RelativeSystematics.eps");