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