Adding Correlation analysis for Leading V0 particle AliLeadingV0Correlation (Sandun...
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliLeadingV0Correlation.h
1
2 //            Class for Leading Charged Track+V0 Correlations Analysis
3
4 #ifndef ALILEADINGV0CORRELATIONH
5 #define ALILEADINGV0CORRELATIONH
6
7 #include "AliAnalysisTask.h"
8 #include "AliUEHist.h"
9 #include "TString.h"
10 #include "AliVParticle.h"
11 #include "AliLog.h"
12 #include "AliPID.h"
13
14 class TList;
15 class TH2;
16 class AliAODEvent;
17 class AliEventPoolManager;
18 class AliEventPool;
19 class AliAnalyseLeadingTrackUE;
20 class AliVParticle;
21 class AliPIDResponse;
22 class AliPID;
23 class AliAODv0;
24
25
26 #ifndef ALIANALYSISTASKSEH
27 #include "AliAnalysisTaskSE.h"
28 #endif
29
30 //---------------------------------------------------------------------------------------
31 class AliLeadingV0Correlation : public AliAnalysisTaskSE {
32 public:
33    AliLeadingV0Correlation();
34    AliLeadingV0Correlation(const char *name);
35    virtual ~AliLeadingV0Correlation();
36
37    virtual void     UserCreateOutputObjects();
38    virtual void     UserExec(Option_t *option);
39    virtual void     Terminate(Option_t *);
40         
41         void SetMaxNEventsInPool(Int_t events){fPoolMaxNEvents=events;}
42         void SetMinNTracksInPool(Int_t tracks){fPoolMinNTracks=tracks;}
43         void SetMinEventsToMix(Int_t events){fMinEventsToMix=events;}
44
45         void SetPoolPVzBinLimits(Int_t Nzvtxbins,const Double_t *ZvtxBins){
46                 fNzVtxBins = Nzvtxbins;
47                 for(int ix = 0;ix<fNzVtxBins+1;ix++){fZvtxBins[ix] = ZvtxBins[ix];}
48         }
49         
50         void SetPoolCentBinLimits(Int_t Ncentbins,const Double_t *CentBins){
51                 fNCentBins = Ncentbins;
52                 for(int ix = 0;ix<fNCentBins+1;ix++){fCentBins[ix] = CentBins[ix];}
53         }
54         
55         void SetCollidingSystem(TString system){fcollidingSys = system;}
56         void SetPrimeryVertexCut(Double_t pvzcut){fpvzcut = pvzcut;}
57         void SetEatCut(Double_t  TrackEtaCut){fTrackEtaCut = TrackEtaCut;}
58         void SetFilterBit(UInt_t  filterBit){fFilterBit = filterBit;}
59         void SetTrigPtBinLimits(Double_t trigPtLow,Double_t trigPtHigh){
60                 ftrigPtLow = trigPtLow;
61                 ftrigPtHigh = trigPtHigh;
62         }
63         void SetAssocPtBinLimits(Double_t assocPtLow,Double_t assocPtHigh){
64                 fassocPtLow = assocPtLow;
65                 fassocPtHigh = assocPtHigh;
66         }
67         void SetMCAnalysis(Bool_t aAnalysisMC){fAnalysisMC=aAnalysisMC;}
68         void SetCutV0Radius(Double_t aCutV0Radius){fCutV0Radius=aCutV0Radius;}
69         void SetCutDCANToP(Double_t aCutDCANegToPV){fCutDCANegToPV=aCutDCANegToPV;}
70         void SetCutDCAPToP(Double_t aCutDCAPosToPV){fCutDCAPosToPV=aCutDCAPosToPV;}
71         void SetCutDCADaughters(Double_t aCutDCAV0Daughters){fCutDCAV0Daughters=aCutDCAV0Daughters;}
72         void SetCutCPA(Double_t aCutV0CosPA){fCutV0CosPA=aCutV0CosPA;}
73         void SetCutArmenterosK0s(Double_t aSpecialArmenterosCutK0s){fSpecialArmenterosCutK0s=aSpecialArmenterosCutK0s;}
74         
75         
76 private:
77         AliVParticle* ParticleWithCuts(TObject* obj, Int_t ipart);
78         TObjArray* CloneAndReduceTrackList(TObjArray* tracks);
79         TObjArray* FindLeadingObjectsK0(TObject *obj);
80         TObjArray* FindLeadingObjectsLambda(TObject *obj);
81         TObjArray* FindLeadingObjectsK0MC(TObject *obj);
82         TObjArray* FindLeadingObjectsLambdaMC(TObject *obj);
83         Int_t  NParticles(TObject* obj);
84         Double_t RangePhi(Double_t DPhi);
85         void FillCorrelations(TObjArray* particles, TObjArray* mixed,TH2F*histo);
86         void QSortTracks(TObjArray &a, Int_t first, Int_t last);
87         Bool_t IsAcseptedPrimaryTrack(const AliAODTrack *itrack);
88         Bool_t IsAcseptedDaughterTrack(const AliAODTrack *itrack);
89         Bool_t IsK0InvMass(const Double_t mass)const;
90         Bool_t IsLambdaInvMass(const Double_t mass) const;
91         Bool_t IsTrackNotFromLambda(const Int_t indexTrack);
92         Bool_t IsTrackNotFromK0(const Int_t indexTrack);
93         Bool_t IsAcseptedV0(const AliAODEvent*aod,const AliAODv0* aodV0, const AliAODTrack* myTrackPos, const AliAODTrack* myTrackNeg);
94         Double_t IsAccepteddEdx(const Double_t mom, const Double_t dEdx, AliPID::EParticleType n);
95
96     
97         AliAODEvent              * fAODEvent;        //! AOD Event
98         AliEventPoolManager      * fPoolMgr;         //! event pool manager
99         AliEventPool             * fPoolK0;          //! Pool for event mixing
100         AliEventPool             * fPoolLambda;      //! Pool for event mixing
101         AliPIDResponse           * fPIDResponse;     //  PID response
102         
103         Int_t fPoolMaxNEvents;                       // set maximum number of events in the pool
104         Int_t fPoolMinNTracks;                       // set minimum number of tracks in the pool
105         Int_t fMinEventsToMix;                       // set the minimum number of events want to mix
106         Int_t fNzVtxBins;                            // number of z vrtx bins
107         Double_t fZvtxBins[100];                     // [fNzVtxBinsDim]
108         Int_t fNCentBins;                            // number of centrality bins
109         Double_t fCentBins[100];                     // [fNCentBinsDim]
110         
111         TString         fcollidingSys;               // "PP" or "PbPb"
112         Double_t        fpvzcut;                     // PVz cut of event
113     Double_t            fTrackEtaCut;                // Eta cut on particles
114     UInt_t              fFilterBit;                  // Select tracks from an specific track cut (default 0xFF all track selected)
115         Double_t        ftrigPtLow;
116         Double_t        ftrigPtHigh;
117         Double_t        fassocPtLow;
118         Double_t        fassocPtHigh;
119         Bool_t           fAnalysisMC;
120         TString         fUsePID;
121         Double_t        fRapidityCut;                // Rapidity cut V0
122         
123         
124         //--- 5 Topological Selections
125         Double_t        fCutV0Radius;
126         Double_t        fCutDCANegToPV;
127         Double_t        fCutDCAPosToPV;
128         Double_t        fCutDCAV0Daughters;
129         Double_t        fCutV0CosPA;
130         Double_t        fSpecialArmenterosCutK0s;
131         
132         //Perform multiplicity selection 
133         // --- (mult estimator in pp, centrality in PbPb)
134         Bool_t fPerformMultiplicityStudy; 
135         Int_t fLoMultBound;
136         Int_t fHiMultBound;
137         
138         
139         TList       * fOutputList;           // Output list
140         // MC histograms
141         TH1F        *fHistMCPrimaryVertexX;       
142         TH1F        *fHistMCPrimaryVertexY;       
143         TH1F        *fHistMCPrimaryVertexZ;  
144         TH2F        *fHistMCtracksProdRadiusK0s;        
145         TH2F        *fHistMCtracksProdRadiusLambda;        
146         TH2F        *fHistMCtracksProdRadiusAntiLambda;
147         TH1F        *fHistMCtracksDecayRadiusK0s;        
148         TH1F        *fHistMCtracksDecayRadiusLambda;        
149         TH1F        *fHistMCtracksDecayRadiusAntiLambda;
150         TH1F        *fHistMCPtAllK0s;
151         TH1F        *fHistMCPtAllLambda;
152         TH1F        *fHistMCPtAllAntiLambda;
153         TH1F        *fHistMCProdRadiusK0s;
154         TH1F        *fHistMCProdRadiusLambda;
155         TH1F        *fHistMCProdRadiusAntiLambda;
156         TH1F        *fHistMCRapK0s;                  
157         TH1F        *fHistMCRapLambda;
158         TH1F        *fHistMCRapAntiLambda;           
159         TH1F        *fHistMCPtK0s;        
160         TH1F        *fHistMCPtLambda;        
161         TH1F        *fHistMCPtAntiLambda;
162         TH1F        *fHistMCPtLambdaFromSigma;        
163         TH1F        *fHistMCPtAntiLambdaFromSigma; 
164         TH1F        *fHistNTimesRecK0s;        
165         TH1F        *fHistNTimesRecLambda;        
166         TH1F        *fHistNTimesRecAntiLambda;        
167         TH2F        *fHistNTimesRecK0sVsPt;        
168         TH2F        *fHistNTimesRecLambdaVsPt;        
169         TH2F        *fHistNTimesRecAntiLambdaVsPt; 
170         TH1F        *fHistMCDaughterTrack;
171         TH2F        *fHistPrimRawPtVsYK0s;
172         TH2F        *fHistPrimRawPtVsYLambda;
173         TH2F        *fHistPrimRawPtVsYAntiLambda;
174         
175         TH1F        *fHistPrimaryVertexX;
176         TH1F        *fHistPrimaryVertexY;
177         TH1F        *fHistPrimaryVertexZ;
178         ///////////////////////////K0s 2D histos: cut vs on fly status/////////////////
179         
180         TH2F        *fHistDcaPosToPrimVertexK0;        
181         TH2F        *fHistDcaNegToPrimVertexK0;        
182         TH2F        *fHistRadiusV0K0;        
183         TH2F        *fHistDecayLengthV0K0;        
184         TH2F        *fHistDcaV0DaughtersK0;        
185         TH2F        *fHistChi2K0;        
186         TH2F        *fHistCosPointAngleK0;        
187         
188         ///////////////////////////K0s 2D histos: cut vs mass//////////////
189         TH2F        *fHistDcaPosToPrimVertexK0vsMassK0;   
190         TH2F        *fHistDcaNegToPrimVertexK0vsMassK0;   
191         TH2F        *fHistRadiusV0K0vsMassK0;             
192         TH2F        *fHistDecayLengthV0K0vsMassK0;        
193         TH2F        *fHistDcaV0DaughtersK0vsMassK0;       
194         TH2F        *fHistCosPointAngleK0vsMassK0;        
195         
196         //////////////////////////Lambda 2D histos: cut vs on fly status////////////////////
197         
198         TH2F        *fHistDcaPosToPrimVertexL;        
199         TH2F        *fHistDcaNegToPrimVertexL;        
200         TH2F        *fHistRadiusV0L;        
201         TH2F        *fHistDecayLengthV0L;        
202         TH2F        *fHistDcaV0DaughtersL;        
203         TH2F        *fHistChi2L;        
204         TH2F        *fHistCosPointAngleL;        
205         TH1F        *fHistcTauL;                   
206         
207         //////////////////////////Lambda 2D histos: cut vs mass////////////////
208         TH2F        *fHistDcaPosToPrimVertexLvsMassL;       
209         TH2F        *fHistDcaNegToPrimVertexLvsMassL;       
210         TH2F        *fHistRadiusV0LvsMassL;                  
211         TH2F        *fHistDecayLengthV0LvsMassL;             
212         TH2F        *fHistDcaV0DaughtersLvsMassL;          
213         TH2F        *fHistCosPointAngleLvsMassL;             
214         TH3F        *fHistCosPointAngleLVsMassVsPtsigL;     
215         TH3F        *fHistCosPointAngleLVsMassVsPtbackL;     
216         //////////////////////////Lambda 2D histos: cut vs on fly status////////////////////
217         
218         TH2F        *fHistDcaPosToPrimVertexAntiL;        
219         TH2F        *fHistDcaNegToPrimVertexAntiL;        
220         TH2F        *fHistRadiusV0AntiL;        
221         TH2F        *fHistDecayLengthV0AntiL;        
222         TH2F        *fHistDcaV0DaughtersAntiL;        
223         TH2F        *fHistChi2AntiL;        
224         TH2F        *fHistCosPointAngleAntiL;           
225         
226         //////////////////////////AntiLambda 2D histos: cut vs mass////////////////
227         TH2F        *fHistDcaPosToPrimVertexAntiLvsMass;       
228         TH2F        *fHistDcaNegToPrimVertexAntiLvsMass;       
229         TH2F        *fHistRadiusV0AntiLvsMass;                  
230         TH2F        *fHistDecayLengthV0AntiLvsMass;             
231         TH2F        *fHistDcaV0DaughtersAntiLvsMass;          
232         TH2F        *fHistCosPointAngleAntiLvsMass;             
233         ////////////////////////////////////////////////////////////////////// 
234         TH1F        *fHistMassK0;        
235         TH1F        *fHistMassLambda;        
236         TH1F        *fHistMassAntiLambda;        
237         TH2F        *fHistMassVsRadiusK0;        
238         TH2F        *fHistMassVsRadiusLambda;        
239         TH2F        *fHistMassVsRadiusAntiLambda;        
240         ////////////////////////////////////////////////////////////////////////////
241         TH2F        *fHistPtVsMassK0;        
242         TH2F        *fHistPtVsMassLambda;        
243         TH2F        *fHistPtVsMassAntiLambda;        
244         TH2F        *fHistArmenterosPodolanski;                       
245         
246         //PID
247         TH1F        *fHistNsigmaPosProtonLambda;         
248         TH1F        *fHistNsigmaNegPionLambda;            
249         TH1F        *fHistNsigmaPosProtonAntiLambda;         
250         TH1F        *fHistNsigmaNegPionAntiLambda;            
251         TH1F        *fHistNsigmaPosPionK0;                 
252         TH1F        *fHistNsigmaNegPionK0;                  
253         
254         // Associated particles histograms       
255         ////////////////////////////////////////////////////////////////////
256         TH1F        *fHistAsMcPtK0;        
257         TH1F        *fHistAsMcPtLambda;        
258         TH1F        *fHistAsMcPtAntiLambda;        
259         /////////////////////////////////////////////////////////////////////       
260         TH1F        *fHistAsMcProdRadiusK0;        
261         TH1F        *fHistAsMcProdRadiusLambda;        
262         TH1F        *fHistAsMcProdRadiusAntiLambda;        
263         TH2F        *fHistAsMcProdRadiusXvsYK0s;        
264         TH2F        *fHistAsMcProdRadiusXvsYLambda;        
265         TH2F        *fHistAsMcProdRadiusXvsYAntiLambda;        
266         TH1F        *fHistPidMcMassK0;        
267         TH1F        *fHistPidMcMassLambda;        
268         TH1F        *fHistPidMcMassAntiLambda;        
269         TH1F        *fHistAsMcMassK0;        
270         TH1F        *fHistAsMcMassLambda;        
271         TH1F        *fHistAsMcMassAntiLambda;        
272         TH2F        *fHistAsMcPtVsMassK0;        
273         TH2F        *fHistAsMcPtVsMassLambda;        
274         TH2F        *fHistAsMcPtVsMassAntiLambda;        
275         TH2F        *fHistAsMcMassVsRadiusK0;        
276         TH2F        *fHistAsMcMassVsRadiusLambda;        
277         TH2F        *fHistAsMcMassVsRadiusAntiLambda;        
278         TH1F        *fHistAsMcResxK0;        
279         TH1F        *fHistAsMcResyK0;        
280         TH1F        *fHistAsMcReszK0;        
281         TH2F        *fHistAsMcResrVsRadiusK0;        
282         TH2F        *fHistAsMcReszVsRadiusK0;        
283         TH1F        *fHistAsMcResxLambda;        
284         TH1F        *fHistAsMcResyLambda;        
285         TH1F        *fHistAsMcReszLambda;        
286         TH2F        *fHistAsMcResrVsRadiusLambda;        
287         TH2F        *fHistAsMcReszVsRadiusLambda;        
288         TH1F        *fHistAsMcResxAntiLambda;        
289         TH1F        *fHistAsMcResyAntiLambda;        
290         TH1F        *fHistAsMcReszAntiLambda;        
291         TH2F        *fHistAsMcResrVsRadiusAntiLambda;        
292         TH2F        *fHistAsMcReszVsRadiusAntiLambda;        
293         TH1F        *fHistAsMcResPtK0;        
294         TH1F        *fHistAsMcResPtLambda;        
295         TH1F        *fHistAsMcResPtAntiLambda;        
296         TH2F        *fHistAsMcResPtVsRapK0;        
297         TH2F        *fHistAsMcResPtVsRapLambda;        
298         TH2F        *fHistAsMcResPtVsRapAntiLambda;        
299         TH2F        *fHistAsMcResPtVsPtK0;        
300         TH2F        *fHistAsMcResPtVsPtLambda;        
301         TH2F        *fHistAsMcResPtVsPtAntiLambda;        
302         TH1F        *fHistAsMcPtLambdaFromSigma;        
303         TH1F        *fHistAsMcPtAntiLambdaFromSigma;        
304         // Associated secondary particles:
305         TH2F        *fHistAsMcSecondaryPtVsRapK0s;        
306         TH2F        *fHistAsMcSecondaryPtVsRapLambda;        
307         TH2F        *fHistAsMcSecondaryPtVsRapAntiLambda;        
308         TH1F        *fHistAsMcSecondaryProdRadiusK0s;        
309         TH1F        *fHistAsMcSecondaryProdRadiusLambda;        
310         TH1F        *fHistAsMcSecondaryProdRadiusAntiLambda;        
311         TH2F        *fHistAsMcSecondaryProdRadiusXvsYK0s;        
312         TH2F        *fHistAsMcSecondaryProdRadiusXvsYLambda;        
313         TH2F        *fHistAsMcSecondaryProdRadiusXvsYAntiLambda;        
314         TH1F        *fHistAsMcSecondaryPtLambdaFromSigma;        
315         TH1F        *fHistAsMcSecondaryPtAntiLambdaFromSigma;        
316         
317         TH2F        * fHistSibK0;
318         TH2F        * fHistMixK0;
319         TH2F        * fHistSibLambda;
320         TH2F        * fHistMixLambda;
321         
322         TH2F        * fHistSibK0MC;
323         TH2F        * fHistMixK0MC;
324         TH2F        * fHistSibLambdaMC;
325         TH2F        * fHistMixLambdaMC;
326         
327         ClassDef(AliLeadingV0Correlation, 1); 
328 };
329 //---------------------------------------------------------------------------------------
330 class AliLeadingBasicParticle : public AliVParticle
331 {
332   public:
333     AliLeadingBasicParticle(Float_t eta, Float_t phi, Float_t pt)
334       : fEta(eta), fPhi(phi), fpT(pt)
335     {
336     }
337     virtual ~AliLeadingBasicParticle() {}
338
339     // kinematics
340     virtual Double_t Px() const { AliFatal("Not implemented"); return 0; }
341     virtual Double_t Py() const { AliFatal("Not implemented"); return 0; }
342     virtual Double_t Pz() const { AliFatal("Not implemented"); return 0; }
343     virtual Double_t Pt() const { return fpT; }
344     virtual Double_t P() const { AliFatal("Not implemented"); return 0; }
345     virtual Bool_t   PxPyPz(Double_t[3]) const { AliFatal("Not implemented"); return 0; }
346
347     virtual Double_t Xv() const { AliFatal("Not implemented"); return 0; }
348     virtual Double_t Yv() const { AliFatal("Not implemented"); return 0; }
349     virtual Double_t Zv() const { AliFatal("Not implemented"); return 0; }
350     virtual Bool_t   XvYvZv(Double_t[3]) const { AliFatal("Not implemented"); return 0; }
351
352     virtual Double_t OneOverPt()  const { AliFatal("Not implemented"); return 0; }
353     virtual Double_t Phi()        const { return fPhi; }
354     virtual Double_t Theta()      const { AliFatal("Not implemented"); return 0; }
355
356
357     virtual Double_t E()          const { AliFatal("Not implemented"); return 0; }
358     virtual Double_t M()          const { AliFatal("Not implemented"); return 0; }
359
360     virtual Double_t Eta()        const { return fEta; }
361     virtual Double_t Y()          const { AliFatal("Not implemented"); return 0; }
362
363     virtual Short_t Charge()      const { AliFatal("Not implemented"); return 0; }
364     virtual Int_t   GetLabel()    const { AliFatal("Not implemented"); return 0; }
365     // PID
366     virtual Int_t   PdgCode()     const { AliFatal("Not implemented"); return 0; }
367     virtual const Double_t *PID() const { AliFatal("Not implemented"); return 0; }
368         
369
370   private:
371     Float_t fEta;      // eta
372     Float_t fPhi;      // phi
373     Float_t fpT;       // pT
374
375     ClassDef( AliLeadingBasicParticle, 1); // class required for event mixing
376 };
377
378 #endif
379