PWGLF/STRANGENESS/Cascades/AliAnalysisTaskStrangenessVsMultiplicity.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
  17 //
  18 // This task is meant to explore the possibility of using a VZERO amplitude
  19 // based multiplicity estimator for proton-proton collisions. For this, two
  20 // main operation methods for this task are foreseen:
  21 //
  22 //  1) (under development) it should act as an auxiliary task and provide a
  23 //     calibrated estimator
  24 //
  25 //  2) "Debug mode" which will also create a ROOT TTree object with event
  26 //     by event info potentially used for exploration / calibration. This
  27 //     includes the following info:
  28 //
  29 //      --- All VZERO Amplitudes (saved as Float_t)
  30 //      --- (optional) time for each channel
  31 //      --- (optional) time width for each channel
  32 //      --- GetReferenceMultiplicity Estimator, |eta|<0.5
  33 //      --- GetReferenceMultiplicity Estimator, |eta|<0.8
  34 //      --- (if MC) True Multiplicity, |eta|<0.5
  35 //      --- (if MC) True Multiplicity,  2.8 < eta < 5.8 (VZEROA region)
  36 //      --- (if MC) True Multiplicity, -3.7 < eta <-1.7 (VZEROC region)
  37 //      --- Run Number
  38 //
  39 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
  40
  41 class TTree;
  42 class TParticle;
  43 class TVector3;
  44
  45 //class AliMCEventHandler;
  46 //class AliMCEvent;
  47 //class AliStack;
  48
  49 class AliESDVertex;
  50 class AliAODVertex;
  51 class AliESDv0;
  52 class AliAODv0;
  53
  54 #include <Riostream.h>
  55 #include "TList.h"
  56 #include "TH1.h"
  57 #include "TH2.h"
  58 #include "TH3.h"
  59 #include "TFile.h"
  60 #include "THnSparse.h"
  61 #include "TVector3.h"
  62 #include "TCanvas.h"
  63 #include "TMath.h"
  64 #include "TLegend.h"
  65 //#include "AliLog.h"
  66
  67 #include "AliESDEvent.h"
  68 #include "AliAODEvent.h"
  69 #include "AliV0vertexer.h"
  70 #include "AliCascadeVertexer.h"
  71 #include "AliESDpid.h"
  72 #include "AliESDtrack.h"
  73 #include "AliESDtrackCuts.h"
  74 #include "AliInputEventHandler.h"
  75 #include "AliAnalysisManager.h"
  76 #include "AliMCEventHandler.h"
  77 #include "AliMCEvent.h"
  78 #include "AliStack.h"
  79 #include "AliCentrality.h"
  80 #include "AliPPVsMultUtils.h"
  81
  82 #include "AliCFContainer.h"
  83 #include "AliMultiplicity.h"
  84 #include "AliAODMCParticle.h"
  85 #include "AliESDcascade.h"
  86 #include "AliAODcascade.h"
  87 #include "AliESDUtils.h"
  88 #include "AliGenEventHeader.h"
  89 #include "AliAnalysisTaskSE.h"
  90 #include "AliAnalysisUtils.h"
  91 #include "AliAnalysisTaskStrangenessVsMultiplicity.h"
  92
  93 using std::cout;
  94 using std::endl;
  95
  96 ClassImp(AliAnalysisTaskStrangenessVsMultiplicity)
  97
  98 AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity()
  99   : AliAnalysisTaskSE(), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0), fPPVsMultUtils(0),
 100   fkSaveV0Tree      ( kFALSE ),
 101   fkSaveCascadeTree ( kTRUE  ),
 102   fkRunVertexers    ( kTRUE  ),
 103   fkSkipEventSelection( kFALSE ),
 104   //---> Variables for fTreeEvent
 105   fAmplitude_V0A   (0),
 106   fAmplitude_V0C   (0),
 107   fAmplitude_V0AEq (0),
 108   fAmplitude_V0CEq (0),
 109   fCentrality_V0A(0),
 110   fCentrality_V0C(0),
 111   fCentrality_V0M(0),
 112   fCentrality_V0AEq(0),
 113   fCentrality_V0CEq(0),
 114   fCentrality_V0MEq(0),
 115   fCustomCentrality_V0M(0),
 116   fCustomCentrality_V0MEq(0),
 117   fRefMultEta5(0),
 118   fRefMultEta8(0),
 119   fRunNumber(0),
 120   fEvSel_HasAtLeastSPDVertex(0),
 121   fEvSel_VtxZCut(0),
 122   fEvSel_IsNotPileup(0),
 123   fEvSel_IsNotPileupInMultBins(0),
 124   fEvSel_HasVtxContributor(0),
 125   fEvSel_nTracklets(0),
 126   fEvSel_nSPDClusters(0),
 127
 128   //---> Variables for fTreeV0
 129         fTreeVariableChi2V0(0),
 130         fTreeVariableDcaV0Daughters(0),
 131         fTreeVariableDcaV0ToPrimVertex(0),
 132         fTreeVariableDcaPosToPrimVertex(0),
 133         fTreeVariableDcaNegToPrimVertex(0),
 134         fTreeVariableV0CosineOfPointingAngle(0),
 135         fTreeVariableV0Radius(0),
 136         fTreeVariablePt(0),
 137         fTreeVariableRapK0Short(0),
 138         fTreeVariableRapLambda(0),
 139         fTreeVariableInvMassK0s(0),
 140         fTreeVariableInvMassLambda(0),
 141         fTreeVariableInvMassAntiLambda(0),
 142         fTreeVariableAlphaV0(0),
 143         fTreeVariablePtArmV0(0),
 144         fTreeVariableNegEta(0),
 145         fTreeVariablePosEta(0),
 146
 147         fTreeVariableNSigmasPosProton(0),
 148         fTreeVariableNSigmasPosPion(0),
 149         fTreeVariableNSigmasNegProton(0),
 150         fTreeVariableNSigmasNegPion(0),
 151
 152         fTreeVariableDistOverTotMom(0),
 153         fTreeVariableLeastNbrCrossedRows(0),
 154         fTreeVariableLeastRatioCrossedRowsOverFindable(0),
 155
 156   fTreeVariableCentV0A(0),
 157   fTreeVariableCentV0C(0),
 158   fTreeVariableCentV0M(0),
 159   fTreeVariableCentV0AEq(0),
 160   fTreeVariableCentV0CEq(0),
 161   fTreeVariableCentV0MEq(0),
 162   fTreeVariableAmpV0A(0),
 163   fTreeVariableAmpV0C(0),
 164   fTreeVariableAmpV0AEq(0),
 165   fTreeVariableAmpV0CEq(0),
 166   fTreeVariableRefMultEta8(0),
 167   fTreeVariableRefMultEta5(0),
 168   fTreeVariableRunNumber(0),
 169   //---> Variables for fTreeCascade
 170   fTreeCascVarCharge(0),
 171         fTreeCascVarMassAsXi(0),
 172         fTreeCascVarMassAsOmega(0),
 173         fTreeCascVarPt(0),
 174         fTreeCascVarRapXi(0),
 175         fTreeCascVarRapOmega(0),
 176         fTreeCascVarNegEta(0),
 177         fTreeCascVarPosEta(0),
 178         fTreeCascVarBachEta(0),
 179         fTreeCascVarDCACascDaughters(0),
 180         fTreeCascVarDCABachToPrimVtx(0),
 181         fTreeCascVarDCAV0Daughters(0),
 182         fTreeCascVarDCAV0ToPrimVtx(0),
 183         fTreeCascVarDCAPosToPrimVtx(0),
 184         fTreeCascVarDCANegToPrimVtx(0),
 185         fTreeCascVarCascCosPointingAngle(0),
 186         fTreeCascVarCascRadius(0),
 187         fTreeCascVarV0Mass(0),
 188         fTreeCascVarV0CosPointingAngle(0),
 189         fTreeCascVarV0CosPointingAngleSpecial(0),
 190         fTreeCascVarV0Radius(0),
 191   fTreeCascVarLeastNbrClusters(0),
 192         fTreeCascVarDistOverTotMom(0),
 193         fTreeCascVarNegNSigmaPion(0),
 194         fTreeCascVarNegNSigmaProton(0),
 195         fTreeCascVarPosNSigmaPion(0),
 196         fTreeCascVarPosNSigmaProton(0),
 197         fTreeCascVarBachNSigmaPion(0),
 198         fTreeCascVarBachNSigmaKaon(0),
 199         fTreeCascVarCentV0A(0),
 200         fTreeCascVarCentV0C(0),
 201         fTreeCascVarCentV0M(0),
 202         fTreeCascVarCentV0AEq(0),
 203         fTreeCascVarCentV0CEq(0),
 204         fTreeCascVarCentV0MEq(0),
 205         fTreeCascVarAmpV0A(0),
 206         fTreeCascVarAmpV0C(0),
 207         fTreeCascVarAmpV0AEq(0),
 208         fTreeCascVarAmpV0CEq(0),
 209   fTreeCascVarRefMultEta8(0),
 210   fTreeCascVarRefMultEta5(0),
 211   fTreeCascVarRunNumber(0),
 212   //Histos
 213   fHistEventCounter(0)
 214 //------------------------------------------------
 215 // Tree Variables
 216 {
 217
 218 }
 219
 220 AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity(const char *name)
 221   : AliAnalysisTaskSE(name), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0), fPPVsMultUtils(0),
 222   fkSaveV0Tree      ( kFALSE ),
 223   fkSaveCascadeTree ( kTRUE  ),
 224   fkRunVertexers    ( kTRUE  ),
 225   fkSkipEventSelection( kFALSE ),
 226   //---> Variables for fTreeEvent
 227   fAmplitude_V0A (0),
 228   fAmplitude_V0C (0),
 229   fAmplitude_V0AEq (0),
 230   fAmplitude_V0CEq (0),
 231   fCentrality_V0A(0),
 232   fCentrality_V0C(0),
 233   fCentrality_V0M(0),
 234   fCentrality_V0AEq(0),
 235   fCentrality_V0CEq(0),
 236   fCentrality_V0MEq(0),
 237   fCustomCentrality_V0M(0),
 238   fCustomCentrality_V0MEq(0),
 239   fRefMultEta5(0),
 240   fRefMultEta8(0),
 241   fRunNumber(0),
 242   fEvSel_HasAtLeastSPDVertex(0),
 243   fEvSel_VtxZCut(0),
 244   fEvSel_IsNotPileup(0),
 245   fEvSel_IsNotPileupInMultBins(0),
 246   fEvSel_HasVtxContributor(0),
 247   fEvSel_nTracklets(0),
 248   fEvSel_nSPDClusters(0),
 249
 250   //---> Variables for fTreeV0
 251         fTreeVariableChi2V0(0),
 252         fTreeVariableDcaV0Daughters(0),
 253         fTreeVariableDcaV0ToPrimVertex(0),
 254         fTreeVariableDcaPosToPrimVertex(0),
 255         fTreeVariableDcaNegToPrimVertex(0),
 256         fTreeVariableV0CosineOfPointingAngle(0),
 257         fTreeVariableV0Radius(0),
 258         fTreeVariablePt(0),
 259         fTreeVariableRapK0Short(0),
 260         fTreeVariableRapLambda(0),
 261         fTreeVariableInvMassK0s(0),
 262         fTreeVariableInvMassLambda(0),
 263         fTreeVariableInvMassAntiLambda(0),
 264         fTreeVariableAlphaV0(0),
 265         fTreeVariablePtArmV0(0),
 266         fTreeVariableNegEta(0),
 267         fTreeVariablePosEta(0),
 268
 269         fTreeVariableNSigmasPosProton(0),
 270         fTreeVariableNSigmasPosPion(0),
 271         fTreeVariableNSigmasNegProton(0),
 272         fTreeVariableNSigmasNegPion(0),
 273
 274         fTreeVariableDistOverTotMom(0),
 275         fTreeVariableLeastNbrCrossedRows(0),
 276         fTreeVariableLeastRatioCrossedRowsOverFindable(0),
 277
 278   fTreeVariableCentV0A(0),
 279   fTreeVariableCentV0C(0),
 280   fTreeVariableCentV0M(0),
 281   fTreeVariableCentV0AEq(0),
 282   fTreeVariableCentV0CEq(0),
 283   fTreeVariableCentV0MEq(0),
 284   fTreeVariableAmpV0A(0),
 285   fTreeVariableAmpV0C(0),
 286   fTreeVariableAmpV0AEq(0),
 287   fTreeVariableAmpV0CEq(0),
 288   fTreeVariableRefMultEta8(0),
 289   fTreeVariableRefMultEta5(0),
 290   fTreeVariableRunNumber(0),
 291   //---> Variables for fTreeCascade
 292   fTreeCascVarCharge(0),
 293         fTreeCascVarMassAsXi(0),
 294         fTreeCascVarMassAsOmega(0),
 295         fTreeCascVarPt(0),
 296         fTreeCascVarRapXi(0),
 297         fTreeCascVarRapOmega(0),
 298         fTreeCascVarNegEta(0),
 299         fTreeCascVarPosEta(0),
 300         fTreeCascVarBachEta(0),
 301         fTreeCascVarDCACascDaughters(0),
 302         fTreeCascVarDCABachToPrimVtx(0),
 303         fTreeCascVarDCAV0Daughters(0),
 304         fTreeCascVarDCAV0ToPrimVtx(0),
 305         fTreeCascVarDCAPosToPrimVtx(0),
 306         fTreeCascVarDCANegToPrimVtx(0),
 307         fTreeCascVarCascCosPointingAngle(0),
 308         fTreeCascVarCascRadius(0),
 309         fTreeCascVarV0Mass(0),
 310         fTreeCascVarV0CosPointingAngle(0),
 311         fTreeCascVarV0CosPointingAngleSpecial(0),
 312         fTreeCascVarV0Radius(0),
 313   fTreeCascVarLeastNbrClusters(0),
 314         fTreeCascVarDistOverTotMom(0),
 315         fTreeCascVarNegNSigmaPion(0),
 316         fTreeCascVarNegNSigmaProton(0),
 317         fTreeCascVarPosNSigmaPion(0),
 318         fTreeCascVarPosNSigmaProton(0),
 319         fTreeCascVarBachNSigmaPion(0),
 320         fTreeCascVarBachNSigmaKaon(0),
 321         fTreeCascVarCentV0A(0),
 322         fTreeCascVarCentV0C(0),
 323         fTreeCascVarCentV0M(0),
 324         fTreeCascVarCentV0AEq(0),
 325         fTreeCascVarCentV0CEq(0),
 326         fTreeCascVarCentV0MEq(0),
 327         fTreeCascVarAmpV0A(0),
 328         fTreeCascVarAmpV0C(0),
 329         fTreeCascVarAmpV0AEq(0),
 330         fTreeCascVarAmpV0CEq(0),
 331   fTreeCascVarRefMultEta8(0),
 332   fTreeCascVarRefMultEta5(0),
 333   fTreeCascVarRunNumber(0),
 334   //Histos
 335   fHistEventCounter(0)
 336 {
 337
 338   //Re-vertex: Will only apply for cascade candidates
 339
 340   fV0VertexerSels[0] =  33.  ;  // max allowed chi2
 341   fV0VertexerSels[1] =   0.02;  // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 342   fV0VertexerSels[2] =   0.02;  // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 343   fV0VertexerSels[4] =   0.95;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 344   fV0VertexerSels[5] =   1.0 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 345   fV0VertexerSels[6] = 200.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 346
 347   fCascadeVertexerSels[0] =  33.   ;  // max allowed chi2 (same as PDC07)
 348   fCascadeVertexerSels[1] =   0.05 ;  // min allowed V0 impact parameter                    (PDC07 : 0.05   / LHC09a4 : 0.025 )
 349   fCascadeVertexerSels[2] =   0.010;  // "window" around the Lambda mass                    (PDC07 : 0.008  / LHC09a4 : 0.010 )
 350   fCascadeVertexerSels[3] =   0.03 ;  // min allowed bachelor's impact parameter            (PDC07 : 0.035  / LHC09a4 : 0.025 )
 351   fCascadeVertexerSels[4] =   2.0  ;  // max allowed DCA between the V0 and the bachelor    (PDC07 : 0.1    / LHC09a4 : 0.2   )
 352   fCascadeVertexerSels[5] =   0.95 ;  // min allowed cosine of the cascade pointing angle   (PDC07 : 0.9985 / LHC09a4 : 0.998 )
 353   fCascadeVertexerSels[6] =   0.4  ;  // min radius of the fiducial volume                  (PDC07 : 0.9    / LHC09a4 : 0.2   )
 354   fCascadeVertexerSels[7] = 100.   ;  // max radius of the fiducial volume                  (PDC07 : 100    / LHC09a4 : 100   )
 355
 356
 357   DefineOutput(1, TList::Class()); // Event Counter Histo
 358   DefineOutput(2, TTree::Class()); // Event Tree
 359   DefineOutput(3, TTree::Class()); // V0 Tree
 360   DefineOutput(4, TTree::Class()); // Cascade Tree
 361 }
 362
 363
 364 AliAnalysisTaskStrangenessVsMultiplicity::~AliAnalysisTaskStrangenessVsMultiplicity()
 365 {
 366 //------------------------------------------------
 367 // DESTRUCTOR
 368 //------------------------------------------------
 369
 370    if (fListHist){
 371       delete fListHist;
 372       fListHist = 0x0;
 373    }
 374    if (fTreeEvent){
 375       delete fTreeEvent;
 376       fTreeEvent = 0x0;
 377    }
 378    if (fTreeV0){
 379       delete fTreeV0;
 380       fTreeV0 = 0x0;
 381    }
 382    if (fTreeCascade){
 383       delete fTreeCascade;
 384       fTreeCascade = 0x0;
 385    }
 386     if (fPPVsMultUtils){
 387         delete fPPVsMultUtils;
 388         fPPVsMultUtils = 0x0;
 389     }
 390 }
 391
 392 //________________________________________________________________________
 393 void AliAnalysisTaskStrangenessVsMultiplicity::UserCreateOutputObjects()
 394 {
 395
 396    OpenFile(2);
 397    // Called once
 398
 399 //------------------------------------------------
 400
 401    fTreeEvent = new TTree("fTreeEvent","Event");
 402
 403 //------------------------------------------------
 404 // fTree Branch definitions - Event by Event info
 405 //------------------------------------------------
 406
 407 //-----------BASIC-INFO---------------------------
 408
 409   //--- VZERO Data (Integrated)
 410   fTreeEvent->Branch("fAmplitude_V0A",&fAmplitude_V0A,"fAmplitude_V0A/F");
 411   fTreeEvent->Branch("fAmplitude_V0C",&fAmplitude_V0C,"fAmplitude_V0C/F");
 412   fTreeEvent->Branch("fAmplitude_V0AEq",&fAmplitude_V0AEq,"fAmplitude_V0AEq/F");
 413   fTreeEvent->Branch("fAmplitude_V0CEq",&fAmplitude_V0CEq,"fAmplitude_V0CEq/F");
 414
 415   //Info from AliCentrality (not necessarily 'centrality' per se)
 416   fTreeEvent->Branch("fCentrality_V0A",&fCentrality_V0A,"fCentrality_V0A/F");
 417   fTreeEvent->Branch("fCentrality_V0C",&fCentrality_V0C,"fCentrality_V0C/F");
 418   fTreeEvent->Branch("fCentrality_V0M",&fCentrality_V0M,"fCentrality_V0M/F");
 419   fTreeEvent->Branch("fCentrality_V0AEq",&fCentrality_V0AEq,"fCentrality_V0AEq/F");
 420   fTreeEvent->Branch("fCentrality_V0CEq",&fCentrality_V0CEq,"fCentrality_V0CEq/F");
 421   fTreeEvent->Branch("fCentrality_V0MEq",&fCentrality_V0MEq,"fCentrality_V0MEq/F");
 422   fTreeEvent->Branch("fCustomCentrality_V0M",&fCustomCentrality_V0M,"fCustomCentrality_V0M/F");
 423   fTreeEvent->Branch("fCustomCentrality_V0MEq",&fCustomCentrality_V0MEq,"fCustomCentrality_V0MEq/F");
 424
 425   //Official GetReferenceMultiplicity
 426   fTreeEvent->Branch("fRefMultEta5",&fRefMultEta5,"fRefMultEta5/I");
 427   fTreeEvent->Branch("fRefMultEta8",&fRefMultEta8,"fRefMultEta8/I");
 428
 429   //Run Number
 430   fTreeEvent->Branch("fRunNumber", &fRunNumber, "fRunNumber/I");
 431
 432   //Booleans for Event Selection
 433   fTreeEvent->Branch("fEvSel_HasAtLeastSPDVertex", &fEvSel_HasAtLeastSPDVertex, "fEvSel_HasAtLeastSPDVertex/O");
 434   fTreeEvent->Branch("fEvSel_VtxZCut", &fEvSel_VtxZCut, "fEvSel_VtxZCut/O");
 435   fTreeEvent->Branch("fEvSel_IsNotPileup", &fEvSel_IsNotPileup, "fEvSel_IsNotPileup/O");
 436   fTreeEvent->Branch("fEvSel_IsNotPileupInMultBins", &fEvSel_IsNotPileupInMultBins, "fEvSel_IsNotPileupInMultBins/O");
 437   fTreeEvent->Branch("fEvSel_HasVtxContributor", &fEvSel_HasVtxContributor, "fEvSel_HasVtxContributor/O");
 438
 439   //Tracklets vs clusters test
 440   fTreeEvent->Branch("fEvSel_nTracklets", &fEvSel_nTracklets, "fEvSel_nTracklets/I");
 441   fTreeEvent->Branch("fEvSel_nSPDClusters", &fEvSel_nSPDClusters, "fEvSel_nSPDClusters/I");
 442
 443   //Create Basic V0 Output Tree
 444   fTreeV0 = new TTree( "fTreeV0", "V0 Candidates");
 445
 446 //------------------------------------------------
 447 // fTreeV0 Branch definitions
 448 //------------------------------------------------
 449
 450 //-----------BASIC-INFO---------------------------
 451   fTreeV0->Branch("fTreeVariableChi2V0",&fTreeVariableChi2V0,"fTreeVariableChi2V0/F");
 452   fTreeV0->Branch("fTreeVariableDcaV0Daughters",&fTreeVariableDcaV0Daughters,"fTreeVariableDcaV0Daughters/F");
 453   fTreeV0->Branch("fTreeVariableDcaV0ToPrimVertex",&fTreeVariableDcaV0ToPrimVertex,"fTreeVariableDcaV0ToPrimVertex/F");
 454   fTreeV0->Branch("fTreeVariableDcaPosToPrimVertex",&fTreeVariableDcaPosToPrimVertex,"fTreeVariableDcaPosToPrimVertex/F");
 455   fTreeV0->Branch("fTreeVariableDcaNegToPrimVertex",&fTreeVariableDcaNegToPrimVertex,"fTreeVariableDcaNegToPrimVertex/F");
 456   fTreeV0->Branch("fTreeVariableV0Radius",&fTreeVariableV0Radius,"fTreeVariableV0Radius/F");
 457   fTreeV0->Branch("fTreeVariablePt",&fTreeVariablePt,"fTreeVariablePt/F");
 458   fTreeV0->Branch("fTreeVariableRapK0Short",&fTreeVariableRapK0Short,"fTreeVariableRapK0Short/F");
 459   fTreeV0->Branch("fTreeVariableRapLambda",&fTreeVariableRapLambda,"fTreeVariableRapLambda/F");
 460   fTreeV0->Branch("fTreeVariableInvMassK0s",&fTreeVariableInvMassK0s,"fTreeVariableInvMassK0s/F");
 461   fTreeV0->Branch("fTreeVariableInvMassLambda",&fTreeVariableInvMassLambda,"fTreeVariableInvMassLambda/F");
 462   fTreeV0->Branch("fTreeVariableInvMassAntiLambda",&fTreeVariableInvMassAntiLambda,"fTreeVariableInvMassAntiLambda/F");
 463   fTreeV0->Branch("fTreeVariableV0CosineOfPointingAngle",&fTreeVariableV0CosineOfPointingAngle,"fTreeVariableV0CosineOfPointingAngle/F");
 464   fTreeV0->Branch("fTreeVariableAlphaV0",&fTreeVariableAlphaV0,"fTreeVariableAlphaV0/F");
 465   fTreeV0->Branch("fTreeVariablePtArmV0",&fTreeVariablePtArmV0,"fTreeVariablePtArmV0/F");
 466   fTreeV0->Branch("fTreeVariableLeastNbrCrossedRows",&fTreeVariableLeastNbrCrossedRows,"fTreeVariableLeastNbrCrossedRows/I");
 467   fTreeV0->Branch("fTreeVariableLeastRatioCrossedRowsOverFindable",&fTreeVariableLeastRatioCrossedRowsOverFindable,"fTreeVariableLeastRatioCrossedRowsOverFindable/F");
 468   fTreeV0->Branch("fTreeVariableDistOverTotMom",&fTreeVariableDistOverTotMom,"fTreeVariableDistOverTotMom/F");
 469   fTreeV0->Branch("fTreeVariableNSigmasPosProton",&fTreeVariableNSigmasPosProton,"fTreeVariableNSigmasPosProton/F");
 470   fTreeV0->Branch("fTreeVariableNSigmasPosPion",&fTreeVariableNSigmasPosPion,"fTreeVariableNSigmasPosPion/F");
 471   fTreeV0->Branch("fTreeVariableNSigmasNegProton",&fTreeVariableNSigmasNegProton,"fTreeVariableNSigmasNegProton/F");
 472   fTreeV0->Branch("fTreeVariableNSigmasNegPion",&fTreeVariableNSigmasNegPion,"fTreeVariableNSigmasNegPion/F");
 473   fTreeV0->Branch("fTreeVariableNegEta",&fTreeVariableNegEta,"fTreeVariableNegEta/F");
 474   fTreeV0->Branch("fTreeVariablePosEta",&fTreeVariablePosEta,"fTreeVariablePosEta/F");
 475 //-----------MULTIPLICITY-INFO--------------------
 476   fTreeV0->Branch("fTreeVariableCentV0A",&fTreeVariableCentV0A,"fTreeVariableCentV0A/F");
 477   fTreeV0->Branch("fTreeVariableCentV0C",&fTreeVariableCentV0C,"fTreeVariableCentV0C/F");
 478   fTreeV0->Branch("fTreeVariableCentV0M",&fTreeVariableCentV0M,"fTreeVariableCentV0M/F");
 479   fTreeV0->Branch("fTreeVariableCentV0AEq",&fTreeVariableCentV0AEq,"fTreeVariableCentV0AEq/F");
 480   fTreeV0->Branch("fTreeVariableCentV0CEq",&fTreeVariableCentV0CEq,"fTreeVariableCentV0CEq/F");
 481   fTreeV0->Branch("fTreeVariableCentV0MEq",&fTreeVariableCentV0MEq,"fTreeVariableCentV0MEq/F");
 482   fTreeV0->Branch("fTreeVariableAmpV0A",&fTreeVariableAmpV0A,"fTreeVariableAmpV0A/F");
 483   fTreeV0->Branch("fTreeVariableAmpV0C",&fTreeVariableAmpV0C,"fTreeVariableAmpV0C/F");
 484   fTreeV0->Branch("fTreeVariableAmpV0AEq",&fTreeVariableAmpV0AEq,"fTreeVariableAmpV0AEq/F");
 485   fTreeV0->Branch("fTreeVariableAmpV0CEq",&fTreeVariableAmpV0CEq,"fTreeVariableAmpV0CEq/F");
 486   fTreeV0->Branch("fTreeVariableRefMultEta8",&fTreeVariableRefMultEta8,"fTreeVariableRefMultEta8/I");
 487   fTreeV0->Branch("fTreeVariableRefMultEta5",&fTreeVariableRefMultEta5,"fTreeVariableRefMultEta5/I");
 488   fTreeV0->Branch("fTreeVariableRunNumber",&fTreeVariableRunNumber,"fTreeVariableRunNumber/I");
 489 //------------------------------------------------
 490
 491   //Create Cascade output tree
 492   fTreeCascade = new TTree("fTreeCascade","CascadeCandidates");
 493
 494 //------------------------------------------------
 495 // fTreeCascade Branch definitions - Cascade Tree
 496 //------------------------------------------------
 497
 498 //-----------BASIC-INFO---------------------------
 499         fTreeCascade->Branch("fTreeCascVarCharge",&fTreeCascVarCharge,"fTreeCascVarCharge/I");
 500   fTreeCascade->Branch("fTreeCascVarMassAsXi",&fTreeCascVarMassAsXi,"fTreeCascVarMassAsXi/F");
 501   fTreeCascade->Branch("fTreeCascVarMassAsOmega",&fTreeCascVarMassAsOmega,"fTreeCascVarMassAsOmega/F");
 502   fTreeCascade->Branch("fTreeCascVarPt",&fTreeCascVarPt,"fTreeCascVarPt/F");
 503   fTreeCascade->Branch("fTreeCascVarRapXi",&fTreeCascVarRapXi,"fTreeCascVarRapXi/F");
 504   fTreeCascade->Branch("fTreeCascVarRapOmega",&fTreeCascVarRapOmega,"fTreeCascVarRapOmega/F");
 505   fTreeCascade->Branch("fTreeCascVarNegEta",&fTreeCascVarNegEta,"fTreeCascVarNegEta/F");
 506   fTreeCascade->Branch("fTreeCascVarPosEta",&fTreeCascVarPosEta,"fTreeCascVarPosEta/F");
 507   fTreeCascade->Branch("fTreeCascVarBachEta",&fTreeCascVarBachEta,"fTreeCascVarBachEta/F");
 508 //-----------INFO-FOR-CUTS------------------------
 509   fTreeCascade->Branch("fTreeCascVarDCACascDaughters",&fTreeCascVarDCACascDaughters,"fTreeCascVarDCACascDaughters/F");
 510   fTreeCascade->Branch("fTreeCascVarDCABachToPrimVtx",&fTreeCascVarDCABachToPrimVtx,"fTreeCascVarDCABachToPrimVtx/F");
 511   fTreeCascade->Branch("fTreeCascVarDCAV0Daughters",&fTreeCascVarDCAV0Daughters,"fTreeCascVarDCAV0Daughters/F");
 512   fTreeCascade->Branch("fTreeCascVarDCAV0ToPrimVtx",&fTreeCascVarDCAV0ToPrimVtx,"fTreeCascVarDCAV0ToPrimVtx/F");
 513   fTreeCascade->Branch("fTreeCascVarDCAPosToPrimVtx",&fTreeCascVarDCAPosToPrimVtx,"fTreeCascVarDCAPosToPrimVtx/F");
 514   fTreeCascade->Branch("fTreeCascVarDCANegToPrimVtx",&fTreeCascVarDCANegToPrimVtx,"fTreeCascVarDCANegToPrimVtx/F");
 515   fTreeCascade->Branch("fTreeCascVarCascCosPointingAngle",&fTreeCascVarCascCosPointingAngle,"fTreeCascVarCascCosPointingAngle/F");
 516   fTreeCascade->Branch("fTreeCascVarCascRadius",&fTreeCascVarCascRadius,"fTreeCascVarCascRadius/F");
 517   fTreeCascade->Branch("fTreeCascVarV0Mass",&fTreeCascVarV0Mass,"fTreeCascVarV0Mass/F");
 518   fTreeCascade->Branch("fTreeCascVarV0CosPointingAngle",&fTreeCascVarV0CosPointingAngle,"fTreeCascVarV0CosPointingAngle/F");
 519   fTreeCascade->Branch("fTreeCascVarV0CosPointingAngleSpecial",&fTreeCascVarV0CosPointingAngleSpecial,"fTreeCascVarV0CosPointingAngleSpecial/F");
 520   fTreeCascade->Branch("fTreeCascVarV0Radius",&fTreeCascVarV0Radius,"fTreeCascVarV0Radius/F");
 521   fTreeCascade->Branch("fTreeCascVarLeastNbrClusters",&fTreeCascVarLeastNbrClusters,"fTreeCascVarLeastNbrClusters/I");
 522 //-----------MULTIPLICITY-INFO--------------------
 523   fTreeCascade->Branch("fTreeCascVarCentV0A",&fTreeCascVarCentV0A,"fTreeCascVarCentV0A/F");
 524   fTreeCascade->Branch("fTreeCascVarCentV0C",&fTreeCascVarCentV0C,"fTreeCascVarCentV0C/F");
 525   fTreeCascade->Branch("fTreeCascVarCentV0M",&fTreeCascVarCentV0M,"fTreeCascVarCentV0M/F");
 526   fTreeCascade->Branch("fTreeCascVarCentV0AEq",&fTreeCascVarCentV0AEq,"fTreeCascVarCentV0AEq/F");
 527   fTreeCascade->Branch("fTreeCascVarCentV0CEq",&fTreeCascVarCentV0CEq,"fTreeCascVarCentV0CEq/F");
 528   fTreeCascade->Branch("fTreeCascVarCentV0MEq",&fTreeCascVarCentV0MEq,"fTreeCascVarCentV0MEq/F");
 529   fTreeCascade->Branch("fTreeCascVarAmpV0A",&fTreeCascVarAmpV0A,"fTreeCascVarAmpV0A/F");
 530   fTreeCascade->Branch("fTreeCascVarAmpV0C",&fTreeCascVarAmpV0C,"fTreeCascVarAmpV0C/F");
 531   fTreeCascade->Branch("fTreeCascVarAmpV0AEq",&fTreeCascVarAmpV0AEq,"fTreeCascVarAmpV0AEq/F");
 532   fTreeCascade->Branch("fTreeCascVarAmpV0CEq",&fTreeCascVarAmpV0CEq,"fTreeCascVarAmpV0CEq/F");
 533   fTreeCascade->Branch("fTreeCascVarRefMultEta8",&fTreeCascVarRefMultEta8,"fTreeCascVarRefMultEta8/I");
 534   fTreeCascade->Branch("fTreeCascVarRefMultEta5",&fTreeCascVarRefMultEta5,"fTreeCascVarRefMultEta5/I");
 535   fTreeCascade->Branch("fTreeCascVarRunNumber",&fTreeCascVarRunNumber,"fTreeCascVarRunNumber/I");
 536 //-----------DECAY-LENGTH-INFO--------------------
 537   fTreeCascade->Branch("fTreeCascVarDistOverTotMom",&fTreeCascVarDistOverTotMom,"fTreeCascVarDistOverTotMom/F");
 538 //------------------------------------------------
 539   fTreeCascade->Branch("fTreeCascVarNegNSigmaPion",&fTreeCascVarNegNSigmaPion,"fTreeCascVarNegNSigmaPion/F");
 540   fTreeCascade->Branch("fTreeCascVarNegNSigmaProton",&fTreeCascVarNegNSigmaProton,"fTreeCascVarNegNSigmaProton/F");
 541   fTreeCascade->Branch("fTreeCascVarPosNSigmaPion",&fTreeCascVarPosNSigmaPion,"fTreeCascVarPosNSigmaPion/F");
 542   fTreeCascade->Branch("fTreeCascVarPosNSigmaProton",&fTreeCascVarPosNSigmaProton,"fTreeCascVarPosNSigmaProton/F");
 543   fTreeCascade->Branch("fTreeCascVarBachNSigmaPion",&fTreeCascVarBachNSigmaPion,"fTreeCascVarBachNSigmaPion/F");
 544   fTreeCascade->Branch("fTreeCascVarBachNSigmaKaon",&fTreeCascVarBachNSigmaKaon,"fTreeCascVarBachNSigmaKaon/F");
 545
 546 //------------------------------------------------
 547 // Particle Identification Setup
 548 //------------------------------------------------
 549
 550    AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 551    AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 552    fPIDResponse = inputHandler->GetPIDResponse();
 553
 554   // Multiplicity
 555   if(! fESDtrackCuts ){
 556     fESDtrackCuts = new AliESDtrackCuts();
 557   }
 558     //Helper
 559     if(! fPPVsMultUtils ){
 560         fPPVsMultUtils = new AliPPVsMultUtils();
 561     }
 562
 563 //------------------------------------------------
 564 // V0 Multiplicity Histograms
 565 //------------------------------------------------
 566
 567    // Create histograms
 568    OpenFile(1);
 569    fListHist = new TList();
 570    fListHist->SetOwner();  // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
 571
 572    if(! fHistEventCounter ) {
 573     //Histogram Output: Event-by-Event
 574     fHistEventCounter = new TH1D( "fHistEventCounter", ";Evt. Sel. Step;Count",5,0,5);
 575     fHistEventCounter->GetXaxis()->SetBinLabel(1, "Processed");
 576     fHistEventCounter->GetXaxis()->SetBinLabel(2, "Phys-Sel");
 577     fHistEventCounter->GetXaxis()->SetBinLabel(3, "Has Vtx");
 578     fHistEventCounter->GetXaxis()->SetBinLabel(4, "Vtx |z|<10cm");
 579     fHistEventCounter->GetXaxis()->SetBinLabel(5, "Isn't Pileup");
 580     fListHist->Add(fHistEventCounter);
 581    }
 582
 583    //List of Histograms: Normal
 584    PostData(1, fListHist);
 585
 586    //TTree Object: Saved to base directory. Should cache to disk while saving.
 587    //(Important to avoid excessive memory usage, particularly when merging)
 588    PostData(2, fTreeEvent);
 589    PostData(3, fTreeV0);
 590    PostData(4, fTreeCascade);
 591
 592 }// end UserCreateOutputObjects
 593
 594
 595 //________________________________________________________________________
 596 void AliAnalysisTaskStrangenessVsMultiplicity::UserExec(Option_t *)
 597 {
 598   // Main loop
 599   // Called for each event
 600
 601    AliESDEvent *lESDevent = 0x0;
 602
 603   //Zero all booleans, etc
 604   fEvSel_HasAtLeastSPDVertex    = kFALSE;
 605   fEvSel_VtxZCut                = kFALSE;
 606   fEvSel_IsNotPileup            = kFALSE;
 607   fEvSel_IsNotPileupInMultBins  = kFALSE;
 608   fEvSel_HasVtxContributor      = kFALSE;
 609
 610   fEvSel_nTracklets   = -1;
 611   fEvSel_nSPDClusters = -1;
 612
 613   // Connect to the InputEvent
 614   // After these lines, we should have an ESD/AOD event + the number of V0s in it.
 615
 616    // Appropriate for ESD analysis!
 617
 618    lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
 619    if (!lESDevent) {
 620       AliWarning("ERROR: lESDevent not available \n");
 621       return;
 622    }
 623
 624   //Get VZERO Information for multiplicity later
 625   AliVVZERO* esdV0 = lESDevent->GetVZEROData();
 626   if (!esdV0) {
 627     AliError("AliVVZERO not available");
 628     return;
 629   }
 630
 631    fRunNumber = lESDevent->GetRunNumber();
 632
 633   Double_t lMagneticField = -10;
 634   lMagneticField = lESDevent->GetMagneticField( );
 635
 636 //------------------------------------------------
 637 // Variable Definition
 638 //------------------------------------------------
 639
 640 //------------------------------------------------
 641 // Physics Selection
 642 //------------------------------------------------
 643
 644   fHistEventCounter->Fill(0.5);
 645
 646   UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
 647   Bool_t isSelected = 0;
 648   isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
 649
 650   //Standard Min-Bias Selection - always do this!
 651   if ( ! isSelected ) {
 652     PostData(1, fListHist);
 653     PostData(2, fTreeEvent);
 654     PostData(3, fTreeV0);
 655     PostData(4, fTreeCascade);
 656     return;
 657   }
 658
 659   fHistEventCounter->Fill(1.5);
 660
 661   //------------------------------------------------
 662   // Primary Vertex Requirements Section:
 663   //  ---> pp: has vertex, |z|<10cm
 664   //------------------------------------------------
 665
 666   //classical Proton-proton like selection
 667   const AliESDVertex *lPrimaryBestESDVtx     = lESDevent->GetPrimaryVertex();
 668   const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
 669   const AliESDVertex *lPrimarySPDVtx         = lESDevent->GetPrimaryVertexSPD();
 670
 671   Double_t lBestPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
 672   lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
 673
 674   //Only accept if Tracking or SPD vertex is fine
 675   if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() && !fkSkipEventSelection ){
 676     AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
 677     PostData(1, fListHist);
 678     PostData(2, fTreeEvent);
 679     PostData(3, fTreeV0);
 680     PostData(4, fTreeCascade);
 681     return;
 682   }
 683
 684   if(! (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus()) ){
 685     //Passed selection!
 686     fEvSel_HasAtLeastSPDVertex = kTRUE;
 687   }
 688
 689   //Has SPD or Tracking Vertex
 690   fHistEventCounter -> Fill(2.5);
 691
 692   //Always do Primary Vertex Selection
 693   if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0 && !fkSkipEventSelection ) {
 694     AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
 695     PostData(1, fListHist);
 696     PostData(2, fTreeEvent);
 697     PostData(3, fTreeV0);
 698     PostData(4, fTreeCascade);
 699     return;
 700   }
 701
 702   if(TMath::Abs(lBestPrimaryVtxPos[2]) <= 10.0 ){
 703     //Passed selection!
 704     fEvSel_VtxZCut = kTRUE;
 705   }
 706
 707   //Fill Event selected counter
 708   fHistEventCounter -> Fill(3.5);
 709
 710   //------------------------------------------------
 711   // Check if this isn't pileup
 712   //------------------------------------------------
 713
 714   if(lESDevent->IsPileupFromSPDInMultBins() && !fkSkipEventSelection ){
 715     // minContributors=3, minZdist=0.8, nSigmaZdist=3., nSigmaDiamXY=2., nSigmaDiamZ=5.
 716     //-> see http://alisoft.cern.ch/viewvc/trunk/STEER/AliESDEvent.h?root=AliRoot&r1=41914&r2=42199&pathrev=42199
 717     AliWarning("Pb / Event tagged as pile-up by SPD... return !");
 718     PostData(1, fListHist);
 719     PostData(2, fTreeEvent);
 720     PostData(3, fTreeV0);
 721     PostData(4, fTreeCascade);
 722     return;
 723   }
 724
 725   if( !lESDevent->IsPileupFromSPD()           ) fEvSel_IsNotPileup           = kTRUE;
 726   if( !lESDevent->IsPileupFromSPDInMultBins() ) fEvSel_IsNotPileupInMultBins = kTRUE;
 727
 728   //Fill Event isn't pileup counter
 729   fHistEventCounter -> Fill(4.5);
 730
 731 //------------------------------------------------
 732 // Multiplicity Information Acquistion
 733 //------------------------------------------------
 734
 735   //Standard GetReferenceMultiplicity Estimator (0.5 and 0.8)
 736   fRefMultEta5 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.5);
 737   fRefMultEta8 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.8);
 738
 739   // VZERO PART
 740   Float_t  multV0A  = 0;            //  multiplicity from V0 reco side A
 741   Float_t  multV0C  = 0;            //  multiplicity from V0 reco side C
 742   Float_t  multV0AEq  = 0;          //  multiplicity from V0 reco side A
 743   Float_t  multV0CEq  = 0;          //  multiplicity from V0 reco side C
 744   Float_t  multV0ACorr  = 0;            //  multiplicity from V0 reco side A
 745   Float_t  multV0CCorr  = 0;            //  multiplicity from V0 reco side C
 746
 747   //Non-Equalized Signal: copy of multV0ACorr and multV0CCorr from AliCentralitySelectionTask
 748   //Getters for uncorrected multiplicity
 749   multV0A=esdV0->GetMTotV0A();
 750   multV0C=esdV0->GetMTotV0C();
 751
 752   //Get Z vertex position of SPD vertex (why not Tracking if available?)
 753   Float_t zvtx = lPrimarySPDVtx->GetZ();
 754
 755   //Acquire Corrected multV0A
 756   multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);
 757   multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);
 758
 759   //Copy to Event Tree for extra information
 760   fAmplitude_V0A = multV0ACorr;
 761   fAmplitude_V0C = multV0CCorr;
 762
 763   // Equalized signals // From AliCentralitySelectionTask // Updated
 764   for(Int_t iCh = 32; iCh < 64; ++iCh) {
 765     Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
 766     multV0AEq += mult;
 767   }
 768   for(Int_t iCh = 0; iCh < 32; ++iCh) {
 769     Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
 770     multV0CEq += mult;
 771   }
 772   fAmplitude_V0AEq = multV0AEq;
 773   fAmplitude_V0CEq = multV0CEq;
 774
 775   fCentrality_V0A   = -100;
 776   fCentrality_V0C   = -100;
 777   fCentrality_V0M   = -100;
 778   fCentrality_V0AEq = -100;
 779   fCentrality_V0CEq = -100;
 780   fCentrality_V0MEq = -100;
 781
 782   //AliCentrality... Check if working?
 783   AliCentrality* centrality;
 784   centrality = lESDevent->GetCentrality();
 785   if ( !(centrality->GetQuality()>1) ){
 786     fCentrality_V0A   = centrality->GetCentralityPercentile( "V0A"   );
 787     fCentrality_V0C   = centrality->GetCentralityPercentile( "V0C"   );
 788     fCentrality_V0M   = centrality->GetCentralityPercentile( "V0M"   );
 789     fCentrality_V0AEq = centrality->GetCentralityPercentile( "V0AEq" );
 790     fCentrality_V0CEq = centrality->GetCentralityPercentile( "V0CEq" );
 791     fCentrality_V0MEq = centrality->GetCentralityPercentile( "V0MEq" );
 792   }
 793
 794     //Cross-check/debug
 795     fCustomCentrality_V0M   = fPPVsMultUtils->GetMultiplicityPercentile(lESDevent,"V0M");
 796     fCustomCentrality_V0MEq = fPPVsMultUtils->GetMultiplicityPercentile(lESDevent,"V0MEq");
 797
 798   //Tracklets vs Clusters Exploratory data
 799   fEvSel_nTracklets     = lESDevent->GetMultiplicity()->GetNumberOfTracklets();
 800   fEvSel_nSPDClusters   = lESDevent->GetNumberOfITSClusters(0) + lESDevent->GetNumberOfITSClusters(1);
 801
 802   //Event-level fill
 803   fTreeEvent->Fill() ;
 804
 805   //STOP HERE if skipping event selections (no point in doing the rest...)
 806   if( fkSkipEventSelection ){
 807     PostData(1, fListHist);
 808     PostData(2, fTreeEvent);
 809     PostData(3, fTreeV0);
 810     PostData(4, fTreeCascade);
 811     return;
 812   }
 813
 814 //------------------------------------------------
 815
 816 //------------------------------------------------
 817 // Fill V0 Tree as needed
 818 //------------------------------------------------
 819
 820 //Variable definition
 821    Int_t    lOnFlyStatus = 0;// nv0sOn = 0, nv0sOff = 0;
 822    Double_t lChi2V0 = 0;
 823    Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
 824    Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
 825    Double_t lV0CosineOfPointingAngle = 0;
 826    Double_t lV0Radius = 0, lPt = 0;
 827    Double_t lRapK0Short = 0, lRapLambda = 0;
 828    Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
 829    Double_t lAlphaV0 = 0, lPtArmV0 = 0;
 830
 831    Double_t fMinV0Pt = 0;
 832    Double_t fMaxV0Pt = 100;
 833
 834    Int_t nv0s = 0;
 835    nv0s = lESDevent->GetNumberOfV0s();
 836
 837    for (Int_t iV0 = 0; iV0 < nv0s; iV0++) //extra-crazy test
 838    {// This is the begining of the V0 loop
 839       AliESDv0 *v0 = ((AliESDEvent*)lESDevent)->GetV0(iV0);
 840       if (!v0) continue;
 841
 842       Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]);
 843
 844       Double_t tV0mom[3];
 845       v0->GetPxPyPz( tV0mom[0],tV0mom[1],tV0mom[2] );
 846       Double_t lV0TotalMomentum = TMath::Sqrt(
 847       tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );
 848
 849       lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);
 850
 851       lPt = v0->Pt();
 852       lRapK0Short = v0->RapK0Short();
 853       lRapLambda  = v0->RapLambda();
 854       if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;
 855
 856       UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex());
 857       UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex());
 858
 859       Double_t lMomPos[3]; v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
 860       Double_t lMomNeg[3]; v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
 861
 862       AliESDtrack *pTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyPos);
 863       AliESDtrack *nTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyNeg);
 864       if (!pTrack || !nTrack) {
 865          Printf("ERROR: Could not retreive one of the daughter track");
 866          continue;
 867       }
 868
 869       //Daughter Eta for Eta selection, afterwards
 870       fTreeVariableNegEta = nTrack->Eta();
 871       fTreeVariablePosEta = pTrack->Eta();
 872
 873       // Filter like-sign V0 (next: add counter and distribution)
 874       if ( pTrack->GetSign() == nTrack->GetSign()){
 875          continue;
 876       }
 877
 878       //________________________________________________________________________
 879       // Track quality cuts
 880       Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);
 881       Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);
 882       fTreeVariableLeastNbrCrossedRows = (Int_t) lPosTrackCrossedRows;
 883       if( lNegTrackCrossedRows < fTreeVariableLeastNbrCrossedRows )
 884          fTreeVariableLeastNbrCrossedRows = (Int_t) lNegTrackCrossedRows;
 885
 886       // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
 887       if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
 888       if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
 889
 890
 891       if ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) ) continue;
 892
 893       //GetKinkIndex condition
 894       if( pTrack->GetKinkIndex(0)>0 || nTrack->GetKinkIndex(0)>0 ) continue;
 895
 896       //Findable clusters > 0 condition
 897       if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;
 898
 899       //Compute ratio Crossed Rows / Findable clusters
 900       //Note: above test avoids division by zero!
 901       Float_t lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF()));
 902       Float_t lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF()));
 903
 904       fTreeVariableLeastRatioCrossedRowsOverFindable = lPosTrackCrossedRowsOverFindable;
 905       if( lNegTrackCrossedRowsOverFindable < fTreeVariableLeastRatioCrossedRowsOverFindable )
 906          fTreeVariableLeastRatioCrossedRowsOverFindable = lNegTrackCrossedRowsOverFindable;
 907
 908       //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here
 909       if ( fTreeVariableLeastRatioCrossedRowsOverFindable < 0.8 ) continue;
 910
 911       //End track Quality Cuts
 912       //________________________________________________________________________
 913
 914       lDcaPosToPrimVertex = TMath::Abs(pTrack->GetD(lBestPrimaryVtxPos[0],
 915                                                         lBestPrimaryVtxPos[1],
 916                                                         lMagneticField) );
 917
 918       lDcaNegToPrimVertex = TMath::Abs(nTrack->GetD(lBestPrimaryVtxPos[0],
 919                                                         lBestPrimaryVtxPos[1],
 920                                                         lMagneticField) );
 921
 922       lOnFlyStatus = v0->GetOnFlyStatus();
 923       lChi2V0 = v0->GetChi2V0();
 924       lDcaV0Daughters = v0->GetDcaV0Daughters();
 925       lDcaV0ToPrimVertex = v0->GetD(lBestPrimaryVtxPos[0],lBestPrimaryVtxPos[1],lBestPrimaryVtxPos[2]);
 926       lV0CosineOfPointingAngle = v0->GetV0CosineOfPointingAngle(lBestPrimaryVtxPos[0],lBestPrimaryVtxPos[1],lBestPrimaryVtxPos[2]);
 927       fTreeVariableV0CosineOfPointingAngle=lV0CosineOfPointingAngle;
 928
 929       // Getting invariant mass infos directly from ESD
 930       v0->ChangeMassHypothesis(310);
 931       lInvMassK0s = v0->GetEffMass();
 932       v0->ChangeMassHypothesis(3122);
 933       lInvMassLambda = v0->GetEffMass();
 934       v0->ChangeMassHypothesis(-3122);
 935       lInvMassAntiLambda = v0->GetEffMass();
 936       lAlphaV0 = v0->AlphaV0();
 937       lPtArmV0 = v0->PtArmV0();
 938
 939       fTreeVariablePt = v0->Pt();
 940       fTreeVariableChi2V0 = lChi2V0;
 941       fTreeVariableDcaV0ToPrimVertex = lDcaV0ToPrimVertex;
 942       fTreeVariableDcaV0Daughters = lDcaV0Daughters;
 943       fTreeVariableV0CosineOfPointingAngle = lV0CosineOfPointingAngle;
 944       fTreeVariableV0Radius = lV0Radius;
 945       fTreeVariableDcaPosToPrimVertex = lDcaPosToPrimVertex;
 946       fTreeVariableDcaNegToPrimVertex = lDcaNegToPrimVertex;
 947       fTreeVariableInvMassK0s = lInvMassK0s;
 948       fTreeVariableInvMassLambda = lInvMassLambda;
 949       fTreeVariableInvMassAntiLambda = lInvMassAntiLambda;
 950       fTreeVariableRapK0Short = lRapK0Short;
 951       fTreeVariableRapLambda = lRapLambda;
 952       fTreeVariableAlphaV0 = lAlphaV0;
 953       fTreeVariablePtArmV0 = lPtArmV0;
 954
 955       //Official means of acquiring N-sigmas
 956       fTreeVariableNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );
 957       fTreeVariableNSigmasPosPion   = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );
 958       fTreeVariableNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );
 959       fTreeVariableNSigmasNegPion   = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );
 960
 961 //This requires an Invariant Mass Hypothesis afterwards
 962       fTreeVariableDistOverTotMom = TMath::Sqrt(
 963                                                 TMath::Power( tDecayVertexV0[0] - lBestPrimaryVtxPos[0] , 2) +
 964                                                 TMath::Power( tDecayVertexV0[1] - lBestPrimaryVtxPos[1] , 2) +
 965                                                 TMath::Power( tDecayVertexV0[2] - lBestPrimaryVtxPos[2] , 2)
 966                                         );
 967       fTreeVariableDistOverTotMom /= (lV0TotalMomentum+1e-10); //avoid division by zero, to be sure
 968
 969       //Copy Multiplicity information
 970       fTreeVariableCentV0A = fCentrality_V0A;
 971       fTreeVariableCentV0C = fCentrality_V0C;
 972       fTreeVariableCentV0M = fCentrality_V0M;
 973       fTreeVariableCentV0AEq = fCentrality_V0AEq;
 974       fTreeVariableCentV0CEq = fCentrality_V0CEq;
 975       fTreeVariableCentV0MEq = fCentrality_V0MEq;
 976       fTreeVariableAmpV0A = fAmplitude_V0A;
 977       fTreeVariableAmpV0C = fAmplitude_V0C;
 978       fTreeVariableAmpV0AEq = fAmplitude_V0AEq;
 979       fTreeVariableAmpV0CEq = fAmplitude_V0CEq;
 980       fTreeVariableRefMultEta8 = fRefMultEta8;
 981       fTreeVariableRefMultEta5 = fRefMultEta5;
 982       fTreeVariableRunNumber = fRunNumber;
 983
 984 //------------------------------------------------
 985 // Fill Tree!
 986 //------------------------------------------------
 987
 988      // The conditionals are meant to decrease excessive
 989      // memory usage!
 990
 991      //First Selection: Reject OnFly
 992      if( lOnFlyStatus == 0 ){
 993        //Second Selection: rough 20-sigma band, parametric.
 994        //K0Short: Enough to parametrize peak broadening with linear function.
 995        Double_t lUpperLimitK0Short = (5.63707e-01) + (1.14979e-02)*fTreeVariablePt;
 996        Double_t lLowerLimitK0Short = (4.30006e-01) - (1.10029e-02)*fTreeVariablePt;
 997        //Lambda: Linear (for higher pt) plus exponential (for low-pt broadening)
 998        //[0]+[1]*x+[2]*TMath::Exp(-[3]*x)
 999        Double_t lUpperLimitLambda = (1.13688e+00) + (5.27838e-03)*fTreeVariablePt + (8.42220e-02)*TMath::Exp(-(3.80595e+00)*fTreeVariablePt);
1000        Double_t lLowerLimitLambda = (1.09501e+00) - (5.23272e-03)*fTreeVariablePt - (7.52690e-02)*TMath::Exp(-(3.46339e+00)*fTreeVariablePt);
1001        //Do Selection
1002        if( (fTreeVariableInvMassLambda    < lUpperLimitLambda  && fTreeVariableInvMassLambda     > lLowerLimitLambda     ) ||
1003           (fTreeVariableInvMassAntiLambda < lUpperLimitLambda  && fTreeVariableInvMassAntiLambda > lLowerLimitLambda     ) ||
1004           (fTreeVariableInvMassK0s        < lUpperLimitK0Short && fTreeVariableInvMassK0s        > lLowerLimitK0Short    ) ){
1005          //Pre-selection in case this is AA...
1006          if ( TMath::Abs(fTreeVariableNegEta)<0.8 && TMath::Abs(fTreeVariablePosEta)<0.8 && fkSaveV0Tree ) fTreeV0->Fill();
1007        }
1008      }
1009    }// This is the end of the V0 loop
1010
1011 //------------------------------------------------
1012 // Fill V0 tree over.
1013 //------------------------------------------------
1014
1015
1016
1017 //------------------------------------------------
1018 // Rerun cascade vertexer!
1019 //------------------------------------------------
1020
1021   if( fkRunVertexers ){
1022     lESDevent->ResetCascades();
1023     lESDevent->ResetV0s();
1024
1025     AliV0vertexer lV0vtxer;
1026     AliCascadeVertexer lCascVtxer;
1027
1028     lV0vtxer.SetDefaultCuts(fV0VertexerSels);
1029     lCascVtxer.SetDefaultCuts(fCascadeVertexerSels);
1030
1031     lV0vtxer.Tracks2V0vertices(lESDevent);
1032     lCascVtxer.V0sTracks2CascadeVertices(lESDevent);
1033   }
1034
1035 //------------------------------------------------
1036 // MAIN CASCADE LOOP STARTS HERE
1037 //------------------------------------------------
1038 // Code Credit: Antonin Maire (thanks^100)
1039 // ---> This is an adaptation
1040
1041   Long_t ncascades = 0;
1042         ncascades = lESDevent->GetNumberOfCascades();
1043
1044   for (Int_t iXi = 0; iXi < ncascades; iXi++){
1045     //------------------------------------------------
1046     // Initializations
1047     //------------------------------------------------
1048           //Double_t lTrkgPrimaryVtxRadius3D = -500.0;
1049           //Double_t lBestPrimaryVtxRadius3D = -500.0;
1050
1051           // - 1st part of initialisation : variables needed to store AliESDCascade data members
1052           Double_t lEffMassXi      = 0. ;
1053           //Double_t lChi2Xi         = -1. ;
1054           Double_t lDcaXiDaughters = -1. ;
1055           Double_t lXiCosineOfPointingAngle = -1. ;
1056           Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 };
1057           Double_t lXiRadius = -1000. ;
1058
1059           // - 2nd part of initialisation : Nbr of clusters within TPC for the 3 daughter cascade tracks
1060           Int_t    lPosTPCClusters    = -1; // For ESD only ...//FIXME : wait for availability in AOD
1061           Int_t    lNegTPCClusters    = -1; // For ESD only ...
1062           Int_t    lBachTPCClusters   = -1; // For ESD only ...
1063
1064           // - 3rd part of initialisation : about V0 part in cascades
1065           Double_t lInvMassLambdaAsCascDghter = 0.;
1066           //Double_t lV0Chi2Xi         = -1. ;
1067           Double_t lDcaV0DaughtersXi = -1.;
1068
1069           Double_t lDcaBachToPrimVertexXi = -1., lDcaV0ToPrimVertexXi = -1.;
1070           Double_t lDcaPosToPrimVertexXi  = -1.;
1071           Double_t lDcaNegToPrimVertexXi  = -1.;
1072           Double_t lV0CosineOfPointingAngleXi = -1. ;
1073           Double_t lV0CosineOfPointingAngleXiSpecial = -1. ;
1074           Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. }; // Position of VO coming from cascade
1075           Double_t lV0RadiusXi = -1000.0;
1076           Double_t lV0quality  = 0.;
1077
1078           // - 4th part of initialisation : Effective masses
1079           Double_t lInvMassXiMinus    = 0.;
1080           Double_t lInvMassXiPlus     = 0.;
1081           Double_t lInvMassOmegaMinus = 0.;
1082           Double_t lInvMassOmegaPlus  = 0.;
1083
1084           // - 6th part of initialisation : extra info for QA
1085           Double_t lXiMomX       = 0. , lXiMomY = 0., lXiMomZ = 0.;
1086           Double_t lXiTransvMom  = 0. ;
1087           //Double_t lXiTransvMomMC= 0. ;
1088           Double_t lXiTotMom     = 0. ;
1089
1090           Double_t lBachMomX       = 0., lBachMomY  = 0., lBachMomZ   = 0.;
1091           //Double_t lBachTransvMom  = 0.;
1092           //Double_t lBachTotMom     = 0.;
1093
1094     fTreeCascVarNegNSigmaPion   = -100;
1095     fTreeCascVarNegNSigmaProton = -100;
1096     fTreeCascVarPosNSigmaPion   = -100;
1097     fTreeCascVarPosNSigmaProton = -100;
1098     fTreeCascVarBachNSigmaPion  = -100;
1099     fTreeCascVarBachNSigmaKaon  = -100;
1100
1101           Short_t  lChargeXi = -2;
1102           //Double_t lV0toXiCosineOfPointingAngle = 0. ;
1103
1104           Double_t lRapXi   = -20.0, lRapOmega = -20.0; //  lEta = -20.0, lTheta = 360., lPhi = 720. ;
1105           //Double_t lAlphaXi = -200., lPtArmXi  = -200.0;
1106
1107     // -------------------------------------
1108     // II.ESD - Calculation Part dedicated to Xi vertices (ESD)
1109
1110           AliESDcascade *xi = lESDevent->GetCascade(iXi);
1111           if (!xi) continue;
1112
1113                 // - II.Step 2 : Assigning the necessary variables for specific AliESDcascade data members (ESD)
1114                 //-------------
1115           lV0quality = 0.;
1116           xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
1117
1118           lEffMassXi                    = xi->GetEffMassXi();
1119           //lChi2Xi                         = xi->GetChi2Xi();
1120           lDcaXiDaughters       = xi->GetDcaXiDaughters();
1121           lXiCosineOfPointingAngle                  = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0],
1122                                                                                  lBestPrimaryVtxPos[1],
1123                                                                                  lBestPrimaryVtxPos[2] );
1124                   // Take care : the best available vertex should be used (like in AliCascadeVertexer)
1125
1126           xi->GetXYZcascade( lPosXi[0],  lPosXi[1], lPosXi[2] );
1127           lXiRadius                     = TMath::Sqrt( lPosXi[0]*lPosXi[0]  +  lPosXi[1]*lPosXi[1] );
1128
1129                 // - II.Step 3 : around the tracks : Bach + V0 (ESD)
1130                 // ~ Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
1131                 //-------------
1132
1133         UInt_t lIdxPosXi        = (UInt_t) TMath::Abs( xi->GetPindex() );
1134         UInt_t lIdxNegXi        = (UInt_t) TMath::Abs( xi->GetNindex() );
1135         UInt_t lBachIdx         = (UInt_t) TMath::Abs( xi->GetBindex() );
1136                 // Care track label can be negative in MC production (linked with the track quality)
1137                 // However = normally, not the case for track index ...
1138
1139           // FIXME : rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
1140           if(lBachIdx == lIdxNegXi) {
1141                   AliWarning("Pb / Idx(Bach. track) = Idx(Neg. track) ... continue!"); continue;
1142           }
1143     if(lBachIdx == lIdxPosXi) {
1144         AliWarning("Pb / Idx(Bach. track) = Idx(Pos. track) ... continue!"); continue;
1145           }
1146
1147           AliESDtrack *pTrackXi         = lESDevent->GetTrack( lIdxPosXi );
1148           AliESDtrack *nTrackXi         = lESDevent->GetTrack( lIdxNegXi );
1149           AliESDtrack *bachTrackXi      = lESDevent->GetTrack( lBachIdx );
1150
1151           if (!pTrackXi || !nTrackXi || !bachTrackXi ) {
1152                   AliWarning("ERROR: Could not retrieve one of the 3 ESD daughter tracks of the cascade ...");
1153                   continue;
1154           }
1155
1156       fTreeCascVarPosEta = pTrackXi->Eta();
1157       fTreeCascVarNegEta = nTrackXi->Eta();
1158       fTreeCascVarBachEta = bachTrackXi->Eta();
1159
1160       Double_t lBMom[3], lNMom[3], lPMom[3];
1161       xi->GetBPxPyPz( lBMom[0], lBMom[1], lBMom[2] );
1162       xi->GetPPxPyPz( lPMom[0], lPMom[1], lPMom[2] );
1163       xi->GetNPxPyPz( lNMom[0], lNMom[1], lNMom[2] );
1164
1165       //fTreeCascVarBachTransMom = TMath::Sqrt( lBMom[0]*lBMom[0] + lBMom[1]*lBMom[1] );
1166       //fTreeCascVarPosTransMom  = TMath::Sqrt( lPMom[0]*lPMom[0] + lPMom[1]*lPMom[1] );
1167       //fTreeCascVarNegTransMom  = TMath::Sqrt( lNMom[0]*lNMom[0] + lNMom[1]*lNMom[1] );
1168
1169     //------------------------------------------------
1170     // TPC dEdx information
1171     //------------------------------------------------
1172     fTreeCascVarNegNSigmaPion   = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kPion   );
1173     fTreeCascVarNegNSigmaProton = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kProton );
1174     fTreeCascVarPosNSigmaPion   = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kPion );
1175     fTreeCascVarPosNSigmaProton = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kProton );
1176     fTreeCascVarBachNSigmaPion  = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kPion );
1177     fTreeCascVarBachNSigmaKaon  = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kKaon );
1178
1179     //------------------------------------------------
1180     // TPC Number of clusters info
1181     // --- modified to save the smallest number
1182     // --- of TPC clusters for the 3 tracks
1183     //------------------------------------------------
1184
1185           lPosTPCClusters   = pTrackXi->GetTPCNcls();
1186           lNegTPCClusters   = nTrackXi->GetTPCNcls();
1187           lBachTPCClusters  = bachTrackXi->GetTPCNcls();
1188
1189     // 1 - Poor quality related to TPCrefit
1190           ULong_t pStatus    = pTrackXi->GetStatus();
1191           ULong_t nStatus    = nTrackXi->GetStatus();
1192           ULong_t bachStatus = bachTrackXi->GetStatus();
1193
1194     //fTreeCascVarkITSRefitBachelor = kTRUE;
1195     //fTreeCascVarkITSRefitNegative = kTRUE;
1196     //fTreeCascVarkITSRefitPositive = kTRUE;
1197
1198     if ((pStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning("Pb / V0 Pos. track has no TPCrefit ... continue!"); continue; }
1199     if ((nStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning("Pb / V0 Neg. track has no TPCrefit ... continue!"); continue; }
1200     if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / Bach.   track has no TPCrefit ... continue!"); continue; }
1201
1202           // 2 - Poor quality related to TPC clusters: lowest cut of 70 clusters
1203     if(lPosTPCClusters  < 70) { AliWarning("Pb / V0 Pos. track has less than 70 TPC clusters ... continue!"); continue; }
1204           if(lNegTPCClusters  < 70) { AliWarning("Pb / V0 Neg. track has less than 70 TPC clusters ... continue!"); continue; }
1205           if(lBachTPCClusters < 70) { AliWarning("Pb / Bach.   track has less than 70 TPC clusters ... continue!"); continue; }
1206           Int_t leastnumberofclusters = 1000;
1207           if( lPosTPCClusters < leastnumberofclusters ) leastnumberofclusters = lPosTPCClusters;
1208           if( lNegTPCClusters < leastnumberofclusters ) leastnumberofclusters = lNegTPCClusters;
1209           if( lBachTPCClusters < leastnumberofclusters ) leastnumberofclusters = lBachTPCClusters;
1210
1211           lInvMassLambdaAsCascDghter    = xi->GetEffMass();
1212           // This value shouldn't change, whatever the working hyp. is : Xi-, Xi+, Omega-, Omega+
1213           lDcaV0DaughtersXi             = xi->GetDcaV0Daughters();
1214           //lV0Chi2Xi                   = xi->GetChi2V0();
1215
1216           lV0CosineOfPointingAngleXi    = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0],
1217                                                                             lBestPrimaryVtxPos[1],
1218                                                                             lBestPrimaryVtxPos[2] );
1219     //Modification: V0 CosPA wrt to Cascade decay vertex
1220           lV0CosineOfPointingAngleXiSpecial     = xi->GetV0CosineOfPointingAngle( lPosXi[0],
1221                                                                             lPosXi[1],
1222                                                                             lPosXi[2] );
1223
1224           lDcaV0ToPrimVertexXi          = xi->GetD( lBestPrimaryVtxPos[0],
1225                                                       lBestPrimaryVtxPos[1],
1226                                                       lBestPrimaryVtxPos[2] );
1227
1228           lDcaBachToPrimVertexXi = TMath::Abs( bachTrackXi->GetD(       lBestPrimaryVtxPos[0],
1229                                                                 lBestPrimaryVtxPos[1],
1230                                                                 lMagneticField  ) );
1231                                           // Note : AliExternalTrackParam::GetD returns an algebraic value ...
1232
1233           xi->GetXYZ( lPosV0Xi[0],  lPosV0Xi[1], lPosV0Xi[2] );
1234           lV0RadiusXi           = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0]  +  lPosV0Xi[1]*lPosV0Xi[1] );
1235
1236           lDcaPosToPrimVertexXi         = TMath::Abs( pTrackXi  ->GetD( lBestPrimaryVtxPos[0],
1237                                                                 lBestPrimaryVtxPos[1],
1238                                                                 lMagneticField  )     );
1239
1240           lDcaNegToPrimVertexXi         = TMath::Abs( nTrackXi  ->GetD( lBestPrimaryVtxPos[0],
1241                                                                 lBestPrimaryVtxPos[1],
1242                                                                 lMagneticField  )     );
1243
1244           // - II.Step 4 : around effective masses (ESD)
1245           // ~ change mass hypotheses to cover all the possibilities :  Xi-/+, Omega -/+
1246
1247           if( bachTrackXi->Charge() < 0 )       {
1248                   lV0quality = 0.;
1249                   xi->ChangeMassHypothesis(lV0quality , 3312);
1250                           // Calculate the effective mass of the Xi- candidate.
1251                           // pdg code 3312 = Xi-
1252                   lInvMassXiMinus = xi->GetEffMassXi();
1253
1254                   lV0quality = 0.;
1255                   xi->ChangeMassHypothesis(lV0quality , 3334);
1256                           // Calculate the effective mass of the Xi- candidate.
1257                           // pdg code 3334 = Omega-
1258                   lInvMassOmegaMinus = xi->GetEffMassXi();
1259
1260                   lV0quality = 0.;
1261                   xi->ChangeMassHypothesis(lV0quality , 3312);  // Back to default hyp.
1262           }// end if negative bachelor
1263
1264
1265           if( bachTrackXi->Charge() >  0 ){
1266                   lV0quality = 0.;
1267                   xi->ChangeMassHypothesis(lV0quality , -3312);
1268                           // Calculate the effective mass of the Xi+ candidate.
1269                           // pdg code -3312 = Xi+
1270                   lInvMassXiPlus = xi->GetEffMassXi();
1271
1272                   lV0quality = 0.;
1273                   xi->ChangeMassHypothesis(lV0quality , -3334);
1274                           // Calculate the effective mass of the Xi+ candidate.
1275                           // pdg code -3334  = Omega+
1276                   lInvMassOmegaPlus = xi->GetEffMassXi();
1277
1278                   lV0quality = 0.;
1279                   xi->ChangeMassHypothesis(lV0quality , -3312);         // Back to "default" hyp.
1280           }// end if positive bachelor
1281                   // - II.Step 6 : extra info for QA (ESD)
1282                   // miscellaneous pieces of info that may help regarding data quality assessment.
1283                   //-------------
1284
1285           xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
1286                   lXiTransvMom          = TMath::Sqrt( lXiMomX*lXiMomX   + lXiMomY*lXiMomY );
1287                   lXiTotMom     = TMath::Sqrt( lXiMomX*lXiMomX   + lXiMomY*lXiMomY   + lXiMomZ*lXiMomZ );
1288
1289           xi->GetBPxPyPz(  lBachMomX,  lBachMomY,  lBachMomZ );
1290                   //lBachTransvMom  = TMath::Sqrt( lBachMomX*lBachMomX   + lBachMomY*lBachMomY );
1291                   //lBachTotMom         = TMath::Sqrt( lBachMomX*lBachMomX   + lBachMomY*lBachMomY  +  lBachMomZ*lBachMomZ  );
1292
1293           lChargeXi = xi->Charge();
1294
1295           //lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
1296
1297           lRapXi    = xi->RapXi();
1298           lRapOmega = xi->RapOmega();
1299           //lEta      = xi->Eta();
1300           //lTheta    = xi->Theta() *180.0/TMath::Pi();
1301           //lPhi      = xi->Phi()   *180.0/TMath::Pi();
1302           //lAlphaXi  = xi->AlphaXi();
1303           //lPtArmXi  = xi->PtArmXi();
1304
1305   //------------------------------------------------
1306   // Set Variables for adding to tree
1307   //------------------------------------------------
1308
1309           fTreeCascVarCharge    = lChargeXi;
1310           if(lInvMassXiMinus!=0)    fTreeCascVarMassAsXi = lInvMassXiMinus;
1311           if(lInvMassXiPlus!=0)     fTreeCascVarMassAsXi = lInvMassXiPlus;
1312           if(lInvMassOmegaMinus!=0) fTreeCascVarMassAsOmega = lInvMassOmegaMinus;
1313           if(lInvMassOmegaPlus!=0)  fTreeCascVarMassAsOmega = lInvMassOmegaPlus;
1314           fTreeCascVarPt = lXiTransvMom;
1315           fTreeCascVarRapXi = lRapXi ;
1316           fTreeCascVarRapOmega = lRapOmega ;
1317           fTreeCascVarDCACascDaughters = lDcaXiDaughters;
1318           fTreeCascVarDCABachToPrimVtx = lDcaBachToPrimVertexXi;
1319           fTreeCascVarDCAV0Daughters = lDcaV0DaughtersXi;
1320           fTreeCascVarDCAV0ToPrimVtx = lDcaV0ToPrimVertexXi;
1321           fTreeCascVarDCAPosToPrimVtx = lDcaPosToPrimVertexXi;
1322           fTreeCascVarDCANegToPrimVtx = lDcaNegToPrimVertexXi;
1323           fTreeCascVarCascCosPointingAngle = lXiCosineOfPointingAngle;
1324           fTreeCascVarCascRadius = lXiRadius;
1325           fTreeCascVarV0Mass = lInvMassLambdaAsCascDghter;
1326           fTreeCascVarV0CosPointingAngle = lV0CosineOfPointingAngleXi;
1327           fTreeCascVarV0CosPointingAngleSpecial = lV0CosineOfPointingAngleXiSpecial;
1328           fTreeCascVarV0Radius = lV0RadiusXi;
1329           fTreeCascVarLeastNbrClusters = leastnumberofclusters;
1330
1331           //Copy Multiplicity information
1332           fTreeCascVarCentV0A = fCentrality_V0A;
1333           fTreeCascVarCentV0C = fCentrality_V0C;
1334           fTreeCascVarCentV0M = fCentrality_V0M;
1335           fTreeCascVarCentV0AEq = fCentrality_V0AEq;
1336           fTreeCascVarCentV0CEq = fCentrality_V0CEq;
1337           fTreeCascVarCentV0MEq = fCentrality_V0MEq;
1338           fTreeCascVarAmpV0A = fAmplitude_V0A;
1339           fTreeCascVarAmpV0C = fAmplitude_V0C;
1340           fTreeCascVarAmpV0AEq = fAmplitude_V0AEq;
1341           fTreeCascVarAmpV0CEq = fAmplitude_V0CEq;
1342           fTreeCascVarRefMultEta8 = fRefMultEta8;
1343           fTreeCascVarRefMultEta5 = fRefMultEta5;
1344           fTreeCascVarRunNumber = fRunNumber;
1345
1346           fTreeCascVarDistOverTotMom = TMath::Sqrt(
1347                                                 TMath::Power( lPosXi[0] - lBestPrimaryVtxPos[0] , 2) +
1348                                                 TMath::Power( lPosXi[1] - lBestPrimaryVtxPos[1] , 2) +
1349                                                 TMath::Power( lPosXi[2] - lBestPrimaryVtxPos[2] , 2)
1350                                         );
1351           fTreeCascVarDistOverTotMom /= (lXiTotMom+1e-13);
1352
1353 //All vars not specified here: specified elsewhere!
1354
1355 //------------------------------------------------
1356 // Fill Tree!
1357 //------------------------------------------------
1358
1359 // The conditional is meant to decrease excessive
1360 // memory usage! Be careful when loosening the
1361 // cut!
1362
1363   //Xi    Mass window: 150MeV wide
1364   //Omega mass window: 150MeV wide
1365
1366   if( fkSaveCascadeTree && ( (fTreeCascVarMassAsXi<1.32+0.075&&fTreeCascVarMassAsXi>1.32-0.075) ||
1367       (fTreeCascVarMassAsOmega<1.68+0.075&&fTreeCascVarMassAsOmega>1.68-0.075) ) ){
1368       fTreeCascade->Fill();
1369   }
1370
1371 //------------------------------------------------
1372 // Fill tree over.
1373 //------------------------------------------------
1374
1375         }// end of the Cascade loop (ESD or AOD)
1376
1377   // Post output data.
1378   PostData(1, fListHist);
1379   PostData(2, fTreeEvent);
1380   PostData(3, fTreeV0);
1381   PostData(4, fTreeCascade);
1382 }
1383
1384 //________________________________________________________________________
1385 void AliAnalysisTaskStrangenessVsMultiplicity::Terminate(Option_t *)
1386 {
1387    // Draw result to the screen
1388    // Called once at the end of the query
1389
1390    TList *cRetrievedList = 0x0;
1391    cRetrievedList = (TList*)GetOutputData(1);
1392    if(!cRetrievedList){
1393       Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : ouput data container list not available\n");
1394       return;
1395    }
1396
1397    fHistEventCounter = dynamic_cast<TH1D*> (  cRetrievedList->FindObject("fHistEventCounter")  );
1398    if (!fHistEventCounter) {
1399       Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : fHistEventCounter not available");
1400       return;
1401    }
1402
1403    TCanvas *canCheck = new TCanvas("AliAnalysisTaskStrangenessVsMultiplicity","V0 Multiplicity",10,10,510,510);
1404    canCheck->cd(1)->SetLogy();
1405
1406    fHistEventCounter->SetMarkerStyle(22);
1407    fHistEventCounter->DrawCopy("E");
1408 }
1409
1410 //----------------------------------------------------------------------------
1411
1412 Double_t AliAnalysisTaskStrangenessVsMultiplicity::MyRapidity(Double_t rE, Double_t rPz) const
1413 {
1414    // Local calculation for rapidity
1415    Double_t ReturnValue = -100;
1416    if( (rE-rPz+1.e-13) != 0 && (rE+rPz) != 0 ){
1417       ReturnValue =  0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
1418    }
1419    return ReturnValue;
1420 }