01b935462cfdefd872e222d2908ffd5dfaa3959c
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / dNdPt / AlidNdPtAnalysisPbPbAOD.h
1 #ifndef AlidNdPtAnalysisPbPbAOD_H
2 #define AlidNdPtAnalysisPbPbAOD_H
3
4
5 //------------------------------------------------------------------------------
6 // AlidNdPtAnalysisPbPbAOD class used for dNdPt analysis in PbPb collision
7 // via AODs 
8 // 
9 // Author: P. Luettig, 15.05.2013
10 // last modified: 08.10.2013
11 //------------------------------------------------------------------------------
12
13 class iostream;
14
15 class TObject;
16 class TFile;
17 class TCint;
18 class THnSparse;
19
20 #include "AliAnalysisTaskSE.h"
21
22
23 #include "TList.h"
24 #include "TFile.h"
25 #include "TH1.h"
26 #include "TH2.h"
27 #include "TH3.h"
28 #include "THnSparse.h"
29 #include "THn.h"
30 #include "TClonesArray.h"
31
32 #include "TParticlePDG.h"
33 #include "TDatabasePDG.h"
34
35 #include "AliLog.h"
36 #include "AliCentrality.h"
37 #include "AliAODEvent.h"
38 #include "AliVEvent.h"
39
40 #include "AliInputEventHandler.h"
41 #include "AliAODInputHandler.h"
42 #include "AliAnalysisManager.h"
43 #include "AliMCEventHandler.h"
44 #include "AliAODMCHeader.h"
45 #include "AliAODMCParticle.h"
46 #include "AliGenHijingEventHeader.h"
47 #include "AliGenPythiaEventHeader.h"
48 #include "AliExternalTrackParam.h"
49 #include "AliESDtrack.h"
50
51 #include "TSystem.h"
52 #include "TROOT.h"
53
54 class AlidNdPtAnalysisPbPbAOD : public AliAnalysisTaskSE {
55   public :
56     enum CheckQuantity { cqCrossedRows = 0, cqNcluster = 1, cqChi = 2, cqLength = 3 };
57     enum KinematicQuantity { kqPt = 0, kqEta = 1, kqPhi = 2 };
58     enum MaxCheckQuantity { cqMax = 4 };
59     enum MaxKinematicQuantity { kqMax = 3 };
60     
61     AlidNdPtAnalysisPbPbAOD(const char *name = "dNdPtPbPbAOD");
62     ~AlidNdPtAnalysisPbPbAOD();
63     
64     virtual void UserCreateOutputObjects();
65     virtual void UserExec(Option_t *option);
66     virtual void Terminate(Option_t *);
67     
68     // Set binning for Histograms (if not set default binning is used)
69     void SetBinsMult(Int_t nbins, Double_t* edges)                      { Printf("[I] Setting Mult Bins"); fMultNbins = nbins; fBinsMult = GetArrayClone(nbins,edges); }
70     void SetBinsPt(Int_t nbins, Double_t* edges)                        { Printf("[I] Setting pT Bins"); fPtNbins = nbins; fBinsPt = GetArrayClone(nbins,edges); }
71     void SetBinsPtCorr(Int_t nbins, Double_t* edges)                    { Printf("[I] Setting pTcorr Bins"); fPtCorrNbins = nbins; fBinsPtCorr = GetArrayClone(nbins,edges); }
72     void SetBinsPtCheck(Int_t nbins, Double_t* edges)                   { Printf("[I] Setting pTcheck Bins"); fPtCheckNbins = nbins; fBinsPtCheck = GetArrayClone(nbins,edges); }
73     void SetBinsEta(Int_t nbins, Double_t* edges)                       { Printf("[I] Setting Eta Bins"); fEtaNbins = nbins; fBinsEta = GetArrayClone(nbins,edges); }
74     void SetBinsEtaCheck(Int_t nbins, Double_t* edges)                  { Printf("[I] Setting EtaCheck Bins"); fEtaCheckNbins = nbins; fBinsEtaCheck = GetArrayClone(nbins,edges); }
75     void SetBinsZv(Int_t nbins, Double_t* edges)                        { Printf("[I] Setting Zv Bins"); fZvNbins = nbins; fBinsZv= GetArrayClone(nbins,edges); }
76     void SetBinsCentrality(Int_t nbins, Double_t* edges)                { Printf("[I] Setting Cent Bins"); fCentralityNbins = nbins; fBinsCentrality = GetArrayClone(nbins,edges); }
77     void SetBinsPhi(Int_t nbins, Double_t* edges)                       { Printf("[I] Setting Phi Bins"); fPhiNbins = nbins; fBinsPhi = GetArrayClone(nbins,edges); }
78     
79     // set event cut variables
80     void SetCutMaxZVertex( Double_t d)                                  { fCutMaxZVertex = d; }
81     Double_t GetCutMaxZVertex()                                         { return fCutMaxZVertex; }
82     
83     // set track kinematic cut parameters
84     void SetCutPtRange(Double_t ptmin, Double_t ptmax)                  { fCutPtMin = ptmin; fCutPtMax = ptmax; }
85     Double_t GetCutPtMin()                                              { return fCutPtMin; }
86     Double_t GetCutPtMax()                                              { return fCutPtMax; }
87     
88     void SetCutEtaRange(Double_t etamin, Double_t etamax)               { fCutEtaMin = etamin; fCutEtaMax = etamax; }
89     Double_t GetCutEtaMin()                                             { return fCutEtaMin; }
90     Double_t GetCutEtaMax()                                             { return fCutEtaMax; }
91     
92     void EnableRelativeCuts()                                           { Printf("[I] Relative Cuts enabled"); fUseRelativeCuts = kTRUE; }
93     Bool_t AreRelativeCutsEnabled()                                     { return fUseRelativeCuts; }
94     
95     // setter and getter track quality cut parameters
96     void SetFilterBit(Int_t b)                                          { fFilterBit = b; };
97     Int_t GetFilterBit()                                                { return fFilterBit; }
98     
99     void SetCutRequireTPCRefit(Bool_t *b)                               { fCutRequireTPCRefit = b; } 
100     Bool_t IsTPCRefitRequired()                                         { return fCutRequireTPCRefit; } 
101     
102     void SetCutRequireITSRefit(Bool_t *b)                               { fCutRequireITSRefit = b; } 
103     Bool_t IsITSRefitRequired()                                         { return fCutRequireITSRefit; } 
104     
105     void SetCutMinNClustersTPC(Double_t d)                              { fCutMinNumberOfClusters = d; }
106     Double_t GetCutMinNClustersTPC()                                    { return fCutMinNumberOfClusters; }
107     
108     void SetCutPercMinNClustersTPC(Double_t d)                          { Printf("[I] Take only %.2f%% tracks with most clusters", d*100.); fCutPercMinNumberOfClusters = d; }
109     Double_t GetCutPercMinNClustersTPC()                                { return fCutPercMinNumberOfClusters; }
110     
111     void SetCutMinNCrossedRowsTPC(Double_t d)                           { fCutMinNumberOfCrossedRows = d; }    
112     Double_t GetCutMinNCrossedRowsTPC()                                 { return fCutMinNumberOfCrossedRows; }
113     
114     void SetCutPercMinNCrossedRowsTPC(Double_t d)                       { Printf("[I] Take only %.2f%% tracks with most crossedRows", d*100.); fCutPercMinNumberOfCrossedRows = d; }    
115     Double_t GetCutPercMinNCrossedRowsTPC()                             { return fCutPercMinNumberOfCrossedRows; }
116     
117     void SetCutMinRatioCrossedRowsOverFindableClustersTPC(Double_t d)   { fCutMinRatioCrossedRowsOverFindableClustersTPC = d; }
118     Double_t GetCutMinRatioCrossedRowsOverFindableClustersTPC()         { return fCutMinRatioCrossedRowsOverFindableClustersTPC; }
119     
120     void SetCutLengthInTPCPtDependent()                                 { fCutLengthInTPCPtDependent = kTRUE; }
121     Bool_t DoCutLengthInTPCPtDependent()                                { return fCutLengthInTPCPtDependent; }
122     
123     void SetPrefactorLengthInTPCPtDependent(Double_t d)                 { fPrefactorLengthInTPCPtDependent = d; }
124     Double_t GetPrefactorLengthInTPCPtDependent()                       { return fPrefactorLengthInTPCPtDependent; }
125      
126     void SetCutMaxChi2PerClusterTPC(Double_t d)                         { fCutMaxChi2PerClusterTPC = d; }
127     void SetCutMaxFractionSharedTPCClusters(Double_t d)                 { fCutMaxFractionSharedTPCClusters = d; }
128     void SetCutMaxDCAToVertexZ(Double_t d)                              { fCutMaxDCAToVertexZ = d; }
129     void SetCutMaxDCAToVertexXY(Double_t d)                             { fCutMaxDCAToVertexXY = d; }
130     void SetCutMaxChi2PerClusterITS(Double_t d)                         { fCutMaxChi2PerClusterITS = d; }
131     void SetCutDCAToVertex2D(Bool_t *b)                                 { fCutDCAToVertex2D = b; } 
132     void SetCutRequireSigmaToVertex(Bool_t *b)                          { fCutRequireSigmaToVertex = b; } 
133     void SetCutMaxDCAToVertexXYPtDep(Double_t d0, Double_t d1, Double_t d2)
134     {
135       fCutMaxDCAToVertexXYPtDepPar0 = d0;
136       fCutMaxDCAToVertexXYPtDepPar1 = d1;
137       fCutMaxDCAToVertexXYPtDepPar2 = d2;
138     }
139     void SetCutAcceptKinkDaughters(Bool_t *b)                           { fCutAcceptKinkDaughters = b; } 
140     void SetCutMaxChi2TPCConstrainedGlobal(Double_t d)                  { fCutMaxChi2TPCConstrainedGlobal = d; }
141        
142     // fill function for cross check histos
143     Bool_t FillDebugHisto(Double_t *dCrossCheckVar, Double_t *dKineVar, Double_t dCentrality, Bool_t bIsAccepted);
144     
145     // getter for DCA
146     Bool_t GetDCA(const AliAODTrack *track, AliAODEvent *evt, Double_t d0z0[2]);
147     
148     THnSparseF *GetHistZvPtEtaCent() const { return fZvPtEtaCent; }
149     TH1F *GetHistEventStatistics() const { return fEventStatistics; }
150     
151     const char * GetParticleName(Int_t pdg);
152     
153     AliGenHijingEventHeader* GetHijingEventHeader(AliAODMCHeader *header);
154     AliGenPythiaEventHeader* GetPythiaEventHeader(AliAODMCHeader *header);
155     
156     
157     Bool_t SetRelativeCuts(AliAODEvent *event);
158     
159     Bool_t IsTrackAccepted(AliAODTrack *tr, Double_t dCentrality, Double_t bMagZ);
160     Bool_t IsMCTrackAccepted(AliAODMCParticle *part);
161     
162     Bool_t IsHijingParticle(const AliAODMCParticle *part, AliGenHijingEventHeader* hijingGenHeader);
163     Bool_t IsPythiaParticle(const AliAODMCParticle *part, AliGenPythiaEventHeader* pythiaGenHeader);
164     
165     static Double_t* GetArrayClone(Int_t n, Double_t* source);
166     
167   private :
168     
169     // Output List
170     TList       *fOutputList;
171     
172     // Histograms
173     TH1F        *fPt; // simple pT histogramm
174     TH1F        *fMCPt; // simple pT truth histogramm
175     THnSparseF  *fZvPtEtaCent; //-> Zv:Pt:Eta:Cent
176     THnSparseF  *fPhiPtEtaCent; //-> Phi:Pt:Eta:Cent
177     THnSparseF  *fPtResptCent; //-> 1/pt:ResolutionPt:Cent
178     THnSparseF  *fMCRecPrimZvPtEtaCent; //-> MC Zv:Pt:Eta:Cent
179     THnSparseF  *fMCGenZvPtEtaCent; //-> MC Zv:Pt:Eta:Cent
180     THnSparseF  *fMCRecSecZvPtEtaCent; //-> MC Zv:Pt:Eta:Cent, only secondaries
181     THnSparseF  *fMCRecPrimPhiPtEtaCent; //-> MC Phi:Pt:Eta:Cent
182     THnSparseF  *fMCGenPhiPtEtaCent; //-> MC Phi:Pt:Eta:Cent
183     THnSparseF  *fMCRecSecPhiPtEtaCent; //-> MC Phi:Pt:Eta:Cent, only secondaries
184     TH1F        *fEventStatistics; // contains statistics of number of events after each cut
185     TH1F        *fEventStatisticsCentrality; // contains number of events vs centrality, events need to have a track in kinematic range
186     TH1F        *fMCEventStatisticsCentrality; // contains MC number of events vs centrality, events need to have a track in kinematic range
187     TH1F        *fAllEventStatisticsCentrality; // contains number of events vs centrality, events need to be triggered
188     TH2F        *fEventStatisticsCentralityTrigger; // contains number of events vs centrality in 1% bins vs trigger
189     THnSparseF  *fZvMultCent; // Zv:Mult:Cent
190     TH1F        *fTriggerStatistics; // contains number of events per trigger
191     TH1F        *fMCTrackPdgCode; // contains statistics of pdg codes of tracks
192     TH1F        *fMCTrackStatusCode; // contains statistics of status codes of tracks
193     TH1F        *fCharge; // charge distribution in data
194     TH1F        *fMCCharge; // charge distribution in MC
195     TH2F        *fMCPdgPt; // PDGvs PT for MC Particles
196     TH1F        *fMCHijingPrim; // number of particles, which are Hijing particles and primaries
197     THnSparseF  *fDCAPtAll; //control histo: DCAz:DCAxy:pT:eta:phi for all reconstructed tracks
198     THnSparseF  *fDCAPtAccepted; //control histo: DCAz:DCAxy:pT:eta:phi for all accepted reco tracks
199     THnSparseF  *fMCDCAPtSecondary; //control histo: DCAz:DCAxy:pT:eta:phi for all accepted reco track, which are secondaries (using MC info)
200     THnSparseF  *fMCDCAPtPrimary; //control histo: DCAz:DCAxy:pT:eta:phi for all accepted reco track, which are primaries (using MC info)
201     THnF        *fCrossCheckAll[4]; //control histo: {CrossedRows,Ncluster,Chi} vs pT,eta,phi,Centrality for all tracks
202     THnF        *fCrossCheckAcc[4]; //control histo: {CrossedRows,Ncluster,Chi} vs pT,eta,phi,Centrality after cuts
203     TH1F        *fCutPercClusters; // control histo: number of clusters, where the relative cut has been set e-by-e
204     TH1F        *fCutPercCrossed; // control histo: number of crossed rows, where the relative cut has been set e-by-e
205     TH2F        *fCrossCheckRowsLength; // control histo: number of crossed rows vs length in TPC
206     TH2F        *fCrossCheckClusterLength; // control histo: number of clusters vs length in TPC
207     TH2F        *fCrossCheckRowsLengthAcc; // control histo: number of crossed rows vs length in TPC for all accepted tracks
208     TH2F        *fCrossCheckClusterLengthAcc; // control histo: number of clusters vs length in TPC for all accepted tracks
209     
210     
211     // global variables
212     Bool_t fIsMonteCarlo;
213     
214     // event cut variables
215     Double_t fCutMaxZVertex;
216     
217     // track kinematic cut variables
218     Double_t fCutPtMin;
219     Double_t fCutPtMax;
220     Double_t fCutEtaMin;
221     Double_t fCutEtaMax;
222     
223     // track quality cut variables
224     Int_t       fFilterBit;
225     Bool_t      fUseRelativeCuts;
226     Bool_t      fCutRequireTPCRefit;
227     Bool_t      fCutRequireITSRefit;
228     Double_t    fCutMinNumberOfClusters;
229     Double_t    fCutPercMinNumberOfClusters;
230     Double_t    fCutMinNumberOfCrossedRows;
231     Double_t    fCutPercMinNumberOfCrossedRows;
232     Double_t    fCutMinRatioCrossedRowsOverFindableClustersTPC;
233     Double_t    fCutMaxChi2PerClusterTPC;
234     Double_t    fCutMaxFractionSharedTPCClusters;
235     Double_t    fCutMaxDCAToVertexZ;
236     Double_t    fCutMaxDCAToVertexXY;
237     Double_t    fCutMaxChi2PerClusterITS;
238     Bool_t      fCutDCAToVertex2D;
239     Bool_t      fCutRequireSigmaToVertex;
240     Double_t    fCutMaxDCAToVertexXYPtDepPar0;
241     Double_t    fCutMaxDCAToVertexXYPtDepPar1;
242     Double_t    fCutMaxDCAToVertexXYPtDepPar2;
243     Bool_t      fCutAcceptKinkDaughters;
244     Double_t    fCutMaxChi2TPCConstrainedGlobal;
245     Bool_t      fCutLengthInTPCPtDependent;
246     Double_t    fPrefactorLengthInTPCPtDependent;
247     
248     //binning for THNsparse
249     Int_t fMultNbins;
250     Int_t fPtNbins;
251     Int_t fPtCorrNbins;
252     Int_t fPtCheckNbins;
253     Int_t fEtaNbins;
254     Int_t fEtaCheckNbins;
255     Int_t fZvNbins;
256     Int_t fCentralityNbins;
257     Int_t fPhiNbins;
258     Double_t* fBinsMult; //[fMultNbins]
259     Double_t* fBinsPt; //[fPtNbins]
260     Double_t* fBinsPtCorr; //[fPtCorrNbins]
261     Double_t* fBinsPtCheck; //[fPtCheckNbins]
262     Double_t* fBinsEta; //[fEtaNbins]
263     Double_t* fBinsEtaCheck; //[fEtaCheckNbins]
264     Double_t* fBinsZv; //[fZvNbins]
265     Double_t* fBinsCentrality; //[fCentralityNbins]
266     Double_t* fBinsPhi; //[fPhiNbins]
267     
268     AlidNdPtAnalysisPbPbAOD(const AlidNdPtAnalysisPbPbAOD&); // not implemented
269     AlidNdPtAnalysisPbPbAOD& operator=(const AlidNdPtAnalysisPbPbAOD&); // not implemented  
270     
271     ClassDef(AlidNdPtAnalysisPbPbAOD,4); // has to be at least 1, otherwise not streamable...
272 };
273
274 #endif