PWGLF/QATasks/AliAnalysisTaskQAV0.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 // QA Task designed to investigate V0 characteristics in any data
  19 //
  20 // Stripped down and adapted from AliAnalysisTaskExtractV0:
  21 // --- smaller, more convenient output for debugging
  22 // --- optional TH3 list of histos to enable a light-weight analysis
  23 // --- This code is under development, so...
  24 //
  25 //  Please Report Any Bugs!
  26 //
  27 //   --- David Dobrigkeit Chinellato
  28 //        (david.chinellato@gmail.com)
  29 //
  30 // +-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
  31
  32 class TTree;
  33 class TParticle;
  34 class TVector3;
  35
  36 //class AliMCEventHandler;
  37 //class AliMCEvent;
  38 //class AliStack;
  39
  40 class AliESDVertex;
  41 class AliAODVertex;
  42 class AliESDv0;
  43 class AliAODv0;
  44
  45 #include <Riostream.h>
  46 #include "TList.h"
  47 #include "TH1.h"
  48 #include "TH2.h"
  49 #include "TH3.h"
  50 #include "THnSparse.h"
  51 #include "TVector3.h"
  52 #include "TCanvas.h"
  53 #include "TMath.h"
  54 #include "TLegend.h"
  55 #include "AliLog.h"
  56 #include "AliCentrality.h"
  57 #include "AliESDEvent.h"
  58 #include "AliAODEvent.h"
  59 #include "AliV0vertexer.h"
  60 #include "AliCascadeVertexer.h"
  61 #include "AliESDpid.h"
  62 #include "AliESDtrack.h"
  63 #include "AliESDtrackCuts.h"
  64 #include "AliInputEventHandler.h"
  65 #include "AliAnalysisManager.h"
  66 #include "AliMCEventHandler.h"
  67
  68 #include "AliCFContainer.h"
  69 #include "AliMultiplicity.h"
  70
  71 #include "AliESDcascade.h"
  72 #include "AliAODcascade.h"
  73 #include "AliESDUtils.h"
  74 #include "AliESDHeader.h"
  75
  76 #include "AliAnalysisUtils.h"
  77 #include "AliAnalysisTaskQAV0.h"
  78
  79 //debugging purposes
  80 #include "TObjectTable.h"
  81
  82 ClassImp(AliAnalysisTaskQAV0)
  83
  84 AliAnalysisTaskQAV0::AliAnalysisTaskQAV0()
  85   : AliAnalysisTaskSE(),
  86   //Output lists
  87   fOutput(0),
  88
  89   //Histos
  90   fHistEvent(0),
  91   fHistTopDCANegToPV(0),
  92   fHistTopDCAPosToPV(0),
  93   fHistTopDCAV0Daughters(0),
  94   fHistTopCosinePA(0),
  95   fHistTopV0Radius(0),
  96   fHistSelectedTopDCANegToPV(0),
  97   fHistSelectedTopDCAPosToPV(0),
  98   fHistSelectedTopDCAV0Daughters(0),
  99   fHistSelectedTopCosinePA(0),
 100   fHistSelectedTopV0Radius(0),
 101
 102   f2dHistInvMassK0Short(0),
 103   f2dHistInvMassLambda(0),
 104   f2dHistInvMassAntiLambda(0),
 105
 106   f2dHistInvMassWithdEdxK0Short(0),
 107   f2dHistInvMassWithdEdxLambda(0),
 108   f2dHistInvMassWithdEdxAntiLambda(0),
 109
 110   f2dHistResponseNegativeAsPion(0),
 111   f2dHistResponseNegativeAsProton(0),
 112   f2dHistResponsePositiveAsPion(0),
 113   f2dHistResponsePositiveAsProton(0),
 114
 115   f2dHistdEdxSignalPionFromLambda(0),
 116   f2dHistdEdxSignalProtonFromLambda(0),
 117   f2dHistResponsePionFromLambda(0),
 118   f2dHistResponseProtonFromLambda(0),
 119
 120
 121   //Task Control / Utils
 122   fPIDResponse(0),
 123   fkRunV0Vertexer ( kFALSE ),
 124   fdEdxCut (3)
 125 {
 126   // Dummy Constructor
 127   for(Int_t iV0selIdx   = 0; iV0selIdx   < 7; iV0selIdx++   ) { fV0Sels          [iV0selIdx   ] = -1.; }
 128 }
 129
 130 AliAnalysisTaskQAV0::AliAnalysisTaskQAV0(const char *name)
 131   : AliAnalysisTaskSE(name),
 132   //Output lists
 133   fOutput(0),
 134
 135   //Histos
 136   fHistEvent(0),
 137   fHistTopDCANegToPV(0),
 138   fHistTopDCAPosToPV(0),
 139   fHistTopDCAV0Daughters(0),
 140   fHistTopCosinePA(0),
 141   fHistTopV0Radius(0),
 142   fHistSelectedTopDCANegToPV(0),
 143   fHistSelectedTopDCAPosToPV(0),
 144   fHistSelectedTopDCAV0Daughters(0),
 145   fHistSelectedTopCosinePA(0),
 146   fHistSelectedTopV0Radius(0),
 147
 148   f2dHistInvMassK0Short(0),
 149   f2dHistInvMassLambda(0),
 150   f2dHistInvMassAntiLambda(0),
 151
 152   f2dHistInvMassWithdEdxK0Short(0),
 153   f2dHistInvMassWithdEdxLambda(0),
 154   f2dHistInvMassWithdEdxAntiLambda(0),
 155
 156   f2dHistResponseNegativeAsPion(0),
 157   f2dHistResponseNegativeAsProton(0),
 158   f2dHistResponsePositiveAsPion(0),
 159   f2dHistResponsePositiveAsProton(0),
 160
 161   f2dHistdEdxSignalPionFromLambda(0),
 162   f2dHistdEdxSignalProtonFromLambda(0),
 163   f2dHistResponsePionFromLambda(0),
 164   f2dHistResponseProtonFromLambda(0),
 165
 166   //Task Control / Utils
 167   fPIDResponse(0),
 168   fkRunV0Vertexer ( kFALSE ),
 169   fdEdxCut (3)
 170 {
 171   // Constructor
 172   //VERTEXER CUTS
 173   // REALLY LOOSE? Be careful when attempting to run over PbPb if fkRunV0Vertexer is set!
 174   fV0VertexerSels[0] =  33.  ;  // max allowed chi2
 175   fV0VertexerSels[1] =   0.02;  // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 176   fV0VertexerSels[2] =   0.02;  // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 177   fV0VertexerSels[3] =   2.0 ;  // max allowed DCA between the daughter tracks       (LHC09a4 : 0.5)
 178   fV0VertexerSels[4] =   0.95;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 179   fV0VertexerSels[5] =   0.5 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 180   fV0VertexerSels[6] = 200.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 181
 182   //SELECTION CUTS
 183   // REALLY LOOSE? Be careful when attempting to run over PbPb if fkRunV0Vertexer is set!
 184   fV0Sels[0] =  33.  ;  // max allowed chi2
 185   fV0Sels[1] =   0.02;  // min allowed impact parameter for the 1st daughter (LHC09a4 : 0.05)
 186   fV0Sels[2] =   0.02;  // min allowed impact parameter for the 2nd daughter (LHC09a4 : 0.05)
 187   fV0Sels[3] =   2.0 ;  // max allowed DCA between the daughter tracks       (LHC09a4 : 0.5)
 188   fV0Sels[4] =   0.95;  // min allowed cosine of V0's pointing angle         (LHC09a4 : 0.99)
 189   fV0Sels[5] =   0.5 ;  // min radius of the fiducial volume                 (LHC09a4 : 0.2)
 190   fV0Sels[6] = 200.  ;  // max radius of the fiducial volume                 (LHC09a4 : 100.0)
 191
 192   // Output slot #0 writes into a TList container (Lambda Histos and fTree)
 193   DefineOutput(1, TList::Class());
 194 }
 195
 196
 197 AliAnalysisTaskQAV0::~AliAnalysisTaskQAV0()
 198 {
 199 //------------------------------------------------
 200 // DESTRUCTOR
 201 //------------------------------------------------
 202
 203    if (fOutput){
 204       delete fOutput;
 205       fOutput = 0x0;
 206    }
 207 }
 208
 209 //________________________________________________________________________
 210 void AliAnalysisTaskQAV0::UserCreateOutputObjects()
 211 {
 212   //Define Output Lists
 213   fOutput = new TList();
 214   fOutput->SetOwner();
 215
 216   //Histogram Output: Event-by-Event
 217   fHistEvent = new TH1D( "fHistEvent", ";Evt. Sel. Step;Count",4,0,4);
 218   fHistEvent->GetXaxis()->SetBinLabel(1, "Processed");
 219   fHistEvent->GetXaxis()->SetBinLabel(2, "Phys-Sel");
 220   fHistEvent->GetXaxis()->SetBinLabel(3, "Has Vtx");
 221   fHistEvent->GetXaxis()->SetBinLabel(4, "Vtx |z|<10cm");
 222   fOutput->Add(fHistEvent);
 223
 224   //Topological Selection Histograms, 1D
 225   fHistTopDCANegToPV        = new TH1D( "fHistTopDCANegToPV",";DCA Neg. Daughter to PV (cm);Counts",200,0,1);
 226   fHistTopDCAPosToPV        = new TH1D( "fHistTopDCAPosToPV",";DCA Pos. Daughter to PV (cm);Counts",200,0,1);
 227   fHistTopDCAV0Daughters    = new TH1D( "fHistTopDCAV0Daughters",";DCA V0 Daughters (#sigma);Counts",200,0,2);
 228   fHistTopCosinePA          = new TH1D( "fHistTopCosinePA",";Cosine of PA;Counts",10010,-1.001,1.001);
 229   fHistTopV0Radius          = new TH1D( "fHistTopV0Radius",";Decay Radius (cm);Counts",200,0.,10);
 230
 231   fOutput->Add( fHistTopDCANegToPV     );
 232   fOutput->Add( fHistTopDCAPosToPV     );
 233   fOutput->Add( fHistTopDCAV0Daughters );
 234   fOutput->Add( fHistTopCosinePA       );
 235   fOutput->Add( fHistTopV0Radius       );
 236
 237   //Zoomed In
 238   fHistSelectedTopDCANegToPV        = new TH1D( "fHistSelectedTopDCANegToPV",";DCA Neg. Daughter to PV (cm);Counts",200,fV0Sels[1],1);
 239   fHistSelectedTopDCAPosToPV        = new TH1D( "fHistSelectedTopDCAPosToPV",";DCA Pos. Daughter to PV (cm);Counts",200,fV0Sels[2],1);
 240   fHistSelectedTopDCAV0Daughters    = new TH1D( "fHistSelectedTopDCAV0Daughters",";DCA V0 Daughters (#sigma);Counts",200,0,fV0Sels[3]);
 241   fHistSelectedTopCosinePA          = new TH1D( "fHistSelectedTopCosinePA",";Cosine of PA;Counts",400,fV0Sels[4],1);
 242   fHistSelectedTopV0Radius          = new TH1D( "fHistSelectedTopV0Radius",";Decay Radius (cm);Counts",200,fV0Sels[5],10);
 243
 244   fOutput->Add( fHistSelectedTopDCANegToPV     );
 245   fOutput->Add( fHistSelectedTopDCAPosToPV     );
 246   fOutput->Add( fHistSelectedTopDCAV0Daughters );
 247   fOutput->Add( fHistSelectedTopCosinePA       );
 248   fOutput->Add( fHistSelectedTopV0Radius       );
 249
 250   //Invariant Mass Plots
 251   f2dHistInvMassK0Short     = new TH2D( "f2dHistInvMassK0Short"   , ";p_{T};M(#pi^{+},#pi^{-})"   ,250,0,25,500,0.25,0.75);
 252   f2dHistInvMassLambda      = new TH2D( "f2dHistInvMassLambda"    , ";p_{T};M(p,#pi^{-})"         ,250,0,25,300,1.07,1.115+0.255);
 253   f2dHistInvMassAntiLambda  = new TH2D( "f2dHistInvMassAntiLambda", ";p_{T};M(#pi^{+},#bar{p})"   ,250,0,25,300,1.07,1.115+0.255);
 254
 255   fOutput->Add( f2dHistInvMassK0Short     );
 256   fOutput->Add( f2dHistInvMassLambda      );
 257   fOutput->Add( f2dHistInvMassAntiLambda  );
 258
 259   //Invariant Mass Plots, with dEdx
 260   f2dHistInvMassWithdEdxK0Short     = new TH2D( "f2dHistInvMassWithdEdxK0Short"   , ";p_{T};M(#pi^{+},#pi^{-})"   ,250,0,25,500,0.25,0.75);
 261   f2dHistInvMassWithdEdxLambda      = new TH2D( "f2dHistInvMassWithdEdxLambda"    , ";p_{T};M(p,#pi^{-})"         ,250,0,25,300,1.07,1.115+0.255);
 262   f2dHistInvMassWithdEdxAntiLambda  = new TH2D( "f2dHistInvMassWithdEdxAntiLambda", ";p_{T};M(#pi^{+},#bar{p})"   ,250,0,25,300,1.07,1.115+0.255);
 263
 264   fOutput->Add( f2dHistInvMassWithdEdxK0Short     );
 265   fOutput->Add( f2dHistInvMassWithdEdxLambda      );
 266   fOutput->Add( f2dHistInvMassWithdEdxAntiLambda  );
 267
 268   //dE/dx QA for main analysis and for calibration check
 269   f2dHistResponseNegativeAsPion    = new TH2D( "f2dHistResponseNegativeAsPion", ";p_{T}^{V0};N#sigma",500,0,5,400,-20,20);
 270   f2dHistResponseNegativeAsProton  = new TH2D( "f2dHistResponseNegativeAsProton", ";p_{T}^{V0};N#sigma",500,0,5,400,-20,20);
 271   f2dHistResponsePositiveAsPion    = new TH2D( "f2dHistResponsePositiveAsPion", ";p_{T}^{V0};N#sigma",500,0,5,400,-20,20);
 272   f2dHistResponsePositiveAsProton  = new TH2D( "f2dHistResponsePositiveAsProton", ";p_{T}^{V0};N#sigma",500,0,5,400,-20,20);
 273
 274   //Clean Signal Check from Lambdas: stricter cuts, raw signal check
 275   f2dHistdEdxSignalPionFromLambda    = new TH2D( "f2dHistdEdxSignalPionFromLambda", ";p_{T}^{V0};TPC Signal",500,0,5,8000,0,800);
 276   f2dHistdEdxSignalProtonFromLambda  = new TH2D( "f2dHistdEdxSignalProtonFromLambda", ";p_{T}^{V0};TPC Signal",500,0,5,8000,0,800);
 277   f2dHistResponsePionFromLambda      = new TH2D( "f2dHistResponsePionFromLambda", ";p_{T}^{V0};N#sigma",500,0,5,400,-20,20);
 278   f2dHistResponseProtonFromLambda    = new TH2D( "f2dHistResponseProtonFromLambda", ";p_{T}^{V0};N#sigma",500,0,5,400,-20,20);
 279
 280
 281   fOutput->Add( f2dHistResponseNegativeAsPion        );
 282   fOutput->Add( f2dHistResponseNegativeAsProton      );
 283   fOutput->Add( f2dHistResponsePositiveAsPion        );
 284   fOutput->Add( f2dHistResponsePositiveAsProton      );
 285
 286   fOutput->Add( f2dHistdEdxSignalPionFromLambda        );
 287   fOutput->Add( f2dHistdEdxSignalProtonFromLambda      );
 288   fOutput->Add( f2dHistResponsePionFromLambda          );
 289   fOutput->Add( f2dHistResponseProtonFromLambda        );
 290
 291 //------------------------------------------------
 292 // Particle Identification Setup
 293 //------------------------------------------------
 294
 295   AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
 296   AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
 297   fPIDResponse = inputHandler->GetPIDResponse();
 298
 299    //Regular output: Histograms
 300    PostData(1, fOutput);
 301 }// end UserCreateOutputObjects
 302
 303
 304 //________________________________________________________________________
 305 void AliAnalysisTaskQAV0::UserExec(Option_t *)
 306 {
 307    // Main loop
 308    // Called for each event
 309    //gObjectTable->Print();
 310    AliESDEvent *lESDevent = 0x0;
 311
 312    //AliAODEvent *lAODevent = 0x0;
 313    Int_t    nV0s                        = -1;
 314
 315    Double_t lTrkgPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
 316    Double_t lBestPrimaryVtxPos[3]          = {-100.0, -100.0, -100.0};
 317    Double_t lMagneticField                 = -10.;
 318
 319    // Connect to the InputEvent
 320    // After these lines, we should have an ESD/AOD event + the number of cascades in it.
 321
 322    // Appropriate for ESD analysis!
 323    lESDevent = dynamic_cast<AliESDEvent*>( InputEvent() );
 324    if (!lESDevent) {
 325       AliWarning("ERROR: lESDevent not available \n");
 326       return;
 327    }
 328
 329   //------------------------------------------------
 330   // Rerun V0 vertexer, if asked for
 331   // --- WARNING: Be careful when using in PbPb
 332   //------------------------------------------------
 333   if( fkRunV0Vertexer ){
 334     lESDevent->ResetV0s();
 335     AliV0vertexer lV0vtxer;
 336     lV0vtxer.SetDefaultCuts(fV0VertexerSels);
 337     lV0vtxer.Tracks2V0vertices(lESDevent);
 338   }
 339
 340   fHistEvent->Fill(0.5);
 341
 342 //------------------------------------------------
 343 // Physics Selection
 344 //------------------------------------------------
 345
 346   // new method
 347   UInt_t maskIsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
 348   Bool_t isSelected = 0;
 349   isSelected = (maskIsSelected & AliVEvent::kMB) == AliVEvent::kMB;
 350
 351   //pA triggering: CINT7
 352   //if( fkSwitchINT7 ) isSelected = (maskIsSelected & AliVEvent::kINT7) == AliVEvent::kINT7;
 353
 354   //Standard Min-Bias Selection
 355   if ( !isSelected ) {
 356     PostData(1, fOutput);
 357     return;
 358   }
 359   fHistEvent->Fill(1.5);
 360
 361 //------------------------------------------------
 362 // After Trigger Selection
 363 //------------------------------------------------
 364
 365
 366         nV0s = lESDevent->GetNumberOfV0s();
 367
 368 //------------------------------------------------
 369 // Getting: Primary Vertex + MagField Info
 370 //------------------------------------------------
 371
 372         const AliESDVertex *lPrimaryTrackingESDVtx = lESDevent->GetPrimaryVertexTracks();
 373         // get the vtx stored in ESD found with tracks
 374         lPrimaryTrackingESDVtx->GetXYZ( lTrkgPrimaryVtxPos );
 375
 376         const AliESDVertex *lPrimaryBestESDVtx = lESDevent->GetPrimaryVertex();
 377   // get the best primary vertex available for the event
 378         // As done in AliCascadeVertexer, we keep the one which is the best one available.
 379         // between : Tracking vertex > SPD vertex > TPC vertex > default SPD vertex
 380         // This one will be used for next calculations (DCA essentially)
 381    lPrimaryBestESDVtx->GetXYZ( lBestPrimaryVtxPos );
 382
 383    Double_t tPrimaryVtxPosition[3];
 384    const AliVVertex *primaryVtx = lESDevent->GetPrimaryVertex();
 385    tPrimaryVtxPosition[0] = primaryVtx->GetX();
 386    tPrimaryVtxPosition[1] = primaryVtx->GetY();
 387    tPrimaryVtxPosition[2] = primaryVtx->GetZ();
 388
 389   //------------------------------------------------
 390   // Primary Vertex Requirements Section:
 391   //  ---> pp and PbPb: Only requires |z|<10cm
 392   //  ---> pPb: all requirements checked at this stage
 393   //------------------------------------------------
 394
 395   //Roberto's PV selection criteria, implemented 17th April 2013
 396
 397   // vertex selection
 398   Bool_t fHasVertex = kFALSE;
 399   const AliESDVertex *vertex = lESDevent->GetPrimaryVertexTracks();
 400   if (vertex->GetNContributors() < 1) {
 401     vertex = lESDevent->GetPrimaryVertexSPD();
 402     if (vertex->GetNContributors() < 1) fHasVertex = kFALSE;
 403     else fHasVertex = kTRUE;
 404     TString vtxTyp = vertex->GetTitle();
 405     Double_t cov[6]={0};
 406     vertex->GetCovarianceMatrix(cov);
 407     Double_t zRes = TMath::Sqrt(cov[5]);
 408     if (vtxTyp.Contains("vertexer:Z") && (zRes>0.25)) fHasVertex = kFALSE;
 409   }
 410   else fHasVertex = kTRUE;
 411
 412   //Is First event in chunk rejection: Still present!
 413   if(fHasVertex == kFALSE) {
 414     AliWarning("Pb / | PV does not satisfy selection criteria!");
 415     PostData(1, fOutput);
 416     return;
 417   }
 418
 419   fHistEvent->Fill(2.5);
 420
 421   //17 April Fix: Always do primary vertex Z selection, after pA vertex selection from Roberto
 422   if(TMath::Abs(lBestPrimaryVtxPos[2]) > 10.0) {
 423     AliWarning("Pb / | Z position of Best Prim Vtx | > 10.0 cm ... return !");
 424     PostData(1, fOutput);
 425     return;
 426   }
 427
 428   fHistEvent->Fill(3.5);
 429
 430   lMagneticField = lESDevent->GetMagneticField( );
 431
 432 //------------------------------------------------
 433 // Only look at events with well-established PV
 434 //------------------------------------------------
 435
 436    //const AliESDVertex *lPrimaryTrackingESDVtxCheck = lESDevent->GetPrimaryVertexTracks();
 437    const AliESDVertex *lPrimarySPDVtx = lESDevent->GetPrimaryVertexSPD();
 438    //if (!lPrimarySPDVtx->GetStatus() && !lPrimaryTrackingESDVtxCheck->GetStatus() ){
 439    //   AliWarning("Pb / No SPD prim. vertex nor prim. Tracking vertex ... return !");
 440    //   PostData(1, fOutput);
 441    //   return;
 442    //}
 443
 444
 445
 446 //------------------------------------------------
 447 // MAIN LAMBDA LOOP STARTS HERE
 448 //------------------------------------------------
 449
 450 //Variable definition
 451   Int_t    lOnFlyStatus = 0;// nv0sOn = 0, nv0sOff = 0;
 452   Double_t lChi2V0 = 0;
 453   Double_t lDcaV0Daughters = 0, lDcaV0ToPrimVertex = 0;
 454   Double_t lDcaPosToPrimVertex = 0, lDcaNegToPrimVertex = 0;
 455   Double_t lV0CosineOfPointingAngle = 0;
 456   Double_t lV0Radius = 0, lPt = 0;
 457   Double_t lRapK0Short = 0, lRapLambda = 0;
 458   Double_t lInvMassK0s = 0, lInvMassLambda = 0, lInvMassAntiLambda = 0;
 459   Double_t lAlphaV0 = 0, lPtArmV0 = 0;
 460   Double_t lNegEta = -100, lPosEta = -100;
 461   Double_t fMinV0Pt = 0;
 462   Double_t fMaxV0Pt = 100;
 463
 464    Int_t nv0s = 0;
 465    nv0s = lESDevent->GetNumberOfV0s();
 466
 467    //for (Int_t iV0 = 0; iV0 < nv0s; iV0++)
 468    for (Int_t iV0 = 0; iV0 < nv0s; iV0++) //extra-crazy test
 469    {// This is the begining of the V0 loop
 470       AliESDv0 *v0 = ((AliESDEvent*)lESDevent)->GetV0(iV0);
 471       if (!v0) continue;
 472
 473       //Only use Offline Candidates for QA
 474       lOnFlyStatus = v0->GetOnFlyStatus();
 475       if( lOnFlyStatus == kTRUE ) continue;
 476
 477       Double_t tDecayVertexV0[3]; v0->GetXYZ(tDecayVertexV0[0],tDecayVertexV0[1],tDecayVertexV0[2]);
 478
 479       Double_t tV0mom[3];
 480       v0->GetPxPyPz( tV0mom[0],tV0mom[1],tV0mom[2] );
 481       Double_t lV0TotalMomentum = TMath::Sqrt(
 482       tV0mom[0]*tV0mom[0]+tV0mom[1]*tV0mom[1]+tV0mom[2]*tV0mom[2] );
 483
 484       lV0Radius = TMath::Sqrt(tDecayVertexV0[0]*tDecayVertexV0[0]+tDecayVertexV0[1]*tDecayVertexV0[1]);
 485
 486       lPt = v0->Pt();
 487       lRapK0Short = v0->RapK0Short();
 488       lRapLambda  = v0->RapLambda();
 489       if ((lPt<fMinV0Pt)||(fMaxV0Pt<lPt)) continue;
 490
 491       UInt_t lKeyPos = (UInt_t)TMath::Abs(v0->GetPindex());
 492       UInt_t lKeyNeg = (UInt_t)TMath::Abs(v0->GetNindex());
 493
 494       Double_t lMomPos[3]; v0->GetPPxPyPz(lMomPos[0],lMomPos[1],lMomPos[2]);
 495       Double_t lMomNeg[3]; v0->GetNPxPyPz(lMomNeg[0],lMomNeg[1],lMomNeg[2]);
 496
 497       AliESDtrack *pTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyPos);
 498       AliESDtrack *nTrack=((AliESDEvent*)lESDevent)->GetTrack(lKeyNeg);
 499       if (!pTrack || !nTrack) {
 500          Printf("ERROR: Could not retreive one of the daughter track");
 501          continue;
 502       }
 503
 504       //Daughter Eta for Eta selection, afterwards
 505       lNegEta = nTrack->Eta();
 506       lPosEta = pTrack->Eta();
 507
 508       // Filter like-sign V0 (next: add counter and distribution)
 509       if ( pTrack->GetSign() == nTrack->GetSign()){
 510          continue;
 511       }
 512
 513       //________________________________________________________________________
 514       // Track quality cuts
 515       Float_t lPosTrackCrossedRows = pTrack->GetTPCClusterInfo(2,1);
 516       Float_t lNegTrackCrossedRows = nTrack->GetTPCClusterInfo(2,1);
 517       Int_t lLeastNbrCrossedRows = (Int_t) lPosTrackCrossedRows;
 518       if( lNegTrackCrossedRows < lLeastNbrCrossedRows )
 519          lLeastNbrCrossedRows = (Int_t) lNegTrackCrossedRows;
 520
 521       // TPC refit condition (done during reconstruction for Offline but not for On-the-fly)
 522       if( !(pTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
 523       if( !(nTrack->GetStatus() & AliESDtrack::kTPCrefit)) continue;
 524
 525       //Get status flags
 526       //lPosTrackStatus = pTrack->GetStatus();
 527       //lNegTrackStatus = nTrack->GetStatus();
 528
 529       if ( ( ( pTrack->GetTPCClusterInfo(2,1) ) < 70 ) || ( ( nTrack->GetTPCClusterInfo(2,1) ) < 70 ) ) continue;
 530
 531       //GetKinkIndex condition
 532       if( pTrack->GetKinkIndex(0)>0 || nTrack->GetKinkIndex(0)>0 ) continue;
 533
 534       //Findable clusters > 0 condition
 535       if( pTrack->GetTPCNclsF()<=0 || nTrack->GetTPCNclsF()<=0 ) continue;
 536
 537       //Compute ratio Crossed Rows / Findable clusters
 538       //Note: above test avoids division by zero!
 539       Float_t lPosTrackCrossedRowsOverFindable = lPosTrackCrossedRows / ((double)(pTrack->GetTPCNclsF()));
 540       Float_t lNegTrackCrossedRowsOverFindable = lNegTrackCrossedRows / ((double)(nTrack->GetTPCNclsF()));
 541
 542       Float_t lLeastRatioCrossedRowsOverFindable = lPosTrackCrossedRowsOverFindable;
 543       if( lNegTrackCrossedRowsOverFindable < lLeastRatioCrossedRowsOverFindable )
 544          lLeastRatioCrossedRowsOverFindable = lNegTrackCrossedRowsOverFindable;
 545
 546       //Lowest Cut Level for Ratio Crossed Rows / Findable = 0.8, set here
 547       if ( lLeastRatioCrossedRowsOverFindable < 0.8 ) continue;
 548
 549       //End track Quality Cuts
 550       //________________________________________________________________________
 551
 552       lDcaPosToPrimVertex = TMath::Abs(pTrack->GetD(tPrimaryVtxPosition[0],
 553                                                         tPrimaryVtxPosition[1],
 554                                                         lMagneticField) );
 555
 556       lDcaNegToPrimVertex = TMath::Abs(nTrack->GetD(tPrimaryVtxPosition[0],
 557                                                         tPrimaryVtxPosition[1],
 558                                                         lMagneticField) );
 559
 560
 561       lChi2V0 = v0->GetChi2V0();
 562       lDcaV0Daughters = v0->GetDcaV0Daughters();
 563       lDcaV0ToPrimVertex = v0->GetD(tPrimaryVtxPosition[0],tPrimaryVtxPosition[1],tPrimaryVtxPosition[2]);
 564       lV0CosineOfPointingAngle = v0->GetV0CosineOfPointingAngle(tPrimaryVtxPosition[0],tPrimaryVtxPosition[1],tPrimaryVtxPosition[2]);
 565       //fTreeVariableV0CosineOfPointingAngle=lV0CosineOfPointingAngle;
 566
 567       // Getting invariant mass infos directly from ESD
 568       v0->ChangeMassHypothesis(310);
 569       lInvMassK0s = v0->GetEffMass();
 570       v0->ChangeMassHypothesis(3122);
 571       lInvMassLambda = v0->GetEffMass();
 572       v0->ChangeMassHypothesis(-3122);
 573       lInvMassAntiLambda = v0->GetEffMass();
 574       lAlphaV0 = v0->AlphaV0();
 575       lPtArmV0 = v0->PtArmV0();
 576
 577       //Official means of acquiring N-sigmas
 578       Float_t lNSigmasPosProton = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kProton );
 579       Float_t lNSigmasPosPion   = fPIDResponse->NumberOfSigmasTPC( pTrack, AliPID::kPion );
 580       Float_t lNSigmasNegProton = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kProton );
 581       Float_t lNSigmasNegPion   = fPIDResponse->NumberOfSigmasTPC( nTrack, AliPID::kPion );
 582
 583 //This requires an Invariant Mass Hypothesis afterwards
 584       Float_t lDistOverTotMom = TMath::Sqrt(
 585                                                 TMath::Power( tDecayVertexV0[0] - lBestPrimaryVtxPos[0] , 2) +
 586                                                 TMath::Power( tDecayVertexV0[1] - lBestPrimaryVtxPos[1] , 2) +
 587                                                 TMath::Power( tDecayVertexV0[2] - lBestPrimaryVtxPos[2] , 2)
 588                                         );
 589       lDistOverTotMom /= (lV0TotalMomentum+1e-10); //avoid division by zero, to be sure
 590
 591 //------------------------------------------------
 592 // Fill Main Output Histograms
 593 //------------------------------------------------
 594
 595   //Topological Variable Checks, One-Dimensional
 596   fHistTopDCANegToPV      -> Fill( lDcaNegToPrimVertex      ) ;
 597   fHistTopDCAPosToPV      -> Fill( lDcaPosToPrimVertex      ) ;
 598   fHistTopDCAV0Daughters  -> Fill( lDcaV0Daughters          ) ;
 599   fHistTopCosinePA        -> Fill( lV0CosineOfPointingAngle ) ;
 600   fHistTopV0Radius        -> Fill( lV0Radius                ) ;
 601
 602
 603   if( lDcaNegToPrimVertex > fV0Sels[1] && lDcaPosToPrimVertex > fV0Sels[2]      &&
 604       lDcaV0Daughters     < fV0Sels[3] && lV0CosineOfPointingAngle > fV0Sels[4] &&
 605       lV0Radius           > fV0Sels[5] && lV0Radius < fV0Sels [6] ){
 606
 607     //Topological Variables zoomed in at selection level (whatever that may be)
 608     //May be slightly redundant if no specific extra configuration was done
 609     fHistSelectedTopDCANegToPV      -> Fill( lDcaNegToPrimVertex      ) ;
 610     fHistSelectedTopDCAPosToPV      -> Fill( lDcaPosToPrimVertex      ) ;
 611     fHistSelectedTopDCAV0Daughters  -> Fill( lDcaV0Daughters          ) ;
 612     fHistSelectedTopCosinePA        -> Fill( lV0CosineOfPointingAngle ) ;
 613     fHistSelectedTopV0Radius        -> Fill( lV0Radius                ) ;
 614
 615     //Specific fV0Sel selection level, but no dEdx applied
 616     f2dHistInvMassK0Short     -> Fill ( lPt , lInvMassK0s        )   ;
 617     f2dHistInvMassLambda      -> Fill ( lPt , lInvMassLambda     )   ;
 618     f2dHistInvMassAntiLambda  -> Fill ( lPt , lInvMassAntiLambda )   ;
 619
 620     //General dE/dx QA
 621     f2dHistResponseNegativeAsPion     -> Fill( lPt, lNSigmasNegPion      );
 622     f2dHistResponseNegativeAsProton   -> Fill( lPt, lNSigmasNegProton    );
 623     f2dHistResponsePositiveAsPion     -> Fill( lPt, lNSigmasPosPion      );
 624     f2dHistResponsePositiveAsProton   -> Fill( lPt, lNSigmasPosProton    );
 625
 626     //Clean Sample From Lambdas
 627     //Very strict cuts to ensure dealing with good Lambdas
 628     if ( lDcaV0Daughters < 1.0 && lV0CosineOfPointingAngle > 0.999 && TMath::Abs( lInvMassK0s - 0.497614 ) > 0.012
 629           && TMath::Abs( lInvMassAntiLambda - 1.115683) > 0.08 && TMath::Abs( lInvMassLambda - 1.115683) < 0.002 ) {
 630
 631       f2dHistdEdxSignalPionFromLambda     -> Fill( lPt, nTrack-> GetTPCsignal() );
 632       f2dHistdEdxSignalProtonFromLambda   -> Fill( lPt, pTrack-> GetTPCsignal() );
 633       f2dHistResponsePionFromLambda     -> Fill( lPt, lNSigmasNegPion   );
 634       f2dHistResponseProtonFromLambda   -> Fill( lPt, lNSigmasPosProton );
 635     }
 636
 637     //Specific fV0Sel selection level, dE/dx applied
 638     if ( TMath::Abs(lNSigmasPosPion)   < fdEdxCut && TMath::Abs(lNSigmasNegPion)   < fdEdxCut )    f2dHistInvMassWithdEdxK0Short       -> Fill ( lPt , lInvMassK0s        )   ;
 639     if ( TMath::Abs(lNSigmasPosProton) < fdEdxCut && TMath::Abs(lNSigmasNegPion)   < fdEdxCut )    f2dHistInvMassWithdEdxLambda        -> Fill ( lPt , lInvMassLambda     )   ;
 640     if ( TMath::Abs(lNSigmasPosPion)   < fdEdxCut && TMath::Abs(lNSigmasNegProton) < fdEdxCut )    f2dHistInvMassWithdEdxAntiLambda    -> Fill ( lPt , lInvMassAntiLambda )   ;
 641
 642   }
 643
 644
 645   //Topological Variable Checks, Two-Dimensional
 646
 647 //------------------------------------------------
 648 // End Filling of main Histograms
 649 //------------------------------------------------
 650
 651   }// This is the end of the V0 loop
 652
 653
 654   // Post output data.
 655   PostData(1, fOutput);
 656
 657 }
 658
 659 //________________________________________________________________________
 660 void AliAnalysisTaskQAV0::Terminate(Option_t *)
 661 {
 662    // Draw result to the screen
 663    // Called once at the end of the query
 664    // This will draw the V0 candidate multiplicity, whose
 665    // number of entries corresponds to the number of triggered events.
 666    TList *cRetrievedList = 0x0;
 667    cRetrievedList = (TList*)GetOutputData(1);
 668    if(!cRetrievedList){
 669       Printf("ERROR - AliAnalysisTaskQAV0 : ouput data container list not available\n");
 670       return;
 671    }
 672    fHistEvent = dynamic_cast<TH1D*> (  cRetrievedList->FindObject("fHistEvent")  );
 673    if (!fHistEvent) {
 674       Printf("ERROR - AliAnalysisTaskQAV0 : fHistEvent not available");
 675       return;
 676    }
 677    TCanvas *canCheck = new TCanvas("AliAnalysisTaskQAV0","V0 Multiplicity",10,10,510,510);
 678    canCheck->cd(1)->SetLogy();
 679    fHistEvent->SetMarkerStyle(22);
 680    fHistEvent->DrawCopy("E");
 681 }