]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGLF/STRANGENESS/LambdaK0/AliAnalysisTaskExtractPerformanceV0.cxx
Add a boolean flag to enable a more lightweight operation
[u/mrichter/AliRoot.git] / PWGLF / STRANGENESS / LambdaK0 / AliAnalysisTaskExtractPerformanceV0.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
17 //
18 // Modified version of AliAnalysisTaskCheckCascade.cxx.
19 // This is a 'hybrid' output version, in that it uses a classic TTree
20 // ROOT object to store the candidates, plus a couple of histograms filled on
21 // a per-event basis for storing variables too numerous to put in a tree. 
22 //
23 // --- Adapted to look for lambdas as well, using code from 
24 //        AliAnalysisTaskCheckPerformanceStrange.cxx
25 //
26 //  --- Algorithm Description 
27 //   1. Loop over primaries in stack to acquire generated charged Xi
28 //   2. Loop over stack to find V0s, fill TH3Fs "PrimRawPt"s for Efficiency
29 //   3. Perform Physics Selection
30 //   4. Perform Primary Vertex |z|<10cm selection
31 //   5. Perform Primary Vertex NoTPCOnly vertexing selection (>0 contrib.)
32 //   6. Perform Pileup Rejection
33 //   7. Analysis Loops: 
34 //    7a. Fill TH3Fs "PrimAnalysisPt" for control purposes only
35 //    7b. Fill TTree object with V0 information, candidates
36 //
37 //  Please Report Any Bugs! 
38 //
39 //   --- David Dobrigkeit Chinellato
40 //        (david.chinellato@gmail.com)
41 //
42 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
43
44 class TTree;
45 class TParticle;
46 class TVector3;
47
48 //class AliMCEventHandler;
49 //class AliMCEvent;
50 //class AliStack;
51
52 class AliESDVertex;
53 class AliAODVertex;
54 class AliESDv0;
55 class AliAODv0;
56
57 #include <Riostream.h>
58 #include "TList.h"
59 #include "TH1.h"
60 #include "TH2.h"
61 #include "TH3.h"
62 #include "TFile.h"
63 #include "THnSparse.h"
64 #include "TVector3.h"
65 #include "TCanvas.h"
66 #include "TMath.h"
67 #include "TLegend.h"
68 #include "AliLog.h"
69
70 #include "AliESDEvent.h"
71 #include "AliAODEvent.h"
72 #include "AliV0vertexer.h"
73 #include "AliCascadeVertexer.h"
74 #include "AliESDpid.h"
75 #include "AliESDtrack.h"
76 #include "AliESDtrackCuts.h"
77 #include "AliInputEventHandler.h"
78 #include "AliAnalysisManager.h"
79 #include "AliMCEventHandler.h"
80 #include "AliMCEvent.h"
81 #include "AliStack.h"
82
83 #include "AliCFContainer.h"
84 #include "AliMultiplicity.h"
85 #include "AliAODMCParticle.h"
86 #include "AliESDcascade.h"
87 #include "AliAODcascade.h"
88 #include "AliESDUtils.h"
89 #include "AliGenEventHeader.h"
90
91 #include "AliAnalysisTaskExtractPerformanceV0.h"
92
93 using std::cout;
94 using std::endl;
95
96 ClassImp(AliAnalysisTaskExtractPerformanceV0)
97
98 AliAnalysisTaskExtractPerformanceV0::AliAnalysisTaskExtractPerformanceV0() 
99   : AliAnalysisTaskSE(), fListHistV0(0), fTree(0), fPIDResponse(0), fESDtrackCuts(0),
100    fkIsNuclear   ( kFALSE ), 
101    fkSwitchINT7  ( kFALSE ),
102    fkUseOnTheFly ( kFALSE ),
103    fkTakeAllTracks ( kFALSE ),
104    fpArapidityShift ( 0.465 ),
105   fCentralityEstimator("V0M"),
106   fkLightWeight  ( kFALSE ),  
107 //------------------------------------------------
108 // Tree Variables 
109
110   fTreeVariablePrimaryStatus(0),
111   fTreeVariablePrimaryStatusMother(0),
112   fTreeVariableChi2V0(0),
113         fTreeVariableDcaV0Daughters(0),
114         fTreeVariableDcaV0ToPrimVertex(0),
115         fTreeVariableDcaPosToPrimVertex(0),
116         fTreeVariableDcaNegToPrimVertex(0),
117         fTreeVariableV0CosineOfPointingAngle(0),
118         fTreeVariableV0Radius(0),
119         fTreeVariablePt(0),
120         fTreeVariablePtMC(0),
121         fTreeVariableRapK0Short(0),
122         fTreeVariableRapLambda(0),
123         fTreeVariableRapMC(0),
124         fTreeVariableInvMassK0s(0),
125         fTreeVariableInvMassLambda(0),
126         fTreeVariableInvMassAntiLambda(0),
127         fTreeVariableAlphaV0(0),
128         fTreeVariablePtArmV0(0),
129         fTreeVariableNegTotMomentum(0),
130         fTreeVariablePosTotMomentum(0),
131         fTreeVariableNegTransvMomentum(0),
132         fTreeVariablePosTransvMomentum(0),
133         fTreeVariableNegTransvMomentumMC(0),
134         fTreeVariablePosTransvMomentumMC(0),
135    
136         fTreeVariableNSigmasPosProton(0),
137         fTreeVariableNSigmasPosPion(0),
138         fTreeVariableNSigmasNegProton(0),
139         fTreeVariableNSigmasNegPion(0),
140
141         fTreeVariablePtMother(0),
142         fTreeVariableV0CreationRadius(0),
143   fTreeVariablePID(0),
144   fTreeVariablePIDPositive(0),
145   fTreeVariablePIDNegative(0),
146   fTreeVariablePIDMother(0),
147   fTreeVariableIndexStatus(0),
148   fTreeVariableIndexStatusMother(0),
149
150   fTreeVariableRunNumber(0),
151   fTreeVariableEventNumber(0),
152
153         fTreeVariableDistOverTotMom(0),
154
155         fTreeVariablePosEta(0),
156         fTreeVariableNegEta(0),
157
158         fTreeVariableVertexZ(0),
159
160   fTreeVariableLeastNbrCrossedRows(0),
161   fTreeVariableLeastRatioCrossedRowsOverFindable(0),
162   fTreeVariableMultiplicity(0),
163   fTreeVariableMultiplicityMC(0),
164
165   fTreeVariableV0x(0),
166   fTreeVariableV0y(0),
167   fTreeVariableV0z(0),
168
169   fTreeVariableV0Px(0),
170   fTreeVariableV0Py(0),
171   fTreeVariableV0Pz(0),
172
173   fTreeVariableMCV0x(0),
174   fTreeVariableMCV0y(0),
175   fTreeVariableMCV0z(0),
176
177   fTreeVariableMCV0Px(0),
178   fTreeVariableMCV0Py(0),
179   fTreeVariableMCV0Pz(0),
180
181   fTreeVariablePVx(0),
182   fTreeVariablePVy(0),
183   fTreeVariablePVz(0),
184
185   fTreeVariableMCPVx(0),
186   fTreeVariableMCPVy(0),
187   fTreeVariableMCPVz(0),
188
189   fTreeVariableIsNonInjected(0),
190
191   fTreeVariableNegTrackStatus(0),
192   fTreeVariablePosTrackStatus(0),
193
194 //------------------------------------------------
195 // HISTOGRAMS
196 // --- Filled on an Event-by-event basis
197 //------------------------------------------------
198    fHistV0MultiplicityBeforeTrigSel(0),
199    fHistV0MultiplicityForTrigEvt(0), 
200    fHistV0MultiplicityForSelEvt(0),
201    fHistV0MultiplicityForSelEvtNoTPCOnly(0),
202    fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup(0),
203    fHistMultiplicityBeforeTrigSel(0),
204    fHistMultiplicityForTrigEvt(0),
205    fHistMultiplicity(0),
206    fHistMultiplicityNoTPCOnly(0),
207    fHistMultiplicityNoTPCOnlyNoPileup(0),
208
209         f2dHistMultiplicityVsTrueBeforeTrigSel(0),
210         f2dHistMultiplicityVsTrueForTrigEvt(0),
211         f2dHistMultiplicityVsTrue(0),
212         f2dHistMultiplicityVsTrueNoTPCOnly(0),
213         f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup(0),
214
215   //Raw Data for Vertex Z position estimator change
216         f2dHistMultiplicityVsVertexZBeforeTrigSel(0),
217         f2dHistMultiplicityVsVertexZForTrigEvt(0),
218         f2dHistMultiplicityVsVertexZ(0),
219         f2dHistMultiplicityVsVertexZNoTPCOnly(0),
220         f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup(0),
221
222   fHistGenVertexZBeforeTrigSel(0),     
223   fHistGenVertexZForTrigEvt(0),
224   fHistGenVertexZ(0),
225   fHistGenVertexZNoTPCOnly(0),
226   fHistGenVertexZNoTPCOnlyNoPileup(0),
227
228 //------------------------------------------------
229 // PARTICLE HISTOGRAMS
230 // --- Filled on a Particle-by-Particle basis
231 //------------------------------------------------
232    f3dHistPrimAnalysisPtVsYVsMultLambda(0),
233    f3dHistPrimAnalysisPtVsYVsMultAntiLambda(0),
234    f3dHistPrimAnalysisPtVsYVsMultK0Short(0),
235    f3dHistPrimAnalysisPtVsYCMSVsMultLambda(0),
236    f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda(0),
237    f3dHistPrimAnalysisPtVsYCMSVsMultK0Short(0),
238    f3dHistPrimRawPtVsYVsMultLambda(0),
239    f3dHistPrimRawPtVsYVsMultAntiLambda(0),
240    f3dHistPrimRawPtVsYVsMultK0Short(0),
241    f3dHistPrimRawPtVsYCMSVsMultLambda(0),
242    f3dHistPrimRawPtVsYCMSVsMultAntiLambda(0),
243    f3dHistPrimRawPtVsYCMSVsMultK0Short(0),
244    f3dHistPrimRawPtVsYVsMultNonInjLambda(0),
245    f3dHistPrimRawPtVsYVsMultNonInjAntiLambda(0),
246    f3dHistPrimRawPtVsYVsMultNonInjK0Short(0),
247    f3dHistPrimRawPtVsYVsMultMCLambda(0),
248    f3dHistPrimRawPtVsYVsMultMCAntiLambda(0),
249    f3dHistPrimRawPtVsYVsMultMCK0Short(0),
250    f3dHistPrimRawPtVsYVsVertexZLambda(0),
251    f3dHistPrimRawPtVsYVsVertexZAntiLambda(0),
252    f3dHistPrimRawPtVsYVsVertexZK0Short(0),
253    f3dHistPrimCloseToPVPtVsYVsMultLambda(0),
254    f3dHistPrimCloseToPVPtVsYVsMultAntiLambda(0),
255    f3dHistPrimCloseToPVPtVsYVsMultK0Short(0),
256    f3dHistPrimRawPtVsYVsDecayLengthLambda(0),
257    f3dHistPrimRawPtVsYVsDecayLengthAntiLambda(0),
258    f3dHistPrimRawPtVsYVsDecayLengthK0Short(0),
259    f3dHistGenPtVsYVsMultXiMinus(0),
260    f3dHistGenPtVsYVsMultXiPlus(0),
261    f3dHistGenPtVsYVsMultOmegaMinus(0),
262    f3dHistGenPtVsYVsMultOmegaPlus(0),
263    f3dHistGenSelectedPtVsYVsMultXiMinus(0),
264    f3dHistGenSelectedPtVsYVsMultXiPlus(0),
265    f3dHistGenSelectedPtVsYVsMultOmegaMinus(0),
266    f3dHistGenSelectedPtVsYVsMultOmegaPlus(0),
267    f3dHistGenPtVsYCMSVsMultXiMinus(0),
268    f3dHistGenPtVsYCMSVsMultXiPlus(0),
269    f3dHistGenPtVsYCMSVsMultOmegaMinus(0),
270    f3dHistGenPtVsYCMSVsMultOmegaPlus(0),
271    f3dHistGenSelectedPtVsYCMSVsMultXiMinus(0),
272    f3dHistGenSelectedPtVsYCMSVsMultXiPlus(0),
273    f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus(0),
274    f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus(0),
275    fHistPVx(0),
276    fHistPVy(0),
277    fHistPVz(0),
278    fHistPVxAnalysis(0),
279    fHistPVyAnalysis(0),
280    fHistPVzAnalysis(0),
281    fHistPVxAnalysisHasHighPtLambda(0),
282    fHistPVyAnalysisHasHighPtLambda(0),
283    fHistPVzAnalysisHasHighPtLambda(0),
284    fHistSwappedV0Counter(0)
285 {
286   // Dummy Constructor
287 }
288
289 AliAnalysisTaskExtractPerformanceV0::AliAnalysisTaskExtractPerformanceV0(const char *name) 
290   : AliAnalysisTaskSE(name), fListHistV0(0), fTree(0), fPIDResponse(0), fESDtrackCuts(0),
291    fkIsNuclear   ( kFALSE ), 
292    fkSwitchINT7  ( kFALSE ),
293    fkUseOnTheFly ( kFALSE ),
294    fkTakeAllTracks ( kFALSE ),
295    fpArapidityShift ( 0.465 ),
296   fCentralityEstimator("V0M"),
297   fkLightWeight  ( kFALSE ),  
298 //------------------------------------------------
299 // Tree Variables 
300
301   fTreeVariablePrimaryStatus(0),
302   fTreeVariablePrimaryStatusMother(0),
303   fTreeVariableChi2V0(0),
304         fTreeVariableDcaV0Daughters(0),
305         fTreeVariableDcaV0ToPrimVertex(0),
306         fTreeVariableDcaPosToPrimVertex(0),
307         fTreeVariableDcaNegToPrimVertex(0),
308         fTreeVariableV0CosineOfPointingAngle(0),
309         fTreeVariableV0Radius(0),
310         fTreeVariablePt(0),
311         fTreeVariablePtMC(0),
312         fTreeVariableRapK0Short(0),
313         fTreeVariableRapLambda(0),
314         fTreeVariableRapMC(0),
315         fTreeVariableInvMassK0s(0),
316         fTreeVariableInvMassLambda(0),
317         fTreeVariableInvMassAntiLambda(0),
318         fTreeVariableAlphaV0(0),
319         fTreeVariablePtArmV0(0),
320         fTreeVariableNegTotMomentum(0),
321         fTreeVariablePosTotMomentum(0),
322         fTreeVariableNegTransvMomentum(0),
323         fTreeVariablePosTransvMomentum(0),
324         fTreeVariableNegTransvMomentumMC(0),
325         fTreeVariablePosTransvMomentumMC(0),
326    
327         fTreeVariableNSigmasPosProton(0),
328         fTreeVariableNSigmasPosPion(0),
329         fTreeVariableNSigmasNegProton(0),
330         fTreeVariableNSigmasNegPion(0),
331
332         fTreeVariablePtMother(0),
333         fTreeVariableV0CreationRadius(0),
334   fTreeVariablePID(0),
335   fTreeVariablePIDPositive(0),
336   fTreeVariablePIDNegative(0),
337   fTreeVariablePIDMother(0),
338   fTreeVariableIndexStatus(0),
339   fTreeVariableIndexStatusMother(0),
340
341   fTreeVariableRunNumber(0),
342   fTreeVariableEventNumber(0),
343
344         fTreeVariableDistOverTotMom(0),
345
346         fTreeVariablePosEta(0),
347         fTreeVariableNegEta(0),
348
349         fTreeVariableVertexZ(0),
350
351   fTreeVariableLeastNbrCrossedRows(0),
352   fTreeVariableLeastRatioCrossedRowsOverFindable(0),
353   fTreeVariableMultiplicity(0),
354   fTreeVariableMultiplicityMC(0),
355
356   fTreeVariableV0x(0),
357   fTreeVariableV0y(0),
358   fTreeVariableV0z(0),
359
360   fTreeVariableV0Px(0),
361   fTreeVariableV0Py(0),
362   fTreeVariableV0Pz(0),
363
364   fTreeVariableMCV0x(0),
365   fTreeVariableMCV0y(0),
366   fTreeVariableMCV0z(0),
367
368   fTreeVariableMCV0Px(0),
369   fTreeVariableMCV0Py(0),
370   fTreeVariableMCV0Pz(0),
371
372   fTreeVariablePVx(0),
373   fTreeVariablePVy(0),
374   fTreeVariablePVz(0),
375
376   fTreeVariableMCPVx(0),
377   fTreeVariableMCPVy(0),
378   fTreeVariableMCPVz(0),
379
380   fTreeVariableIsNonInjected(0),
381
382   fTreeVariableNegTrackStatus(0),
383   fTreeVariablePosTrackStatus(0),
384
385 //------------------------------------------------
386 // HISTOGRAMS
387 // --- Filled on an Event-by-event basis
388 //------------------------------------------------
389    fHistV0MultiplicityBeforeTrigSel(0),
390    fHistV0MultiplicityForTrigEvt(0), 
391    fHistV0MultiplicityForSelEvt(0),
392    fHistV0MultiplicityForSelEvtNoTPCOnly(0),
393    fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup(0),
394    fHistMultiplicityBeforeTrigSel(0),
395    fHistMultiplicityForTrigEvt(0),
396    fHistMultiplicity(0),
397    fHistMultiplicityNoTPCOnly(0),
398    fHistMultiplicityNoTPCOnlyNoPileup(0),
399
400         f2dHistMultiplicityVsTrueBeforeTrigSel(0),
401         f2dHistMultiplicityVsTrueForTrigEvt(0),
402         f2dHistMultiplicityVsTrue(0),
403         f2dHistMultiplicityVsTrueNoTPCOnly(0),
404         f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup(0),
405
406   //Raw Data for Vertex Z position estimator change
407         f2dHistMultiplicityVsVertexZBeforeTrigSel(0),
408         f2dHistMultiplicityVsVertexZForTrigEvt(0),
409         f2dHistMultiplicityVsVertexZ(0),
410         f2dHistMultiplicityVsVertexZNoTPCOnly(0),
411         f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup(0),
412
413   fHistGenVertexZBeforeTrigSel(0),     
414   fHistGenVertexZForTrigEvt(0),
415   fHistGenVertexZ(0),
416   fHistGenVertexZNoTPCOnly(0),
417   fHistGenVertexZNoTPCOnlyNoPileup(0),
418
419 //------------------------------------------------
420 // PARTICLE HISTOGRAMS
421 // --- Filled on a Particle-by-Particle basis
422 //------------------------------------------------
423    f3dHistPrimAnalysisPtVsYVsMultLambda(0),
424    f3dHistPrimAnalysisPtVsYVsMultAntiLambda(0),
425    f3dHistPrimAnalysisPtVsYVsMultK0Short(0),
426    f3dHistPrimAnalysisPtVsYCMSVsMultLambda(0),
427    f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda(0),
428    f3dHistPrimAnalysisPtVsYCMSVsMultK0Short(0),
429    f3dHistPrimRawPtVsYVsMultLambda(0),
430    f3dHistPrimRawPtVsYVsMultAntiLambda(0),
431    f3dHistPrimRawPtVsYVsMultK0Short(0),
432    f3dHistPrimRawPtVsYCMSVsMultLambda(0),
433    f3dHistPrimRawPtVsYCMSVsMultAntiLambda(0),
434    f3dHistPrimRawPtVsYCMSVsMultK0Short(0),
435    f3dHistPrimRawPtVsYVsMultNonInjLambda(0),
436    f3dHistPrimRawPtVsYVsMultNonInjAntiLambda(0),
437    f3dHistPrimRawPtVsYVsMultNonInjK0Short(0),
438    f3dHistPrimRawPtVsYVsMultMCLambda(0),
439    f3dHistPrimRawPtVsYVsMultMCAntiLambda(0),
440    f3dHistPrimRawPtVsYVsMultMCK0Short(0),
441    f3dHistPrimRawPtVsYVsVertexZLambda(0),
442    f3dHistPrimRawPtVsYVsVertexZAntiLambda(0),
443    f3dHistPrimRawPtVsYVsVertexZK0Short(0),
444    f3dHistPrimCloseToPVPtVsYVsMultLambda(0),
445    f3dHistPrimCloseToPVPtVsYVsMultAntiLambda(0),
446    f3dHistPrimCloseToPVPtVsYVsMultK0Short(0),
447    f3dHistPrimRawPtVsYVsDecayLengthLambda(0),
448    f3dHistPrimRawPtVsYVsDecayLengthAntiLambda(0),
449    f3dHistPrimRawPtVsYVsDecayLengthK0Short(0),
450    f3dHistGenPtVsYVsMultXiMinus(0),
451    f3dHistGenPtVsYVsMultXiPlus(0),
452    f3dHistGenPtVsYVsMultOmegaMinus(0),
453    f3dHistGenPtVsYVsMultOmegaPlus(0),
454    f3dHistGenSelectedPtVsYVsMultXiMinus(0),
455    f3dHistGenSelectedPtVsYVsMultXiPlus(0),
456    f3dHistGenSelectedPtVsYVsMultOmegaMinus(0),
457    f3dHistGenSelectedPtVsYVsMultOmegaPlus(0),
458    f3dHistGenPtVsYCMSVsMultXiMinus(0),
459    f3dHistGenPtVsYCMSVsMultXiPlus(0),
460    f3dHistGenPtVsYCMSVsMultOmegaMinus(0),
461    f3dHistGenPtVsYCMSVsMultOmegaPlus(0),
462    f3dHistGenSelectedPtVsYCMSVsMultXiMinus(0),
463    f3dHistGenSelectedPtVsYCMSVsMultXiPlus(0),
464    f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus(0),
465    f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus(0),
466    fHistPVx(0),
467    fHistPVy(0),
468    fHistPVz(0),
469    fHistPVxAnalysis(0),
470    fHistPVyAnalysis(0),
471    fHistPVzAnalysis(0),
472    fHistPVxAnalysisHasHighPtLambda(0),
473    fHistPVyAnalysisHasHighPtLambda(0),
474    fHistPVzAnalysisHasHighPtLambda(0),
475    fHistSwappedV0Counter(0)
476 {
477    // Constructor
478    // Output slot #0 writes into a TList container (Cascade)
479    DefineOutput(1, TList::Class());
480    DefineOutput(2, TTree::Class());
481 }
482
483
484 AliAnalysisTaskExtractPerformanceV0::~AliAnalysisTaskExtractPerformanceV0()
485 {
486 //------------------------------------------------
487 // DESTRUCTOR
488 //------------------------------------------------
489
490    if (fListHistV0){
491       delete fListHistV0;
492       fListHistV0 = 0x0;
493    }
494    if (fTree){
495       delete fTree;
496       fTree = 0x0;
497    }
498     //cleanup esd track cuts object too...
499    if (fESDtrackCuts){
500     delete fESDtrackCuts;
501     fESDtrackCuts = 0x0; 
502   }
503 }
504
505 //________________________________________________________________________
506 void AliAnalysisTaskExtractPerformanceV0::UserCreateOutputObjects()
507 {
508
509    OpenFile(2); 
510    // Called once
511
512 //------------------------------------------------
513
514    fTree = new TTree("fTree","V0Candidates");
515
516 //------------------------------------------------
517 // fTree Branch definitions - V0 Tree
518 //------------------------------------------------
519
520 //-----------BASIC-INFO---------------------------
521 /* 1*/   fTree->Branch("fTreeVariablePrimaryStatus",&fTreeVariablePrimaryStatus,"fTreeVariablePrimaryStatus/I");        
522 /* 2*/   fTree->Branch("fTreeVariablePrimaryStatusMother",&fTreeVariablePrimaryStatusMother,"fTreeVariablePrimaryStatusMother/I");      
523 /* 2*/   fTree->Branch("fTreeVariableChi2V0",&fTreeVariableChi2V0,"Chi2V0/F");
524 /* 3*/   fTree->Branch("fTreeVariableDcaV0Daughters",&fTreeVariableDcaV0Daughters,"fTreeVariableDcaV0Daughters/F");
525 /* 4*/   fTree->Branch("fTreeVariableDcaPosToPrimVertex",&fTreeVariableDcaPosToPrimVertex,"fTreeVariableDcaPosToPrimVertex/F");
526 /* 5*/   fTree->Branch("fTreeVariableDcaNegToPrimVertex",&fTreeVariableDcaNegToPrimVertex,"fTreeVariableDcaNegToPrimVertex/F");
527 /* 6*/   fTree->Branch("fTreeVariableV0Radius",&fTreeVariableV0Radius,"fTreeVariableV0Radius/F");
528 /* 7*/   fTree->Branch("fTreeVariablePt",&fTreeVariablePt,"fTreeVariablePt/F");
529 /* 7*/   fTree->Branch("fTreeVariablePtMC",&fTreeVariablePtMC,"fTreeVariablePtMC/F");
530 /* 8*/   fTree->Branch("fTreeVariableRapK0Short",&fTreeVariableRapK0Short,"fTreeVariableRapK0Short/F");
531 /* 9*/   fTree->Branch("fTreeVariableRapLambda",&fTreeVariableRapLambda,"fTreeVariableRapLambda/F");
532 /*10*/   fTree->Branch("fTreeVariableRapMC",&fTreeVariableRapMC,"fTreeVariableRapMC/F");
533 /*11*/   fTree->Branch("fTreeVariableInvMassK0s",&fTreeVariableInvMassK0s,"fTreeVariableInvMassK0s/F");
534 /*12*/   fTree->Branch("fTreeVariableInvMassLambda",&fTreeVariableInvMassLambda,"fTreeVariableInvMassLambda/F");
535 /*13*/   fTree->Branch("fTreeVariableInvMassAntiLambda",&fTreeVariableInvMassAntiLambda,"fTreeVariableInvMassAntiLambda/F");
536 /*14*/   fTree->Branch("fTreeVariableAlphaV0",&fTreeVariableAlphaV0,"fTreeVariableAlphaV0/F");
537 /*15*/   fTree->Branch("fTreeVariablePtArmV0",&fTreeVariablePtArmV0,"fTreeVariablePtArmV0/F");
538 /*16*/   fTree->Branch("fTreeVariableNegTransvMomentum",&fTreeVariableNegTransvMomentum,"fTreeVariableNegTransvMomentum/F");
539 /*17*/   fTree->Branch("fTreeVariablePosTransvMomentum",&fTreeVariablePosTransvMomentum,"fTreeVariablePosTransvMomentum/F");
540 /*18*/   fTree->Branch("fTreeVariableNegTransvMomentumMC",&fTreeVariableNegTransvMomentumMC,"fTreeVariableNegTransvMomentumMC/F");
541 /*19*/   fTree->Branch("fTreeVariablePosTransvMomentumMC",&fTreeVariablePosTransvMomentumMC,"fTreeVariablePosTransvMomentumMC/F");
542 /*20*/   fTree->Branch("fTreeVariableLeastNbrCrossedRows",&fTreeVariableLeastNbrCrossedRows,"fTreeVariableLeastNbrCrossedRows/I");
543 /*21*/   fTree->Branch("fTreeVariableLeastRatioCrossedRowsOverFindable",&fTreeVariableLeastRatioCrossedRowsOverFindable,"fTreeVariableLeastRatioCrossedRowsOverFindable/F");
544 /*22*/   fTree->Branch("fTreeVariablePID",&fTreeVariablePID,"fTreeVariablePID/I");
545 /*23*/   fTree->Branch("fTreeVariablePIDPositive",&fTreeVariablePIDPositive,"fTreeVariablePIDPositive/I");
546 /*24*/   fTree->Branch("fTreeVariablePIDNegative",&fTreeVariablePIDNegative,"fTreeVariablePIDNegative/I");
547 /*25*/   fTree->Branch("fTreeVariablePIDMother",&fTreeVariablePIDMother,"fTreeVariablePIDMother/I");
548 /*26*/   fTree->Branch("fTreeVariablePtXiMother",&fTreeVariablePtMother,"fTreeVariablePtMother/F");
549 /*27*/   fTree->Branch("fTreeVariableV0CosineOfPointingAngle",&fTreeVariableV0CosineOfPointingAngle,"fTreeVariableV0CosineOfPointingAngle/F");
550 //-----------MULTIPLICITY-INFO--------------------
551 /*28*/   fTree->Branch("fTreeVariableMultiplicity",&fTreeVariableMultiplicity,"fTreeVariableMultiplicity/I");
552 /*28*/   fTree->Branch("fTreeVariableMultiplicityMC",&fTreeVariableMultiplicityMC,"fTreeVariableMultiplicityMC/I");
553 //------------------------------------------------
554 /*29*/   fTree->Branch("fTreeVariableDistOverTotMom",&fTreeVariableDistOverTotMom,"fTreeVariableDistOverTotMom/F");
555 /*30*/   fTree->Branch("fTreeVariableNSigmasPosProton",&fTreeVariableNSigmasPosProton,"fTreeVariableNSigmasPosProton/F");
556 /*31*/   fTree->Branch("fTreeVariableNSigmasPosPion",&fTreeVariableNSigmasPosPion,"fTreeVariableNSigmasPosPion/F");
557 /*32*/   fTree->Branch("fTreeVariableNSigmasNegProton",&fTreeVariableNSigmasNegProton,"fTreeVariableNSigmasNegProton/F");
558 /*33*/   fTree->Branch("fTreeVariableNSigmasNegPion",&fTreeVariableNSigmasNegPion,"fTreeVariableNSigmasNegPion/F");
559 //------------------------------------------------
560 /*34*/   fTree->Branch("fTreeVariableNegEta",&fTreeVariableNegEta,"fTreeVariableNegEta/F");
561 /*35*/   fTree->Branch("fTreeVariablePosEta",&fTreeVariablePosEta,"fTreeVariablePosEta/F");
562 /*36*/   fTree->Branch("fTreeVariableV0CreationRadius",&fTreeVariableV0CreationRadius,"fTreeVariableV0CreationRadius/F");
563   
564     if ( fkLightWeight == kFALSE ){ // these are debugging branches!
565 /*37*/   fTree->Branch("fTreeVariableIndexStatus",&fTreeVariableIndexStatus,"fTreeVariableIndexStatus/I");
566 /*38*/   fTree->Branch("fTreeVariableIndexStatusMother",&fTreeVariableIndexStatusMother,"fTreeVariableIndexStatusMother/I");
567     }
568
569 /*39*/   fTree->Branch("fTreeVariableRunNumber",&fTreeVariableRunNumber,"fTreeVariableRunNumber/I");
570 /*40*/   fTree->Branch("fTreeVariableEventNumber",&fTreeVariableEventNumber,"fTreeVariableEventNumber/l");
571
572     if ( fkLightWeight == kFALSE ){ // these are debugging branches!
573 /*34*/   fTree->Branch("fTreeVariableVertexZ",&fTreeVariableVertexZ,"fTreeVariableVertexZ/F");
574     }
575       
576     if ( fkLightWeight == kFALSE ){ // these are debugging branches!
577 //-----------FOR CTAU DEBUGGING: Full Phase Space + Decay Position Info 
578         fTree->Branch("fTreeVariableV0x",&fTreeVariableV0x,"fTreeVariableV0x/F");
579         fTree->Branch("fTreeVariableV0y",&fTreeVariableV0y,"fTreeVariableV0y/F");
580         fTree->Branch("fTreeVariableV0z",&fTreeVariableV0z,"fTreeVariableV0z/F");
581
582         fTree->Branch("fTreeVariableV0Px",&fTreeVariableV0Px,"fTreeVariableV0Px/F");
583         fTree->Branch("fTreeVariableV0Py",&fTreeVariableV0Py,"fTreeVariableV0Py/F");
584         fTree->Branch("fTreeVariableV0Pz",&fTreeVariableV0Pz,"fTreeVariableV0Pz/F");
585
586 //-----------FOR CTAU DEBUGGING: Full Phase Space + Decay Position Info, perfect info from MC
587         fTree->Branch("fTreeVariableMCV0x",&fTreeVariableMCV0x,"fTreeVariableMCV0x/F");
588         fTree->Branch("fTreeVariableMCV0y",&fTreeVariableMCV0y,"fTreeVariableMCV0y/F");
589         fTree->Branch("fTreeVariableMCV0z",&fTreeVariableMCV0z,"fTreeVariableMCV0z/F");
590
591         fTree->Branch("fTreeVariableMCV0Px",&fTreeVariableMCV0Px,"fTreeVariableMCV0Px/F");
592         fTree->Branch("fTreeVariableMCV0Py",&fTreeVariableMCV0Py,"fTreeVariableMCV0Py/F");
593         fTree->Branch("fTreeVariableMCV0Pz",&fTreeVariableMCV0Pz,"fTreeVariableMCV0Pz/F");
594
595 //-----------FOR CTAU DEBUGGING: Primary vertex info 
596         fTree->Branch("fTreeVariablePVx",&fTreeVariablePVx,"fTreeVariablePVx/F");
597         fTree->Branch("fTreeVariablePVy",&fTreeVariablePVy,"fTreeVariablePVy/F");
598         fTree->Branch("fTreeVariablePVz",&fTreeVariablePVz,"fTreeVariablePVz/F");
599
600         fTree->Branch("fTreeVariableMCPVx",&fTreeVariableMCPVx,"fTreeVariableMCPVx/F");
601         fTree->Branch("fTreeVariableMCPVy",&fTreeVariableMCPVy,"fTreeVariableMCPVy/F");
602         fTree->Branch("fTreeVariableMCPVz",&fTreeVariableMCPVz,"fTreeVariableMCPVz/F");
603     }
604
605         fTree->Branch("fTreeVariableIsNonInjected",&fTreeVariableIsNonInjected,"fTreeVariableIsNonInjected/O"); //O for bOOlean...
606   
607   if ( fkLightWeight == kFALSE ){ // these are debugging branches!
608         fTree->Branch("fTreeVariableNegTrackStatus",&fTreeVariableNegTrackStatus,"fTreeVariableNegTrackStatus/l");
609         fTree->Branch("fTreeVariablePosTrackStatus",&fTreeVariablePosTrackStatus,"fTreeVariablePosTrackStatus/l");
610   }
611   
612 //------------------------------------------------
613 // Particle Identification Setup
614 //------------------------------------------------
615
616    AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
617    AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
618    fPIDResponse = inputHandler->GetPIDResponse();
619
620   // Multiplicity 
621
622     if(! fESDtrackCuts ){
623           fESDtrackCuts = new AliESDtrackCuts();
624     }
625
626 //------------------------------------------------
627 // V0 Multiplicity Histograms
628 //------------------------------------------------
629
630    // Create histograms
631    OpenFile(1);
632    fListHistV0 = new TList();
633    fListHistV0->SetOwner();  // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
634
635
636    if(! fHistV0MultiplicityBeforeTrigSel) {
637       fHistV0MultiplicityBeforeTrigSel = new TH1F("fHistV0MultiplicityBeforeTrigSel", 
638          "V0s per event (before Trig. Sel.);Nbr of V0s/Evt;Events", 
639          25, 0, 25);
640       fListHistV0->Add(fHistV0MultiplicityBeforeTrigSel);
641    }
642            
643    if(! fHistV0MultiplicityForTrigEvt) {
644       fHistV0MultiplicityForTrigEvt = new TH1F("fHistV0MultiplicityForTrigEvt", 
645          "V0s per event (for triggered evt);Nbr of V0s/Evt;Events", 
646          25, 0, 25);
647       fListHistV0->Add(fHistV0MultiplicityForTrigEvt);
648    }
649
650    if(! fHistV0MultiplicityForSelEvt) {
651       fHistV0MultiplicityForSelEvt = new TH1F("fHistV0MultiplicityForSelEvt", 
652          "V0s per event;Nbr of V0s/Evt;Events", 
653          25, 0, 25);
654       fListHistV0->Add(fHistV0MultiplicityForSelEvt);
655    }
656
657    if(! fHistV0MultiplicityForSelEvtNoTPCOnly) {
658       fHistV0MultiplicityForSelEvtNoTPCOnly = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnly", 
659          "V0s per event;Nbr of V0s/Evt;Events", 
660          25, 0, 25);
661       fListHistV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnly);
662    }
663    if(! fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup) {
664       fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup = new TH1F("fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup", 
665          "V0s per event;Nbr of V0s/Evt;Events", 
666          25, 0, 25);
667       fListHistV0->Add(fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup);
668    }
669
670 //------------------------------------------------
671 // Track Multiplicity Histograms
672 //------------------------------------------------
673
674    if(! fHistMultiplicityBeforeTrigSel) {
675       fHistMultiplicityBeforeTrigSel = new TH1F("fHistMultiplicityBeforeTrigSel", 
676          "Tracks per event;Nbr of Tracks;Events", 
677          200, 0, 200);          
678       fListHistV0->Add(fHistMultiplicityBeforeTrigSel);
679    }
680    if(! fHistMultiplicityForTrigEvt) {
681       fHistMultiplicityForTrigEvt = new TH1F("fHistMultiplicityForTrigEvt", 
682          "Tracks per event;Nbr of Tracks;Events", 
683          200, 0, 200);          
684       fListHistV0->Add(fHistMultiplicityForTrigEvt);
685    }
686    if(! fHistMultiplicity) {
687       fHistMultiplicity = new TH1F("fHistMultiplicity", 
688          "Tracks per event;Nbr of Tracks;Events", 
689          200, 0, 200);          
690       fListHistV0->Add(fHistMultiplicity);
691    }
692    if(! fHistMultiplicityNoTPCOnly) {
693       fHistMultiplicityNoTPCOnly = new TH1F("fHistMultiplicityNoTPCOnly", 
694          "Tracks per event;Nbr of Tracks;Events", 
695          200, 0, 200);          
696       fListHistV0->Add(fHistMultiplicityNoTPCOnly);
697    }
698    if(! fHistMultiplicityNoTPCOnlyNoPileup) {
699       fHistMultiplicityNoTPCOnlyNoPileup = new TH1F("fHistMultiplicityNoTPCOnlyNoPileup", 
700          "Tracks per event;Nbr of Tracks;Events", 
701          200, 0, 200);          
702       fListHistV0->Add(fHistMultiplicityNoTPCOnlyNoPileup);
703    }
704
705   //Raw Data for J/Psi paper Technique
706         //TH2F    *f2dHistMultiplicityVsTrueBeforeTrigSel;              //! multiplicity distribution    
707         //TH2F    *f2dHistMultiplicityVsTrueForTrigEvt;                         //! multiplicity distribution
708         //TH2F    *f2dHistMultiplicityVsTrue;                                                   //! multiplicity distribution
709         //TH2F    *f2dHistMultiplicityVsTrueNoTPCOnly;                          //! multiplicity distribution
710         //TH2F    *f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup;                  //! multiplicity distribution
711
712    if(! f2dHistMultiplicityVsTrueBeforeTrigSel) {
713       f2dHistMultiplicityVsTrueBeforeTrigSel = new TH2F("f2dHistMultiplicityVsTrueBeforeTrigSel", 
714          "Tracks per event", 200, 0, 200, 200, 0, 200);                 
715       fListHistV0->Add(f2dHistMultiplicityVsTrueBeforeTrigSel);
716    }
717    if(! f2dHistMultiplicityVsTrueForTrigEvt) {
718       f2dHistMultiplicityVsTrueForTrigEvt = new TH2F("f2dHistMultiplicityVsTrueForTrigEvt", 
719          "Tracks per event", 200, 0, 200, 200, 0, 200);                 
720       fListHistV0->Add(f2dHistMultiplicityVsTrueForTrigEvt);
721    }
722    if(! f2dHistMultiplicityVsTrue) {
723       f2dHistMultiplicityVsTrue = new TH2F("f2dHistMultiplicityVsTrue", 
724          "Tracks per event", 200, 0, 200, 200, 0, 200);                 
725       fListHistV0->Add(f2dHistMultiplicityVsTrue);
726    }
727    if(! f2dHistMultiplicityVsTrueNoTPCOnly) {
728       f2dHistMultiplicityVsTrueNoTPCOnly = new TH2F("f2dHistMultiplicityVsTrueNoTPCOnly", 
729          "Tracks per event", 200, 0, 200, 200, 0, 200);                 
730       fListHistV0->Add(f2dHistMultiplicityVsTrueNoTPCOnly);
731    }
732    if(! f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup) {
733       f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup = new TH2F("f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup", 
734          "Tracks per event", 200, 0, 200, 200, 0, 200);                 
735       fListHistV0->Add(f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup);
736    }
737
738
739   //Raw Data for Vertex Z position estimator change
740         //TH2F    *f2dHistMultiplicityVsVertexZBeforeTrigSel;           //! multiplicity distribution    
741         //TH2F    *f2dHistMultiplicityVsVertexZForTrigEvt;                      //! multiplicity distribution
742         //TH2F    *f2dHistMultiplicityVsVertexZ;                                                //! multiplicity distribution
743         //TH2F    *f2dHistMultiplicityVsVertexZNoTPCOnly;                               //! multiplicity distribution
744         //TH2F    *f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup;                       //! multiplicity distribution
745
746    if(! f2dHistMultiplicityVsVertexZBeforeTrigSel) {
747       f2dHistMultiplicityVsVertexZBeforeTrigSel = new TH2F("f2dHistMultiplicityVsVertexZBeforeTrigSel", 
748          "Tracks per event", 200, 0, 200,400, -20, 20);                 
749       fListHistV0->Add(f2dHistMultiplicityVsVertexZBeforeTrigSel);
750    }
751    if(! f2dHistMultiplicityVsVertexZForTrigEvt) {
752       f2dHistMultiplicityVsVertexZForTrigEvt = new TH2F("f2dHistMultiplicityVsVertexZForTrigEvt", 
753          "Tracks per event", 200, 0, 200, 400, -20, 20);                
754       fListHistV0->Add(f2dHistMultiplicityVsVertexZForTrigEvt);
755    }
756    if(! f2dHistMultiplicityVsVertexZ) {
757       f2dHistMultiplicityVsVertexZ = new TH2F("f2dHistMultiplicityVsVertexZ", 
758          "Tracks per event", 200, 0, 200, 400, -20, 20);                
759       fListHistV0->Add(f2dHistMultiplicityVsVertexZ);
760    }
761    if(! f2dHistMultiplicityVsVertexZNoTPCOnly) {
762       f2dHistMultiplicityVsVertexZNoTPCOnly = new TH2F("f2dHistMultiplicityVsVertexZNoTPCOnly", 
763          "Tracks per event", 200, 0, 200, 400, -20, 20);                
764       fListHistV0->Add(f2dHistMultiplicityVsVertexZNoTPCOnly);
765    }
766    if(! f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup) {
767       f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup = new TH2F("f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup", 
768          "Tracks per event", 200, 0, 200, 400, -20, 20);                
769       fListHistV0->Add(f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup);
770    }
771
772 //Generated PVz test
773
774 //   TH1F      *fHistGenVertexZBeforeTrigSel;     //! multiplicity distribution      
775 //   TH1F      *fHistGenVertexZForTrigEvt;        //! multiplicity distribution
776 //   TH1F      *fHistGenVertexZ;                  //! multiplicity distribution
777 //   TH1F      *fHistGenVertexZNoTPCOnly;         //! multiplicity distribution
778 //   TH1F      *fHistGenVertexZNoTPCOnlyNoPileup; //! multiplicity distribution
779
780    if(! fHistGenVertexZBeforeTrigSel) {
781          fHistGenVertexZBeforeTrigSel = new TH1F("fHistGenVertexZBeforeTrigSel", 
782             "PV z position;Nbr of Evts;z", 
783             400, -20, 20);       
784       fListHistV0->Add(fHistGenVertexZBeforeTrigSel);
785    }
786    if(! fHistGenVertexZForTrigEvt) {
787          fHistGenVertexZForTrigEvt = new TH1F("fHistGenVertexZForTrigEvt", 
788             "PV z position;Nbr of Evts;z", 
789             400, -20, 20);       
790       fListHistV0->Add(fHistGenVertexZForTrigEvt);
791    }
792    if(! fHistGenVertexZ) {
793          fHistGenVertexZ = new TH1F("fHistGenVertexZ", 
794             "PV z position;Nbr of Evts;z", 
795             400, -20, 20);       
796       fListHistV0->Add(fHistGenVertexZ);
797    }
798    if(! fHistGenVertexZNoTPCOnly) {
799          fHistGenVertexZNoTPCOnly = new TH1F("fHistGenVertexZNoTPCOnly", 
800             "PV z position;Nbr of Evts;z", 
801             400, -20, 20);       
802       fListHistV0->Add(fHistGenVertexZNoTPCOnly);
803    }
804    if(! fHistGenVertexZNoTPCOnlyNoPileup) {
805          fHistGenVertexZNoTPCOnlyNoPileup = new TH1F("fHistGenVertexZNoTPCOnlyNoPileup", 
806             "PV z position;Nbr of Evts;z", 
807             400, -20, 20);       
808       fListHistV0->Add(fHistGenVertexZNoTPCOnlyNoPileup);
809    }
810
811
812 //------------------------------------------------
813 // Generated Particle Histograms
814 //------------------------------------------------
815
816    Int_t lCustomNBins = 200; 
817    Double_t lCustomPtUpperLimit = 20; 
818    Int_t lCustomNBinsMultiplicity = 100;
819
820 //----------------------------------
821 // Raw Generated (Pre-physics-selection)
822 //----------------------------------
823
824 //--- 3D Histo (Pt, Y, Multiplicity)  
825
826    if(! f3dHistPrimRawPtVsYVsMultLambda) {
827       f3dHistPrimRawPtVsYVsMultLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
828       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultLambda);
829    }
830    if(! f3dHistPrimRawPtVsYVsMultAntiLambda) {
831       f3dHistPrimRawPtVsYVsMultAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
832       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultAntiLambda);
833    }
834    if(! f3dHistPrimRawPtVsYVsMultK0Short) {
835       f3dHistPrimRawPtVsYVsMultK0Short = new TH3F( "f3dHistPrimRawPtVsYVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
836       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultK0Short);
837    }
838
839    if(! f3dHistPrimRawPtVsYCMSVsMultLambda) {
840       f3dHistPrimRawPtVsYCMSVsMultLambda = new TH3F( "f3dHistPrimRawPtVsYCMSVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
841       fListHistV0->Add(f3dHistPrimRawPtVsYCMSVsMultLambda);
842    }
843    if(! f3dHistPrimRawPtVsYCMSVsMultAntiLambda) {
844       f3dHistPrimRawPtVsYCMSVsMultAntiLambda = new TH3F( "f3dHistPrimRawPtVsYCMSVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
845       fListHistV0->Add(f3dHistPrimRawPtVsYCMSVsMultAntiLambda);
846    }
847    if(! f3dHistPrimRawPtVsYCMSVsMultK0Short) {
848       f3dHistPrimRawPtVsYCMSVsMultK0Short = new TH3F( "f3dHistPrimRawPtVsYCMSVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
849       fListHistV0->Add(f3dHistPrimRawPtVsYCMSVsMultK0Short);
850    }
851
852 //---> Non-injected particles
853
854    if(! f3dHistPrimRawPtVsYVsMultNonInjLambda) {
855       f3dHistPrimRawPtVsYVsMultNonInjLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultNonInjLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
856       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultNonInjLambda);
857    }
858    if(! f3dHistPrimRawPtVsYVsMultNonInjAntiLambda) {
859       f3dHistPrimRawPtVsYVsMultNonInjAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultNonInjAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
860       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultNonInjAntiLambda);
861    }
862    if(! f3dHistPrimRawPtVsYVsMultNonInjK0Short) {
863       f3dHistPrimRawPtVsYVsMultNonInjK0Short = new TH3F( "f3dHistPrimRawPtVsYVsMultNonInjK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
864       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultNonInjK0Short);
865    }
866
867 //--- 3D Histo (Pt, Y, MultiplicityMC)  
868
869    if(! f3dHistPrimRawPtVsYVsMultMCLambda) {
870       f3dHistPrimRawPtVsYVsMultMCLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultMCLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; MultMC", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
871       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultMCLambda);
872    }
873    if(! f3dHistPrimRawPtVsYVsMultMCAntiLambda) {
874       f3dHistPrimRawPtVsYVsMultMCAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsMultMCAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; MultMC", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
875       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultMCAntiLambda);
876    }
877    if(! f3dHistPrimRawPtVsYVsMultMCK0Short) {
878       f3dHistPrimRawPtVsYVsMultMCK0Short = new TH3F( "f3dHistPrimRawPtVsYVsMultMCK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; MultMC", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
879       fListHistV0->Add(f3dHistPrimRawPtVsYVsMultMCK0Short);
880    }
881
882 //--- 3D Histo (Pt, Y, VertexZ)  
883
884    if(! f3dHistPrimRawPtVsYVsVertexZLambda) {
885       f3dHistPrimRawPtVsYVsVertexZLambda = new TH3F( "f3dHistPrimRawPtVsYVsVertexZLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs VertexZiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; VertexZ", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,40,-10,10);
886       fListHistV0->Add(f3dHistPrimRawPtVsYVsVertexZLambda);
887    }
888    if(! f3dHistPrimRawPtVsYVsVertexZAntiLambda) {
889       f3dHistPrimRawPtVsYVsVertexZAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsVertexZAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs VertexZiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; VertexZ", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,40,-10,10);
890       fListHistV0->Add(f3dHistPrimRawPtVsYVsVertexZAntiLambda);
891    }
892    if(! f3dHistPrimRawPtVsYVsVertexZK0Short) {
893       f3dHistPrimRawPtVsYVsVertexZK0Short = new TH3F( "f3dHistPrimRawPtVsYVsVertexZK0Short", "Pt_{K0S} Vs Y_{K0S} Vs VertexZiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; VertexZ", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,40,-10,10);
894       fListHistV0->Add(f3dHistPrimRawPtVsYVsVertexZK0Short);
895    }
896
897 //--- 3D Histo (Pt, Y, Multiplicity), close to PV criterion
898
899    if(! f3dHistPrimCloseToPVPtVsYVsMultLambda) {
900       f3dHistPrimCloseToPVPtVsYVsMultLambda = new TH3F( "f3dHistPrimCloseToPVPtVsYVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
901       fListHistV0->Add(f3dHistPrimCloseToPVPtVsYVsMultLambda);
902    }
903    if(! f3dHistPrimCloseToPVPtVsYVsMultAntiLambda) {
904       f3dHistPrimCloseToPVPtVsYVsMultAntiLambda = new TH3F( "f3dHistPrimCloseToPVPtVsYVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
905       fListHistV0->Add(f3dHistPrimCloseToPVPtVsYVsMultAntiLambda);
906    }
907    if(! f3dHistPrimCloseToPVPtVsYVsMultK0Short) {
908       f3dHistPrimCloseToPVPtVsYVsMultK0Short = new TH3F( "f3dHistPrimCloseToPVPtVsYVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
909       fListHistV0->Add(f3dHistPrimCloseToPVPtVsYVsMultK0Short);
910    }
911
912
913 //--- 3D Histo (Pt, Y, Proper Decay Length)
914
915    if(! f3dHistPrimRawPtVsYVsDecayLengthLambda) {
916       f3dHistPrimRawPtVsYVsDecayLengthLambda = new TH3F( "f3dHistPrimRawPtVsYVsDecayLengthLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs DecayLength; Pt_{lambda} (GeV/c); Y_{#Lambda} ; DecayLength", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,200,0,50);
917       fListHistV0->Add(f3dHistPrimRawPtVsYVsDecayLengthLambda);
918    }
919    if(! f3dHistPrimRawPtVsYVsDecayLengthAntiLambda) {
920       f3dHistPrimRawPtVsYVsDecayLengthAntiLambda = new TH3F( "f3dHistPrimRawPtVsYVsDecayLengthAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs DecayLength; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; DecayLength", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,200,0,50);
921       fListHistV0->Add(f3dHistPrimRawPtVsYVsDecayLengthAntiLambda);
922    }
923    if(! f3dHistPrimRawPtVsYVsDecayLengthK0Short) {
924       f3dHistPrimRawPtVsYVsDecayLengthK0Short = new TH3F( "f3dHistPrimRawPtVsYVsDecayLengthK0Short", "Pt_{K0S} Vs Y_{K0S} Vs DecayLength; Pt_{K0S} (GeV/c); Y_{K0S} ; DecayLength", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,200,0,50);
925       fListHistV0->Add(f3dHistPrimRawPtVsYVsDecayLengthK0Short);
926    }
927
928 //--------------------------------------------------------------------------------------
929 //--- 3D Histo (Pt, Y, Multiplicity) for generated XiMinus/Plus, all generated
930
931    if(! f3dHistGenPtVsYVsMultXiMinus) {
932       f3dHistGenPtVsYVsMultXiMinus = new TH3F( "f3dHistGenPtVsYVsMultXiMinus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
933       fListHistV0->Add(f3dHistGenPtVsYVsMultXiMinus);
934    }
935    if(! f3dHistGenPtVsYVsMultXiPlus) {
936       f3dHistGenPtVsYVsMultXiPlus = new TH3F( "f3dHistGenPtVsYVsMultXiPlus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
937       fListHistV0->Add(f3dHistGenPtVsYVsMultXiPlus);
938    }
939 //--- 3D Histo (Pt, Y, Multiplicity) for generated OmegaMinus/Plus
940
941    if(! f3dHistGenPtVsYVsMultOmegaMinus) {
942       f3dHistGenPtVsYVsMultOmegaMinus = new TH3F( "f3dHistGenPtVsYVsMultOmegaMinus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
943       fListHistV0->Add(f3dHistGenPtVsYVsMultOmegaMinus);
944    }
945    if(! f3dHistGenPtVsYVsMultOmegaPlus) {
946       f3dHistGenPtVsYVsMultOmegaPlus = new TH3F( "f3dHistGenPtVsYVsMultOmegaPlus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
947       fListHistV0->Add(f3dHistGenPtVsYVsMultOmegaPlus);
948    }
949
950 //CASCADEs, Y CMS 
951
952    if(! f3dHistGenPtVsYCMSVsMultXiMinus) {
953       f3dHistGenPtVsYCMSVsMultXiMinus = new TH3F( "f3dHistGenPtVsYCMSVsMultXiMinus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
954       fListHistV0->Add(f3dHistGenPtVsYCMSVsMultXiMinus);
955    }
956    if(! f3dHistGenPtVsYCMSVsMultXiPlus) {
957       f3dHistGenPtVsYCMSVsMultXiPlus = new TH3F( "f3dHistGenPtVsYCMSVsMultXiPlus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
958       fListHistV0->Add(f3dHistGenPtVsYCMSVsMultXiPlus);
959    }
960 //--- 3D Histo (Pt, Y, Multiplicity) for generated OmegaMinus/Plus
961
962    if(! f3dHistGenPtVsYCMSVsMultOmegaMinus) {
963       f3dHistGenPtVsYCMSVsMultOmegaMinus = new TH3F( "f3dHistGenPtVsYCMSVsMultOmegaMinus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
964       fListHistV0->Add(f3dHistGenPtVsYCMSVsMultOmegaMinus);
965    }
966    if(! f3dHistGenPtVsYCMSVsMultOmegaPlus) {
967       f3dHistGenPtVsYCMSVsMultOmegaPlus = new TH3F( "f3dHistGenPtVsYCMSVsMultOmegaPlus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
968       fListHistV0->Add(f3dHistGenPtVsYCMSVsMultOmegaPlus);
969    }
970
971
972
973 //--------------------------------------------------------------------------------------
974 //--- 3D Histo (Pt, Y, Multiplicity) for generated XiMinus/Plus, at selected analysis evts
975
976    if(! f3dHistGenSelectedPtVsYVsMultXiMinus) {
977       f3dHistGenSelectedPtVsYVsMultXiMinus = new TH3F( "f3dHistGenSelectedPtVsYVsMultXiMinus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
978       fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultXiMinus);
979    }
980    if(! f3dHistGenSelectedPtVsYVsMultXiPlus) {
981       f3dHistGenSelectedPtVsYVsMultXiPlus = new TH3F( "f3dHistGenSelectedPtVsYVsMultXiPlus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
982       fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultXiPlus);
983    }
984 //--- 3D Histo (Pt, Y, Multiplicity) for generated OmegaMinus/Plus
985
986    if(! f3dHistGenSelectedPtVsYVsMultOmegaMinus) {
987       f3dHistGenSelectedPtVsYVsMultOmegaMinus = new TH3F( "f3dHistGenSelectedPtVsYVsMultOmegaMinus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
988       fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultOmegaMinus);
989    }
990    if(! f3dHistGenSelectedPtVsYVsMultOmegaPlus) {
991       f3dHistGenSelectedPtVsYVsMultOmegaPlus = new TH3F( "f3dHistGenSelectedPtVsYVsMultOmegaPlus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
992       fListHistV0->Add(f3dHistGenSelectedPtVsYVsMultOmegaPlus);
993    }
994
995 //CASCADES, analysis level, y CMS
996
997 //--------------------------------------------------------------------------------------
998 //--- 3D Histo (Pt, Y, Multiplicity) for generated XiMinus/Plus, at selected analysis evts
999
1000    if(! f3dHistGenSelectedPtVsYCMSVsMultXiMinus) {
1001       f3dHistGenSelectedPtVsYCMSVsMultXiMinus = new TH3F( "f3dHistGenSelectedPtVsYCMSVsMultXiMinus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1002       fListHistV0->Add(f3dHistGenSelectedPtVsYCMSVsMultXiMinus);
1003    }
1004    if(! f3dHistGenSelectedPtVsYCMSVsMultXiPlus) {
1005       f3dHistGenSelectedPtVsYCMSVsMultXiPlus = new TH3F( "f3dHistGenSelectedPtVsYCMSVsMultXiPlus", "Pt_{#Xi} Vs Y_{#Xi} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Xi} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1006       fListHistV0->Add(f3dHistGenSelectedPtVsYCMSVsMultXiPlus);
1007    }
1008 //--- 3D Histo (Pt, Y, Multiplicity) for generated OmegaMinus/Plus
1009
1010    if(! f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus) {
1011       f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus = new TH3F( "f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1012       fListHistV0->Add(f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus);
1013    }
1014    if(! f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus) {
1015       f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus = new TH3F( "f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus", "Pt_{#Omega} Vs Y_{#Omega} Vs Multiplicity; Pt_{cascade} (GeV/c); Y_{#Omega} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1016       fListHistV0->Add(f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus);
1017    }
1018
1019
1020 //----------------------------------
1021 // Histos at analysis level 
1022 //----------------------------------
1023
1024    if(! f3dHistPrimAnalysisPtVsYVsMultLambda) {
1025       f3dHistPrimAnalysisPtVsYVsMultLambda = new TH3F( "f3dHistPrimAnalysisPtVsYVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1026       fListHistV0->Add(f3dHistPrimAnalysisPtVsYVsMultLambda);
1027    }
1028    if(! f3dHistPrimAnalysisPtVsYVsMultAntiLambda) {
1029       f3dHistPrimAnalysisPtVsYVsMultAntiLambda = new TH3F( "f3dHistPrimAnalysisPtVsYVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1030       fListHistV0->Add(f3dHistPrimAnalysisPtVsYVsMultAntiLambda);
1031    }
1032    if(! f3dHistPrimAnalysisPtVsYVsMultK0Short) {
1033       f3dHistPrimAnalysisPtVsYVsMultK0Short = new TH3F( "f3dHistPrimAnalysisPtVsYVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1034       fListHistV0->Add(f3dHistPrimAnalysisPtVsYVsMultK0Short);
1035    }
1036
1037    if(! f3dHistPrimAnalysisPtVsYCMSVsMultLambda) {
1038       f3dHistPrimAnalysisPtVsYCMSVsMultLambda = new TH3F( "f3dHistPrimAnalysisPtVsYCMSVsMultLambda", "Pt_{lambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{lambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1039       fListHistV0->Add(f3dHistPrimAnalysisPtVsYCMSVsMultLambda);
1040    }
1041    if(! f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda) {
1042       f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda = new TH3F( "f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda", "Pt_{antilambda} Vs Y_{#Lambda} Vs Multiplicity; Pt_{antilambda} (GeV/c); Y_{#Lambda} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1043       fListHistV0->Add(f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda);
1044    }
1045    if(! f3dHistPrimAnalysisPtVsYCMSVsMultK0Short) {
1046       f3dHistPrimAnalysisPtVsYCMSVsMultK0Short = new TH3F( "f3dHistPrimAnalysisPtVsYCMSVsMultK0Short", "Pt_{K0S} Vs Y_{K0S} Vs Multiplicity; Pt_{K0S} (GeV/c); Y_{K0S} ; Mult", lCustomNBins, 0., lCustomPtUpperLimit, 48, -1.2,1.2,lCustomNBinsMultiplicity,0,lCustomNBinsMultiplicity);
1047       fListHistV0->Add(f3dHistPrimAnalysisPtVsYCMSVsMultK0Short);
1048    }
1049
1050 //----------------------------------
1051 // Primary Vertex Position Histos
1052 //----------------------------------
1053
1054    if(! fHistPVx) {
1055          fHistPVx = new TH1F("fHistPVx", 
1056             "PV x position;Nbr of Evts;x", 
1057             2000, -0.5, 0.5);       
1058       fListHistV0->Add(fHistPVx);
1059    }
1060    if(! fHistPVy) {
1061          fHistPVy = new TH1F("fHistPVy", 
1062             "PV y position;Nbr of Evts;y", 
1063             2000, -0.5, 0.5);       
1064       fListHistV0->Add(fHistPVy);
1065    }
1066    if(! fHistPVz) {
1067          fHistPVz = new TH1F("fHistPVz", 
1068             "PV z position;Nbr of Evts;z", 
1069             400, -20, 20);       
1070       fListHistV0->Add(fHistPVz);
1071    }
1072
1073    if(! fHistPVxAnalysis) {
1074          fHistPVxAnalysis = new TH1F("fHistPVxAnalysis", 
1075             "PV x position;Nbr of Evts;x", 
1076             2000, -0.5, 0.5);       
1077       fListHistV0->Add(fHistPVxAnalysis);
1078    }
1079    if(! fHistPVyAnalysis) {
1080          fHistPVyAnalysis = new TH1F("fHistPVyAnalysis", 
1081             "PV y position;Nbr of Evts;y", 
1082             2000, -0.5, 0.5);       
1083       fListHistV0->Add(fHistPVyAnalysis);
1084    }
1085    if(! fHistPVzAnalysis) {
1086          fHistPVzAnalysis = new TH1F("fHistPVzAnalysis", 
1087             "PV z position;Nbr of Evts;z", 
1088             400, -20, 20);       
1089       fListHistV0->Add(fHistPVzAnalysis);
1090    }
1091
1092    if(! fHistPVxAnalysisHasHighPtLambda) {
1093          fHistPVxAnalysisHasHighPtLambda = new TH1F("fHistPVxAnalysisHasHighPtLambda", 
1094             "PV x position;Nbr of Evts;x", 
1095             2000, -0.5, 0.5);       
1096       fListHistV0->Add(fHistPVxAnalysisHasHighPtLambda);
1097    }
1098    if(! fHistPVyAnalysisHasHighPtLambda) {
1099          fHistPVyAnalysisHasHighPtLambda = new TH1F("fHistPVyAnalysisHasHighPtLambda", 
1100             "PV y position;Nbr of Evts;y", 
1101             2000, -0.5, 0.5);       
1102       fListHistV0->Add(fHistPVyAnalysisHasHighPtLambda);
1103    }
1104    if(! fHistPVzAnalysisHasHighPtLambda) {
1105          fHistPVzAnalysisHasHighPtLambda = new TH1F("fHistPVzAnalysisHasHighPtLambda", 
1106             "PV z position;Nbr of Evts;z", 
1107             400, -20, 20);       
1108       fListHistV0->Add(fHistPVzAnalysisHasHighPtLambda);
1109    }
1110    if(! fHistSwappedV0Counter) {
1111       fHistSwappedV0Counter = new TH1F("fHistSwappedV0Counter", 
1112          "Swap or not histo;Swapped (1) or not (0); count", 
1113          2, 0, 2);              
1114       fListHistV0->Add(fHistSwappedV0Counter);
1115    }
1116
1117    //List of Histograms: Normal
1118    PostData(1, fListHistV0);
1119
1120    //TTree Object: Saved to base directory. Should cache to disk while saving. 
1121    //(Important to avoid excessive memory usage, particularly when merging)
1122    PostData(2, fTree);
1123
1124 }// end UserCreateOutputObjects
1125
1126
1127 //________________________________________________________________________
1128 void AliAnalysisTaskExtractPerformanceV0::UserExec(Option_t *) 
1129 {
1130   // Main loop
1131   // Called for each event
1132
1133    AliESDEvent *lESDevent = 0x0;
1134    AliMCEvent  *lMCevent  = 0x0; 
1135    AliStack    *lMCstack  = 0x0; 
1136
1137    Int_t    lNumberOfV0s                      = -1;
1138    Double_t lTrkgPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
1139    Double_t lBestPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
1140    Double_t lMagneticField                 = -10.;
1141    
1142   // Connect to the InputEvent   
1143   // After these lines, we should have an ESD/AOD event + the number of V0s in it.
1144
1145    // Appropriate for ESD analysis! 
1146       
1147    lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
1148    if (!lESDevent) {
1149       AliWarning("ERROR: lESDevent not available \n");
1150       return;
1151    }
1152         
1153    fTreeVariableRunNumber = lESDevent->GetRunNumber();
1154    fTreeVariableEventNumber =  
1155     ( ( ((ULong64_t)lESDevent->GetPeriodNumber() ) << 36 ) |
1156       ( ((ULong64_t)lESDevent->GetOrbitNumber () ) << 12 ) |
1157         ((ULong64_t)lESDevent->GetBunchCrossNumber() )  );
1158
1159    lMCevent = MCEvent();
1160    if (!lMCevent) {
1161       Printf("ERROR: Could not retrieve MC event \n");
1162       cout << "Name of the file with pb :" <<  fInputHandler->GetTree()->GetCurrentFile()->GetName() << endl;   
1163       return;
1164    }
1165
1166    lMCstack = lMCevent->Stack();
1167    if (!lMCstack) {
1168       Printf("ERROR: Could not retrieve MC stack \n");
1169       cout << "Name of the file with pb :" <<  fInputHandler->GetTree()->GetCurrentFile()->GetName() << endl;
1170       return;
1171    }
1172    TArrayF mcPrimaryVtx;
1173    AliGenEventHeader* mcHeader=lMCevent->GenEventHeader();
1174    if(!mcHeader) return;
1175    mcHeader->PrimaryVertex(mcPrimaryVtx);
1176         
1177 //------------------------------------------------
1178 // Multiplicity Information Acquistion
1179 //------------------------------------------------
1180
1181    //REVISED multiplicity estimator after 'multiplicity day' (2011)
1182    Int_t lMultiplicity = -100; 
1183
1184    //testing purposes
1185    if(fkIsNuclear == kFALSE) lMultiplicity = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.5);
1186
1187    //---> If this is a nuclear collision, then go nuclear on "multiplicity" variable...
1188    //---> Warning: Experimental
1189    if(fkIsNuclear == kTRUE){ 
1190       AliCentrality* centrality;
1191       centrality = lESDevent->GetCentrality();
1192       lMultiplicity = ( ( Int_t ) ( centrality->GetCentralityPercentile(   fCentralityEstimator.Data() ) ) );
1193       if (centrality->GetQuality()>1) {
1194         PostData(1, fListHistV0);
1195         PostData(2, fTree);
1196         return;
1197       }
1198    }
1199   
1200    //Set variable for filling tree afterwards!
1201    //---> pp case......: GetReferenceMultiplicity
1202    //---> Pb-Pb case...: Centrality by V0M
1203    fTreeVariableMultiplicity = lMultiplicity;
1204
1205    fHistV0MultiplicityBeforeTrigSel->Fill ( lESDevent->GetNumberOfV0s() );
1206    fHistMultiplicityBeforeTrigSel->Fill ( lMultiplicity );
1207         
1208 //------------------------------------------------
1209 // MC Information Acquistion
1210 //------------------------------------------------
1211
1212    Int_t iNumberOfPrimaries = -1;
1213    iNumberOfPrimaries = lMCstack->GetNprimary();
1214    if(iNumberOfPrimaries < 1) return; 
1215    Bool_t lHasHighPtLambda = kFALSE;
1216
1217 //------------------------------------------------
1218 // Variable Definition
1219 //------------------------------------------------
1220
1221    Int_t lNbMCPrimary        = 0;
1222
1223    Int_t lPdgcodeCurrentPart = 0;
1224    Double_t lRapCurrentPart  = 0;
1225    Double_t lPtCurrentPart   = 0;
1226   
1227    //Int_t lComeFromSigma      = 0;
1228
1229    // current mc particle 's mother
1230    //Int_t iCurrentMother  = 0;
1231    lNbMCPrimary = lMCstack->GetNprimary();
1232
1233 //------------------------------------------------
1234 // Pre-Physics Selection
1235 //------------------------------------------------
1236
1237 //----- Loop on primary Xi, Omega --------------------------------------------------------------
1238    for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < lNbMCPrimary; iCurrentLabelStack++) 
1239    {// This is the begining of the loop on primaries
1240       
1241       TParticle* lCurrentParticlePrimary = 0x0; 
1242       lCurrentParticlePrimary = lMCstack->Particle( iCurrentLabelStack );
1243       if(!lCurrentParticlePrimary){
1244          Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
1245          continue;
1246       }
1247       if ( TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3312 || TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3334 ) { 
1248          Double_t lRapXiMCPrimary = -100;
1249          if( (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) != 0 ) { 
1250            if ( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) !=0 ){
1251              lRapXiMCPrimary = 0.5*TMath::Log( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) );
1252            }
1253          }
1254
1255          //=================================================================================
1256          // Xi Histograms
1257          if( lCurrentParticlePrimary->GetPdgCode() == 3312 ){ 
1258             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1259             f3dHistGenPtVsYVsMultXiMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1260             f3dHistGenPtVsYCMSVsMultXiMinus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1261          }
1262          if( lCurrentParticlePrimary->GetPdgCode() == -3312 ){ 
1263             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1264             f3dHistGenPtVsYVsMultXiPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1265             f3dHistGenPtVsYCMSVsMultXiPlus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1266          }
1267          // Omega Histograms
1268          if( lCurrentParticlePrimary->GetPdgCode() == 3334 ){ 
1269             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1270             f3dHistGenPtVsYVsMultOmegaMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1271             f3dHistGenPtVsYCMSVsMultOmegaMinus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1272          }
1273          if( lCurrentParticlePrimary->GetPdgCode() == -3334 ){ 
1274             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1275             f3dHistGenPtVsYVsMultOmegaPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1276             f3dHistGenPtVsYCMSVsMultOmegaPlus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1277          }
1278       } 
1279    }
1280 //----- End Loop on primary Xi, Omega ----------------------------------------------------------
1281
1282 //--------- GENERATED NUMBER OF CHARGED PARTICLES
1283 // ---> Set Variables to Zero again
1284 // ---> Variable Definition
1285
1286   Long_t lNumberOfCharged = 0; 
1287
1288 //----- Loop on Stack ----------------------------------------------------------------
1289    for (Int_t iCurrentLabelStack = 0;  iCurrentLabelStack < (lMCstack->GetNtrack()); iCurrentLabelStack++) 
1290    {// This is the begining of the loop on tracks
1291       TParticle* particleOne = lMCstack->Particle(iCurrentLabelStack);
1292       if(!particleOne) continue;
1293       if(!particleOne->GetPDG()) continue;
1294       Double_t lThisCharge = particleOne->GetPDG()->Charge()/3.;
1295       if(TMath::Abs(lThisCharge)<0.001) continue;
1296       if(! (lMCstack->IsPhysicalPrimary(iCurrentLabelStack)) ) continue;
1297       
1298       Double_t gpt = particleOne -> Pt();
1299       Double_t geta = particleOne -> Eta(); 
1300
1301       if( TMath::Abs(geta) < 0.5) lNumberOfCharged++; 
1302    }//End of loop on tracks
1303 //----- End Loop on Stack ------------------------------------------------------------
1304
1305    //Double_t lpArapidityShift = 0.465;
1306    Bool_t lStackNatural = kTRUE;
1307 //----- Loop on Lambda, K0Short ----------------------------------------------------------------
1308    for (Int_t iCurrentLabelStack = 0;  iCurrentLabelStack < (lMCstack->GetNtrack()); iCurrentLabelStack++) 
1309    {// This is the begining of the loop on tracks
1310       
1311       TParticle* lCurrentParticleForLambdaCheck = 0x0; 
1312       lCurrentParticleForLambdaCheck = lMCstack->Particle( iCurrentLabelStack );
1313       if(!lCurrentParticleForLambdaCheck){
1314          Printf("V0s loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
1315          continue;
1316       }
1317
1318       //=================================================================================
1319       //Single-Strange checks
1320       // Keep only K0s, Lambda and AntiLambda:
1321       lPdgcodeCurrentPart = lCurrentParticleForLambdaCheck->GetPdgCode();             
1322
1323       if ( (lCurrentParticleForLambdaCheck->GetPdgCode() == 310   ) ||
1324            (lCurrentParticleForLambdaCheck->GetPdgCode() == 3122  ) ||
1325            (lCurrentParticleForLambdaCheck->GetPdgCode() == -3122 ) )
1326            {
1327          lRapCurrentPart   = MyRapidity(lCurrentParticleForLambdaCheck->Energy(),lCurrentParticleForLambdaCheck->Pz());
1328          lPtCurrentPart    = lCurrentParticleForLambdaCheck->Pt();
1329
1330           //Use Close to PV for filling CloseToPV histograms!
1331          Double_t dx, dy, dz; 
1332
1333          dx = ( (mcPrimaryVtx.At(0)) - (lCurrentParticleForLambdaCheck->Vx()) ); 
1334          dy = ( (mcPrimaryVtx.At(1)) - (lCurrentParticleForLambdaCheck->Vy()) );
1335          dz = ( (mcPrimaryVtx.At(2)) - (lCurrentParticleForLambdaCheck->Vz()) );
1336          Double_t lDistToPV = TMath::Sqrt(dx*dx + dy*dy + dz*dz);
1337          if( lDistToPV <= 0.001){ 
1338            if( lPdgcodeCurrentPart == 3122 ){
1339               f3dHistPrimCloseToPVPtVsYVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1340            }
1341            if( lPdgcodeCurrentPart == -3122 ){
1342               f3dHistPrimCloseToPVPtVsYVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1343            }
1344            if( lPdgcodeCurrentPart == 310 ){
1345               f3dHistPrimCloseToPVPtVsYVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1346            }
1347          }
1348
1349          //Use Physical Primaries only for filling PrimRaw Histograms!
1350          if ( lMCstack->IsPhysicalPrimary(iCurrentLabelStack)!=kTRUE ) continue;
1351
1352           lStackNatural = lMCevent->IsFromBGEvent(iCurrentLabelStack); //Is it? 
1353           if (!lStackNatural){
1354             if (!(lCurrentParticleForLambdaCheck->GetFirstMother()<0)) 
1355               {lStackNatural = kTRUE;} // because there are primaries (ALICE definition) not produced in the collision
1356           }
1357
1358          if( lPdgcodeCurrentPart == 3122 ){
1359             f3dHistPrimRawPtVsYVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1360             f3dHistPrimRawPtVsYCMSVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart+fpArapidityShift, lMultiplicity);
1361             if(lStackNatural){f3dHistPrimRawPtVsYVsMultNonInjLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);}
1362             f3dHistPrimRawPtVsYVsMultMCLambda->Fill(lPtCurrentPart, lRapCurrentPart, lNumberOfCharged);
1363             f3dHistPrimRawPtVsYVsVertexZLambda->Fill(lPtCurrentPart, lRapCurrentPart, mcPrimaryVtx.At(2));
1364             if( TMath::Abs( lCurrentParticleForLambdaCheck->Eta() )<1.2 && lPtCurrentPart>2 ){
1365                lHasHighPtLambda = kTRUE; //Keep track of events with Lambda within |eta|<1.2 and pt>2
1366             }
1367          }
1368          if( lPdgcodeCurrentPart == -3122 ){
1369             f3dHistPrimRawPtVsYVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1370             f3dHistPrimRawPtVsYCMSVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart+fpArapidityShift, lMultiplicity);
1371             if(lStackNatural){f3dHistPrimRawPtVsYVsMultNonInjAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);}
1372             f3dHistPrimRawPtVsYVsMultMCAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lNumberOfCharged);
1373             f3dHistPrimRawPtVsYVsVertexZAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, mcPrimaryVtx.At(2));
1374          }
1375          if( lPdgcodeCurrentPart == 310 ){
1376             f3dHistPrimRawPtVsYVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1377             f3dHistPrimRawPtVsYCMSVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart+fpArapidityShift, lMultiplicity);
1378             if(lStackNatural){f3dHistPrimRawPtVsYVsMultNonInjK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);}
1379             f3dHistPrimRawPtVsYVsMultMCK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lNumberOfCharged);
1380             f3dHistPrimRawPtVsYVsVertexZK0Short->Fill(lPtCurrentPart, lRapCurrentPart, mcPrimaryVtx.At(2));
1381          }
1382          //Decay Length Acquisition=====================================================
1383          Double_t decaylength = -1; 
1384          Double_t lV0Mass = -1; 
1385           
1386          if( !(lCurrentParticleForLambdaCheck->GetDaughter(0) < 0) ) {
1387             TParticle* lDght0ofV0 = lMCstack->Particle(  lCurrentParticleForLambdaCheck->GetDaughter(0) ); //get first daughter
1388             if(lDght0ofV0){ // skip if not defined. 
1389                decaylength = TMath::Sqrt(
1390                                         TMath::Power( lCurrentParticleForLambdaCheck->Vx() - lDght0ofV0->Vx() , 2) + 
1391                                         TMath::Power( lCurrentParticleForLambdaCheck->Vy() - lDght0ofV0->Vy() , 2) + 
1392                                         TMath::Power( lCurrentParticleForLambdaCheck->Vz() - lDght0ofV0->Vz() , 2)
1393                );
1394                //Need to correct for relativitity! Involves multiplying by mass and dividing by momentum. 
1395                if(TMath::Abs( lPdgcodeCurrentPart ) == 3122 ) { lV0Mass = 1.115683; }
1396                if(TMath::Abs( lPdgcodeCurrentPart ) == 310 ) { lV0Mass = 0.497614; }
1397                if( lCurrentParticleForLambdaCheck->P() + 1e-10 != 0 ) decaylength = ( lV0Mass * decaylength ) / ( lCurrentParticleForLambdaCheck->P() + 1e-10 );
1398                if( lCurrentParticleForLambdaCheck->P() + 1e-10 == 0 ) decaylength = 1e+5;
1399             }
1400          }
1401          if( lPdgcodeCurrentPart == 3122) f3dHistPrimRawPtVsYVsDecayLengthLambda ->Fill( lPtCurrentPart, lRapCurrentPart , decaylength ); 
1402          if( lPdgcodeCurrentPart == -3122) f3dHistPrimRawPtVsYVsDecayLengthAntiLambda ->Fill( lPtCurrentPart, lRapCurrentPart , decaylength ); 
1403          if( lPdgcodeCurrentPart == 310) f3dHistPrimRawPtVsYVsDecayLengthK0Short ->Fill( lPtCurrentPart, lRapCurrentPart , decaylength ); 
1404       }
1405    }//End of loop on tracks
1406 //----- End Loop on Lambda, K0Short ------------------------------------------------------------
1407
1408
1409    f2dHistMultiplicityVsTrueBeforeTrigSel->Fill ( lMultiplicity , lNumberOfCharged );
1410
1411     fTreeVariableMultiplicityMC = lNumberOfCharged;
1412
1413    fHistGenVertexZBeforeTrigSel->Fill( (mcPrimaryVtx.At(2)) );
1414
1415    lPdgcodeCurrentPart = 0;
1416    lRapCurrentPart  = 0;
1417    lPtCurrentPart   = 0;
1418
1419 //------------------------------------------------
1420 // Physics Selection
1421 //------------------------------------------------
1422
1423    UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1424    Bool_t isSelected = 0;
1425    isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
1426
1427    //pA triggering: CINT7
1428    if( fkSwitchINT7 ) isSelected = (maskIsSelected & AliVEvent::kINT7) == AliVEvent::kINT7;
1429
1430    //Standard Min-Bias Selection
1431    if ( ! isSelected ) { 
1432       PostData(1, fListHistV0);
1433       PostData(2, fTree);
1434       return;
1435    }
1436
1437    f2dHistMultiplicityVsTrueForTrigEvt->Fill ( lMultiplicity , lNumberOfCharged );
1438    fHistGenVertexZForTrigEvt->Fill( mcPrimaryVtx.At(2) );
1439 //------------------------------------------------
1440 // After Trigger Selection
1441 //------------------------------------------------
1442
1443    lNumberOfV0s          = lESDevent->GetNumberOfV0s();
1444   
1445    //Set variable for filling tree afterwards!
1446    fHistV0MultiplicityForTrigEvt->Fill(lNumberOfV0s);
1447    fHistMultiplicityForTrigEvt->Fill ( lMultiplicity );
1448
1449 //------------------------------------------------
1450 // Getting: Primary Vertex + MagField Info
1451 //------------------------------------------------
1452
1453    const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
1454    // get the vtx stored in ESD found with tracks
1455    lPrimaryTrackingESDVtx->GetXYZ( lTrkgPrimaryVtxPos );
1456         
1457    const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();      
1458    // get the best primary vertex available for the event
1459    // As done in AliCascadeVertexer, we keep the one which is the best one available.
1460    // between : Tracking vertex > SPD vertex > TPC vertex > default SPD vertex
1461    // This one will be used for next calculations (DCA essentially)
1462    lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
1463
1464    Double_t lPrimaryVtxPosition[3];
1465    const AliVVertex *primaryVtx = lESDevent->GetPrimaryVertex();
1466    lPrimaryVtxPosition[0] = primaryVtx->GetX();
1467    lPrimaryVtxPosition[1] = primaryVtx->GetY();
1468    lPrimaryVtxPosition[2] = primaryVtx->GetZ();
1469    fHistPVx->Fill( lPrimaryVtxPosition[0] );
1470    fHistPVy->Fill( lPrimaryVtxPosition[1] );
1471    fHistPVz->Fill( lPrimaryVtxPosition[2] );
1472
1473    f2dHistMultiplicityVsVertexZForTrigEvt->Fill( lMultiplicity, lPrimaryVtxPosition[2] );
1474
1475 //------------------------------------------------
1476 // Primary Vertex Z position: SKIP
1477 //------------------------------------------------
1478
1479    if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0 ) { 
1480       AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !"); 
1481       PostData(1, fListHistV0);
1482       PostData(2, fTree);
1483       return; 
1484    }
1485
1486    f2dHistMultiplicityVsVertexZ->Fill( lMultiplicity, lPrimaryVtxPosition[2] );
1487
1488    lMagneticField = lESDevent->GetMagneticField( );
1489    fHistV0MultiplicityForSelEvt ->Fill( lNumberOfV0s );
1490    fHistMultiplicity->Fill(lMultiplicity);
1491    f2dHistMultiplicityVsTrue->Fill ( lMultiplicity , lNumberOfCharged );
1492    fHistGenVertexZ->Fill( (mcPrimaryVtx.At(2)) );
1493 //------------------------------------------------
1494 // SKIP: Events with well-established PVtx
1495 //------------------------------------------------
1496         
1497    const AliESDVertex *lPrimaryTrackingESDVtxCheck = lESDevent->GetPrimaryVertexTracks();
1498    const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
1499    if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtxCheck->GetStatus() ){
1500       AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
1501       PostData(1, fListHistV0);
1502       PostData(2, fTree);
1503       return;
1504    }
1505
1506    f2dHistMultiplicityVsVertexZNoTPCOnly->Fill( lMultiplicity, lPrimaryVtxPosition[2] );
1507    fHistV0MultiplicityForSelEvtNoTPCOnly ->Fill( lNumberOfV0s );
1508    fHistMultiplicityNoTPCOnly->Fill(lMultiplicity);
1509    f2dHistMultiplicityVsTrueNoTPCOnly->Fill ( lMultiplicity , lNumberOfCharged );
1510    fHistGenVertexZNoTPCOnly->Fill( (mcPrimaryVtx.At(2)) );
1511 //------------------------------------------------
1512 // Pileup Rejection Studies
1513 //------------------------------------------------
1514
1515    // FIXME : quality selection regarding pile-up rejection 
1516    if(lESDevent->IsPileupFromSPD() && !fkIsNuclear ){// minContributors=3, minZdist=0.8, nSigmaZdist=3., nSigmaDiamXY=2., nSigmaDiamZ=5.  -> see http://alisoft.cern.ch/viewvc/trunk/STEER/AliESDEvent.h?root=AliRoot&r1=41914&r2=42199&pathrev=42199
1517       AliWarning("Pb / This is tagged as Pileup from SPD... return !");
1518       PostData(1, fListHistV0);
1519       PostData(2, fTree);
1520       return;
1521    }
1522    f2dHistMultiplicityVsVertexZNoTPCOnlyNoPileup->Fill( lMultiplicity, lPrimaryVtxPosition[2] );
1523    fHistV0MultiplicityForSelEvtNoTPCOnlyNoPileup ->Fill( lNumberOfV0s );
1524    fHistMultiplicityNoTPCOnlyNoPileup->Fill(lMultiplicity);
1525    f2dHistMultiplicityVsTrueNoTPCOnlyNoPileup->Fill ( lMultiplicity , lNumberOfCharged );
1526    fHistGenVertexZNoTPCOnlyNoPileup->Fill( (mcPrimaryVtx.At(2)) );
1527    //Do control histograms without the IsFromVertexerZ events, but consider them in analysis...
1528    if( ! (lESDevent->GetPrimaryVertex()->IsFromVertexerZ() )     ){ 
1529       fHistPVxAnalysis->Fill( lPrimaryVtxPosition[0] );
1530       fHistPVyAnalysis->Fill( lPrimaryVtxPosition[1] );
1531       fHistPVzAnalysis->Fill( lPrimaryVtxPosition[2] );
1532       if ( lHasHighPtLambda == kTRUE ){ 
1533          fHistPVxAnalysisHasHighPtLambda->Fill( lPrimaryVtxPosition[0] );
1534          fHistPVyAnalysisHasHighPtLambda->Fill( lPrimaryVtxPosition[1] );
1535          fHistPVzAnalysisHasHighPtLambda->Fill( lPrimaryVtxPosition[2] );
1536       }
1537    }
1538
1539   fTreeVariableVertexZ = lPrimaryVtxPosition[2];
1540
1541   fTreeVariablePVx = lPrimaryVtxPosition[0];
1542   fTreeVariablePVy = lPrimaryVtxPosition[1];
1543   fTreeVariablePVz = lPrimaryVtxPosition[2];
1544
1545   fTreeVariableMCPVx = (mcPrimaryVtx.At(0));
1546   fTreeVariableMCPVy = (mcPrimaryVtx.At(1));
1547   fTreeVariableMCPVz = (mcPrimaryVtx.At(2));
1548
1549 //------------------------------------------------
1550 // stack loop starts here
1551 //------------------------------------------------
1552
1553 //---> Loop over ALL PARTICLES
1554  
1555    for (Int_t iMc = 0; iMc < (lMCstack->GetNtrack()); iMc++) {  
1556       TParticle *p0 = lMCstack->Particle(iMc); 
1557       if (!p0) {
1558          //Printf("ERROR: particle with label %d not found in lMCstack (mc loop)", iMc);
1559          continue;
1560       }
1561       lPdgcodeCurrentPart = p0->GetPdgCode();
1562
1563       // Keep only K0s, Lambda and AntiLambda:
1564       if ( (lPdgcodeCurrentPart != 310 ) && (lPdgcodeCurrentPart != 3122 ) && (lPdgcodeCurrentPart != -3122 ) ) continue;
1565         
1566       lRapCurrentPart   = MyRapidity(p0->Energy(),p0->Pz());
1567       lPtCurrentPart    = p0->Pt();
1568
1569         //Use Physical Primaries only for filling PrimRaw Histograms!
1570       if ( lMCstack->IsPhysicalPrimary(iMc)!=kTRUE ) continue;
1571
1572       if( lPdgcodeCurrentPart == 3122 ){
1573          f3dHistPrimAnalysisPtVsYVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1574          f3dHistPrimAnalysisPtVsYCMSVsMultLambda->Fill(lPtCurrentPart, lRapCurrentPart+fpArapidityShift, lMultiplicity);
1575       }
1576       if( lPdgcodeCurrentPart == -3122 ){
1577          f3dHistPrimAnalysisPtVsYVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1578          f3dHistPrimAnalysisPtVsYCMSVsMultAntiLambda->Fill(lPtCurrentPart, lRapCurrentPart+fpArapidityShift, lMultiplicity);
1579       }
1580       if( lPdgcodeCurrentPart == 310 ){
1581          f3dHistPrimAnalysisPtVsYVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart, lMultiplicity);
1582          f3dHistPrimAnalysisPtVsYCMSVsMultK0Short->Fill(lPtCurrentPart, lRapCurrentPart+fpArapidityShift, lMultiplicity);
1583       }
1584    }
1585
1586 //----- Loop on primary Xi, Omega --------------------------------------------------------------
1587    for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < lNbMCPrimary; iCurrentLabelStack++) 
1588    {// This is the begining of the loop on primaries
1589       
1590       TParticle* lCurrentParticlePrimary = 0x0; 
1591       lCurrentParticlePrimary = lMCstack->Particle( iCurrentLabelStack );
1592       if(!lCurrentParticlePrimary){
1593          Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
1594          continue;
1595       }
1596       if ( TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3312 || TMath::Abs(lCurrentParticlePrimary->GetPdgCode()) == 3334 ) { 
1597          Double_t lRapXiMCPrimary = -100;
1598          if( (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) != 0 ) { 
1599            if ( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) !=0 ){
1600              lRapXiMCPrimary = 0.5*TMath::Log( (lCurrentParticlePrimary->Energy() + lCurrentParticlePrimary->Pz()) / (lCurrentParticlePrimary->Energy() - lCurrentParticlePrimary->Pz() +1.e-13) );
1601            }
1602          }
1603
1604          //=================================================================================
1605          // Xi Histograms
1606          if( lCurrentParticlePrimary->GetPdgCode() == 3312 ){ 
1607             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1608             f3dHistGenSelectedPtVsYVsMultXiMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1609             f3dHistGenSelectedPtVsYCMSVsMultXiMinus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1610          }
1611          if( lCurrentParticlePrimary->GetPdgCode() == -3312 ){ 
1612             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1613             f3dHistGenSelectedPtVsYVsMultXiPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1614             f3dHistGenSelectedPtVsYCMSVsMultXiPlus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1615          }
1616          // Omega Histograms
1617          if( lCurrentParticlePrimary->GetPdgCode() == 3334 ){ 
1618             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1619             f3dHistGenSelectedPtVsYVsMultOmegaMinus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1620             f3dHistGenSelectedPtVsYCMSVsMultOmegaMinus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1621          }
1622          if( lCurrentParticlePrimary->GetPdgCode() == -3334 ){ 
1623             lPtCurrentPart    = lCurrentParticlePrimary->Pt();
1624             f3dHistGenSelectedPtVsYVsMultOmegaPlus->Fill(lPtCurrentPart, lRapXiMCPrimary, lMultiplicity);
1625             f3dHistGenSelectedPtVsYCMSVsMultOmegaPlus->Fill(lPtCurrentPart, lRapXiMCPrimary+fpArapidityShift, lMultiplicity);
1626          }
1627       } 
1628    }
1629 //----- End Loop on primary Xi, Omega ----------------------------------------------------------
1630
1631 //------------------------------------------------
1632 // MAIN LAMBDA LOOP STARTS HERE
1633 //------------------------------------------------
1634
1635    //Variable definition
1636    Int_t    lOnFlyStatus = 0;
1637    Double_t lChi2V0 = 0;
1638    Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
1639    Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
1640    Double_t lV0CosineOfPointingAngle = 0;
1641    Double_t lV0Radius = 0, lPt = 0;
1642    Double_t lRapK0Short = 0, lRapLambda = 0;
1643    Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
1644    Double_t lAlphaV0 = 0, lPtArmV0 = 0;
1645    Double_t fMinV0Pt = 0; 
1646    Double_t fMaxV0Pt = 100; 
1647
1648    Int_t nv0s = 0;
1649    nv0s = lESDevent->GetNumberOfV0s();
1650    
1651    for (Int_t iV0 = 0; iV0 < nv0s; iV0++) 
1652         {// This is the begining of the V0 loop
1653       AliESDv0 *v0 = ((AliESDEvent*)lESDevent)->GetV0(iV0);
1654       if (!v0) continue;
1655
1656       //---> Fix On-the-Fly candidates
1657       if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ){
1658         fHistSwappedV0Counter -> Fill( 1 );
1659       }else{
1660         fHistSwappedV0Counter -> Fill( 0 ); 
1661       }
1662       if ( fkUseOnTheFly ) CheckChargeV0(v0); 
1663
1664
1665       Double_t tV0mom[3];
1666       v0->GetPxPyPz( tV0mom[0],tV0mom[1],tV0mom[2] ); 
1667       Double_t lV0TotalMomentum = TMath::Sqrt(
1668          tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );
1669
1670       Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]); 
1671       lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);
1672       lPt = v0->Pt();
1673       lRapK0Short = v0->RapK0Short();
1674       lRapLambda  = v0->RapLambda();
1675
1676       //Set Variables for later filling
1677       fTreeVariableV0x = tDecayVertexV0[0];
1678       fTreeVariableV0y = tDecayVertexV0[1];
1679       fTreeVariableV0z = tDecayVertexV0[2];
1680
1681       //Set Variables for later filling
1682       fTreeVariableV0Px = tV0mom[0];
1683       fTreeVariableV0Py = tV0mom[1];
1684       fTreeVariableV0Pz = tV0mom[2];
1685
1686       if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;
1687
1688       UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex());
1689       UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex());
1690
1691       Double_t lMomPos[3]; v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
1692       Double_t lMomNeg[3]; v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
1693
1694       AliESDtrack *pTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyPos);
1695       AliESDtrack *nTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyNeg);
1696       if (!pTrack || !nTrack) {
1697          Printf("ERROR: Could not retreive one of the daughter track");
1698          continue;
1699       }
1700
1701       fTreeVariableNegEta = nTrack->Eta();
1702       fTreeVariablePosEta = pTrack->Eta();
1703
1704       // Filter like-sign V0 (next: add counter and distribution)
1705       if ( pTrack->GetSign() == nTrack->GetSign()){
1706          continue;
1707       } 
1708
1709       //________________________________________________________________________
1710       // Track quality cuts 
1711       Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);
1712       Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);
1713       fTreeVariableLeastNbrCrossedRows = (Int_t) lPosTrackCrossedRows;
1714       if( lNegTrackCrossedRows < fTreeVariableLeastNbrCrossedRows )
1715          fTreeVariableLeastNbrCrossedRows = (Int_t) lNegTrackCrossedRows;
1716
1717       // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
1718       if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;      
1719       if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
1720
1721       //Get status flags
1722       fTreeVariablePosTrackStatus = pTrack->GetStatus();
1723       fTreeVariableNegTrackStatus = nTrack->GetStatus();
1724     
1725       if ( ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) )&&(fkTakeAllTracks==kFALSE) ) continue;
1726         
1727       //GetKinkIndex condition
1728       if( pTrack->GetKinkIndex(0)>0 || nTrack->GetKinkIndex(0)>0 ) continue;
1729
1730       //Findable clusters > 0 condition
1731       if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;
1732
1733       //Compute ratio Crossed Rows / Findable clusters
1734       //Note: above test avoids division by zero! 
1735       Float_t lPosTrackCrossedRowsOverFindable = -1;
1736       Float_t lNegTrackCrossedRowsOverFindable = -1;
1737       if ( ((double)(pTrack->GetTPCNclsF()) ) != 0 ) lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF())); 
1738       if ( ((double)(nTrack->GetTPCNclsF()) ) != 0 ) lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF())); 
1739
1740       fTreeVariableLeastRatioCrossedRowsOverFindable = lPosTrackCrossedRowsOverFindable;
1741       if( lNegTrackCrossedRowsOverFindable < fTreeVariableLeastRatioCrossedRowsOverFindable )
1742          fTreeVariableLeastRatioCrossedRowsOverFindable = lNegTrackCrossedRowsOverFindable;
1743
1744       //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here
1745       if ( (fTreeVariableLeastRatioCrossedRowsOverFindable < 0.8)&&(fkTakeAllTracks==kFALSE) ) continue;
1746
1747       //End track Quality Cuts
1748       //________________________________________________________________________
1749
1750       lDcaPosToPrimVertex = TMath::Abs(pTrack->GetD(lPrimaryVtxPosition[0],
1751                                                         lPrimaryVtxPosition[1],
1752                                                         lMagneticField) );
1753
1754       lDcaNegToPrimVertex = TMath::Abs(nTrack->GetD(lPrimaryVtxPosition[0],
1755                                                         lPrimaryVtxPosition[1],
1756                                                         lMagneticField) );
1757
1758       lOnFlyStatus = v0->GetOnFlyStatus();
1759       lChi2V0 = v0->GetChi2V0();
1760       lDcaV0Daughters = v0->GetDcaV0Daughters();
1761       lDcaV0ToPrimVertex = v0->GetD(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lPrimaryVtxPosition[2]);
1762       lV0CosineOfPointingAngle = v0->GetV0CosineOfPointingAngle(lPrimaryVtxPosition[0],lPrimaryVtxPosition[1],lPrimaryVtxPosition[2]);
1763       fTreeVariableV0CosineOfPointingAngle=lV0CosineOfPointingAngle;
1764
1765       // Getting invariant mass infos directly from ESD
1766       v0->ChangeMassHypothesis(310);
1767       lInvMassK0s = v0->GetEffMass();
1768       v0->ChangeMassHypothesis(3122);
1769       lInvMassLambda = v0->GetEffMass();
1770       v0->ChangeMassHypothesis(-3122);
1771       lInvMassAntiLambda = v0->GetEffMass();
1772       lAlphaV0 = v0->AlphaV0();
1773       lPtArmV0 = v0->PtArmV0();
1774
1775       //fTreeVariableOnFlyStatus = lOnFlyStatus;
1776       //fHistV0OnFlyStatus->Fill(lOnFlyStatus);
1777
1778 //===============================================
1779 // Monte Carlo Association starts here
1780 //===============================================
1781
1782       //---> Set Everything to "I don't know" before starting
1783
1784       fTreeVariablePIDPositive = 0;
1785       fTreeVariablePIDNegative = 0;
1786
1787       fTreeVariableIndexStatus = 0;
1788       fTreeVariableIndexStatusMother = 0;
1789
1790       fTreeVariablePtMother = -1;
1791       fTreeVariablePtMC = -1;
1792       fTreeVariableRapMC = -100;
1793
1794       fTreeVariablePID = -1; 
1795       fTreeVariablePIDMother = -1;
1796
1797       fTreeVariablePrimaryStatus = 0; 
1798       fTreeVariablePrimaryStatusMother = 0; 
1799       fTreeVariableV0CreationRadius = -1;
1800
1801       Int_t lblPosV0Dghter = (Int_t) TMath::Abs( pTrack->GetLabel() );  
1802       Int_t lblNegV0Dghter = (Int_t) TMath::Abs( nTrack->GetLabel() );
1803                 
1804       TParticle* mcPosV0Dghter = lMCstack->Particle( lblPosV0Dghter );
1805       TParticle* mcNegV0Dghter = lMCstack->Particle( lblNegV0Dghter );
1806             
1807       fTreeVariablePosTransvMomentumMC = mcPosV0Dghter->Pt();
1808       fTreeVariableNegTransvMomentumMC = mcNegV0Dghter->Pt();
1809
1810       Int_t lPIDPositive = mcPosV0Dghter -> GetPdgCode();
1811       Int_t lPIDNegative = mcNegV0Dghter -> GetPdgCode();
1812
1813       fTreeVariablePIDPositive = lPIDPositive;
1814       fTreeVariablePIDNegative = lPIDNegative;
1815
1816       Int_t lblMotherPosV0Dghter = mcPosV0Dghter->GetFirstMother() ; 
1817       Int_t lblMotherNegV0Dghter = mcNegV0Dghter->GetFirstMother();
1818
1819       if( lblMotherPosV0Dghter == lblMotherNegV0Dghter && lblMotherPosV0Dghter > -1 ){
1820          //either label is fine, they're equal at this stage
1821          TParticle* pThisV0 = lMCstack->Particle( lblMotherPosV0Dghter ); 
1822          //Set tree variables
1823          fTreeVariablePID   = pThisV0->GetPdgCode(); //PDG Code
1824          fTreeVariablePtMC  = pThisV0->Pt(); //Perfect Pt
1825
1826           fTreeVariableIsNonInjected = lMCevent->IsFromBGEvent(lblMotherPosV0Dghter); //Is it? 
1827           if (!fTreeVariableIsNonInjected){
1828             if (!(pThisV0->GetFirstMother()<0)) 
1829               {fTreeVariableIsNonInjected = kTRUE;} // because there are primaries (ALICE definition) not produced in the collision
1830           }
1831
1832          //Set Variables for later filling
1833          //Be careful: Vx, Vy, Vz: Creation vertex. So decay position is the 
1834          //Creation vertex of any one of the daughters!
1835          fTreeVariableMCV0x = mcPosV0Dghter->Vx();
1836          fTreeVariableMCV0y = mcPosV0Dghter->Vy();
1837          fTreeVariableMCV0z = mcPosV0Dghter->Vz();
1838
1839          //Set Variables for later filling
1840          fTreeVariableMCV0Px = pThisV0->Px();
1841          fTreeVariableMCV0Py = pThisV0->Py();
1842          fTreeVariableMCV0Pz = pThisV0->Pz();
1843
1844          //Only Interested if it's a Lambda, AntiLambda or K0s 
1845          //Avoid the Junction Bug! PYTHIA has particles with Px=Py=Pz=E=0 occasionally, 
1846          //having particle code 88 (unrecognized by PDG), for documentation purposes.
1847          //Even ROOT's TParticle::Y() is not prepared to deal with that exception!
1848          //Note that TParticle::Pt() is immune (that would just return 0)...
1849          //Though granted that that should be extremely rare in this precise condition...
1850          if( TMath::Abs(fTreeVariablePID) == 3122 || fTreeVariablePID==310 ){
1851             fTreeVariableRapMC = pThisV0->Y(); //Perfect Y
1852          }
1853          fTreeVariableV0CreationRadius = TMath::Sqrt(
1854           TMath::Power(  ( (mcPrimaryVtx.At(0)) - (pThisV0->Vx()) ) , 2) + 
1855           TMath::Power(  ( (mcPrimaryVtx.At(1)) - (pThisV0->Vy()) ) , 2) + 
1856           TMath::Power(  ( (mcPrimaryVtx.At(2)) - (pThisV0->Vz()) ) , 2) 
1857          );
1858          if( lblMotherPosV0Dghter  < lNbMCPrimary ) fTreeVariableIndexStatus = 1; //looks primary
1859          if( lblMotherPosV0Dghter >= lNbMCPrimary ) fTreeVariableIndexStatus = 2; //looks secondary
1860          if( lMCstack->IsPhysicalPrimary       (lblMotherPosV0Dghter) ) fTreeVariablePrimaryStatus = 1; //Is Primary!
1861          if( lMCstack->IsSecondaryFromWeakDecay(lblMotherPosV0Dghter) ) fTreeVariablePrimaryStatus = 2; //Weak Decay!
1862          if( lMCstack->IsSecondaryFromMaterial (lblMotherPosV0Dghter) ) fTreeVariablePrimaryStatus = 3; //Material Int!
1863          
1864          //Now we try to acquire the V0 parent particle, if possible
1865          Int_t lblThisV0Parent = pThisV0->GetFirstMother();
1866          if ( lblThisV0Parent > -1 ){ //if it has a parent, get it and store specs
1867             TParticle* pThisV0Parent = lMCstack->Particle( lblThisV0Parent );
1868             fTreeVariablePIDMother   = pThisV0Parent->GetPdgCode(); //V0 Mother PDG
1869             fTreeVariablePtMother    = pThisV0Parent->Pt();         //V0 Mother Pt
1870             //Primary Status for the V0 Mother particle 
1871             if( lblThisV0Parent  < lNbMCPrimary ) fTreeVariableIndexStatusMother = 1; //looks primary
1872             if( lblThisV0Parent >= lNbMCPrimary ) fTreeVariableIndexStatusMother = 2; //looks secondary
1873             if( lMCstack->IsPhysicalPrimary       (lblThisV0Parent) ) fTreeVariablePrimaryStatusMother = 1; //Is Primary!
1874             if( lMCstack->IsSecondaryFromWeakDecay(lblThisV0Parent) ) fTreeVariablePrimaryStatusMother = 2; //Weak Decay!
1875             if( lMCstack->IsSecondaryFromMaterial (lblThisV0Parent) ) fTreeVariablePrimaryStatusMother = 3; //Material Int!
1876          }
1877       }
1878
1879       fTreeVariablePt = v0->Pt();
1880       fTreeVariableChi2V0 = lChi2V0; 
1881       fTreeVariableDcaV0ToPrimVertex = lDcaV0ToPrimVertex;
1882       fTreeVariableDcaV0Daughters = lDcaV0Daughters;
1883       fTreeVariableV0CosineOfPointingAngle = lV0CosineOfPointingAngle; 
1884       fTreeVariableV0Radius = lV0Radius;
1885       fTreeVariableDcaPosToPrimVertex = lDcaPosToPrimVertex;
1886       fTreeVariableDcaNegToPrimVertex = lDcaNegToPrimVertex;
1887       fTreeVariableInvMassK0s = lInvMassK0s;
1888       fTreeVariableInvMassLambda = lInvMassLambda;
1889       fTreeVariableInvMassAntiLambda = lInvMassAntiLambda;
1890       fTreeVariableRapK0Short = lRapK0Short;
1891
1892       fTreeVariableRapLambda = lRapLambda;
1893       fTreeVariableAlphaV0 = lAlphaV0;
1894       fTreeVariablePtArmV0 = lPtArmV0;
1895
1896       //Official means of acquiring N-sigmas 
1897       fTreeVariableNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );
1898       fTreeVariableNSigmasPosPion   = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );
1899       fTreeVariableNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );
1900       fTreeVariableNSigmasNegPion   = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );
1901
1902 //tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]
1903       Double_t lDistanceTravelled = TMath::Sqrt(
1904                                                 TMath::Power( tDecayVertexV0[0] - lBestPrimaryVtxPos[0] , 2) +
1905                                                 TMath::Power( tDecayVertexV0[1] - lBestPrimaryVtxPos[1] , 2) +
1906                                                 TMath::Power( tDecayVertexV0[2] - lBestPrimaryVtxPos[2] , 2)
1907                                         );
1908       fTreeVariableDistOverTotMom = 1e+5;
1909       if( lV0TotalMomentum + 1e-10 != 0 ) fTreeVariableDistOverTotMom = lDistanceTravelled / (lV0TotalMomentum + 1e-10); //avoid division by zero, to be sure
1910
1911       Double_t lMomentumPosTemp[3];
1912       pTrack->GetPxPyPz(lMomentumPosTemp);
1913       Double_t lPtPosTemporary = sqrt(pow(lMomentumPosTemp[0],2) + pow(lMomentumPosTemp[1],2));
1914
1915       Double_t lMomentumNegTemp[3];
1916       nTrack->GetPxPyPz(lMomentumNegTemp);
1917       Double_t lPtNegTemporary = sqrt(pow(lMomentumNegTemp[0],2) + pow(lMomentumNegTemp[1],2));
1918
1919       fTreeVariablePosTransvMomentum = lPtPosTemporary;
1920       fTreeVariableNegTransvMomentum = lPtNegTemporary;
1921
1922
1923 //------------------------------------------------
1924 // Fill Tree! 
1925 //------------------------------------------------
1926
1927       // The conditionals are meant to decrease excessive
1928       // memory usage! 
1929
1930       //Modified version: Keep only OnFlyStatus == 0
1931       //Keep only if included in a parametric InvMass Region 20 sigmas away from peak
1932
1933       //First Selection: Reject OnFly
1934       if( (lOnFlyStatus == 0 && fkUseOnTheFly == kFALSE) || (lOnFlyStatus != 0 && fkUseOnTheFly == kTRUE ) ){
1935          //Second Selection: rough 20-sigma band, parametric. 
1936          //K0Short: Enough to parametrize peak broadening with linear function.    
1937          Double_t lUpperLimitK0Short = (5.63707e-01) + (1.14979e-02)*fTreeVariablePt; 
1938          Double_t lLowerLimitK0Short = (4.30006e-01) - (1.10029e-02)*fTreeVariablePt;
1939          //Lambda: Linear (for higher pt) plus exponential (for low-pt broadening)
1940          //[0]+[1]*x+[2]*TMath::Exp(-[3]*x)
1941          Double_t lUpperLimitLambda = (1.13688e+00) + (5.27838e-03)*fTreeVariablePt + (8.42220e-02)*TMath::Exp(-(3.80595e+00)*fTreeVariablePt); 
1942          Double_t lLowerLimitLambda = (1.09501e+00) - (5.23272e-03)*fTreeVariablePt - (7.52690e-02)*TMath::Exp(-(3.46339e+00)*fTreeVariablePt);
1943          //Do Selection      
1944          if( (fTreeVariableInvMassLambda     < lUpperLimitLambda  && fTreeVariableInvMassLambda     > lLowerLimitLambda     ) || 
1945              (fTreeVariableInvMassAntiLambda < lUpperLimitLambda  && fTreeVariableInvMassAntiLambda > lLowerLimitLambda     ) || 
1946              (fTreeVariableInvMassK0s        < lUpperLimitK0Short && fTreeVariableInvMassK0s        > lLowerLimitK0Short    ) ){
1947              //Pre-selection in case this is AA...
1948              if( fkIsNuclear == kFALSE ) fTree->Fill();
1949              if( fkIsNuclear == kTRUE){ 
1950              //If this is a nuclear collision___________________
1951              // ... pre-filter with daughter eta selection only (not TPC)
1952                if ( TMath::Abs(fTreeVariableNegEta)<0.8 && TMath::Abs(fTreeVariablePosEta)<0.8 ) fTree->Fill();
1953              }//end nuclear_____________________________________
1954          }
1955       }
1956
1957 //------------------------------------------------
1958 // Fill tree over.
1959 //------------------------------------------------
1960
1961
1962    }// This is the end of the V0 loop
1963
1964 //------------------------------------------------
1965
1966    // Post output data.
1967    PostData(1, fListHistV0);
1968    PostData(2, fTree);
1969 }
1970
1971 //________________________________________________________________________
1972 void AliAnalysisTaskExtractPerformanceV0::Terminate(Option_t *)
1973 {
1974    // Draw result to the screen
1975    // Called once at the end of the query
1976
1977    TList *cRetrievedList = 0x0;
1978    cRetrievedList = (TList*)GetOutputData(1);
1979    if(!cRetrievedList){
1980       Printf("ERROR - AliAnalysisTaskExtractV0 : ouput data container list not available\n");
1981       return;
1982    }    
1983         
1984    fHistV0MultiplicityForTrigEvt = dynamic_cast<TH1F*> (  cRetrievedList->FindObject("fHistV0MultiplicityForTrigEvt")  );
1985    if (!fHistV0MultiplicityForTrigEvt) {
1986       Printf("ERROR - AliAnalysisTaskExtractV0 : fHistV0MultiplicityForTrigEvt not available");
1987       return;
1988    }
1989   
1990    TCanvas *canCheck = new TCanvas("AliAnalysisTaskExtractV0","V0 Multiplicity",10,10,510,510);
1991    canCheck->cd(1)->SetLogy();
1992
1993    fHistV0MultiplicityForTrigEvt->SetMarkerStyle(22);
1994    fHistV0MultiplicityForTrigEvt->DrawCopy("E");
1995 }
1996
1997 //----------------------------------------------------------------------------
1998
1999 Double_t AliAnalysisTaskExtractPerformanceV0::MyRapidity(Double_t rE, Double_t rPz) const
2000 {
2001    // Local calculation for rapidity
2002    Double_t ReturnValue = -100;
2003    if( (rE-rPz+1.e-13) != 0 && (rE+rPz) != 0 ){ 
2004       ReturnValue =  0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
2005    }
2006    return ReturnValue;
2007
2008
2009 //________________________________________________________________________
2010 void AliAnalysisTaskExtractPerformanceV0::CheckChargeV0(AliESDv0 *v0)
2011 {
2012    // This function checks charge of negative and positive daughter tracks. 
2013    // If incorrectly defined (onfly vertexer), swaps out. 
2014    if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ){
2015       //V0 daughter track swapping is required! Note: everything is swapped here... P->N, N->P
2016       Long_t lCorrectNidx = v0->GetPindex();
2017       Long_t lCorrectPidx = v0->GetNindex();
2018       Double32_t        lCorrectNmom[3];
2019       Double32_t        lCorrectPmom[3];
2020       v0->GetPPxPyPz( lCorrectNmom[0], lCorrectNmom[1], lCorrectNmom[2] );
2021       v0->GetNPxPyPz( lCorrectPmom[0], lCorrectPmom[1], lCorrectPmom[2] );
2022
2023       AliExternalTrackParam     lCorrectParamN(
2024         v0->GetParamP()->GetX() , 
2025         v0->GetParamP()->GetAlpha() , 
2026         v0->GetParamP()->GetParameter() , 
2027         v0->GetParamP()->GetCovariance() 
2028       );
2029       AliExternalTrackParam     lCorrectParamP(
2030         v0->GetParamN()->GetX() , 
2031         v0->GetParamN()->GetAlpha() , 
2032         v0->GetParamN()->GetParameter() , 
2033         v0->GetParamN()->GetCovariance() 
2034       );
2035       lCorrectParamN.SetMostProbablePt( v0->GetParamP()->GetMostProbablePt() );
2036       lCorrectParamP.SetMostProbablePt( v0->GetParamN()->GetMostProbablePt() );
2037
2038       //Get Variables___________________________________________________
2039       Double_t lDcaV0Daughters = v0 -> GetDcaV0Daughters();
2040       Double_t lCosPALocal     = v0 -> GetV0CosineOfPointingAngle(); 
2041       Bool_t lOnFlyStatusLocal = v0 -> GetOnFlyStatus();
2042
2043       //Create Replacement Object_______________________________________
2044       AliESDv0 *v0correct = new AliESDv0(lCorrectParamN,lCorrectNidx,lCorrectParamP,lCorrectPidx);
2045       v0correct->SetDcaV0Daughters          ( lDcaV0Daughters   );
2046       v0correct->SetV0CosineOfPointingAngle ( lCosPALocal       );
2047       v0correct->ChangeMassHypothesis       ( kK0Short          );
2048       v0correct->SetOnFlyStatus             ( lOnFlyStatusLocal );
2049
2050       //Reverse Cluster info..._________________________________________
2051       v0correct->SetClusters( v0->GetClusters( 1 ), v0->GetClusters ( 0 ) );
2052
2053       *v0 = *v0correct;
2054       //Proper cleanup..._______________________________________________
2055       v0correct->Delete();
2056       v0correct = 0x0;
2057
2058       //Just another cross-check and output_____________________________
2059       if( v0->GetParamN()->Charge() > 0 && v0->GetParamP()->Charge() < 0 ) {
2060         AliWarning("Found Swapped Charges, tried to correct but something FAILED!");
2061       }else{
2062         //AliWarning("Found Swapped Charges and fixed.");
2063       }
2064       //________________________________________________________________
2065    }else{
2066       //Don't touch it! ---
2067       //Printf("Ah, nice. Charges are already ordered...");
2068    }
2069    return;
2070
2071