]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/JetTasks/AliAnalysisTaskPWG4PidDetEx.cxx
Init with values, to satisfy CINT compilation
[u/mrichter/AliRoot.git] / PWG4 / JetTasks / AliAnalysisTaskPWG4PidDetEx.cxx
CommitLineData
de6f8090 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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#include <TROOT.h>
17#include <TSystem.h>
18#include <TChain.h>
19#include <TTree.h>
20#include <TFile.h>
21#include <TList.h>
22#include <TStyle.h>
23#include <TCanvas.h>
24#include <TMath.h>
25#include <TH1.h>
26#include <TH2.h>
27#include <TF1.h>
28#include <TProfile.h>
29#include <TLegend.h>
30#include <TDatabasePDG.h>
31#include <TParticle.h>
32
33#include "AliAnalysisManager.h"
34#include "AliAnalysisFilter.h"
35#include "AliESDInputHandler.h"
36#include "AliESDEvent.h"
37#include "AliESDVertex.h"
38#include "AliMCEventHandler.h"
39#include "AliMCEvent.h"
40#include "AliStack.h"
41#include "AliAODInputHandler.h"
42#include "AliAODEvent.h"
43#include "AliAODVertex.h"
44#include "AliAODMCParticle.h"
45#include "AliLog.h"
46
47#include "AliAnalysisTaskPWG4PidDetEx.h"
48
49// STL includes
50#include <iostream>
51using namespace std;
52
53//
54// Analysis class example for PID using detector signals.
55// Works on ESD and AOD; has a flag for MC analysis
56//
57// Alexandru.Dobrin@hep.lu.se
58// Peter.Christiansen@hep.lu.se
59//
60
61ClassImp(AliAnalysisTaskPWG4PidDetEx)
62
63const Double_t AliAnalysisTaskPWG4PidDetEx::fgkCtau = 370; // distance for kaon decay
792df7ad 64const Double_t AliAnalysisTaskPWG4PidDetEx::fgkPionMass = 1.39570000000000000e-01;
65const Double_t AliAnalysisTaskPWG4PidDetEx::fgkKaonMass = 4.93599999999999983e-01;
66const Double_t AliAnalysisTaskPWG4PidDetEx::fgkProtonMass = 9.38270000000000048e-01;
67const Double_t AliAnalysisTaskPWG4PidDetEx::fgkC = 2.99792458e-2;
68
69
70
71
de6f8090 72
73//_____________________________________________________________________________
74AliAnalysisTaskPWG4PidDetEx::AliAnalysisTaskPWG4PidDetEx():
75 AliAnalysisTaskSE(),
76 fESD(0),
77 fAOD(0),
78 fListOfHists(0),
79 fEtaCut(0.9),
80 fPtCut(0.5),
81 fXbins(20),
82 fXmin(0.),
83 fTOFCutP(2.),
84 fTrackFilter(0x0),
85 fAnalysisMC(kTRUE),
86 fAnalysisType("ESD"),
87 fTriggerMode(kMB2),
88 fEvents(0), fEffTot(0), fEffPID(0), fAccP(0), fAccPt(0), fKaonDecayCorr(0),
89 fdNdPt(0), fMassAll(0),
90 fdNdPtPion(0), fMassPion(0), fdEdxTPCPion(0), fbgTPCPion(0),
91 fdNdPtKaon(0), fMassKaon(0), fdEdxTPCKaon(0), fbgTPCKaon(0),
92 fdNdPtProton(0), fMassProton(0), fdEdxTPCProton(0), fbgTPCProton(0),
93 fdNdPtMC(0), fdNdPtMCPion(0), fdNdPtMCKaon(0), fdNdPtMCProton(0)
94{
95 // Default constructor (should not be used)
96}
97
98//______________________________________________________________________________
99AliAnalysisTaskPWG4PidDetEx::AliAnalysisTaskPWG4PidDetEx(const char *name):
100 AliAnalysisTaskSE(name),
101 fESD(0),
102 fAOD(0),
103 fListOfHists(0),
104 fEtaCut(0.9),
105 fPtCut(0.5),
106 fXbins(20),
107 fXmin(0.),
108 fTOFCutP(2.),
109 fTrackFilter(0x0),
110 fAnalysisMC(kTRUE),
111 fAnalysisType("ESD"),
112 fTriggerMode(kMB2),
113 fEvents(0), fEffTot(0), fEffPID(0), fAccP(0), fAccPt(0), fKaonDecayCorr(0),
114 fdNdPt(0), fMassAll(0),
115 fdNdPtPion(0), fMassPion(0), fdEdxTPCPion(0), fbgTPCPion(0),
116 fdNdPtKaon(0), fMassKaon(0), fdEdxTPCKaon(0), fbgTPCKaon(0),
117 fdNdPtProton(0), fMassProton(0), fdEdxTPCProton(0), fbgTPCProton(0),
118 fdNdPtMC(0), fdNdPtMCPion(0), fdNdPtMCKaon(0), fdNdPtMCProton(0)
119{
de6f8090 120 // Output slot #0 writes into a TList
18179af0 121 DefineOutput(1, TList::Class());
de6f8090 122
123}
124
125//_____________________________________________________________________________
126AliAnalysisTaskPWG4PidDetEx::~AliAnalysisTaskPWG4PidDetEx()
127{
128 // Destructor
129 // histograms are in the output list and deleted when the output
130 // list is deleted by the TSelector dtor
131 if (fListOfHists) {
132 delete fListOfHists;
133 fListOfHists = 0;
134 }
135}
136
18179af0 137// //______________________________________________________________________________
138// void AliAnalysisTaskPWG4PidDetEx::ConnectInputData(Option_t *)
139// {
140// // Connect AOD here
141// // Called once
142// if (fDebug > 1) AliInfo("ConnectInputData() \n");
de6f8090 143
18179af0 144// TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
145// if (!tree) {
146// Printf("ERROR: Could not read chain from input slot 0");
147// }
148// else {
149// if(fAnalysisType == "ESD") {
150// AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
151// if (!esdH) {
152// Printf("ERROR: Could not get ESDInputHandler");
153// } else
154// fESD = esdH->GetEvent();
155// }
156// else if(fAnalysisType == "AOD") {
157// AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
158// if (!aodH) {
159// Printf("ERROR: Could not get AODInputHandler");
160// } else
161// fAOD = aodH->GetEvent();
162// }
163// }
164
165// }
de6f8090 166
167//______________________________________________________________________________
18179af0 168void AliAnalysisTaskPWG4PidDetEx::UserCreateOutputObjects()
de6f8090 169{
18179af0 170 OpenFile(1);
de6f8090 171 fListOfHists = new TList();
172
173 fEvents = new TH1I("fEvents","Number of analyzed events; Events; Counts", 1, 0, 1);
174 fListOfHists->Add(fEvents);
175
176 fEffTot = new TH1F ("fEffTot","Efficiency; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
177 fEffTot->Sumw2();
178 fListOfHists->Add(fEffTot);
179
180 fEffPID = new TH1F ("fEffPID","Efficiency; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
181 fEffPID->Sumw2();
182 fListOfHists->Add(fEffPID);
183
184 fAccP = new TH1F ("fAccP","Acceptance; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
185 fAccP->Sumw2();
186 fListOfHists->Add(fAccP);
187
188 fAccPt = new TH1F ("fAccPt","Acceptance; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
189 fAccPt->Sumw2();
190 fListOfHists->Add(fAccPt);
191
192 fKaonDecayCorr = new TProfile("fKaonDecayCorr","Kaon decay correction vs p_{T}; p_{T} [GeV/c]; Kaon decay correction", fXbins, fXmin, fTOFCutP);
193 fListOfHists->Add(fKaonDecayCorr);
194
195 fdNdPt = new TH1F("fdNdPt","p_{T} distribution; p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", 100 , 0, 10);
196 fdNdPt->Sumw2();
197 fListOfHists->Add(fdNdPt);
198
199 fMassAll = new TH2F("fMassAll","Mass^{2} vs momentum (ToF); p [GeV/c]; M^{2} [GeV^{2}/c^{4}]", 100, -5, 5, 100, -0.2, 1.2);
200 fListOfHists->Add(fMassAll);
201
202 //Pion
203 fdNdPtPion = new TH1F("fdNdPtPion","p_{T} distribution (Pions); p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", fXbins, fXmin, fTOFCutP);
204 fdNdPtPion->Sumw2();
205 fListOfHists->Add(fdNdPtPion);
206
207 fMassPion = new TH2F("fMassPion","Mass squared for pions after cuts vs pmean; p [GeV/c]; Mass^{2} [GeV^{2}/c^{4}]", 100, -5, 5, 100, -0.2, 1.2);
208 fListOfHists->Add(fMassPion);
209
210 fdEdxTPCPion = new TH2F("fdEdxTPCPion","Energy loss vs momentum; p [GeV/c]; dE/dx [a. u.]", 40, -2, 2, 100, 0, 200);
211 fListOfHists->Add(fdEdxTPCPion);
212
213 fbgTPCPion = new TH2F("fbgTPCPion", "Energy loss vs #beta#gamma (Pions);#beta#gamma; dE/dx [a. u.]", 100, 0, 15, 100, 0, 200);
214 fListOfHists->Add(fbgTPCPion);
215
216 //Kaon
217 fdNdPtKaon = new TH1F("fdNdPtKaon","p_{T} distribution (Kaons);p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", fXbins, fXmin, fTOFCutP);
218 fdNdPtKaon->Sumw2();
219 fListOfHists->Add(fdNdPtKaon);
220
221 fMassKaon = new TH2F("fMassKaon","Mass squared for kaons after cuts vs pmean; p [GeV/c]; Mass^{2} [GeV^{2}/c^{4}]", 100, -5, 5, 100, -0.2, 1.2);
222 fListOfHists->Add(fMassKaon);
223
224 fdEdxTPCKaon = new TH2F("fdEdxTPCKaon","Energy loss vs momentum; p [GeV/c]; dE/dx [a. u.]", 40, -2, 2, 100, 0, 200);
225 fListOfHists->Add(fdEdxTPCKaon);
226
227 fbgTPCKaon = new TH2F("fbgTPCKaon", "Energy loss vs #beta#gamma (Kaons);#beta#gamma; dE/dx [a. u.]", 100, 0, 15, 100, 0, 200);
228 fListOfHists->Add(fbgTPCKaon);
229
230 //Proton
231 fdNdPtProton = new TH1F("fdNdPtProton","p_{T} distribution (Protons);p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", fXbins, fXmin, fTOFCutP);
232 fdNdPtProton->Sumw2();
233 fListOfHists->Add(fdNdPtProton);
234
235 fMassProton = new TH2F("fMassProton","Mass squared for protons after cuts vs pmean; p [GeV/c]; Mass^{2} [GeV^{2}/c^{4}]", 100, -5, 5, 100, -0.2, 1.2);
236 fListOfHists->Add(fMassProton);
237
238 fdEdxTPCProton = new TH2F("fdEdxTPCProton","Energy loss vs momentum; p [GeV/c]; dE/dx [a. u.]", 40, -2, 2, 100, 0, 200);
239 fListOfHists->Add(fdEdxTPCProton);
240
241 fbgTPCProton = new TH2F("fbgTPCProton", "Energy loss vs #beta#gamma (Protons);#beta#gamma; dE/dx [a. u.]", 100, 0, 15, 100, 0, 200);
242 fListOfHists->Add(fbgTPCProton);
243
244
245 //MC
246 if(fAnalysisMC){
247 fdNdPtMC = new TH1F("fdNdPtMC","p_{T} distribution;p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", 100, 0, 10);
248 fdNdPtMC->SetLineColor(2);
249 fdNdPtMC->Sumw2();
250 fListOfHists->Add(fdNdPtMC);
251
252 fdNdPtMCPion = new TH1F("fdNdPtMCPion","p_{T} distribution;p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", fXbins, fXmin, fTOFCutP);
253 fdNdPtMCPion->SetLineColor(2);
254 fdNdPtMCPion->Sumw2();
255 fListOfHists->Add(fdNdPtMCPion);
256
257 fdNdPtMCKaon = new TH1F("fdNdPtMCKaon","p_{T} distribution;p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", fXbins, fXmin, fTOFCutP);
258 fdNdPtMCKaon->SetLineColor(2);
259 fdNdPtMCKaon->Sumw2();
260 fListOfHists->Add(fdNdPtMCKaon);
261
262 fdNdPtMCProton = new TH1F("fdNdPtMCProton","p_{T} distribution;p_{T} [GeV/c]; #frac{1}{2#pip_{T}} #frac{1}{N_{ev}} #frac{dN}{dp_{T}}", fXbins, fXmin, fTOFCutP);
263 fdNdPtMCProton->SetLineColor(2);
264 fdNdPtMCProton->Sumw2();
265 fListOfHists->Add(fdNdPtMCProton);
266 }
267
268}
269
270//______________________________________________________________________________
18179af0 271void AliAnalysisTaskPWG4PidDetEx::UserExec(Option_t *)
de6f8090 272{
18179af0 273
274 // Create histograms
275 // Called once
276 if (fAnalysisType == "AOD") {
277 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
278 if(!fAOD){
279 Printf("%s:%d AODEvent not found in Input Manager",(char*)__FILE__,__LINE__);
280 return;
281 }
282 else{
283 // assume that the AOD is in the general output...
284 // fAOD = AODEvent();
285 // if(!fAOD){
286 // Printf("%s:%d AODEvent not found in the Output",(char*)__FILE__,__LINE__);
287 }
288 }
289 else if (fAnalysisType == "ESD"){
290 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
291 if(!fESD){
292 Printf("%s:%d ESDEvent not found in Input Manager",(char*)__FILE__,__LINE__);
293 this->Dump();
294 return;
295 }
296 }
297
de6f8090 298 // Main loop
299 // Called for each event
300 if (fDebug > 1) AliInfo("Exec() \n" );
301
302 if(fAnalysisType == "ESD") {
303 if (!fESD) {
304 Printf("ERROR: fESD not available");
305 return;
306 }
307 if(IsEventTriggered(fESD, fTriggerMode)) {
308 Printf("PWG4 PID ESD analysis");
309 AnalyzeESD(fESD);
310 }//triggered event
311 }//ESD analysis
312
313 else if(fAnalysisType == "AOD") {
314 if (!fAOD) {
315 Printf("ERROR: fAOD not available");
316 return;
317 }
318 if(IsEventTriggered(fAOD, fTriggerMode)) {
319 Printf("PWG4 PID AOD analysis");
320 AnalyzeAOD(fAOD);
321 }//triggered event
322 }//AOD analysis
323
324 // Post output data.
18179af0 325 PostData(1, fListOfHists);
de6f8090 326}
327
328//________________________________________________________________________
329void AliAnalysisTaskPWG4PidDetEx::AnalyzeESD(AliESDEvent* esd)
330{
331 // Get vertex
332 const AliESDVertex* primaryVertex = esd->GetPrimaryVertex();
333 const Double_t vertexResZ = primaryVertex->GetZRes();
334
335 // Select only events with a reconstructed vertex
336 if(vertexResZ<5.0) {
337 const Int_t nESDTracks = esd->GetNumberOfTracks();
338 for(Int_t i = 0; i < nESDTracks; i++) {
339 AliESDtrack* trackESD = esd->GetTrack(i);
340
341 //Cuts
342 UInt_t selectInfo = 0;
343 if (fTrackFilter) {
344 selectInfo = fTrackFilter->IsSelected(trackESD);
345 if (!selectInfo) continue;
346 }
347
348 if ((trackESD->Pt() < fPtCut) || (TMath::Abs(trackESD->Eta()) > fEtaCut ))
349 continue;
350
351 //Corrections
352 fEffTot->Fill(trackESD->Pt());
353 if (trackESD->GetTOFsignal() !=0)
354 fEffPID->Fill(trackESD->Pt());
355
356 if (trackESD->P() < fTOFCutP) fAccP->Fill(trackESD->Pt());
357 if (trackESD->Pt() < fTOFCutP) fAccPt->Fill(trackESD->Pt());
358
359 //Analysis
360 fdNdPt->Fill(trackESD->Pt(), 1.0/trackESD->Pt());
361
362 //TOF
363 if ((trackESD->GetIntegratedLength() == 0) || (trackESD->GetTOFsignal() == 0))
364 continue;
365
366 fMassAll->Fill(trackESD->Charge()*trackESD->P(), MassSquared(trackESD));
367
368 if ((MassSquared(trackESD) < 0.15) && (MassSquared(trackESD) > -0.15) && (trackESD->P() < fTOFCutP)){
369 fdNdPtPion->Fill(trackESD->Pt(), 1.0/trackESD->Pt());
370 fMassPion->Fill(trackESD->Charge()*trackESD->P(), MassSquared(trackESD));
371 fdEdxTPCPion->Fill(trackESD->Charge()*trackESD->P(),trackESD->GetTPCsignal());
372 fbgTPCPion->Fill(trackESD->P()/fgkPionMass,trackESD->GetTPCsignal());
373 }
374
375 if ((MassSquared(trackESD) > 0.15) && (MassSquared(trackESD) < 0.45) && (trackESD->P() < fTOFCutP)){
376 fdNdPtKaon->Fill(trackESD->Pt(), 1.0/trackESD->Pt());
377 fMassKaon->Fill(trackESD->Charge()*trackESD->P(), MassSquared(trackESD));
378 fdEdxTPCKaon->Fill(trackESD->Charge()*trackESD->P(), trackESD->GetTPCsignal());
379 fbgTPCKaon->Fill(trackESD->P()/fgkKaonMass, trackESD->GetTPCsignal());
380 //Kaon decay correction
381 fKaonDecayCorr->Fill(trackESD->Pt(), TMath::Exp(KaonDecay(trackESD)));
382 }
383
384 if ((MassSquared(trackESD) > 0.75) && (MassSquared(trackESD) < 1.05) && (trackESD->P() < fTOFCutP)){
385 fdNdPtProton->Fill(trackESD->Pt(), 1.0/trackESD->Pt());
386 fMassProton->Fill(trackESD->Charge()*trackESD->P(), MassSquared(trackESD));
387 fdEdxTPCProton->Fill(trackESD->Charge()*trackESD->P(),trackESD->GetTPCsignal());
388 fbgTPCProton->Fill(trackESD->P()/fgkProtonMass,trackESD->GetTPCsignal());
389 }
390 }//ESD track loop
391
392 //MC
393 if(fAnalysisMC){
394 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
395 if (!mcHandler) {
396 Printf("ERROR: Could not retrieve MC event handler");
397 return;
398 }
399
400 AliMCEvent* mcEvent = mcHandler->MCEvent();
401 if (!mcEvent) {
402 Printf("ERROR: Could not retrieve MC event");
403 return;
404 }
405
406 AliStack* mcStack = mcEvent->Stack();
407 if (!mcStack) {
408 Printf("ERROR: Could not retrieve MC stack");
409 return;
410 }
411
412 const Int_t nTracksMC = mcStack->GetNtrack();
413 for (Int_t iTracks = 0; iTracks < nTracksMC; iTracks++) {
414 //Cuts
415 if(!(mcStack->IsPhysicalPrimary(iTracks)))
416 continue;
417
418 TParticle* mcTrack = mcStack->Particle(iTracks);
419
420 Double_t charge = mcTrack->GetPDG()->Charge();
421 if (charge == 0)
422 continue;
423
424 if ((mcTrack->Pt() < fPtCut) || (TMath::Abs(mcTrack->Eta()) > fEtaCut ))
425 continue;
426
427 //Analysis
428 fdNdPtMC->Fill(mcTrack->Pt(), 1.0/mcTrack->Pt());
429
430 if ((mcTrack->GetPdgCode() == 211) || (mcTrack->GetPdgCode() == -211))
431 fdNdPtMCPion->Fill(mcTrack->Pt(),1.0/mcTrack->Pt());
432
433 if ((mcTrack->GetPdgCode() == 321) || (mcTrack->GetPdgCode() == -321))
434 fdNdPtMCKaon->Fill(mcTrack->Pt(),1.0/mcTrack->Pt());
435
436 if ((mcTrack->GetPdgCode() == 2212) || (mcTrack->GetPdgCode() == -2212))
437 fdNdPtMCProton->Fill(mcTrack->Pt(),1.0/mcTrack->Pt());
438 }//MC track loop
439 }//if MC
440 fEvents->Fill(0);
441 }//if vertex
442
443}
444
445//________________________________________________________________________
446void AliAnalysisTaskPWG4PidDetEx::AnalyzeAOD(AliAODEvent* aod)
447{
448 // Get vertex
449 AliAODVertex* primaryVertex = aod->GetPrimaryVertex();
450 const Double_t vertexResZ = TMath::Sqrt(primaryVertex->RotatedCovMatrixZZ());
451
452 // Select only events with a reconstructed vertex
453 if (vertexResZ < 5) {
454 //AOD
455 Int_t nGoodTracks = aod->GetNumberOfTracks();
456
457 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) {
458 AliAODTrack* trackAOD = aod->GetTrack(iTracks);
459 //Cuts
460 if (!(trackAOD->IsPrimaryCandidate()))
461 continue;
462
463 if ((trackAOD->Pt() < fPtCut) || (TMath::Abs(trackAOD->Eta()) > fEtaCut ))
464 continue;
465
466 //Corrections
467 fEffTot->Fill(trackAOD->Pt());
468 if (trackAOD->GetDetPid()->GetTOFsignal() !=0 )
469 fEffPID->Fill(trackAOD->Pt());
470
471 if (trackAOD->P() < fTOFCutP) fAccP->Fill(trackAOD->Pt());
472 if (trackAOD->Pt() < fTOFCutP) fAccPt->Fill(trackAOD->Pt());
473
474 //Analysis
475 fdNdPt->Fill(trackAOD->Pt(), 1.0/trackAOD->Pt());
476
477 //TOF
478 if ((IntegratedLength(trackAOD) == 0) || (trackAOD->GetDetPid()->GetTOFsignal() == 0))
479 continue;
480
481 fMassAll->Fill(trackAOD->Charge()*trackAOD->P(), MassSquared(trackAOD));
482
483 if ((MassSquared(trackAOD) < 0.15) && (MassSquared(trackAOD) > -0.15) && (trackAOD->P() < fTOFCutP)){
484 fdNdPtPion->Fill(trackAOD->Pt(), 1.0/trackAOD->Pt());
485 fMassPion->Fill(trackAOD->Charge()*trackAOD->P(), MassSquared(trackAOD));
486 fdEdxTPCPion->Fill(trackAOD->Charge()*trackAOD->P(),trackAOD->GetDetPid()->GetTPCsignal());
487 fbgTPCPion->Fill(trackAOD->P()/fgkPionMass,trackAOD->GetDetPid()->GetTPCsignal());
488 }
489
490 if ((MassSquared(trackAOD) > 0.15) && (MassSquared(trackAOD) < 0.45) && (trackAOD->P() < fTOFCutP)){
491 fdNdPtKaon->Fill(trackAOD->Pt(), 1.0/trackAOD->Pt());
492 fMassKaon->Fill(trackAOD->Charge()*trackAOD->P(), MassSquared(trackAOD));
493 fdEdxTPCKaon->Fill(trackAOD->Charge()*trackAOD->P(),trackAOD->GetDetPid()->GetTPCsignal());
494 fbgTPCKaon->Fill(trackAOD->P()/fgkKaonMass,trackAOD->GetDetPid()->GetTPCsignal());
495 //Kaon decay correction
496 fKaonDecayCorr->Fill(trackAOD->Pt(), TMath::Exp(KaonDecay(trackAOD)));
497 }
498
499 if ((MassSquared(trackAOD) > 0.75) && (MassSquared(trackAOD) < 1.05) && (trackAOD->P() < fTOFCutP)){
500 fdNdPtProton->Fill(trackAOD->Pt(), 1.0/trackAOD->Pt());
501 fMassProton->Fill(trackAOD->Charge()*trackAOD->P(), MassSquared(trackAOD));
502 fdEdxTPCProton->Fill(trackAOD->Charge()*trackAOD->P(),trackAOD->GetDetPid()->GetTPCsignal());
503 fbgTPCProton->Fill(trackAOD->P()/fgkProtonMass,trackAOD->GetDetPid()->GetTPCsignal());
504 }
505 }//AOD track loop
506
507 //MC
508 if(fAnalysisMC){
509 TClonesArray* farray = (TClonesArray*)aod->FindListObject("mcparticles");
510 Int_t ntrks = farray->GetEntries();
511 for(Int_t i =0 ; i < ntrks; i++){
512 AliAODMCParticle* trk = (AliAODMCParticle*)farray->At(i);
513 //Cuts
514 if (!(trk->IsPhysicalPrimary()))
515 continue;
516
517 if (trk->Charge() == 0)
518 continue;
519
520 if ((trk->Pt() < fPtCut) || (TMath::Abs(trk->Eta()) > fEtaCut ))
521 continue;
522
523 //Analysis
524 fdNdPtMC->Fill(trk->Pt(), 1.0/trk->Pt());
525
526 if ((trk->GetPdgCode() == 211) || (trk->GetPdgCode() == -211))
527 fdNdPtMCPion->Fill(trk->Pt(),1.0/trk->Pt());
528
529 if ((trk->GetPdgCode() == 321) || (trk->GetPdgCode() == -321))
530 fdNdPtMCKaon->Fill(trk->Pt(),1.0/trk->Pt());
531
532 if ((trk->GetPdgCode() == 2212) || (trk->GetPdgCode() == -2212))
533 fdNdPtMCProton->Fill(trk->Pt(),1.0/trk->Pt());
534 }//MC track loop
535 }//if MC
536 fEvents->Fill(0);
537 }//if vertex resolution
538
539}
540
541//________________________________________________________________________
542Bool_t AliAnalysisTaskPWG4PidDetEx::IsEventTriggered(AliVEvent* ev, TriggerMode trigger)
543{
544 //adapted from PWG2 (AliAnalysisTaskProton)
545
546 ULong64_t triggerMask = ev->GetTriggerMask();
547
548 // definitions from p-p.cfg
549 ULong64_t spdFO = (1 << 14);
550 ULong64_t v0left = (1 << 11);
551 ULong64_t v0right = (1 << 12);
552
553 switch (trigger) {
554 case kMB1:
555 if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
556 return kTRUE;
557 break;
558
559 case kMB2:
560 if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
561 return kTRUE;
562 break;
563
564 case kSPDFASTOR:
565 if (triggerMask & spdFO)
566 return kTRUE;
567 break;
568
569 }//switch
570
571 return kFALSE;
572}
573
574//_____________________________________________________________________________
575Double_t AliAnalysisTaskPWG4PidDetEx::IntegratedLength(AliVTrack* track) const
576{
577 Double_t intTime [5];
578 for (Int_t i = 0; i < 5; i++) intTime[i] = -100.;
579 Double_t timeElectron = 0, intLength = 0;
580
581 AliAODTrack* trackAOD = (AliAODTrack*)track;
582 trackAOD->GetDetPid()->GetIntegratedTimes(intTime);
583 timeElectron = intTime[0];
792df7ad 584 intLength = fgkC*timeElectron;
de6f8090 585
586 return intLength;
587}
588
589//_____________________________________________________________________________
590Double_t AliAnalysisTaskPWG4PidDetEx::MassSquared(AliVTrack* track) const
591{
592 Double_t beta = -10, mass = -10;
593
594 if(fAnalysisType == "ESD"){
595 AliESDtrack* trackESD = (AliESDtrack*)track;
792df7ad 596 beta = trackESD->GetIntegratedLength()/trackESD->GetTOFsignal()/fgkC;
de6f8090 597 mass = trackESD->P()*trackESD->P()*(1./(beta*beta) - 1.0);
598 }
599
600 if(fAnalysisType == "AOD"){
601 AliAODTrack* trackAOD = (AliAODTrack*)track;
792df7ad 602 beta =IntegratedLength(trackAOD)/trackAOD->GetDetPid()->GetTOFsignal()/fgkC;
de6f8090 603 mass = trackAOD->P()*trackAOD->P()*(1./(beta*beta) - 1.0);
604 }
605
606 return mass;
607}
608
609//_____________________________________________________________________________
610Double_t AliAnalysisTaskPWG4PidDetEx::KaonDecay(AliVTrack* track) const
611{
612 Double_t decay = -10;
613
614 if(fAnalysisType == "ESD"){
615 AliESDtrack* trackESD = (AliESDtrack*)track;
616 decay = trackESD->GetIntegratedLength()*fgkKaonMass/fgkCtau/trackESD->P();
617 }
618
619 if (fAnalysisType == "AOD"){
620 AliAODTrack* trackAOD = (AliAODTrack*)track;
621 decay = IntegratedLength(trackAOD)*fgkKaonMass/fgkCtau/trackAOD->P();
622 }
623
624 return decay;
625}
626
627//_____________________________________________________________________________
628void AliAnalysisTaskPWG4PidDetEx::Terminate(Option_t *)
629{
630 // Terminate loop
631 if (fDebug > 1) Printf("Terminate()");
632
633 //
634 // The followig code has now been moved to drawPid.C
635 //
636
637// fListOfHists = dynamic_cast<TList*> (GetOutputData(0));
638// if (!fListOfHists) {
639// Printf("ERROR: fListOfHists not available");
640// return;
641// }
642
643// TH1I* hevents = dynamic_cast<TH1I*> (fListOfHists->FindObject("fEvents"));
644// Int_t nEvents = (Int_t) hevents->GetBinContent(1);
645
646// const Float_t normalization = 2.0*fEtaCut*2.0*TMath::Pi();
647
648// //-----------------
649// TCanvas* c1 = new TCanvas("c1", "c1");
650// c1->cd();
651
652// TH2F* hMassAll = dynamic_cast<TH2F*> (fListOfHists->FindObject("fMassAll"));
653// hMassAll->SetLineColor(1);
654// hMassAll->SetMarkerColor(1);
655// hMassAll->DrawCopy();
656
657// TH2F* hMassPion = dynamic_cast<TH2F*> (fListOfHists->FindObject("fMassPion"));
658// hMassPion->SetLineColor(2);
659// hMassPion->SetMarkerColor(2);
660// hMassPion->SetMarkerStyle(7);
661// hMassPion->DrawCopy("same");
662
663// TH2F* hMassKaon = dynamic_cast<TH2F*> (fListOfHists->FindObject("fMassKaon"));
664// hMassKaon->SetLineColor(3);
665// hMassKaon->SetMarkerColor(3);
666// hMassKaon->SetMarkerStyle(7);
667// hMassKaon->DrawCopy("same");
668
669// TH2F* hMassProton = dynamic_cast<TH2F*> (fListOfHists->FindObject("fMassProton"));
670// hMassProton->SetLineColor(4);
671// hMassProton->SetMarkerColor(4);
672// hMassProton->SetMarkerStyle(7);
673// hMassProton->DrawCopy("same");
674
675// TLegend* legend1 = new TLegend(0.8, 0.8, 1, 1);
676// legend1->SetBorderSize(0);
677// legend1->SetFillColor(0);
678// legend1->AddEntry(hMassAll, "All","LP");
679// legend1->AddEntry(hMassPion, "Pions","LP");
680// legend1->AddEntry(hMassKaon, "Kaons","LP");
681// legend1->AddEntry(hMassProton, "Protons","LP");
682// legend1->Draw();
683
684// c1->Update();
685
686// //-----------------
687// TCanvas* c2 = new TCanvas("c2", "c2");
688// c2->cd();
689
690// TH2F* hdEdxTPCPion = dynamic_cast<TH2F*> (fListOfHists->FindObject("fdEdxTPCPion"));
691// hdEdxTPCPion->SetTitle("dE/dx vs p (TPC)");
692// hdEdxTPCPion->SetLineColor(2);
693// hdEdxTPCPion->SetMarkerColor(2);
694// hdEdxTPCPion->SetMarkerStyle(7);
695// hdEdxTPCPion->DrawCopy();
696
697// TH2F* hdEdxTPCKaon = dynamic_cast<TH2F*> (fListOfHists->FindObject("fdEdxTPCKaon"));
698// hdEdxTPCKaon->SetLineColor(3);
699// hdEdxTPCKaon->SetMarkerColor(3);
700// hdEdxTPCKaon->SetMarkerStyle(7);
701// hdEdxTPCKaon->DrawCopy("same");
702
703// TH2F* hdEdxTPCProton = dynamic_cast<TH2F*> (fListOfHists->FindObject("fdEdxTPCProton"));
704// hdEdxTPCProton->SetLineColor(4);
705// hdEdxTPCProton->SetMarkerColor(4);
706// hdEdxTPCProton->SetMarkerStyle(7);
707// hdEdxTPCProton->DrawCopy("same");
708
709// TLegend* legend2 = new TLegend(0.66, 0.66, 0.88, 0.88);
710// legend2->SetBorderSize(0);
711// legend2->SetFillColor(0);
712// legend2->AddEntry(hdEdxTPCPion, "Pions","LP");
713// legend2->AddEntry(hdEdxTPCKaon, "Kaons","LP");
714// legend2->AddEntry(hdEdxTPCProton, "Protons","LP");
715// legend2->Draw();
716
717// c2->Update();
718
719// //-----------------
720// TCanvas* c3 = new TCanvas("c3", "c3");
721// c3->cd();
722
723// TH2F* hdEdxTPCbgPion = dynamic_cast<TH2F*> (fListOfHists->FindObject("fbgTPCPion"));
724// hdEdxTPCbgPion->SetTitle("dE/dx vs #beta#gamma (TPC)");
725// hdEdxTPCbgPion->SetLineColor(2);
726// hdEdxTPCbgPion->SetMarkerColor(2);
727// hdEdxTPCbgPion->SetMarkerStyle(7);
728// hdEdxTPCbgPion->DrawCopy();
729
730// TH2F* hdEdxTPCbgKaon = dynamic_cast<TH2F*> (fListOfHists->FindObject("fbgTPCKaon"));
731// hdEdxTPCbgKaon->SetLineColor(3);
732// hdEdxTPCbgKaon->SetMarkerColor(3);
733// hdEdxTPCbgKaon->SetMarkerStyle(7);
734// hdEdxTPCbgKaon->DrawCopy("same");
735
736// TH2F* hdEdxTPCbgProton = dynamic_cast<TH2F*> (fListOfHists->FindObject("fbgTPCProton"));
737// hdEdxTPCbgProton->SetLineColor(4);
738// hdEdxTPCbgProton->SetMarkerColor(4);
739// hdEdxTPCbgProton->SetMarkerStyle(7);
740// hdEdxTPCbgProton->DrawCopy("same");
741
742// TLegend* legend3 = new TLegend(0.66, 0.66, 0.88, 0.88);
743// legend3->SetBorderSize(0);
744// legend3->SetFillColor(0);
745// legend3->AddEntry(hdEdxTPCbgPion, "Pions","LP");
746// legend3->AddEntry(hdEdxTPCbgKaon, "Kaons","LP");
747// legend3->AddEntry(hdEdxTPCbgProton, "Protons","LP");
748// legend3->Draw();
749
750// c3->Update();
751
752// //-----------------
753// TCanvas* c4 = new TCanvas("c4", "c4", 100, 100, 500, 900);
754// c4->Divide(1,3);
755
756// c4->cd(1);
757// TH1F* hAccepatncePt = dynamic_cast<TH1F*> (fListOfHists->FindObject("fAccPt"));
758// TH1F* hAcceptanceP = dynamic_cast<TH1F*> (fListOfHists->FindObject("fAccP"));
759// hAccepatncePt->Divide(hAcceptanceP);
760// hAccepatncePt->SetTitle("Acceptance correction");
761// hAccepatncePt->GetYaxis()->SetTitle("Acceptance correction");
762// hAccepatncePt->DrawCopy();
763
764// c4->cd(2);
765// TH1F* hEfficiencyPID = dynamic_cast<TH1F*> (fListOfHists->FindObject("fEffPID"));
766// TH1F* hEfficiencyTot = dynamic_cast<TH1F*> (fListOfHists->FindObject("fEffTot"));
767// hEfficiencyPID->Divide(hEfficiencyTot);
768// hEfficiencyPID->SetTitle("Efficiency correction");
769// hEfficiencyPID->GetYaxis()->SetTitle("Efficiency correction");
770// hEfficiencyPID->DrawCopy();
771
772// c4->cd(3);
773// TProfile* hKDecayCorr = dynamic_cast<TProfile*> (fListOfHists->FindObject("fKaonDecayCorr"));
774// hKDecayCorr->SetTitle("Kaon decay correction");
775// hKDecayCorr->GetYaxis()->SetTitle("Kaon decay correction");
776// hKDecayCorr->GetXaxis()->SetTitle(" p_{T} [GeV/c]");
777// hKDecayCorr->DrawCopy();
778
779// c4->Update();
780
781// //---------------------
782// TCanvas* c5 = new TCanvas("c5", "c5", 100, 100, 600, 900);
783// c5->Divide(1,2);
784
785// c5->cd(1);
786// TH1F* hPtPionNoCorr = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtPion"));
787// hPtPionNoCorr->Scale(1.0/nEvents/normalization/hPtPionNoCorr->GetBinWidth(1));
788// hPtPionNoCorr->SetTitle("p_{T} distribution (no corrections)");
789// hPtPionNoCorr->SetLineColor(2);
790// hPtPionNoCorr->GetYaxis()->SetRangeUser(1E-3,1E1);
791// hPtPionNoCorr->DrawCopy();
792
793// TH1F* hPtKaonNoCorr = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtKaon"));
794// hPtKaonNoCorr->Scale(1.0/nEvents/normalization/hPtKaonNoCorr->GetBinWidth(1));
795// hPtKaonNoCorr->SetLineColor(3);
796// hPtKaonNoCorr->GetYaxis()->SetRangeUser(1E-3,1E1);
797// hPtKaonNoCorr->DrawCopy("same");
798
799// TH1F* hPtProtonNoCorr = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtProton"));
800// hPtProtonNoCorr->Scale(1.0/nEvents/normalization/hPtProtonNoCorr->GetBinWidth(1));
801// hPtProtonNoCorr->SetLineColor(4);
802// hPtProtonNoCorr->GetYaxis()->SetRangeUser(1E-3,1E1);
803// hPtProtonNoCorr->DrawCopy("same");
804
805// TH1F* hPt = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPt"));
806// hPt->Scale(1.0/nEvents/normalization/hPt->GetBinWidth(1));
807// hPt->GetYaxis()->SetRangeUser(1E-3,1E1);
808// hPt->DrawCopy("same");
809
810// TLegend* legend4 = new TLegend(0.63, 0.63, 0.88, 0.88);
811// legend4->SetBorderSize(0);
812// legend4->SetFillColor(0);
813// legend4->AddEntry(hPt, "p_{T} dist (all)","L");
814// legend4->AddEntry(hPtPionNoCorr, "p_{T} dist (Pions)","L");
815// legend4->AddEntry(hPtKaonNoCorr, "p_{T} dist (Kaons)","L");
816// legend4->AddEntry(hPtProtonNoCorr, "p_{T} dist (Protons)","L");
817// legend4->Draw();
818
819// gPad->SetLogy();
820
821
822// c5->cd(2);
823// TH1F* hPtPionCorr = static_cast<TH1F*>(hPtPionNoCorr->Clone());
824// hPtPionCorr->SetTitle("p_{T} distribution (with corrections)");
825// hPtPionCorr->Multiply(hAccepatncePt);
826// hPtPionCorr->Divide(hEfficiencyPID);
827// hPtPionCorr->GetYaxis()->SetRangeUser(1E-2,1E1);
828// hPtPionCorr->DrawCopy();
829
830// TH1F* hPtKaonCorr = static_cast<TH1F*>(hPtKaonNoCorr->Clone());
831// hPtKaonCorr->Multiply(hAccepatncePt);
832// hPtKaonCorr->Divide(hEfficiencyPID);
833// hPtKaonCorr->Multiply(hKDecayCorr);
834// hPtKaonCorr->GetYaxis()->SetRangeUser(1E-2,1E1);
835// hPtKaonCorr->DrawCopy("same");
836
837// TH1F* hPtProtonCorr = static_cast<TH1F*>(hPtProtonNoCorr->Clone());
838// hPtProtonCorr->Multiply(hAccepatncePt);
839// hPtProtonCorr->Divide(hEfficiencyPID);
840// hPtProtonCorr->GetYaxis()->SetRangeUser(1E-2,1E1);
841// hPtProtonCorr->DrawCopy("same");
842
843// hPt->GetYaxis()->SetRangeUser(1E-2,1E1);
844// hPt->DrawCopy("same");
845
846// TLegend* legend5 = new TLegend(0.63, 0.63, 0.88, 0.88);
847// legend5->SetBorderSize(0);
848// legend5->SetFillColor(0);
849// legend5->AddEntry(hPt, "p_{T} dist (all)","L");
850// legend5->AddEntry(hPtPionCorr, "p_{T} dist (Pions)","L");
851// legend5->AddEntry(hPtKaonCorr, "p_{T} dist (Kaons)","L");
852// legend5->AddEntry(hPtProtonCorr, "p_{T} dist (Protons)","L");
853// legend5->Draw();
854
855// gPad->SetLogy();
856
857// c5->Update();
858
859// //-----------------
860// if (fAnalysisMC){
861// TCanvas* c6 = new TCanvas("c6", "c6", 100, 100, 1200, 800);
862// c6->Divide(2,2);
863
864// c6->cd(1);
865// TH1F* hPt_clone = static_cast<TH1F*>(hPt->Clone());
866// hPt_clone->GetYaxis()->SetRangeUser(1E-5,1E1);
867// hPt_clone->SetTitle("p_{T} distribution (all)");
868// hPt_clone->SetLineColor(1);
869// hPt_clone->DrawCopy();
870
871// TH1F* hPtMC = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtMC"));
872// hPtMC->Scale(1.0/nEvents/normalization/hPtMC->GetBinWidth(1));
873// hPtMC->GetYaxis()->SetRangeUser(1E-5,1E1);
874// hPtMC->DrawCopy("same");
875
876// TLegend* legend6 = new TLegend(0.57, 0.57, 0.87, 0.87);
877// legend6->SetBorderSize(0);
878// legend6->SetFillColor(0);
879// legend6->AddEntry(hPt_clone, "p_{T} dist (Rec)", "L");
880// legend6->AddEntry(hPtMC, "p_{T} dist (MC)", "L");
881// legend6->Draw();
882
883// gPad->SetLogy();
884
885// c6->cd(2);
886// TH1F* hPtPion_clone = static_cast<TH1F*>(hPtPionCorr->Clone());
887// hPtPion_clone->GetYaxis()->SetRangeUser(1E-2,1E1);
888// hPtPion_clone->SetTitle("p_{T} distribution (Pions)");
889// hPtPion_clone->SetLineColor(1);
890// hPtPion_clone->DrawCopy();
891
892// TH1F* hPtMCPion = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtMCPion"));
893// hPtMCPion->Scale(1.0/nEvents/normalization/hPtMCPion->GetBinWidth(1));
894// hPtMCPion->GetYaxis()->SetRangeUser(1E-2,1E1);
895// hPtMCPion->DrawCopy("same");
896
897// TLegend* legend7 = new TLegend(0.57, 0.57, 0.87, 0.87);
898// legend7->SetBorderSize(0);
899// legend7->SetFillColor(0);
900// legend7->AddEntry(hPtPion_clone, "p_{T} dist (Rec)", "L");
901// legend7->AddEntry(hPtMCPion, "p_{T} dist (MC)", "L");
902// legend7->Draw();
903
904// gPad->SetLogy();
905
906// c6->cd(3);
907// TH1F* hPtKaon_clone = static_cast<TH1F*>(hPtKaonCorr->Clone());
908// hPtKaon_clone->GetYaxis()->SetRangeUser(1E-2,1E0);
909// hPtKaon_clone->SetLineColor(1);
910// hPtKaon_clone->DrawCopy();
911
912// TH1F* hPtMCKaon = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtMCKaon"));
913// hPtMCKaon->Scale(1.0/nEvents/normalization/hPtMCKaon->GetBinWidth(1));
914// hPtMCKaon->GetYaxis()->SetRangeUser(1E-2,1E0);
915// hPtMCKaon->DrawCopy("same");
916
917// TLegend* legend8 = new TLegend(0.57, 0.57, 0.87, 0.87);
918// legend8->SetBorderSize(0);
919// legend8->SetFillColor(0);
920// legend8->AddEntry(hPtKaon_clone, "p_{T} dist (Rec)", "L");
921// legend8->AddEntry(hPtMCKaon, "p_{T} dist (MC)", "L");
922// legend8->Draw();
923
924// gPad->SetLogy();
925
926// c6->cd(4);
927// TH1F* hPtProton_clone = static_cast<TH1F*>(hPtProtonCorr->Clone());
928// hPtProton_clone->GetYaxis()->SetRangeUser(1E-2,1E-1);
929// hPtProton_clone->SetLineColor(1);
930// hPtProton_clone->DrawCopy();
931
932// TH1F* hPtMCProton = dynamic_cast<TH1F*> (fListOfHists->FindObject("fdNdPtMCProton"));
933// hPtMCProton->Scale(1.0/nEvents/normalization/hPtMCProton->GetBinWidth(1));
934// hPtMCProton->GetYaxis()->SetRangeUser(1E-2,1E-1);
935// hPtMCProton->DrawCopy("same");
936
937// TLegend* legend9 = new TLegend(0.2, 0.25, 0.5, 0.55);
938// legend9->SetBorderSize(0);
939// legend9->SetFillColor(0);
940// legend9->AddEntry(hPtProton_clone, "p_{T} dist (Rec)", "L");
941// legend9->AddEntry(hPtMCProton, "p_{T} dist (MC)", "L");
942// legend9->Draw();
943
944// gPad->SetLogy();
945
946// c6->Update();
947// }
948
949}