First preliminary tasks for PWGLF production QA (cont'd)
[u/mrichter/AliRoot.git] / PWGLF / QATasks / AliQAProdMultistrange.cxx
1 /**************************************************************************
2  *  Authors : Antonin Maire, Boris Hippolyte                              *
3  * Contributors are mentioned in the code where appropriate.              *
4  *                                                                        *
5  * Permission to use, copy, modify and distribute this software and its   *
6  * documentation strictly for non-commercial purposes is hereby granted   *
7  * without fee, provided that the above copyright notice appears in all   *
8  * copies and that both the copyright notice and this permission notice   *
9  * appear in the supporting documentation. The authors make no claims     *
10  * about the suitability of this software for any purpose. It is          *
11  * provided "as is" without express or implied warranty.                  *
12  **************************************************************************/
13
14 //-----------------------------------------------------------------
15 //            AliQAProdMultistrange class
16 //
17 //            Origin AliAnalysisTaskCheckCascade which has four roles :
18 //              1. QAing the Cascades from ESD and AOD
19 //                 Origin:  AliAnalysisTaskESDCheckV0 by Boris Hippolyte Nov2007, hippolyt@in2p3.fr
20 //              2. Prepare the plots which stand as raw material for yield extraction (wi/wo PID)
21 //              3. Supply an AliCFContainer meant to define the optimised topological selections
22 //              4. Rough azimuthal correlation study (Eta, Phi)
23 //              Adapted to Cascade : A.Maire Mar2008, antonin.maire@ires.in2p3.fr
24 //              Modified :           A.Maire Mar2010 
25 //
26 //              Adapted to PbPb analysis: M. Nicassio, maria.nicassio@ba.infn.it
27 //               Feb-August2011
28 //                - Physics selection moved to the run.C macro
29 //                - Centrality selection added (+ setters) 
30 //                - setters added (vertex range)
31 //                - histo added and histo/container binning changed 
32 //                - protection in the destructor for CAF usage          
33 //                - AliWarning disabled
34 //                - automatic settings for PID
35 //               September2011
36 //                - proper time histos/container added (V0 and Cascades)
37 //               November2011
38 //                - re-run V0's and cascade's vertexers (SetCuts instead of SetDefaultCuts!!)
39 //               Genuary2012 
40 //                - AOD analysis part completed 
41 //               March2012
42 //                - min number of TPC clusters for track selection as a parameter       
43 //               July2012
44 //                - cut on min pt for daughter tracks as a parameter (+control histos)
45 //               August2012 
46 //                - cut on pseudorapidity for daughter tracks as a parameter (+control histos for Xi-)
47 //-----------------------------------------------------------------
48
49 class TTree;
50 class TParticle;
51 class TVector3;
52
53 class AliESDVertex;
54 class AliAODVertex;
55 class AliESDv0;
56 class AliAODv0;
57
58 #include <Riostream.h>
59 #include "THnSparse.h"
60 #include "TVector3.h"
61 #include "TMath.h"
62
63 #include "AliLog.h"
64 #include "AliCentrality.h"
65 #include "AliESDEvent.h"
66 #include "AliAODEvent.h"
67 #include "AliESDtrackCuts.h"
68 #include "AliPIDResponse.h"
69
70 #include "AliInputEventHandler.h"
71 #include "AliAnalysisManager.h"
72 #include "AliESDInputHandler.h" 
73 #include "AliAODInputHandler.h"
74 #include "AliCFContainer.h"
75 #include "AliMultiplicity.h"
76
77 #include "AliESDcascade.h"
78 #include "AliAODcascade.h"
79 #include "AliAODTrack.h"
80
81 #include "AliQAProdMultistrange.h"
82
83 ClassImp(AliQAProdMultistrange)
84
85
86
87 //________________________________________________________________________
88 AliQAProdMultistrange::AliQAProdMultistrange() 
89   : AliAnalysisTaskSE(), 
90     fAnalysisType               ("ESD"), 
91     fESDtrackCuts               (0),
92     fCollidingSystem            ("PbPb"),
93     fPIDResponse                (0),
94     fkSDDSelectionOn            (kTRUE),
95     fkQualityCutZprimVtxPos     (kTRUE),
96     fkQualityCutNoTPConlyPrimVtx(kTRUE),
97     fkQualityCutTPCrefit        (kTRUE),
98     fkQualityCutnTPCcls         (kTRUE),
99     fkQualityCutPileup          (kTRUE),
100     fwithSDD                    (kTRUE),
101     fMinnTPCcls                 (0),  
102     fCentrLowLim                (0),
103     fCentrUpLim                 (0),
104     fCentrEstimator             (0),
105     fkUseCleaning               (0),
106     fVtxRange                   (0),
107     fMinPtCutOnDaughterTracks   (0),
108     fEtaCutOnDaughterTracks     (0),
109
110     
111     fCFContCascadeCuts(0)
112     
113
114 {
115   // Dummy Constructor
116 }
117
118
119 //________________________________________________________________________
120 AliQAProdMultistrange::AliQAProdMultistrange(const char *name) 
121   : AliAnalysisTaskSE(name),
122     fAnalysisType               ("ESD"), 
123     fESDtrackCuts               (0),
124     fCollidingSystem            ("PbPb"),
125     fPIDResponse                (0),
126     fkSDDSelectionOn            (kTRUE),
127     fkQualityCutZprimVtxPos     (kTRUE),
128     fkQualityCutNoTPConlyPrimVtx(kTRUE),
129     fkQualityCutTPCrefit        (kTRUE),
130     fkQualityCutnTPCcls         (kTRUE),
131     fkQualityCutPileup          (kTRUE),
132     fwithSDD                    (kTRUE),
133     fMinnTPCcls                 (0),
134     fCentrLowLim                (0),
135     fCentrUpLim                 (0),
136     fCentrEstimator             (0),
137     fkUseCleaning               (0),
138     fVtxRange                   (0),
139     fMinPtCutOnDaughterTracks   (0),
140     fEtaCutOnDaughterTracks     (0),
141
142
143     fCFContCascadeCuts(0)
144     
145
146 {
147   // Constructor
148   // Output slot #0 writes into a TList container (Cascade)
149   DefineOutput(1, AliCFContainer::Class());
150
151   AliLog::SetClassDebugLevel("AliQAProdMultistrange",1); // MN this should (?) enable only AliFatal
152 }
153
154
155 AliQAProdMultistrange::~AliQAProdMultistrange() {
156   //
157   // Destructor
158   //
159   // For all TH1, 2, 3 HnSparse and CFContainer are in the fListCascade TList.
160   // They will be deleted when fListCascade is deleted by the TSelector dtor
161   // Because of TList::SetOwner() ...
162   if (fCFContCascadeCuts && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) { delete fCFContCascadeCuts;     fCFContCascadeCuts = 0x0;  }
163   if (fESDtrackCuts)         { delete fESDtrackCuts;        fESDtrackCuts = 0x0; }
164 }
165
166
167
168 //________________________________________________________________________
169 void AliQAProdMultistrange::UserCreateOutputObjects() {
170   // Create histograms
171   // Called once
172
173  //-----------------------------------------------
174  // Particle Identification Setup (new PID object)
175  //-----------------------------------------------
176  AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
177  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
178  fPIDResponse = inputHandler->GetPIDResponse();
179
180
181  // Only used to get the number of primary reconstructed tracks
182  if (fAnalysisType == "ESD" && (! fESDtrackCuts )){
183    fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
184    //Printf("CheckCascade - ESDtrackCuts set up to 2010 std ITS-TPC cuts...");
185  }
186
187
188  //---------------------------------------------------
189  // Define the container for the topological variables
190  //---------------------------------------------------
191   if(! fCFContCascadeCuts) {
192       // Container meant to store all the relevant distributions corresponding to the cut variables.
193       // NB: overflow/underflow of variables on which we want to cut later should be 0!!! 
194       const Int_t  lNbSteps      =  4 ;
195       const Int_t  lNbVariables  =  21 ;
196       //Array for the number of bins in each dimension :
197       Int_t lNbBinsPerVar[lNbVariables] = {0};
198       lNbBinsPerVar[0]  = 25;     //DcaCascDaughters             :  [0.0,2.4,3.0]       -> Rec.Cut = 2.0;
199       lNbBinsPerVar[1]  = 25;     //DcaBachToPrimVertex          :  [0.0,0.24,100.0]    -> Rec.Cut = 0.01; 
200       lNbBinsPerVar[2]  = 30;     //CascCosineOfPointingAngle    :  [0.97,1.0]          -> Rec.Cut = 0.98;
201       lNbBinsPerVar[3]  = 40;     //CascRadius                   :  [0.0,3.9,1000.0]    -> Rec.Cut = 0.2;
202       lNbBinsPerVar[4]  = 30;     //InvMassLambdaAsCascDghter    :  [1.1,1.3]           -> Rec.Cut = 0.008;
203       lNbBinsPerVar[5]  = 20;     //DcaV0Daughters               :  [0.0,2.0]           -> Rec.Cut = 1.5;
204       lNbBinsPerVar[6]  = 201;    //V0CosineOfPointingAngleToPV  :  [0.89,1.0]          -> Rec.Cut = 0.9;
205       lNbBinsPerVar[7]  = 40;     //V0Radius                     :  [0.0,3.9,1000.0]    -> Rec.Cut = 0.2;
206       lNbBinsPerVar[8]  = 40;     //DcaV0ToPrimVertex            :  [0.0,0.39,110.0]    -> Rec.Cut = 0.01;  
207       lNbBinsPerVar[9]  = 25;     //DcaPosToPrimVertex           :  [0.0,0.24,100.0]    -> Rec.Cut = 0.05;
208       lNbBinsPerVar[10] = 25;     //DcaNegToPrimVertex           :  [0.0,0.24,100.0]    -> Rec.Cut = 0.05
209       lNbBinsPerVar[11] = 150;    //InvMassXi                    :   2-MeV/c2 bins
210       lNbBinsPerVar[12] = 120;    //InvMassOmega                 :   2-MeV/c2 bins
211       lNbBinsPerVar[13] = 100;    //XiTransvMom                  :  [0.0,10.0]
212       lNbBinsPerVar[14] = 110;    //Y(Xi)                        :   0.02 in rapidity units
213       lNbBinsPerVar[15] = 110;    //Y(Omega)                     :   0.02 in rapidity units
214       lNbBinsPerVar[16] = 112;    //Proper lenght of cascade       
215       lNbBinsPerVar[17] = 112;    //Proper lenght of V0
216       lNbBinsPerVar[18] = 201;    //V0CosineOfPointingAngleToXiV
217       lNbBinsPerVar[19] = 11;     //Centrality
218       lNbBinsPerVar[20] = 100;    //ESD track multiplicity
219       //define the container
220       fCFContCascadeCuts = new AliCFContainer("fCFContCascadeCuts","Container for Cascade cuts", lNbSteps, lNbVariables, lNbBinsPerVar );
221       //Setting the bin limits 
222       //0 -  DcaXiDaughters
223       Double_t *lBinLim0  = new Double_t[ lNbBinsPerVar[0] + 1 ];
224          for(Int_t i=0; i< lNbBinsPerVar[0]; i++) lBinLim0[i] = (Double_t)0.0 + (2.4 - 0.0)/(lNbBinsPerVar[0] - 1) * (Double_t)i;
225          lBinLim0[ lNbBinsPerVar[0] ] = 3.0;
226       fCFContCascadeCuts -> SetBinLimits(0, lBinLim0);
227       delete [] lBinLim0;
228       //1 - DcaToPrimVertexXi
229       Double_t *lBinLim1  = new Double_t[ lNbBinsPerVar[1] + 1 ];
230          for(Int_t i=0; i<lNbBinsPerVar[1]; i++) lBinLim1[i] = (Double_t)0.0 + (0.24  - 0.0)/(lNbBinsPerVar[1] - 1) * (Double_t)i;
231          lBinLim1[ lNbBinsPerVar[1] ] = 100.0;
232       fCFContCascadeCuts -> SetBinLimits(1, lBinLim1);
233       delete [] lBinLim1;
234       //2 - CascCosineOfPointingAngle 
235       fCFContCascadeCuts->SetBinLimits(2, 0.97, 1.);
236       //3 - CascRadius
237       Double_t *lBinLim3  = new Double_t[ lNbBinsPerVar[3]+1 ];
238          for(Int_t i=0; i< lNbBinsPerVar[3]; i++)   lBinLim3[i]  = (Double_t)0.0   + (3.9  - 0.0 )/(lNbBinsPerVar[3] - 1)  * (Double_t)i ;
239          lBinLim3[ lNbBinsPerVar[3] ] = 1000.0;
240       fCFContCascadeCuts -> SetBinLimits(3,  lBinLim3 );
241       delete [] lBinLim3;
242       //4 - InvMassLambdaAsCascDghter
243       fCFContCascadeCuts->SetBinLimits(4, 1.1, 1.13);
244       //5 - DcaV0Daughters
245       fCFContCascadeCuts -> SetBinLimits(5, 0., 2.);
246       //6 - V0CosineOfPointingAngleToPV
247       fCFContCascadeCuts -> SetBinLimits(6, 0.8, 1.001);
248       //7 - V0Radius
249       Double_t *lBinLim7 = new Double_t[ lNbBinsPerVar[7] + 1];
250          for(Int_t i=0; i< lNbBinsPerVar[7];i++) lBinLim7[i] = (Double_t)0.0 + (3.9 - 0.0)/(lNbBinsPerVar[7] - 1) * (Double_t)i;
251          lBinLim7[ lNbBinsPerVar[7] ] = 1000.0;
252       fCFContCascadeCuts -> SetBinLimits(7, lBinLim7);
253       delete [] lBinLim7;
254       //8 - DcaV0ToPrimVertex
255       Double_t *lBinLim8  = new Double_t[ lNbBinsPerVar[8]+1 ];
256          for(Int_t i=0; i< lNbBinsPerVar[8];i++)   lBinLim8[i]  = (Double_t)0.0   + (0.39  - 0.0 )/(lNbBinsPerVar[8]-1)  * (Double_t)i ;
257          lBinLim8[ lNbBinsPerVar[8]  ] = 100.0;
258       fCFContCascadeCuts -> SetBinLimits(8,  lBinLim8 );
259       delete [] lBinLim8;
260       //9 - DcaPosToPrimVertex
261       Double_t *lBinLim9  = new Double_t[ lNbBinsPerVar[9]+1 ];
262          for(Int_t i=0; i< lNbBinsPerVar[9];i++)   lBinLim9[i]  = (Double_t)0.0   + (0.24  - 0.0 )/(lNbBinsPerVar[9]-1)  * (Double_t)i ;
263          lBinLim9[ lNbBinsPerVar[9]  ] = 100.0;
264       fCFContCascadeCuts -> SetBinLimits(9,  lBinLim9 );
265       delete [] lBinLim9;
266       //10 - DcaNegToPrimVertex
267       Double_t *lBinLim10  = new Double_t[ lNbBinsPerVar[10]+1 ];
268          for(Int_t i=0; i< lNbBinsPerVar[10];i++)   lBinLim10[i]  = (Double_t)0.0   + (0.24  - 0.0 )/(lNbBinsPerVar[10]-1)  * (Double_t)i ;
269          lBinLim10[ lNbBinsPerVar[10]  ] = 100.0;
270       fCFContCascadeCuts -> SetBinLimits(10,  lBinLim10 );     
271       delete [] lBinLim10;
272       //11 - InvMassXi
273       fCFContCascadeCuts->SetBinLimits(11, 1.25, 1.40);
274       //12 - InvMassOmega
275       fCFContCascadeCuts->SetBinLimits(12, 1.62, 1.74);
276       //13 - XiTransvMom
277       fCFContCascadeCuts->SetBinLimits(13, 0.0, 10.0);
278       //14 - Y(Xi)
279       fCFContCascadeCuts->SetBinLimits(14, -1.1, 1.1);
280       //15 - Y(Omega)
281       fCFContCascadeCuts->SetBinLimits(15, -1.1, 1.1);
282       //16 - Proper time of cascade
283       Double_t *lBinLim16  = new Double_t[ lNbBinsPerVar[16]+1 ];
284          for(Int_t i=0; i< lNbBinsPerVar[16];i++) lBinLim16[i] = (Double_t) -1. + (110. + 1.0 ) / (lNbBinsPerVar[16] - 1) * (Double_t) i;
285          lBinLim16[ lNbBinsPerVar[16] ] = 2000.0;
286       fCFContCascadeCuts->SetBinLimits(16, lBinLim16);
287       //17 - Proper time of V0
288       fCFContCascadeCuts->SetBinLimits(17, lBinLim16);
289       //18 - V0CosineOfPointingAngleToXiV
290       fCFContCascadeCuts -> SetBinLimits(18, 0.8, 1.001);
291       //19
292       Double_t *lBinLim19  = new Double_t[ lNbBinsPerVar[19]+1 ];
293          for(Int_t i=3; i< lNbBinsPerVar[19]+1;i++)   lBinLim19[i]  = (Double_t)(i-1)*10.;
294          lBinLim19[0] = 0.0; 
295          lBinLim19[1] = 5.0; 
296          lBinLim19[2] = 10.0;
297       fCFContCascadeCuts->SetBinLimits(19,  lBinLim19 );     
298       delete [] lBinLim19;
299       //20
300       fCFContCascadeCuts->SetBinLimits(20, 0.0, 6000.0);
301       // Setting the number of steps : one for each cascade species (Xi-, Xi+ and Omega-, Omega+)
302       fCFContCascadeCuts->SetStepTitle(0, "#Xi^{-} candidates");
303       fCFContCascadeCuts->SetStepTitle(1, "#bar{#Xi}^{+} candidates");
304       fCFContCascadeCuts->SetStepTitle(2, "#Omega^{-} candidates");
305       fCFContCascadeCuts->SetStepTitle(3, "#bar{#Omega}^{+} candidates");
306       // Setting the variable title, per axis
307       fCFContCascadeCuts->SetVarTitle(0,  "Dca(cascade daughters) (cm)");
308       fCFContCascadeCuts->SetVarTitle(1,  "ImpactParamToPV(bachelor) (cm)");
309       fCFContCascadeCuts->SetVarTitle(2,  "cos(cascade PA)");
310       fCFContCascadeCuts->SetVarTitle(3,  "R_{2d}(cascade decay) (cm)");
311       fCFContCascadeCuts->SetVarTitle(4,  "M_{#Lambda}(as casc dghter) (GeV/c^{2})");
312       fCFContCascadeCuts->SetVarTitle(5,  "Dca(V0 daughters) in Xi (cm)");
313       fCFContCascadeCuts->SetVarTitle(6,  "cos(V0 PA) in cascade to PV");
314       fCFContCascadeCuts->SetVarTitle(7,  "R_{2d}(V0 decay) (cm)");
315       fCFContCascadeCuts->SetVarTitle(8,  "ImpactParamToPV(V0) (cm)");
316       fCFContCascadeCuts->SetVarTitle(9,  "ImpactParamToPV(Pos) (cm)");
317       fCFContCascadeCuts->SetVarTitle(10, "ImpactParamToPV(Neg) (cm)");
318       fCFContCascadeCuts->SetVarTitle(11, "Inv. Mass(Xi) (GeV/c^{2})");
319       fCFContCascadeCuts->SetVarTitle(12, "Inv. Mass(Omega) (GeV/c^{2})");
320       fCFContCascadeCuts->SetVarTitle(13, "pt(cascade) (GeV/c)");
321       fCFContCascadeCuts->SetVarTitle(14, "Y(Xi)");
322       fCFContCascadeCuts->SetVarTitle(15, "Y(Omega)");
323       fCFContCascadeCuts->SetVarTitle(16, "mL/p (cascade) (cm)");
324       fCFContCascadeCuts->SetVarTitle(17, "mL/p (V0) (cm)");
325       fCFContCascadeCuts->SetVarTitle(18,  "cos(V0 PA) in cascade to XiV");
326       fCFContCascadeCuts->SetVarTitle(19, "Centrality");
327       fCFContCascadeCuts->SetVarTitle(20, "ESD track multiplicity");
328   }
329
330 PostData(1, fCFContCascadeCuts);
331
332 }// end UserCreateOutputObjects
333
334
335 //________________________________________________________________________
336 void AliQAProdMultistrange::UserExec(Option_t *) {
337
338   //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
339   // Main loop (called for each event)
340   //%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
341
342   //-----------------------
343   //Define ESD/AOD handlers 
344   AliESDEvent *lESDevent = 0x0;
345   AliAODEvent *lAODevent = 0x0;
346
347   //---------------------
348   //Check the PIDresponse
349   if(!fPIDResponse) {
350        AliError("Cannot get pid response");
351        return;
352   }
353
354   //__________________________________________________
355   // After these lines we should have an ESD/AOD event
356
357   //---------------------------------------------------------
358   //Load the InputEvent and check, before any event selection   
359   //---------------------------------------------------------
360   Float_t  lPrimaryTrackMultiplicity = -1.;
361   AliCentrality* centrality;
362   if (fAnalysisType == "ESD") {
363       lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
364       if (!lESDevent) {
365           AliWarning("ERROR: lESDevent not available \n");
366           return;
367       }
368       if (fCollidingSystem == "PbPb") lPrimaryTrackMultiplicity = fESDtrackCuts->CountAcceptedTracks(lESDevent);
369       if (fCollidingSystem == "PbPb") centrality = lESDevent->GetCentrality();
370       
371   } else if (fAnalysisType == "AOD") {
372       lAODevent = dynamic_cast<AliAODEvent*>( InputEvent() );
373       if (!lAODevent) {
374           AliWarning("ERROR: lAODevent not available \n");
375           return;
376       }
377       if (fCollidingSystem == "PbPb") {
378           lPrimaryTrackMultiplicity = 0;
379           Int_t    nTrackMultiplicity = (InputEvent())->GetNumberOfTracks();
380           for (Int_t itrack = 0; itrack < nTrackMultiplicity; itrack++) {
381                AliAODTrack* track = lAODevent->GetTrack(itrack);
382                if (track->TestFilterBit(AliAODTrack::kTrkGlobalNoDCA)) lPrimaryTrackMultiplicity++; 
383           }
384       }
385       if (fCollidingSystem == "PbPb") centrality = lAODevent->GetCentrality();
386   } else {
387     Printf("Analysis type (ESD or AOD) not specified \n");
388     return;
389   }
390
391   //-----------------------------------------
392   // Centrality selection for PbPb collisions
393   //-----------------------------------------
394   Float_t lcentrality = 0.;
395   if (fCollidingSystem == "PbPb") { 
396        if (fkUseCleaning) lcentrality = centrality->GetCentralityPercentile(fCentrEstimator.Data());
397        else {
398            lcentrality = centrality->GetCentralityPercentileUnchecked(fCentrEstimator.Data());
399            if (centrality->GetQuality()>1) {
400                PostData(1, fCFContCascadeCuts);
401                return;
402            }
403        }
404        if (lcentrality<fCentrLowLim||lcentrality>=fCentrUpLim) { 
405            PostData(1, fCFContCascadeCuts);
406            return;
407        }
408   } else if (fCollidingSystem == "pp") lcentrality = 0.;
409
410
411   //----------------------------------------
412   // SDD selection for pp@2.76TeV collisions
413   //----------------------------------------
414   if (fCollidingSystem == "pp") {
415       if (fAnalysisType == "ESD") {
416           if (fkSDDSelectionOn) {
417               TString trcl = lESDevent->GetFiredTriggerClasses();
418               if      (fwithSDD) { if(!(trcl.Contains("ALLNOTRD"))) { PostData(1, fCFContCascadeCuts); return; } }
419               else if (!fwithSDD){ if((trcl.Contains("ALLNOTRD")))  { PostData(1, fCFContCascadeCuts); return; } }
420           }
421       } else if (fAnalysisType == "AOD") {
422           if (fkSDDSelectionOn) {
423               TString trcl = lAODevent->GetFiredTriggerClasses();
424               if      (fwithSDD)  { if(!(trcl.Contains("ALLNOTRD"))) { PostData(1, fCFContCascadeCuts); return; } }
425               else if (!fwithSDD) { if((trcl.Contains("ALLNOTRD")))  { PostData(1, fCFContCascadeCuts); return; } }
426           }
427       }
428   }
429
430   //--------------------------------------------
431   // Physics selection for pp@2.76TeV collisions
432   //--------------------------------------------
433   // - moved to the runGrid for the PbPb collisions
434   if (fCollidingSystem == "pp") {
435       if (fAnalysisType == "ESD") {
436           UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
437           Bool_t isSelected = 0;
438           isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
439           if(! isSelected){ PostData(1, fCFContCascadeCuts); return; }
440       } else if (fAnalysisType == "AOD") {
441           UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
442           Bool_t isSelected = 0;
443           isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
444           if(! isSelected){ PostData(1, fCFContCascadeCuts); return; }
445       }
446   }
447
448   //------------------------------
449   // Well-established PV selection
450   //------------------------------
451   Double_t lBestPrimaryVtxPos[3] = {-100.0, -100.0, -100.0}; 
452   Double_t lMagneticField        = -10.;
453   if (fAnalysisType == "ESD") {
454        const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();   
455        const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();  
456        if (fkQualityCutNoTPConlyPrimVtx) {
457            const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
458            if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() ){
459                AliWarning(" No SPD prim. vertex nor prim. Tracking vertex ... return !");
460                PostData(1, fCFContCascadeCuts);
461                return;
462            }
463        }
464        lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
465        lMagneticField = lESDevent->GetMagneticField( );
466   } else if (fAnalysisType == "AOD") {
467        const AliAODVertex *lPrimaryBestAODVtx = lAODevent->GetPrimaryVertex();
468        if (!lPrimaryBestAODVtx){
469            AliWarning("No prim. vertex in AOD... return!");
470            PostData(1, fCFContCascadeCuts);
471            return;
472        }
473        lPrimaryBestAODVtx->GetXYZ( lBestPrimaryVtxPos );
474        lMagneticField = lAODevent->GetMagneticField();  
475   }
476
477   //------------------------------------------
478   // Pilup selection for pp@2.76TeV collisions
479   //------------------------------------------
480   if (fCollidingSystem == "pp") { 
481       if (fAnalysisType == "ESD") {
482           if (fkQualityCutPileup) { if(lESDevent->IsPileupFromSPD()){ PostData(1, fCFContCascadeCuts); return; } }
483       } else if (fAnalysisType == "AOD") {
484           if (fkQualityCutPileup) { if(lAODevent->IsPileupFromSPD()){ PostData(1, fCFContCascadeCuts); return; } }
485       }
486   }
487
488   //----------------------------
489   // Vertex Z position selection
490   //----------------------------
491   if (fkQualityCutZprimVtxPos) {
492       if (TMath::Abs(lBestPrimaryVtxPos[2]) > fVtxRange ) {
493           PostData(1, fCFContCascadeCuts);
494           return;
495       }
496   }
497
498
499
500   //////////////////////////////
501   // CASCADE RECONSTRUCTION PART
502   //////////////////////////////
503
504   //%%%%%%%%%%%%%
505   // Cascade loop
506   Int_t ncascades = 0;
507   if      (fAnalysisType == "ESD") ncascades = lESDevent->GetNumberOfCascades();
508   else if (fAnalysisType == "AOD") ncascades = lAODevent->GetNumberOfCascades();
509
510   for (Int_t iXi = 0; iXi < ncascades; iXi++) {// This is the begining of the Cascade loop (ESD or AOD)
511            
512     // -------------------------------------
513     // - Initialisation of the local variables that will be needed for ESD/AOD
514     // -- Container variables (1st round)
515     Double_t lDcaXiDaughters              = -1. ;                   //[Container]
516     Double_t lXiCosineOfPointingAngle     = -1. ;                   //[Container]
517     Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 };             //Useful to define other variables: radius fid. vol., ctau, etc. for cascade
518     Double_t lXiRadius                    = -1000. ;                //[Container]
519     UShort_t lPosTPCClusters              = -1;                     //To check the quality of the tracks. For ESD only ...
520     UShort_t lNegTPCClusters              = -1;                     //To check the quality of the tracks. For ESD only ...
521     UShort_t lBachTPCClusters             = -1;                     //To check the quality of the tracks. For ESD only ...
522     Double_t lInvMassLambdaAsCascDghter   = 0.;                     //[Container]
523     Double_t lDcaV0DaughtersXi            = -1.;                    //[Container]
524     Double_t lDcaBachToPrimVertexXi       = -1.;                    //[Container]
525     Double_t lDcaV0ToPrimVertexXi         = -1.;                    //[Container]
526     Double_t lDcaPosToPrimVertexXi        = -1.;                    //[Container]
527     Double_t lDcaNegToPrimVertexXi        = -1.;                    //[Container]
528     Double_t lV0CosineOfPointingAngle     = -1.;                    //[Container]
529     Double_t lV0toXiCosineOfPointingAngle = -1.;                    //[Container] 
530     Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. };             //Useful to define other variables: radius fid. vol., ctau, etc. for VO 
531     Double_t lV0RadiusXi                  = -1000.0;                //[Container]
532     Double_t lV0quality                   = 0.;                     //  ??
533     Double_t lInvMassXiMinus              = 0.;                     //[Container]
534     Double_t lInvMassXiPlus               = 0.;                     //[Container]
535     Double_t lInvMassOmegaMinus           = 0.;                     //[Container]
536     Double_t lInvMassOmegaPlus            = 0.;                     //[Container]
537     // -- PID treatment
538     Bool_t   lIsBachelorKaonForTPC = kFALSE; 
539     Bool_t   lIsBachelorPionForTPC = kFALSE; 
540     Bool_t   lIsNegPionForTPC      = kFALSE; 
541     Bool_t   lIsPosPionForTPC      = kFALSE; 
542     Bool_t   lIsNegProtonForTPC    = kFALSE; 
543     Bool_t   lIsPosProtonForTPC    = kFALSE; 
544     // -- More container variables and quality checks
545     Double_t lXiMomX           = 0.;                               //Useful to define other variables: lXiTransvMom, lXiTotMom
546     Double_t lXiMomY           = 0.;                               //Useful to define other variables: lXiTransvMom, lXiTotMom
547     Double_t lXiMomZ           = 0.;                               //Useful to define other variables: lXiTransvMom, lXiTotMom
548     Double_t lXiTransvMom      = 0.;                               //[Container]
549     Double_t lXiTotMom         = 0.;                               //Useful to define other variables: cTau
550     Double_t lV0PMomX          = 0.;                               //Useful to define other variables: lV0TotMom, lpTrackTransvMom
551     Double_t lV0PMomY          = 0.;                               //Useful to define other variables: lV0TotMom, lpTrackTransvMom
552     Double_t lV0PMomZ          = 0.;                               //Useful to define other variables: lV0TotMom, lpTrackTransvMom
553     Double_t lV0NMomX          = 0.;                               //Useful to define other variables: lV0TotMom, lnTrackTransvMom
554     Double_t lV0NMomY          = 0.;                               //Useful to define other variables: lV0TotMom, lnTrackTransvMom
555     Double_t lV0NMomZ          = 0.;                               //Useful to define other variables: lV0TotMom, lnTrackTransvMom
556     Double_t lV0TotMom         = 0.;                               //Useful to define other variables: lctauV0
557     Double_t lBachMomX         = 0.;                               //Useful to define other variables: lBachTransvMom
558     Double_t lBachMomY         = 0.;                               //Useful to define other variables: lBachTransvMom
559     Double_t lBachMomZ         = 0.;                               //Useful to define other variables: lBachTransvMom
560     Double_t lBachTransvMom    = 0.;                               //Selection on the min bachelor pT
561     Double_t lpTrackTransvMom  = 0.;                               //Selection on the min bachelor pT
562     Double_t lnTrackTransvMom  = 0.;                               //Selection on the min bachelor pT
563     Short_t  lChargeXi         = -2;                               //Useful to select the particles based on the charge
564     Double_t lRapXi            = -20.0;                            //[Container]
565     Double_t lRapOmega         = -20.0;                            //[Container]
566     Float_t  etaBach           = 0.;                               //Selection on the eta range
567     Float_t  etaPos            = 0.;                               //Selection on the eta range
568     Float_t  etaNeg            = 0.;                               //Selection on the eta range
569     // --  variables for the AliCFContainer dedicated to cascade cut optmisiation: ESD and AOD 
570     if (fAnalysisType == "ESD") { 
571   
572           // -------------------------------------
573           // - Load the cascades from the handler 
574           AliESDcascade *xi = lESDevent->GetCascade(iXi);
575           if (!xi) continue;
576         
577           // ---------------------------------------------------------------------------
578           // - Assigning the necessary variables for specific AliESDcascade data members        
579           lV0quality = 0.;
580           xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
581           lDcaXiDaughters          = xi->GetDcaXiDaughters();
582           lXiCosineOfPointingAngle = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
583                                        // Take care : the best available vertex should be used (like in AliCascadeVertexer)
584           xi->GetXYZcascade( lPosXi[0],  lPosXi[1], lPosXi[2] ); 
585           lXiRadius             = TMath::Sqrt( lPosXi[0]*lPosXi[0]  +  lPosXi[1]*lPosXi[1] );
586                 
587           // -------------------------------------------------------------------------------------------------------------------------------
588           // - Around the tracks : Bach + V0 (ESD). Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
589           UInt_t lIdxPosXi      = (UInt_t) TMath::Abs( xi->GetPindex() );
590           UInt_t lIdxNegXi      = (UInt_t) TMath::Abs( xi->GetNindex() );
591           UInt_t lBachIdx       = (UInt_t) TMath::Abs( xi->GetBindex() );
592                                     // Care track label can be negative in MC production (linked with the track quality)
593                                     // However = normally, not the case for track index ...
594           // - Rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
595           if (lBachIdx == lIdxNegXi) continue;    
596           if (lBachIdx == lIdxPosXi) continue; 
597           // - Get the track for the daughters
598           AliESDtrack *pTrackXi         = lESDevent->GetTrack( lIdxPosXi );
599           AliESDtrack *nTrackXi         = lESDevent->GetTrack( lIdxNegXi );
600           AliESDtrack *bachTrackXi      = lESDevent->GetTrack( lBachIdx );
601           if (!pTrackXi || !nTrackXi || !bachTrackXi )  continue;
602           // - Get the TPCnumber of cluster for the daughters
603           lPosTPCClusters   = pTrackXi->GetTPCNcls();
604           lNegTPCClusters   = nTrackXi->GetTPCNcls();
605           lBachTPCClusters  = bachTrackXi->GetTPCNcls();
606       
607           // ------------------------------------
608           // - Rejection of a poor quality tracks
609           if (fkQualityCutTPCrefit) {
610                 // 1 - Poor quality related to TPCrefit
611                 ULong_t pStatus    = pTrackXi->GetStatus();
612                 ULong_t nStatus    = nTrackXi->GetStatus();
613                 ULong_t bachStatus = bachTrackXi->GetStatus();
614                 if ((pStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning(" V0 Pos. track has no TPCrefit ... continue!"); continue; }
615                 if ((nStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning(" V0 Neg. track has no TPCrefit ... continue!"); continue; }
616                 if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning(" Bach.   track has no TPCrefit ... continue!"); continue; }
617           }
618           if (fkQualityCutnTPCcls) {
619                 // 2 - Poor quality related to TPC clusters
620                 if (lPosTPCClusters  < fMinnTPCcls) { AliWarning(" V0 Pos. track has less than minn TPC clusters ... continue!"); continue; }
621                 if (lNegTPCClusters  < fMinnTPCcls) { AliWarning(" V0 Neg. track has less than minn TPC clusters ... continue!"); continue; }
622                 if (lBachTPCClusters < fMinnTPCcls) { AliWarning(" Bach.   track has less than minn TPC clusters ... continue!"); continue; }
623           }
624
625           // ------------------------------
626           etaPos  = pTrackXi->Eta();             
627           etaNeg  = nTrackXi->Eta();
628           etaBach = bachTrackXi->Eta();
629           lInvMassLambdaAsCascDghter = xi->GetEffMass();
630           lDcaV0DaughtersXi          = xi->GetDcaV0Daughters(); 
631           lV0CosineOfPointingAngle   = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] );
632           lDcaV0ToPrimVertexXi       = xi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lBestPrimaryVtxPos[2] ); 
633           lDcaBachToPrimVertexXi     = TMath::Abs( bachTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1],     lMagneticField) ); 
634           xi->GetXYZ( lPosV0Xi[0],  lPosV0Xi[1], lPosV0Xi[2] ); 
635           lV0RadiusXi                = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0] + lPosV0Xi[1]*lPosV0Xi[1] );
636           lDcaPosToPrimVertexXi      = TMath::Abs( pTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField ) ); 
637           lDcaNegToPrimVertexXi      = TMath::Abs( nTrackXi->GetD( lBestPrimaryVtxPos[0], lBestPrimaryVtxPos[1], lMagneticField ) ); 
638         
639            //----------------------------------------------------------------------------------------------------       
640            // - Around effective masses. Change mass hypotheses to cover all the possibilities:  Xi-/+, Omega -/+
641            if (bachTrackXi->Charge() < 0) {
642                 // Calculate the effective mass of the Xi- candidate. pdg code 3312 = Xi-
643                 lV0quality = 0.;
644                 xi->ChangeMassHypothesis(lV0quality , 3312);    
645                 lInvMassXiMinus = xi->GetEffMassXi();
646                 // Calculate the effective mass of the Xi- candidate. pdg code 3334 = Omega-
647                 lV0quality = 0.;
648                 xi->ChangeMassHypothesis(lV0quality , 3334);    
649                 lInvMassOmegaMinus = xi->GetEffMassXi();
650                 // Back to default hyp.                 
651                 lV0quality = 0.;
652                 xi->ChangeMassHypothesis(lV0quality , 3312);
653            }// end if negative bachelor
654            if ( bachTrackXi->Charge() >  0 ) {
655                 // Calculate the effective mass of the Xi+ candidate. pdg code -3312 = Xi+
656                 lV0quality = 0.;
657                 xi->ChangeMassHypothesis(lV0quality , -3312);   
658                 lInvMassXiPlus = xi->GetEffMassXi();
659                 // Calculate the effective mass of the Xi+ candidate. pdg code -3334  = Omega+
660                 lV0quality = 0.;
661                 xi->ChangeMassHypothesis(lV0quality , -3334);   
662                 lInvMassOmegaPlus = xi->GetEffMassXi();
663                 // Back to "default" hyp.
664                 lV0quality = 0.;
665                 xi->ChangeMassHypothesis(lV0quality , -3312); 
666            }// end if positive bachelor
667
668            // ----------------------------------------------    
669            // - TPC PID : 3-sigma bands on Bethe-Bloch curve
670            // Bachelor
671            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
672            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
673            // Negative V0 daughter
674            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kPion   )) < 4) lIsNegPionForTPC   = kTRUE;
675            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
676            // Positive V0 daughter
677            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kPion   )) < 4) lIsPosPionForTPC   = kTRUE;
678            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
679         
680            // ------------------------------
681            // - Miscellaneous pieces of info that may help regarding data quality assessment.
682            xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
683            lXiTransvMom = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY );
684            lXiTotMom    = TMath::Sqrt( lXiMomX*lXiMomX + lXiMomY*lXiMomY + lXiMomZ*lXiMomZ );
685            xi->GetNPxPyPz(lV0NMomX,lV0NMomY,lV0NMomZ);
686            xi->GetPPxPyPz(lV0PMomX,lV0PMomY,lV0PMomZ);
687            lV0TotMom = TMath::Sqrt(TMath::Power(lV0NMomX+lV0PMomX,2)+TMath::Power(lV0NMomY+lV0PMomY,2)+TMath::Power(lV0NMomZ+lV0PMomZ,2));      
688            xi->GetBPxPyPz( lBachMomX, lBachMomY, lBachMomZ );
689            lBachTransvMom  = TMath::Sqrt( lBachMomX*lBachMomX + lBachMomY*lBachMomY );
690            lnTrackTransvMom = TMath::Sqrt( lV0NMomX*lV0NMomX + lV0NMomY*lV0NMomY );
691            lpTrackTransvMom = TMath::Sqrt( lV0PMomX*lV0PMomX + lV0PMomY*lV0PMomY );
692            lChargeXi = xi->Charge();
693            lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
694            lRapXi    = xi->RapXi();
695            lRapOmega = xi->RapOmega();
696         
697     } else if (fAnalysisType == "AOD") {
698
699            // -------------------------------------
700            // - Load the cascades from the handler      
701            const AliAODcascade *xi = lAODevent->GetCascade(iXi);
702            if (!xi) continue;
703                 
704            //----------------------------------------------------------------------------        
705            // - Assigning the necessary variables for specific AliESDcascade data members  
706            lDcaXiDaughters              = xi->DcaXiDaughters();
707            lXiCosineOfPointingAngle     = xi->CosPointingAngleXi( lBestPrimaryVtxPos[0], 
708                                                           lBestPrimaryVtxPos[1], 
709                                                           lBestPrimaryVtxPos[2] );
710            lPosXi[0] = xi->DecayVertexXiX();
711            lPosXi[1] = xi->DecayVertexXiY();
712            lPosXi[2] = xi->DecayVertexXiZ();
713            lXiRadius = TMath::Sqrt( lPosXi[0]*lPosXi[0]  +  lPosXi[1]*lPosXi[1] );      
714
715            //-------------------------------------------------------------------------------------------------------------------------------
716            // - Around the tracks: Bach + V0 (ESD). Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
717            AliAODTrack *pTrackXi    = dynamic_cast<AliAODTrack*>( xi->GetDaughter(0) );
718            AliAODTrack *nTrackXi    = dynamic_cast<AliAODTrack*>( xi->GetDaughter(1) );
719            AliAODTrack *bachTrackXi = dynamic_cast<AliAODTrack*>( xi->GetDecayVertexXi()->GetDaughter(0) );
720            if (!pTrackXi || !nTrackXi || !bachTrackXi ) continue;
721            UInt_t lIdxPosXi  = (UInt_t) TMath::Abs( pTrackXi->GetID() );  
722            UInt_t lIdxNegXi  = (UInt_t) TMath::Abs( nTrackXi->GetID() );
723            UInt_t lBachIdx   = (UInt_t) TMath::Abs( bachTrackXi->GetID() );
724                                 // Care track label can be negative in MC production (linked with the track quality)
725                                 // However = normally, not the case for track index ...
726
727            // - Rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
728            if (lBachIdx == lIdxNegXi) continue; 
729            if (lBachIdx == lIdxPosXi) continue;
730            // - Get the TPCnumber of cluster for the daughters
731            lPosTPCClusters   = pTrackXi->GetTPCNcls(); // FIXME: Is this ok? or something like in LambdaK0PbPb task AOD?
732            lNegTPCClusters   = nTrackXi->GetTPCNcls();
733            lBachTPCClusters  = bachTrackXi->GetTPCNcls();
734
735            // ------------------------------------
736            // - Rejection of a poor quality tracks
737            if (fkQualityCutTPCrefit) {
738                 // - Poor quality related to TPCrefit
739                 if (!(pTrackXi->IsOn(AliAODTrack::kTPCrefit)))    { AliWarning(" V0 Pos. track has no TPCrefit ... continue!"); continue; }
740                 if (!(nTrackXi->IsOn(AliAODTrack::kTPCrefit)))    { AliWarning(" V0 Neg. track has no TPCrefit ... continue!"); continue; }
741                 if (!(bachTrackXi->IsOn(AliAODTrack::kTPCrefit))) { AliWarning(" Bach.   track has no TPCrefit ... continue!"); continue; }
742            }
743            if (fkQualityCutnTPCcls) {
744                 // - Poor quality related to TPC clusters
745                 if (lPosTPCClusters  < fMinnTPCcls) continue; 
746                 if (lNegTPCClusters  < fMinnTPCcls) continue; 
747                 if (lBachTPCClusters < fMinnTPCcls) continue; 
748            }
749   
750            // ------------------------------------------------------------------------------------------------------------------------------
751            // - Around the tracks: Bach + V0 (AOD). Necessary variables for AODcascade data members coming from the AODv0 part (inheritance)
752            etaPos  = pTrackXi->Eta();
753            etaNeg  = nTrackXi->Eta();
754            etaBach = bachTrackXi->Eta();
755            lChargeXi                    = xi->ChargeXi();
756            if ( lChargeXi < 0)  lInvMassLambdaAsCascDghter      = xi->MassLambda();
757            else                 lInvMassLambdaAsCascDghter      = xi->MassAntiLambda();
758            lDcaV0DaughtersXi            = xi->DcaV0Daughters(); 
759            lDcaV0ToPrimVertexXi                 = xi->DcaV0ToPrimVertex();
760            lDcaBachToPrimVertexXi               = xi->DcaBachToPrimVertex(); 
761            lPosV0Xi[0] = xi->DecayVertexV0X();
762            lPosV0Xi[1] = xi->DecayVertexV0Y();
763            lPosV0Xi[2] = xi->DecayVertexV0Z(); 
764            lV0RadiusXi  = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0]  +  lPosV0Xi[1]*lPosV0Xi[1] );
765            lV0CosineOfPointingAngle     = xi->CosPointingAngle( lBestPrimaryVtxPos ); 
766            lDcaPosToPrimVertexXi        = xi->DcaPosToPrimVertex(); 
767            lDcaNegToPrimVertexXi        = xi->DcaNegToPrimVertex(); 
768
769            // ---------------------------------------------------------------------------------------------------       
770            // - Around effective masses. Change mass hypotheses to cover all the possibilities:  Xi-/+, Omega -/+
771            if ( lChargeXi < 0 )         lInvMassXiMinus         = xi->MassXi();
772            if ( lChargeXi > 0 )         lInvMassXiPlus          = xi->MassXi();
773            if ( lChargeXi < 0 )         lInvMassOmegaMinus      = xi->MassOmega();
774            if ( lChargeXi > 0 )         lInvMassOmegaPlus       = xi->MassOmega();
775
776            // ----------------------------------------------
777            // - TPC PID : 3-sigma bands on Bethe-Bloch curve
778            // Bachelor
779            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
780            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
781            // Negative V0 daughter
782            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kPion   )) < 4) lIsNegPionForTPC   = kTRUE;
783            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
784            // Positive V0 daughter
785            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kPion   )) < 4) lIsPosPionForTPC   = kTRUE;
786            if (TMath::Abs(fPIDResponse->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
787
788            //---------------------------------
789            // - Extra info for QA (AOD)
790            // Miscellaneous pieces of info that may help regarding data quality assessment.
791            // Cascade transverse and total momentum     
792            lXiMomX = xi->MomXiX();
793            lXiMomY = xi->MomXiY();
794            lXiMomZ = xi->MomXiZ();
795            lXiTransvMom         = TMath::Sqrt( lXiMomX*lXiMomX   + lXiMomY*lXiMomY );
796            lXiTotMom    = TMath::Sqrt( lXiMomX*lXiMomX   + lXiMomY*lXiMomY   + lXiMomZ*lXiMomZ );
797            Double_t lV0MomX = xi->MomV0X();
798            Double_t lV0MomY = xi->MomV0Y();
799            Double_t lV0MomZ = xi->MomV0Z();
800            lV0TotMom = TMath::Sqrt(TMath::Power(lV0MomX,2)+TMath::Power(lV0MomY,2)+TMath::Power(lV0MomZ,2));
801            lBachMomX = xi->MomBachX();
802            lBachMomY = xi->MomBachY();
803            lBachMomZ = xi->MomBachZ();          
804            lBachTransvMom  = TMath::Sqrt( lBachMomX*lBachMomX   + lBachMomY*lBachMomY );
805            lV0NMomX = xi->MomNegX();
806            lV0NMomY = xi->MomNegY();
807            lV0PMomX = xi->MomPosX();
808            lV0PMomY = xi->MomPosY(); 
809            lnTrackTransvMom = TMath::Sqrt( lV0NMomX*lV0NMomX   + lV0NMomY*lV0NMomY );
810            lpTrackTransvMom = TMath::Sqrt( lV0PMomX*lV0PMomX   + lV0PMomY*lV0PMomY );
811            lV0toXiCosineOfPointingAngle = xi->CosPointingAngle( xi->GetDecayVertexXi() );
812            lRapXi    = xi->RapXi();
813            lRapOmega = xi->RapOmega();
814
815     }// end of AOD treatment
816
817     //---------------------------------------
818     // Cut on pt of the three daughter tracks
819     if (lBachTransvMom<fMinPtCutOnDaughterTracks) continue;
820     if (lpTrackTransvMom<fMinPtCutOnDaughterTracks) continue;
821     if (lnTrackTransvMom<fMinPtCutOnDaughterTracks) continue;
822  
823     //---------------------------------------------------
824     // Cut on pseudorapidity of the three daughter tracks
825     if (TMath::Abs(etaBach)>fEtaCutOnDaughterTracks) continue;
826     if (TMath::Abs(etaPos)>fEtaCutOnDaughterTracks) continue;
827     if (TMath::Abs(etaNeg)>fEtaCutOnDaughterTracks) continue;
828
829     //----------------------------------
830     // Calculate proper time for cascade
831     Double_t cascadeMass = 0.;
832     if ( ( (lChargeXi<0) && lIsBachelorPionForTPC && lIsPosProtonForTPC && lIsNegPionForTPC ) ||
833          ( (lChargeXi>0) && lIsBachelorPionForTPC && lIsNegProtonForTPC && lIsPosPionForTPC )  ) cascadeMass = 1.321;
834     if ( ( (lChargeXi<0) && lIsBachelorKaonForTPC   && lIsPosProtonForTPC    && lIsNegPionForTPC ) ||
835          ( (lChargeXi>0) && lIsBachelorKaonForTPC   && lIsNegProtonForTPC    && lIsPosPionForTPC )  ) cascadeMass = 1.672; 
836     Double_t lctau =  TMath::Sqrt(TMath::Power((lPosXi[0]-lBestPrimaryVtxPos[0]),2)+TMath::Power((lPosXi[1]-lBestPrimaryVtxPos[1]),2)+TMath::Power(( lPosXi[2]-lBestPrimaryVtxPos[2]),2));
837     if (lXiTotMom!=0)         lctau = lctau*cascadeMass/lXiTotMom;
838     else lctau = -1.;
839     // Calculate proper time for Lambda (reconstructed)
840     Float_t lambdaMass = 1.115683; // PDG mass
841     Float_t distV0Xi =  TMath::Sqrt(TMath::Power((lPosV0Xi[0]-lPosXi[0]),2)+TMath::Power((lPosV0Xi[1]-lPosXi[1]),2)+TMath::Power((lPosV0Xi[2]-lPosXi[2]),2));
842     Float_t lctauV0 = -1.;
843     if (lV0TotMom!=0) lctauV0 = distV0Xi*lambdaMass/lV0TotMom;
844
845                 
846     // Fill the AliCFContainer (optimisation of topological selections)
847     Double_t lContainerCutVars[21] = {0.0};
848     lContainerCutVars[0]  = lDcaXiDaughters;
849     lContainerCutVars[1]  = lDcaBachToPrimVertexXi;
850     lContainerCutVars[2]  = lXiCosineOfPointingAngle;
851     lContainerCutVars[3]  = lXiRadius;
852     lContainerCutVars[4]  = lInvMassLambdaAsCascDghter;
853     lContainerCutVars[5]  = lDcaV0DaughtersXi;
854     lContainerCutVars[6]  = lV0CosineOfPointingAngle;
855     lContainerCutVars[7]  = lV0RadiusXi;
856     lContainerCutVars[8]  = lDcaV0ToPrimVertexXi;
857     lContainerCutVars[9]  = lDcaPosToPrimVertexXi;
858     lContainerCutVars[10] = lDcaNegToPrimVertexXi;
859     lContainerCutVars[13] = lXiTransvMom;
860     lContainerCutVars[16] = lctau;
861     lContainerCutVars[17] = lctauV0;
862     lContainerCutVars[18] = lV0toXiCosineOfPointingAngle;
863     lContainerCutVars[19] = lcentrality;
864     lContainerCutVars[20] = lPrimaryTrackMultiplicity;
865     if ( lChargeXi < 0 ) {
866          lContainerCutVars[11] = lInvMassXiMinus;
867          lContainerCutVars[12] = lInvMassOmegaMinus;
868          lContainerCutVars[14] = lRapXi;
869          lContainerCutVars[15] = -1.;
870          if (lIsBachelorPionForTPC && lIsPosProtonForTPC && lIsNegPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,0); // for Xi-
871          lContainerCutVars[11] = lInvMassXiMinus;
872          lContainerCutVars[12] = lInvMassOmegaMinus;
873          lContainerCutVars[14] = -1.;
874          lContainerCutVars[15] = lRapOmega;
875          if (lIsBachelorKaonForTPC && lIsPosProtonForTPC && lIsNegPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,2); // for Omega-
876     } else {
877          lContainerCutVars[11] = lInvMassXiPlus;
878          lContainerCutVars[12] = lInvMassOmegaPlus;
879          lContainerCutVars[14] = lRapXi;
880          lContainerCutVars[15] = -1.;
881          if (lIsBachelorPionForTPC && lIsNegProtonForTPC && lIsPosPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,1); // for Xi+
882          lContainerCutVars[11] = lInvMassXiPlus;
883          lContainerCutVars[12] = lInvMassOmegaPlus;
884          lContainerCutVars[14] = -1.;
885          lContainerCutVars[15] = lRapOmega;
886          if (lIsBachelorKaonForTPC && lIsNegProtonForTPC && lIsPosPionForTPC) fCFContCascadeCuts->Fill(lContainerCutVars,3); // for Omega+ 
887     }
888
889
890   }// end of the Cascade loop (ESD or AOD)
891     
892   
893   // Post output data.
894   PostData(1, fCFContCascadeCuts); 
895 }
896
897 //________________________________________________________________________
898 void AliQAProdMultistrange::Terminate(Option_t *) 
899 {
900
901 }