PWG3/vertexingHF/AliAnalysisTaskSELambdac.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /////////////////////////////////////////////////////////////
  17 //
  18 // AliAnalysisTaskSE for the extraction of signal(e.g D+) of heavy flavor
  19 // decay candidates with the MC truth.
  20 // Authors: Renu Bala, bala@to.infn.it
  21 // F. Prino, prino@to.infn.it
  22 // G. Ortona, ortona@to.infn.it
  23 /////////////////////////////////////////////////////////////
  24
  25 #include <TClonesArray.h>
  26 #include <TNtuple.h>
  27 #include <TCanvas.h>
  28 #include <TList.h>
  29 #include <TString.h>
  30 #include <TH1F.h>
  31 #include <TH2F.h>
  32 #include <TDatabasePDG.h>
  33
  34 #include "AliAnalysisManager.h"
  35 #include "AliAODHandler.h"
  36 #include "AliAODEvent.h"
  37 #include "AliAODVertex.h"
  38 #include "AliAODTrack.h"
  39 #include "AliAODMCHeader.h"
  40 #include "AliAODMCParticle.h"
  41 #include "AliAODRecoDecayHF3Prong.h"
  42 #include "AliAnalysisVertexingHF.h"
  43 #include "AliAnalysisTaskSE.h"
  44 #include "AliAnalysisTaskSELambdac.h"
  45 #include "AliKFParticle.h"
  46 #include "AliAODPidHF.h"
  47 #include "AliRDHFCutsLctopKpi.h"
  48 #include "AliRDHFCuts.h"
  49
  50 ClassImp(AliAnalysisTaskSELambdac)
  51
  52
  53 //________________________________________________________________________
  54 AliAnalysisTaskSELambdac::AliAnalysisTaskSELambdac():
  55 AliAnalysisTaskSE(),
  56 fOutput(0),
  57 fHistNEvents(0),
  58 fWellIDProt3sig(0),
  59 fWellIDProt2sig(0),
  60 fWellIDProt3sigSe(0),
  61 fWellIDProt2sigSe(0),
  62 fWellIDProt3sigRe(0),
  63 fWellIDProt2sigRe(0),
  64 fWellIDKaon3sig(0),
  65 fWellIDKaon3sigSe(0),
  66 fWellIDKaon3sigRe(0),
  67 fWellIDKaon2sig(0),
  68 fWellIDKaon2sigSe(0),
  69 fWellIDKaon2sigRe(0),
  70 fRealProt(0),
  71 fRealKaon(0),
  72 fFakeProt3sig(0),
  73 fFakeProt2sig(0),
  74 fFakeProt3sigSe(0),
  75 fFakeProt2sigSe(0),
  76 fFakeProt3sigRe(0),
  77 fFakeProt2sigRe(0),
  78 fFakeKaon3sig(0),
  79 fFakeKaon3sigSe(0),
  80 fFakeKaon3sigRe(0),
  81 fFakeKaon2sig(0),
  82 fFakeKaon2sigSe(0),
  83 fFakeKaon2sigRe(0),
  84 fTPCSignal3Sigma(0),
  85 fTPCSignal3SigmaReK(0),
  86 fTPCSignal3SigmaSedK(0),
  87 fTPCSignal3SigmaRep(0),
  88 fTPCSignal3SigmaSedp(0),
  89 fTPCSignal2Sigma(0),
  90 fTPCSignal2SigmaReK(0),
  91 fTPCSignal2SigmaSedK(0),
  92 fTPCSignal2SigmaRep(0),
  93 fTPCSignal2SigmaSedp(0),
  94 fNtupleLambdac(0),
  95 fUpmasslimit(2.486),
  96 fLowmasslimit(2.086),
  97 fNPtBins(0),
  98 fRDCutsAnalysis(0),
  99 fRDCutsProduction(0),
 100 fListCuts(0),
 101 fFillNtuple(kFALSE),
 102 fReadMC(kFALSE),
 103 fMCPid(kFALSE),
 104 fRealPid(kFALSE),
 105 fResPid(kTRUE),
 106 fUseKF(kFALSE),
 107 fVHF(0)
 108 {
 109    // Default constructor
 110 }
 111
 112 //________________________________________________________________________
 113 AliAnalysisTaskSELambdac::AliAnalysisTaskSELambdac(const char *name,Bool_t fillNtuple,AliRDHFCutsLctopKpi *lccutsana,AliRDHFCutsLctopKpi *lccutsprod):
 114 AliAnalysisTaskSE(name),
 115 fOutput(0),
 116 fHistNEvents(0),
 117 fWellIDProt3sig(0),
 118 fWellIDProt2sig(0),
 119 fWellIDProt3sigSe(0),
 120 fWellIDProt2sigSe(0),
 121 fWellIDProt3sigRe(0),
 122 fWellIDProt2sigRe(0),
 123 fWellIDKaon3sig(0),
 124 fWellIDKaon3sigSe(0),
 125 fWellIDKaon3sigRe(0),
 126 fWellIDKaon2sig(0),
 127 fWellIDKaon2sigSe(0),
 128 fWellIDKaon2sigRe(0),
 129 fRealProt(0),
 130 fRealKaon(0),
 131 fFakeProt3sig(0),
 132 fFakeProt2sig(0),
 133 fFakeProt3sigSe(0),
 134 fFakeProt2sigSe(0),
 135 fFakeProt3sigRe(0),
 136 fFakeProt2sigRe(0),
 137 fFakeKaon3sig(0),
 138 fFakeKaon3sigSe(0),
 139 fFakeKaon3sigRe(0),
 140 fFakeKaon2sig(0),
 141 fFakeKaon2sigSe(0),
 142 fFakeKaon2sigRe(0),
 143 fTPCSignal3Sigma(0),
 144 fTPCSignal3SigmaReK(0),
 145 fTPCSignal3SigmaSedK(0),
 146 fTPCSignal3SigmaRep(0),
 147 fTPCSignal3SigmaSedp(0),
 148 fTPCSignal2Sigma(0),
 149 fTPCSignal2SigmaReK(0),
 150 fTPCSignal2SigmaSedK(0),
 151 fTPCSignal2SigmaRep(0),
 152 fTPCSignal2SigmaSedp(0),
 153 fNtupleLambdac(0),
 154 fUpmasslimit(2.486),
 155 fLowmasslimit(2.086),
 156 fNPtBins(0),
 157 fRDCutsAnalysis(lccutsana),
 158 fRDCutsProduction(lccutsprod),
 159 fListCuts(0),
 160 fFillNtuple(fillNtuple),
 161 fReadMC(kFALSE),
 162 fMCPid(kFALSE),
 163 fRealPid(kTRUE),
 164 fResPid(kFALSE),
 165 fUseKF(kFALSE),
 166 fVHF(0)
 167 {
 168   //Double_t ptlim[5]={0.,2.,3.,5,9999999.};
 169   //SetPtBinLimit(5, ptlim);
 170    SetPtBinLimit(fRDCutsAnalysis->GetNPtBins()+1,fRDCutsAnalysis->GetPtBinLimits());
 171   // Default constructor
 172    // Output slot #1 writes into a TList container
 173   DefineOutput(1,TList::Class());  //My private output
 174   DefineOutput(2,TList::Class());
 175
 176   if(fFillNtuple){
 177     // Output slot #2 writes into a TNtuple container
 178     DefineOutput(3,TNtuple::Class());  //My private output
 179   }
 180 }
 181
 182 //________________________________________________________________________
 183 AliAnalysisTaskSELambdac::~AliAnalysisTaskSELambdac()
 184 {
 185   // Destructor
 186   if (fOutput) {
 187     delete fOutput;
 188     fOutput = 0;
 189   }
 190   if (fVHF) {
 191     delete fVHF;
 192     fVHF = 0;
 193   }
 194
 195  if(fRDCutsAnalysis){
 196     delete fRDCutsAnalysis;
 197     fRDCutsAnalysis = 0;
 198  }
 199  if(fRDCutsProduction){
 200     delete fRDCutsProduction;
 201     fRDCutsProduction = 0;
 202  }
 203
 204  if (fListCuts) {
 205    delete fListCuts;
 206    fListCuts = 0;
 207  }
 208 }
 209 //_________________________________________________________________
 210 void  AliAnalysisTaskSELambdac::SetMassLimits(Float_t range){
 211   fUpmasslimit = 2.286+range;
 212   fLowmasslimit = 2.286-range;
 213 }
 214 //_________________________________________________________________
 215 void  AliAnalysisTaskSELambdac::SetMassLimits(Float_t lowlimit, Float_t uplimit){
 216   if(uplimit>lowlimit)
 217     {
 218       fUpmasslimit = lowlimit;
 219       fLowmasslimit = uplimit;
 220     }
 221 }
 222
 223
 224 //________________________________________________________________________
 225 void AliAnalysisTaskSELambdac::SetPtBinLimit(Int_t n, Float_t* lim){
 226   // define pt bins for analysis
 227   if(n>kMaxPtBins){
 228     printf("Max. number of Pt bins = %d\n",kMaxPtBins);
 229     fNPtBins=kMaxPtBins;
 230     fArrayBinLimits[0]=0.;
 231     fArrayBinLimits[1]=2.;
 232     fArrayBinLimits[2]=3.;
 233     fArrayBinLimits[3]=5.;
 234     for(Int_t i=4; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
 235   }else{
 236     fNPtBins=n-1;
 237     fArrayBinLimits[0]=lim[0];
 238     for(Int_t i=1; i<fNPtBins+1; i++)
 239       if(lim[i]>fArrayBinLimits[i-1]){
 240         fArrayBinLimits[i]=lim[i];
 241       }
 242       else {
 243         fArrayBinLimits[i]=fArrayBinLimits[i-1];
 244       }
 245     for(Int_t i=fNPtBins; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
 246   }
 247   if(fDebug > 1){
 248     printf("Number of Pt bins = %d\n",fNPtBins);
 249     for(Int_t i=0; i<fNPtBins+1; i++) printf(" Bin%d = %8.2f-%8.2f\n",i,fArrayBinLimits[i],fArrayBinLimits[i+1]);
 250   }
 251 }
 252 //_________________________________________________________________
 253 Double_t  AliAnalysisTaskSELambdac::GetPtBinLimit(Int_t ibin){
 254   if(ibin>fNPtBins)return -1;
 255   return fArrayBinLimits[ibin];
 256 }
 257
 258 //_________________________________________________________________
 259 void AliAnalysisTaskSELambdac::Init()
 260 {
 261   // Initialization
 262
 263   if(fDebug > 1) printf("AnalysisTaskSELambdac::Init() \n");
 264
 265   fListCuts=new TList();
 266   AliRDHFCutsLctopKpi *production = new AliRDHFCutsLctopKpi();
 267     production=fRDCutsProduction;
 268
 269   AliRDHFCutsLctopKpi *analysis = new AliRDHFCutsLctopKpi();
 270   analysis=fRDCutsAnalysis;
 271
 272    fListCuts->Add(analysis);
 273    fListCuts->Add(production);
 274    PostData(2,fListCuts);
 275   return;
 276 }
 277
 278 //________________________________________________________________________
 279 void AliAnalysisTaskSELambdac::UserCreateOutputObjects()
 280 {
 281   // Create the output container
 282   //
 283   if(fDebug > 1) printf("AnalysisTaskSELambdac::UserCreateOutputObjects() \n");
 284
 285   // Several histograms are more conveniently managed in a TList
 286   fOutput = new TList();
 287   fOutput->SetOwner();
 288   fOutput->SetName("OutputHistos");
 289
 290   TString hisname;
 291   Int_t index=0;
 292   Int_t indexLS=0;
 293   for(Int_t i=0;i<fNPtBins;i++){
 294
 295     index=GetHistoIndex(i);
 296     indexLS=GetLSHistoIndex(i);
 297
 298     hisname.Form("hMassPt%d",i);
 299     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 300     fMassHist[index]->Sumw2();
 301     hisname.Form("hCosPAllPt%d",i);
 302     fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
 303     fCosPHist[index]->Sumw2();
 304     hisname.Form("hDLenAllPt%d",i);
 305     fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
 306     fDLenHist[index]->Sumw2();
 307     hisname.Form("hSumd02AllPt%d",i);
 308     fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
 309     fSumd02Hist[index]->Sumw2();
 310     hisname.Form("hSigVertAllPt%d",i);
 311     fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 312     fSigVertHist[index]->Sumw2();
 313     hisname.Form("hPtMaxAllPt%d",i);
 314     fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
 315     fPtMaxHist[index]->Sumw2();
 316
 317     hisname.Form("hDCAAllPt%d",i);
 318     fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 319     fDCAHist[index]->Sumw2();
 320
 321
 322
 323     hisname.Form("hMassPt%dTC",i);
 324     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 325     fMassHistTC[index]->Sumw2();
 326
 327
 328
 329
 330
 331     hisname.Form("hCosPAllPt%dLS",i);
 332     fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
 333     fCosPHistLS[index]->Sumw2();
 334     hisname.Form("hDLenAllPt%dLS",i);
 335     fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
 336     fDLenHistLS[index]->Sumw2();
 337     hisname.Form("hSumd02AllPt%dLS",i);
 338     fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
 339     fSumd02HistLS[index]->Sumw2();
 340     hisname.Form("hSigVertAllPt%dLS",i);
 341     fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 342     fSigVertHistLS[index]->Sumw2();
 343     hisname.Form("hPtMaxAllPt%dLS",i);
 344     fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
 345     fPtMaxHistLS[index]->Sumw2();
 346
 347     hisname.Form("hDCAAllPt%dLS",i);
 348     fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 349     fDCAHistLS[index]->Sumw2();
 350
 351     hisname.Form("hLSPt%dLC",i);
 352     fMassHistLS[indexLS] = new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 353     fMassHistLS[indexLS]->Sumw2();
 354
 355     hisname.Form("hLSPt%dTC",i);
 356     fMassHistLSTC[indexLS] = new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 357     fMassHistLSTC[indexLS]->Sumw2();
 358
 359
 360
 361     index=GetSignalHistoIndex(i);
 362     indexLS++;
 363     hisname.Form("hSigPt%d",i);
 364     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 365     fMassHist[index]->Sumw2();
 366     hisname.Form("hCosPSigPt%d",i);
 367     fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
 368     fCosPHist[index]->Sumw2();
 369     hisname.Form("hDLenSigPt%d",i);
 370     fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
 371     fDLenHist[index]->Sumw2();
 372     hisname.Form("hSumd02SigPt%d",i);
 373     fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
 374     fSumd02Hist[index]->Sumw2();
 375     hisname.Form("hSigVertSigPt%d",i);
 376     fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 377     fSigVertHist[index]->Sumw2();
 378     hisname.Form("hPtMaxSigPt%d",i);
 379     fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
 380     fPtMaxHist[index]->Sumw2();
 381
 382     hisname.Form("hDCASigPt%d",i);
 383     fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 384     fDCAHist[index]->Sumw2();
 385
 386
 387     hisname.Form("hSigPt%dTC",i);
 388     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 389     fMassHistTC[index]->Sumw2();
 390
 391     hisname.Form("hLSPt%dLCnw",i);
 392     fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 393     fMassHistLS[indexLS]->Sumw2();
 394     hisname.Form("hLSPt%dTCnw",i);
 395     fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 396     fMassHistLSTC[indexLS]->Sumw2();
 397
 398
 399
 400     hisname.Form("hCosPSigPt%dLS",i);
 401     fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
 402     fCosPHistLS[index]->Sumw2();
 403     hisname.Form("hDLenSigPt%dLS",i);
 404     fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
 405     fDLenHistLS[index]->Sumw2();
 406     hisname.Form("hSumd02SigPt%dLS",i);
 407     fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
 408     fSumd02HistLS[index]->Sumw2();
 409     hisname.Form("hSigVertSigPt%dLS",i);
 410     fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 411     fSigVertHistLS[index]->Sumw2();
 412     hisname.Form("hPtMaxSigPt%dLS",i);
 413     fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
 414     fPtMaxHistLS[index]->Sumw2();
 415
 416     hisname.Form("hDCASigPt%dLS",i);
 417     fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 418     fDCAHistLS[index]->Sumw2();
 419
 420
 421
 422     index=GetBackgroundHistoIndex(i);
 423     indexLS++;
 424     hisname.Form("hBkgPt%d",i);
 425     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 426     fMassHist[index]->Sumw2();
 427     hisname.Form("hCosPBkgPt%d",i);
 428     fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
 429     fCosPHist[index]->Sumw2();
 430     hisname.Form("hDLenBkgPt%d",i);
 431     fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
 432     fDLenHist[index]->Sumw2();
 433     hisname.Form("hSumd02BkgPt%d",i);
 434     fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
 435     fSumd02Hist[index]->Sumw2();
 436     hisname.Form("hSigVertBkgPt%d",i);
 437     fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 438     fSigVertHist[index]->Sumw2();
 439     hisname.Form("hPtMaxBkgPt%d",i);
 440     fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
 441     fPtMaxHist[index]->Sumw2();
 442
 443     hisname.Form("hDCABkgPt%d",i);
 444     fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 445     fDCAHist[index]->Sumw2();
 446
 447
 448     hisname.Form("hBkgPt%dTC",i);
 449     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 450     fMassHistTC[index]->Sumw2();
 451
 452     hisname.Form("hLSPt%dLCntrip",i);
 453     fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 454     fMassHistLS[indexLS]->Sumw2();
 455     hisname.Form("hLSPt%dTCntrip",i);
 456     fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 457     fMassHistLSTC[indexLS]->Sumw2();
 458
 459
 460     hisname.Form("hCosPBkgPt%dLS",i);
 461     fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
 462     fCosPHistLS[index]->Sumw2();
 463     hisname.Form("hDLenBkgPt%dLS",i);
 464     fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
 465     fDLenHistLS[index]->Sumw2();
 466     hisname.Form("hSumd02BkgPt%dLS",i);
 467     fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
 468     fSumd02HistLS[index]->Sumw2();
 469     hisname.Form("hSigVertBkgPt%dLS",i);
 470     fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 471     fSigVertHistLS[index]->Sumw2();
 472     hisname.Form("hPtMaxBkgPt%dLS",i);
 473     fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
 474     fPtMaxHistLS[index]->Sumw2();
 475     hisname.Form("hDCABkgPt%dLS",i);
 476     fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
 477     fDCAHistLS[index]->Sumw2();
 478
 479
 480     indexLS++;
 481     hisname.Form("hLSPt%dLCntripsinglecut",i);
 482     fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 483     fMassHistLS[indexLS]->Sumw2();
 484     hisname.Form("hLSPt%dTCntripsinglecut",i);
 485     fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 486     fMassHistLSTC[indexLS]->Sumw2();
 487
 488     indexLS++;
 489     hisname.Form("hLSPt%dLCspc",i);
 490     fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 491     fMassHistLS[indexLS]->Sumw2();
 492     hisname.Form("hLSPt%dTCspc",i);
 493     fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
 494     fMassHistLSTC[indexLS]->Sumw2();
 495   }
 496
 497   fRealProt=new TH1F("fRealProt","fRealProt",100,0, 10);
 498   fRealKaon=new TH1F("fRealKaon","fRealKaon",100,0, 10);
 499   fWellIDProt3sig=new TH1F("fWellIDProt3sig","fWellIDProt3sig",100,0, 10);
 500   fWellIDProt2sig=new TH1F("fWellIDProt2sig","fWellIDProt2sig",100,0, 10);
 501   fWellIDProt3sigSe=new TH1F("fWellIDProt3sigSe","fWellIDProt3sigSe",100,0, 10);
 502   fWellIDProt2sigSe=new TH1F("fWellIDProt2sigSe","fWellIDProt2sigSe",100,0, 10);
 503   fWellIDProt3sigRe=new TH1F("fWellIDProt3sigRe","fWellIDProt3sigRe",100,0, 10);
 504   fWellIDProt2sigRe=new TH1F("fWellIDProt2sigRe","fWellIDProt2sigRe",100,0, 10);
 505   fWellIDKaon3sig=new TH1F("fWellIDKaon3sig","fWellIDKaon3sig",100,0, 10);
 506   fWellIDKaon2sig=new TH1F("fWellIDKaon2sig","fWellIDKaon2sig",100,0, 10);
 507   fWellIDKaon3sigSe=new TH1F("fWellIDKaon3sigSe","fWellIDKaon3sigSe",100,0, 10);
 508   fWellIDKaon2sigSe=new TH1F("fWellIDKaon2sigSe","fWellIDKaon2sigSe",100,0, 10);
 509   fWellIDKaon3sigRe=new TH1F("fWellIDKaon3sigRe","fWellIDKaon3sigRe",100,0, 10);
 510   fWellIDKaon2sigRe=new TH1F("fWellIDKaon2sigRe","fWellIDKaon2sigRe",100,0, 10);
 511   fTPCSignal3Sigma=new TH2F("fTPCSignal3Sigma","fTPCSignal3Sigma",100,0, 10, 100, 0, 100);
 512   fTPCSignal3SigmaReK=new TH2F("fTPCSignal3SigmaReK","fTPCSignal3SigmaReK",100,0, 10, 100, 0, 100);
 513   fTPCSignal3SigmaSedK=new TH2F("fTPCSignal3SigmaSedK","fTPCSignal3SigmaSedK",100,0, 10, 100, 0, 100);
 514   fTPCSignal3SigmaRep=new TH2F("fTPCSignal3SigmaRep","fTPCSignal3SigmaRep",100,0, 10, 100, 0, 100);
 515   fTPCSignal3SigmaSedp=new TH2F("fTPCSignal3SigmaSedp","fTPCSignal3SigmaSedp",100,0, 10, 100, 0, 100);
 516   fTPCSignal2Sigma=new TH2F("fTPCSignal2Sigma","fTPCSignal2Sigma",100,0, 10, 100, 0, 100);
 517   fTPCSignal2SigmaReK=new TH2F("fTPCSignal2SigmaReK","fTPCSignal2SigmaReK",100,0, 10, 100, 0, 100);
 518   fTPCSignal2SigmaSedK=new TH2F("fTPCSignal2SigmaSedK","fTPCSignal2SigmaSedK",100,0, 10, 100, 0, 100);
 519   fTPCSignal2SigmaRep=new TH2F("fTPCSignal2SigmaRep","fTPCSignal2SigmaRep",100,0, 10, 100, 0, 100);
 520   fTPCSignal2SigmaSedp=new TH2F("fTPCSignal2SigmaSedp","fTPCSignal2SigmaSedp",100,0, 10, 100, 0, 100);
 521   fFakeProt3sig=new TH1F("fFakeProt3sig","fFakeProt3sig",100,0, 10);
 522   fFakeProt2sig=new TH1F("fFakeProt2sig","fFakeProt2sig",100,0, 10);
 523   fFakeProt3sigSe=new TH1F("fFakeProt3sigSe","fFakeProt3sigSe",100,0, 10);
 524   fFakeProt2sigSe=new TH1F("fFakeProt2sigSe","fFakeProt2sigSe",100,0, 10);
 525   fFakeProt3sigRe=new TH1F("fFakeProt3sigRe","fFakeProt3sigRe",100,0, 10);
 526   fFakeProt2sigRe=new TH1F("fFakeProt2sigRe","fFakeProt2sigRe",100,0, 10);
 527   fFakeKaon3sig=new TH1F("fFakeKaon3sig","fFakeKaon3sig",100,0, 10);
 528   fFakeKaon2sig=new TH1F("fFakeKaon2sig","fFakeKaon2sig",100,0, 10);
 529   fFakeKaon3sigSe=new TH1F("fFakeKaon3sigSe","fFakeKaon3sigSe",100,0, 10);
 530   fFakeKaon2sigSe=new TH1F("fFakeKaon2sigSe","fFakeKaon2sigSe",100,0, 10);
 531   fFakeKaon3sigRe=new TH1F("fFakeKaon3sigRe","fFakeKaon3sigRe",100,0, 10);
 532   fFakeKaon2sigRe=new TH1F("fFakeKaon2sigRe","fFakeKaon2sigRe",100,0, 10);
 533
 534   for(Int_t i=0; i<3*fNPtBins; i++){
 535     fOutput->Add(fMassHist[i]);
 536     fOutput->Add(fMassHistTC[i]);
 537     fOutput->Add(fCosPHist[i]);
 538     fOutput->Add(fDLenHist[i]);
 539     fOutput->Add(fSumd02Hist[i]);
 540     fOutput->Add(fSigVertHist[i]);
 541     fOutput->Add(fPtMaxHist[i]);
 542     fOutput->Add(fDCAHist[i]);
 543   }
 544
 545     fOutput->Add(fTPCSignal3Sigma);
 546     fOutput->Add(fTPCSignal3SigmaReK);
 547     fOutput->Add(fTPCSignal3SigmaSedK);
 548     fOutput->Add(fTPCSignal3SigmaRep);
 549     fOutput->Add(fTPCSignal3SigmaSedp);
 550     fOutput->Add(fTPCSignal2Sigma);
 551     fOutput->Add(fTPCSignal2SigmaReK);
 552     fOutput->Add(fTPCSignal2SigmaSedK);
 553     fOutput->Add(fTPCSignal2SigmaRep);
 554     fOutput->Add(fTPCSignal2SigmaSedp);
 555     fOutput->Add(fWellIDProt3sig);
 556     fOutput->Add(fWellIDProt2sig);
 557     fOutput->Add(fWellIDProt3sigSe);
 558     fOutput->Add(fWellIDProt2sigSe);
 559     fOutput->Add(fWellIDProt3sigRe);
 560     fOutput->Add(fWellIDProt2sigRe);
 561     fOutput->Add(fWellIDKaon3sig);
 562     fOutput->Add(fWellIDKaon2sig);
 563     fOutput->Add(fWellIDKaon3sigSe);
 564     fOutput->Add(fWellIDKaon2sigSe);
 565     fOutput->Add(fWellIDKaon3sigRe);
 566     fOutput->Add(fWellIDKaon2sigRe);
 567
 568     fOutput->Add(fRealProt);
 569     fOutput->Add(fRealKaon);
 570     fOutput->Add(fFakeProt3sig);
 571     fOutput->Add(fFakeProt2sig);
 572     fOutput->Add(fFakeProt3sigSe);
 573     fOutput->Add(fFakeProt2sigSe);
 574     fOutput->Add(fFakeProt3sigRe);
 575     fOutput->Add(fFakeProt2sigRe);
 576     fOutput->Add(fFakeKaon3sig);
 577     fOutput->Add(fFakeKaon2sig);
 578     fOutput->Add(fFakeKaon3sigSe);
 579     fOutput->Add(fFakeKaon2sigSe);
 580     fOutput->Add(fFakeKaon3sigRe);
 581     fOutput->Add(fFakeKaon2sigRe);
 582   fHistNEvents = new TH1F("fHistNEvents", "Number of processed events; ; Events",3,-1.5,1.5);
 583   fHistNEvents->Sumw2();
 584   fHistNEvents->SetMinimum(0);
 585   fOutput->Add(fHistNEvents);
 586
 587
 588
 589   if(fFillNtuple){
 590     //OpenFile(3); // 2 is the slot number of the ntuple
 591
 592     fNtupleLambdac = new TNtuple("fNtupleLambdac","D +","pdg:Px:Py:Pz:PtTrue:VxTrue:VyTrue:VzTrue:Ptpi:PtK:Ptpi2:PtRec:PointingAngle:DecLeng:VxRec:VyRec:VzRec:InvMass:sigvert:d0Pi:d0K:d0Pi2:dca:d0square");
 593
 594   }
 595
 596   return;
 597 }
 598
 599 //________________________________________________________________________
 600 void AliAnalysisTaskSELambdac::UserExec(Option_t */*option*/)
 601 {
 602   // Execute analysis for current event:
 603   // heavy flavor candidates association to MC truth
 604
 605    AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
 606    //tmp!
 607  fHistNEvents->Fill(0); // count event
 608   // Post the data already here
 609   PostData(1,fOutput);
 610
 611   TClonesArray *array3Prong = 0;
 612   TClonesArray *arrayLikeSign =0;
 613   if(!aod && AODEvent() && IsStandardAOD()) {
 614     // In case there is an AOD handler writing a standard AOD, use the AOD
 615     // event in memory rather than the input (ESD) event.
 616     aod = dynamic_cast<AliAODEvent*> (AODEvent());
 617     // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
 618     // have to taken from the AOD event hold by the AliAODExtension
 619     AliAODHandler* aodHandler = (AliAODHandler*)
 620       ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
 621     if(aodHandler->GetExtensions()) {
 622       AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
 623       AliAODEvent *aodFromExt = ext->GetAOD();
 624       array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
 625       arrayLikeSign=(TClonesArray*)aodFromExt->GetList()->FindObject("LikeSign3Prong");
 626     }
 627   } else {
 628     array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
 629     arrayLikeSign=(TClonesArray*)aod->GetList()->FindObject("LikeSign3Prong");
 630   }
 631
 632   if(!array3Prong) {
 633     printf("AliAnalysisTaskSELambdac::UserExec: Charm3Prong branch not found!\n");
 634     return;
 635   }
 636   if(!arrayLikeSign) {
 637     printf("AliAnalysisTaskSELambdac::UserExec: LikeSign3Prong branch not found!\n");
 638   //  return;
 639   }
 640
 641
 642   TClonesArray *arrayMC=0;
 643   AliAODMCHeader *mcHeader=0;
 644
 645   // AOD primary vertex
 646   AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
 647
 648   // load MC particles
 649   if(fReadMC){
 650
 651     arrayMC =  (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 652     if(!arrayMC) {
 653       printf("AliAnalysisTaskSELambdac::UserExec: MC particles branch not found!\n");
 654       return;
 655     }
 656
 657   // load MC header
 658     mcHeader =  (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 659     if(!mcHeader) {
 660     printf("AliAnalysisTaskSELambdac::UserExec: MC header branch not found!\n");
 661     return;
 662     }
 663   }
 664
 665   Int_t n3Prong = array3Prong->GetEntriesFast();
 666
 667
 668   Int_t nOS=0;
 669   Int_t index;
 670   for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
 671     AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
 672
 673
 674     Bool_t unsetvtx=kFALSE;
 675     if(!d->GetOwnPrimaryVtx()){
 676       d->SetOwnPrimaryVtx(vtx1);
 677       unsetvtx=kTRUE;
 678     }
 679
 680     if(fRDCutsProduction->IsSelected(d,AliRDHFCuts::kCandidate)) {
 681       Int_t iPtBin = -1;
 682       Double_t ptCand = d->Pt();
 683
 684       for(Int_t ibin=0;ibin<fNPtBins&&iPtBin<0&&ptCand>fArrayBinLimits[0]&&ptCand<fArrayBinLimits[fNPtBins];ibin++){
 685         if(ptCand<fArrayBinLimits[ibin+1])iPtBin=ibin;
 686       }
 687
 688       Bool_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kCandidate);
 689       Int_t labDp=-1;
 690       Float_t deltaPx=0.;
 691       Float_t deltaPy=0.;
 692       Float_t deltaPz=0.;
 693       Float_t truePt=0.;
 694       Float_t xDecay=0.;
 695       Float_t yDecay=0.;
 696       Float_t zDecay=0.;
 697       Float_t pdgCode=-2;
 698       if(fReadMC){
 699         labDp = MatchToMCLambdac(d,arrayMC);
 700         //
 701         if(labDp>0){
 702           AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
 703           AliAODMCParticle *dg0 = (AliAODMCParticle*)arrayMC->At(partDp->GetDaughter(0));
 704           deltaPx=partDp->Px()-d->Px();
 705           deltaPy=partDp->Py()-d->Py();
 706           deltaPz=partDp->Pz()-d->Pz();
 707           truePt=partDp->Pt();
 708           xDecay=dg0->Xv();
 709           yDecay=dg0->Yv();
 710           zDecay=dg0->Zv();
 711           pdgCode=TMath::Abs(partDp->GetPdgCode());
 712         }else{
 713           pdgCode=-1;
 714         }
 715       }
 716
 717       Double_t invMasspKpi=-1.;
 718       Double_t invMasspiKp=-1.;
 719       Int_t pdgs[3]={0,0,0};
 720       Double_t field=aod->GetMagneticField();
 721       if(fReadMC && fMCPid){
 722
 723        if(IspKpiMC(d,arrayMC,pdgs)) {
 724         invMasspKpi=d->InvMassLcpKpi();
 725         if(fUseKF){
 726          pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
 727          if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
 728         }
 729        }
 730        if(IspiKpMC(d,arrayMC)) {
 731         invMasspiKp=d->InvMassLcpiKp();
 732         if(fUseKF){
 733          pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
 734          if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
 735         }
 736        }
 737       }
 738       if(fRealPid){
 739        if(IspKpiReal(d)) {
 740         invMasspKpi=d->InvMassLcpKpi();
 741         if(fUseKF){
 742          pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
 743          if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
 744         }
 745        }
 746        if(fReadMC) IspKpiMC(d,arrayMC,pdgs);
 747        if(IspiKpReal(d,pdgs)) {
 748         invMasspiKp=d->InvMassLcpiKp();
 749         if(fUseKF){
 750          pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
 751          if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
 752         }
 753        }
 754       }
 755       if(fResPid){
 756        if(IspKpiResonant(d,field)) {
 757         invMasspKpi=d->InvMassLcpKpi();
 758         if(fUseKF){
 759          pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
 760          if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
 761         }
 762        }
 763        if(IspiKpResonant(d,field)) {
 764         invMasspiKp=d->InvMassLcpiKp();
 765         if(fUseKF){
 766          pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
 767          if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
 768         }
 769        }
 770       }
 771       if(!fResPid && !fRealPid && !fMCPid){
 772        invMasspiKp=d->InvMassLcpiKp();
 773        if(fUseKF){
 774          pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
 775          if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
 776         }
 777        invMasspKpi=d->InvMassLcpKpi();
 778        if(fUseKF){
 779         pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
 780         if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
 781        }
 782       }
 783
 784       if(invMasspiKp<0. && invMasspKpi<0.) continue;
 785
 786
 787       Float_t tmp[24];
 788       if(fFillNtuple){
 789         tmp[0]=pdgCode;
 790         tmp[1]=deltaPx;
 791         tmp[2]=deltaPy;
 792         tmp[3]=deltaPz;
 793         tmp[4]=truePt;
 794         tmp[5]=xDecay;
 795         tmp[6]=yDecay;
 796         tmp[7]=zDecay;
 797         tmp[8]=d->PtProng(0);
 798         tmp[9]=d->PtProng(1);
 799         tmp[10]=d->PtProng(2);
 800         tmp[11]=d->Pt();
 801         tmp[12]=d->CosPointingAngle();
 802         tmp[13]=d->DecayLength();
 803         tmp[14]=d->Xv();
 804         tmp[15]=d->Yv();
 805         tmp[16]=d->Zv();
 806         if(invMasspiKp>0.) tmp[17]=invMasspiKp;
 807         if(invMasspKpi>0.) tmp[17]=invMasspKpi;
 808         tmp[18]=d->GetSigmaVert();
 809         tmp[19]=d->Getd0Prong(0);
 810         tmp[20]=d->Getd0Prong(1);
 811         tmp[21]=d->Getd0Prong(2);
 812         tmp[22]=d->GetDCA();
 813         tmp[23]=d->Prodd0d0();
 814         fNtupleLambdac->Fill(tmp);
 815         PostData(3,fNtupleLambdac);
 816       }
 817       Double_t dlen=d->DecayLength();
 818       Double_t cosp=d->CosPointingAngle();
 819       Double_t sumD02=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
 820       Double_t dca=d->GetDCA();
 821 Double_t ptmax=0;
 822       for(Int_t i=0;i<3;i++){
 823         if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
 824       }
 825       if(iPtBin>=0){
 826
 827         index=GetHistoIndex(iPtBin);
 828         if(invMasspiKp>0. && invMasspKpi>0.){
 829         if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp,0.5);
 830         if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi,0.5);
 831         }else{
 832          if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
 833          if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
 834         }
 835         fCosPHist[index]->Fill(cosp);
 836         fDLenHist[index]->Fill(dlen);
 837         fSumd02Hist[index]->Fill(sumD02);
 838         fPtMaxHist[index]->Fill(ptmax);
 839         fDCAHist[index]->Fill(dca);
 840
 841         if(passTightCuts){
 842          if(invMasspiKp>0. && invMasspKpi>0.){
 843           if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp,0.5);
 844           if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi,0.5);
 845          }else{
 846           if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
 847           if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
 848          }
 849         }
 850
 851         if(fReadMC){
 852           if(labDp>0) {
 853             index=GetSignalHistoIndex(iPtBin);
 854             if(invMasspiKp>0. && invMasspKpi>0.){
 855              if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp,0.5);
 856              if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi,0.5);
 857             }else{
 858              if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
 859              if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
 860             }
 861             fCosPHist[index]->Fill(cosp);
 862             fDLenHist[index]->Fill(dlen);
 863             fSumd02Hist[index]->Fill(sumD02);
 864             fPtMaxHist[index]->Fill(ptmax);
 865             fDCAHist[index]->Fill(dca);
 866             if(passTightCuts){
 867              if(invMasspiKp>0. && invMasspKpi>0.){
 868               if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp,0.5);
 869               if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi,0.5);
 870              }else{
 871               if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
 872               if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
 873              }
 874             }
 875
 876           }else{
 877             index=GetBackgroundHistoIndex(iPtBin);
 878             if(invMasspiKp>0. && invMasspKpi>0.){
 879              fMassHist[index]->Fill(invMasspiKp,0.5);
 880              fMassHist[index]->Fill(invMasspKpi,0.5);
 881             }else{
 882              if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
 883              if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
 884            }
 885             fCosPHist[index]->Fill(cosp);
 886             fDLenHist[index]->Fill(dlen);
 887             fSumd02Hist[index]->Fill(sumD02);
 888             fPtMaxHist[index]->Fill(ptmax);
 889             fDCAHist[index]->Fill(dca);
 890             if(passTightCuts){
 891             if(invMasspiKp>0. && invMasspKpi>0.){
 892                fMassHistTC[index]->Fill(invMasspiKp,0.5);
 893                fMassHistTC[index]->Fill(invMasspKpi,0.5);
 894               }else{
 895                if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
 896                if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
 897               }
 898             }
 899           }
 900         }
 901       }
 902       /*
 903       //start OS analysis
 904       if(labDp<0)fHistOSbkg->Fill(d->InvMassDplus());
 905       fHistOS->Fill(d->InvMassDplus());
 906       */
 907       nOS++;
 908     }
 909     if(unsetvtx) d->UnsetOwnPrimaryVtx();
 910   }
 911
 912   PostData(1,fOutput);
 913   return;
 914 }
 915
 916
 917
 918 //________________________________________________________________________
 919 void AliAnalysisTaskSELambdac::Terminate(Option_t */*option*/)
 920 {
 921   // Terminate analysis
 922   //
 923   if(fDebug > 1) printf("AnalysisTaskSELambdac: Terminate() \n");
 924
 925   fOutput = dynamic_cast<TList*> (GetOutputData(1));
 926   if (!fOutput) {
 927     printf("ERROR: fOutput not available\n");
 928     return;
 929   }
 930  fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
 931
 932  TString hisname;
 933  Int_t index=0;
 934
 935
 936  Int_t indexLS=0;
 937  for(Int_t i=0;i<fNPtBins;i++){
 938     index=GetHistoIndex(i);
 939     hisname.Form("hMassPt%d",i);
 940     fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 941      hisname.Form("hCosPAllPt%d",i);
 942     fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 943      hisname.Form("hDLenAllPt%d",i);
 944     fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 945      hisname.Form("hSumd02AllPt%d",i);
 946      fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 947      hisname.Form("hSigVertAllPt%d",i);
 948      fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 949      hisname.Form("hPtMaxAllPt%d",i);
 950      fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 951      hisname.Form("hDCAAllPt%d",i);
 952      fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 953      hisname.Form("hMassPt%dTC",i);
 954     fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 955
 956     index=GetSignalHistoIndex(i);
 957     hisname.Form("hSigPt%d",i);
 958     fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 959     hisname.Form("hCosPSigPt%d",i);
 960     fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 961     hisname.Form("hDLenSigPt%d",i);
 962     fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 963     hisname.Form("hSumd02SigPt%d",i);
 964     fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 965     hisname.Form("hSigVertSigPt%d",i);
 966     fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 967     hisname.Form("hPtMaxSigPt%d",i);
 968     fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 969     hisname.Form("hDCASigPt%d",i);
 970     fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 971
 972     hisname.Form("hSigPt%dTC",i);
 973     fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 974
 975     index=GetBackgroundHistoIndex(i);
 976     hisname.Form("hBkgPt%d",i);
 977     fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 978     hisname.Form("hCosPBkgPt%d",i);
 979     fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 980     hisname.Form("hDLenBkgPt%d",i);
 981     fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 982     hisname.Form("hSumd02BkgPt%d",i);
 983     fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 984     hisname.Form("hSigVertBkgPt%d",i);
 985     fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 986     hisname.Form("hPtMaxBkgPt%d",i);
 987     fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 988     hisname.Form("hDCABkgPt%d",i);
 989     fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 990     hisname.Form("hBkgPt%dTC",i);
 991     fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
 992
 993  }
 994     fWellIDProt3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt3sig));
 995     fWellIDProt2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt2sig));
 996     fWellIDProt3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt3sigSe));
 997     fWellIDProt2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt2sigSe));
 998     fWellIDProt2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt2sigRe));
 999     fWellIDProt3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDProt3sigRe));
