1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 /* $Id: AliAnalysisTaskPIDqa.cxx 43811 2010-09-23 14:13:31Z wiechula $ */
19 #include <TObjArray.h>
27 #include <AliAnalysisManager.h>
28 #include <AliInputEventHandler.h>
29 #include <AliVEventHandler.h>
30 #include <AliVEvent.h>
31 #include <AliVParticle.h>
32 #include <AliVTrack.h>
35 #include <AliPIDResponse.h>
36 #include <AliITSPIDResponse.h>
37 #include <AliTPCPIDResponse.h>
38 #include <AliTRDPIDResponse.h>
39 #include <AliTOFPIDResponse.h>
41 #include <AliESDEvent.h>
43 #include "AliAnalysisTaskPIDqa.h"
46 ClassImp(AliAnalysisTaskPIDqa)
48 //______________________________________________________________________________
49 AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa():
55 fListQAitsPureSA(0x0),
67 //______________________________________________________________________________
68 AliAnalysisTaskPIDqa::AliAnalysisTaskPIDqa(const char* name):
69 AliAnalysisTaskSE(name),
74 fListQAitsPureSA(0x0),
82 // Default constructor
84 DefineInput(0,TChain::Class());
85 DefineOutput(1,TList::Class());
88 //______________________________________________________________________________
89 AliAnalysisTaskPIDqa::~AliAnalysisTaskPIDqa()
97 //______________________________________________________________________________
98 void AliAnalysisTaskPIDqa::UserCreateOutputObjects()
101 // Create the output QA objects
104 AliLog::SetClassDebugLevel("AliAnalysisTaskPIDqa",10);
107 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
108 AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(man->GetInputEventHandler());
109 if (!inputHandler) AliFatal("Input handler needed");
111 //pid response object
112 fPIDResponse=inputHandler->GetPIDResponse();
113 if (!fPIDResponse) AliError("PIDResponse object was not created");
119 fListQAits=new TList;
120 fListQAits->SetOwner();
121 fListQAits->SetName("ITS");
123 fListQAitsSA=new TList;
124 fListQAitsSA->SetOwner();
125 fListQAitsSA->SetName("ITS_SA");
127 fListQAitsPureSA=new TList;
128 fListQAitsPureSA->SetOwner();
129 fListQAitsPureSA->SetName("ITS_PureSA");
131 fListQAtpc=new TList;
132 fListQAtpc->SetOwner();
133 fListQAtpc->SetName("TPC");
135 fListQAtrd=new TList;
136 fListQAtrd->SetOwner();
137 fListQAtrd->SetName("TRD");
139 fListQAtof=new TList;
140 fListQAtof->SetOwner();
141 fListQAtof->SetName("TOF");
143 fListQAemcal=new TList;
144 fListQAemcal->SetOwner();
145 fListQAemcal->SetName("EMCAL");
147 fListQAtpctof=new TList;
148 fListQAtpctof->SetOwner();
149 fListQAtpctof->SetName("TPC_TOF");
151 fListQA->Add(fListQAits);
152 fListQA->Add(fListQAitsSA);
153 fListQA->Add(fListQAitsPureSA);
154 fListQA->Add(fListQAtpc);
155 fListQA->Add(fListQAtrd);
156 fListQA->Add(fListQAtof);
157 fListQA->Add(fListQAemcal);
158 fListQA->Add(fListQAtpctof);
171 //______________________________________________________________________________
172 void AliAnalysisTaskPIDqa::UserExec(Option_t */*option*/)
175 // Setup the PID response functions and fill the QA histograms
178 AliVEvent *event=InputEvent();
179 if (!event||!fPIDResponse) return;
192 //______________________________________________________________________________
193 void AliAnalysisTaskPIDqa::FillITSqa()
196 // Fill PID qa histograms for the ITS
199 AliVEvent *event=InputEvent();
201 Int_t ntracks=event->GetNumberOfTracks();
202 for(Int_t itrack = 0; itrack < ntracks; itrack++){
203 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
204 ULong_t status=track->GetStatus();
205 // not that nice. status bits not in virtual interface
206 // ITS refit + ITS pid selection
207 if (!( ( (status & AliVTrack::kITSrefit)==AliVTrack::kITSrefit ) ||
208 ! ( (status & AliVTrack::kITSpid )==AliVTrack::kITSpid ) )) continue;
209 Double_t mom=track->P();
211 TList *theList = 0x0;
212 if(( (status & AliVTrack::kTPCin)==AliVTrack::kTPCin )){
216 if(!( (status & AliVTrack::kITSpureSA)==AliVTrack::kITSpureSA )){
217 //ITS Standalone tracks
218 theList=fListQAitsSA;
220 //ITS Pure Standalone tracks
221 theList=fListQAitsPureSA;
226 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
227 TH2 *h=(TH2*)theList->At(ispecie);
229 Double_t nSigma=fPIDResponse->NumberOfSigmasITS(track, (AliPID::EParticleType)ispecie);
232 TH2 *h=(TH2*)theList->At(AliPID::kSPECIES);
234 Double_t sig=track->GetITSsignal();
240 //______________________________________________________________________________
241 void AliAnalysisTaskPIDqa::FillTPCqa()
244 // Fill PID qa histograms for the TPC
247 AliVEvent *event=InputEvent();
249 Int_t ntracks=event->GetNumberOfTracks();
250 for(Int_t itrack = 0; itrack < ntracks; itrack++){
251 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
256 ULong_t status=track->GetStatus();
257 // not that nice. status bits not in virtual interface
258 // TPC refit + ITS refit + TPC pid
259 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
260 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
261 !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) ) continue;
263 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
264 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
265 if (track->GetTPCNclsF()>0) {
266 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
269 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
271 Double_t mom=track->GetTPCmomentum();
273 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
274 TH2 *h=(TH2*)fListQAtpc->At(ispecie);
276 Double_t nSigma=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
280 TH2 *h=(TH2*)fListQAtpc->At(AliPID::kSPECIES);
282 Double_t sig=track->GetTPCsignal();
288 //______________________________________________________________________________
289 void AliAnalysisTaskPIDqa::FillTRDqa()
292 // Fill PID qa histograms for the TRD
294 AliVEvent *event=InputEvent();
295 Int_t ntracks = event->GetNumberOfTracks();
296 for(Int_t itrack = 0; itrack < ntracks; itrack++){
297 AliVTrack *track = (AliVTrack *)event->GetTrack(itrack);
302 ULong_t status=track->GetStatus();
303 // not that nice. status bits not in virtual interface
304 // TPC refit + ITS refit + TPC pid + TRD out
305 if (!( (status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
306 !( (status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
307 !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) ||
308 !( (status & AliVTrack::kTRDout ) == AliVTrack::kTRDout )) continue;
310 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
311 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
312 if (track->GetTPCNclsF()>0) {
313 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
316 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
318 Double_t likelihoods[AliPID::kSPECIES];
319 if(fPIDResponse->ComputeTRDProbability(track, AliPID::kSPECIES, likelihoods) != AliPIDResponse::kDetPidOk) continue;
320 Int_t ntracklets = 0;
321 Double_t momentum = -1.;
322 for(Int_t itl = 0; itl < 6; itl++)
323 if(track->GetTRDmomentum(itl) > 0.){
325 if(momentum < 0) momentum = track->GetTRDmomentum(itl);
327 for(Int_t ispecie = 0; ispecie < AliPID::kSPECIES; ispecie++){
328 TH2F *hLike = (TH2F *)fListQAtrd->At(ntracklets*AliPID::kSPECIES+ispecie);
329 if (hLike) hLike->Fill(momentum,likelihoods[ispecie]);
334 //______________________________________________________________________________
335 void AliAnalysisTaskPIDqa::FillTOFqa()
337 AliVEvent *event=InputEvent();
341 Int_t ntracks=event->GetNumberOfTracks();
342 Int_t tracksAtTof = 0;
343 for(Int_t itrack = 0; itrack < ntracks; itrack++){
344 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
349 ULong_t status=track->GetStatus();
350 // TPC refit + ITS refit +
351 // TOF out + TOFpid +
353 // (we don't use kTOFmismatch because it depends on TPC....)
354 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
355 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
356 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
357 !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
358 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
360 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
361 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
362 if (track->GetTPCNclsF()>0) {
363 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
366 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
370 Double_t mom=track->P();
372 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
373 TH2 *h=(TH2*)fListQAtof->At(ispecie);
375 Double_t nSigma=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
379 TH2 *h=(TH2*)fListQAtof->FindObject("hSigP_TOF");
381 Double_t sig=track->GetTOFsignal();
385 Double_t mask = (Double_t)fPIDResponse->GetTOFResponse().GetStartTimeMask(mom) + 0.5;
386 ((TH1F*)fListQAtof->FindObject("hStartTimeMask_TOF"))->Fill(mask);
388 if (mom >= 1.0 && mom <= 2.0 ) {
389 Double_t nsigma= fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)AliPID::kKaon);
391 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Kaon_T0-Fill"))->Fill(nsigma);
392 } else if (mask == 1) {
393 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Kaon_T0-TOF"))->Fill(nsigma);
394 } else if ( (mask == 2) || (mask == 4) || (mask == 6) ) {
395 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Kaon_T0-T0"))->Fill(nsigma);
397 ((TH1F*)fListQAtof->FindObject("hNsigma_TOF_Kaon_T0-Best"))->Fill(nsigma);
401 Double_t res = (Double_t)fPIDResponse->GetTOFResponse().GetStartTimeRes(mom);
402 ((TH1F*)fListQAtof->FindObject("hStartTimeRes_TOF"))->Fill(res);
404 AliESDEvent *esd = dynamic_cast<AliESDEvent *>(event);
406 Double_t startTime = esd->GetT0TOF(0);
407 if (startTime < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeAC_T0"))->Fill(startTime);
409 startTime = esd->GetT0TOF(1);
410 if (startTime < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeA_T0"))->Fill(startTime);
411 startTime = esd->GetT0TOF(2);
412 if (startTime < 90000) ((TH1F*)fListQAtof->FindObject("hStartTimeC_T0"))->Fill(startTime);
416 if (tracksAtTof > 0) {
417 ((TH1F* )fListQAtof->FindObject("hnTracksAt_TOF"))->Fill(tracksAtTof);
418 Int_t mask = fPIDResponse->GetTOFResponse().GetStartTimeMask(5.);
419 if (mask & 0x1) ((TH1F*)fListQAtof->FindObject("hT0MakerEff"))->Fill(tracksAtTof);
425 //______________________________________________________________________________
426 void AliAnalysisTaskPIDqa::FillEMCALqa()
429 // Fill PID qa histograms for the EMCAL
432 AliVEvent *event=InputEvent();
434 Int_t ntracks=event->GetNumberOfTracks();
435 for(Int_t itrack = 0; itrack < ntracks; itrack++){
436 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
441 ULong_t status=track->GetStatus();
442 // not that nice. status bits not in virtual interface
443 // TPC refit + ITS refit +
444 // TOF out + TOFpid +
446 if (!( (status & AliVTrack::kEMCALmatch) == AliVTrack::kEMCALmatch) ) continue;
448 Double_t pt=track->Pt();
450 //EMCAL nSigma (only for electrons at the moment)
451 TH2 *h=(TH2*)fListQAemcal->At(0);
453 Double_t nSigma=fPIDResponse->NumberOfSigmasEMCAL(track, (AliPID::EParticleType)0);
456 //EMCAL signal (E/p vs. pT)
457 h=(TH2*)fListQAemcal->At(1);
460 Int_t nMatchClus = track->GetEMCALcluster();
461 Double_t mom = track->P();
466 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
470 // matched cluster is EMCAL
471 if(matchedClus->IsEMCAL()){
473 Double_t fClsE = matchedClus->E();
482 Printf("status status = AliVTrack::kEMCALmatch, BUT no matched cluster!");
488 //______________________________________________________________________________
489 void AliAnalysisTaskPIDqa::FillTPCTOFqa()
492 // Fill PID qa histograms for the TOF
493 // Here also the TPC histograms after TOF selection are filled
496 AliVEvent *event=InputEvent();
498 Int_t ntracks=event->GetNumberOfTracks();
499 for(Int_t itrack = 0; itrack < ntracks; itrack++){
500 AliVTrack *track=(AliVTrack*)event->GetTrack(itrack);
505 ULong_t status=track->GetStatus();
506 // not that nice. status bits not in virtual interface
507 // TPC refit + ITS refit +
508 // TOF out + TOFpid +
510 if (!((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) ||
511 !((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) ||
512 !( (status & AliVTrack::kTPCpid ) == AliVTrack::kTPCpid ) ||
513 !((status & AliVTrack::kTOFout ) == AliVTrack::kTOFout ) ||
514 !((status & AliVTrack::kTOFpid ) == AliVTrack::kTOFpid ) ||
515 !((status & AliVTrack::kTIME ) == AliVTrack::kTIME ) ) continue;
517 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
518 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
519 if (track->GetTPCNclsF()>0) {
520 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC/track->GetTPCNclsF();
523 if ( nCrossedRowsTPC<70 || ratioCrossedRowsOverFindableClustersTPC<.8 ) continue;
526 Double_t mom=track->P();
527 Double_t momTPC=track->GetTPCmomentum();
529 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
531 Double_t nSigmaTOF=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
532 Double_t nSigmaTPC=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
535 TH2 *h=(TH2*)fListQAtpctof->At(ispecie);
536 if (h && TMath::Abs(nSigmaTOF)<3.) h->Fill(momTPC,nSigmaTPC);
539 h=(TH2*)fListQAtpctof->At(ispecie+AliPID::kSPECIES);
540 if (h && TMath::Abs(nSigmaTPC)<3.) h->Fill(mom,nSigmaTOF);
542 //EMCAL after TOF and TPC cut
543 h=(TH2*)fListQAtpctof->At(ispecie+2*AliPID::kSPECIES);
544 if (h && TMath::Abs(nSigmaTOF)<3. && TMath::Abs(nSigmaTPC)<3. ){
546 Int_t nMatchClus = track->GetEMCALcluster();
547 Double_t pt = track->Pt();
552 AliVCluster *matchedClus = (AliVCluster*)event->GetCaloCluster(nMatchClus);
556 // matched cluster is EMCAL
557 if(matchedClus->IsEMCAL()){
559 Double_t fClsE = matchedClus->E();
573 //______________________________________________________________________________
574 void AliAnalysisTaskPIDqa::SetupITSqa()
577 // Create the ITS qa objects
580 TVectorD *vX=MakeLogBinning(200,.1,30);
583 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
584 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_ITS_%s",AliPID::ParticleName(ispecie)),
585 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
586 vX->GetNrows()-1,vX->GetMatrixArray(),
588 fListQAits->Add(hNsigmaP);
590 TH2F *hSig = new TH2F("hSigP_ITS",
591 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
592 vX->GetNrows()-1,vX->GetMatrixArray(),
594 fListQAits->Add(hSig);
596 //ITS Standalone tracks
597 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
598 TH2F *hNsigmaPSA = new TH2F(Form("hNsigmaP_ITSSA_%s",AliPID::ParticleName(ispecie)),
599 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
600 vX->GetNrows()-1,vX->GetMatrixArray(),
602 fListQAitsSA->Add(hNsigmaPSA);
604 TH2F *hSigSA = new TH2F("hSigP_ITSSA",
605 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
606 vX->GetNrows()-1,vX->GetMatrixArray(),
608 fListQAitsSA->Add(hSigSA);
610 //ITS Pure Standalone tracks
611 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
612 TH2F *hNsigmaPPureSA = new TH2F(Form("hNsigmaP_ITSPureSA_%s",AliPID::ParticleName(ispecie)),
613 Form("ITS n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
614 vX->GetNrows()-1,vX->GetMatrixArray(),
616 fListQAitsPureSA->Add(hNsigmaPPureSA);
618 TH2F *hSigPureSA = new TH2F("hSigP_ITSPureSA",
619 "ITS signal vs. p;p [GeV]; ITS signal [arb. units]",
620 vX->GetNrows()-1,vX->GetMatrixArray(),
622 fListQAitsPureSA->Add(hSigPureSA);
627 //______________________________________________________________________________
628 void AliAnalysisTaskPIDqa::SetupTPCqa()
631 // Create the TPC qa objects
634 TVectorD *vX=MakeLogBinning(200,.1,30);
636 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
637 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_%s",AliPID::ParticleName(ispecie)),
638 Form("TPC n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
639 vX->GetNrows()-1,vX->GetMatrixArray(),
641 fListQAtpc->Add(hNsigmaP);
645 TH2F *hSig = new TH2F("hSigP_TPC",
646 "TPC signal vs. p;p [GeV]; TPC signal [arb. units]",
647 vX->GetNrows()-1,vX->GetMatrixArray(),
649 fListQAtpc->Add(hSig);
654 //______________________________________________________________________________
655 void AliAnalysisTaskPIDqa::SetupTRDqa()
658 // Create the TRD qa objects
660 TVectorD *vX=MakeLogBinning(200,.1,30);
661 for(Int_t itl = 0; itl < 6; ++itl){
662 for(Int_t ispecie = 0; ispecie < AliPID::kSPECIES; ispecie++){
663 TH2F *hLikeP = new TH2F(Form("hLikeP_TRD_%dtls_%s", itl, AliPID::ParticleName(ispecie)),
664 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)),
665 vX->GetNrows()-1, vX->GetMatrixArray(),
667 fListQAtrd->Add(hLikeP);
673 //______________________________________________________________________________
674 void AliAnalysisTaskPIDqa::SetupTOFqa()
677 // Create the TOF qa objects
680 TVectorD *vX=MakeLogBinning(200,.1,30);
682 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
683 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_%s",AliPID::ParticleName(ispecie)),
684 Form("TOF n#sigma %s vs. p;p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
685 vX->GetNrows()-1,vX->GetMatrixArray(),
687 fListQAtof->Add(hNsigmaP);
690 // for Kaons PID we differentiate on Time Zero
691 TH1F *hnSigT0Fill = new TH1F("hNsigma_TOF_Kaon_T0-Fill","TOF n#sigma (Kaon) T0-FILL [1-2. GeV/c]",200,-10,10);
692 fListQAtof->Add(hnSigT0Fill);
693 TH1F *hnSigT0T0 = new TH1F("hNsigma_TOF_Kaon_T0-T0","TOF n#sigma (Kaon) T0-T0 [1-2. GeV/c]",200,-10,10);
694 fListQAtof->Add(hnSigT0T0);
695 TH1F *hnSigT0TOF = new TH1F("hNsigma_TOF_Kaon_T0-TOF","TOF n#sigma (Kaon) T0-TOF [1.-2. GeV/c]",200,-10,10);
696 fListQAtof->Add(hnSigT0TOF);
697 TH1F *hnSigT0Best = new TH1F("hNsigma_TOF_Kaon_T0-Best","TOF n#sigma (Kaon) T0-Best [1-2. GeV/c]",200,-10,10);
698 fListQAtof->Add(hnSigT0Best);
700 TH2F *hSig = new TH2F("hSigP_TOF",
701 "TOF signal vs. p;p [GeV]; TOF signal [arb. units]",
702 vX->GetNrows()-1,vX->GetMatrixArray(),
707 fListQAtof->Add(hSig);
709 TH1F *hStartTimeMaskTOF = new TH1F("hStartTimeMask_TOF","StartTime mask",8,0,8);
710 fListQAtof->Add(hStartTimeMaskTOF);
711 TH1F *hStartTimeResTOF = new TH1F("hStartTimeRes_TOF","StartTime resolution [ps]",100,0,500);
712 fListQAtof->Add(hStartTimeResTOF);
714 TH1F *hnTracksAtTOF = new TH1F("hnTracksAt_TOF","Matched tracks at TOF",20,0,20);
715 fListQAtof->Add(hnTracksAtTOF);
716 TH1F *hT0MakerEff = new TH1F("hT0MakerEff","Events with T0-TOF vs nTracks",20,0,20);
717 fListQAtof->Add(hT0MakerEff);
719 // this in principle should stay on a T0 PID QA, but are just the data prepared for TOF use
720 TH1F *hStartTimeAT0 = new TH1F("hStartTimeA_T0","StartTime from T0A [ps]",1000,-1000,1000);
721 fListQAtof->Add(hStartTimeAT0);
722 TH1F *hStartTimeCT0 = new TH1F("hStartTimeC_T0","StartTime from T0C [ps]",1000,-1000,1000);
723 fListQAtof->Add(hStartTimeCT0);
724 TH1F *hStartTimeACT0 = new TH1F("hStartTimeAC_T0","StartTime from T0AC [ps]",1000,-1000,1000);;
725 fListQAtof->Add(hStartTimeACT0);
728 //______________________________________________________________________________
729 void AliAnalysisTaskPIDqa::SetupEMCALqa()
732 // Create the EMCAL qa objects
735 TVectorD *vX=MakeLogBinning(200,.1,30);
737 TH2F *hNsigmaPt = new TH2F(Form("hNsigmaPt_EMCAL_%s",AliPID::ParticleName(0)),
738 Form("EMCAL n#sigma %s vs. p_{T};p_{T} [GeV]; n#sigma",AliPID::ParticleName(0)),
739 vX->GetNrows()-1,vX->GetMatrixArray(),
741 fListQAemcal->Add(hNsigmaPt);
743 TH2F *hSigPt = new TH2F("hSigPt_EMCAL",
744 "EMCAL signal (E/p) vs. p_{T};p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",
745 vX->GetNrows()-1,vX->GetMatrixArray(),
747 fListQAemcal->Add(hSigPt);
752 //______________________________________________________________________________
753 void AliAnalysisTaskPIDqa::SetupTPCTOFqa()
756 // Create the qa objects for TPC + TOF combination
759 TVectorD *vX=MakeLogBinning(200,.1,30);
761 //TPC signals after TOF cut
762 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
763 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TPC_TOF_%s",AliPID::ParticleName(ispecie)),
764 Form("TPC n#sigma %s vs. p (after TOF 3#sigma cut);p_{TPC} [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
765 vX->GetNrows()-1,vX->GetMatrixArray(),
767 fListQAtpctof->Add(hNsigmaP);
770 //TOF signals after TPC cut
771 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
772 TH2F *hNsigmaP = new TH2F(Form("hNsigmaP_TOF_TPC_%s",AliPID::ParticleName(ispecie)),
773 Form("TOF n#sigma %s vs. p (after TPC n#sigma cut);p [GeV]; n#sigma",AliPID::ParticleName(ispecie)),
774 vX->GetNrows()-1,vX->GetMatrixArray(),
776 fListQAtpctof->Add(hNsigmaP);
779 //EMCAL signal after TOF and TPC cut
780 for (Int_t ispecie=0; ispecie<AliPID::kSPECIES; ++ispecie){
781 TH2F *heopPt = new TH2F(Form("heopPt_TOF_TPC_%s",AliPID::ParticleName(ispecie)),
782 Form("EMCAL signal (E/p) %s vs. p_{T};p_{T} [GeV]; EMCAL signal (E/p) [arb. units]",AliPID::ParticleName(ispecie)),
783 vX->GetNrows()-1,vX->GetMatrixArray(),
785 fListQAtpctof->Add(heopPt);
791 //______________________________________________________________________________
792 TVectorD* AliAnalysisTaskPIDqa::MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
795 // Make logarithmic binning
796 // the user has to delete the array afterwards!!!
800 if (xmin<1e-20 || xmax<1e-20){
801 AliError("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
802 return MakeLinBinning(nbinsX, xmin, xmax);
809 TVectorD *binLim=new TVectorD(nbinsX+1);
812 Double_t expMax=TMath::Log(last/first);
813 for (Int_t i=0; i<nbinsX+1; ++i){
814 (*binLim)[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
819 //______________________________________________________________________________
820 TVectorD* AliAnalysisTaskPIDqa::MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
823 // Make linear binning
824 // the user has to delete the array afterwards!!!
831 TVectorD *binLim=new TVectorD(nbinsX+1);
834 Double_t binWidth=(last-first)/nbinsX;
835 for (Int_t i=0; i<nbinsX+1; ++i){
836 (*binLim)[i]=first+binWidth*(Double_t)i;
841 //_____________________________________________________________________________
842 TVectorD* AliAnalysisTaskPIDqa::MakeArbitraryBinning(const char* bins)
845 // Make arbitrary binning, bins separated by a ','
847 TString limits(bins);
848 if (limits.IsNull()){
849 AliError("Bin Limit string is empty, cannot add the variable");
853 TObjArray *arr=limits.Tokenize(",");
854 Int_t nLimits=arr->GetEntries();
856 AliError("Need at leas 2 bin limits, cannot add the variable");
861 TVectorD *binLimits=new TVectorD(nLimits);
862 for (Int_t iLim=0; iLim<nLimits; ++iLim){
863 (*binLimits)[iLim]=(static_cast<TObjString*>(arr->At(iLim)))->GetString().Atof();