1 /**************************************************************************
2 * Copyright(c) 1998-1999, 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 **************************************************************************/
17 // Drawing nice pictures containing
19 // - Signal/Background
21 // More Post Analysis code will be added in time
23 // Autor: Markus Fasel
32 #include "AliCFContainer.h"
33 #include "AliCFEffGrid.h"
35 #include "AliHFEcontainer.h"
36 #include "AliHFEcuts.h"
37 #include "AliHFEpostAnalysis.h"
39 ClassImp(AliHFEpostAnalysis)
41 //____________________________________________________________
42 AliHFEpostAnalysis::AliHFEpostAnalysis():
46 fEfficiencyContainer(NULL),
47 fPIDperformance(NULL),
48 fSignalToBackgroundMC(NULL)
51 // Default Constructor
55 //____________________________________________________________
56 AliHFEpostAnalysis::AliHFEpostAnalysis(const AliHFEpostAnalysis &ref):
58 fResults(ref.fResults),
59 fAnalysisObjects(ref.fAnalysisObjects),
60 fEfficiencyContainer(ref.fEfficiencyContainer),
61 fPIDperformance(ref.fPIDperformance),
62 fSignalToBackgroundMC(ref.fSignalToBackgroundMC)
69 //____________________________________________________________
70 AliHFEpostAnalysis& AliHFEpostAnalysis::operator=(const AliHFEpostAnalysis &ref){
72 // Assignment Operator
74 if(this == &ref) return *this;
75 TObject::operator=(ref);
76 fResults = ref.fResults;
77 fAnalysisObjects = ref.fAnalysisObjects;
78 fPIDperformance = ref.fPIDperformance;
79 fSignalToBackgroundMC = ref.fSignalToBackgroundMC;
84 //____________________________________________________________
85 AliHFEpostAnalysis::~AliHFEpostAnalysis(){
87 // Do not delete objects where we are not owner
89 if(fResults) delete fResults;
92 //____________________________________________________________
93 Int_t AliHFEpostAnalysis::SetTaskQA(const TList *input){
95 // Publish the results to the post analysis
98 fPIDperformance = dynamic_cast<THnSparseF *>(input->FindObject("PIDperformance"));
100 AliError("Histogram fPIDperformance not found in the List of Outputs");
102 SETBIT(fAnalysisObjects, kPIDperf);
105 fSignalToBackgroundMC = dynamic_cast<THnSparseF *>(input->FindObject("SignalToBackgroundMC"));
106 if(!fSignalToBackgroundMC){
107 AliError("Histogram fSignalToBackgroundMC not found in the list of outputs");
109 SETBIT(fAnalysisObjects, kSigBackg);
112 AliInfo(Form("Found %d analysis objects", nFound));
116 //____________________________________________________________
117 void AliHFEpostAnalysis::StoreOutput(const char *filename){
119 // Save the results produced in a rootfile
122 TFile *outfile = new TFile(filename, "RECREATE");
124 fResults->Write("HFEresults", TObject::kSingleKey);
130 //____________________________________________________________
131 void AliHFEpostAnalysis::DrawMCSignal2Background(){
133 // Draw the MC signal/background plots
135 if(!fSignalToBackgroundMC) return;
137 // First Select everything within the first ITS Layer
138 fSignalToBackgroundMC->GetAxis(5)->SetRange(2,2);
139 TH1 *hEff[3], *hSB[3];
140 // Select for different charge
141 hEff[0] = CreateHistoSignalToBackgroundMC(0, 0);
142 hEff[1] = CreateHistoSignalToBackgroundMC(0, 1);
143 hEff[2] = CreateHistoSignalToBackgroundMC(0, 2);
145 hSB[0] = CreateHistoSignalToBackgroundMC(1, 0);
146 hSB[1] = CreateHistoSignalToBackgroundMC(1, 1);
147 hSB[2] = CreateHistoSignalToBackgroundMC(1, 2);
150 fSignalToBackgroundMC->GetAxis(5)->SetRange(0, fSignalToBackgroundMC->GetAxis(4)->GetNbins());
153 TCanvas *cMCSB = new TCanvas("cMCSB", "MC Sig/Backg studies", 800, 400);
156 TH1 **sample[2] = {&hEff[0], &hSB[0]};
157 const char *chargename[3] = {"All Charges", "Negative Charge", "Positive Charge"};
158 for(Int_t isample = 0; isample < 2; isample++){
159 leg = new TLegend(0.7, 0.1, 0.89, 0.3);
160 leg->SetBorderSize(1);
161 leg->SetFillColor(kWhite);
162 cMCSB->cd(isample + 1);
163 for(Int_t icharge = 0; icharge < 3; icharge++){
164 sample[isample][icharge]->Draw(icharge > 0? "psame" : "p");
165 leg->AddEntry(sample[isample][icharge], chargename[icharge], "p");
172 //____________________________________________________________
173 TH1 *AliHFEpostAnalysis::CreateHistoSignalToBackgroundMC(Int_t mode, Int_t charge){
175 // Make Efficiency / SB histograms for different charges
177 TH1 *hNom = NULL, *hDenom = NULL;
178 if(charge) fSignalToBackgroundMC->GetAxis(3)->SetRange(charge, charge);
179 // For Signal electrons we project axis 4 to everything > 0
180 fSignalToBackgroundMC->GetAxis(4)->SetRange(2,3);
181 hNom = fSignalToBackgroundMC->Projection(0);
182 hNom->SetName("nom");
183 fSignalToBackgroundMC->GetAxis(4)->SetRange(0, fSignalToBackgroundMC->GetAxis(4)->GetLast() + 1);
184 if(mode == 1) fSignalToBackgroundMC->GetAxis(4)->SetRange(1,1);
185 hDenom = fSignalToBackgroundMC->Projection(0);
186 hDenom->SetName("denom");
187 if(mode == 1) fSignalToBackgroundMC->GetAxis(4)->SetRange(0, fSignalToBackgroundMC->GetAxis(4)->GetLast() + 1);
188 if(charge) fSignalToBackgroundMC->GetAxis(3)->SetRange(0, fSignalToBackgroundMC->GetAxis(3)->GetLast() + 1);
190 TH1 *hEff = dynamic_cast<TH1D *>(hNom->Clone());
192 TString hname, cname;
193 hname = mode ? "sigToBack" : "sigEff";
194 Color_t mycolor = kBlack;
196 case 0: mycolor = kBlue; cname = "All"; break;
197 case 1: mycolor = kBlack; cname = "Neg"; break;
198 case 2: mycolor = kRed; cname ="Pos"; break;
202 hEff->SetName(hname);
203 hEff->SetTitle(mode ? "Signal/Background" : "Signal/(Signal+Background)");
204 hEff->Divide(hDenom);
207 hEff->GetXaxis()->SetTitle("p_{T} / GeV/c");
208 hEff->GetYaxis()->SetTitle("Efficiency");
209 hEff->SetStats(kFALSE);
210 hEff->SetLineColor(kBlack);
211 hEff->SetLineWidth(1);
212 hEff->SetMarkerStyle(22);
213 hEff->SetMarkerColor(mycolor);
216 delete hNom; delete hDenom;
220 //____________________________________________________________
221 void AliHFEpostAnalysis::DrawEfficiency(){
223 // Draw the Efficiency
225 // + InAcceptance / Generated
226 // + Signal / Generated
227 // + Selected / Generated
228 // + Selected / InAcceptance (Reconstructible)
230 TCanvas *cEff = new TCanvas("cEff", "Efficiency", 800, 600);
232 if(!fEfficiencyContainer) return;
233 AliCFContainer *tracks = fEfficiencyContainer->MakeMergedCFContainer("trackContCombined", "MC + Rec(reco) Track Information", "MCTrackCont:recTrackContReco");
234 AliCFEffGrid *effCalc = new AliCFEffGrid("effCalc", "Efficiency Calculation Grid", *tracks);
235 effCalc->CalculateEfficiency(AliHFEcuts::kStepMCInAcceptance, AliHFEcuts::kStepMCGenerated);
236 TH1 *effReconstructibleP = effCalc->Project(0);
237 effReconstructibleP->SetName("effReconstructibleP");
238 effReconstructibleP->SetTitle("Efficiency of reconstructible tracks");
239 effReconstructibleP->GetXaxis()->SetTitle("p_{T} / GeV/c");
240 effReconstructibleP->GetYaxis()->SetTitle("Efficiency");
241 effReconstructibleP->GetYaxis()->SetRangeUser(0.,1.);
242 effReconstructibleP->SetMarkerStyle(22);
243 effReconstructibleP->SetMarkerColor(kBlue);
244 effReconstructibleP->SetLineColor(kBlack);
245 effReconstructibleP->SetStats(kFALSE);
247 effReconstructibleP->Draw("e");
248 effCalc->CalculateEfficiency(AliHFEcuts::kStepMCGeneratedZOutNoPileUpCentralityFine, AliHFEcuts::kStepMCGenerated);
249 TH1 *effSignal = effCalc->Project(0);
250 effSignal->SetName("effSignal");
251 effSignal->SetTitle("Efficiency of Signal Electrons");
252 effSignal->GetXaxis()->SetTitle("p_{T} / GeV/c");
253 effSignal->GetYaxis()->SetTitle("Efficiency");
254 effSignal->GetYaxis()->SetRangeUser(0., 1.);
255 effSignal->SetMarkerStyle(22);
256 effSignal->SetMarkerColor(kBlue);
257 effSignal->SetLineColor(kBlack);
258 effSignal->SetStats(kFALSE);
260 effSignal->Draw("e");
261 effCalc->CalculateEfficiency(tracks->GetNStep() - 1, AliHFEcuts::kStepMCGenerated);
262 TH1 *effPIDP = effCalc->Project(0);
263 effPIDP->SetName("effPIDP");
264 effPIDP->SetTitle("Efficiency of selected tracks");
265 effPIDP->GetXaxis()->SetTitle("p_{T} / GeV/c");
266 effPIDP->GetYaxis()->SetTitle("Efficiency");
267 effPIDP->GetYaxis()->SetRangeUser(0.,1.);
268 effPIDP->SetMarkerStyle(22);
269 effPIDP->SetMarkerColor(kBlue);
270 effPIDP->SetLineColor(kBlack);
271 effPIDP->SetStats(kFALSE);
274 effCalc->CalculateEfficiency(tracks->GetNStep() - 1, AliHFEcuts::kStepMCInAcceptance);
275 TH1 *effPIDAcc = effCalc->Project(0);
276 effPIDAcc->SetName("effPIDAcc");
277 effPIDAcc->SetTitle("Efficiency of selected tracks in acceptance");
278 effPIDAcc->GetXaxis()->SetTitle("p_{T} / GeV/c");
279 effPIDAcc->GetYaxis()->SetTitle("Efficiency");
280 effPIDAcc->GetYaxis()->SetRangeUser(0.,1.);
281 effPIDAcc->SetMarkerStyle(22);
282 effPIDAcc->SetMarkerColor(kBlue);
283 effPIDAcc->SetLineColor(kBlack);
284 effPIDAcc->SetStats(kFALSE);
286 effPIDAcc->Draw("e");
290 //____________________________________________________________
291 void AliHFEpostAnalysis::DrawPIDperformance(){
293 // Plotting Ratio histograms
294 // + All electrons / all candidates (Purity for Electrons)
295 // + All signal electrons / all electrons (Purity for signals)
298 if(!fPIDperformance) return;
300 TH1 *electronPurity[3], *signalPurity[3], *fakeContamination[3];
301 electronPurity[0] = CreateHistoPIDperformance(0, 0);
302 electronPurity[1] = CreateHistoPIDperformance(0, 1);
303 electronPurity[2] = CreateHistoPIDperformance(0, 2);
305 signalPurity[0] = CreateHistoPIDperformance(1, 0);
306 signalPurity[1] = CreateHistoPIDperformance(1, 1);
307 signalPurity[2] = CreateHistoPIDperformance(1, 2);
309 fakeContamination[0] = CreateHistoPIDperformance(2, 0);
310 fakeContamination[1] = CreateHistoPIDperformance(2, 1);
311 fakeContamination[2] = CreateHistoPIDperformance(2, 2);
314 TCanvas *cRatios = new TCanvas("cRatios", "Ratio Plots", 800, 600);
315 const char *chargename[3] = {"All Charges", "Negative Charge", "Positive Charge"};
316 cRatios->Divide(2,2);
317 TH1 **sample[3] = {&electronPurity[0], &signalPurity[0], &fakeContamination[0]};
319 for(Int_t isample = 0; isample < 3; isample++){
320 cRatios->cd(isample + 1);
321 leg = new TLegend(0.7, 0.1, 0.89, 0.3);
322 leg->SetBorderSize(1);
323 leg->SetFillColor(kWhite);
324 for(Int_t icharge = 0; icharge < 3; icharge++){
325 leg->AddEntry(sample[isample][icharge], chargename[icharge], "p");
326 sample[isample][icharge]->Draw(icharge > 0 ? "esame" : "e");
333 //____________________________________________________________
334 TH1 *AliHFEpostAnalysis::CreateHistoPIDperformance(Int_t mode, Int_t charge){
336 // Make Histograms for PID performance plots
338 fPIDperformance->GetAxis(4)->SetRange(0, fPIDperformance->GetAxis(4)->GetNbins()+1);
339 fPIDperformance->GetAxis(3)->SetRange(0, fPIDperformance->GetAxis(3)->GetNbins() + 1);
341 TH1 *hNom = NULL, *hDenom = NULL;
342 TString hname, htitle, cname;
343 Color_t mycolor = kBlack;
344 if(charge) fPIDperformance->GetAxis(3)->SetRange(charge, charge);
345 // Normalisation by all candidates - no restriction in axis 4 - only for mode == 1
346 if(mode == 1) fPIDperformance->GetAxis(4)->SetRange(2,3);
347 hDenom = fPIDperformance->Projection(0);
349 hDenom->SetName("hDenom");
350 if(mode == 1) fPIDperformance->GetAxis(4)->SetRange(0, fPIDperformance->GetAxis(4)->GetNbins() + 1);
351 // Nominator need a restriction in the 4th axis
353 case 0: // Electron purity
354 fPIDperformance->GetAxis(4)->SetRange(2,3);
355 hname = "electronPurity";
356 htitle = "Electron Purity";
358 case 1: // Signal purity
359 fPIDperformance->GetAxis(4)->SetRange(3,3); // here signal not divided into charm and beauty
360 hname = "signalPurity";
361 htitle = "Signal Purity";
363 case 2: // Fake contamination
364 fPIDperformance->GetAxis(4)->SetRange(1,1);
365 hname = "fakeContamination";
366 htitle = "Contamination of misidentified hadrons";
385 hNom = fPIDperformance->Projection(0);
387 hNom->SetName("hNom");
389 fPIDperformance->GetAxis(4)->SetRange(0, fPIDperformance->GetAxis(4)->GetNbins()+1);
390 if(charge) fPIDperformance->GetAxis(3)->SetRange(0, fPIDperformance->GetAxis(3)->GetNbins() + 1);
392 // Create Efficiency histogram
393 TH1 *hEff = dynamic_cast<TH1D *>(hNom->Clone());
395 hEff->SetName(hname.Data());
396 hEff->SetTitle(htitle.Data());
397 hEff->Divide(hDenom);
399 hEff->GetXaxis()->SetTitle("p_{T} / GeV/c");
400 hEff->GetYaxis()->SetTitle(mode < 2 ? "Purity / %" : "Contamination / %");
401 hEff->GetYaxis()->SetRangeUser(0., 100.);
402 hEff->SetStats(kFALSE);
403 hEff->SetLineColor(kBlack);
404 hEff->SetLineWidth(1);
405 hEff->SetMarkerColor(mycolor);
406 hEff->SetMarkerStyle(22);
407 delete hNom; delete hDenom;
412 //____________________________________________________________
413 void AliHFEpostAnalysis::DrawCutEfficiency(Bool_t MC, Int_t source){
415 // Calculate efficiency for each cut step
416 // Starting from MC steps
418 TCanvas *output = new TCanvas("effCut", "Cut Step efficiency", 800, 600);
420 TLegend *leg = new TLegend(0.6, 0.7, 0.89, 0.89);
421 leg->SetHeader("Cut Step:");
422 leg->SetBorderSize(0);
423 leg->SetFillColor(kWhite);
424 leg->SetFillStyle(0);
426 AliCFContainer *tracks = fEfficiencyContainer->MakeMergedCFContainer("mergedTracks", "Container for MC and reconstructed Track information", "MCTrackCont:recTrackContReco");
427 Int_t nStepMC = fEfficiencyContainer->GetCFContainer("MCTrackCont")->GetNStep();
428 AliDebug(1, Form("Number of MC Cut Steps: %d", nStepMC));
429 const Int_t markerStart = 24;
430 const Int_t colorStart = 1;
434 if(source > -1 && source < 4){
435 AliInfo(Form("Setting source to %d", source));
436 for(Int_t istep = 0; istep < tracks->GetNStep(); istep++) tracks->GetAxis(4, istep)->SetRange(source + 1, source + 1);
439 AliCFEffGrid effcalc("cutEfficiency", "Cut step efficiency calculation", *tracks);
440 Bool_t first = kTRUE;
441 // Calculate efficiency for MC Steps
442 Int_t histcounter = 0;
445 // Draw plots related to MC values
447 effcalc.CalculateEfficiency(AliHFEcuts::kStepMCInAcceptance, AliHFEcuts::kStepMCGenerated);
448 hTemp = effcalc.Project(0);
449 hTemp->SetName(Form("hEff%d", histcounter));
450 hTemp->SetTitle("Cut Step Efficiency");
451 hTemp->SetMarkerColor(colorStart + histcounter);
452 hTemp->SetMarkerStyle(markerStart + histcounter);
453 hTemp->GetXaxis()->SetTitle("p / GeV/c");
454 hTemp->GetYaxis()->SetTitle("Efficiency");
455 hTemp->GetYaxis()->SetRangeUser(0., 2.);
456 hTemp->SetStats(kFALSE);
458 leg->AddEntry(hTemp, tracks->GetStepTitle(AliHFEcuts::kStepMCInAcceptance), "p");
460 effcalc.CalculateEfficiency(nStepMC, AliHFEcuts::kStepMCGenerated);
461 hTemp = effcalc.Project(0);
462 hTemp->SetName("hEff2");
463 hTemp->SetTitle("Cut Step Efficiency");
464 hTemp->SetMarkerColor(colorStart + histcounter);
465 hTemp->SetMarkerStyle(markerStart + histcounter);
466 hTemp->GetXaxis()->SetTitle("p / GeV/c");
467 hTemp->GetYaxis()->SetTitle("Efficiency");
468 hTemp->GetYaxis()->SetRangeUser(0., 2.);
469 hTemp->SetStats(kFALSE);
470 hTemp->Draw("epsame");
471 leg->AddEntry(hTemp, tracks->GetStepTitle(nStepMC), "p");
475 for(Int_t istep = nStepMC+1; istep < tracks->GetNStep(); istep++){
476 effcalc.CalculateEfficiency(istep, istep - 1);
477 hTemp = effcalc.Project(0);
478 hTemp->SetName(Form("hEff%d", istep));
479 hTemp->SetTitle("Cut Step Efficiency");
480 hTemp->SetMarkerColor(colorStart + histcounter);
481 hTemp->SetMarkerStyle(markerStart + histcounter);
482 hTemp->SetStats(kFALSE);
483 hTemp->Draw(first ? "ep" : "epsame");
484 hTemp->GetXaxis()->SetTitle("P / GeV/c");
485 hTemp->GetYaxis()->SetTitle("Efficiency");
486 hTemp->GetYaxis()->SetRangeUser(0., 2.);
487 leg->AddEntry(hTemp, tracks->GetStepTitle(istep), "p");