PWGHF/vertexingHF/AliAnalysisTaskSELc2V0bachelor.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2009, 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  * appeuear 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 /* $Id$ */
  17
  18 //
  19 //
  20 //               Base class for Lc2V0 Analysis
  21 //
  22 //
  23 //  The Lc spectra study is done 2D histograms:
  24 //   cascade_invariantMass VS cascade_pT
  25 //
  26 //  Cuts have been centralized in AliRDHFCutsLctoV0 class
  27 //
  28 //-------------------------------------------------------------------------
  29 //
  30 //                 Authors: A.De Caro(a,b), P. Pagano(b)
  31 //  (a) Centro 'E.Fermi' - Roma
  32 //  (b) INFN and University of Salerno
  33 //
  34 //  Contatcs: decaro@sa.infn.it
  35 //            paola.pagano@sa.infn.it
  36 //-------------------------------------------------------------------------
  37
  38 #include <TSystem.h>
  39 #include <TParticle.h>
  40 #include <TParticlePDG.h>
  41 #include <TH1F.h>
  42 #include <TH1F.h>
  43 #include <TH2F.h>
  44 #include <TTree.h>
  45 #include "TROOT.h"
  46 #include <TDatabasePDG.h>
  47 #include <AliAnalysisDataSlot.h>
  48 #include <AliAnalysisDataContainer.h>
  49 #include "AliStack.h"
  50 #include "AliMCEvent.h"
  51 #include "AliAnalysisManager.h"
  52 #include "AliAODMCHeader.h"
  53 #include "AliAODHandler.h"
  54 #include "AliLog.h"
  55 #include "AliExternalTrackParam.h"
  56 #include "AliAODVertex.h"
  57 #include "AliAODRecoDecay.h"
  58 #include "AliAODRecoDecayHF.h"
  59 #include "AliAODRecoCascadeHF.h"
  60 #include "AliAnalysisVertexingHF.h"
  61 #include "AliESDtrack.h"
  62 #include "AliAODTrack.h"
  63 #include "AliAODv0.h"
  64 #include "AliAODMCParticle.h"
  65 #include "AliAnalysisTaskSE.h"
  66 #include "AliAnalysisTaskSELc2V0bachelor.h"
  67 #include "AliNormalizationCounter.h"
  68 #include "AliAODPidHF.h"
  69 #include "AliInputEventHandler.h"
  70 #include "AliESDtrackCuts.h"
  71 #include "AliNeutralTrackParam.h"
  72
  73 using std::cout;
  74 using std::endl;
  75
  76 ClassImp(AliAnalysisTaskSELc2V0bachelor)
  77
  78 //__________________________________________________________________________
  79 AliAnalysisTaskSELc2V0bachelor::AliAnalysisTaskSELc2V0bachelor() : AliAnalysisTaskSE(),
  80   fUseMCInfo(kFALSE),
  81   fOutput(0),
  82   fOutputAll(0),
  83   fOutputPIDBach(0),
  84   fCEvents(0),
  85   fIsK0SAnalysis(kFALSE),
  86   fCounter(0),
  87   fAnalCuts(0),
  88   fUseOnTheFlyV0(kFALSE),
  89   fIsEventSelected(kFALSE),
  90   fWriteVariableTree(kFALSE),
  91   fVariablesTree(0),
  92   fCandidateVariables(),
  93   fVtx1(0),
  94   fBzkG(0),
  95   fAdditionalChecks(kFALSE)
  96 {
  97   //
  98   // Default ctor
  99   //
 100 }
 101 //___________________________________________________________________________
 102 AliAnalysisTaskSELc2V0bachelor::AliAnalysisTaskSELc2V0bachelor(const Char_t* name,
 103                                                                AliRDHFCutsLctoV0* analCuts, Bool_t useOnTheFly,
 104                                                                Bool_t writeVariableTree, Bool_t additionalChecks) :
 105   AliAnalysisTaskSE(name),
 106   fUseMCInfo(kFALSE),
 107   fOutput(0),
 108   fOutputAll(0),
 109   fOutputPIDBach(0),
 110   fCEvents(0),
 111   fIsK0SAnalysis(kFALSE),
 112   fCounter(0),
 113   fAnalCuts(analCuts),
 114   fUseOnTheFlyV0(useOnTheFly),
 115   fIsEventSelected(kFALSE),
 116   fWriteVariableTree(writeVariableTree),
 117   fVariablesTree(0),
 118   fCandidateVariables(),
 119   fVtx1(0),
 120   fBzkG(0),
 121   fAdditionalChecks(additionalChecks)
 122 {
 123   //
 124   // Constructor. Initialization of Inputs and Outputs
 125   //
 126   Info("AliAnalysisTaskSELc2V0bachelor","Calling Constructor");
 127
 128   DefineOutput(1,TList::Class());  //conters
 129   DefineOutput(2,AliNormalizationCounter::Class());
 130   DefineOutput(3,AliRDHFCutsLctoV0::Class());
 131   if (!writeVariableTree) {
 132     DefineOutput(4,TList::Class());  //All Entries output
 133     DefineOutput(5,TList::Class());  //3sigma PID output
 134   } else {
 135     // Output slot #4 keeps a tree of the candidate variables after track selection
 136     DefineOutput(4,TTree::Class());  //My private output
 137   }
 138
 139 }
 140
 141 //___________________________________________________________________________
 142 AliAnalysisTaskSELc2V0bachelor::~AliAnalysisTaskSELc2V0bachelor() {
 143   //
 144   // destructor
 145   //
 146   Info("~AliAnalysisTaskSELc2V0bachelor","Calling Destructor");
 147
 148   if (fOutput) {
 149     delete fOutput;
 150     fOutput = 0;
 151   }
 152
 153   if (fOutputAll) {
 154     delete fOutputAll;
 155     fOutputAll = 0;
 156   }
 157
 158   if (fOutputPIDBach) {
 159     delete fOutputPIDBach;
 160     fOutputPIDBach = 0;
 161   }
 162
 163   if (fCounter) {
 164     delete fCounter;
 165     fCounter = 0;
 166   }
 167
 168   if (fAnalCuts) {
 169     delete fAnalCuts;
 170     fAnalCuts = 0;
 171   }
 172
 173   if (fVariablesTree) {
 174     delete fVariablesTree;
 175     fVariablesTree = 0;
 176   }
 177
 178 }
 179 //_________________________________________________
 180 void AliAnalysisTaskSELc2V0bachelor::Init() {
 181   //
 182   // Initialization
 183   //
 184
 185   fIsEventSelected=kFALSE;
 186
 187   if (fDebug > 1) AliInfo("Init");
 188
 189   PostData(3,fAnalCuts);
 190
 191   return;
 192 }
 193
 194 //_________________________________________________
 195 void AliAnalysisTaskSELc2V0bachelor::UserExec(Option_t *)
 196 {
 197   // user exec
 198   if (!fInputEvent) {
 199     AliError("NO EVENT FOUND!");
 200     return;
 201   }
 202
 203   AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
 204   TClonesArray *arrayLctopKos=0;
 205
 206   if (!aodEvent && AODEvent() && IsStandardAOD()) {
 207     // In case there is an AOD handler writing a standard AOD, use the AOD
 208     // event in memory rather than the input (ESD) event.
 209     aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
 210     // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
 211     // have to taken from the AOD event hold by the AliAODExtension
 212     AliAODHandler* aodHandler = (AliAODHandler*)
 213       ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
 214
 215     if (aodHandler->GetExtensions()) {
 216       AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
 217       AliAODEvent *aodFromExt = ext->GetAOD();
 218       arrayLctopKos=(TClonesArray*)aodFromExt->GetList()->FindObject("CascadesHF");
 219     }
 220   } else {
 221     arrayLctopKos=(TClonesArray*)aodEvent->GetList()->FindObject("CascadesHF");
 222   }
 223
 224   fCEvents->Fill(1);
 225
 226   if (fUseMCInfo)
 227     fAnalCuts->SetTriggerClass("");
 228
 229   // AOD primary vertex
 230   fVtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex();
 231   if (!fVtx1) return;
 232
 233   fIsEventSelected = fAnalCuts->IsEventSelected(aodEvent); // better to initialize before CheckEventSelection call
 234
 235   CheckEventSelection(aodEvent);
 236
 237
 238   // fix for temporary bug in ESDfilter
 239   fBzkG = (Double_t)aodEvent->GetMagneticField();
 240   if (TMath::Abs(fBzkG)<0.001) return;
 241   fCEvents->Fill(2);
 242
 243   if (!arrayLctopKos) {
 244     AliInfo("Could not find array of HF cascades, skipping the event");
 245     return;
 246   } else {
 247     if (arrayLctopKos->GetEntriesFast()) {
 248       AliInfo(Form("Found %d cascades",arrayLctopKos->GetEntriesFast()));
 249     }
 250   }
 251   fCEvents->Fill(3);
 252
 253   // mc analysis
 254   TClonesArray *mcArray = 0;
 255   AliAODMCHeader *mcHeader=0;
 256
 257   if (fUseMCInfo) {
 258     // MC array need for maching
 259     mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
 260     if (!mcArray) {
 261       AliError("Could not find Monte-Carlo in AOD");
 262       return;
 263     }
 264     fCEvents->Fill(4); // in case of MC events
 265
 266     // load MC header
 267     mcHeader = (AliAODMCHeader*)aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 268     if (!mcHeader) {
 269       AliError("AliAnalysisTaskSELc2V0bachelor::UserExec: MC header branch not found!\n");
 270       return;
 271     }
 272     fCEvents->Fill(5); // in case of MC events
 273
 274     Double_t zMCVertex = mcHeader->GetVtxZ();
 275     if (TMath::Abs(zMCVertex) > fAnalCuts->GetMaxVtxZ()) {
 276       AliDebug(2,Form("Event rejected: abs(zVtxMC)=%f > fAnalCuts->GetMaxVtxZ()=%f",zMCVertex,fAnalCuts->GetMaxVtxZ()));
 277       return;
 278     } else {
 279       fCEvents->Fill(17); // in case of MC events
 280     }
 281   }
 282
 283   fCounter->StoreEvent(aodEvent,fAnalCuts,fUseMCInfo); // it is very important that it stays BEFORE any other event selection
 284
 285   if (fVtx1->GetNContributors()>0) // this check is done in IsEventSelected
 286     fCEvents->Fill(6);
 287
 288   if ( !fIsEventSelected ) return; // don't take into account not selected events
 289   fCEvents->Fill(7);
 290
 291   Int_t nSelectedAnal = 0;
 292   if (fIsK0SAnalysis) {
 293     MakeAnalysisForLc2prK0S(arrayLctopKos,mcArray, nSelectedAnal, fAnalCuts);
 294
 295     if (nSelectedAnal)
 296       CheckEventSelectionWithCandidates(aodEvent);
 297
 298   }
 299
 300   fCounter->StoreCandidates(aodEvent,nSelectedAnal,kTRUE);
 301   fCounter->StoreCandidates(aodEvent,nSelectedAnal,kFALSE);
 302
 303   PostData(1,fOutput);
 304   PostData(2,fCounter);
 305   if (!fWriteVariableTree) {
 306     PostData(4,fOutputAll);
 307     PostData(5,fOutputPIDBach);
 308   } else {
 309     PostData(4,fVariablesTree);
 310   }
 311
 312   fIsEventSelected=kFALSE;
 313
 314   return;
 315 }
 316
 317 //________________________________________ terminate ___________________________
 318 void AliAnalysisTaskSELc2V0bachelor::Terminate(Option_t*)
 319 {
 320   // The Terminate() function is the last function to be called during
 321   // a query. It always runs on the client, it can be used to present
 322   // the results graphically or save the results to file.
 323
 324   //AliInfo("Terminate","");
 325   AliAnalysisTaskSE::Terminate();
 326
 327   fOutput = dynamic_cast<TList*> (GetOutputData(1));
 328   if (!fOutput) {
 329     AliError("fOutput not available");
 330     return;
 331   }
 332
 333   //fCEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fCEvents"));
 334   if (!fWriteVariableTree) {
 335     fOutputAll = dynamic_cast<TList*> (GetOutputData(4));
 336     if (!fOutputAll) {
 337       AliError("fOutputAll not available");
 338       return;
 339     }
 340
 341     fOutputPIDBach = dynamic_cast<TList*> (GetOutputData(5));
 342     if (!fOutputPIDBach) {
 343       AliError("fOutputPIDBach not available");
 344       return;
 345     }
 346   }
 347
 348   return;
 349 }
 350 //___________________________________________________________________________
 351 void AliAnalysisTaskSELc2V0bachelor::UserCreateOutputObjects() {
 352   // output
 353   AliInfo(Form("CreateOutputObjects of task %s\n", GetName()));
 354
 355   fOutput = new TList();
 356   fOutput->SetOwner();
 357   fOutput->SetName("chist0");
 358   DefineGeneralHistograms(); // define general histograms
 359   PostData(1,fOutput);
 360
 361   fCounter = new AliNormalizationCounter("NormalizationCounter");
 362   fCounter->Init();
 363   PostData(2,fCounter);
 364
 365   if (!fWriteVariableTree) {
 366
 367     fOutputAll = new TList();
 368     fOutputAll->SetOwner();
 369     fOutputAll->SetName("listAll");
 370
 371     fOutputPIDBach = new TList();
 372     fOutputPIDBach->SetOwner();
 373     fOutputPIDBach->SetName("listPIDBach");
 374
 375     DefineAnalysisHistograms(); // define analysis histograms
 376
 377     PostData(4,fOutputAll);
 378     PostData(5,fOutputPIDBach);
 379   }
 380   else {
 381     DefineTreeVariables();
 382     PostData(4,fVariablesTree);
 383   }
 384
 385   return;
 386 }
 387
 388 //-------------------------------------------------------------------------------
 389 void AliAnalysisTaskSELc2V0bachelor::MakeAnalysisForLc2prK0S(TClonesArray *arrayLctopKos,
 390                                                              TClonesArray *mcArray,
 391                                                              Int_t &nSelectedAnal,
 392                                                              AliRDHFCutsLctoV0 *cutsAnal)
 393 {
 394
 395   // make the analysis
 396
 397   Int_t pdgCand = 4122;
 398   Int_t pdgDgLctoV0bachelor[2]={2212,310}; // always 1st bachelor, 2nd V0
 399   Int_t pdgDgV0toDaughters[2]={211,211};
 400
 401   // loop over cascades to search for candidates Lc->p+K0S
 402   Int_t nCascades= arrayLctopKos->GetEntriesFast();
 403   if (nCascades==0) {
 404     AliInfo("Could not find cascades, skipping the event");
 405     return;
 406   }
 407
 408   for (Int_t iLctopK0S = 0; iLctopK0S<nCascades; iLctopK0S++) {
 409
 410     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(0);
 411
 412     // Lc candidates and K0S from Lc
 413     AliAODRecoCascadeHF* lcK0Spr = dynamic_cast<AliAODRecoCascadeHF*>(arrayLctopKos->At(iLctopK0S));
 414     if (!lcK0Spr) {
 415       AliDebug(2,Form("Cascade %d doens't exist, skipping",iLctopK0S));
 416       continue;
 417     }
 418
 419     Bool_t unsetvtx=kFALSE;
 420     if (!lcK0Spr->GetOwnPrimaryVtx()) {
 421       lcK0Spr->SetOwnPrimaryVtx(fVtx1);
 422       unsetvtx=kTRUE;
 423     }
 424
 425     if (!lcK0Spr->GetSecondaryVtx()) {
 426       AliInfo("No secondary vertex"); // it will be done in AliRDHFCutsLctoV0::IsSelected
 427       continue;
 428     }
 429
 430     if (lcK0Spr->GetNDaughters()!=2) {
 431       AliDebug(2,Form("Cascade %d has not 2 daughters (nDaughters=%d)",iLctopK0S,lcK0Spr->GetNDaughters())); // it will be done in AliRDHFCutsLctoV0::IsSelected
 432       continue;
 433     }
 434
 435     AliAODv0 * v0part = dynamic_cast<AliAODv0*>(lcK0Spr->Getv0());
 436     AliAODTrack * bachPart = dynamic_cast<AliAODTrack*>(lcK0Spr->GetBachelor());
 437     if (!v0part || !bachPart) {
 438       AliDebug(2,Form("Cascade %d has no V0 or no bachelor object",iLctopK0S)); // it will be done in AliRDHFCutsLctoV0::IsSelected
 439       continue;
 440     }
 441
 442     if (!v0part->GetSecondaryVtx()) {
 443       AliDebug(2,Form("No secondary vertex for V0 by cascade %d",iLctopK0S)); // it will be done in AliRDHFCutsLctoV0::IsSelected
 444       continue;
 445     }
 446
 447     if (v0part->GetNDaughters()!=2) {
 448       AliDebug(2,Form("current V0 has not 2 daughters (onTheFly=%d, nDaughters=%d)",v0part->GetOnFlyStatus(),v0part->GetNDaughters())); // it will be done in AliRDHFCutsLctoV0::IsSelected
 449       continue;
 450     }
 451
 452     AliAODTrack * v0Pos = dynamic_cast<AliAODTrack*>(lcK0Spr->Getv0PositiveTrack());
 453     AliAODTrack * v0Neg = dynamic_cast<AliAODTrack*>(lcK0Spr->Getv0NegativeTrack());
 454     if (!v0Neg || !v0Pos) {
 455       AliDebug(2,Form("V0 by cascade %d has no V0positive of V0negative object",iLctopK0S)); // it will be done in AliRDHFCutsLctoV0::IsSelected
 456       continue;
 457     }
 458
 459     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(1);
 460
 461     if (v0Pos->Charge() == v0Neg->Charge()) continue;
 462
 463     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(2);
 464
 465     Int_t isLc = 0;
 466
 467     if (fUseMCInfo) {
 468
 469       Int_t pdgCode=-2;
 470
 471       // find associated MC particle for Lc -> p+K0 and K0S->pi+pi
 472       Int_t mcLabel = lcK0Spr->MatchToMC(pdgCand,pdgDgLctoV0bachelor[1],pdgDgLctoV0bachelor,pdgDgV0toDaughters,mcArray,kTRUE);
 473       if (mcLabel!=-1) {
 474         AliDebug(2,Form(" ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~cascade number %d (total cascade number = %d)", iLctopK0S,nCascades));
 475
 476         AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel));
 477         if (partLc) {
 478           pdgCode = partLc->GetPdgCode();
 479           if (pdgCode<0) AliDebug(2,Form(" ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ MClabel=%d ~~~~~~~~~~ pdgCode=%d", mcLabel, pdgCode));
 480           pdgCode = TMath::Abs(pdgCode);
 481           isLc = 1;
 482         }
 483       } else {
 484         AliDebug(2,Form("No MC candidate (cascade number %d -total cascade number = %d -)", iLctopK0S,nCascades));
 485         pdgCode=-1;
 486       }
 487     }
 488
 489     FillLc2pK0Sspectrum(lcK0Spr, isLc,
 490                         nSelectedAnal, cutsAnal,
 491                         mcArray);
 492
 493     if (unsetvtx) lcK0Spr->UnsetOwnPrimaryVtx();
 494
 495   }
 496
 497   AliDebug(2, Form("Found %d Reco particles that are Lc!!", nSelectedAnal));
 498
 499   return;
 500 }
 501
 502 //________________________________________________________________________
 503 void AliAnalysisTaskSELc2V0bachelor::FillLc2pK0Sspectrum(AliAODRecoCascadeHF *part,
 504                                                          Int_t isLc,
 505                                                          Int_t &nSelectedAnal,
 506                                                          AliRDHFCutsLctoV0 *cutsAnal,
 507                                                          TClonesArray *mcArray)
 508 {
 509   //
 510   // Fill histos for Lc -> K0S+proton
 511   //
 512
 513   TString fillthis="";
 514
 515   AliAODTrack *bachelor = (AliAODTrack*)part->GetBachelor();
 516   Double_t momBach  = bachelor->P();
 517
 518   AliAODv0 * v0part = (AliAODv0*)part->Getv0();
 519   Bool_t onFlyV0 = v0part->GetOnFlyStatus(); // on-the-flight V0s
 520
 521   Bool_t areCutsUsingPID = cutsAnal->GetIsUsePID();
 522   cutsAnal->SetUsePID(kFALSE);
 523   Bool_t isInCascadeWindow = (((cutsAnal->IsSelectedSingleCut(part,AliRDHFCuts::kCandidate,0))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)); // cut on Lc->p+K0S invMass
 524   cutsAnal->SetUsePID(areCutsUsingPID);
 525
 526   //if ( !( !onFlyV0 || (onFlyV0 && fUseOnTheFlyV0) ) ) return;
 527   if ( onFlyV0 && !fUseOnTheFlyV0 ) return;
 528
 529   if (fAdditionalChecks) CheckCandidatesAtDifferentLevels(part,cutsAnal);
 530
 531   if ( !(cutsAnal->IsInFiducialAcceptance(part->Pt(),part->Y(4122))) ) return;
 532
 533   if ( !( ( (cutsAnal->IsSelected(part,AliRDHFCuts::kTracks))&(AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr) ) ) return;
 534
 535   if ( ( ( (cutsAnal->IsSelected(part,AliRDHFCuts::kAll))&(AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr) ) ) nSelectedAnal++;
 536
 537   // Fill candidate variable Tree (track selection, V0 invMass selection)
 538   if ( fWriteVariableTree ) {
 539     Double_t invmassK0S = v0part->MassK0Short();
 540     Double_t mk0sPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
 541     if ( !onFlyV0 && isInCascadeWindow && part->CosV0PointingAngle()>0.99 && TMath::Abs(invmassK0S-mk0sPDG)<=0.05)
 542       FillTheTree(part,cutsAnal,mcArray,isLc);
 543     return;
 544   }
 545
 546   cutsAnal->SetUsePID(kFALSE);
 547   Bool_t isCandidateSelectedCuts = (((cutsAnal->IsSelected(part,AliRDHFCuts::kCandidate))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)); // kinematic/topological cuts
 548   cutsAnal->SetUsePID(areCutsUsingPID);
 549   Bool_t isBachelorID = (((cutsAnal->IsSelected(part,AliRDHFCuts::kPID))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)); // ID x bachelor
 550
 551   //if (((cutsAnal->IsSelected(part,AliRDHFCuts::kCandidate))&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr) {
 552   if (((cutsAnal->IsSelected(part,AliRDHFCuts::kAll))&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr) {
 553     if (fUseMCInfo && isLc && !fWriteVariableTree) {
 554       Int_t pdgCand1 = 4122;
 555       Int_t pdgDgLctoV0bachelor1[2]={2212,310};
 556       Int_t pdgDgV0toDaughters1[2]={211,211};
 557       Int_t mcLabel1=part->MatchToMC(pdgCand1,pdgDgLctoV0bachelor1[1],pdgDgLctoV0bachelor1,pdgDgV0toDaughters1,mcArray,kTRUE);
 558       AliDebug(2,Form(" Found true MC candidate: Lc->pK0S(%d) - onTheFly=%1d",mcLabel1,onFlyV0));
 559     }
 560   }
 561
 562   Double_t nSigmaTPCpr=-999.;
 563   cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,4,nSigmaTPCpr);
 564   Double_t nSigmaTOFpr=-999.;
 565   cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,4,nSigmaTOFpr);
 566
 567   Double_t nSigmaTPCpi=-999.;
 568   cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,2,nSigmaTPCpi);
 569   Double_t nSigmaTOFpi=-999.;
 570   cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,2,nSigmaTOFpi);
 571
 572   Double_t nSigmaTPCka=-999.;
 573   cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,3,nSigmaTPCka);
 574   Double_t nSigmaTOFka=-999.;
 575   cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,3,nSigmaTOFka);
 576
 577   //if (onFlyV0 && fUseOnTheFlyV0) {
 578   if (onFlyV0) {
 579
 580     if (isCandidateSelectedCuts)
 581       FillAnalysisHistograms(part,isBachelorID,"");
 582
 583   }
 584   //else if (!onFlyV0) {
 585   else {
 586
 587     if (isCandidateSelectedCuts) {
 588
 589       fillthis="histoprotonBachSigmaVspTOF";
 590       ((TH2F*)(fOutput->FindObject(fillthis)))->Fill(momBach,nSigmaTOFpr);
 591       fillthis="histoprotonBachSigmaVspTPC";
 592       ((TH2F*)(fOutput->FindObject(fillthis)))->Fill(momBach,nSigmaTPCpr);
 593
 594       FillAnalysisHistograms(part,isBachelorID,"Offline");
 595
 596     }
 597
 598   }
 599
 600
 601   if (fUseMCInfo) {
 602     if (isLc==1) {
 603       //if (onFlyV0 && fUseOnTheFlyV0) {
 604       if (onFlyV0) {
 605
 606         if (isCandidateSelectedCuts)
 607           FillAnalysisHistograms(part,isBachelorID,"Sgn");
 608
 609       }
 610       //else if (!onFlyV0) {
 611       else {
 612
 613         if (isCandidateSelectedCuts) {
 614
 615           fillthis="histoprotonBachSigmaVspTOFsgn";
 616           ((TH2F*)(fOutput->FindObject(fillthis)))->Fill(momBach,nSigmaTOFpr);
 617           fillthis="histoprotonBachSigmaVspTPCsgn";
 618           ((TH2F*)(fOutput->FindObject(fillthis)))->Fill(momBach,nSigmaTPCpr);
 619
 620           FillAnalysisHistograms(part,isBachelorID,"OfflineSgn");
 621
 622         }
 623
 624       }
 625
 626     }// sgn
 627     else { // bkg
 628       //if (onFlyV0 && fUseOnTheFlyV0) {
 629       if (onFlyV0) {
 630
 631         if (isCandidateSelectedCuts)
 632           FillAnalysisHistograms(part,isBachelorID,"Bkg");
 633
 634       }
 635       //else if (!onFlyV0) {
 636       else {
 637
 638         if (isCandidateSelectedCuts) {
 639
 640           fillthis="histoprotonBachSigmaVspTOFbkg";
 641           ((TH2F*)(fOutput->FindObject(fillthis)))->Fill(momBach,nSigmaTOFpr);
 642           fillthis="histoprotonBachSigmaVspTPCbkg";
 643           ((TH2F*)(fOutput->FindObject(fillthis)))->Fill(momBach,nSigmaTPCpr);
 644
 645           FillAnalysisHistograms(part,isBachelorID,"OfflineBkg");
 646         }
 647
 648       }
 649
 650     }
 651   } // if fUseMCInfo
 652
 653   return;
 654 }
 655
 656 //----------------------------------------------------
 657 void AliAnalysisTaskSELc2V0bachelor::DefineK0SHistos()
 658 {
 659
 660   TString nameMass=" ", nameSgn=" ", nameBkg=" ";
 661
 662   Double_t mLcPDG  = TDatabasePDG::Instance()->GetParticle(4122)->Mass();
 663   Double_t mK0SPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
 664
 665   if (fUseOnTheFlyV0) {
 666
 667     // V0 invariant masses (on-the-fly)
 668     nameMass="histK0SMass";
 669     TH2F* spectrumK0SMass = new TH2F(nameMass.Data(),"K^{0}_{S} invariant mass VS p_{T}; M(#pi^{+}#pi^{-}) [GeV/c^{2}]; p_{T}(K^{0}_{S}) [GeV/c]; Entries",
 670                                     1000,mK0SPDG-0.050,mK0SPDG+0.050,175,0.,35.);
 671
 672     // Lc invariant masses (x K0S on-the-fly)
 673     nameMass="histLcMassByK0S";
 674     TH2F* spectrumLcMassByK0S = new TH2F(nameMass.Data(),"#Lambda_{c} invariant mass (by K^{0}_{S}) vs p_{T} ; m_{inv}(p-K^{0}_{S}) [GeV/c^{2}]; p_{T}(#Lambda_{c}) [GeV/c]",
 675                                          1000,mLcPDG-0.250,mLcPDG+0.250,175,0.,35.);
 676
 677     nameMass="histpK0Svsp";
 678     TH2F* momentumDistributionK0Svsp = new TH2F(nameMass.Data(),"#Lambda_{c}: p(K^{0}_{S}) vs p(p);  p_{p}; p_{K^{0}_{S}}  ",
 679                                                 175,0.,35.,175,0.,35.);
 680
 681     nameMass="histArmPodK0S";
 682     TH2F* armenterosPodK0S = new TH2F(nameMass.Data(),"K^{0}_{S}: Armenteros-Podolanski distribution; #frac{p_{L}^{+}-p_{L}^{-}}{p_{L}^{+}+p_{L}^{-}}; p_{T}^{+} [GeV/c]",
 683                                       200,-1.,1.,300,0.,0.3);
 684
 685     nameMass="histDCAtoPVvspK0S";
 686     TH2F *dcatoPVvspK0S = new TH2F(nameMass.Data(),"K^{0}_{S}: DCA to Primary Vertex vs K^{0}_{S} momentum ; p(K^{0}_{S}) [GeV/c]; DCA to Primary Vertex [n#sigmas]; Entries",
 687                                    175,0.,35.,50,0.,5.);
 688
 689     nameMass="histK0ScosPAwrtPVvspK0S";
 690     TH2F *cosPAwrtPVvspK0S = new TH2F(nameMass.Data(),"K^{0}_{S}: cosine of pointing angle wrt primary vertex vs K^{0}_{S} momentum ; p(K^{0}_{S}) [GeV/c]; cosine; Entries",
 691                                       175,0.,35.,100,0.99,1.);
 692
 693     TH2F* allspectrumK0SMass = (TH2F*)spectrumK0SMass->Clone();
 694     TH2F* allspectrumLcMassByK0S    = (TH2F*)spectrumLcMassByK0S->Clone();
 695     TH2F* allmomentumDistributionK0Svsp= (TH2F*)momentumDistributionK0Svsp->Clone();
 696     TH2F* alldcatoPVvspK0S=(TH2F*)dcatoPVvspK0S->Clone();
 697     TH2F* allcosV0PAwrtPVvspK0S=(TH2F*)cosPAwrtPVvspK0S->Clone();
 698
 699     TH2F* pidBachspectrumK0SMass = (TH2F*)spectrumK0SMass->Clone();
 700     TH2F* pidBachspectrumLcMassByK0S    = (TH2F*)spectrumLcMassByK0S->Clone();
 701     TH2F* pidBachmomentumDistributionK0Svsp= (TH2F*)momentumDistributionK0Svsp->Clone();
 702     TH2F* pidBachdcatoPVvspK0S=(TH2F*)dcatoPVvspK0S->Clone();
 703     TH2F* pidBachcosV0PAwrtPVvspK0S=(TH2F*)cosPAwrtPVvspK0S->Clone();
 704
 705     TH2F* allArmenterosPodK0S = (TH2F*)armenterosPodK0S->Clone();
 706     TH2F* pidBachArmenterosPodK0S = (TH2F*)armenterosPodK0S->Clone();
 707
 708     fOutputAll->Add(allspectrumK0SMass);
 709     fOutputAll->Add(allspectrumLcMassByK0S);
 710     fOutputAll->Add(allmomentumDistributionK0Svsp);
 711     fOutputAll->Add(allArmenterosPodK0S);
 712     fOutputAll->Add(alldcatoPVvspK0S);
 713     fOutputAll->Add(allcosV0PAwrtPVvspK0S);
 714
 715     fOutputPIDBach->Add(pidBachspectrumK0SMass);
 716     fOutputPIDBach->Add(pidBachspectrumLcMassByK0S);
 717     fOutputPIDBach->Add(pidBachmomentumDistributionK0Svsp);
 718     fOutputPIDBach->Add(pidBachArmenterosPodK0S);
 719     fOutputPIDBach->Add(pidBachdcatoPVvspK0S);
 720     fOutputPIDBach->Add(pidBachcosV0PAwrtPVvspK0S);
 721
 722   }
 723
 724   // V0 invariant masses (offline)
 725   nameMass="histK0SMassOffline";
 726   TH2F* spectrumK0SMassOffline = new TH2F(nameMass.Data(),"K^{0}_{S} invariant mass VS p_{T}; M(#pi^{+}#pi^{-}) [GeV/c^{2}]; p_{T}(K^{0}_{S}) [GeV/c]; Entries",
 727                                          1000,mK0SPDG-0.050,mK0SPDG+0.050,175,0.,35.);
 728
 729   // Lc invariant masses (x K0S offline)
 730   nameMass="histLcMassByK0SOffline";
 731   TH2F* spectrumLcMassOfflineByK0S = new TH2F(nameMass.Data(),"#Lambda_{c} invariant mass (by K^{0}_{S}) vs p_{T}; M(K^{0}_{S}p) [GeV/c^{2}]; p_{T} [GeV/c]",
 732                                               1000,mLcPDG-0.250,mLcPDG+0.250,175,0.,35.);
 733
 734   nameMass="histpK0SvspOffline";
 735   TH2F* momentumDistributionK0SvspOffline = new TH2F(nameMass.Data(),"#Lambda_{c}: p(K^{0}_{S}) vs p(p) - Offline ;  p_{p} [GeV/c]; p_{K^{0}_{S}} [GeV/c]",
 736                                                      175,0.,35.,175,0.,35.);
 737
 738   nameMass="histArmPodK0SOffline";
 739   TH2F* armenterosPodK0SOff = new TH2F(nameMass.Data(),"K^{0}_{S}  Armenteros-Podolanski distribution - Offline; #frac{p_{L}^{+}-p_{L}^{-}}{p_{L}^{+}+p_{L}^{-}}; p_{T}^{+} [GeV/c]",
 740                                        200,-1.,1.,300,0.,0.3);
 741
 742   nameMass="histDCAtoPVvspK0SOffline";
 743   TH2F *dcatoPVvspK0SOffline = new TH2F(nameMass.Data(),"K^{0}_{S}: DCA to Primary Vertex vs  K^{0}_{S} invariant mass - Offline; p(K^{0}_{S}) [GeV/c]; DCA to Primary Vertex [n#sigmas]; Entries",
 744                                         175,0.,35.,50,0.,5.);
 745
 746   nameMass="histK0ScosPAwrtPVvspK0SOffline";
 747   TH2F *cosPAwrtPVvspK0SOffline = new TH2F(nameMass.Data(),"K^{0}_{S}: cosine of pointing angle wrt primary vertex vs K^{0}_{S} momentum - Offline; p(K^{0}_{S}) [GeV/c]; cosine; Entries",
 748                                            175,0.,35.,100,0.99,1.);
 749
 750
 751
 752   TH2F* allspectrumK0SMassOffline = (TH2F*)spectrumK0SMassOffline->Clone();
 753   TH2F* allspectrumLcMassOfflineByK0S    = (TH2F*)spectrumLcMassOfflineByK0S->Clone();
 754   TH2F* allmomentumDistributionK0SvspOffline= (TH2F*)momentumDistributionK0SvspOffline->Clone();
 755   TH2F* alldcatoPVvspK0SOffline=(TH2F*)dcatoPVvspK0SOffline->Clone();
 756   TH2F* allcosPAwrtPVvspK0SOffline=(TH2F*)cosPAwrtPVvspK0SOffline->Clone();
 757   TH2F* allArmenterosPodK0SOff = (TH2F*)armenterosPodK0SOff->Clone();
 758
 759   TH2F* pidBachspectrumK0SMassOffline = (TH2F*)spectrumK0SMassOffline->Clone();
 760   TH2F* pidBachspectrumLcMassOfflineByK0S    = (TH2F*)spectrumLcMassOfflineByK0S->Clone();
 761   TH2F* pidBachmomentumDistributionK0SvspOffline= (TH2F*)momentumDistributionK0SvspOffline->Clone();
 762   TH2F* pidBachdcatoPVvspK0SOffline=(TH2F*)dcatoPVvspK0SOffline->Clone();
 763   TH2F* pidBachcosPAwrtPVvspK0SOffline=(TH2F*)cosPAwrtPVvspK0SOffline->Clone();
 764   TH2F* pidBachArmenterosPodK0SOff = (TH2F*)armenterosPodK0SOff->Clone();
 765
 766
 767   fOutputAll->Add(allspectrumK0SMassOffline);
 768   fOutputAll->Add(allspectrumLcMassOfflineByK0S);
 769   fOutputAll->Add(allmomentumDistributionK0SvspOffline);
 770   fOutputAll->Add(allArmenterosPodK0SOff);
 771   fOutputAll->Add(alldcatoPVvspK0SOffline);
 772   fOutputAll->Add(allcosPAwrtPVvspK0SOffline);
 773
 774   fOutputPIDBach->Add(pidBachspectrumK0SMassOffline);
 775   fOutputPIDBach->Add(pidBachspectrumLcMassOfflineByK0S);
 776   fOutputPIDBach->Add(pidBachmomentumDistributionK0SvspOffline);
 777   fOutputPIDBach->Add(pidBachArmenterosPodK0SOff);
 778   fOutputPIDBach->Add(pidBachdcatoPVvspK0SOffline);
 779   fOutputPIDBach->Add(pidBachcosPAwrtPVvspK0SOffline);
 780
 781   if (fUseMCInfo) {
 782
 783     if (fUseOnTheFlyV0) {
 784
 785       nameSgn="histK0SMassSgn";
 786       nameBkg="histK0SMassBkg";
 787       TH2F* spectrumK0SMassSgn = new TH2F(nameSgn.Data(), "K^{0}_{S} - sgn: invariant mass VS p_{T} - MC; m_{inv}(#pi^{+}#pi^{-}) [GeV/c^{2}]; p_{T}(K^{0}_{S}) [GeV/c]; Entries",
 788                                           1000,mK0SPDG-0.050,mK0SPDG+0.050,175,0.,35.);
 789       TH2F* spectrumK0SMassBkg = new TH2F(nameBkg.Data(), "K^{0}_{S} - bkg: invariant mass VS p_{T} - MC; m_{inv}(#pi^{+}#pi^{-}) [GeV/c^{2}]; p_{T}(K^{0}_{S}) [GeV/c]; Entries",
 790                                           1000,mK0SPDG-0.050,mK0SPDG+0.050,175,0.,35.);
 791
 792       TH2F* allspectrumK0SMassSgn = (TH2F*)spectrumK0SMassSgn->Clone();
 793       TH2F* allspectrumK0SMassBkg = (TH2F*) spectrumK0SMassBkg->Clone();
 794       TH2F* pidBachspectrumK0SMassSgn = (TH2F*)spectrumK0SMassSgn->Clone();
 795       TH2F* pidBachspectrumK0SMassBkg = (TH2F*) spectrumK0SMassBkg->Clone();
 796       fOutputAll->Add(allspectrumK0SMassSgn);
 797       fOutputAll->Add(allspectrumK0SMassBkg);
 798       fOutputPIDBach->Add(pidBachspectrumK0SMassSgn);
 799       fOutputPIDBach->Add(pidBachspectrumK0SMassBkg);
 800
 801
 802       nameSgn="histLcMassByK0SSgn";
 803       nameBkg="histLcMassByK0SBkg";
 804       TH2F* spectrumLcMassByK0SSgn = new TH2F(nameSgn.Data(), "#Lambda_{c} - sgn: invariant mass (by K^{0}_{S}) vs p_{T}  - MC; m_{inv}(p-K^{0}_{S}) [GeV/c^{2}];  p_{T}",
 805                                               1000,mLcPDG-0.250,mLcPDG+0.250,175,0.,35.);
 806       TH2F* spectrumLcMassByK0SBkg = new TH2F(nameBkg.Data(), "#Lambda_{c} - bkg: invariant mass (by K^{0}_{S}) vs p_{T}  - MC; m_{inv}(p-K^{0}_{S}) [GeV/c^{2}]; p_{T}",
 807                                               1000,mLcPDG-0.250,mLcPDG+0.250,175,0.,35.);
 808
 809       TH2F* allspectrumLcMassByK0SSgn = (TH2F*)spectrumLcMassByK0SSgn->Clone();
 810       TH2F* allspectrumLcMassByK0SBkg = (TH2F*) spectrumLcMassByK0SBkg->Clone();
 811       TH2F* pidBachspectrumLcMassByK0SSgn = (TH2F*)spectrumLcMassByK0SSgn->Clone();
 812       TH2F* pidBachspectrumLcMassByK0SBkg = (TH2F*) spectrumLcMassByK0SBkg->Clone();
 813       fOutputAll->Add(allspectrumLcMassByK0SSgn);
 814       fOutputAll->Add(allspectrumLcMassByK0SBkg);
 815       fOutputPIDBach->Add(pidBachspectrumLcMassByK0SSgn);
 816       fOutputPIDBach->Add(pidBachspectrumLcMassByK0SBkg);
 817
 818       nameSgn="histpK0SvspSgn";
 819       nameBkg="histpK0SvspBkg";
 820       TH2F* momentumDistributionK0SvspSgn= new TH2F(nameSgn.Data(),"#Lambda_{c} - sgn: K^{0}_{S} vs p Total Momentum Distribution - MC; p_{p}; p_{K^{0}_{S}}",
 821                                                     175,0.,35.,175,0.,35.);
 822       TH2F* momentumDistributionK0SvspBkg= new TH2F(nameBkg.Data(),"#Lambda_{c} - bkg: K^{0}_{S} vs p Total Momentum Distribution - MC; p_{p}; p_{K^{0}_{S}}",
 823                                                     175,0.,35.,175,0.,35.);
 824
 825       TH2F* allmomentumDistributionK0SvspSgn= (TH2F*)momentumDistributionK0SvspSgn->Clone();
 826       TH2F* allmomentumDistributionK0SvspBkg= (TH2F*)momentumDistributionK0SvspBkg->Clone();
 827       TH2F* pidBachmomentumDistributionK0SvspSgn= (TH2F*)momentumDistributionK0SvspSgn->Clone();
 828       TH2F* pidBachmomentumDistributionK0SvspBkg= (TH2F*)momentumDistributionK0SvspBkg->Clone();
 829       fOutputAll->Add(allmomentumDistributionK0SvspSgn);
 830       fOutputAll->Add(allmomentumDistributionK0SvspBkg);
 831       fOutputPIDBach->Add(pidBachmomentumDistributionK0SvspSgn);
 832       fOutputPIDBach->Add(pidBachmomentumDistributionK0SvspBkg);
 833
 834
 835       // armenteros-podolanski plots K0S
 836       nameSgn="histArmPodK0SSgn";
 837       nameBkg="histArmPodK0SBkg";
 838       TH2F* armenterosPodK0SSgn = new TH2F(nameSgn.Data(),"K^{0}_{S}  Armenteros-Podolanski distribution (sgn); #frac{p_{L}^{+}-p_{L}^{-}}{p_{L}^{+}+p_{L}^{-}}; p_{T}^{+} [GeV/c]",
 839                                            200,-1.,1.,300,0.,0.3);
 840       TH2F* armenterosPodK0SBkg = new TH2F(nameBkg.Data(),"K^{0}_{S}  Armenteros-Podolanski distribution (bkg); #frac{p_{L}^{+}-p_{L}^{-}}{p_{L}^{+}+p_{L}^{-}}; p_{T}^{+} [GeV/c]",
 841                                            200,-1.,1.,300,0.,0.3);
 842
 843       TH2F* allArmenterosPodK0SSgn = (TH2F*)armenterosPodK0SSgn->Clone();
 844       TH2F* allArmenterosPodK0SBkg = (TH2F*)armenterosPodK0SBkg->Clone();
 845       TH2F* pidBachArmenterosPodK0SSgn = (TH2F*)armenterosPodK0SSgn->Clone();
 846       TH2F* pidBachArmenterosPodK0SBkg = (TH2F*)armenterosPodK0SBkg->Clone();
 847
 848       fOutputAll->Add(allArmenterosPodK0SSgn);
 849       fOutputAll->Add(allArmenterosPodK0SBkg);
 850       fOutputPIDBach->Add(pidBachArmenterosPodK0SSgn);
 851       fOutputPIDBach->Add(pidBachArmenterosPodK0SBkg);
 852
 853
 854       nameSgn="histDCAtoPVvspK0SSgn";
 855       nameBkg="histDCAtoPVvspK0SBkg";
 856       TH2F *dcatoPVvspK0SSgn=new TH2F(nameSgn.Data(),"K^{0}_{S} - sgn: DCA to Primary Vertex vs  K^{0}_{S} invariant mass (sgn); p(K^{0}_{S}) [GeV/c]; DCA to Primary Vertex [n#sigmas]; Entries",175,0.,35.,50,0.,5.);
 857       TH2F *dcatoPVvspK0SBkg=new TH2F(nameBkg.Data(),"K^{0}_{S} - bkg: DCA to Primary Vertex vs  K^{0}_{S} invariant mass (bkg); p(K^{0}_{S}) [GeV/c]; DCA to Primary Vertex [n#sigmas]; Entries",175,0.,35.,50,0.,5.);
 858
 859       TH2F* alldcatoPVvspK0SSgn= (TH2F*)dcatoPVvspK0SSgn->Clone();
 860       TH2F* alldcatoPVvspK0SBkg= (TH2F*)dcatoPVvspK0SBkg->Clone();
 861       TH2F* pidBachdcatoPVvspK0SSgn= (TH2F*)dcatoPVvspK0SSgn->Clone();
 862       TH2F* pidBachdcatoPVvspK0SBkg= (TH2F*)dcatoPVvspK0SBkg->Clone();
 863
 864       fOutputAll->Add(alldcatoPVvspK0SSgn);
 865       fOutputPIDBach->Add(pidBachdcatoPVvspK0SSgn);
 866       fOutputAll->Add(alldcatoPVvspK0SBkg);
 867       fOutputPIDBach->Add(pidBachdcatoPVvspK0SBkg);
 868
 869     }
 870
 871
 872     nameSgn="histK0SMassOfflineSgn";
 873     nameBkg="histK0SMassOfflineBkg";
 874     TH2F* spectrumK0SMassOfflineSgn = new TH2F(nameSgn.Data(), "K^{0}_{S} - sgn: invariant mass VS p_{T} - MC; m_{inv}(#pi^{+}-#pi^{-}) [GeV/c^{2}]; p_{T}(K^{0}_{S}) [GeV/c]; Entries",
 875                                               1000,mK0SPDG-0.050,mK0SPDG+0.050,175,0.,35.);
 876     TH2F* spectrumK0SMassOfflineBkg = new TH2F(nameBkg.Data(), "K^{0}_{S} - bkg: invariant mass VS p_{T} - MC; m_{inv}(#pi^{+}-#pi^{-}) [GeV/c^{2}]; p_{T}(K^{0}_{S}) [GeV/c]; Entries",
 877                                               1000,mK0SPDG-0.050,mK0SPDG+0.050,175,0.,35.);
 878
 879     TH2F* allspectrumK0SMassOfflineSgn = (TH2F*)spectrumK0SMassOfflineSgn->Clone();
 880     TH2F* allspectrumK0SMassOfflineBkg = (TH2F*) spectrumK0SMassOfflineBkg->Clone();
 881     fOutputAll->Add(allspectrumK0SMassOfflineSgn);
 882     fOutputAll->Add(allspectrumK0SMassOfflineBkg);
 883
 884
 885     TH2F* pidBachspectrumK0SMassOfflineSgn = (TH2F*)spectrumK0SMassOfflineSgn->Clone();
 886     TH2F* pidBachspectrumK0SMassOfflineBkg = (TH2F*) spectrumK0SMassOfflineBkg->Clone();
 887     fOutputPIDBach->Add(pidBachspectrumK0SMassOfflineSgn);
 888     fOutputPIDBach->Add(pidBachspectrumK0SMassOfflineBkg);
 889
 890
 891     nameSgn="histLcMassByK0SOfflineSgn";
 892     nameBkg="histLcMassByK0SOfflineBkg";
 893     TH2F* spectrumLcMassOfflineByK0SSgn = new TH2F(nameSgn.Data(), "#Lambda_{c} - sgn: invariant mass (by K^{0}_{S})  vs p_{T} - MC; M(#Lambda_{c}) [GeV/c^{2}]; p_{T}",
 894                                                    1000,mLcPDG-0.250,mLcPDG+0.250,175,0.,35.);
 895     TH2F* spectrumLcMassOfflineByK0SBkg = new TH2F(nameBkg.Data(), "#Lambda_{c} - bkg: invariant mass (by K^{0}_{S})  vs p_{T} - MC; M(#Lambda_{c}) [GeV/c^{2}]; p_{T}",
 896                                                    1000,mLcPDG-0.250,mLcPDG+0.250,175,0.,35.);
 897
 898
 899     TH2F* allspectrumLcMassOfflineByK0SSgn = (TH2F*)spectrumLcMassOfflineByK0SSgn->Clone();
 900     TH2F* allspectrumLcMassOfflineByK0SBkg = (TH2F*) spectrumLcMassOfflineByK0SBkg->Clone();
 901     TH2F* pidBachspectrumLcMassOfflineByK0SSgn = (TH2F*)spectrumLcMassOfflineByK0SSgn->Clone();
 902     TH2F* pidBachspectrumLcMassOfflineByK0SBkg = (TH2F*) spectrumLcMassOfflineByK0SBkg->Clone();
 903     fOutputAll->Add(allspectrumLcMassOfflineByK0SSgn);
 904     fOutputAll->Add(allspectrumLcMassOfflineByK0SBkg);
 905     fOutputPIDBach->Add(pidBachspectrumLcMassOfflineByK0SSgn);
 906     fOutputPIDBach->Add(pidBachspectrumLcMassOfflineByK0SBkg);
 907
 908
 909     nameSgn="histpK0SvspOfflineSgn";
 910     nameBkg="histpK0SvspOfflineBkg";
 911     TH2F* momentumDistributionK0SvspOfflineSgn= new TH2F(nameSgn.Data(),"#Lambda_{c} - sgn: K^{0}_{S} vs p Total Momentum Distribution - Offline  - MC; p_{p}; p_{K^{0}_{S}}",
 912                                                          175,0.,35.,175,0.,35.);
 913     TH2F* momentumDistributionK0SvspOfflineBkg= new TH2F(nameBkg.Data(),"#Lambda_{c} - bkg: K^{0}_{S} vs p Total Momentum Distribution - Offline  - MC; p_{p}; p_{K^{0}_{S}}",
 914                                                          175,0.,35.,175,0.,35.);
 915
 916
 917     TH2F* allmomentumDistributionK0SvspOfflineSgn= (TH2F*)momentumDistributionK0SvspOfflineSgn->Clone();
 918     TH2F* allmomentumDistributionK0SvspOfflineBkg= (TH2F*)momentumDistributionK0SvspOfflineBkg->Clone();
 919     TH2F* pidBachmomentumDistributionK0SvspOfflineSgn= (TH2F*)momentumDistributionK0SvspOfflineSgn->Clone();
 920     TH2F* pidBachmomentumDistributionK0SvspOfflineBkg= (TH2F*)momentumDistributionK0SvspOfflineBkg->Clone();
 921     fOutputAll->Add(allmomentumDistributionK0SvspOfflineSgn);
 922     fOutputAll->Add(allmomentumDistributionK0SvspOfflineBkg);
 923     fOutputPIDBach->Add(pidBachmomentumDistributionK0SvspOfflineSgn);
 924     fOutputPIDBach->Add(pidBachmomentumDistributionK0SvspOfflineBkg);
 925
 926
 927
 928
 929
 930     // armenteros-podolanski plots K0S (offline)
 931     nameSgn="histArmPodK0SOfflineSgn";
 932     nameBkg="histArmPodK0SOfflineBkg";
 933     TH2F* armenterosPodK0SOffSgn = new TH2F(nameSgn.Data(),"K^{0}_{S}  Armenteros-Podolanski distribution (sgn) -offline-; #frac{p_{L}^{+}-p_{L}^{-}}{p_{L}^{+}+p_{L}^{-}}; p_{T}^{+} [GeV/c]",
 934                                             200,-1.,1.,300,0.,0.3);
 935     TH2F* armenterosPodK0SOffBkg = new TH2F(nameBkg.Data(),"K^{0}_{S}  Armenteros-Podolanski distribution (bkg) -offline-; #frac{p_{L}^{+}-p_{L}^{-}}{p_{L}^{+}+p_{L}^{-}}; p_{T}^{+} [GeV/c]",
 936                                             200,-1.,1.,300,0.,0.3);
 937
 938
 939     TH2F* allArmenterosPodK0SOffSgn = (TH2F*)armenterosPodK0SOffSgn->Clone();
 940     TH2F* allArmenterosPodK0SOffBkg = (TH2F*)armenterosPodK0SOffBkg->Clone();
 941     TH2F* pidBachArmenterosPodK0SOffSgn = (TH2F*)armenterosPodK0SOffSgn->Clone();
 942     TH2F* pidBachArmenterosPodK0SOffBkg = (TH2F*)armenterosPodK0SOffBkg->Clone();
 943
 944
 945     fOutputAll->Add(allArmenterosPodK0SOffSgn);
 946     fOutputAll->Add(allArmenterosPodK0SOffBkg);
 947     fOutputPIDBach->Add(pidBachArmenterosPodK0SOffSgn);
 948     fOutputPIDBach->Add(pidBachArmenterosPodK0SOffBkg);
 949
 950
 951     nameSgn="histDCAtoPVvspK0SOfflineSgn";
 952     nameBkg="histDCAtoPVvspK0SOfflineBkg";
 953     TH2F *dcatoPVvspK0SOfflineSgn=new TH2F(nameSgn.Data(),"K^{0}_{S} -offline - (sgn): DCA to Primary Vertex vs  K^{0}_{S} invariant mass; p(K^{0}_{S}) [GeV/c]; DCA to Primary Vertex [n#sigmas]; Entries",
 954                                            175,0.,35.,50,0.,5.);
 955     TH2F *dcatoPVvspK0SOfflineBkg=new TH2F(nameBkg.Data(),"K^{0}_{S} -offline - (bkg): DCA to Primary Vertex vs  K^{0}_{S} invariant mass; p(K^{0}_{S}) [GeV/c]; DCA to Primary Vertex [n#sigmas]; Entries",
 956                                            175,0.,35.,50,0.,5.);
 957
 958
 959     TH2F* alldcatoPVvspK0SOfflineSgn= (TH2F*)dcatoPVvspK0SOfflineSgn->Clone();
 960     TH2F* pidBachdcatoPVvspK0SOfflineSgn= (TH2F*)dcatoPVvspK0SOfflineSgn->Clone();
 961     TH2F* alldcatoPVvspK0SOfflineBkg= (TH2F*)dcatoPVvspK0SOfflineBkg->Clone();
 962     TH2F* pidBachdcatoPVvspK0SOfflineBkg= (TH2F*)dcatoPVvspK0SOfflineBkg->Clone();
 963
 964
 965
 966     fOutputAll->Add(alldcatoPVvspK0SOfflineSgn);
 967     fOutputPIDBach->Add(pidBachdcatoPVvspK0SOfflineSgn);
 968     fOutputAll->Add(alldcatoPVvspK0SOfflineBkg);
 969     fOutputPIDBach->Add(pidBachdcatoPVvspK0SOfflineBkg);
 970
 971   }
 972
 973   return;
 974 }
 975
 976 //---------------------------
 977 void AliAnalysisTaskSELc2V0bachelor::CheckEventSelection(AliAODEvent *aodEvent) {
 978   //
 979   // To fill control histograms
 980   //
 981
 982   TClonesArray *arrayLctopKos=0;
 983   if (!aodEvent){
 984     if(AODEvent() && IsStandardAOD()) {
 985       // In case there is an AOD handler writing a standard AOD, use the AOD
 986       // event in memory rather than the input (ESD) event.
 987       aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
 988       // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
 989       // have to taken from the AOD event hold by the AliAODExtension
 990       AliAODHandler* aodHandler = (AliAODHandler*)
 991         ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
 992
 993       if (aodHandler->GetExtensions()) {
 994         AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
 995         AliAODEvent *aodFromExt = ext->GetAOD();
 996         arrayLctopKos=(TClonesArray*)aodFromExt->GetList()->FindObject("CascadesHF");
 997       }
 998     }
 999   } else {
1000     arrayLctopKos=(TClonesArray*)aodEvent->GetList()->FindObject("CascadesHF");
1001   }
1002
1003   Float_t zVertex = fVtx1->GetZ();
1004   TString titleVtx=fVtx1->GetTitle();
1005
1006   if (TMath::Abs(fBzkG)>=0.001) {
1007
1008     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ2")))->Fill(zVertex);
1009
1010     if (arrayLctopKos) {
1011
1012       if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ3")))->Fill(zVertex);
1013
1014       // mc analysis
1015       TClonesArray *mcArray = 0;
1016       AliAODMCHeader *mcHeader=0;
1017
1018       if (fUseMCInfo) {
1019         // MC array need for maching
1020         mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1021
1022         if (mcArray) {
1023           // load MC header
1024           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ4")))->Fill(zVertex);
1025           mcHeader = (AliAODMCHeader*)aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName());
1026
1027           if (mcHeader && fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ5")))->Fill(zVertex);
1028
1029           // check on MC Lc Daughter
1030           if (fAdditionalChecks) {
1031             for (Int_t iii=0; iii<mcArray->GetEntries(); iii++)
1032               SearchLcDaughter(mcArray,iii);
1033           }
1034
1035         }
1036
1037       }
1038
1039       if (fVtx1->GetNContributors()>0) {
1040         if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ6")))->Fill(zVertex);
1041
1042         TString firedTriggerClasses = aodEvent->GetFiredTriggerClasses(); // trigger class
1043         ULong64_t fTriggerMask=AliVEvent::kAnyINT;
1044         Bool_t check1 = kFALSE;
1045         if ( !fUseMCInfo && // don't do for MC...
1046              (aodEvent->GetRunNumber()<136851 || aodEvent->GetRunNumber()>139517) ) { // ...and for PbPb 2010 data
1047           if ( !(firedTriggerClasses.Contains("CINT1")) ) {
1048             AliInfo(Form(" ======================== firedTriggerClasses.Data() = %s",firedTriggerClasses.Data()));
1049             fCEvents->Fill(8);
1050             if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ8")))->Fill(zVertex);
1051             check1 = kTRUE;
1052           }
1053         }
1054
1055         Bool_t isSelectedAAA = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerMask);
1056         if (!isSelectedAAA) {
1057           fCEvents->Fill(9);
1058           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ9")))->Fill(zVertex);
1059         }
1060
1061         if (!isSelectedAAA || check1) {
1062           fCEvents->Fill(16);
1063           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ16")))->Fill(zVertex);
1064         }
1065
1066         fTriggerMask=AliVEvent::kAny;
1067         Bool_t isSelectedBBB = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerMask);
1068         if (!isSelectedBBB) {
1069           fCEvents->Fill(10);
1070           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ10")))->Fill(zVertex);
1071         }
1072
1073         if (titleVtx.Contains("Z")) {
1074           fCEvents->Fill(11);
1075           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ11")))->Fill(zVertex);
1076         }
1077         else if (titleVtx.Contains("3D")) {
1078           fCEvents->Fill(12);
1079           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ12")))->Fill(zVertex);
1080         } else {
1081           fCEvents->Fill(13);
1082           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ13")))->Fill(zVertex);
1083         }
1084
1085         if (TMath::Abs(zVertex)<=fAnalCuts->GetMaxVtxZ()) {
1086           fCEvents->Fill(14);
1087           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ14")))->Fill(zVertex);
1088         }
1089
1090         if ( fIsEventSelected ) {
1091           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ7")))->Fill(zVertex);
1092         } else {
1093           fCEvents->Fill(15);
1094           if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ15")))->Fill(zVertex);
1095         }
1096
1097       } // nContributors>=1
1098     } // analysisArray exists
1099   } // magnetic field exists
1100
1101   return;
1102 }
1103
1104 //-----------------
1105 void AliAnalysisTaskSELc2V0bachelor::CheckEventSelectionWithCandidates(AliAODEvent *aodEvent) {
1106   //
1107   // To fill control histograms
1108   //
1109
1110   Float_t zVertex = fVtx1->GetZ();
1111   TString titleVtx=fVtx1->GetTitle();
1112   TString firedTriggerClasses = aodEvent->GetFiredTriggerClasses(); // trigger class
1113   ULong64_t fTriggerMask=AliVEvent::kAnyINT;
1114
1115   ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(6);
1116   if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ6a")))->Fill(zVertex);
1117
1118   Bool_t check1a = kFALSE;
1119   if ( !fUseMCInfo && // don't do for MC...
1120        (aodEvent->GetRunNumber()<136851 || aodEvent->GetRunNumber()>139517) ) { // ...and for PbPb 2010 data
1121     if ( !(firedTriggerClasses.Contains("CINT1")) ) {
1122       AliInfo(Form(" ======================== firedTriggerClasses.Data() = %s",firedTriggerClasses.Data()));
1123       ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(8);
1124       if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ8a")))->Fill(zVertex);
1125       check1a = kTRUE;
1126     }
1127   }
1128
1129   Bool_t isSelectedAAAa = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerMask);
1130   if (!isSelectedAAAa) {
1131     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(9);
1132     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ9a")))->Fill(zVertex);
1133   }
1134
1135   if (!isSelectedAAAa || check1a) {
1136     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(16);
1137     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ16a")))->Fill(zVertex);
1138   }
1139
1140   fTriggerMask=AliVEvent::kAny;
1141   Bool_t isSelectedBBBa = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerMask);
1142   if (!isSelectedBBBa) {
1143     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(10);
1144     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ10a")))->Fill(zVertex);
1145   }
1146
1147   if (titleVtx.Contains("Z")) {
1148     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(11);
1149     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ11a")))->Fill(zVertex);
1150   }
1151   else if (titleVtx.Contains("3D")) {
1152     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(12);
1153     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ12a")))->Fill(zVertex);
1154   } else {
1155     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(13);
1156     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ13a")))->Fill(zVertex);
1157   }
1158
1159   if (TMath::Abs(zVertex)<=fAnalCuts->GetMaxVtxZ()) {
1160     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(14);
1161     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ14a")))->Fill(zVertex);
1162   }
1163
1164   if ( fIsEventSelected ) {
1165     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(7);
1166     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ7a")))->Fill(zVertex);
1167   } else {
1168     ((TH1F*)(fOutput->FindObject("hEventsWithCandidates")))->Fill(15);
1169     if (fAdditionalChecks) ((TH1F*)(fOutput->FindObject("hZ15a")))->Fill(zVertex);
1170   }
1171
1172   return;
1173 }
1174
1175 //-----------------------------
1176 Int_t AliAnalysisTaskSELc2V0bachelor::MatchToMC(AliAODRecoCascadeHF *lc2bacV0,
1177                                                 Int_t *pdgDgLc2bacV0, Int_t *pdgDgV0,
1178                                                 TClonesArray *mcArray) {
1179   //
1180   // This is now implemented in AliAODRecoCascadeHF
1181   //
1182
1183   // bachelor
1184   AliAODTrack *bachelor = (AliAODTrack*)lc2bacV0->GetBachelor();
1185   if (!bachelor) return -1;
1186   Int_t labBachelor = TMath::Abs(bachelor->GetLabel());
1187   if (labBachelor<0) return -1;
1188   AliAODMCParticle *partBachelor = (AliAODMCParticle*)mcArray->At(labBachelor);
1189   if (!partBachelor) return -1;
1190   if (TMath::Abs(partBachelor->GetPdgCode())!=pdgDgLc2bacV0[0]) return -1;
1191
1192   Int_t labBacMother = partBachelor->GetMother();
1193   if (labBacMother<0) return -1;
1194   AliAODMCParticle *partBacMother = (AliAODMCParticle*)mcArray->At(labBacMother);
1195   if (!partBacMother) return -1;
1196   if (TMath::Abs(partBacMother->GetPdgCode())!=4122) return -1;
1197
1198   // V0
1199   AliAODTrack *posV0Daugh = (AliAODTrack*)lc2bacV0->Getv0PositiveTrack();
1200   AliAODTrack *negV0Daugh = (AliAODTrack*)lc2bacV0->Getv0NegativeTrack();
1201   if (!posV0Daugh || !negV0Daugh) return -1;
1202
1203   Int_t labV0pos = TMath::Abs(posV0Daugh->GetLabel());
1204   Int_t labV0neg = TMath::Abs(negV0Daugh->GetLabel());
1205   if (labV0pos<0 || labV0neg<0) return -1;
1206
1207   AliAODMCParticle *partV0pos = (AliAODMCParticle*)mcArray->At(labV0neg);
1208   AliAODMCParticle *partV0neg = (AliAODMCParticle*)mcArray->At(labV0pos);
1209   if (!partV0pos || !partV0neg) return -1;
1210
1211   if ( ! ( (TMath::Abs(partV0pos->GetPdgCode())==pdgDgV0[0] &&
1212             TMath::Abs(partV0neg->GetPdgCode())==pdgDgV0[1]) ||
1213            (TMath::Abs(partV0pos->GetPdgCode())==pdgDgV0[1] &&
1214             TMath::Abs(partV0neg->GetPdgCode())==pdgDgV0[0]) ) ) return -1;
1215   Int_t labV0posMother = partV0pos->GetMother();
1216   Int_t labV0negMother = partV0neg->GetMother();
1217
1218   if (labV0posMother<0 || labV0negMother<0) return -1;
1219   if (labV0posMother!=labV0negMother) return -1;
1220
1221   AliAODMCParticle *motherV0 = (AliAODMCParticle*)mcArray->At(labV0posMother);
1222   if (!motherV0) return-1;
1223
1224   if (TMath::Abs(motherV0->GetPdgCode())!=pdgDgLc2bacV0[1]) return -1;
1225   Int_t labV0mother = motherV0->GetMother();
1226   if (labV0mother<0) return -1;
1227   AliAODMCParticle *gMotherV0 = (AliAODMCParticle*)mcArray->At(labV0mother);
1228   if (!gMotherV0) return-1;
1229
1230   if ( !(pdgDgLc2bacV0[1]==310 && TMath::Abs(gMotherV0->GetPdgCode())==311) &&
1231        !(pdgDgLc2bacV0[1]==3122 && TMath::Abs(motherV0->GetPdgCode())==3122) ) return -1;
1232
1233   if ( (pdgDgLc2bacV0[1]==310 && TMath::Abs(gMotherV0->GetPdgCode())==311) ) {
1234     Int_t labV0GMother = gMotherV0->GetMother();
1235     if (labV0GMother<0) return -1;
1236     AliAODMCParticle *ggMotherV0 = (AliAODMCParticle*)mcArray->At(labV0GMother);
1237     if (!ggMotherV0) return-1;
1238
1239     if (TMath::Abs(ggMotherV0->GetPdgCode())!=4122) return -1;
1240     gMotherV0 = (AliAODMCParticle*)ggMotherV0;
1241     labV0mother=labV0GMother;
1242   }
1243   else if (pdgDgLc2bacV0[1]==3122 && TMath::Abs(motherV0->GetPdgCode())==3122) {
1244     if (TMath::Abs(gMotherV0->GetPdgCode())!=4122) return -1;
1245   }
1246
1247   if (labBacMother!=labV0mother) {
1248     return -1;
1249   }
1250
1251   return labBacMother;
1252
1253 }
1254
1255 //________________________________________________________________
1256 Int_t AliAnalysisTaskSELc2V0bachelor::SearchLcDaughter(TClonesArray *arrayMC, Int_t iii) {
1257   //
1258   // This is to check Lc dinasty
1259   //
1260
1261   Int_t indexToBeReturned=-999;
1262
1263   Int_t pdgLc=4122;
1264   Int_t pdgLambda=3122;
1265   Int_t pdgV0=310;
1266   Int_t pdgK0=311;
1267   Int_t pdgBachelor=2212;
1268   Int_t pdgBachelorPi=211;
1269
1270   TString fillthis="";
1271   fillthis="histMcStatLc";
1272
1273   AliAODMCParticle *searchLc = dynamic_cast<AliAODMCParticle*>(arrayMC->At(iii));
1274   if(!searchLc) return -999;
1275   if (TMath::Abs(searchLc->GetPdgCode()) != pdgLc) return -999;
1276
1277   ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(0);
1278   indexToBeReturned = 0;
1279
1280   ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*1);
1281   indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*1;
1282
1283   Int_t nDaughLc = searchLc->GetNDaughters();
1284   if (nDaughLc!=2) {
1285     ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*10);
1286     indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*10;
1287     return indexToBeReturned;
1288   }
1289
1290   Int_t index1=searchLc->GetDaughter(0);
1291   Int_t index2=searchLc->GetDaughter(1);
1292   if (index1<=0 || index2<=0) {
1293     return -999;
1294   }
1295
1296   AliAODMCParticle *daugh1 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index1));
1297   AliAODMCParticle *daugh2 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index2));
1298   if (!daugh1 || !daugh2) return -999;
1299
1300   Int_t daughPdg1 = TMath::Abs(daugh1->GetPdgCode());
1301   Int_t daughPdg2 = TMath::Abs(daugh2->GetPdgCode());
1302   if ( !( (daughPdg1==pdgBachelor && daughPdg2==pdgK0) ||
1303           (daughPdg2==pdgBachelor && daughPdg1==pdgK0) ||
1304           (daughPdg1==pdgLambda && daughPdg2==pdgBachelorPi) ||
1305           (daughPdg2==pdgLambda && daughPdg1==pdgBachelorPi) ) ) {
1306     ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*10);
1307     indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*10;
1308     return indexToBeReturned;
1309   }
1310
1311   if (daughPdg1==pdgK0 || daughPdg1==pdgLambda) {
1312     index1=searchLc->GetDaughter(1);
1313     index2=searchLc->GetDaughter(0);
1314   }
1315   daugh1 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index1));
1316   daugh2 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index2));
1317   if (!daugh1 || !daugh2) return -999;
1318
1319   daughPdg1=TMath::Abs(daugh1->GetPdgCode());
1320   daughPdg2=TMath::Abs(daugh2->GetPdgCode());
1321
1322   if ( daughPdg1==pdgBachelor && daughPdg2==pdgK0 ) { // Lc+ -> p K0bar
1323
1324     ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*2);
1325     indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*2;
1326
1327     Int_t nDaughK0 = daugh2->GetNDaughters();
1328     if (nDaughK0!=1) return -999;
1329
1330     Int_t indexK0daugh=daugh2->GetDaughter(0);
1331     if (indexK0daugh<=0) return -999;
1332
1333     AliAODMCParticle *daughK0 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(indexK0daugh));
1334     if (!daughK0) return -999;
1335
1336     Int_t daughK0Pdg=TMath::Abs(daughK0->GetPdgCode());
1337     if (daughK0Pdg!=pdgV0) {
1338       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*4); // K0L
1339       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*4;
1340       return indexToBeReturned;
1341     }
1342     ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*3); // K0S
1343     indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*3;
1344
1345     Int_t nDaughK0S = daughK0->GetNDaughters();
1346     if (nDaughK0S!=2) {
1347       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*5); // other decays for K0S
1348       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*5;
1349       return indexToBeReturned;
1350     }
1351
1352     index1=daughK0->GetDaughter(0);
1353     index2=daughK0->GetDaughter(1);
1354     if(index1<=0 || index2<=0) {
1355       return -999;
1356     }
1357
1358     AliAODMCParticle *daughK0S1 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index1));
1359     AliAODMCParticle *daughK0S2 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index2));
1360     if (!daughK0S1 || !daughK0S2) return -999;
1361
1362     Int_t daughK0S1pdg=TMath::Abs(daughK0S1->GetPdgCode());
1363     Int_t daughK0S2pdg=TMath::Abs(daughK0S2->GetPdgCode());
1364
1365     if ( daughK0S1pdg==211 && daughK0S2pdg==211 ) {
1366       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*6); // K0S -> pi+ pi-
1367       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*6;
1368     } else {
1369       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*5); // other decays for K0S
1370       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*5;
1371     }
1372
1373   } //if (daughPdg1==pdgBachelor && daughPdg2==pdgK0)
1374   else if ( daughPdg1==pdgBachelorPi && daughPdg2==pdgLambda ) { // Lc+ -> pi+ Lambda
1375
1376     ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*7);
1377     indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*7;
1378
1379     Int_t nDaughL = daugh2->GetNDaughters();
1380     if (nDaughL!=2) {
1381       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*8);
1382       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*8;
1383       return indexToBeReturned;
1384     }
1385
1386     index1=daugh2->GetDaughter(0);
1387     index2=daugh2->GetDaughter(1);
1388     if(index1<=0 || index2<=0) {
1389       return -999;
1390     }
1391
1392     AliAODMCParticle *daughL1 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index1));
1393     AliAODMCParticle *daughL2 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(index2));
1394     if (!daughL1 || !daughL2) return -999;
1395
1396     Int_t daughL1pdg=TMath::Abs(daughL1->GetPdgCode());
1397     Int_t daughL2pdg=TMath::Abs(daughL2->GetPdgCode());
1398     if ( (daughL1pdg==211 && daughL2pdg==2212) ||
1399          (daughL2pdg==211 && daughL1pdg==2212) ) {
1400       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*9);
1401       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*9;
1402     } else {
1403       ((TH1F*)(fOutput->FindObject(fillthis)))->Fill(TMath::Nint(searchLc->Charge()/3.)*8);
1404       indexToBeReturned = TMath::Nint(searchLc->Charge()/3.)*8;
1405     }
1406
1407   } //else if (daughPdg1==pdgBachelorPi && daughPdg2==pdgLambda)
1408
1409   return indexToBeReturned;
1410 }
1411
1412 //________________________________________________________________
1413 void AliAnalysisTaskSELc2V0bachelor::FillArmPodDistribution(AliAODv0 *vZero,
1414                                                             TString histoTitle,
1415                                                             TList *histoList) {
1416   //
1417   // This is to fill Armenteros Podolanski plots
1418   //
1419
1420   Double_t alpha = vZero->AlphaV0();
1421   Double_t qT    = vZero->PtArmV0();
1422
1423   ((TH2F*)(histoList->FindObject(histoTitle)))->Fill(alpha,qT);
1424
1425 }
1426
1427 //-------------------------------------------------------------------------------
1428 void AliAnalysisTaskSELc2V0bachelor::CheckCandidatesAtDifferentLevels(AliAODRecoCascadeHF *part, AliRDHFCutsLctoV0* cutsAnal) {
1429   //
1430   // This is to check candidates at different levels
1431   //
1432
1433   Bool_t areCutsUsingPID = cutsAnal->GetIsUsePID();
1434
1435   AliAODv0 * v0part = (AliAODv0*)part->Getv0();
1436   Bool_t onFlyV0 = v0part->GetOnFlyStatus(); // on-the-flight V0s
1437
1438   AliAODTrack *bachelor = (AliAODTrack*)part->GetBachelor();
1439
1440   if ( !onFlyV0 )
1441     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(3); // it counts number of candidates coming from offline V0s
1442
1443   if ( cutsAnal->IsInFiducialAcceptance(part->Pt(),part->Y(4122)) )
1444     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(4);
1445   if ( (((cutsAnal->IsSelected(part,AliRDHFCuts::kTracks))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)) )
1446     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(5);
1447   cutsAnal->SetUsePID(kFALSE);
1448   if ( (((cutsAnal->IsSelected(part,AliRDHFCuts::kCandidate))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)) )
1449     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(6);
1450   cutsAnal->SetUsePID(areCutsUsingPID);
1451   if ( (((cutsAnal->IsSelected(part,AliRDHFCuts::kPID))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)) )
1452     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(7);
1453   if ( (((cutsAnal->IsSelected(part,AliRDHFCuts::kCandidate))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)) )
1454     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(8);
1455   if ( (((cutsAnal->IsSelected(part,AliRDHFCuts::kAll))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)) )
1456     ((TH1F*)(fOutput->FindObject("hCandidateSelection")))->Fill(9);
1457
1458   if ( cutsAnal->IsInFiducialAcceptance(part->Pt(),part->Y(4122)) ) {
1459
1460     if ( ( (cutsAnal->IsSelected(part,AliRDHFCuts::kTracks))&(AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr) ) {
1461
1462       Int_t aaa = cutsAnal->IsSelected(part,AliRDHFCuts::kTracks);
1463       if ( (aaa&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr ) {
1464         if ( ( (aaa&AliRDHFCutsLctoV0::kLcToLpi)==AliRDHFCutsLctoV0::kLcToLpi && bachelor->Charge()==-1)  ||
1465              ( (aaa&AliRDHFCutsLctoV0::kLcToLBarpi)==AliRDHFCutsLctoV0::kLcToLBarpi && bachelor->Charge()==+1) )
1466           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates1")))->Fill( -aaa );
1467         else
1468           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates1")))->Fill( aaa );
1469       }
1470
1471       cutsAnal->SetUsePID(kFALSE);
1472       aaa = cutsAnal->IsSelected(part,AliRDHFCuts::kCandidate);
1473       if ((aaa&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr) {
1474         if ( ( (aaa&AliRDHFCutsLctoV0::kLcToLpi)==AliRDHFCutsLctoV0::kLcToLpi && bachelor->Charge()==-1) ||
1475              ( (aaa&AliRDHFCutsLctoV0::kLcToLBarpi)==AliRDHFCutsLctoV0::kLcToLBarpi && bachelor->Charge()==+1) )
1476           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates2")))->Fill( -aaa );
1477         else
1478           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates2")))->Fill( aaa );
1479       }
1480       cutsAnal->SetUsePID(areCutsUsingPID);
1481
1482       aaa = cutsAnal->IsSelected(part,AliRDHFCuts::kPID);
1483       if ((aaa&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr) {
1484         if ( ( (aaa&AliRDHFCutsLctoV0::kLcToLpi)==AliRDHFCutsLctoV0::kLcToLpi && bachelor->Charge()==-1) ||
1485              ( (aaa&AliRDHFCutsLctoV0::kLcToLBarpi)==AliRDHFCutsLctoV0::kLcToLBarpi && bachelor->Charge()==+1) )
1486           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates3")))->Fill( -aaa );
1487         else
1488           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates3")))->Fill( aaa );
1489       }
1490
1491       aaa = cutsAnal->IsSelected(part,AliRDHFCuts::kAll);
1492       if ((aaa&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr) {
1493         if ( ( (aaa&AliRDHFCutsLctoV0::kLcToLpi)==AliRDHFCutsLctoV0::kLcToLpi && bachelor->Charge()==-1) ||
1494              ( (aaa&AliRDHFCutsLctoV0::kLcToLBarpi)==AliRDHFCutsLctoV0::kLcToLBarpi && bachelor->Charge()==+1) )
1495           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates4")))->Fill( -aaa );
1496         else
1497           ((TH1F*)(fOutput->FindObject("hSwitchOnCandidates4")))->Fill( aaa );
1498       }
1499
1500     }
1501   }
1502
1503   return;
1504 }
1505
1506 //-------------------------------------------------------------------------------
1507 void AliAnalysisTaskSELc2V0bachelor::FillTheTree(AliAODRecoCascadeHF *part, AliRDHFCutsLctoV0 *cutsAnal, TClonesArray *mcArray, Int_t isLc) {
1508   //
1509   // This is to fill tree
1510   //
1511
1512   Double_t mk0sPDG = TDatabasePDG::Instance()->GetParticle(310)->Mass();
1513   Double_t mLPDG   = TDatabasePDG::Instance()->GetParticle(3122)->Mass();
1514
1515   Double_t invmassLc = part->InvMassLctoK0sP();
1516   Double_t invmassLc2Lpi = part->InvMassLctoLambdaPi();
1517
1518   AliAODTrack *bachelor = (AliAODTrack*)part->GetBachelor();
1519
1520   AliAODv0 * v0part = (AliAODv0*)part->Getv0();
1521   Double_t dcaV0ptp = v0part->GetDCA();
1522   Double_t invmassK0S = v0part->MassK0Short();
1523   Double_t invmassLambda = v0part->MassLambda();
1524   Double_t invmassLambdaBar = v0part->MassAntiLambda();
1525
1526   Int_t isLc2LBarpi=0, isLc2Lpi=0;
1527   Int_t mcLabel = -1;
1528   Int_t isDp2K0Spi=0, isDs2K0SK=0;
1529   Int_t mcLabel2 = -1;
1530   Int_t mcLabel3 = -1;
1531   Int_t isKstar12K0Spi=0, isKstar22K0Spi=0;
1532   Int_t mcLabel4 = -1;
1533   Int_t mcLabel5 = -1;
1534   Double_t ptCandByMC = 0.;//fmcPartCandidate->Pt();
1535   Double_t yCandByMC  = 0.;//fmcPartCandidate->Y() ;
1536   if (fUseMCInfo) {
1537     if (isLc) {
1538       Int_t pdgCand0 = 4122;
1539       Int_t pdgDgLctoV0bachelor0[2]={2212,310};
1540       Int_t pdgDgV0toDaughters0[2]={211,211};
1541       Int_t mcLabelLc2pK0S = part->MatchToMC(pdgCand0,pdgDgLctoV0bachelor0[1],pdgDgLctoV0bachelor0,pdgDgV0toDaughters0,mcArray,kTRUE);
1542       AliAODMCParticle *lambdaCpartMC = (AliAODMCParticle*)mcArray->At(mcLabelLc2pK0S);
1543       if (lambdaCpartMC) {
1544         ptCandByMC = lambdaCpartMC->Pt();
1545         yCandByMC  = lambdaCpartMC->Y() ;
1546       }
1547     }
1548
1549     Int_t pdgCand = 4122;
1550     Int_t pdgDgLctoV0bachelor[2]={211,3122};
1551     Int_t pdgDgV0toDaughters[2]={2212,211};
1552     mcLabel = part->MatchToMC(pdgCand,pdgDgLctoV0bachelor[1],pdgDgLctoV0bachelor,pdgDgV0toDaughters,mcArray,kTRUE);
1553     if (mcLabel!=-1) {
1554       if (bachelor->Charge()==-1) isLc2LBarpi=1;
1555       if (bachelor->Charge()==+1) isLc2Lpi=1;
1556       AliAODMCParticle *lambdaCpartMC = (AliAODMCParticle*)mcArray->At(mcLabel);
1557       if (lambdaCpartMC) {
1558         ptCandByMC = lambdaCpartMC->Pt();
1559         yCandByMC  = lambdaCpartMC->Y() ;
1560       }
1561     }
1562
1563     Int_t pdgCand2 = 411; // D+ -> pi+ K0S
1564     Int_t pdgCand3 = 431; // Ds+ -> K+ K0S
1565     Int_t pdgDgCand2[2]={211,310};
1566     Int_t pdgDgCand3[2]={321,310};
1567     pdgDgV0toDaughters[0]=211;
1568     pdgDgV0toDaughters[1]=211;
1569     mcLabel2 = part->MatchToMC(pdgCand2,pdgDgCand2[1],pdgDgCand2,pdgDgV0toDaughters,mcArray,kTRUE);
1570     mcLabel3 = part->MatchToMC(pdgCand3,pdgDgCand3[1],pdgDgCand3,pdgDgV0toDaughters,mcArray,kTRUE);
1571     if (mcLabel2!=-1) {
1572       isDp2K0Spi=1;
1573       AliAODMCParticle *lambdaCpartMC = (AliAODMCParticle*)mcArray->At(mcLabel2);
1574       if (lambdaCpartMC) {
1575         ptCandByMC = lambdaCpartMC->Pt();
1576         yCandByMC  = lambdaCpartMC->Y() ;
1577       }
1578     }
1579     if (mcLabel3!=-1) {
1580       isDs2K0SK=1;
1581       AliAODMCParticle *lambdaCpartMC = (AliAODMCParticle*)mcArray->At(mcLabel3);
1582       if (lambdaCpartMC) {
1583         ptCandByMC = lambdaCpartMC->Pt();
1584         yCandByMC  = lambdaCpartMC->Y() ;
1585       }
1586     }
1587
1588     Int_t pdgCand4 = 313; // K*(892)+ -> pi+ K0S
1589     Int_t pdgCand5 = 325; // K*(1430)+ -> pi+ K0S
1590     Int_t pdgDgCand4[2]={211,310};
1591     Int_t pdgDgCand5[2]={211,310};
1592     pdgDgV0toDaughters[0]=211;
1593     pdgDgV0toDaughters[1]=211;
1594     mcLabel4 = part->MatchToMC(pdgCand4,pdgDgCand4[1],pdgDgCand4,pdgDgV0toDaughters,mcArray,kTRUE);
1595     mcLabel5 = part->MatchToMC(pdgCand5,pdgDgCand5[1],pdgDgCand5,pdgDgV0toDaughters,mcArray,kTRUE);
1596     if (mcLabel4!=-1) {
1597       isKstar12K0Spi=1;
1598       AliAODMCParticle *lambdaCpartMC = (AliAODMCParticle*)mcArray->At(mcLabel4);
1599       if (lambdaCpartMC) {
1600         ptCandByMC = lambdaCpartMC->Pt();
1601         yCandByMC  = lambdaCpartMC->Y() ;
1602       }
1603     }
1604     if (mcLabel5!=-1) {
1605       isKstar22K0Spi=1;
1606       AliAODMCParticle *lambdaCpartMC = (AliAODMCParticle*)mcArray->At(mcLabel5);
1607       if (lambdaCpartMC) {
1608         ptCandByMC = lambdaCpartMC->Pt();
1609         yCandByMC  = lambdaCpartMC->Y() ;
1610       }
1611     }
1612   }
1613
1614   Int_t isLcByMC = isLc+isLc2LBarpi*2+isLc2Lpi*4+isDp2K0Spi*8+isDs2K0SK*16+isKstar12K0Spi*32+isKstar22K0Spi*64;
1615
1616   Bool_t isMCparticleInFiducialAcceptance = kTRUE;
1617   if (isLc || isLc2LBarpi || isLc2Lpi || isDp2K0Spi || isDs2K0SK || isKstar12K0Spi || isKstar22K0Spi) {
1618     isMCparticleInFiducialAcceptance = cutsAnal->IsInFiducialAcceptance(ptCandByMC,yCandByMC);
1619   }
1620
1621   Int_t isK0S = 0;
1622   Int_t isLambda = 0;
1623   Int_t isLambdaBar = 0;
1624   Int_t isGamma = 0;
1625   if (fUseMCInfo) {
1626     Int_t pdgDg2prong[2] = {211, 211};
1627     Int_t labelK0S = v0part->MatchToMC(310,mcArray,2,pdgDg2prong);
1628     if (labelK0S>=0) isK0S = 1;
1629
1630     pdgDg2prong[0] = 211;
1631     pdgDg2prong[1] = 2212;
1632     Int_t lambdaLabel = v0part->MatchToMC(3122,mcArray,2,pdgDg2prong);
1633     if (lambdaLabel>=0) {
1634       AliAODMCParticle *lambdaTrack = (AliAODMCParticle*)mcArray->At(lambdaLabel);
1635       if (lambdaTrack->GetPdgCode()==3122) isLambda = 1;
1636       else if (lambdaTrack->GetPdgCode()==-3122) isLambdaBar = 1;
1637     }
1638
1639     pdgDg2prong[0] = 11;
1640     pdgDg2prong[1] = 11;
1641     Int_t gammaLabel = v0part->MatchToMC(22,mcArray,2,pdgDg2prong);
1642     if (gammaLabel>=0) {
1643       AliAODMCParticle *gammaTrack = (AliAODMCParticle*)mcArray->At(gammaLabel);
1644       if (gammaTrack->GetPdgCode()==22) isGamma = 1;
1645     }
1646   }
1647
1648   Int_t isV0ByMC = isK0S+isLambdaBar*2+isLambda*4+isGamma*8;
1649
1650   Int_t isBachelorSelected = (bachelor->TestFilterMask(BIT(4)))*1 + (!(bachelor->TestFilterMask(BIT(4))))*2;
1651   isBachelorSelected += (bachelor->GetLabel()<0)*4 + (bachelor->GetLabel()>=0)*8;
1652   if ( ( !(bachelor->HasPointOnITSLayer(0)) && !(bachelor->HasPointOnITSLayer(1)) ) )
1653     isBachelorSelected += 16;
1654   else {
1655     if ( bachelor->HasPointOnITSLayer(0) && !(bachelor->HasPointOnITSLayer(1)) )
1656       isBachelorSelected += 32;
1657     else if ( !(bachelor->HasPointOnITSLayer(0)) && bachelor->HasPointOnITSLayer(1) )
1658       isBachelorSelected += 64;
1659     else
1660       isBachelorSelected += 128;
1661   }
1662
1663   AliAODTrack *v0pos = (AliAODTrack*)part->Getv0PositiveTrack();
1664   AliAODTrack *v0neg = (AliAODTrack*)part->Getv0NegativeTrack();
1665
1666   Int_t areV0daughtersSelected = (v0pos->TestFilterMask(BIT(4)))*1 + (!(v0pos->TestFilterMask(BIT(4))))*2;
1667   areV0daughtersSelected += (v0pos->GetLabel()<0)*4 + (v0pos->GetLabel()>=0)*8;
1668   areV0daughtersSelected += (v0neg->TestFilterMask(BIT(4)))*16 + (!(v0neg->TestFilterMask(BIT(4))))*32;
1669   areV0daughtersSelected += (v0neg->GetLabel()<0)*64 + (v0neg->GetLabel()>=0)*128;
1670
1671   Double_t nSigmaITSpr=-999.;
1672   cutsAnal->GetPidHF()->GetnSigmaITS(bachelor,4,nSigmaITSpr);
1673   Double_t nSigmaTPCpr=-999.;
1674   cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,4,nSigmaTPCpr);
1675   Double_t nSigmaTOFpr=-999.;
1676   cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,4,nSigmaTOFpr);
1677
1678   Double_t nSigmaITSpi=-999.;
1679   cutsAnal->GetPidHF()->GetnSigmaITS(bachelor,2,nSigmaITSpi);
1680   Double_t nSigmaTPCpi=-999.;
1681   cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,2,nSigmaTPCpi);
1682   Double_t nSigmaTOFpi=-999.;
1683   cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,2,nSigmaTOFpi);
1684
1685   Double_t nSigmaITSka=-999.;
1686   cutsAnal->GetPidHF()->GetnSigmaITS(bachelor,3,nSigmaITSka);
1687   Double_t nSigmaTPCka=-999.;
1688   cutsAnal->GetPidHF()->GetnSigmaTPC(bachelor,3,nSigmaTPCka);
1689   Double_t nSigmaTOFka=-999.;
1690   cutsAnal->GetPidHF()->GetnSigmaTOF(bachelor,3,nSigmaTOFka);
1691
1692
1693   Int_t flagToCheckCandidate = 1*(TMath::Abs(invmassK0S-mk0sPDG)<=0.050);
1694   flagToCheckCandidate+=2*((TMath::Abs(invmassLambdaBar-mLPDG)<=0.050) && (bachelor->Charge()==-1));
1695   flagToCheckCandidate+=4*((TMath::Abs(invmassLambda-mLPDG)<=0.050) && (bachelor->Charge()==+1));
1696   flagToCheckCandidate+=8*((TMath::Abs(invmassLambdaBar-mLPDG)<=0.050) && (bachelor->Charge()==+1));
1697   flagToCheckCandidate+=16*((TMath::Abs(invmassLambda-mLPDG)<=0.050) && (bachelor->Charge()==-1));
1698
1699   /*
1700   Bool_t areCutsUsingPID = cutsAnal->GetIsUsePID();
1701   cutsAnal->SetUsePID(kFALSE);
1702   Int_t aaa = cutsAnal->IsSelected(part,AliRDHFCuts::kCandidate);
1703   if ( (aaa&AliRDHFCutsLctoV0::kLcToK0Spr)==AliRDHFCutsLctoV0::kLcToK0Spr ) {
1704     if ( aaa==AliRDHFCutsLctoV0::kLcToK0Spr ) {
1705       flagToCheckCandidate = aaa; // Lc->K0S+p OK
1706     } else {
1707       if ( (aaa&AliRDHFCutsLctoV0::kLcToLpi)==AliRDHFCutsLctoV0::kLcToLpi ) {
1708         if (bachelor->Charge()==+1)
1709           flagToCheckCandidate = aaa; // Lc->Lambda+pi+
1710         else if (bachelor->Charge()==-1)
1711           flagToCheckCandidate =-aaa; // Lambda+pi- AS Lc->K0S+p candidate
1712       }
1713       if ( (aaa&AliRDHFCutsLctoV0::kLcToLBarpi)==AliRDHFCutsLctoV0::kLcToLBarpi ) {
1714         if (bachelor->Charge()==-1)
1715           flagToCheckCandidate = aaa; // Lc->LambdaBar+pi-
1716         else if (bachelor->Charge()==+1)
1717           flagToCheckCandidate =-aaa; // LambdaBar+pi+ AS Lc->K0S+p candidate
1718       }
1719     }
1720   } else {
1721     //if ( aaa==AliRDHFCutsLctoV0::kLcToK0Spr ) {
1722     //flagToCheckCandidate = -10-(AliRDHFCutsLctoV0::kLcToK0Spr); // NEVER
1723     //} else {
1724       if ( (aaa&AliRDHFCutsLctoV0::kLcToLpi)==AliRDHFCutsLctoV0::kLcToLpi ) {
1725         if (bachelor->Charge()==+1)
1726           flagToCheckCandidate = aaa; // Lc->Lambda+pi+ OK
1727         else if (bachelor->Charge()==-1)
1728           flagToCheckCandidate =-aaa; // Lambda+pi- AS Lc->Lambda+pi+ candidate
1729       }
1730       if ( (aaa&AliRDHFCutsLctoV0::kLcToLBarpi)==AliRDHFCutsLctoV0::kLcToLBarpi ) {
1731         if (bachelor->Charge()==-1)
1732           flagToCheckCandidate = aaa; // Lc->LambdaBar+pi- OK
1733         else if (bachelor->Charge()==+1)
1734           flagToCheckCandidate =-aaa; // LambdaBar+pi+ AS Lc->LambdaBar+pi- candidate
1735       }
1736       //}
1737   }
1738   cutsAnal->SetUsePID(areCutsUsingPID);
1739   */
1740
1741   fCandidateVariables[ 0] = fUseMCInfo+isLcByMC; // 0: real data; 1: bkg; 2: Lc->K0S+p; 3: Lc->LambdaBar+pbar; 5: Lc->Lambda+p
1742   fCandidateVariables[ 1] = fUseMCInfo+isV0ByMC; // 0: real data; 1: bkg; 2: K0S->pi+pi; 3: LambdaBar->pbar+pi+; 5: Lambda->p+pi-
1743   fCandidateVariables[ 2] = isBachelorSelected;
1744   fCandidateVariables[ 3] = areV0daughtersSelected;
1745   fCandidateVariables[ 4] = flagToCheckCandidate;
1746   fCandidateVariables[ 5] = invmassLc;
1747   fCandidateVariables[ 6] = invmassLc2Lpi;
1748   fCandidateVariables[ 7] = part->InvMass2Prongs(0,1,211,310); // D+ -> pi+ K0S
1749   fCandidateVariables[ 8] = part->InvMass2Prongs(0,1,321,310); // D+S -> K+ K0S
1750   fCandidateVariables[ 9] = invmassK0S;
1751   fCandidateVariables[10] = invmassLambda;
1752   fCandidateVariables[11] = invmassLambdaBar;
1753   fCandidateVariables[12] = v0part->InvMass2Prongs(0,1,11,11);
1754   fCandidateVariables[13] = part->GetDCA();
1755   fCandidateVariables[14] = dcaV0ptp;
1756   fCandidateVariables[15] = part->Getd0Prong(0);
1757   fCandidateVariables[16] = part->Getd0Prong(1);
1758   fCandidateVariables[17] = v0part->Getd0Prong(0);
1759   fCandidateVariables[18] = v0part->Getd0Prong(1);
1760   fCandidateVariables[19] = part->CosPointingAngle();
1761   fCandidateVariables[20] = part->CosV0PointingAngle();
1762   fCandidateVariables[21] = v0part->RadiusSecVtx();
1763   fCandidateVariables[22] = nSigmaITSpr;
1764   fCandidateVariables[23] = nSigmaITSpi;
1765   fCandidateVariables[24] = nSigmaITSka;
1766   fCandidateVariables[25] = nSigmaTPCpr;
1767   fCandidateVariables[26] = nSigmaTPCpi;
1768   fCandidateVariables[27] = nSigmaTPCka;
1769   fCandidateVariables[28] = nSigmaTOFpr;
1770   fCandidateVariables[29] = nSigmaTOFpi;
1771   fCandidateVariables[30] = nSigmaTOFka;
1772   fCandidateVariables[31] = part->Y(4122);
1773   fCandidateVariables[32] = bachelor->Eta();
1774   fCandidateVariables[33] = v0pos->Eta();
1775   fCandidateVariables[34] = v0neg->Eta();
1776   fCandidateVariables[35] = part->P();
1777   fCandidateVariables[36] = part->Pt();
1778   fCandidateVariables[37] = v0part->P();
1779   fCandidateVariables[38] = v0part->Pt();
1780   fCandidateVariables[39] = bachelor->P();
1781   fCandidateVariables[40] = bachelor->Pt();
1782   fCandidateVariables[41] = v0pos->P();
1783   fCandidateVariables[42] = v0pos->Pt();
1784   fCandidateVariables[43] = v0neg->P();
1785   fCandidateVariables[44] = v0neg->Pt();
1786   fCandidateVariables[45] = part->DecayLength();
1787   fCandidateVariables[46] = part->DecayLengthV0();
1788   fCandidateVariables[47] = part->CosPointingAngleXY();
1789   fCandidateVariables[48] = part->CosV0PointingAngleXY();
1790   fCandidateVariables[49] = part->DecayLengthXY();
1791   fCandidateVariables[50] = part->DecayLengthXYV0();
1792   fCandidateVariables[51] = part->NormalizedDecayLength();
1793   fCandidateVariables[52] = part->NormalizedV0DecayLength();
1794   fCandidateVariables[53] = part->NormalizedDecayLengthXY();
1795   fCandidateVariables[54] = part->NormalizedV0DecayLengthXY();
1796   Double_t xVtxLc=0, yVtxLc=0, zVtxLc=0;
1797   Double_t xLcMC=0,yLcMC=0,zLcMC=0;
1798   Double_t pxVtxBachelor=0, pyVtxBachelor=0, pzVtxBachelor=0;
1799   Double_t dcaForLc = PropagateToDCA(v0part,bachelor,fBzkG, xVtxLc, yVtxLc, zVtxLc, pxVtxBachelor, pyVtxBachelor, pzVtxBachelor);
1800   if (isLc) {
1801     Int_t pdgCand0 = 4122;
1802     Int_t pdgDgLctoV0bachelor0[2]={2212,310};
1803     Int_t pdgDgV0toDaughters0[2]={211,211};
1804     Int_t mcLabel0 = part->MatchToMC(pdgCand0,pdgDgLctoV0bachelor0[1],pdgDgLctoV0bachelor0,pdgDgV0toDaughters0,mcArray,kTRUE);
1805     AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel0));
1806     if(partLc){
1807       AliAODMCParticle *partLcDaug0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(partLc->GetDaughter(0)));
1808       if(partLcDaug0){
1809         xLcMC=partLcDaug0->Xv(), yLcMC=partLcDaug0->Yv(), zLcMC=partLcDaug0->Zv();
1810       }
1811     }
1812   } else if (isLc2LBarpi || isLc2Lpi) {
1813     AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel));
1814     AliAODMCParticle *partLcDaug0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(partLc->GetDaughter(0)));
1815     xLcMC=partLcDaug0->Xv(), yLcMC=partLcDaug0->Yv(), zLcMC=partLcDaug0->Zv();
1816   } else if (isDp2K0Spi) {
1817     AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel2));
1818     AliAODMCParticle *partLcDaug0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(partLc->GetDaughter(0)));
1819     xLcMC=partLcDaug0->Xv(), yLcMC=partLcDaug0->Yv(), zLcMC=partLcDaug0->Zv();
1820   } else if (isDs2K0SK) {
1821     AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel3));
1822     AliAODMCParticle *partLcDaug0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(partLc->GetDaughter(0)));
1823     xLcMC=partLcDaug0->Xv(), yLcMC=partLcDaug0->Yv(), zLcMC=partLcDaug0->Zv();
1824   } else if (isKstar12K0Spi) {
1825     AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel4));
1826     AliAODMCParticle *partLcDaug0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(partLc->GetDaughter(0)));
1827     xLcMC=partLcDaug0->Xv(), yLcMC=partLcDaug0->Yv(), zLcMC=partLcDaug0->Zv();
1828   } else if (isKstar22K0Spi) {
1829     AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(mcLabel5));
1830     AliAODMCParticle *partLcDaug0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(partLc->GetDaughter(0)));
1831     xLcMC=partLcDaug0->Xv(), yLcMC=partLcDaug0->Yv(), zLcMC=partLcDaug0->Zv();
1832   }
1833   fCandidateVariables[55]=dcaForLc;
1834   fCandidateVariables[56]=part->GetSecVtxX();
1835   fCandidateVariables[57]=part->GetSecVtxY();
1836   fCandidateVariables[58]=part->GetSecVtxZ();
1837   fCandidateVariables[59]=xVtxLc;
1838   fCandidateVariables[60]=yVtxLc;
1839   fCandidateVariables[61]=zVtxLc;
1840   fCandidateVariables[62]=xLcMC;
1841   fCandidateVariables[63]=yLcMC;
1842   fCandidateVariables[64]=zLcMC;
1843   fCandidateVariables[65]=bachelor->Px();
1844   fCandidateVariables[66]=bachelor->Py();
1845   fCandidateVariables[67]=pxVtxBachelor;
1846   fCandidateVariables[68]=pyVtxBachelor;
1847   fCandidateVariables[69]=v0part->Px();
1848   fCandidateVariables[70]=v0part->Py();
1849   fCandidateVariables[71]=fVtx1->GetX();
1850   fCandidateVariables[72]=fVtx1->GetY();
1851   fCandidateVariables[73]=fVtx1->GetZ();
1852   fCandidateVariables[74]=part->CosThetaStar(0,4122,2212,310);
1853   fCandidateVariables[75]=part->CosThetaStar(1,4122,2212,310);
1854   fCandidateVariables[76]=v0part->Eta();
1855   fCandidateVariables[77]=v0part->Y(310);
1856   fCandidateVariables[78]=pzVtxBachelor;
1857   fCandidateVariables[79]=v0part->Pz();
1858   fCandidateVariables[80]=bachelor->Charge();
1859   fCandidateVariables[81]=isMCparticleInFiducialAcceptance;
1860   if (fUseMCInfo) {
1861     fCandidateVariables[82]=0;
1862     if (bachelor->GetLabel()!=-1) {
1863       AliAODMCParticle *partBachelor = dynamic_cast<AliAODMCParticle*>(mcArray->At(TMath::Abs(bachelor->GetLabel())));
1864       if(partBachelor) fCandidateVariables[82]=partBachelor->GetPdgCode();
1865     }
1866   } else {
1867     fCandidateVariables[82]=-1;
1868   }
1869   fCandidateVariables[83] = part->InvMass2Prongs(0,1,211,310); // Kstar( 892)+ -> pi+K0S
1870   fCandidateVariables[84] = part->InvMass2Prongs(0,1,321,310); // Kstar(1430)+ -> pi+K0S
1871
1872   fCandidateVariables[85]=0;
1873   fCandidateVariables[86]=0;
1874   fCandidateVariables[87]=0;
1875   if (fUseMCInfo) {
1876     if (bachelor->GetLabel()!=-1 &&
1877         v0pos->GetLabel()!=-1 &&
1878         v0neg->GetLabel()!=-1) {
1879       const Int_t ndg=3;
1880       Int_t dgLabels[ndg]={TMath::Abs(bachelor->GetLabel()),
1881                            TMath::Abs(v0pos->GetLabel()),
1882                            TMath::Abs(v0neg->GetLabel())};
1883       Int_t ndgCk=0;
1884       Int_t *pdgDg=0;
1885       Int_t absLabelMother=-1;
1886       fCandidateVariables[85]=SearchForCommonMother(mcArray,
1887                                                     dgLabels,ndg,ndgCk,pdgDg,absLabelMother);
1888       AliAODMCParticle *part1 = dynamic_cast<AliAODMCParticle*>(mcArray->At(TMath::Abs(v0pos->GetLabel())));
1889       AliAODMCParticle *part2 = dynamic_cast<AliAODMCParticle*>(mcArray->At(TMath::Abs(v0neg->GetLabel())));
1890       if(!part1 || !part2) {
1891         AliDebug(2,"Daughter particles not found\n");
1892         return;
1893       }
1894       fCandidateVariables[86]=part1->GetPdgCode();
1895       fCandidateVariables[87]=part2->GetPdgCode();
1896     }
1897   }
1898
1899   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
1900   AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
1901   AliPIDResponse *pidResponse=inputHandler->GetPIDResponse();
1902   fCandidateVariables[88]=pidResponse->GetTOFResponse().GetStartTimeMask(bachelor->P());
1903
1904   //fCandidateVariables[65] = bachelor->Px();
1905   //fCandidateVariables[66] = bachelor->Py();
1906   //fCandidateVariables[67] = bachelor->Pz();
1907   //fCandidateVariables[68] = v0pos->Px();
1908   //fCandidateVariables[69] = v0pos->Py();
1909   //fCandidateVariables[70] = v0pos->Pz();
1910   //fCandidateVariables[71] = v0neg->Px();
1911   //fCandidateVariables[72] = v0neg->Py();
1912   //fCandidateVariables[73] = v0neg->Pz();
1913   //fCandidateVariables[74] = part->PxProng(0);
1914   //fCandidateVariables[75] = part->PyProng(0);
1915   //fCandidateVariables[76] = part->PzProng(0);
1916   //fCandidateVariables[77] = part->PxProng(1);
1917   //fCandidateVariables[78] = part->PyProng(1);
1918   //fCandidateVariables[79] = part->PzProng(1);
1919   //fCandidateVariables[80] = v0part->PxProng(0);
1920   //fCandidateVariables[81] = v0part->PyProng(0);
1921   //fCandidateVariables[82] = v0part->PzProng(0);
1922   //fCandidateVariables[83] = v0part->PxProng(1);
1923   //fCandidateVariables[84] = v0part->PyProng(1);
1924   //fCandidateVariables[85] = v0part->PzProng(1);
1925   //fCandidateVariables[86] = part->QtProng(0);
1926   //fCandidateVariables[87] = part->Alpha();
1927
1928   fVariablesTree->Fill();
1929
1930   return;
1931 }
1932
1933 //-------------------------------------------------------------------------------
1934 void AliAnalysisTaskSELc2V0bachelor::DefineTreeVariables() {
1935   //
1936   // This is to define tree variables
1937   //
1938
1939   const char* nameoutput = GetOutputSlot(4)->GetContainer()->GetName();
1940   fVariablesTree = new TTree(nameoutput,"Candidates variables tree");
1941   Int_t nVar = 89;
1942   fCandidateVariables = new Float_t [nVar];
1943   TString * fCandidateVariableNames = new TString[nVar];
1944   fCandidateVariableNames[ 0]="isLcByMC";
1945   fCandidateVariableNames[ 1]="isV0ByMC";
1946   fCandidateVariableNames[ 2]="flagToCheckBachelor";
1947   fCandidateVariableNames[ 3]="flagToCheckV0daughters";
1948   fCandidateVariableNames[ 4]="flagToCheckCandidate";
1949   fCandidateVariableNames[ 5]="massLc2K0Sp";
1950   fCandidateVariableNames[ 6]="massLc2Lambdapi";
1951   fCandidateVariableNames[ 7]="massD2K0Spi"; // D+ -> pi+ K0S
1952   fCandidateVariableNames[ 8]="massDS2K0SK"; // D+S -> K+ K0S
1953   fCandidateVariableNames[ 9]="massK0S";
1954   fCandidateVariableNames[10]="massLambda";
1955   fCandidateVariableNames[11]="massLambdaBar";
1956   fCandidateVariableNames[12]="massGamma";
1957   fCandidateVariableNames[13]="dcaLcptp"; // DCA Lc prong-to-prong
1958   fCandidateVariableNames[14]="dcaV0ptp";
1959   fCandidateVariableNames[15]="tImpParBach";
1960   fCandidateVariableNames[16]="tImpParV0";
1961   fCandidateVariableNames[17]="dcaV0postoPV";
1962   fCandidateVariableNames[18]="dcaV0negtoPV";
1963   fCandidateVariableNames[19]="cosPALc";
1964   fCandidateVariableNames[20]="cosPAK0S";
1965   fCandidateVariableNames[21]="rhoV0";
1966   fCandidateVariableNames[22]="nSigmaITSpr";
1967   fCandidateVariableNames[23]="nSigmaITSpi";
1968   fCandidateVariableNames[24]="nSigmaITSka";
1969   fCandidateVariableNames[25]="nSigmaTPCpr";
1970   fCandidateVariableNames[26]="nSigmaTPCpi";
1971   fCandidateVariableNames[27]="nSigmaTPCka";
1972   fCandidateVariableNames[28]="nSigmaTOFpr";
1973   fCandidateVariableNames[29]="nSigmaTOFpi";
1974   fCandidateVariableNames[30]="nSigmaTOFka";
1975   fCandidateVariableNames[31]="yLc";
1976   fCandidateVariableNames[32]="etaBach"; // etaBachelor
1977   fCandidateVariableNames[33]="etaV0pos"; // etaV0pos
1978   fCandidateVariableNames[34]="etaV0neg"; // etaV0neg
1979   fCandidateVariableNames[35]="LcP"; // @ DCA
1980   fCandidateVariableNames[36]="LcPt"; // @ DCA
1981   fCandidateVariableNames[37]="v0P"; // @ V0 DCA
1982   fCandidateVariableNames[38]="v0Pt"; // @ V0 DCA
1983   fCandidateVariableNames[39]="bachelorP"; // @ prim vtx
1984   fCandidateVariableNames[40]="bachelorPt"; // @ prim vtx
1985   fCandidateVariableNames[41]="V0positiveP"; // @ prim vtx
1986   fCandidateVariableNames[42]="V0positivePt"; // @ prim vtx
1987   fCandidateVariableNames[43]="V0negativeP"; // @ prim vtx
1988   fCandidateVariableNames[44]="V0negativePt"; // @ prim vtx
1989   fCandidateVariableNames[45]="decayLengthLc";
1990   fCandidateVariableNames[46]="decayLengthV0";
1991   fCandidateVariableNames[47]="cosPALcXY"; // cosPA XY x Lc
1992   fCandidateVariableNames[48]="cosPAV0XY"; // cosPA XY x V0
1993   fCandidateVariableNames[49]="decayLengthLcXY"; // decay length XY x Lc
1994   fCandidateVariableNames[50]="decayLengthV0XY"; // decay length XY x V0
1995   fCandidateVariableNames[51]="normalizedDecayLengthLc"; // normalized decay length x Lc
1996   fCandidateVariableNames[52]="normalizedDecayLengthV0"; // normalized decay length x V0
1997   fCandidateVariableNames[53]="normalizedDecayLengthXYLc"; // normalized decay length XY x Lc
1998   fCandidateVariableNames[54]="normalizedDecayLengthXYV0"; // normalized decay length XY x V0
1999   fCandidateVariableNames[55]="newLcDCA";
2000   fCandidateVariableNames[56]="xVtxLcBad";
2001   fCandidateVariableNames[57]="yVtxLcBad";
2002   fCandidateVariableNames[58]="zVtxLcBad";
2003   fCandidateVariableNames[59]="xVtxLcGood";
2004   fCandidateVariableNames[60]="yVtxLcGood";
2005   fCandidateVariableNames[61]="zVtxLcGood";
2006   fCandidateVariableNames[62]="xVtxLcMC";
2007   fCandidateVariableNames[63]="yVtxLcMC";
2008   fCandidateVariableNames[64]="zVtxLcMC";
2009   fCandidateVariableNames[65]="pxVtxBachelorBad";
2010   fCandidateVariableNames[66]="pyVtxBachelorBad";
2011   fCandidateVariableNames[67]="pxVtxBachelorGood";
2012   fCandidateVariableNames[68]="pyVtxBachelorGood";
2013   fCandidateVariableNames[69]="pxVtxV0";
2014   fCandidateVariableNames[70]="pyVtxV0";
2015   fCandidateVariableNames[71]="xPvtx";
2016   fCandidateVariableNames[72]="yPvtx";
2017   fCandidateVariableNames[73]="zPvtx";
2018   fCandidateVariableNames[74]="cosThetaStarBachelor";
2019   fCandidateVariableNames[75]="cosThetaStarV0";
2020   fCandidateVariableNames[76]="etaV0";
2021   fCandidateVariableNames[77]="yV0";
2022   fCandidateVariableNames[78]="pzVtxBachelorGood";
2023   fCandidateVariableNames[79]="pzVtxV0";
2024   fCandidateVariableNames[80]="bachelorCharge";
2025   fCandidateVariableNames[81]="isMCparticleInFiducialAcceptance";
2026   fCandidateVariableNames[82]="pdgBachelor"; // pdg MC bachelor
2027   fCandidateVariableNames[83]="massKstar12K0Spi"; // Kstar( 892)+ -> pi+ K0S
2028   fCandidateVariableNames[84]="massKstar22K0Spi"; // Kstar(1430)+ -> pi+ K0S
2029   fCandidateVariableNames[85]="pdgCandidate"; // pdg MC candidate recovered via new method
2030   fCandidateVariableNames[86]="pdgV0pos"; // pdg MC V0 positive
2031   fCandidateVariableNames[87]="pdgV0neg"; // pdg MC V0 negative
2032   fCandidateVariableNames[88]="startTimeMask"; // start time mask
2033
2034   //fCandidateVariableNames[65]="bachelorPx";
2035   //fCandidateVariableNames[66]="bachelorPy";
2036   //fCandidateVariableNames[67]="bachelorPz";
2037   //fCandidateVariableNames[68]="V0positivePx";
2038   //fCandidateVariableNames[69]="V0positivePy";
2039   //fCandidateVariableNames[70]="V0positivePz";
2040   //fCandidateVariableNames[71]="V0negativePx";
2041   //fCandidateVariableNames[72]="V0negativePy";
2042   //fCandidateVariableNames[73]="V0negativePz";
2043   //fCandidateVariableNames[74]="bachelorPxDCA";
2044   //fCandidateVariableNames[75]="bachelorPyDCA";
2045   //fCandidateVariableNames[76]="bachelorPzDCA";
2046   //fCandidateVariableNames[77]="v0PxDCA";
2047   //fCandidateVariableNames[78]="v0PyDCA";
2048   //fCandidateVariableNames[79]="v0PzDCA";
2049   //fCandidateVariableNames[80]="V0positivePxDCA";
2050   //fCandidateVariableNames[81]="V0positivePyDCA";
2051   //fCandidateVariableNames[82]="V0positivePzDCA";
2052   //fCandidateVariableNames[83]="V0negativePxDCA";
2053   //fCandidateVariableNames[84]="V0negativePyDCA";
2054   //fCandidateVariableNames[85]="V0negativePzDCA";
2055   //fCandidateVariableNames[86]="qtLc";
2056   //fCandidateVariableNames[87]="alphaLc";
2057
2058   for (Int_t ivar=0; ivar<nVar; ivar++) {
2059     fVariablesTree->Branch(fCandidateVariableNames[ivar].Data(),&fCandidateVariables[ivar],Form("%s/f",fCandidateVariableNames[ivar].Data()));
2060   }
2061
2062   return;
2063 }
2064
2065 //__________________________________________________________________________
2066 void  AliAnalysisTaskSELc2V0bachelor::DefineGeneralHistograms() {
2067   //
2068   // This is to define general histograms
2069   //
2070
2071   fCEvents = new TH1F("fCEvents","conter",18,0,18);
2072   fCEvents->SetStats(kTRUE);
2073   fCEvents->GetXaxis()->SetBinLabel(1,"X1");
2074   fCEvents->GetXaxis()->SetBinLabel(2,"Analyzed events");
2075   fCEvents->GetXaxis()->SetBinLabel(3,"AliAODVertex exists");
2076   fCEvents->GetXaxis()->SetBinLabel(4,"CascadesHF exists");
2077   fCEvents->GetXaxis()->SetBinLabel(5,"MCarray exists");
2078   fCEvents->GetXaxis()->SetBinLabel(6,"MCheader exists");
2079   fCEvents->GetXaxis()->SetBinLabel(7,"GetNContributors()>0");
2080   fCEvents->GetXaxis()->SetBinLabel(8,"IsEventSelected");
2081   fCEvents->GetXaxis()->SetBinLabel(9,"triggerClass!=CINT1");
2082   fCEvents->GetXaxis()->SetBinLabel(10,"triggerMask!=kAnyINT");
2083   fCEvents->GetXaxis()->SetBinLabel(11,"triggerMask!=kAny");
2084   fCEvents->GetXaxis()->SetBinLabel(12,"vtxTitle.Contains(Z)");
2085   fCEvents->GetXaxis()->SetBinLabel(13,"vtxTitle.Contains(3D)");
2086   fCEvents->GetXaxis()->SetBinLabel(14,"vtxTitle.Doesn'tContain(Z-3D)");
2087   fCEvents->GetXaxis()->SetBinLabel(15,Form("zVtx<=%2.0fcm",fAnalCuts->GetMaxVtxZ()));
2088   fCEvents->GetXaxis()->SetBinLabel(16,"!IsEventSelected");
2089   fCEvents->GetXaxis()->SetBinLabel(17,"triggerMask!=kAnyINT || triggerClass!=CINT1");
2090   fCEvents->GetXaxis()->SetBinLabel(18,Form("zVtxMC<=%2.0fcm",fAnalCuts->GetMaxVtxZ()));
2091   //fCEvents->GetXaxis()->SetTitle("");
2092   fCEvents->GetYaxis()->SetTitle("counts");
2093
2094   fOutput->Add(fCEvents);
2095   TString fillthis="";
2096
2097   if (fUseMCInfo && fAdditionalChecks) {
2098     fillthis="histMcStatLc";
2099     TH1F* mcStatisticLc = new TH1F(fillthis.Data(),"#Lambda_{c} generated and their decays",21,-10.5,10.5);
2100     fOutput->Add(mcStatisticLc);
2101   }
2102
2103   //fillthis="histopionV0SigmaVspTOF";
2104   //TH2F *hpionV0SigmaVspTOF=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2105   fillthis="histoprotonBachSigmaVspTOF";
2106   TH2F *hprotonBachSigmaVspTOF=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2107
2108   //fOutput->Add(hpionV0SigmaVspTOF);
2109   fOutput->Add(hprotonBachSigmaVspTOF);
2110
2111   //fillthis="histopionV0SigmaVspTPC";
2112   //TH2F *hpionV0SigmaVspTPC=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2113   fillthis="histoprotonBachSigmaVspTPC";
2114   TH2F *hprotonBachSigmaVspTPC=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2115
2116   //fOutput->Add(hpionV0SigmaVspTPC);
2117   fOutput->Add(hprotonBachSigmaVspTPC);
2118
2119   if (fUseMCInfo) {
2120
2121     //fillthis="histopionV0SigmaVspTOFsgn";
2122     //TH2F *hpionV0SigmaVspTOFsgn=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2123     fillthis="histoprotonBachSigmaVspTOFsgn";
2124     TH2F *hprotonBachSigmaVspTOFsgn=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2125
2126     //fOutput->Add(hpionV0SigmaVspTOFsgn);
2127     fOutput->Add(hprotonBachSigmaVspTOFsgn);
2128
2129     //fillthis="histopionV0SigmaVspTPCsgn";
2130     //TH2F *hpionV0SigmaVspTPCsgn=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2131     fillthis="histoprotonBachSigmaVspTPCsgn";
2132     TH2F *hprotonBachSigmaVspTPCsgn=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2133
2134     //fOutput->Add(hpionV0SigmaVspTPCsgn);
2135     fOutput->Add(hprotonBachSigmaVspTPCsgn);
2136
2137
2138     //fillthis="histopionV0SigmaVspTOFbkg";
2139     //TH2F *hpionV0SigmaVspTOFbkg=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2140     fillthis="histoprotonBachSigmaVspTOFbkg";
2141     TH2F *hprotonBachSigmaVspTOFbkg=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2142
2143     //fOutput->Add(hpionV0SigmaVspTOFbkg);
2144     fOutput->Add(hprotonBachSigmaVspTOFbkg);
2145
2146     //fillthis="histopionV0SigmaVspTPCbkg";
2147     //TH2F *hpionV0SigmaVspTPCbkg=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2148     fillthis="histoprotonBachSigmaVspTPCbkg";
2149     TH2F *hprotonBachSigmaVspTPCbkg=new TH2F(fillthis.Data(),fillthis.Data(),300,0.,30.,100,-5.,5.);
2150
2151     //fOutput->Add(hpionV0SigmaVspTPCbkg);
2152     fOutput->Add(hprotonBachSigmaVspTPCbkg);
2153
2154   }
2155
2156   if (fAdditionalChecks) {
2157
2158     TH1F *hZ2 = new TH1F("hZ2","",100,-50.,50.);
2159     fOutput->Add(hZ2);
2160     TH1F *hZ3 = new TH1F("hZ3","",100,-50.,50.);
2161     fOutput->Add(hZ3);
2162     TH1F *hZ4 = new TH1F("hZ4","",100,-50.,50.);
2163     fOutput->Add(hZ4);
2164     TH1F *hZ5 = new TH1F("hZ5","",100,-50.,50.);
2165     fOutput->Add(hZ5);
2166     TH1F *hZ6 = new TH1F("hZ6","",100,-50.,50.);
2167     fOutput->Add(hZ6);
2168     TH1F *hZ7 = new TH1F("hZ7","",100,-50.,50.);
2169     fOutput->Add(hZ7);
2170     TH1F *hZ8 = new TH1F("hZ8","",100,-50.,50.);
2171     fOutput->Add(hZ8);
2172     TH1F *hZ9 = new TH1F("hZ9","",100,-50.,50.);
2173     fOutput->Add(hZ9);
2174     TH1F *hZ10 = new TH1F("hZ10","",100,-50.,50.);
2175     fOutput->Add(hZ10);
2176     TH1F *hZ11 = new TH1F("hZ11","",100,-50.,50.);
2177     fOutput->Add(hZ11);
2178     TH1F *hZ12 = new TH1F("hZ12","",100,-50.,50.);
2179     fOutput->Add(hZ12);
2180     TH1F *hZ13 = new TH1F("hZ13","",100,-50.,50.);
2181     fOutput->Add(hZ13);
2182     TH1F *hZ14 = new TH1F("hZ14","",100,-50.,50.);
2183     fOutput->Add(hZ14);
2184     TH1F *hZ15 = new TH1F("hZ15","",100,-50.,50.);
2185     fOutput->Add(hZ15);
2186     TH1F *hZ16 = new TH1F("hZ16","",100,-50.,50.);
2187     fOutput->Add(hZ16);
2188   }
2189
2190   TH1F *hCandidateSelection = new TH1F("hCandidateSelection","",10,-0.5,9.5);
2191   hCandidateSelection->GetXaxis()->SetBinLabel(1,"IsEventSelected");
2192   hCandidateSelection->GetXaxis()->SetBinLabel(2,"IsSecondaryVtx");
2193   hCandidateSelection->GetXaxis()->SetBinLabel(3,"V0toPosNeg");
2194   hCandidateSelection->GetXaxis()->SetBinLabel(4,"offlineV0");
2195   hCandidateSelection->GetXaxis()->SetBinLabel(5,"isInFiducialAcceptance");
2196   hCandidateSelection->GetXaxis()->SetBinLabel(6,"analCuts::kTracks");
2197   hCandidateSelection->GetXaxis()->SetBinLabel(7,"analCuts::kCandidateNoPID");
2198   hCandidateSelection->GetXaxis()->SetBinLabel(8,"analCuts::kPID");
2199   hCandidateSelection->GetXaxis()->SetBinLabel(9,"analCuts::kCandidateWithPID");
2200   hCandidateSelection->GetXaxis()->SetBinLabel(10,"analCuts::kAll");
2201   fOutput->Add(hCandidateSelection);
2202
2203   TH1F *hEventsWithCandidates = new TH1F("hEventsWithCandidates","conter",11,5.5,16.5);
2204   hEventsWithCandidates->GetXaxis()->SetBinLabel(1,"GetNContributors()>0");
2205   hEventsWithCandidates->GetXaxis()->SetBinLabel(2,"IsEventSelected");
2206   hEventsWithCandidates->GetXaxis()->SetBinLabel(3,"triggerClass!=CINT1");
2207   hEventsWithCandidates->GetXaxis()->SetBinLabel(4,"triggerMask!=kAnyINT");
2208   hEventsWithCandidates->GetXaxis()->SetBinLabel(5,"triggerMask!=kAny");
2209   hEventsWithCandidates->GetXaxis()->SetBinLabel(6,"vtxTitle.Contains(Z)");
2210   hEventsWithCandidates->GetXaxis()->SetBinLabel(7,"vtxTitle.Contains(3D)");
2211   hEventsWithCandidates->GetXaxis()->SetBinLabel(8,"vtxTitle.Doesn'tContain(Z-3D)");
2212   hEventsWithCandidates->GetXaxis()->SetBinLabel(9,Form("zVtx<=%2.0fcm",fAnalCuts->GetMaxVtxZ()));
2213   hEventsWithCandidates->GetXaxis()->SetBinLabel(10,"!IsEventSelected");
2214   hEventsWithCandidates->GetXaxis()->SetBinLabel(11,"triggerMask!=kAnyINT || triggerClass!=CINT1");
2215   fOutput->Add(hEventsWithCandidates);
2216
2217   if (fAdditionalChecks) {
2218
2219     TH1F *hZ6a = new TH1F("hZ6a","",100,-50.,50.);
2220     fOutput->Add(hZ6a);
2221     TH1F *hZ7a = new TH1F("hZ7a","",100,-50.,50.);
2222     fOutput->Add(hZ7a);
2223     TH1F *hZ8a = new TH1F("hZ8a","",100,-50.,50.);
2224     fOutput->Add(hZ8a);
2225     TH1F *hZ9a = new TH1F("hZ9a","",100,-50.,50.);
2226     fOutput->Add(hZ9a);
2227     TH1F *hZ10a = new TH1F("hZ10a","",100,-50.,50.);
2228     fOutput->Add(hZ10a);
2229     TH1F *hZ11a = new TH1F("hZ11a","",100,-50.,50.);
2230     fOutput->Add(hZ11a);
2231     TH1F *hZ12a = new TH1F("hZ12a","",100,-50.,50.);
2232     fOutput->Add(hZ12a);
2233     TH1F *hZ13a = new TH1F("hZ13a","",100,-50.,50.);
2234     fOutput->Add(hZ13a);
2235     TH1F *hZ14a = new TH1F("hZ14a","",100,-50.,50.);
2236     fOutput->Add(hZ14a);
2237     TH1F *hZ15a = new TH1F("hZ15a","",100,-50.,50.);
2238     fOutput->Add(hZ15a);
2239     TH1F *hZ16a = new TH1F("hZ16a","",100,-50.,50.);
2240     fOutput->Add(hZ16a);
2241   }
2242
2243   TH1F *hSwitchOnCandidates1 = new TH1F("hSwitchOnCandidates1","",15,-7.5,7.5);
2244   fOutput->Add(hSwitchOnCandidates1);
2245   TH1F *hSwitchOnCandidates2 = new TH1F("hSwitchOnCandidates2","",15,-7.5,7.5);
2246   fOutput->Add(hSwitchOnCandidates2);
2247   TH1F *hSwitchOnCandidates3 = new TH1F("hSwitchOnCandidates3","",15,-7.5,7.5);
2248   fOutput->Add(hSwitchOnCandidates3);
2249   TH1F *hSwitchOnCandidates4 = new TH1F("hSwitchOnCandidates4","",15,-7.5,7.5);
2250   fOutput->Add(hSwitchOnCandidates4);
2251
2252   return;
2253 }
2254
2255 //________________________________________________________________________
2256 void  AliAnalysisTaskSELc2V0bachelor::DefineAnalysisHistograms() {
2257   //
2258   // This is to define analysis histograms
2259   //
2260
2261   if (fIsK0SAnalysis) DefineK0SHistos();// hK0S histos declarations
2262
2263   return;
2264 }
2265
2266 //________________________________________________________________________
2267 void  AliAnalysisTaskSELc2V0bachelor::FillAnalysisHistograms(AliAODRecoCascadeHF *part, Bool_t isBachelorID, TString appendthis) {
2268   //
2269   // This is to fill analysis histograms
2270   //
2271
2272   TString fillthis="";
2273
2274   Double_t invmassLc = part->InvMassLctoK0sP();
2275   Double_t lambdacpt = part->Pt();
2276
2277   AliAODTrack *bachelor = (AliAODTrack*)part->GetBachelor();
2278   Double_t momBach  = bachelor->P();
2279
2280   AliAODv0 * v0part = (AliAODv0*)part->Getv0();
2281   Double_t momK0S = v0part->P();
2282   Double_t ptK0S = v0part->Pt();
2283   Double_t dcaV0ptp = v0part->GetDCA();
2284   Double_t invmassK0S = v0part->MassK0Short();
2285
2286     fillthis="histK0SMass"+appendthis;
2287     //    cout << fillthis << endl;
2288     ((TH2F*)(fOutputAll->FindObject(fillthis)))->Fill(invmassK0S,ptK0S);
2289     if (isBachelorID)  ((TH2F*)(fOutputPIDBach->FindObject(fillthis)))->Fill(invmassK0S,ptK0S);
2290
2291     fillthis="histpK0Svsp"+appendthis;
2292     //    cout << fillthis << endl;
2293     ((TH2F*)(fOutputAll->FindObject(fillthis)))->Fill(momBach,momK0S);
2294     if (isBachelorID)  ((TH2F*)(fOutputPIDBach->FindObject(fillthis)))->Fill(momBach,momK0S);
2295
2296     fillthis="histDCAtoPVvspK0S"+appendthis;
2297     //    cout << fillthis << endl;
2298     ((TH2F*)(fOutputAll->FindObject(fillthis)))->Fill(momK0S,dcaV0ptp);
2299     if (isBachelorID)  ((TH2F*)(fOutputPIDBach->FindObject(fillthis)))->Fill(momK0S,dcaV0ptp);
2300
2301     fillthis="histArmPodK0S"+appendthis;
2302     //    cout << fillthis << endl;
2303     FillArmPodDistribution(v0part,fillthis,fOutputAll);
2304     if (isBachelorID) FillArmPodDistribution(v0part,fillthis,fOutputPIDBach);
2305
2306     fillthis="histLcMassByK0S"+appendthis;
2307     //    cout << fillthis << endl;
2308     ((TH2F*)(fOutputAll->FindObject(fillthis)))->Fill(invmassLc,lambdacpt);
2309     if (isBachelorID)((TH2F*)(fOutputPIDBach->FindObject(fillthis)))->Fill(invmassLc,lambdacpt);
2310
2311     return;
2312 }
2313 //---------------------------
2314 Double_t AliAnalysisTaskSELc2V0bachelor::PropagateToDCA(AliAODv0 *v, AliAODTrack *bachelor, Double_t b,
2315                                                         Double_t &xVtxLc, Double_t &yVtxLc, Double_t &zVtxLc,
2316                                                         Double_t &pxVtxBachelor, Double_t &pyVtxBachelor, Double_t &pzVtxBachelor) {
2317   //--------------------------------------------------------------------
2318   // This function returns the DCA between the V0 and the track
2319   // This is a copy of AliCascadeVertexer::PropagateToDCA(...) method
2320   //--------------------------------------------------------------------
2321
2322   // Get AliExternalTrackParam out of the AliAODTracks
2323   Double_t xyz[3], pxpypz[3], cv[21]; Short_t sign;
2324   bachelor->PxPyPz(pxpypz);
2325   bachelor->XvYvZv(xyz);
2326   bachelor->GetCovarianceXYZPxPyPz(cv);
2327   sign=bachelor->Charge();
2328   AliExternalTrackParam *t = new AliExternalTrackParam(xyz,pxpypz,cv,sign);
2329
2330   Double_t alpha=t->GetAlpha(), cs1=TMath::Cos(alpha), sn1=TMath::Sin(alpha);
2331   //Double_t alpha = GetAlpha(xyz,pxpypz), cs1=TMath::Cos(alpha), sn1=TMath::Sin(alpha);
2332
2333   // position and momentum of bachelor
2334   Double_t x1=xyz[0], y1=xyz[1], z1=xyz[2];
2335   Double_t px1=pxpypz[0], py1=pxpypz[1], pz1=pxpypz[2];
2336
2337   // position and momentum of V0
2338   Double_t x2=v->DecayVertexV0X(),
2339     y2=v->DecayVertexV0Y(),
2340     z2=v->DecayVertexV0Z();
2341   Double_t px2=v->Px(),
2342     py2=v->Py(),
2343     pz2=v->Pz();
2344
2345   /*
2346   AliAODTrack *trackP = (AliAODTrack*) v->GetDaughter(0);
2347   //Double_t pxpypzP[3]; trackP->PxPyPz(pxpypzP);
2348   //Double_t xyzP[3]; trackP->XvYvZv(xyzP);
2349   Double_t cvP[21]; trackP->GetCovarianceXYZPxPyPz(cvP);
2350   //Short_t signP=trackP->Charge();
2351   //AliExternalTrackParam *tP = new AliExternalTrackParam(xyzP,pxpypzP,cvP,signP);
2352
2353   // Get AliExternalTrackParam out of the AliAODTrack
2354   AliAODTrack *trackN = (AliAODTrack*) v->GetDaughter(1);
2355   //Double_t pxpypzN[3]; trackN->PxPyPz(pxpypzN);
2356   //Double_t xyzN[3]; trackN->XvYvZv(xyzN);
2357   Double_t cvN[21]; trackN->GetCovarianceXYZPxPyPz(cvN);
2358   //Short_t signN=trackN->Charge();
2359   //AliExternalTrackParam *tN = new AliExternalTrackParam(xyzN,pxpypzN,cvN,signN);
2360
2361   Double_t xyzV0[3]={x2,y2,z2};
2362   Double_t pxpypzV0[3]={px2,py2,pz2};
2363   Double_t cvV0[21]; for (Int_t ii=0; ii<21; ii++) cvV0[ii]=cvP[ii]+cvN[ii];
2364   AliNeutralTrackParam *trackV0 = new AliNeutralTrackParam(xyzV0,pxpypzV0,cvV0,0);
2365   */
2366
2367   // calculation dca
2368   Double_t dd= Det(x2-x1,y2-y1,z2-z1,px1,py1,pz1,px2,py2,pz2);
2369   Double_t ax= Det(py1,pz1,py2,pz2);
2370   Double_t ay=-Det(px1,pz1,px2,pz2);
2371   Double_t az= Det(px1,py1,px2,py2);
2372
2373   Double_t dca=TMath::Abs(dd)/TMath::Sqrt(ax*ax + ay*ay + az*az);
2374
2375   // bachelor point @ the DCA
2376   Double_t t1 = Det(x2-x1,y2-y1,z2-z1,px2,py2,pz2,ax,ay,az)/
2377     Det(px1,py1,pz1,px2,py2,pz2,ax,ay,az);
2378   x1 += px1*t1; y1 += py1*t1; z1 += pz1*t1;
2379
2380   //propagate track to the point of DCA
2381   Double_t rho1=x1*cs1 + y1*sn1;
2382   if (!t->PropagateTo(rho1,b)) {
2383     Error("PropagateToDCA","Propagation failed !");
2384     delete t; t=NULL;
2385     return 1.e+33;
2386   }
2387
2388   Double_t pBachelorDCA[3]; t->GetPxPyPz(pBachelorDCA);
2389   pxVtxBachelor=pBachelorDCA[0], pyVtxBachelor=pBachelorDCA[1], pzVtxBachelor=pBachelorDCA[2];
2390
2391   delete t; t=NULL;
2392
2393   // V0 point @ the DCA
2394   Double_t t2 = Det(x1-x2,y1-y2,z1-z2,px1,py1,pz1,ax,ay,az)/
2395     Det(px2,py2,pz2,px1,py1,pz1,ax,ay,az);
2396   x2 += px2*t2; y2 += py2*t2; z2 += pz2*t2;
2397
2398
2399   // Lc decay vtx
2400   xVtxLc = 0.5*(x1+x2);
2401   yVtxLc = 0.5*(y1+y2);
2402   zVtxLc = 0.5*(z1+z2);
2403
2404   return dca;
2405
2406 }
2407
2408 //---------------------------
2409 Double_t AliAnalysisTaskSELc2V0bachelor::GetAlpha(Double_t xyz[3],Double_t pxpypz[3])
2410 {
2411   //
2412   // To estimate alpha according to what done in the AliExternalTrackParam::Set(...) method
2413   //
2414
2415   Double_t alpha = 0.;
2416
2417   const double kSafe = 1e-5;
2418   Double_t radPos2 = xyz[0]*xyz[0]+xyz[1]*xyz[1];
2419   Double_t radMax  = 45.; // approximately ITS outer radius
2420   if (radPos2 < radMax*radMax) { // inside the ITS
2421     alpha = TMath::ATan2(pxpypz[1],pxpypz[0]);
2422   } else { // outside the ITS
2423     Float_t phiPos = TMath::Pi()+TMath::ATan2(-xyz[1], -xyz[0]);
2424      alpha =
2425        TMath::DegToRad()*(20*((((Int_t)(phiPos*TMath::RadToDeg()))/20))+10);
2426   }
2427
2428   Double_t cs=TMath::Cos(alpha), sn=TMath::Sin(alpha);
2429   // protection: avoid alpha being too close to 0 or +-pi/2
2430   if (TMath::Abs(sn)<2*kSafe) {
2431     if (alpha>0) alpha += alpha< TMath::Pi()/2. ?  2*kSafe : -2*kSafe;
2432     else         alpha += alpha>-TMath::Pi()/2. ? -2*kSafe :  2*kSafe;
2433     cs=TMath::Cos(alpha);
2434     sn=TMath::Sin(alpha);
2435   }
2436   else if (TMath::Abs(cs)<2*kSafe) {
2437     if (alpha>0) alpha += alpha> TMath::Pi()/2. ? 2*kSafe : -2*kSafe;
2438     else         alpha += alpha>-TMath::Pi()/2. ? 2*kSafe : -2*kSafe;
2439     cs=TMath::Cos(alpha);
2440     sn=TMath::Sin(alpha);
2441   }
2442
2443
2444   return alpha;
2445 }
2446
2447 //---------------------------
2448 Double_t AliAnalysisTaskSELc2V0bachelor::Det(Double_t a00, Double_t a01,
2449                                              Double_t a10, Double_t a11) const {
2450   //--------------------------------------------------------------------
2451   // This function calculates locally a 2x2 determinant.
2452   // This is a copy of the AliCascadeVertexer::Det(...) method
2453   //--------------------------------------------------------------------
2454   return a00*a11 - a01*a10;
2455 }
2456
2457 //---------------------------
2458 Double_t AliAnalysisTaskSELc2V0bachelor::Det(Double_t a00,Double_t a01,Double_t a02,
2459                                              Double_t a10,Double_t a11,Double_t a12,
2460                                              Double_t a20,Double_t a21,Double_t a22) const {
2461   //--------------------------------------------------------------------
2462   // This function calculates locally a 3x3 determinant
2463   // This is a copy of the AliCascadeVertexer::Det(...) method
2464   //--------------------------------------------------------------------
2465   return  a00*Det(a11,a12,a21,a22)-a01*Det(a10,a12,a20,a22)+a02*Det(a10,a11,a20,a21);
2466 }
2467
2468 //----------------------------------------------------------------------------
2469 Int_t AliAnalysisTaskSELc2V0bachelor::MatchToMClabelC(AliAODRecoCascadeHF *candidate,
2470                                                       TClonesArray *mcArray)
2471 {
2472   //
2473   // Check if this candidate is matched to a MC signal  Lc -> p K0S + X
2474   // If no, return -1
2475   // If yes, return label (>=0) of the AliAODMCParticle
2476   //
2477
2478   AliAODv0 *theV0 = dynamic_cast<AliAODv0*>(candidate->Getv0()); // the V0
2479   AliVTrack *trk = dynamic_cast<AliVTrack*>(candidate->GetBachelor()); // the bachelor
2480   if (!trk || !theV0) return -1;
2481
2482   if (trk->GetLabel()==-1) return -1;
2483   Int_t bachLabels = TMath::Abs(trk->GetLabel());
2484   AliAODMCParticle*bachelorMC = dynamic_cast<AliAODMCParticle*>(mcArray->At(bachLabels));
2485   if (!bachelorMC) return -1;
2486   if (TMath::Abs(bachelorMC->GetPdgCode())!=2212) return -1;
2487   Int_t indexMotherBach = bachelorMC->GetMother();
2488   if (indexMotherBach==-1) return -1;
2489
2490   Int_t pdgDg2prong[2] = {211,211};
2491   Int_t lab2Prong = theV0->MatchToMC(310,mcArray,2,pdgDg2prong); // the V0
2492   if(lab2Prong<0) return -1;
2493   AliAODMCParticle*partK0S = dynamic_cast<AliAODMCParticle*>(mcArray->At(lab2Prong));
2494   if (!partK0S) return -1;
2495   Int_t indexMotherK0S = partK0S->GetMother();
2496   if (indexMotherK0S==-1) return -1;
2497   AliAODMCParticle*partK0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(indexMotherK0S));
2498   if (!partK0) return -1;
2499   Int_t indexMotherK0 = partK0->GetMother();
2500   if (indexMotherK0==-1) return -1;
2501
2502   if (indexMotherBach!=indexMotherK0) return -1; // p e K0S sono fratelli
2503
2504   AliAODMCParticle*partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(indexMotherK0));
2505   if (!partLc) return -1;
2506   Int_t ndg2 = partLc->GetDaughter(1)-partLc->GetDaughter(0)+1;
2507   if (ndg2==2) return -1;
2508
2509   TString stringaCheck = Form(">>>>>>>> %d -> ",partLc->GetPdgCode());
2510   for(Int_t ii=0; ii<ndg2; ii++) {
2511     AliAODMCParticle* partDau=(AliAODMCParticle*)(mcArray->At(partLc->GetDaughter(0)+ii));
2512     stringaCheck.Append(Form("  %d",partDau->GetPdgCode()));
2513   }
2514   printf("%s \n",stringaCheck.Data());
2515
2516   return indexMotherBach;
2517
2518 }
2519 //-----------------------------------------------------------------------------
2520 Int_t AliAnalysisTaskSELc2V0bachelor::SearchForCommonMother(TClonesArray *mcArray,
2521                                                             Int_t dgLabels[10],Int_t ndg,
2522                                                             Int_t &ndgCk, Int_t *pdgDg, Int_t &absLabelMother) const
2523 {
2524   //
2525   // Check if this candidate is matched to a MC signal
2526   // If no, return 0
2527   // If yes, return pdgCode of particle
2528   //
2529
2530   Int_t lab=-1,labMother=-1,pdgMother=0;
2531   AliAODMCParticle *part=0;
2532   AliAODMCParticle *mother=0;
2533
2534   // loop on daughter labels
2535   TArrayI **labelMother = new TArrayI*[ndg];
2536   for(Int_t i=0; i<ndg; i++) labelMother[i] = new TArrayI(0);
2537   for(Int_t i=0; i<ndg; i++) {
2538     lab = TMath::Abs(dgLabels[i]);
2539     if(lab<0) {
2540       AliDebug(2,Form("daughter with negative label %d\n",lab));
2541       delete [] labelMother;
2542       return 0;
2543     }
2544     part = (AliAODMCParticle*)mcArray->At(lab);
2545     if(!part) {
2546       AliDebug(2,"no MC particle\n");
2547       delete [] labelMother;
2548       return 0;
2549     }
2550
2551     mother = part;
2552     while(mother->GetMother()>=0) {
2553       labMother=mother->GetMother();
2554       mother = (AliAODMCParticle*)mcArray->At(labMother);
2555       if(!mother) {
2556         AliDebug(2,"no MC mother particle\n");
2557         break;
2558       }
2559       pdgMother = TMath::Abs(mother->GetPdgCode());
2560       if (pdgMother<10 || (pdgMother>18 && pdgMother<111)) {
2561         break;
2562       }
2563       labelMother[i]->Set(labelMother[i]->GetSize()+1);
2564       labelMother[i]->AddAt(labMother,labelMother[i]->GetSize()-1);
2565     }
2566
2567   } // end loop on daughters
2568
2569   for(Int_t i=0; i<ndg; i++) {
2570     AliAODMCParticle*part0 = (AliAODMCParticle*)mcArray->At(TMath::Abs(dgLabels[i]));
2571     AliInfo(Form("part[%d]->GetLabel()=%d(%d) | ",i,dgLabels[i],part0->GetPdgCode()));
2572     AliInfo(Form("labelMother[%d] = ",i));
2573     for (Int_t jj=0;jj<labelMother[i]->GetSize(); jj++)
2574       AliInfo(Form("%d, ",labelMother[i]->At(jj)));
2575     AliInfo("\n");
2576   }
2577
2578   Int_t pdgToBeReturned=0;
2579   Bool_t found=kFALSE;
2580   for (Int_t ii=0;ii<labelMother[0]->GetSize(); ii++) {
2581     for (Int_t jj=0;jj<labelMother[1]->GetSize(); jj++) {
2582       for (Int_t kk=0;kk<labelMother[2]->GetSize(); kk++) {
2583         if (labelMother[0]->At(ii)==labelMother[1]->At(jj) &&
2584             labelMother[1]->At(jj)==labelMother[2]->At(kk) &&
2585             labelMother[0]->At(ii)!=0 && labelMother[0]->At(ii)!=1 && !found) {
2586           mother = (AliAODMCParticle*)mcArray->At(labelMother[0]->At(ii));
2587           pdgToBeReturned=mother->GetPdgCode();
2588           absLabelMother=labelMother[0]->At(ii);
2589           AliDebug(2,Form("FOUND label for the mother of this candidate: %d (PDG=%d)\n",labelMother[0]->At(ii),pdgToBeReturned));
2590           mother->Print();
2591           found = kTRUE;
2592           ndgCk=3;
2593           pdgDg = new Int_t[ndgCk];
2594           for (Int_t aa=0; aa<ndgCk; aa++) {
2595             AliAODMCParticle *partMC = (AliAODMCParticle*)mcArray->At(dgLabels[aa]);
2596             pdgDg[aa]=partMC->GetPdgCode();
2597           }
2598           break;
2599         }
2600       }
2601     }
2602   }
2603
2604   delete [] labelMother;
2605   if(pdgDg) delete [] pdgDg;
2606
2607   return pdgToBeReturned;
2608
2609 }