1000     fWellIDKaon3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon3sig));
1001     fWellIDKaon2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon2sig));
1002     fWellIDKaon3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon3sigSe));
1003     fWellIDKaon2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon2sigSe));
1004     fWellIDKaon2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon2sigRe));
1005     fWellIDKaon3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fWellIDKaon3sigRe));
1006     fTPCSignal3Sigma=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3Sigma));
1007     fTPCSignal3SigmaReK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaReK));
1008     fTPCSignal3SigmaSedK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaSedK));
1009     fTPCSignal3SigmaRep=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaRep));
1010     fTPCSignal3SigmaSedp=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaSedp));
1011     fTPCSignal2Sigma=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2Sigma));
1012     fTPCSignal2SigmaReK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaReK));
1013     fTPCSignal2SigmaSedK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaSedK));
1014     fTPCSignal2SigmaSedp=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaSedp));
1015     fTPCSignal2SigmaRep=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaRep));
1016     fFakeProt3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt3sig));
1017     fFakeProt2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt2sig));
1018     fFakeProt3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt3sigSe));
1019     fFakeProt2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt2sigSe));
1020     fFakeProt2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt2sigRe));
1021     fFakeProt3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeProt3sigRe));
1022     fFakeKaon3sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon3sig));
1023     fFakeKaon2sig=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon2sig));
1024     fFakeKaon3sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon3sigSe));
1025     fFakeKaon2sigSe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon2sigSe));
1026     fFakeKaon2sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon2sigRe));
1027     fFakeKaon3sigRe=dynamic_cast<TH1F*>(fOutput->FindObject(fFakeKaon3sigRe));
1028
1029     fRealProt=dynamic_cast<TH1F*>(fOutput->FindObject(fRealProt));
1030     fRealKaon=dynamic_cast<TH1F*>(fOutput->FindObject(fRealKaon));
1031
1032   if(fFillNtuple){
1033     fNtupleLambdac = dynamic_cast<TNtuple*>(GetOutputData(3));
1034   }
1035
1036   TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
1037   //c1->cd();
1038   //TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
1039  fRealProt->Draw();
1040
1041  return;
1042 }
1043
1044 //________________________________________________________________________
1045 Int_t AliAnalysisTaskSELambdac::MatchToMCLambdac(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
1046
1047   Int_t lambdacLab[3]={0,0,0};
1048   Int_t pdgs[3]={0,0,0};
1049   for(Int_t i=0;i<3;i++){
1050    AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
1051    Int_t lab=daugh->GetLabel();
1052    if(lab<0) return 0;
1053    AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab);
1054    if(!part) continue;
1055    pdgs[i]=part->GetPdgCode();
1056    Int_t partPdgcode = TMath::Abs(part->GetPdgCode());
1057    if(partPdgcode==211 || partPdgcode==321 || partPdgcode==2212){
1058         Int_t motherLabel=part->GetMother();
1059              if(motherLabel<0) return 0;
1060    AliAODMCParticle *motherPart = (AliAODMCParticle*)arrayMC->At(motherLabel);
1061    if(!motherPart) continue;
1062         Int_t motherPdg = TMath::Abs(motherPart->GetPdgCode());
1063    if(motherPdg==4122) {
1064         if(GetLambdacDaugh(motherPart,arrayMC)){lambdacLab[i]=motherLabel;continue;}
1065               }
1066    if(motherPdg==313 || motherPdg==2224 || motherPdg==3124){
1067          Int_t GmotherLabel=motherPart->GetMother();
1068                if(GmotherLabel<0) return 0;
1069          AliAODMCParticle *GmotherPart = (AliAODMCParticle*)arrayMC->At(GmotherLabel);
1070          if(!GmotherPart) continue;
1071                Int_t GmotherPdg = TMath::Abs(GmotherPart->GetPdgCode());
1072          if(GmotherPdg==4122) {
1073                 if(GetLambdacDaugh(GmotherPart,arrayMC)) {lambdacLab[i]=GmotherLabel;continue;}
1074               }
1075         }
1076      }
1077   }
1078
1079  if(lambdacLab[0]==lambdacLab[1] && lambdacLab[1]==lambdacLab[2]) {return lambdacLab[0];}
1080  return 0;
1081
1082 }
1083 //------------------------
1084 Bool_t AliAnalysisTaskSELambdac::GetLambdacDaugh(AliAODMCParticle *part,TClonesArray *arrayMC){
1085
1086  Int_t numberOfLambdac=0;
1087  if(TMath::Abs(part->GetPdgCode())!=4122) return kFALSE;
1088  Int_t daugh_tmp[2];
1089  daugh_tmp[0]=part->GetDaughter(0);
1090  daugh_tmp[1]=part->GetDaughter(1);
1091  Int_t nDaugh = (Int_t)part->GetNDaughters();
1092  if(nDaugh<2) return kFALSE;
1093  if(nDaugh>3) return kFALSE;
1094  AliAODMCParticle* pdaugh1 = (AliAODMCParticle*)arrayMC->At(part->GetDaughter(0));
1095    if(!pdaugh1) {return kFALSE;}
1096  Int_t number1 = TMath::Abs(pdaugh1->GetPdgCode());
1097  AliAODMCParticle* pdaugh2 = (AliAODMCParticle*)arrayMC->At(part->GetDaughter(1));
1098    if(!pdaugh2) {return kFALSE;}
1099   Int_t number2 = TMath::Abs(pdaugh2->GetPdgCode());
1100
1101   if(nDaugh==3){
1102    Int_t thirdDaugh=part->GetDaughter(1)-1;
1103             AliAODMCParticle* pdaugh3 = (AliAODMCParticle*)arrayMC->At(thirdDaugh);
1104   Int_t number3 = TMath::Abs(pdaugh3->GetPdgCode());
1105    if((number1==321 && number2==211 && number3==2212) || (number1==211 && number2==321 && number3==2212) || (number1==211 && number2==2212 && number3==321) || (number1==321 && number2==2212 && number3==211) || (number1==2212 && number2==321 && number3==211) || (number1==2212 && number2==211 && number3==321)) {
1106    numberOfLambdac++;
1107    }
1108   }
1109
1110  if(nDaugh==2){
1111
1112          //Lambda in Lambda pi
1113 //  if((number1==211 && number2==3122)|| (number1==3122 && number2==211))return kFALSE;
1114
1115   //Lambda resonant
1116
1117   //Lambda -> p K*0
1118   //
1119   Int_t nfiglieK=0;
1120
1121    if((number1==2212 && number2==313)){
1122      nfiglieK=pdaugh2->GetNDaughters();
1123      if(nfiglieK!=2) return kFALSE;
1124      AliAODMCParticle* pdaughK1 = (AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
1125      AliAODMCParticle* pdaughK2 = (AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
1126      if(!pdaughK1) return kFALSE;
1127      if(!pdaughK2) return kFALSE;
1128      if((TMath::Abs(pdaughK1->GetPdgCode())==211 && TMath::Abs(pdaughK2->GetPdgCode())==321) || (TMath::Abs(pdaughK1->GetPdgCode())==321 && TMath::Abs(pdaughK2->GetPdgCode())==211)) numberOfLambdac++;
1129     }
1130
1131    if((number1==313 && number2==2212)){
1132     nfiglieK=pdaugh1->GetNDaughters();
1133     if(nfiglieK!=2) return kFALSE;
1134     AliAODMCParticle* pdaughK1 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
1135     AliAODMCParticle* pdaughK2 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
1136     if(!pdaughK1) return kFALSE;
1137     if(!pdaughK2) return kFALSE;
1138      if((TMath::Abs(pdaughK1->GetPdgCode())==211 && TMath::Abs(pdaughK2->GetPdgCode())==321) || (TMath::Abs(pdaughK1->GetPdgCode())==321 && TMath::Abs(pdaughK2->GetPdgCode())==211)) numberOfLambdac++;
1139    }
1140
1141    //Lambda -> Delta++ k
1142    Int_t nfiglieDelta=0;
1143    if(number1==321 && number2==2224){
1144     nfiglieDelta=pdaugh2->GetNDaughters();
1145     if(nfiglieDelta!=2) return kFALSE;
1146     AliAODMCParticle *pdaughD1=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
1147     AliAODMCParticle *pdaughD2=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
1148     if(!pdaughD1) return kFALSE;
1149     if(!pdaughD2) return kFALSE;
1150     if((TMath::Abs(pdaughD1->GetPdgCode())==211 && TMath::Abs(pdaughD2->GetPdgCode())==2212) || (TMath::Abs(pdaughD1->GetPdgCode())==2212 && TMath::Abs(pdaughD2->GetPdgCode())==211)) numberOfLambdac++;
1151    }
1152    if(number1==2224 && number2==321){
1153     nfiglieDelta=pdaugh1->GetNDaughters();
1154     if(nfiglieDelta!=2) return kFALSE;
1155     AliAODMCParticle* pdaughD1 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
1156     AliAODMCParticle* pdaughD2 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
1157     if(!pdaughD1) return kFALSE;
1158     if(!pdaughD2) return kFALSE;
1159     if((TMath::Abs(pdaughD1->GetPdgCode())==211 && TMath::Abs(pdaughD2->GetPdgCode())==2212) || (TMath::Abs(pdaughD1->GetPdgCode())==2212 && TMath::Abs(pdaughD2->GetPdgCode())==211)) numberOfLambdac++;
1160    }
1161
1162
1163    //Lambdac -> Lambda(1520) pi
1164    Int_t nfiglieLa=0;
1165    if(number1==3124 && number2==211){
1166     nfiglieLa=pdaugh1->GetNDaughters();
1167     if(nfiglieLa!=2) return kFALSE;
1168     AliAODMCParticle *pdaughL1=(AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
1169     AliAODMCParticle *pdaughL2=(AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
1170     if(!pdaughL1) return kFALSE;
1171     if(!pdaughL2) return kFALSE;
1172     if((TMath::Abs(pdaughL1->GetPdgCode())==321 && TMath::Abs(pdaughL2->GetPdgCode())==2212) || (TMath::Abs(pdaughL1->GetPdgCode())==2212 && TMath::Abs(pdaughL2->GetPdgCode())==321)) numberOfLambdac++;
1173    }
1174    if(number1==211 && number2==3124){
1175     nfiglieLa=pdaugh2->GetNDaughters();
1176     if(nfiglieLa!=2) return kFALSE;
1177     AliAODMCParticle *pdaughL1=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
1178     AliAODMCParticle *pdaughL2=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
1179     if(!pdaughL1) return kFALSE;
1180     if(!pdaughL2) return kFALSE;
1181     if((TMath::Abs(pdaughL1->GetPdgCode())==321 && TMath::Abs(pdaughL2->GetPdgCode())==2212) || (TMath::Abs(pdaughL1->GetPdgCode())==2212 && TMath::Abs(pdaughL2->GetPdgCode())==321)) numberOfLambdac++;
1182
1183    }
1184  }
1185
1186  if(numberOfLambdac>0) {return kTRUE;}
1187          return kFALSE;
1188 }
1189 //-----------------------------
1190 Bool_t AliAnalysisTaskSELambdac::IspKpiMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC,Int_t *pdgs){
1191
1192    Int_t lab[3]={0,0,0};//,pdgs[3]={0,0,0};
1193    for(Int_t i=0;i<3;i++){
1194     AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
1195     lab[i]=daugh->GetLabel();
1196     if(lab[i]<0) return kFALSE;
1197     AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab[i]);
1198     if(!part) return kFALSE;
1199     pdgs[i]=TMath::Abs(part->GetPdgCode());
1200    }
1201
1202    if(pdgs[0]==2212 && pdgs[1]==321 && pdgs[2]==211) return kTRUE;
1203
1204    return kFALSE;
1205 }
1206 //-----------------------------
1207 Bool_t AliAnalysisTaskSELambdac::IspiKpMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
1208
1209    Int_t lab[3]={0,0,0},pdgs[3]={0,0,0};
1210    for(Int_t i=0;i<3;i++){
1211     AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
1212     lab[i]=daugh->GetLabel();
1213     if(lab[i]<0) return kFALSE;
1214     AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab[i]);
1215     if(!part) return kFALSE;
1216     pdgs[i]=TMath::Abs(part->GetPdgCode());
1217    }
1218
1219    if(pdgs[2]==2212 && pdgs[1]==321 && pdgs[0]==211) {return kTRUE;}
1220
1221    return kFALSE;
1222 }
1223 //--------------------------------------
1224 Bool_t AliAnalysisTaskSELambdac::IspiKpReal(AliAODRecoDecayHF3Prong *d,Int_t *pdgs){
1225
1226   AliAODPidHF* pidObj=new AliAODPidHF();
1227   Bool_t type[2]={kFALSE,kFALSE};
1228   for(Int_t i=0;i<3;i++){
1229 //   Bool_t pid[3]={kFALSE,kFALSE,kFALSE};
1230   AliAODTrack *track=(AliAODTrack*)d->GetDaughter(i);
1231    //test asymmetric TPC PID
1232      Double_t plim[2]={0.6,2.};
1233      Double_t sigmas[5]={2.,1.,2.,3.,0.};
1234      pidObj->SetPLimit(plim);
1235      pidObj->SetSigma(sigmas);
1236      if(i==1) type[i]=pidObj->MatchTPCTOF(track,2,3,kFALSE);
1237      if(i==2) type[i]=pidObj->MatchTPCTOF(track,2,4,kFALSE);
1238      //pidinfo=pidObj->MatchTPCTOF(track,3,3,kFALSE);
1239
1240   }
1241
1242   if(type[0] && type[1]) {delete pidObj;return kTRUE;}
1243   delete pidObj;
1244  return kFALSE;
1245 }
1246 //--------------------------------------
1247 Bool_t AliAnalysisTaskSELambdac::IspKpiReal(AliAODRecoDecayHF3Prong *d){
1248
1249   AliAODPidHF* pidObj=new AliAODPidHF();
1250   Bool_t type[2]={kFALSE,kFALSE};
1251   for(Int_t i=0;i<2;i++){
1252    //Bool_t pid[3]={kFALSE,kFALSE,kFALSE};
1253    AliAODTrack *track=(AliAODTrack*)d->GetDaughter(i);
1254      Double_t plim[2]={0.6,2.};
1255      Double_t sigmas[5]={2.,1.,2.,3.,0.};
1256      pidObj->SetPLimit(plim);
1257      pidObj->SetSigma(sigmas);
1258      if(i==1) type[i]=pidObj->MatchTPCTOF(track,2,3,kFALSE);
1259      if(i==0) type[i]=pidObj->MatchTPCTOF(track,2,4,kFALSE);
1260   }
1261   if(type[0] && type[1]) {delete pidObj;return kTRUE;}
1262   delete pidObj;
1263  return kFALSE;
1264 }
1265 //------------------------------------
1266 Bool_t AliAnalysisTaskSELambdac::VertexingKF(AliAODRecoDecayHF3Prong *d,Int_t *pdgs,Double_t field){
1267
1268    Int_t iprongs[3]={0,1,2};
1269    Double_t mass[2]={0.,0.};
1270  //topological constr
1271    AliKFParticle *Lambdac=d->ApplyVertexingKF(iprongs,3,pdgs,kTRUE,field,mass);
1272    if(!Lambdac) return kFALSE;
1273   Double_t probTot=TMath::Prob(Lambdac->GetChi2(),Lambdac->GetNDF());
1274   if(probTot<fCutsKF[0]) return kFALSE;
1275
1276   if(d->PtProng(0)<1.5 && d->Getd0Prong(0)>0.02) return kFALSE;
1277     if(d->PtProng(1)<1.5 && d->Getd0Prong(1)>0.02) return kFALSE;
1278       if(d->PtProng(2)<1.5 && d->Getd0Prong(2)>0.02) return kFALSE;
1279  //mass constr for K*
1280    Int_t ipRes[2];
1281    Int_t pdgres[2];
1282    mass[0]=0.8961;mass[1]=0.03;
1283    if(TMath::Abs(pdgs[0])==211){
1284     ipRes[0]=0;ipRes[1]=1;
1285     pdgres[0]=pdgs[0];pdgres[1]=321;
1286    }
1287    if(TMath::Abs(pdgs[2])==211){
1288     ipRes[0]=2;ipRes[1]=1;
1289     pdgres[0]=pdgs[2];pdgres[1]=321;
1290    }
1291    AliKFParticle *KappaStar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1292
1293   Double_t probKstar=TMath::Prob(KappaStar->GetChi2(),KappaStar->GetNDF());
1294   if(probKstar>fCutsKF[1]) {
1295     AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
1296     AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
1297     AliKFParticle prong1(*esdProng1,pdgres[0]);
1298     AliKFParticle prong2(*esdProng2,pdgres[1]);
1299     if(KappaStar->GetPt()<fCutsKF[2] && prong1.GetAngle(prong2)>fCutsKF[3]) return kFALSE;
1300   }
1301  //mass constr for Lambda
1302    mass[0]=1.520;mass[1]=0.005;
1303    if(TMath::Abs(pdgs[0])==2212){
1304     ipRes[0]=0;ipRes[1]=1;
1305     pdgres[0]=pdgs[0];pdgres[1]=pdgs[1];
1306    }
1307    if(TMath::Abs(pdgs[2])==2212){
1308     ipRes[0]=2;ipRes[1]=1;
1309     pdgres[0]=pdgs[2];pdgres[1]=pdgs[1];
1310    }
1311    AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1312   Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
1313   if(probLa>fCutsKF[4]) {
1314     AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
1315     AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
1316     AliKFParticle prong1(*esdProng1,pdgres[0]);
1317     AliKFParticle prong2(*esdProng2,pdgres[1]);
1318     if(Lambda1520->GetPt()<fCutsKF[5] && prong1.GetAngle(prong2)>fCutsKF[6]) return kFALSE;
1319   }
1320  //mass constr for Delta
1321    mass[0]=1.2;mass[1]=0.15;
1322    ipRes[0]=0;ipRes[1]=2;
1323    pdgres[0]=pdgs[0];pdgres[2]=pdgs[2];
1324    AliKFParticle *Delta=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1325   Double_t probDelta=TMath::Prob(Delta->GetChi2(),Delta->GetNDF());
1326   if(probDelta>fCutsKF[7]) {
1327     AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
1328     AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
1329     AliKFParticle prong1(*esdProng1,pdgres[0]);
1330     AliKFParticle prong2(*esdProng2,pdgres[1]);
1331     if(Delta->GetPt()<fCutsKF[8] && prong1.GetAngle(prong2)>fCutsKF[9]) return kFALSE;
1332   }
1333  return kTRUE;
1334 }
1335 //-------------------------------------
1336 Bool_t AliAnalysisTaskSELambdac::IspiKpResonant(AliAODRecoDecayHF3Prong *d,Double_t field){
1337
1338  //if lambda* -> pk
1339         Double_t mass[2]={1.520,0.005};
1340         Int_t ipRes[2]={1,2};
1341         Int_t pdgres[2]={321,2212};
1342         AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1343         Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
1344         if(probLa>0.1) return kTRUE;
1345  //K* -> kpi
1346 //        mass[0]=0.8961;mass[1]=0.03;
1347 //        ipRes[0]=0;ipRes[1]=1;
1348 //        pdgres[0]=211;pdgres[1]=321;
1349 //        AliKFParticle *Kstar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1350 //      Double_t probKa=TMath::Prob(Kstar->GetChi2(),Kstar->GetNDF());
1351 //      if(probKa>0.1) return kTRUE;
1352
1353  return kFALSE;
1354
1355 }
1356 //-------------------------------------
1357 Bool_t AliAnalysisTaskSELambdac::IspKpiResonant(AliAODRecoDecayHF3Prong *d,Double_t field){
1358
1359  //if lambda* -> pk
1360         Double_t mass[2]={1.520,0.005};
1361         Int_t ipRes[2]={0,1};
1362         Int_t pdgres[2]={2212,321};
1363         AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1364         Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
1365         if(probLa>0.1) return kTRUE;
1366  //K* -> kpi
1367 //        mass[0]=0.8961;mass[1]=0.03;
1368 //        ipRes[0]=1;ipRes[1]=2;
1369 //        pdgres[1]=211;pdgres[0]=321;
1370 //        AliKFParticle *Kstar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
1371 //      Double_t probKa=TMath::Prob(Kstar->GetChi2(),Kstar->GetNDF());
1372 //      if(probKa>0.1) return kTRUE;
1373
1374  return kFALSE;
1375
1376 }