591824176e69c0a1b71ac8bb1d9983b7082badb9
[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     //Dummy Constructor
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     DefineInput(0,TChain::Class());
175     DefineOutput(1,TList::Class());
176 }
177
178 //_____________________________________________________________________
179 AliAnalysisTaskPIDconfig::~AliAnalysisTaskPIDconfig()
180 {
181     //Destructor
182     
183  //   delete fPID;
184   //  delete fPIDqa;
185  //   delete fTrackCuts;
186    // delete fSparseSpecies;
187     //delete []fvalueSpecies;
188
189
190 }
191 //______________________________________________________________________
192 void AliAnalysisTaskPIDconfig::UserCreateOutputObjects()
193 {
194     //
195     // Create the output QA objects
196     //
197     
198     AliLog::SetClassDebugLevel("AliAnalysisTaskPIDconfig",10);
199     
200     //input hander
201     AliAnalysisManager *mgr=AliAnalysisManager::GetAnalysisManager();
202     AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(mgr->GetInputEventHandler());
203     if (!inputHandler) AliFatal("Input handler needed");
204     
205     //pid response object
206     fPIDResponse=inputHandler->GetPIDResponse();
207     if (!fPIDResponse) AliError("PIDResponse object was not created");
208     
209     //
210     fListQA=new TList;
211     fListQA->SetOwner();
212     
213     fListQAtpctof=new TList;
214     fListQAtpctof->SetOwner();
215     fListQAtpctof->SetName("PID_TPC_TOF");
216     
217     fListQAInfo=new TList;
218     fListQAInfo->SetOwner();
219     fListQAInfo->SetName("Event_Track_Info");
220     
221     fListQA->Add(fListQAtpctof);
222     fListQA->Add(fListQAInfo);
223     
224     SetupTPCTOFqa();
225     SetupEventInfo();
226     
227     PostData(1,fListQA);
228 }
229 //______________________________________________________________________
230 void AliAnalysisTaskPIDconfig::UserExec(Option_t*){
231     //Main loop
232     //Called for each event
233     
234     // create pointer to event
235     fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
236     fESD = dynamic_cast<AliESDEvent*>(InputEvent());
237     
238     Int_t ntracks=fAOD->GetNumberOfTracks();
239
240
241     if(!(fESD || fAOD)){
242         printf("ERROR: fESD & fAOD not available\n");
243         return;
244     }
245     fVevent = dynamic_cast<AliVEvent*>(InputEvent());
246     if (!fVevent) {
247         printf("ERROR: fVevent not available\n");
248         return;
249     }
250     
251     Bool_t pass = kFALSE;
252     
253     CheckCentrality(fVevent,pass);
254     
255     if(!pass){ return;}
256     
257     const AliVVertex *pVtx = fVevent->GetPrimaryVertex();
258     
259     Double_t pVtxZ = -999;
260     pVtxZ = pVtx->GetZ();
261     
262     if(TMath::Abs(pVtxZ)>10) return;
263
264     TH1F *hNoEvents = (TH1F*)fListQAInfo->At(1);
265     TH1F *histpVtxZ = (TH1F*)fListQAInfo->At(2);
266
267     if(hNoEvents) hNoEvents->Fill(0);
268     if(histpVtxZ) histpVtxZ->Fill(pVtxZ);
269
270     if(ntracks<2) return;
271     
272    // if(!pass) return;
273     
274     TH2F *HistTPCvsGlobalMultBeforeOutliers = (TH2F*)fListQAInfo->At(3);
275
276     TH2F *HistTPCvsGlobalMultAfterOutliers = (TH2F*)fListQAInfo->At(4);
277
278     
279     Float_t multTPC(0.); // tpc mult estimate
280     Float_t multGlobal(0.); // global multiplicity
281     
282     const Int_t nGoodTracks = fVevent->GetNumberOfTracks();
283     if(!fData2011) { // cut on outliers
284       
285         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill tpc mult
286             AliAODTrack* AODtrack =dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
287             if (!AODtrack) continue;
288             if (!(AODtrack->TestFilterBit(1))) continue;
289             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;
290             multTPC++;
291         }//track loop
292
293         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill global mult
294             AliAODTrack *AODtrack=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
295             if (!AODtrack) continue;
296             if (!(AODtrack->TestFilterBit(16))) continue;
297             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;
298             Double_t b[2] = {-99., -99.};
299             Double_t bCov[3] = {-99., -99., -99.};
300             AliAODTrack copy(*AODtrack);
301             if (!(copy.PropagateToDCA(fVevent->GetPrimaryVertex(), fVevent->GetMagneticField(), 100., b, bCov))) continue;
302             if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
303             multGlobal++;
304         } //track loop
305         
306         HistTPCvsGlobalMultBeforeOutliers->Fill(multGlobal,multTPC);
307
308         if(multTPC < (-40.3+1.22*multGlobal) || multTPC > (32.1+1.59*multGlobal)){ pass = kFALSE;}
309         
310         if(!pass) return;
311         HistTPCvsGlobalMultAfterOutliers->Fill(multGlobal,multTPC);
312
313     }
314     
315     
316     if(fData2011) { // cut on outliers
317         //Float_t multTPC(0.); // tpc mult estimate
318         //Float_t multGlob(0.); // global multiplicity
319         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill tpc mult
320             AliAODTrack *AODtrack=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
321             if (!AODtrack) continue;
322             if (!(AODtrack->TestFilterBit(1))) continue;
323             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;
324             multTPC++;
325         }
326         for(Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++) { // fill global mult
327             AliAODTrack *AODtrack=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
328             if (!AODtrack) continue;
329             if (!(AODtrack->TestFilterBit(16))) continue;
330             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;
331             Double_t b[2] = {-99., -99.};
332             Double_t bCov[3] = {-99., -99., -99.};
333             AliAODTrack copy(*AODtrack);
334             if (!(copy.PropagateToDCA(fVevent->GetPrimaryVertex(), fVevent->GetMagneticField(), 100., b, bCov))) continue;
335             if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
336             multGlobal++;
337             
338         } //track loop
339         
340         HistTPCvsGlobalMultBeforeOutliers->Fill(multGlobal,multTPC);
341
342         if(multTPC < (-36.73 + 1.48*multGlobal) || multTPC > (62.87 + 1.78*multGlobal)){pass = kFALSE;}
343         
344         if(!pass) return;
345         HistTPCvsGlobalMultAfterOutliers->Fill(multGlobal,multTPC);
346
347     }
348  
349
350     
351      for(Int_t itrack = 0; itrack < ntracks; itrack++){
352
353         AliAODTrack *track=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(itrack));
354         if(!track) continue;
355         
356         Float_t dcaXY = track->DCA();
357         Float_t dcaZ  = track->ZAtDCA();
358          
359         TH2F* HistDCAbefore =(TH2F*)fListQAInfo->At(5);
360         HistDCAbefore->Fill(dcaZ,dcaXY);
361         
362         Double_t p = -999, pTPC = -999, pT = -999, phi = -999, eta = -999, dEdx =-999;
363         Double_t length = -999., beta =-999, tofTime = -999., tof = -999.;
364         Double_t c = TMath::C()*1.E-9;// m/ns
365
366         //cout<<"track->GetFilterMap()= "<<track->GetFilterMap()<<endl;
367         if(!track->TestFilterBit(fFilterBit)) continue;
368
369         //Float_t dcaXY = -999, dcaZ = -999;
370         p=track->P();
371         pTPC=track->GetTPCmomentum();
372         pT=track->Pt();
373         phi=track->Phi();
374         eta=track->Eta();
375         dEdx=track->GetDetPid()->GetTPCsignal();
376             
377         if ( (track->IsOn(AliAODTrack::kTOFin)) &&
378         (track->IsOn(AliAODTrack::kTIME)) &&  (track->IsOn(AliAODTrack::kTOFout))) {
379       //  if ( (track->IsOn(AliAODTrack::kTOFin)) &&
380       //              (track->IsOn(AliAODTrack::kTOFout))  ) {
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          
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             Bool_t pWithinRange = kFALSE;
430             Int_t pRange = -999;
431             TCutG *cut[3][10];
432             if(fPIDcuts){
433                 TGraph *ContourCut[3][10];
434                 Double_t plow[10] = {0.2,0.5,1,1.5,2,2.5,3,3.5,4,4.5};
435                 Double_t phigh[10] = {0.5,1,1.5,2,2.5,3,3.5,4,4.5,5};
436                // TString species[3] = {pion,kaon,proton};
437                 for(int i=0;i<3;i++){
438                     for(int j=0;j<10;j++){
439                         if(p>plow[j] && p<phigh[j]){
440                             pWithinRange = kTRUE;
441                             pRange = j;
442                             TList *Scontours = (TList*)fContourCutList->At(i);
443                             TList *Pcontours = (TList*)Scontours->At(j);
444                             if (!Pcontours || !Scontours) return;
445                             ContourCut[i][j] = (TGraph*)Pcontours->First();
446                             cut[i][j] = new TCutG("cut",ContourCut[i][j]->GetN(),ContourCut[i][j]->GetX(),ContourCut[i][j]->GetY());
447                             //ContourCut[i][j] = (TGraph *)fContourCutList->At(10*i+j);
448                         }
449                     }
450                 }
451             }
452
453             for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
454                 //TOF nSigma
455                 Double_t nSigmaTOF=fPIDResponse->NumberOfSigmasTOF(track, (AliPID::EParticleType)ispecie);
456                 Double_t nSigmaTPC=fPIDResponse->NumberOfSigmasTPC(track, (AliPID::EParticleType)ispecie);
457                 if(fPIDcuts && ispecie>1 && ispecie<5 && pWithinRange){// for pions, kaons and protons only
458                     if(cut[ispecie-2][pRange]->IsInside(nSigmaTOF,nSigmaTPC)){
459                         pass = kTRUE;
460                     }
461                     else{
462                         pass = kFALSE;
463                         continue;
464                     }
465                 }
466                 
467                 //TPC and TOF cuts, TPC TOF nsigma vs. momentum
468                 if(pass){
469                     TH3 *hist1 = (TH3*)fListQAtpctof->At(ispecie);
470                     if (hist1){
471                         hist1->Fill(nSigmaTPC,nSigmaTOF,p);}
472                 
473                     TH3 *hist2 = (TH3*)fListQAtpctof->At(ispecie+AliPID::kSPECIESC);
474                     if (hist2){
475                         hist2->Fill(nSigmaTPC,nSigmaTOF,pT);}
476                 }
477                 
478             }
479          
480
481      }//track loop
482     
483     TH2F *HistTPCvsGlobalMultAfter = (TH2F*) fListQAInfo->At(13);
484     HistTPCvsGlobalMultAfter->Fill(multGlobal,multTPC);
485     
486 }
487 //_________________________________________
488 void AliAnalysisTaskPIDconfig::CheckCentrality(AliVEvent* event, Bool_t &centralitypass)
489 {
490     // Check if event is within the set centrality range. Falls back to V0 centrality determination if no method is set
491     if (!fUseCentrality) AliFatal("No centrality method set! FATAL ERROR!");
492     Double_t centvalue = event->GetCentrality()->GetCentralityPercentile(fCentralityEstimator);
493      //cout << "Centrality evaluated-------------------------: " << centvalue <<endl;
494     if ((centvalue >= fCentralityPercentileMin) && (centvalue < fCentralityPercentileMax))
495     {
496         TH1F *hCentralityPass = (TH1F*)fListQAInfo->At(0);
497         hCentralityPass->Fill(centvalue);
498         //cout << "--------------Fill pass-------------------------"<<endl;
499         centralitypass = kTRUE;
500     }
501     
502 }
503 //______________________________________________________________________________
504 void AliAnalysisTaskPIDconfig::SetupTPCTOFqa()
505 {
506     //
507     // Create the qa objects for TPC + TOF combination
508     
509
510     //TPC and TOF signal vs. momentum
511     for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
512         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);
513         fListQAtpctof->Add(fhistNsigmaP);
514     }
515     //TPC and TOF signal vs. transverse momentum
516     for (Int_t ispecie=0; ispecie<AliPID::kSPECIESC; ++ispecie){
517         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);
518         fListQAtpctof->Add(fhistNsigmaPt);
519     }
520
521 }
522 //______________________________________________________________________________
523 void AliAnalysisTaskPIDconfig::SetupEventInfo()
524 {
525     //event and track info
526     
527     fhistCentralityPass = new TH1F("fcentralityPass","centralityPass", 100,0,100);
528     fListQAInfo->Add(fhistCentralityPass);
529
530     fNoEvents = new TH1F("number of events","no. of events",1,0,1);
531     fListQAInfo->Add(fNoEvents);
532
533     fpVtxZ = new TH1F("pVtxZ","pVtxZ",100,-20,20);
534     fListQAInfo->Add(fpVtxZ);
535     
536     fTPCvsGlobalMultBeforeOutliers = new TH2F("TPC vs. Global Multiplicity Before","TPC vs. Global Multiplicity Before",500,0,6000,500,0,6000);
537     fListQAInfo->Add(fTPCvsGlobalMultBeforeOutliers);
538     
539     fTPCvsGlobalMultAfterOutliers = new TH2F("TPC vs. Global Multiplicity After outliers","TPC vs. Global Multiplicity After outliers",500,0,6000,500,0,6000);
540     fListQAInfo->Add(fTPCvsGlobalMultAfterOutliers);
541     
542     fhistDCABefore = new TH2F("DCA xy vs z (before)","DCA before",200,0,10,200,0,10);
543     fListQAInfo->Add(fhistDCABefore);
544     
545     fHistBetavsPTOFbeforePID = new TH2F("momentum vs beta before PID","momentum vs beta before PID",1000,-10.,10.,1000,0,1.2);
546     fListQAInfo->Add(fHistBetavsPTOFbeforePID);
547     
548     fHistdEdxVsPTPCbeforePID = new TH2F("momentum vs dEdx before PID","momentum vs dEdx before PID",1000,-10.,10.,1000,0,1000);
549     fListQAInfo->Add(fHistdEdxVsPTPCbeforePID);
550     
551     fhistPhiDistBefore = new TH1F("Phi Distribution Before Cuts","Phi Distribution Before Cuts",200,0,6.4);
552     fListQAInfo->Add(fhistPhiDistBefore);
553
554     fhistEtaDistBefore = new TH1F("Eta Distribution Before Cuts","Eta Distribution Before Cuts",200,-10,10);
555     fListQAInfo->Add(fhistEtaDistBefore);
556     
557     fhistDCAAfter = new TH2F("DCA xy vs z (after)","DCA after",200,0,10,200,0,10);
558     fListQAInfo->Add(fhistDCAAfter);
559     
560     fhistPhiDistAfter = new TH1F("Phi Distribution After Cuts","Phi Distribution After Cuts",200,0,6.4);
561     fListQAInfo->Add(fhistPhiDistAfter);
562     
563     fhistEtaDistAfter = new TH1F("Eta Distribution After Cuts","Eta Distribution After Cuts",200,-10,10);
564     fListQAInfo->Add(fhistEtaDistAfter);
565     
566     fTPCvsGlobalMultAfter = new TH2F("TPC vs. Global Multiplicity After","TPC vs. Global Multiplicity After",500,0,6000,500,0,6000);
567     fListQAInfo->Add(fTPCvsGlobalMultAfter);
568     
569 //    fHistBetavsPTOFafterPID = new TH2F("momentum vs beta after PID","momentum vs beta after PID",1000,-10.,10.,1000,0,1.2);
570 //    fListQAInfo->Add(fHistBetavsPTOFafterPID);
571     
572 //    fHistdEdxVsPTPCafterPID = new TH2F("momentum vs dEdx after PID","momentum vs dEdx after PID",1000,-10.,10.,1000,0,1000);
573 //    fListQAInfo->Add(fHistdEdxVsPTPCafterPID);
574
575 }
576
577