AliPhysics  b43479f (b43479f)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DrawEfficiency.C
Go to the documentation of this file.
1 
2 // macro to calculate the pt RecAnCuts/MCAcc
3 // efficiency, from a CF container
4 // Author: C. Zampolli
5 
6 // channel could be:
7 // D0: D0 -> Kpi, from old task
8 // D0_New: D0 -> Kpi, from new CF common implementation
9 // Dplus_New: D+ -> Kpipi, from new CF common implementation
10 //
11 // eff = sum of the efficiency indexes to compute
12 // MCAcc_over_MCLimAcc = 0x001
13 // RecPPR_over_MCAcc = 0x002
14 // RecPID_over_MCAcc = 0x004
15 // all = 0x007
16 //
17 // selection = D origin
18 // 0 --> from c only
19 // 1 --> from b only
20 // 2 --> from both c and b
21 
22 #include <Riostream.h>
23 
24 extern TRandom *gRandom;
25 extern TBenchmark *gBenchmark;
26 extern TSystem *gSystem;
27 
28 void DrawEfficiency(const char* channel, Int_t selection = 0, Int_t ieff = 7){
29 
30  gROOT->SetStyle("Plain");
31  gStyle->SetPalette(1);
32  gStyle->SetOptStat(0);
33  gStyle->SetPalette(1);
34  gStyle->SetCanvasColor(0);
35  gStyle->SetFrameFillColor(0);
36  gStyle->SetOptTitle(0);
37 
38  gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include");
39  gSystem->Load("libANALYSIS");
40  gSystem->Load("libANALYSISalice");
41  gSystem->Load("$ALICE_ROOT/CORRFW/libCORRFW") ;
42  gSystem->Load("libPWGHFbase");
43  gSystem->Load("libPWGHFvertexingHF");
44 
45  Int_t mcAcc_over_mcLimAcc = 0x001;
46  Int_t recPPR_over_mcAcc = 0x002;
47  Int_t recPID_over_mcAcc = 0x004;
48 
49  // pt index
50  Int_t ipt =0;
51 
52  Int_t stepNum = -1;
53  Int_t stepDen = -1;
54 
55  // Read the container from file
56  TFile* f = new TFile("AnalysisResults.root");
57  TString directoryName;
58  TString containerName;
59  TString cutName;
60  TString outfileName;
61 
62  if (channel == "D0") {
63  if (selection == 0){
64  directoryName = "PWG3_D2H_CFtaskD0toKpi";
65  containerName = "CFHFccontainer0";
66  cutName = "Cuts";
67  outfileName = "fileEff_D0_from_c.root";
68  }
69  else if (selection == 1){
70  directoryName = "PWG3_D2H_CFtaskD0toKpiKeepD0fromBOnly";
71  containerName = "CFHFccontainer0D0fromB";
72  cutName = "Cuts";
73  outfileName = "fileEff_D0_from_b.root";
74  }
75  else if (selection == 2){
76  directoryName = "PWG3_D2H_CFtaskD0toKpiKeepD0fromB";
77  containerName = "CFHFccontainer0allD0";
78  cutName = "Cuts";
79  outfileName = "fileEff_D0_from_c_and_b.root";
80  }
81  else {
82  Printf("not a valid selection, return");
83  return;
84  }
85  }
86  else if (channel == "D0_New"){
87  if (selection == 0){
88  directoryName = "PWG3_D2H_CFtaskD0toKpi_NEW";
89  containerName = "CFHFccontainer0_New";
90  cutName = "Cuts_New";
91  //directoryName = "PWG3_D2H_CFtaskD0toKpi_CommonFramework";
92  //containerName = "CFHFccontainer0_CommonFramework";
93  //cutName = "Cuts_CommonFramework";
94  outfileName = "fileEff_D0_CommonFramework_from_c.root";
95  }
96  else if (selection == 1){
97  directoryName = "PWG3_D2H_CFtaskD0toKpiKeepDfromBOnly";
98  containerName = "CFHFccontainer0DfromB_New";
99  cutName = "Cuts_New";
100  //directoryName = "PWG3_D2H_CFtaskD0toKpiKeepDfromBOnly_CommonFramework";
101  //containerName = "CFHFccontainer0DfromB_CommonFramework";
102  //cutName = "Cuts_CommonFramework";
103  outfileName = "fileEff_D0_CommonFramework_from_b.root";
104  }
105  else if (selection == 2){
106  directoryName = "PWG3_D2H_CFtaskD0toKpiKeepDfromB_NEW";
107  containerName = "CFHFccontainer0allD_New";
108  cutName = "Cuts_New";
109  //directoryName = "PWG3_D2H_CFtaskD0toKpiKeepDfromB_CommonFramework";
110  //containerName = "CFHFccontainer0allD_CommonFramework";
111  //cutName = "Cuts_CommonFramework";
112  outfileName = "fileEff_D0_CommonFramework_from_c_and_b.root";
113  }
114  else {
115  Printf("not a valid selection, return");
116  return;
117  }
118  }
119  else if (channel == "Dplus_New"){
120  if (selection == 0){
121  directoryName = "PWG3_D2H_CFtaskDplustoKpipi_NEW";
122  containerName = "CFHFccontainer0_New_3Prong";
123  cutName = "Cuts_3Prong";
124  //directoryName = "PWG3_D2H_CFtaskDplustoKpipi_CommonFramework";
125  //containerName = "CFHFccontainer0_3Prong_CommonFramework";
126  //cutName = "Cuts_3Prong_CommonFramework";
127  outfileName = "fileEff_Dplus_CommonFramework_from_c.root";
128  }
129  else if (selection == 1){
130  directoryName = "PWG3_D2H_CFtaskDplustoKpipiKeepDfromBOnly";
131  containerName = "CFHFccontainer0DfromB_New_3Prong";
132  cutName = "Cuts_3Prong";
133  //directoryName = "PWG3_D2H_CFtaskDplustoKpipiKeepDfromBOnly_CommonFramework";
134  //containerName = "CFHFccontainer0DfromB_3Prong_CommonFramework";
135  //cutName = "Cuts_3Prong_CommonFramework";
136  outfileName = "fileEff_Dplus_CommonFramework_from_b.root";
137  }
138  else if (selection == 2){
139  directoryName = "PWG3_D2H_CFtaskDplustoKpipiKeepDfromB_NEW";
140  containerName = "CFHFccontainer0allD_New_3Prong";
141  cutName = "Cuts_3Prong";
142  //directoryName = "PWG3_D2H_CFtaskDplustoKpipiKeepDfromB_CommonFramework";
143  //containerName = "CFHFccontainer0allD_3Prong_CommonFramework";
144  //cutName = "Cuts_3Prong_CommonFramework";
145  outfileName = "fileEff_Dplus_CommonFramework_from_c_and_b.root";
146  }
147  else {
148  Printf("not a valid selection, return");
149  return;
150  }
151  }
152  else {
153  Printf("not a valid channel, return");
154  return;
155  }
156 
157  Printf("Opening file Analysisresults.root");
158  Printf("Reading Directory \"%s\"",directoryName.Data());
159  Printf("Getting CF Container \"%s\"",containerName.Data());
160  Printf("Getting Cut Object \"%s\"",cutName.Data());
161 
162 
163  TDirectoryFile* d = (TDirectoryFile*)f->Get(directoryName.Data());
164  if (!d){
165  Printf("Directory does not exist! Check directory name (%s) in file AnalysisResults.root - returning...", directoryName.Data());
166  return;
167  }
168  AliCFContainer *data = (AliCFContainer*) (d->Get(containerName.Data()));
169  AliRDHFCuts *cutsRDHF = (AliRDHFCuts*)(d->Get(cutName.Data()));
170 
171  if (!data){
172  Printf("Container does not exist! Check container name (%s) in directory %s - returning...", containerName.Data(), directoryName.Data());
173  return;
174  }
175 
176  TFile* fileEff = new TFile(outfileName.Data(), "RECREATE");
177  TString plotDir(Form("EffPlots/%s",channel));
178  gSystem->Exec(Form("mkdir EffPlots"));
179  gSystem->Exec(Form("mkdir %s",plotDir.Data()));
180 
181  //construct the efficiency grid from the data container
182  AliCFEffGrid *eff = new AliCFEffGrid("eff"," The efficiency",*data);
183 
184  TCanvas *ceffpt = new TCanvas("ceffpt","Efficiency vs pt",50,50,550,550);
185  ceffpt->cd();
186  ceffpt->SetLeftMargin(0.15);
187  ceffpt->SetRightMargin(0.05);
188  TH1D *hpteffCF = 0x0; //the efficiency vs pt
189 
190  if (ieff & mcAcc_over_mcLimAcc){
191  AliCFEffGrid *eff = new AliCFEffGrid("eff"," The efficiency",*data);
192  stepDen = (Int_t)(AliCFTaskVertexingHF::kStepGeneratedLimAcc);
193  stepNum = (Int_t)(AliCFTaskVertexingHF::kStepAcceptance);
194  printf("Calculating efficiency for mcAcc_over_mcLimAcc: stepDen = %d, stepNum = %d\n",stepDen,stepNum);
195  eff->CalculateEfficiency(stepNum,stepDen); //eff= step1/step0
196 
197  //canvas
198  ceffpt->cd();
199 
200  //The efficiency along the variables
201  hpteffCF = eff->Project(ipt);
202  SetHistoEff(hpteffCF,8,20,"mcAcc_over_mcLimAcc");
203  hpteffCF->Draw("hist");
204  hpteffCF->Draw("err same");
205  fileEff->cd();
206  hpteffCF->Write("hpteffCF_mcAcc_over_mcLimAcc");
207 
208  // printing png files
209  ceffpt->Print(Form("%s/effpt_mcAcc_over_mcLimAcc.png", plotDir.Data()));
210  ceffpt->Print(Form("%s/effpt_mcAcc_over_mcLimAcc.eps", plotDir.Data()));
211  ceffpt->Print(Form("%s/effpt_mcAcc_over_mcLimAcc.gif", plotDir.Data()));
212  delete eff;
213  eff = 0x0;
214  }
215 
216  if (ieff & recPPR_over_mcAcc){
217  AliCFEffGrid *eff = new AliCFEffGrid("eff"," The efficiency",*data);
218  stepDen = (Int_t)(AliCFTaskVertexingHF::kStepAcceptance);
219  stepNum = (Int_t)(AliCFTaskVertexingHF::kStepRecoPPR);
220  printf("Calculating efficiency for RecPPR_over_mcAcc: stepDen = %d, stepNum = %d\n",stepDen,stepNum);
221  eff->CalculateEfficiency(stepNum,stepDen); //eff= step1/step0
222 
223  //canvas
224  ceffpt->cd();
225 
226  //The efficiency along the variables
227  hpteffCF = eff->Project(ipt);
228  SetHistoEff(hpteffCF,8,20, "recAnCuts_over_mcAcc");
229  hpteffCF->Draw("hist");
230  hpteffCF->Draw("err same");
231  fileEff->cd();
232  hpteffCF->Write("hpteffCF_RecAnCut_over_mcAcc");
233 
234  // printing png files
235  ceffpt->Print(Form("%s/effpt_RecAnCut_over_mcAcc.png", plotDir.Data()));
236  ceffpt->Print(Form("%s/effpt_RecAnCut_over_mcAcc.eps", plotDir.Data()));
237  ceffpt->Print(Form("%s/effpt_RecAnCut_over_mcAcc.gif", plotDir.Data()));
238  delete eff;
239  eff = 0x0;
240  }
241 
242  if (ieff & recPID_over_mcAcc){
243  AliCFEffGrid *eff = new AliCFEffGrid("eff"," The efficiency",*data);
244  stepDen = (Int_t)(AliCFTaskVertexingHF::kStepAcceptance);
245  stepNum = (Int_t)(AliCFTaskVertexingHF::kStepRecoPID);
246  printf("Calculating efficiency for RecPID_over_mcAcc: stepDen = %d, stepNum = %d\n",stepDen,stepNum);
247  eff->CalculateEfficiency(stepNum,stepDen); //eff= step1/step0
248 
249  //canvas
250  ceffpt->cd();
251 
252  //The efficiency along the variables
253  hpteffCF = eff->Project(ipt);
254  SetHistoEff(hpteffCF,8,20,"recPID_over_mcAcc");
255  hpteffCF->Draw("hist");
256  hpteffCF->Draw("err same");
257  fileEff->cd();
258  hpteffCF->Write("hpteffCF_RecPID_over_mcAcc");
259 
260  // printing png files
261  ceffpt->Print(Form("%s/effpt_RecPID_over_mcAcc.png", plotDir.Data()));
262  ceffpt->Print(Form("%s/effpt_RecPID_over_mcAcc.eps", plotDir.Data()));
263  ceffpt->Print(Form("%s/effpt_RecPID_over_mcAcc.gif", plotDir.Data()));
264  delete eff;
265  eff = 0x0;
266  }
267 
268  cutsRDHF->Write("Cuts");
269 
270  // writing single distributions
271  TH1D *hMCAccpt = data->ShowProjection(ipt, AliCFTaskVertexingHF::kStepAcceptance);
272  SetHistoDistribution(hMCAccpt, 1, 20);
273  hMCAccpt->Draw();
274  TH1D *hMCLimAccpt = data->ShowProjection(ipt, AliCFHeavyFlavourTaskMultiVarMultiStep::kStepGeneratedLimAcc);
275  SetHistoDistribution(hMCLimAccpt, 4, 20);
276  TH1D *hRecoAnCutspt = data->ShowProjection(ipt, AliCFTaskVertexingHF::kStepRecoPPR);
277  SetHistoDistribution(hRecoAnCutspt, 8, 20);
278  TH1D *hRecoPIDpt = data->ShowProjection(ipt, AliCFTaskVertexingHF::kStepRecoPID);
279  SetHistoDistribution(hRecoPIDpt, 6, 20);
280  hMCAccpt->Write("hMCAccpt");
281  hMCLimAccpt->Write("hMCLimAccpt");
282  hRecoAnCutspt->Write("hRecoAnCutspt");
283  hRecoPIDpt->Write("hRecoPIDpt");
284 
285  // fileEff->Close(); // commented out to see the canvas on the screen....
286 
287 }
288 
289 void SetHistoEff(TH1D* h, Int_t color, Int_t style, const char* effType){
290 
291  h->SetLineColor(color);
292  h->SetLineWidth(3);
293  h->SetMarkerStyle(style);
294  h->SetMarkerColor(color);
295  h->SetMarkerSize(1.2);
296  h->GetYaxis()->SetTitleOffset(1.5);
297  h->GetXaxis()->SetTitleOffset(1.2);
298  h->GetXaxis()->SetTitle("p_{t} [GeV/c]");
299  h->GetYaxis()->SetTitle(Form("%s, Efficiency",effType));
300  return;
301 }
302 void SetHistoDistribution(TH1D* h, Int_t color, Int_t style){
303 
304  h->SetLineColor(color);
305  h->SetLineWidth(3);
306  h->SetMarkerStyle(style);
307  h->SetMarkerColor(color);
308  h->SetMarkerSize(1.2);
309  h->GetXaxis()->SetTitleOffset(1.2);
310  h->GetXaxis()->SetTitle("p_{t} [GeV/c]");
311  return;
312 }
color
Int_t color[]
Definition: DrawAnaCaloTrackQA.C:21
SetHistoDistribution
void SetHistoDistribution(TH1D *h, Int_t color, Int_t style)
Definition: DrawEfficiency.C:302
gSystem
TSystem * gSystem
gRandom
TRandom * gRandom
SetHistoEff
void SetHistoEff(TH1D *h, Int_t color, Int_t style, const char *effType)
Definition: DrawEfficiency.C:289
DrawEfficiency
void DrawEfficiency(const char *channel, Int_t selection=0, Int_t ieff=7)
Definition: DrawEfficiency.C:28
gBenchmark
TBenchmark * gBenchmark