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 pures 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 pure. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
18 // AliAnalysisTaskPIDflowQA:
21 //origin: Marek Chojnacki, Marek.Chojnacki@cern.ch
22 //modified: Mikolaj Krzewicki, Mikolaj.Krzewicki@cern.ch
26 #include "AliAnalysisTaskPIDflowQA.h"
27 #include "AliAnalysisManager.h"
28 #include "AliESDEvent.h"
30 #include "AliMCEvent.h"
35 #include "AliVEvent.h"
37 #include "AliCDBManager.h"
38 #include "AliFlowEventCuts.h"
39 #include "AliFlowTrackCuts.h"
40 #include "AliVEventHandler.h"
41 #include "AliInputEventHandler.h"
45 ClassImp( AliAnalysisTaskPIDflowQA)
47 //________________________________________________________________________
48 AliAnalysisTaskPIDflowQA:: AliAnalysisTaskPIDflowQA():
49 AliAnalysisTaskSE("AliAnalysisTaskPIDflowQA"),
55 fUseDebugFile(kFALSE),
79 fTOFbetaAfterElectronsCuts(NULL),
80 fTOFbetaAfterPionCuts(NULL),
81 fTOFbetaAfterKaonCuts(NULL),
82 fTOFbetaAfterProtonCuts(NULL),
83 fTPCsignalAfterPionCuts(NULL),
84 fTPCsignalAfterKaonCuts(NULL),
85 fTPCsignalAfterProtonCuts(NULL),
86 fTOFbetaAfterElectronsCuts1(NULL),
87 fTOFbetaAfterPionCuts1(NULL),
88 fTOFbetaAfterKaonCuts1(NULL),
89 fTOFbetaAfterProtonCuts1(NULL),
90 fTPCsignalAfterPionCuts1(NULL),
91 fTPCsignalAfterKaonCuts1(NULL),
92 fTPCsignalAfterProtonCuts1(NULL),
94 fTOFbetaPiafter(NULL),
97 fTPCsignalPiafter(NULL),
98 fTPCsignalKafter(NULL),
99 fTPCsignalPafter(NULL),
100 fTOFyieldSelEmcE(NULL),
101 fTOFyieldSelPimcE(NULL),
102 fTOFyieldSelKmcE(NULL),
103 fTOFyieldSelPmcE(NULL),
104 fTOFyieldSelEmcM(NULL),
105 fTOFyieldSelPimcM(NULL),
106 fTOFyieldSelKmcM(NULL),
107 fTOFyieldSelPmcM(NULL),
108 fTOFyieldSelEmcPi(NULL),
109 fTOFyieldSelPimcPi(NULL),
110 fTOFyieldSelKmcPi(NULL),
111 fTOFyieldSelPmcPi(NULL),
112 fTOFyieldSelEmcK(NULL),
113 fTOFyieldSelPimcK(NULL),
114 fTOFyieldSelKmcK(NULL),
115 fTOFyieldSelPmcK(NULL),
116 fTOFyieldSelEmcP(NULL),
117 fTOFyieldSelPimcP(NULL),
118 fTOFyieldSelKmcP(NULL),
119 fTOFyieldSelPmcP(NULL),
120 fTPCyieldSelEmcE(NULL),
121 fTPCyieldSelPimcE(NULL),
122 fTPCyieldSelKmcE(NULL),
123 fTPCyieldSelPmcE(NULL),
124 fTPCyieldSelEmcM(NULL),
125 fTPCyieldSelPimcM(NULL),
126 fTPCyieldSelKmcM(NULL),
127 fTPCyieldSelPmcM(NULL),
128 fTPCyieldSelEmcPi(NULL),
129 fTPCyieldSelPimcPi(NULL),
130 fTPCyieldSelKmcPi(NULL),
131 fTPCyieldSelPmcPi(NULL),
132 fTPCyieldSelEmcK(NULL),
133 fTPCyieldSelPimcK(NULL),
134 fTPCyieldSelKmcK(NULL),
135 fTPCyieldSelPmcK(NULL),
136 fTPCyieldSelEmcP(NULL),
137 fTPCyieldSelPimcP(NULL),
138 fTPCyieldSelKmcP(NULL),
139 fTPCyieldSelPmcP(NULL),
140 fTPCdedxAfterTOFpidPions(NULL),
141 fTPCdedxAfterTOFpidKaons(NULL),
142 fTPCdedxAfterTOFpidProtons(NULL),
145 fTPCvsGlobalMult(NULL),
146 fStandardGlobalCuts(NULL),
147 fStandardTPCCuts(NULL),
148 fCutsTOFElectrons(NULL),
151 fCutsTOFProtons(NULL),
152 fCutsTPCElectrons(NULL),
155 fCutsTPCProtons(NULL),
161 //________________________________________________________________________
162 AliAnalysisTaskPIDflowQA:: AliAnalysisTaskPIDflowQA(const char *name):
163 AliAnalysisTaskSE(name),
169 fUseDebugFile(kFALSE),
175 fTPCsignalPimc(NULL),
193 fTOFbetaAfterElectronsCuts(NULL),
194 fTOFbetaAfterPionCuts(NULL),
195 fTOFbetaAfterKaonCuts(NULL),
196 fTOFbetaAfterProtonCuts(NULL),
197 fTPCsignalAfterPionCuts(NULL),
198 fTPCsignalAfterKaonCuts(NULL),
199 fTPCsignalAfterProtonCuts(NULL),
200 fTOFbetaAfterElectronsCuts1(NULL),
201 fTOFbetaAfterPionCuts1(NULL),
202 fTOFbetaAfterKaonCuts1(NULL),
203 fTOFbetaAfterProtonCuts1(NULL),
204 fTPCsignalAfterPionCuts1(NULL),
205 fTPCsignalAfterKaonCuts1(NULL),
206 fTPCsignalAfterProtonCuts1(NULL),
207 fTOFbetaEafter(NULL),
208 fTOFbetaPiafter(NULL),
209 fTOFbetaKafter(NULL),
210 fTOFbetaPafter(NULL),
211 fTPCsignalPiafter(NULL),
212 fTPCsignalKafter(NULL),
213 fTPCsignalPafter(NULL),
214 fTOFyieldSelEmcE(NULL),
215 fTOFyieldSelPimcE(NULL),
216 fTOFyieldSelKmcE(NULL),
217 fTOFyieldSelPmcE(NULL),
218 fTOFyieldSelEmcM(NULL),
219 fTOFyieldSelPimcM(NULL),
220 fTOFyieldSelKmcM(NULL),
221 fTOFyieldSelPmcM(NULL),
222 fTOFyieldSelEmcPi(NULL),
223 fTOFyieldSelPimcPi(NULL),
224 fTOFyieldSelKmcPi(NULL),
225 fTOFyieldSelPmcPi(NULL),
226 fTOFyieldSelEmcK(NULL),
227 fTOFyieldSelPimcK(NULL),
228 fTOFyieldSelKmcK(NULL),
229 fTOFyieldSelPmcK(NULL),
230 fTOFyieldSelEmcP(NULL),
231 fTOFyieldSelPimcP(NULL),
232 fTOFyieldSelKmcP(NULL),
233 fTOFyieldSelPmcP(NULL),
234 fTPCyieldSelEmcE(NULL),
235 fTPCyieldSelPimcE(NULL),
236 fTPCyieldSelKmcE(NULL),
237 fTPCyieldSelPmcE(NULL),
238 fTPCyieldSelEmcM(NULL),
239 fTPCyieldSelPimcM(NULL),
240 fTPCyieldSelKmcM(NULL),
241 fTPCyieldSelPmcM(NULL),
242 fTPCyieldSelEmcPi(NULL),
243 fTPCyieldSelPimcPi(NULL),
244 fTPCyieldSelKmcPi(NULL),
245 fTPCyieldSelPmcPi(NULL),
246 fTPCyieldSelEmcK(NULL),
247 fTPCyieldSelPimcK(NULL),
248 fTPCyieldSelKmcK(NULL),
249 fTPCyieldSelPmcK(NULL),
250 fTPCyieldSelEmcP(NULL),
251 fTPCyieldSelPimcP(NULL),
252 fTPCyieldSelKmcP(NULL),
253 fTPCyieldSelPmcP(NULL),
254 fTPCdedxAfterTOFpidPions(NULL),
255 fTPCdedxAfterTOFpidKaons(NULL),
256 fTPCdedxAfterTOFpidProtons(NULL),
259 fTPCvsGlobalMult(NULL),
260 fStandardGlobalCuts(NULL),
261 fStandardTPCCuts(NULL),
262 fCutsTOFElectrons(NULL),
265 fCutsTOFProtons(NULL),
266 fCutsTPCElectrons(NULL),
269 fCutsTPCProtons(NULL),
273 fESDpid=new AliESDpid();
276 fESDpid->GetTPCResponse().SetBetheBlochParameters(0.0283086,
282 //fESDpid->GetTPCResponse().SetBetheBlochParameters(1.28949/50.,
284 // TMath::Exp(-3.21763e+01),
288 DefineOutput(1, TList::Class());
291 //________________________________________________________________________
292 void AliAnalysisTaskPIDflowQA::UserCreateOutputObjects()
294 //UserCreateOutputObject
295 if (fOutputList) fOutputList->Delete();
297 fOutputList=new TList();
298 fOutputList->SetOwner(kTRUE);
302 const Int_t npredec=50;
303 Double_t tabx[ndec*npredec+1];
304 for (Int_t i=0; i<ndec; i++)
306 for (Int_t j=0; j<npredec; j++)
308 tabx[npredec*i+j]=TMath::Power(10,((Double_t)i)+((Double_t)startvalue)+((Double_t)j)/((Double_t)npredec));
311 tabx[ndec*npredec]=TMath::Power(10,ndec+startvalue);
314 Double_t binsPtDummy[kPtBins+1];
316 for(int i=1; i<=kPtBins+1; i++)
318 if(binsPtDummy[i-1]+0.05<1.01)
319 binsPtDummy[i]=binsPtDummy[i-1]+0.05;
321 binsPtDummy[i]=binsPtDummy[i-1]+0.1;
325 Double_t binsPDummy[kPBins+1];
327 for(int i=1; i<=kPBins+1; i++)
329 if(binsPDummy[i-1]+0.05<1.01)
330 binsPDummy[i]=binsPDummy[i-1]+0.05;
332 binsPDummy[i]=binsPDummy[i-1]+0.1;
335 fTPCsignal=new TH2F("fTPCsignal",";p [GeV/c];dEdx",kPBins,binsPDummy,500,0,500);
336 fOutputList->Add(fTPCsignal);
337 fTPCdedxAfterTOFpidPions=new TH2F("fTPCsignalAfterTOFpions",";p [GeV/c];dEdx",kPBins,binsPDummy,500,0,500);
338 fOutputList->Add(fTPCdedxAfterTOFpidPions);
339 fTPCdedxAfterTOFpidKaons=new TH2F("fTPCsignalAfterTOFkaons",";p [GeV/c];dEdx",kPBins,binsPDummy,500,0,500);
340 fOutputList->Add(fTPCdedxAfterTOFpidKaons);
341 fTPCdedxAfterTOFpidProtons=new TH2F("fTPCsignalAfterTOFprotons",";p [GeV/c];dEdx",kPBins,binsPDummy,500,0,500);
342 fOutputList->Add(fTPCdedxAfterTOFpidProtons);
343 fTPCsignalAfterPionCuts=new TH2F("fTPCsignalAfterPionCuts",";p [GeV/c];dE/dx",kPBins,binsPDummy,500, 0., 500.);//
344 fTPCsignalAfterKaonCuts=new TH2F("fTPCsignalAfterKaonCuts",";p [GeV/c];dE/dx",kPBins,binsPDummy,500, 0., 500.);//
345 fTPCsignalAfterProtonCuts=new TH2F("fTPCsignalAfterProtonCuts",";p [GeV/c];dE/dx",kPBins,binsPDummy,500, 0., 500.);//
346 fOutputList->Add(fTPCsignalAfterPionCuts);
347 fOutputList->Add(fTPCsignalAfterKaonCuts);
348 fOutputList->Add(fTPCsignalAfterProtonCuts);
349 fTPCsignalAfterPionCuts1=new TH2F("fTPCsignalAfterPionCuts1",";p [GeV/c];dE/dx",kPBins,binsPDummy,500, 0., 500.);//
350 fTPCsignalAfterKaonCuts1=new TH2F("fTPCsignalAfterKaonCuts1",";p [GeV/c];dE/dx",kPBins,binsPDummy,500, 0., 500.);//
351 fTPCsignalAfterProtonCuts1=new TH2F("fTPCsignalAfterProtonCuts1",";p [GeV/c];dE/dx",kPBins,binsPDummy,500, 0., 500.);//
352 fOutputList->Add(fTPCsignalAfterPionCuts1);
353 fOutputList->Add(fTPCsignalAfterKaonCuts1);
354 fOutputList->Add(fTPCsignalAfterProtonCuts1);
356 fTPCsignalPi=new TH2F("fTPCsignalPi",";p [GeV/c];signal",kPBins,binsPDummy,300,-2,2);//TPC PID signal as function of p for pi+
357 fTPCsignalK=new TH2F("fTPCsignalK",";p [GeV/c];signal",kPBins,binsPDummy,300,-2,2);//TPC PID signal as function of p for K+
358 fTPCsignalP=new TH2F("fTPCsignalP",";p [GeV/c];signal",kPBins,binsPDummy,300,-2,2);//TPC PID signal as function of p for p
359 fOutputList->Add(fTPCsignalPi);
360 fOutputList->Add(fTPCsignalK);
361 fOutputList->Add(fTPCsignalP);
362 fTPCsignalPiafter=new TH2F("fTPCsignalPiafter",";p [GeV/c];(dE/dx-dE/dx_{#pi})/dE/dx_{#pi}",kPBins,binsPDummy,300, -2., 2.);//
363 fTPCsignalKafter=new TH2F("fTPCsignalKafter",";p [GeV/c];(dE/dx-dE/dx_{K})/dE/dx_{K}",kPBins,binsPDummy,300, -2., 2.);//
364 fTPCsignalPafter=new TH2F("fTPCsignalPafter",";p [GeV/c];(dE/dx-dE/dx_{p})/dE/dx_{p}",kPBins,binsPDummy,300, -2., 2.);//
365 fOutputList->Add(fTPCsignalPiafter);
366 fOutputList->Add(fTPCsignalKafter);
367 fOutputList->Add(fTPCsignalPafter);
371 fTPCsignalPimc=new TH2F("fTPCsignalPionsMC",";p [GeV/c];signal",kPBins,binsPDummy,600,-2,2);//TPC PID signal as function of pt for pi+
372 fTPCsignalKmc=new TH2F("fTPCsignalKaonsMC",";p [GeV/c];signal",kPBins,binsPDummy,600,-2,2);//TPC PID signal as function of pt for K+
373 fTPCsignalPmc=new TH2F("fTPCsignalProtonsMC",";p [GeV/c];signal",kPBins,binsPDummy,600,-2,2);//TPC PID signal as function of pt for p
374 fOutputList->Add(fTPCsignalPimc);
375 fOutputList->Add(fTPCsignalKmc);
376 fOutputList->Add(fTPCsignalPmc);
379 fTOFtime=new TH2F("fTOFtime",";p[GeV/c];#time",kPBins,binsPDummy,1000, 0.4, 1.1);//
380 fOutputList->Add(fTOFtime);
381 fTOFtimeE=new TH2F("fTOFtimeE",";p [GeV/c];#time-#time_{#pi}",kPBins,binsPDummy,500, -0.25, 0.25);//
382 fTOFtimePi=new TH2F("fTOFtimePi",";p [GeV/c];#time-#time_{#pi}",kPBins,binsPDummy,500, -0.25, 0.25);//
383 fTOFtimeK=new TH2F("fTOFtimeK",";p [GeV/c];#time-#time_{K}",kPBins,binsPDummy,500, -0.25, 0.25);//
384 fTOFtimeP=new TH2F("fTOFtimeP",";p [GeV/c];#time-#time_{p}",kPBins,binsPDummy,500, -0.25, 0.25);//
385 //fOutputList->Add(fTOFtimeE);
386 //fOutputList->Add(fTOFtimePi);
387 //fOutputList->Add(fTOFtimeK);
388 //fOutputList->Add(fTOFtimeP);
390 fTOFbeta=new TH2F("fTOFbeta",";p[GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
391 fOutputList->Add(fTOFbeta);
392 fTOFbetaE=new TH2F("fTOFbetaE",";p [GeV/c];#beta-#beta_{#pi}",kPBins,binsPDummy,500, -0.25, 0.25);//
393 fTOFbetaPi=new TH2F("fTOFbetaPi",";p [GeV/c];#beta-#beta_{#pi}",kPBins,binsPDummy,500, -0.25, 0.25);//
394 fTOFbetaK=new TH2F("fTOFbetaK",";p [GeV/c];#beta-#beta_{K}",kPBins,binsPDummy,500, -0.25, 0.25);//
395 fTOFbetaP=new TH2F("fTOFbetaP",";p [GeV/c];#beta-#beta_{p}",kPBins,binsPDummy,500, -0.25, 0.25);//
396 fOutputList->Add(fTOFbetaE);
397 fOutputList->Add(fTOFbetaPi);
398 fOutputList->Add(fTOFbetaK);
399 fOutputList->Add(fTOFbetaP);
401 fTOFinvbeta=new TH2F("fTOFinvbeta",";p[GeV/c];1/#beta",kPBins,binsPDummy,1000, 0.90, 2.5);//
402 fOutputList->Add(fTOFinvbeta);
403 fTOFinvbetaE=new TH2F("fTOFinvbetaE",";p [GeV/c];1/#beta-1/#beta_{#pi}",kPBins,binsPDummy,600, -0.3, 0.3);//
404 fTOFinvbetaPi=new TH2F("fTOFinvbetaPi",";p [GeV/c];1/#beta-1/#beta_{#pi}",kPBins,binsPDummy,600, -0.3, 0.3);//
405 fTOFinvbetaK=new TH2F("fTOFinvbetaK",";p [GeV/c];1/#beta-1/#beta_{K}",kPBins,binsPDummy,600, -0.3, 0.3);//
406 fTOFinvbetaP=new TH2F("fTOFinvbetaP",";p [GeV/c];1/#beta-1/#beta_{p}",kPBins,binsPDummy,600, -0.3, 0.3);//
407 fOutputList->Add(fTOFinvbetaE);
408 fOutputList->Add(fTOFinvbetaPi);
409 fOutputList->Add(fTOFinvbetaK);
410 fOutputList->Add(fTOFinvbetaP);
412 fTOFbetaAfterElectronsCuts=new TH2F("fTOFbetaAfterElectronsCuts",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
413 fTOFbetaAfterPionCuts=new TH2F("fTOFbetaAfterPionCuts",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
414 fTOFbetaAfterKaonCuts=new TH2F("fTOFbetaAfterKaonCuts",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
415 fTOFbetaAfterProtonCuts=new TH2F("fTOFbetaAfterProtonCuts",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
416 fOutputList->Add(fTOFbetaAfterElectronsCuts);
417 fOutputList->Add(fTOFbetaAfterPionCuts);
418 fOutputList->Add(fTOFbetaAfterKaonCuts);
419 fOutputList->Add(fTOFbetaAfterProtonCuts);
420 fTOFbetaAfterElectronsCuts1=new TH2F("fTOFbetaAfterElectronsCuts1",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
421 fTOFbetaAfterPionCuts1=new TH2F("fTOFbetaAfterPionCuts1",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
422 fTOFbetaAfterKaonCuts1=new TH2F("fTOFbetaAfterKaonCuts1",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
423 fTOFbetaAfterProtonCuts1=new TH2F("fTOFbetaAfterProtonCuts1",";p [GeV/c];#beta",kPBins,binsPDummy,1000, 0.4, 1.1);//
424 fOutputList->Add(fTOFbetaAfterElectronsCuts1);
425 fOutputList->Add(fTOFbetaAfterPionCuts1);
426 fOutputList->Add(fTOFbetaAfterKaonCuts1);
427 fOutputList->Add(fTOFbetaAfterProtonCuts1);
429 fTOFbetaEafter=new TH2F("fTOFbetaEafter",";p [GeV/c];#beta-#beta_{#pi}",kPBins,binsPDummy,500, -0.25, 0.25 );//
430 fTOFbetaPiafter=new TH2F("fTOFbetaPiafter",";p [GeV/c];#beta-#beta_{#pi}",kPBins,binsPDummy,500, -0.25, 0.25 );//
431 fTOFbetaKafter=new TH2F("fTOFbetaKafter",";p [GeV/c];#beta-#beta_{K}",kPBins,binsPDummy,500, -0.25, 0.25 );//
432 fTOFbetaPafter=new TH2F("fTOFbetaPafter",";p [GeV/c];#beta-#beta_{p}",kPBins,binsPDummy,500, -0.25, 0.25 );//
433 fOutputList->Add(fTOFbetaEafter);
434 fOutputList->Add(fTOFbetaPiafter);
435 fOutputList->Add(fTOFbetaKafter);
436 fOutputList->Add(fTOFbetaPafter);
440 fTOFyieldSelEmcE = new TH1F("fTOFyieldSelEmcE",";p [Gev/c];",kPBins,binsPDummy);
441 fTOFyieldSelPimcE = new TH1F("fTOFyieldSelPimcE",";p [Gev/c];",kPBins,binsPDummy);
442 fTOFyieldSelKmcE = new TH1F("fTOFyieldSelKmcE",";p [Gev/c];",kPBins,binsPDummy);
443 fTOFyieldSelPmcE = new TH1F("fTOFyieldSelPmcE",";p [Gev/c];",kPBins,binsPDummy);
444 fOutputList->Add(fTOFyieldSelEmcE);
445 fOutputList->Add(fTOFyieldSelPimcE);
446 fOutputList->Add(fTOFyieldSelKmcE);
447 fOutputList->Add(fTOFyieldSelPmcE);
448 fTOFyieldSelEmcM = new TH1F("fTOFyieldSelEmcM",";p [Gev/c];",kPBins,binsPDummy);
449 fTOFyieldSelPimcM = new TH1F("fTOFyieldSelPimcM",";p [Gev/c];",kPBins,binsPDummy);
450 fTOFyieldSelKmcM = new TH1F("fTOFyieldSelKmcM",";p [Gev/c];",kPBins,binsPDummy);
451 fTOFyieldSelPmcM = new TH1F("fTOFyieldSelPmcM",";p [Gev/c];",kPBins,binsPDummy);
452 fOutputList->Add(fTOFyieldSelEmcM);
453 fOutputList->Add(fTOFyieldSelPimcM);
454 fOutputList->Add(fTOFyieldSelKmcM);
455 fOutputList->Add(fTOFyieldSelPmcM);
456 fTOFyieldSelEmcPi = new TH1F("fTOFyieldSelEmcPi",";p [Gev/c];",kPBins,binsPDummy);
457 fTOFyieldSelPimcPi = new TH1F("fTOFyieldSelPimcPi",";p [Gev/c];",kPBins,binsPDummy);
458 fTOFyieldSelKmcPi = new TH1F("fTOFyieldSelKmcPi",";p [Gev/c];",kPBins,binsPDummy);
459 fTOFyieldSelPmcPi = new TH1F("fTOFyieldSelPmcPi",";p [Gev/c];",kPBins,binsPDummy);
460 fOutputList->Add(fTOFyieldSelEmcPi);
461 fOutputList->Add(fTOFyieldSelPimcPi);
462 fOutputList->Add(fTOFyieldSelKmcPi);
463 fOutputList->Add(fTOFyieldSelPmcPi);
464 fTOFyieldSelEmcK = new TH1F("fTOFyieldSelEmcK",";p [Gev/c];",kPBins,binsPDummy);
465 fTOFyieldSelPimcK = new TH1F("fTOFyieldSelPimcK",";p [Gev/c];",kPBins,binsPDummy);
466 fTOFyieldSelKmcK = new TH1F("fTOFyieldSelKmcK",";p [Gev/c];",kPBins,binsPDummy);
467 fTOFyieldSelPmcK = new TH1F("fTOFyieldSelPmcK",";p [Gev/c];",kPBins,binsPDummy);
468 fOutputList->Add(fTOFyieldSelEmcK);
469 fOutputList->Add(fTOFyieldSelPimcK);
470 fOutputList->Add(fTOFyieldSelKmcK);
471 fOutputList->Add(fTOFyieldSelPmcK);
472 fTOFyieldSelEmcP = new TH1F("fTOFyieldSelEmcP",";p [Gev/c];",kPBins,binsPDummy);
473 fTOFyieldSelPimcP = new TH1F("fTOFyieldSelPimcP",";p [Gev/c];",kPBins,binsPDummy);
474 fTOFyieldSelKmcP = new TH1F("fTOFyieldSelKmcP",";p [Gev/c];",kPBins,binsPDummy);
475 fTOFyieldSelPmcP = new TH1F("fTOFyieldSelPmcP",";p [Gev/c];",kPBins,binsPDummy);
476 fOutputList->Add(fTOFyieldSelEmcP);
477 fOutputList->Add(fTOFyieldSelPimcP);
478 fOutputList->Add(fTOFyieldSelKmcP);
479 fOutputList->Add(fTOFyieldSelPmcP);
480 fTPCyieldSelEmcE = new TH1F("fTPCyieldSelEmcE",";p [Gev/c];",kPBins,binsPDummy);
481 fTPCyieldSelPimcE = new TH1F("fTPCyieldSelPimcE",";p [Gev/c];",kPBins,binsPDummy);
482 fTPCyieldSelKmcE = new TH1F("fTPCyieldSelKmcE",";p [Gev/c];",kPBins,binsPDummy);
483 fTPCyieldSelPmcE = new TH1F("fTPCyieldSelPmcE",";p [Gev/c];",kPBins,binsPDummy);
484 fOutputList->Add(fTPCyieldSelEmcE);
485 fOutputList->Add(fTPCyieldSelPimcE);
486 fOutputList->Add(fTPCyieldSelKmcE);
487 fOutputList->Add(fTPCyieldSelPmcE);
488 fTPCyieldSelEmcM = new TH1F("fTPCyieldSelEmcM",";p [Gev/c];",kPBins,binsPDummy);
489 fTPCyieldSelPimcM = new TH1F("fTPCyieldSelPimcM",";p [Gev/c];",kPBins,binsPDummy);
490 fTPCyieldSelKmcM = new TH1F("fTPCyieldSelKmcM",";p [Gev/c];",kPBins,binsPDummy);
491 fTPCyieldSelPmcM = new TH1F("fTPCyieldSelPmcM",";p [Gev/c];",kPBins,binsPDummy);
492 fOutputList->Add(fTPCyieldSelEmcM);
493 fOutputList->Add(fTPCyieldSelPimcM);
494 fOutputList->Add(fTPCyieldSelKmcM);
495 fOutputList->Add(fTPCyieldSelPmcM);
496 fTPCyieldSelEmcPi = new TH1F("fTPCyieldSelEmcPi",";p [Gev/c];",kPBins,binsPDummy);
497 fTPCyieldSelPimcPi = new TH1F("fTPCyieldSelPimcPi",";p [Gev/c];",kPBins,binsPDummy);
498 fTPCyieldSelKmcPi = new TH1F("fTPCyieldSelKmcPi",";p [Gev/c];",kPBins,binsPDummy);
499 fTPCyieldSelPmcPi = new TH1F("fTPCyieldSelPmcPi",";p [Gev/c];",kPBins,binsPDummy);
500 fOutputList->Add(fTPCyieldSelEmcPi);
501 fOutputList->Add(fTPCyieldSelPimcPi);
502 fOutputList->Add(fTPCyieldSelKmcPi);
503 fOutputList->Add(fTPCyieldSelPmcPi);
504 fTPCyieldSelEmcK = new TH1F("fTPCyieldSelEmcK",";p [Gev/c];",kPBins,binsPDummy);
505 fTPCyieldSelPimcK = new TH1F("fTPCyieldSelPimcK",";p [Gev/c];",kPBins,binsPDummy);
506 fTPCyieldSelKmcK = new TH1F("fTPCyieldSelKmcK",";p [Gev/c];",kPBins,binsPDummy);
507 fTPCyieldSelPmcK = new TH1F("fTPCyieldSelPmcK",";p [Gev/c];",kPBins,binsPDummy);
508 fOutputList->Add(fTPCyieldSelEmcK);
509 fOutputList->Add(fTPCyieldSelPimcK);
510 fOutputList->Add(fTPCyieldSelKmcK);
511 fOutputList->Add(fTPCyieldSelPmcK);
512 fTPCyieldSelEmcP = new TH1F("fTPCyieldSelEmcP",";p [Gev/c];",kPBins,binsPDummy);
513 fTPCyieldSelPimcP = new TH1F("fTPCyieldSelPimcP",";p [Gev/c];",kPBins,binsPDummy);
514 fTPCyieldSelKmcP = new TH1F("fTPCyieldSelKmcP",";p [Gev/c];",kPBins,binsPDummy);
515 fTPCyieldSelPmcP = new TH1F("fTPCyieldSelPmcP",";p [Gev/c];",kPBins,binsPDummy);
516 fOutputList->Add(fTPCyieldSelEmcP);
517 fOutputList->Add(fTPCyieldSelPimcP);
518 fOutputList->Add(fTPCyieldSelKmcP);
519 fOutputList->Add(fTPCyieldSelPmcP);
522 fPvsPt=new TH2F("fPvsPt","p vs p_{t};p [GeV/c];p_{t} [GeV/c]",kPBins,binsPDummy,kPtBins,binsPtDummy);
523 fOutputList->Add(fPvsPt);
525 fMeanPvsP = new TProfile("fMeanPvsP","Mean P vs P;p [Gev/c];<p> [GeV/c]",kPBins,binsPDummy);
526 fOutputList->Add(fMeanPvsP);
528 fTPCvsGlobalMult = new TH2F("fTPCvsGlobalMult","TPC only vs Global track multiplicity;global;TPC only",500,0,2500,500,0,3500);
529 fOutputList->Add(fTPCvsGlobalMult);
531 fStandardGlobalCuts = AliFlowTrackCuts::GetStandardGlobalTrackCuts2010();
532 fStandardTPCCuts = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
534 fCutsTOFElectrons = new AliFlowTrackCuts("tof electron cuts");
535 fCutsTOFElectrons->SetPID(AliPID::kElectron, AliFlowTrackCuts::kTOFbeta);
536 fCutsTOFElectrons->SetQA();
537 fCutsTOFPions = new AliFlowTrackCuts("tof pion cuts");
538 fCutsTOFPions->SetPID(AliPID::kPion, AliFlowTrackCuts::kTOFbeta);
539 fCutsTOFPions->SetQA();
540 fCutsTOFKaons = new AliFlowTrackCuts("tof kaon cuts");
541 fCutsTOFKaons->SetPID(AliPID::kKaon, AliFlowTrackCuts::kTOFbeta);
542 fCutsTOFKaons->SetQA();
543 fCutsTOFProtons = new AliFlowTrackCuts("tof proton cuts");
544 fCutsTOFProtons->SetPID(AliPID::kProton, AliFlowTrackCuts::kTOFbeta);
545 fCutsTOFProtons->SetQA();
546 fCutsTPCElectrons = new AliFlowTrackCuts("tpc electron cuts");
547 fCutsTPCElectrons->SetPID(AliPID::kElectron, AliFlowTrackCuts::kTPCdedx);
548 fCutsTPCElectrons->SetQA();
549 fCutsTPCPions = new AliFlowTrackCuts("tpc pion cuts");
550 fCutsTPCPions->SetPID(AliPID::kPion, AliFlowTrackCuts::kTPCdedx);
551 fCutsTPCPions->SetQA();
552 fCutsTPCKaons = new AliFlowTrackCuts("tpc kaon cuts");
553 fCutsTPCKaons->SetPID(AliPID::kKaon, AliFlowTrackCuts::kTPCdedx);
554 fCutsTPCKaons->SetQA();
555 fCutsTPCProtons = new AliFlowTrackCuts("tpc proton cuts");
556 fCutsTPCProtons->SetPID(AliPID::kProton, AliFlowTrackCuts::kTPCdedx);
557 fCutsTPCProtons->SetQA();
559 //fOutputList->Add(fESDpid);
562 h = static_cast<TH1*>(fCutsTOFPions->GetQA()->At(0));
563 h->SetName(Form("pion direct %s",h->GetName()));
564 h = static_cast<TH1*>(fCutsTOFKaons->GetQA()->At(0));
565 h->SetName(Form("kaon direct %s",h->GetName()));
566 h = static_cast<TH1*>(fCutsTOFProtons->GetQA()->At(0));
567 h->SetName(Form("proton direct %s",h->GetName()));
568 fOutputList->Add(fCutsTOFPions->GetQA()->At(0));
569 fOutputList->Add(fCutsTOFKaons->GetQA()->At(0));
570 fOutputList->Add(fCutsTOFProtons->GetQA()->At(0));
572 h = static_cast<TH1*>(fCutsTOFPions->GetQA()->At(1));
573 h->SetName(Form("pion direct %s",h->GetName()));
574 h = static_cast<TH1*>(fCutsTOFKaons->GetQA()->At(1));
575 h->SetName(Form("kaon direct %s",h->GetName()));
576 h = static_cast<TH1*>(fCutsTOFProtons->GetQA()->At(1));
577 h->SetName(Form("proton direct %s",h->GetName()));
578 fOutputList->Add(fCutsTOFPions->GetQA()->At(1));
579 fOutputList->Add(fCutsTOFKaons->GetQA()->At(1));
580 fOutputList->Add(fCutsTOFProtons->GetQA()->At(1));
582 h = static_cast<TH1*>(fCutsTOFPions->GetQA()->At(2));
583 h->SetName(Form("pion direct %s",h->GetName()));
584 h = static_cast<TH1*>(fCutsTOFKaons->GetQA()->At(2));
585 h->SetName(Form("kaon direct %s",h->GetName()));
586 h = static_cast<TH1*>(fCutsTOFProtons->GetQA()->At(2));
587 h->SetName(Form("proton direct %s",h->GetName()));
588 fOutputList->Add(fCutsTOFPions->GetQA()->At(2));
589 fOutputList->Add(fCutsTOFKaons->GetQA()->At(2));
590 fOutputList->Add(fCutsTOFProtons->GetQA()->At(2));
592 h = static_cast<TH1*>(fCutsTOFPions->GetQA()->At(3));
593 h->SetName(Form("pion direct %s",h->GetName()));
594 h = static_cast<TH1*>(fCutsTOFKaons->GetQA()->At(3));
595 h->SetName(Form("kaon direct %s",h->GetName()));
596 h = static_cast<TH1*>(fCutsTOFProtons->GetQA()->At(3));
597 h->SetName(Form("proton direct %s",h->GetName()));
598 fOutputList->Add(fCutsTOFPions->GetQA()->At(3));
599 fOutputList->Add(fCutsTOFKaons->GetQA()->At(3));
600 fOutputList->Add(fCutsTOFProtons->GetQA()->At(3));
602 if (fUseDebugFile) fFile = fopen("debug.txt","w");
604 PostData(1, fOutputList);
607 //________________________________________________________________________
608 void AliAnalysisTaskPIDflowQA::UserExec(Option_t *)
610 fESD = dynamic_cast<AliESDEvent*> (InputEvent());
613 //do the calibration bit
614 fESDpid->SetTOFResponse(fESD,AliESDpid::kTOF_T0); // to use T0-TOF
615 fESDpid->MakePID(fESD,kFALSE);
617 if(!fCuts || !fEventCuts)
619 Printf("No CUTS Defined.........\n");
620 PostData(1, fOutputList);
624 if (!(fEventCuts->IsSelected(fESD)))
630 AliMCEvent* mcEvent=NULL;
633 mcEvent = (AliMCEvent*) MCEvent();
634 Printf("MC particles: %d", mcEvent->GetNumberOfTracks());
635 stack = mcEvent->Stack();
638 Printf("There are %d tracks in this event", fESD->GetNumberOfTracks());
639 Int_t nTracks=fESD->GetNumberOfTracks();
641 AliESDtrack *trackESD=0;
643 for(int tr1=0; tr1<nTracks; tr1++)
645 trackESD=fESD->GetTrack(tr1);
646 if (!trackESD) continue;
648 Double_t p=trackESD->GetP();
649 Double_t pt=trackESD->Pt();
651 if(!(fCuts->IsSelected(trackESD))) continue;
654 if(fMC) label=trackESD->GetLabel();
659 TParticle* particle2 = stack->Particle(TMath::Abs(label));
660 pdgcode=particle2->GetPdgCode();
664 fMeanPvsP->Fill(p,p);
666 pidTPC(trackESD,pdgcode);
667 pidTOF(trackESD,pdgcode);
670 //check the correlation between the global and TPConly number of tracks
671 fStandardGlobalCuts->SetEvent(fESD);
672 fStandardTPCCuts->SetEvent(fESD);
673 Int_t multGlobal = fStandardGlobalCuts->Count();
674 Int_t multTPC = fStandardTPCCuts->Count();
675 fTPCvsGlobalMult->Fill(multGlobal,multTPC);
679 const AliESDVertex* pvtx = fESD->GetPrimaryVertex();
680 const AliESDVertex* tpcvtx = fESD->GetPrimaryVertexTPC();
681 const AliESDVertex* spdvtx = fESD->GetPrimaryVertexSPD();
682 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
683 AliVEventHandler* handler = mgr->GetInputEventHandler();
684 TTree* tree = handler->GetTree();
685 TFile* file = tree->GetCurrentFile();
686 if (multTPC>(23+1.216*multGlobal) || multTPC<(-20+1.087*multGlobal))
688 fprintf(fFile, "%i %i %s %i\n",multTPC,multGlobal,file->GetName(),fESD->GetEventNumberInFile());
689 fprintf(fFile, " primary vertex: x: %.2f, y: %.2f, z: %.2f, n: %i\n", pvtx->GetX(), pvtx->GetY(), pvtx->GetZ(), pvtx->GetNContributors());
690 fprintf(fFile, " SPD vertex: x: %.2f, y: %.2f, z: %.2f, n: %i\n", spdvtx->GetX(), spdvtx->GetY(), spdvtx->GetZ(), spdvtx->GetNContributors());
691 fprintf(fFile, " TPC vertex: x: %.2f, y: %.2f, z: %.2f, n: %i\n", tpcvtx->GetX(), tpcvtx->GetY(), tpcvtx->GetZ(), tpcvtx->GetNContributors());
696 //________________________________________________________________________
697 void AliAnalysisTaskPIDflowQA::Terminate(Option_t *)
705 Printf("AliAnalysisTaskPIDflowQA: end of Terminate");
709 //________________________________________________________________________
710 void AliAnalysisTaskPIDflowQA::pidTPC(AliESDtrack* t, Int_t pdgcode)
713 const AliExternalTrackParam* innerParam = t->GetInnerParam();
714 if (!innerParam) return;
715 Double_t pinTPCglobal=innerParam->GetP();
716 Double_t tpcSignal =t ->GetTPCsignal();
717 Float_t p=innerParam->P();
718 Float_t pt=innerParam->Pt();
719 Float_t sigPion = fESDpid->GetTPCResponse().GetExpectedSignal(pinTPCglobal, AliPID::kPion);
720 Float_t sigKaon = fESDpid->GetTPCResponse().GetExpectedSignal(pinTPCglobal, AliPID::kKaon);
721 Float_t sigProton = fESDpid->GetTPCResponse().GetExpectedSignal(pinTPCglobal, AliPID::kProton);
722 if(!(sigPion>0.0&&sigKaon>0.0&&sigProton>0.0))
725 fTPCsignal->Fill(pinTPCglobal,tpcSignal);
727 fTPCsignalPi->Fill(p,(tpcSignal-sigPion)/sigPion);
728 fTPCsignalK->Fill(p,(tpcSignal-sigKaon)/sigKaon);
729 fTPCsignalP->Fill(p,(tpcSignal-sigProton)/sigProton);
731 if (fCutsTPCPions->IsSelected(t))
733 fTPCsignalAfterPionCuts->Fill(p,tpcSignal);
734 fTPCsignalPiafter->Fill(p,(tpcSignal-sigPion)/sigPion);
737 if (TMath::Abs(pdgcode)==11)
739 fTPCyieldSelPimcE->Fill(pt);
741 if (TMath::Abs(pdgcode)==13)
743 fTPCyieldSelPimcM->Fill(pt);
745 if (TMath::Abs(pdgcode)==211)
747 fTPCyieldSelPimcPi->Fill(pt);
749 else if(TMath::Abs(pdgcode)==321)
751 fTPCyieldSelPimcK->Fill(pt);
753 else if (TMath::Abs(pdgcode)==2212)
755 fTPCyieldSelPimcP->Fill(pt);
759 if (fCutsTPCKaons->IsSelected(t))
761 fTPCsignalAfterKaonCuts->Fill(p,tpcSignal);
762 fTPCsignalKafter->Fill(p,(tpcSignal-sigKaon)/sigKaon);
765 if (TMath::Abs(pdgcode)==11)
767 fTPCyieldSelKmcE->Fill(pt);
769 if (TMath::Abs(pdgcode)==13)
771 fTPCyieldSelKmcM->Fill(pt);
773 if (TMath::Abs(pdgcode)==211)
775 fTPCyieldSelKmcPi->Fill(pt);
777 else if(TMath::Abs(pdgcode)==321)
779 fTPCyieldSelKmcK->Fill(pt);
781 else if (TMath::Abs(pdgcode)==2212)
783 fTPCyieldSelKmcP->Fill(pt);
787 if (fCutsTPCProtons->IsSelected(t))
789 fTPCsignalAfterProtonCuts->Fill(p,tpcSignal);
790 fTPCsignalPafter->Fill(p,(tpcSignal-sigProton)/sigProton);
793 if (TMath::Abs(pdgcode)==11)
795 fTPCyieldSelPmcE->Fill(pt);
797 if (TMath::Abs(pdgcode)==13)
799 fTPCyieldSelPmcM->Fill(pt);
801 if (TMath::Abs(pdgcode)==211)
803 fTPCyieldSelPmcPi->Fill(pt);
805 else if(TMath::Abs(pdgcode)==321)
807 fTPCyieldSelPmcK->Fill(pt);
809 else if (TMath::Abs(pdgcode)==2212)
811 fTPCyieldSelPmcP->Fill(pt);
818 if(TMath::Abs(pdgcode)==211)
819 fTPCsignalPimc->Fill(p,(tpcSignal-sigPion)/sigPion);
820 else if(TMath::Abs(pdgcode)==321)
821 fTPCsignalKmc->Fill(p,(tpcSignal-sigKaon)/sigKaon);
822 else if (TMath::Abs(pdgcode)==2212)
823 fTPCsignalPmc->Fill(p,(tpcSignal-sigProton)/sigProton);
827 //______________________________________________________________________________
828 void AliAnalysisTaskPIDflowQA::pidTOF(AliESDtrack* track, Int_t pdgcode)
831 Bool_t goodtrack = (track) &&
832 (track->GetStatus() & AliESDtrack::kTOFpid) &&
833 (track->GetTOFsignal() > 12000) &&
834 (track->GetTOFsignal() < 100000) &&
835 (track->GetIntegratedLength() > 365) &&
836 !(track->GetStatus() & AliESDtrack::kTOFmismatch);
838 if (!goodtrack) return;
840 const AliExternalTrackParam* innerParam = track->GetInnerParam();
841 if (!innerParam) return;
842 Double_t pinTPCglobal=innerParam->GetP();
843 Double_t tpcSignal =track->GetTPCsignal();
844 const Float_t c = 2.99792457999999984e-02;
845 Float_t p = track->GetP();
846 Float_t pt = track->Pt();
847 Float_t l = track->GetIntegratedLength();
848 Float_t trackT0 = fESDpid->GetTOFResponse().GetStartTime(p);
849 Float_t timeTOF = track->GetTOFsignal()- trackT0;
850 Double_t integratedTimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
851 track->GetIntegratedTimes(integratedTimes);
854 Float_t beta = l/timeTOF/c;
855 Float_t betaHypothesis[5] = {0.0,0.0,0.0,0.0,0.0};
856 Float_t betadiff[5] = {0.0,0.0,0.0,0.0,0.0};
857 for (Int_t i=0;i<5;i++)
859 betaHypothesis[i] = l/integratedTimes[i]/c;
860 betadiff[i] = beta-betaHypothesis[i];
864 Float_t invbeta = 1/beta;
865 Float_t invbetaHypothesis[5] = {0.0,0.0,0.0,0.0,0.0};
866 Float_t invbetadiff[5] = {0.0,0.0,0.0,0.0,0.0};
867 for (Int_t i=0;i<5;i++)
869 invbetaHypothesis[i] = 1/betaHypothesis[i];
870 invbetadiff[i] = invbeta-invbetaHypothesis[i];
873 /////////////simple cuts
874 Bool_t isPion = ( (betadiff[2]<0.015) && (betadiff[2]>-0.015) &&
875 (betadiff[3]>0.025) &&
876 (betadiff[4]>0.03) );
878 Bool_t isKaon = ( (betadiff[3]<0.015) && (betadiff[3]>-0.015) &&
879 (betadiff[2]<-0.03) &&
880 (betadiff[4]>0.03) );
882 Bool_t isProton = ( (betadiff[4]<0.015) && (betadiff[4]>-0.015) &&
883 (betadiff[3]<-0.025) &&
884 (betadiff[2]<-0.025) );
886 if (isPion) fTOFbetaAfterPionCuts1->Fill(p,beta);
887 if (isKaon) fTOFbetaAfterKaonCuts1->Fill(p,beta);
888 if (isProton) fTOFbetaAfterProtonCuts1->Fill(p,beta);
891 fTOFbeta->Fill(p,beta);
892 fTOFbetaE->Fill(p,betadiff[0]);
893 fTOFbetaPi->Fill(p,betadiff[2]);
894 fTOFbetaK->Fill(p,betadiff[3]);
895 fTOFbetaP->Fill(p,betadiff[4]);
897 fTOFinvbeta->Fill(p,invbeta);
898 fTOFinvbetaE->Fill(p,invbetadiff[0]);
899 fTOFinvbetaPi->Fill(p,invbetadiff[2]);
900 fTOFinvbetaK->Fill(p,invbetadiff[3]);
901 fTOFinvbetaP->Fill(p,invbetadiff[4]);
903 if (fCutsTOFElectrons->PassesTOFbetaCut(track))
905 fTOFbetaAfterElectronsCuts->Fill(p,beta);
906 fTOFbetaEafter->Fill(p,beta-betaHypothesis[0]);
909 if (TMath::Abs(pdgcode)==11)
911 fTOFyieldSelEmcE->Fill(pt);
913 if (TMath::Abs(pdgcode)==13)
915 fTOFyieldSelEmcM->Fill(pt);
917 if (TMath::Abs(pdgcode)==211)
919 fTOFyieldSelEmcPi->Fill(pt);
921 else if(TMath::Abs(pdgcode)==321)
923 fTOFyieldSelEmcK->Fill(pt);
925 else if (TMath::Abs(pdgcode)==2212)
927 fTOFyieldSelEmcP->Fill(pt);
931 if (fCutsTOFPions->PassesTOFbetaCut(track))
933 fTPCdedxAfterTOFpidPions->Fill(pinTPCglobal,tpcSignal);
934 fTOFbetaAfterPionCuts->Fill(p,beta);
935 fTOFbetaPiafter->Fill(p,beta-betaHypothesis[2]);
938 if (TMath::Abs(pdgcode)==11)
940 fTOFyieldSelPimcE->Fill(pt);
942 if (TMath::Abs(pdgcode)==13)
944 fTOFyieldSelPimcM->Fill(pt);
946 if (TMath::Abs(pdgcode)==211)
948 fTOFyieldSelPimcPi->Fill(pt);
950 else if(TMath::Abs(pdgcode)==321)
952 fTOFyieldSelPimcK->Fill(pt);
954 else if (TMath::Abs(pdgcode)==2212)
956 fTOFyieldSelPimcP->Fill(pt);
960 if (fCutsTOFKaons->PassesTOFbetaCut(track))
962 fTPCdedxAfterTOFpidKaons->Fill(pinTPCglobal,tpcSignal);
963 fTOFbetaAfterKaonCuts->Fill(p,beta);
964 fTOFbetaKafter->Fill(p,beta-betaHypothesis[3]);
967 if (TMath::Abs(pdgcode)==11)
969 fTOFyieldSelKmcE->Fill(pt);
971 if (TMath::Abs(pdgcode)==13)
973 fTOFyieldSelKmcM->Fill(pt);
975 if (TMath::Abs(pdgcode)==211)
977 fTOFyieldSelKmcPi->Fill(pt);
979 else if(TMath::Abs(pdgcode)==321)
981 fTOFyieldSelKmcK->Fill(pt);
983 else if (TMath::Abs(pdgcode)==2212)
985 fTOFyieldSelKmcP->Fill(pt);
989 if (fCutsTOFProtons->PassesTOFbetaCut(track))
991 fTPCdedxAfterTOFpidProtons->Fill(pinTPCglobal,tpcSignal);
992 fTOFbetaAfterProtonCuts->Fill(p,beta);
993 fTOFbetaPafter->Fill(p,beta-betaHypothesis[4]);
996 if (TMath::Abs(pdgcode)==11)
998 fTOFyieldSelPmcE->Fill(pt);
1000 if (TMath::Abs(pdgcode)==13)
1002 fTOFyieldSelPmcM->Fill(pt);
1004 if (TMath::Abs(pdgcode)==211)
1006 fTOFyieldSelPmcPi->Fill(pt);
1008 else if(TMath::Abs(pdgcode)==321)
1010 fTOFyieldSelPmcK->Fill(pt);
1012 else if (TMath::Abs(pdgcode)==2212)
1014 fTOFyieldSelPmcP->Fill(pt);
1021 //______________________________________________________________________________
1022 Float_t AliAnalysisTaskPIDflowQA::Beta(Float_t m, Float_t p)
1024 //get theoretical beta
1025 return TMath::Sqrt(1. / (1. + m * m / (p * p)));