Changes to compile with Root6 on macosx64
[u/mrichter/AliRoot.git] / PWGCF / FEMTOSCOPY / PLamAnalysis / AliAnalysisTaskProtonLambda.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 #if !defined(__CINT__) || defined(__MAKECINT__)
16
17 #include <TH1F.h>
18 #include <TH2F.h>
19 #include <TH3F.h>
20
21 #include <TAxis.h>
22 // #include "TObjArray.h"
23
24 #include <AliAnalysisTask.h>
25 #include <AliAnalysisManager.h>
26
27 #include <AliAODEvent.h>
28 #include <AliAODVertex.h>
29 #include <AliAODv0.h>
30 #include <AliAODInputHandler.h>
31
32 #include "AliAnalysisTaskProtonLambda.h"
33 #include <AliCentrality.h>
34 //#include "AliAODpid.h"
35 #include <AliPID.h>
36 #include <AliPIDResponse.h>
37 // #include <../STEER/STEER/AliV0.h>
38 #include <AliExternalTrackParam.h>
39 //#include <AliAODTrack.h>
40 //#include <AliESDtrack.h>
41
42 //#include "EventCollection.h"
43
44 // Task to study femtoscopic proton-lambda correlations
45 // Author: Hans Beck
46
47 ClassImp(AliAnalysisTaskProtonLambda)
48 //ClassImp(AliAnalysisTaskProtonLambda::GlobalTrackInfo)
49 //ClassImp(AliAnalysisTaskProtonLambda::GTIContainer)
50
51 #endif  
52 //________________________________________________________________________
53 AliAnalysisTaskProtonLambda::AliAnalysisTaskProtonLambda() 
54   : AliAnalysisTaskSE(),
55   fkUseOnTheFly(kTRUE),
56   fkAbsZvertexCut(10.0),
57   fkCentCut(20.0),
58   fkLamMass(1.115683),
59   fkProMass(0.9382720),
60   fkPioMass(0.13957018),
61   fPIDResponse(0), 
62   fTpcResponse(0),
63   fFemtoBuffer(0),
64   fAOD(0), fPrimaryVtx(0), fOutputList(0), fOutputPrimaries(0),
65   fOutput2Part(0),
66   fGTI(0),fTrackBuffSize(18000),
67   fHistGoodEvent(0),
68   // fHistPrimaryVertexPosXY(0), fHistPrimaryVertexPosZ(0),        
69   // fHistTrackMultiplicity(0),    
70   // fHistShareV0pos(0),fHistShareV0neg(0),
71   // fHistPosTofBeforeCut(0), fHistPosTofAfterCut(0),           
72   // fHistNegTofBeforeCut(0), fHistNegTofAfterCut(0),           
73   // fHistPosTpcBeforeCut(0), fHistPosTpcAfterCut(0),            
74   // fHistNegTpcBeforeCut(0), fHistNegTpcAfterCut(0),            
75   // fHistGoodV0(0), fHistCorrectSigns(0),              
76   // fHistDcaPosToPrimVertex(0), fHistDcaNegToPrimVertex(0),        
77   // fHistDcaPosToPrimVertexZoom(0), fHistDcaNegToPrimVertexZoom(0),  
78   // fHistRadiusV0(0), fHistDecayLengthV0(0), fHistDcaV0Daughters(0),          
79   // fHistChi2(0), fHistCosPointAngle(0), fHistCosPointAngleZoom(0),
80   fHistSideBandOffLam(0), fHistSideBandOffALam(0), fHistTPCNclsPosOffLam(0),       
81   fHistTPCNclsNegOffLam(0), fHistTPCNclsPosOffALam(0), fHistTPCNclsNegOffALam(0),      
82   // fHistPosNsigmaTpcOffLam(0), fHistPosNsigmaTpcOffALam(0), fHistNegNsigmaTpcOffLam(0),
83   // fHistNegNsigmaTpcOffALam(0), fHistUseTofOffLam(0), fHistUseTofOffALam(0),
84   // fHistDcaPosOffLam(0), fHistDcaPosOffALam(0), fHistDcaNegOffLam(0),           
85   // fHistDcaNegOffALam(0), fHistDcaV0DaughtersOffLam(0), fHistDcaV0DaughtersOffALam(0),  
86   // fHistCosPointLamOff(0), fHistCosPointALamOff(0), fHistCosPointLamZoomOff(0),     
87   // fHistCosPointALamZoomOff(0), fHistV0RadiusLamOff(0), fHistV0RadiusALamOff(0),        
88   // fHistV0DecayLengthLamOff(0), fHistV0DecayLengthALamOff(0), fHistDcaV0PriVertexLamOff(0),     
89   // fHistDcaV0PriVertexALamOff(0),
90   fHistMassLambdaOff(0), fHistMassAntiLambdaOff(0),        
91   // fHistPtVsMassLambdaOff(0), fHistPtVsMassAntiLambdaOff(0),
92   fHistYPtMassLamOff(0), fHistYPtMassALamOff(0),
93   // fHistPtVsYLambdaOff(0), fHistPtVsYAntiLambdaOff(0),       
94   fHistSideBandOnLam(0), fHistSideBandOnALam(0),
95   // fHistLikeSignOnLam(0), fHistLikeSignOnALam(0),         
96   fHistTPCNclsPosOnLam(0), fHistTPCNclsNegOnLam(0), fHistTPCNclsPosOnALam(0),fHistTPCNclsNegOnALam(0),     
97   // fHistPosNsigmaTpcOnLam(0), fHistPosNsigmaTpcOnALam(0), fHistNegNsigmaTpcOnLam(0), fHistNegNsigmaTpcOnALam(0),        
98   // fHistUseTofOnLam(0),fHistUseTofOnALam(0),fHistDcaPosOnLam(0),fHistDcaPosOnALam(0),fHistDcaNegOnLam(0),               
99   // fHistDcaNegOnALam(0),fHistDcaV0DaughtersOnLam(0),fHistDcaV0DaughtersOnALam(0),fHistCosPointLamOn(0),             
100   // fHistCosPointALamOn(0),fHistCosPointLamZoomOn(0),fHistCosPointALamZoomOn(0),fHistV0RadiusLamOn(0),             
101   // fHistV0RadiusALamOn(0),fHistV0DecayLengthLamOn(0),fHistV0DecayLengthALamOn(0),fHistDcaV0PriVertexLamOn(0),       
102   // fHistDcaV0PriVertexALamOn(0),
103   // fHistChi2TPCPosLamOn(0),  fHistChi2TPCPosALamOn(0),  fHistChi2TPCNegLamOn(0),  fHistChi2TPCNegALamOn(0),
104   // fHistMinvTPConlyLamOn(0),  fHistMinvTPConlyALamOn(0),
105   fHistMassLambdaOn(0),fHistMassAntiLambdaOn(0),
106   // fHistPtVsMassLambdaOn(0), fHistPtVsMassAntiLambdaOn(0),
107   fHistYPtMassLamOn(0),fHistYPtMassALamOn(0),
108   // fHistPtVsYLambdaOn(0), fHistPtVsYAntiLambdaOn(0),
109   // fHistMomDiffLam(0),fHistMomDiffALam(0),fHistMomDiffBgLam(0),fHistMomDiffBgALam(0),
110   // fHistMomDiffWoSPDLam(0),fHistMomDiffWoSPDALam(0),fHistMomDiffWoSPDBgLam(0),fHistMomDiffWoSPDBgALam(0),
111   fPriHistShare(0),
112   // fPriHistPosNsigmaTof(0),
113   fPriHistPosNsigmaTofVsP(0),fPriHistPosNsigmaTofVsPt(0),     
114   // fPriHistNegNsigmaTof(0),
115   fPriHistNegNsigmaTofVsP(0),fPriHistNegNsigmaTofVsPt(0),fPriHistTOFsignalPosVsP(0),      
116   fPriHistTOFsignalPosVsPt(0),fPriHistTOFsignalNegVsP(0),fPriHistTOFsignalNegVsPt(0),fPriHistHybridTOFsigPosWoTPC(0), 
117   fPriHistHybridTOFsigPosTPCok(0),fPriHistHybridTOFsigNegWoTPC(0),fPriHistHybridTOFsigNegTPCok(0),
118   // fPriHistHasTofPos(0),          
119   fPriHistTPCsignalPos(0),
120   // fPriHistNsigmaTPCPos(0), fPriHistTPCsignalTOFcutPos(0),fPriHistNsigmaTPCTOFcutPos(0),   
121   fPriHistTPCsignalLowPPos(0),fPriHistTPCsignalMedPPos(0),fPriHistTPCsignalHigPPos(0),
122   // fPriHistHasTofNeg(0),         
123   fPriHistTPCsignalNeg(0),
124   // fPriHistNsigmaTPCNeg(0),fPriHistTPCsignalTOFcutNeg(0),fPriHistNsigmaTPCTOFcutNeg(0),    
125   fPriHistTPCsignalLowPNeg(0),fPriHistTPCsignalMedPNeg(0),fPriHistTPCsignalHigPNeg(0),
126   fPriHistDCAxyYPtPro(0),fPriHistDCAxyYPtAPro(0),
127   // f2HistLamLamMeanMinDistProReal(0),
128   // f2HistLamLamMeanMinDistPioReal(0),f2HistLamProMeanMinDistProReal(0),f2HistALamALamMeanMinDistAProReal(0), 
129   // f2HistALamALamMeanMinDistPioReal(0),f2HistALamAProMeanMinDistAProReal(0),
130   // f2HistSftLamLamMeanMinDistProReal(0),
131   // f2HistSftLamLamMeanMinDistPioReal(0),f2HistSftLamProMeanMinDistProReal(0),f2HistSftALamALamMeanMinDistAProReal(0), 
132   // f2HistSftALamALamMeanMinDistPioReal(0),f2HistSftALamAProMeanMinDistAProReal(0),
133   // f2HistSftIrocLamLamMeanMinDistProReal(0),
134   // f2HistSftIrocLamLamMeanMinDistPioReal(0),f2HistSftIrocLamProMeanMinDistProReal(0),f2HistSftIrocALamALamMeanMinDistAProReal(0), 
135   // f2HistSftIrocALamALamMeanMinDistPioReal(0),f2HistSftIrocALamAProMeanMinDistAProReal(0),
136   // f2HistSftOrocLamLamMeanMinDistProReal(0),
137   // f2HistSftOrocLamLamMeanMinDistPioReal(0),f2HistSftOrocLamProMeanMinDistProReal(0),f2HistSftOrocALamALamMeanMinDistAProReal(0), 
138   // f2HistSftOrocALamALamMeanMinDistPioReal(0),f2HistSftOrocALamAProMeanMinDistAProReal(0),
139   // f2HistMtLamLamReal(0), 
140   f2HistMtLamProReal(0), 
141   // f2HistMtALamALamReal(0), 
142   f2HistMtALamAProReal(0),
143   // f2HistMtLowQLamLamReal(0), 
144   f2HistMtLowQLamProReal(0), 
145   // f2HistMtLowQALamALamReal(0), 
146   f2HistMtLowQALamAProReal(0),
147   LamProReal(0),ALamAProReal(0),
148   // f3HistLamLamQinvReal(0),               
149   // f3HistALamALamQinvReal(0),f3HistLamLamMinvReal(0),               
150   // f3HistLamProMinvReal(0),f3HistALamALamMinvReal(0),f3HistALamAProMinvReal(0),
151   // f2HistBgLamBgLamMeanMinDistProReal(0),f2HistBgLamBgLamMeanMinDistPioReal(0),
152   // f2HistBgLamProMeanMinDistProReal(0),f2HistBgALamBgALamMeanMinDistAProReal(0),
153   // f2HistBgALamBgALamMeanMinDistPioReal(0),f2HistBgALamAProMeanMinDistAProReal(0),
154   // f2HistSftBgLamBgLamMeanMinDistProReal(0),f2HistSftBgLamBgLamMeanMinDistPioReal(0),
155   // f2HistSftBgLamProMeanMinDistProReal(0),f2HistSftBgALamBgALamMeanMinDistAProReal(0),
156   // f2HistSftBgALamBgALamMeanMinDistPioReal(0),f2HistSftBgALamAProMeanMinDistAProReal(0),
157   // f2HistSftIrocBgLamBgLamMeanMinDistProReal(0),f2HistSftIrocBgLamBgLamMeanMinDistPioReal(0),
158   // f2HistSftIrocBgLamProMeanMinDistProReal(0),f2HistSftIrocBgALamBgALamMeanMinDistAProReal(0),
159   // f2HistSftIrocBgALamBgALamMeanMinDistPioReal(0),f2HistSftIrocBgALamAProMeanMinDistAProReal(0),
160   // f2HistSftOrocBgLamBgLamMeanMinDistProReal(0),f2HistSftOrocBgLamBgLamMeanMinDistPioReal(0),
161   // f2HistSftOrocBgLamProMeanMinDistProReal(0),f2HistSftOrocBgALamBgALamMeanMinDistAProReal(0),
162   // f2HistSftOrocBgALamBgALamMeanMinDistPioReal(0),f2HistSftOrocBgALamAProMeanMinDistAProReal(0),
163   BgLamProReal(0),BgALamAProReal(0),
164   // f3HistBgLamBgLamQinvReal(0),             
165   // f3HistBgALamBgALamQinvReal(0),
166   // f2HistLamLamMeanMinDistProMixed(0),f2HistLamLamMeanMinDistPioMixed(0),
167   // f2HistLamProMeanMinDistProMixed(0),f2HistALamALamMeanMinDistAProMixed(0),   
168   // f2HistALamALamMeanMinDistPioMixed(0),f2HistALamAProMeanMinDistAProMixed(0),
169   // f2HistSftLamLamMeanMinDistProMixed(0),f2HistSftLamLamMeanMinDistPioMixed(0),
170   // f2HistSftLamProMeanMinDistProMixed(0),f2HistSftALamALamMeanMinDistAProMixed(0),   
171   // f2HistSftALamALamMeanMinDistPioMixed(0),f2HistSftALamAProMeanMinDistAProMixed(0),
172   // f2HistSftIrocLamLamMeanMinDistProMixed(0),f2HistSftIrocLamLamMeanMinDistPioMixed(0),
173   // f2HistSftIrocLamProMeanMinDistProMixed(0),f2HistSftIrocALamALamMeanMinDistAProMixed(0),   
174   // f2HistSftIrocALamALamMeanMinDistPioMixed(0),f2HistSftIrocALamAProMeanMinDistAProMixed(0),
175   // f2HistSftOrocLamLamMeanMinDistProMixed(0),f2HistSftOrocLamLamMeanMinDistPioMixed(0),
176   // f2HistSftOrocLamProMeanMinDistProMixed(0),f2HistSftOrocALamALamMeanMinDistAProMixed(0),   
177   // f2HistSftOrocALamALamMeanMinDistPioMixed(0),f2HistSftOrocALamAProMeanMinDistAProMixed(0),
178   LamProMixed(0),ALamAProMixed(0),
179   // f3HistLamLamQinvMixed(0),                
180   // f3HistALamALamQinvMixed(0),f3HistLamLamMinvMixed(0),                
181   // f3HistLamProMinvMixed(0),f3HistALamALamMinvMixed(0),f3HistALamAProMinvMixed(0),
182   // f2HistBgLamBgLamMeanMinDistProMixed(0),f2HistBgLamBgLamMeanMinDistPioMixed(0),
183   // f2HistBgLamProMeanMinDistProMixed(0),f2HistBgALamBgALamMeanMinDistAProMixed(0),
184   // f2HistBgALamBgALamMeanMinDistPioMixed(0),f2HistBgALamAProMeanMinDistAProMixed(0),
185   // f2HistSftBgLamBgLamMeanMinDistProMixed(0),f2HistSftBgLamBgLamMeanMinDistPioMixed(0),
186   // f2HistSftBgLamProMeanMinDistProMixed(0),f2HistSftBgALamBgALamMeanMinDistAProMixed(0),
187   // f2HistSftBgALamBgALamMeanMinDistPioMixed(0),f2HistSftBgALamAProMeanMinDistAProMixed(0),
188   // f2HistSftIrocBgLamBgLamMeanMinDistProMixed(0),f2HistSftIrocBgLamBgLamMeanMinDistPioMixed(0),
189   // f2HistSftIrocBgLamProMeanMinDistProMixed(0),f2HistSftIrocBgALamBgALamMeanMinDistAProMixed(0),
190   // f2HistSftIrocBgALamBgALamMeanMinDistPioMixed(0),f2HistSftIrocBgALamAProMeanMinDistAProMixed(0),
191   // f2HistSftOrocBgLamBgLamMeanMinDistProMixed(0),f2HistSftOrocBgLamBgLamMeanMinDistPioMixed(0),
192   // f2HistSftOrocBgLamProMeanMinDistProMixed(0),f2HistSftOrocBgALamBgALamMeanMinDistAProMixed(0),
193   // f2HistSftOrocBgALamBgALamMeanMinDistPioMixed(0),f2HistSftOrocBgALamAProMeanMinDistAProMixed(0),
194   BgLamProMixed(0),BgALamAProMixed(0)
195   // f3HistBgLamBgLamQinvMixed(0),             
196   // f3HistBgALamBgALamQinvMixed(0)
197 {
198   // Dummy constructor
199   fPrimaryVtxPosition[0]=0;
200   fPrimaryVtxPosition[1]=0;
201   fPrimaryVtxPosition[2]=0;
202 }
203 //________________________________________________________________________
204 AliAnalysisTaskProtonLambda::AliAnalysisTaskProtonLambda(const char *name) 
205   : AliAnalysisTaskSE(name),
206   fkUseOnTheFly(kTRUE),
207   fkAbsZvertexCut(10.0),
208   fkCentCut(20.0),
209   fkLamMass(1.115683),
210   fkProMass(0.9382720),
211   fkPioMass(0.13957018),
212   fPIDResponse(0), 
213   fTpcResponse(0),
214   fFemtoBuffer(0),
215   fAOD(0), fPrimaryVtx(0), fOutputList(0), fOutputPrimaries(0),
216   fOutput2Part(0),
217   fGTI(0),fTrackBuffSize(18000),
218   fHistGoodEvent(0),
219   // fHistPrimaryVertexPosXY(0), fHistPrimaryVertexPosZ(0),        
220   // fHistTrackMultiplicity(0),    
221   // fHistShareV0pos(0),fHistShareV0neg(0),
222   // fHistPosTofBeforeCut(0), fHistPosTofAfterCut(0),           
223   // fHistNegTofBeforeCut(0), fHistNegTofAfterCut(0),           
224   // fHistPosTpcBeforeCut(0), fHistPosTpcAfterCut(0),            
225   // fHistNegTpcBeforeCut(0), fHistNegTpcAfterCut(0),            
226   // fHistGoodV0(0), fHistCorrectSigns(0),              
227   // fHistDcaPosToPrimVertex(0), fHistDcaNegToPrimVertex(0),        
228   // fHistDcaPosToPrimVertexZoom(0), fHistDcaNegToPrimVertexZoom(0),  
229   // fHistRadiusV0(0), fHistDecayLengthV0(0), fHistDcaV0Daughters(0),          
230   // fHistChi2(0), fHistCosPointAngle(0), fHistCosPointAngleZoom(0),
231   fHistSideBandOffLam(0), fHistSideBandOffALam(0), fHistTPCNclsPosOffLam(0),       
232   fHistTPCNclsNegOffLam(0), fHistTPCNclsPosOffALam(0), fHistTPCNclsNegOffALam(0),      
233   // fHistPosNsigmaTpcOffLam(0), fHistPosNsigmaTpcOffALam(0), fHistNegNsigmaTpcOffLam(0),
234   // fHistNegNsigmaTpcOffALam(0), fHistUseTofOffLam(0), fHistUseTofOffALam(0),
235   // fHistDcaPosOffLam(0), fHistDcaPosOffALam(0), fHistDcaNegOffLam(0),           
236   // fHistDcaNegOffALam(0), fHistDcaV0DaughtersOffLam(0), fHistDcaV0DaughtersOffALam(0),  
237   // fHistCosPointLamOff(0), fHistCosPointALamOff(0), fHistCosPointLamZoomOff(0),     
238   // fHistCosPointALamZoomOff(0), fHistV0RadiusLamOff(0), fHistV0RadiusALamOff(0),        
239   // fHistV0DecayLengthLamOff(0), fHistV0DecayLengthALamOff(0), fHistDcaV0PriVertexLamOff(0),     
240   // fHistDcaV0PriVertexALamOff(0),
241   fHistMassLambdaOff(0), fHistMassAntiLambdaOff(0),        
242   // fHistPtVsMassLambdaOff(0), fHistPtVsMassAntiLambdaOff(0),
243   fHistYPtMassLamOff(0), fHistYPtMassALamOff(0),
244   // fHistPtVsYLambdaOff(0), fHistPtVsYAntiLambdaOff(0),       
245   fHistSideBandOnLam(0), fHistSideBandOnALam(0),
246   // fHistLikeSignOnLam(0), fHistLikeSignOnALam(0),         
247   fHistTPCNclsPosOnLam(0), fHistTPCNclsNegOnLam(0), fHistTPCNclsPosOnALam(0),fHistTPCNclsNegOnALam(0),     
248   // fHistPosNsigmaTpcOnLam(0), fHistPosNsigmaTpcOnALam(0), fHistNegNsigmaTpcOnLam(0), fHistNegNsigmaTpcOnALam(0),        
249   // fHistUseTofOnLam(0),fHistUseTofOnALam(0),fHistDcaPosOnLam(0),fHistDcaPosOnALam(0),fHistDcaNegOnLam(0),               
250   // fHistDcaNegOnALam(0),fHistDcaV0DaughtersOnLam(0),fHistDcaV0DaughtersOnALam(0),fHistCosPointLamOn(0),             
251   // fHistCosPointALamOn(0),fHistCosPointLamZoomOn(0),fHistCosPointALamZoomOn(0),fHistV0RadiusLamOn(0),             
252   // fHistV0RadiusALamOn(0),fHistV0DecayLengthLamOn(0),fHistV0DecayLengthALamOn(0),fHistDcaV0PriVertexLamOn(0),       
253   // fHistDcaV0PriVertexALamOn(0),
254   // fHistChi2TPCPosLamOn(0),  fHistChi2TPCPosALamOn(0),  fHistChi2TPCNegLamOn(0),  fHistChi2TPCNegALamOn(0),
255   // fHistMinvTPConlyLamOn(0),  fHistMinvTPConlyALamOn(0),
256   fHistMassLambdaOn(0),fHistMassAntiLambdaOn(0),
257   // fHistPtVsMassLambdaOn(0), fHistPtVsMassAntiLambdaOn(0),
258   fHistYPtMassLamOn(0),fHistYPtMassALamOn(0),
259   // fHistPtVsYLambdaOn(0), fHistPtVsYAntiLambdaOn(0),
260   // fHistMomDiffLam(0),fHistMomDiffALam(0),fHistMomDiffBgLam(0),fHistMomDiffBgALam(0),
261   // fHistMomDiffWoSPDLam(0),fHistMomDiffWoSPDALam(0),fHistMomDiffWoSPDBgLam(0),fHistMomDiffWoSPDBgALam(0),
262   fPriHistShare(0),
263   // fPriHistPosNsigmaTof(0),
264   fPriHistPosNsigmaTofVsP(0),fPriHistPosNsigmaTofVsPt(0),     
265   // fPriHistNegNsigmaTof(0),
266   fPriHistNegNsigmaTofVsP(0),fPriHistNegNsigmaTofVsPt(0),fPriHistTOFsignalPosVsP(0),      
267   fPriHistTOFsignalPosVsPt(0),fPriHistTOFsignalNegVsP(0),fPriHistTOFsignalNegVsPt(0),fPriHistHybridTOFsigPosWoTPC(0), 
268   fPriHistHybridTOFsigPosTPCok(0),fPriHistHybridTOFsigNegWoTPC(0),fPriHistHybridTOFsigNegTPCok(0),
269   // fPriHistHasTofPos(0),          
270   fPriHistTPCsignalPos(0),
271   // fPriHistNsigmaTPCPos(0), fPriHistTPCsignalTOFcutPos(0),fPriHistNsigmaTPCTOFcutPos(0),   
272   fPriHistTPCsignalLowPPos(0),fPriHistTPCsignalMedPPos(0),fPriHistTPCsignalHigPPos(0),
273   // fPriHistHasTofNeg(0),         
274   fPriHistTPCsignalNeg(0),
275   // fPriHistNsigmaTPCNeg(0),fPriHistTPCsignalTOFcutNeg(0),fPriHistNsigmaTPCTOFcutNeg(0),    
276   fPriHistTPCsignalLowPNeg(0),fPriHistTPCsignalMedPNeg(0),fPriHistTPCsignalHigPNeg(0),
277   fPriHistDCAxyYPtPro(0),fPriHistDCAxyYPtAPro(0),
278   // f2HistLamLamMeanMinDistProReal(0),
279   // f2HistLamLamMeanMinDistPioReal(0),f2HistLamProMeanMinDistProReal(0),f2HistALamALamMeanMinDistAProReal(0), 
280   // f2HistALamALamMeanMinDistPioReal(0),f2HistALamAProMeanMinDistAProReal(0),
281   // f2HistSftLamLamMeanMinDistProReal(0),
282   // f2HistSftLamLamMeanMinDistPioReal(0),f2HistSftLamProMeanMinDistProReal(0),f2HistSftALamALamMeanMinDistAProReal(0), 
283   // f2HistSftALamALamMeanMinDistPioReal(0),f2HistSftALamAProMeanMinDistAProReal(0),
284   // f2HistSftIrocLamLamMeanMinDistProReal(0),
285   // f2HistSftIrocLamLamMeanMinDistPioReal(0),f2HistSftIrocLamProMeanMinDistProReal(0),f2HistSftIrocALamALamMeanMinDistAProReal(0), 
286   // f2HistSftIrocALamALamMeanMinDistPioReal(0),f2HistSftIrocALamAProMeanMinDistAProReal(0),
287   // f2HistSftOrocLamLamMeanMinDistProReal(0),
288   // f2HistSftOrocLamLamMeanMinDistPioReal(0),f2HistSftOrocLamProMeanMinDistProReal(0),f2HistSftOrocALamALamMeanMinDistAProReal(0), 
289   // f2HistSftOrocALamALamMeanMinDistPioReal(0),f2HistSftOrocALamAProMeanMinDistAProReal(0),
290   // f2HistMtLamLamReal(0), 
291   f2HistMtLamProReal(0), 
292   // f2HistMtALamALamReal(0), 
293   f2HistMtALamAProReal(0),
294   // f2HistMtLowQLamLamReal(0), 
295   f2HistMtLowQLamProReal(0), 
296   // f2HistMtLowQALamALamReal(0), 
297   f2HistMtLowQALamAProReal(0),
298   LamProReal(0),ALamAProReal(0),
299   // f3HistLamLamQinvReal(0),               
300   // f3HistALamALamQinvReal(0),f3HistLamLamMinvReal(0),               
301   // f3HistLamProMinvReal(0),f3HistALamALamMinvReal(0),f3HistALamAProMinvReal(0),
302   // f2HistBgLamBgLamMeanMinDistProReal(0),f2HistBgLamBgLamMeanMinDistPioReal(0),
303   // f2HistBgLamProMeanMinDistProReal(0),f2HistBgALamBgALamMeanMinDistAProReal(0),
304   // f2HistBgALamBgALamMeanMinDistPioReal(0),f2HistBgALamAProMeanMinDistAProReal(0),
305   // f2HistSftBgLamBgLamMeanMinDistProReal(0),f2HistSftBgLamBgLamMeanMinDistPioReal(0),
306   // f2HistSftBgLamProMeanMinDistProReal(0),f2HistSftBgALamBgALamMeanMinDistAProReal(0),
307   // f2HistSftBgALamBgALamMeanMinDistPioReal(0),f2HistSftBgALamAProMeanMinDistAProReal(0),
308   // f2HistSftIrocBgLamBgLamMeanMinDistProReal(0),f2HistSftIrocBgLamBgLamMeanMinDistPioReal(0),
309   // f2HistSftIrocBgLamProMeanMinDistProReal(0),f2HistSftIrocBgALamBgALamMeanMinDistAProReal(0),
310   // f2HistSftIrocBgALamBgALamMeanMinDistPioReal(0),f2HistSftIrocBgALamAProMeanMinDistAProReal(0),
311   // f2HistSftOrocBgLamBgLamMeanMinDistProReal(0),f2HistSftOrocBgLamBgLamMeanMinDistPioReal(0),
312   // f2HistSftOrocBgLamProMeanMinDistProReal(0),f2HistSftOrocBgALamBgALamMeanMinDistAProReal(0),
313   // f2HistSftOrocBgALamBgALamMeanMinDistPioReal(0),f2HistSftOrocBgALamAProMeanMinDistAProReal(0),
314   BgLamProReal(0),BgALamAProReal(0),
315   // f3HistBgLamBgLamQinvReal(0),             
316   // f3HistBgALamBgALamQinvReal(0),
317   // f2HistLamLamMeanMinDistProMixed(0),f2HistLamLamMeanMinDistPioMixed(0),
318   // f2HistLamProMeanMinDistProMixed(0),f2HistALamALamMeanMinDistAProMixed(0),   
319   // f2HistALamALamMeanMinDistPioMixed(0),f2HistALamAProMeanMinDistAProMixed(0),
320   // f2HistSftLamLamMeanMinDistProMixed(0),f2HistSftLamLamMeanMinDistPioMixed(0),
321   // f2HistSftLamProMeanMinDistProMixed(0),f2HistSftALamALamMeanMinDistAProMixed(0),   
322   // f2HistSftALamALamMeanMinDistPioMixed(0),f2HistSftALamAProMeanMinDistAProMixed(0),
323   // f2HistSftIrocLamLamMeanMinDistProMixed(0),f2HistSftIrocLamLamMeanMinDistPioMixed(0),
324   // f2HistSftIrocLamProMeanMinDistProMixed(0),f2HistSftIrocALamALamMeanMinDistAProMixed(0),   
325   // f2HistSftIrocALamALamMeanMinDistPioMixed(0),f2HistSftIrocALamAProMeanMinDistAProMixed(0),
326   // f2HistSftOrocLamLamMeanMinDistProMixed(0),f2HistSftOrocLamLamMeanMinDistPioMixed(0),
327   // f2HistSftOrocLamProMeanMinDistProMixed(0),f2HistSftOrocALamALamMeanMinDistAProMixed(0),   
328   // f2HistSftOrocALamALamMeanMinDistPioMixed(0),f2HistSftOrocALamAProMeanMinDistAProMixed(0),
329   LamProMixed(0),ALamAProMixed(0),
330   // f3HistLamLamQinvMixed(0),                
331   // f3HistALamALamQinvMixed(0),f3HistLamLamMinvMixed(0),                
332   // f3HistLamProMinvMixed(0),f3HistALamALamMinvMixed(0),f3HistALamAProMinvMixed(0),
333   // f2HistBgLamBgLamMeanMinDistProMixed(0),f2HistBgLamBgLamMeanMinDistPioMixed(0),
334   // f2HistBgLamProMeanMinDistProMixed(0),f2HistBgALamBgALamMeanMinDistAProMixed(0),
335   // f2HistBgALamBgALamMeanMinDistPioMixed(0),f2HistBgALamAProMeanMinDistAProMixed(0),
336   // f2HistSftBgLamBgLamMeanMinDistProMixed(0),f2HistSftBgLamBgLamMeanMinDistPioMixed(0),
337   // f2HistSftBgLamProMeanMinDistProMixed(0),f2HistSftBgALamBgALamMeanMinDistAProMixed(0),
338   // f2HistSftBgALamBgALamMeanMinDistPioMixed(0),f2HistSftBgALamAProMeanMinDistAProMixed(0),
339   // f2HistSftIrocBgLamBgLamMeanMinDistProMixed(0),f2HistSftIrocBgLamBgLamMeanMinDistPioMixed(0),
340   // f2HistSftIrocBgLamProMeanMinDistProMixed(0),f2HistSftIrocBgALamBgALamMeanMinDistAProMixed(0),
341   // f2HistSftIrocBgALamBgALamMeanMinDistPioMixed(0),f2HistSftIrocBgALamAProMeanMinDistAProMixed(0),
342   // f2HistSftOrocBgLamBgLamMeanMinDistProMixed(0),f2HistSftOrocBgLamBgLamMeanMinDistPioMixed(0),
343   // f2HistSftOrocBgLamProMeanMinDistProMixed(0),f2HistSftOrocBgALamBgALamMeanMinDistAProMixed(0),
344   // f2HistSftOrocBgALamBgALamMeanMinDistPioMixed(0),f2HistSftOrocBgALamAProMeanMinDistAProMixed(0),
345   BgLamProMixed(0),BgALamAProMixed(0)
346   // f3HistBgLamBgLamQinvMixed(0),             
347   // f3HistBgALamBgALamQinvMixed(0)
348 {
349   // Constructor
350   fPrimaryVtxPosition[0]=0;
351   fPrimaryVtxPosition[1]=0;
352   fPrimaryVtxPosition[2]=0;
353
354   // Define output slots only here
355   // Output slot #1 writes into a TList container
356   DefineOutput(1, TList::Class());
357   DefineOutput(2, TList::Class());
358   DefineOutput(3, TList::Class());
359 }
360 //________________________________________________________________________
361 AliAnalysisTaskProtonLambda::~AliAnalysisTaskProtonLambda() {
362   // Destructor, go through the data member and delete them
363
364   // fPIDResponse is just a pointer to the pid response task,
365   // we don't create it so we don't delete it. It comes from 
366   // the AliInputEventHandler
367
368   if (fTpcResponse){
369     delete fTpcResponse;
370     fTpcResponse=0;
371   }
372   if(fFemtoBuffer){
373     delete fFemtoBuffer;
374     fFemtoBuffer=0;
375   }
376   // fAOD also just comes from a function of the AliAnalysisTaskSE
377   // fPrimaryVtx comes from the fAOD
378
379   // The lists containing the histograms
380   if (fOutputList){
381     fOutputList->Delete();
382     delete fOutputList;
383     fOutputList=0;
384   }
385   if (fOutputPrimaries){
386     fOutputPrimaries->Delete();
387     delete fOutputPrimaries;
388     fOutputPrimaries=0;
389   }
390   if (fOutput2Part){
391     fOutput2Part->Delete();
392     delete fOutput2Part;
393     fOutput2Part=0;
394   }
395
396   // Array, note the [] with the delete
397   if (fGTI)
398     delete[] fGTI;
399   fGTI=0;
400
401 }
402 //________________________________________________________________________
403 void AliAnalysisTaskProtonLambda::UserCreateOutputObjects()
404 {
405   // Create histograms and other objects and variables
406   // Called once
407
408   // Get the PID response object
409   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
410   if(!man){AliError("Couldn't get the analysis manager!");}
411   AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
412   if(!inputHandler){AliError("Couldn't get the input handler!");}
413   fPIDResponse = inputHandler->GetPIDResponse();
414   if(!fPIDResponse){AliError("Couldn't get the PID response task!");}
415
416   // Create dE/dx spectra cut. use it by calling
417   // fTpcResponse->GetExpectedSignal(mom, AliPID::kProton)
418   fTpcResponse = new AliTPCPIDResponse();
419   Double_t AlephParameters[5];
420   // They are only valid for data, see $ALICE_ROOT/PWG2/SPECTRA/AliProtonAnalysisBase.cxx
421   // for monte carlo parameters
422   AlephParameters[0] = 0.0283086;
423   AlephParameters[1] = 2.63394e+01;
424   AlephParameters[2] = 5.04114e-11;
425   AlephParameters[3] = 2.12543e+00;
426   AlephParameters[4] = 4.88663e+00;
427   fTpcResponse->SetBetheBlochParameters(AlephParameters[0],AlephParameters[1],AlephParameters[2],AlephParameters[3],AlephParameters[4]);
428   
429   // Create the buffer for event mixing
430   // Standard values are
431   //  fkZvertexBins(10),
432   //  fkCentBins(10),
433   //  fkMixBuff(5),
434   //  fkPriTrackLim(100),
435   //  fkV0Lim(50),
436   //  fFemtoBuffer = new FemtoBuffer(10,10,5,100,50,fkAbsZvertexCut,fkCentCut);
437   fFemtoBuffer = new FemtoBuffer(4,2,7,100,50,fkAbsZvertexCut,fkCentCut);
438
439   // In AODs, TPC only tracks don't have the pid information stored.
440   // Also, the TPC only tracks don't have any resolution in the DCAxy
441   // to distinguish between primaries and secondaries so we need the
442   // corresponding global track for every TPC only track. The way to do 
443   // this is to just store the pointer to the global track for every id.
444   fGTI = new AliAODTrack *[fTrackBuffSize]; // Array of pointers 
445
446   // Create the output list
447   fOutputList = new TList();
448   fOutputList->SetOwner();
449   fOutputPrimaries = new TList();
450   fOutputPrimaries->SetOwner();
451   fOutput2Part = new TList();
452   fOutput2Part->SetOwner();
453
454   // Invariant mass binning for lambdas
455   const Int_t nMinvBins = 140;
456   const Float_t minvLowEdge=1.060683, minvHiEdge=1.200683;
457
458   // Control hist for event cuts
459   fHistGoodEvent = new TH1F("h1GoodEvent","No of events passing the cuts.",10,-.5,9.5);
460   fOutputList->Add(fHistGoodEvent);
461
462   // Primary Vertex:
463   // fHistPrimaryVertexPosXY       = new TH2F("h2PrimaryVertexPosXY", "Primary Vertex Position XY;Primary Vertex Position X (cm);Primary Vertex Position Y (cm)",100,-0.5,0.5,100,-0.5,0.5);
464   // fOutputList->Add(fHistPrimaryVertexPosXY);
465   // fHistPrimaryVertexPosZ       = new TH1F("h1PrimaryVertexPosZ", "Primary Vertex Position Z;Primary Vertex Position Z (cm);Events",200,-12.0,12.0);
466   // fOutputList->Add(fHistPrimaryVertexPosZ);
467
468   // // Multiplicity
469   // fHistTrackMultiplicity = new TH1F("fHistTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 200, 0, 40000);
470   // fOutputList->Add(fHistTrackMultiplicity);
471
472   // //
473   // //  V0 histograms
474   // //
475   // // Shared clusters
476   // fHistShareV0pos = new TH1F("h1ShareV0pos","Shared clusters pos V0 daughters;#shared clusters;counts"
477   //                         ,160,0,160);
478   // fOutputList->Add(fHistShareV0pos);
479   // fHistShareV0neg = new TH1F("h1ShareV0neg","Shared clusters neg V0 daughters;#shared clusters;counts"
480   //                         ,160,0,160);
481   // fOutputList->Add(fHistShareV0neg);
482   
483   // // PID via TPC dE/dx
484   // fHistPosTpcBeforeCut = new TH2F ("h2PosTpcBeforeCut","TPC signal (pos daughter) before cut;momentum (GeV/c);TPC signal",40,0,4,100,0,400);
485   // fOutputList->Add(fHistPosTpcBeforeCut);
486   // fHistPosTpcAfterCut  = new TH2F ("h2PosTpcAfterCut","TPC signal (pos daughter) after cut;momentum (GeV/c);TPC signal",40,0,4,100,0,400);
487   // fOutputList->Add(fHistPosTpcAfterCut);
488   // fHistNegTpcBeforeCut = new TH2F ("h2NegTpcBeforeCut","TPC signal (neg daughter) before cut;momentum (GeV/c);TPC signal",40,0,4,100,0,400);
489   // fOutputList->Add(fHistNegTpcBeforeCut);
490   // fHistNegTpcAfterCut  = new TH2F ("h2NegTpcAfterCut","TPC signal (neg daughter) after cut;momentum (GeV/c);TPC signal",40,0,4,100,0,400);
491   // fOutputList->Add(fHistNegTpcAfterCut);
492
493   // // Histograms comparing offline and on-the-fly
494   // fHistGoodV0                  = new TH2F("h2GoodV0","0: all, 1: two daughters, 2: like-sign, 3: 80 clusters4: tpcrefit;id;Status",10,-.5,9.5,2,-.5,1.5);
495   // fOutputList->Add(fHistGoodV0);
496   // fHistCorrectSigns            = new TH2F ("h2CorrectSigns","0: correct, 1: swapped, 2: like-sign;sign;Status",3,-.5,2.5,2,-.5,1.5);
497   // fOutputList->Add(fHistCorrectSigns);
498   // fHistDcaPosToPrimVertex      = new TH2F("h2DcaPosToPrimVertex", "Positive V0 daughter;dca(cm);Status",500,0,5,2,-0.5,1.5);
499   // fOutputList->Add(fHistDcaPosToPrimVertex);
500   // fHistDcaNegToPrimVertex      = new TH2F("h2DcaNegToPrimVertex", "Negative V0 daughter;dca(cm);Status",500,0,5,2,-0.5,1.5);
501   // fOutputList->Add(fHistDcaNegToPrimVertex);
502   // fHistDcaPosToPrimVertexZoom  = new TH2F("h2DcaPosToPrimVertexZoom", "Positive V0 daughter;dca(cm);Status",100,0,0.1,2,-0.5,1.5);
503   // fOutputList->Add(fHistDcaPosToPrimVertexZoom);
504   // fHistDcaNegToPrimVertexZoom  = new TH2F("h2DcaNegToPrimVertexZoom", "Negative V0 daughter;dca(cm);Status",100,0,0.1,2,-0.5,1.5);
505   // fOutputList->Add(fHistDcaNegToPrimVertexZoom);
506   // fHistRadiusV0                = new TH2F("h2RadiusV0", "Radius;Radius(cm);Status",1000,0,100,2,-0.5,1.5);
507   // fOutputList->Add(fHistRadiusV0);
508   // fHistDecayLengthV0           = new TH2F("h2DecayLengthV0", "V0s decay Length;decay length(cm);Status", 200, 0, 100,2,-0.5,1.5);
509   // fOutputList->Add(fHistDecayLengthV0);
510   // fHistDcaV0Daughters          = new TH2F("h2DcaV0Daughters", "DCA between daughters;dca(cm);Status", 160, 0, 4,2,-0.5,1.5);
511   // fOutputList->Add(fHistDcaV0Daughters);
512   // fHistChi2                    = new TH2F("h2Chi2", "V0s chi2;chi2;Status", 12, 0, 1.2,2,-0.5,1.5);
513   // fOutputList->Add(fHistChi2);
514   // fHistCosPointAngle           = new TH2F("h2CosPointAngle", "Cosine of V0's pointing angle", 100,0,1,2,-0.5,1.5);
515   // fOutputList->Add(fHistCosPointAngle);
516   // fHistCosPointAngleZoom       = new TH2F("h2CosPointAngleZoom", "Cosine of V0's pointing angle", 100,0.9,1,2,-0.5,1.5);
517   // fOutputList->Add(fHistCosPointAngleZoom);
518
519   //
520   // V0 offline distributons
521   //
522   
523   // Invariant mass distribution for the side band background
524   fHistSideBandOffLam = new TH1F ("h1SideBandOffLam","m_{inv}(#Lambda) w/o any cuts;m_{inv}(#Lambda)",nMinvBins,minvLowEdge,minvHiEdge);
525   fOutputList->Add(fHistSideBandOffLam);
526   fHistSideBandOffALam  = new TH1F ("h1SideBandOffALam","m_{inv}(#bar{#Lambda}) w/o any cuts;m_{inv}(#bar{#Lambda})",nMinvBins,minvLowEdge,minvHiEdge);
527   fOutputList->Add(fHistSideBandOffALam);
528
529   // Number of TPC clusters
530   fHistTPCNclsPosOffLam = new TH2F ("h2MassLamOffTPCNclsPos","m_{inv}(#Lambda) vs NTPCcls(pos);NTPCcls(pos);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
531   fOutputList->Add(fHistTPCNclsPosOffLam);
532   fHistTPCNclsNegOffLam = new TH2F ("h2MassLamOffTPCNclsNeg","m_{inv}(#Lambda) vs NTPCcls(neg);NTPCcls(neg);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
533   fOutputList->Add(fHistTPCNclsNegOffLam);
534   fHistTPCNclsPosOffALam = new TH2F ("h2MassALamOffTPCNclsPos","m_{inv}(#bar{#Lambda}) vs NTPCcls(pos);NTPCcls(pos);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
535   fOutputList->Add(fHistTPCNclsPosOffALam);
536   fHistTPCNclsNegOffALam = new TH2F ("h2MassALamOffTPCNclsNeg","m_{inv}(#bar{#Lambda}) vs NTPCcls(neg);NTPCcls(neg);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
537   fOutputList->Add(fHistTPCNclsNegOffALam);
538
539   // fHistPosNsigmaTpcOffLam  = new TH2F ("h2PosNsigmaTpcOffLam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);
540   // fOutputList->Add(fHistPosNsigmaTpcOffLam);
541   // fHistPosNsigmaTpcOffALam = new TH2F ("h2PosNsigmaTpcOffALam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);        
542   // fOutputList->Add(fHistPosNsigmaTpcOffALam);
543   // fHistNegNsigmaTpcOffLam  = new TH2F ("h2NegNsigmaTpcOffLam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);        
544   // fOutputList->Add(fHistNegNsigmaTpcOffLam);
545   // fHistNegNsigmaTpcOffALam = new TH2F ("h2NegNsigmaTpcOffALam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);        
546   // fOutputList->Add(fHistNegNsigmaTpcOffALam);
547   // fHistUseTofOffLam   = new TH2F ("h2UseTofOffLam","0: no tof or within 5sigma, 1: tof more than 5sigma;m_{inv};TOF",2,-.5,1.5,nMinvBins,minvLowEdge,minvHiEdge);        
548   // fOutputList->Add(fHistUseTofOffLam);
549   // fHistUseTofOffALam  = new TH2F ("h2UseTofOffALam","0: no tof or within 5sigma, 1: tof more than 5sigma;m_{inv};TOF",2,-.5,1.5,nMinvBins,minvLowEdge,minvHiEdge);        
550   // fOutputList->Add(fHistUseTofOffALam);
551
552   // // DCA of daughters to primary vertex
553   // fHistDcaPosOffLam = new TH2F ("h2DcaPosOffLam","m_{inv}(#Lambda) vs dca pos daughter;dca (cm);m_{inv}(p#pi^{-})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
554   // fOutputList->Add(fHistDcaPosOffLam);
555   // fHistDcaPosOffALam = new TH2F ("h2DcaPosOffALam","m_{inv}(#bar{#Lambda}) vs dca pos daughter;dca (cm);m_{inv}(#bar{p}#pi^{+})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
556   // fOutputList->Add(fHistDcaPosOffALam);
557   // fHistDcaNegOffLam = new TH2F ("h2DcaNegOffLam","m_{inv}(#Lambda) vs dca neg daughter;dca (cm);m_{inv}(p#pi^{-})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
558   // fOutputList->Add(fHistDcaNegOffLam);
559   // fHistDcaNegOffALam = new TH2F ("h2DcaNegOffALam","m_{inv}(#bar{#Lambda}) vs dca neg daughter;dca (cm);m_{inv}(#bar{p}#pi^{+})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
560   // fOutputList->Add(fHistDcaNegOffALam);
561
562   // // DCA of V0 daughters
563   // fHistDcaV0DaughtersOffLam = new TH2F ("h2DcaLamDaughtersOff","DCA of #Lambda daughters vs minv;dca(cm);minv",20,0,2,nMinvBins,minvLowEdge,minvHiEdge);
564   // fOutputList->Add(fHistDcaV0DaughtersOffLam);
565   // fHistDcaV0DaughtersOffALam = new TH2F ("h2DcaALamDaughtersOff","DCA of #bar{#Lambda} daughters vs minv;dca(cm);minv",20,0,2,nMinvBins,minvLowEdge,minvHiEdge);
566   // fOutputList->Add(fHistDcaV0DaughtersOffALam);
567
568   // // Cosine of pointing angle
569   // fHistCosPointLamOff = new TH2F ("h2CosPointLamOff","m_{inv}(#Lambda) vs cos(pointing angle);cos(pointing angle);m_{inv}(#Lambda)",10,0.99,1.0,nMinvBins,minvLowEdge,minvHiEdge);
570   // fOutputList->Add(fHistCosPointLamOff);
571   // fHistCosPointALamOff = new TH2F ("h2CosPointALamOff","m_{inv}(#bar{#Lambda}) vs cos(pointing angle);cos(pointing angle);m_{inv}(#bar{#Lambda})",10,0.99,1.0,nMinvBins,minvLowEdge,minvHiEdge);
572   // fOutputList->Add(fHistCosPointALamOff);
573   // fHistCosPointLamZoomOff = new TH2F ("h2CosPointLamZoomOff","m_{inv}(#Lambda) vs cos(pointing angle);cos(pointing angle);m_{inv}(#Lambda)",10,0.999,1.0,nMinvBins,minvLowEdge,minvHiEdge);
574   // fOutputList->Add(fHistCosPointLamZoomOff);
575   // fHistCosPointALamZoomOff = new TH2F ("h2CosPointALamZoomOff","m_{inv}(#bar{#Lambda}) vs cos(pointing angle);cos(pointing angle);m_{inv}(#bar{#Lambda})",10,0.999,1.0,nMinvBins,minvLowEdge,minvHiEdge);
576   // fOutputList->Add(fHistCosPointALamZoomOff);
577
578   // // Radius of V0 vertex position
579   // fHistV0RadiusLamOff = new TH2F ("h2V0RadiusLamOff","m_{inv}(#Lambda) vs V0 radius of V0 vertex;radius(cm);m_{inv}",20,0,10,nMinvBins,minvLowEdge,minvHiEdge);
580   // fOutputList->Add(fHistV0RadiusLamOff);
581   // fHistV0RadiusALamOff = new TH2F ("h2V0RadiusALamOff","m_{inv}(#bar{#Lambda}) vs V0 radius of V0 vertex;radius(cm);m_{inv}",20,0,10,nMinvBins,minvLowEdge,minvHiEdge);
582   // fOutputList->Add(fHistV0RadiusALamOff);
583
584   // // Decay length of V0
585   // fHistV0DecayLengthLamOff = new TH2F ("h2V0DecayLengthLamOff","m_{inv}(#Lambda) vs decay length of V0;decay length (cm);m_{inv}(#Lambda)",100,0,20,nMinvBins,minvLowEdge,minvHiEdge);
586   // fOutputList->Add(fHistV0DecayLengthLamOff);
587   // fHistV0DecayLengthALamOff = new TH2F ("h2V0DecayLengthALamOff","m_{inv}(#bar{#Lambda}) vs decay length of V0;decay length (cm);m_{inv}(#bar{#Lambda})",100,0,20,nMinvBins,minvLowEdge,minvHiEdge);
588   // fOutputList->Add(fHistV0DecayLengthALamOff);
589
590   // // DCA of primary vertex and V0
591   // fHistDcaV0PriVertexLamOff = new TH2F ("h2DcaV0PriVertexLamOff","m_{inv}(#Lambda) vs dca (V0 - prim. vertex);dca(cm);m_{inv}(#Lambda)",200,0,20,nMinvBins,minvLowEdge,minvHiEdge);
592   // fOutputList->Add(fHistDcaV0PriVertexLamOff);
593   // fHistDcaV0PriVertexALamOff = new TH2F ("h2DcaV0PriVertexALamOff","m_{inv}(#bar{#Lambda}) vs dca (V0 - prim. vertex);dca(cm);m_{inv}(#bar{#Lambda})",200,0,20,nMinvBins,minvLowEdge,minvHiEdge);
594   // fOutputList->Add(fHistDcaV0PriVertexALamOff);
595
596   // Invariant mass, invariant mass vs pt and y-pt
597   fHistMassLambdaOff            = new TH1F("h1MassLambdaOff", "#Lambda^{0} Offline candidates;M(p#pi^{-}) (GeV/c^{2});Counts", nMinvBins, minvLowEdge, minvHiEdge);
598   fOutputList->Add(fHistMassLambdaOff);
599   fHistMassAntiLambdaOff          = new TH1F("h1MassAntiLambdaOff", "#bar{#Lambda}^{0} Offline candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", nMinvBins, minvLowEdge, minvHiEdge);
600   fOutputList->Add(fHistMassAntiLambdaOff);
601   // fHistPtVsMassLambdaOff         = new TH2F("h2PtVsMassLambdaOff","#Lambda^{0} Offline candidates;p_{t} (GeV/c);M(p#pi^{-}) (GeV/c^{2})",100,0,10,nMinvBins, minvLowEdge, minvHiEdge);
602   // fOutputList->Add(fHistPtVsMassLambdaOff);
603   // fHistPtVsMassAntiLambdaOff     = new TH2F("h2PtVsMassAntiLambdaOff","#bar{#Lambda}^{0} Offline candidates;p_{t} (GeV/c);M(#bar{p}#pi^{+}) (GeV/c^{2})",100,0,10,nMinvBins, minvLowEdge, minvHiEdge);
604   // fOutputList->Add(fHistPtVsMassAntiLambdaOff);
605   // fHistPtVsYLambdaOff          = new TH2F("h2PtVsYLambdaOff", "#Lambda^{0} Offline candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
606   // fOutputList->Add(fHistPtVsYLambdaOff);
607   // fHistPtVsYAntiLambdaOff      = new TH2F("h2PtVsYAntiLambdaOff", "#bar{#Lambda}^{0} Offline candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
608   // fOutputList->Add(fHistPtVsYAntiLambdaOff);
609
610     // 3d y pt mass
611   fHistYPtMassLamOff = new TH3F ("h3YPtMassLamOff","m_{inv}(#Lambda) vs y and pt;y;pt;mass",30,-1.5,1.5,30,0,15,nMinvBins,minvLowEdge,minvHiEdge);
612   fOutputList->Add(fHistYPtMassLamOff);
613   fHistYPtMassALamOff = new TH3F ("h3YPtMassALamOff","m_{inv}(#bar{#Lambda}) vs y and pt;y;pt;mass",30,-1.5,1.5,30,0,15,nMinvBins,minvLowEdge,minvHiEdge);
614   fOutputList->Add(fHistYPtMassALamOff);
615
616   //
617   // V0 on-the-fly distributons
618   //
619   
620   // Invariant mass distribution for the side band background
621   fHistSideBandOnLam = new TH1F ("h1SideBandOnLam","m_{inv}(#Lambda) w/o any cuts;m_{inv}(#Lambda)",nMinvBins,minvLowEdge,minvHiEdge);
622   fOutputList->Add(fHistSideBandOnLam);
623   fHistSideBandOnALam  = new TH1F ("h1SideBandOnALam","m_{inv}(#bar{#Lambda}) w/o any cuts;m_{inv}(#bar{#Lambda})",nMinvBins,minvLowEdge,minvHiEdge);
624   fOutputList->Add(fHistSideBandOnALam);
625   
626   // // Like-sign
627   // fHistLikeSignOnLam = new TH2F ("h2MassLamOnLikeSign"," 0 = ok, 1 = swapped, 2 = like sign;sign;m_{inv} p#pi^{-}",3,-.5,2.5,nMinvBins,minvLowEdge,minvHiEdge);
628   // fOutputList->Add(fHistLikeSignOnLam);
629   // fHistLikeSignOnALam = new TH2F ("h2MassALamOnLikeSign"," 0 = ok, 1 = swapped, 2= like sign;sign;m_{inv} #bar{p}#pi^{+}",3,-.5,2.5,nMinvBins,minvLowEdge,minvHiEdge);
630   // fOutputList->Add(fHistLikeSignOnALam);
631
632   // Number of TPC clusters
633   fHistTPCNclsPosOnLam = new TH2F ("h2MassLamOnTPCNclsPos","m_{inv}(#Lambda) vs NTPCcls(pos);NTPCcls(pos);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
634   fOutputList->Add(fHistTPCNclsPosOnLam);
635   fHistTPCNclsNegOnLam = new TH2F ("h2MassLamOnTPCNclsNeg","m_{inv}(#Lambda) vs NTPCcls(neg);NTPCcls(neg);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
636   fOutputList->Add(fHistTPCNclsNegOnLam);
637   fHistTPCNclsPosOnALam = new TH2F ("h2MassALamOnTPCNclsPos","m_{inv}(#bar{#Lambda}) vs NTPCcls(pos);NTPCcls(pos);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
638   fOutputList->Add(fHistTPCNclsPosOnALam);
639   fHistTPCNclsNegOnALam = new TH2F ("h2MassALamOnTPCNclsNeg","m_{inv}(#bar{#Lambda}) vs NTPCcls(neg);NTPCcls(neg);minv",18,0,180,nMinvBins,minvLowEdge,minvHiEdge);
640   fOutputList->Add(fHistTPCNclsNegOnALam);
641
642 //  fHistPosNsigmaTpcOnLam  = new TH2F ("h2PosNsigmaTpcOnLam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);
643 //   fOutputList->Add(fHistPosNsigmaTpcOnLam);
644 //   fHistPosNsigmaTpcOnALam = new TH2F ("h2PosNsigmaTpcOnALam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);        
645 //   fOutputList->Add(fHistPosNsigmaTpcOnALam);
646 //   fHistNegNsigmaTpcOnLam  = new TH2F ("h2NegNsigmaTpcOnLam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);        
647 //   fOutputList->Add(fHistNegNsigmaTpcOnLam);
648 //   fHistNegNsigmaTpcOnALam = new TH2F ("h2NegNsigmaTpcOnALam","minv vs nsigma tpc pos;Nsigma TPC;m_{inv}",50,0,5,nMinvBins,minvLowEdge,minvHiEdge);        
649 //   fOutputList->Add(fHistNegNsigmaTpcOnALam);
650 //   fHistUseTofOnLam   = new TH2F ("h2UseTofOnLam","0: no tof or within 5sigma, 1: tof more than 5sigma;m_{inv};TOF",2,-.5,1.5,nMinvBins,minvLowEdge,minvHiEdge);        
651 //   fOutputList->Add(fHistUseTofOnLam);
652 //   fHistUseTofOnALam  = new TH2F ("h2UseTofOnALam","0: no tof or within 5sigma, 1: tof more than 5sigma;m_{inv};TOF",2,-.5,1.5,nMinvBins,minvLowEdge,minvHiEdge);        
653 //   fOutputList->Add(fHistUseTofOnALam);
654
655 //   // DCA of daughters to primary vertex
656 //   fHistDcaPosOnLam = new TH2F ("h2DcaPosOnLam","m_{inv}(#Lambda) vs dca pos daughter;dca (cm);m_{inv}(p#pi^{-})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
657 //   fOutputList->Add(fHistDcaPosOnLam);
658 //   fHistDcaPosOnALam = new TH2F ("h2DcaPosOnALam","m_{inv}(#bar{#Lambda}) vs dca pos daughter;dca (cm);m_{inv}(#bar{p}#pi^{+})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
659 //   fOutputList->Add(fHistDcaPosOnALam);
660 //   fHistDcaNegOnLam = new TH2F ("h2DcaNegOnLam","m_{inv}(#Lambda) vs dca neg daughter;dca (cm);m_{inv}(p#pi^{-})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
661 //   fOutputList->Add(fHistDcaNegOnLam);
662 //   fHistDcaNegOnALam = new TH2F ("h2DcaNegOnALam","m_{inv}(#bar{#Lambda}) vs dca neg daughter;dca (cm);m_{inv}(#bar{p}#pi^{+})",50,0,0.5,nMinvBins,minvLowEdge,minvHiEdge);
663 //   fOutputList->Add(fHistDcaNegOnALam);
664
665 //   // DCA of V0 daughters
666 //   fHistDcaV0DaughtersOnLam = new TH2F ("h2DcaLamDaughtersOn","DCA of #Lambda daughters vs minv;dca(cm);minv",20,0,2,nMinvBins,minvLowEdge,minvHiEdge);
667 //   fOutputList->Add(fHistDcaV0DaughtersOnLam);
668 //   fHistDcaV0DaughtersOnALam = new TH2F ("h2DcaALamDaughtersOn","DCA of #bar{#Lambda} daughters vs minv;dca(cm);minv",20,0,2,nMinvBins,minvLowEdge,minvHiEdge);
669 //   fOutputList->Add(fHistDcaV0DaughtersOnALam);
670
671 // // Cosine of pointing angle
672 //   fHistCosPointLamOn = new TH2F ("h2CosPointLamOn","m_{inv}(#Lambda) vs cos(pointing angle);cos(pointing angle);m_{inv}(#Lambda)",10,0.99,1.0,nMinvBins,minvLowEdge,minvHiEdge);
673 //   fOutputList->Add(fHistCosPointLamOn);
674 //   fHistCosPointALamOn = new TH2F ("h2CosPointALamOn","m_{inv}(#bar{#Lambda}) vs cos(pointing angle);cos(pointing angle);m_{inv}(#bar{#Lambda})",10,0.99,1.0,nMinvBins,minvLowEdge,minvHiEdge);
675 //   fOutputList->Add(fHistCosPointALamOn);
676 //   fHistCosPointLamZoomOn = new TH2F ("h2CosPointLamZoomOn","m_{inv}(#Lambda) vs cos(pointing angle);cos(pointing angle);m_{inv}(#Lambda)",10,0.999,1.0,nMinvBins,minvLowEdge,minvHiEdge);
677 //   fOutputList->Add(fHistCosPointLamZoomOn);
678 //   fHistCosPointALamZoomOn = new TH2F ("h2CosPointALamZoomOn","m_{inv}(#bar{#Lambda}) vs cos(pointing angle);cos(pointing angle);m_{inv}(#bar{#Lambda})",10,0.999,1.0,nMinvBins,minvLowEdge,minvHiEdge);
679 //   fOutputList->Add(fHistCosPointALamZoomOn);
680
681 //   // Radius of V0 vertex position
682 //   fHistV0RadiusLamOn = new TH2F ("h2V0RadiusLamOn","m_{inv}(#Lambda) vs V0 radius of V0 vertex;radius(cm);m_{inv}",20,0,10,nMinvBins,minvLowEdge,minvHiEdge);
683 //   fOutputList->Add(fHistV0RadiusLamOn);
684 //   fHistV0RadiusALamOn = new TH2F ("h2V0RadiusALamOn","m_{inv}(#bar{#Lambda}) vs V0 radius of V0 vertex;radius(cm);m_{inv}",20,0,10,nMinvBins,minvLowEdge,minvHiEdge);
685 //   fOutputList->Add(fHistV0RadiusALamOn);
686
687 //   // Decay length of V0
688 //   fHistV0DecayLengthLamOn = new TH2F ("h2V0DecayLengthLamOn","m_{inv}(#Lambda) vs decay length of V0;decay length (cm);m_{inv}(#Lambda)",100,0,20,nMinvBins,minvLowEdge,minvHiEdge);
689 //   fOutputList->Add(fHistV0DecayLengthLamOn);
690 //   fHistV0DecayLengthALamOn = new TH2F ("h2V0DecayLengthALamOn","m_{inv}(#bar{#Lambda}) vs decay length of V0;decay length (cm);m_{inv}(#bar{#Lambda})",100,0,20,nMinvBins,minvLowEdge,minvHiEdge);
691 //   fOutputList->Add(fHistV0DecayLengthALamOn);
692
693 //   // DCA of primary vertex and V0
694 //   fHistDcaV0PriVertexLamOn = new TH2F ("h2DcaV0PriVertexLamOn","m_{inv}(#Lambda) vs dca (V0 - prim. vertex);dca(cm);m_{inv}(#Lambda)",200,0,20,nMinvBins,minvLowEdge,minvHiEdge);
695 //   fOutputList->Add(fHistDcaV0PriVertexLamOn);
696 //   fHistDcaV0PriVertexALamOn = new TH2F ("h2DcaV0PriVertexALamOn","m_{inv}(#bar{#Lambda}) vs dca (V0 - prim. vertex);dca(cm);m_{inv}(#bar{#Lambda})",200,0,20,nMinvBins,minvLowEdge,minvHiEdge);
697 //   fOutputList->Add(fHistDcaV0PriVertexALamOn);
698
699 //   // Chi2 of TPC clusters
700 //   fHistChi2TPCPosLamOn = new TH2F ("h2Chi2TPCPosLamOn","m_{inv}(#Lambda) vs TPC #Chi^{2} / ndf pos daughter;TPC #Chi^{2}/ndf pos daughter;m_{inv}(#Lambda)",100,0.,10.,nMinvBins,minvLowEdge,minvHiEdge);
701 //   fOutputList->Add(fHistChi2TPCPosLamOn);
702 //   fHistChi2TPCPosALamOn = new TH2F ("h2Chi2TPCPosALamOn","m_{inv}(#bar{#Lambda}) vs TPC #Chi^{2} / ndf pos daughter;TPC #Chi^{2}/ndf pos daughter;m_{inv}(#bar{#Lambda})",100,0.,10.,nMinvBins,minvLowEdge,minvHiEdge);
703 //   fOutputList->Add(fHistChi2TPCPosALamOn);
704 //   fHistChi2TPCNegLamOn = new TH2F ("h2Chi2TPCNegLamOn","m_{inv}(#Lambda) vs TPC #Chi^{2} / ndf neg daughter;TPC #Chi^{2}/ndf neg daughter;m_{inv}(#Lambda)",100,0.,10.,nMinvBins,minvLowEdge,minvHiEdge);
705 //   fOutputList->Add(fHistChi2TPCNegLamOn);
706 //   fHistChi2TPCNegALamOn = new TH2F ("h2Chi2TPCNegALamOn","m_{inv}(#bar{#Lambda}) vs TPC #Chi^{2} / ndf neg daughter;TPC #Chi^{2}/ndf neg daughter;m_{inv}(#bar{#Lambda})",100,0.,10.,nMinvBins,minvLowEdge,minvHiEdge);
707 //   fOutputList->Add(fHistChi2TPCNegALamOn);
708 //   // Invariant mass with TPC only tracks
709 //   fHistMinvTPConlyLamOn = new TH1F ("h1MinvTPConlyLamOn","m_{inv}(#Lambda) using TPC only daughters;m_{inv}(p#pi^{-})[GeV/c^{2}]",nMinvBins,minvLowEdge,minvHiEdge);
710 //   fOutputList->Add(fHistMinvTPConlyLamOn);
711 //   fHistMinvTPConlyALamOn = new TH1F ("h1MinvTPConlyALamOn","m_{inv}(#bar{#Lambda}) using TPC only daughters;m_{inv}(#bar{p}#pi^{+})[GeV/c^{2}]",nMinvBins,minvLowEdge,minvHiEdge);
712 //   fOutputList->Add(fHistMinvTPConlyALamOn);
713
714   // Invariant mass, invariant mass vs pt and y-pt
715   fHistMassLambdaOn            = new TH1F("h1MassLambdaOn", "#Lambda^{0} Online candidates;M(p#pi^{-}) (GeV/c^{2});Counts", nMinvBins, minvLowEdge, minvHiEdge);
716   fOutputList->Add(fHistMassLambdaOn);
717   fHistMassAntiLambdaOn          = new TH1F("h1MassAntiLambdaOn", "#bar{#Lambda}^{0} Online candidates;M(#bar{p}#pi^{+}) (GeV/c^{2});Counts", nMinvBins, minvLowEdge, minvHiEdge);
718   fOutputList->Add(fHistMassAntiLambdaOn);
719   // fHistPtVsMassLambdaOn         = new TH2F("h2PtVsMassLambdaOn","#Lambda^{0} Online candidates;p_{t} (GeV/c);M(p#pi^{-}) (GeV/c^{2})",100,0,10,nMinvBins, minvLowEdge, minvHiEdge);
720   // fOutputList->Add(fHistPtVsMassLambdaOn);
721   // fHistPtVsMassAntiLambdaOn     = new TH2F("h2PtVsMassAntiLambdaOn","#bar{#Lambda}^{0} Online candidates;p_{t} (GeV/c);M(#bar{p}#pi^{+}) (GeV/c^{2})",100,0,10,nMinvBins, minvLowEdge, minvHiEdge);
722   // fOutputList->Add(fHistPtVsMassAntiLambdaOn);
723   // fHistPtVsYLambdaOn          = new TH2F("h2PtVsYLambdaOn", "#Lambda^{0} Online candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
724   // fOutputList->Add(fHistPtVsYLambdaOn);
725   // fHistPtVsYAntiLambdaOn      = new TH2F("h2PtVsYAntiLambdaOn", "#bar{#Lambda}^{0} Online candidates;p_{t} (GeV/c);rapidity",30,0,15,30,-1.5,1.5);
726   // fOutputList->Add(fHistPtVsYAntiLambdaOn);
727
728   // 3d y pt mass
729   fHistYPtMassLamOn = new TH3F ("h3YPtMassLamOn","m_{inv}(#Lambda) vs y and pt;y;pt;mass",30,-1.5,1.5,30,0,15,nMinvBins,minvLowEdge,minvHiEdge);
730   fOutputList->Add(fHistYPtMassLamOn);
731   fHistYPtMassALamOn = new TH3F ("h3YPtMassALamOn","m_{inv}(#bar{#Lambda}) vs y and pt;y;pt;mass",30,-1.5,1.5,30,0,15,nMinvBins,minvLowEdge,minvHiEdge);
732   fOutputList->Add(fHistYPtMassALamOn);
733
734   // // Momentum difference of standard (on-the-fly/offline) V0 and TPC only V0
735   // Int_t nPBins=200; Float_t AbsPRange=1.;
736   // fHistMomDiffLam = new TH3F ("h3MomDiffLam","momentum difference #DeltaP standard V0 / TPConly V0 #Lambda;#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
737   //                          ,nPBins,-AbsPRange,AbsPRange
738   //                          ,nPBins,-AbsPRange,AbsPRange
739   //                          ,nPBins,-AbsPRange,AbsPRange);
740   // fOutputList->Add(fHistMomDiffLam);
741   // fHistMomDiffALam = new TH3F ("h3MomDiffALam","momentum difference #DeltaP standard V0 / TPConly V0 #bar{#Lamdba};#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
742   //                           ,nPBins,-AbsPRange,AbsPRange
743   //                           ,nPBins,-AbsPRange,AbsPRange
744   //                           ,nPBins,-AbsPRange,AbsPRange);
745   // fOutputList->Add(fHistMomDiffALam);
746   // fHistMomDiffBgLam = new TH3F ("h3MomDiffBgLam","momentum difference #DeltaP standard V0 / TPConly V0 Bg#Lambda;#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
747   //                            ,nPBins,-AbsPRange,AbsPRange
748   //                            ,nPBins,-AbsPRange,AbsPRange
749   //                            ,nPBins,-AbsPRange,AbsPRange);
750   // fOutputList->Add(fHistMomDiffBgLam);
751   // fHistMomDiffBgALam = new TH3F ("h3MomDiffBgALam","momentum difference #DeltaP standard V0 / TPConly V0 Bg#bar{#Lambda};#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
752   //                             ,nPBins,-AbsPRange,AbsPRange
753   //                             ,nPBins,-AbsPRange,AbsPRange
754   //                             ,nPBins,-AbsPRange,AbsPRange);
755   // fOutputList->Add(fHistMomDiffBgALam);
756
757   // // Same momentum difference with rejecting tracks with SPD hits
758   // fHistMomDiffWoSPDLam = new TH3F ("h3MomDiffWoSPDLam","momentum difference #DeltaP standard V0 / TPConly V0 #Lambda;#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
759   //                          ,nPBins,-AbsPRange,AbsPRange
760   //                          ,nPBins,-AbsPRange,AbsPRange
761   //                          ,nPBins,-AbsPRange,AbsPRange);
762   // fOutputList->Add(fHistMomDiffWoSPDLam);
763   // fHistMomDiffWoSPDALam = new TH3F ("h3MomDiffWoSPDALam","momentum difference #DeltaP standard V0 / TPConly V0 #bar{#Lamdba};#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
764   //                           ,nPBins,-AbsPRange,AbsPRange
765   //                           ,nPBins,-AbsPRange,AbsPRange
766   //                           ,nPBins,-AbsPRange,AbsPRange);
767   // fOutputList->Add(fHistMomDiffWoSPDALam);
768   // fHistMomDiffWoSPDBgLam = new TH3F ("h3MomDiffWoSPDBgLam","momentum difference #DeltaP standard V0 / TPConly V0 Bg#Lambda;#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
769   //                            ,nPBins,-AbsPRange,AbsPRange
770   //                            ,nPBins,-AbsPRange,AbsPRange
771   //                            ,nPBins,-AbsPRange,AbsPRange);
772   // fOutputList->Add(fHistMomDiffWoSPDBgLam);
773   // fHistMomDiffWoSPDBgALam = new TH3F ("h3MomDiffWoSPDBgALam","momentum difference #DeltaP standard V0 / TPConly V0 Bg#bar{#Lambda};#DeltaP_{x}[GeV/c];#DeltaP_{y}[GeV/c];#DeltaP_{z}[GeV/c]"
774   //                             ,nPBins,-AbsPRange,AbsPRange
775   //                             ,nPBins,-AbsPRange,AbsPRange
776   //                             ,nPBins,-AbsPRange,AbsPRange);
777   // fOutputList->Add(fHistMomDiffWoSPDBgALam);
778
779   //
780   // Distributions for the primaries 
781   //
782   // Shared clusters
783   fPriHistShare = new TH1F ("h1PriShare","Shared clusters, primaries;#shared clusters;counts",
784                             160,0,160);
785   fOutputPrimaries->Add(fPriHistShare);
786
787   // Nsigma TOF distribution when forcing TOF
788   // fPriHistPosNsigmaTof = new TH1F ("h1PosNsigmaTof","Nsigma distribution for positives;n_{#sigma,TOF}(proton);counts",200,-50,50);
789   // fOutputPrimaries->Add(fPriHistPosNsigmaTof);
790   fPriHistPosNsigmaTofVsP = new TH2F ("h2PosNsigmaTofVsP","Nsigma distribution for positives;total momentum [GeV/c];n_{#sigma,TOF}(proton);counts",20,0,5,200,-50,50);
791   fOutputPrimaries->Add(fPriHistPosNsigmaTofVsP);
792   fPriHistPosNsigmaTofVsPt   = new TH2F ("h2PosNsigmaTofVsPt","Nsigma distribution for positives;transverse momentum [GeV/c];n_{#sigma,TOF}(proton);counts",20,0,5,200,-50,50);
793   fOutputPrimaries->Add(fPriHistPosNsigmaTofVsPt);
794
795   // fPriHistNegNsigmaTof = new TH1F ("h1NegNsigmaTof","Nsigma distribution for negatives;n_{#sigma,TOF}(anti-proton);counts",200,-50,50);
796   // fOutputPrimaries->Add(fPriHistNegNsigmaTof);
797   fPriHistNegNsigmaTofVsP = new TH2F ("h2NegNsigmaTofVsP","Nsigma distribution for negatives;total momentum [GeV/c];n_{#sigma,TOF}(anti-proton);counts",20,0,5,200,-50,50);
798   fOutputPrimaries->Add(fPriHistNegNsigmaTofVsP);
799   fPriHistNegNsigmaTofVsPt   = new TH2F ("h2NegNsigmaTofVsPt","Nsigma distribution for negatives;transverse momentum [GeV/c];n_{#sigma,TOF}(anti-proton);counts",20,0,5,200,-50,50);
800   fOutputPrimaries->Add(fPriHistNegNsigmaTofVsPt);
801   fPriHistTOFsignalPosVsP = new TH2F ("h2TOFsignalPosVsP","tof signal vs p (positives);p [GeV/c];t_{meas} - t_{0} - t_{expected} [ps]",20,0.0,5.0,120,-10000.0,5000.0);
802   fOutputPrimaries->Add(fPriHistTOFsignalPosVsP);
803   fPriHistTOFsignalPosVsPt = new TH2F ("h2TOFsignalPosVsPt","tof signal vs pt (positives);pt [GeV/c];t_{meas} - t_{0} - t_{expected} [ps]",20,0.0,5.0,120,-10000.0,5000.0);
804   fOutputPrimaries->Add(fPriHistTOFsignalPosVsPt);
805   fPriHistTOFsignalNegVsP = new TH2F ("h2TOFsignalNegVsP","tof signal vs p (negatives);p [GeV/c];t_{meas} - t_{0} - t_{expected} [ps]",20,0.0,5.0,120,-10000.0,5000.0);
806   fOutputPrimaries->Add(fPriHistTOFsignalNegVsP);
807   fPriHistTOFsignalNegVsPt = new TH2F ("h2TOFsignalNegVsPt","tof signal vs pt (negatives);pt [GeV/c];t_{meas} - t_{0} - t_{expected} [ps]",20,0.0,5.0,120,-10000.0,5000.0);
808   fOutputPrimaries->Add(fPriHistTOFsignalNegVsPt);
809   // Hybrid analysis
810   fPriHistHybridTOFsigPosWoTPC = new TH1F ("h1HybridTOFsigPosWoTPC","tof signal pos (p=.75-1.0GeV) w/o dedx sel.;t_{meas} - t_{0} - t_{expected} [ps]",120,-10000.0,5000.0);
811   fOutputPrimaries->Add(fPriHistHybridTOFsigPosWoTPC);
812   fPriHistHybridTOFsigPosTPCok = new TH1F ("h1HybridTOFsigPosTPCok","tof signal pos (p=.75-1.0GeV) with dedx sel.;t_{meas} - t_{0} - t_{expected} [ps]",120,-10000.0,5000.0);
813   fOutputPrimaries->Add(fPriHistHybridTOFsigPosTPCok);
814   fPriHistHybridTOFsigNegWoTPC = new TH1F ("h1HybridTOFsigNegWoTPC","tof signal neg (p=.75-1.0GeV) w/o dedx sel.;t_{meas} - t_{0} - t_{expected} [ps]",120,-10000.0,5000.0);
815   fOutputPrimaries->Add(fPriHistHybridTOFsigNegWoTPC);
816   fPriHistHybridTOFsigNegTPCok = new TH1F ("h1HybridTOFsigNegTPCok","tof signal neg (p=.75-1.0GeV) with dedx sel.;t_{meas} - t_{0} - t_{expected} [ps]",120,-10000.0,5000.0);
817   fOutputPrimaries->Add(fPriHistHybridTOFsigNegTPCok);
818   // dEdx analysis
819   // fPriHistHasTofPos = new TH1F ("h1HasTofPos","Positives: 0 = no TOF, 1 = TOFpid bit there",2,-.5,1.5);
820   // fOutputPrimaries->Add(fPriHistHasTofPos);
821   fPriHistTPCsignalPos = new TH2F ("h2TPCsignalPos","TPC signal for positives;p_{tot};dEdx",40,0,4,100,0,400);
822   fOutputPrimaries->Add(fPriHistTPCsignalPos);
823   // fPriHistNsigmaTPCPos = new TH2F ("h2NsigmaTPCPos","Nsigma TPC for positives;p_{tot};N_{#sigma}",40,0,4,100,-5.0,5.0);
824   // fOutputPrimaries->Add(fPriHistNsigmaTPCPos);
825   // fPriHistTPCsignalTOFcutPos = new TH2F ("h2TPCsignalTOFcutPos","TPC signal for positives using a +/- 10sigma TOF cut;p_{tot};dEdx",40,0.0,4.0,100,0.0,400.0);
826   // fOutputPrimaries->Add(fPriHistTPCsignalTOFcutPos);
827   // fPriHistNsigmaTPCTOFcutPos = new TH2F ("h2NsigmaTPCTOFcutPos","Nsigma TPC for positives using a +/- 10sigma TOF cut;p_{tot};N_{#sigma}",40,0.0,4.0,100,-5.0,5.0);
828   // fOutputPrimaries->Add(fPriHistNsigmaTPCTOFcutPos);
829
830   // fPriHistHasTofNeg = new TH1F ("h1HasTofNeg","Negatives: 0 = no TOF, 1 = TOFpid bit there",2,-.5,1.5);
831   // fOutputPrimaries->Add(fPriHistHasTofNeg);
832   fPriHistTPCsignalNeg = new TH2F ("h2TPCsignalNeg","TPC signal for negatives;p_{tot};dEdx",40,0.0,4.0,100,0.0,400.0);
833   fOutputPrimaries->Add(fPriHistTPCsignalNeg);
834   // fPriHistNsigmaTPCNeg = new TH2F ("h2NsigmaTPCNeg","Nsigma TPC for negatives;p_{tot};N_{#sigma}",40,0.0,4.0,100,-5.0,5.0);
835   // fOutputPrimaries->Add(fPriHistNsigmaTPCNeg);
836   // fPriHistTPCsignalTOFcutNeg = new TH2F ("h2TPCsignalTOFcutNeg","TPC signal for negatives using a +/- 10sigma TOF cut;p_{tot};dEdx",40,0.0,4.0,100,0.0,400.0);
837   // fOutputPrimaries->Add(fPriHistTPCsignalTOFcutNeg);
838   // fPriHistNsigmaTPCTOFcutNeg = new TH2F ("h2NsigmaTPCTOFcutNeg","Nsigma TPC for negatives using a +/- 10sigma TOF cut;p_{tot};N_{#sigma}",40,0.0,4.0,100,-5.0,5.0);
839   // fOutputPrimaries->Add(fPriHistNsigmaTPCTOFcutNeg);  
840
841   fPriHistTPCsignalLowPPos = new TH2F ("h2TPCsignalLowPPos","dEdx for low momenta, positives",20,0.1,0.3,3000,0,3000);
842   fOutputPrimaries->Add(fPriHistTPCsignalLowPPos);
843   fPriHistTPCsignalMedPPos = new TH2F ("h2TPCsignalMedPPos","dEdx for medium momenta, positives",60,0.3,0.9,500,0,500);
844   fOutputPrimaries->Add(fPriHistTPCsignalMedPPos);
845   fPriHistTPCsignalHigPPos = new TH2F ("h2TPCsignalHigPPos","dEdx for high momenta, positives",100,0.9,1.9,120,0,120);
846   fOutputPrimaries->Add(fPriHistTPCsignalHigPPos);
847   fPriHistTPCsignalLowPNeg = new TH2F ("h2TPCsignalLowPNeg","dEdx for low momenta, negatives",20,0.1,0.3,3000,0,3000);
848   fOutputPrimaries->Add(fPriHistTPCsignalLowPNeg);
849   fPriHistTPCsignalMedPNeg = new TH2F ("h2TPCsignalMedPNeg","dEdx for medium momenta, negatives",60,0.3,0.9,500,0,500);
850   fOutputPrimaries->Add(fPriHistTPCsignalMedPNeg);
851   fPriHistTPCsignalHigPNeg = new TH2F ("h2TPCsignalHigPNeg","dEdx for high momenta, negatives",100,0.9,1.9,120,0,120);
852   fOutputPrimaries->Add(fPriHistTPCsignalHigPNeg);
853   
854   //  Common for all protons
855
856   // DCA xy distribution to determine primaries, secondaries from weak decay and secondaries from material
857   fPriHistDCAxyYPtPro = new TH3F ("h3DCAxyYPtPro","DCAxy vs (y,pt) protons",100,-3.,3.,30,-1.5,1.5,14,0.,3.5);
858   fOutputPrimaries->Add(fPriHistDCAxyYPtPro);
859   fPriHistDCAxyYPtAPro = new TH3F ("h3DCAxyYPtAPro","DCAxy vs (y,pt) anti-protons",100,-3.,3.,30,-1.5,1.5,14,0.,3.5);
860   fOutputPrimaries->Add(fPriHistDCAxyYPtAPro);
861
862   //  2 particle histograms fOutput2Part
863   // Common binning for TTR
864   Int_t nDistBins=200;
865   Float_t distLow=0.,distHig=20.;
866   //  Two-track resolution: real events
867   // f2HistLamLamMeanMinDistProReal = new TH2F ("h2LamLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
868   // fOutput2Part->Add(f2HistLamLamMeanMinDistProReal);
869   // f2HistLamLamMeanMinDistPioReal = new TH2F ("h2LamLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
870   // fOutput2Part->Add(f2HistLamLamMeanMinDistPioReal);
871   // f2HistLamProMeanMinDistProReal = new TH2F ("h2LamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
872   // fOutput2Part->Add(f2HistLamProMeanMinDistProReal);
873   // f2HistALamALamMeanMinDistAProReal = new TH2F ("h2ALamALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
874   // fOutput2Part->Add(f2HistALamALamMeanMinDistAProReal);
875   // f2HistALamALamMeanMinDistPioReal = new TH2F ("h2ALamALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
876   // fOutput2Part->Add(f2HistALamALamMeanMinDistPioReal);
877   // f2HistALamAProMeanMinDistAProReal = new TH2F ("h2ALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
878   // fOutput2Part->Add(f2HistALamAProMeanMinDistAProReal);
879
880   // f2HistSftLamLamMeanMinDistProReal = new TH2F ("h2SftLamLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
881   // fOutput2Part->Add(f2HistSftLamLamMeanMinDistProReal);
882   // f2HistSftLamLamMeanMinDistPioReal = new TH2F ("h2SftLamLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
883   // fOutput2Part->Add(f2HistSftLamLamMeanMinDistPioReal);
884   // f2HistSftLamProMeanMinDistProReal = new TH2F ("h2SftLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
885   // fOutput2Part->Add(f2HistSftLamProMeanMinDistProReal);
886   // f2HistSftALamALamMeanMinDistAProReal = new TH2F ("h2SftALamALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
887   // fOutput2Part->Add(f2HistSftALamALamMeanMinDistAProReal);
888   // f2HistSftALamALamMeanMinDistPioReal = new TH2F ("h2SftALamALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
889   // fOutput2Part->Add(f2HistSftALamALamMeanMinDistPioReal);
890   // f2HistSftALamAProMeanMinDistAProReal = new TH2F ("h2SftALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
891   // fOutput2Part->Add(f2HistSftALamAProMeanMinDistAProReal);
892
893   // f2HistSftIrocLamLamMeanMinDistProReal = new TH2F ("h2SftIrocLamLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
894   // fOutput2Part->Add(f2HistSftIrocLamLamMeanMinDistProReal);
895   // f2HistSftIrocLamLamMeanMinDistPioReal = new TH2F ("h2SftIrocLamLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
896   // fOutput2Part->Add(f2HistSftIrocLamLamMeanMinDistPioReal);
897   // f2HistSftIrocLamProMeanMinDistProReal = new TH2F ("h2SftIrocLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
898   // fOutput2Part->Add(f2HistSftIrocLamProMeanMinDistProReal);
899   // f2HistSftIrocALamALamMeanMinDistAProReal = new TH2F ("h2SftIrocALamALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
900   // fOutput2Part->Add(f2HistSftIrocALamALamMeanMinDistAProReal);
901   // f2HistSftIrocALamALamMeanMinDistPioReal = new TH2F ("h2SftIrocALamALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
902   // fOutput2Part->Add(f2HistSftIrocALamALamMeanMinDistPioReal);
903   // f2HistSftIrocALamAProMeanMinDistAProReal = new TH2F ("h2SftIrocALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
904   // fOutput2Part->Add(f2HistSftIrocALamAProMeanMinDistAProReal);
905
906   // f2HistSftOrocLamLamMeanMinDistProReal = new TH2F ("h2SftOrocLamLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
907   // fOutput2Part->Add(f2HistSftOrocLamLamMeanMinDistProReal);
908   // f2HistSftOrocLamLamMeanMinDistPioReal = new TH2F ("h2SftOrocLamLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
909   // fOutput2Part->Add(f2HistSftOrocLamLamMeanMinDistPioReal);
910   // f2HistSftOrocLamProMeanMinDistProReal = new TH2F ("h2SftOrocLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
911   // fOutput2Part->Add(f2HistSftOrocLamProMeanMinDistProReal);
912   // f2HistSftOrocALamALamMeanMinDistAProReal = new TH2F ("h2SftOrocALamALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
913   // fOutput2Part->Add(f2HistSftOrocALamALamMeanMinDistAProReal);
914   // f2HistSftOrocALamALamMeanMinDistPioReal = new TH2F ("h2SftOrocALamALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
915   // fOutput2Part->Add(f2HistSftOrocALamALamMeanMinDistPioReal);
916   // f2HistSftOrocALamAProMeanMinDistAProReal = new TH2F ("h2SftOrocALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
917   // fOutput2Part->Add(f2HistSftOrocALamAProMeanMinDistAProReal);
918   
919   // Mt of the pairs
920   Int_t nMtBins=25;
921   Float_t mtLow=1.0,mtHig=3.5;
922   // f2HistMtLamLamReal = new TH1F("h1MtLamLamReal"
923   //                            ,"m_{t}(#Lambda #Lambda);m{t} [GeV];counts"
924   //                            ,nMtBins,mtLow,mtHig);
925   // fOutput2Part->Add(f2HistMtLamLamReal);
926   f2HistMtLamProReal = new TH1F("h1MtLamProReal"
927                                 ,"m_{t}(p #Lambda);m{t} [GeV];counts"
928                                 ,nMtBins,mtLow,mtHig);
929   fOutput2Part->Add(f2HistMtLamProReal);
930   // f2HistMtALamALamReal = new TH1F("h1MtALamALamReal"
931   //                              ,"m_{t}(#bar{#Lambda} #bar{#Lambda});m{t} [GeV];counts"
932   //                              ,nMtBins,mtLow,mtHig);
933   // fOutput2Part->Add(f2HistMtALamALamReal);
934   f2HistMtALamAProReal  = new TH1F("h1MtALamAProReal"
935                                    ,"m_{t}(#bar{#Lambda} #bar{#Lambda});m{t} [GeV];counts"
936                                    ,nMtBins,mtLow,mtHig);
937   fOutput2Part->Add(f2HistMtALamAProReal);
938   // The same only filling for low q pairs
939   // f2HistMtLowQLamLamReal = new TH1F("h1MtLowQLamLamReal"
940   //                            ,"m_{t}(#Lambda #Lambda);m{t} [GeV];counts"
941   //                            ,nMtBins,mtLow,mtHig);
942   // fOutput2Part->Add(f2HistMtLowQLamLamReal);
943   f2HistMtLowQLamProReal = new TH1F("h1MtLowQLamProReal"
944                                 ,"m_{t}(p #Lambda);m{t} [GeV];counts"
945                                 ,nMtBins,mtLow,mtHig);
946   fOutput2Part->Add(f2HistMtLowQLamProReal);
947   // f2HistMtLowQALamALamReal = new TH1F("h1MtLowQALamALamReal"
948   //                              ,"m_{t}(#bar{#Lambda} #bar{#Lambda});m{t} [GeV];counts"
949   //                              ,nMtBins,mtLow,mtHig);
950   // fOutput2Part->Add(f2HistMtLowQALamALamReal);
951   f2HistMtLowQALamAProReal  = new TH1F("h1MtLowQALamAProReal"
952                                    ,"m_{t}(#bar{#Lambda} #bar{#Lambda});m{t} [GeV];counts"
953                                    ,nMtBins,mtLow,mtHig);
954   fOutput2Part->Add(f2HistMtLowQALamAProReal);
955
956   // Common qinv binning
957   Int_t nQinvBins = 400; // also for minv
958   Float_t QinvLow = 0.0;
959   Float_t QinvHig = 2.5;
960
961   // Sept'12 Use a THnSparse for (Bg)(A)Lam(A)Pro with 4 dimensions:
962   // qinv (ppri,ppri), mean dist (ppri,ppri), min dist(ppri,ppri)
963   // qinv (lam pro)
964   Int_t HnSpBins[4]={nQinvBins,nDistBins,nDistBins,nQinvBins};
965   Double_t HnSpMin[4]={QinvLow,distLow,distLow,QinvLow};
966   Double_t HnSpMax[4]={QinvHig,distHig,distHig,QinvHig};
967   LamProReal = new THnSparseF("HnSp4LamProReal","lamProRealQppMeanMinQlamp"
968                               ,4,HnSpBins,HnSpMin,HnSpMax);
969   fOutput2Part->Add(LamProReal);
970   ALamAProReal = new THnSparseF("HnSp4ALamAProReal","alamAProRealQppMeanMinQlamp"
971                               ,4,HnSpBins,HnSpMin,HnSpMax);
972   fOutput2Part->Add(ALamAProReal);
973
974   // Qinv: real events
975   // Since March 4th 2012 do corr. fcts vs distances
976   // f3HistLamLamQinvReal = new TH3F ("h3LamLamQinvReal", "Qinv LamLam;q_{inv} [GeV/c];min(d) p [cm];min(d) #pi [cm]"
977   //                               ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
978   // fOutput2Part->Add(f3HistLamLamQinvReal);
979   // f3HistALamALamQinvReal = new TH3F ("h3ALamALamQinvReal", "Qinv ALamALam;q_{inv} [GeV/c];min(d) #bar{p} [cm];min (d) #pi [cm]"
980   //                                 ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
981   // fOutput2Part->Add(f3HistALamALamQinvReal);
982   // // minv (H2 dibaryon??)
983   // f3HistLamLamMinvReal = new TH3F ("h3LamLamMinvReal", "Minv LamLam;q_{inv} [GeV/c];min(d) p [cm];min(d) #pi [cm]"
984   //                               ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
985   // fOutput2Part->Add(f3HistLamLamMinvReal);
986   // f3HistLamProMinvReal = new TH3F ("h3LamProMinvReal", "Minv LamPro;q_{inv} [GeV/c];<d> p [cm];min(d) p [cm]"
987   //                               ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
988   // fOutput2Part->Add(f3HistLamProMinvReal);
989   // f3HistALamALamMinvReal = new TH3F ("h3ALamALamMinvReal", "Minv ALamALam;q_{inv} [GeV/c];min(d) #bar{p} [cm];min(d) #pi [cm]"
990   //                                 ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
991   // fOutput2Part->Add(f3HistALamALamMinvReal);
992   // f3HistALamAProMinvReal = new TH3F ("h3ALamAProMinvReal", "Minv ALamAPro;q_{inv} [GeV/c];<d> #bar{p} [cm];min(d) #bar{p} [cm]"
993   //                                 ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
994   // fOutput2Part->Add(f3HistALamAProMinvReal);
995   
996   
997   // // Two-track resolution: mixed events
998   // f2HistLamLamMeanMinDistProMixed = new TH2F ("h2LamLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
999   // fOutput2Part->Add(f2HistLamLamMeanMinDistProMixed);
1000   // f2HistLamLamMeanMinDistPioMixed = new TH2F ("h2LamLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1001   // fOutput2Part->Add(f2HistLamLamMeanMinDistPioMixed);
1002   // f2HistLamProMeanMinDistProMixed = new TH2F ("h2LamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1003   // fOutput2Part->Add(f2HistLamProMeanMinDistProMixed);
1004   // f2HistALamALamMeanMinDistAProMixed = new TH2F ("h2ALamALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1005   // fOutput2Part->Add(f2HistALamALamMeanMinDistAProMixed);
1006   // f2HistALamALamMeanMinDistPioMixed = new TH2F ("h2ALamALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1007   // fOutput2Part->Add(f2HistALamALamMeanMinDistPioMixed);
1008   // f2HistALamAProMeanMinDistAProMixed = new TH2F ("h2ALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1009   // fOutput2Part->Add(f2HistALamAProMeanMinDistAProMixed);
1010
1011   // f2HistSftLamLamMeanMinDistProMixed = new TH2F ("h2SftLamLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1012   // fOutput2Part->Add(f2HistSftLamLamMeanMinDistProMixed);
1013   // f2HistSftLamLamMeanMinDistPioMixed = new TH2F ("h2SftLamLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1014   // fOutput2Part->Add(f2HistSftLamLamMeanMinDistPioMixed);
1015   // f2HistSftLamProMeanMinDistProMixed = new TH2F ("h2SftLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1016   // fOutput2Part->Add(f2HistSftLamProMeanMinDistProMixed);
1017   // f2HistSftALamALamMeanMinDistAProMixed = new TH2F ("h2SftALamALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1018   // fOutput2Part->Add(f2HistSftALamALamMeanMinDistAProMixed);
1019   // f2HistSftALamALamMeanMinDistPioMixed = new TH2F ("h2SftALamALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1020   // fOutput2Part->Add(f2HistSftALamALamMeanMinDistPioMixed);
1021   // f2HistSftALamAProMeanMinDistAProMixed = new TH2F ("h2SftALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1022   // fOutput2Part->Add(f2HistSftALamAProMeanMinDistAProMixed);
1023
1024   // f2HistSftIrocLamLamMeanMinDistProMixed = new TH2F ("h2SftIrocLamLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1025   // fOutput2Part->Add(f2HistSftIrocLamLamMeanMinDistProMixed);
1026   // f2HistSftIrocLamLamMeanMinDistPioMixed = new TH2F ("h2SftIrocLamLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1027   // fOutput2Part->Add(f2HistSftIrocLamLamMeanMinDistPioMixed);
1028   // f2HistSftIrocLamProMeanMinDistProMixed = new TH2F ("h2SftIrocLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1029   // fOutput2Part->Add(f2HistSftIrocLamProMeanMinDistProMixed);
1030   // f2HistSftIrocALamALamMeanMinDistAProMixed = new TH2F ("h2SftIrocALamALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1031   // fOutput2Part->Add(f2HistSftIrocALamALamMeanMinDistAProMixed);
1032   // f2HistSftIrocALamALamMeanMinDistPioMixed = new TH2F ("h2SftIrocALamALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1033   // fOutput2Part->Add(f2HistSftIrocALamALamMeanMinDistPioMixed);
1034   // f2HistSftIrocALamAProMeanMinDistAProMixed = new TH2F ("h2SftIrocALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1035   // fOutput2Part->Add(f2HistSftIrocALamAProMeanMinDistAProMixed);
1036
1037   // f2HistSftOrocLamLamMeanMinDistProMixed = new TH2F ("h2SftOrocLamLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1038   // fOutput2Part->Add(f2HistSftOrocLamLamMeanMinDistProMixed);
1039   // f2HistSftOrocLamLamMeanMinDistPioMixed = new TH2F ("h2SftOrocLamLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1040   // fOutput2Part->Add(f2HistSftOrocLamLamMeanMinDistPioMixed);
1041   // f2HistSftOrocLamProMeanMinDistProMixed = new TH2F ("h2SftOrocLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1042   // fOutput2Part->Add(f2HistSftOrocLamProMeanMinDistProMixed);
1043   // f2HistSftOrocALamALamMeanMinDistAProMixed = new TH2F ("h2SftOrocALamALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1044   // fOutput2Part->Add(f2HistSftOrocALamALamMeanMinDistAProMixed);
1045   // f2HistSftOrocALamALamMeanMinDistPioMixed = new TH2F ("h2SftOrocALamALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1046   // fOutput2Part->Add(f2HistSftOrocALamALamMeanMinDistPioMixed);
1047   // f2HistSftOrocALamAProMeanMinDistAProMixed = new TH2F ("h2SftOrocALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1048   // fOutput2Part->Add(f2HistSftOrocALamAProMeanMinDistAProMixed);
1049
1050   // Sept'12 Use a THnSparse for (Bg)(A)Lam(A)Pro with 4 dimensions:
1051   // qinv (ppri,ppri), mean dist (ppri,ppri), min dist(ppri,ppri)
1052   // qinv (lam pro)
1053   LamProMixed = new THnSparseF("HnSp4LamProMixed","lamProMixedQppMeanMinQlamp"
1054                               ,4,HnSpBins,HnSpMin,HnSpMax);
1055   fOutput2Part->Add(LamProMixed);
1056   ALamAProMixed = new THnSparseF("HnSp4ALamAProMixed","alamAProMixedQppMeanMinQlamp"
1057                               ,4,HnSpBins,HnSpMin,HnSpMax);
1058   fOutput2Part->Add(ALamAProMixed);
1059
1060   // // Qinv: mixed events
1061   // f3HistLamLamQinvMixed = new TH3F ("h3LamLamQinvMixed", "Qinv LamLam;q_{inv} [GeV/c];min(d) p [cm];min(d) #pi [cm]"
1062   //                                ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
1063   // fOutput2Part->Add(f3HistLamLamQinvMixed);
1064   // f3HistALamALamQinvMixed = new TH3F ("h3ALamALamQinvMixed", "Qinv ALamALam;q_{inv} [GeV/c];min(d) #bar{p} [cm];min(d) #pi [cm]"
1065   //                                  ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
1066   // fOutput2Part->Add(f3HistALamALamQinvMixed);
1067   // // minv (H2 di-baryon??)
1068   // f3HistLamLamMinvMixed = new TH3F ("h3LamLamMinvMixed", "Minv LamLam;q_{inv} [GeV/c];min(d) p [cm];min(d) #pi [cm]"
1069   //                                ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
1070   // fOutput2Part->Add(f3HistLamLamMinvMixed);
1071   // f3HistLamProMinvMixed = new TH3F ("h3LamProMinvMixed", "Minv LamPro;q_{inv} [GeV/c];<d> p [cm];min(d) p [cm]"
1072   //                                ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
1073   // fOutput2Part->Add(f3HistLamProMinvMixed);
1074   // f3HistALamALamMinvMixed = new TH3F ("h3ALamALamMinvMixed", "Minv ALamALam;q_{inv} [GeV/c];min(d) #bar{p} [cm];min(d) #pi [cm]"
1075   //                                  ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
1076   // fOutput2Part->Add(f3HistALamALamMinvMixed);
1077   // f3HistALamAProMinvMixed = new TH3F ("h3ALamAProMinvMixed", "Minv ALamAPro;q_{inv} [GeV/c];<d> #bar{p} [cm];min(d) #bar{p} [cm]"
1078   //                                  ,nQinvBins,2.0,3.0,100,0.,10.,100,0.,10.);
1079   // fOutput2Part->Add(f3HistALamAProMinvMixed);
1080
1081   // Same for Background (anti-)lambdas
1082
1083   // // Two-track resolution: real events
1084   // f2HistBgLamBgLamMeanMinDistProReal = new TH2F ("h2BgLamBgLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1085   // fOutput2Part->Add(f2HistBgLamBgLamMeanMinDistProReal);
1086   // f2HistBgLamBgLamMeanMinDistPioReal = new TH2F ("h2BgLamBgLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1087   // fOutput2Part->Add(f2HistBgLamBgLamMeanMinDistPioReal);
1088   // f2HistBgLamProMeanMinDistProReal = new TH2F ("h2BgLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1089   // fOutput2Part->Add(f2HistBgLamProMeanMinDistProReal);
1090   // f2HistBgALamBgALamMeanMinDistAProReal = new TH2F ("h2BgALamBgALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1091   // fOutput2Part->Add(f2HistBgALamBgALamMeanMinDistAProReal);
1092   // f2HistBgALamBgALamMeanMinDistPioReal = new TH2F ("h2BgALamBgALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1093   // fOutput2Part->Add(f2HistBgALamBgALamMeanMinDistPioReal);
1094   // f2HistBgALamAProMeanMinDistAProReal = new TH2F ("h2BgALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1095   // fOutput2Part->Add(f2HistBgALamAProMeanMinDistAProReal);
1096
1097   // f2HistSftBgLamBgLamMeanMinDistProReal = new TH2F ("h2SftBgLamBgLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1098   // fOutput2Part->Add(f2HistSftBgLamBgLamMeanMinDistProReal);
1099   // f2HistSftBgLamBgLamMeanMinDistPioReal = new TH2F ("h2SftBgLamBgLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1100   // fOutput2Part->Add(f2HistSftBgLamBgLamMeanMinDistPioReal);
1101   // f2HistSftBgLamProMeanMinDistProReal = new TH2F ("h2SftBgLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1102   // fOutput2Part->Add(f2HistSftBgLamProMeanMinDistProReal);
1103   // f2HistSftBgALamBgALamMeanMinDistAProReal = new TH2F ("h2SftBgALamBgALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1104   // fOutput2Part->Add(f2HistSftBgALamBgALamMeanMinDistAProReal);
1105   // f2HistSftBgALamBgALamMeanMinDistPioReal = new TH2F ("h2SftBgALamBgALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1106   // fOutput2Part->Add(f2HistSftBgALamBgALamMeanMinDistPioReal);
1107   // f2HistSftBgALamAProMeanMinDistAProReal = new TH2F ("h2SftBgALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1108   // fOutput2Part->Add(f2HistSftBgALamAProMeanMinDistAProReal);
1109
1110   // f2HistSftIrocBgLamBgLamMeanMinDistProReal = new TH2F ("h2SftIrocBgLamBgLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1111   // fOutput2Part->Add(f2HistSftIrocBgLamBgLamMeanMinDistProReal);
1112   // f2HistSftIrocBgLamBgLamMeanMinDistPioReal = new TH2F ("h2SftIrocBgLamBgLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1113   // fOutput2Part->Add(f2HistSftIrocBgLamBgLamMeanMinDistPioReal);
1114   // f2HistSftIrocBgLamProMeanMinDistProReal = new TH2F ("h2SftIrocBgLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1115   // fOutput2Part->Add(f2HistSftIrocBgLamProMeanMinDistProReal);
1116   // f2HistSftIrocBgALamBgALamMeanMinDistAProReal = new TH2F ("h2SftIrocBgALamBgALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1117   // fOutput2Part->Add(f2HistSftIrocBgALamBgALamMeanMinDistAProReal);
1118   // f2HistSftIrocBgALamBgALamMeanMinDistPioReal = new TH2F ("h2SftIrocBgALamBgALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1119   // fOutput2Part->Add(f2HistSftIrocBgALamBgALamMeanMinDistPioReal);
1120   // f2HistSftIrocBgALamAProMeanMinDistAProReal = new TH2F ("h2SftIrocBgALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1121   // fOutput2Part->Add(f2HistSftIrocBgALamAProMeanMinDistAProReal);
1122
1123   // f2HistSftOrocBgLamBgLamMeanMinDistProReal = new TH2F ("h2SftOrocBgLamBgLamMeanMinDistProReal","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1124   // fOutput2Part->Add(f2HistSftOrocBgLamBgLamMeanMinDistProReal);
1125   // f2HistSftOrocBgLamBgLamMeanMinDistPioReal = new TH2F ("h2SftOrocBgLamBgLamMeanMinDistPioReal","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1126   // fOutput2Part->Add(f2HistSftOrocBgLamBgLamMeanMinDistPioReal);
1127   // f2HistSftOrocBgLamProMeanMinDistProReal = new TH2F ("h2SftOrocBgLamProMeanMinDistProReal","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1128   // fOutput2Part->Add(f2HistSftOrocBgLamProMeanMinDistProReal);
1129   // f2HistSftOrocBgALamBgALamMeanMinDistAProReal = new TH2F ("h2SftOrocBgALamBgALamMeanMinDistAProReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1130   // fOutput2Part->Add(f2HistSftOrocBgALamBgALamMeanMinDistAProReal);
1131   // f2HistSftOrocBgALamBgALamMeanMinDistPioReal = new TH2F ("h2SftOrocBgALamBgALamMeanMinDistPioReal","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1132   // fOutput2Part->Add(f2HistSftOrocBgALamBgALamMeanMinDistPioReal);
1133   // f2HistSftOrocBgALamAProMeanMinDistAProReal = new TH2F ("h2SftOrocBgALamAProMeanMinDistProReal","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1134   // fOutput2Part->Add(f2HistSftOrocBgALamAProMeanMinDistAProReal);
1135
1136   // Sept'12 Use a THnSparse for (Bg)(A)Lam(A)Pro with 4 dimensions:
1137   // qinv (ppri,ppri), mean dist (ppri,ppri), min dist(ppri,ppri)
1138   // qinv (lam pro)
1139   BgLamProReal = new THnSparseF("HnSp4BgLamProReal","lamProRealQppMeanMinQlamp"
1140                               ,4,HnSpBins,HnSpMin,HnSpMax);
1141   fOutput2Part->Add(BgLamProReal);
1142   BgALamAProReal = new THnSparseF("HnSp4BgALamAProReal","alamAProRealQppMeanMinQlamp"
1143                               ,4,HnSpBins,HnSpMin,HnSpMax);
1144   fOutput2Part->Add(BgALamAProReal);
1145
1146   // Qinv: real events
1147   // f3HistBgLamBgLamQinvReal = new TH3F ("h3BgLamBgLamQinvReal", "Qinv BgLamBgLam;q_{inv} [GeV/c];min(d) p [cm];min(d) #pi [cm]"
1148   //                                   ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
1149   // fOutput2Part->Add(f3HistBgLamBgLamQinvReal);
1150   // f3HistBgALamBgALamQinvReal = new TH3F ("h3BgALamBgALamQinvReal", "Qinv BgALamBgALam;q_{inv} [GeV/c];min(d) #bar{p} [cm];min(d) #pi [cm]"
1151   //                                     ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
1152   // fOutput2Part->Add(f3HistBgALamBgALamQinvReal);
1153   
1154   // // Two-track resolution: mixed events
1155   // f2HistBgLamBgLamMeanMinDistProMixed = new TH2F ("h2BgLamBgLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1156   // fOutput2Part->Add(f2HistBgLamBgLamMeanMinDistProMixed);
1157   // f2HistBgLamBgLamMeanMinDistPioMixed = new TH2F ("h2BgLamBgLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1158   // fOutput2Part->Add(f2HistBgLamBgLamMeanMinDistPioMixed);
1159   // f2HistBgLamProMeanMinDistProMixed = new TH2F ("h2BgLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1160   // fOutput2Part->Add(f2HistBgLamProMeanMinDistProMixed);
1161   // f2HistBgALamBgALamMeanMinDistAProMixed = new TH2F ("h2BgALamBgALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1162   // fOutput2Part->Add(f2HistBgALamBgALamMeanMinDistAProMixed);
1163   // f2HistBgALamBgALamMeanMinDistPioMixed = new TH2F ("h2BgALamBgALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1164   // fOutput2Part->Add(f2HistBgALamBgALamMeanMinDistPioMixed);
1165   // f2HistBgALamAProMeanMinDistAProMixed = new TH2F ("h2BgALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1166   // fOutput2Part->Add(f2HistBgALamAProMeanMinDistAProMixed);
1167
1168   // f2HistSftBgLamBgLamMeanMinDistProMixed = new TH2F ("h2SftBgLamBgLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1169   // fOutput2Part->Add(f2HistSftBgLamBgLamMeanMinDistProMixed);
1170   // f2HistSftBgLamBgLamMeanMinDistPioMixed = new TH2F ("h2SftBgLamBgLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1171   // fOutput2Part->Add(f2HistSftBgLamBgLamMeanMinDistPioMixed);
1172   // f2HistSftBgLamProMeanMinDistProMixed = new TH2F ("h2SftBgLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1173   // fOutput2Part->Add(f2HistSftBgLamProMeanMinDistProMixed);
1174   // f2HistSftBgALamBgALamMeanMinDistAProMixed = new TH2F ("h2SftBgALamBgALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1175   // fOutput2Part->Add(f2HistSftBgALamBgALamMeanMinDistAProMixed);
1176   // f2HistSftBgALamBgALamMeanMinDistPioMixed = new TH2F ("h2SftBgALamBgALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1177   // fOutput2Part->Add(f2HistSftBgALamBgALamMeanMinDistPioMixed);
1178   // f2HistSftBgALamAProMeanMinDistAProMixed = new TH2F ("h2SftBgALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1179   // fOutput2Part->Add(f2HistSftBgALamAProMeanMinDistAProMixed);
1180
1181   // f2HistSftIrocBgLamBgLamMeanMinDistProMixed = new TH2F ("h2SftIrocBgLamBgLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1182   // fOutput2Part->Add(f2HistSftIrocBgLamBgLamMeanMinDistProMixed);
1183   // f2HistSftIrocBgLamBgLamMeanMinDistPioMixed = new TH2F ("h2SftIrocBgLamBgLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1184   // fOutput2Part->Add(f2HistSftIrocBgLamBgLamMeanMinDistPioMixed);
1185   // f2HistSftIrocBgLamProMeanMinDistProMixed = new TH2F ("h2SftIrocBgLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1186   // fOutput2Part->Add(f2HistSftIrocBgLamProMeanMinDistProMixed);
1187   // f2HistSftIrocBgALamBgALamMeanMinDistAProMixed = new TH2F ("h2SftIrocBgALamBgALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1188   // fOutput2Part->Add(f2HistSftIrocBgALamBgALamMeanMinDistAProMixed);
1189   // f2HistSftIrocBgALamBgALamMeanMinDistPioMixed = new TH2F ("h2SftIrocBgALamBgALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1190   // fOutput2Part->Add(f2HistSftIrocBgALamBgALamMeanMinDistPioMixed);
1191   // f2HistSftIrocBgALamAProMeanMinDistAProMixed = new TH2F ("h2SftIrocBgALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1192   // fOutput2Part->Add(f2HistSftIrocBgALamAProMeanMinDistAProMixed);
1193
1194   // f2HistSftOrocBgLamBgLamMeanMinDistProMixed = new TH2F ("h2SftOrocBgLamBgLamMeanMinDistProMixed","#Lambda#Lambda Mean vs min dist of decay protons;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1195   // fOutput2Part->Add(f2HistSftOrocBgLamBgLamMeanMinDistProMixed);
1196   // f2HistSftOrocBgLamBgLamMeanMinDistPioMixed = new TH2F ("h2SftOrocBgLamBgLamMeanMinDistPioMixed","#Lambda#Lambda Mean vs min dist of decay pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1197   // fOutput2Part->Add(f2HistSftOrocBgLamBgLamMeanMinDistPioMixed);
1198   // f2HistSftOrocBgLamProMeanMinDistProMixed = new TH2F ("h2SftOrocBgLamProMeanMinDistProMixed","p#Lambda Mean vs min dist of pri. p - dec p;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1199   // fOutput2Part->Add(f2HistSftOrocBgLamProMeanMinDistProMixed);
1200   // f2HistSftOrocBgALamBgALamMeanMinDistAProMixed = new TH2F ("h2SftOrocBgALamBgALamMeanMinDistAProMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist decay #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig); 
1201   // fOutput2Part->Add(f2HistSftOrocBgALamBgALamMeanMinDistAProMixed);
1202   // f2HistSftOrocBgALamBgALamMeanMinDistPioMixed = new TH2F ("h2SftOrocBgALamBgALamMeanMinDistPioMixed","#bar{#Lambda}#bar{#Lambda} Mean vs min dist of dec. pions;mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1203   // fOutput2Part->Add(f2HistSftOrocBgALamBgALamMeanMinDistPioMixed);
1204   // f2HistSftOrocBgALamAProMeanMinDistAProMixed = new TH2F ("h2SftOrocBgALamAProMeanMinDistProMixed","#bar{p}#bar{#Lambda} Mean vs min dist of pri. #bar{p} - dec #bar{p};mean dist. [cm];min dist [cm]",nDistBins,distLow,distHig,nDistBins,distLow,distHig);
1205   // fOutput2Part->Add(f2HistSftOrocBgALamAProMeanMinDistAProMixed);
1206
1207
1208   // Sept'12 Use a THnSparse for (Bg)(A)Lam(A)Pro with 4 dimensions:
1209   // qinv (ppri,ppri), mean dist (ppri,ppri), min dist(ppri,ppri)
1210   // qinv (lam pro)
1211   BgLamProMixed = new THnSparseF("HnSp4BgLamProMixed","lamProMixedQppMeanMinQlamp"
1212                               ,4,HnSpBins,HnSpMin,HnSpMax);
1213   fOutput2Part->Add(BgLamProMixed);
1214   BgALamAProMixed = new THnSparseF("HnSp4BgALamAProMixed","alamAProMixedQppMeanMinQlamp"
1215                               ,4,HnSpBins,HnSpMin,HnSpMax);
1216   fOutput2Part->Add(BgALamAProMixed);
1217
1218   // Qinv: mixed events
1219   // f3HistBgLamBgLamQinvMixed = new TH3F ("h3BgLamBgLamQinvMixed", "Qinv BgLamBgLam;q_{inv} [GeV/c];min(d) p [cm];min(d) #pi [cm]"
1220   //                                    ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
1221   // fOutput2Part->Add(f3HistBgLamBgLamQinvMixed);
1222   // f3HistBgALamBgALamQinvMixed = new TH3F ("h3BgALamBgALamQinvMixed", "Qinv BgALamBgALam;q_{inv} [GeV/c];min(d) #bar{p} [cm];min(d) #pi [cm]"
1223   //                                      ,nQinvBins,QinvLow,QinvHig,100,0.,10.,100,0.,10.);
1224   // fOutput2Part->Add(f3HistBgALamBgALamQinvMixed);
1225   
1226   // Post the data
1227   PostData(1, fOutputList);
1228   PostData(2, fOutputPrimaries);
1229   PostData(3, fOutput2Part);
1230
1231 }
1232
1233 //________________________________________________________________________
1234 void AliAnalysisTaskProtonLambda::UserExec(Option_t *) 
1235 {
1236   // Main loop
1237   // Called for each event
1238   
1239   // Fill a control histogram
1240   fHistGoodEvent->Fill(0.0);
1241
1242   // Get the event
1243   fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
1244   if (!fAOD) {
1245     printf("ERROR: fAOD not available\n");
1246     return;
1247   }
1248
1249   // Fill a control histogram
1250   fHistGoodEvent->Fill(1.0);  
1251
1252   // Get the centrality selection
1253   AliCentrality *centrality=NULL;
1254   centrality = fAOD->GetCentrality();
1255   if (!centrality) {
1256     printf ("ERROR: couldn't get the AliCentrality\n");
1257     return;
1258   }
1259   
1260   // Fill a control histogram
1261   fHistGoodEvent->Fill(2.0);  
1262
1263   // Check the fQuality flag of the centrality task
1264   // for details see
1265   // https://twiki.cern.ch/twiki/bin/viewauth/ALICE/CentStudies#How_we_determine_centrality
1266   if (centrality->GetQuality()){
1267     return;
1268   }
1269
1270   // Fill a control histogram
1271   fHistGoodEvent->Fill(3.0);  
1272
1273   // Analyze only 20% most central events using multiplicity in V0 detector (standard)
1274   Float_t centralityPercentile = centrality->GetCentralityPercentileUnchecked("V0M");
1275   if ( centralityPercentile > fkCentCut){
1276     return;
1277   }
1278
1279   // Fill a control histogram
1280   fHistGoodEvent->Fill(4.0);  
1281
1282   // Primary vertex, GetPrimaryVertex() returns the "best" reconstructed vertex
1283   fPrimaryVtx = fAOD->GetPrimaryVertex();
1284   if (!fPrimaryVtx){
1285     printf ("ERROR: no primary vertex\n");
1286     return;
1287   }
1288
1289   // Fill a control histogram
1290   fHistGoodEvent->Fill(5.0);  
1291   fPrimaryVtx->GetXYZ(fPrimaryVtxPosition);
1292   // fHistPrimaryVertexPosXY->Fill(fPrimaryVtxPosition[0],fPrimaryVtxPosition[1]);
1293   // fHistPrimaryVertexPosZ->Fill(fPrimaryVtxPosition[2]);
1294   
1295   // Zvertex cut, probably done anyhow in centrality task
1296   if (TMath::Abs(fPrimaryVtxPosition[2]) > fkAbsZvertexCut)
1297     return;
1298   
1299   // Fill a control histogram
1300   fHistGoodEvent->Fill(6.0);
1301
1302   // Multiplicity
1303   if (!(fAOD->GetNumberOfTracks())) {
1304     return;
1305   }
1306   
1307   // Fill a control histogram
1308   fHistGoodEvent->Fill(7.0);
1309
1310   // fHistTrackMultiplicity->Fill(fAOD->GetNumberOfTracks());
1311
1312   // Set up the event buffer to store this event
1313   fFemtoBuffer->ShiftAndAdd(fAOD);
1314
1315   // // Debugging: print number of stored tracks in the event
1316   // for(UChar_t i=0;i<fFemtoBuffer->GetMixBuffSize();i++)
1317   //   printf("iMix: %u, NPro %u, NAPro %u, NLam %u, NALam %u"
1318   //       "NBgLam %u, NBgALam %u\n"
1319   //       ,i
1320   //       ,fFemtoBuffer->GetEvt(i)->GetNPro()
1321   //       ,fFemtoBuffer->GetEvt(i)->GetNAPro()
1322   //       ,fFemtoBuffer->GetEvt(i)->GetNLam()
1323   //       ,fFemtoBuffer->GetEvt(i)->GetNALam()
1324   //       ,fFemtoBuffer->GetEvt(i)->GetNBgLam()
1325   //       ,fFemtoBuffer->GetEvt(i)->GetNBgALam()
1326   //       );
1327   // printf("\n");
1328
1329
1330   // Reset the reference array to the global tracks..
1331   ResetGlobalTrackReference();
1332   // ..and set it
1333   AliAODTrack *track=NULL;
1334   for (Int_t iTrack=0;iTrack<fAOD->GetNumberOfTracks();iTrack++){
1335     track = fAOD->GetTrack(iTrack);
1336     if (!track) continue;
1337     
1338     // Store the reference of the global tracks
1339     StoreGlobalTrackReference(track);
1340   }
1341   
1342   // V0 loop
1343   const Int_t nV0s = fAOD->GetNumberOfV0s();
1344   AliAODv0 *v0=NULL;
1345   AliAODTrack *pTrack=NULL;
1346   AliAODTrack *nTrack=NULL;
1347   for (Int_t iV0 = 0; iV0 < nV0s; iV0++) {
1348     v0 = fAOD->GetV0(iV0);
1349
1350     // Skip if V0 is not there
1351     if((!v0))
1352       continue;
1353
1354     // Check that the array fGTI isn't too small
1355     // for the track ids
1356     if(v0->GetPosID() >= fTrackBuffSize||
1357        v0->GetNegID() >= fTrackBuffSize)
1358       continue;
1359
1360     // This is AODs: find the track for given id:
1361     pTrack=fGTI[v0->GetPosID()];
1362     nTrack=fGTI[v0->GetNegID()];
1363         
1364     // Skip if one of the daughter is not there
1365     if ((!pTrack) || (!nTrack)) continue;
1366
1367     // Famous crossed rows / findable clusters cut,
1368     // rejects split tracks very well
1369     // (Don't do it for the V0s as we require 80 clusters 
1370     // and reject shared clusters)
1371     //    if( (!acceptTrack(pTrack)) || (!acceptTrack(nTrack)) )
1372     //      continue;
1373
1374     // Reject tracks with shared clusters
1375     if(!GoodTPCFitMapSharedMap(pTrack,nTrack))
1376       continue;
1377
1378     // Analysis done seperately for offline and on-the-fly
1379     if (!(v0->GetOnFlyStatus()))
1380       ProcessOffline(v0, pTrack, nTrack);
1381     else
1382       ProcessOnTheFly(v0, pTrack, nTrack);
1383
1384     // V0s get added to the mixed events in the 'Process..' fcts
1385     
1386   } // End of V0 loop
1387   
1388
1389   // Loop over primary tracks
1390   for (Int_t iTrack=0;iTrack<fAOD->GetNumberOfTracks();iTrack++){
1391     track = fAOD->GetTrack(iTrack);
1392     if (!track) continue;
1393     
1394     if(!track->TestFilterBit(128))
1395       continue;
1396
1397     // Famous crossed rows / findable clusters cut,
1398     // rejects split tracks very well
1399     if(!acceptTrack(track))
1400       continue;
1401
1402     // Reject tracks with shared clusters
1403     if(!GoodTPCFitMapSharedMap(track))
1404       continue;
1405
1406     // Check that the array fGTI isn't too small
1407     // for the track id
1408     if(-track->GetID()-1 >= fTrackBuffSize)
1409       continue;
1410
1411     // Without a corresponding global track it's useless
1412     if(!fGTI[-track->GetID()-1]){
1413       printf ("No global info! iTrack %d, ID %d\n",iTrack,track->GetID());
1414       continue;
1415     }
1416
1417     // Visualization of TPC dE/dx
1418     FillDedxHist(track);
1419
1420     // Depending on momentum choose pid method
1421     if (track->P() < 0.75){
1422        ProcessTPC(track);
1423     }
1424     else if (track->P() < 1.0){
1425        ProcessHybrid(track);
1426     }
1427     else if (track->P() < 3.25){
1428       ProcessTOF(track);
1429     }
1430
1431     
1432     // Tracks get added to the mixed events in the 'Process..' fcts
1433
1434   } // End of loop over primary tracks
1435
1436   // Track cuts do not allow for split tracks
1437
1438   //
1439   // TODO: Use Adam's shared cluster cut!
1440   //
1441
1442
1443   // Cleaning procedure for lambdas & lambdas, lambdas & protons,
1444   // anti-lambdas & anti-lambdas, anti-lambdas & protons + (anti-)lambda background
1445   CleaningProcedure();
1446
1447   // Process real events
1448   ProcessReal();
1449   ProcessRealBackground();
1450   
1451   // Process mixed events
1452   ProcessMixed();
1453   ProcessMixedBackground();
1454
1455   // Post output data.
1456   PostData(1, fOutputList);
1457   PostData(2, fOutputPrimaries);
1458   PostData(3, fOutput2Part);
1459
1460 }
1461 //________________________________________________________________________
1462 void AliAnalysisTaskProtonLambda::ProcessOffline(AliAODv0 *v0, AliAODTrack *pTrack, AliAODTrack *nTrack) 
1463 {
1464
1465   // For clarity in code: Fill some hists with on-the-fly status
1466   //  const Float_t kOnTheFlyStat = 0.0;
1467
1468   // All cuts are checked with invariant mass histograms
1469   //  v0->ChangeMassHypothesis(3122);
1470   Float_t minvLam = v0->MassLambda();
1471   //  v0->ChangeMassHypothesis(-3122);
1472   Float_t minvALam = v0->MassAntiLambda();
1473   // Cosine as local variable as this is some computation
1474   const Float_t lCosPoint = v0->CosPointingAngle(fPrimaryVtxPosition);
1475
1476   // Also calculate a V0 momentum with TPC only daughters
1477   //  Double_t TPConlyV0Mom[3], TPConlyV0MinvLam=0, TPConlyV0MinvALam=0;
1478   //  getTPConlyV0Info(pTrack, nTrack,
1479   //               TPConlyV0Mom, TPConlyV0MinvLam, TPConlyV0MinvALam);
1480
1481   // Fill a minv hist w/o any cuts. Select background from the sideband
1482   fHistSideBandOffLam->Fill(minvLam);
1483   fHistSideBandOffALam->Fill(minvALam);
1484   // Fill the event buffer w/ background
1485   if (!fkUseOnTheFly){
1486     if ( TMath::Abs(minvLam - fkLamMass) > 0.015 &&
1487          TMath::Abs(minvLam - fkLamMass) < 0.035 ){
1488       fFemtoBuffer->GetEvt(0)->AddBgLam(v0, pTrack, nTrack);
1489     }
1490     if ( TMath::Abs(minvALam - fkLamMass) > 0.015 &&
1491          TMath::Abs(minvALam - fkLamMass) < 0.035 ){
1492       fFemtoBuffer->GetEvt(0)->AddBgALam(v0, pTrack, nTrack);
1493     }
1494   }
1495
1496   // Control histogram: fill all v0s
1497   // fHistGoodV0->Fill(0.0,kOnTheFlyStat);
1498   // fHistGoodV0->Fill(1.0,kOnTheFlyStat);
1499
1500   // Require 80 TPC clusters for both pos and neg daughter
1501   fHistTPCNclsPosOffLam->Fill(pTrack->GetTPCNcls(),minvLam);
1502   fHistTPCNclsNegOffLam->Fill(nTrack->GetTPCNcls(),minvLam);
1503   fHistTPCNclsPosOffALam->Fill(pTrack->GetTPCNcls(),minvALam);
1504   fHistTPCNclsNegOffALam->Fill(nTrack->GetTPCNcls(),minvALam);
1505
1506   if ( ( (pTrack->GetTPCNcls()) < 80 ) || ( (nTrack->GetTPCNcls()) < 80 ) ) 
1507     return;
1508   //  fHistGoodV0->Fill(2.0,kOnTheFlyStat);
1509
1510   // Require a maximum dca of the daughters of 0.6cm
1511   // fHistDcaV0DaughtersOffLam->Fill(v0->DcaV0Daughters(),minvLam);
1512   // fHistDcaV0DaughtersOffALam->Fill(v0->DcaV0Daughters(),minvALam);
1513   // fHistDcaV0Daughters->Fill(v0->DcaV0Daughters(),kOnTheFlyStat);
1514   if (v0->DcaV0Daughters() > 0.6)
1515     return;
1516   //  fHistGoodV0->Fill(3.0,kOnTheFlyStat);
1517   
1518   // Force TPC PID to be present
1519   if (!(pTrack->GetStatus() & AliVTrack::kTPCpid) ||
1520       !(nTrack->GetStatus() & AliVTrack::kTPCpid))
1521     return;
1522   //  fHistGoodV0->Fill(4.0,kOnTheFlyStat);
1523
1524   // Visualize TPC signal before performing selection
1525   // fHistPosTpcBeforeCut->Fill(pTrack->P(),pTrack->GetTPCsignal());
1526   // fHistNegTpcBeforeCut->Fill(nTrack->P(),nTrack->GetTPCsignal());
1527   // The Nsigma distribution for TPC dE/dx
1528   // fHistPosNsigmaTpcOffLam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kProton)),minvLam);
1529   // fHistPosNsigmaTpcOffALam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kPion)),minvALam);
1530   // fHistNegNsigmaTpcOffLam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kPion)),minvLam);
1531   // fHistNegNsigmaTpcOffALam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kProton)),minvALam);
1532   // Perform cut on TPC dE/dx
1533   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kProton)) > 3.4)
1534     minvLam=0.0;
1535   // else 
1536   //   fHistPosTpcAfterCut->Fill(pTrack->P(),pTrack->GetTPCsignal());
1537   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kPion)) > 4.4)
1538     minvLam=0.0;
1539   // else
1540   //   fHistNegTpcAfterCut->Fill(nTrack->P(),nTrack->GetTPCsignal());
1541   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kPion)) > 4.2)
1542     minvALam=0.0;
1543   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kProton)) > 3.4)
1544     minvALam=0.0;
1545
1546   // Don't use a tof cut for pions
1547
1548   // Check whether to use a 5sigma tof cut or none for protons
1549   // if (pTrack->GetStatus() & AliVTrack::kTOFpid){
1550   //   if (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(pTrack, AliPID::kProton)) > 5.0)
1551   //     fHistUseTofOffLam->Fill(1.0,minvLam);
1552   //   else
1553   //     fHistUseTofOffLam->Fill(0.0,minvLam);
1554   // }
1555   // else
1556   //   fHistUseTofOffLam->Fill(0.0,minvLam);
1557   // Check whether to use a 5sigma tof cut or none for anti-protons
1558   // if (nTrack->GetStatus() & AliVTrack::kTOFpid){
1559   //   if (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(nTrack, AliPID::kProton)) > 5.0)
1560   //     fHistUseTofOffALam->Fill(1.0,minvALam);
1561   //   else
1562   //     fHistUseTofOffALam->Fill(0.0,minvALam);
1563   // }
1564   // else
1565   //   fHistUseTofOffALam->Fill(0.0,minvALam);
1566
1567   // Don't use a TOF cut for offline
1568   
1569   // Don't need to check for sign of pairs as this is always
1570   // correct for offline finder
1571
1572   // Don't need to check for TPC refit as it is required
1573   // by the offline finder itself
1574
1575   //
1576   // Require a minimum distance between daughters and primary vertex
1577   //
1578   // Fill histograms with the distributions before cutting
1579   // fHistDcaPosOffLam->Fill(v0->DcaPosToPrimVertex(),minvLam);
1580   // fHistDcaPosOffALam->Fill(v0->DcaPosToPrimVertex(),minvALam);
1581   // fHistDcaNegOffLam->Fill(v0->DcaNegToPrimVertex(),minvLam);
1582   // fHistDcaNegOffALam->Fill(v0->DcaNegToPrimVertex(),minvALam);
1583   
1584   // fHistDcaPosToPrimVertex->Fill(v0->DcaPosToPrimVertex(),kOnTheFlyStat);
1585   // fHistDcaNegToPrimVertex->Fill(v0->DcaNegToPrimVertex(),kOnTheFlyStat);
1586   // fHistDcaPosToPrimVertexZoom->Fill(v0->DcaPosToPrimVertex(),kOnTheFlyStat);
1587   // fHistDcaNegToPrimVertexZoom->Fill(v0->DcaNegToPrimVertex(),kOnTheFlyStat);
1588   
1589   // Do the cut
1590   if (v0->DcaPosToPrimVertex() < 0.1)
1591     minvLam=0.0;
1592   if (v0->DcaPosToPrimVertex() < 0.3)
1593     minvALam=0.0;
1594   if (v0->DcaNegToPrimVertex() < 0.1)
1595     minvALam=0.0;
1596   if (v0->DcaNegToPrimVertex() < 0.3)
1597     minvLam=0.0;
1598
1599   // Cosine of pointing angle. Computed at the beginning.
1600   // Fill historgrams before cutting
1601   // fHistCosPointLamOff->Fill(lCosPoint,minvLam);
1602   // fHistCosPointALamOff->Fill(lCosPoint,minvALam);
1603   // fHistCosPointLamZoomOff->Fill(lCosPoint,minvLam);
1604   // fHistCosPointALamZoomOff->Fill(lCosPoint,minvALam);
1605   
1606   // fHistCosPointAngle->Fill(lCosPoint,kOnTheFlyStat);
1607   // fHistCosPointAngleZoom->Fill(lCosPoint,kOnTheFlyStat);
1608   
1609   // Do the cut in cos (pointing angle) 
1610   // (note the difference 0.9996 for offline and 0.9999 for on-the-fly)
1611   if (lCosPoint < 0.9996)
1612     return;
1613   
1614   // fHistGoodV0->Fill(7.0,kOnTheFlyStat);
1615   
1616   // Fill some histograms with cut variables
1617   // fHistChi2->Fill(v0->Chi2V0(),kOnTheFlyStat);
1618   
1619   // Idea to cut on the radius
1620   // fHistRadiusV0->Fill(v0->RadiusV0(),kOnTheFlyStat);
1621   // fHistV0RadiusLamOff->Fill(v0->RadiusV0(),minvLam);
1622   // fHistV0RadiusALamOff->Fill(v0->RadiusV0(),minvALam);
1623
1624   // Idea to cut on the decay length
1625   // fHistDecayLengthV0->Fill(v0->DecayLengthV0(fPrimaryVtxPosition),kOnTheFlyStat);
1626   // fHistV0DecayLengthLamOff->Fill(v0->DecayLengthV0(fPrimaryVtxPosition),minvLam);
1627   // fHistV0DecayLengthALamOff->Fill(v0->DecayLengthV0(fPrimaryVtxPosition),minvALam);
1628  
1629   // Idea to cut on DCA of V0 and primay vertex
1630   // fHistDcaV0PriVertexLamOff->Fill(v0->DcaV0ToPrimVertex(),minvLam);
1631   // fHistDcaV0PriVertexALamOff->Fill(v0->DcaV0ToPrimVertex(),minvALam);
1632      
1633   // Fill some invariant mass distributions
1634   fHistMassLambdaOff->Fill(minvLam);
1635   fHistMassAntiLambdaOff->Fill(minvALam);
1636   // fHistPtVsMassLambdaOff->Fill(v0->Pt(),minvLam);
1637   // fHistPtVsMassAntiLambdaOff->Fill(v0->Pt(),minvALam);
1638
1639   // 3d histogram: rapidity, pt and mass
1640   fHistYPtMassLamOff->Fill(v0->Y(3122),v0->Pt(),minvLam);
1641   fHistYPtMassALamOff->Fill(v0->Y(-3122),v0->Pt(),minvALam);
1642   
1643   // Invariant mass cut lambda :: fill a y-pt hist
1644   // if ( TMath::Abs(minvLam - fkLamMass) < 0.01 ){
1645   //   fHistPtVsYLambdaOff->Fill(v0->Pt(),v0->Y(3122));
1646   // }
1647   // // Invariant mass cut anti-lambda :: fill a y-pt hist
1648   // if ( TMath::Abs(minvALam - fkLamMass) < 0.01 ){
1649   //   fHistPtVsYAntiLambdaOff->Fill(v0->Pt(),v0->Y(-3122));
1650   // }
1651
1652   // Fill the mixed events when offline V0 finder is used
1653   if (!fkUseOnTheFly){
1654     // Highest significance for minv +/- 4 MeV
1655     if ( TMath::Abs(minvLam - fkLamMass) < 0.004 ){
1656       fFemtoBuffer->GetEvt(0)->AddLam(v0, pTrack, nTrack);
1657     }
1658     if ( TMath::Abs(minvALam - fkLamMass) < 0.004 ){
1659       fFemtoBuffer->GetEvt(0)->AddALam(v0, pTrack, nTrack);
1660     }
1661   }
1662 }   
1663 //________________________________________________________________________
1664 void AliAnalysisTaskProtonLambda::ProcessOnTheFly(AliAODv0 *v0, AliAODTrack *pTrack, AliAODTrack *nTrack) 
1665 {
1666   // For clarity in code: Fill some hists with on-the-fly status
1667   //  const Float_t kOnTheFlyStat = 1.0;
1668
1669   // All cuts are checked with invariant mass histograms
1670   Float_t minvLam = v0->MassLambda();
1671   Float_t minvALam = v0->MassAntiLambda();
1672   const Float_t lCosPoint = v0->CosPointingAngle(fPrimaryVtxPosition);
1673
1674   // Control histogram: fill all v0s
1675   //  fHistGoodV0->Fill(0.0,kOnTheFlyStat);
1676   // Control hist: after require two daughter tracks
1677   //  fHistGoodV0->Fill(1.0,kOnTheFlyStat);
1678   // Check the right sign of the tracks (mainly on-the-fly)
1679   if (pTrack->Charge() > 0 && nTrack->Charge() < 0){
1680     // Correct assignment
1681     // fHistCorrectSigns->Fill(0.0,kOnTheFlyStat);
1682
1683     // fHistLikeSignOnLam->Fill(0.0,minvLam);
1684     // fHistLikeSignOnALam->Fill(0.0,minvALam);    
1685   }
1686   else if (pTrack->Charge() < 0 && nTrack->Charge() > 0){
1687     // Swapped sings
1688     //    fHistCorrectSigns->Fill(1.0,kOnTheFlyStat);
1689
1690     pTrack = fGTI[v0->GetNegID()];
1691     nTrack = fGTI[v0->GetPosID()];
1692     
1693
1694     // See http://savannah.cern.ch/bugs/?90749
1695     // For AODs it depends on with which root version 
1696     // the AODs got produced.
1697
1698     // See above: swapping mass assignment
1699     minvLam = v0->MassAntiLambda();
1700     minvALam = v0->MassLambda();
1701
1702     //    fHistLikeSignOnLam->Fill(1.0,minvLam);
1703     //    fHistLikeSignOnALam->Fill(1.0,minvALam);    
1704   }
1705   else {
1706     // Like sign pairs
1707     //    fHistCorrectSigns->Fill(2.0,kOnTheFlyStat);
1708     
1709     //    fHistLikeSignOnLam->Fill(2.0,minvLam);
1710     //    fHistLikeSignOnALam->Fill(2.0,minvALam);    
1711
1712     // Don't use like sign-pairs
1713     return;
1714   }
1715   //  fHistGoodV0->Fill(2.0,kOnTheFlyStat);
1716
1717   // V0 momentum
1718   Double_t V0Mom[3];
1719   v0->PxPyPz(V0Mom);
1720   // Also calculate a V0 momentum with TPC only daughters
1721   //  Double_t TPConlyV0Mom[3], TPConlyV0MinvLam=0, TPConlyV0MinvALam=0;
1722   //  getTPConlyV0Info(pTrack, nTrack,
1723   //               TPConlyV0Mom, TPConlyV0MinvLam, TPConlyV0MinvALam);
1724
1725    // Fill a minv hist w/o any cuts. Select background from the sideband
1726   fHistSideBandOnLam->Fill(minvLam);
1727   fHistSideBandOnALam->Fill(minvALam);
1728   // Fill the event buffer w/ background
1729   if (fkUseOnTheFly){
1730     // Select side band aka background lambdas
1731     if (TMath::Abs(minvLam - fkLamMass) > 0.015 &&
1732         TMath::Abs(minvLam - fkLamMass) < 0.035 ){
1733       
1734       fFemtoBuffer->GetEvt(0)->AddBgLam(v0, pTrack, nTrack);
1735       // Momentum difference of standard V0 / TPC only V0
1736       //      fHistMomDiffBgLam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1737       //                            V0Mom[1] - TPConlyV0Mom[1],
1738       //                            V0Mom[2] - TPConlyV0Mom[2]);
1739       // Same excluding V0s with daughters with SPD hits
1740       //      if( !(pTrack->HasPointOnITSLayer(0) || pTrack->HasPointOnITSLayer(1) ||
1741       //            nTrack->HasPointOnITSLayer(0) || nTrack->HasPointOnITSLayer(1) )){
1742       // No SPD hits
1743       //        fHistMomDiffWoSPDBgLam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1744       //                                   V0Mom[1] - TPConlyV0Mom[1],
1745       //                                   V0Mom[2] - TPConlyV0Mom[2]);
1746
1747       //    }
1748     } // End of background lambdas
1749     // Select side band aka background anti-lambdas
1750     if ( TMath::Abs(minvALam - fkLamMass) > 0.015 &&
1751          TMath::Abs(minvALam - fkLamMass) < 0.035 ){
1752
1753       fFemtoBuffer->GetEvt(0)->AddBgALam(v0, pTrack, nTrack);
1754       // Momentum difference of standard V0 / TPC only V0
1755       //      fHistMomDiffBgALam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1756       //                               V0Mom[1] - TPConlyV0Mom[1],
1757       //                               V0Mom[2] - TPConlyV0Mom[2]);
1758       // Same excluding V0s with daughters with SPD hits
1759       //      if( !(pTrack->HasPointOnITSLayer(0) || pTrack->HasPointOnITSLayer(1) ||
1760       //            nTrack->HasPointOnITSLayer(0) || nTrack->HasPointOnITSLayer(1) )){
1761       // No SPD hits
1762       //        fHistMomDiffWoSPDBgALam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1763       //                              V0Mom[1] - TPConlyV0Mom[1],
1764       //                                      V0Mom[2] - TPConlyV0Mom[2]);
1765       //      } // No SPD hits
1766     } // End of background anti-lambda
1767   } // End of if use on-the-fly finder
1768
1769   //
1770   // Require 80 TPC clusters for both daughters
1771   //
1772   // There's a lambda signal for 0-9 clusters of the proton 
1773   // as it's for 110-120?!
1774   // There was a bug in the finding of the global track, since 
1775   // fixing it, offline is fine (and the problem looks less 
1776   // severe for on-the-fly). Still there is a problem here. 
1777   // There are tracks with 0 clusters. This is not the case
1778   // for the offline finder. The speculation would be that 
1779   // 1-9 clusters are treated correctly also here, it's just
1780   // the 0 cluster tracks. Should be a filter issue: on-the-fly
1781   // finds a V0, stores the daughter but info doesn't get written.
1782   if(pTrack->GetTPCNcls()){
1783     // More than zero clusters
1784     fHistTPCNclsPosOnLam->Fill(pTrack->GetTPCNcls(),minvLam);
1785     fHistTPCNclsPosOnALam->Fill(pTrack->GetTPCNcls(),minvALam);
1786   }
1787   else {
1788     // Zero clusters, fill the underflow to distinguish
1789     fHistTPCNclsPosOnLam->Fill(-1,minvLam);
1790     fHistTPCNclsPosOnALam->Fill(-1,minvALam);
1791   }
1792   if(nTrack->GetTPCNcls()){
1793     // More than zero clusters
1794     fHistTPCNclsNegOnLam->Fill(nTrack->GetTPCNcls(),minvLam);
1795     fHistTPCNclsNegOnALam->Fill(nTrack->GetTPCNcls(),minvALam);
1796   }
1797   else {
1798     // Zero clusters, fill the underflow to distinguish
1799     fHistTPCNclsNegOnLam->Fill(-1,minvLam);
1800     fHistTPCNclsNegOnALam->Fill(-1,minvALam);
1801   }
1802   
1803   // Do the cut on the TPC clusters, 0 OR at least 80
1804   if ( ( pTrack->GetTPCNcls() < 80 && pTrack->GetTPCNcls() ) ||
1805        ( nTrack->GetTPCNcls() < 80 && nTrack->GetTPCNcls() ) ) 
1806     return;
1807   //  fHistGoodV0->Fill(3.0,kOnTheFlyStat);
1808
1809   // Require a maximum dca of the daughters of 0.2cm
1810   // fHistDcaV0DaughtersOnLam->Fill(v0->DcaV0Daughters(),minvLam);
1811   // fHistDcaV0DaughtersOnALam->Fill(v0->DcaV0Daughters(),minvALam);
1812   // fHistDcaV0Daughters->Fill(v0->DcaV0Daughters(),kOnTheFlyStat);
1813   if (v0->DcaV0Daughters() > 0.2)
1814     return;
1815   //  fHistGoodV0->Fill(4.0,kOnTheFlyStat);
1816   
1817   // Require cosine of pointing angle bigger than 0.9999
1818   // fHistCosPointAngle->Fill(lCosPoint,kOnTheFlyStat);
1819   // fHistCosPointAngleZoom->Fill(lCosPoint,kOnTheFlyStat);
1820   // fHistCosPointLamOn->Fill(lCosPoint,minvLam);
1821   // fHistCosPointALamOn->Fill(lCosPoint,minvALam);
1822   // fHistCosPointLamZoomOn->Fill(lCosPoint,minvLam);
1823   // fHistCosPointALamZoomOn->Fill(lCosPoint,minvALam);
1824   if (lCosPoint<0.9999)
1825     return;
1826   //  fHistGoodV0->Fill(5.0,kOnTheFlyStat);
1827   // Force TPC PID to be present
1828   if (!(pTrack->GetStatus() & AliVTrack::kTPCpid) ||
1829       !(nTrack->GetStatus() & AliVTrack::kTPCpid)) {
1830     // No TPC pid present for this track
1831     return;
1832   }
1833   //  fHistGoodV0->Fill(6.0,kOnTheFlyStat);
1834   // Visualize TPC signal before performing selection
1835   // fHistPosTpcBeforeCut->Fill(pTrack->P(),pTrack->GetTPCsignal());
1836   // fHistNegTpcBeforeCut->Fill(nTrack->P(),nTrack->GetTPCsignal());
1837   // // The Nsigma distribution for TPC dE/dx
1838   // fHistPosNsigmaTpcOnLam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kProton)),minvLam);
1839   // fHistPosNsigmaTpcOnALam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kPion)),minvALam);
1840   // fHistNegNsigmaTpcOnLam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kPion)),minvLam);
1841   // fHistNegNsigmaTpcOnALam->Fill(TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kProton)),minvALam);
1842
1843   // Perform cut on TPC dE/dx
1844   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kProton)) > 3.7)
1845     minvLam=0.0;
1846   // else 
1847   //   fHistPosTpcAfterCut->Fill(pTrack->P(),pTrack->GetTPCsignal());
1848   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kPion)) > 3.8)
1849     minvLam=0.0;
1850   // else
1851   //   fHistNegTpcAfterCut->Fill(nTrack->P(),nTrack->GetTPCsignal());
1852   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(pTrack, AliPID::kPion)) > 4.2)
1853     minvALam=0.0;
1854   if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC(nTrack, AliPID::kProton)) > 3.9)
1855     minvALam=0.0;
1856
1857   // Don't use a tof cut for pions
1858
1859   // Check whether to use a 5sigma tof cut or none for protons
1860   // if (pTrack->GetStatus() & AliVTrack::kTOFpid){
1861   //   if (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(pTrack, AliPID::kProton)) > 5.0)
1862   //     fHistUseTofOnLam->Fill(1.0,minvLam);
1863   //   else
1864   //     fHistUseTofOnLam->Fill(0.0,minvLam);
1865   // }
1866   // else
1867   //   fHistUseTofOnLam->Fill(0.0,minvLam);
1868   // // Check whether to use a 5sigma tof cut or none for anti-protons
1869   // if (nTrack->GetStatus() & AliVTrack::kTOFpid){
1870   //   if (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(nTrack, AliPID::kProton)) > 5.0)
1871   //     fHistUseTofOnALam->Fill(1.0,minvALam);
1872   //   else
1873   //     fHistUseTofOnALam->Fill(0.0,minvALam);
1874   // }
1875   // else
1876   //   fHistUseTofOnALam->Fill(0.0,minvALam);
1877
1878   // Reject (anti-)protons with more than 5sigma TOF
1879   if (nTrack->GetStatus() & AliVTrack::kTOFpid){
1880     if (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(nTrack, AliPID::kProton)) > 5.0)
1881       minvALam=0.0;
1882   }
1883   if (pTrack->GetStatus() & AliVTrack::kTOFpid){
1884     if (TMath::Abs(fPIDResponse->NumberOfSigmasTOF(pTrack, AliPID::kProton)) > 5.0)
1885       minvLam=0.0;
1886   }
1887   
1888   // Don't require TPC refit. You would kill nearly your whole signal  
1889
1890   // Distance between daughters and primary vertex
1891   // fHistDcaPosToPrimVertex->Fill(v0->DcaPosToPrimVertex(),kOnTheFlyStat);
1892   // fHistDcaNegToPrimVertex->Fill(v0->DcaNegToPrimVertex(),kOnTheFlyStat);
1893   // fHistDcaPosToPrimVertexZoom->Fill(v0->DcaPosToPrimVertex(),kOnTheFlyStat);
1894   // fHistDcaNegToPrimVertexZoom->Fill(v0->DcaNegToPrimVertex(),kOnTheFlyStat);
1895   // fHistDcaPosOnLam->Fill(v0->DcaPosToPrimVertex(),minvLam);
1896   // fHistDcaPosOnALam->Fill(v0->DcaPosToPrimVertex(),minvALam);
1897   // fHistDcaNegOnLam->Fill(v0->DcaNegToPrimVertex(),minvLam);
1898   // fHistDcaNegOnALam->Fill(v0->DcaNegToPrimVertex(),minvALam);
1899   // Require at least 0.02 cm distance from the primary vertex for the (anti-)protons
1900   if (v0->DcaPosToPrimVertex() < 0.02)
1901     minvLam=0.0;
1902   if (v0->DcaNegToPrimVertex() < 0.02)
1903     minvALam=0.0;
1904   // Require at least 0.05 cm distance from the primary vertex for the pions
1905   if (v0->DcaPosToPrimVertex() < 0.05)
1906     minvALam=0.0;
1907   if (v0->DcaNegToPrimVertex() < 0.05)
1908     minvLam=0.0;
1909   
1910   // Fill some histograms with cut variables
1911   //  fHistChi2->Fill(v0->Chi2V0(),kOnTheFlyStat);
1912   
1913   
1914   // Idea to cut on the radius
1915   // fHistRadiusV0->Fill(v0->RadiusV0(),kOnTheFlyStat);
1916   // fHistV0RadiusLamOn->Fill(v0->RadiusV0(),minvLam);
1917   // fHistV0RadiusALamOn->Fill(v0->RadiusV0(),minvALam);
1918   
1919   // Idea to cut on the decay length
1920   // fHistDecayLengthV0->Fill(v0->DecayLengthV0(fPrimaryVtxPosition),kOnTheFlyStat);
1921   // fHistV0DecayLengthLamOn->Fill(v0->DecayLengthV0(fPrimaryVtxPosition),minvLam);
1922   // fHistV0DecayLengthALamOn->Fill(v0->DecayLengthV0(fPrimaryVtxPosition),minvALam);
1923
1924   // Idea to cut on DCA of V0 and primay vertex
1925   // fHistDcaV0PriVertexLamOn->Fill(v0->DcaV0ToPrimVertex(),minvLam);
1926   // fHistDcaV0PriVertexALamOn->Fill(v0->DcaV0ToPrimVertex(),minvALam);
1927
1928   // TPC Chi2 / number of degrees of freedom
1929   // A cut on at least 80 clusters is already done before,
1930   // no concern to divide by zero
1931   // fHistChi2TPCPosLamOn->Fill(pTrack->Chi2perNDF(),minvLam);
1932   // fHistChi2TPCPosALamOn->Fill(pTrack->Chi2perNDF(),minvALam);
1933   // fHistChi2TPCNegLamOn->Fill(nTrack->Chi2perNDF(),minvLam);
1934   // fHistChi2TPCNegALamOn->Fill(nTrack->Chi2perNDF(),minvALam);
1935   // Don't cut like Chi2/ndf < 4! One might throw away the tracks
1936   // with Chi2/ndf roughly one as they are good primaries
1937
1938   // Fill some invariant mass distributions
1939   fHistMassLambdaOn->Fill(minvLam);
1940   fHistMassAntiLambdaOn->Fill(minvALam);
1941   // fHistPtVsMassLambdaOn->Fill(v0->Pt(),minvLam);
1942   // fHistPtVsMassAntiLambdaOn->Fill(v0->Pt(),minvALam);
1943
1944   // TPC only invariant mass distributions
1945   //  if(minvLam > .1){
1946     // Lambda is good
1947     //    fHistMinvTPConlyLamOn->Fill(TPConlyV0MinvLam);
1948   //  }
1949   //  if (minvALam > .1){
1950     // Anti-lambda is good
1951     //    fHistMinvTPConlyALamOn->Fill(TPConlyV0MinvALam);
1952   //  }
1953   
1954   // 3d histogram: rapidity, pt and mass
1955   fHistYPtMassLamOn->Fill(v0->Y(3122),v0->Pt(),minvLam);
1956   fHistYPtMassALamOn->Fill(v0->Y(-3122),v0->Pt(),minvALam);
1957   
1958   // // Invariant mass cut lambda :: fill a y-pt hists
1959   // if ( TMath::Abs(minvLam - fkLamMass) < 0.01 ){
1960   //   fHistPtVsYLambdaOn->Fill(v0->Pt(),v0->Y(3122));
1961   // }
1962   // // Invariant mass cut anti-lambda :: fill a y-pt hists
1963   // if ( TMath::Abs(minvALam - fkLamMass) < 0.01 ){
1964   //   fHistPtVsYAntiLambdaOn->Fill(v0->Pt(),v0->Y(-3122));
1965   // }
1966   
1967   // Fill the mixed events when on-the-fly V0 finder is used
1968   if (fkUseOnTheFly){
1969
1970     // Highest significance for minv +/- 4 MeV
1971     if ( TMath::Abs(minvLam - fkLamMass) < 0.004 ){
1972       fFemtoBuffer->GetEvt(0)->AddLam(v0, pTrack, nTrack);
1973       // Momentum difference of standard V0 / TPC only V0
1974       //      fHistMomDiffLam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1975       //                    V0Mom[1] - TPConlyV0Mom[1],
1976       //                    V0Mom[2] - TPConlyV0Mom[2]);
1977       // Same excluding V0s with daughters with SPD hits
1978       //      if( !(pTrack->HasPointOnITSLayer(0) || pTrack->HasPointOnITSLayer(1) ||
1979       //            nTrack->HasPointOnITSLayer(0) || nTrack->HasPointOnITSLayer(1) )){
1980       //        // No SPD hits
1981       //        fHistMomDiffWoSPDLam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1982       //                                   V0Mom[1] - TPConlyV0Mom[1],
1983       //                                   V0Mom[2] - TPConlyV0Mom[2]);
1984       // } // No SPD hits
1985     } // Good lambda
1986     if ( TMath::Abs(minvALam - fkLamMass) < 0.004 ) {
1987       fFemtoBuffer->GetEvt(0)->AddALam(v0, pTrack, nTrack);
1988       // Momentum difference of standard V0 / TPC only V0
1989       //      fHistMomDiffALam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1990       //                            V0Mom[1] - TPConlyV0Mom[1],
1991       //                            V0Mom[2] - TPConlyV0Mom[2]);
1992       // Same excluding V0s with daughters with SPD hits
1993       //      if( !(pTrack->HasPointOnITSLayer(0) || pTrack->HasPointOnITSLayer(1) ||
1994       //            nTrack->HasPointOnITSLayer(0) || nTrack->HasPointOnITSLayer(1) )){
1995         // No SPD hits
1996       //        fHistMomDiffWoSPDALam->Fill(V0Mom[0] - TPConlyV0Mom[0], 
1997       //                                   V0Mom[1] - TPConlyV0Mom[1],
1998       //                                   V0Mom[2] - TPConlyV0Mom[2]);
1999       //    } // No SPD hits
2000     } // Good anti-lambda
2001   } // Use on-the-fly finder for Femto analysis
2002 } // ProcessOnTheFly
2003 //________________________________________________________________________
2004 void AliAnalysisTaskProtonLambda::ProcessTOF(AliAODTrack* track) 
2005 {
2006   // Request the kTOFpid bit. There are tracks with kTOFout and wihthout kTOFpid,
2007   // but these tracks have a bad TOF signal.
2008   if(!((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTOFpid))
2009     return;
2010
2011   // TOF signal corrected for expected time and (if neccessary) for start time
2012   Float_t corrTOFsig = GetCorrectedTOFSignal(track);
2013
2014   // Distinguish between charges
2015   if (track->Charge() > 0){
2016     // Simple Nsigma TOF distribution
2017     //    fPriHistPosNsigmaTof->Fill(fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton));
2018     // Nsigma TOF in bins of total momentum
2019     fPriHistPosNsigmaTofVsP->Fill(track->P(),fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton));
2020     // Nsigma TOF in bins of transverse momentum
2021     fPriHistPosNsigmaTofVsPt->Fill(track->Pt(),fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton));
2022     
2023     // Try the tof signal instead of nsigma
2024     fPriHistTOFsignalPosVsP->Fill(track->P(), corrTOFsig);
2025     fPriHistTOFsignalPosVsPt->Fill(track->Pt(), corrTOFsig);
2026     
2027   }
2028   else if (track->Charge() < 0){
2029     // Simple Nsigma TOF distribution
2030     //    fPriHistNegNsigmaTof->Fill(fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton));
2031     // Nsigma TOF in bins of total momentum
2032     fPriHistNegNsigmaTofVsP->Fill(track->P(),fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton));
2033     // Nsigma TOF in bins of transverse momentum
2034     fPriHistNegNsigmaTofVsPt->Fill(track->Pt(),fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton));
2035     
2036     // Try the tof signal instead of nsigma
2037     fPriHistTOFsignalNegVsP->Fill(track->P(), corrTOFsig);
2038     fPriHistTOFsignalNegVsPt->Fill(track->Pt(), corrTOFsig);
2039   }
2040   
2041   // Final judging: simple first idea. 
2042   // min -800 up to 2 GeV and 0 up to 3.25GeV
2043   if (track->P() < 2.0){
2044     //    if (corrTOFsig > -800.0){
2045     // In AODs, the resolution is better, do -500 (AODs) 
2046     // instead of -800 (ESDs)
2047     if (corrTOFsig > -500.0) {
2048       // Create additional TPC only constrained tp pri. vtx track parameters
2049       //        constrainTrack(track);
2050       if (track->Charge()>0){
2051         // Cut .1 cm on DCAxy and fill a histogram
2052         if(goodDCA(track)){
2053           // Add to the femto event
2054           fFemtoBuffer->GetEvt(0)->AddPro(track);
2055         }
2056       }
2057       else{
2058         // Cut .1 cm on DCAxy and fill a histogram
2059         if(goodDCA(track)){
2060           // Add to the femto event
2061           fFemtoBuffer->GetEvt(0)->AddAPro(track);
2062         }
2063       }
2064     }
2065   }
2066   else if (track->P() < 3.25){
2067     if (corrTOFsig > 0){
2068       // Create additional TPC only constrained tp pri. vtx track parameters
2069       //        constrainTrack(track);
2070       if (track->Charge()>0){
2071         // Cut .1 cm on DCAxy and fill a histogram
2072         if(goodDCA(track)){
2073           // Add to the femto event
2074           fFemtoBuffer->GetEvt(0)->AddPro(track);
2075         }
2076       }
2077       else{
2078         // Cut .1 cm on DCAxy and fill a histogram
2079         if(goodDCA(track)){
2080           // Add to the femto event
2081           fFemtoBuffer->GetEvt(0)->AddAPro(track);
2082         }
2083       }
2084     }
2085   }
2086 } // End of void ProcessTOF
2087 //________________________________________________________________________
2088 void AliAnalysisTaskProtonLambda::ProcessTPC(AliAODTrack* track){
2089
2090   // Require the TPCpid bit
2091   if (!((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTPCpid))
2092     return;
2093     
2094   // In contrast to ESDs one doesn't check for AliESDtrack::kTOFpid
2095   // but for AliVTrack::kTOFout?? 
2096   // Check how many particles have TOFout bit
2097   // if (track->Charge() > 0){
2098   //   if ((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTOFpid)
2099   //     fPriHistHasTofPos->Fill(1.0);
2100   //   else 
2101   //     fPriHistHasTofPos->Fill(0.0);
2102   // }
2103   // else{
2104   //   if ((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTOFpid)
2105   //     fPriHistHasTofNeg->Fill(1.0);
2106   //   else 
2107   //     fPriHistHasTofNeg->Fill(0.0);
2108   // }
2109
2110   // For all plots <dE/dx> vs p one should use
2111   // the momentum at the inner wall of the TPC.
2112
2113   // Use a TOF cut and fill the same dE/dx histograms
2114   // Bool_t acceptedTOF=kFALSE;
2115   // if ((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTOFpid){
2116   //   if (fPIDResponse->NumberOfSigmasTOF((fGTI[-track->GetID()-1]), AliPID::kProton) > -10.0)
2117   //     acceptedTOF=kTRUE;
2118   // }
2119   // if (acceptedTOF){
2120   //     if (track->Charge() > 0){
2121   //    fPriHistTPCsignalTOFcutPos->Fill((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2122   //                                     (fGTI[-track->GetID()-1])->GetTPCsignal());
2123   //    fPriHistNsigmaTPCTOFcutPos->Fill((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2124   //                                     fPIDResponse->NumberOfSigmasTPC((fGTI[-track->GetID()-1]), AliPID::kProton));
2125   //     }
2126   //     else{
2127   //    fPriHistTPCsignalTOFcutNeg->Fill((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2128   //                                     (fGTI[-track->GetID()-1])->GetTPCsignal());
2129   //    fPriHistNsigmaTPCTOFcutNeg->Fill((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2130   //                                     fPIDResponse->NumberOfSigmasTPC((fGTI[-track->GetID()-1]), AliPID::kProton));
2131   //     }
2132   // }  
2133     
2134   // A first idea of a cut: use the spectra cut.
2135   // (should perhaps change for momenta ~ 0.75 GeV)
2136   if ( ((fGTI[-track->GetID()-1])->GetTPCsignal() > 
2137         fTpcResponse->GetExpectedSignal((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2138                                         AliPID::kProton))
2139       // New since Sept 10th 2012: Also use a cut to reject deuterons.
2140       // I checked: The cut is good!
2141        && ((fGTI[-track->GetID()-1])->GetTPCsignal() <
2142            2.0*fTpcResponse->GetExpectedSignal((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2143                                       AliPID::kProton))
2144                         ) {
2145     // Distinguish between charges
2146     if (track->Charge()>0){
2147       // Cut .1 cm on DCAxy and fill a histogram
2148       if(goodDCA(track)){
2149         // Add to the femto event
2150         fFemtoBuffer->GetEvt(0)->AddPro(track);
2151       }
2152     }
2153     else{
2154       // Cut .1 cm on DCAxy and fill a histogram
2155       if(goodDCA(track)){
2156         // Add to the femto event
2157         fFemtoBuffer->GetEvt(0)->AddAPro(track);
2158       }
2159     }
2160   }
2161 } // End of void ProcessTPC
2162 //________________________________________________________________________
2163 void AliAnalysisTaskProtonLambda::ProcessHybrid(AliAODTrack *track){
2164   
2165   // Intermediate momentum: use dEdx for a pre-selection
2166   // and do the pid with tof
2167   
2168   // Boolean for extra! tpc pid cuts
2169   Bool_t acceptTPC = kTRUE;
2170
2171   // Require the TPCpid bit
2172   if (!((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTPCpid))
2173     acceptTPC = kFALSE;
2174  
2175   // Pre-selection cut with TPC, don't return immediately to be able
2176   // to visualize the effect
2177   if (acceptTPC){
2178     // Do a mild dEdx cut
2179     if ((fGTI[-track->GetID()-1])->GetTPCsignal() < 
2180         fTpcResponse->GetExpectedSignal((fGTI[-track->GetID()-1])->GetTPCmomentum(),
2181                                         AliPID::kElectron))
2182       acceptTPC = kFALSE;
2183   }
2184     
2185   // Ask for TOF pid flag and fill
2186   if (!((fGTI[-track->GetID()-1])->GetStatus() & AliVTrack::kTOFpid))
2187     return;
2188   
2189   // The corrected TOF signal
2190   Double_t corrTOFsig = GetCorrectedTOFSignal(track);
2191   
2192   // Distinguish between charges
2193   if (track->Charge() > 0) {
2194     // Fill the tof signal w/o dedx pre-selection
2195     fPriHistHybridTOFsigPosWoTPC->Fill(corrTOFsig);
2196     // Do the pre-selection
2197     if (acceptTPC){
2198       fPriHistHybridTOFsigPosTPCok->Fill(corrTOFsig);
2199
2200       // Do the tof cut
2201       // Sept '12: also include an upper cut
2202       if ( (corrTOFsig > -1000.0) && (corrTOFsig < 1250.) ){
2203         // Create additional TPC only constrained to pri. vtx track parameters
2204         //      constrainTrack(track);
2205         // Cut .1 cm on DCAxy and fill a histogram
2206         if(goodDCA(track)){
2207           // Add to the femto event
2208           fFemtoBuffer->GetEvt(0)->AddPro(track);
2209         }
2210       }
2211     }
2212   }
2213   else {
2214     // Fill the tof signal w/o dedx pre-selection
2215     fPriHistHybridTOFsigNegWoTPC->Fill(corrTOFsig);
2216     // Do the pre-selection
2217     if (acceptTPC){
2218       fPriHistHybridTOFsigNegTPCok->Fill(corrTOFsig);
2219       
2220       // Do the tof cut
2221       // Sept '12: also include an upper cut
2222       if ( (corrTOFsig > -1000.0) && (corrTOFsig < 1250.) ){
2223         // Create additional TPC only constrained to pri. vtx track parameters
2224         //      constrainTrack(track);
2225         // Cut .1 cm on DCAxy and fill a histogram
2226         if(goodDCA(track)){
2227           // add to the femto event
2228           fFemtoBuffer->GetEvt(0)->AddAPro(track);
2229         }
2230       }
2231     }
2232   }
2233 } // End of ProcessHybrid
2234 //________________________________________________________________________
2235 void AliAnalysisTaskProtonLambda::CleaningProcedure() {
2236   // fFemtoBuffer->GetEvt(0) pointer must be set
2237   // Checks that no tracks are shared between Lam & Lam, Lam & Pro, ALam & ALam, ALam & APro
2238
2239   // printf ("Cleaning procedure. Lam: %d, ALam: %d, Pro: %d, APro:%d\n"
2240   //      ,fFemtoBuffer->GetEvt(0)->GetNLam(),fFemtoBuffer->GetEvt(0)->GetNALam(),fFemtoBuffer->GetEvt(0)->GetNPro(),fFemtoBuffer->GetEvt(0)->GetNAPro());
2241
2242   //
2243   // Check for lambdas..
2244   //
2245   for (Int_t i=0;i<fFemtoBuffer->GetEvt(0)->GetNLam();i++) {
2246     if (!fFemtoBuffer->GetEvt(0)->fLamTracks[i].UseIt())
2247       continue;
2248     // Unique track id's for first V0
2249     Int_t posId1 = fFemtoBuffer->GetEvt(0)->fLamTracks[i].fPosDaughter.fID;
2250     Int_t negId1 = fFemtoBuffer->GetEvt(0)->fLamTracks[i].fNegDaughter.fID;
2251
2252     // .. & lambdas
2253     for (Int_t j=i+1;j<fFemtoBuffer->GetEvt(0)->GetNLam();j++){
2254       if (!fFemtoBuffer->GetEvt(0)->fLamTracks[j].UseIt())
2255         continue;
2256       // Unique track id's for second V0
2257       Int_t posId2 = fFemtoBuffer->GetEvt(0)->fLamTracks[j].fPosDaughter.fID;
2258       Int_t negId2 = fFemtoBuffer->GetEvt(0)->fLamTracks[j].fPosDaughter.fID;
2259       
2260       // If V0s share a track remove one
2261       if (posId1 == posId2 || negId1 == negId2){
2262
2263         // printf ("shared track lamlam! posId1: %d, posId2: %d, negId1: %d, negId2: %d\n",
2264         //      posId1, posId2, negId1, negId2);
2265         
2266         // Use a criterion to select best V0
2267         if (fFemtoBuffer->GetEvt(0)->fLamTracks[i].fCosPoint > fFemtoBuffer->GetEvt(0)->fLamTracks[j].fCosPoint){
2268           fFemtoBuffer->GetEvt(0)->fLamTracks[j].SetBadFlag();
2269         }
2270         else{
2271           fFemtoBuffer->GetEvt(0)->fLamTracks[i].SetBadFlag();
2272         }
2273       }
2274       
2275     } // Scnd V0 loop
2276
2277     if (!fFemtoBuffer->GetEvt(0)->fLamTracks[i].UseIt())
2278       continue;
2279
2280     // .. & protons
2281     for (Int_t j=0;j<fFemtoBuffer->GetEvt(0)->GetNPro();j++){
2282       if (!fFemtoBuffer->GetEvt(0)->fProTracks[j].UseIt())
2283         continue;
2284       // Unique track id's for second V0
2285       Int_t posId2 = fFemtoBuffer->GetEvt(0)->fProTracks[j].fID;
2286       
2287       // If V0 and proton share a track
2288       if (posId1 == posId2){
2289         //      printf ("shared track lam p! id:%d\n",posId1);
2290         
2291         // Remove the proton
2292         fFemtoBuffer->GetEvt(0)->fProTracks[j].SetBadFlag();
2293       }
2294       
2295     } // Proton loop
2296
2297   } // First V0 loop
2298
2299   //
2300   // Check for anti-lambdas..
2301   //
2302   for (Int_t i=0;i<fFemtoBuffer->GetEvt(0)->GetNALam();i++){
2303     if (!fFemtoBuffer->GetEvt(0)->fALamTracks[i].UseIt())
2304       continue;
2305     // Unique track id's for first V0
2306     Int_t posId1 = fFemtoBuffer->GetEvt(0)->fALamTracks[i].fPosDaughter.fID;
2307     Int_t negId1 = fFemtoBuffer->GetEvt(0)->fALamTracks[i].fNegDaughter.fID;
2308
2309     // .. & anti-lambdas
2310     for (Int_t j=i+1;j<fFemtoBuffer->GetEvt(0)->GetNALam();j++){
2311       if (!fFemtoBuffer->GetEvt(0)->fALamTracks[j].UseIt())
2312         continue;
2313       // Unique track id's for second V0
2314       Int_t posId2 = fFemtoBuffer->GetEvt(0)->fALamTracks[j].fPosDaughter.fID;
2315       Int_t negId2 = fFemtoBuffer->GetEvt(0)->fALamTracks[j].fNegDaughter.fID;
2316       
2317       // If V0s share a track remove one
2318       if (posId1 == posId2 || negId1 == negId2){
2319         
2320         // printf ("shared track ALamALam! posId1: %d, posId2: %d, negId1: %d, negId2: %d\n",
2321         //      posId1, posId2, negId1, negId2);
2322
2323         // Use a criterion to select best V0
2324         if (fFemtoBuffer->GetEvt(0)->fALamTracks[i].fCosPoint > fFemtoBuffer->GetEvt(0)->fALamTracks[j].fCosPoint){
2325           fFemtoBuffer->GetEvt(0)->fALamTracks[j].SetBadFlag();
2326         }
2327         else{
2328           fFemtoBuffer->GetEvt(0)->fALamTracks[i].SetBadFlag();
2329         }
2330       }
2331       
2332     } // Scnd anti-V0 loop
2333
2334     if (!fFemtoBuffer->GetEvt(0)->fALamTracks[i].UseIt())
2335       continue;
2336     
2337     // .. & anti-protons
2338     for (Int_t j=0;j<fFemtoBuffer->GetEvt(0)->GetNAPro();j++){
2339       if (!fFemtoBuffer->GetEvt(0)->fAProTracks[j].UseIt())
2340         continue;
2341       // Unique track id's for second V0
2342       Int_t negId2 = fFemtoBuffer->GetEvt(0)->fAProTracks[j].fID;
2343       
2344       // If V0 and proton share a track
2345       if (negId1 == negId2){
2346         //      printf ("shared track alam ap! id:%d\n",posId1);
2347
2348         // Remove the proton
2349         fFemtoBuffer->GetEvt(0)->fAProTracks[j].SetBadFlag();
2350       }
2351       
2352     } // Anti-proton loop
2353
2354   } // First anti-V0 loop
2355
2356   //
2357   // Do the same with the side band background.
2358   // Discard background when sharing track with primary proton.
2359   //
2360
2361    for (Int_t i=0;i<fFemtoBuffer->GetEvt(0)->GetNBgLam();i++){
2362     if (!fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].UseIt())
2363       continue;
2364     // Unique track id's for first V0
2365     Int_t posId1 = fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].fPosDaughter.fID;
2366     Int_t negId1 = fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].fNegDaughter.fID;
2367
2368     // .. & lambdas
2369     for (Int_t j=i+1;j<fFemtoBuffer->GetEvt(0)->GetNBgLam();j++){
2370       if (!fFemtoBuffer->GetEvt(0)->fBgLamTracks[j].UseIt())
2371         continue;
2372       // Unique track id's for second V0
2373       Int_t posId2 = fFemtoBuffer->GetEvt(0)->fBgLamTracks[j].fPosDaughter.fID;
2374       Int_t negId2 = fFemtoBuffer->GetEvt(0)->fBgLamTracks[j].fNegDaughter.fID;
2375       
2376       // If V0s share a track remove one
2377       if (posId1 == posId2 || negId1 == negId2){
2378
2379         // printf ("shared track bglambglam! posId1: %d, posId2: %d, negId1: %d, negId2: %d\n",
2380         //      posId1, posId2, negId1, negId2);
2381         
2382         // Use a criterion to select best V0
2383         if (fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].fCosPoint > fFemtoBuffer->GetEvt(0)->fBgLamTracks[j].fCosPoint){
2384           fFemtoBuffer->GetEvt(0)->fBgLamTracks[j].SetBadFlag();
2385         }
2386         else{
2387           fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].SetBadFlag();
2388         }
2389       }
2390       
2391     } // Scnd V0 loop
2392
2393     if (!fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].UseIt())
2394       continue;
2395
2396     // .. & protons
2397     for (Int_t j=0;j<fFemtoBuffer->GetEvt(0)->GetNPro();j++) {
2398       if (!fFemtoBuffer->GetEvt(0)->fProTracks[j].UseIt())
2399         continue;
2400       if (!fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].UseIt())
2401         continue;
2402
2403       // Unique track id's for second V0
2404       Int_t posId2 = fFemtoBuffer->GetEvt(0)->fProTracks[j].fID;
2405       
2406       // If V0 and proton share a track
2407       if (posId1 == posId2){
2408         //      printf ("shared track bglam p! id:%d\n",posId1);
2409         // Remove the background lambda
2410         fFemtoBuffer->GetEvt(0)->fBgLamTracks[i].SetBadFlag();
2411       }
2412       
2413     } // Proton loop
2414
2415   } // First V0 loop
2416
2417   //
2418   // Check for anti-lambdas..
2419   //
2420   for (Int_t i=0;i<fFemtoBuffer->GetEvt(0)->GetNBgALam();i++){
2421     if (!fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].UseIt())
2422       continue;
2423     // Unique track id's for first V0
2424     Int_t posId1 = fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].fPosDaughter.fID;
2425     Int_t negId1 = fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].fNegDaughter.fID;
2426
2427     // .. & anti-lambdas
2428     for (Int_t j=i+1;j<fFemtoBuffer->GetEvt(0)->GetNBgALam();j++){
2429       if (!fFemtoBuffer->GetEvt(0)->fBgALamTracks[j].UseIt())
2430         continue;
2431       // Unique track id's for second V0
2432       Int_t posId2 = fFemtoBuffer->GetEvt(0)->fBgALamTracks[j].fPosDaughter.fID;
2433       Int_t negId2 = fFemtoBuffer->GetEvt(0)->fBgALamTracks[j].fNegDaughter.fID;
2434       
2435       // If V0s share a track remove one
2436       if (posId1 == posId2 || negId1 == negId2){
2437         
2438         // printf ("shared track BgALamBgALam! posId1: %d, posId2: %d, negId1: %d, negId2: %d\n",
2439         //      posId1, posId2, negId1, negId2);
2440
2441         // Use a criterion to select best V0
2442         if (fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].fCosPoint > fFemtoBuffer->GetEvt(0)->fBgALamTracks[j].fCosPoint){
2443           fFemtoBuffer->GetEvt(0)->fBgALamTracks[j].SetBadFlag();
2444         }
2445         else{
2446           fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].SetBadFlag();
2447         }
2448       }
2449       
2450     } // Scnd anti-V0 loop
2451
2452     if (!fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].UseIt())
2453       continue;
2454     
2455     // .. & anti-protons
2456     for (Int_t j=0;j<fFemtoBuffer->GetEvt(0)->GetNAPro();j++){
2457       if (!fFemtoBuffer->GetEvt(0)->fAProTracks[j].UseIt())
2458         continue;
2459       if (!fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].UseIt())
2460         continue;
2461       
2462       // Unique track id's for second V0
2463       Int_t negId2 = fFemtoBuffer->GetEvt(0)->fAProTracks[j].fID;
2464       
2465       // If V0 and proton share a track
2466       if (negId1 == negId2){
2467         //      printf ("shared track bgalam ap! id:%d\n",posId1);
2468         // Remove the background anti-lambda
2469         fFemtoBuffer->GetEvt(0)->fBgALamTracks[i].SetBadFlag();
2470       }
2471       
2472     } // Anti-proton loop
2473
2474   } // First anti-V0 loop
2475
2476   
2477 }
2478 //________________________________________________________________________
2479 void AliAnalysisTaskProtonLambda::ProcessReal() {
2480   // Process real events
2481   
2482   // // Count the number of pairs before TTR cut
2483   // Int_t nLamLamPairsWoTTR = 0,nLamProPairsWoTTR=0,nALamALamPairsWoTTR=0,nALamAProPairsWoTTR=0;
2484   // // and with TTR cut
2485   // Int_t nLamLamPairs = 0,nLamProPairs=0,nALamALamPairs=0,nALamAProPairs=0;
2486
2487   // Declare numbers to speed up the code
2488   Int_t iLam,//iLam2,
2489     iRadius,iPro,iALam,
2490     //iALam2,
2491     iAPro,nMeasPro,
2492     //nMeasPio,
2493     nMeasAPro;
2494   //  Int_t nMeasProIroc,nMeasPioIroc,nMeasAProIroc,nMeasProOroc,nMeasPioOroc,nMeasAProOroc;
2495   // Float_t distPro,distPio,minDistPro,meanDistPro,minDistPio,meanDistPio,
2496   //   distAPro,minDistAPro,meanDistAPro;
2497   Float_t distSftPro,//distSftPio,
2498     minDistSftPro,meanDistSftPro,//minDistSftPio,meanDistSftPio,
2499     distSftAPro,minDistSftAPro,meanDistSftAPro;
2500   // Float_t minDistSftIrocPro,meanDistSftIrocPro,minDistSftIrocPio,meanDistSftIrocPio,
2501   //   minDistSftIrocAPro,meanDistSftIrocAPro;
2502   // Float_t minDistSftOrocPro,meanDistSftOrocPro,minDistSftOrocPio,meanDistSftOrocPio,
2503   //   minDistSftOrocAPro,meanDistSftOrocAPro;
2504   
2505   // printf("Real event, NLam: %d, NPro %d, NALam %d, NAPro %d\n",
2506   //     fFemtoBuffer->GetEvt(0)->GetNLam(),
2507   //     fFemtoBuffer->GetEvt(0)->GetNPro(),
2508   //     fFemtoBuffer->GetEvt(0)->GetNALam(),
2509   //     fFemtoBuffer->GetEvt(0)->GetNAPro()
2510   //     );
2511   // Lambda loop
2512   for (iLam = 0; iLam < fFemtoBuffer->GetEvt(0)->GetNLam(); iLam++){
2513
2514     // Skip if unUseIt() entry
2515     if (!fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].UseIt())
2516       continue;
2517     
2518     // // Second lambda loop
2519     // for (iLam2 = iLam+1; iLam2 < fFemtoBuffer->GetEvt(0)->GetNLam(); iLam2++){
2520
2521     //   // Skip if unUseIt() entry
2522     //   if (!fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2].UseIt())
2523     //  continue;
2524
2525     //   // Count the number of lam-lam pairs
2526     //   nLamLamPairsWoTTR++;
2527
2528     //   // Reset the distances for each pair
2529     //   minDistPro=999.0;meanDistPro=0.0;minDistPio=999.0;meanDistPio=0.0;
2530     //   minDistSftPro=999.0;meanDistSftPro=0.0;minDistSftPio=999.0;meanDistSftPio=0.0;
2531     //   minDistSftIrocPro=999.0;meanDistSftIrocPro=0.0;minDistSftIrocPio=999.0;meanDistSftIrocPio=0.0;
2532     //   minDistSftOrocPro=999.0;meanDistSftOrocPro=0.0;minDistSftOrocPio=999.0;meanDistSftOrocPio=0.0;
2533     //   // Reset the number of measurements for the mean
2534     //   nMeasPro=0;nMeasPio=0;nMeasProIroc=0;nMeasPioIroc=0;
2535     //   nMeasProOroc=0;nMeasPioOroc=0;
2536
2537     //   // Check for two-track resolution
2538     //   for (iRadius=0;iRadius<9;iRadius++){
2539     //  // Get the spatial distance at each radius
2540     //  distPro = calcDist(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fPosDaughter.fXglobal[iRadius],fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2].fPosDaughter.fXglobal[iRadius]);
2541     //  distPio = calcDist(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fNegDaughter.fXglobal[iRadius],fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2].fNegDaughter.fXglobal[iRadius]);
2542     //  // Shifted distances
2543     //  distSftPro = calcDist(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fPosDaughter.fXshifted[iRadius],fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2].fPosDaughter.fXshifted[iRadius]);//,fFemtoBuffer->GetEvt(0),fFemtoBuffer->GetEvt(0));
2544     //  distSftPio = calcDist(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fNegDaughter.fXshifted[iRadius],fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2].fNegDaughter.fXshifted[iRadius]);//,fFemtoBuffer->GetEvt(0),fFemtoBuffer->GetEvt(0));
2545
2546     //  // calcDist returns -2.0 if one shouldn't use this distance, see calcDist function
2547     //  if (distPro > -1.0) {
2548     //    // Minimum distance
2549     //    if (distPro < minDistPro)
2550     //      minDistPro = distPro;
2551     //    if (distSftPro < minDistSftPro)
2552     //      minDistSftPro = distSftPro;
2553     //    // Mean distance
2554     //    meanDistPro+=distPro;
2555     //    meanDistSftPro+=distSftPro;
2556     //    nMeasPro++;
2557         
2558     //    // IROC
2559     //    if (iRadius<3){
2560     //      if (distSftPro < minDistSftIrocPro)
2561     //        minDistSftIrocPro = distSftPro;
2562     //      meanDistSftIrocPro+=distSftPro;
2563     //      nMeasProIroc++;
2564     //    }
2565     //    // OROC
2566     //    else {
2567     //      if (distSftPro < minDistSftOrocPro)
2568     //        minDistSftOrocPro = distSftPro;
2569     //      meanDistSftOrocPro+=distSftPro;
2570     //      nMeasProOroc++;         
2571     //    }
2572     //  }
2573     //  if (distPio > -1.0){
2574     //    // Minimum distance
2575     //    if (distPio < minDistPio)
2576     //      minDistPio = distPio;
2577     //    if (distSftPio < minDistSftPio)
2578     //      minDistSftPio = distSftPio;
2579     //    // Mean distance
2580     //    meanDistPio+=distPio;
2581     //    meanDistSftPio+=distSftPio;
2582     //    nMeasPio++; 
2583           
2584     //    // IROC
2585     //    if (iRadius<3){
2586     //      if (distSftPio < minDistSftIrocPio)
2587     //        minDistSftIrocPio = distSftPio;
2588     //      meanDistSftIrocPio+=distSftPio;
2589     //      nMeasPioIroc++;
2590     //    }
2591     //    // OROC
2592     //    else {
2593     //      if (distSftPio < minDistSftOrocPio)
2594     //        minDistSftOrocPio = distSftPio;
2595     //      meanDistSftOrocPio+=distSftPio;
2596     //      nMeasPioOroc++;         
2597     //    }
2598
2599     //  }
2600
2601     //   } // Loop over iRadius
2602       
2603     //   // Require at least one measurement
2604     //   if ( (!nMeasPio) || (!nMeasPro) )
2605     //  continue;
2606
2607     //   // Divide by the number of measurements to get the mean
2608     //   meanDistPro /= (Float_t)nMeasPro;
2609     //   meanDistPio /= (Float_t)nMeasPio;
2610     //   meanDistSftPro /= (Float_t)nMeasPro;
2611     //   meanDistSftPio /= (Float_t)nMeasPio;
2612       
2613     //   // Fill the two track resolution histograms
2614     //   f2HistLamLamMeanMinDistProReal->Fill(meanDistPro,minDistPro);
2615     //   f2HistLamLamMeanMinDistPioReal->Fill(meanDistPio,minDistPio);
2616
2617     //   f2HistSftLamLamMeanMinDistProReal->Fill(meanDistSftPro,minDistSftPro);
2618     //   f2HistSftLamLamMeanMinDistPioReal->Fill(meanDistSftPio,minDistSftPio);
2619
2620     //   // Fill IROC / OROC histograms only with at least one measurement
2621     //   if (nMeasProIroc){
2622     //  meanDistSftIrocPro /= (Float_t)nMeasProIroc;
2623     //  f2HistSftIrocLamLamMeanMinDistProReal->Fill(meanDistSftIrocPro,minDistSftIrocPro);
2624     //   }
2625     //   if (nMeasPioIroc){
2626     //  meanDistSftIrocPio /= (Float_t)nMeasPioIroc;
2627     //  f2HistSftIrocLamLamMeanMinDistPioReal->Fill(meanDistSftIrocPio,minDistSftIrocPio);
2628     //   }
2629     //   if (nMeasProOroc){
2630     //  meanDistSftOrocPro /= (Float_t)nMeasProOroc;
2631     //  f2HistSftOrocLamLamMeanMinDistProReal->Fill(meanDistSftOrocPro,minDistSftOrocPro);
2632     //   }
2633     //   if (nMeasPioOroc){
2634     //  meanDistSftOrocPio /= (Float_t)nMeasPioOroc;
2635     //  f2HistSftOrocLamLamMeanMinDistPioReal->Fill(meanDistSftOrocPio,minDistSftOrocPio);
2636     //   }
2637
2638     //   // Do a cut (value needs to be refined)
2639     //   //      if ( meanDistSftPro < 2.0 || meanDistSftPio < 2.0 )
2640     //   //     continue;
2641       
2642     //   // Count the number of pairs
2643     //   nLamLamPairs++;
2644
2645     //   // Mt of the pair
2646     //   f2HistMtLamLamReal->Fill(mt(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2647     //                                    fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2]));
2648
2649     //   // Fill the qinv, minv histogram
2650     //   f3HistLamLamQinvReal->Fill(Qinv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam], fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2]),minDistSftPro,minDistSftPio);
2651     //   f3HistLamLamMinvReal->Fill(Minv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam], fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2]),minDistSftPro,minDistSftPio);
2652       
2653     //   // Mt of the pair fr low q pairs only
2654     //   if(Qinv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2655     //        fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2]) < .2)
2656     //  f2HistMtLowQLamLamReal->Fill(mt(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2657     //                                  fFemtoBuffer->GetEvt(0)->fLamTracks[iLam2]));
2658
2659     // } // Second lambda loop
2660
2661     // Proton loop
2662     for (iPro=0;iPro<fFemtoBuffer->GetEvt(0)->GetNPro();iPro++){
2663
2664       // Skip if unUseIt() entry
2665       if (!fFemtoBuffer->GetEvt(0)->fProTracks[iPro].UseIt())
2666         continue;
2667
2668       //      printf(".");
2669       //      nLamProPairsWoTTR++;
2670
2671       // Reset the distances for each pair
2672       // minDistPro=999.0;meanDistPro=0.0;
2673       minDistSftPro=999.0;meanDistSftPro=0.0;
2674       //minDistSftIrocPro=999.0;meanDistSftIrocPro=0.0;minDistSftOrocPro=999.0;meanDistSftOrocPro=0.0;
2675       // Reset the number of measurements for the mean
2676       nMeasPro=0;//nMeasProIroc=0;nMeasProOroc=0;
2677       
2678       // Check for two-track resolution
2679       for (iRadius=0;iRadius<9;iRadius++){
2680         // Get the spatial distance at each radius
2681         //      distPro = calcDist(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fPosDaughter.fXglobal[iRadius],fFemtoBuffer->GetEvt(0)->fProTracks[iPro].fXglobal[iRadius]);
2682         distSftPro = calcDist(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fPosDaughter.fXshifted[iRadius],fFemtoBuffer->GetEvt(0)->fProTracks[iPro].fXshifted[iRadius]);//,fFemtoBuffer->GetEvt(0),fFemtoBuffer->GetEvt(0));
2683
2684         // calcDist returns -2.0 if one shouldn't use this distance, see calcDist function
2685         if (distSftPro > -1.0){
2686           // Minimum distance
2687           // if (distPro < minDistPro)
2688           //   minDistPro = distPro;
2689           if (distSftPro < minDistSftPro)
2690             minDistSftPro = distSftPro;
2691           // Mean distance
2692           //      meanDistPro+=distPro;
2693           meanDistSftPro+=distSftPro;
2694           nMeasPro++;
2695           
2696           // // IROC                                   
2697           // if (iRadius<3){
2698           //   if (distSftPro < minDistSftIrocPro)
2699           //     minDistSftIrocPro = distSftPro;
2700           //   meanDistSftIrocPro+=distSftPro;
2701           //   nMeasProIroc++;
2702           // }
2703           // // OROC                                   
2704           // else {
2705           //   if (distSftPro < minDistSftOrocPro)
2706           //     minDistSftOrocPro = distSftPro;
2707           //   meanDistSftOrocPro+=distSftPro;
2708           //   nMeasProOroc++;
2709           // }
2710
2711
2712         }
2713       } // Loop over iRadius
2714
2715       // Require at least one measurement
2716       if ( !nMeasPro )
2717         continue;
2718
2719       // Divide by the number of measurements to get the mean
2720       //      meanDistPro /= (Float_t)nMeasPro;
2721       meanDistSftPro /= (Float_t)nMeasPro;
2722       
2723       // Fill the two track resolution histogram
2724       // f2HistLamProMeanMinDistProReal->Fill(meanDistPro,minDistPro);
2725       // f2HistSftLamProMeanMinDistProReal->Fill(meanDistSftPro,minDistSftPro);
2726
2727       // // Fill IROC / OROC histograms only with at least one measurement                                
2728       // if (nMeasProIroc){
2729       //   meanDistSftIrocPro /= (Float_t)nMeasProIroc;
2730       //   f2HistSftIrocLamProMeanMinDistProReal->Fill(meanDistSftIrocPro,minDistSftIrocPro);
2731       // }
2732       // if (nMeasProOroc){
2733       //   meanDistSftOrocPro /= (Float_t)nMeasProOroc;
2734       //   f2HistSftOrocLamProMeanMinDistProReal->Fill(meanDistSftOrocPro,minDistSftOrocPro);
2735       // }
2736         
2737       // Do a cut (value needs to be refined)
2738       //      if ( meanDistSftPro < 2.0 )
2739       //        continue;
2740
2741       // Look at possible residual correlations of the daughters
2742       // vs mean dist of protons
2743       //      f2HistLamPosDProQinvReal->Fill(QinvProPro(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fPosDaughter, fFemtoBuffer->GetEvt(0)->fProTracks[iPro]),meanDistSftPro);
2744       //      f2HistLamNegDProQinvReal->Fill(QinvPioPro(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fNegDaughter, fFemtoBuffer->GetEvt(0)->fProTracks[iPro]),meanDistSftPro);
2745   
2746
2747       // Count the number of pairs
2748       //      nLamProPairs++;
2749
2750       // Mt of the pair
2751       f2HistMtLamProReal->Fill(mt(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2752                                   fFemtoBuffer->GetEvt(0)->fProTracks[iPro]));
2753
2754       // THnSparse with qinvpropro, mean dist propro, min dist propro, qinv lampro
2755       Double_t x[4]={
2756         QinvProPro(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam].fPosDaughter,
2757                    fFemtoBuffer->GetEvt(0)->fProTracks[iPro]),
2758         meanDistSftPro,
2759         minDistSftPro,
2760         Qinv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2761              fFemtoBuffer->GetEvt(0)->fProTracks[iPro])
2762       };
2763       LamProReal->Fill(x);
2764
2765       // Fill the qinv histogram, using TPC only momentum for primary protons
2766       //      f3HistLamProQinvReal->Fill(Qinv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam], fFemtoBuffer->GetEvt(0)->fProTracks[iPro]),meanDistSftPro,minDistSftPro); // Using THnSparse since sept '12
2767       //      f3HistLamProMinvReal->Fill(Minv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam], fFemtoBuffer->GetEvt(0)->fProTracks[iPro]),meanDistSftPro,minDistSftPro);
2768       // and using TPC only plus primary vertex constraint momentum for primary protons
2769       //      f3HistLamProQinvConstrReal->Fill(QinvConstr(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam], fFemtoBuffer->GetEvt(0)->fProTracks[iPro]),meanDistSftPro,minDistSftPro);
2770
2771       // Mt of the pair for low q pairs only
2772       if(Qinv(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2773               fFemtoBuffer->GetEvt(0)->fProTracks[iPro]) < .2)
2774         f2HistMtLowQLamProReal->Fill(mt(fFemtoBuffer->GetEvt(0)->fLamTracks[iLam],
2775                                         fFemtoBuffer->GetEvt(0)->fProTracks[iPro]));
2776       
2777     }// Proton loop
2778
2779   }// Lambda loop
2780
2781
2782   
2783   // Anti-lambda loop
2784   for (iALam = 0; iALam < fFemtoBuffer->GetEvt(0)->GetNALam(); iALam++){
2785
2786     // Skip if unUseIt() entry
2787     if (!fFemtoBuffer->GetEvt(0)->fALamTracks[iALam].UseIt())
2788       continue;
2789
2790     // // Second anti-lambda loop
2791     // for (iALam2 = iALam+1; iALam2 < fFemtoBuffer->GetEvt(0)->GetNALam(); iALam2++){
2792
2793     //   // Skip if unUseIt() entry
2794     //   if (!fFemtoBuffer->GetEvt(0)->fALamTracks[iALam2].UseIt())
2795     //  continue;
2796
2797     //   // Count the number of antilam-antilam pairs
2798     //   nALamALamPairsWoTTR++;
2799
2800     //   // Reset the distances for each pair
2801     //   minDistAPro=999.0;meanDistAPro=0.0;minDistPio=999.0;meanDistPio=0.0;
2802     //   minDistSftAPro=999.0;meanDistSftAPro=0.0;minDistSftPio=999.0;meanDistSftPio=0.0;
2803     //   minDistSftIrocAPro=999.0;meanDistSftIrocAPro=0.0;minDistSftIrocPio=999.0;meanDistSftIrocPio=0.0;
2804     //   minDistSftOrocAPro=999.0;meanDistSftOrocAPro=0.0;minDistSftOrocPio=999.0;meanDistSftOrocPio=0.0;
2805     //   // Reset the number of measurements for the mean
2806     //   nMeasAPro=0;nMeasPio=0;nMeasAProIroc=0;nMeasPioIroc=0;nMeasAProOroc=0;nMeasPioOroc=0;
2807
2808     //   // Check for two-track resolution
2809     //   for (iRadius=0;iRadius<9;iRadius++){
2810     //  // Get the spatial distance at each radius
2811     //  distPio = calcDist(fFemtoBuffer->GetEvt(0)->fALamTracks[iALam].fPosDaughter.fXglobal[iRadius],fFemtoBuffer->GetEvt(0)->fALamTracks[iALam2].fPosDaughter.fXglobal[iRadius]);
2812     //  distAPro = calcDist(fFemtoBuffer->GetEvt(0)->fALamTracks[iALam].fNegDaughter.fXglobal[iRadius],fFemtoBuffer->GetEvt(0)->fALamTracks[iALam2].fNegDaughter.fXglobal[iRadius]);
2813     //  // Shifted distances
2814     //  distSftPio = calcDist(fFemtoBuffer->GetEvt(0)->fALamTracks[iALam].fPosDaughter.fXshifted[iRadius],fFemtoBuffer->GetEvt(0)->fALamTracks[iALam2].fPosDaughter.fXshifted[iRadius]);//,fFemtoBuffer->GetEvt(0),fFemtoBuffer->GetEvt(0));
2815     //  distSftAPro = calcDist(fFemtoBuffer->GetEvt(0)->fALamTracks[iALam].fNegDaughter.fXshifted[iRadius],fFemtoBuffer->GetEvt(0)->fALamTracks[iALam2].fNegDaughter.fXshifted[iRadius]);//,fFemtoBuffer->GetEvt(0),fFemtoBuffer->GetEvt(0));
2816
2817     //  // calcDist returns -2.0 if one shouldn't use this distance, see calcDist function
2818     //     if (distAPro > -1.0){
2819     //       // Minimum distance
2820     //       if (distAPro < minDistAPro)
2821     //         minDistAPro = distAPro;
2822     //    if (distSftAPro < minDistSftAPro)
2823     //         minDistSftAPro = distSftAPro;
2824     //       // Mean distance
2825     //       meanDistAPro+=distAPro;
2826     //       meanDistSftAPro+=distSftAPro;
2827     //       nMeasAPro++;
2828
2829     //    // IROC                                                      
2830     //       if (iRadius<3){
2831     //         if (distSftAPro < minDistSftIrocAPro)
2832     //           minDistSftIrocAPro = distSftAPro;
2833     //         meanDistSftIrocAPro+=distSftAPro;
2834     //         nMeasAProIroc++;
2835     //       }
2836     //       // OROC                                                         
2837     //       else {
2838     //         if (distSftAPro < minDistSftOrocAPro)
2839     //           minDistSftOrocAPro = distSftAPro;
2840     //         meanDistSftOrocAPro+=distSftAPro;
2841     //         nMeasAProOroc++;
2842     //       }
2843
2844     //     }
2845     //     if (distPio > -1.0){
2846     //       // Minimum distance
2847     //       if (distPio < minDistPio)
2848     //         minDistPio = distPio;
2849     //       if (distSftPio < minDistSftPio)
2850     //         minDistSftPio = distSftPio;
2851     //       // Mean distance
2852     //       meanDistPio+=distPio;
2853     //       meanDistSftPio+=distSftPio;
2854     //       nMeasPio++; 
2855
2856     //    // IROC                                                         
2857     //       if (iRadius<3){
2858     //         if (distSftPio < minDistSftIrocPio)
2859     //           minDistSftIrocPio = distSftPio;
2860     //         meanDistSftIrocPio+=distSftPio;
2861     //         nMeasPioIroc++;
2862     //       }
2863     //       // OROC                                                         
2864     //       else {
2865     //         if (distSftPio < minDistSftOrocPio)
2866     //           minDistSftOrocPio = distSftPio;
2867     //         meanDistSftOrocPio+=distSftPio;
2868     //         nMeasPioOroc++;
2869     //       }
2870     //     }
2871     //   } // Loop over iRadius
2872       
2873     //   // Require at least one measurement
2874     //   if ( (!nMeasPio) || (!nMeasAPro) )
2875     //     continue;
2876
2877     //   // Divide by the number of measurements to get the mean
2878     //   meanDistAPro /= (Float_t)nMeasAPro;
2879     //   meanDistPio /= (Float_t)nMeasPio;
2880     //   meanDistSftAPro /= (Float_t)nMeasAPro;
2881     //   meanDistSftPio /= (Float_t)nMeasPio;
2882
2883     //   // Fill the two track resolution histograms
2884     //   f2HistALamALamMeanMinDistAProReal->Fill(meanDistAPro,minDistAPro);
2885     //   f2HistALamALamMeanMinDistPioReal->Fill(meanDistPio,minDistPio);
2886
2887     //   f2HistSftALamALamMeanMinDistAProReal->Fill(meanDistSftAPro,minDistSftAPro);