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