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