]> git.uio.no Git - u/mrichter/AliRoot.git/blame - ANALYSIS/AliAnalysisTaskPIDqa.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
ATO-80
[u/mrichter/AliRoot.git] / ANALYSIS / AliAnalysisTaskPIDqa.cxx
CommitLineData
5a9dc560 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>
ce33212f 40#include <AliTPCdEdxInfo.h>
5a9dc560 41
42#include <AliESDEvent.h>
43#include <AliAODEvent.h>
44#include <AliESDv0.h>
45#include <AliAODv0.h>
46#include <AliESDv0KineCuts.h>
ce33212f 47#include <AliESDtrackCuts.h>
48
49#include <AliMCEvent.h>
5a9dc560 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),
ce33212f 70fListQAtpcBasic(0x0),
71fListQAtpcMCtruth(0x0),
72fListQAtpcHybrid(0x0),
73fListQAtpcOROChigh(0x0),
74fListQAtpcV0(0x0),
5a9dc560 75fListQAtrd(0x0),
6f3523b5 76fListQAtrdNsig(0x0),
c2f39c0f 77fListQAtrdNsigTPCTOF(0x0),
5a9dc560 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),
ce33212f 106fListQAtpcBasic(0x0),
107fListQAtpcMCtruth(0x0),
108fListQAtpcHybrid(0x0),
109fListQAtpcOROChigh(0x0),
110fListQAtpcV0(0x0),
5a9dc560 111fListQAtrd(0x0),
6f3523b5 112fListQAtrdNsig(0x0),
c2f39c0f 113fListQAtrdNsigTPCTOF(0x0),
5a9dc560 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");
ce33212f 192
5a9dc560 193 fListQAtrd=new TList;
194 fListQAtrd->SetOwner();
195 fListQAtrd->SetName("TRD");
6f3523b5 196
197 fListQAtrdNsig=new TList;
198 fListQAtrdNsig->SetOwner();
199 fListQAtrdNsig->SetName("TRDnSigma");
5a9dc560 200
c2f39c0f 201 fListQAtrdNsigTPCTOF=new TList;
202 fListQAtrdNsigTPCTOF->SetOwner();
203 fListQAtrdNsigTPCTOF->SetName("TRDnSigma_TPCTOF");
204
5a9dc560 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();
ce33212f 252// SetupTPCqa(kFALSE, kTRUE, kFALSE);
5a9dc560 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
ce33212f 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
5a9dc560 759//______________________________________________________________________________
760void AliAnalysisTaskPIDqa::FillTPCqa()
761{
762 //
763 // Fill PID qa histograms for the TPC
764 //
ce33212f 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
5a9dc560 779 AliVEvent *event=InputEvent();
ce33212f 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 mult=fAODevent->GetHeader()->GetTPConlyRefMultiplicity();
838 }
839
840 /*if (mult < 0) {
841 printf("Reference multiplicity not available \n");
842 //return;
843 }*/
844
845 //printf("The multiplicity is = %i ",mult);
846
847
848 // -- Begin: track loop --
5a9dc560 849 Int_t ntracks=event->GetNumberOfTracks();
850 for(Int_t itrack = 0; itrack < ntracks; itrack++){
851 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
ce33212f 852
5a9dc560 853 //
854 //basic track cuts
855 //
856 ULong_t status=track->GetStatus();
857 // not that nice. status bits not in virtual interface
858 // TPC refit + ITS refit + TPC pid
859 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
860 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
861
862 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
863 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
864 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
865 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
866 if (track->GetTPCNclsF()>0) {
867 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
868 }
ce33212f 869
5a9dc560 870 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
ce33212f 871
872 eta=track->Eta();
873 if ( TMath::Abs(eta)>0.9 ) continue;
874
875 //nacc++; // counter for accepted multiplicity
876
877 // the default ("basic") scenario
878 if (scBasic == 1) {
879 FillTPCHistogramsNsigma(fListQAtpcBasic,0,track,mult);
880 FillTPCHistogramsSignal(fListQAtpcBasic,0,track,mult);
5a9dc560 881 }
5a9dc560 882
ce33212f 883 // only MC truth identified particles
884 if (scMCtruth == 1) {
885 FillTPCHistogramsNsigma(fListQAtpcMCtruth,1,track,mult);
5a9dc560 886 }
5a9dc560 887
ce33212f 888 // the "hybrid" scenario (only for LHC11h)
889 if (scHybrid == 1) {
890 FillTPCHistogramsNsigma(fListQAtpcHybrid,2,track,mult);
891 }
892
893 // the "OROC high" scenario (only for LHC11h)
894 if (scOROChigh == 1) {
895 FillTPCHistogramsNsigma(fListQAtpcOROChigh,3,track,mult);
5a9dc560 896 }
ce33212f 897
898 } // -- End: track loop --
899
900
901 // -- Begin: track loops for V0 candidates --
902 if (scV0 == 1) {
903
904 // - Begin: track loop for electrons from V0 -
905 for(Int_t itrack = 0; itrack < fV0electrons->GetEntries(); itrack++){
906 AliVTrack *track=(AliVTrack*)fV0electrons->At(itrack);
907
908 //
909 //basic track cuts
910 //
911 ULong_t status=track->GetStatus();
912 // not that nice. status bits not in virtual interface
913 // TPC refit + ITS refit + TPC pid
914 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
915 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
916
917 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
918 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
919 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
920 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
921 if (track->GetTPCNclsF()>0) {
922 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
923 }
924
925 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
926
927 eta=track->Eta();
928 if ( TMath::Abs(eta)>0.9 ) continue;
929
930 // fill histograms for V0 candidates
931 FillTPCHistogramsNsigma(fListQAtpcV0,40,track,mult);
932
933 } // - End: track loop for electrons from V0 -
934
935
936 // - Begin: track loop for pions from V0 -
937 for(Int_t itrack = 0; itrack < fV0pions->GetEntries(); itrack++){
938 AliVTrack *track=(AliVTrack*)fV0pions->At(itrack);
939
940 //
941 //basic track cuts
942 //
943 ULong_t status=track->GetStatus();
944 // not that nice. status bits not in virtual interface
945 // TPC refit + ITS refit + TPC pid
946 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
947 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
948
949 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
950 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
951 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
952 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
953 if (track->GetTPCNclsF()>0) {
954 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
955 }
956
957 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
958
959 eta=track->Eta();
960 if ( TMath::Abs(eta)>0.9 ) continue;
961
962 // fill histograms for V0 candidates
963 FillTPCHistogramsNsigma(fListQAtpcV0,42,track,mult);
964
965 } // - End: track loop for pions from V0 -
966
967
968 // - Begin: track loop for kaons from V0 -
969 for(Int_t itrack = 0; itrack < fV0kaons->GetEntries(); itrack++){
970 AliVTrack *track=(AliVTrack*)fV0kaons->At(itrack);
971
972 //
973 //basic track cuts
974 //
975 ULong_t status=track->GetStatus();
976 // not that nice. status bits not in virtual interface
977 // TPC refit + ITS refit + TPC pid
978 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
979 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
980
981 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
982 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
983 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
984 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
985 if (track->GetTPCNclsF()>0) {
986 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
987 }
988
989 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
990
991 eta=track->Eta();
992 if ( TMath::Abs(eta)>0.9 ) continue;
993
994 // fill histograms for V0 candidates
995 FillTPCHistogramsNsigma(fListQAtpcV0,43,track,mult);
996
997 } // - End: track loop for kaons from V0 -
998
999
1000 // - Begin: track loop for protons from V0 -
1001 for(Int_t itrack = 0; itrack < fV0protons->GetEntries(); itrack++){
1002 AliVTrack *track=(AliVTrack*)fV0protons->At(itrack);
1003
1004 //
1005 //basic track cuts
1006 //
1007 ULong_t status=track->GetStatus();
1008 // not that nice. status bits not in virtual interface
1009 // TPC refit + ITS refit + TPC pid
1010 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1011 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
1012
1013 // The TPC pid cut removes the light nuclei (>5 sigma from proton line)
1014 //|| !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid )
1015 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1016 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1017 if (track->GetTPCNclsF()>0) {
1018 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1019 }
1020
1021 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1022
1023 eta=track->Eta();
1024 if ( TMath::Abs(eta)>0.9 ) continue;
1025
1026 // fill histograms for V0 candidates
1027 FillTPCHistogramsNsigma(fListQAtpcV0,44,track,mult);
1028
1029 } // - End: track loop for protons from V0 -
1030
1031 } // -- End: track loops for V0 candidates --
1032
1033
1034 // Multiplicity distribution
1035 TH1 *hm=(TH1*)fListQAtpc->At(scCounter);
1036 if (hm) {
1037 hm->Fill(mult);
5a9dc560 1038 }
ce33212f 1039
1040 //printf("\nAccepted multiplicity = %i \n --- END of event --- \n",nacc);
1041
5a9dc560 1042}
1043
1044//______________________________________________________________________________
1045void AliAnalysisTaskPIDqa::FillTRDqa()
1046{
1047 //
1048 // Fill PID qa histograms for the TRD
1049 //
1050 AliVEvent *event=InputEvent();
1051 Int_t ntracks = event->GetNumberOfTracks();
1052 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1053 AliVTrack *track = (AliVTrack *)event->GetTrack(itrack);
1054
1055 //
1056 //basic track cuts
1057 //
1058 ULong_t status=track->GetStatus();
1059 // not that nice. status bits not in virtual interface
1060 // TPC refit + ITS refit + TPC pid + TRD out
1061 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1062 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1063// !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) || //removes light nuclei. So it is out for the moment
1064 !( (status & AliVTrack::kTRDout ) == AliVTrack::kTRDout )) continue;
1065
1066 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1067 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1068 if (track->GetTPCNclsF()>0) {
1069 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1070 }
1071
1072 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1073
1074 Double_t likelihoods[AliPID::kSPECIES];
1075 if(fPIDResponse->ComputeTRDProbability(track, AliPID::kSPECIES, likelihoods) != AliPIDResponse::kDetPidOk) continue;
1076 Int_t ntracklets = 0;
1077 Double_t momentum = -1.;
6f3523b5 1078 for(Int_t itl = 0; itl < 6; itl++) {
1079 if(track->GetTRDmomentum(itl) > 0.) {
5a9dc560 1080 ntracklets++;
1081 if(momentum < 0) momentum = track->GetTRDmomentum(itl);
6f3523b5 1082 }
1083 }
1084
5a9dc560 1085 for(Int_t ispecie = 0; ispecie < AliPID::kSPECIES; ispecie++){
1086 TH2F *hLike = (TH2F *)fListQAtrd->At(ntracklets*AliPID::kSPECIES+ispecie);
1087 if (hLike) hLike->Fill(momentum,likelihoods[ispecie]);
1088 }
6f3523b5 1089
1090 //=== nSigma and signal ===
1091 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1092 TH2 *h=(TH2*)fListQAtrdNsig->At(ispecie);
c2f39c0f 1093 TH2 *hTPCTOF=(TH2*)fListQAtrdNsigTPCTOF->At(ispecie);
1094 if (!h || !hTPCTOF) continue;
1095 Float_t nSigmaTPC=fPIDResponse->NumberOfSigmas(AliPIDResponse::kTPC, track, (AliPID::EParticleType)ispecie);
1096 Float_t nSigmaTRD=fPIDResponse->NumberOfSigmas(AliPIDResponse::kTRD, track, (AliPID::EParticleType)ispecie);
1097 Float_t nSigmaTOF=fPIDResponse->NumberOfSigmas(AliPIDResponse::kTOF, track, (AliPID::EParticleType)ispecie);
1098 h->Fill(momentum,nSigmaTRD);
1099
1100 if (TMath::Abs(nSigmaTPC)<3 && TMath::Abs(nSigmaTOF)<3) {
1101 hTPCTOF->Fill(momentum,nSigmaTRD);
1102 }
6f3523b5 1103 }
1104
1105 TH2 *h=(TH2*)fListQAtrdNsig->Last();
1106
1107 if (h) {
1108 Double_t sig=track->GetTRDsignal();
1109 h->Fill(momentum,sig);
1110 }
1111
5a9dc560 1112 }
1113}
1114
1115//______________________________________________________________________________
1116void AliAnalysisTaskPIDqa::FillTOFqa()
1117{
1118 //
1119 // Fill TOF information
1120 //
1121 AliVEvent *event=InputEvent();
1122
1123 Int_t ntracks=event->GetNumberOfTracks();
1124 Int_t tracksAtTof = 0;
1125 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1126 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1127
1128 //
1129 //basic track cuts
1130 //
1131 ULong_t status=track->GetStatus();
1132 // TPC refit + ITS refit +
1133 // TOF out + kTIME
1134 // kTIME
1135 // (we don't use kTOFmismatch because it depends on TPC and kTOFpid because it prevents light nuclei
1136 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1137 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1138 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
1139 // !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
1140 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
1141
1142 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1143 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1144 if (track->GetTPCNclsF()>0) {
1145 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1146 }
1147
1148 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1149
1150 tracksAtTof++;
1151
1152 Double_t mom=track->P();
1153
1154 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1155 TH2 *h=(TH2*)fListQAtof->At(ispecie);
1156 if (!h) continue;
1157 Double_t nSigma=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
1158 h->Fill(mom,nSigma);
1159 }
1160
1161 TH2 *h=(TH2*)fListQAtof->FindObject("hSigP_TOF");
1162 if (h) {
1163 Double_t sig=track->GetTOFsignal()/1000.;
1164 h->Fill(mom,sig);
1165 }
1166
1167 Int_t mask = fPIDResponse->GetTOFResponse().GetStartTimeMask(mom);
1168 ((TH1F*)fListQAtof->FindObject("hStartTimeMask_TOF"))->Fill((Double_t)(mask+0.5));
1169
1170 if (mom >= 0.75 && mom <= 1.25 ) {
1171 Double_t nsigma= fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)AliPID::kPion);
1172 if (mask == 0) {
1173 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-Fill"))->Fill(nsigma);
1174 } else if (mask == 1) {
1175 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-TOF"))->Fill(nsigma);
1176 } else if ( (mask == 2) || (mask == 4) || (mask == 6) ) {
1177 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-T0"))->Fill(nsigma);
1178 } else {
1179 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Pion_T0-Best"))->Fill(nsigma);
1180 }
1181 if (mask & 0x1) { //at least TOF-T0 present
1182 Double_t delta=0;
1183 (void)fPIDResponse->GetSignalDelta((AliPIDResponse::EDetector)AliPIDResponse::kTOF,track,(AliPID::EParticleType)AliPID::kPion,delta);
1184 ((TH1F*)fListQAtof->FindObject("hDelta_TOF_Pion"))->Fill(delta);
1185 }
1186 }
1187
1188 Double_t res = (Double_t)fPIDResponse->GetTOFResponse().GetStartTimeRes(mom);
1189 ((TH1F*)fListQAtof->FindObject("hStartTimeRes_TOF"))->Fill(res);
1190
1191 Double_t startTimeT0 = event->GetT0TOF(0);
1192 if (startTimeT0 < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeAC_T0"))->Fill(startTimeT0);
1193 else {
1194 startTimeT0 = event->GetT0TOF(1);
1195 if (startTimeT0 < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeA_T0"))->Fill(startTimeT0);
1196 startTimeT0 = event->GetT0TOF(2);
1197 if (startTimeT0 < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeC_T0"))->Fill(startTimeT0);
1198 }
1199 }
1200 if (tracksAtTof > 0) {
1201 ((TH1F* )fListQAtof->FindObject("hnTracksAt_TOF"))->Fill(tracksAtTof);
1202 Int_t mask = fPIDResponse->GetTOFResponse().GetStartTimeMask(5.);
1203 if (mask & 0x1) ((TH1F*)fListQAtof->FindObject("hT0MakerEff"))->Fill(tracksAtTof);
1204 }
1205}
1206
1207//______________________________________________________________________________
1208void AliAnalysisTaskPIDqa::FillT0qa()
1209{
1210 //
1211 // Fill TOF information
1212 //
1213 AliVEvent *event=InputEvent();
1214
1215 Int_t ntracks=event->GetNumberOfTracks();
1216
1217 Int_t tracksAtT0 = 0;
1218
1219 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1220 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1221
1222 //
1223 //basic track cuts
1224 //
1225 ULong_t status=track->GetStatus();
1226 // TPC refit + ITS refit +
1227 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1228 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
1229 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1230 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1231 if (track->GetTPCNclsF()>0) {
1232 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1233 }
1234 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1235
1236 tracksAtT0++;
1237 }
1238
1239 Bool_t t0A = kFALSE;
1240 Bool_t t0C = kFALSE;
1241 Bool_t t0And = kFALSE;
1242 Double_t startTimeT0 = event->GetT0TOF(0); // AND
1243 if (startTimeT0 < 90000) {
1244 t0And = kTRUE;
1245 ((TH1F*)fListQAt0->FindObject("hStartTimeAC_T0"))->Fill(startTimeT0);
1246 }
1247 startTimeT0 = event->GetT0TOF(1); // T0A
1248 if (startTimeT0 < 90000) {
1249 t0A = kTRUE;
1250 ((TH1F*)fListQAt0->FindObject("hStartTimeA_T0"))->Fill(startTimeT0);
1251
1252 }
1253 startTimeT0 = event->GetT0TOF(2); // T0C
1254 if (startTimeT0 < 90000) {
1255 t0C = kTRUE;
1256 ((TH1F*)fListQAt0->FindObject("hStartTimeC_T0"))->Fill(startTimeT0);
1257 }
1258
1259 ((TH1F* )fListQAt0->FindObject("hnTracksAt_T0"))->Fill(tracksAtT0);
1260 if (t0A) ((TH1F*)fListQAt0->FindObject("hT0AEff"))->Fill(tracksAtT0);
1261 if (t0C) ((TH1F*)fListQAt0->FindObject("hT0CEff"))->Fill(tracksAtT0);
1262 if (t0And) ((TH1F*)fListQAt0->FindObject("hT0AndEff"))->Fill(tracksAtT0);
1263 if (t0A || t0C) ((TH1F*)fListQAt0->FindObject("hT0OrEff"))->Fill(tracksAtT0);
1264}
1265
1266
1267//______________________________________________________________________________
1268void AliAnalysisTaskPIDqa::FillEMCALqa()
1269{
1270 //
1271 // Fill PID qa histograms for the EMCAL
1272 //
1273
1274 AliVEvent *event=InputEvent();
1275
1276 Int_t ntracks=event->GetNumberOfTracks();
1277 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1278 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1279
1280 //
1281 //basic track cuts
1282 //
1283 ULong_t status=track->GetStatus();
1284 // not that nice. status bits not in virtual interface
1285 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1286
1287 Double_t pt=track->Pt();
1288
1289 //EMCAL nSigma (only for electrons at the moment)
1290 TH2 *h=(TH2*)fListQAemcal->At(0);
1291 if (!h) continue;
1292 Double_t nSigma=fPIDResponse->NumberOfSigmasEMCAL(track, (AliPID::EParticleType)0);
1293 h->Fill(pt,nSigma);
1294
1295 }
1296
1297 //EMCAL signal (E/p vs. pT) for electrons from V0
1298 for(Int_t itrack = 0; itrack < fV0electrons->GetEntries(); itrack++){
1299 AliVTrack *track=(AliVTrack*)fV0electrons->At(itrack);
1300
1301 //
1302 //basic track cuts
1303 //
1304 ULong_t status=track->GetStatus();
1305 // not that nice. status bits not in virtual interface
1306 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1307
1308 Double_t pt=track->Pt();
1309
1310 TH2 *h=(TH2*)fListQAemcal->At(1);
1311 if (h) {
1312
1313 Int_t nMatchClus = track->GetEMCALcluster();
1314 Double_t mom = track->P();
1315 Double_t eop = -1.;
1316
1317 if(nMatchClus > -1){
1318
1319 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1320
1321 if(matchedClus){
1322
1323 // matched cluster is EMCAL
1324 if(matchedClus->IsEMCAL()){
1325
1326 Double_t fClsE = matchedClus->E();
1327 eop = fClsE/mom;
1328
1329 h->Fill(pt,eop);
1330
1331 }
1332 }
1333 }
1334 }
1335 }
1336
1337 //EMCAL signal (E/p vs. pT) for pions from V0
1338 for(Int_t itrack = 0; itrack < fV0pions->GetEntries(); itrack++){
1339 AliVTrack *track=(AliVTrack*)fV0pions->At(itrack);
1340
1341 //
1342 //basic track cuts
1343 //
1344 ULong_t status=track->GetStatus();
1345 // not that nice. status bits not in virtual interface
1346 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1347
1348 Double_t pt=track->Pt();
1349
1350 TH2 *h=(TH2*)fListQAemcal->At(2);
1351 if (h) {
1352
1353 Int_t nMatchClus = track->GetEMCALcluster();
1354 Double_t mom = track->P();
1355 Double_t eop = -1.;
1356
1357 if(nMatchClus > -1){
1358
1359 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1360
1361 if(matchedClus){
1362
1363 // matched cluster is EMCAL
1364 if(matchedClus->IsEMCAL()){
1365
1366 Double_t fClsE = matchedClus->E();
1367 eop = fClsE/mom;
1368
1369 h->Fill(pt,eop);
1370
1371 }
1372 }
1373 }
1374 }
1375 }
1376
1377 //EMCAL signal (E/p vs. pT) for protons from V0
1378 for(Int_t itrack = 0; itrack < fV0protons->GetEntries(); itrack++){
1379 AliVTrack *track=(AliVTrack*)fV0protons->At(itrack);
1380
1381 //
1382 //basic track cuts
1383 //
1384 ULong_t status=track->GetStatus();
1385 // not that nice. status bits not in virtual interface
1386 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
1387
1388 Double_t pt=track->Pt();
1389
1390 TH2 *hP=(TH2*)fListQAemcal->At(3);
1391 TH2 *hAP=(TH2*)fListQAemcal->At(4);
1392 if (hP && hAP) {
1393
1394 Int_t nMatchClus = track->GetEMCALcluster();
1395 Double_t mom = track->P();
1396 Int_t charge = track->Charge();
1397 Double_t eop = -1.;
1398
1399 if(nMatchClus > -1){
1400
1401 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1402
1403 if(matchedClus){
1404
1405 // matched cluster is EMCAL
1406 if(matchedClus->IsEMCAL()){
1407
1408 Double_t fClsE = matchedClus->E();
1409 eop = fClsE/mom;
1410
1411 if(charge > 0) hP->Fill(pt,eop);
1412 else if(charge < 0) hAP->Fill(pt,eop);
1413
1414 }
1415 }
1416 }
1417 }
1418 }
1419
1420}
1421
1422
1423//______________________________________________________________________________
1424void AliAnalysisTaskPIDqa::FillHMPIDqa()
1425{
1426 //
1427 // Fill PID qa histograms for the HMPID
1428 //
1429
1430 AliVEvent *event=InputEvent();
1431
1432 Int_t ntracks=event->GetNumberOfTracks();
1433 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1434 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1435
1436 //
1437 //basic track cuts
1438 //
1439 const ULong_t status=track->GetStatus();
1440 // not that nice. status bits not in virtual interface
1441 // TPC refit + ITS refit +
1442 // TOF out + TOFpid +
1443 // kTIME
1444 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1445 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ) continue;
1446
1447 const Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1448 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1449 if (track->GetTPCNclsF()>0) {
1450 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1451 }
1452
1453 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1454
1455 const Double_t mom = track->P();
1456 const Double_t ckovAngle = track->GetHMPIDsignal();
1457
1458 Int_t nhists=0;
1459 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1460 if (ispecie==AliPID::kElectron || ispecie==AliPID::kMuon) continue;
1461 TH2 *h=(TH2*)fListQAhmpid->At(nhists);
1462 if (!h) {++nhists; continue;}
1463 const Double_t nSigma=fPIDResponse->NumberOfSigmasHMPID(track, (AliPID::EParticleType)ispecie);
1464 h->Fill(mom,nSigma);
1465 ++nhists;
1466 }
1467
1468 TH1F *hThetavsMom = (TH1F*)fListQAhmpid->At(AliPID::kSPECIESC);
1469
1470 if (hThetavsMom) hThetavsMom->Fill(mom,ckovAngle);
1471
1472 }
1473}
1474//______________________________________________________________________________
1475void AliAnalysisTaskPIDqa::FillTOFHMPIDqa()
1476{
1477 //
1478 // Fill PID qa histograms for the HMPID
1479 //
1480
1481 AliVEvent *event=InputEvent();
1482
1483 Int_t ntracks=event->GetNumberOfTracks();
1484 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1485 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1486
1487 //
1488 //basic track cuts
1489 //
1490 ULong_t status=track->GetStatus();
1491 // not that nice. status bits not in virtual interface
1492 // TPC refit + ITS refit +
1493 // TOF out + TOFpid +
1494 // kTIME
1495 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1496 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1497 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
1498 !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
1499 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
1500
1501 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1502 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1503 if (track->GetTPCNclsF()>0) {
1504 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1505 }
1506
1507 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1508
1509 Double_t mom = track->P();
1510 Double_t ckovAngle = track->GetHMPIDsignal();
1511
1512 Double_t nSigmaTOF[3];
1513 TH1F *h[3];
1514
1515 for (Int_t ispecie=2; ispecie<5; ++ispecie){
1516 //TOF nSigma
1517 nSigmaTOF[ispecie-2]=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
1518 h[ispecie-2] = (TH1F*)fListQAtofhmpid->At(ispecie-2);}
1519
1520 if(TMath::Abs(nSigmaTOF[0])<2) h[0]->Fill(mom,ckovAngle);
1521
1522 if(TMath::Abs(nSigmaTOF[1])<2 && TMath::Abs(nSigmaTOF[0])>3) h[1]->Fill(mom,ckovAngle);
1523
1524 if(TMath::Abs(nSigmaTOF[2])<2 && TMath::Abs(nSigmaTOF[1])>3 && TMath::Abs(nSigmaTOF[0])>3) h[2]->Fill(mom,ckovAngle);
1525
1526 }
1527
1528}
1529
1530//______________________________________________________________________________
1531void AliAnalysisTaskPIDqa::FillTPCTOFqa()
1532{
1533 //
1534 // Fill PID qa histograms for the TOF
1535 // Here also the TPC histograms after TOF selection are filled
1536 //
1537
1538 AliVEvent *event=InputEvent();
1539
1540 Int_t ntracks=event->GetNumberOfTracks();
1541 for(Int_t itrack = 0; itrack < ntracks; itrack++){
1542 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
1543
1544 //
1545 //basic track cuts
1546 //
1547 ULong_t status=track->GetStatus();
1548 // not that nice. status bits not in virtual interface
1549 // TPC refit + ITS refit +
1550 // TOF out + TOFpid +
1551 // kTIME
1552 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
1553 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
1554// !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) || //removes light nuclei, so it is out for the moment
1555 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
1556 !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
1557 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
1558
1559 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
1560 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
1561 if (track->GetTPCNclsF()>0) {
1562 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
1563 }
1564
1565 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
1566
1567
1568 Double_t mom=track->P();
1569 Double_t momTPC=track->GetTPCmomentum();
1570
1571 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1572 //TOF nSigma
1573 Double_t nSigmaTOF=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
1574 Double_t nSigmaTPC=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
1575
1576 //TPC after TOF cut
1577 TH2 *h=(TH2*)fListQAtpctof->At(ispecie);
1578 if (h && TMath::Abs(nSigmaTOF)<3.) h->Fill(momTPC,nSigmaTPC);
1579
1580 //TOF after TPC cut
1581 h=(TH2*)fListQAtpctof->At(ispecie+AliPID::kSPECIESC);
1582 if (h && TMath::Abs(nSigmaTPC)<3.) h->Fill(mom,nSigmaTOF);
1583
1584 //EMCAL after TOF and TPC cut
1585 h=(TH2*)fListQAtpctof->At(ispecie+2*AliPID::kSPECIESC);
1586 if (h && TMath::Abs(nSigmaTOF)<3. && TMath::Abs(nSigmaTPC)<3. ){
1587
1588 Int_t nMatchClus = track->GetEMCALcluster();
1589 Double_t pt = track->Pt();
1590 Double_t eop = -1.;
1591
1592 if(nMatchClus > -1){
1593
1594 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
1595
1596 if(matchedClus){
1597
1598 // matched cluster is EMCAL
1599 if(matchedClus->IsEMCAL()){
1600
1601 Double_t fClsE = matchedClus->E();
1602 eop = fClsE/mom;
1603
1604 h->Fill(pt,eop);
1605
1606
1607 }
1608 }
1609 }
1610 }
1611 }
1612 }
1613}
1614
1615//_____________________________________________________________________________
1616void AliAnalysisTaskPIDqa::FillQAinfo()
1617{
1618 //
1619 // Fill the QA information
1620 //
1621
1622
1623 //TPC QA info
1624 TObjArray *arrTPC=static_cast<TObjArray*>(fListQAinfo->At(0));
1625 if (fPIDResponse && arrTPC){
1626 AliTPCPIDResponse &tpcResp=fPIDResponse->GetTPCResponse();
1627 // fill spline names
1628 if (!arrTPC->UncheckedAt(0)){
1629
1630 TObjArray *arrTPCsplineNames=new TObjArray(AliPID::kSPECIESC);
1631 arrTPCsplineNames->SetOwner();
1632 arrTPCsplineNames->SetName("TPC_spline_names");
1633 arrTPC->AddAt(arrTPCsplineNames,0);
1634
1635 for (Int_t iresp=0; iresp<AliPID::kSPECIESC; ++iresp){
1636 const TObject *o=tpcResp.GetResponseFunction((AliPID::EParticleType)iresp);
1637 if (!o) continue;
1638 arrTPCsplineNames->Add(new TObjString(Form("%02d: %s",iresp, o->GetName())));
1639 }
1640 }
1641
1642 // tpc response config
1643 if (!arrTPC->UncheckedAt(1)){
1644
1645 TObjArray *arrTPCconfigInfo=new TObjArray;
1646 arrTPCconfigInfo->SetOwner();
1647 arrTPCconfigInfo->SetName("TPC_config_info");
1648 arrTPC->AddAt(arrTPCconfigInfo,1);
1649
1650 TObjString *ostr=0x0;
1651 ostr=new TObjString;
1652 ostr->String().Form("Eta Corr map: %s", tpcResp.GetEtaCorrMap()?tpcResp.GetEtaCorrMap()->GetName():"none");
1653 arrTPCconfigInfo->Add(ostr);
1654
1655 ostr=new TObjString;
1656 ostr->String().Form("Sigma Par map: %s", tpcResp.GetSigmaPar1Map()?tpcResp.GetSigmaPar1Map()->GetName():"none");
1657 arrTPCconfigInfo->Add(ostr);
1658
1659 ostr=new TObjString;
1660 ostr->String().Form("MIP: %.2f", tpcResp.GetMIP());
1661 arrTPCconfigInfo->Add(ostr);
1662
1663 ostr=new TObjString;
1664 ostr->String().Form("Res: Def %.3g (%.3g) : AllHigh %.3g (%.3g) : OROC high %.3g (%.3g)",
1665 tpcResp.GetRes0(AliTPCPIDResponse::kDefault), tpcResp.GetResN2(AliTPCPIDResponse::kDefault),
1666 tpcResp.GetRes0(AliTPCPIDResponse::kALLhigh), tpcResp.GetResN2(AliTPCPIDResponse::kALLhigh),
1667 tpcResp.GetRes0(AliTPCPIDResponse::kOROChigh), tpcResp.GetResN2(AliTPCPIDResponse::kOROChigh)
1668 );
1669 arrTPCconfigInfo->Add(ostr);
1670 }
1671 }
1672}
1673
1674//______________________________________________________________________________
1675void AliAnalysisTaskPIDqa::SetupITSqa()
1676{
1677 //
1678 // Create the ITS qa objects
1679 //
1680
1681 TVectorD *vX=MakeLogBinning(200,.1,30);
1682
1683 //ITS+TPC tracks
1684 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1685 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_ITS_%s",AliPID::ParticleName(ispecie)),
1686 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
1687 vX->GetNrows()-1,vX->GetMatrixArray(),
1688 200,-10,10);
1689 fListQAits->Add(hNsigmaP);
1690 }
1691 TH2F *hSig = new TH2F("hSigP_ITS",
1692 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
1693 vX->GetNrows()-1,vX->GetMatrixArray(),
1694 300,0,300);
1695 fListQAits->Add(hSig);
1696
1697 //ITS Standalone tracks
1698 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1699 TH2F *hNsigmaPSA = new TH2F(Form("hNsigmaP_ITSSA_%s",AliPID::ParticleName(ispecie)),
1700 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
1701 vX->GetNrows()-1,vX->GetMatrixArray(),
1702 200,-10,10);
1703 fListQAitsSA->Add(hNsigmaPSA);
1704 }
1705 TH2F *hSigSA = new TH2F("hSigP_ITSSA",
1706 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
1707 vX->GetNrows()-1,vX->GetMatrixArray(),
1708 300,0,300);
1709 fListQAitsSA->Add(hSigSA);
1710
1711 //ITS Pure Standalone tracks
1712 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
1713 TH2F *hNsigmaPPureSA = new TH2F(Form("hNsigmaP_ITSPureSA_%s",AliPID::ParticleName(ispecie)),
1714 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
1715 vX->GetNrows()-1,vX->GetMatrixArray(),
1716 200,-10,10);
1717 fListQAitsPureSA->Add(hNsigmaPPureSA);
1718 }
1719 TH2F *hSigPureSA = new TH2F("hSigP_ITSPureSA",
1720 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
1721 vX->GetNrows()-1,vX->GetMatrixArray(),
1722 300,0,300);
1723 fListQAitsPureSA->Add(hSigPureSA);
1724
1725 delete vX;
1726}
1727
ce33212f 1728//_____________________________________________________________________________
1729void AliAnalysisTaskPIDqa::AddTPCHistogramsSignal(TList *sublist, const char *scenario)
5a9dc560 1730{
1731 //
ce33212f 1732 // Create the TPC qa objects: create histograms for the TPC signal for different settings
5a9dc560 1733 //
ce33212f 1734
5a9dc560 1735 TVectorD *vX=MakeLogBinning(200,.1,30);
ce33212f 1736 Int_t nBinsMult = 38;
1737 Double_t xBinsMult[39] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100,
1738 120, 140, 160, 180, 200,
1739 300, 400, 500, 600, 700, 800, 900, 1000,
1740 1200, 1400, 1600, 1800, 2000,
1741 2200, 2400, 2600, 2800, 3000,
1742 3200, 3400, 3600, 3800, 4000
1743 };
1744 const Int_t binsEta=110;
1745 Float_t etaMin=-1.1;
1746 Float_t etaMax=1.1;
1747
1748 char signal[4][12]={"std","IROC","OROCmedium","OROClong"};
1749
1750
1751 // TPC signal vs. p for all particles (standard, IROC, OROCmedium, OROClong)
1752 for (Int_t iSig=0; iSig<4; iSig++) {
1753 TH2F *hSigP = new TH2F(Form("hSigP_TPC_%s_%s",signal[iSig],scenario),
1754 Form("TPC_%s signal (%s) vs. p;p [GeV]; TPC signal [arb. units]",scenario,signal[iSig]),
1755 vX->GetNrows()-1,vX->GetMatrixArray(),
1756 300,0,300);
1757 sublist->Add(hSigP);
1758 }
1759
1760 // MIP pions: TPC signal vs. eta
1761 for (Int_t iSig=0; iSig<4; iSig++) {
1762 TH2F *hSigEtaMIPpi = new TH2F(Form("hSigEta_TPC_%s_%s_MIPpi",signal[iSig],scenario),
1763 Form("TPC_%s signal (%s) MIPpi vs. eta;#eta;TPC signal [arb. units]",scenario,signal[iSig]),
1764 binsEta,etaMin,etaMax,
1765 300,0,300);
1766 sublist->Add(hSigEtaMIPpi);
1767 }
1768
1769 // MIP pions: TPC signal vs. multiplicity
1770 for (Int_t iSig=0; iSig<4; iSig++) {
1771 TH2F *hSigMultMPIpi = new TH2F(Form("hSigMult_TPC_%s_%s_MIPpi",signal[iSig],scenario),
1772 Form("TPC_%s signal (%s) MIPpi vs. mult;multiplicity;TPC signal [arb. units]",scenario,signal[iSig]),
1773 nBinsMult,xBinsMult,
1774 300,0,300);
1775 sublist->Add(hSigMultMPIpi);
1776 }
1777
1778 // Electrons: TPC signal vs. eta
1779 for (Int_t iSig=0; iSig<4; iSig++) {
1780 TH2F *hSigEtaEle = new TH2F(Form("hSigEta_TPC_%s_%s_Ele",signal[iSig],scenario),
1781 Form("TPC_%s signal (%s) electrons vs. eta;#eta;TPC signal [arb. units]",scenario,signal[iSig]),
1782 binsEta,etaMin,etaMax,
1783 300,0,300);
1784 sublist->Add(hSigEtaEle);
1785 }
1786
1787 // Electrons: TPC signal vs. multiplicity
1788 for (Int_t iSig=0; iSig<4; iSig++) {
1789 TH2F *hSigMultEle = new TH2F(Form("hSigMult_TPC_%s_%s_Ele",signal[iSig],scenario),
1790 Form("TPC_%s signal (%s) electrons vs. mult;multiplicity;TPC signal [arb. units]",scenario,signal[iSig]),
1791 nBinsMult,xBinsMult,
1792 300,0,300);
1793 sublist->Add(hSigMultEle);
1794 }
1795
1796 delete vX;
1797
1798}
1799
1800//_____________________________________________________________________________
1801void AliAnalysisTaskPIDqa::AddTPCHistogramsNsigma(TList *sublist, const char *scenario, Int_t scnumber)
1802{
1803 //
1804 // Create the TPC qa objects: create histograms for TPC Nsigma for different settings
1805 //
1806
1807 TVectorD *vX=MakeLogBinning(200,.1,30.);
1808 Int_t nBinsMult = 38;
1809 Double_t xBinsMult[39] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100,
1810 120, 140, 160, 180, 200,
1811 300, 400, 500, 600, 700, 800, 900, 1000,
1812 1200, 1400, 1600, 1800, 2000,
1813 2200, 2400, 2600, 2800, 3000,
1814 3200, 3400, 3600, 3800, 4000
1815 };
1816 const Int_t binsEta=110;
1817 Float_t etaMin=-1.1;
1818 Float_t etaMax=1.1;
1819
1820 Int_t nSpecies=0;
5a9dc560 1821
ce33212f 1822 if (scnumber == 4) nSpecies=(Int_t)AliPID::kSPECIES;
1823 else nSpecies=(Int_t)AliPID::kSPECIESC;
1824
1825 // Nsigma vs. p for different particle species
1826 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
1827 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_%s_%s",scenario,AliPID::ParticleName(ispecie)),
1828 Form("TPC_%s n#sigma %s vs. p;p [GeV]; n#sigma",scenario,AliPID::ParticleName(ispecie)),
5a9dc560 1829 vX->GetNrows()-1,vX->GetMatrixArray(),
1830 200,-10,10);
ce33212f 1831 sublist->Add(hNsigmaP);
5a9dc560 1832 }
1833
ce33212f 1834 // Nsigma vs. eta for different particle species (only for some scenarios)
1835 if ( scnumber == 1 || scnumber == 4 ) {
1836 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
1837 TH2F *hNsigmaEta = new TH2F(Form("hNsigmaEta_TPC_%s_%s",scenario,AliPID::ParticleName(ispecie)),
1838 Form("TPC_%s n#sigma %s vs. eta;#eta; n#sigma",scenario,AliPID::ParticleName(ispecie)),
1839 binsEta,etaMin,etaMax,
5a9dc560 1840 200,-10,10);
ce33212f 1841 sublist->Add(hNsigmaEta);
1842 }
5a9dc560 1843 }
ce33212f 1844
1845 // Nsigma vs. multiplicity for different particle species (only for some scenarios)
1846 if ( scnumber == 1 || scnumber == 4 ) {
1847 for (Int_t ispecie=0; ispecie<nSpecies; ++ispecie){
1848 TH2F *hNsigmaMult = new TH2F(Form("hNsigmaMult_TPC_%s_%s",scenario,AliPID::ParticleName(ispecie)),
1849 Form("TPC_%s n#sigma %s vs. mult;multiplicity; n#sigma",scenario,AliPID::ParticleName(ispecie)),
1850 nBinsMult,xBinsMult,
5a9dc560 1851 200,-10,10);
ce33212f 1852 sublist->Add(hNsigmaMult);
1853 }
5a9dc560 1854 }
ce33212f 1855
1856 // - Beginn: Adding histograms for MIP pions and electrons (only for some scenarios) -
1857 if ( scnumber == 0 || scnumber == 2 || scnumber == 3 ) {
1858
1859 // MIP pions: Nsigma vs. eta
1860 TH2F *hNsigmaEtaMIPpi = new TH2F(Form("hNsigmaEta_TPC_%s_MIPpi",scenario),
1861 Form("TPC_%s n#sigma MIPpi vs. eta;#eta; n#sigma",scenario),
1862 binsEta,etaMin,etaMax,
1863 200,-10,10);
1864 sublist->Add(hNsigmaEtaMIPpi);
1865
1866 // MIP pions: Nsigma vs. multiplicity
1867 TH2F *hNsigmaMultMIPpi = new TH2F(Form("hNsigmaMult_TPC_%s_MIPpi",scenario),
1868 Form("TPC_%s n#sigma MIPpi vs. mult;multiplicity; n#sigma",scenario),
1869 nBinsMult,xBinsMult,
1870 200,-10,10);
1871 sublist->Add(hNsigmaMultMIPpi);
1872
1873 // Electrons: Nsigma vs. eta
1874 TH2F *hNsigmaEtaEle = new TH2F(Form("hNsigmaEta_TPC_%s_Ele",scenario),
1875 Form("TPC_%s n#sigma electrons vs. eta;#eta; n#sigma",scenario),
1876 binsEta,etaMin,etaMax,
1877 200,-10,10);
1878 sublist->Add(hNsigmaEtaEle);
1879
1880 // Electrons: Nsigma vs. multiplicity
1881 TH2F *hNsigmaMultEle = new TH2F(Form("hNsigmaMult_TPC_%s_Ele",scenario),
1882 Form("TPC_%s n#sigma electrons vs. mult;multiplicity; n#sigma",scenario),
1883 nBinsMult,xBinsMult,
1884 200,-10,10);
1885 sublist->Add(hNsigmaMultEle);
1886 } // - End: Adding histograms for MIP pions and electrons
1887
1888 delete vX;
1889
1890}
1891
1892//______________________________________________________________________________
1893void AliAnalysisTaskPIDqa::SetupTPCqa(Bool_t fillMC, Bool_t fill11h, Bool_t fillV0)
1894{
1895 //
1896 // Create the TPC qa objects
1897 //
5a9dc560 1898
ce33212f 1899 // Set up the multiplicity binning
1900 Int_t nBinsMult = 38;
1901 Double_t xBinsMult[39] = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100,
1902 120, 140, 160, 180, 200,
1903 300, 400, 500, 600, 700, 800, 900, 1000,
1904 1200, 1400, 1600, 1800, 2000,
1905 2200, 2400, 2600, 2800, 3000,
1906 3200, 3400, 3600, 3800, 4000
1907 };
1908
1909
1910 // Create TPC sublists for different scenarios
1911 // corresponding to available information,
1912 // e.g. MC or not, special settings for LHC11h
1913
1914 // basic/default scenario, used always
1915 fListQAtpcBasic=new TList;
1916 fListQAtpcBasic->SetOwner();
1917 fListQAtpcBasic->SetName("TPCBasic");
1918 fListQAtpc->Add(fListQAtpcBasic);
1919
1920 // MC truth scenario: use only MC truth identified particles
1921 // only available for MC
1922 if (fillMC == kTRUE) {
1923 fListQAtpcMCtruth=new TList;
1924 fListQAtpcMCtruth->SetOwner();
1925 fListQAtpcMCtruth->SetName("TPCMCtruth");
1926 fListQAtpc->Add(fListQAtpcMCtruth);
1927 }
5a9dc560 1928
ce33212f 1929 // Hybrid and OROChigh scenarios,
1930 // special settings only available for PbPb LHC11h data
1931 if (fill11h == kTRUE) {
1932 fListQAtpcHybrid=new TList;
1933 fListQAtpcHybrid->SetOwner();
1934 fListQAtpcHybrid->SetName("TPCHybrid");
1935 fListQAtpc->Add(fListQAtpcHybrid);
5a9dc560 1936
ce33212f 1937 fListQAtpcOROChigh=new TList;
1938 fListQAtpcOROChigh->SetOwner();
1939 fListQAtpcOROChigh->SetName("TPCOROChigh");
1940 fListQAtpc->Add(fListQAtpcOROChigh);
1941 }
1942
1943 // scenario only for V0s,
1944 // only available for ESDs
1945 if (fillV0 == kTRUE) {
1946 fListQAtpcV0=new TList;
1947 fListQAtpcV0->SetOwner();
1948 fListQAtpcV0->SetName("TPCV0");
1949 fListQAtpc->Add(fListQAtpcV0);
1950 }
1951
1952
1953 // the default ("basic") scenario
1954 AddTPCHistogramsNsigma(fListQAtpcBasic,"Basic",0);
1955 AddTPCHistogramsSignal(fListQAtpcBasic,"Basic");
1956
1957 // only MC truth identified particles
1958 if (fillMC) {
1959 AddTPCHistogramsNsigma(fListQAtpcMCtruth,"MCtruth",1);
1960 }
1961
1962 // the "hybrid" scenario (only for period LHC11h)
1963 if (fill11h) {
1964 AddTPCHistogramsNsigma(fListQAtpcHybrid,"Hybrid",2);
1965 }
1966
1967 // the "OROC high" scenario (only for period LHC11h)
1968 if (fill11h) {
1969 AddTPCHistogramsNsigma(fListQAtpcOROChigh,"OROChigh",3);
1970 }
1971
1972 // only for V0s
1973 if (fillV0) {
1974 AddTPCHistogramsNsigma(fListQAtpcV0,"V0",4);
1975 }
1976
1977
1978 // Multiplicity distribution --- as check
1979 TH1F *hMult = new TH1F("hMult_TPC",
1980 "Multiplicity distribution;multiplicity;counts",
1981 nBinsMult,xBinsMult);
1982 fListQAtpc->Add(hMult);
5a9dc560 1983
5a9dc560 1984}
1985
1986//______________________________________________________________________________
1987void AliAnalysisTaskPIDqa::SetupTRDqa()
1988{
1989 //
1990 // Create the TRD qa objects
1991 //
1992 TVectorD *vX=MakeLogBinning(200,.1,30);
1993 for(Int_t itl = 0; itl < 6; ++itl){
1994 for(Int_t ispecie = 0; ispecie < AliPID::kSPECIES; ispecie++){
1995 TH2F *hLikeP = new TH2F(Form("hLikeP_TRD_%dtls_%s", itl, AliPID::ParticleName(ispecie)),
1996 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)),
1997 vX->GetNrows()-1, vX->GetMatrixArray(),
1998 100, 0., 1.);
1999 fListQAtrd->Add(hLikeP);
2000 }
2001 }
6f3523b5 2002
2003 // === nSigma Values and signal ===
2004 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2005 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TRD_%s",AliPID::ParticleName(ispecie)),
2006 Form("TRD n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2007 vX->GetNrows()-1,vX->GetMatrixArray(),
2008 100,-10,10);
2009 fListQAtrdNsig->Add(hNsigmaP);
2010 }
2011
2012 TH2F *hSig = new TH2F("hSigP_TRD",
2013 "TRD signal vs. p;p [GeV]; TRD signal [arb. units]",
2014 vX->GetNrows()-1,vX->GetMatrixArray(),
2015 100,0,100);
2016 fListQAtrdNsig->Add(hSig);
2017
2018 fListQAtrd->Add(fListQAtrdNsig);
c2f39c0f 2019
2020 // === Same after 3 sigma in TPC and TOF
2021 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2022 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TRD_TPCTOF_%s",AliPID::ParticleName(ispecie)),
2023 Form("TRD n#sigma %s vs. p after 3#sigma cut in TPC&TOF;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2024 vX->GetNrows()-1,vX->GetMatrixArray(),
2025 100,-10,10);
2026 fListQAtrdNsigTPCTOF->Add(hNsigmaP);
2027 }
2028
2029 fListQAtrd->Add(fListQAtrdNsigTPCTOF);
6f3523b5 2030
5a9dc560 2031 delete vX;
2032}
2033
2034//______________________________________________________________________________
2035void AliAnalysisTaskPIDqa::SetupTOFqa()
2036{
2037 //
2038 // Create the TOF qa objects
2039 //
2040
2041 TVectorD *vX=MakeLogBinning(200,.1,30);
2042
2043 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2044 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_%s",AliPID::ParticleName(ispecie)),
2045 Form("TOF n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2046 vX->GetNrows()-1,vX->GetMatrixArray(),
2047 200,-10,10);
2048 fListQAtof->Add(hNsigmaP);
2049 }
2050
2051 TH1F *hnSigT0Fill = new TH1F("hNsigma_TOF_Pion_T0-Fill","TOF n#sigma (Pion) T0-FILL [0.75-1.25. GeV/c]",200,-10,10);
2052 fListQAtof->Add(hnSigT0Fill);
2053 TH1F *hnSigT0T0 = new TH1F("hNsigma_TOF_Pion_T0-T0","TOF n#sigma (Pion) T0-T0 [0.75-1.25 GeV/c]",200,-10,10);
2054 fListQAtof->Add(hnSigT0T0);
2055 TH1F *hnSigT0TOF = new TH1F("hNsigma_TOF_Pion_T0-TOF","TOF n#sigma (Pion) T0-TOF [0.75-1.25 GeV/c]",200,-10,10);
2056 fListQAtof->Add(hnSigT0TOF);
2057 TH1F *hnSigT0Best = new TH1F("hNsigma_TOF_Pion_T0-Best","TOF n#sigma (Pion) T0-Best [0.75-1.25 GeV/c]",200,-10,10);
2058 fListQAtof->Add(hnSigT0Best);
2059 TH1F *hnDeltaPi = new TH1F("hDelta_TOF_Pion","DeltaT (Pion) [0.75-1.25 GeV/c]",50,-500,500);
2060 fListQAtof->Add(hnDeltaPi);
2061
2062 TH2F *hSig = new TH2F("hSigP_TOF",
2063 "TOF signal vs. p;p [GeV]; TOF signal [ns]",
2064 vX->GetNrows()-1,vX->GetMatrixArray(),
2065 300,0,30);
2066
2067 delete vX;
2068
2069 fListQAtof->Add(hSig);
2070
2071 TH1F *hStartTimeMaskTOF = new TH1F("hStartTimeMask_TOF","StartTime mask",8,0,8);
2072 fListQAtof->Add(hStartTimeMaskTOF);
2073 TH1F *hStartTimeResTOF = new TH1F("hStartTimeRes_TOF","StartTime resolution [ps]",100,0,500);
2074 fListQAtof->Add(hStartTimeResTOF);
2075
2076 TH1F *hnTracksAtTOF = new TH1F("hnTracksAt_TOF","Matched tracks at TOF",100,0,100);
2077 fListQAtof->Add(hnTracksAtTOF);
2078 TH1F *hT0MakerEff = new TH1F("hT0MakerEff","Events with T0-TOF vs nTracks",100,0,100);
2079 fListQAtof->Add(hT0MakerEff);
2080
2081 // this in principle should stay on a T0 PID QA, but are just the data prepared for TOF use
2082 TH1F *hStartTimeAT0 = new TH1F("hStartTimeA_T0","StartTime from T0A [ps]",1000,-1000,1000);
2083 fListQAtof->Add(hStartTimeAT0);
2084 TH1F *hStartTimeCT0 = new TH1F("hStartTimeC_T0","StartTime from T0C [ps]",1000,-1000,1000);
2085 fListQAtof->Add(hStartTimeCT0);
2086 TH1F *hStartTimeACT0 = new TH1F("hStartTimeAC_T0","StartTime from T0AC [ps]",1000,-1000,1000);;
2087 fListQAtof->Add(hStartTimeACT0);
2088}
2089
2090
2091//______________________________________________________________________________
2092void AliAnalysisTaskPIDqa::SetupT0qa()
2093{
2094 //
2095 // Create the T0 qa objects
2096 //
2097
2098 // these are similar to plots inside TOFqa, but these are for all events
2099 TH1F *hStartTimeAT0 = new TH1F("hStartTimeA_T0","StartTime from T0A [ps]",1000,-1000,1000);
2100 fListQAt0->Add(hStartTimeAT0);
2101 TH1F *hStartTimeCT0 = new TH1F("hStartTimeC_T0","StartTime from T0C [ps]",1000,-1000,1000);
2102 fListQAt0->Add(hStartTimeCT0);
2103 TH1F *hStartTimeACT0 = new TH1F("hStartTimeAC_T0","StartTime from T0AC [ps]",1000,-1000,1000);;
2104 fListQAt0->Add(hStartTimeACT0);
2105
2106 TH1F *hnTracksAtT0 = new TH1F("hnTracksAt_T0","Tracks for events selected for T0",100,0,100);
2107 fListQAt0->Add(hnTracksAtT0);
2108 TH1F *hT0AEff = new TH1F("hT0AEff","Events with T0A vs nTracks",100,0,100);
2109 fListQAt0->Add(hT0AEff);
2110 TH1F *hT0CEff = new TH1F("hT0CEff","Events with T0C vs nTracks",100,0,100);
2111 fListQAt0->Add(hT0CEff);
2112 TH1F *hT0AndEff = new TH1F("hT0AndEff","Events with T0AC (AND) vs nTracks",100,0,100);
2113 fListQAt0->Add(hT0AndEff);
2114 TH1F *hT0OrEff = new TH1F("hT0OrEff","Events with T0AC (OR) vs nTracks",100,0,100);
2115 fListQAt0->Add(hT0OrEff);
2116
2117
2118}
2119
2120//______________________________________________________________________________
2121void AliAnalysisTaskPIDqa::SetupEMCALqa()
2122{
2123 //
2124 // Create the EMCAL qa objects
2125 //
2126
2127 TVectorD *vX=MakeLogBinning(200,.1,30);
2128
2129 TH2F *hNsigmaPt = new TH2F(Form("hNsigmaPt_EMCAL_%s",AliPID::ParticleName(0)),
2130 Form("EMCAL n#sigma %s vs. p_{T};p_{T} [GeV]; n#sigma",AliPID::ParticleName(0)),
2131 vX->GetNrows()-1,vX->GetMatrixArray(),
2132 200,-10,10);
2133 fListQAemcal->Add(hNsigmaPt);
2134
2135
2136 TH2F *hSigPtEle = new TH2F("hSigPt_EMCAL_Ele",
2137 "EMCAL signal (E/p) vs. p_{T} for electrons;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2138 vX->GetNrows()-1,vX->GetMatrixArray(),
2139 200,0,2);
2140 fListQAemcal->Add(hSigPtEle);
2141
2142 TH2F *hSigPtPions = new TH2F("hSigPt_EMCAL_Pions",
2143 "EMCAL signal (E/p) vs. p_{T} for pions;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2144 vX->GetNrows()-1,vX->GetMatrixArray(),
2145 200,0,2);
2146 fListQAemcal->Add(hSigPtPions);
2147
2148 TH2F *hSigPtProtons = new TH2F("hSigPt_EMCAL_Protons",
2149 "EMCAL signal (E/p) vs. p_{T} for protons;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2150 vX->GetNrows()-1,vX->GetMatrixArray(),
2151 200,0,2);
2152 fListQAemcal->Add(hSigPtProtons);
2153
2154 TH2F *hSigPtAntiProtons = new TH2F("hSigPt_EMCAL_Antiprotons",
2155 "EMCAL signal (E/p) vs. p_{T} for antiprotons;p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
2156 vX->GetNrows()-1,vX->GetMatrixArray(),
2157 200,0,2);
2158 fListQAemcal->Add(hSigPtAntiProtons);
2159
2160 delete vX;
2161}
2162
2163//______________________________________________________________________________
2164void AliAnalysisTaskPIDqa::SetupHMPIDqa()
2165{
2166 //
2167 // Create the HMPID qa objects
2168 //
2169
2170 TVectorD *vX=MakeLogBinning(200,.1,30);
2171
2172 // nSigmas
2173 Int_t nhists=0;
2174 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
2175 if (ispecie==AliPID::kElectron || ispecie==AliPID::kMuon) continue;
2176 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_HMPID_%s",AliPID::ParticleName(ispecie)),
2177 Form("HMPID n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2178 vX->GetNrows()-1,vX->GetMatrixArray(),
2179 200,-10,10);
2180 fListQAhmpid->AddAt(hNsigmaP, nhists);
2181 ++nhists;
2182 }
2183
2184 // cherenkov angle
2185 TH2F *hCkovAnglevsMom = new TH2F("hCkovAnglevsMom", "Cherenkov angle vs momentum",
2186 vX->GetNrows()-1,vX->GetMatrixArray(),
2187 500,0,1);
2188 fListQAhmpid->AddAt(hCkovAnglevsMom,nhists);
2189
2190 delete vX;
2191}
2192
2193//______________________________________________________________________________
2194void AliAnalysisTaskPIDqa::SetupTOFHMPIDqa()
2195{
2196 //
2197 // Create the HMPID qa objects
2198 //
2199
2200 TH2F *hCkovAnglevsMomPion = new TH2F("hCkovAnglevsMom_pion", "Cherenkov angle vs momentum for pions",500,0,5.,500,0,1);
2201 fListQAtofhmpid->Add(hCkovAnglevsMomPion);
2202
2203 TH2F *hCkovAnglevsMomKaon = new TH2F("hCkovAnglevsMom_kaon", "Cherenkov angle vs momentum for kaons",500,0,5.,500,0,1);
2204 fListQAtofhmpid->Add(hCkovAnglevsMomKaon);
2205
2206 TH2F *hCkovAnglevsMomProton = new TH2F("hCkovAnglevsMom_proton","Cherenkov angle vs momentum for protons",500,0,5.,500,0,1);
2207 fListQAtofhmpid->Add(hCkovAnglevsMomProton);
2208
2209
2210}
2211
2212//______________________________________________________________________________
2213void AliAnalysisTaskPIDqa::SetupTPCTOFqa()
2214{
2215 //
2216 // Create the qa objects for TPC + TOF combination
2217 //
2218
2219 TVectorD *vX=MakeLogBinning(200,.1,30);
2220
2221 //TPC signals after TOF cut
2222 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2223 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_TOF_%s",AliPID::ParticleName(ispecie)),
2224 Form("TPC n#sigma %s vs. p (after TOF 3#sigma cut);p_{TPC} [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2225 vX->GetNrows()-1,vX->GetMatrixArray(),
2226 200,-10,10);
2227 fListQAtpctof->Add(hNsigmaP);
2228 }
2229
2230 //TOF signals after TPC cut
2231 for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
2232 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_TPC_%s",AliPID::ParticleName(ispecie)),
2233 Form("TOF n#sigma %s vs. p (after TPC n#sigma cut);p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
2234 vX->GetNrows()-1,vX->GetMatrixArray(),
2235 200,-10,10);
2236 fListQAtpctof->Add(hNsigmaP);
2237 }
2238
2239 //EMCAL signal after TOF and TPC cut
2240 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
2241 TH2F *heopPt = new TH2F(Form("heopPt_TOF_TPC_%s",AliPID::ParticleName(ispecie)),
2242 Form("EMCAL signal (E/p) %s vs. p_{T};p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",AliPID::ParticleName(ispecie)),
2243 vX->GetNrows()-1,vX->GetMatrixArray(),
2244 200,0,2);
2245 fListQAtpctof->Add(heopPt);
2246 }
2247
2248 delete vX;
2249}
2250//______________________________________________________________________________
2251void AliAnalysisTaskPIDqa::SetupV0qa()
2252{
2253 //
2254 // Create the qa objects for V0 Kine cuts
2255 //
2256
2257 TH2F *hArmenteros = new TH2F("hArmenteros", "Armenteros plot",200,-1.,1.,200,0.,0.4);
2258 fListQAV0->Add(hArmenteros);
2259
2260}
2261
2262//_____________________________________________________________________________
2263void AliAnalysisTaskPIDqa::SetupQAinfo(){
2264 //
2265 // Setup the info of QA objects
2266 //
2267
2268 TObjArray *arr=new TObjArray;
2269 arr->SetName("TPC_info");
2270 fListQAinfo->Add(arr);
2271}
2272
2273//______________________________________________________________________________
2274TVectorD* AliAnalysisTaskPIDqa::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
2275{
2276 //
2277 // Make logarithmic binning
2278 // the user has to delete the array afterwards!!!
2279 //
2280
2281 //check limits
2282 if (xmin<1e-20 || xmax<1e-20){
2283 AliError("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
2284 return MakeLinBinning(nbinsX, xmin, xmax);
2285 }
2286 if (xmax<xmin){
2287 Double_t tmp=xmin;
2288 xmin=xmax;
2289 xmax=tmp;
2290 }
2291 TVectorD *binLim=new TVectorD(nbinsX+1);
2292 Double_t first=xmin;
2293 Double_t last=xmax;
2294 Double_t expMax=TMath::Log(last/first);
2295 for (Int_t i=0; i<nbinsX+1; ++i){
2296 (*binLim)[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
2297 }
2298 return binLim;
2299}
2300
2301//______________________________________________________________________________
2302TVectorD* AliAnalysisTaskPIDqa::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
2303{
2304 //
2305 // Make linear binning
2306 // the user has to delete the array afterwards!!!
2307 //
2308 if (xmax<xmin){
2309 Double_t tmp=xmin;
2310 xmin=xmax;
2311 xmax=tmp;
2312 }
2313 TVectorD *binLim=new TVectorD(nbinsX+1);
2314 Double_t first=xmin;
2315 Double_t last=xmax;
2316 Double_t binWidth=(last-first)/nbinsX;
2317 for (Int_t i=0; i<nbinsX+1; ++i){
2318 (*binLim)[i]=first+binWidth*(Double_t)i;
2319 }
2320 return binLim;
2321}
2322
2323//_____________________________________________________________________________
2324TVectorD* AliAnalysisTaskPIDqa::MakeArbitraryBinning(const char* bins)
2325{
2326 //
2327 // Make arbitrary binning, bins separated by a ','
2328 //
2329 TString limits(bins);
2330 if (limits.IsNull()){
2331 AliError("Bin Limit string is empty, cannot add the variable");
2332 return 0x0;
2333 }
2334
2335 TObjArray *arr=limits.Tokenize(",");
2336 Int_t nLimits=arr->GetEntries();
2337 if (nLimits<2){
2338 AliError("Need at leas 2 bin limits, cannot add the variable");
2339 delete arr;
2340 return 0x0;
2341 }
2342
2343 TVectorD *binLimits=new TVectorD(nLimits);
2344 for (Int_t iLim=0; iLim<nLimits; ++iLim){
2345 (*binLimits)[iLim]=(static_cast<TObjString*>(arr->At(iLim)))->GetString().Atof();
2346 }
2347
2348 delete arr;
2349 return binLimits;
2350}
2351