update to AOD analysis
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / PiKaPr / TestAOD / AliSpectraAODTrackCuts.cxx
1
2 /**************************************************************************
3  * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
4  *                                                                        *
5  * Author: The ALICE Off-line Project.                                    *
6  * Contributors are mentioned in the code where appropriate.              *
7  *                                                                        *
8  * Permission to use, copy, modify and distribute this software and its   *
9  * documentation strictly for non-commercial purposes is hereby granted   *
10  * without fee, provided that the above copyright notice appears in all   *
11  * copies and that both the copyright notice and this permission notice   *
12  * appear in the supporting documentation. The authors make no claims     *
13  * about the suitability of this software for any purpose. It is          *
14  * provided "as is" without express or implied warranty.                  *
15  **************************************************************************/
16
17 //-----------------------------------------------------------------
18 //         AliSpectraAODTrackCuts class
19 //-----------------------------------------------------------------
20
21 #include "TChain.h"
22 #include "TTree.h"
23 #include "TLegend.h"
24 #include "TH1F.h"
25 #include "TH2F.h"
26 #include "TCanvas.h"
27 #include "AliAnalysisTask.h"
28 #include "AliAnalysisManager.h"
29 #include "AliAODTrack.h"
30 #include "AliAODMCParticle.h"
31 #include "AliAODEvent.h"
32 #include "AliAODInputHandler.h"
33 #include "AliAnalysisTaskESDfilter.h"
34 #include "AliAnalysisDataContainer.h"
35 #include "AliSpectraAODTrackCuts.h"
36 #include "AliSpectraAODHistoManager.h"
37 #include <iostream>
38
39 using namespace std;
40
41 ClassImp(AliSpectraAODTrackCuts)
42
43
44 AliSpectraAODTrackCuts::AliSpectraAODTrackCuts(const char *name) : TNamed(name, "AOD Track Cuts"), fIsSelected(0), fTrackBits(0), fEtaCut(0), fDCACut(0), fPCut(0), fPtCut(0),fQvecCutMin(0),fQvecCutMax(0), fHistoCuts(0), fTrack(0)
45
46 {
47    // Constructor
48    fHistoCuts = new TH1I("fTrkCuts", "Track Cuts", kNTrkCuts, -0.5, kNTrkCuts - 0.5);
49    fEtaCut = 100000.0; // default value of eta cut ~ no cut
50    fDCACut = 100000.0; // default value of dca cut ~ no cut
51    fPCut = 100000.0; // default value of p cut ~ no cut
52    fPtCut = 100000.0; // default value of pt cut ~ no cut 
53    fQvecCutMin = -100000.0; // default value of qvec cut ~ no cut 
54    fQvecCutMax = 100000.0; // default value of qvec cut ~ no cut 
55    
56 }
57
58 //_______________________________________________________
59 Bool_t AliSpectraAODTrackCuts::IsSelected(AliAODTrack * track)
60 {
61 // Returns true if Track Cuts are selected and applied
62    if (!track)
63    {
64       printf("ERROR: Could not receive track");
65       return kFALSE;
66    }
67    fTrack = track;
68    fIsSelected = (CheckTrackType() && CheckEtaCut() && CheckDCACut() && CheckPCut() && CheckPtCut() && CheckTOFMatching());
69    return fIsSelected ;
70 }
71 //_________________________________________________________
72
73 Bool_t AliSpectraAODTrackCuts::CheckTrackType()
74 {
75    // Check track cuts
76    if (fTrack->TestFilterBit(fTrackBits)) return kTRUE;
77    fHistoCuts->Fill(kTrkBit);
78    return kFALSE;
79 }
80 //________________________________________________________
81 Bool_t AliSpectraAODTrackCuts::CheckEtaCut()
82 {
83    // Check eta cut
84    if (fTrack->Eta() < fEtaCut && fTrack->Eta() > - fEtaCut) return kTRUE;
85    fHistoCuts->Fill(kTrkEta);
86    return kFALSE;
87 }
88 //_______________________________________________________
89 Bool_t AliSpectraAODTrackCuts::CheckDCACut()
90 {
91    // Check DCA cut
92    if (fTrack->DCA() < fDCACut) return kTRUE;
93    fHistoCuts->Fill(kTrkDCA);
94    return kFALSE;
95 }
96 //________________________________________________________
97 Bool_t AliSpectraAODTrackCuts::CheckPCut()
98 {
99    // Check P cut
100    if (fTrack->P() < fPCut) return kTRUE;
101    fHistoCuts->Fill(kTrkP);
102    return kFALSE;
103 }
104 //_______________________________________________________
105 Bool_t AliSpectraAODTrackCuts::CheckPtCut()
106 {
107     // check Pt cut
108 //    if ((fTrack->Pt() < fPtCut) && (fTrack->Pt() > 0.3 )) return kTRUE;
109    if (fTrack->Pt() < fPtCut) return kTRUE;
110     fHistoCuts->Fill(kTrkPt);
111     return kFALSE;
112 }
113
114 //_______________________________________________________
115 Bool_t AliSpectraAODTrackCuts::CheckTOFMatching()
116 {
117   // check Pt cut
118   //    if ((fTrack->Pt() < fPtCut) && (fTrack->Pt() > 0.3 )) return kTRUE;
119   
120   if (fTrack->Pt() < fPtCutTOFMatching) return kTRUE;
121   else{
122     UInt_t status; 
123     status=fTrack->GetStatus();
124     if((status&AliAODTrack::kTOFout)==0 || (status&AliAODTrack::kTIME)==0) return kFALSE; //tof matching and PID
125     return kTRUE;
126   }
127 }
128 //_______________________________________________________
129 void AliSpectraAODTrackCuts::PrintCuts() const
130 {
131   // Print cuts
132     cout << "Track Cuts" << endl;
133     cout << " > TrackBit\t" << fTrackBits << endl;
134     cout << " > Eta cut\t" << fEtaCut << endl;
135     cout << " > DCA cut\t" << fDCACut << endl;
136     cout << " > P cut\t" << fPCut << endl;
137     cout << " > Pt cut \t" << fPtCut << endl;
138     cout << " > Q vactor Min \t" << fQvecCutMin << endl;
139     cout << " > Q vactor Max \t" << fQvecCutMax << endl;
140 }
141 //_______________________________________________________
142 void AliSpectraAODTrackCuts::SetTrackType(UInt_t bit)
143 {
144    // Set the type of track to be used. The argument should be the bit number. The mask is produced automatically.
145    fTrackBits = (0x1 << (bit - 1));
146 }
147 //_______________________________________________________
148
149 Long64_t AliSpectraAODTrackCuts::Merge(TCollection* list)
150 {
151   // Merge a list of AliSpectraAODTrackCuts objects with this.
152   // Returns the number of merged objects (including this).
153
154   //  AliInfo("Merging");
155
156   if (!list)
157     return 0;
158
159   if (list->IsEmpty())
160     return 1;
161
162   TIterator* iter = list->MakeIterator();
163   TObject* obj;
164
165   // collections of all histograms
166   TList collections;
167
168   Int_t count = 0;
169
170   while ((obj = iter->Next())) {
171     AliSpectraAODTrackCuts* entry = dynamic_cast<AliSpectraAODTrackCuts*> (obj);
172     if (entry == 0) 
173       continue;
174
175     TH1I * histo = entry->GetHistoCuts();      
176     collections.Add(histo);
177     count++;
178   }
179   
180   fHistoCuts->Merge(&collections);
181   
182   delete iter;
183
184   return count+1;
185 }
186