6e03e3711893d70e9760dae5ba35f6eb0acffe26
[u/mrichter/AliRoot.git] / PWG / FLOW / Tasks / AliAnalysisTaskPIDconfig.cxx
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: AliAnalysisTaskPIDconfig.cxx 43811 2014-10-11 Naghmeh Mohammadi $ */
17
18 #include "TChain.h"
19 #include "TTree.h"
20 #include "TList.h"
21 #include "TMath.h"
22 #include "TObjArray.h"
23 #include "TCanvas.h"
24 #include "TGraphErrors.h"
25 #include "TString.h"
26 #include "TFile.h"
27 #include "TH1F.h"
28 #include "TH2F.h"
29 #include "TH3F.h" 
30 #include "TH2D.h"                  
31 #include "TH3D.h"
32 #include "TArrayF.h"
33 #include "TF1.h"
34 #include "TROOT.h"
35 #include "stdio.h"
36 #include "TCutG.h"
37
38
39 #include "AliTHn.h"
40 #include "AliLog.h"
41 #include "AliAnalysisManager.h"
42 #include "AliESDEvent.h"
43 #include "AliAODInputHandler.h"
44 #include "AliAODEvent.h"
45 #include "AliAODTrack.h"
46 #include "AliAODInputHandler.h"
47 #include "AliCollisionGeometry.h"
48 #include "AliGenEventHeader.h"
49 #include "AliAnalysisUtils.h"
50 #include "AliPIDCombined.h"        
51 #include "AliAnalysisTask.h"
52 #include "AliAODHandler.h"
53 #include <AliInputEventHandler.h>
54 #include <AliVEventHandler.h>
55 #include <AliVParticle.h>
56 #include <AliVTrack.h>
57 #include <AliTPCPIDResponse.h>
58 #include <AliTOFPIDResponse.h>
59 #include "AliAnalysisTaskPIDconfig.h"
60 #include "AliAnalysisTaskSE.h"
61 #include "AliAODPid.h"
62 #include "AliPhysicsSelection.h"
63 #include "AliCentralitySelectionTask.h"
64 #include "AliCentrality.h"
65 #include "AliKFParticle.h"
66 #include "AliKFVertex.h"
67 #include "AliPID.h"
68 #include "AliPIDResponse.h"
69 #include "AliCFContainer.h"
70 #include "AliCFManager.h"
71 #include "AliVEvent.h"
72 #include "AliAODVZERO.h"
73
74
75 ClassImp(AliAnalysisTaskPIDconfig)
76 //ClassImp()
77 //___________________________________________________________________
78 AliAnalysisTaskPIDconfig::AliAnalysisTaskPIDconfig():
79 AliAnalysisTaskSE(),
80 fVevent(0),
81 fESD(0),
82 fAOD(0),
83 fPIDResponse(0),
84 fTriggerSelection(0),
85 fCentralityPercentileMin(0.),
86 fCentralityPercentileMax(5.),
87 fFilterBit(128),
88 fDCAxyCut(-1),
89 fDCAzCut(-1),
90 fData2011(kFALSE),
91 fTriggerMB(kTRUE),
92 fTriggerCentral(kFALSE),
93 fUseCentrality(kTRUE),
94 fCutTPCmultiplicityOutliersAOD(kFALSE),
95 fPIDcuts(kFALSE),
96 fCentralityEstimator("V0M"),
97 fContourCutList(0),
98 fListQA(0x0),
99 fListQAtpctof(0x0),
100 fListQAInfo(0x0),
101 fhistCentralityPass(0),
102 fNoEvents(0),
103 fpVtxZ(0),
104 fhistDCABefore(0),
105 fhistDCAAfter(0),
106 fhistPhiDistBefore(0),
107 fhistPhiDistAfter(0),
108 fhistEtaDistBefore(0),
109 fhistEtaDistAfter(0),
110 fTPCvsGlobalMultBeforeOutliers(0),
111 fTPCvsGlobalMultAfterOutliers(0),
112 fTPCvsGlobalMultAfter(0),
113 fHistBetavsPTOFbeforePID(0),
114 fHistdEdxVsPTPCbeforePID(0),
115 fHistBetavsPTOFafterPID(0),
116 fHistdEdxVsPTPCafterPID(0),
117 fhistNsigmaP(0),
118 fhistNsigmaPt(0)
119 //fSparseSpecies(0),
120 //fvalueSpecies(0)
121 {
122     for(int i=0;i<3;i++){for(int j=0;j<10;j++){fContourCut[i][j]= NULL;}}
123 }
124
125
126 //___________________________________________________________________
127
128 AliAnalysisTaskPIDconfig::AliAnalysisTaskPIDconfig(const char *name):
129 AliAnalysisTaskSE(name),
130 fVevent(0),
131 fESD(0),
132 fAOD(0),
133 fPIDResponse(0),
134 fTriggerSelection(0),
135 fCentralityPercentileMin(0.),
136 fCentralityPercentileMax(5.),
137 fFilterBit(1),
138 fDCAxyCut(-1),
139 fDCAzCut(-1),
140 fData2011(kFALSE),
141 fTriggerMB(kTRUE),
142 fTriggerCentral(kFALSE),
143 fUseCentrality(kTRUE),
144 fCutTPCmultiplicityOutliersAOD(kFALSE),
145 fPIDcuts(kFALSE),
146 fCentralityEstimator("V0M"),
147 fContourCutList(0),
148 fListQA(0x0),
149 fListQAtpctof(0x0),
150 fListQAInfo(0x0),
151 fhistCentralityPass(0),
152 fNoEvents(0),
153 fpVtxZ(0),
154 fhistDCABefore(0),
155 fhistDCAAfter(0),
156 fhistPhiDistBefore(0),
157 fhistPhiDistAfter(0),
158 fhistEtaDistBefore(0),
159 fhistEtaDistAfter(0),
160 fTPCvsGlobalMultBeforeOutliers(0),
161 fTPCvsGlobalMultAfterOutliers(0),
162 fTPCvsGlobalMultAfter(0),
163 fHistBetavsPTOFbeforePID(0),
164 fHistdEdxVsPTPCbeforePID(0),
165 fHistBetavsPTOFafterPID(0),
166 fHistdEdxVsPTPCafterPID(0),
167 fhistNsigmaP(0),
168 fhistNsigmaPt(0)
169 //fSparseSpecies(0),
170 //fvalueSpecies(0)
171 {
172     //fvalueSpecies = new Double_t[9];
173     //Default Constructor
174     //fContourCut[3][10]=NULL;
175     DefineInput(0,TChain::Class());
176     DefineOutput(1,TList::Class());
177 }
178
179 //_____________________________________________________________________
180 AliAnalysisTaskPIDconfig::~AliAnalysisTaskPIDconfig()
181 {
182     //Destructor
183     
184  //   delete fPID;
185   //  delete fPIDqa;
186  //   delete fTrackCuts;
187    // delete fSparseSpecies;
188     //delete []fvalueSpecies;
189
190
191 }
192 //______________________________________________________________________
193 void AliAnalysisTaskPIDconfig::UserCreateOutputObjects()
194 {
195     //
196     // Create the output QA objects
197     //
198     
199     AliLog::SetClassDebugLevel("AliAnalysisTaskPIDconfig",10);
200     
201     //input hander
202     AliAnalysisManager *mgr=AliAnalysisManager::GetAnalysisManager();
203     AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(mgr->GetInputEventHandler());
204     if (!inputHandler) AliFatal("Input handler needed");
205     
206     //pid response object
207     fPIDResponse=inputHandler->GetPIDResponse();
208     if (!fPIDResponse) AliError("PIDResponse object was not created");
209     
210     if(fPIDcuts){ GetPIDContours();}
211     //
212     fListQA=new TList;
213     fListQA->SetOwner();
214     
215     fListQAtpctof=new TList;
216     fListQAtpctof->SetOwner();
217     fListQAtpctof->SetName("PID_TPC_TOF");
218     
219     fListQAInfo=new TList;
220     fListQAInfo->SetOwner();
221     fListQAInfo->SetName("Event_Track_Info");
222     
223     fListQA->Add(fListQAtpctof);
224     fListQA->Add(fListQAInfo);
225     
226     SetupTPCTOFqa();
227     SetupEventInfo();
228     
229     PostData(1,fListQA);
230 }
231 //______________________________________________________________________
232 void AliAnalysisTaskPIDconfig::UserExec(Option_t*){
233     //Main loop
234     //Called for each event
235     
236     // create pointer to event
237     fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
238     fESD = dynamic_cast<AliESDEvent*>(InputEvent());
239     
240     Int_t ntracks=fAOD->GetNumberOfTracks();
241
242
243     if(!(fESD || fAOD)){
244         printf("ERROR: fESD & fAOD not available\n");
245         return;
246     }
247     fVevent = dynamic_cast<AliVEvent*>(InputEvent());
248     if (!fVevent) {
249         printf("ERROR: fVevent not available\n");
250         return;
251     }
252     
253     Bool_t pass = kFALSE;
254     CheckCentrality(fVevent,pass);
255     
256     if(!pass){ return;}
257     
258     const AliVVertex *pVtx = fVevent->GetPrimaryVertex();
259     
260     Double_t pVtxZ = -999;
261     pVtxZ = pVtx->GetZ();
262     
263     if(TMath::Abs(pVtxZ)>10) return;
264
265     TH1F *hNoEvents = (TH1F*)fListQAInfo->At(1);
266     TH1F *histpVtxZ = (TH1F*)fListQAInfo->At(2);
267
268     if(hNoEvents) hNoEvents->Fill(0);
269     if(histpVtxZ) histpVtxZ->Fill(pVtxZ);
270
271     if(ntracks<2) return;
272     
273    // if(!pass) return;
274     
275     TH2F *HistTPCvsGlobalMultBeforeOutliers = (TH2F*)fListQAInfo->At(3);
276
277     TH2F *HistTPCvsGlobalMultAfterOutliers = (TH2F*)fListQAInfo->At(4);
278
279     
280     Float_t multTPC(0.); // tpc mult estimate
281     Float_t multGlobal(0.); // global multiplicity
282     
283     const Int_t nGoodTracks = fVevent->GetNumberOfTracks();
284     if(!fData2011) { // cut on outliers
285       
286         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill tpc mult
287             AliAODTrack* AODtrack =dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
288             if (!AODtrack) continue;
289             if (!(AODtrack->TestFilterBit(1))) continue;
290             if ((AODtrack->Pt() < .2) || (AODtrack->Pt() > 5.0) || (TMath::Abs(AODtrack->Eta()) > .8) || (AODtrack->GetTPCNcls() < 70)  || (AODtrack->GetDetPid()->GetTPCsignal() < 10.0) || (AODtrack->Chi2perNDF() < 0.2)) continue;
291             multTPC++;
292         }//track loop
293
294         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill global mult
295             AliAODTrack *AODtrack=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
296             if (!AODtrack) continue;
297             if (!(AODtrack->TestFilterBit(16))) continue;
298             if ((AODtrack->Pt() < .2) || (AODtrack->Pt() > 5.0) || (TMath::Abs(AODtrack->Eta()) > .8) || (AODtrack->GetTPCNcls() < 70) || (AODtrack->GetDetPid()->GetTPCsignal() < 10.0) || (AODtrack->Chi2perNDF() < 0.1)) continue;
299             Double_t b[2] = {-99., -99.};
300             Double_t bCov[3] = {-99., -99., -99.};
301             AliAODTrack copy(*AODtrack);
302             if (!(copy.PropagateToDCA(fVevent->GetPrimaryVertex(), fVevent->GetMagneticField(), 100., b, bCov))) continue;
303             if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
304             multGlobal++;
305         } //track loop
306         
307         HistTPCvsGlobalMultBeforeOutliers->Fill(multGlobal,multTPC);
308
309         if(multTPC < (-40.3+1.22*multGlobal) || multTPC > (32.1+1.59*multGlobal)){ pass = kFALSE;}
310         
311         if(!pass) return;
312         HistTPCvsGlobalMultAfterOutliers->Fill(multGlobal,multTPC);
313
314     }
315     
316     
317     if(fData2011) { // cut on outliers
318         //Float_t multTPC(0.); // tpc mult estimate
319         //Float_t multGlob(0.); // global multiplicity
320         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill tpc mult
321             AliAODTrack *AODtrack=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
322             if (!AODtrack) continue;
323             if (!(AODtrack->TestFilterBit(1))) continue;
324             if ((AODtrack->Pt() < .2) || (AODtrack->Pt() > 5.0) || (TMath::Abs(AODtrack->Eta()) > .8) || (AODtrack->GetTPCNcls() < 70)  || (AODtrack->GetDetPid()->GetTPCsignal() < 10.0) || (AODtrack->Chi2perNDF() < 0.2)) continue;
325             multTPC++;
326         }
327         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill global mult
328             AliAODTrack *AODtrack=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
329             if (!AODtrack) continue;
330             if (!(AODtrack->TestFilterBit(16))) continue;
331             if ((AODtrack->Pt() < .2) || (AODtrack->Pt() > 5.0) || (TMath::Abs(AODtrack->Eta()) > .8) || (AODtrack->GetTPCNcls() < 70) || (AODtrack->GetDetPid()->GetTPCsignal() < 10.0) || (AODtrack->Chi2perNDF() < 0.1)) continue;
332             Double_t b[2] = {-99., -99.};
333             Double_t bCov[3] = {-99., -99., -99.};
334             AliAODTrack copy(*AODtrack);
335             if (!(copy.PropagateToDCA(fVevent->GetPrimaryVertex(), fVevent->GetMagneticField(), 100., b, bCov))) continue;
336             if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
337             multGlobal++;
338             
339         } //track loop
340         
341         HistTPCvsGlobalMultBeforeOutliers->Fill(multGlobal,multTPC);
342
343         if(multTPC < (-36.73 + 1.48*multGlobal) || multTPC > (62.87 + 1.78*multGlobal)){pass = kFALSE;}
344         
345         if(!pass) return;
346         HistTPCvsGlobalMultAfterOutliers->Fill(multGlobal,multTPC);
347
348     }
349     
350      for(Int_t itrack = 0; itrack < ntracks; itrack++){
351
352         AliAODTrack *track=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(itrack));
353         if(!track) continue;
354         
355         Float_t dcaXY = track->DCA();
356         Float_t dcaZ  = track->ZAtDCA();
357          
358         TH2F* HistDCAbefore =(TH2F*)fListQAInfo->At(5);
359         HistDCAbefore->Fill(dcaZ,dcaXY);
360         
361         Double_t p = -999, /*pTPC = -999,*/ pT = -999, phi = -999, eta = -999, dEdx =-999;
362         Double_t length = -999., beta =-999, tofTime = -999., tof = -999.;
363         Double_t c = TMath::C()*1.E-9;// m/ns
364
365         //cout<<"track->GetFilterMap()= "<<track->GetFilterMap()<<endl;
366         if(!track->TestFilterBit(fFilterBit)) continue;
367
368         //Float_t dcaXY = -999, dcaZ = -999;
369         p=track->P();
370         //pTPC=track->GetTPCmomentum();
371         pT=track->Pt();
372         phi=track->Phi();
373         eta=track->Eta();
374         dEdx=track->GetDetPid()->GetTPCsignal();
375          
376         Float_t probMis = fPIDResponse->GetTOFMismatchProbability(track);
377         if (probMis < 0.01) { //if u want to reduce mismatch using also TPC
378
379             //if ( (track->IsOn(AliAODTrack::kTOFin)) && (track->IsOn(AliAODTrack::kTIME)) &&  (track->IsOn(AliAODTrack::kTOFout))) {
380             if ( (track->IsOn(AliAODTrack::kITSin)) && (track->IsOn(AliAODTrack::kTOFpid)) ) {
381
382                 tofTime = track->GetTOFsignal();//in ps
383                 length = track->GetIntegratedLength();
384                 
385                 tof = tofTime*1E-3; // ns
386                 //cout<<"tof = "<<tof<<endl;
387                 if (tof <= 0)continue;
388                 //cout<<"length = "<<length<<endl;
389                 if (length <= 0) continue;
390                 
391                 length = length*0.01; // in meters
392                 tof = tof*c;
393                 beta = length/tof;
394                 
395                 TH2F *HistBetavsPTOFbeforePID = (TH2F*)fListQAInfo->At(6);
396                 HistBetavsPTOFbeforePID ->Fill(track->P()*track->Charge(),beta);
397             }//TOF signal
398
399             TH2F *HistdEdxVsPTPCbeforePID = (TH2F*)fListQAInfo->At(7);
400             HistdEdxVsPTPCbeforePID -> Fill(p*track->Charge(),dEdx); //TPC signal
401         }//probMis
402          
403         //QA plot
404         TH1F *HistPhiDistBefore = (TH1F*)fListQAInfo->At(8);
405         HistPhiDistBefore->Fill(phi);
406         //
407         TH1F *HistEtaDistBefore = (TH1F*)fListQAInfo->At(9);
408         HistEtaDistBefore->Fill(eta);
409          
410
411         if(pT<0.1) continue;
412         if(TMath::Abs(eta)>0.8) continue;
413         
414         Int_t TPCNcls = track->GetTPCNcls();
415             
416         if(TPCNcls<70 || dEdx<10) continue;
417         
418         // fill QA histograms
419
420         TH2F* HistDCAAfter =(TH2F*)fListQAInfo->At(10);
421         HistDCAAfter->Fill(dcaZ,dcaXY);
422         
423         TH1F *HistPhiDistAfter = (TH1F*)fListQAInfo->At(11);
424         HistPhiDistAfter->Fill(phi);
425                 
426         TH1F *HistEtaDistAfter = (TH1F*)fListQAInfo->At(12);
427         HistEtaDistAfter->Fill(eta);
428          
429          
430         Bool_t pWithinRange = kFALSE;
431         Int_t pRange = -999;
432         Double_t pBins[11] = {0.2,0.5,1,1.5,2,2.5,3,3.5,4,4.5,5};
433         for(int i=0;i<10;i++){
434             if(p>pBins[i] && p<pBins[i+1]){
435                 //cout<<"Inside if(p>pBins[i] && p<pBins[i+1])"<<endl;
436                 pWithinRange = kTRUE;
437                 pRange = i;
438             }
439         }
440         for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
441             //TOF nSigma
442             Double_t nSigmaTOF=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
443             Double_t nSigmaTPC=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
444                 
445             if(fPIDcuts && ispecie>1 && ispecie<5 && pWithinRange){// for pions, kaons and protons only
446
447
448                 if(fContourCut[ispecie-2][pRange]){cout<<"4) fContourCut exists"<<endl;}
449                     
450                 if(fContourCut[ispecie-2][pRange]->IsInside(nSigmaTOF,nSigmaTPC)){
451                     pass = kTRUE;
452                 }
453                 else{
454                     pass = kFALSE;
455                     continue;
456                 }
457             }
458                 
459             //TPC and TOF cuts, TPC TOF nsigma vs. momentum
460             if(pass){
461                 TH3 *hist1 = (TH3*)fListQAtpctof->At(ispecie);
462                 if (hist1){
463                     hist1->Fill(nSigmaTPC,nSigmaTOF,p);}
464                 
465                 TH3 *hist2 = (TH3*)fListQAtpctof->At(ispecie+AliPID::kSPECIESC);
466                 if (hist2){
467                     hist2->Fill(nSigmaTPC,nSigmaTOF,pT);}
468             }
469                 
470         }
471          
472
473     }//track loop
474     
475     TH2F *HistTPCvsGlobalMultAfter = (TH2F*) fListQAInfo->At(13);
476     HistTPCvsGlobalMultAfter->Fill(multGlobal,multTPC);
477     
478 }
479 //_________________________________________
480 void AliAnalysisTaskPIDconfig::CheckCentrality(AliVEvent* event, Bool_t &centralitypass)
481 {
482     // Check if event is within the set centrality range. Falls back to V0 centrality determination if no method is set
483     if (!fUseCentrality) AliFatal("No centrality method set! FATAL ERROR!");
484     Double_t centvalue = event->GetCentrality()->GetCentralityPercentile(fCentralityEstimator);
485      //cout << "Centrality evaluated-------------------------: " << centvalue <<endl;
486     if ((centvalue >= fCentralityPercentileMin) && (centvalue < fCentralityPercentileMax))
487     {
488         TH1F *hCentralityPass = (TH1F*)fListQAInfo->At(0);
489         hCentralityPass->Fill(centvalue);
490         //cout << "--------------Fill pass-------------------------"<<endl;
491         centralitypass = kTRUE;
492     }
493     
494 }
495 //______________________________________________________________________________
496 void AliAnalysisTaskPIDconfig::GetPIDContours()
497 {
498     if(fContourCutList){cout<<"+++++++++++++++++The contour file has been retrieved+++++++++++++++++++"<<endl;}
499     
500     TGraph *CutGraph[3][10];
501     //TCutG *cut[3][10];
502     Double_t phigh=0 ,plow=0;
503     Double_t pBinning[11] = {0,0.5,1,1.5,2,2.5,3,3.5,4,4.5,5};
504     TString species[3] = {"pion","kaon","proton"};
505     for(int i=0;i<3;i++){
506         TList *Species_contours = (TList*)fContourCutList->Get(species[i]);
507         if(Species_contours){cout<<"Species_contours exists"<<endl;}
508         for(int j=0;j<10;j++){
509             phigh = pBinning[j+1]*10;
510             plow = pBinning[j]*10;
511             TString Graph_Name = "contourlines_";
512             Graph_Name += species[i];
513             Graph_Name += Form("%.f%.f-%i%icent",plow,phigh,fCentralityPercentileMin,fCentralityPercentileMax);
514             cout<<Graph_Name<<endl;
515             CutGraph[i][j] = (TGraph*)Species_contours->FindObject(Graph_Name);
516             if(!CutGraph[i][j]){cout<<"Contour Graph does not exist"<<endl; continue;}
517             
518             fContourCut[i][j] = new TCutG(Graph_Name.Data(),CutGraph[i][j]->GetN(),CutGraph[i][j]->GetX(),CutGraph[i][j]->GetY());
519             if(!fContourCut[i][j]){cout<<"i,j "<<i<<" "<<j<<endl;}
520
521         }
522     }
523 }
524 //______________________________________________________________________________
525 void AliAnalysisTaskPIDconfig::SetupTPCTOFqa()
526 {
527     //
528     // Create the qa objects for TPC + TOF combination
529     
530
531     //TPC and TOF signal vs. momentum
532     for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
533         fhistNsigmaP = new TH3F(Form("NsigmaP_TPC_TOF_%s",AliPID::ParticleName(ispecie)),Form("TPC n#sigma vs. TOF n#sigma %s vs. p ;TPC n#sigma;TOF n#sigma;p [GeV]",AliPID::ParticleName(ispecie)),200,-20,20,200,-20,20,60,0.1,6);
534         fListQAtpctof->Add(fhistNsigmaP);
535     }
536     //TPC and TOF signal vs. transverse momentum
537     for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
538         fhistNsigmaPt = new TH3F(Form("NsigmaPt_TPC_TOF_%s",AliPID::ParticleName(ispecie)),Form("TPC n#sigma vs. TOF n#sigma %s vs. Pt ;TPC n#sigma;TOF n#sigma;pT [GeV]",AliPID::ParticleName(ispecie)),200,-20,20,200,-20,20,60,0.1,6);
539         fListQAtpctof->Add(fhistNsigmaPt);
540     }
541
542 }
543 //______________________________________________________________________________
544 void AliAnalysisTaskPIDconfig::SetupEventInfo()
545 {
546     //event and track info
547     
548     fhistCentralityPass = new TH1F("fcentralityPass","centralityPass", 100,0,100);
549     fListQAInfo->Add(fhistCentralityPass);
550
551     fNoEvents = new TH1F("number of events","no. of events",1,0,1);
552     fListQAInfo->Add(fNoEvents);
553
554     fpVtxZ = new TH1F("pVtxZ","pVtxZ",100,-20,20);
555     fListQAInfo->Add(fpVtxZ);
556     
557     fTPCvsGlobalMultBeforeOutliers = new TH2F("TPC vs. Global Multiplicity Before","TPC vs. Global Multiplicity Before",500,0,6000,500,0,6000);
558     fListQAInfo->Add(fTPCvsGlobalMultBeforeOutliers);
559     
560     fTPCvsGlobalMultAfterOutliers = new TH2F("TPC vs. Global Multiplicity After outliers","TPC vs. Global Multiplicity After outliers",500,0,6000,500,0,6000);
561     fListQAInfo->Add(fTPCvsGlobalMultAfterOutliers);
562     
563     fhistDCABefore = new TH2F("DCA xy vs z (before)","DCA before",200,0,10,200,0,10);
564     fListQAInfo->Add(fhistDCABefore);
565     
566     fHistBetavsPTOFbeforePID = new TH2F("momentum vs beta before PID","momentum vs beta before PID",1000,-10.,10.,1000,0,1.2);
567     fListQAInfo->Add(fHistBetavsPTOFbeforePID);
568     
569     fHistdEdxVsPTPCbeforePID = new TH2F("momentum vs dEdx before PID","momentum vs dEdx before PID",1000,-10.,10.,1000,0,1000);
570     fListQAInfo->Add(fHistdEdxVsPTPCbeforePID);
571     
572     fhistPhiDistBefore = new TH1F("Phi Distribution Before Cuts","Phi Distribution Before Cuts",200,0,6.4);
573     fListQAInfo->Add(fhistPhiDistBefore);
574
575     fhistEtaDistBefore = new TH1F("Eta Distribution Before Cuts","Eta Distribution Before Cuts",200,-10,10);
576     fListQAInfo->Add(fhistEtaDistBefore);
577     
578     fhistDCAAfter = new TH2F("DCA xy vs z (after)","DCA after",200,0,10,200,0,10);
579     fListQAInfo->Add(fhistDCAAfter);
580     
581     fhistPhiDistAfter = new TH1F("Phi Distribution After Cuts","Phi Distribution After Cuts",200,0,6.4);
582     fListQAInfo->Add(fhistPhiDistAfter);
583     
584     fhistEtaDistAfter = new TH1F("Eta Distribution After Cuts","Eta Distribution After Cuts",200,-10,10);
585     fListQAInfo->Add(fhistEtaDistAfter);
586     
587     fTPCvsGlobalMultAfter = new TH2F("TPC vs. Global Multiplicity After","TPC vs. Global Multiplicity After",500,0,6000,500,0,6000);
588     fListQAInfo->Add(fTPCvsGlobalMultAfter);
589     
590 //    fHistBetavsPTOFafterPID = new TH2F("momentum vs beta after PID","momentum vs beta after PID",1000,-10.,10.,1000,0,1.2);
591 //    fListQAInfo->Add(fHistBetavsPTOFafterPID);
592     
593 //    fHistdEdxVsPTPCafterPID = new TH2F("momentum vs dEdx after PID","momentum vs dEdx after PID",1000,-10.,10.,1000,0,1000);
594 //    fListQAInfo->Add(fHistdEdxVsPTPCafterPID);
595
596 }
597
598