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