]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - ANALYSIS/AliAnalysisTaskPIDqa.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
adding cross check histograms
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskPIDqa.cxx
... / ...
CommitLineData
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/* $Id: AliAnalysisTaskPIDqa.cxx 43811 2010-09-23 14:13:31Z wiechula $ */
17#include <TList.h>
18#include <TVectorD.h>
19#include <TObjArray.h>
20#include <TH2.h>
21#include <TFile.h>
22#include <TPRegexp.h>
23#include <TChain.h>
24#include <TF1.h>
25#include <TSpline.h>
26
27#include <AliAnalysisManager.h>
28#include <AliInputEventHandler.h>
29#include <AliVEventHandler.h>
30#include <AliVEvent.h>
31#include <AliVParticle.h>
32#include <AliVTrack.h>
33#include <AliLog.h>
34#include <AliPID.h>
35#include <AliPIDResponse.h>
36#include <AliITSPIDResponse.h>
37#include <AliTPCPIDResponse.h>
38#include <AliTRDPIDResponse.h>
39#include <AliTOFPIDResponse.h>
40#include <AliTPCdEdxInfo.h>
41
42#include <AliESDEvent.h>
43#include <AliAODEvent.h>
44#include <AliESDv0.h>
45#include <AliAODv0.h>
46#include <AliESDv0KineCuts.h>
47#include <AliESDtrackCuts.h>
48
49#include <AliMCEvent.h>
50
51#include "AliAnalysisTaskPIDqa.h"
52
53
54ClassImp(AliAnalysisTaskPIDqa)
55
56//______________________________________________________________________________
57AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa():
58AliAnalysisTaskSE(),
59fPIDResponse(0x0),
60fV0cuts(0x0),
61fV0electrons(0x0),
62fV0pions(0x0),
63fV0kaons(0x0),
64fV0protons(0x0),
65fListQA(0x0),
66fListQAits(0x0),
67fListQAitsSA(0x0),
68fListQAitsPureSA(0x0),
69fListQAtpc(0x0),
70fListQAtpcBasic(0x0),
71fListQAtpcMCtruth(0x0),
72fListQAtpcHybrid(0x0),
73fListQAtpcOROChigh(0x0),
74fListQAtpcV0(0x0),
75fListQAtrd(0x0),
76fListQAtrdNsig(0x0),
77fListQAtrdNsigTPCTOF(0x0),
78fListQAtof(0x0),
79fListQAt0(0x0),
80fListQAemcal(0x0),
81fListQAhmpid(0x0),
82fListQAtofhmpid(0x0),
83fListQAtpctof(0x0),
84fListQAV0(0x0),
85fListQAinfo(0x0)
86{
87 //
88 // Dummy constructor
89 //
90}
91
92//______________________________________________________________________________
93AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa(const char* name):
94AliAnalysisTaskSE(name),
95fPIDResponse(0x0),
96fV0cuts(0x0),
97fV0electrons(0x0),
98fV0pions(0x0),
99fV0kaons(0x0),
100fV0protons(0x0),
101fListQA(0x0),
102fListQAits(0x0),
103fListQAitsSA(0x0),
104fListQAitsPureSA(0x0),
105fListQAtpc(0x0),
106fListQAtpcBasic(0x0),
107fListQAtpcMCtruth(0x0),
108fListQAtpcHybrid(0x0),
109fListQAtpcOROChigh(0x0),
110fListQAtpcV0(0x0),
111fListQAtrd(0x0),
112fListQAtrdNsig(0x0),
113fListQAtrdNsigTPCTOF(0x0),
114fListQAtof(0x0),
115fListQAt0(0x0),
116fListQAemcal(0x0),
117fListQAhmpid(0x0),
118fListQAtofhmpid(0x0),
119fListQAtpctof(0x0),
120fListQAV0(0x0),
121fListQAinfo(0x0)
122{
123 //
124 // Default constructor
125 //
126 DefineInput(0,TChain::Class());
127 DefineOutput(1,TList::Class());
128}
129
130//______________________________________________________________________________
131AliAnalysisTaskPIDqa::~AliAnalysisTaskPIDqa()
132{
133 //
134 // Destructor
135 //
136
137 delete fV0cuts;
138 delete fV0electrons;
139 delete fV0pions;
140 delete fV0kaons;
141 delete fV0protons;
142
143 if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fListQA;
144}
145
146//______________________________________________________________________________
147void AliAnalysisTaskPIDqa::UserCreateOutputObjects()
148{
149 //
150 // Create the output QA objects
151 //
152
153 AliLog::SetClassDebugLevel("AliAnalysisTaskPIDqa",10);
154
155 //input hander
156 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
157 AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
158 if (!inputHandler) AliFatal("Input handler needed");
159
160 //pid response object
161 fPIDResponse=inputHandler->GetPIDResponse();
162 if (!fPIDResponse) AliError("PIDResponse object was not created");
163
164 // V0 Kine cuts
165 fV0cuts = new AliESDv0KineCuts;
166
167 // V0 PID Obj arrays
168 fV0electrons = new TObjArray;
169 fV0pions = new TObjArray;
170 fV0kaons = new TObjArray;
171 fV0protons = new TObjArray;
172
173 //
174 fListQA=new TList;
175 fListQA->SetOwner();
176
177 fListQAits=new TList;
178 fListQAits->SetOwner();
179 fListQAits->SetName("ITS");
180
181 fListQAitsSA=new TList;
182 fListQAitsSA->SetOwner();
183 fListQAitsSA->SetName("ITS_SA");
184
185 fListQAitsPureSA=new TList;
186 fListQAitsPureSA->SetOwner();
187 fListQAitsPureSA->SetName("ITS_PureSA");
188
189 fListQAtpc=new TList;
190 fListQAtpc->SetOwner();
191 fListQAtpc->SetName("TPC");
192
193 fListQAtrd=new TList;
194 fListQAtrd->SetOwner();
195 fListQAtrd->SetName("TRD");
196
197 fListQAtrdNsig=new TList;
198 fListQAtrdNsig->SetOwner();
199 fListQAtrdNsig->SetName("TRDnSigma");
200
201 fListQAtrdNsigTPCTOF=new TList;
202 fListQAtrdNsigTPCTOF->SetOwner();
203 fListQAtrdNsigTPCTOF->SetName("TRDnSigma_TPCTOF");
204
205 fListQAtof=new TList;
206 fListQAtof->SetOwner();
207 fListQAtof->SetName("TOF");
208
209 fListQAt0=new TList;
210 fListQAt0->SetOwner();
211 fListQAt0->SetName("T0");
212
213 fListQAemcal=new TList;
214 fListQAemcal->SetOwner();
215 fListQAemcal->SetName("EMCAL");
216
217 fListQAhmpid=new TList;
218 fListQAhmpid->SetOwner();
219 fListQAhmpid->SetName("HMPID");
220
221 fListQAtpctof=new TList;
222 fListQAtpctof->SetOwner();
223 fListQAtpctof->SetName("TPC_TOF");
224
225 fListQAtofhmpid=new TList;
226 fListQAtofhmpid->SetOwner();
227 fListQAtofhmpid->SetName("TOF_HMPID");
228
229 fListQAV0=new TList;
230 fListQAV0->SetOwner();
231 fListQAV0->SetName("V0decay");
232
233 fListQAinfo=new TList;
234 fListQAinfo->SetOwner();
235 fListQAinfo->SetName("QAinfo");
236
237 fListQA->Add(fListQAits);
238 fListQA->Add(fListQAitsSA);
239 fListQA->Add(fListQAitsPureSA);
240 fListQA->Add(fListQAtpc);
241 fListQA->Add(fListQAtrd);
242 fListQA->Add(fListQAtof);
243 fListQA->Add(fListQAt0);
244 fListQA->Add(fListQAemcal);
245 fListQA->Add(fListQAhmpid);
246 fListQA->Add(fListQAtpctof);
247 fListQA->Add(fListQAtofhmpid);
248 fListQA->Add(fListQAV0);
249 fListQA->Add(fListQAinfo);
250
251 SetupITSqa();
252// SetupTPCqa(kFALSE, kTRUE, kFALSE);
253 SetupTRDqa();
254 SetupTOFqa();
255 SetupT0qa();
256 SetupEMCALqa();
257 SetupHMPIDqa();
258 SetupTPCTOFqa();
259 SetupTOFHMPIDqa();
260 SetupV0qa();
261 SetupQAinfo();
262
263 PostData(1,fListQA);
264}
265
266
267//______________________________________________________________________________
268void AliAnalysisTaskPIDqa::UserExec(Option_t */*option*/)
269{
270 //
271 // Setup the PID response functions and fill the QA histograms
272 //
273
274 AliVEvent *event=InputEvent();
275 if (!event||!fPIDResponse) return;
276
277 // Start with the V0 task (only possible for ESDs?)
278 FillV0PIDlist();
279
280 FillITSqa();
281 FillTPCqa();
282 FillTRDqa();
283 FillTOFqa();
284 FillEMCALqa();
285 FillHMPIDqa();
286 FillT0qa();
287
288 //combined detector QA
289 FillTPCTOFqa();
290 FillTOFHMPIDqa();
291
292 // Clear the V0 PID arrays
293 ClearV0PIDlist();
294
295 //QA info
296 FillQAinfo();
297
298 PostData(1,fListQA);
299}
300
301//______________________________________________________________________________
302void AliAnalysisTaskPIDqa::FillV0PIDlist(){
303
304 //
305 // Fill the PID object arrays holding the pointers to identified particle tracks
306 //
307
308 // Dynamic cast to ESD events (DO NOTHING for AOD events)
309 AliESDEvent *event = dynamic_cast<AliESDEvent *>(InputEvent());
310 if ( !event ) return;
311
312 if(TString(event->GetBeamType())=="Pb-Pb" || TString(event->GetBeamType())=="A-A"){
313 fV0cuts->SetMode(AliESDv0KineCuts::kPurity,AliESDv0KineCuts::kPbPb);
314 }
315 else{
316 fV0cuts->SetMode(AliESDv0KineCuts::kPurity,AliESDv0KineCuts::kPP);
317 }
318
319 // V0 selection
320 // set event
321 fV0cuts->SetEvent(event);
322
323 // loop over V0 particles
324 for(Int_t iv0=0; iv0<event->GetNumberOfV0s();iv0++){
325
326 AliESDv0 *v0 = (AliESDv0 *) event->GetV0(iv0);
327
328 if(!v0) continue;
329 if(v0->GetOnFlyStatus()) continue;
330
331 // Get the particle selection
332 Bool_t foundV0 = kFALSE;
333 Int_t pdgV0, pdgP, pdgN;
334
335 foundV0 = fV0cuts->ProcessV0(v0, pdgV0, pdgP, pdgN);
336 if(!foundV0) continue;
337
338 Int_t iTrackP = v0->GetPindex(); // positive track
339 Int_t iTrackN = v0->GetNindex(); // negative track
340
341 // v0 Armenteros plot (QA)
342 Float_t armVar[2] = {0.0,0.0};
343 fV0cuts->Armenteros(v0, armVar);
344
345 TH2 *h=(TH2*)fListQAV0->At(0);
346 if (!h) continue;
347 h->Fill(armVar[0],armVar[1]);
348
349 // fill the Object arrays
350 // positive particles
351 if( pdgP == -11){
352 fV0electrons->Add((AliVTrack*)event->GetTrack(iTrackP));
353 }
354 else if( pdgP == 211){
355 fV0pions->Add((AliVTrack*)event->GetTrack(iTrackP));
356 }
357 else if( pdgP == 321){
358 fV0kaons->Add((AliVTrack*)event->GetTrack(iTrackP));
359 }
360 else if( pdgP == 2212){
361 fV0protons->Add((AliVTrack*)event->GetTrack(iTrackP));
362 }
363
364 // negative particles
365 if( pdgN == 11){
366 fV0electrons->Add((AliVTrack*)event->GetTrack(iTrackN));
367 }
368 else if( pdgN == -211){
369 fV0pions->Add((AliVTrack*)event->GetTrack(iTrackN));
370 }
371 else if( pdgN == -321){
372 fV0kaons->Add((AliVTrack*)event->GetTrack(iTrackN));
373 }
374 else if( pdgN == -2212){
375 fV0protons->Add((AliVTrack*)event->GetTrack(iTrackN));
376 }
377
378
379 }
380}
381//______________________________________________________________________________
382void AliAnalysisTaskPIDqa::ClearV0PIDlist(){
383
384 //
385 // Clear the PID object arrays
386 //
387
388 fV0electrons->Clear();
389 fV0pions->Clear();
390 fV0kaons->Clear();
391 fV0protons->Clear();
392
393}
394//______________________________________________________________________________
395void AliAnalysisTaskPIDqa::FillITSqa()
396{
397 //
398 // Fill PID qa histograms for the ITS
399 //
400
401 AliVEvent *event=InputEvent();
402
403 Int_t ntracks=event->GetNumberOfTracks();
404 for(Int_t itrack = 0; itrack < ntracks; itrack++){
405 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
406 ULong_t status=track->GetStatus();
407 // not that nice. status bits not in virtual interface
408 // ITS refit + ITS pid selection
409 if (!( ( (status & AliVTrack::kITSrefit)==AliVTrack::kITSrefit ) ||
410 ! ( (status & AliVTrack::kITSpid )==AliVTrack::kITSpid ) )) continue;
411 Double_t mom=track->P();
412
413 TList *theList = 0x0;
414 if(( (status & AliVTrack::kTPCin)==AliVTrack::kTPCin )){
415 //ITS+TPC tracks
416 theList=fListQAits;
417 }else{
418 if(!( (status & AliVTrack::kITSpureSA)==AliVTrack::kITSpureSA )){
419 //ITS Standalone tracks
420 theList=fListQAitsSA;
421 }else{
422 //ITS Pure Standalone tracks
423 theList=fListQAitsPureSA;
424 }
425 }
426
427
428 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
429 TH2 *h=(TH2*)theList->At(ispecie);
430 if (!h) continue;
431 Double_t nSigma=fPIDResponse->NumberOfSigmasITS(track, (AliPID::EParticleType)ispecie);
432 h->Fill(mom,nSigma);
433 }
434 TH2 *h=(TH2*)theList->At(AliPID::kSPECIESC);
435 if (h) {
436 Double_t sig=track->GetITSsignal();
437 h->Fill(mom,sig);
438 }
439 }
440}
441
442
443//______________________________________________________________________________
444void AliAnalysisTaskPIDqa::FillTPCHistogramsSignal(TList *sublist, Int_t scenario, AliVTrack *track, Int_t mult)
445{
446 //
447 // Fill PID qa histograms for the TPC: Fill the histograms for the TPC signal for different settings
448 //
449
450 AliMCEvent *eventMC=MCEvent(); // MC event for MC truth PID
451
452 Double_t mom=0.; // track momentum
453 Double_t eta=0.; // track eta
454 Double_t sig=0.; // TPC dE/dx signal
455 Double_t sigStd=0.; // TPC dE/dx signal (standard = all ROCs)
456 Double_t sigIROC=0.; // TPC dE/dx signal (IROC)
457 Double_t sigOROCmedium=0.; // TPC dE/dx signal (OROCmedium)
458 Double_t sigOROClong=0.; // TPC dE/dx signal (OROClong)
459 Double_t eleLineDist=0.; // difference between TPC signal and electron expectation
460 Int_t trackLabel=0; // label of the AliVTrack to identify the corresponding MCtrack
461 Int_t pdgCode=0; // pdgcode of MC track for MC truth scenario
462 Int_t pdgCodeAbs=0; // absolute value of pdgcode to get both particles and antiparticles
463 Int_t iSigMax=1; // number of TPC signals (std = 1, set automatically higher if available)
464 Int_t nSpecies=0; // number of particle species under study
465 Int_t count=0; // counter for the number of plot sets for all species (i.e. nsigma vs. p, eta and mult)
466
467 mom=track->GetTPCmomentum();
468 eta=track->Eta();
469 sigStd=track->GetTPCsignal();
470
471 eleLineDist=sigStd-fPIDResponse->GetTPCResponse().GetExpectedSignal(track,AliPID::kElectron);
472
473 // Get number of particle species (less for V0 candidates = scenarios 40-44)
474 if (scenario > 39) nSpecies=(Int_t)AliPID::kSPECIES;
475 else nSpecies=(Int_t)AliPID::kSPECIESC;
476
477 // Set number of plot sets for all species
478 // (i.e. only nsigma vs. p => count=1; also vs. eta and mult => count=3)
479 if ( scenario == 1 || scenario > 39) count=3;
480 else count=1;
481
482 // Get MC track ( --> can be deleted if TPC signal is NOT filled for scenario=1 (MC truth)
483 if (eventMC) {
484 trackLabel=TMath::Abs(track->GetLabel());
485 AliVTrack *mcTrack=(AliVTrack*)eventMC->GetTrack(trackLabel);
486 pdgCode=mcTrack->PdgCode();
487 pdgCodeAbs=TMath::Abs(pdgCode);
488 }
489
490 // Get TPC dE/dx info and different TPC signals (IROC, OROCmedium, OROClong)
491 AliTPCdEdxInfo* fTPCdEdxInfo = 0x0;
492 fTPCdEdxInfo = track->GetTPCdEdxInfo();
493
494 if (fTPCdEdxInfo) {
495 sigIROC=fTPCdEdxInfo->GetTPCsignalShortPad();
496 sigOROCmedium=fTPCdEdxInfo->GetTPCsignalMediumPad();
497 sigOROClong=fTPCdEdxInfo->GetTPCsignalLongPad();
498 iSigMax=4;
499
500 //printf("mom = %.3f sigStd = %.3f sigIROC = %.3f sigOROCmedium = %.3f sigOROClong = %.3f \n",mom,sigStd,sigIROC,sigOROCmedium,sigOROClong);
501 }
502
503
504 // TPC signal for all particles vs. momentum (standard, IROC, OROCmedium, OROClong)
505 TH2 *h1std=(TH2*)sublist->At(count*nSpecies+4);
506 if (h1std) {
507 h1std->Fill(mom,sigStd);
508 }
509
510 TH2 *h1iroc=(TH2*)sublist->At(count*nSpecies+5);
511 if ( h1iroc && sigIROC ) {
512 h1iroc->Fill(mom,sigIROC);
513 }
514
515 TH2 *h1orocm=(TH2*)sublist->At(count*nSpecies+6);
516 if (h1orocm && sigOROCmedium ) {
517 h1orocm->Fill(mom,sigOROCmedium);
518 }
519
520 TH2 *h1orocl=(TH2*)sublist->At(count*nSpecies+7);
521 if ( h1orocl && sigOROClong ) {
522 h1orocl->Fill(mom,sigOROClong);
523 }
524
525
526 // - Beginn: MIP pions: TPC signal vs. eta, TPC signal vs. mult -
527 if (mom>0.45 && mom<0.5 && sigStd>40 && sigStd<60) {
528
529 Bool_t isPionMC=kTRUE;
530
531 if (scenario == 1) {
532 if ( pdgCodeAbs != 211 && pdgCodeAbs != 111 ) isPionMC=kFALSE;
533 }
534
535 // MIP pions: TPC signal vs. eta (standard, IROC, OROCmedium, OROClong)
536 for (Int_t iSig=0; iSig<iSigMax; iSig++) {
537 if (iSig==0) sig=sigStd;
538 else if (iSig==1) sig=sigIROC;
539 else if (iSig==2) sig=sigOROCmedium;
540 else if (iSig==3) sig=sigOROClong;
541
542 TH2 *h2=(TH2*)sublist->At(count*nSpecies+8+iSig);
543 if ( h2 && isPionMC ) {
544 h2->Fill(eta,sig);
545 }
546 }
547
548 // MIP pions: TPC signal vs. mult (standard, IROC, OROCmedium, OROClong)
549 for (Int_t iSig=0; iSig<iSigMax; iSig++) {
550 if (iSig==0) sig=sigStd;
551 else if (iSig==1) sig=sigIROC;
552 else if (iSig==2) sig=sigOROCmedium;
553 else if (iSig==3) sig=sigOROClong;
554
555 TH2 *h3=(TH2*)sublist->At(count*nSpecies+12+iSig);
556 if ( h3 && isPionMC && mult > 0 ) {
557 h3->Fill(mult,sig);
558 }
559 }
560 } // - End: MIP pions -
561
562 // - Beginn: Electrons: TPC signal vs. eta, TPC signal vs. mult -
563 if (mom>0.32 && mom<0.38 && eleLineDist>-10. && eleLineDist<15.) {
564
565 Bool_t isElectronMC=kTRUE;
566
567 if (scenario == 1) {
568 if ( pdgCodeAbs != 11 ) isElectronMC=kFALSE;
569 }
570
571 // Electrons: TPC signal vs. eta (standard, IROC, OROCmedium, OROClong)
572 for (Int_t iSig=0; iSig<iSigMax; iSig++) {
573 if (iSig==0) sig=sigStd;
574 else if (iSig==1) sig=sigIROC;
575 else if (iSig==2) sig=sigOROCmedium;
576 else if (iSig==3) sig=sigOROClong;
577
578 TH2 *h4=(TH2*)sublist->At(count*nSpecies+16+iSig);
579 if ( h4 && isElectronMC ) {
580 h4->Fill(eta,sig);
581 }
582 }
583
584 // Electrons: TPC signal vs. mult (standard, IROC, OROCmedium, OROClong)
585 for (Int_t iSig=0; iSig<iSigMax; iSig++) {
586 if (iSig==0) sig=sigStd;
587 else if (iSig==1) sig=sigIROC;
588 else if (iSig==2) sig=sigOROCmedium;
589 else if (iSig==3) sig=sigOROClong;
590
591 TH2 *h5=(TH2*)sublist->At(count*nSpecies+20+iSig);
592 if ( h5 && isElectronMC && mult > 0 ) {
593 h5->Fill(mult,sig);
594 }
595 }
596 } // - End: Electrons -
597
598}
599
600//______________________________________________________________________________
601void AliAnalysisTaskPIDqa::FillTPCHistogramsNsigma(TList *sublist, Int_t scenario, AliVTrack *track, Int_t mult)
602{
603 //
604 // Fill PID qa histograms for the TPC: Fill the histograms for TPC Nsigma for different settings
605 //
606
607 AliMCEvent *eventMC=MCEvent(); // MC event for MC truth PID
608
609 Double_t mom=0.; // track momentum
610 Double_t eta=0.; // track eta
611 Double_t nSigma=0.; // number of sigmas wrt. expected signal
612 Double_t sig=0.; // TPC dE/dx signal
613 Double_t eleLineDist=0.; // difference between TPC signal and electron expectation
614 Int_t trackLabel=0; // label of the AliVTrack to identify the corresponding MCtrack
615 Int_t pdgCode=0; // pdgcode of MC track for MC truth scenario
616 Int_t pdgCodeAbs=0; // absolute value of pdgcode to get both particles and antiparticles
617 Int_t nSpecies=0; // number of particle species under study
618 Int_t count=0; // counter for the number of plot sets for all species (i.e. vs. p, eta and mult)
619
620 mom=track->GetTPCmomentum();
621 eta=track->Eta();
622 sig=track->GetTPCsignal();
623
624 eleLineDist=sig-fPIDResponse->GetTPCResponse().GetExpectedSignal(track,AliPID::kElectron);
625
626 // Get number of particle species (less for V0 candidates = scenarios 40-44)
627 if (scenario > 39) nSpecies=(Int_t)AliPID::kSPECIES;
628 else nSpecies=(Int_t)AliPID::kSPECIESC;
629
630 // Set number of plot sets for all species
631 // (i.e. only vs. p => count=1; also vs. eta and mult => count=3)
632 if ( scenario == 1 || scenario > 39 ) count=3;
633 else count=1;
634
635 // Get MC track
636 if (eventMC) {
637 trackLabel=TMath::Abs(track->GetLabel());
638 AliVTrack *mcTrack=(AliVTrack*)eventMC->GetTrack(trackLabel);
639 pdgCode=mcTrack->PdgCode();
640 pdgCodeAbs=TMath::Abs(pdgCode);
641 }
642
643
644 // - Beginn: Nsigma vs. p, vs. eta and vs. multiplicity for different particle species -
645 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
646
647 TH2 *h=(TH2*)sublist->At(ispecie);
648 if (!h) continue;
649
650 if (scenario == 1) {
651 if ( ispecie == 0 && pdgCodeAbs != 11 ) continue; // Electron
652 if ( ispecie == 1 && pdgCodeAbs != 13 ) continue; // Muon
653 if ( ispecie == 2 && pdgCodeAbs != 211 && pdgCodeAbs!=111 ) continue; // Pion
654 if ( ispecie == 3 && pdgCodeAbs != 321 && pdgCodeAbs!=311 ) continue; // Kaon
655 if ( ispecie == 4 && pdgCodeAbs != 2212 ) continue; // Proton
656 if ( ispecie == 5 && pdgCodeAbs != 1000010020 ) continue; // Deuteron
657 if ( ispecie == 6 && pdgCodeAbs != 1000010030 ) continue; // Triton
658 if ( ispecie == 7 && pdgCodeAbs != 1000020030 ) continue; // Helium-3
659 if ( ispecie == 8 && pdgCodeAbs != 1000020040 ) continue; // Alpha
660 }
661 else if (scenario > 39) {
662 if ( ispecie == 0 && scenario != 40 ) continue; // Electron
663 if ( ispecie == 1 ) continue; // Muon
664 if ( ispecie == 2 && scenario != 42 ) continue; // Pion
665 if ( ispecie == 3 && scenario != 43 ) continue; // Kaon
666 if ( ispecie == 4 && scenario != 44 ) continue; // Proton
667 }
668
669 if (scenario == 2) {
670 nSigma=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie, AliTPCPIDResponse::kdEdxHybrid);
671 }
672 else if (scenario == 3) {
673 nSigma=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie, AliTPCPIDResponse::kdEdxOROC);
674 }
675 else {
676 nSigma=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
677 }
678
679 h->Fill(mom,nSigma);
680
681 if (count == 3) {
682 TH2 *hEta=(TH2*)sublist->At(ispecie+nSpecies);
683 TH2 *hMult=(TH2*)sublist->At(ispecie+2*nSpecies);
684
685 if ( hEta ) hEta->Fill(eta,nSigma);
686 if ( hMult && mult > 0 ) hMult->Fill(mult,nSigma);
687 }
688 } // - End: different particle species -
689
690
691 // -- Beginn: Fill histograms for MIP pions and electrons (only for some scenarios) --
692 if ( scenario == 0 || scenario == 2 || scenario == 3 ) {
693
694 // - Beginn: MIP pions: Nsigma vs. eta, Nsigma vs. mult -
695 if (mom>0.45 && mom<0.5 && sig>40 && sig<60) {
696
697 Bool_t isPionMC=kTRUE;
698
699 TH2 *h1=(TH2*)sublist->At(count*nSpecies);
700 if (h1) {
701 if (scenario == 1) {
702 if ( pdgCodeAbs != 211 && pdgCodeAbs != 111 ) isPionMC=kFALSE;
703 if (isPionMC) {
704 nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion);
705 }
706 }
707 else if (scenario == 2) {
708 nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion, AliTPCPIDResponse::kdEdxHybrid);
709 }
710 else if (scenario == 3) {
711 nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion, AliTPCPIDResponse::kdEdxOROC);
712 }
713 else nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion);
714
715 if (isPionMC) h1->Fill(eta,nSigma);
716 }
717
718 TH2 *h2m=(TH2*)sublist->At(count*nSpecies+1);
719 if ( h2m && isPionMC && mult > 0 ) {
720 h2m->Fill(mult,nSigma);
721 }
722
723 } // - End: MIP pions -
724
725 // - Beginn: Electrons: Nsigma vs. eta, Nsigma vs. mult -
726 if (mom>0.32 && mom<0.38 && eleLineDist>-10. && eleLineDist<15.) {
727
728 Bool_t isElectronMC=kTRUE;
729
730 TH2 *h3=(TH2*)sublist->At(count*nSpecies+2);
731 if (h3) {
732 if (scenario == 1) {
733 if ( pdgCodeAbs != 11 ) isElectronMC=kFALSE;
734 if (isElectronMC) {
735 nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
736 }
737 }
738 if (scenario == 2) {
739 nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kElectron, AliTPCPIDResponse::kdEdxHybrid);
740 }
741 else if (scenario == 3) {
742 nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kElectron, AliTPCPIDResponse::kdEdxOROC);
743 }
744 else nSigma=fPIDResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
745
746 if (isElectronMC) h3->Fill(eta,nSigma);
747 }
748
749 TH2 *h4m=(TH2*)sublist->At(count*nSpecies+3);
750 if ( h4m && isElectronMC && mult > 0 ) {
751 h4m->Fill(mult,nSigma);
752 }
753
754 } // - End: Electrons -
755 } // -- End: Fill histograms for MIP pions and electrons --
756
757}
758
759//______________________________________________________________________________
760void AliAnalysisTaskPIDqa::FillTPCqa()
761{
762 //
763 // Fill PID qa histograms for the TPC
764 //
765
766 // switches for the different scenarios
767 Bool_t scBasic=1; // default/basic
768 Bool_t scMCtruth=1; // for MC truth tracks
769 Bool_t scHybrid=1; // for hybrid PID (only LHC11h)
770 Bool_t scOROChigh=1; // only OROC signal (only LHC11h)
771 Bool_t scV0=1; // for V0 candidates (only for ESDs available)
772 Int_t scCounter=0; // counter of scenarios, used for the histograms at the end of FillTPCqa
773
774 // input handler
775 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
776 AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
777 if (!inputHandler) AliFatal("Input handler needed");
778
779 AliVEvent *event=InputEvent();
780
781 // ESD or AOD event needed to get reference multiplicity (not in AliVEvent)
782 AliAODEvent *fAODevent = 0x0; // AOD event
783 AliESDEvent *fESDevent = 0x0; // ESD event
784 AliESDtrackCuts *esdTrackCuts = 0x0; // ESD track Cuts (ref mult is in AliESDtrackCuts)
785
786 Double_t eta=0.; // track eta
787 Int_t mult=0; // event multiplicity (TPConlyRefMult)
788 //Int_t nacc=0; // counter for accepted multiplicity
789
790 // Check for MC
791 scMCtruth=(MCEvent()!=0x0);
792
793 // Check if period is data LHC11h by checking if
794 // the splines for ALLhigh have been set by AliPIDResponse
795 AliTPCPIDResponse &tpcResp=fPIDResponse->GetTPCResponse();
796 if (tpcResp.GetResponseFunction(AliPID::kPion, AliTPCPIDResponse::kALLhigh)==0x0) {
797 scHybrid = kFALSE;
798 scOROChigh = kFALSE;
799 }
800
801 // Check if "ESD" or "AOD" and get the corresponding event and the beam type (or centrality)
802 TString analysisType = inputHandler->GetDataType(); // can be "ESD" or "AOD"
803 if (analysisType == "ESD") {
804 fESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
805 esdTrackCuts = new AliESDtrackCuts("esdTrackCuts");
806 //printf("\n--- New event - event type = ESD \n");
807 }
808 else if (analysisType == "AOD") {
809 fAODevent = dynamic_cast<AliAODEvent*>( InputEvent() );
810 //printf("\n--- New event - event type = AOD \n");
811
812 // disable V0 scenario, because V0s are not available for AODs in this current implementation
813 scV0=0;
814 }
815
816 // Check if Basic list is already created
817 // If not: Go to SetupTPCqa and creat lists and histograms
818 if(!fListQAtpcBasic) {
819 //printf("\n--- No list QA TPC Basic found -> go to SetupTPCqa! ---\n");
820 SetupTPCqa(scMCtruth, scHybrid, scV0);
821 }
822
823 // Get the number of scenarios by counting those, which are switched on
824 if (scBasic) scCounter++;
825 if (scMCtruth) scCounter++;
826 if (scHybrid) scCounter++;
827 if (scOROChigh) scCounter++;
828 if (scV0) scCounter++;
829
830 // Get reference multiplicity for ESDs
831 if ( analysisType == "ESD" && esdTrackCuts ) {
832 mult=esdTrackCuts->GetReferenceMultiplicity(fESDevent,kTRUE);
833 }
834
835 // Get reference multiplicity for AODs
836 if ( analysisType == "AOD" && fAODevent ) {
837 AliAODHeader * header=dynamic_cast<AliAODHeader*>(fAODevent->GetHeader());
838 if(!header) AliFatal("Not a standard AOD");
839 mult=header->GetTPConlyRefMultiplicity();
840 }
841
842 /*if (mult < 0) {
843 printf("Reference multiplicity not available \n");
844 //return;
845 }*/
846
847 //printf("The multiplicity is = %i ",mult);
848
849
850 // -- Begin: track loop --
851 Int_t ntracks=event->GetNumberOfTracks();
852 for(Int_t itrack = 0; itrack < ntracks; itrack++){
853 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
854
855 //
856 //basic track cuts
857 //
858 ULong_t status=track->GetStatus();
859 // not that nice. status bits not in virtual interface
860 // TPC refit + ITS refit + TPC pid
861 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
862 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
863
864 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
865 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
866 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
867 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
868 if (track->GetTPCNclsF()>0) {
869 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
870 }
871
872 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
873
874 eta=track->Eta();
875 if ( TMath::Abs(eta)>0.9 ) continue;
876
877 //nacc++; // counter for accepted multiplicity
878
879 // the default ("basic") scenario
880 if (scBasic == 1) {
881 FillTPCHistogramsNsigma(fListQAtpcBasic,0,track,mult);
882 FillTPCHistogramsSignal(fListQAtpcBasic,0,track,mult);
883 }
884
885 // only MC truth identified particles
886 if (scMCtruth == 1) {
887 FillTPCHistogramsNsigma(fListQAtpcMCtruth,1,track,mult);
888 }
889
890 // the "hybrid" scenario (only for LHC11h)
891 if (scHybrid == 1) {
892 FillTPCHistogramsNsigma(fListQAtpcHybrid,2,track,mult);
893 }
894
895 // the "OROC high" scenario (only for LHC11h)
896 if (scOROChigh == 1) {
897 FillTPCHistogramsNsigma(fListQAtpcOROChigh,3,track,mult);
898 }
899
900 } // -- End: track loop --
901
902
903 // -- Begin: track loops for V0 candidates --
904 if (scV0 == 1) {
905
906 // - Begin: track loop for electrons from V0 -
907 for(Int_t itrack = 0; itrack < fV0electrons->GetEntries(); itrack++){
908 AliVTrack *track=(AliVTrack*)fV0electrons->At(itrack);
909
910 //
911 //basic track cuts
912 //
913 ULong_t status=track->GetStatus();
914 // not that nice. status bits not in virtual interface
915 // TPC refit + ITS refit + TPC pid
916 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
917 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
918
919 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
920 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
921 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
922 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
923 if (track->GetTPCNclsF()>0) {
924 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
925 }
926
927 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
928
929 eta=track->Eta();
930 if ( TMath::Abs(eta)>0.9 ) continue;
931
932 // fill histograms for V0 candidates
933 FillTPCHistogramsNsigma(fListQAtpcV0,40,track,mult);
934
935 } // - End: track loop for electrons from V0 -
936
937
938 // - Begin: track loop for pions from V0 -
939 for(Int_t itrack = 0; itrack < fV0pions->GetEntries(); itrack++){
940 AliVTrack *track=(AliVTrack*)fV0pions->At(itrack);
941
942 //
943 //basic track cuts
944 //
945 ULong_t status=track->GetStatus();
946 // not that nice. status bits not in virtual interface
947 // TPC refit + ITS refit + TPC pid
948 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
949 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
950
951 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
952 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
953 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
954 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
955 if (track->GetTPCNclsF()>0) {
956 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
957 }
958
959 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
960
961 eta=track->Eta();
962 if ( TMath::Abs(eta)>0.9 ) continue;
963
964 // fill histograms for V0 candidates
965 FillTPCHistogramsNsigma(fListQAtpcV0,42,track,mult);
966
967 } // - End: track loop for pions from V0 -
968
969
970 // - Begin: track loop for kaons from V0 -
971 for(Int_t itrack = 0; itrack < fV0kaons->GetEntries(); itrack++){
972 AliVTrack *track=(AliVTrack*)fV0kaons->At(itrack);
973
974 //
975 //basic track cuts
976 //
977 ULong_t status=track->GetStatus();
978 // not that nice. status bits not in virtual interface
979 // TPC refit + ITS refit + TPC pid
980 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
981 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
982
983 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
984 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
985 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
986 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
987 if (track->GetTPCNclsF()>0) {
988 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
989 }
990
991 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
992
993 eta=track->Eta();
994 if ( TMath::Abs(eta)>0.9 ) continue;
995
996 // fill histograms for V0 candidates
997 FillTPCHistogramsNsigma(fListQAtpcV0,43,track,mult);
998
999 } // - End: track loop for kaons from V0 -
1000
1001
1002 // - Begin: track loop for protons from V0 -
1003 for(Int_t itrack = 0; itrack < fV0protons->GetEntries(); itrack++){
1004 AliVTrack *track=(AliVTrack*)fV0protons->At(itrack);
1005
1006 //
1007 //basic track cuts
1008 //
1009 ULong_t status=track->GetStatus();
1010 // not that nice. status bits not in virtual interface
1011 // TPC refit + ITS refit + TPC pid
1012 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1013 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
1014
1015 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
1016 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
1017 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1018 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1019 if (track->GetTPCNclsF()>0) {
1020 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1021 }
1022
1023 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1024
1025 eta=track->Eta();
1026 if ( TMath::Abs(eta)>0.9 ) continue;
1027
1028 // fill histograms for V0 candidates
1029 FillTPCHistogramsNsigma(fListQAtpcV0,44,track,mult);
1030
1031 } // - End: track loop for protons from V0 -
1032
1033 } // -- End: track loops for V0 candidates --
1034
1035
1036 // Multiplicity distribution
1037 TH1 *hm=(TH1*)fListQAtpc->At(scCounter);
1038 if (hm) {
1039 hm->Fill(mult);
1040 }
1041
1042 //printf("\nAccepted multiplicity = %i \n --- END of event --- \n",nacc);
1043
1044}
1045
1046//______________________________________________________________________________
1047void AliAnalysisTaskPIDqa::FillTRDqa()
1048{
1049 //
1050 // Fill PID qa histograms for the TRD
1051 //
1052 AliVEvent *event=InputEvent();
1053 Int_t ntracks = event->GetNumberOfTracks();
1054 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1055 AliVTrack *track = (AliVTrack *)event->GetTrack(itrack);
1056
1057 //
1058 //basic track cuts
1059 //
1060 ULong_t status=track->GetStatus();
1061 // not that nice. status bits not in virtual interface
1062 // TPC refit + ITS refit + TPC pid + TRD out
1063 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1064 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1065// !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) || //removes light nuclei. So it is out for the moment
1066 !( (status & AliVTrack::kTRDout ) == AliVTrack::kTRDout )) continue;
1067
1068 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1069 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1070 if (track->GetTPCNclsF()>0) {
1071 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1072 }
1073
1074 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1075
1076 Double_t likelihoods[AliPID::kSPECIES];
1077 if(fPIDResponse->ComputeTRDProbability(track, AliPID::kSPECIES, likelihoods) != AliPIDResponse::kDetPidOk) continue;
1078 Int_t ntracklets = 0;
1079 Double_t momentum = -1.;
1080 for(Int_t itl = 0; itl < 6; itl++) {
1081 if(track->GetTRDmomentum(itl) > 0.) {
1082 ntracklets++;
1083 if(momentum < 0) momentum = track->GetTRDmomentum(itl);
1084 }
1085 }
1086
1087 for(Int_t ispecie = 0; ispecie < AliPID::kSPECIES; ispecie++){
1088 TH2F *hLike = (TH2F *)fListQAtrd->At(ntracklets*AliPID::kSPECIES+ispecie);
1089 if (hLike) hLike->Fill(momentum,likelihoods[ispecie]);
1090 }
1091
1092 //=== nSigma and signal ===
1093 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1094 TH2 *h=(TH2*)fListQAtrdNsig->At(ispecie);
1095 TH2 *hTPCTOF=(TH2*)fListQAtrdNsigTPCTOF->At(ispecie);
1096 if (!h || !hTPCTOF) continue;
1097 Float_t nSigmaTPC=fPIDResponse->NumberOfSigmas(AliPIDResponse::kTPC, track, (AliPID::EParticleType)ispecie);
1098 Float_t nSigmaTRD=fPIDResponse->NumberOfSigmas(AliPIDResponse::kTRD, track, (AliPID::EParticleType)ispecie);
1099 Float_t nSigmaTOF=fPIDResponse->NumberOfSigmas(AliPIDResponse::kTOF, track, (AliPID::EParticleType)ispecie);
1100 h->Fill(momentum,nSigmaTRD);
1101
1102 if (TMath::Abs(nSigmaTPC)<3 && TMath::Abs(nSigmaTOF)<3) {
1103 hTPCTOF->Fill(momentum,nSigmaTRD);
1104 }
1105 }
1106
1107 TH2 *h=(TH2*)fListQAtrdNsig->Last();
1108
1109 if (h) {
1110 Double_t sig=track->GetTRDsignal();
1111 h->Fill(momentum,sig);
1112 }
1113
1114 }
1115}
1116
1117//______________________________________________________________________________
1118void AliAnalysisTaskPIDqa::FillTOFqa()
1119{
1120 //
1121 // Fill TOF information
1122 //
1123 AliVEvent *event=InputEvent();
1124
1125 Int_t ntracks=event->GetNumberOfTracks();
1126 Int_t tracksAtTof = 0;
1127 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1128 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1129
1130 //
1131 //basic track cuts
1132 //
1133 ULong_t status=track->GetStatus();
1134 // TPC refit + ITS refit +
1135 // TOF out + kTIME
1136 // kTIME
1137 // (we don't use kTOFmismatch because it depends on TPC and kTOFpid because it prevents light nuclei
1138 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1139 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1140 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
1141 // !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
1142 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
1143
1144 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1145 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1146 if (track->GetTPCNclsF()>0) {
1147 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1148 }
1149
1150 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1151
1152 tracksAtTof++;
1153
1154 Double_t mom=track->P();
1155
1156 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1157 TH2 *h=(TH2*)fListQAtof->At(ispecie);
1158 if (!h) continue;
1159 Double_t nSigma=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
1160 h->Fill(mom,nSigma);
1161 }
1162
1163 TH2 *h=(TH2*)fListQAtof->FindObject("hSigP_TOF");
1164 if (h) {
1165 Double_t sig=track->GetTOFsignal()/1000.;
1166 h->Fill(mom,sig);
1167 }
1168
1169 Int_t mask = fPIDResponse->GetTOFResponse().GetStartTimeMask(mom);
1170 ((TH1F*)fListQAtof->FindObject("hStartTimeMask_TOF"))->Fill((Double_t)(mask+0.5));
1171
1172 if (mom >= 0.75 && mom <= 1.25 ) {
1173 Double_t nsigma= fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)AliPID::kPion);
1174 if (mask == 0) {
1175 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-Fill"))->Fill(nsigma);
1176 } else if (mask == 1) {
1177 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-TOF"))->Fill(nsigma);
1178 } else if ( (mask == 2) || (mask == 4) || (mask == 6) ) {
1179 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-T0"))->Fill(nsigma);
1180 } else {
1181 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-Best"))->Fill(nsigma);
1182 }
1183 if (mask & 0x1) { //at least TOF-T0 present
1184 Double_t delta=0;
1185 (void)fPIDResponse->GetSignalDelta((AliPIDResponse::EDetector)AliPIDResponse::kTOF,track,(AliPID::EParticleType)AliPID::kPion,delta);
1186 ((TH1F*)fListQAtof->FindObject("hDelta_TOF_Pion"))->Fill(delta);
1187 }
1188 }
1189
1190 Double_t res = (Double_t)fPIDResponse->GetTOFResponse().GetStartTimeRes(mom);
1191 ((TH1F*)fListQAtof->FindObject("hStartTimeRes_TOF"))->Fill(res);
1192
1193 Double_t startTimeT0 = event->GetT0TOF(0);
1194 if (startTimeT0 < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeAC_T0"))->Fill(startTimeT0);
1195 else {
1196 startTimeT0 = event->GetT0TOF(1);
1197 if (startTimeT0 < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeA_T0"))->Fill(startTimeT0);
1198 startTimeT0 = event->GetT0TOF(2);
1199 if (startTimeT0 < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeC_T0"))->Fill(startTimeT0);
1200 }
1201 }
1202 if (tracksAtTof > 0) {
1203 ((TH1F* )fListQAtof->FindObject("hnTracksAt_TOF"))->Fill(tracksAtTof);
1204 Int_t mask = fPIDResponse->GetTOFResponse().GetStartTimeMask(5.);
1205 if (mask & 0x1) ((TH1F*)fListQAtof->FindObject("hT0MakerEff"))->Fill(tracksAtTof);
1206 }
1207}
1208
1209//______________________________________________________________________________
1210void AliAnalysisTaskPIDqa::FillT0qa()
1211{
1212 //
1213 // Fill TOF information
1214 //
1215 AliVEvent *event=InputEvent();
1216
1217 Int_t ntracks=event->GetNumberOfTracks();
1218
1219 Int_t tracksAtT0 = 0;
1220
1221 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1222 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1223
1224 //
1225 //basic track cuts
1226 //
1227 ULong_t status=track->GetStatus();
1228 // TPC refit + ITS refit +
1229 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1230 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
1231 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1232 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1233 if (track->GetTPCNclsF()>0) {
1234 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1235 }
1236 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1237
1238 tracksAtT0++;
1239 }
1240
1241 Bool_t t0A = kFALSE;
1242 Bool_t t0C = kFALSE;
1243 Bool_t t0And = kFALSE;
1244 Double_t startTimeT0 = event->GetT0TOF(0); // AND
1245 if (startTimeT0 < 90000) {
1246 t0And = kTRUE;
1247 ((TH1F*)fListQAt0->FindObject("hStartTimeAC_T0"))->Fill(startTimeT0);
1248 }
1249 startTimeT0 = event->GetT0TOF(1); // T0A
1250 if (startTimeT0 < 90000) {
1251 t0A = kTRUE;
1252 ((TH1F*)fListQAt0->FindObject("hStartTimeA_T0"))->Fill(startTimeT0);
1253
1254 }
1255 startTimeT0 = event->GetT0TOF(2); // T0C
1256 if (startTimeT0 < 90000) {
1257 t0C = kTRUE;
1258 ((TH1F*)fListQAt0->FindObject("hStartTimeC_T0"))->Fill(startTimeT0);
1259 }
1260
1261 ((TH1F* )fListQAt0->FindObject("hnTracksAt_T0"))->Fill(tracksAtT0);
1262 if (t0A) ((TH1F*)fListQAt0->FindObject("hT0AEff"))->Fill(tracksAtT0);
1263 if (t0C) ((TH1F*)fListQAt0->FindObject("hT0CEff"))->Fill(tracksAtT0);
1264 if (t0And) ((TH1F*)fListQAt0->FindObject("hT0AndEff"))->Fill(tracksAtT0);
1265 if (t0A || t0C) ((TH1F*)fListQAt0->FindObject("hT0OrEff"))->Fill(tracksAtT0);
1266}
1267
1268
1269//______________________________________________________________________________
1270void AliAnalysisTaskPIDqa::FillEMCALqa()
1271{
1272 //
1273 // Fill PID qa histograms for the EMCAL
1274 //
1275
1276 AliVEvent *event=InputEvent();
1277
1278 Int_t ntracks=event->GetNumberOfTracks();
1279 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1280 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1281
1282 //
1283 //basic track cuts
1284 //
1285 ULong_t status=track->GetStatus();
1286 // not that nice. status bits not in virtual interface
1287 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1288
1289 Double_t pt=track->Pt();
1290
1291 //EMCAL nSigma (only for electrons at the moment)
1292 TH2 *h=(TH2*)fListQAemcal->At(0);
1293 if (!h) continue;
1294 Double_t nSigma=fPIDResponse->NumberOfSigmasEMCAL(track, (AliPID::EParticleType)0);
1295 h->Fill(pt,nSigma);
1296
1297 }
1298
1299 //EMCAL signal (E/p vs. pT) for electrons from V0
1300 for(Int_t itrack = 0; itrack < fV0electrons->GetEntries(); itrack++){
1301 AliVTrack *track=(AliVTrack*)fV0electrons->At(itrack);
1302
1303 //
1304 //basic track cuts
1305 //
1306 ULong_t status=track->GetStatus();
1307 // not that nice. status bits not in virtual interface
1308 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1309
1310 Double_t pt=track->Pt();
1311
1312 TH2 *h=(TH2*)fListQAemcal->At(1);
1313 if (h) {
1314
1315 Int_t nMatchClus = track->GetEMCALcluster();
1316 Double_t mom = track->P();
1317 Double_t eop = -1.;
1318
1319 if(nMatchClus > -1){
1320
1321 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1322
1323 if(matchedClus){
1324
1325 // matched cluster is EMCAL
1326 if(matchedClus->IsEMCAL()){
1327
1328 Double_t fClsE = matchedClus->E();
1329 eop = fClsE/mom;
1330
1331 h->Fill(pt,eop);
1332
1333 }
1334 }
1335 }
1336 }
1337 }
1338
1339 //EMCAL signal (E/p vs. pT) for pions from V0
1340 for(Int_t itrack = 0; itrack < fV0pions->GetEntries(); itrack++){
1341 AliVTrack *track=(AliVTrack*)fV0pions->At(itrack);
1342
1343 //
1344 //basic track cuts
1345 //
1346 ULong_t status=track->GetStatus();
1347 // not that nice. status bits not in virtual interface
1348 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1349
1350 Double_t pt=track->Pt();
1351
1352 TH2 *h=(TH2*)fListQAemcal->At(2);
1353 if (h) {
1354
1355 Int_t nMatchClus = track->GetEMCALcluster();
1356 Double_t mom = track->P();
1357 Double_t eop = -1.;
1358
1359 if(nMatchClus > -1){
1360
1361 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1362
1363 if(matchedClus){
1364
1365 // matched cluster is EMCAL
1366 if(matchedClus->IsEMCAL()){
1367
1368 Double_t fClsE = matchedClus->E();
1369 eop = fClsE/mom;
1370
1371 h->Fill(pt,eop);
1372
1373 }
1374 }
1375 }
1376 }
1377 }
1378
1379 //EMCAL signal (E/p vs. pT) for protons from V0
1380 for(Int_t itrack = 0; itrack < fV0protons->GetEntries(); itrack++){
1381 AliVTrack *track=(AliVTrack*)fV0protons->At(itrack);
1382
1383 //
1384 //basic track cuts
1385 //
1386 ULong_t status=track->GetStatus();
1387 // not that nice. status bits not in virtual interface
1388 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1389
1390 Double_t pt=track->Pt();
1391
1392 TH2 *hP=(TH2*)fListQAemcal->At(3);
1393 TH2 *hAP=(TH2*)fListQAemcal->At(4);
1394 if (hP && hAP) {
1395
1396 Int_t nMatchClus = track->GetEMCALcluster();
1397 Double_t mom = track->P();
1398 Int_t charge = track->Charge();
1399 Double_t eop = -1.;
1400
1401 if(nMatchClus > -1){
1402
1403 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1404
1405 if(matchedClus){
1406
1407 // matched cluster is EMCAL
1408 if(matchedClus->IsEMCAL()){
1409
1410 Double_t fClsE = matchedClus->E();
1411 eop = fClsE/mom;
1412
1413 if(charge > 0) hP->Fill(pt,eop);
1414 else if(charge < 0) hAP->Fill(pt,eop);
1415
1416 }
1417 }
1418 }
1419 }
1420 }
1421
1422}
1423
1424
1425//______________________________________________________________________________
1426void AliAnalysisTaskPIDqa::FillHMPIDqa()
1427{
1428 //
1429 // Fill PID qa histograms for the HMPID
1430 //
1431
1432 AliVEvent *event=InputEvent();
1433
1434 Int_t ntracks=event->GetNumberOfTracks();
1435 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1436 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1437
1438 //
1439 //basic track cuts
1440 //
1441 const ULong_t status=track->GetStatus();
1442 // not that nice. status bits not in virtual interface
1443 // TPC refit + ITS refit +
1444 // TOF out + TOFpid +
1445 // kTIME
1446 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1447 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
1448
1449 const Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1450 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1451 if (track->GetTPCNclsF()>0) {
1452 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1453 }
1454
1455 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1456
1457 const Double_t mom = track->P();
1458 const Double_t ckovAngle = track->GetHMPIDsignal();
1459
1460 Int_t nhists=0;
1461 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1462 if (ispecie==AliPID::kElectron || ispecie==AliPID::kMuon) continue;
1463 TH2 *h=(TH2*)fListQAhmpid->At(nhists);
1464 if (!h) {++nhists; continue;}
1465 const Double_t nSigma=fPIDResponse->NumberOfSigmasHMPID(track, (AliPID::EParticleType)ispecie);
1466 h->Fill(mom,nSigma);
1467 ++nhists;
1468 }
1469
1470 TH1F *hThetavsMom = (TH1F*)fListQAhmpid->At(AliPID::kSPECIESC);
1471
1472 if (hThetavsMom) hThetavsMom->Fill(mom,ckovAngle);
1473
1474 }
1475}
1476//______________________________________________________________________________
1477void AliAnalysisTaskPIDqa::FillTOFHMPIDqa()
1478{
1479 //
1480 // Fill PID qa histograms for the HMPID
1481 //
1482
1483 AliVEvent *event=InputEvent();
1484
1485 Int_t ntracks=event->GetNumberOfTracks();
1486 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1487 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1488
1489 //
1490 //basic track cuts
1491 //
1492 ULong_t status=track->GetStatus();
1493 // not that nice. status bits not in virtual interface
1494 // TPC refit + ITS refit +
1495 // TOF out + TOFpid +
1496 // kTIME
1497 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1498 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1499 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
1500 !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
1501 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
1502
1503 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1504 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1505 if (track->GetTPCNclsF()>0) {
1506 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1507 }
1508
1509 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1510
1511 Double_t mom = track->P();
1512 Double_t ckovAngle = track->GetHMPIDsignal();
1513
1514 Double_t nSigmaTOF[3];
1515 TH1F *h[3];
1516
1517 for (Int_t ispecie=2; ispecie<5; ++ispecie){
1518 //TOF nSigma
1519 nSigmaTOF[ispecie-2]=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
1520 h[ispecie-2] = (TH1F*)fListQAtofhmpid->At(ispecie-2);}
1521
1522 if(TMath::Abs(nSigmaTOF[0])<2) h[0]->Fill(mom,ckovAngle);
1523
1524 if(TMath::Abs(nSigmaTOF[1])<2 && TMath::Abs(nSigmaTOF[0])>3) h[1]->Fill(mom,ckovAngle);
1525
1526 if(TMath::Abs(nSigmaTOF[2])<2 && TMath::Abs(nSigmaTOF[1])>3 && TMath::Abs(nSigmaTOF[0])>3) h[2]->Fill(mom,ckovAngle);
1527
1528 }
1529
1530}
1531
1532//______________________________________________________________________________
1533void AliAnalysisTaskPIDqa::FillTPCTOFqa()
1534{
1535 //
1536 // Fill PID qa histograms for the TOF
1537 // Here also the TPC histograms after TOF selection are filled
1538 //
1539
1540 AliVEvent *event=InputEvent();
1541
1542 Int_t ntracks=event->GetNumberOfTracks();
1543 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1544 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1545
1546 //
1547 //basic track cuts
1548 //
1549 ULong_t status=track->GetStatus();
1550 // not that nice. status bits not in virtual interface
1551 // TPC refit + ITS refit +
1552 // TOF out + TOFpid +
1553 // kTIME
1554 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1555 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1556// !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) || //removes light nuclei, so it is out for the moment
1557 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
1558 !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
1559 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
1560
1561 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1562 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1563 if (track->GetTPCNclsF()>0) {
1564 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1565 }
1566
1567 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1568
1569
1570 Double_t mom=track->P();
1571 Double_t momTPC=track->GetTPCmomentum();
1572
1573 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1574 //TOF nSigma
1575 Double_t nSigmaTOF=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
1576 Double_t nSigmaTPC=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
1577
1578 //TPC after TOF cut
1579 TH2 *h=(TH2*)fListQAtpctof->At(ispecie);
1580 if (h && TMath::Abs(nSigmaTOF)<3.) h->Fill(momTPC,nSigmaTPC);
1581
1582 //TOF after TPC cut
1583 h=(TH2*)fListQAtpctof->At(ispecie+AliPID::kSPECIESC);
1584 if (h && TMath::Abs(nSigmaTPC)<3.) h->Fill(mom,nSigmaTOF);
1585
1586 //EMCAL after TOF and TPC cut
1587 h=(TH2*)fListQAtpctof->At(ispecie+2*AliPID::kSPECIESC);
1588 if (h && TMath::Abs(nSigmaTOF)<3. && TMath::Abs(nSigmaTPC)<3. ){
1589
1590 Int_t nMatchClus = track->GetEMCALcluster();
1591 Double_t pt = track->Pt();
1592 Double_t eop = -1.;
1593
1594 if(nMatchClus > -1){
1595
1596 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1597
1598 if(matchedClus){
1599
1600 // matched cluster is EMCAL
1601 if(matchedClus->IsEMCAL()){
1602
1603 Double_t fClsE = matchedClus->E();
1604 eop = fClsE/mom;
1605
1606 h->Fill(pt,eop);
1607
1608
1609 }
1610 }
1611 }
1612 }
1613 }
1614 }
1615}
1616
1617//_____________________________________________________________________________
1618void AliAnalysisTaskPIDqa::FillQAinfo()
1619{
1620 //
1621 // Fill the QA information
1622 //
1623
1624
1625 //TPC QA info
1626 TObjArray *arrTPC=static_cast<TObjArray*>(fListQAinfo->At(0));
1627 if (fPIDResponse && arrTPC){
1628 AliTPCPIDResponse &tpcResp=fPIDResponse->GetTPCResponse();
1629 // fill spline names
1630 if (!arrTPC->UncheckedAt(0)){
1631
1632 TObjArray *arrTPCsplineNames=new TObjArray(AliPID::kSPECIESC);
1633 arrTPCsplineNames->SetOwner();
1634 arrTPCsplineNames->SetName("TPC_spline_names");
1635 arrTPC->AddAt(arrTPCsplineNames,0);
1636
1637 for (Int_t iresp=0; iresp<AliPID::kSPECIESC; ++iresp){
1638 const TObject *o=tpcResp.GetResponseFunction((AliPID::EParticleType)iresp);
1639 if (!o) continue;
1640 arrTPCsplineNames->Add(new TObjString(Form("%02d: %s",iresp, o->GetName())));
1641 }
1642 }
1643
1644 // tpc response config
1645 if (!arrTPC->UncheckedAt(1)){
1646
1647 TObjArray *arrTPCconfigInfo=new TObjArray;
1648 arrTPCconfigInfo->SetOwner();
1649 arrTPCconfigInfo->SetName("TPC_config_info");
1650 arrTPC->AddAt(arrTPCconfigInfo,1);
1651
1652 TObjString *ostr=0x0;
1653 ostr=new TObjString;
1654 ostr->String().Form("Eta Corr map: %s", tpcResp.GetEtaCorrMap()?tpcResp.GetEtaCorrMap()->GetName():"none");
1655 arrTPCconfigInfo->Add(ostr);
1656
1657 ostr=new TObjString;
1658 ostr->String().Form("Sigma Par map: %s", tpcResp.GetSigmaPar1Map()?tpcResp.GetSigmaPar1Map()->GetName():"none");
1659 arrTPCconfigInfo->Add(ostr);
1660
1661 ostr=new TObjString;
1662 ostr->String().Form("MIP: %.2f", tpcResp.GetMIP());
1663 arrTPCconfigInfo->Add(ostr);
1664
1665 ostr=new TObjString;
1666 ostr->String().Form("Res: Def %.3g (%.3g) : AllHigh %.3g (%.3g) : OROC high %.3g (%.3g)",
1667 tpcResp.GetRes0(AliTPCPIDResponse::kDefault), tpcResp.GetResN2(AliTPCPIDResponse::kDefault),
1668 tpcResp.GetRes0(AliTPCPIDResponse::kALLhigh), tpcResp.GetResN2(AliTPCPIDResponse::kALLhigh),
1669 tpcResp.GetRes0(AliTPCPIDResponse::kOROChigh), tpcResp.GetResN2(AliTPCPIDResponse::kOROChigh)
1670 );
1671 arrTPCconfigInfo->Add(ostr);
1672 }
1673 }
1674}
1675
1676//______________________________________________________________________________
1677void AliAnalysisTaskPIDqa::SetupITSqa()
1678{
1679 //
1680 // Create the ITS qa objects
1681 //
1682
1683 TVectorD *vX=MakeLogBinning(200,.1,30);
1684
1685 //ITS+TPC tracks
1686 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1687 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_ITS_%s",AliPID::ParticleName(ispecie)),
1688 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
1689 vX->GetNrows()-1,vX->GetMatrixArray(),
1690 200,-10,10);
1691 fListQAits->Add(hNsigmaP);
1692 }
1693 TH2F *hSig = new TH2F("hSigP_ITS",
1694 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
1695 vX->GetNrows()-1,vX->GetMatrixArray(),
1696 300,0,300);
1697 fListQAits->Add(hSig);
1698
1699 //ITS Standalone tracks
1700 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1701 TH2F *hNsigmaPSA = new TH2F(Form("hNsigmaP_ITSSA_%s",AliPID::ParticleName(ispecie)),
1702 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
1703 vX->GetNrows()-1,vX->GetMatrixArray(),
1704 200,-10,10);
1705 fListQAitsSA->Add(hNsigmaPSA);
1706 }
1707 TH2F *hSigSA = new TH2F("hSigP_ITSSA",
1708 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
1709 vX->GetNrows()-1,vX->GetMatrixArray(),
1710 300,0,300);
1711 fListQAitsSA->Add(hSigSA);
1712
1713 //ITS Pure Standalone tracks
1714 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1715 TH2F *hNsigmaPPureSA = new TH2F(Form("hNsigmaP_ITSPureSA_%s",AliPID::ParticleName(ispecie)),
1716 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
1717 vX->GetNrows()-1,vX->GetMatrixArray(),
1718 200,-10,10);
1719 fListQAitsPureSA->Add(hNsigmaPPureSA);
1720 }
1721 TH2F *hSigPureSA = new TH2F("hSigP_ITSPureSA",
1722 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
1723 vX->GetNrows()-1,vX->GetMatrixArray(),
1724 300,0,300);
1725 fListQAitsPureSA->Add(hSigPureSA);
1726
1727 delete vX;
1728}
1729
1730//_____________________________________________________________________________
1731void AliAnalysisTaskPIDqa::AddTPCHistogramsSignal(TList *sublist, const char *scenario)
1732{
1733 //
1734 // Create the TPC qa objects: create histograms for the TPC signal for different settings
1735 //
1736
1737 TVectorD *vX=MakeLogBinning(200,.1,30);
1738 Int_t nBinsMult = 38;
1739 Double_t xBinsMult[39] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100,
1740 120, 140, 160, 180, 200,
1741 300, 400, 500, 600, 700, 800, 900, 1000,
1742 1200, 1400, 1600, 1800, 2000,
1743 2200, 2400, 2600, 2800, 3000,
1744 3200, 3400, 3600, 3800, 4000
1745 };
1746 const Int_t binsEta=110;
1747 Float_t etaMin=-1.1;
1748 Float_t etaMax=1.1;
1749
1750 char signal[4][12]={"std","IROC","OROCmedium","OROClong"};
1751
1752
1753 // TPC signal vs. p for all particles (standard, IROC, OROCmedium, OROClong)
1754 for (Int_t iSig=0; iSig<4; iSig++) {
1755 TH2F *hSigP = new TH2F(Form("hSigP_TPC_%s_%s",signal[iSig],scenario),
1756 Form("TPC_%s signal (%s) vs. p;p [GeV]; TPC signal [arb. units]",scenario,signal[iSig]),
1757 vX->GetNrows()-1,vX->GetMatrixArray(),
1758 300,0,300);
1759 sublist->Add(hSigP);
1760 }
1761
1762 // MIP pions: TPC signal vs. eta
1763 for (Int_t iSig=0; iSig<4; iSig++) {
1764 TH2F *hSigEtaMIPpi = new TH2F(Form("hSigEta_TPC_%s_%s_MIPpi",signal[iSig],scenario),
1765 Form("TPC_%s signal (%s) MIPpi vs. eta;#eta;TPC signal [arb. units]",scenario,signal[iSig]),
1766 binsEta,etaMin,etaMax,
1767 300,0,300);
1768 sublist->Add(hSigEtaMIPpi);
1769 }
1770
1771 // MIP pions: TPC signal vs. multiplicity
1772 for (Int_t iSig=0; iSig<4; iSig++) {
1773 TH2F *hSigMultMPIpi = new TH2F(Form("hSigMult_TPC_%s_%s_MIPpi",signal[iSig],scenario),
1774 Form("TPC_%s signal (%s) MIPpi vs. mult;multiplicity;TPC signal [arb. units]",scenario,signal[iSig]),
1775 nBinsMult,xBinsMult,
1776 300,0,300);
1777 sublist->Add(hSigMultMPIpi);
1778 }
1779
1780 // Electrons: TPC signal vs. eta
1781 for (Int_t iSig=0; iSig<4; iSig++) {
1782 TH2F *hSigEtaEle = new TH2F(Form("hSigEta_TPC_%s_%s_Ele",signal[iSig],scenario),
1783 Form("TPC_%s signal (%s) electrons vs. eta;#eta;TPC signal [arb. units]",scenario,signal[iSig]),
1784 binsEta,etaMin,etaMax,
1785 300,0,300);
1786 sublist->Add(hSigEtaEle);
1787 }
1788
1789 // Electrons: TPC signal vs. multiplicity
1790 for (Int_t iSig=0; iSig<4; iSig++) {
1791 TH2F *hSigMultEle = new TH2F(Form("hSigMult_TPC_%s_%s_Ele",signal[iSig],scenario),
1792 Form("TPC_%s signal (%s) electrons vs. mult;multiplicity;TPC signal [arb. units]",scenario,signal[iSig]),
1793 nBinsMult,xBinsMult,
1794 300,0,300);
1795 sublist->Add(hSigMultEle);
1796 }
1797
1798 delete vX;
1799
1800}
1801
1802//_____________________________________________________________________________
1803void AliAnalysisTaskPIDqa::AddTPCHistogramsNsigma(TList *sublist, const char *scenario, Int_t scnumber)
1804{
1805 //
1806 // Create the TPC qa objects: create histograms for TPC Nsigma for different settings
1807 //
1808
1809 TVectorD *vX=MakeLogBinning(200,.1,30.);
1810 Int_t nBinsMult = 38;
1811 Double_t xBinsMult[39] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100,
1812 120, 140, 160, 180, 200,
1813 300, 400, 500, 600, 700, 800, 900, 1000,
1814 1200, 1400, 1600, 1800, 2000,
1815 2200, 2400, 2600, 2800, 3000,
1816 3200, 3400, 3600, 3800, 4000
1817 };
1818 const Int_t binsEta=110;
1819 Float_t etaMin=-1.1;
1820 Float_t etaMax=1.1;
1821
1822 Int_t nSpecies=0;
1823
1824 if (scnumber == 4) nSpecies=(Int_t)AliPID::kSPECIES;
1825 else nSpecies=(Int_t)AliPID::kSPECIESC;
1826
1827 // Nsigma vs. p for different particle species
1828 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
1829 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_%s_%s",scenario,AliPID::ParticleName(ispecie)),
1830 Form("TPC_%s n#sigma %s vs. p;p [GeV]; n#sigma",scenario,AliPID::ParticleName(ispecie)),
1831 vX->GetNrows()-1,vX->GetMatrixArray(),
1832 200,-10,10);
1833 sublist->Add(hNsigmaP);
1834 }
1835
1836 // Nsigma vs. eta for different particle species (only for some scenarios)
1837 if ( scnumber == 1 || scnumber == 4 ) {
1838 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
1839 TH2F *hNsigmaEta = new TH2F(Form("hNsigmaEta_TPC_%s_%s",scenario,AliPID::ParticleName(ispecie)),
1840 Form("TPC_%s n#sigma %s vs. eta;#eta; n#sigma",scenario,AliPID::ParticleName(ispecie)),
1841 binsEta,etaMin,etaMax,
1842 200,-10,10);
1843 sublist->Add(hNsigmaEta);
1844 }
1845 }
1846
1847 // Nsigma vs. multiplicity for different particle species (only for some scenarios)
1848 if ( scnumber == 1 || scnumber == 4 ) {
1849 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
1850 TH2F *hNsigmaMult = new TH2F(Form("hNsigmaMult_TPC_%s_%s",scenario,AliPID::ParticleName(ispecie)),
1851 Form("TPC_%s n#sigma %s vs. mult;multiplicity; n#sigma",scenario,AliPID::ParticleName(ispecie)),
1852 nBinsMult,xBinsMult,
1853 200,-10,10);
1854 sublist->Add(hNsigmaMult);
1855 }
1856 }
1857
1858 // - Beginn: Adding histograms for MIP pions and electrons (only for some scenarios) -
1859 if ( scnumber == 0 || scnumber == 2 || scnumber == 3 ) {
1860
1861 // MIP pions: Nsigma vs. eta
1862 TH2F *hNsigmaEtaMIPpi = new TH2F(Form("hNsigmaEta_TPC_%s_MIPpi",scenario),
1863 Form("TPC_%s n#sigma MIPpi vs. eta;#eta; n#sigma",scenario),
1864 binsEta,etaMin,etaMax,
1865 200,-10,10);
1866 sublist->Add(hNsigmaEtaMIPpi);
1867
1868 // MIP pions: Nsigma vs. multiplicity
1869 TH2F *hNsigmaMultMIPpi = new TH2F(Form("hNsigmaMult_TPC_%s_MIPpi",scenario),
1870 Form("TPC_%s n#sigma MIPpi vs. mult;multiplicity; n#sigma",scenario),
1871 nBinsMult,xBinsMult,
1872 200,-10,10);
1873 sublist->Add(hNsigmaMultMIPpi);
1874
1875 // Electrons: Nsigma vs. eta
1876 TH2F *hNsigmaEtaEle = new TH2F(Form("hNsigmaEta_TPC_%s_Ele",scenario),
1877 Form("TPC_%s n#sigma electrons vs. eta;#eta; n#sigma",scenario),
1878 binsEta,etaMin,etaMax,
1879 200,-10,10);
1880 sublist->Add(hNsigmaEtaEle);
1881
1882 // Electrons: Nsigma vs. multiplicity
1883 TH2F *hNsigmaMultEle = new TH2F(Form("hNsigmaMult_TPC_%s_Ele",scenario),
1884 Form("TPC_%s n#sigma electrons vs. mult;multiplicity; n#sigma",scenario),
1885 nBinsMult,xBinsMult,
1886 200,-10,10);
1887 sublist->Add(hNsigmaMultEle);
1888 } // - End: Adding histograms for MIP pions and electrons
1889
1890 delete vX;
1891
1892}
1893
1894//______________________________________________________________________________
1895void AliAnalysisTaskPIDqa::SetupTPCqa(Bool_t fillMC, Bool_t fill11h, Bool_t fillV0)
1896{
1897 //
1898 // Create the TPC qa objects
1899 //
1900
1901 // Set up the multiplicity binning
1902 Int_t nBinsMult = 38;
1903 Double_t xBinsMult[39] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100,
1904 120, 140, 160, 180, 200,
1905 300, 400, 500, 600, 700, 800, 900, 1000,
1906 1200, 1400, 1600, 1800, 2000,
1907 2200, 2400, 2600, 2800, 3000,
1908 3200, 3400, 3600, 3800, 4000
1909 };
1910
1911
1912 // Create TPC sublists for different scenarios
1913 // corresponding to available information,
1914 // e.g. MC or not, special settings for LHC11h
1915
1916 // basic/default scenario, used always
1917 fListQAtpcBasic=new TList;
1918 fListQAtpcBasic->SetOwner();
1919 fListQAtpcBasic->SetName("TPCBasic");
1920 fListQAtpc->Add(fListQAtpcBasic);
1921
1922 // MC truth scenario: use only MC truth identified particles
1923 // only available for MC
1924 if (fillMC == kTRUE) {
1925 fListQAtpcMCtruth=new TList;
1926 fListQAtpcMCtruth->SetOwner();
1927 fListQAtpcMCtruth->SetName("TPCMCtruth");
1928 fListQAtpc->Add(fListQAtpcMCtruth);
1929 }
1930
1931 // Hybrid and OROChigh scenarios,
1932 // special settings only available for PbPb LHC11h data
1933 if (fill11h == kTRUE) {
1934 fListQAtpcHybrid=new TList;
1935 fListQAtpcHybrid->SetOwner();
1936 fListQAtpcHybrid->SetName("TPCHybrid");
1937 fListQAtpc->Add(fListQAtpcHybrid);
1938
1939 fListQAtpcOROChigh=new TList;
1940 fListQAtpcOROChigh->SetOwner();
1941 fListQAtpcOROChigh->SetName("TPCOROChigh");
1942 fListQAtpc->Add(fListQAtpcOROChigh);
1943 }
1944
1945 // scenario only for V0s,
1946 // only available for ESDs
1947 if (fillV0 == kTRUE) {
1948 fListQAtpcV0=new TList;
1949 fListQAtpcV0->SetOwner();
1950 fListQAtpcV0->SetName("TPCV0");
1951 fListQAtpc->Add(fListQAtpcV0);
1952 }
1953
1954
1955 // the default ("basic") scenario
1956 AddTPCHistogramsNsigma(fListQAtpcBasic,"Basic",0);
1957 AddTPCHistogramsSignal(fListQAtpcBasic,"Basic");
1958
1959 // only MC truth identified particles
1960 if (fillMC) {
1961 AddTPCHistogramsNsigma(fListQAtpcMCtruth,"MCtruth",1);
1962 }
1963
1964 // the "hybrid" scenario (only for period LHC11h)
1965 if (fill11h) {
1966 AddTPCHistogramsNsigma(fListQAtpcHybrid,"Hybrid",2);
1967 }
1968
1969 // the "OROC high" scenario (only for period LHC11h)
1970 if (fill11h) {
1971 AddTPCHistogramsNsigma(fListQAtpcOROChigh,"OROChigh",3);
1972 }
1973
1974 // only for V0s
1975 if (fillV0) {
1976 AddTPCHistogramsNsigma(fListQAtpcV0,"V0",4);
1977 }
1978
1979
1980 // Multiplicity distribution --- as check
1981 TH1F *hMult = new TH1F("hMult_TPC",
1982 "Multiplicity distribution;multiplicity;counts",
1983 nBinsMult,xBinsMult);
1984 fListQAtpc->Add(hMult);
1985
1986}
1987
1988//______________________________________________________________________________
1989void AliAnalysisTaskPIDqa::SetupTRDqa()
1990{
1991 //
1992 // Create the TRD qa objects
1993 //
1994 TVectorD *vX=MakeLogBinning(200,.1,30);
1995 for(Int_t itl = 0; itl < 6; ++itl){
1996 for(Int_t ispecie = 0; ispecie < AliPID::kSPECIES; ispecie++){
1997 TH2F *hLikeP = new TH2F(Form("hLikeP_TRD_%dtls_%s", itl, AliPID::ParticleName(ispecie)),
1998 Form("TRD Likelihood to be %s %s for tracks having %d %s; p (GeV/c); TRD %s Likelihood", ispecie == 0 ? "an" : "a", AliPID::ParticleName(ispecie), itl+1, itl == 0 ? "tracklet" : "tracklets", AliPID::ParticleName(ispecie)),
1999 vX->GetNrows()-1, vX->GetMatrixArray(),
2000 100, 0., 1.);
2001 fListQAtrd->Add(hLikeP);
2002 }
2003 }
2004
2005 // === nSigma Values and signal ===
2006 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2007 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TRD_%s",AliPID::ParticleName(ispecie)),
2008 Form("TRD n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2009 vX->GetNrows()-1,vX->GetMatrixArray(),
2010 100,-10,10);
2011 fListQAtrdNsig->Add(hNsigmaP);
2012 }
2013
2014 TH2F *hSig = new TH2F("hSigP_TRD",
2015 "TRD signal vs. p;p [GeV]; TRD signal [arb. units]",
2016 vX->GetNrows()-1,vX->GetMatrixArray(),
2017 100,0,100);
2018 fListQAtrdNsig->Add(hSig);
2019
2020 fListQAtrd->Add(fListQAtrdNsig);
2021
2022 // === Same after 3 sigma in TPC and TOF
2023 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2024 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TRD_TPCTOF_%s",AliPID::ParticleName(ispecie)),
2025 Form("TRD n#sigma %s vs. p after 3#sigma cut in TPC&TOF;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2026 vX->GetNrows()-1,vX->GetMatrixArray(),
2027 100,-10,10);
2028 fListQAtrdNsigTPCTOF->Add(hNsigmaP);
2029 }
2030
2031 fListQAtrd->Add(fListQAtrdNsigTPCTOF);
2032
2033 delete vX;
2034}
2035
2036//______________________________________________________________________________
2037void AliAnalysisTaskPIDqa::SetupTOFqa()
2038{
2039 //
2040 // Create the TOF qa objects
2041 //
2042
2043 TVectorD *vX=MakeLogBinning(200,.1,30);
2044
2045 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2046 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_%s",AliPID::ParticleName(ispecie)),
2047 Form("TOF n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2048 vX->GetNrows()-1,vX->GetMatrixArray(),
2049 200,-10,10);
2050 fListQAtof->Add(hNsigmaP);
2051 }
2052
2053 TH1F *hnSigT0Fill = new TH1F("hNsigma_TOF_Pion_T0-Fill","TOF n#sigma (Pion) T0-FILL [0.75-1.25. GeV/c]",200,-10,10);
2054 fListQAtof->Add(hnSigT0Fill);
2055 TH1F *hnSigT0T0 = new TH1F("hNsigma_TOF_Pion_T0-T0","TOF n#sigma (Pion) T0-T0 [0.75-1.25 GeV/c]",200,-10,10);
2056 fListQAtof->Add(hnSigT0T0);
2057 TH1F *hnSigT0TOF = new TH1F("hNsigma_TOF_Pion_T0-TOF","TOF n#sigma (Pion) T0-TOF [0.75-1.25 GeV/c]",200,-10,10);
2058 fListQAtof->Add(hnSigT0TOF);
2059 TH1F *hnSigT0Best = new TH1F("hNsigma_TOF_Pion_T0-Best","TOF n#sigma (Pion) T0-Best [0.75-1.25 GeV/c]",200,-10,10);
2060 fListQAtof->Add(hnSigT0Best);
2061 TH1F *hnDeltaPi = new TH1F("hDelta_TOF_Pion","DeltaT (Pion) [0.75-1.25 GeV/c]",50,-500,500);
2062 fListQAtof->Add(hnDeltaPi);
2063
2064 TH2F *hSig = new TH2F("hSigP_TOF",
2065 "TOF signal vs. p;p [GeV]; TOF signal [ns]",
2066 vX->GetNrows()-1,vX->GetMatrixArray(),
2067 300,0,30);
2068
2069 delete vX;
2070
2071 fListQAtof->Add(hSig);
2072
2073 TH1F *hStartTimeMaskTOF = new TH1F("hStartTimeMask_TOF","StartTime mask",8,0,8);
2074 fListQAtof->Add(hStartTimeMaskTOF);
2075 TH1F *hStartTimeResTOF = new TH1F("hStartTimeRes_TOF","StartTime resolution [ps]",100,0,500);
2076 fListQAtof->Add(hStartTimeResTOF);
2077
2078 TH1F *hnTracksAtTOF = new TH1F("hnTracksAt_TOF","Matched tracks at TOF",100,0,100);
2079 fListQAtof->Add(hnTracksAtTOF);
2080 TH1F *hT0MakerEff = new TH1F("hT0MakerEff","Events with T0-TOF vs nTracks",100,0,100);
2081 fListQAtof->Add(hT0MakerEff);
2082
2083 // this in principle should stay on a T0 PID QA, but are just the data prepared for TOF use
2084 TH1F *hStartTimeAT0 = new TH1F("hStartTimeA_T0","StartTime from T0A [ps]",1000,-1000,1000);
2085 fListQAtof->Add(hStartTimeAT0);
2086 TH1F *hStartTimeCT0 = new TH1F("hStartTimeC_T0","StartTime from T0C [ps]",1000,-1000,1000);
2087 fListQAtof->Add(hStartTimeCT0);
2088 TH1F *hStartTimeACT0 = new TH1F("hStartTimeAC_T0","StartTime from T0AC [ps]",1000,-1000,1000);;
2089 fListQAtof->Add(hStartTimeACT0);
2090}
2091
2092
2093//______________________________________________________________________________
2094void AliAnalysisTaskPIDqa::SetupT0qa()
2095{
2096 //
2097 // Create the T0 qa objects
2098 //
2099
2100 // these are similar to plots inside TOFqa, but these are for all events
2101 TH1F *hStartTimeAT0 = new TH1F("hStartTimeA_T0","StartTime from T0A [ps]",1000,-1000,1000);
2102 fListQAt0->Add(hStartTimeAT0);
2103 TH1F *hStartTimeCT0 = new TH1F("hStartTimeC_T0","StartTime from T0C [ps]",1000,-1000,1000);
2104 fListQAt0->Add(hStartTimeCT0);
2105 TH1F *hStartTimeACT0 = new TH1F("hStartTimeAC_T0","StartTime from T0AC [ps]",1000,-1000,1000);;
2106 fListQAt0->Add(hStartTimeACT0);
2107
2108 TH1F *hnTracksAtT0 = new TH1F("hnTracksAt_T0","Tracks for events selected for T0",100,0,100);
2109 fListQAt0->Add(hnTracksAtT0);
2110 TH1F *hT0AEff = new TH1F("hT0AEff","Events with T0A vs nTracks",100,0,100);
2111 fListQAt0->Add(hT0AEff);
2112 TH1F *hT0CEff = new TH1F("hT0CEff","Events with T0C vs nTracks",100,0,100);
2113 fListQAt0->Add(hT0CEff);
2114 TH1F *hT0AndEff = new TH1F("hT0AndEff","Events with T0AC (AND) vs nTracks",100,0,100);
2115 fListQAt0->Add(hT0AndEff);
2116 TH1F *hT0OrEff = new TH1F("hT0OrEff","Events with T0AC (OR) vs nTracks",100,0,100);
2117 fListQAt0->Add(hT0OrEff);
2118
2119
2120}
2121
2122//______________________________________________________________________________
2123void AliAnalysisTaskPIDqa::SetupEMCALqa()
2124{
2125 //
2126 // Create the EMCAL qa objects
2127 //
2128
2129 TVectorD *vX=MakeLogBinning(200,.1,30);
2130
2131 TH2F *hNsigmaPt = new TH2F(Form("hNsigmaPt_EMCAL_%s",AliPID::ParticleName(0)),
2132 Form("EMCAL n#sigma %s vs. p_{T};p_{T} [GeV]; n#sigma",AliPID::ParticleName(0)),
2133 vX->GetNrows()-1,vX->GetMatrixArray(),
2134 200,-10,10);
2135 fListQAemcal->Add(hNsigmaPt);
2136
2137
2138 TH2F *hSigPtEle = new TH2F("hSigPt_EMCAL_Ele",
2139 "EMCAL signal (E/p) vs. p_{T} for electrons;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2140 vX->GetNrows()-1,vX->GetMatrixArray(),
2141 200,0,2);
2142 fListQAemcal->Add(hSigPtEle);
2143
2144 TH2F *hSigPtPions = new TH2F("hSigPt_EMCAL_Pions",
2145 "EMCAL signal (E/p) vs. p_{T} for pions;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2146 vX->GetNrows()-1,vX->GetMatrixArray(),
2147 200,0,2);
2148 fListQAemcal->Add(hSigPtPions);
2149
2150 TH2F *hSigPtProtons = new TH2F("hSigPt_EMCAL_Protons",
2151 "EMCAL signal (E/p) vs. p_{T} for protons;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2152 vX->GetNrows()-1,vX->GetMatrixArray(),
2153 200,0,2);
2154 fListQAemcal->Add(hSigPtProtons);
2155
2156 TH2F *hSigPtAntiProtons = new TH2F("hSigPt_EMCAL_Antiprotons",
2157 "EMCAL signal (E/p) vs. p_{T} for antiprotons;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2158 vX->GetNrows()-1,vX->GetMatrixArray(),
2159 200,0,2);
2160 fListQAemcal->Add(hSigPtAntiProtons);
2161
2162 delete vX;
2163}
2164
2165//______________________________________________________________________________
2166void AliAnalysisTaskPIDqa::SetupHMPIDqa()
2167{
2168 //
2169 // Create the HMPID qa objects
2170 //
2171
2172 TVectorD *vX=MakeLogBinning(200,.1,30);
2173
2174 // nSigmas
2175 Int_t nhists=0;
2176 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
2177 if (ispecie==AliPID::kElectron || ispecie==AliPID::kMuon) continue;
2178 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_HMPID_%s",AliPID::ParticleName(ispecie)),
2179 Form("HMPID n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2180 vX->GetNrows()-1,vX->GetMatrixArray(),
2181 200,-10,10);
2182 fListQAhmpid->AddAt(hNsigmaP, nhists);
2183 ++nhists;
2184 }
2185
2186 // cherenkov angle
2187 TH2F *hCkovAnglevsMom = new TH2F("hCkovAnglevsMom", "Cherenkov angle vs momentum",
2188 vX->GetNrows()-1,vX->GetMatrixArray(),
2189 500,0,1);
2190 fListQAhmpid->AddAt(hCkovAnglevsMom,nhists);
2191
2192 delete vX;
2193}
2194
2195//______________________________________________________________________________
2196void AliAnalysisTaskPIDqa::SetupTOFHMPIDqa()
2197{
2198 //
2199 // Create the HMPID qa objects
2200 //
2201
2202 TH2F *hCkovAnglevsMomPion = new TH2F("hCkovAnglevsMom_pion", "Cherenkov angle vs momentum for pions",500,0,5.,500,0,1);
2203 fListQAtofhmpid->Add(hCkovAnglevsMomPion);
2204
2205 TH2F *hCkovAnglevsMomKaon = new TH2F("hCkovAnglevsMom_kaon", "Cherenkov angle vs momentum for kaons",500,0,5.,500,0,1);
2206 fListQAtofhmpid->Add(hCkovAnglevsMomKaon);
2207
2208 TH2F *hCkovAnglevsMomProton = new TH2F("hCkovAnglevsMom_proton","Cherenkov angle vs momentum for protons",500,0,5.,500,0,1);
2209 fListQAtofhmpid->Add(hCkovAnglevsMomProton);
2210
2211
2212}
2213
2214//______________________________________________________________________________
2215void AliAnalysisTaskPIDqa::SetupTPCTOFqa()
2216{
2217 //
2218 // Create the qa objects for TPC + TOF combination
2219 //
2220
2221 TVectorD *vX=MakeLogBinning(200,.1,30);
2222
2223 //TPC signals after TOF cut
2224 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2225 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_TOF_%s",AliPID::ParticleName(ispecie)),
2226 Form("TPC n#sigma %s vs. p (after TOF 3#sigma cut);p_{TPC} [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2227 vX->GetNrows()-1,vX->GetMatrixArray(),
2228 200,-10,10);
2229 fListQAtpctof->Add(hNsigmaP);
2230 }
2231
2232 //TOF signals after TPC cut
2233 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2234 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_TPC_%s",AliPID::ParticleName(ispecie)),
2235 Form("TOF n#sigma %s vs. p (after TPC n#sigma cut);p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2236 vX->GetNrows()-1,vX->GetMatrixArray(),
2237 200,-10,10);
2238 fListQAtpctof->Add(hNsigmaP);
2239 }
2240
2241 //EMCAL signal after TOF and TPC cut
2242 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
2243 TH2F *heopPt = new TH2F(Form("heopPt_TOF_TPC_%s",AliPID::ParticleName(ispecie)),
2244 Form("EMCAL signal (E/p) %s vs. p_{T};p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",AliPID::ParticleName(ispecie)),
2245 vX->GetNrows()-1,vX->GetMatrixArray(),
2246 200,0,2);
2247 fListQAtpctof->Add(heopPt);
2248 }
2249
2250 delete vX;
2251}
2252//______________________________________________________________________________
2253void AliAnalysisTaskPIDqa::SetupV0qa()
2254{
2255 //
2256 // Create the qa objects for V0 Kine cuts
2257 //
2258
2259 TH2F *hArmenteros = new TH2F("hArmenteros", "Armenteros plot",200,-1.,1.,200,0.,0.4);
2260 fListQAV0->Add(hArmenteros);
2261
2262}
2263
2264//_____________________________________________________________________________
2265void AliAnalysisTaskPIDqa::SetupQAinfo(){
2266 //
2267 // Setup the info of QA objects
2268 //
2269
2270 TObjArray *arr=new TObjArray;
2271 arr->SetName("TPC_info");
2272 fListQAinfo->Add(arr);
2273}
2274
2275//______________________________________________________________________________
2276TVectorD* AliAnalysisTaskPIDqa::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
2277{
2278 //
2279 // Make logarithmic binning
2280 // the user has to delete the array afterwards!!!
2281 //
2282
2283 //check limits
2284 if (xmin<1e-20 || xmax<1e-20){
2285 AliError("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
2286 return MakeLinBinning(nbinsX, xmin, xmax);
2287 }
2288 if (xmax<xmin){
2289 Double_t tmp=xmin;
2290 xmin=xmax;
2291 xmax=tmp;
2292 }
2293 TVectorD *binLim=new TVectorD(nbinsX+1);
2294 Double_t first=xmin;
2295 Double_t last=xmax;
2296 Double_t expMax=TMath::Log(last/first);
2297 for (Int_t i=0; i<nbinsX+1; ++i){
2298 (*binLim)[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
2299 }
2300 return binLim;
2301}
2302
2303//______________________________________________________________________________
2304TVectorD* AliAnalysisTaskPIDqa::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
2305{
2306 //
2307 // Make linear binning
2308 // the user has to delete the array afterwards!!!
2309 //
2310 if (xmax<xmin){
2311 Double_t tmp=xmin;
2312 xmin=xmax;
2313 xmax=tmp;
2314 }
2315 TVectorD *binLim=new TVectorD(nbinsX+1);
2316 Double_t first=xmin;
2317 Double_t last=xmax;
2318 Double_t binWidth=(last-first)/nbinsX;
2319 for (Int_t i=0; i<nbinsX+1; ++i){
2320 (*binLim)[i]=first+binWidth*(Double_t)i;
2321 }
2322 return binLim;
2323}
2324
2325//_____________________________________________________________________________
2326TVectorD* AliAnalysisTaskPIDqa::MakeArbitraryBinning(const char* bins)
2327{
2328 //
2329 // Make arbitrary binning, bins separated by a ','
2330 //
2331 TString limits(bins);
2332 if (limits.IsNull()){
2333 AliError("Bin Limit string is empty, cannot add the variable");
2334 return 0x0;
2335 }
2336
2337 TObjArray *arr=limits.Tokenize(",");
2338 Int_t nLimits=arr->GetEntries();
2339 if (nLimits<2){
2340 AliError("Need at leas 2 bin limits, cannot add the variable");
2341 delete arr;
2342 return 0x0;
2343 }
2344
2345 TVectorD *binLimits=new TVectorD(nLimits);
2346 for (Int_t iLim=0; iLim<nLimits; ++iLim){
2347 (*binLimits)[iLim]=(static_cast<TObjString*>(arr->At(iLim)))->GetString().Atof();
2348 }
2349
2350 delete arr;
2351 return binLimits;
2352}
2353