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 **************************************************************************/
30 #include <TDatabasePDG.h>
31 #include <TParticle.h>
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"
41 #include "AliAODInputHandler.h"
42 #include "AliAODEvent.h"
43 #include "AliAODVertex.h"
44 #include "AliAODMCParticle.h"
47 #include "AliAnalysisTaskPWG4PidDetEx.h"
54 // Analysis class example for PID using detector signals.
55 // Works on ESD and AOD; has a flag for MC analysis
57 // Alexandru.Dobrin@hep.lu.se
58 // Peter.Christiansen@hep.lu.se
61 ClassImp(AliAnalysisTaskPWG4PidDetEx)
63 const Double_t AliAnalysisTaskPWG4PidDetEx::fgkCtau = 370; // distance for kaon decay
64 const Double_t AliAnalysisTaskPWG4PidDetEx::fgkPionMass = 1.39570000000000000e-01;
65 const Double_t AliAnalysisTaskPWG4PidDetEx::fgkKaonMass = 4.93599999999999983e-01;
66 const Double_t AliAnalysisTaskPWG4PidDetEx::fgkProtonMass = 9.38270000000000048e-01;
67 const Double_t AliAnalysisTaskPWG4PidDetEx::fgkC = 2.99792458e-2;
73 //_____________________________________________________________________________
74 AliAnalysisTaskPWG4PidDetEx::AliAnalysisTaskPWG4PidDetEx():
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)
95 // Default constructor (should not be used)
98 //______________________________________________________________________________
99 AliAnalysisTaskPWG4PidDetEx::AliAnalysisTaskPWG4PidDetEx(const char *name):
100 AliAnalysisTaskSE(name),
111 fAnalysisType("ESD"),
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)
120 // Output slot #0 writes into a TList
121 DefineOutput(1, TList::Class());
125 //_____________________________________________________________________________
126 AliAnalysisTaskPWG4PidDetEx::~AliAnalysisTaskPWG4PidDetEx()
131 // //______________________________________________________________________________
132 // void AliAnalysisTaskPWG4PidDetEx::ConnectInputData(Option_t *)
134 // // Connect AOD here
136 // if (fDebug > 1) AliInfo("ConnectInputData() \n");
138 // TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
140 // Printf("ERROR: Could not read chain from input slot 0");
143 // if(fAnalysisType == "ESD") {
144 // AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
146 // Printf("ERROR: Could not get ESDInputHandler");
148 // fESD = esdH->GetEvent();
150 // else if(fAnalysisType == "AOD") {
151 // AliAODInputHandler *aodH = dynamic_cast<AliAODInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
153 // Printf("ERROR: Could not get AODInputHandler");
155 // fAOD = aodH->GetEvent();
161 //______________________________________________________________________________
162 void AliAnalysisTaskPWG4PidDetEx::UserCreateOutputObjects()
165 fListOfHists = new TList();
167 fEvents = new TH1I("fEvents","Number of analyzed events; Events; Counts", 1, 0, 1);
168 fListOfHists->Add(fEvents);
170 fEffTot = new TH1F ("fEffTot","Efficiency; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
172 fListOfHists->Add(fEffTot);
174 fEffPID = new TH1F ("fEffPID","Efficiency; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
176 fListOfHists->Add(fEffPID);
178 fAccP = new TH1F ("fAccP","Acceptance; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
180 fListOfHists->Add(fAccP);
182 fAccPt = new TH1F ("fAccPt","Acceptance; p_{T} [GeV/c]; Counts", fXbins, fXmin, fTOFCutP);
184 fListOfHists->Add(fAccPt);
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);
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);
191 fListOfHists->Add(fdNdPt);
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);
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);
199 fListOfHists->Add(fdNdPtPion);
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);
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);
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);
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);
213 fListOfHists->Add(fdNdPtKaon);
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);
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);
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);
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);
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);
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);
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);
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);
244 fListOfHists->Add(fdNdPtMC);
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);
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);
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);
264 //______________________________________________________________________________
265 void AliAnalysisTaskPWG4PidDetEx::UserExec(Option_t *)
270 if (fAnalysisType == "AOD") {
271 fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
273 Printf("%s:%d AODEvent not found in Input Manager",(char*)__FILE__,__LINE__);
277 // assume that the AOD is in the general output...
278 // fAOD = AODEvent();
280 // Printf("%s:%d AODEvent not found in the Output",(char*)__FILE__,__LINE__);
283 else if (fAnalysisType == "ESD"){
284 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
286 Printf("%s:%d ESDEvent not found in Input Manager",(char*)__FILE__,__LINE__);
293 // Called for each event
294 if (fDebug > 1) AliInfo("Exec() \n" );
296 if(fAnalysisType == "ESD") {
298 Printf("ERROR: fESD not available");
301 if(IsEventTriggered(fESD, fTriggerMode)) {
302 Printf("PWG4 PID ESD analysis");
307 else if(fAnalysisType == "AOD") {
309 Printf("ERROR: fAOD not available");
312 if(IsEventTriggered(fAOD, fTriggerMode)) {
313 Printf("PWG4 PID AOD analysis");
319 PostData(1, fListOfHists);
322 //________________________________________________________________________
323 void AliAnalysisTaskPWG4PidDetEx::AnalyzeESD(AliESDEvent* esd)
326 const AliESDVertex* primaryVertex = esd->GetPrimaryVertex();
327 const Double_t vertexResZ = primaryVertex->GetZRes();
329 // Select only events with a reconstructed vertex
331 const Int_t nESDTracks = esd->GetNumberOfTracks();
332 for(Int_t i = 0; i < nESDTracks; i++) {
333 AliESDtrack* trackESD = esd->GetTrack(i);
336 UInt_t selectInfo = 0;
338 selectInfo = fTrackFilter->IsSelected(trackESD);
339 if (!selectInfo) continue;
342 if ((trackESD->Pt() < fPtCut) || (TMath::Abs(trackESD->Eta()) > fEtaCut ))
346 fEffTot->Fill(trackESD->Pt());
347 if (trackESD->GetTOFsignal() !=0)
348 fEffPID->Fill(trackESD->Pt());
350 if (trackESD->P() < fTOFCutP) fAccP->Fill(trackESD->Pt());
351 if (trackESD->Pt() < fTOFCutP) fAccPt->Fill(trackESD->Pt());
354 fdNdPt->Fill(trackESD->Pt(), 1.0/trackESD->Pt());
357 if ((trackESD->GetIntegratedLength() == 0) || (trackESD->GetTOFsignal() == 0))
360 fMassAll->Fill(trackESD->Charge()*trackESD->P(), MassSquared(trackESD));
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());
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)));
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());
388 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
390 Printf("ERROR: Could not retrieve MC event handler");
394 AliMCEvent* mcEvent = mcHandler->MCEvent();
396 Printf("ERROR: Could not retrieve MC event");
400 AliStack* mcStack = mcEvent->Stack();
402 Printf("ERROR: Could not retrieve MC stack");
406 const Int_t nTracksMC = mcStack->GetNtrack();
407 for (Int_t iTracks = 0; iTracks < nTracksMC; iTracks++) {
409 if(!(mcStack->IsPhysicalPrimary(iTracks)))
412 TParticle* mcTrack = mcStack->Particle(iTracks);
414 Double_t charge = mcTrack->GetPDG()->Charge();
418 if ((mcTrack->Pt() < fPtCut) || (TMath::Abs(mcTrack->Eta()) > fEtaCut ))
422 fdNdPtMC->Fill(mcTrack->Pt(), 1.0/mcTrack->Pt());
424 if ((mcTrack->GetPdgCode() == 211) || (mcTrack->GetPdgCode() == -211))
425 fdNdPtMCPion->Fill(mcTrack->Pt(),1.0/mcTrack->Pt());
427 if ((mcTrack->GetPdgCode() == 321) || (mcTrack->GetPdgCode() == -321))
428 fdNdPtMCKaon->Fill(mcTrack->Pt(),1.0/mcTrack->Pt());
430 if ((mcTrack->GetPdgCode() == 2212) || (mcTrack->GetPdgCode() == -2212))
431 fdNdPtMCProton->Fill(mcTrack->Pt(),1.0/mcTrack->Pt());
439 //________________________________________________________________________
440 void AliAnalysisTaskPWG4PidDetEx::AnalyzeAOD(AliAODEvent* aod)
443 AliAODVertex* primaryVertex = aod->GetPrimaryVertex();
444 const Double_t vertexResZ = TMath::Sqrt(primaryVertex->RotatedCovMatrixZZ());
446 // Select only events with a reconstructed vertex
447 if (vertexResZ < 5) {
449 Int_t nGoodTracks = aod->GetNumberOfTracks();
451 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) {
452 AliAODTrack* trackAOD = aod->GetTrack(iTracks);
454 if (!(trackAOD->IsPrimaryCandidate()))
457 if ((trackAOD->Pt() < fPtCut) || (TMath::Abs(trackAOD->Eta()) > fEtaCut ))
461 fEffTot->Fill(trackAOD->Pt());
462 if (trackAOD->GetDetPid()->GetTOFsignal() !=0 )
463 fEffPID->Fill(trackAOD->Pt());
465 if (trackAOD->P() < fTOFCutP) fAccP->Fill(trackAOD->Pt());
466 if (trackAOD->Pt() < fTOFCutP) fAccPt->Fill(trackAOD->Pt());
469 fdNdPt->Fill(trackAOD->Pt(), 1.0/trackAOD->Pt());
472 if ((IntegratedLength(trackAOD) == 0) || (trackAOD->GetDetPid()->GetTOFsignal() == 0))
475 fMassAll->Fill(trackAOD->Charge()*trackAOD->P(), MassSquared(trackAOD));
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());
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)));
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());
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);
508 if (!(trk->IsPhysicalPrimary()))
511 if (trk->Charge() == 0)
514 if ((trk->Pt() < fPtCut) || (TMath::Abs(trk->Eta()) > fEtaCut ))
518 fdNdPtMC->Fill(trk->Pt(), 1.0/trk->Pt());
520 if ((trk->GetPdgCode() == 211) || (trk->GetPdgCode() == -211))
521 fdNdPtMCPion->Fill(trk->Pt(),1.0/trk->Pt());
523 if ((trk->GetPdgCode() == 321) || (trk->GetPdgCode() == -321))
524 fdNdPtMCKaon->Fill(trk->Pt(),1.0/trk->Pt());
526 if ((trk->GetPdgCode() == 2212) || (trk->GetPdgCode() == -2212))
527 fdNdPtMCProton->Fill(trk->Pt(),1.0/trk->Pt());
531 }//if vertex resolution
535 //________________________________________________________________________
536 Bool_t AliAnalysisTaskPWG4PidDetEx::IsEventTriggered(AliVEvent* ev, TriggerMode trigger)
538 //adapted from PWG2 (AliAnalysisTaskProton)
540 ULong64_t triggerMask = ev->GetTriggerMask();
542 // definitions from p-p.cfg
543 ULong64_t spdFO = (1 << 14);
544 ULong64_t v0left = (1 << 11);
545 ULong64_t v0right = (1 << 12);
549 if (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)))
554 if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
559 if (triggerMask & spdFO)
568 //_____________________________________________________________________________
569 Double_t AliAnalysisTaskPWG4PidDetEx::IntegratedLength(AliVTrack* track) const
571 Double_t intTime [5];
572 for (Int_t i = 0; i < 5; i++) intTime[i] = -100.;
573 Double_t timeElectron = 0, intLength = 0;
575 AliAODTrack* trackAOD = (AliAODTrack*)track;
576 trackAOD->GetDetPid()->GetIntegratedTimes(intTime);
577 timeElectron = intTime[0];
578 intLength = fgkC*timeElectron;
583 //_____________________________________________________________________________
584 Double_t AliAnalysisTaskPWG4PidDetEx::MassSquared(AliVTrack* track) const
586 Double_t beta = -10, mass = -10;
588 if(fAnalysisType == "ESD"){
589 AliESDtrack* trackESD = (AliESDtrack*)track;
590 beta = trackESD->GetIntegratedLength()/trackESD->GetTOFsignal()/fgkC;
591 mass = trackESD->P()*trackESD->P()*(1./(beta*beta) - 1.0);
594 if(fAnalysisType == "AOD"){
595 AliAODTrack* trackAOD = (AliAODTrack*)track;
596 beta =IntegratedLength(trackAOD)/trackAOD->GetDetPid()->GetTOFsignal()/fgkC;
597 mass = trackAOD->P()*trackAOD->P()*(1./(beta*beta) - 1.0);
603 //_____________________________________________________________________________
604 Double_t AliAnalysisTaskPWG4PidDetEx::KaonDecay(AliVTrack* track) const
606 Double_t decay = -10;
608 if(fAnalysisType == "ESD"){
609 AliESDtrack* trackESD = (AliESDtrack*)track;
610 decay = trackESD->GetIntegratedLength()*fgkKaonMass/fgkCtau/trackESD->P();
613 if (fAnalysisType == "AOD"){
614 AliAODTrack* trackAOD = (AliAODTrack*)track;
615 decay = IntegratedLength(trackAOD)*fgkKaonMass/fgkCtau/trackAOD->P();
621 //_____________________________________________________________________________
622 void AliAnalysisTaskPWG4PidDetEx::Terminate(Option_t *)
625 if (fDebug > 1) Printf("Terminate()");
628 // The followig code has now been moved to drawPid.C
631 // fListOfHists = dynamic_cast<TList*> (GetOutputData(0));
632 // if (!fListOfHists) {
633 // Printf("ERROR: fListOfHists not available");
637 // TH1I* hevents = dynamic_cast<TH1I*> (fListOfHists->FindObject("fEvents"));
638 // Int_t nEvents = (Int_t) hevents->GetBinContent(1);
640 // const Float_t normalization = 2.0*fEtaCut*2.0*TMath::Pi();
642 // //-----------------
643 // TCanvas* c1 = new TCanvas("c1", "c1");
646 // TH2F* hMassAll = dynamic_cast<TH2F*> (fListOfHists->FindObject("fMassAll"));
647 // hMassAll->SetLineColor(1);
648 // hMassAll->SetMarkerColor(1);
649 // hMassAll->DrawCopy();
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");
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");
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");
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");
680 // //-----------------
681 // TCanvas* c2 = new TCanvas("c2", "c2");
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();
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");
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");
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");
713 // //-----------------
714 // TCanvas* c3 = new TCanvas("c3", "c3");
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();
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");
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");
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");
746 // //-----------------
747 // TCanvas* c4 = new TCanvas("c4", "c4", 100, 100, 500, 900);
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();
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();
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();
775 // //---------------------
776 // TCanvas* c5 = new TCanvas("c5", "c5", 100, 100, 600, 900);
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();
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");
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");
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");
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");
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();
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");
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");
837 // hPt->GetYaxis()->SetRangeUser(1E-2,1E1);
838 // hPt->DrawCopy("same");
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");
853 // //-----------------
855 // TCanvas* c6 = new TCanvas("c6", "c6", 100, 100, 1200, 800);
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();
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");
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");
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();
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");
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");
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();
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");
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");
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();
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");
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");