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