PWG2/SPECTRA/AliAnalysisTaskCheckPerformanceCascade.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 //               AliAnalysisTaskCheckPerformanceCascade class
  16 //            This task is for a performance study of cascade identification.
  17 //            It works with MC info and ESD.
  18 //              Use with AOD tree = under development
  19 //            Origin   : A.Maire Mar2009, antonin.maire@ires.in2p3.fr
  20 //            Modified : A.Maire Nov2010, antonin.maire@ires.in2p3.fr
  21 //-----------------------------------------------------------------
  22
  23
  24 #include <Riostream.h>
  25
  26 #include "TList.h"
  27 #include "TFile.h"
  28 #include "TH1F.h"
  29 #include "TH2F.h"
  30 #include "TVector3.h"
  31 #include "TCanvas.h"
  32 #include "TParticle.h"
  33 #include "TMath.h"
  34
  35 #include "AliLog.h"
  36 #include "AliHeader.h"
  37 #include "AliMCEvent.h"
  38 #include "AliStack.h"
  39 #include "AliMultiplicity.h"
  40 #include "AliInputEventHandler.h"
  41 #include "AliAnalysisManager.h"
  42
  43 #include "AliCFContainer.h"
  44 #include "AliESDpid.h"
  45 #include "AliESDtrackCuts.h"
  46 //   #include "AliV0vertexer.h"
  47 //   #include "AliCascadeVertexer.h"
  48 #include "AliESDEvent.h"
  49 #include "AliESDcascade.h"
  50
  51 #include "AliAODEvent.h"
  52
  53 #include "AliAnalysisTaskCheckPerformanceCascade.h"
  54
  55 ClassImp(AliAnalysisTaskCheckPerformanceCascade)
  56
  57
  58
  59      //_____Dummy constructor________________________________________________________________
  60 AliAnalysisTaskCheckPerformanceCascade::AliAnalysisTaskCheckPerformanceCascade()
  61 : AliAnalysisTaskSE(), // <- take care to AliAnalysisTask( empty )
  62   fDebugCascade(0), fAnalysisType("ESD"), fTriggerMaskType("kMB"), fCollidingSystems(0), fESDpid(0), fESDtrackCuts(0), /*fPaveTextBookKeeping(0),*/
  63     fkRerunV0CascVertexers         (0),
  64     fkQualityCutZprimVtxPos        (kTRUE),
  65     fkRejectEventPileUp            (kTRUE),
  66     fkQualityCutNoTPConlyPrimVtx   (kTRUE),
  67     fkQualityCutTPCrefit           (kTRUE),
  68     fkQualityCut80TPCcls           (kTRUE),
  69     fkIsDataRecoWith1PadTPCCluster (kTRUE),
  70     fkExtraSelections              (0),
  71
  72         // - Cascade part initialisation
  73 fListHistCascade(0),
  74     fHistMCTrackMultiplicity(0),
  75        // - Resolution of the multiplicity estimator
  76     f2dHistRecoPrimTrckMultVsMCMult(0), f2dHistRecoEstimateMultVsMCMult(0),
  77
  78     fHistEtaGenProton(0),
  79     fHistEtaGenAntiProton(0),
  80
  81    // Xi-
  82    fHistEtaGenCascXiMinus(0),
  83    f2dHistGenPtVsGenYGenXiMinus(0),
  84
  85     fHistThetaGenCascXiMinus(0),
  86     f2dHistGenPtVsGenYFdblXiMinus(0),
  87
  88     fHistThetaLambdaXiMinus(0),
  89     fHistThetaBachXiMinus(0),
  90
  91     fHistThetaMesDghterXiMinus(0),
  92     fHistThetaBarDghterXiMinus(0),
  93
  94     fHistPtBachXiMinus(0),
  95     fHistPtMesDghterXiMinus(0),
  96     fHistPtBarDghterXiMinus(0),
  97
  98
  99    // Xi+
 100    fHistEtaGenCascXiPlus(0),
 101    f2dHistGenPtVsGenYGenXiPlus(0),
 102
 103     fHistThetaGenCascXiPlus(0),
 104     f2dHistGenPtVsGenYFdblXiPlus(0),
 105
 106     fHistThetaLambdaXiPlus(0),
 107     fHistThetaBachXiPlus(0),
 108
 109     fHistThetaMesDghterXiPlus(0),
 110     fHistThetaBarDghterXiPlus(0),
 111
 112     fHistPtBachXiPlus(0),
 113     fHistPtMesDghterXiPlus(0),
 114     fHistPtBarDghterXiPlus(0),
 115
 116    // Omega-
 117    fHistEtaGenCascOmegaMinus(0),
 118    f2dHistGenPtVsGenYGenOmegaMinus(0),
 119
 120     fHistThetaGenCascOmegaMinus(0),
 121     f2dHistGenPtVsGenYFdblOmegaMinus(0),
 122
 123     fHistThetaLambdaOmegaMinus(0),
 124     fHistThetaBachOmegaMinus(0),
 125
 126     fHistThetaMesDghterOmegaMinus(0),
 127     fHistThetaBarDghterOmegaMinus(0),
 128
 129     fHistPtBachOmegaMinus(0),
 130     fHistPtMesDghterOmegaMinus(0),
 131     fHistPtBarDghterOmegaMinus(0),
 132
 133
 134    // Omega+
 135    fHistEtaGenCascOmegaPlus(0),
 136    f2dHistGenPtVsGenYGenOmegaPlus(0),
 137
 138     fHistThetaGenCascOmegaPlus(0),
 139     f2dHistGenPtVsGenYFdblOmegaPlus(0),
 140
 141     fHistThetaLambdaOmegaPlus(0),
 142     fHistThetaBachOmegaPlus(0),
 143
 144     fHistThetaMesDghterOmegaPlus(0),
 145     fHistThetaBarDghterOmegaPlus(0),
 146
 147     fHistPtBachOmegaPlus(0),
 148     fHistPtMesDghterOmegaPlus(0),
 149     fHistPtBarDghterOmegaPlus(0),
 150
 151 // Part 2 - Association to MC
 152
 153     fHistMassXiMinus(0),
 154     fHistMassXiPlus(0),
 155     fHistMassOmegaMinus(0),
 156     fHistMassOmegaPlus(0),
 157
 158         // - Effective mass histos with combined PID
 159     fHistMassWithCombPIDXiMinus(0), fHistMassWithCombPIDXiPlus(0),
 160     fHistMassWithCombPIDOmegaMinus(0), fHistMassWithCombPIDOmegaPlus(0),
 161
 162         // - PID Probability versus MC Pt(bachelor track)
 163     f2dHistPIDprobaKaonVsMCPtBach(0), f2dHistPIDprobaPionVsMCPtBach(0),
 164
 165         // - Effective mass histos with perfect MC PID on the bachelor
 166     fHistMassWithMcPIDXiMinus(0), fHistMassWithMcPIDXiPlus(0),
 167     fHistMassWithMcPIDOmegaMinus(0), fHistMassWithMcPIDOmegaPlus(0),
 168
 169
 170         // - Effective mass histos for the cascade candidates associated with MC
 171     fHistAsMCMassXiMinus(0),
 172     fHistAsMCMassXiPlus(0),
 173     fHistAsMCMassOmegaMinus(0),
 174     fHistAsMCMassOmegaPlus(0),
 175
 176         // - Generated Pt Vs generated y, for the cascade candidates associated with MC + Info Comb. PID
 177     f2dHistAsMCandCombPIDGenPtVsGenYXiMinus(0),
 178     f2dHistAsMCandCombPIDGenPtVsGenYXiPlus(0),
 179     f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus(0),
 180     f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus(0),
 181
 182         // - Generated Pt Vs generated y, for the cascade candidates associated with MC
 183     f2dHistAsMCGenPtVsGenYXiMinus(0),
 184     f2dHistAsMCGenPtVsGenYXiPlus(0),
 185     f2dHistAsMCGenPtVsGenYOmegaMinus(0),
 186     f2dHistAsMCGenPtVsGenYOmegaPlus(0),
 187
 188         // - Generated Eta of the the cascade candidates associated with MC
 189     fHistAsMCGenEtaXiMinus(0),
 190     fHistAsMCGenEtaXiPlus(0),
 191     fHistAsMCGenEtaOmegaMinus(0),
 192     fHistAsMCGenEtaOmegaPlus(0),
 193
 194         // - Resolution in Pt as function of generated Pt
 195     f2dHistAsMCResPtXiMinus(0),
 196     f2dHistAsMCResPtXiPlus(0),
 197     f2dHistAsMCResPtOmegaMinus(0),
 198     f2dHistAsMCResPtOmegaPlus(0),
 199
 200         // - Resolution in R(2D) as function of generated R
 201     f2dHistAsMCResRXiMinus(0),
 202     f2dHistAsMCResRXiPlus(0),
 203     f2dHistAsMCResROmegaMinus(0),
 204     f2dHistAsMCResROmegaPlus(0),
 205
 206         // - Resolution in phi as function of generated Pt
 207     f2dHistAsMCResPhiXiMinus(0),
 208     f2dHistAsMCResPhiXiPlus(0),
 209     f2dHistAsMCResPhiOmegaMinus(0),
 210     f2dHistAsMCResPhiOmegaPlus(0),
 211
 212         // - Correlation in Pt between the cascade and its (anti)proton daughter
 213     f2dHistAsMCPtProtonVsPtXiMinus(0),
 214     f2dHistAsMCPtAntiProtonVsPtXiPlus(0),
 215     f2dHistAsMCPtProtonVsPtOmegaMinus(0),
 216     f2dHistAsMCPtAntiProtonVsPtOmegaPlus(0),
 217
 218     fCFContCascadePIDAsXiMinus(0),
 219     fCFContCascadePIDAsXiPlus(0),
 220     fCFContCascadePIDAsOmegaMinus(0),
 221     fCFContCascadePIDAsOmegaPlus(0),
 222
 223     fCFContAsCascadeCuts(0)
 224
 225 {
 226 // Dummy constructor
 227         for(Int_t iAlephIdx   = 0; iAlephIdx   < 5; iAlephIdx++   ) { fAlephParameters [iAlephIdx]    = -1.; }
 228         for(Int_t iV0selIdx   = 0; iV0selIdx   < 7; iV0selIdx++   ) { fV0Sels          [iV0selIdx   ] = -1.; }
 229         for(Int_t iCascSelIdx = 0; iCascSelIdx < 8; iCascSelIdx++ ) { fCascSels        [iCascSelIdx ] = -1.; }
 230 }
 231
 232
 233
 234
 235 //_____Non-default Constructor________________________________________________________________
 236 AliAnalysisTaskCheckPerformanceCascade::AliAnalysisTaskCheckPerformanceCascade(const char *name)
 237   : AliAnalysisTaskSE(name),
 238     fDebugCascade(0), fAnalysisType("ESD"), fTriggerMaskType("kMB"), fCollidingSystems(0), fESDpid(0), fESDtrackCuts(0), /*fPaveTextBookKeeping(0),*/
 239     fkRerunV0CascVertexers         (0),
 240     fkQualityCutZprimVtxPos        (kTRUE),
 241     fkRejectEventPileUp            (kTRUE),
 242     fkQualityCutNoTPConlyPrimVtx   (kTRUE),
 243     fkQualityCutTPCrefit           (kTRUE),
 244     fkQualityCut80TPCcls           (kTRUE),
 245     fkIsDataRecoWith1PadTPCCluster (kTRUE),
 246     fkExtraSelections              (0),
 247
 248         // - Cascade part initialisation
 249 fListHistCascade(0),
 250     fHistMCTrackMultiplicity(0),
 251        // - Resolution of the multiplicity estimator
 252     f2dHistRecoPrimTrckMultVsMCMult(0), f2dHistRecoEstimateMultVsMCMult(0),
 253
 254     fHistEtaGenProton(0),
 255     fHistEtaGenAntiProton(0),
 256
 257 // Xi-
 258    fHistEtaGenCascXiMinus(0),
 259    f2dHistGenPtVsGenYGenXiMinus(0),
 260
 261     fHistThetaGenCascXiMinus(0),
 262     f2dHistGenPtVsGenYFdblXiMinus(0),
 263
 264     fHistThetaLambdaXiMinus(0),
 265     fHistThetaBachXiMinus(0),
 266
 267     fHistThetaMesDghterXiMinus(0),
 268     fHistThetaBarDghterXiMinus(0),
 269
 270     fHistPtBachXiMinus(0),
 271     fHistPtMesDghterXiMinus(0),
 272     fHistPtBarDghterXiMinus(0),
 273
 274
 275    // Xi+
 276    fHistEtaGenCascXiPlus(0),
 277    f2dHistGenPtVsGenYGenXiPlus(0),
 278
 279     fHistThetaGenCascXiPlus(0),
 280     f2dHistGenPtVsGenYFdblXiPlus(0),
 281
 282     fHistThetaLambdaXiPlus(0),
 283     fHistThetaBachXiPlus(0),
 284
 285     fHistThetaMesDghterXiPlus(0),
 286     fHistThetaBarDghterXiPlus(0),
 287
 288     fHistPtBachXiPlus(0),
 289     fHistPtMesDghterXiPlus(0),
 290     fHistPtBarDghterXiPlus(0),
 291
 292    // Omega-
 293    fHistEtaGenCascOmegaMinus(0),
 294    f2dHistGenPtVsGenYGenOmegaMinus(0),
 295
 296     fHistThetaGenCascOmegaMinus(0),
 297     f2dHistGenPtVsGenYFdblOmegaMinus(0),
 298
 299     fHistThetaLambdaOmegaMinus(0),
 300     fHistThetaBachOmegaMinus(0),
 301
 302     fHistThetaMesDghterOmegaMinus(0),
 303     fHistThetaBarDghterOmegaMinus(0),
 304
 305     fHistPtBachOmegaMinus(0),
 306     fHistPtMesDghterOmegaMinus(0),
 307     fHistPtBarDghterOmegaMinus(0),
 308
 309
 310    // Omega+
 311    fHistEtaGenCascOmegaPlus(0),
 312    f2dHistGenPtVsGenYGenOmegaPlus(0),
 313
 314     fHistThetaGenCascOmegaPlus(0),
 315     f2dHistGenPtVsGenYFdblOmegaPlus(0),
 316
 317     fHistThetaLambdaOmegaPlus(0),
 318     fHistThetaBachOmegaPlus(0),
 319
 320     fHistThetaMesDghterOmegaPlus(0),
 321     fHistThetaBarDghterOmegaPlus(0),
 322
 323     fHistPtBachOmegaPlus(0),
 324     fHistPtMesDghterOmegaPlus(0),
 325     fHistPtBarDghterOmegaPlus(0),
 326
 327 // Part 2 - Association to MC
 328
 329     fHistMassXiMinus(0),
 330     fHistMassXiPlus(0),
 331     fHistMassOmegaMinus(0),
 332     fHistMassOmegaPlus(0),
 333
 334         // - Effective mass histos with combined PID
 335     fHistMassWithCombPIDXiMinus(0), fHistMassWithCombPIDXiPlus(0),
 336     fHistMassWithCombPIDOmegaMinus(0), fHistMassWithCombPIDOmegaPlus(0),
 337
 338         // - PID Probability versus MC Pt(bachelor track)
 339     f2dHistPIDprobaKaonVsMCPtBach(0), f2dHistPIDprobaPionVsMCPtBach(0),
 340
 341         // - Effective mass histos with perfect MC PID on the bachelor
 342     fHistMassWithMcPIDXiMinus(0), fHistMassWithMcPIDXiPlus(0),
 343     fHistMassWithMcPIDOmegaMinus(0), fHistMassWithMcPIDOmegaPlus(0),
 344
 345
 346         // - Effective mass histos for the cascade candidates associated with MC
 347     fHistAsMCMassXiMinus(0),
 348     fHistAsMCMassXiPlus(0),
 349     fHistAsMCMassOmegaMinus(0),
 350     fHistAsMCMassOmegaPlus(0),
 351
 352         // - Generated Pt Vs generated y, for the cascade candidates associated with MC + Info Comb. PID
 353     f2dHistAsMCandCombPIDGenPtVsGenYXiMinus(0),
 354     f2dHistAsMCandCombPIDGenPtVsGenYXiPlus(0),
 355     f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus(0),
 356     f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus(0),
 357
 358         // - Generated Pt Vs generated y, for the cascade candidates associated with MC
 359     f2dHistAsMCGenPtVsGenYXiMinus(0),
 360     f2dHistAsMCGenPtVsGenYXiPlus(0),
 361     f2dHistAsMCGenPtVsGenYOmegaMinus(0),
 362     f2dHistAsMCGenPtVsGenYOmegaPlus(0),
 363
 364         // - Generated Eta of the the cascade candidates associated with MC
 365     fHistAsMCGenEtaXiMinus(0),
 366     fHistAsMCGenEtaXiPlus(0),
 367     fHistAsMCGenEtaOmegaMinus(0),
 368     fHistAsMCGenEtaOmegaPlus(0),
 369
 370         // - Resolution in Pt as function of generated Pt
 371     f2dHistAsMCResPtXiMinus(0),
 372     f2dHistAsMCResPtXiPlus(0),
 373     f2dHistAsMCResPtOmegaMinus(0),
 374     f2dHistAsMCResPtOmegaPlus(0),
 375
 376         // - Resolution in R(2D) as function of generated R
 377     f2dHistAsMCResRXiMinus(0),
 378     f2dHistAsMCResRXiPlus(0),
 379     f2dHistAsMCResROmegaMinus(0),
 380     f2dHistAsMCResROmegaPlus(0),
 381
 382        // - Resolution in phi as function of generated Pt
 383     f2dHistAsMCResPhiXiMinus(0),
 384     f2dHistAsMCResPhiXiPlus(0),
 385     f2dHistAsMCResPhiOmegaMinus(0),
 386     f2dHistAsMCResPhiOmegaPlus(0),
 387
 388         // - Correlation in Pt between the cascade and its (anti)proton daughter
 389     f2dHistAsMCPtProtonVsPtXiMinus(0),
 390     f2dHistAsMCPtAntiProtonVsPtXiPlus(0),
 391     f2dHistAsMCPtProtonVsPtOmegaMinus(0),
 392     f2dHistAsMCPtAntiProtonVsPtOmegaPlus(0),
 393
 394     fCFContCascadePIDAsXiMinus(0),
 395     fCFContCascadePIDAsXiPlus(0),
 396     fCFContCascadePIDAsOmegaMinus(0),
 397     fCFContCascadePIDAsOmegaPlus(0),
 398
 399     fCFContAsCascadeCuts(0)
 400
 401 {
 402   // Constructor
 403
 404   // Define input and output slots here
 405   // Input slot #0 works with a TChain
 406   // Output slot #1 writes into a TList container (cascade)
 407
 408         for(Int_t iAlephIdx   = 0; iAlephIdx   < 5; iAlephIdx++   ) { fAlephParameters [iAlephIdx]    = -1.; }
 409
 410         // New Loose : 1st step for the 7 TeV pp analysis
 411
 412         fV0Sels[0] =  33.  ;  // max allowed chi2
 413         fV0Sels[1] =   0.02;  // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 414         fV0Sels[2] =   0.02;  // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 415         fV0Sels[3] =   2.0 ;  // max allowed DCA between the daughter tracks       (LHC09a4 : 0.5)
 416         fV0Sels[4] =   0.95;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 417         fV0Sels[5] =   1.0 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 418         fV0Sels[6] = 100.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 419
 420         fCascSels[0] =  33.   ;  // max allowed chi2 (same as PDC07)
 421         fCascSels[1] =   0.05 ;  // min allowed V0 impact parameter                    (PDC07 : 0.05   / LHC09a4 : 0.025 )
 422         fCascSels[2] =   0.010;  // "window" around the Lambda mass                    (PDC07 : 0.008  / LHC09a4 : 0.010 )
 423         fCascSels[3] =   0.03 ;  // min allowed bachelor's impact parameter            (PDC07 : 0.035  / LHC09a4 : 0.025 )
 424         fCascSels[4] =   2.0  ;  // max allowed DCA between the V0 and the bachelor    (PDC07 : 0.1    / LHC09a4 : 0.2   )
 425         fCascSels[5] =   0.95 ;  // min allowed cosine of the cascade pointing angle   (PDC07 : 0.9985 / LHC09a4 : 0.998 )
 426         fCascSels[6] =   0.4  ;  // min radius of the fiducial volume                  (PDC07 : 0.9    / LHC09a4 : 0.2   )
 427         fCascSels[7] = 100.   ;  // max radius of the fiducial volume                  (PDC07 : 100    / LHC09a4 : 100   )
 428
 429
 430         // Hyper Loose "à la 900 GeV 2009 data", with lower cosine of pointing angle for Xi (0.95 down to 0.82) = 900 GeV paper
 431         /*
 432         fV0Sels[0] =  33.  ;  // max allowed chi2
 433         fV0Sels[1] =   0.001; // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 434         fV0Sels[2] =   0.001; // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 435         fV0Sels[3] =   5.0 ;  // max allowed DCA between the daughter tracks       (LHC09a4 : 0.5)
 436         fV0Sels[4] =   0.0 ;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 437         fV0Sels[5] =   0.1 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 438         fV0Sels[6] = 100.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 439
 440         fCascSels[0] =  33.   ;  // max allowed chi2 (same as PDC07)
 441         fCascSels[1] =   0.001;  // min allowed V0 impact parameter                    (PDC07 : 0.05   / LHC09a4 : 0.025 )
 442         fCascSels[2] =   0.008;  // "window" around the Lambda mass                    (PDC07 : 0.008  / LHC09a4 : 0.010 )
 443         fCascSels[3] =   0.001;  // min allowed bachelor's impact parameter            (PDC07 : 0.035  / LHC09a4 : 0.025 )
 444         fCascSels[4] =   5.0  ;  // max allowed DCA between the V0 and the bachelor    (PDC07 : 0.1    / LHC09a4 : 0.2   )
 445         fCascSels[5] =   0.82 ;  //FIXME min allowed cosine of the cascade pointing angle   (PDC07 : 0.9985 / LHC09a4 : 0.998 )
 446         fCascSels[6] =   0.1  ;  // min radius of the fiducial volume                  (PDC07 : 0.9    / LHC09a4 : 0.2   )
 447         fCascSels[7] = 100.   ;  // max radius of the fiducial volume                  (PDC07 : 100    / LHC09a4 : 100   )
 448         */
 449
 450         //New default vtxR (http://alisoft.cern.ch/viewvc?view=rev&root=AliRoot&revision=40955, 5 May 2010)
 451         /*
 452         fV0Sels[0] =  33.  ;  // max allowed chi2
 453         fV0Sels[1] =   0.05;  // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 454         fV0Sels[2] =   0.05;  // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 455         fV0Sels[3] =   1.5 ;  // max allowed DCA between the daughter tracks       (LHC09a4 : 0.5)
 456         fV0Sels[4] =   0.9 ;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 457         fV0Sels[5] =   0.2 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 458         fV0Sels[6] = 100.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 459
 460         fCascSels[0] =  33.   ;  // max allowed chi2 (same as PDC07)
 461         fCascSels[1] =   0.01 ;  // min allowed V0 impact parameter                    (PDC07 : 0.05   / LHC09a4 : 0.025 )
 462         fCascSels[2] =   0.008;  // "window" around the Lambda mass                    (PDC07 : 0.008  / LHC09a4 : 0.010 )
 463         fCascSels[3] =   0.01 ;  // min allowed bachelor's impact parameter            (PDC07 : 0.035  / LHC09a4 : 0.025 )
 464         fCascSels[4] =   2.0  ;  // max allowed DCA between the V0 and the bachelor    (PDC07 : 0.1    / LHC09a4 : 0.2   )
 465         fCascSels[5] =   0.98 ;  // min allowed cosine of the cascade pointing angle   (PDC07 : 0.9985 / LHC09a4 : 0.998 )
 466         fCascSels[6] =   0.2  ;  // min radius of the fiducial volume                  (PDC07 : 0.9    / LHC09a4 : 0.2   )
 467         fCascSels[7] = 100.   ;  // max radius of the fiducial volume                  (PDC07 : 100    / LHC09a4 : 100   )
 468         */
 469
 470   DefineOutput(1, TList::Class());
 471
 472 }
 473
 474
 475 AliAnalysisTaskCheckPerformanceCascade::~AliAnalysisTaskCheckPerformanceCascade()
 476 {
 477   //
 478   // Destructor
 479   //
 480
 481   // For all TH1, 2, 3 HnSparse and CFContainer are in the fListCascade TList.
 482   // They will be deleted when fListCascade is deleted by the TSelector dtor
 483   // Because of TList::SetOwner()
 484
 485   if (fListHistCascade)      { delete fListHistCascade;     fListHistCascade = 0x0;  }
 486   if (fESDpid)               { delete fESDpid;              fESDpid = 0x0;} // fESDpid is not stored into the TList
 487   if (fESDtrackCuts)         { delete fESDtrackCuts;        fESDtrackCuts = 0x0; }
 488   /*if (fPaveTextBookKeeping)  { delete fPaveTextBookKeeping; fPaveTextBookKeeping = 0x0; } // fPaveTextBookKeeping is not stored into the TList*/
 489 }
 490
 491
 492 //________________________________________________________________________
 493 void AliAnalysisTaskCheckPerformanceCascade::UserCreateOutputObjects()
 494 {
 495   // Create histograms
 496   // Called once
 497
 498
 499    // Option for AliLog
 500         AliLog::SetGlobalLogLevel(AliLog::kError);
 501         // to suppress the extensive info prompted by a run with MC
 502
 503    // Definition of the output datamembers
 504    fListHistCascade = new TList();
 505    fListHistCascade->SetOwner(); // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
 506
 507 if(! fESDpid){
 508
 509          if(fkIsDataRecoWith1PadTPCCluster){
 510                 // Home made parameterization for LHC10f6a production = p+p 7 TeV
 511                 fAlephParameters[0] = 0.04;
 512                 fAlephParameters[1] = 17.5;
 513                 fAlephParameters[2] = 3.4e-09;
 514                 fAlephParameters[3] = 2.15;
 515                 fAlephParameters[4] = 3.91720e+00;
 516
 517                 // Home made parameterization for LHC10e13 production = p+p 900 GeV/c
 518          }
 519          else {
 520                 // Reasonable parameters extracted for p-p simulation (LHC09a4) - A.Kalweit
 521                  // fAlephParameters[0] = 4.23232575531564326e+00/50;//50*0.76176e-1;
 522                  // fAlephParameters[1] = 8.68482806165147636e+00;//10.632;
 523                  // fAlephParameters[2] = 1.34000000000000005e-05;//0.13279e-4;
 524                  // fAlephParameters[3] = 2.30445734159456084e+00;//1.8631;
 525                  // fAlephParameters[4] = 2.25624744086878559e+00;//1.9479;
 526
 527                 // Param for LHC09d10 prod - A.Kalweit
 528                 fAlephParameters[0] = 2.15898e+00/50.;
 529                 fAlephParameters[1] = 1.75295e+01;
 530                 fAlephParameters[2] = 3.40030e-09;
 531                 fAlephParameters[3] = 1.96178e+00;
 532                 fAlephParameters[4] = 3.91720e+00;
 533          }
 534          Printf("CheckPerfCascade - Check Aleph Param in case of MC Data   (fAlephParameters[3] = %f) (To be compared with : 2.15 for 1-pad-cluster prod. / 1.96178 otherwise)\n",  fAlephParameters[3]);
 535
 536   fESDpid = new AliESDpid();
 537   fESDpid->GetTPCResponse().SetBetheBlochParameters( fAlephParameters[0],
 538                                                      fAlephParameters[1],
 539                                                      fAlephParameters[2],
 540                                                      fAlephParameters[3],
 541                                                      fAlephParameters[4] );
 542 }
 543
 544
 545 if(! fESDtrackCuts ){
 546       fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE); // Std definition of primary (see kTRUE argument) tracks for 2010
 547       fESDtrackCuts->SetEtaRange(-0.8,+0.8);
 548       fESDtrackCuts->SetPtRange(0.15, 1e10);
 549       Printf("CheckCascade - ESDtrackCuts set up to 2010 std ITS-TPC cuts...");
 550 }
 551
 552
 553 /*
 554 if( !fPaveTextBookKeeping){
 555         fPaveTextBookKeeping = new TPaveText(0.1, 0.1, 0.9, 0.9,"NDC");
 556         fPaveTextBookKeeping->SetName("fPaveTextBookKeeping");
 557         fPaveTextBookKeeping->SetBorderSize(0);
 558         fPaveTextBookKeeping->SetTextAlign(12);
 559         fPaveTextBookKeeping->SetFillColor(kWhite);
 560         fPaveTextBookKeeping->SetTextFont(42);        // regular Arial or Helvetica,
 561         fPaveTextBookKeeping->SetTextColor(kGray+3);
 562
 563
 564         fPaveTextBookKeeping->AddText( "Task CHECK PERFORMANCE CASCADE analysis" );
 565         fPaveTextBookKeeping->AddText("- - - - - - - - - - - ");
 566         fPaveTextBookKeeping->AddText( Form("AnalysisType : %s ", fAnalysisType.Data() ));
 567         if(!fCollidingSystems)  fPaveTextBookKeeping->AddText("Colliding system : p-p collisions ");
 568         else                    fPaveTextBookKeeping->AddText("Colliding system : A-A collisions ");
 569
 570         fPaveTextBookKeeping->AddText("- - - - - - - - - - - ");
 571
 572         if(fkRerunV0CascVertexers){
 573                 fPaveTextBookKeeping->AddText("A.1. With V0 vertexer : ");
 574                 fPaveTextBookKeeping->AddText( Form("  - V0 #chi^{2} _________________ <  %.3f ",               fV0Sels[0] ));
 575                 fPaveTextBookKeeping->AddText( Form("  - DCA(prim. Vtx/ 1^{st} daughter) ___ >  %.3f     cm ",  fV0Sels[1] ));
 576                 fPaveTextBookKeeping->AddText( Form("  - DCA(prim. Vtx/ 2^{nd} daughter) __  >  %.3f     cm",   fV0Sels[2] ));
 577                 fPaveTextBookKeeping->AddText( Form("  - DCA between V0 daughters ___ <  %.3f      cm",         fV0Sels[3] ));
 578                 fPaveTextBookKeeping->AddText( Form("  - cos(V0 pointing angle) _______ >  %.3f ",              fV0Sels[4] ));
 579                 fPaveTextBookKeeping->AddText( Form("  - R_{transv}(V0 decay) ________ >  %.3f             cm", fV0Sels[5] ));
 580                 fPaveTextBookKeeping->AddText( Form("  - R_{transv}(V0 decay) ________ <  %.3f         cm",     fV0Sels[6] ));
 581
 582                 fPaveTextBookKeeping->AddText(" ");
 583
 584                 fPaveTextBookKeeping->AddText("A.2. With Casc. vertexer : ");
 585                 fPaveTextBookKeeping->AddText( Form("  - Casc. #chi^{2} ______________  <  %.3f ",                               fCascSels[0] ));
 586                 fPaveTextBookKeeping->AddText( Form("  - DCA(prim. Vtx/ V0) _________ >  %.3f    cm",                            fCascSels[1] ));
 587                 fPaveTextBookKeeping->AddText( Form("  - | M_{#Lambda}(reco) - M_{#Lambda}(pdg) | _______ <  %.3f    GeV/c^{2}", fCascSels[2] ));
 588                 fPaveTextBookKeeping->AddText( Form("  - DCA(prim. Vtx/ Bach) _______ >  %.3f    cm",                            fCascSels[3] ));
 589                 fPaveTextBookKeeping->AddText( Form("  - DCA between Bach/ #Lambda ______ <  %.3f    cm",                        fCascSels[4] ));
 590                 fPaveTextBookKeeping->AddText( Form("  - cos(Casc. pointing angle) ____ >  %.3f ",                               fCascSels[5] ));
 591                 fPaveTextBookKeeping->AddText( Form("  - R_{transv}(Casc. decay) ______ >  %.3f         cm",                     fCascSels[6] ));
 592                 fPaveTextBookKeeping->AddText( Form("  - R_{transv}(Casc. decay) ______ <  %.3f     cm",                         fCascSels[7] ));
 593         }
 594         else{   fPaveTextBookKeeping->AddText("A. No rerunning of the V0/Casc. vertexers ... See std cuts in (AliRoot+Rec.C) used for this prod. cycle");}
 595
 596         fPaveTextBookKeeping->AddText("- - - - - - - - - - - ");
 597
 598         if(fkQualityCutZprimVtxPos)      fPaveTextBookKeeping->AddText("B. Quality Cut(prim. Vtx z-Pos)    = ON  ");
 599         else                             fPaveTextBookKeeping->AddText("B. Quality Cut(prim. Vtx z-Pos)    = Off ");
 600         if(fkQualityCutNoTPConlyPrimVtx) fPaveTextBookKeeping->AddText("C. Quality Cut(No TPC prim. vtx) = ON  ");
 601         else                             fPaveTextBookKeeping->AddText("C. Quality Cut(No TPC prim. vtx) = Off ");
 602         if(fkQualityCutTPCrefit)         fPaveTextBookKeeping->AddText("D. Quality Cut(TPCrefit)               = ON  ");
 603         else                             fPaveTextBookKeeping->AddText("D. Quality Cut(TPCrefit)               = Off ");
 604         if(fkQualityCut80TPCcls)         fPaveTextBookKeeping->AddText("E. Quality Cut(80 TPC clusters)   = ON  ");
 605         else                             fPaveTextBookKeeping->AddText("E. Quality Cut(80 TPC clusters)   = Off ");
 606         if(fkExtraSelections)            fPaveTextBookKeeping->AddText("F. Extra Analysis Selections         = ON  ");
 607         else                             fPaveTextBookKeeping->AddText("F. Extra Analysis Selections         = Off ");
 608
 609         fPaveTextBookKeeping->AddText("- - - - - - - - - - - ");
 610
 611         fPaveTextBookKeeping->AddText("G. TPC Aleph Param : ");
 612         fPaveTextBookKeeping->AddText( Form("   - fAlephParam [0] =  %.5g", fAlephParameters[0] ));
 613         fPaveTextBookKeeping->AddText( Form("   - fAlephParam [1] =  %.5g", fAlephParameters[1] ));
 614         fPaveTextBookKeeping->AddText( Form("   - fAlephParam [2] =  %.5g", fAlephParameters[2] ));
 615         fPaveTextBookKeeping->AddText( Form("   - fAlephParam [3] =  %.5g", fAlephParameters[3] ));
 616         fPaveTextBookKeeping->AddText( Form("   - fAlephParam [4] =  %.5g", fAlephParameters[4] ));
 617
 618         fListHistCascade->Add(fPaveTextBookKeeping);
 619 }
 620 */
 621
 622   // - General
 623
 624   if (!fHistMCTrackMultiplicity) {
 625      fHistMCTrackMultiplicity = new TH1F("fHistMCTrackMultiplicity", "MC Track Multiplicity;Number of MC tracks;Events", 100, 0, 500);
 626    //  fHistMCTrackMultiplicity = new TH1F("fHistMCTrackMultiplicity", "Multiplicity distribution;Number of tracks;Events", 200, 0, 40000); //HERE
 627   fListHistCascade->Add(fHistMCTrackMultiplicity);
 628   }
 629
 630     // - Resolution of the multiplicity estimator
 631   if(! f2dHistRecoPrimTrckMultVsMCMult){
 632        f2dHistRecoPrimTrckMultVsMCMult = new TH2F("f2dHistRecoPrimTrckMultVsMCMult", "Resolution of the multiplicity estimator (prim. tracks in |#eta| < 0.8); Reco Multiplicity (prim. tracks); MC multiplicity (gen. part. in |#eta| < 0.8)", 120, 0., 120., 300, 0., 300.);
 633        fListHistCascade->Add(f2dHistRecoPrimTrckMultVsMCMult);
 634   }
 635
 636   if(! f2dHistRecoEstimateMultVsMCMult){
 637        f2dHistRecoEstimateMultVsMCMult = new TH2F("f2dHistRecoEstimateMultVsMCMult", "Resolution of the multiplicity estimator (EstimateMult. in |#eta| < 1.0); Reco Multiplicity (tr(ITS-TPC)+ITSsa+tracklets); MC multiplicity (gen. part. in |#eta| < 1.0)", 160, 0., 160., 300, 0., 300.);
 638        fListHistCascade->Add(f2dHistRecoEstimateMultVsMCMult);
 639   }
 640
 641   if (!fHistEtaGenProton) {
 642      fHistEtaGenProton = new TH1F("fHistEtaGenProton", "#eta of any gen. p^{+};#eta;Number of prim. protons", 200, -10, 10);
 643      fListHistCascade->Add(fHistEtaGenProton);
 644   }
 645
 646   if (!fHistEtaGenAntiProton) {
 647      fHistEtaGenAntiProton = new TH1F("fHistEtaGenAntiProton", "#eta of any gen. #bar{p}^{-};#eta;Number of prim. #bar{p}", 200, -10, 10);
 648      fListHistCascade->Add(fHistEtaGenAntiProton);
 649   }
 650
 651
 652
 653
 654
 655
 656
 657   //--------
 658   // I - Xi-
 659   // - Pseudo-Rapidity distribution
 660   if (!fHistEtaGenCascXiMinus) {
 661      fHistEtaGenCascXiMinus = new TH1F("fHistEtaGenCascXiMinus", "#eta of any gen. #Xi^{-};#eta;Number of Casc", 200, -10, 10);
 662      fListHistCascade->Add(fHistEtaGenCascXiMinus);
 663   }
 664
 665   if (!f2dHistGenPtVsGenYGenXiMinus) {
 666      f2dHistGenPtVsGenYGenXiMinus = new TH2F("f2dHistGenPtVsGenYGenXiMinus", "MC P_{t} Vs MC Y of Gen #Xi^{-} ;Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 667      fListHistCascade->Add(f2dHistGenPtVsGenYGenXiMinus);
 668   }
 669
 670
 671   // - Info at the generation level of multi-strange particle
 672
 673   if (!fHistThetaGenCascXiMinus) {
 674      fHistThetaGenCascXiMinus = new TH1F("fHistThetaGenCascXiMinus", "#theta of gen. #Xi^{-};#theta;Number of Casc.", 200, -10, 190);
 675      fListHistCascade->Add(fHistThetaGenCascXiMinus);
 676   }
 677
 678   if (!f2dHistGenPtVsGenYFdblXiMinus) {
 679      f2dHistGenPtVsGenYFdblXiMinus = new TH2F("f2dHistGenPtVsGenYFdblXiMinus", "MC P_{t} Vs MC Y of findable Gen #Xi^{-}; Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 680      fListHistCascade->Add(f2dHistGenPtVsGenYFdblXiMinus);
 681   }
 682
 683                 // - Theta distribution the daughters (control plots)
 684
 685   if (!fHistThetaLambdaXiMinus) {
 686      fHistThetaLambdaXiMinus = new TH1F("fHistThetaLambdaXiMinus", "#theta of gen. #Lambda (Xi dghter);#theta_{#Lambda};Number of #Lambda^0", 200, -10, 190);
 687      fListHistCascade->Add(fHistThetaLambdaXiMinus);
 688   }
 689
 690   if (!fHistThetaBachXiMinus) {
 691      fHistThetaBachXiMinus = new TH1F("fHistThetaBachXiMinus", "#theta of gen. Bach.;#theta_{Bach};Number of Bach.", 200, -10, 190);
 692      fListHistCascade->Add(fHistThetaBachXiMinus);
 693   }
 694
 695   if (!fHistThetaMesDghterXiMinus) {
 696      fHistThetaMesDghterXiMinus = new TH1F("fHistThetaMesDghterXiMinus", "#theta of gen. Meson #Lambda dghter;#theta_{MesDght};Number of Mes.", 200, -10, 190);
 697      fListHistCascade->Add(fHistThetaMesDghterXiMinus);
 698   }
 699
 700   if (!fHistThetaBarDghterXiMinus) {
 701      fHistThetaBarDghterXiMinus = new TH1F("fHistThetaBarDghterXiMinus", "#theta of gen. Baryon #Lambda dghter;#theta_{BarDght};Number of Bar.", 200, -10, 190);
 702      fListHistCascade->Add(fHistThetaBarDghterXiMinus);
 703   }
 704
 705                 // - Pt distribution (control plots)
 706
 707   if (!fHistPtBachXiMinus) {
 708      fHistPtBachXiMinus = new TH1F("fHistPtBachXiMinus", "p_{t} of gen. Bach.;pt_{Bach};Number of Bach.", 200, 0, 10);
 709      fListHistCascade->Add(fHistPtBachXiMinus);
 710   }
 711
 712   if (!fHistPtMesDghterXiMinus) {
 713      fHistPtMesDghterXiMinus = new TH1F("fHistPtMesDghterXiMinus", "p_{t} of gen. Meson #Lambda dghter;pt_{MesDght};Number of Mes.", 200, 0, 10);
 714      fListHistCascade->Add(fHistPtMesDghterXiMinus);
 715   }
 716
 717   if (!fHistPtBarDghterXiMinus) {
 718      fHistPtBarDghterXiMinus = new TH1F("fHistPtBarDghterXiMinus", "p_{t} of gen. Baryon #Lambda dghter;pt_{BarDght};Number of Bar.", 200, 0, 10);
 719      fListHistCascade->Add(fHistPtBarDghterXiMinus);
 720   }
 721
 722
 723
 724   //--------
 725   // II - Xi+
 726   // - Pseudo-Rapidity distribution
 727   if (!fHistEtaGenCascXiPlus) {
 728      fHistEtaGenCascXiPlus = new TH1F("fHistEtaGenCascXiPlus", "#eta of any gen. #bar{#Xi}^{+};#eta;Number of Casc", 200, -10, 10);
 729      fListHistCascade->Add(fHistEtaGenCascXiPlus);
 730   }
 731
 732   if (!f2dHistGenPtVsGenYGenXiPlus) {
 733      f2dHistGenPtVsGenYGenXiPlus = new TH2F("f2dHistGenPtVsGenYGenXiPlus", "MC P_{t} Vs MC Y of Gen #bar{#Xi}^{+} ;Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 734      fListHistCascade->Add(f2dHistGenPtVsGenYGenXiPlus);
 735   }
 736
 737
 738   // - Info at the generation level of multi-strange particle
 739
 740   if (!fHistThetaGenCascXiPlus) {
 741      fHistThetaGenCascXiPlus = new TH1F("fHistThetaGenCascXiPlus", "#theta of gen. #bar{#Xi}^{+};#theta;Number of Casc.", 200, -10, 190);
 742      fListHistCascade->Add(fHistThetaGenCascXiPlus);
 743   }
 744
 745   if (!f2dHistGenPtVsGenYFdblXiPlus) {
 746      f2dHistGenPtVsGenYFdblXiPlus = new TH2F("f2dHistGenPtVsGenYFdblXiPlus", "MC P_{t} Vs MC Y of findable Gen #bar{#Xi}^{+} ;Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 747      fListHistCascade->Add(f2dHistGenPtVsGenYFdblXiPlus);
 748   }
 749
 750                 // - Theta distribution the daughters (control plots)
 751
 752   if (!fHistThetaLambdaXiPlus) {
 753      fHistThetaLambdaXiPlus = new TH1F("fHistThetaLambdaXiPlus", "#theta of gen. #Lambda (Xi dghter);#theta_{#Lambda};Number of #Lambda", 200, -10, 190);
 754      fListHistCascade->Add(fHistThetaLambdaXiPlus);
 755   }
 756
 757   if (!fHistThetaBachXiPlus) {
 758      fHistThetaBachXiPlus = new TH1F("fHistThetaBachXiPlus", "#theta of gen. Bach.;#theta_{Bach};Number of Bach.", 200, -10, 190);
 759      fListHistCascade->Add(fHistThetaBachXiPlus);
 760   }
 761
 762   if (!fHistThetaMesDghterXiPlus) {
 763      fHistThetaMesDghterXiPlus = new TH1F("fHistThetaMesDghterXiPlus", "#theta of gen. Meson #Lambda dghter;#theta_{MesDght};Number of Mes.", 200, -10, 190);
 764      fListHistCascade->Add(fHistThetaMesDghterXiPlus);
 765   }
 766
 767   if (!fHistThetaBarDghterXiPlus) {
 768      fHistThetaBarDghterXiPlus = new TH1F("fHistThetaBarDghterXiPlus", "#theta of gen. Baryon #Lambda dghter;#theta_{BarDght};Number of Bar.", 200, -10, 190);
 769      fListHistCascade->Add(fHistThetaBarDghterXiPlus);
 770   }
 771
 772                 // - Pt distribution (control plots)
 773
 774   if (!fHistPtBachXiPlus) {
 775      fHistPtBachXiPlus = new TH1F("fHistPtBachXiPlus", "p_{t} of gen. Bach.;pt_{Bach};Number of Bach.", 200, 0, 10);
 776      fListHistCascade->Add(fHistPtBachXiPlus);
 777   }
 778
 779   if (!fHistPtMesDghterXiPlus) {
 780      fHistPtMesDghterXiPlus = new TH1F("fHistPtMesDghterXiPlus", "p_{t} of gen. Meson #Lambda dghter);pt_{MesDght};Number of Mes.", 200, 0, 10);
 781      fListHistCascade->Add(fHistPtMesDghterXiPlus);
 782   }
 783
 784   if (!fHistPtBarDghterXiPlus) {
 785      fHistPtBarDghterXiPlus = new TH1F("fHistPtBarDghterXiPlus", "p_{t} of gen. Baryon #Lambda dghter);pt_{BarDght};Number of Bar.", 200, 0, 10);
 786      fListHistCascade->Add(fHistPtBarDghterXiPlus);
 787   }
 788
 789
 790   //---------
 791   // III - Omega-
 792   // - Pseudo-Rapidity distribution
 793   if (!fHistEtaGenCascOmegaMinus) {
 794      fHistEtaGenCascOmegaMinus = new TH1F("fHistEtaGenCascOmegaMinus", "#eta of any gen. #Omega^{-};#eta;Number of Casc", 200, -10, 10);
 795      fListHistCascade->Add(fHistEtaGenCascOmegaMinus);
 796   }
 797
 798   if (!f2dHistGenPtVsGenYGenOmegaMinus) {
 799      f2dHistGenPtVsGenYGenOmegaMinus = new TH2F("f2dHistGenPtVsGenYGenOmegaMinus", "MC P_{t} Vs MC Y of Gen #Omega^{-} ;Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 800      fListHistCascade->Add(f2dHistGenPtVsGenYGenOmegaMinus);
 801   }
 802
 803
 804   // - Info at the generation level of multi-strange particle
 805
 806   if (!fHistThetaGenCascOmegaMinus) {
 807      fHistThetaGenCascOmegaMinus = new TH1F("fHistThetaGenCascOmegaMinus", "#theta of gen. #Omega^{-};#theta;Number of Casc.", 200, -10, 190);
 808      fListHistCascade->Add(fHistThetaGenCascOmegaMinus);
 809   }
 810
 811   if (!f2dHistGenPtVsGenYFdblOmegaMinus) {
 812      f2dHistGenPtVsGenYFdblOmegaMinus = new TH2F("f2dHistGenPtVsGenYFdblOmegaMinus", "MC P_{t} Vs MC Y of findable Gen #Omega^{-}; Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 813      fListHistCascade->Add(f2dHistGenPtVsGenYFdblOmegaMinus);
 814   }
 815
 816                 // - Theta distribution the daughters (control plots)
 817
 818   if (!fHistThetaLambdaOmegaMinus) {
 819      fHistThetaLambdaOmegaMinus = new TH1F("fHistThetaLambdaOmegaMinus", "#theta of gen. #Lambda (Omega dghter);#theta_{#Lambda};Number of #Lambda", 200, -10, 190);
 820      fListHistCascade->Add(fHistThetaLambdaOmegaMinus);
 821   }
 822
 823   if (!fHistThetaBachOmegaMinus) {
 824      fHistThetaBachOmegaMinus = new TH1F("fHistThetaBachOmegaMinus", "#theta of gen. Bach.;#theta_{Bach};Number of Bach.", 200, -10, 190);
 825      fListHistCascade->Add(fHistThetaBachOmegaMinus);
 826   }
 827
 828   if (!fHistThetaMesDghterOmegaMinus) {
 829      fHistThetaMesDghterOmegaMinus = new TH1F("fHistThetaMesDghterOmegaMinus", "#theta of gen. Meson #Lambda dghter;#theta_{MesDght};Number of Mes.", 200, -10, 190);
 830      fListHistCascade->Add(fHistThetaMesDghterOmegaMinus);
 831   }
 832
 833   if (!fHistThetaBarDghterOmegaMinus) {
 834      fHistThetaBarDghterOmegaMinus = new TH1F("fHistThetaBarDghterOmegaMinus", "#theta of gen. Baryon #Lambda dghter;#theta_{BarDght};Number of Bar.", 200, -10, 190);
 835      fListHistCascade->Add(fHistThetaBarDghterOmegaMinus);
 836   }
 837
 838                 // - Pt distribution (control plots)
 839
 840   if (!fHistPtBachOmegaMinus) {
 841      fHistPtBachOmegaMinus = new TH1F("fHistPtBachOmegaMinus", "p_{t} of gen. Bach.;pt_{Bach};Number of Bach.", 200, 0, 10);
 842      fListHistCascade->Add(fHistPtBachOmegaMinus);
 843   }
 844
 845   if (!fHistPtMesDghterOmegaMinus) {
 846      fHistPtMesDghterOmegaMinus = new TH1F("fHistPtMesDghterOmegaMinus", "p_{t} of gen. Meson #Lambda dghter);pt_{MesDght};Number of Mes.", 200, 0, 10);
 847      fListHistCascade->Add(fHistPtMesDghterOmegaMinus);
 848   }
 849
 850   if (!fHistPtBarDghterOmegaMinus) {
 851      fHistPtBarDghterOmegaMinus = new TH1F("fHistPtBarDghterOmegaMinus", "p_{t} of gen. Baryon #Lambda dghter);pt_{BarDght};Number of Bar.", 200, 0, 10);
 852      fListHistCascade->Add(fHistPtBarDghterOmegaMinus);
 853   }
 854
 855
 856   //---------
 857   // IV - Omega+
 858   // - Pseudo-Rapidity distribution
 859   if (!fHistEtaGenCascOmegaPlus) {
 860      fHistEtaGenCascOmegaPlus = new TH1F("fHistEtaGenCascOmegaPlus", "#eta of any gen. #bar{#Omega}^{+};#eta;Number of Casc", 200, -10, 10);
 861      fListHistCascade->Add(fHistEtaGenCascOmegaPlus);
 862   }
 863
 864   if (!f2dHistGenPtVsGenYGenOmegaPlus) {
 865      f2dHistGenPtVsGenYGenOmegaPlus = new TH2F("f2dHistGenPtVsGenYGenOmegaPlus", "MC P_{t} Vs MC Y of Gen #bar{#Omega}^{+} ;Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 866      fListHistCascade->Add(f2dHistGenPtVsGenYGenOmegaPlus);
 867   }
 868
 869
 870
 871   // - Info at the generation level of multi-strange particle
 872
 873   if (!fHistThetaGenCascOmegaPlus) {
 874      fHistThetaGenCascOmegaPlus = new TH1F("fHistThetaGenCascOmegaPlus", "#theta of gen. #bar{#Omega}^{+};#theta;Number of Casc.", 200, -10, 190);
 875      fListHistCascade->Add(fHistThetaGenCascOmegaPlus);
 876   }
 877
 878   if (!f2dHistGenPtVsGenYFdblOmegaPlus) {
 879      f2dHistGenPtVsGenYFdblOmegaPlus = new TH2F("f2dHistGenPtVsGenYFdblOmegaPlus", "MC P_{t} Vs MC Y of findable Gen #bar{#Omega}^{+}; Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
 880      fListHistCascade->Add(f2dHistGenPtVsGenYFdblOmegaPlus);
 881   }
 882
 883
 884                 // - Theta distribution the daughters (control plots)
 885
 886   if (!fHistThetaLambdaOmegaPlus) {
 887      fHistThetaLambdaOmegaPlus = new TH1F("fHistThetaLambdaOmegaPlus", "#theta of gen. #Lambda (Omega dghter);#theta_{#Lambda};Number of #Lambda", 200, -10, 190);
 888      fListHistCascade->Add(fHistThetaLambdaOmegaPlus);
 889   }
 890
 891   if (!fHistThetaBachOmegaPlus) {
 892      fHistThetaBachOmegaPlus = new TH1F("fHistThetaBachOmegaPlus", "#theta of gen. Bach.;#theta_{Bach};Number of Bach.", 200, -10, 190);
 893      fListHistCascade->Add(fHistThetaBachOmegaPlus);
 894   }
 895
 896   if (!fHistThetaMesDghterOmegaPlus) {
 897      fHistThetaMesDghterOmegaPlus = new TH1F("fHistThetaMesDghterOmegaPlus", "#theta of gen. Meson #Lambda dghter;#theta_{MesDght};Number of Mes.", 200, -10, 190);
 898      fListHistCascade->Add(fHistThetaMesDghterOmegaPlus);
 899   }
 900
 901   if (!fHistThetaBarDghterOmegaPlus) {
 902      fHistThetaBarDghterOmegaPlus = new TH1F("fHistThetaBarDghterOmegaPlus", "#theta of gen. Baryon #Lambda dghter;#theta_{BarDght};Number of Bar.", 200, -10, 190);
 903      fListHistCascade->Add(fHistThetaBarDghterOmegaPlus);
 904   }
 905
 906                 // - Pt distribution (control plots)
 907
 908   if (!fHistPtBachOmegaPlus) {
 909      fHistPtBachOmegaPlus = new TH1F("fHistPtBachOmegaPlus", "p_{t} of gen. Bach.;pt_{Bach};Number of Bach.", 200, 0, 10);
 910      fListHistCascade->Add(fHistPtBachOmegaPlus);
 911   }
 912
 913   if (!fHistPtMesDghterOmegaPlus) {
 914      fHistPtMesDghterOmegaPlus = new TH1F("fHistPtMesDghterOmegaPlus", "p_{t} of gen. Meson #Lambda dghter);pt_{MesDght};Number of Mes.", 200, 0, 10);
 915      fListHistCascade->Add(fHistPtMesDghterOmegaPlus);
 916   }
 917
 918   if (!fHistPtBarDghterOmegaPlus) {
 919      fHistPtBarDghterOmegaPlus = new TH1F("fHistPtBarDghterOmegaPlus", "p_{t} of gen. Baryon #Lambda dghter);pt_{BarDght};Number of Bar.", 200, 0, 10);
 920      fListHistCascade->Add(fHistPtBarDghterOmegaPlus);
 921   }
 922
 923
 924 //--------------------------------------------------------------------------------
 925 // Part 2 - Any reconstructed cascades + reconstructed cascades associated with MC
 926
 927                 // - Effective mass histos for cascades candidates.
 928
 929   if (! fHistMassXiMinus) {
 930           fHistMassXiMinus = new TH1F("fHistMassXiMinus","#Xi^{-} candidates;M( #Lambda , #pi^{-} ) (GeV/c^{2});Counts", 400,1.2,2.0);
 931           fListHistCascade->Add(fHistMassXiMinus);
 932   }
 933
 934   if (! fHistMassXiPlus) {
 935           fHistMassXiPlus = new TH1F("fHistMassXiPlus","#bar{#Xi}^{+} candidates;M( #bar{#Lambda}^{0} , #pi^{+} ) (GeV/c^{2});Counts",400,1.2,2.0);
 936           fListHistCascade->Add(fHistMassXiPlus);
 937   }
 938
 939   if (! fHistMassOmegaMinus) {
 940           fHistMassOmegaMinus = new TH1F("fHistMassOmegaMinus","#Omega^{-} candidates;M( #Lambda , K^{-} ) (GeV/c^{2});Counts", 500,1.5,2.5);
 941           fListHistCascade->Add(fHistMassOmegaMinus);
 942   }
 943
 944   if (! fHistMassOmegaPlus) {
 945           fHistMassOmegaPlus = new TH1F("fHistMassOmegaPlus","#bar{#Omega}^{+} candidates;M( #bar{#Lambda}^{0} , K^{+} ) (GeV/c^{2});Counts", 500,1.5,2.5);
 946           fListHistCascade->Add(fHistMassOmegaPlus);
 947   }
 948
 949
 950
 951                 // - Effective mass histos with combined PID
 952
 953   if (! fHistMassWithCombPIDXiMinus) {
 954     fHistMassWithCombPIDXiMinus = new TH1F("fHistMassWithCombPIDXiMinus","#Xi^{-} candidates, with Bach. comb. PID;M( #Lambda , #pi^{-} ) (GeV/c^{2});Counts", 400,1.2,2.0);
 955     fListHistCascade->Add(fHistMassWithCombPIDXiMinus);
 956   }
 957
 958   if (! fHistMassWithCombPIDXiPlus) {
 959     fHistMassWithCombPIDXiPlus = new TH1F("fHistMassWithCombPIDXiPlus","#bar{#Xi}^{+} candidates, with Bach. comb. PID;M( #bar{#Lambda}^{0} , #pi^{+} ) (GeV/c^{2});Counts",400,1.2,2.0);
 960     fListHistCascade->Add(fHistMassWithCombPIDXiPlus);
 961   }
 962
 963   if (! fHistMassWithCombPIDOmegaMinus) {
 964         fHistMassWithCombPIDOmegaMinus = new TH1F("fHistMassWithCombPIDOmegaMinus","#Omega^{-} candidates, with Bach. comb. PID;M( #Lambda , K^{-} ) (GeV/c^{2});Counts", 500,1.5,2.5);
 965     fListHistCascade->Add(fHistMassWithCombPIDOmegaMinus);
 966   }
 967
 968   if (! fHistMassWithCombPIDOmegaPlus) {
 969         fHistMassWithCombPIDOmegaPlus = new TH1F("fHistMassWithCombPIDOmegaPlus","#bar{#Omega}^{+} candidates, with Bach. comb. PID;M( #bar{#Lambda}^{0} , K^{+} ) (GeV/c^{2});Counts", 500,1.5,2.5);
 970     fListHistCascade->Add(fHistMassWithCombPIDOmegaPlus);
 971   }
 972
 973                 // - PID Probability versus MC Pt(bachelor track)
 974   if(! f2dHistPIDprobaKaonVsMCPtBach ){
 975         f2dHistPIDprobaKaonVsMCPtBach  = new TH2F( "f2dHistPIDprobaKaonVsMCPtBach" , "Comb. PID proba to be K^{#pm} Vs MC Bach. Pt ; Pt_{MC}(Bach.) (GeV/c); Comb. PID Proba (Bach. = K^{#pm})", 100, 0.0, 5.0, 110, 0.0, 1.10 );
 976         fListHistCascade->Add(f2dHistPIDprobaKaonVsMCPtBach);
 977   }
 978
 979   if(! f2dHistPIDprobaPionVsMCPtBach ){
 980         f2dHistPIDprobaPionVsMCPtBach  = new TH2F( "f2dHistPIDprobaPionVsMCPtBach" , "Comb. PID proba to be #pi^{#pm} Vs MC Bach. Pt ; Pt_{MC}(Bach.) (GeV/c); Comb. PID Proba (Bach. = #pi^{#pm})", 100, 0.0, 5.0, 110, 0.0, 1.10 );
 981         fListHistCascade->Add(f2dHistPIDprobaPionVsMCPtBach);
 982   }
 983
 984
 985                 // - Effective mass histos with perfect MC PID on the bachelor
 986
 987   if (! fHistMassWithMcPIDXiMinus) {
 988     fHistMassWithMcPIDXiMinus = new TH1F("fHistMassWithMcPIDXiMinus","#Xi^{-} candidates, with Bach. MC PID;M( #Lambda , #pi^{-} ) (GeV/c^{2});Counts", 400,1.2,2.0);
 989     fListHistCascade->Add(fHistMassWithMcPIDXiMinus);
 990   }
 991
 992   if (! fHistMassWithMcPIDXiPlus) {
 993     fHistMassWithMcPIDXiPlus = new TH1F("fHistMassWithMcPIDXiPlus","#bar{#Xi}^{+} candidates, with Bach. MC PID;M( #bar{#Lambda}^{0} , #pi^{+} ) (GeV/c^{2});Counts",400,1.2,2.0);
 994     fListHistCascade->Add(fHistMassWithMcPIDXiPlus);
 995   }
 996
 997   if (! fHistMassWithMcPIDOmegaMinus) {
 998         fHistMassWithMcPIDOmegaMinus = new TH1F("fHistMassWithMcPIDOmegaMinus","#Omega^{-} candidates, with Bach. MC PID;M( #Lambda , K^{-} ) (GeV/c^{2});Counts", 500,1.5,2.5);
 999     fListHistCascade->Add(fHistMassWithMcPIDOmegaMinus);
1000   }
1001
1002   if (! fHistMassWithMcPIDOmegaPlus) {
1003         fHistMassWithMcPIDOmegaPlus = new TH1F("fHistMassWithMcPIDOmegaPlus","#bar{#Omega}^{+} candidates, with Bach. MC PID;M( #bar{#Lambda}^{0} , K^{+} ) (GeV/c^{2});Counts", 500,1.5,2.5);
1004     fListHistCascade->Add(fHistMassWithMcPIDOmegaPlus);
1005   }
1006
1007
1008                 // - Effective mass histos for cascades candidates ASSOCIATED with MC.
1009
1010   if (! fHistAsMCMassXiMinus) {
1011           fHistAsMCMassXiMinus = new TH1F("fHistAsMCMassXiMinus","#Xi^{-} candidates associated to MC;M( #Lambda , #pi^{-} ) (GeV/c^{2});Counts", 400,1.2,2.0);
1012           fListHistCascade->Add(fHistAsMCMassXiMinus);
1013   }
1014
1015   if (! fHistAsMCMassXiPlus) {
1016           fHistAsMCMassXiPlus = new TH1F("fHistAsMCMassXiPlus","#bar{#Xi}^{+} candidates associated to MC;M( #bar{#Lambda}^{0} , #pi^{+} ) (GeV/c^{2});Counts",400,1.2,2.0);
1017           fListHistCascade->Add(fHistAsMCMassXiPlus);
1018   }
1019
1020   if (! fHistAsMCMassOmegaMinus) {
1021           fHistAsMCMassOmegaMinus = new TH1F("fHistAsMCMassOmegaMinus","#Omega^{-} candidates associated to MC;M( #Lambda , K^{-} ) (GeV/c^{2});Counts", 500,1.5,2.5);
1022           fListHistCascade->Add(fHistAsMCMassOmegaMinus);
1023   }
1024
1025   if (! fHistAsMCMassOmegaPlus) {
1026           fHistAsMCMassOmegaPlus = new TH1F("fHistAsMCMassOmegaPlus","#bar{#Omega}^{+} candidates associated to MC;M( #bar{#Lambda}^{0} , K^{+} ) (GeV/c^{2});Counts", 500,1.5,2.5);
1027           fListHistCascade->Add(fHistAsMCMassOmegaPlus);
1028   }
1029
1030
1031                 // -  Generated Pt Vs generated Y of the cascade candidates associated with MC
1032                 //     + having the proper maximum proba of combined PID for the bachelor
1033
1034   if (!f2dHistAsMCandCombPIDGenPtVsGenYXiMinus) {
1035      f2dHistAsMCandCombPIDGenPtVsGenYXiMinus = new TH2F("f2dHistAsMCandCombPIDGenPtVsGenYXiMinus", "MC P_{t} Vs MC Y of #Xi^{-} (associated+Bach.PID); Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
1036      fListHistCascade->Add(f2dHistAsMCandCombPIDGenPtVsGenYXiMinus);
1037   }
1038
1039   if (!f2dHistAsMCandCombPIDGenPtVsGenYXiPlus) {
1040      f2dHistAsMCandCombPIDGenPtVsGenYXiPlus = new TH2F("f2dHistAsMCandCombPIDGenPtVsGenYXiPlus", "MC P_{t} Vs MC Y of #bar{#Xi}^{+} (associated+Bach.PID); Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
1041      fListHistCascade->Add(f2dHistAsMCandCombPIDGenPtVsGenYXiPlus);
1042   }
1043
1044   if (!f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus) {
1045      f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus = new TH2F("f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus", "MC P_{t} Vs MC Y of #Omega^{-} (associated+Bach.PID); Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
1046      fListHistCascade->Add(f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus);
1047   }
1048
1049   if (!f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus) {
1050      f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus = new TH2F("f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus", "MC P_{t} Vs MC Y of #bar{#Omega}^{+} (associated+Bach.PID); Pt_{MC} (GeV/c); Y_{MC}", 200, 0., 10., 220, -1.1, 1.1);
1051      fListHistCascade->Add(f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus);
1052   }
1053
1054
1055                 // - Generated Pt Vs Generated Y, for the cascade candidates associated with MC
1056
1057   if (!f2dHistAsMCGenPtVsGenYXiMinus) {
1058         f2dHistAsMCGenPtVsGenYXiMinus = new TH2F("f2dHistAsMCGenPtVsGenYXiMinus", "MC P_{t} Vs MC Y of gen. #Xi^{-} (associated);Pt_{MC} (GeV/c); Rapidity, Y_{MC}",200, 0., 10., 220, -1.1, 1.1);
1059           fListHistCascade->Add(f2dHistAsMCGenPtVsGenYXiMinus );
1060   }
1061
1062   if (!f2dHistAsMCGenPtVsGenYXiPlus) {
1063           f2dHistAsMCGenPtVsGenYXiPlus = new TH2F("f2dHistAsMCGenPtVsGenYXiPlus", "MC P_{t} Vs MC Y of gen. #bar{#Xi}^{+} (associated);Pt_{MC} (GeV/c); Rapidity, Y_{MC}",200, 0., 10., 220, -1.1, 1.1);
1064           fListHistCascade->Add(f2dHistAsMCGenPtVsGenYXiPlus );
1065   }
1066
1067   if (!f2dHistAsMCGenPtVsGenYOmegaMinus) {
1068           f2dHistAsMCGenPtVsGenYOmegaMinus = new TH2F("f2dHistAsMCGenPtVsGenYOmegaMinus", "MC P_{t} Vs MC Y of gen. #Omega^{-} (associated);Pt_{MC} (GeV/c); Rapidity, Y_{MC}",200, 0., 10., 220, -1.1, 1.1);
1069           fListHistCascade->Add(f2dHistAsMCGenPtVsGenYOmegaMinus );
1070   }
1071
1072   if (!f2dHistAsMCGenPtVsGenYOmegaPlus) {
1073           f2dHistAsMCGenPtVsGenYOmegaPlus = new TH2F("f2dHistAsMCGenPtVsGenYOmegaPlus", "MC P_{t} Vs MC Y of gen. #bar{#Omega}^{+} (associated);Pt_{MC} (GeV/c); Rapidity, Y_{MC}",200, 0., 10., 220, -1.1, 1.1);
1074           fListHistCascade->Add(f2dHistAsMCGenPtVsGenYOmegaPlus );
1075   }
1076
1077
1078                 // - Generated Eta of the the cascade candidates associated with MC
1079   if (!fHistAsMCGenEtaXiMinus) {
1080           fHistAsMCGenEtaXiMinus = new TH1F("fHistAsMCGenEtaXiMinus", "#eta of gen. #Xi^{-} (associated);#eta;Number of Casc", 100, -5, 5);
1081           fListHistCascade->Add( fHistAsMCGenEtaXiMinus );
1082   }
1083
1084   if (!fHistAsMCGenEtaXiPlus) {
1085           fHistAsMCGenEtaXiPlus = new TH1F("fHistAsMCGenEtaXiPlus", "#eta of gen. #bar{#Xi}^{+} (associated);#eta;Number of Casc", 100, -5, 5);
1086           fListHistCascade->Add( fHistAsMCGenEtaXiPlus );
1087   }
1088
1089   if (!fHistAsMCGenEtaOmegaMinus) {
1090           fHistAsMCGenEtaOmegaMinus = new TH1F("fHistAsMCGenEtaOmegaMinus", "#eta of gen. #Omega^{-} (associated);#eta;Number of Casc", 100, -5, 5);
1091           fListHistCascade->Add( fHistAsMCGenEtaOmegaMinus );
1092   }
1093
1094   if (!fHistAsMCGenEtaOmegaPlus) {
1095           fHistAsMCGenEtaOmegaPlus = new TH1F("fHistAsMCGenEtaOmegaPlus", "#eta of gen. #bar{#Omega}^{+} (associated);#eta;Number of Casc", 100, -5, 5);
1096           fListHistCascade->Add( fHistAsMCGenEtaOmegaPlus );
1097   }
1098
1099
1100
1101                 // - Resolution in Pt as function of generated Pt
1102
1103   if(! f2dHistAsMCResPtXiMinus) {
1104           f2dHistAsMCResPtXiMinus = new TH2F( "f2dHistAsMCResPtXiMinus", "Resolution in Pt reconstruction for #Xi^{-}; Pt_{MC} (GeV/c); (Pt_{reco} - Pt_{MC}) / Pt_{MC}", 200, 0., 10., 200, -0.1, 0.1);
1105           fListHistCascade->Add(f2dHistAsMCResPtXiMinus);
1106   }
1107
1108   if(! f2dHistAsMCResPtXiPlus) {
1109           f2dHistAsMCResPtXiPlus = new TH2F( "f2dHistAsMCResPtXiPlus", "Resolution in Pt reconstruction for #bar{#Xi}^{+}; Pt_{MC} (GeV/c); (Pt_{reco} - Pt_{MC}) / Pt_{MC}", 200, 0., 10., 200, -0.1, 0.1);
1110           fListHistCascade->Add(f2dHistAsMCResPtXiPlus);
1111   }
1112
1113   if(! f2dHistAsMCResPtOmegaMinus) {
1114           f2dHistAsMCResPtOmegaMinus = new TH2F( "f2dHistAsMCResPtOmegaMinus", "Resolution in Pt reconstruction for #Omega^{-}; Pt_{MC} (GeV/c); (Pt_{reco} - Pt_{MC}) / Pt_{MC}", 200, 0., 10., 200, -0.1, 0.1);
1115           fListHistCascade->Add(f2dHistAsMCResPtOmegaMinus);
1116   }
1117
1118   if(! f2dHistAsMCResPtOmegaPlus) {
1119           f2dHistAsMCResPtOmegaPlus = new TH2F( "f2dHistAsMCResPtOmegaPlus", "Resolution in Pt reconstruction for #bar{#Omega}^{+}; Pt_{MC} (GeV/c); (Pt_{reco} - Pt_{MC}) / Pt_{MC}", 200, 0., 10., 200, -0.1, 0.1);
1120           fListHistCascade->Add(f2dHistAsMCResPtOmegaPlus);
1121   }
1122
1123                 // - Resolution in R(2D) as function of generated R
1124
1125   if(! f2dHistAsMCResRXiMinus) {
1126           f2dHistAsMCResRXiMinus = new TH2F( "f2dHistAsMCResRXiMinus", "Resolution in transv. position for #Xi^{-}; R_{MC} (cm); (R_{reco} - R_{MC}) / R_{MC}", 450, 0., 45.0, 240, -0.3, 0.3);
1127           fListHistCascade->Add(f2dHistAsMCResRXiMinus);
1128   }
1129
1130   if(! f2dHistAsMCResRXiPlus) {
1131           f2dHistAsMCResRXiPlus = new TH2F( "f2dHistAsMCResRXiPlus", "Resolution in transv. position for #bar{#Xi}^{+}; R_{MC} (cm); (R_{reco} - R_{MC}) / R_{MC}", 450, 0., 45.0, 240, -0.3, 0.3);
1132           fListHistCascade->Add(f2dHistAsMCResRXiPlus);
1133   }
1134
1135   if(! f2dHistAsMCResROmegaMinus) {
1136           f2dHistAsMCResROmegaMinus = new TH2F( "f2dHistAsMCResROmegaMinus", "Resolution in transv. position for #Omega^{-}; R_{MC} (cm); (R_{reco} - R_{MC}) / R_{MC}", 450, 0., 45.0, 240, -0.3, 0.3);
1137           fListHistCascade->Add(f2dHistAsMCResROmegaMinus);
1138   }
1139
1140   if(! f2dHistAsMCResROmegaPlus) {
1141           f2dHistAsMCResROmegaPlus = new TH2F( "f2dHistAsMCResROmegaPlus", "Resolution in transv. position for #bar{#Omega}^{+}; R_{MC} (cm); (R_{reco} - R_{MC}) / R_{MC}", 450, 0., 45.0, 240, -0.3, 0.3);
1142           fListHistCascade->Add(f2dHistAsMCResROmegaPlus);
1143   }
1144
1145                 // - Resolution in phi as function of generated Pt
1146
1147   if(! f2dHistAsMCResPhiXiMinus) {
1148           f2dHistAsMCResPhiXiMinus = new TH2F( "f2dHistAsMCResPhiXiMinus", "Resolution in #phi for #Xi^{-}; Pt_{MC} (GeV/c); #phi(MC) - #phi(reco)   (deg)", 200, 0., 10., 60, -30., 30.);
1149           fListHistCascade->Add(f2dHistAsMCResPhiXiMinus);
1150   }
1151
1152   if(! f2dHistAsMCResPhiXiPlus) {
1153           f2dHistAsMCResPhiXiPlus = new TH2F( "f2dHistAsMCResPhiXiPlus", "Resolution in #phi for #bar{#Xi}^{+}; Pt_{MC} (GeV/c); #phi(MC) - #phi(reco)   (deg)", 200, 0., 10., 60, -30., 30.);
1154           fListHistCascade->Add(f2dHistAsMCResPhiXiPlus);
1155   }
1156
1157   if(! f2dHistAsMCResPhiOmegaMinus) {
1158           f2dHistAsMCResPhiOmegaMinus = new TH2F( "f2dHistAsMCResPhiOmegaMinus", "Resolution in #phi for #Omega^{-}; Pt_{MC} (GeV/c); #phi(MC) - #phi(reco)   (deg)", 200, 0., 10., 60, -30., 30.);
1159           fListHistCascade->Add(f2dHistAsMCResPhiOmegaMinus);
1160   }
1161
1162   if(! f2dHistAsMCResPhiOmegaPlus) {
1163           f2dHistAsMCResPhiOmegaPlus = new TH2F( "f2dHistAsMCResPhiOmegaPlus", "Resolution in #phi for #bar{#Omega}^{+}; Pt_{MC} (GeV/c); #phi(MC) - #phi(reco)   (deg)", 200, 0., 10., 60, -30., 30.);
1164           fListHistCascade->Add(f2dHistAsMCResPhiOmegaPlus);
1165   }
1166
1167
1168                 // - Correlation in Pt between the cascade and its (anti)proton daughter
1169   if(! f2dHistAsMCPtProtonVsPtXiMinus) {
1170         f2dHistAsMCPtProtonVsPtXiMinus = new TH2F( "f2dHistAsMCPtProtonVsPtXiMinus", "Correlation Pt(p) Vs Pt(#Xi^{-}), associated to MC; Pt_{MC}(p) (GeV/c); Pt_{MC}(#Xi^{-}) (GeV/c)", 100, 0., 10., 100, 0., 10.);
1171         fListHistCascade->Add(f2dHistAsMCPtProtonVsPtXiMinus);
1172   }
1173
1174   if(! f2dHistAsMCPtAntiProtonVsPtXiPlus) {
1175         f2dHistAsMCPtAntiProtonVsPtXiPlus = new TH2F( "f2dHistAsMCPtAntiProtonVsPtXiPlus", "Correlation Pt(#bar{p}) Vs Pt(#bar{#Xi}^{+}), associated to MC; Pt_{MC}(#bar{p}) (GeV/c); Pt_{MC}(#bar{#Xi}^{+}) (GeV/c)", 100, 0., 10., 100, 0., 10.);
1176         fListHistCascade->Add(f2dHistAsMCPtAntiProtonVsPtXiPlus);
1177   }
1178
1179   if(! f2dHistAsMCPtProtonVsPtOmegaMinus) {
1180         f2dHistAsMCPtProtonVsPtOmegaMinus = new TH2F( "f2dHistAsMCPtProtonVsPtOmegaMinus", "Correlation Pt(p) Vs Pt(#Omega^{-}), associated to MC; Pt_{MC}(p) (GeV/c); Pt_{MC}(#Omega^{-}) (GeV/c)", 100, 0., 10., 100, 0., 10.);
1181         fListHistCascade->Add(f2dHistAsMCPtProtonVsPtOmegaMinus);
1182   }
1183
1184
1185   if(! f2dHistAsMCPtAntiProtonVsPtOmegaPlus) {
1186         f2dHistAsMCPtAntiProtonVsPtOmegaPlus = new TH2F( "f2dHistAsMCPtAntiProtonVsPtOmegaPlus", "Correlation Pt(#bar{p}) Vs Pt(#bar{#Omega}^{+}), associated to MC; Pt_{MC}(#bar{p}) (GeV/c); Pt_{MC}(#bar{#Omega}^{+}) (GeV/c)", 100, 0., 10., 100, 0., 10.);
1187         fListHistCascade->Add(f2dHistAsMCPtAntiProtonVsPtOmegaPlus);
1188   }
1189
1190
1191
1192                 // - PID container
1193 if(! fCFContCascadePIDAsXiMinus)  {
1194   const Int_t  lNbSteps      =  7 ;
1195   const Int_t  lNbVariables  =  4 ;
1196
1197   //array for the number of bins in each dimension :
1198   Int_t lNbBinsPerVar[4] = {0};
1199   lNbBinsPerVar[0] = 200;
1200   lNbBinsPerVar[1] = 75;
1201   lNbBinsPerVar[2] = 44;
1202   lNbBinsPerVar[3] = 250;
1203
1204
1205   fCFContCascadePIDAsXiMinus = new AliCFContainer("fCFContCascadePIDAsXiMinus","Pt_{cascade} Vs M_{#Xi^{-} candidates} Vs Y_{#Xi}", lNbSteps, lNbVariables, lNbBinsPerVar );
1206
1207   //setting the bin limits (valid  for v4-18-10-AN)
1208   fCFContCascadePIDAsXiMinus->SetBinLimits(0,   0.0  ,  10.0  );        // Pt(Cascade)
1209   fCFContCascadePIDAsXiMinus->SetBinLimits(1,   1.25 ,   1.40 );        // Xi Effective mass
1210   fCFContCascadePIDAsXiMinus->SetBinLimits(2,  -1.1  ,   1.1  );        // Rapidity
1211   if(fCollidingSystems)
1212         fCFContCascadePIDAsXiMinus->SetBinLimits(3, 0.0, 20000.0  );    // nTrackPrimaryMultiplicity
1213   else
1214         fCFContCascadePIDAsXiMinus->SetBinLimits(3, 0.0, 250.0  );     // nTrackPrimaryMultiplicity
1215
1216   // Setting the step title : one per PID case
1217   fCFContCascadePIDAsXiMinus->SetStepTitle(0, "No PID");
1218   fCFContCascadePIDAsXiMinus->SetStepTitle(1, "TPC PID / 4-#sigma cut on Bachelor track");
1219   fCFContCascadePIDAsXiMinus->SetStepTitle(2, "TPC PID / 4-#sigma cut on Bachelor+Baryon tracks");
1220   fCFContCascadePIDAsXiMinus->SetStepTitle(3, "TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks");
1221   fCFContCascadePIDAsXiMinus->SetStepTitle(4, "Comb. PID / Bachelor");
1222   fCFContCascadePIDAsXiMinus->SetStepTitle(5, "Comb. PID / Bachelor+Baryon");
1223   fCFContCascadePIDAsXiMinus->SetStepTitle(6, "Comb. PID / Bachelor+Baryon+Meson");
1224
1225   // Setting the variable title, per axis
1226   fCFContCascadePIDAsXiMinus->SetVarTitle(0, "Pt_{cascade} (GeV/c)");
1227   fCFContCascadePIDAsXiMinus->SetVarTitle(1, "M( #Lambda , #pi^{-} ) (GeV/c^{2})");
1228   fCFContCascadePIDAsXiMinus->SetVarTitle(2, "Y_{#Xi}");
1229   fCFContCascadePIDAsXiMinus->SetVarTitle(3, "Primary Track Multiplicity");
1230
1231   fListHistCascade->Add(fCFContCascadePIDAsXiMinus);
1232
1233 }
1234
1235 if(! fCFContCascadePIDAsXiPlus)  {
1236   const Int_t  lNbSteps      =  7 ;
1237   const Int_t  lNbVariables  =  4 ;
1238
1239   //array for the number of bins in each dimension :
1240   Int_t lNbBinsPerVar[4] = {0};
1241   lNbBinsPerVar[0] = 200;
1242   lNbBinsPerVar[1] = 75;
1243   lNbBinsPerVar[2] = 44;
1244   lNbBinsPerVar[3] = 250;
1245
1246
1247   fCFContCascadePIDAsXiPlus = new AliCFContainer("fCFContCascadePIDAsXiPlus","Pt_{cascade} Vs M_{#bar{#Xi}^{+} candidates} Vs Y_{#Xi}", lNbSteps, lNbVariables, lNbBinsPerVar );
1248
1249
1250   //setting the bin limits (valid  for v4-18-10-AN)
1251   fCFContCascadePIDAsXiPlus->SetBinLimits(0,   0.0  ,  10.0  );         // Pt(Cascade)
1252   fCFContCascadePIDAsXiPlus->SetBinLimits(1,   1.25 ,   1.40 );         // Xi Effective mass
1253   fCFContCascadePIDAsXiPlus->SetBinLimits(2,  -1.1  ,   1.1  );         // Rapidity
1254   if(fCollidingSystems)
1255         fCFContCascadePIDAsXiPlus->SetBinLimits(3, 0.0, 20000.0  );    // nTrackPrimaryMultiplicity
1256   else
1257         fCFContCascadePIDAsXiPlus->SetBinLimits(3, 0.0, 250.0  );     // nTrackPrimaryMultiplicity
1258
1259   // Setting the step title : one per PID case
1260   fCFContCascadePIDAsXiPlus->SetStepTitle(0, "No PID");
1261   fCFContCascadePIDAsXiPlus->SetStepTitle(1, "TPC PID / 4-#sigma cut on Bachelor track");
1262   fCFContCascadePIDAsXiPlus->SetStepTitle(2, "TPC PID / 4-#sigma cut on Bachelor+Baryon tracks");
1263   fCFContCascadePIDAsXiPlus->SetStepTitle(3, "TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks");
1264   fCFContCascadePIDAsXiPlus->SetStepTitle(4, "Comb. PID / Bachelor");
1265   fCFContCascadePIDAsXiPlus->SetStepTitle(5, "Comb. PID / Bachelor+Baryon");
1266   fCFContCascadePIDAsXiPlus->SetStepTitle(6, "Comb. PID / Bachelor+Baryon+Meson");
1267
1268   // Setting the variable title, per axis
1269   fCFContCascadePIDAsXiPlus->SetVarTitle(0, "Pt_{cascade} (GeV/c)");
1270   fCFContCascadePIDAsXiPlus->SetVarTitle(1, "M( #Lambda , #pi^{+} ) (GeV/c^{2})");
1271   fCFContCascadePIDAsXiPlus->SetVarTitle(2, "Y_{#Xi}");
1272   fCFContCascadePIDAsXiPlus->SetVarTitle(3, "Primary Track Multiplicity");
1273
1274   fListHistCascade->Add(fCFContCascadePIDAsXiPlus);
1275
1276 }
1277
1278
1279 if(! fCFContCascadePIDAsOmegaMinus)  {
1280   const Int_t  lNbSteps      =  7 ;
1281   const Int_t  lNbVariables  =  4 ;
1282
1283   //array for the number of bins in each dimension :
1284   Int_t lNbBinsPerVar[4] = {0};
1285   lNbBinsPerVar[0] = 200;
1286   lNbBinsPerVar[1] = 60;
1287   lNbBinsPerVar[2] = 44;
1288   lNbBinsPerVar[3] = 250;
1289
1290
1291   fCFContCascadePIDAsOmegaMinus = new AliCFContainer("fCFContCascadePIDAsOmegaMinus","Pt_{cascade} Vs M_{#Omega^{-} candidates} Vs Y_{#Omega}", lNbSteps, lNbVariables, lNbBinsPerVar );
1292
1293
1294   //setting the bin limits (valid  for v4-18-10-AN)
1295   fCFContCascadePIDAsOmegaMinus->SetBinLimits(0,   0.0  ,  10.0  );     // Pt(Cascade)
1296   fCFContCascadePIDAsOmegaMinus->SetBinLimits(1,   1.62 ,   1.74 );     // Omega Effective mass
1297   fCFContCascadePIDAsOmegaMinus->SetBinLimits(2,  -1.1  ,   1.1  );     // Rapidity
1298   if(fCollidingSystems)
1299         fCFContCascadePIDAsOmegaMinus->SetBinLimits(3, 0.0, 20000.0  );    // nTrackPrimaryMultiplicity
1300   else
1301         fCFContCascadePIDAsOmegaMinus->SetBinLimits(3, 0.0, 250.0  );     // nTrackPrimaryMultiplicity
1302
1303   // Setting the step title : one per PID case
1304   fCFContCascadePIDAsOmegaMinus->SetStepTitle(0, "No PID");
1305   fCFContCascadePIDAsOmegaMinus->SetStepTitle(1, "TPC PID / 4-#sigma cut on Bachelor track");
1306   fCFContCascadePIDAsOmegaMinus->SetStepTitle(2, "TPC PID / 4-#sigma cut on Bachelor+Baryon tracks");
1307   fCFContCascadePIDAsOmegaMinus->SetStepTitle(3, "TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks");
1308   fCFContCascadePIDAsOmegaMinus->SetStepTitle(4, "Comb. PID / Bachelor");
1309   fCFContCascadePIDAsOmegaMinus->SetStepTitle(5, "Comb. PID / Bachelor+Baryon");
1310   fCFContCascadePIDAsOmegaMinus->SetStepTitle(6, "Comb. PID / Bachelor+Baryon+Meson");
1311
1312   // Setting the variable title, per axis
1313   fCFContCascadePIDAsOmegaMinus->SetVarTitle(0, "Pt_{cascade} (GeV/c)");
1314   fCFContCascadePIDAsOmegaMinus->SetVarTitle(1, "M( #Lambda , K^{-} ) (GeV/c^{2})");
1315   fCFContCascadePIDAsOmegaMinus->SetVarTitle(2, "Y_{#Omega}");
1316   fCFContCascadePIDAsOmegaMinus->SetVarTitle(3, "Primary Track Multiplicity");
1317
1318   fListHistCascade->Add(fCFContCascadePIDAsOmegaMinus);
1319
1320 }
1321
1322 if(! fCFContCascadePIDAsOmegaPlus)  {
1323   const Int_t  lNbSteps      =  7 ;
1324   const Int_t  lNbVariables  =  4 ;
1325
1326   //array for the number of bins in each dimension :
1327   Int_t lNbBinsPerVar[4]= {0};
1328   lNbBinsPerVar[0] = 200;
1329   lNbBinsPerVar[1] = 60;
1330   lNbBinsPerVar[2] = 44;
1331   lNbBinsPerVar[3] = 250;
1332
1333
1334   fCFContCascadePIDAsOmegaPlus = new AliCFContainer("fCFContCascadePIDAsOmegaPlus","Pt_{cascade} Vs M_{#bar{#Omega}^{+} candidates} Vs Y_{#Omega}", lNbSteps, lNbVariables, lNbBinsPerVar );
1335
1336
1337   //setting the bin limits (valid  for v4-18-10-AN)
1338   fCFContCascadePIDAsOmegaPlus->SetBinLimits(0,   0.0  ,  10.0  );      // Pt(Cascade)
1339   fCFContCascadePIDAsOmegaPlus->SetBinLimits(1,   1.62 ,   1.74 );      // Omega Effective mass
1340   fCFContCascadePIDAsOmegaPlus->SetBinLimits(2,  -1.1  ,   1.1  );      // Rapidity
1341   if(fCollidingSystems)
1342         fCFContCascadePIDAsOmegaPlus->SetBinLimits(3, 0.0, 20000.0  );    // nTrackPrimaryMultiplicity
1343   else
1344         fCFContCascadePIDAsOmegaPlus->SetBinLimits(3, 0.0, 250.0  );     // nTrackPrimaryMultiplicity
1345
1346   // Setting the step title : one per PID case
1347   fCFContCascadePIDAsOmegaPlus->SetStepTitle(0, "No PID");
1348   fCFContCascadePIDAsOmegaPlus->SetStepTitle(1, "TPC PID / 4-#sigma cut on Bachelor track");
1349   fCFContCascadePIDAsOmegaPlus->SetStepTitle(2, "TPC PID / 4-#sigma cut on Bachelor+Baryon tracks");
1350   fCFContCascadePIDAsOmegaPlus->SetStepTitle(3, "TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks");
1351   fCFContCascadePIDAsOmegaPlus->SetStepTitle(4, "Comb. PID / Bachelor");
1352   fCFContCascadePIDAsOmegaPlus->SetStepTitle(5, "Comb. PID / Bachelor+Baryon");
1353   fCFContCascadePIDAsOmegaPlus->SetStepTitle(6, "Comb. PID / Bachelor+Baryon+Meson");
1354
1355   // Setting the variable title, per axis
1356   fCFContCascadePIDAsOmegaPlus->SetVarTitle(0, "Pt_{cascade} (GeV/c)");
1357   fCFContCascadePIDAsOmegaPlus->SetVarTitle(1, "M( #Lambda , K^{+} ) (GeV/c^{2})");
1358   fCFContCascadePIDAsOmegaPlus->SetVarTitle(2, "Y_{#Omega}");
1359   fCFContCascadePIDAsOmegaPlus->SetVarTitle(3, "Primary Track Multiplicity");
1360
1361   fListHistCascade->Add(fCFContCascadePIDAsOmegaPlus);
1362
1363 }
1364
1365 // Part 3 : Towards the optimisation of topological selections -------
1366 if(! fCFContAsCascadeCuts){
1367
1368         // Container meant to store all the relevant distributions corresponding to the cut variables.
1369         // So far, 20 variables have been identified.
1370         // The following will be done in quite a brut force way ...
1371         // FIXME Improvement expected later (before Pb-Pb data at least)
1372         //          - Define a user binning to have less bins in each dimension
1373         //          - boolean for enabling/disbaling this CFContainer
1374   const Int_t  lNbSteps      =  4 ;
1375   const Int_t  lNbVariables  =  20 ;
1376
1377   //array for the number of bins in each dimension :
1378   Int_t lNbBinsPerVar[20] = {0};
1379   lNbBinsPerVar[0]  = 25;
1380   lNbBinsPerVar[1]  = 25;
1381   lNbBinsPerVar[2]  = 20;
1382   lNbBinsPerVar[3]  = 40;
1383   lNbBinsPerVar[4]  = 30;
1384   lNbBinsPerVar[5]  = 25;
1385
1386   lNbBinsPerVar[6]  = 20;
1387   lNbBinsPerVar[7]  = 40;
1388   lNbBinsPerVar[8]  = 40;
1389   lNbBinsPerVar[9]  = 25;
1390   lNbBinsPerVar[10] = 25;
1391
1392   lNbBinsPerVar[11] = 75; // 2-MeV/c2 bins
1393   lNbBinsPerVar[12] = 60; // 2-MeV/c2 bins
1394
1395   lNbBinsPerVar[13] = 100;
1396
1397   lNbBinsPerVar[14] = 44; // 0.05 in rapidity units
1398   lNbBinsPerVar[15] = 44; // 0.05 in rapidity units
1399
1400   lNbBinsPerVar[16] = 20;
1401
1402   lNbBinsPerVar[17] = 50;
1403   lNbBinsPerVar[18] = 100;
1404   lNbBinsPerVar[19] = 24;
1405
1406
1407   fCFContAsCascadeCuts = new AliCFContainer("fCFContAsCascadeCuts","Cut Container for Asso. Cascades", lNbSteps, lNbVariables, lNbBinsPerVar );
1408
1409   //0
1410   Double_t *lBinLim0  = new Double_t[ lNbBinsPerVar[0]+1 ];
1411         for(Int_t i=0; i< lNbBinsPerVar[0];i++)  lBinLim0[i]  = (Double_t)0.0   + (4.8  - 0.0 )/(lNbBinsPerVar[0]-1)  * (Double_t)i ;
1412         lBinLim0[ lNbBinsPerVar[0]  ] = 20.0;
1413   fCFContAsCascadeCuts -> SetBinLimits(0,  lBinLim0 );            // DcaXiDaughters : 0.0 to 5.0
1414   delete [] lBinLim0;
1415   //1
1416   Double_t *lBinLim1  = new Double_t[ lNbBinsPerVar[1]+1 ];
1417         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 ;
1418         lBinLim1[ lNbBinsPerVar[1]  ] = 100.0;
1419   fCFContAsCascadeCuts -> SetBinLimits(1,  lBinLim1 );            // DcaBachToPrimVertexXi : 0.0 to 0.25
1420   delete [] lBinLim1;
1421   //2
1422   Double_t *lBinLim2  = new Double_t[ lNbBinsPerVar[2]+1 ];
1423         for(Int_t i=1; i< lNbBinsPerVar[2]+1;i++)   lBinLim2[i]  = (Double_t)0.81   + (1.0  - 0.81 )/(lNbBinsPerVar[2]-1)  * (Double_t) (i-1) ;
1424         lBinLim2[0] = 0.0;
1425   fCFContAsCascadeCuts -> SetBinLimits(2,  lBinLim2 );            // XiCosineOfPointingAngle : 0.80 to 1.0
1426   delete [] lBinLim2;
1427   //3
1428   Double_t *lBinLim3  = new Double_t[ lNbBinsPerVar[3]+1 ];
1429         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 ;
1430         lBinLim3[ lNbBinsPerVar[3]  ] = 110.0;
1431   fCFContAsCascadeCuts -> SetBinLimits(3,  lBinLim3 );            // XiRadius : 0.0 to 4.0
1432   delete [] lBinLim3;
1433   //4
1434   fCFContAsCascadeCuts->SetBinLimits(4,    1.1  ,  1.13 );        // InvMassLambdaAsCascDghter
1435   //5
1436   Double_t *lBinLim5  = new Double_t[ lNbBinsPerVar[5]+1 ];
1437         for(Int_t i=0; i< lNbBinsPerVar[5];i++)   lBinLim5[i]  = (Double_t)0.0   + (4.8  - 0.0 )/(lNbBinsPerVar[5]-1)  * (Double_t)i ;
1438         lBinLim5[ lNbBinsPerVar[5]  ] = 20.0;
1439   fCFContAsCascadeCuts -> SetBinLimits(5,  lBinLim5 );            // DcaV0DaughtersXi : 0.0 to 5.0
1440   delete [] lBinLim5;
1441
1442
1443   //6
1444   Double_t *lBinLim6  = new Double_t[ lNbBinsPerVar[6]+1 ];
1445         for(Int_t i=1; i< lNbBinsPerVar[6]+1 ;i++)   lBinLim6[i]  = (Double_t)0.81   + (1.0  - 0.81 )/(lNbBinsPerVar[6]-1)  * (Double_t) (i-1) ;
1446         lBinLim6[0] = 0.0;
1447   fCFContAsCascadeCuts -> SetBinLimits(6,  lBinLim6 );            // V0CosineOfPointingAngleXi : 0.80 to 1.0
1448   delete [] lBinLim6;
1449   //7
1450   Double_t *lBinLim7  = new Double_t[ lNbBinsPerVar[7]+1 ];
1451         for(Int_t i=0; i< lNbBinsPerVar[7];i++)   lBinLim7[i]  = (Double_t)0.0   + (7.8  - 0.0 )/(lNbBinsPerVar[7]-1)  * (Double_t)i ;
1452         lBinLim7[ lNbBinsPerVar[7]  ] = 100.0;
1453   fCFContAsCascadeCuts -> SetBinLimits(7,  lBinLim7 );            // V0RadiusXi : 0.0 to 8.0
1454   delete [] lBinLim7;
1455   //8
1456   Double_t *lBinLim8  = new Double_t[ lNbBinsPerVar[8]+1 ];
1457         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 ;
1458         lBinLim8[ lNbBinsPerVar[8]  ] = 100.0;
1459   fCFContAsCascadeCuts -> SetBinLimits(8,  lBinLim8 );            // DcaV0ToPrimVertexXi : 0.0 to 0.4
1460   delete [] lBinLim8;
1461   //9
1462   Double_t *lBinLim9  = new Double_t[ lNbBinsPerVar[9]+1 ];
1463         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 ;
1464         lBinLim9[ lNbBinsPerVar[9]  ] = 100.0;
1465   fCFContAsCascadeCuts -> SetBinLimits(9,  lBinLim9 );            // DcaPosToPrimVertexXi : 0.0 to 0.25
1466   delete [] lBinLim9;
1467   //10
1468   Double_t *lBinLim10  = new Double_t[ lNbBinsPerVar[10]+1 ];
1469         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 ;
1470         lBinLim10[ lNbBinsPerVar[10]  ] = 100.0;
1471   fCFContAsCascadeCuts -> SetBinLimits(10,  lBinLim10 );            // DcaPosToPrimVertexXi : 0.0 to 0.25
1472   delete [] lBinLim10;
1473
1474
1475   //11
1476   fCFContAsCascadeCuts->SetBinLimits(11,   1.25 ,  1.40 );      // InvMassXi
1477   fCFContAsCascadeCuts->SetBinLimits(12,   1.62  , 1.74 );      // InvMassOmega
1478   fCFContAsCascadeCuts->SetBinLimits(13,   0.0  , 10.0  );      // XiTransvMom
1479   fCFContAsCascadeCuts->SetBinLimits(14,  -1.1  ,  1.1  );      // Y(Xi)
1480   fCFContAsCascadeCuts->SetBinLimits(15,  -1.1  ,  1.1  );      // Y(Omega)
1481   fCFContAsCascadeCuts->SetBinLimits(16, -10.0  , 10.0  );      // BestPrimaryVtxPosZ
1482   if(fCollidingSystems){
1483         fCFContAsCascadeCuts->SetBinLimits(17,   0.0, 10000.0  );    // nTrackPrimaryMultiplicity
1484         fCFContAsCascadeCuts->SetBinLimits(18,   0.0, 10000.0  );    // nITSandTPCtracksAndSPDtracklets
1485   }
1486   else{
1487         //fCFContAsCascadeCuts->SetBinLimits(17,   0.0, 250.0  );     // nTrackPrimaryMultiplicity
1488         Double_t *lBinLim17  = new Double_t[ lNbBinsPerVar[17]+1 ];
1489                 lBinLim17[0] = 0;       lBinLim17[10] = 10;      lBinLim17[20] = 24;     lBinLim17[30] = 45;    lBinLim17[40] = 95;    lBinLim17[50] = 250;
1490                 lBinLim17[1] = 1;       lBinLim17[11] = 11;      lBinLim17[21] = 25;     lBinLim17[31] = 50;    lBinLim17[41] = 100;
1491                 lBinLim17[2] = 2;       lBinLim17[12] = 13;      lBinLim17[22] = 27;     lBinLim17[32] = 55;    lBinLim17[42] = 105;
1492                 lBinLim17[3] = 3;       lBinLim17[13] = 14;      lBinLim17[23] = 30;     lBinLim17[33] = 60;    lBinLim17[43] = 110;
1493                 lBinLim17[4] = 4;       lBinLim17[14] = 15;      lBinLim17[24] = 31;     lBinLim17[34] = 65;    lBinLim17[44] = 115;
1494                 lBinLim17[5] = 5;       lBinLim17[15] = 16;      lBinLim17[25] = 32;     lBinLim17[35] = 70;    lBinLim17[45] = 120;
1495                 lBinLim17[6] = 6;       lBinLim17[16] = 20;      lBinLim17[26] = 33;     lBinLim17[36] = 75;    lBinLim17[46] = 125;
1496                 lBinLim17[7] = 7;       lBinLim17[17] = 21;      lBinLim17[27] = 34;     lBinLim17[37] = 80;    lBinLim17[47] = 130;
1497                 lBinLim17[8] = 8;       lBinLim17[18] = 22;      lBinLim17[28] = 35;     lBinLim17[38] = 85;    lBinLim17[48] = 135;
1498                 lBinLim17[9] = 9;       lBinLim17[19] = 23;      lBinLim17[29] = 40;     lBinLim17[39] = 90;    lBinLim17[49] = 140;
1499
1500                 fCFContAsCascadeCuts -> SetBinLimits(17,  lBinLim17 );            // nTrackPrimaryMultiplicity : 0 to 250
1501         delete [] lBinLim17;
1502
1503
1504         fCFContAsCascadeCuts->SetBinLimits(18,   0.0, 200.0  );     // nITSandTPCtracksAndSPDtracklets
1505   }
1506   fCFContAsCascadeCuts->SetBinLimits(19,  68.0  ,164.0  );      // BachTPCClusters
1507
1508
1509   // Regular binning definition (valid for v4-18-10-AN on)
1510   /*
1511   //setting the bin limits
1512   fCFContAsCascadeCuts->SetBinLimits(0,    0.0  ,  2.5  );      // DcaXiDaughters
1513   fCFContAsCascadeCuts->SetBinLimits(1,    0.0  ,  0.25 );      // DcaBachToPrimVertexXi
1514   fCFContAsCascadeCuts->SetBinLimits(2,    0.99 ,  1.0  );      // XiCosineOfPointingAngle
1515   fCFContAsCascadeCuts->SetBinLimits(3,    0.0  ,  4.0  );      // XiRadius
1516   fCFContAsCascadeCuts->SetBinLimits(4,    1.1  ,  1.15 );      // InvMassLambdaAsCascDghter
1517   fCFContAsCascadeCuts->SetBinLimits(5,    0.0  ,  1.0  );      // DcaV0DaughtersXi
1518   fCFContAsCascadeCuts->SetBinLimits(6,    0.98 ,  1.0  );      // V0CosineOfPointingAngleXi
1519   fCFContAsCascadeCuts->SetBinLimits(7,    0.0  , 20.0  );      // V0RadiusXi
1520   fCFContAsCascadeCuts->SetBinLimits(8,    0.0  ,  1.0  );      // DcaV0ToPrimVertexXi
1521   fCFContAsCascadeCuts->SetBinLimits(9,    0.0  ,  0.25 );      // DcaPosToPrimVertexXi
1522   fCFContAsCascadeCuts->SetBinLimits(10,   0.0  ,  0.25 );      // DcaNegToPrimVertexXi
1523   fCFContAsCascadeCuts->SetBinLimits(11,   1.25 ,  1.40 );      // InvMassXi
1524   fCFContAsCascadeCuts->SetBinLimits(12,   1.62 ,  1.74 );      // InvMassOmega
1525   fCFContAsCascadeCuts->SetBinLimits(13,   0.0  , 10.0  );      // pt_MC(Xi)
1526   fCFContAsCascadeCuts->SetBinLimits(14,  -1.1  ,  1.1  );      // Y_MC(Xi)
1527   fCFContAsCascadeCuts->SetBinLimits(15,  -1.1  ,  1.1  );      // Y_MC(Omega)
1528   fCFContAsCascadeCuts->SetBinLimits(16, -10.0  , 10.0  );      // BestPrimaryVtxPosZ
1529   if(fCollidingSystems){
1530         fCFContAsCascadeCuts->SetBinLimits(17,   0.0, 10000.0  );    // nTrackPrimaryMultiplicity
1531         fCFContAsCascadeCuts->SetBinLimits(18,   0.0, 10000.0  );    // nITSandTPCtracksAndSPDtracklets
1532   }
1533   else{
1534         fCFContAsCascadeCuts->SetBinLimits(17,   0.0, 250.0  );     // nTrackPrimaryMultiplicity
1535         fCFContAsCascadeCuts->SetBinLimits(18,   0.0, 200.0  );     // nITSandTPCtracksAndSPDtracklets
1536   }
1537   fCFContAsCascadeCuts->SetBinLimits(19,  25.0  ,165.0  );      // BachTPCClusters
1538   */
1539
1540   // Setting the number of steps : one for each cascade species (Xi-, Xi+ and Omega-, Omega+)
1541   fCFContAsCascadeCuts->SetStepTitle(0, "#Xi^{-} candidates associated to MC");
1542   fCFContAsCascadeCuts->SetStepTitle(1, "#bar{#Xi}^{+} candidates associated to MC");
1543   fCFContAsCascadeCuts->SetStepTitle(2, "#Omega^{-} candidates associated to MC");
1544   fCFContAsCascadeCuts->SetStepTitle(3, "#bar{#Omega}^{+} candidates associated to MC");
1545
1546   // Setting the variable title, per axis
1547   // fCFContAsCascadeCuts->SetVarTitle(40,  "Chi2Xi");
1548   fCFContAsCascadeCuts->SetVarTitle(0,  "Dca(XiDaughters) (cm)");
1549   fCFContAsCascadeCuts->SetVarTitle(1,  "Dca(Bach/PrimVertex) (cm)");
1550   fCFContAsCascadeCuts->SetVarTitle(2,  "cos(Xi pointing angle)");
1551   fCFContAsCascadeCuts->SetVarTitle(3,  "R_{2d}(Xi decay) (cm)");
1552   fCFContAsCascadeCuts->SetVarTitle(4,  "M_{#Lambda}(As Casc Dghter) (GeV/c^{2})");
1553    // fCFContAsCascadeCuts->SetVarTitle(40,  "V0Chi2Xi");
1554   fCFContAsCascadeCuts->SetVarTitle(5,  "Dca(V0 Daughters Xi) (cm)");
1555
1556   fCFContAsCascadeCuts->SetVarTitle(6,  "cos(V0 pointing Angle) in Casc");
1557   fCFContAsCascadeCuts->SetVarTitle(7,  "R_{2d}(V0 decay) (cm)");
1558   fCFContAsCascadeCuts->SetVarTitle(8,  "Dca(V0/PrimVertex) (cm)");
1559   fCFContAsCascadeCuts->SetVarTitle(9,  "Dca(Pos/PrimVertex) (cm)");
1560   fCFContAsCascadeCuts->SetVarTitle(10, "Dca(Neg/PrimVertex) (cm)");
1561
1562   fCFContAsCascadeCuts->SetVarTitle(11, "Inv. Mass(Xi) (GeV/c^{2})");
1563   fCFContAsCascadeCuts->SetVarTitle(12, "Inv. Mass(Omega) (GeV/c^{2})");
1564
1565   fCFContAsCascadeCuts->SetVarTitle(13, "Pt_{MC}(Casc.) (GeV/c)");
1566   //fCFContAsCascadeCuts->SetVarTitle(40, "V0toXiCosineOfPointingAngle");
1567
1568   fCFContAsCascadeCuts->SetVarTitle(14, "Y_{MC}(Xi)");
1569   fCFContAsCascadeCuts->SetVarTitle(15, "Y_{MC}(Omega)");
1570
1571   fCFContAsCascadeCuts->SetVarTitle(16, "Z-position(BestPrimVtx) (cm)");
1572
1573   fCFContAsCascadeCuts->SetVarTitle(17, "Primary Track Multiplicity");
1574   fCFContAsCascadeCuts->SetVarTitle(18, "(ITS+TPC tracks + SPD tracklets) Multiplicity");
1575   fCFContAsCascadeCuts->SetVarTitle(19, "Bach.TPC Clusters");
1576
1577   fListHistCascade->Add(fCFContAsCascadeCuts);
1578 }
1579
1580
1581 PostData(1, fListHistCascade);
1582 }// end CreateOutputObjects
1583
1584
1585
1586
1587
1588
1589 //________________________________________________________________________
1590 void AliAnalysisTaskCheckPerformanceCascade::UserExec(Option_t *)
1591 {
1592
1593   // Main loop
1594   // Called for each event
1595
1596         AliESDEvent *lESDevent = 0x0;
1597         AliAODEvent *lAODevent = 0x0;
1598         AliMCEvent  *lMCevent  = 0x0;
1599         AliStack    *lMCstack  = 0x0;
1600         Int_t ncascades = -1;
1601
1602
1603   // Connect to the InputEvent
1604   // After these lines, we should have an ESD/AOD event + the number of cascades in it.
1605
1606         if(fAnalysisType == "ESD"){
1607                 lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
1608                 if (!lESDevent) {
1609                         Printf("ERROR: lESDevent not available \n");
1610                         cout << "Name of the file with pb :" <<  CurrentFileName() << endl;  // or AliAnalysisTaskSE::CurrentFileName()
1611                         return;
1612                 }
1613         }
1614
1615         else if(fAnalysisType == "AOD"){
1616                 lAODevent = dynamic_cast<AliAODEvent*>( InputEvent() );
1617                 if (!lAODevent) {
1618                         Printf("ERROR: lAODevent not available \n");
1619                         cout << "Name of the file with pb :" <<  CurrentFileName() << endl;
1620                         return;
1621                 }
1622         }
1623
1624
1625         lMCevent = MCEvent();
1626         if (!lMCevent) {
1627                 Printf("ERROR: Could not retrieve MC event \n");
1628                 cout << "Name of the file with pb :" <<  CurrentFileName() << endl;
1629                 return;
1630         }
1631
1632         lMCstack = lMCevent->Stack();
1633         if (!lMCstack) {
1634                 Printf("ERROR: Could not retrieve MC stack \n");
1635                 cout << "Name of the file with pb :" <<  CurrentFileName() << endl;
1636                 return;
1637
1638         }
1639
1640
1641    // Temporary way : AOD awareness of the code to be developed  FIXME
1642    if(fAnalysisType == "AOD") return;
1643
1644
1645
1646
1647   //-------------------------------------------------
1648   // 0 - Trigger managment + global event selection
1649   // NOTE : Check the availability of the proper trigger
1650
1651    // Note : Presuppose the presence of AliPhysicsSelectionTask
1652
1653         UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1654         Bool_t isSelected = 0;
1655         if(     fTriggerMaskType == "kMB")           isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
1656         else if(fTriggerMaskType == "kHighMult")     isSelected = (maskIsSelected & AliVEvent::kHighMult) == AliVEvent::kHighMult;
1657         else                                         isSelected = 1; // default = select anyway (use case = run without Phys Selection task)
1658
1659         if ( ! isSelected ) {
1660                 PostData(1, fListHistCascade);
1661                 return;
1662         }
1663         //else Printf("Event selected ... \n");
1664
1665
1666   //-------------------------------------------------
1667   // 1 - Cascade vertexer (ESD)
1668         if(fkRerunV0CascVertexers){ // FIXME : relaunch V0 and Cascade vertexers
1669                 if(fAnalysisType == "ESD" ){
1670 //                         lESDevent->ResetCascades();
1671 //                         lESDevent->ResetV0s();
1672 //
1673 //                         AliV0vertexer lV0vtxer;
1674 //                         AliCascadeVertexer lCascVtxer;
1675 //
1676 //                         lV0vtxer.SetDefaultCuts(fV0Sels);
1677 //                         lCascVtxer.SetDefaultCuts(fCascSels);
1678 //
1679 //                         lV0vtxer.Tracks2V0vertices(lESDevent);
1680 //                         lCascVtxer.V0sTracks2CascadeVertices(lESDevent);
1681                 }
1682         }
1683
1684
1685   //------------------------------------------------
1686   // 2 - Preparing the general info about of the event = prim. Vtx + magnetic field (ESD)
1687
1688
1689 //      if(fAnalysisType == "ESD" ){
1690
1691         // Magnetic field
1692                 const Double_t lMagneticField = lESDevent->GetMagneticField( );
1693
1694         // Prim vertex
1695                 const AliESDVertex *lPrimaryTrackingVtx = lESDevent->GetPrimaryVertexTracks();  // get the vtx stored in ESD found with tracks
1696                 const AliESDVertex *lPrimarySPDVtx      = lESDevent->GetPrimaryVertexSPD();     // get the vtx stored in ESD found with SPD tracklets
1697
1698                 const AliESDVertex *lPrimaryBestVtx     = lESDevent->GetPrimaryVertex();
1699                         // get the best primary vertex available for the event
1700                         // As done in AliCascadeVertexer, we keep the one which is the best one available.
1701                         // between : Tracking vertex > SPD vertex > TPC vertex > default SPD vertex
1702                         Double_t lBestPrimaryVtxPos[3]   = {-100.0, -100.0, -100.0};
1703                 lPrimaryBestVtx->GetXYZ( lBestPrimaryVtxPos );
1704
1705                 // FIXME : quality cut on the z-position of the prim vertex.
1706                 if(fkQualityCutZprimVtxPos) {
1707                         if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0 ) {
1708                                 AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
1709                                 PostData(1, fListHistCascade);
1710                                 return;
1711                         }
1712                 }
1713                 // FIXME : quality selection regarding pile-up rejection
1714                 if(fkRejectEventPileUp) {
1715                         if(lESDevent->IsPileupFromSPDInMultBins() ){// 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
1716                                 AliWarning("Pb / Event tagged as pile-up by SPD... return !");
1717                                 PostData(1, fListHistCascade);
1718                                 return;
1719                         }
1720                 }
1721                 // FIXME : remove TPC-only primary vertex : retain only events with tracking + SPD vertex
1722                 if(fkQualityCutNoTPConlyPrimVtx) {
1723                         if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingVtx->GetStatus() ){
1724                                 AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
1725                                 PostData(1, fListHistCascade);
1726                                 return;
1727                         }
1728                 }
1729 //      }// if ESD
1730
1731
1732   //    cout << "Name of the accessed file :" <<  fInputHandler->GetTree()->GetCurrentFile()->GetName() << endl;
1733
1734   //    cout << "Tree characteristics ..." << endl;
1735   //    fInputHandler->GetTree()->Print("toponly");
1736   //    fInputHandler->GetTree()->GetBranch("PrimaryVertex")->Print();
1737   //    fInputHandler->GetTree()->GetBranch("SPDVertex")->Print();
1738
1739
1740
1741   // ---------------------------------------------------------------
1742   // - Initialisation of the part dedicated to cascade vertices
1743
1744   if(fAnalysisType == "ESD")            ncascades = lESDevent->GetNumberOfCascades();
1745   else if(fAnalysisType == "AOD")       ncascades = lAODevent->GetNumberOfCascades();
1746
1747
1748   Int_t   nNumberOfMCPrimaries       = -1;
1749   Int_t   nMCPrimariesInEtaBelow0p8  =  0;
1750   Int_t   nMCPrimariesInEtaBelow1p0  =  0;
1751
1752   Int_t   nTrackPrimaryMultiplicity        = -1;
1753   Int_t   nSPDTracklets                    =  0; // AliESDEvent::EstimateMultiplicity will re-initialise the value to 0
1754   Int_t   nITSandTPCtracksAndSPDtracklets  =  0; // AliESDEvent::EstimateMultiplicity will re-initialise the value to 0
1755   Int_t   nTracksITSSApure                 =  0; // AliESDEvent::EstimateMultiplicity will re-initialise the value to 0
1756
1757         nNumberOfMCPrimaries       = lMCstack->GetNprimary();
1758         if(nNumberOfMCPrimaries < 1) return;
1759
1760         nTrackPrimaryMultiplicity  = fESDtrackCuts->CountAcceptedTracks(lESDevent);
1761         //EstimateMultiplicity(Int_t &tracklets, Int_t &trITSTPC, Int_t &trITSSApure, Double_t eta, Bool_t useDCAFlag,Bool_t useV0Flag)
1762         lESDevent->EstimateMultiplicity( nSPDTracklets, nITSandTPCtracksAndSPDtracklets, nTracksITSSApure, 1.0, kTRUE, kTRUE);
1763
1764
1765         fHistMCTrackMultiplicity->Fill( nNumberOfMCPrimaries );
1766
1767 //_____________________________________________________________________________
1768 // Part 1 - Loop over the MC primaries
1769
1770     for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < nNumberOfMCPrimaries; iCurrentLabelStack++)
1771     {// This is the begining of the loop on primaries
1772
1773         TParticle* lCurrentParticle = 0x0;
1774                    lCurrentParticle = lMCstack->Particle( iCurrentLabelStack );
1775         if(!lCurrentParticle){
1776                 Printf("MC Primary loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
1777                 continue;
1778         }
1779
1780         Double_t lEtaCurrentParticle = TMath::Abs( lCurrentParticle->Eta() );
1781         if( lEtaCurrentParticle < 1.0 ){
1782                 nMCPrimariesInEtaBelow1p0++;
1783                 if( lEtaCurrentParticle < 0.8 ) nMCPrimariesInEtaBelow0p8++;
1784         }
1785     }
1786
1787     f2dHistRecoPrimTrckMultVsMCMult->Fill( nTrackPrimaryMultiplicity,       nMCPrimariesInEtaBelow0p8 );
1788     f2dHistRecoEstimateMultVsMCMult->Fill( nITSandTPCtracksAndSPDtracklets, nMCPrimariesInEtaBelow1p0 );
1789
1790
1791    // For proton
1792    /*
1793    for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < nNumberOfMCPrimaries; iCurrentLabelStack++)
1794         {// This is the begining of the loop on primaries, for protons
1795
1796         TParticle* lCurrentParticle = lMCstack->Particle( iCurrentLabelStack );
1797         if(!lCurrentParticle){
1798                 Printf("Proton loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
1799                 continue;
1800
1801         }
1802
1803         if( lCurrentParticle->GetPdgCode() == 2212 )
1804                 fHistEtaGenProton->Fill( lCurrentParticle->Eta() );
1805
1806         if( lCurrentParticle->GetPdgCode() == -2212 )
1807                 fHistEtaGenAntiProton->Fill( lCurrentParticle->Eta() );
1808         }// end loop over primary proton
1809    */
1810
1811
1812
1813 //_____________________________________________________________________________
1814 // Part 2 - Loop over the different types of GENERATED cascades (Xi-+, Omega-+)
1815
1816         // - Initialisation of useful local variables
1817
1818         Int_t lPdgCodeCasc            = 0;
1819         Int_t lPdgCodeBach            = 0;
1820         Int_t lPdgCodeLambda          = 0;
1821         Int_t lPdgCodeDghtMesV0       = 0;
1822         Int_t lPdgCodeDghtBarV0       = 0;
1823
1824
1825         TH1F *lHistEtaGenCasc         = 0;
1826         TH2F *l2dHistGenPtVsGenYGen   = 0;
1827
1828         TH1F *lHistThetaGenCasc       = 0;
1829         TH2F *l2dHistGenPtVsGenYFdbl  = 0;
1830         TH1F *lHistThetaLambda        = 0;
1831         TH1F *lHistThetaBach          = 0;
1832         TH1F *lHistThetaBarDghter     = 0;
1833         TH1F *lHistThetaMesDghter     = 0;
1834         TH1F *lHistPtBach             = 0;
1835         TH1F *lHistPtBarDghter        = 0;
1836         TH1F *lHistPtMesDghter        = 0;
1837
1838
1839 for(Int_t iCascType = 1; iCascType < 5; iCascType++)
1840 {
1841
1842 switch (iCascType)
1843   {
1844     case 1: // Xi-
1845          lPdgCodeCasc       =   3312;  //Xi-
1846          lPdgCodeBach       =   -211;  //Pi-
1847          lPdgCodeLambda     =   3122;  //Lambda0
1848          lPdgCodeDghtMesV0  =   -211;  //Pi-
1849          lPdgCodeDghtBarV0  =   2212;  //Proton
1850
1851                 // any Xi-
1852          lHistEtaGenCasc        = fHistEtaGenCascXiMinus;
1853          l2dHistGenPtVsGenYGen  = f2dHistGenPtVsGenYGenXiMinus;
1854
1855                 // cascades generated within acceptance (cut in pt + theta)
1856          lHistThetaGenCasc      = fHistThetaGenCascXiMinus;
1857          l2dHistGenPtVsGenYFdbl = f2dHistGenPtVsGenYFdblXiMinus;
1858          lHistThetaLambda       = fHistThetaLambdaXiMinus;
1859          lHistThetaBach         = fHistThetaBachXiMinus;
1860          lHistThetaBarDghter    = fHistThetaBarDghterXiMinus;
1861          lHistThetaMesDghter    = fHistThetaMesDghterXiMinus;
1862          lHistPtBach            = fHistPtBachXiMinus;
1863          lHistPtBarDghter       = fHistPtBarDghterXiMinus;
1864          lHistPtMesDghter       = fHistPtMesDghterXiMinus;
1865         break;
1866
1867     case 2: // Xi+
1868          lPdgCodeCasc        =  -3312;  //Xi+
1869          lPdgCodeBach        =    211;  //Pi+
1870          lPdgCodeLambda      =  -3122;  //AntiLambda0
1871          lPdgCodeDghtMesV0   =    211;  //Pi+
1872          lPdgCodeDghtBarV0   =  -2212;  //AntiProton
1873
1874                 // any Xi+
1875          lHistEtaGenCasc        = fHistEtaGenCascXiPlus;
1876          l2dHistGenPtVsGenYGen  = f2dHistGenPtVsGenYGenXiPlus;
1877
1878                 // cascades generated within acceptance (cut in pt + theta)
1879          lHistThetaGenCasc      = fHistThetaGenCascXiPlus;
1880          l2dHistGenPtVsGenYFdbl = f2dHistGenPtVsGenYFdblXiPlus;
1881          lHistThetaLambda       = fHistThetaLambdaXiPlus;
1882          lHistThetaBach         = fHistThetaBachXiPlus;
1883          lHistThetaBarDghter    = fHistThetaBarDghterXiPlus;
1884          lHistThetaMesDghter    = fHistThetaMesDghterXiPlus;
1885          lHistPtBach            = fHistPtBachXiPlus;
1886          lHistPtBarDghter       = fHistPtBarDghterXiPlus;
1887          lHistPtMesDghter       = fHistPtMesDghterXiPlus;
1888         break;
1889
1890     case 3: // Omega-
1891          lPdgCodeCasc       =   3334;  //Omega-
1892          lPdgCodeBach       =   -321;  //K-
1893          lPdgCodeLambda     =   3122;  //Lambda0
1894          lPdgCodeDghtMesV0  =   -211;  //Pi-
1895          lPdgCodeDghtBarV0  =   2212;  //Proton
1896
1897                 // any Omega-
1898          lHistEtaGenCasc        = fHistEtaGenCascOmegaMinus;
1899          l2dHistGenPtVsGenYGen  = f2dHistGenPtVsGenYGenOmegaMinus;
1900
1901                 // cascades generated within acceptance (cut in pt + theta)
1902          lHistThetaGenCasc      = fHistThetaGenCascOmegaMinus;
1903          l2dHistGenPtVsGenYFdbl = f2dHistGenPtVsGenYFdblOmegaMinus;
1904          lHistThetaLambda       = fHistThetaLambdaOmegaMinus;
1905          lHistThetaBach         = fHistThetaBachOmegaMinus;
1906          lHistThetaBarDghter    = fHistThetaBarDghterOmegaMinus;
1907          lHistThetaMesDghter    = fHistThetaMesDghterOmegaMinus;
1908          lHistPtBach            = fHistPtBachOmegaMinus;
1909          lHistPtBarDghter       = fHistPtBarDghterOmegaMinus;
1910          lHistPtMesDghter       = fHistPtMesDghterOmegaMinus;
1911         break;
1912
1913     case 4:  // Omega+
1914          lPdgCodeCasc       =  -3334;  //Omega+
1915          lPdgCodeBach       =    321;  //K+
1916          lPdgCodeLambda     =  -3122;  //AntiLambda0
1917          lPdgCodeDghtMesV0  =    211;  //Pi+
1918          lPdgCodeDghtBarV0  =  -2212;  //AntiProton
1919
1920                 // any Omega+
1921          lHistEtaGenCasc        = fHistEtaGenCascOmegaPlus;
1922          l2dHistGenPtVsGenYGen  = f2dHistGenPtVsGenYGenOmegaPlus;
1923
1924                 // cascades generated within acceptance (cut in pt + theta)
1925          lHistThetaGenCasc      = fHistThetaGenCascOmegaPlus;
1926          l2dHistGenPtVsGenYFdbl = f2dHistGenPtVsGenYFdblOmegaPlus;
1927          lHistThetaLambda       = fHistThetaLambdaOmegaPlus;
1928          lHistThetaBach         = fHistThetaBachOmegaPlus;
1929          lHistThetaBarDghter    = fHistThetaBarDghterOmegaPlus;
1930          lHistThetaMesDghter    = fHistThetaMesDghterOmegaPlus;
1931          lHistPtBach            = fHistPtBachOmegaPlus;
1932          lHistPtBarDghter       = fHistPtBarDghterOmegaPlus;
1933          lHistPtMesDghter       = fHistPtMesDghterOmegaPlus;
1934         break;
1935
1936   }// end switch cascade
1937
1938
1939    for (Int_t iCurrentLabelStack = 0; iCurrentLabelStack < nNumberOfMCPrimaries; iCurrentLabelStack++)
1940     {// This is the begining of the loop on primaries
1941
1942         TParticle* lCurrentParticle = 0x0;
1943                    lCurrentParticle = lMCstack->Particle( iCurrentLabelStack );
1944         if(!lCurrentParticle){
1945                 Printf("Cascade loop %d - MC TParticle pointer to current stack particle = 0x0 ! Skip ...\n", iCurrentLabelStack );
1946                 continue;
1947
1948         }
1949
1950         if( lCurrentParticle->GetPdgCode() == lPdgCodeCasc ){  // Here !
1951                 //cout << "Xi- within loop " << iCurrentLabelStack << "/ " << nNumberOfMCPrimaries << endl;
1952
1953                 // -  Xi level ... _____________________________________________________________
1954                 TParticle* xiMC = 0x0;
1955                            xiMC = lCurrentParticle;
1956                 if(!xiMC){
1957                         Printf("MC TParticle pointer to Cascade = 0x0 ! Skip ...");
1958                         continue;
1959
1960                 }
1961
1962                 // Fill the first histos : = any generated Xi, not necessarily within the acceptance
1963                 Double_t lRapXiMC = 0.5*TMath::Log((xiMC->Energy() + xiMC->Pz()) / (xiMC->Energy() - xiMC->Pz() +1.e-13));
1964
1965                 lHistEtaGenCasc         ->Fill( xiMC->Eta() );
1966                 l2dHistGenPtVsGenYGen   ->Fill( xiMC->Pt(), lRapXiMC  );
1967
1968
1969
1970                 // Check the emission of particle stays within the acceptance of the detector (cut in theta)
1971                 if( xiMC->Theta() < TMath::Pi()/4.0  ||    xiMC->Theta() > 3.0*TMath::Pi()/4.0 ) continue;
1972                 if( xiMC->GetNDaughters() != 2) continue;
1973                 if( xiMC->GetDaughter(0) < 0 )  continue;
1974                 if( xiMC->GetDaughter(1) < 0 )  continue;
1975
1976                         TParticle* lDght0ofXi = lMCstack->Particle(  xiMC->GetDaughter(0) );
1977                         TParticle* lDght1ofXi = lMCstack->Particle(  xiMC->GetDaughter(1) );
1978
1979                 TParticle* lLambda = 0;
1980                 TParticle* lBach   = 0;
1981
1982                 // Xi - Case 1
1983                         if(     lDght0ofXi->GetPdgCode() == lPdgCodeLambda   &&  // Here !
1984                                 lDght1ofXi->GetPdgCode() == lPdgCodeBach ){      // Here !
1985
1986                                 lLambda = lDght0ofXi;
1987                                 lBach   = lDght1ofXi;
1988                         }// end if dghter 0 = Lambda and    dghter 1 = Pi-
1989
1990                 // Xi - Case 2
1991                         else if( lDght0ofXi->GetPdgCode() == lPdgCodeBach  &&      // Here !
1992                                  lDght1ofXi->GetPdgCode() == lPdgCodeLambda ){     // Here !
1993
1994                                 lBach   = lDght0ofXi;
1995                                 lLambda = lDght1ofXi;
1996                         }//  end if dghter 0 = Pi-  and   dghter 1 = Lambda
1997
1998                 // V0 otherwise - Case 3
1999                         else continue;
2000
2001                         // Check the emission of particle stays within the acceptance of the detector (cut in pt + theta)
2002                         if( lLambda->Theta() < TMath::Pi()/4.0  ||    lLambda->Theta() > 3.0*TMath::Pi()/4.0 ) continue;
2003                         if( lBach->Theta() < TMath::Pi()/4.0    ||    lBach->Theta() > 3.0*TMath::Pi()/4.0 )   continue;
2004
2005                         if( lBach->Pt() < 0.150 ) continue; //FIXME : maybe tuned for Xi but not for K- from Omega ...
2006
2007
2008
2009                 // -  V0 level ... _____________________________________________________________
2010                 TParticle* lDghtBarV0 = 0;
2011                 TParticle* lDghtMesV0 = 0;
2012
2013                 if( lLambda->GetNDaughters() != 2 )  continue;
2014                 if( lLambda->GetDaughter(0) < 0 )    continue;
2015                 if( lLambda->GetDaughter(1) < 0 )    continue;
2016
2017
2018                 TParticle* lDght0ofLambda = lMCstack->Particle(  lLambda->GetDaughter(0) );
2019                 TParticle* lDght1ofLambda = lMCstack->Particle(  lLambda->GetDaughter(1) );
2020
2021                 // V0 - Case 1
2022                         if(     lDght0ofLambda->GetPdgCode() == lPdgCodeDghtBarV0 &&    // Here !
2023                                 lDght1ofLambda->GetPdgCode() == lPdgCodeDghtMesV0 ){    // Here !
2024
2025                                 lDghtBarV0 = lDght0ofLambda;
2026                                 lDghtMesV0 = lDght1ofLambda;
2027                         }// end if dghter 0 = Proton  and   dghter 1 = Pi-
2028
2029                 // V0 - Case 2
2030                         else if( lDght0ofLambda->GetPdgCode() == lPdgCodeDghtMesV0  &&     // Here !
2031                                  lDght1ofLambda->GetPdgCode() == lPdgCodeDghtBarV0 ){      // Here !
2032
2033                                 lDghtMesV0 = lDght0ofLambda;
2034                                 lDghtBarV0 = lDght1ofLambda;
2035                         }//  end if dghter 0 = Pi-  and   dghter 1 = proton
2036
2037                 // V0 otherwise - Case 3
2038                         else continue;
2039
2040
2041                         // Check the emission of particle stays within the acceptance of the detector
2042                         if( lDghtBarV0->Theta() < TMath::Pi()/4.0  ||  lDghtBarV0->Theta() > 3.0*TMath::Pi()/4.0 ) continue;
2043                         if( lDghtMesV0->Theta() < TMath::Pi()/4.0  ||  lDghtMesV0->Theta() > 3.0*TMath::Pi()/4.0 ) continue;
2044
2045                         if( lDghtBarV0->Pt() < 0.250 ) continue;
2046                         if( lDghtMesV0->Pt() < 0.150 ) continue;
2047
2048
2049
2050                 // - Just to know which file is currently open : locate the file containing Xi
2051                 //cout << "Name of the file containing generated Xi :" <<  fInputHandler->GetTree()->GetCurrentFile()->GetName()
2052                 //                                                   <<  endl;
2053
2054                 Double_t lRadToDeg = 180.0/TMath::Pi();
2055
2056                 // - Filling histos ... _________________________________________________________________
2057                         lHistThetaGenCasc       ->Fill( lRadToDeg * xiMC->Theta()  );
2058                         l2dHistGenPtVsGenYFdbl  ->Fill( xiMC->Pt(), lRapXiMC );
2059
2060                         // - Fill theta histos for Lambda and Bach
2061                         lHistThetaLambda        ->Fill( lRadToDeg * lLambda->Theta() );
2062                         lHistThetaBach          ->Fill( lRadToDeg *   lBach->Theta() );
2063
2064                         // - Fill theta histos for V0 daughters
2065                         lHistThetaBarDghter     ->Fill( lRadToDeg * lDghtBarV0->Theta() );
2066                         lHistThetaMesDghter     ->Fill( lRadToDeg * lDghtMesV0->Theta() );
2067
2068                         // - Fill pt histos.
2069                         lHistPtBach             ->Fill(      lBach->Pt() );
2070                         lHistPtBarDghter        ->Fill( lDghtBarV0->Pt() );
2071                         lHistPtMesDghter        ->Fill( lDghtMesV0->Pt() );
2072
2073                 }// end if current particle = Xi-
2074
2075      }// This is the end of the loop on primaries
2076
2077 // - Re-initialisation of the local TH1F pointers
2078 lHistEtaGenCasc         = 0x0;
2079 l2dHistGenPtVsGenYGen   = 0x0;
2080
2081 lHistThetaGenCasc       = 0x0;
2082 l2dHistGenPtVsGenYFdbl  = 0x0;
2083 lHistThetaLambda        = 0x0;
2084 lHistThetaBach          = 0x0;
2085 lHistThetaBarDghter     = 0x0;
2086 lHistThetaMesDghter     = 0x0;
2087 lHistPtBach             = 0x0;
2088 lHistPtBarDghter        = 0x0;
2089 lHistPtMesDghter        = 0x0;
2090
2091 } // end of loop over the different types of cascades (Xi-+, Omega-+)
2092
2093
2094
2095 //__________________________________________________________________________
2096 // Part 3 - Loop over the reconstructed candidates
2097
2098
2099 for (Int_t iXi = 0; iXi < ncascades; iXi++)
2100 {// This is the begining of the Cascade loop
2101
2102         AliESDcascade *xiESD = lESDevent->GetCascade(iXi);
2103         if (!xiESD) continue;
2104
2105         // - Step II.1 : Connection to daughter tracks of the current cascade
2106         //-------------
2107
2108                 UInt_t lIdxPosXi        = (UInt_t) TMath::Abs( xiESD->GetPindex() );
2109                 UInt_t lIdxNegXi        = (UInt_t) TMath::Abs( xiESD->GetNindex() );
2110                 UInt_t lBachIdx         = (UInt_t) TMath::Abs( xiESD->GetBindex() );
2111                 // abs value not needed ; the index should always be positive (!= label ...)
2112
2113
2114         // FIXME : rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
2115         if(lBachIdx == lIdxNegXi) {
2116                 AliWarning("Pb / Idx(Bach. track) = Idx(Neg. track) ... continue!"); continue;
2117         }
2118         if(lBachIdx == lIdxPosXi) {
2119                 AliWarning("Pb / Idx(Bach. track) = Idx(Pos. track) ... continue!"); continue;
2120         }
2121
2122         AliESDtrack *pTrackXi           = lESDevent->GetTrack( lIdxPosXi );
2123         AliESDtrack *nTrackXi           = lESDevent->GetTrack( lIdxNegXi );
2124         AliESDtrack *bachTrackXi        = lESDevent->GetTrack( lBachIdx  );
2125         if (!pTrackXi || !nTrackXi || !bachTrackXi ) {
2126                 Printf("ERROR: Could not retrieve one of the 3 daughter tracks of the cascade ...");
2127                 continue;
2128         }
2129
2130         Int_t lPosTPCClusters   = pTrackXi->GetTPCNcls();
2131         Int_t lNegTPCClusters   = nTrackXi->GetTPCNcls();
2132         Int_t lBachTPCClusters  = bachTrackXi->GetTPCNcls();
2133
2134                 // FIXME : rejection of a poor quality tracks
2135         if(fkQualityCutTPCrefit){
2136                 // 1 - Poor quality related to TPCrefit
2137                 ULong_t pStatus    = pTrackXi->GetStatus();
2138                 ULong_t nStatus    = nTrackXi->GetStatus();
2139                 ULong_t bachStatus = bachTrackXi->GetStatus();
2140                 if ((pStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning("Pb / V0 Pos. track has no TPCrefit ... continue!"); continue; }
2141                 if ((nStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning("Pb / V0 Neg. track has no TPCrefit ... continue!"); continue; }
2142                 if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / Bach.   track has no TPCrefit ... continue!"); continue; }
2143         }
2144         if(fkQualityCut80TPCcls){
2145                 // 2 - Poor quality related to TPC clusters
2146                 if(lPosTPCClusters  < 80) { AliWarning("Pb / V0 Pos. track has less than 80 TPC clusters ... continue!"); continue; }
2147                 if(lNegTPCClusters  < 80) { AliWarning("Pb / V0 Neg. track has less than 80 TPC clusters ... continue!"); continue; }
2148                 if(lBachTPCClusters < 80) { AliWarning("Pb / Bach.   track has less than 80 TPC clusters ... continue!"); continue; }
2149         }
2150
2151         // - Step II.2 : Info over reconstructed cascades
2152         //-------------
2153
2154         Double_t lInvMassXiMinus    = 0.;
2155         Double_t lInvMassXiPlus     = 0.;
2156         Double_t lInvMassOmegaMinus = 0.;
2157         Double_t lInvMassOmegaPlus  = 0.;
2158
2159         Double_t lV0quality = 0.;
2160
2161         if( bachTrackXi->Charge() < 0 ) {
2162                 lV0quality = 0.;
2163                 xiESD->ChangeMassHypothesis(lV0quality , 3312);
2164                         // Calculate the effective mass of the Xi- candidate.
2165                         // pdg code 3312 = Xi-
2166                 lInvMassXiMinus = xiESD->GetEffMassXi();
2167
2168                 lV0quality = 0.;
2169                 xiESD->ChangeMassHypothesis(lV0quality , 3334);
2170                         // Calculate the effective mass of the Xi- candidate.
2171                         // pdg code 3334 = Omega-
2172                 lInvMassOmegaMinus = xiESD->GetEffMassXi();
2173
2174                 lV0quality = 0.;
2175                 xiESD->ChangeMassHypothesis(lV0quality , 3312);         // Back to default hyp.
2176
2177         }
2178
2179         if( bachTrackXi->Charge() >  0 ){
2180                 lV0quality = 0.;
2181                 xiESD->ChangeMassHypothesis(lV0quality , -3312);
2182                         // Calculate the effective mass of the Xi+ candidate.
2183                         // pdg code -3312 = Xi+
2184                 lInvMassXiPlus = xiESD->GetEffMassXi();
2185
2186                 lV0quality = 0.;
2187                 xiESD->ChangeMassHypothesis(lV0quality , -3334);
2188                         // Calculate the effective mass of the Xi+ candidate.
2189                         // pdg code -3334  = Omega+
2190                 lInvMassOmegaPlus = xiESD->GetEffMassXi();
2191
2192                 lV0quality = 0.;
2193                 xiESD->ChangeMassHypothesis(lV0quality , -3312);        // Back to "default" hyp.
2194         }
2195
2196         Double_t lChargeXi = xiESD->Charge();
2197
2198         if( lChargeXi < 0 )     fHistMassXiMinus        ->Fill( lInvMassXiMinus );
2199         if( lChargeXi > 0 )     fHistMassXiPlus         ->Fill( lInvMassXiPlus );
2200         if( lChargeXi < 0 )     fHistMassOmegaMinus     ->Fill( lInvMassOmegaMinus );
2201         if( lChargeXi > 0 )     fHistMassOmegaPlus      ->Fill( lInvMassOmegaPlus );
2202
2203
2204         // - Step II.3 : PID info
2205         //-------------
2206
2207
2208         // 3.1 - PID Information
2209
2210         Bool_t   lIsPosInXiProton      = kFALSE;
2211         Bool_t   lIsPosInXiPion        = kFALSE;
2212         Bool_t   lIsPosInOmegaProton   = kFALSE;
2213         Bool_t   lIsPosInOmegaPion     = kFALSE;
2214
2215         Bool_t   lIsNegInXiProton      = kFALSE;
2216         Bool_t   lIsNegInXiPion        = kFALSE;
2217         Bool_t   lIsNegInOmegaProton   = kFALSE;
2218         Bool_t   lIsNegInOmegaPion     = kFALSE;
2219
2220         Bool_t   lIsBachelorKaon       = kFALSE;
2221         Bool_t   lIsBachelorPion       = kFALSE;
2222
2223         Bool_t   lIsBachelorKaonForTPC = kFALSE; // For ESD only ...//FIXME : wait for availability in AOD
2224         Bool_t   lIsBachelorPionForTPC = kFALSE; // For ESD only ...
2225         Bool_t   lIsNegPionForTPC      = kFALSE; // For ESD only ...
2226         Bool_t   lIsPosPionForTPC      = kFALSE; // For ESD only ...
2227         Bool_t   lIsNegProtonForTPC    = kFALSE; // For ESD only ...
2228         Bool_t   lIsPosProtonForTPC    = kFALSE; // For ESD only ...
2229
2230         // 3.1.A - Combined PID
2231         // Reasonable guess for the priors for the cascade track sample (e-, mu, pi, K, p)
2232         Double_t lPriorsGuessXi[10]    = {0, 0, 2, 0, 1,  0,0,0,0,0};
2233         Double_t lPriorsGuessOmega[10] = {0, 0, 1, 1, 1,  0,0,0,0,0};
2234
2235         // Combined VO-positive-daughter PID
2236         AliPID pPidXi;         pPidXi.SetPriors(    lPriorsGuessXi    , kTRUE); // kTRUE = for charged particle PID
2237         AliPID pPidOmega;      pPidOmega.SetPriors( lPriorsGuessOmega , kTRUE); // kTRUE = for charged particle PID
2238
2239         if( pTrackXi->IsOn(AliESDtrack::kESDpid) ){  // Combined PID exists
2240                 Double_t r[10] = {0.}; pTrackXi->GetESDpid(r);
2241                 pPidXi.SetProbabilities(r);
2242                 pPidOmega.SetProbabilities(r);
2243
2244                 // Check if the V0 positive track is a proton (case for Xi-)
2245                 Double_t pproton = pPidXi.GetProbability(AliPID::kProton);
2246                 if (pproton > pPidXi.GetProbability(AliPID::kElectron) &&
2247                     pproton > pPidXi.GetProbability(AliPID::kMuon)     &&
2248                     pproton > pPidXi.GetProbability(AliPID::kPion)     &&
2249                     pproton > pPidXi.GetProbability(AliPID::kKaon)     )     lIsPosInXiProton = kTRUE;
2250
2251                 // Check if the V0 positive track is a pi+ (case for Xi+)
2252                 Double_t ppion = pPidXi.GetProbability(AliPID::kPion);
2253                 if (ppion > pPidXi.GetProbability(AliPID::kElectron) &&
2254                     ppion > pPidXi.GetProbability(AliPID::kMuon)     &&
2255                     ppion > pPidXi.GetProbability(AliPID::kKaon)     &&
2256                     ppion > pPidXi.GetProbability(AliPID::kProton)   )     lIsPosInXiPion = kTRUE;
2257
2258
2259                 // Check if the V0 positive track is a proton (case for Omega-)
2260                 pproton = 0.;
2261                     pproton = pPidOmega.GetProbability(AliPID::kProton);
2262                 if (pproton > pPidOmega.GetProbability(AliPID::kElectron) &&
2263                     pproton > pPidOmega.GetProbability(AliPID::kMuon)     &&
2264                     pproton > pPidOmega.GetProbability(AliPID::kPion)     &&
2265                     pproton > pPidOmega.GetProbability(AliPID::kKaon)     )  lIsPosInOmegaProton = kTRUE;
2266
2267                 // Check if the V0 positive track is a pi+ (case for Omega+)
2268                 ppion = 0.;
2269                     ppion = pPidOmega.GetProbability(AliPID::kPion);
2270                 if (ppion > pPidOmega.GetProbability(AliPID::kElectron) &&
2271                     ppion > pPidOmega.GetProbability(AliPID::kMuon)     &&
2272                     ppion > pPidOmega.GetProbability(AliPID::kKaon)     &&
2273                     ppion > pPidOmega.GetProbability(AliPID::kProton)   )    lIsPosInOmegaPion = kTRUE;
2274
2275         }// end if V0 positive track with existing combined PID
2276
2277
2278         // Combined VO-negative-daughter PID
2279         AliPID nPidXi;         nPidXi.SetPriors(    lPriorsGuessXi    , kTRUE); // kTRUE = for charged particle PID
2280         AliPID nPidOmega;      nPidOmega.SetPriors( lPriorsGuessOmega , kTRUE); // kTRUE = for charged particle PID
2281
2282         if( nTrackXi->IsOn(AliESDtrack::kESDpid) ){  // Combined PID exists
2283                 Double_t r[10] = {0.}; nTrackXi->GetESDpid(r);
2284                 nPidXi.SetProbabilities(r);
2285                 nPidOmega.SetProbabilities(r);
2286
2287                 // Check if the V0 negative track is a pi- (case for Xi-)
2288                 Double_t ppion = nPidXi.GetProbability(AliPID::kPion);
2289                 if (ppion > nPidXi.GetProbability(AliPID::kElectron) &&
2290                     ppion > nPidXi.GetProbability(AliPID::kMuon)     &&
2291                     ppion > nPidXi.GetProbability(AliPID::kKaon)     &&
2292                     ppion > nPidXi.GetProbability(AliPID::kProton)   )     lIsNegInXiPion = kTRUE;
2293
2294                 // Check if the V0 negative track is an anti-proton (case for Xi+)
2295                 Double_t pproton = nPidXi.GetProbability(AliPID::kProton);
2296                 if (pproton > nPidXi.GetProbability(AliPID::kElectron) &&
2297                     pproton > nPidXi.GetProbability(AliPID::kMuon)     &&
2298                     pproton > nPidXi.GetProbability(AliPID::kPion)     &&
2299                     pproton > nPidXi.GetProbability(AliPID::kKaon)     )     lIsNegInXiProton = kTRUE;
2300
2301                 // Check if the V0 negative track is a pi- (case for Omega-)
2302                 ppion = 0.;
2303                     ppion = nPidOmega.GetProbability(AliPID::kPion);
2304                 if (ppion > nPidOmega.GetProbability(AliPID::kElectron) &&
2305                     ppion > nPidOmega.GetProbability(AliPID::kMuon)     &&
2306                     ppion > nPidOmega.GetProbability(AliPID::kKaon)     &&
2307                     ppion > nPidOmega.GetProbability(AliPID::kProton)   )    lIsNegInOmegaPion = kTRUE;
2308
2309                 // Check if the V0 negative track is an anti-proton (case for Omega+)
2310                 pproton = 0.;
2311                     pproton = nPidOmega.GetProbability(AliPID::kProton);
2312                 if (pproton > nPidOmega.GetProbability(AliPID::kElectron) &&
2313                     pproton > nPidOmega.GetProbability(AliPID::kMuon)     &&
2314                     pproton > nPidOmega.GetProbability(AliPID::kPion)     &&
2315                     pproton > nPidOmega.GetProbability(AliPID::kKaon)     )  lIsNegInOmegaProton = kTRUE;
2316
2317         }// end if V0 negative track with existing combined PID
2318
2319
2320         // Combined bachelor PID
2321         AliPID bachPidXi;      bachPidXi.SetPriors(    lPriorsGuessXi    , kTRUE); // kTRUE = for charged particle PID
2322         AliPID bachPidOmega;   bachPidOmega.SetPriors( lPriorsGuessOmega , kTRUE); // kTRUE = for charged particle PID
2323
2324         Double_t ppionBach = 0.0, pkaonBach = 0.0;
2325
2326         if( bachTrackXi->IsOn(AliESDtrack::kESDpid) ){  // Combined PID exists
2327                 Double_t r[10] = {0.}; bachTrackXi->GetESDpid(r);
2328                 bachPidXi.SetProbabilities(r);
2329                 bachPidOmega.SetProbabilities(r);
2330                 // Check if the bachelor track is a pion
2331                     ppionBach = bachPidXi.GetProbability(AliPID::kPion);
2332                 if (ppionBach > bachPidXi.GetProbability(AliPID::kElectron) &&
2333                     ppionBach > bachPidXi.GetProbability(AliPID::kMuon)     &&
2334                     ppionBach > bachPidXi.GetProbability(AliPID::kKaon)     &&
2335                     ppionBach > bachPidXi.GetProbability(AliPID::kProton)   )     lIsBachelorPion = kTRUE;
2336                 // Check if the bachelor track is a kaon
2337                     pkaonBach = bachPidOmega.GetProbability(AliPID::kKaon);
2338                 if (pkaonBach > bachPidOmega.GetProbability(AliPID::kElectron) &&
2339                     pkaonBach > bachPidOmega.GetProbability(AliPID::kMuon)     &&
2340                     pkaonBach > bachPidOmega.GetProbability(AliPID::kPion)     &&
2341                     pkaonBach > bachPidOmega.GetProbability(AliPID::kProton)   )  lIsBachelorKaon = kTRUE;
2342         }// end if bachelor track with existing combined PID
2343
2344
2345         // 3.1.B - TPC PID : 4-sigma bands on Bethe-Bloch curve
2346
2347         // Bachelor
2348         if (TMath::Abs(fESDpid->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 4) lIsBachelorKaonForTPC = kTRUE;
2349         if (TMath::Abs(fESDpid->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 4) lIsBachelorPionForTPC = kTRUE;
2350
2351         // Negative V0 daughter
2352         if (TMath::Abs(fESDpid->NumberOfSigmasTPC( nTrackXi,AliPID::kPion   )) < 4) lIsNegPionForTPC   = kTRUE;
2353         if (TMath::Abs(fESDpid->NumberOfSigmasTPC( nTrackXi,AliPID::kProton )) < 4) lIsNegProtonForTPC = kTRUE;
2354
2355         // Positive V0 daughter
2356         if (TMath::Abs(fESDpid->NumberOfSigmasTPC( pTrackXi,AliPID::kPion   )) < 4) lIsPosPionForTPC   = kTRUE;
2357         if (TMath::Abs(fESDpid->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 4) lIsPosProtonForTPC = kTRUE;
2358
2359         /*
2360         const AliExternalTrackParam *pInnerWallTrackXi    = pTrackXi    ->GetInnerParam(); // Do not use GetTPCInnerWall
2361         const AliExternalTrackParam *nInnerWallTrackXi    = nTrackXi    ->GetInnerParam();
2362         const AliExternalTrackParam *bachInnerWallTrackXi = bachTrackXi ->GetInnerParam();
2363         if(pInnerWallTrackXi && nInnerWallTrackXi && bachInnerWallTrackXi ){
2364
2365                 Double_t pMomInnerWall    = pInnerWallTrackXi   ->GetP();
2366                 Double_t nMomInnerWall    = nInnerWallTrackXi   ->GetP();
2367                 Double_t bachMomInnerWall = bachInnerWallTrackXi->GetP();
2368
2369                 // Bachelor
2370                 if (TMath::Abs(fESDpid->NumberOfSigmasTPC( bachTrackXi,AliPID::kPion)) < 3)                              lIsBachelorPionForTPC = kTRUE;
2371                 if (bachMomInnerWall < 0.350  && TMath::Abs(fESDpid->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 5) lIsBachelorKaonForTPC = kTRUE;
2372                 if (bachMomInnerWall > 0.350  && TMath::Abs(fESDpid->NumberOfSigmasTPC( bachTrackXi,AliPID::kKaon)) < 3) lIsBachelorKaonForTPC = kTRUE;
2373
2374                 // Negative V0 daughter
2375                 if (TMath::Abs(fESDpid->NumberOfSigmasTPC( nTrackXi,AliPID::kPion   )) < 3  )                           lIsNegPionForTPC   = kTRUE;
2376                 if (nMomInnerWall < 0.6  && TMath::Abs(fESDpid->NumberOfSigmasTPC( nTrackXi,AliPID::kProton ) ) < 5 )   lIsNegProtonForTPC = kTRUE;
2377                 if (nMomInnerWall > 0.6  && TMath::Abs(fESDpid->NumberOfSigmasTPC( nTrackXi,AliPID::kProton ) ) < 3 )   lIsNegProtonForTPC = kTRUE;
2378
2379                 // Positive V0 daughter
2380                 if (TMath::Abs(fESDpid->NumberOfSigmasTPC( pTrackXi,AliPID::kPion   )) < 3 )                            lIsPosPionForTPC   = kTRUE;
2381                 if (pMomInnerWall < 0.6  && TMath::Abs(fESDpid->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 5)     lIsPosProtonForTPC = kTRUE;
2382                 if (pMomInnerWall > 0.6  && TMath::Abs(fESDpid->NumberOfSigmasTPC( pTrackXi,AliPID::kProton )) < 3)     lIsPosProtonForTPC = kTRUE;
2383
2384         }
2385         */
2386
2387         // Combined PID TH1s
2388         if( lChargeXi < 0 && lIsBachelorPion )    fHistMassWithCombPIDXiMinus     ->Fill( lInvMassXiMinus    );
2389         if( lChargeXi > 0 && lIsBachelorPion )    fHistMassWithCombPIDXiPlus      ->Fill( lInvMassXiPlus     );
2390         if( lChargeXi < 0 && lIsBachelorKaon )    fHistMassWithCombPIDOmegaMinus  ->Fill( lInvMassOmegaMinus );
2391         if( lChargeXi > 0 && lIsBachelorKaon )    fHistMassWithCombPIDOmegaPlus   ->Fill( lInvMassOmegaPlus  );
2392
2393
2394         // 3.2 - PID proba Vs Pt(Bach)
2395         Int_t      lblBachForPID  = (Int_t) TMath::Abs( bachTrackXi->GetLabel() );
2396         TParticle* mcBachForPID   = lMCstack->Particle( lblBachForPID );
2397         Double_t   lmcPtBach      = mcBachForPID->Pt();
2398
2399         if(lIsBachelorPion)   f2dHistPIDprobaPionVsMCPtBach->Fill( lmcPtBach, ppionBach );
2400         if(lIsBachelorKaon)   f2dHistPIDprobaKaonVsMCPtBach->Fill( lmcPtBach, pkaonBach );
2401
2402
2403         // 3.3 - MC perfect PID
2404         Bool_t   lIsBachelorMCPiMinus  = kFALSE;
2405         Bool_t   lIsBachelorMCPiPlus   = kFALSE;
2406         Bool_t   lIsBachelorMCKMinus   = kFALSE;
2407         Bool_t   lIsBachelorMCKPlus    = kFALSE;
2408
2409         if( mcBachForPID->GetPdgCode() == -211) lIsBachelorMCPiMinus = kTRUE;
2410         if( mcBachForPID->GetPdgCode() ==  211) lIsBachelorMCPiPlus  = kTRUE;
2411         if( mcBachForPID->GetPdgCode() == -321) lIsBachelorMCKMinus  = kTRUE;
2412         if( mcBachForPID->GetPdgCode() ==  321) lIsBachelorMCKPlus   = kTRUE;
2413
2414         if( lChargeXi < 0 && lIsBachelorMCPiMinus )    fHistMassWithMcPIDXiMinus     ->Fill( lInvMassXiMinus );
2415         if( lChargeXi > 0 && lIsBachelorMCPiPlus  )    fHistMassWithMcPIDXiPlus      ->Fill( lInvMassXiPlus );
2416         if( lChargeXi < 0 && lIsBachelorMCKMinus  )    fHistMassWithMcPIDOmegaMinus  ->Fill( lInvMassOmegaMinus );
2417         if( lChargeXi > 0 && lIsBachelorMCKPlus   )    fHistMassWithMcPIDOmegaPlus   ->Fill( lInvMassOmegaPlus );
2418
2419
2420
2421
2422         // - Step II.4 : MC association (care : lots of "continue;" below this line)
2423         //-------------
2424
2425         Bool_t lAssoXiMinus    = kFALSE;
2426         Bool_t lAssoXiPlus     = kFALSE;
2427         Bool_t lAssoOmegaMinus = kFALSE;
2428         Bool_t lAssoOmegaPlus  = kFALSE;
2429
2430
2431         if(fDebug > 5)
2432                 cout    << "MC EventNumber : " << lMCevent->Header()->GetEvent()
2433                         << " / MC event Number in Run : " << lMCevent->Header()->GetEventNrInRun() << endl;
2434
2435         // - Step 4.1 : level of the V0 daughters
2436
2437         Int_t lblPosV0Dghter = (Int_t) TMath::Abs( pTrackXi->GetLabel() );
2438                 // Abs value = needed ! question of quality track association ...
2439         Int_t lblNegV0Dghter = (Int_t) TMath::Abs( nTrackXi->GetLabel() );
2440
2441         TParticle* mcPosV0Dghter = lMCstack->Particle( lblPosV0Dghter );
2442         TParticle* mcNegV0Dghter = lMCstack->Particle( lblNegV0Dghter );
2443
2444
2445         // - Step 4.2 : level of the Xi daughters
2446
2447         Int_t lblMotherPosV0Dghter = mcPosV0Dghter->GetFirstMother() ;
2448         Int_t lblMotherNegV0Dghter = mcNegV0Dghter->GetFirstMother();
2449
2450                 if( lblMotherPosV0Dghter != lblMotherNegV0Dghter) continue; // same mother
2451                 if( lblMotherPosV0Dghter < 0 ) continue; // mother != primary (!= -1)
2452                 if( lblMotherNegV0Dghter < 0 ) continue;
2453
2454
2455                 // mothers = Lambda candidate ... a priori
2456
2457         TParticle* mcMotherPosV0Dghter = lMCstack->Particle( lblMotherPosV0Dghter );
2458         TParticle* mcMotherNegV0Dghter = lMCstack->Particle( lblMotherNegV0Dghter );
2459
2460         Int_t      lblBach  = (Int_t) TMath::Abs( bachTrackXi->GetLabel() );
2461         TParticle* mcBach   = lMCstack->Particle( lblBach );
2462
2463
2464         // - Step 4.3 : level of Xi candidate
2465
2466         Int_t lblGdMotherPosV0Dghter =   mcMotherPosV0Dghter->GetFirstMother() ;
2467         Int_t lblGdMotherNegV0Dghter =   mcMotherNegV0Dghter->GetFirstMother() ;
2468
2469                 if( lblGdMotherPosV0Dghter != lblGdMotherNegV0Dghter ) continue;
2470                 if( lblGdMotherPosV0Dghter < 0 ) continue; // primary lambda ...
2471                 if( lblGdMotherNegV0Dghter < 0 ) continue; // primary lambda ...
2472
2473
2474                 // Gd mothers = Xi candidate ... a priori
2475
2476         TParticle* mcGdMotherPosV0Dghter = lMCstack->Particle( lblGdMotherPosV0Dghter );
2477         TParticle* mcGdMotherNegV0Dghter = lMCstack->Particle( lblGdMotherNegV0Dghter );
2478
2479         Int_t lblMotherBach = (Int_t) TMath::Abs( mcBach->GetFirstMother()  );
2480
2481                 if( lblMotherBach != lblGdMotherPosV0Dghter ) continue; //same mother for bach and V0 daughters
2482
2483         TParticle* mcMotherBach = lMCstack->Particle( lblMotherBach );
2484
2485
2486         // - Step 4.4 : Manage boolean for association
2487
2488         if( mcMotherBach                ->GetPdgCode() ==   3312 &&
2489             mcGdMotherPosV0Dghter       ->GetPdgCode() ==   3312 &&
2490             mcGdMotherNegV0Dghter       ->GetPdgCode() ==   3312)       lAssoXiMinus = kTRUE;
2491
2492         else if( mcMotherBach           ->GetPdgCode() ==  -3312 &&
2493             mcGdMotherPosV0Dghter       ->GetPdgCode() ==  -3312 &&
2494             mcGdMotherNegV0Dghter       ->GetPdgCode() ==  -3312)       lAssoXiPlus = kTRUE;
2495
2496         else if( mcMotherBach           ->GetPdgCode() ==   3334 &&
2497             mcGdMotherPosV0Dghter       ->GetPdgCode() ==   3334 &&
2498             mcGdMotherNegV0Dghter       ->GetPdgCode() ==   3334)       lAssoOmegaMinus = kTRUE;
2499
2500         else if( mcMotherBach           ->GetPdgCode() ==  -3334 &&
2501             mcGdMotherPosV0Dghter       ->GetPdgCode() ==  -3334 &&
2502             mcGdMotherNegV0Dghter       ->GetPdgCode() ==  -3334)       lAssoOmegaPlus = kTRUE;
2503
2504
2505
2506         if(!lAssoXiMinus && !lAssoXiPlus && !lAssoOmegaMinus && !lAssoOmegaPlus) continue; // no association, skip the rest of the code
2507         // If a proper association  exists ...
2508
2509         if(fDebug > 4){
2510                 cout << "XiMinus    = " << lAssoXiMinus    << endl;
2511                 cout << "XiPlus     = " << lAssoXiPlus     << endl;
2512                 cout << "OmegaMinus = " << lAssoOmegaMinus << endl;
2513                 cout << "OmegaPlus  = " << lAssoOmegaPlus  << endl
2514                      << "----"          << endl;
2515         }
2516
2517
2518         if(fDebug > 5){
2519                 cout << endl;
2520                 cout << "- V0 daughters - " << endl;
2521                 cout << "     + V0 Pos. / Label : " << lblPosV0Dghter
2522                 << " - Pdg Code : " << mcPosV0Dghter->GetTitle() << endl;
2523                 cout << "     - V0 Neg. / Label : " << lblNegV0Dghter
2524                 << " - Pdg Code : " << mcNegV0Dghter->GetTitle() << endl;
2525
2526                 cout << "- Xi daughters - " << endl;
2527                 cout << "     + V0 Pos. mother / Label : " << lblMotherPosV0Dghter
2528                 << " - Pdg Code : " << mcMotherPosV0Dghter->GetTitle() << endl;
2529                 cout << "     - V0 Neg. mother / Label : " << lblMotherNegV0Dghter
2530                 << " - Pdg Code : " << mcMotherNegV0Dghter->GetTitle() << endl;
2531
2532                 cout << "     --  Bach. / Label :" << lblBach
2533                 << " -  Pdg Code : " << mcBach->GetTitle() << endl;
2534
2535                 cout << "- Xi candidate -" << endl;
2536                 cout << "    +  V0 Pos. Gd Mother / Label : " << lblGdMotherPosV0Dghter
2537                 << " - Pdg Code : " << mcGdMotherPosV0Dghter->GetTitle() << endl;
2538                 cout << "    -  V0 Neg. Gd Mother / Label : "  << lblGdMotherNegV0Dghter
2539                 << " - Pdg Code : "<< mcGdMotherNegV0Dghter->GetTitle() << endl;
2540
2541                 cout << "    --  Mother Bach. / Label : " << lblMotherBach
2542                 << " - Pdg Code    : " << mcMotherBach->GetTitle() << endl;
2543                 cout << endl;
2544         }
2545
2546
2547         // - Step 5 : Plots around the cascade candidates associated with MC
2548         //-------------
2549
2550         Double_t lmcPt             = mcMotherBach->Pt();
2551         Double_t lmcRapCasc        = 0.5*TMath::Log( (mcMotherBach->Energy() + mcMotherBach->Pz()) /
2552                                                      (mcMotherBach->Energy() - mcMotherBach->Pz() +1.e-13) );
2553         Double_t lmcPtBaryon       = (mcNegV0Dghter->GetPdgCode() == -2212 ) ? mcNegV0Dghter->Pt() : mcPosV0Dghter->Pt(); // spot the baryon daughter : pbar or p
2554
2555
2556         Double_t lmcEta            = mcMotherBach->Eta();
2557         Double_t lmcTransvRadius   = mcBach->R(); // to get the decay point of Xi, = the production vertex of Bachelor ...
2558
2559         TVector3 lmcTVect3Mom( mcMotherBach->Px(), mcMotherBach->Py(), mcMotherBach->Pz() );
2560
2561         Double_t lrecoPt           = xiESD->Pt();
2562         Double_t lrecoTransvRadius = TMath::Sqrt( xiESD->Xv() * xiESD->Xv() + xiESD->Yv() * xiESD->Yv() );
2563
2564         TVector3 lrecoTVect3Mom( xiESD->Px(), xiESD->Py(), xiESD->Pz()  );
2565         Double_t lDeltaPhiMcReco   = lmcTVect3Mom.DeltaPhi( lrecoTVect3Mom ) * 180.0/TMath::Pi();
2566
2567
2568         // - Histos for the cascade candidates associated with MC
2569
2570         if( lChargeXi < 0 && lAssoXiMinus){
2571                 fHistAsMCMassXiMinus          ->Fill( lInvMassXiMinus  );
2572                 if(lIsBachelorPion)     f2dHistAsMCandCombPIDGenPtVsGenYXiMinus->Fill( lmcPt, lmcRapCasc );
2573                 f2dHistAsMCGenPtVsGenYXiMinus ->Fill( lmcPt, lmcRapCasc);
2574                 fHistAsMCGenEtaXiMinus        ->Fill( lmcEta           );
2575                 f2dHistAsMCResPtXiMinus       ->Fill( lmcPt,           (lrecoPt - lmcPt)/ lmcPt );
2576                 f2dHistAsMCResRXiMinus        ->Fill( lmcTransvRadius, (lrecoTransvRadius - lmcTransvRadius)/ lmcTransvRadius    );
2577                 f2dHistAsMCResPhiXiMinus      ->Fill( lmcPt, lDeltaPhiMcReco );
2578                 f2dHistAsMCPtProtonVsPtXiMinus->Fill( lmcPtBaryon, lmcPt );
2579         }
2580
2581         else if( lChargeXi > 0 && lAssoXiPlus){
2582                 fHistAsMCMassXiPlus                     ->Fill( lInvMassXiPlus   );
2583                 if(lIsBachelorPion)     f2dHistAsMCandCombPIDGenPtVsGenYXiPlus->Fill( lmcPt, lmcRapCasc );
2584                 f2dHistAsMCGenPtVsGenYXiPlus            ->Fill( lmcPt, lmcRapCasc);
2585                 fHistAsMCGenEtaXiPlus                   ->Fill( lmcEta           );
2586                 f2dHistAsMCResPtXiPlus                  ->Fill( lmcPt,           (lrecoPt - lmcPt)/ lmcPt );
2587                 f2dHistAsMCResRXiPlus                   ->Fill( lmcTransvRadius, (lrecoTransvRadius - lmcTransvRadius)/ lmcTransvRadius    );
2588                 f2dHistAsMCResPhiXiPlus                 ->Fill( lmcPt, lDeltaPhiMcReco );
2589                 f2dHistAsMCPtAntiProtonVsPtXiPlus       ->Fill( lmcPtBaryon, lmcPt );
2590
2591         }
2592
2593         else if( lChargeXi < 0 && lAssoOmegaMinus){
2594                 fHistAsMCMassOmegaMinus          ->Fill( lInvMassOmegaMinus );
2595                 if(lIsBachelorKaon)     f2dHistAsMCandCombPIDGenPtVsGenYOmegaMinus->Fill( lmcPt, lmcRapCasc );
2596                 f2dHistAsMCGenPtVsGenYOmegaMinus ->Fill( lmcPt, lmcRapCasc  );
2597                 fHistAsMCGenEtaOmegaMinus        ->Fill( lmcEta             );
2598                 f2dHistAsMCResPtOmegaMinus       ->Fill( lmcPt,           (lrecoPt - lmcPt)/ lmcPt );
2599                 f2dHistAsMCResROmegaMinus        ->Fill( lmcTransvRadius, (lrecoTransvRadius - lmcTransvRadius)/ lmcTransvRadius    );
2600                 f2dHistAsMCResPhiOmegaMinus      ->Fill( lmcPt, lDeltaPhiMcReco );
2601                 f2dHistAsMCPtProtonVsPtOmegaMinus->Fill( lmcPtBaryon, lmcPt );
2602         }
2603
2604         else if( lChargeXi > 0 && lAssoOmegaPlus){
2605                 fHistAsMCMassOmegaPlus                  ->Fill( lInvMassOmegaPlus );
2606                 if(lIsBachelorKaon)     f2dHistAsMCandCombPIDGenPtVsGenYOmegaPlus->Fill( lmcPt, lmcRapCasc );
2607                 f2dHistAsMCGenPtVsGenYOmegaPlus         ->Fill( lmcPt, lmcRapCasc );
2608                 fHistAsMCGenEtaOmegaPlus                ->Fill( lmcEta            );
2609                 f2dHistAsMCResPtOmegaPlus               ->Fill( lmcPt,           (lrecoPt - lmcPt)/ lmcPt );
2610                 f2dHistAsMCResROmegaPlus                ->Fill( lmcTransvRadius, (lrecoTransvRadius - lmcTransvRadius)/ lmcTransvRadius    );
2611                 f2dHistAsMCResPhiOmegaPlus              ->Fill( lmcPt, lDeltaPhiMcReco );
2612                 f2dHistAsMCPtAntiProtonVsPtOmegaPlus    ->Fill( lmcPtBaryon, lmcPt );
2613         }
2614
2615
2616         // - Step 6 : Containers = Cascade cuts + PID
2617         //-------------
2618
2619         // Double_t lChi2Xi         = -1. ;
2620         Double_t lDcaXiDaughters            = -1. ;
2621         Double_t lDcaBachToPrimVertexXi     = -1. ;
2622         Double_t lXiCosineOfPointingAngle   = -1. ;
2623         Double_t lPosXi[3]                  = { -1000.0, -1000.0, -1000.0 };
2624         Double_t lXiRadius                  = -1000. ;
2625
2626         Double_t lInvMassLambdaAsCascDghter = 0.;
2627         Double_t lDcaV0DaughtersXi          = -1.;
2628         // Double_t lV0Chi2Xi               = -1. ;
2629         Double_t lV0CosineOfPointingAngleXi = -1.;
2630         Double_t lPosV0Xi[3]                = { -1000. , -1000., -1000. }; // Position of VO coming from cascade
2631         Double_t lV0RadiusXi                = -1000.;
2632         Double_t lDcaV0ToPrimVertexXi       = -1.;
2633         Double_t lDcaPosToPrimVertexXi      = -1.;
2634         Double_t lDcaNegToPrimVertexXi      = -1.;
2635
2636         //Int_t    nTrackWithTPCrefitMultiplicity  =  0;
2637
2638
2639         // 6.2 - Definition of the needed variables
2640
2641         //lChi2Xi                       = xiESD->GetChi2Xi();
2642         lDcaXiDaughters                 = xiESD->GetDcaXiDaughters();
2643         lDcaBachToPrimVertexXi          = TMath::Abs( bachTrackXi->GetD( lBestPrimaryVtxPos[0],
2644                                                                          lBestPrimaryVtxPos[1],
2645                                                                          lMagneticField  ) );
2646                                         // NOTE : AliExternalTrackParam::GetD returns an algebraic value
2647         lXiCosineOfPointingAngle        = xiESD->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0],
2648                                                                                   lBestPrimaryVtxPos[1],
2649                                                                                   lBestPrimaryVtxPos[2] );
2650                                         // Take care : the best available vertex should be used (like in AliCascadeVertexer)
2651                 xiESD->GetXYZcascade( lPosXi[0],  lPosXi[1], lPosXi[2] );
2652         lXiRadius                       = TMath::Sqrt( lPosXi[0]*lPosXi[0]  +  lPosXi[1]*lPosXi[1] );
2653         lInvMassLambdaAsCascDghter      = xiESD->GetEffMass();
2654                                         // This value shouldn't change, whatever the working hyp. is : Xi-, Xi+, Omega-, Omega+
2655         lDcaV0DaughtersXi               = xiESD->GetDcaV0Daughters();
2656         // lV0Chi2Xi                    = xiESD->GetChi2V0();
2657         lV0CosineOfPointingAngleXi      = xiESD->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0],
2658                                                                              lBestPrimaryVtxPos[1],
2659                                                                              lBestPrimaryVtxPos[2] );
2660                 xiESD->GetXYZ( lPosV0Xi[0],  lPosV0Xi[1], lPosV0Xi[2] );
2661         lV0RadiusXi                     = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0]  +  lPosV0Xi[1]*lPosV0Xi[1] );
2662
2663         lDcaV0ToPrimVertexXi            = xiESD->GetD( lBestPrimaryVtxPos[0],
2664                                                 lBestPrimaryVtxPos[1],
2665                                                 lBestPrimaryVtxPos[2] );
2666
2667         lDcaPosToPrimVertexXi           = TMath::Abs( pTrackXi   ->GetD( lBestPrimaryVtxPos[0],
2668                                                                          lBestPrimaryVtxPos[1],
2669                                                                          lMagneticField  )     );
2670
2671         lDcaNegToPrimVertexXi           = TMath::Abs( nTrackXi   ->GetD( lBestPrimaryVtxPos[0],
2672                                                                          lBestPrimaryVtxPos[1],
2673                                                                          lMagneticField  )     );
2674
2675
2676
2677         //nTrackWithTPCrefitMultiplicity  = DoESDTrackWithTPCrefitMultiplicity(lESDevent);  // FIXME : variable which is not used anymore at the moment ...
2678                                                                                           //    ->  keep it while the task is still under development.
2679
2680
2681
2682         // 6.3 - Filling the AliCFContainer (optimisation of topological selections + systematics)
2683         Double_t lContainerCutVars[20] = {0.0};
2684
2685         lContainerCutVars[0]  = lDcaXiDaughters;
2686         lContainerCutVars[1]  = lDcaBachToPrimVertexXi;
2687         lContainerCutVars[2]  = lXiCosineOfPointingAngle;
2688         lContainerCutVars[3]  = lXiRadius;
2689         lContainerCutVars[4]  = lInvMassLambdaAsCascDghter;
2690         lContainerCutVars[5]  = lDcaV0DaughtersXi;
2691         lContainerCutVars[6]  = lV0CosineOfPointingAngleXi;
2692         lContainerCutVars[7]  = lV0RadiusXi;
2693         lContainerCutVars[8]  = lDcaV0ToPrimVertexXi;
2694         lContainerCutVars[9]  = lDcaPosToPrimVertexXi;
2695         lContainerCutVars[10] = lDcaNegToPrimVertexXi;
2696
2697         lContainerCutVars[13] = lmcPt;
2698
2699         lContainerCutVars[16] = lBestPrimaryVtxPos[2];
2700         lContainerCutVars[17] = nTrackPrimaryMultiplicity;       // FIXME : nTrackPrimaryMultiplicity not checked for AOD ...
2701         lContainerCutVars[18] = nITSandTPCtracksAndSPDtracklets; // FIXME : nITSandTPCtracksAndSPDtracklets is not available for AOD ...
2702         lContainerCutVars[19] = lBachTPCClusters;                // FIXME : BachTPCClusters          is not available for AOD ...
2703
2704         // All cases should be covered below
2705         if( lChargeXi < 0 && lAssoXiMinus    ) {
2706                 lContainerCutVars[11] = lInvMassXiMinus;
2707                 lContainerCutVars[12] = 1.63;
2708                 lContainerCutVars[14] = lmcRapCasc;
2709                 lContainerCutVars[15] = -1.;
2710                         if( lIsBachelorPionForTPC   && lIsPosProtonForTPC    && lIsNegPionForTPC )
2711                                 fCFContAsCascadeCuts->Fill(lContainerCutVars,0); // for Xi-
2712         }
2713         if( lChargeXi > 0 && lAssoXiPlus     ){
2714                 lContainerCutVars[11] = lInvMassXiPlus;
2715                 lContainerCutVars[12] = 1.26;
2716                 lContainerCutVars[14] = lmcRapCasc;
2717                 lContainerCutVars[15] = -1.;
2718                         if( lIsBachelorPionForTPC   && lIsNegProtonForTPC    && lIsPosPionForTPC )
2719                                 fCFContAsCascadeCuts->Fill(lContainerCutVars,1); // for Xi+
2720         }
2721         if( lChargeXi < 0 && lAssoOmegaMinus )  {
2722                 lContainerCutVars[11] = 1.63;
2723                 lContainerCutVars[12] = lInvMassOmegaMinus;
2724                 lContainerCutVars[14] = -1.;
2725                 lContainerCutVars[15] = lmcRapCasc;
2726                         if( lIsBachelorKaonForTPC   && lIsPosProtonForTPC    && lIsNegPionForTPC  && (TMath::Abs( lInvMassXiMinus-1.3217 ) > 0.008) )
2727                                 fCFContAsCascadeCuts->Fill(lContainerCutVars,2); // for Omega-
2728         }
2729         if( lChargeXi > 0 && lAssoOmegaPlus  ){
2730                 lContainerCutVars[11] = 1.26;
2731                 lContainerCutVars[12] = lInvMassOmegaPlus;
2732                 lContainerCutVars[14] = -1.;
2733                 lContainerCutVars[15] = lmcRapCasc;
2734                         if( lIsBachelorKaonForTPC   && lIsNegProtonForTPC    && lIsPosPionForTPC && (TMath::Abs( lInvMassXiPlus-1.3217 ) > 0.008) )
2735                                 fCFContAsCascadeCuts->Fill(lContainerCutVars,3); // for Omega+
2736         }
2737
2738
2739         // 6.4 - Filling the AliCFContainers related to PID
2740
2741         Double_t lContainerPIDVars[4] = {0.0};
2742
2743
2744         // Xi Minus
2745         if( lChargeXi < 0 && lAssoXiMinus ) {
2746                 lContainerPIDVars[0] = lmcPt              ;
2747                 lContainerPIDVars[1] = lInvMassXiMinus    ;
2748                 lContainerPIDVars[2] = lmcRapCasc         ;
2749                 lContainerPIDVars[3] = nTrackPrimaryMultiplicity ;  // FIXME : nTrackPrimaryMultiplicity is not checked for AOD ...
2750
2751                 // No PID
2752                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 0); // No PID
2753                 // TPC PID
2754                 if( lIsBachelorPionForTPC  )
2755                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 1); // TPC PID / 4-#sigma cut on Bachelor track
2756
2757                 if( lIsBachelorPionForTPC &&
2758                     lIsPosProtonForTPC     )
2759                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 2); // TPC PID / 4-#sigma cut on Bachelor+Baryon tracks
2760
2761                 if( lIsBachelorPionForTPC &&
2762                     lIsPosProtonForTPC    &&
2763                     lIsNegPionForTPC       )
2764                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 3); // TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks
2765
2766                 // Combined PID
2767                 if( lIsBachelorPion        )
2768                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 4); // Comb. PID / Bachelor
2769
2770                 if( lIsBachelorPion       &&
2771                     lIsPosInXiProton    )
2772                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 5); // Comb. PID / Bachelor+Baryon
2773
2774                 if(lIsBachelorPion     &&
2775                    lIsPosInXiProton &&
2776                    lIsNegInXiPion    )
2777                         fCFContCascadePIDAsXiMinus->Fill(lContainerPIDVars, 6); // Comb. PID / Bachelor+Baryon+Meson
2778         }
2779
2780         lContainerPIDVars[0] = 0.; lContainerPIDVars[1] = 0.; lContainerPIDVars[2] = 0.; lContainerPIDVars[3] = 0.;
2781
2782         // Xi Plus
2783         if( lChargeXi > 0 && lAssoXiPlus ) {
2784                 lContainerPIDVars[0] = lmcPt              ;
2785                 lContainerPIDVars[1] = lInvMassXiPlus     ;
2786                 lContainerPIDVars[2] = lmcRapCasc           ;
2787                 lContainerPIDVars[3] = nTrackPrimaryMultiplicity ;    // FIXME : nTrackPrimaryMultiplicity is not checked for AOD ...
2788
2789                 // No PID
2790                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 0); // No PID
2791                 // TPC PID
2792                 if( lIsBachelorPionForTPC  )
2793                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 1); // TPC PID / 4-#sigma cut on Bachelor track
2794
2795                 if( lIsBachelorPionForTPC &&
2796                     lIsNegProtonForTPC     )
2797                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 2); // TPC PID / 4-#sigma cut on Bachelor+Baryon tracks
2798
2799                 if( lIsBachelorPionForTPC &&
2800                     lIsNegProtonForTPC    &&
2801                     lIsPosPionForTPC       )
2802                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 3); // TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks
2803
2804                 // Combined PID
2805                 if( lIsBachelorPion        )
2806                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 4); // Comb. PID / Bachelor
2807
2808                 if( lIsBachelorPion       &&
2809                     lIsNegInXiProton    )
2810                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 5); // Comb. PID / Bachelor+Baryon
2811
2812                 if(lIsBachelorPion     &&
2813                    lIsNegInXiProton &&
2814                    lIsPosInXiPion    )
2815                         fCFContCascadePIDAsXiPlus->Fill(lContainerPIDVars, 6); // Comb. PID / Bachelor+Baryon+Meson
2816         }
2817
2818         lContainerPIDVars[0] = 0.; lContainerPIDVars[1] = 0.; lContainerPIDVars[2] = 0.; lContainerPIDVars[3] = 0.;
2819
2820         // Omega Minus
2821         if( lChargeXi < 0 && lAssoOmegaMinus ) {
2822                 lContainerPIDVars[0] = lmcPt              ;
2823                 lContainerPIDVars[1] = lInvMassOmegaMinus ;
2824                 lContainerPIDVars[2] = lmcRapCasc         ;
2825                 lContainerPIDVars[3] = nTrackPrimaryMultiplicity ;  // FIXME : nTrackPrimaryMultiplicity is not checked for AOD ...
2826
2827                 // No PID
2828                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 0); // No PID
2829                 // TPC PID
2830                 if( lIsBachelorKaonForTPC  )
2831                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 1); // TPC PID / 4-#sigma cut on Bachelor track
2832
2833                 if( lIsBachelorKaonForTPC &&
2834                     lIsPosProtonForTPC     )
2835                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 2); // TPC PID / 4-#sigma cut on Bachelor+Baryon tracks
2836
2837                 if( lIsBachelorKaonForTPC &&
2838                     lIsPosProtonForTPC    &&
2839                     lIsNegPionForTPC       )
2840                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 3); // TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks
2841
2842                 // Combined PID
2843                 if( lIsBachelorKaon        )
2844                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 4); // Comb. PID / Bachelor
2845
2846                 if( lIsBachelorKaon       &&
2847                     lIsPosInOmegaProton    )
2848                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 5); // Comb. PID / Bachelor+Baryon
2849
2850                 if(lIsBachelorKaon     &&
2851                    lIsPosInOmegaProton &&
2852                    lIsNegInOmegaPion    )
2853                         fCFContCascadePIDAsOmegaMinus->Fill(lContainerPIDVars, 6); // Comb. PID / Bachelor+Baryon+Meson
2854         }
2855
2856         lContainerPIDVars[0] = 0.; lContainerPIDVars[1] = 0.; lContainerPIDVars[2] = 0.; lContainerPIDVars[3] = 0.;
2857
2858         // Omega Plus
2859         if( lChargeXi > 0 && lAssoOmegaPlus) {
2860                 lContainerPIDVars[0] = lmcPt              ;
2861                 lContainerPIDVars[1] = lInvMassOmegaPlus  ;
2862                 lContainerPIDVars[2] = lmcRapCasc         ;
2863                 lContainerPIDVars[3] = nTrackPrimaryMultiplicity ;  // FIXME : nTrackPrimaryMultiplicity is not checked for AOD ...
2864
2865                 // No PID
2866                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 0); // No PID
2867                 // TPC PID
2868                 if( lIsBachelorKaonForTPC  )
2869                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 1); // TPC PID / 4-#sigma cut on Bachelor track
2870
2871                 if( lIsBachelorKaonForTPC &&
2872                     lIsNegProtonForTPC     )
2873                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 2); // TPC PID / 4-#sigma cut on Bachelor+Baryon tracks
2874
2875                 if( lIsBachelorKaonForTPC &&
2876                     lIsNegProtonForTPC    &&
2877                     lIsPosPionForTPC       )
2878                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 3); // TPC PID / 4-#sigma cut on Bachelor+Baryon+Meson tracks
2879
2880                 // Combined PID
2881                 if( lIsBachelorKaon        )
2882                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 4); // Comb. PID / Bachelor
2883
2884                 if( lIsBachelorKaon       &&
2885                     lIsNegInOmegaProton    )
2886                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 5); // Comb. PID / Bachelor+Baryon
2887
2888                 if(lIsBachelorKaon     &&
2889                    lIsNegInOmegaProton &&
2890                    lIsPosInOmegaPion    )
2891                         fCFContCascadePIDAsOmegaPlus->Fill(lContainerPIDVars, 6); // Comb. PID / Bachelor+Baryon+Meson
2892         }
2893
2894
2895
2896 }// End of loop over reconstructed cascades
2897
2898
2899
2900
2901   // Post output data.
2902  PostData(1, fListHistCascade);
2903 }
2904
2905
2906
2907 Int_t AliAnalysisTaskCheckPerformanceCascade::DoESDTrackWithTPCrefitMultiplicity(const AliESDEvent *lESDevent)
2908 {
2909     // Checking the number of tracks with TPCrefit for each event
2910     // Needed for a rough assessment of the event multiplicity
2911
2912         Int_t nTrackWithTPCrefitMultiplicity = 0;
2913         for(Int_t iTrackIdx = 0; iTrackIdx < (InputEvent())->GetNumberOfTracks(); iTrackIdx++){
2914                 AliESDtrack *esdTrack   = 0x0;
2915                              esdTrack   = lESDevent->GetTrack( iTrackIdx );
2916                 if (!esdTrack) { AliWarning("Pb / Could not retrieve one track within the track loop for TPCrefit check ..."); continue; }
2917
2918                 ULong_t lTrackStatus    = esdTrack->GetStatus();
2919                             if ((lTrackStatus&AliESDtrack::kTPCrefit)    == 0) continue;
2920                             else nTrackWithTPCrefitMultiplicity++;
2921                     // FIXME :
2922                     // The goal here is to get a better assessment of the event multiplicity.
2923                     // (InputEvent())->GetNumberOfTracks() takes into account ITS std alone tracks + global tracks
2924                     // This may introduce a bias. Hence the number of TPC refit tracks.
2925                     // Note : the event multiplicity = analysis on its own... See Jacek's or Jan Fiete's analysis on dN/d(pt) and dN/d(eta)
2926
2927         }// end loop over all event tracks
2928         return  nTrackWithTPCrefitMultiplicity;
2929 }
2930
2931
2932
2933
2934
2935
2936
2937 //________________________________________________________________________
2938 void AliAnalysisTaskCheckPerformanceCascade::Terminate(Option_t *)
2939 {
2940   // Draw result to the screen
2941   // Called once at the end of the query
2942
2943   TList *cRetrievedList = 0x0;
2944   cRetrievedList = (TList*)GetOutputData(1);
2945   if(!cRetrievedList){
2946         Printf("ERROR - AliAnalysisTaskCheckPerformanceCascade : ouput data container list not available\n");
2947         return;
2948   }
2949
2950   fHistMCTrackMultiplicity = dynamic_cast<TH1F*> (  cRetrievedList->FindObject("fHistMCTrackMultiplicity")  );
2951   if (!fHistMCTrackMultiplicity) {
2952     Printf("ERROR - AliAnalysisTaskCheckPerformanceCascade : fHistMCTrackMultiplicity not available");
2953     return;
2954   }
2955
2956
2957   TCanvas *canCheckPerformanceCascade = new TCanvas("AliAnalysisTaskCheckPerformanceCascade","Multiplicity",10,10,510,510);
2958   canCheckPerformanceCascade->cd(1)->SetLogy();
2959
2960   fHistMCTrackMultiplicity->SetMarkerStyle(22);
2961   fHistMCTrackMultiplicity->DrawCopy("E");
2962
2963 }