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
  81 #include "AliCFContainer.h"
  82 #include "AliMultiplicity.h"
  83 #include "AliAODMCParticle.h"
  84 #include "AliESDcascade.h"
  85 #include "AliAODcascade.h"
  86 #include "AliESDUtils.h"
  87 #include "AliGenEventHeader.h"
  88 #include "AliAnalysisTaskSE.h"
  89 #include "AliAnalysisUtils.h"
  90 #include "AliAnalysisTaskStrangenessVsMultiplicity.h"
  91
  92 using std::cout;
  93 using std::endl;
  94
  95 ClassImp(AliAnalysisTaskStrangenessVsMultiplicity)
  96
  97 AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity()
  98   : AliAnalysisTaskSE(), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0),
  99   fkSaveV0Tree      ( kFALSE ),
 100   fkSaveCascadeTree ( kTRUE  ),
 101   fkRunVertexers    ( kTRUE  ),
 102   //---> Variables for fTreeEvent
 103   fAmplitude_V0A   (0),
 104   fAmplitude_V0C   (0),
 105   fAmplitude_V0AEq (0),
 106   fAmplitude_V0CEq (0),
 107   fCentrality_V0A(0),
 108   fCentrality_V0C(0),
 109   fCentrality_V0M(0),
 110   fCentrality_V0AEq(0),
 111   fCentrality_V0CEq(0),
 112   fCentrality_V0MEq(0),
 113   fRefMultEta5(0),
 114   fRefMultEta8(0),
 115   fRunNumber(0),
 116   //---> Variables for fTreeCascade
 117   fTreeCascVarCharge(0),
 118         fTreeCascVarMassAsXi(0),
 119         fTreeCascVarMassAsOmega(0),
 120         fTreeCascVarPt(0),
 121         fTreeCascVarRapXi(0),
 122         fTreeCascVarRapOmega(0),
 123         fTreeCascVarNegEta(0),
 124         fTreeCascVarPosEta(0),
 125         fTreeCascVarBachEta(0),
 126         fTreeCascVarDCACascDaughters(0),
 127         fTreeCascVarDCABachToPrimVtx(0),
 128         fTreeCascVarDCAV0Daughters(0),
 129         fTreeCascVarDCAV0ToPrimVtx(0),
 130         fTreeCascVarDCAPosToPrimVtx(0),
 131         fTreeCascVarDCANegToPrimVtx(0),
 132         fTreeCascVarCascCosPointingAngle(0),
 133         fTreeCascVarCascRadius(0),
 134         fTreeCascVarV0Mass(0),
 135         fTreeCascVarV0CosPointingAngle(0),
 136         fTreeCascVarV0CosPointingAngleSpecial(0),
 137         fTreeCascVarV0Radius(0),
 138   fTreeCascVarLeastNbrClusters(0),
 139         fTreeCascVarDistOverTotMom(0),
 140         fTreeCascVarNegNSigmaPion(0),
 141         fTreeCascVarNegNSigmaProton(0),
 142         fTreeCascVarPosNSigmaPion(0),
 143         fTreeCascVarPosNSigmaProton(0),
 144         fTreeCascVarBachNSigmaPion(0),
 145         fTreeCascVarBachNSigmaKaon(0),
 146         fTreeCascVarCentV0A(0),
 147         fTreeCascVarCentV0C(0),
 148         fTreeCascVarCentV0M(0),
 149         fTreeCascVarCentV0AEq(0),
 150         fTreeCascVarCentV0CEq(0),
 151         fTreeCascVarCentV0MEq(0),
 152         fTreeCascVarAmpV0A(0),
 153         fTreeCascVarAmpV0C(0),
 154         fTreeCascVarAmpV0AEq(0),
 155         fTreeCascVarAmpV0CEq(0),
 156   fTreeCascVarRefMultEta8(0),
 157   fTreeCascVarRefMultEta5(0),
 158   fTreeCascVarRunNumber(0),
 159   //Histos
 160   fHistEventCounter(0)
 161 //------------------------------------------------
 162 // Tree Variables
 163 {
 164
 165 }
 166
 167 AliAnalysisTaskStrangenessVsMultiplicity::AliAnalysisTaskStrangenessVsMultiplicity(const char *name)
 168   : AliAnalysisTaskSE(name), fListHist(0), fTreeEvent(0), fTreeV0(0), fTreeCascade(0), fPIDResponse(0), fESDtrackCuts(0),
 169   fkSaveV0Tree      ( kFALSE ),
 170   fkSaveCascadeTree ( kTRUE  ),
 171   fkRunVertexers    ( kTRUE  ),
 172   //---> Variables for fTreeEvent
 173   fAmplitude_V0A (0),
 174   fAmplitude_V0C (0),
 175   fAmplitude_V0AEq (0),
 176   fAmplitude_V0CEq (0),
 177   fCentrality_V0A(0),
 178   fCentrality_V0C(0),
 179   fCentrality_V0M(0),
 180   fCentrality_V0AEq(0),
 181   fCentrality_V0CEq(0),
 182   fCentrality_V0MEq(0),
 183   fRefMultEta5(0),
 184   fRefMultEta8(0),
 185   fRunNumber(0),
 186   //---> Variables for fTreeCascade
 187   fTreeCascVarCharge(0),
 188         fTreeCascVarMassAsXi(0),
 189         fTreeCascVarMassAsOmega(0),
 190         fTreeCascVarPt(0),
 191         fTreeCascVarRapXi(0),
 192         fTreeCascVarRapOmega(0),
 193         fTreeCascVarNegEta(0),
 194         fTreeCascVarPosEta(0),
 195         fTreeCascVarBachEta(0),
 196         fTreeCascVarDCACascDaughters(0),
 197         fTreeCascVarDCABachToPrimVtx(0),
 198         fTreeCascVarDCAV0Daughters(0),
 199         fTreeCascVarDCAV0ToPrimVtx(0),
 200         fTreeCascVarDCAPosToPrimVtx(0),
 201         fTreeCascVarDCANegToPrimVtx(0),
 202         fTreeCascVarCascCosPointingAngle(0),
 203         fTreeCascVarCascRadius(0),
 204         fTreeCascVarV0Mass(0),
 205         fTreeCascVarV0CosPointingAngle(0),
 206         fTreeCascVarV0CosPointingAngleSpecial(0),
 207         fTreeCascVarV0Radius(0),
 208   fTreeCascVarLeastNbrClusters(0),
 209         fTreeCascVarDistOverTotMom(0),
 210         fTreeCascVarNegNSigmaPion(0),
 211         fTreeCascVarNegNSigmaProton(0),
 212         fTreeCascVarPosNSigmaPion(0),
 213         fTreeCascVarPosNSigmaProton(0),
 214         fTreeCascVarBachNSigmaPion(0),
 215         fTreeCascVarBachNSigmaKaon(0),
 216         fTreeCascVarCentV0A(0),
 217         fTreeCascVarCentV0C(0),
 218         fTreeCascVarCentV0M(0),
 219         fTreeCascVarCentV0AEq(0),
 220         fTreeCascVarCentV0CEq(0),
 221         fTreeCascVarCentV0MEq(0),
 222         fTreeCascVarAmpV0A(0),
 223         fTreeCascVarAmpV0C(0),
 224         fTreeCascVarAmpV0AEq(0),
 225         fTreeCascVarAmpV0CEq(0),
 226   fTreeCascVarRefMultEta8(0),
 227   fTreeCascVarRefMultEta5(0),
 228   fTreeCascVarRunNumber(0),
 229   //Histos
 230   fHistEventCounter(0)
 231 {
 232
 233   //Re-vertex: Will only apply for cascade candidates
 234
 235   fV0VertexerSels[0] =  33.  ;  // max allowed chi2
 236   fV0VertexerSels[1] =   0.02;  // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 237   fV0VertexerSels[2] =   0.02;  // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 238   fV0VertexerSels[4] =   0.95;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 239   fV0VertexerSels[5] =   1.0 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 240   fV0VertexerSels[6] = 200.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 241
 242   fCascadeVertexerSels[0] =  33.   ;  // max allowed chi2 (same as PDC07)
 243   fCascadeVertexerSels[1] =   0.05 ;  // min allowed V0 impact parameter                    (PDC07 : 0.05   / LHC09a4 : 0.025 )
 244   fCascadeVertexerSels[2] =   0.010;  // "window" around the Lambda mass                    (PDC07 : 0.008  / LHC09a4 : 0.010 )
 245   fCascadeVertexerSels[3] =   0.03 ;  // min allowed bachelor's impact parameter            (PDC07 : 0.035  / LHC09a4 : 0.025 )
 246   fCascadeVertexerSels[4] =   2.0  ;  // max allowed DCA between the V0 and the bachelor    (PDC07 : 0.1    / LHC09a4 : 0.2   )
 247   fCascadeVertexerSels[5] =   0.95 ;  // min allowed cosine of the cascade pointing angle   (PDC07 : 0.9985 / LHC09a4 : 0.998 )
 248   fCascadeVertexerSels[6] =   0.4  ;  // min radius of the fiducial volume                  (PDC07 : 0.9    / LHC09a4 : 0.2   )
 249   fCascadeVertexerSels[7] = 100.   ;  // max radius of the fiducial volume                  (PDC07 : 100    / LHC09a4 : 100   )
 250
 251
 252   DefineOutput(1, TList::Class()); // Event Counter Histo
 253   DefineOutput(2, TTree::Class()); // Event Tree
 254   DefineOutput(3, TTree::Class()); // V0 Tree
 255   DefineOutput(4, TTree::Class()); // Cascade Tree
 256 }
 257
 258
 259 AliAnalysisTaskStrangenessVsMultiplicity::~AliAnalysisTaskStrangenessVsMultiplicity()
 260 {
 261 //------------------------------------------------
 262 // DESTRUCTOR
 263 //------------------------------------------------
 264
 265    if (fListHist){
 266       delete fListHist;
 267       fListHist = 0x0;
 268    }
 269    if (fTreeEvent){
 270       delete fTreeEvent;
 271       fTreeEvent = 0x0;
 272    }
 273    if (fTreeV0){
 274       delete fTreeV0;
 275       fTreeV0 = 0x0;
 276    }
 277    if (fTreeCascade){
 278       delete fTreeCascade;
 279       fTreeCascade = 0x0;
 280    }
 281 }
 282
 283 //________________________________________________________________________
 284 void AliAnalysisTaskStrangenessVsMultiplicity::UserCreateOutputObjects()
 285 {
 286
 287    OpenFile(2);
 288    // Called once
 289
 290 //------------------------------------------------
 291
 292    fTreeEvent = new TTree("fTreeEvent","Event");
 293
 294 //------------------------------------------------
 295 // fTree Branch definitions - Event by Event info
 296 //------------------------------------------------
 297
 298 //-----------BASIC-INFO---------------------------
 299
 300   //--- VZERO Data (Integrated)
 301   fTreeEvent->Branch("fAmplitude_V0A",&fAmplitude_V0A,"fAmplitude_V0A/F");
 302   fTreeEvent->Branch("fAmplitude_V0C",&fAmplitude_V0C,"fAmplitude_V0C/F");
 303   fTreeEvent->Branch("fAmplitude_V0AEq",&fAmplitude_V0AEq,"fAmplitude_V0AEq/F");
 304   fTreeEvent->Branch("fAmplitude_V0CEq",&fAmplitude_V0CEq,"fAmplitude_V0CEq/F");
 305
 306   //Info from AliCentrality (not necessarily 'centrality' per se)
 307   fTreeEvent->Branch("fCentrality_V0A",&fCentrality_V0A,"fCentrality_V0A/F");
 308   fTreeEvent->Branch("fCentrality_V0C",&fCentrality_V0C,"fCentrality_V0C/F");
 309   fTreeEvent->Branch("fCentrality_V0M",&fCentrality_V0A,"fCentrality_V0M/F");
 310   fTreeEvent->Branch("fCentrality_V0AEq",&fCentrality_V0AEq,"fCentrality_V0AEq/F");
 311   fTreeEvent->Branch("fCentrality_V0CEq",&fCentrality_V0CEq,"fCentrality_V0CEq/F");
 312   fTreeEvent->Branch("fCentrality_V0MEq",&fCentrality_V0AEq,"fCentrality_V0MEq/F");
 313
 314   //Official GetReferenceMultiplicity
 315   fTreeEvent->Branch("fRefMultEta5",&fRefMultEta5,"fRefMultEta5/I");
 316   fTreeEvent->Branch("fRefMultEta8",&fRefMultEta8,"fRefMultEta8/I");
 317
 318   //Run Number
 319   fTreeEvent->Branch("fRunNumber", &fRunNumber, "fRunNumber/I");
 320
 321   //Create Basic V0 Output Tree
 322   fTreeV0 = new TTree( "fTreeV0", "V0 Candidates");
 323
 324   //Create Cascade output tree
 325   fTreeCascade = new TTree("fTreeCascade","CascadeCandidates");
 326
 327 //------------------------------------------------
 328 // fTreeCascade Branch definitions - Cascade Tree
 329 //------------------------------------------------
 330
 331 //-----------BASIC-INFO---------------------------
 332         fTreeCascade->Branch("fTreeCascVarCharge",&fTreeCascVarCharge,"fTreeCascVarCharge/I");
 333   fTreeCascade->Branch("fTreeCascVarMassAsXi",&fTreeCascVarMassAsXi,"fTreeCascVarMassAsXi/F");
 334   fTreeCascade->Branch("fTreeCascVarMassAsOmega",&fTreeCascVarMassAsOmega,"fTreeCascVarMassAsOmega/F");
 335   fTreeCascade->Branch("fTreeCascVarPt",&fTreeCascVarPt,"fTreeCascVarPt/F");
 336   fTreeCascade->Branch("fTreeCascVarRapXi",&fTreeCascVarRapXi,"fTreeCascVarRapXi/F");
 337   fTreeCascade->Branch("fTreeCascVarRapOmega",&fTreeCascVarRapOmega,"fTreeCascVarRapOmega/F");
 338   fTreeCascade->Branch("fTreeCascVarNegEta",&fTreeCascVarNegEta,"fTreeCascVarNegEta/F");
 339   fTreeCascade->Branch("fTreeCascVarPosEta",&fTreeCascVarPosEta,"fTreeCascVarPosEta/F");
 340   fTreeCascade->Branch("fTreeCascVarBachEta",&fTreeCascVarBachEta,"fTreeCascVarBachEta/F");
 341 //-----------INFO-FOR-CUTS------------------------
 342   fTreeCascade->Branch("fTreeCascVarDCACascDaughters",&fTreeCascVarDCACascDaughters,"fTreeCascVarDCACascDaughters/F");
 343   fTreeCascade->Branch("fTreeCascVarDCABachToPrimVtx",&fTreeCascVarDCABachToPrimVtx,"fTreeCascVarDCABachToPrimVtx/F");
 344   fTreeCascade->Branch("fTreeCascVarDCAV0Daughters",&fTreeCascVarDCAV0Daughters,"fTreeCascVarDCAV0Daughters/F");
 345   fTreeCascade->Branch("fTreeCascVarDCAV0ToPrimVtx",&fTreeCascVarDCAV0ToPrimVtx,"fTreeCascVarDCAV0ToPrimVtx/F");
 346   fTreeCascade->Branch("fTreeCascVarDCAPosToPrimVtx",&fTreeCascVarDCAPosToPrimVtx,"fTreeCascVarDCAPosToPrimVtx/F");
 347   fTreeCascade->Branch("fTreeCascVarDCANegToPrimVtx",&fTreeCascVarDCANegToPrimVtx,"fTreeCascVarDCANegToPrimVtx/F");
 348   fTreeCascade->Branch("fTreeCascVarCascCosPointingAngle",&fTreeCascVarCascCosPointingAngle,"fTreeCascVarCascCosPointingAngle/F");
 349   fTreeCascade->Branch("fTreeCascVarCascRadius",&fTreeCascVarCascRadius,"fTreeCascVarCascRadius/F");
 350   fTreeCascade->Branch("fTreeCascVarV0Mass",&fTreeCascVarV0Mass,"fTreeCascVarV0Mass/F");
 351   fTreeCascade->Branch("fTreeCascVarV0CosPointingAngle",&fTreeCascVarV0CosPointingAngle,"fTreeCascVarV0CosPointingAngle/F");
 352   fTreeCascade->Branch("fTreeCascVarV0CosPointingAngleSpecial",&fTreeCascVarV0CosPointingAngleSpecial,"fTreeCascVarV0CosPointingAngleSpecial/F");
 353   fTreeCascade->Branch("fTreeCascVarV0Radius",&fTreeCascVarV0Radius,"fTreeCascVarV0Radius/F");
 354   fTreeCascade->Branch("fTreeCascVarLeastNbrClusters",&fTreeCascVarLeastNbrClusters,"fTreeCascVarLeastNbrClusters/I");
 355 //-----------MULTIPLICITY-INFO--------------------
 356   fTreeCascade->Branch("fTreeCascVarCentV0A",&fTreeCascVarCentV0A,"fTreeCascVarCentV0A/F");
 357   fTreeCascade->Branch("fTreeCascVarCentV0C",&fTreeCascVarCentV0C,"fTreeCascVarCentV0C/F");
 358   fTreeCascade->Branch("fTreeCascVarCentV0M",&fTreeCascVarCentV0M,"fTreeCascVarCentV0M/F");
 359   fTreeCascade->Branch("fTreeCascVarCentV0AEq",&fTreeCascVarCentV0AEq,"fTreeCascVarCentV0AEq/F");
 360   fTreeCascade->Branch("fTreeCascVarCentV0CEq",&fTreeCascVarCentV0CEq,"fTreeCascVarCentV0CEq/F");
 361   fTreeCascade->Branch("fTreeCascVarCentV0MEq",&fTreeCascVarCentV0MEq,"fTreeCascVarCentV0MEq/F");
 362   fTreeCascade->Branch("fTreeCascVarAmpV0A",&fTreeCascVarAmpV0A,"fTreeCascVarAmpV0A/F");
 363   fTreeCascade->Branch("fTreeCascVarAmpV0C",&fTreeCascVarAmpV0C,"fTreeCascVarAmpV0C/F");
 364   fTreeCascade->Branch("fTreeCascVarAmpV0AEq",&fTreeCascVarAmpV0AEq,"fTreeCascVarAmpV0AEq/F");
 365   fTreeCascade->Branch("fTreeCascVarAmpV0CEq",&fTreeCascVarAmpV0CEq,"fTreeCascVarAmpV0CEq/F");
 366   fTreeCascade->Branch("fTreeCascVarRefMultEta8",&fTreeCascVarRefMultEta8,"fTreeCascVarRefMultEta8/I");
 367   fTreeCascade->Branch("fTreeCascVarRefMultEta5",&fTreeCascVarRefMultEta5,"fTreeCascVarRefMultEta5/I");
 368   fTreeCascade->Branch("fTreeCascVarRunNumber",&fTreeCascVarRunNumber,"fTreeCascVarRunNumber/I");
 369 //-----------DECAY-LENGTH-INFO--------------------
 370   fTreeCascade->Branch("fTreeCascVarDistOverTotMom",&fTreeCascVarDistOverTotMom,"fTreeCascVarDistOverTotMom/F");
 371 //------------------------------------------------
 372   fTreeCascade->Branch("fTreeCascVarNegNSigmaPion",&fTreeCascVarNegNSigmaPion,"fTreeCascVarNegNSigmaPion/F");
 373   fTreeCascade->Branch("fTreeCascVarNegNSigmaProton",&fTreeCascVarNegNSigmaProton,"fTreeCascVarNegNSigmaProton/F");
 374   fTreeCascade->Branch("fTreeCascVarPosNSigmaPion",&fTreeCascVarPosNSigmaPion,"fTreeCascVarPosNSigmaPion/F");
 375   fTreeCascade->Branch("fTreeCascVarPosNSigmaProton",&fTreeCascVarPosNSigmaProton,"fTreeCascVarPosNSigmaProton/F");
 376   fTreeCascade->Branch("fTreeCascVarBachNSigmaPion",&fTreeCascVarBachNSigmaPion,"fTreeCascVarBachNSigmaPion/F");
 377   fTreeCascade->Branch("fTreeCascVarBachNSigmaKaon",&fTreeCascVarBachNSigmaKaon,"fTreeCascVarBachNSigmaKaon/F");
 378
 379 //------------------------------------------------
 380 // Particle Identification Setup
 381 //------------------------------------------------
 382
 383    AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 384    AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 385    fPIDResponse = inputHandler->GetPIDResponse();
 386
 387   // Multiplicity
 388   if(! fESDtrackCuts ){
 389     fESDtrackCuts = new AliESDtrackCuts();
 390   }
 391
 392 //------------------------------------------------
 393 // V0 Multiplicity Histograms
 394 //------------------------------------------------
 395
 396    // Create histograms
 397    OpenFile(1);
 398    fListHist = new TList();
 399    fListHist->SetOwner();  // See http://root.cern.ch/root/html/TCollection.html#TCollection:SetOwner
 400
 401    if(! fHistEventCounter ) {
 402     //Histogram Output: Event-by-Event
 403     fHistEventCounter = new TH1D( "fHistEventCounter", ";Evt. Sel. Step;Count",5,0,5);
 404     fHistEventCounter->GetXaxis()->SetBinLabel(1, "Processed");
 405     fHistEventCounter->GetXaxis()->SetBinLabel(2, "Phys-Sel");
 406     fHistEventCounter->GetXaxis()->SetBinLabel(3, "Has Vtx");
 407     fHistEventCounter->GetXaxis()->SetBinLabel(4, "Vtx |z|<10cm");
 408     fHistEventCounter->GetXaxis()->SetBinLabel(5, "Isn't Pileup");
 409     fListHist->Add(fHistEventCounter);
 410    }
 411
 412    //List of Histograms: Normal
 413    PostData(1, fListHist);
 414
 415    //TTree Object: Saved to base directory. Should cache to disk while saving.
 416    //(Important to avoid excessive memory usage, particularly when merging)
 417    PostData(2, fTreeEvent);
 418    PostData(3, fTreeV0);
 419    PostData(4, fTreeCascade);
 420
 421 }// end UserCreateOutputObjects
 422
 423
 424 //________________________________________________________________________
 425 void AliAnalysisTaskStrangenessVsMultiplicity::UserExec(Option_t *)
 426 {
 427   // Main loop
 428   // Called for each event
 429
 430    AliESDEvent *lESDevent = 0x0;
 431
 432   // Connect to the InputEvent
 433   // After these lines, we should have an ESD/AOD event + the number of V0s in it.
 434
 435    // Appropriate for ESD analysis!
 436
 437    lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
 438    if (!lESDevent) {
 439       AliWarning("ERROR: lESDevent not available \n");
 440       return;
 441    }
 442
 443   //Get VZERO Information for multiplicity later
 444   AliVVZERO* esdV0 = lESDevent->GetVZEROData();
 445   if (!esdV0) {
 446     AliError("AliVVZERO not available");
 447     return;
 448   }
 449
 450    fRunNumber = lESDevent->GetRunNumber();
 451
 452   Double_t lMagneticField = -10;
 453   lMagneticField = lESDevent->GetMagneticField( );
 454
 455 //------------------------------------------------
 456 // Variable Definition
 457 //------------------------------------------------
 458
 459 //------------------------------------------------
 460 // Physics Selection
 461 //------------------------------------------------
 462
 463   fHistEventCounter->Fill(0.5);
 464
 465   UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
 466   Bool_t isSelected = 0;
 467   isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
 468
 469   //Standard Min-Bias Selection
 470   if ( ! isSelected ) {
 471     PostData(1, fListHist);
 472     PostData(2, fTreeEvent);
 473     PostData(3, fTreeV0);
 474     PostData(4, fTreeCascade);
 475     return;
 476   }
 477
 478   fHistEventCounter->Fill(1.5);
 479
 480   //------------------------------------------------
 481   // Primary Vertex Requirements Section:
 482   //  ---> pp: has vertex, |z|<10cm
 483   //------------------------------------------------
 484
 485   //classical Proton-proton like selection
 486   const AliESDVertex *lPrimaryBestESDVtx     = lESDevent->GetPrimaryVertex();
 487   const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
 488   const AliESDVertex *lPrimarySPDVtx         = lESDevent->GetPrimaryVertexSPD();
 489
 490   Double_t lBestPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
 491   lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
 492
 493   //Only accept if Tracking or SPD vertex is fine
 494   if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtx->GetStatus() ){
 495     AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
 496     PostData(1, fListHist);
 497     PostData(2, fTreeEvent);
 498     PostData(3, fTreeV0);
 499     PostData(4, fTreeCascade);
 500     return;
 501   }
 502
 503   //Has SPD or Tracking Vertex
 504   fHistEventCounter -> Fill(2.5);
 505
 506   //Always do Primary Vertex Selection
 507   if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0) {
 508     AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
 509     PostData(1, fListHist);
 510     PostData(2, fTreeEvent);
 511     PostData(3, fTreeV0);
 512     PostData(4, fTreeCascade);
 513     return;
 514   }
 515
 516   //Fill Event selected counter
 517   fHistEventCounter -> Fill(3.5);
 518
 519   //------------------------------------------------
 520   // Check if this isn't pileup
 521   //------------------------------------------------
 522
 523   if(lESDevent->IsPileupFromSPDInMultBins() ){
 524     // minContributors=3, minZdist=0.8, nSigmaZdist=3., nSigmaDiamXY=2., nSigmaDiamZ=5.
 525     //-> see http://alisoft.cern.ch/viewvc/trunk/STEER/AliESDEvent.h?root=AliRoot&r1=41914&r2=42199&pathrev=42199
 526     AliWarning("Pb / Event tagged as pile-up by SPD... return !");
 527     PostData(1, fListHist);
 528     PostData(2, fTreeEvent);
 529     PostData(3, fTreeV0);
 530     PostData(4, fTreeCascade);
 531     return;
 532   }
 533   //Fill Event isn't pileup counter
 534   fHistEventCounter -> Fill(4.5);
 535
 536 //------------------------------------------------
 537 // Multiplicity Information Acquistion
 538 //------------------------------------------------
 539
 540   //Standard GetReferenceMultiplicity Estimator (0.5 and 0.8)
 541   fRefMultEta5 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.5);
 542   fRefMultEta8 = fESDtrackCuts->GetReferenceMultiplicity(lESDevent, AliESDtrackCuts::kTrackletsITSTPC,0.8);
 543
 544   // VZERO PART
 545   Float_t  multV0A  = 0;            //  multiplicity from V0 reco side A
 546   Float_t  multV0C  = 0;            //  multiplicity from V0 reco side C
 547   Float_t  multV0AEq  = 0;          //  multiplicity from V0 reco side A
 548   Float_t  multV0CEq  = 0;          //  multiplicity from V0 reco side C
 549   Float_t  multV0ACorr  = 0;            //  multiplicity from V0 reco side A
 550   Float_t  multV0CCorr  = 0;            //  multiplicity from V0 reco side C
 551
 552   //Non-Equalized Signal: copy of multV0ACorr and multV0CCorr from AliCentralitySelectionTask
 553   //Getters for uncorrected multiplicity
 554   multV0A=esdV0->GetMTotV0A();
 555   multV0C=esdV0->GetMTotV0C();
 556
 557   //Get Z vertex position of SPD vertex (why not Tracking if available?)
 558   Float_t zvtx = lPrimarySPDVtx->GetZ();
 559
 560   //Acquire Corrected multV0A
 561   multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);
 562   multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);
 563
 564   //Copy to Event Tree for extra information
 565   fAmplitude_V0A = multV0ACorr;
 566   fAmplitude_V0C = multV0CCorr;
 567
 568   // Equalized signals // From AliCentralitySelectionTask
 569   for(Int_t iCh = 4; iCh < 7; ++iCh) {
 570     Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
 571     multV0AEq += mult;
 572   }
 573   for(Int_t iCh = 0; iCh < 3; ++iCh) {
 574     Double_t mult = lESDevent->GetVZEROEqMultiplicity(iCh);
 575     multV0CEq += mult;
 576   }
 577   fAmplitude_V0AEq = multV0AEq;
 578   fAmplitude_V0CEq = multV0CEq;
 579
 580   fCentrality_V0A   = -100;
 581   fCentrality_V0C   = -100;
 582   fCentrality_V0M   = -100;
 583   fCentrality_V0AEq = -100;
 584   fCentrality_V0CEq = -100;
 585   fCentrality_V0MEq = -100;
 586
 587   //AliCentrality... Check if working?
 588   AliCentrality* centrality;
 589   centrality = lESDevent->GetCentrality();
 590   if ( !(centrality->GetQuality()>1) ){
 591     fCentrality_V0A   = centrality->GetCentralityPercentile( "V0A"   );
 592     fCentrality_V0C   = centrality->GetCentralityPercentile( "V0C"   );
 593     fCentrality_V0M   = centrality->GetCentralityPercentile( "V0M"   );
 594     fCentrality_V0AEq = centrality->GetCentralityPercentile( "V0AEq" );
 595     fCentrality_V0CEq = centrality->GetCentralityPercentile( "V0CEq" );
 596     fCentrality_V0MEq = centrality->GetCentralityPercentile( "V0MEq" );
 597   }
 598
 599   //Event-level fill
 600   fTreeEvent->Fill() ;
 601
 602 //------------------------------------------------
 603
 604 //------------------------------------------------
 605 // Rerun cascade vertexer!
 606 //------------------------------------------------
 607
 608   if( fkRunVertexers ){
 609     lESDevent->ResetCascades();
 610     lESDevent->ResetV0s();
 611
 612     AliV0vertexer lV0vtxer;
 613     AliCascadeVertexer lCascVtxer;
 614
 615     lV0vtxer.SetDefaultCuts(fV0VertexerSels);
 616     lCascVtxer.SetDefaultCuts(fCascadeVertexerSels);
 617
 618     lV0vtxer.Tracks2V0vertices(lESDevent);
 619     lCascVtxer.V0sTracks2CascadeVertices(lESDevent);
 620   }
 621
 622 //------------------------------------------------
 623 // MAIN CASCADE LOOP STARTS HERE
 624 //------------------------------------------------
 625 // Code Credit: Antonin Maire (thanks^100)
 626 // ---> This is an adaptation
 627
 628   Long_t ncascades = 0;
 629         ncascades = lESDevent->GetNumberOfCascades();
 630
 631   for (Int_t iXi = 0; iXi < ncascades; iXi++){
 632     //------------------------------------------------
 633     // Initializations
 634     //------------------------------------------------
 635           //Double_t lTrkgPrimaryVtxRadius3D = -500.0;
 636           //Double_t lBestPrimaryVtxRadius3D = -500.0;
 637
 638           // - 1st part of initialisation : variables needed to store AliESDCascade data members
 639           Double_t lEffMassXi      = 0. ;
 640           //Double_t lChi2Xi         = -1. ;
 641           Double_t lDcaXiDaughters = -1. ;
 642           Double_t lXiCosineOfPointingAngle = -1. ;
 643           Double_t lPosXi[3] = { -1000.0, -1000.0, -1000.0 };
 644           Double_t lXiRadius = -1000. ;
 645
 646           // - 2nd part of initialisation : Nbr of clusters within TPC for the 3 daughter cascade tracks
 647           Int_t    lPosTPCClusters    = -1; // For ESD only ...//FIXME : wait for availability in AOD
 648           Int_t    lNegTPCClusters    = -1; // For ESD only ...
 649           Int_t    lBachTPCClusters   = -1; // For ESD only ...
 650
 651           // - 3rd part of initialisation : about V0 part in cascades
 652           Double_t lInvMassLambdaAsCascDghter = 0.;
 653           //Double_t lV0Chi2Xi         = -1. ;
 654           Double_t lDcaV0DaughtersXi = -1.;
 655
 656           Double_t lDcaBachToPrimVertexXi = -1., lDcaV0ToPrimVertexXi = -1.;
 657           Double_t lDcaPosToPrimVertexXi  = -1.;
 658           Double_t lDcaNegToPrimVertexXi  = -1.;
 659           Double_t lV0CosineOfPointingAngleXi = -1. ;
 660           Double_t lV0CosineOfPointingAngleXiSpecial = -1. ;
 661           Double_t lPosV0Xi[3] = { -1000. , -1000., -1000. }; // Position of VO coming from cascade
 662           Double_t lV0RadiusXi = -1000.0;
 663           Double_t lV0quality  = 0.;
 664
 665           // - 4th part of initialisation : Effective masses
 666           Double_t lInvMassXiMinus    = 0.;
 667           Double_t lInvMassXiPlus     = 0.;
 668           Double_t lInvMassOmegaMinus = 0.;
 669           Double_t lInvMassOmegaPlus  = 0.;
 670
 671           // - 6th part of initialisation : extra info for QA
 672           Double_t lXiMomX       = 0. , lXiMomY = 0., lXiMomZ = 0.;
 673           Double_t lXiTransvMom  = 0. ;
 674           //Double_t lXiTransvMomMC= 0. ;
 675           Double_t lXiTotMom     = 0. ;
 676
 677           Double_t lBachMomX       = 0., lBachMomY  = 0., lBachMomZ   = 0.;
 678           //Double_t lBachTransvMom  = 0.;
 679           //Double_t lBachTotMom     = 0.;
 680
 681     fTreeCascVarNegNSigmaPion   = -100;
 682     fTreeCascVarNegNSigmaProton = -100;
 683     fTreeCascVarPosNSigmaPion   = -100;
 684     fTreeCascVarPosNSigmaProton = -100;
 685     fTreeCascVarBachNSigmaPion  = -100;
 686     fTreeCascVarBachNSigmaKaon  = -100;
 687
 688           Short_t  lChargeXi = -2;
 689           //Double_t lV0toXiCosineOfPointingAngle = 0. ;
 690
 691           Double_t lRapXi   = -20.0, lRapOmega = -20.0; //  lEta = -20.0, lTheta = 360., lPhi = 720. ;
 692           //Double_t lAlphaXi = -200., lPtArmXi  = -200.0;
 693
 694     // -------------------------------------
 695     // II.ESD - Calculation Part dedicated to Xi vertices (ESD)
 696
 697           AliESDcascade *xi = lESDevent->GetCascade(iXi);
 698           if (!xi) continue;
 699
 700                 // - II.Step 2 : Assigning the necessary variables for specific AliESDcascade data members (ESD)
 701                 //-------------
 702           lV0quality = 0.;
 703           xi->ChangeMassHypothesis(lV0quality , 3312); // default working hypothesis : cascade = Xi- decay
 704
 705           lEffMassXi                    = xi->GetEffMassXi();
 706           //lChi2Xi                         = xi->GetChi2Xi();
 707           lDcaXiDaughters       = xi->GetDcaXiDaughters();
 708           lXiCosineOfPointingAngle                  = xi->GetCascadeCosineOfPointingAngle( lBestPrimaryVtxPos[0],
 709                                                                                  lBestPrimaryVtxPos[1],
 710                                                                                  lBestPrimaryVtxPos[2] );
 711                   // Take care : the best available vertex should be used (like in AliCascadeVertexer)
 712
 713           xi->GetXYZcascade( lPosXi[0],  lPosXi[1], lPosXi[2] );
 714           lXiRadius                     = TMath::Sqrt( lPosXi[0]*lPosXi[0]  +  lPosXi[1]*lPosXi[1] );
 715
 716                 // - II.Step 3 : around the tracks : Bach + V0 (ESD)
 717                 // ~ Necessary variables for ESDcascade data members coming from the ESDv0 part (inheritance)
 718                 //-------------
 719
 720         UInt_t lIdxPosXi        = (UInt_t) TMath::Abs( xi->GetPindex() );
 721         UInt_t lIdxNegXi        = (UInt_t) TMath::Abs( xi->GetNindex() );
 722         UInt_t lBachIdx         = (UInt_t) TMath::Abs( xi->GetBindex() );
 723                 // Care track label can be negative in MC production (linked with the track quality)
 724                 // However = normally, not the case for track index ...
 725
 726           // FIXME : rejection of a double use of a daughter track (nothing but just a crosscheck of what is done in the cascade vertexer)
 727           if(lBachIdx == lIdxNegXi) {
 728                   AliWarning("Pb / Idx(Bach. track) = Idx(Neg. track) ... continue!"); continue;
 729           }
 730     if(lBachIdx == lIdxPosXi) {
 731         AliWarning("Pb / Idx(Bach. track) = Idx(Pos. track) ... continue!"); continue;
 732           }
 733
 734           AliESDtrack *pTrackXi         = lESDevent->GetTrack( lIdxPosXi );
 735           AliESDtrack *nTrackXi         = lESDevent->GetTrack( lIdxNegXi );
 736           AliESDtrack *bachTrackXi      = lESDevent->GetTrack( lBachIdx );
 737
 738           if (!pTrackXi || !nTrackXi || !bachTrackXi ) {
 739                   AliWarning("ERROR: Could not retrieve one of the 3 ESD daughter tracks of the cascade ...");
 740                   continue;
 741           }
 742
 743       fTreeCascVarPosEta = pTrackXi->Eta();
 744       fTreeCascVarNegEta = nTrackXi->Eta();
 745       fTreeCascVarBachEta = bachTrackXi->Eta();
 746
 747       Double_t lBMom[3], lNMom[3], lPMom[3];
 748       xi->GetBPxPyPz( lBMom[0], lBMom[1], lBMom[2] );
 749       xi->GetPPxPyPz( lPMom[0], lPMom[1], lPMom[2] );
 750       xi->GetNPxPyPz( lNMom[0], lNMom[1], lNMom[2] );
 751
 752       //fTreeCascVarBachTransMom = TMath::Sqrt( lBMom[0]*lBMom[0] + lBMom[1]*lBMom[1] );
 753       //fTreeCascVarPosTransMom  = TMath::Sqrt( lPMom[0]*lPMom[0] + lPMom[1]*lPMom[1] );
 754       //fTreeCascVarNegTransMom  = TMath::Sqrt( lNMom[0]*lNMom[0] + lNMom[1]*lNMom[1] );
 755
 756     //------------------------------------------------
 757     // TPC dEdx information
 758     //------------------------------------------------
 759     fTreeCascVarNegNSigmaPion   = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kPion   );
 760     fTreeCascVarNegNSigmaProton = fPIDResponse->NumberOfSigmasTPC( nTrackXi, AliPID::kProton );
 761     fTreeCascVarPosNSigmaPion   = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kPion );
 762     fTreeCascVarPosNSigmaProton = fPIDResponse->NumberOfSigmasTPC( pTrackXi, AliPID::kProton );
 763     fTreeCascVarBachNSigmaPion  = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kPion );
 764     fTreeCascVarBachNSigmaKaon  = fPIDResponse->NumberOfSigmasTPC( bachTrackXi, AliPID::kKaon );
 765
 766     //------------------------------------------------
 767     // TPC Number of clusters info
 768     // --- modified to save the smallest number
 769     // --- of TPC clusters for the 3 tracks
 770     //------------------------------------------------
 771
 772           lPosTPCClusters   = pTrackXi->GetTPCNcls();
 773           lNegTPCClusters   = nTrackXi->GetTPCNcls();
 774           lBachTPCClusters  = bachTrackXi->GetTPCNcls();
 775
 776     // 1 - Poor quality related to TPCrefit
 777           ULong_t pStatus    = pTrackXi->GetStatus();
 778           ULong_t nStatus    = nTrackXi->GetStatus();
 779           ULong_t bachStatus = bachTrackXi->GetStatus();
 780
 781     //fTreeCascVarkITSRefitBachelor = kTRUE;
 782     //fTreeCascVarkITSRefitNegative = kTRUE;
 783     //fTreeCascVarkITSRefitPositive = kTRUE;
 784
 785     if ((pStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning("Pb / V0 Pos. track has no TPCrefit ... continue!"); continue; }
 786     if ((nStatus&AliESDtrack::kTPCrefit)    == 0) { AliWarning("Pb / V0 Neg. track has no TPCrefit ... continue!"); continue; }
 787     if ((bachStatus&AliESDtrack::kTPCrefit) == 0) { AliWarning("Pb / Bach.   track has no TPCrefit ... continue!"); continue; }
 788
 789     //Extra Debug Information: booleans for ITS refit
 790     //if ((pStatus&AliESDtrack::kITSrefit)    == 0) { fTreeCascVarkITSRefitPositive = kFALSE; }
 791     //if ((nStatus&AliESDtrack::kITSrefit)    == 0) { fTreeCascVarkITSRefitNegative = kFALSE; }
 792     //if ((bachStatus&AliESDtrack::kITSrefit) == 0) { fTreeCascVarkITSRefitBachelor = kFALSE; }
 793
 794           // 2 - Poor quality related to TPC clusters: lowest cut of 70 clusters
 795     if(lPosTPCClusters  < 70) { AliWarning("Pb / V0 Pos. track has less than 70 TPC clusters ... continue!"); continue; }
 796           if(lNegTPCClusters  < 70) { AliWarning("Pb / V0 Neg. track has less than 70 TPC clusters ... continue!"); continue; }
 797           if(lBachTPCClusters < 70) { AliWarning("Pb / Bach.   track has less than 70 TPC clusters ... continue!"); continue; }
 798           Int_t leastnumberofclusters = 1000;
 799           if( lPosTPCClusters < leastnumberofclusters ) leastnumberofclusters = lPosTPCClusters;
 800           if( lNegTPCClusters < leastnumberofclusters ) leastnumberofclusters = lNegTPCClusters;
 801           if( lBachTPCClusters < leastnumberofclusters ) leastnumberofclusters = lBachTPCClusters;
 802
 803           lInvMassLambdaAsCascDghter    = xi->GetEffMass();
 804           // This value shouldn't change, whatever the working hyp. is : Xi-, Xi+, Omega-, Omega+
 805           lDcaV0DaughtersXi             = xi->GetDcaV0Daughters();
 806           //lV0Chi2Xi                   = xi->GetChi2V0();
 807
 808           lV0CosineOfPointingAngleXi    = xi->GetV0CosineOfPointingAngle( lBestPrimaryVtxPos[0],
 809                                                                             lBestPrimaryVtxPos[1],
 810                                                                             lBestPrimaryVtxPos[2] );
 811     //Modification: V0 CosPA wrt to Cascade decay vertex
 812           lV0CosineOfPointingAngleXiSpecial     = xi->GetV0CosineOfPointingAngle( lPosXi[0],
 813                                                                             lPosXi[1],
 814                                                                             lPosXi[2] );
 815
 816           lDcaV0ToPrimVertexXi          = xi->GetD( lBestPrimaryVtxPos[0],
 817                                                       lBestPrimaryVtxPos[1],
 818                                                       lBestPrimaryVtxPos[2] );
 819
 820           lDcaBachToPrimVertexXi = TMath::Abs( bachTrackXi->GetD(       lBestPrimaryVtxPos[0],
 821                                                                 lBestPrimaryVtxPos[1],
 822                                                                 lMagneticField  ) );
 823                                           // Note : AliExternalTrackParam::GetD returns an algebraic value ...
 824
 825           xi->GetXYZ( lPosV0Xi[0],  lPosV0Xi[1], lPosV0Xi[2] );
 826           lV0RadiusXi           = TMath::Sqrt( lPosV0Xi[0]*lPosV0Xi[0]  +  lPosV0Xi[1]*lPosV0Xi[1] );
 827
 828           lDcaPosToPrimVertexXi         = TMath::Abs( pTrackXi  ->GetD( lBestPrimaryVtxPos[0],
 829                                                                 lBestPrimaryVtxPos[1],
 830                                                                 lMagneticField  )     );
 831
 832           lDcaNegToPrimVertexXi         = TMath::Abs( nTrackXi  ->GetD( lBestPrimaryVtxPos[0],
 833                                                                 lBestPrimaryVtxPos[1],
 834                                                                 lMagneticField  )     );
 835
 836           // - II.Step 4 : around effective masses (ESD)
 837           // ~ change mass hypotheses to cover all the possibilities :  Xi-/+, Omega -/+
 838
 839           if( bachTrackXi->Charge() < 0 )       {
 840                   lV0quality = 0.;
 841                   xi->ChangeMassHypothesis(lV0quality , 3312);
 842                           // Calculate the effective mass of the Xi- candidate.
 843                           // pdg code 3312 = Xi-
 844                   lInvMassXiMinus = xi->GetEffMassXi();
 845
 846                   lV0quality = 0.;
 847                   xi->ChangeMassHypothesis(lV0quality , 3334);
 848                           // Calculate the effective mass of the Xi- candidate.
 849                           // pdg code 3334 = Omega-
 850                   lInvMassOmegaMinus = xi->GetEffMassXi();
 851
 852                   lV0quality = 0.;
 853                   xi->ChangeMassHypothesis(lV0quality , 3312);  // Back to default hyp.
 854           }// end if negative bachelor
 855
 856
 857           if( bachTrackXi->Charge() >  0 ){
 858                   lV0quality = 0.;
 859                   xi->ChangeMassHypothesis(lV0quality , -3312);
 860                           // Calculate the effective mass of the Xi+ candidate.
 861                           // pdg code -3312 = Xi+
 862                   lInvMassXiPlus = xi->GetEffMassXi();
 863
 864                   lV0quality = 0.;
 865                   xi->ChangeMassHypothesis(lV0quality , -3334);
 866                           // Calculate the effective mass of the Xi+ candidate.
 867                           // pdg code -3334  = Omega+
 868                   lInvMassOmegaPlus = xi->GetEffMassXi();
 869
 870                   lV0quality = 0.;
 871                   xi->ChangeMassHypothesis(lV0quality , -3312);         // Back to "default" hyp.
 872           }// end if positive bachelor
 873                   // - II.Step 6 : extra info for QA (ESD)
 874                   // miscellaneous pieces of info that may help regarding data quality assessment.
 875                   //-------------
 876
 877           xi->GetPxPyPz( lXiMomX, lXiMomY, lXiMomZ );
 878                   lXiTransvMom          = TMath::Sqrt( lXiMomX*lXiMomX   + lXiMomY*lXiMomY );
 879                   lXiTotMom     = TMath::Sqrt( lXiMomX*lXiMomX   + lXiMomY*lXiMomY   + lXiMomZ*lXiMomZ );
 880
 881           xi->GetBPxPyPz(  lBachMomX,  lBachMomY,  lBachMomZ );
 882                   //lBachTransvMom  = TMath::Sqrt( lBachMomX*lBachMomX   + lBachMomY*lBachMomY );
 883                   //lBachTotMom         = TMath::Sqrt( lBachMomX*lBachMomX   + lBachMomY*lBachMomY  +  lBachMomZ*lBachMomZ  );
 884
 885           lChargeXi = xi->Charge();
 886
 887           //lV0toXiCosineOfPointingAngle = xi->GetV0CosineOfPointingAngle( lPosXi[0], lPosXi[1], lPosXi[2] );
 888
 889           lRapXi    = xi->RapXi();
 890           lRapOmega = xi->RapOmega();
 891           //lEta      = xi->Eta();
 892           //lTheta    = xi->Theta() *180.0/TMath::Pi();
 893           //lPhi      = xi->Phi()   *180.0/TMath::Pi();
 894           //lAlphaXi  = xi->AlphaXi();
 895           //lPtArmXi  = xi->PtArmXi();
 896
 897   //------------------------------------------------
 898   // Set Variables for adding to tree
 899   //------------------------------------------------
 900
 901           fTreeCascVarCharge    = lChargeXi;
 902           if(lInvMassXiMinus!=0)    fTreeCascVarMassAsXi = lInvMassXiMinus;
 903           if(lInvMassXiPlus!=0)     fTreeCascVarMassAsXi = lInvMassXiPlus;
 904           if(lInvMassOmegaMinus!=0) fTreeCascVarMassAsOmega = lInvMassOmegaMinus;
 905           if(lInvMassOmegaPlus!=0)  fTreeCascVarMassAsOmega = lInvMassOmegaPlus;
 906           fTreeCascVarPt = lXiTransvMom;
 907           fTreeCascVarRapXi = lRapXi ;
 908           fTreeCascVarRapOmega = lRapOmega ;
 909           fTreeCascVarDCACascDaughters = lDcaXiDaughters;
 910           fTreeCascVarDCABachToPrimVtx = lDcaBachToPrimVertexXi;
 911           fTreeCascVarDCAV0Daughters = lDcaV0DaughtersXi;
 912           fTreeCascVarDCAV0ToPrimVtx = lDcaV0ToPrimVertexXi;
 913           fTreeCascVarDCAPosToPrimVtx = lDcaPosToPrimVertexXi;
 914           fTreeCascVarDCANegToPrimVtx = lDcaNegToPrimVertexXi;
 915           fTreeCascVarCascCosPointingAngle = lXiCosineOfPointingAngle;
 916           fTreeCascVarCascRadius = lXiRadius;
 917           fTreeCascVarV0Mass = lInvMassLambdaAsCascDghter;
 918           fTreeCascVarV0CosPointingAngle = lV0CosineOfPointingAngleXi;
 919           fTreeCascVarV0CosPointingAngleSpecial = lV0CosineOfPointingAngleXiSpecial;
 920           fTreeCascVarV0Radius = lV0RadiusXi;
 921           fTreeCascVarLeastNbrClusters = leastnumberofclusters;
 922
 923           //Copy Multiplicity information
 924           fTreeCascVarCentV0A = fCentrality_V0A;
 925           fTreeCascVarCentV0C = fCentrality_V0C;
 926           fTreeCascVarCentV0M = fCentrality_V0M;
 927           fTreeCascVarCentV0AEq = fCentrality_V0AEq;
 928           fTreeCascVarCentV0CEq = fCentrality_V0CEq;
 929           fTreeCascVarCentV0MEq = fCentrality_V0MEq;
 930           fTreeCascVarAmpV0A = fAmplitude_V0A;
 931           fTreeCascVarAmpV0C = fAmplitude_V0C;
 932           fTreeCascVarAmpV0AEq = fAmplitude_V0AEq;
 933           fTreeCascVarAmpV0CEq = fAmplitude_V0CEq;
 934           fTreeCascVarRefMultEta8 = fRefMultEta8;
 935           fTreeCascVarRefMultEta5 = fRefMultEta5;
 936           fTreeCascVarRunNumber = fRunNumber;
 937
 938           fTreeCascVarDistOverTotMom = TMath::Sqrt(
 939                                                 TMath::Power( lPosXi[0] - lBestPrimaryVtxPos[0] , 2) +
 940                                                 TMath::Power( lPosXi[1] - lBestPrimaryVtxPos[1] , 2) +
 941                                                 TMath::Power( lPosXi[2] - lBestPrimaryVtxPos[2] , 2)
 942                                         );
 943           fTreeCascVarDistOverTotMom /= (lXiTotMom+1e-13);
 944
 945 //All vars not specified here: specified elsewhere!
 946
 947 //------------------------------------------------
 948 // Fill Tree!
 949 //------------------------------------------------
 950
 951 // The conditional is meant to decrease excessive
 952 // memory usage! Be careful when loosening the
 953 // cut!
 954
 955   //Xi    Mass window: 150MeV wide
 956   //Omega mass window: 150MeV wide
 957
 958   if( fkSaveCascadeTree && ( (fTreeCascVarMassAsXi<1.32+0.075&&fTreeCascVarMassAsXi>1.32-0.075) ||
 959       (fTreeCascVarMassAsOmega<1.68+0.075&&fTreeCascVarMassAsOmega>1.68-0.075) ) ){
 960       fTreeCascade->Fill();
 961   }
 962
 963 //------------------------------------------------
 964 // Fill tree over.
 965 //------------------------------------------------
 966
 967         }// end of the Cascade loop (ESD or AOD)
 968
 969   // Post output data.
 970   PostData(1, fListHist);
 971   PostData(2, fTreeEvent);
 972   PostData(3, fTreeV0);
 973   PostData(4, fTreeCascade);
 974 }
 975
 976 //________________________________________________________________________
 977 void AliAnalysisTaskStrangenessVsMultiplicity::Terminate(Option_t *)
 978 {
 979    // Draw result to the screen
 980    // Called once at the end of the query
 981
 982    TList *cRetrievedList = 0x0;
 983    cRetrievedList = (TList*)GetOutputData(1);
 984    if(!cRetrievedList){
 985       Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : ouput data container list not available\n");
 986       return;
 987    }
 988
 989    fHistEventCounter = dynamic_cast<TH1D*> (  cRetrievedList->FindObject("fHistEventCounter")  );
 990    if (!fHistEventCounter) {
 991       Printf("ERROR - AliAnalysisTaskStrangenessVsMultiplicity : fHistEventCounter not available");
 992       return;
 993    }
 994
 995    TCanvas *canCheck = new TCanvas("AliAnalysisTaskStrangenessVsMultiplicity","V0 Multiplicity",10,10,510,510);
 996    canCheck->cd(1)->SetLogy();
 997
 998    fHistEventCounter->SetMarkerStyle(22);
 999    fHistEventCounter->DrawCopy("E");
1000 }
1001
1002 //----------------------------------------------------------------------------
1003
1004 Double_t AliAnalysisTaskStrangenessVsMultiplicity::MyRapidity(Double_t rE, Double_t rPz) const
1005 {
1006    // Local calculation for rapidity
1007    Double_t ReturnValue = -100;
1008    if( (rE-rPz+1.e-13) != 0 && (rE+rPz) != 0 ){
1009       ReturnValue =  0.5*TMath::Log((rE+rPz)/(rE-rPz+1.e-13));
1010    }
1011    return ReturnValue;
1012 }