New task for Lc->pKpi (Rossella)
authordainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 24 May 2010 13:08:24 +0000 (13:08 +0000)
committerdainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 24 May 2010 13:08:24 +0000 (13:08 +0000)
PWG3/PWG3vertexingHFLinkDef.h
PWG3/libPWG3vertexingHF.pkg
PWG3/vertexingHF/AddTaskLambdac.C [new file with mode: 0644]
PWG3/vertexingHF/AliAnalysisTaskSELambdac.cxx [new file with mode: 0644]
PWG3/vertexingHF/AliAnalysisTaskSELambdac.h [new file with mode: 0644]

index 15bde53..f2c0091 100644 (file)
@@ -27,6 +27,7 @@
 #pragma link C++ class AliAnalysisTaskSECompareHF+;
 #pragma link C++ class AliAnalysisTaskSEDplus+;
 #pragma link C++ class AliAnalysisTaskSEDs+;
+#pragma link C++ class AliAnalysisTaskSELambdac+;
 #pragma link C++ class AliAnalysisTaskSED0Mass+;
 #pragma link C++ class AliAnalysisTaskSECharmFraction+;
 #pragma link C++ class AliCFHeavyFlavourTaskMultiVarMultiStep+;
index e8b8999..17a14e8 100644 (file)
@@ -21,6 +21,7 @@ SRCS:=   vertexingHF/AliAODRecoDecayHF.cxx \
        vertexingHF/AliAnalysisTaskSECompareHF.cxx \
        vertexingHF/AliAnalysisTaskSEDplus.cxx \
        vertexingHF/AliAnalysisTaskSEDs.cxx \
+       vertexingHF/AliAnalysisTaskSELambdac.cxx \
        vertexingHF/AliAnalysisTaskSED0Mass.cxx \
        vertexingHF/AliAnalysisTaskSECharmFraction.cxx \
        vertexingHF/AliCFHeavyFlavourTaskMultiVarMultiStep.cxx \
diff --git a/PWG3/vertexingHF/AddTaskLambdac.C b/PWG3/vertexingHF/AddTaskLambdac.C
new file mode 100644 (file)
index 0000000..30150d3
--- /dev/null
@@ -0,0 +1,150 @@
+AliAnalysisTaskSELambdac *AddTaskLambdac(const char *name,Bool_t storeNtuple,Bool_t readMC,Bool_t MCPid,Bool_t realPid,Bool_t resPid,Bool_t useKF)
+{
+  //==============================================================================                                                      
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+    ::Error("AddTaskLambdac", "No analysis manager to connect to.");
+    return NULL;
+  }
+
+  //pt bins
+  const Int_t nptbins=4;
+    Float_t* ptbins;
+    ptbins=new Float_t[nptbins+1];
+    ptbins[0]=0.;
+    ptbins[1]=2.;
+    ptbins[2]=3.;
+    ptbins[3]=5.;
+    ptbins[4]=99999.;
+    const Int_t nvars=12;
+
+  Float_t** prodcutsval;
+  prodcutsval=new Float_t*[nvars];
+  for(Int_t ic=0;ic<nvars;ic++){prodcutsval[ic]=new Float_t[nptbins];}
+  for(Int_t ipt=0;ipt<nptbins;ipt++){
+   prodcutsval[0][ipt]=0.2;
+   prodcutsval[1][ipt]=0.4;
+   prodcutsval[2][ipt]=0.4;
+   prodcutsval[3][ipt]=0.;
+   prodcutsval[4][ipt]=0.;
+   prodcutsval[5][ipt]=0.01;
+   prodcutsval[6][ipt]=0.06;
+   prodcutsval[7][ipt]=0.02;
+   prodcutsval[8][ipt]=0.;
+   prodcutsval[9][ipt]=0.85;
+   prodcutsval[10][ipt]=0.;
+   prodcutsval[11][ipt]=0.1;
+  }
+
+  Float_t** anacutsval;
+  anacutsval=new Float_t*[nvars];
+  for(Int_t ic=0;ic<nvars;ic++){anacutsval[ic]=new Float_t[nptbins];}
+  Int_t ic=0;
+  for(Int_t ipt=0;ipt<nptbins;ipt++){
+       anacutsval[ic][ipt]=0.2;
+          }
+  Int_t ic=1;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.8;
+            }
+  Int_t ic=2;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.4;
+            }
+  Int_t ic=3;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=4;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=5;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=6;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.1;
+            }
+  Int_t ic=7;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=8;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=9;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=10;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.;
+            }
+  Int_t ic=11;
+     for(Int_t ipt=0;ipt<nptbins;ipt++){
+          anacutsval[ic][ipt]=0.1;
+            }
+
+  anacutsval[1][0]=0.6;
+
+  anacutsval[2][2]=0.5;
+  anacutsval[2][3]=0.5;
+
+  anacutsval[9][2]=0.9;
+  anacutsval[9][3]=0.9;
+
+
+  AliRDHFCutsLctopKpi *prodcuts = new AliRDHFCutsLctopKpi();
+    prodcuts->SetPtBins(nptbins+1,ptbins);
+      prodcuts->SetCuts(nvars,nptbins,prodcutsval);
+  //Analysis cuts
+  AliRDHFCutsLctopKpi *analysiscuts = new AliRDHFCutsLctopKpi();
+  ((AliRDHFCuts*)analysiscuts)->SetPtBins(nptbins+1,ptbins);
+  ((AliRDHFCuts*)analysiscuts)->SetCuts(nvars,nptbins,anacutsval);
+  // Aanalysis task                                                                                                                     
+  AliAnalysisTaskSELambdac *lambdacTask = new AliAnalysisTaskSELambdac("LambdacAnalysis",storeNtuple,analysiscuts,prodcuts);
+  lambdacTask->SetReadMC(readMC);
+  if(MCPid) lambdacTask->SetMCPid();
+  if(resPid) lambdacTask->SetResonantPid();
+  if(realPid) lambdacTask->SetRealPid();
+  lambdacTask->SetDebugLevel(0);
+  if(useKF) {
+    lambdacTask->SetUseKF();
+     Float_t cuts[10]={0.1,0.1,1.5,0.5,0.1,1.5,0.5,0.1,1.5,0.5};
+      lambdacTask->SetCutsKF(cuts);
+       }
+  mgr->AddTask(lambdacTask);
+
+  //                                                                                                                                    
+  // Create containers for input/output                                                                                                 
+  TString outputfile = AliAnalysisManager::GetCommonFileName();
+  outputfile += ":PWG3_D2H_InvMassLambdac";
+  AliAnalysisDataContainer *cinputLambdac = mgr->CreateContainer("cinputLambdac",TChain::Class(),
+                                                          AliAnalysisManager::kInputContainer);
+  AliAnalysisDataContainer *coutputLambdacCuts = mgr->CreateContainer("coutputLambdacCuts",TList::Class(),
+                                                           AliAnalysisManager::kOutputContainer,
+                                                               outputfile.Data());
+  AliAnalysisDataContainer *coutputLambdac = mgr->CreateContainer("coutputLambdac",TList::Class(),
+                                                           AliAnalysisManager::kOutputContainer,
+                                                               outputfile.Data());
+  if(storeNtuple){
+    AliAnalysisDataContainer *coutputLambdac2 = mgr->CreateContainer("coutputLambdac2",TNtuple::Class(),
+                                                           AliAnalysisManager::kOutputContainer,
+                                                                "InvMassLambdac_nt1.root");
+
+    coutputLambdac2->SetSpecialOutput();
+  }
+
+  mgr->ConnectInput(lambdacTask,0,mgr->GetCommonInputContainer());
+
+  mgr->ConnectOutput(lambdacTask,1,coutputLambdac);
+  mgr->ConnectOutput(lambdacTask,2,coutputLambdacCuts);
+  
+  if(storeNtuple){
+    mgr->ConnectOutput(lambdacTask,3,coutputLambdac2);
+  }
+  return lambdacTask;
+}
diff --git a/PWG3/vertexingHF/AliAnalysisTaskSELambdac.cxx b/PWG3/vertexingHF/AliAnalysisTaskSELambdac.cxx
new file mode 100644 (file)
index 0000000..467db92
--- /dev/null
@@ -0,0 +1,1266 @@
+/**************************************************************************
+ * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/////////////////////////////////////////////////////////////
+//
+// AliAnalysisTaskSE for the extraction of signal(e.g D+) of heavy flavor
+// decay candidates with the MC truth.
+// Authors: Renu Bala, bala@to.infn.it
+// F. Prino, prino@to.infn.it
+// G. Ortona, ortona@to.infn.it
+/////////////////////////////////////////////////////////////
+
+#include <TClonesArray.h>
+#include <TNtuple.h>
+#include <TCanvas.h>
+#include <TList.h>
+#include <TString.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TDatabasePDG.h>
+
+#include "AliAnalysisManager.h"
+#include "AliAODHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODVertex.h"
+#include "AliAODTrack.h"
+#include "AliAODMCHeader.h"
+#include "AliAODMCParticle.h"
+#include "AliAODRecoDecayHF3Prong.h"
+#include "AliAnalysisVertexingHF.h"
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisTaskSELambdac.h"
+#include "AliKFParticle.h"
+#include "AliAODPidHF.h"
+#include "AliRDHFCutsLctopKpi.h"
+#include "AliRDHFCuts.h"
+
+ClassImp(AliAnalysisTaskSELambdac)
+
+
+//________________________________________________________________________
+AliAnalysisTaskSELambdac::AliAnalysisTaskSELambdac():
+AliAnalysisTaskSE(),
+fOutput(0), 
+fHistNEvents(0),
+fTPCSignal3Sigma(0),
+fTPCSignal3SigmaReK(0),
+fTPCSignal3SigmaRep(0),
+fTPCSignal2Sigma(0),
+fTPCSignal2SigmaReK(0),
+fTPCSignal2SigmaRep(0),
+fNtupleLambdac(0),
+fUpmasslimit(2.486),
+fLowmasslimit(2.086),
+fNPtBins(0),
+fRDCutsAnalysis(0),
+fRDCutsProduction(0),
+fListCuts(0),
+fFillNtuple(kFALSE),
+fReadMC(kFALSE),
+fMCPid(kFALSE),
+fRealPid(kFALSE),
+fResPid(kTRUE),
+fUseKF(kFALSE),
+fVHF(0)
+{
+   // Default constructor
+}
+
+//________________________________________________________________________
+AliAnalysisTaskSELambdac::AliAnalysisTaskSELambdac(const char *name,Bool_t fillNtuple,AliRDHFCutsLctopKpi *lccutsana,AliRDHFCutsLctopKpi *lccutsprod):
+AliAnalysisTaskSE(name),
+fOutput(0),
+fHistNEvents(0),
+fTPCSignal3Sigma(0),
+fTPCSignal3SigmaReK(0),
+fTPCSignal3SigmaRep(0),
+fTPCSignal2Sigma(0),
+fTPCSignal2SigmaReK(0),
+fTPCSignal2SigmaRep(0),
+fNtupleLambdac(0),
+fUpmasslimit(2.486),
+fLowmasslimit(2.086),
+fNPtBins(0),
+fRDCutsAnalysis(lccutsana),
+fRDCutsProduction(lccutsprod),
+fListCuts(0),
+fFillNtuple(fillNtuple),
+fReadMC(kFALSE),
+fMCPid(kFALSE),
+fRealPid(kTRUE),
+fResPid(kFALSE),
+fUseKF(kFALSE),
+fVHF(0)
+{
+  //Double_t ptlim[5]={0.,2.,3.,5,9999999.};
+  //SetPtBinLimit(5, ptlim);
+   SetPtBinLimit(fRDCutsAnalysis->GetNPtBins()+1,fRDCutsAnalysis->GetPtBinLimits());
+  // Default constructor
+   // Output slot #1 writes into a TList container
+  DefineOutput(1,TList::Class());  //My private output
+  DefineOutput(2,TList::Class());
+
+  if(fFillNtuple){
+    // Output slot #2 writes into a TNtuple container
+    DefineOutput(3,TNtuple::Class());  //My private output
+  }
+}
+
+//________________________________________________________________________
+AliAnalysisTaskSELambdac::~AliAnalysisTaskSELambdac()
+{
+  // Destructor
+  if (fOutput) {
+    delete fOutput;
+    fOutput = 0;
+  }
+  if (fVHF) {
+    delete fVHF;
+    fVHF = 0;
+  }
+  
+ if(fRDCutsAnalysis){
+    delete fRDCutsAnalysis;
+    fRDCutsAnalysis = 0;
+ }
+ if(fRDCutsProduction){
+    delete fRDCutsProduction;
+    fRDCutsProduction = 0;
+ }
+
+ if (fListCuts) {
+   delete fListCuts;
+   fListCuts = 0;
+ }
+}  
+//_________________________________________________________________
+void  AliAnalysisTaskSELambdac::SetMassLimits(Float_t range){
+  fUpmasslimit = 2.286+range;
+  fLowmasslimit = 2.286-range;
+}
+//_________________________________________________________________
+void  AliAnalysisTaskSELambdac::SetMassLimits(Float_t lowlimit, Float_t uplimit){
+  if(uplimit>lowlimit)
+    {
+      fUpmasslimit = lowlimit;
+      fLowmasslimit = uplimit;
+    }
+}
+
+
+//________________________________________________________________________
+void AliAnalysisTaskSELambdac::SetPtBinLimit(Int_t n, Float_t* lim){
+  // define pt bins for analysis
+  if(n>kMaxPtBins){
+    printf("Max. number of Pt bins = %d\n",kMaxPtBins);
+    fNPtBins=kMaxPtBins;
+    fArrayBinLimits[0]=0.;
+    fArrayBinLimits[1]=2.;
+    fArrayBinLimits[2]=3.;
+    fArrayBinLimits[3]=5.;
+    for(Int_t i=4; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
+  }else{
+    fNPtBins=n-1;
+    fArrayBinLimits[0]=lim[0];
+    for(Int_t i=1; i<fNPtBins+1; i++) 
+      if(lim[i]>fArrayBinLimits[i-1]){
+       fArrayBinLimits[i]=lim[i];
+      }
+      else {
+       fArrayBinLimits[i]=fArrayBinLimits[i-1];
+      }
+    for(Int_t i=fNPtBins; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
+  }
+  if(fDebug > 1){
+    printf("Number of Pt bins = %d\n",fNPtBins);
+    for(Int_t i=0; i<fNPtBins+1; i++) printf(" Bin%d = %8.2f-%8.2f\n",i,fArrayBinLimits[i],fArrayBinLimits[i+1]);    
+  }
+}
+//_________________________________________________________________
+Double_t  AliAnalysisTaskSELambdac::GetPtBinLimit(Int_t ibin){
+  if(ibin>fNPtBins)return -1;
+  return fArrayBinLimits[ibin];
+} 
+
+//_________________________________________________________________
+void AliAnalysisTaskSELambdac::Init()
+{
+  // Initialization
+
+  if(fDebug > 1) printf("AnalysisTaskSELambdac::Init() \n");
+
+  fListCuts=new TList();
+  AliRDHFCutsLctopKpi *production = new AliRDHFCutsLctopKpi();
+    production=fRDCutsProduction;
+
+  AliRDHFCutsLctopKpi *analysis = new AliRDHFCutsLctopKpi();
+  analysis=fRDCutsAnalysis;
+
+   fListCuts->Add(analysis);
+   fListCuts->Add(production);
+   PostData(2,fListCuts);
+  return;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSELambdac::UserCreateOutputObjects()
+{
+  // Create the output container
+  //
+  if(fDebug > 1) printf("AnalysisTaskSELambdac::UserCreateOutputObjects() \n");
+
+  // Several histograms are more conveniently managed in a TList
+  fOutput = new TList();
+  fOutput->SetOwner();
+  fOutput->SetName("OutputHistos");
+
+  TString hisname;
+  Int_t index=0;
+  Int_t indexLS=0;
+  for(Int_t i=0;i<fNPtBins;i++){
+
+    index=GetHistoIndex(i);
+    indexLS=GetLSHistoIndex(i);
+
+    hisname.Form("hMassPt%d",i);
+    fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHist[index]->Sumw2();
+    hisname.Form("hCosPAllPt%d",i);
+    fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
+    fCosPHist[index]->Sumw2();
+    hisname.Form("hDLenAllPt%d",i);
+    fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
+    fDLenHist[index]->Sumw2();
+    hisname.Form("hSumd02AllPt%d",i);
+    fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
+    fSumd02Hist[index]->Sumw2();
+    hisname.Form("hSigVertAllPt%d",i);
+    fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fSigVertHist[index]->Sumw2();
+    hisname.Form("hPtMaxAllPt%d",i);
+    fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
+    fPtMaxHist[index]->Sumw2();
+
+    hisname.Form("hDCAAllPt%d",i);
+    fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fDCAHist[index]->Sumw2();
+
+
+
+    hisname.Form("hMassPt%dTC",i);
+    fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistTC[index]->Sumw2();
+
+
+
+
+    
+    hisname.Form("hCosPAllPt%dLS",i);
+    fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
+    fCosPHistLS[index]->Sumw2();
+    hisname.Form("hDLenAllPt%dLS",i);
+    fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
+    fDLenHistLS[index]->Sumw2();
+    hisname.Form("hSumd02AllPt%dLS",i);
+    fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
+    fSumd02HistLS[index]->Sumw2();
+    hisname.Form("hSigVertAllPt%dLS",i);
+    fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fSigVertHistLS[index]->Sumw2();
+    hisname.Form("hPtMaxAllPt%dLS",i);
+    fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
+    fPtMaxHistLS[index]->Sumw2();
+    
+    hisname.Form("hDCAAllPt%dLS",i);
+    fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fDCAHistLS[index]->Sumw2();
+    
+    hisname.Form("hLSPt%dLC",i);
+    fMassHistLS[indexLS] = new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLS[indexLS]->Sumw2();
+    
+    hisname.Form("hLSPt%dTC",i);
+    fMassHistLSTC[indexLS] = new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLSTC[indexLS]->Sumw2();
+
+
+    
+    index=GetSignalHistoIndex(i);    
+    indexLS++;
+    hisname.Form("hSigPt%d",i);
+    fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHist[index]->Sumw2();
+    hisname.Form("hCosPSigPt%d",i);
+    fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
+    fCosPHist[index]->Sumw2();
+    hisname.Form("hDLenSigPt%d",i);
+    fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
+    fDLenHist[index]->Sumw2();
+    hisname.Form("hSumd02SigPt%d",i);
+    fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
+    fSumd02Hist[index]->Sumw2();
+    hisname.Form("hSigVertSigPt%d",i);
+    fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fSigVertHist[index]->Sumw2();
+    hisname.Form("hPtMaxSigPt%d",i);
+    fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
+    fPtMaxHist[index]->Sumw2();    
+
+    hisname.Form("hDCASigPt%d",i);
+    fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fDCAHist[index]->Sumw2();    
+
+
+    hisname.Form("hSigPt%dTC",i);
+    fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistTC[index]->Sumw2();
+
+    hisname.Form("hLSPt%dLCnw",i);
+    fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLS[indexLS]->Sumw2();
+    hisname.Form("hLSPt%dTCnw",i);
+    fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLSTC[indexLS]->Sumw2();
+
+
+    
+    hisname.Form("hCosPSigPt%dLS",i);
+    fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
+    fCosPHistLS[index]->Sumw2();
+    hisname.Form("hDLenSigPt%dLS",i);
+    fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
+    fDLenHistLS[index]->Sumw2();
+    hisname.Form("hSumd02SigPt%dLS",i);
+    fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
+    fSumd02HistLS[index]->Sumw2();
+    hisname.Form("hSigVertSigPt%dLS",i);
+    fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fSigVertHistLS[index]->Sumw2();
+    hisname.Form("hPtMaxSigPt%dLS",i);
+    fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
+    fPtMaxHistLS[index]->Sumw2();
+
+    hisname.Form("hDCASigPt%dLS",i);
+    fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fDCAHistLS[index]->Sumw2();
+    
+
+
+    index=GetBackgroundHistoIndex(i); 
+    indexLS++;
+    hisname.Form("hBkgPt%d",i);
+    fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHist[index]->Sumw2();
+    hisname.Form("hCosPBkgPt%d",i);
+    fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
+    fCosPHist[index]->Sumw2();
+    hisname.Form("hDLenBkgPt%d",i);
+    fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
+    fDLenHist[index]->Sumw2();
+    hisname.Form("hSumd02BkgPt%d",i);
+    fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
+    fSumd02Hist[index]->Sumw2();
+    hisname.Form("hSigVertBkgPt%d",i);
+    fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fSigVertHist[index]->Sumw2();
+    hisname.Form("hPtMaxBkgPt%d",i);
+    fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
+    fPtMaxHist[index]->Sumw2();
+
+    hisname.Form("hDCABkgPt%d",i);
+    fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fDCAHist[index]->Sumw2();
+
+
+    hisname.Form("hBkgPt%dTC",i);
+    fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistTC[index]->Sumw2();
+
+    hisname.Form("hLSPt%dLCntrip",i);
+    fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLS[indexLS]->Sumw2();
+    hisname.Form("hLSPt%dTCntrip",i);
+    fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLSTC[indexLS]->Sumw2();
+
+    
+    hisname.Form("hCosPBkgPt%dLS",i);
+    fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
+    fCosPHistLS[index]->Sumw2();
+    hisname.Form("hDLenBkgPt%dLS",i);
+    fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
+    fDLenHistLS[index]->Sumw2();
+    hisname.Form("hSumd02BkgPt%dLS",i);
+    fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
+    fSumd02HistLS[index]->Sumw2();
+    hisname.Form("hSigVertBkgPt%dLS",i);
+    fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fSigVertHistLS[index]->Sumw2();
+    hisname.Form("hPtMaxBkgPt%dLS",i);
+    fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,5.);
+    fPtMaxHistLS[index]->Sumw2();
+    hisname.Form("hDCABkgPt%dLS",i);
+    fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
+    fDCAHistLS[index]->Sumw2();
+    
+
+    indexLS++;
+    hisname.Form("hLSPt%dLCntripsinglecut",i);
+    fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLS[indexLS]->Sumw2();
+    hisname.Form("hLSPt%dTCntripsinglecut",i);
+    fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLSTC[indexLS]->Sumw2();
+
+    indexLS++;
+    hisname.Form("hLSPt%dLCspc",i);
+    fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLS[indexLS]->Sumw2();
+    hisname.Form("hLSPt%dTCspc",i);
+    fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),100,fLowmasslimit,fUpmasslimit);
+    fMassHistLSTC[indexLS]->Sumw2();
+  }
+
+  fTPCSignal3Sigma=new TH2F("fTPCSignal3Sigma","fTPCSignal3Sigma",100,0, 10, 100, 0, 100);
+  fTPCSignal3SigmaReK=new TH2F("fTPCSignal3SigmaReK","fTPCSignal3SigmaRe",100,0, 10, 100, 0, 100);
+  fTPCSignal3SigmaRep=new TH2F("fTPCSignal3SigmaRep","fTPCSignal3SigmaRe",100,0, 10, 100, 0, 100);
+  fTPCSignal2Sigma=new TH2F("fTPCSignal2Sigma","fTPCSignal2Sigma",100,0, 10, 100, 0, 100);
+  fTPCSignal2SigmaReK=new TH2F("fTPCSignal2SigmaReK","fTPCSignal2SigmaRe",100,0, 10, 100, 0, 100);
+  fTPCSignal2SigmaRep=new TH2F("fTPCSignal2SigmaRep","fTPCSignal2SigmaRe",100,0, 10, 100, 0, 100);
+  
+  for(Int_t i=0; i<3*fNPtBins; i++){
+    fOutput->Add(fMassHist[i]);
+    fOutput->Add(fMassHistTC[i]);
+    fOutput->Add(fCosPHist[i]);
+    fOutput->Add(fDLenHist[i]);
+    fOutput->Add(fSumd02Hist[i]);
+    fOutput->Add(fSigVertHist[i]);
+    fOutput->Add(fPtMaxHist[i]);
+    fOutput->Add(fDCAHist[i]);
+  }
+
+    fOutput->Add(fTPCSignal3Sigma);
+    fOutput->Add(fTPCSignal3SigmaReK);
+    fOutput->Add(fTPCSignal3SigmaRep);
+    fOutput->Add(fTPCSignal2Sigma);
+    fOutput->Add(fTPCSignal2SigmaReK);
+    fOutput->Add(fTPCSignal2SigmaRep);
+  
+
+  fHistNEvents = new TH1F("fHistNEvents", "Number of processed events; ; Events",3,-1.5,1.5);
+  fHistNEvents->Sumw2();
+  fHistNEvents->SetMinimum(0);
+  fOutput->Add(fHistNEvents);
+  
+  
+
+  if(fFillNtuple){
+    //OpenFile(3); // 2 is the slot number of the ntuple
+   
+    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");  
+    
+  }
+  
+  return;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSELambdac::UserExec(Option_t */*option*/)
+{
+  // Execute analysis for current event:
+  // heavy flavor candidates association to MC truth
+
+   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
+   //tmp!
+ fHistNEvents->Fill(0); // count event
+  // Post the data already here
+  PostData(1,fOutput);
+
+  TClonesArray *array3Prong = 0;
+  TClonesArray *arrayLikeSign =0;
+  if(!aod && AODEvent() && IsStandardAOD()) {
+    // In case there is an AOD handler writing a standard AOD, use the AOD 
+    // event in memory rather than the input (ESD) event.    
+    aod = dynamic_cast<AliAODEvent*> (AODEvent());
+    // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
+    // have to taken from the AOD event hold by the AliAODExtension
+    AliAODHandler* aodHandler = (AliAODHandler*) 
+      ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
+    if(aodHandler->GetExtensions()) {
+      AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
+      AliAODEvent *aodFromExt = ext->GetAOD();
+      array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
+      arrayLikeSign=(TClonesArray*)aodFromExt->GetList()->FindObject("LikeSign3Prong");
+    }
+  } else {
+    array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
+    arrayLikeSign=(TClonesArray*)aod->GetList()->FindObject("LikeSign3Prong");
+  }
+
+  if(!array3Prong) {
+    printf("AliAnalysisTaskSELambdac::UserExec: Charm3Prong branch not found!\n");
+    return;
+  }
+  if(!arrayLikeSign) {
+    printf("AliAnalysisTaskSELambdac::UserExec: LikeSign3Prong branch not found!\n");
+  //  return;
+  }
+
+  TClonesArray *arrayMC=0;
+  AliAODMCHeader *mcHeader=0;
+
+  // AOD primary vertex
+  AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
+  
+  // load MC particles
+  if(fReadMC){
+    
+    arrayMC =  (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
+    if(!arrayMC) {
+      printf("AliAnalysisTaskSELambdac::UserExec: MC particles branch not found!\n");
+      return;
+    }
+    
+  // load MC header
+    mcHeader =  (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
+    if(!mcHeader) {
+    printf("AliAnalysisTaskSELambdac::UserExec: MC header branch not found!\n");
+    return;
+    }
+  }
+  
+  Int_t n3Prong = array3Prong->GetEntriesFast();
+  
+  
+  Int_t nOS=0;
+  Int_t index;
+  for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
+    AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
+    
+    
+    Bool_t unsetvtx=kFALSE;
+    if(!d->GetOwnPrimaryVtx()){
+      d->SetOwnPrimaryVtx(vtx1);
+      unsetvtx=kTRUE;
+    }
+
+    if(fRDCutsProduction->IsSelected(d,AliRDHFCuts::kCandidate)) {
+      Int_t iPtBin = -1;
+      Double_t ptCand = d->Pt();
+      
+      for(Int_t ibin=0;ibin<fNPtBins&&iPtBin<0&&ptCand>fArrayBinLimits[0]&&ptCand<fArrayBinLimits[fNPtBins];ibin++){
+       if(ptCand<fArrayBinLimits[ibin+1])iPtBin=ibin;
+      }
+
+      Bool_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kCandidate);
+      Int_t labDp=-1;
+      Float_t deltaPx=0.;
+      Float_t deltaPy=0.;
+      Float_t deltaPz=0.;
+      Float_t truePt=0.;
+      Float_t xDecay=0.;
+      Float_t yDecay=0.;
+      Float_t zDecay=0.;
+      Float_t pdgCode=-2;
+      if(fReadMC){
+       labDp = MatchToMCLambdac(d,arrayMC);   
+       //
+       if(labDp>0){
+         AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
+         AliAODMCParticle *dg0 = (AliAODMCParticle*)arrayMC->At(partDp->GetDaughter(0));  
+         deltaPx=partDp->Px()-d->Px();
+         deltaPy=partDp->Py()-d->Py();
+         deltaPz=partDp->Pz()-d->Pz();
+         truePt=partDp->Pt();
+         xDecay=dg0->Xv();       
+         yDecay=dg0->Yv();       
+         zDecay=dg0->Zv();
+         pdgCode=TMath::Abs(partDp->GetPdgCode());
+       }else{
+         pdgCode=-1;
+       }
+      }
+
+      Double_t invMasspKpi=-1.;
+      Double_t invMasspiKp=-1.;
+      Int_t pdgs[3]={0,0,0};
+      Double_t field=aod->GetMagneticField();
+      if(fReadMC && fMCPid){
+       
+       if(IspKpiMC(d,arrayMC)) {
+        invMasspKpi=d->InvMassLcpKpi();
+       if(fUseKF){
+        pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
+        if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
+       }
+       }
+       IspiKpReal(d);
+       IspiKpResonant(d,field);
+       if(IspiKpMC(d,arrayMC)) {
+        invMasspiKp=d->InvMassLcpiKp();
+       if(fUseKF){
+        pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
+        if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
+       }
+       }
+      }
+      if(fRealPid){
+       if(IspKpiReal(d)) {
+        invMasspKpi=d->InvMassLcpKpi();
+       if(fUseKF){
+        pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
+        if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
+       }
+       }
+       if(IspiKpReal(d)) {
+        invMasspiKp=d->InvMassLcpiKp();
+       if(fUseKF){
+        pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
+        if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
+       }
+       }
+      }
+      if(fResPid){
+       if(IspKpiResonant(d,field)) {
+        invMasspKpi=d->InvMassLcpKpi();
+       if(fUseKF){
+        pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
+        if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
+       }
+       }
+       if(IspiKpResonant(d,field)) {
+        invMasspiKp=d->InvMassLcpiKp();
+       if(fUseKF){
+        pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
+        if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
+       }
+       }
+      }
+      if(!fResPid && !fRealPid && !fMCPid){
+       invMasspiKp=d->InvMassLcpiKp(); 
+       if(fUseKF){
+         pdgs[0]=211;pdgs[1]=321;pdgs[2]=2212;
+         if(!VertexingKF(d,pdgs,field)) invMasspiKp=-1.;
+        }
+       invMasspKpi=d->InvMassLcpKpi();
+       if(fUseKF){
+        pdgs[0]=2212;pdgs[1]=321;pdgs[2]=211;
+        if(!VertexingKF(d,pdgs,field)) invMasspKpi=-1.;
+       }
+      }
+
+      if(invMasspiKp<0. && invMasspKpi<0.) continue;
+
+
+      Float_t tmp[24];
+      if(fFillNtuple){           
+       tmp[0]=pdgCode;
+       tmp[1]=deltaPx;
+       tmp[2]=deltaPy;
+       tmp[3]=deltaPz;
+       tmp[4]=truePt;
+       tmp[5]=xDecay;    
+       tmp[6]=yDecay;    
+       tmp[7]=zDecay;    
+       tmp[8]=d->PtProng(0);
+       tmp[9]=d->PtProng(1);
+       tmp[10]=d->PtProng(2);
+       tmp[11]=d->Pt();
+       tmp[12]=d->CosPointingAngle();
+       tmp[13]=d->DecayLength();
+       tmp[14]=d->Xv();
+       tmp[15]=d->Yv();
+       tmp[16]=d->Zv();
+       if(invMasspiKp>0.) tmp[17]=invMasspiKp;
+       if(invMasspKpi>0.) tmp[17]=invMasspKpi;
+       tmp[18]=d->GetSigmaVert();
+       tmp[19]=d->Getd0Prong(0);
+       tmp[20]=d->Getd0Prong(1);
+       tmp[21]=d->Getd0Prong(2);
+       tmp[22]=d->GetDCA();
+       tmp[23]=d->Prodd0d0(); 
+       fNtupleLambdac->Fill(tmp);
+       PostData(3,fNtupleLambdac);
+      }
+      Double_t dlen=d->DecayLength();
+      Double_t cosp=d->CosPointingAngle();
+      Double_t sumD02=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
+      Double_t dca=d->GetDCA();      
+Double_t ptmax=0;
+      for(Int_t i=0;i<3;i++){
+       if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
+      }
+      if(iPtBin>=0){
+      
+       index=GetHistoIndex(iPtBin);
+       if(invMasspiKp>0. && invMasspKpi>0.){
+       if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp,0.5);
+       if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi,0.5);
+       }else{
+         if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
+         if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
+        }
+       fCosPHist[index]->Fill(cosp);
+       fDLenHist[index]->Fill(dlen);
+       fSumd02Hist[index]->Fill(sumD02);
+       fPtMaxHist[index]->Fill(ptmax);
+       fDCAHist[index]->Fill(dca);
+       
+       if(passTightCuts){
+        if(invMasspiKp>0. && invMasspKpi>0.){
+         if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp,0.5);
+         if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi,0.5);
+        }else{
+          if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
+          if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
+         }
+       }
+       
+       if(fReadMC){
+         if(labDp>0) {
+           index=GetSignalHistoIndex(iPtBin);
+           if(invMasspiKp>0. && invMasspKpi>0.){
+            if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp,0.5);
+            if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi,0.5);
+            }else{
+             if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
+             if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
+            }
+           fCosPHist[index]->Fill(cosp);
+           fDLenHist[index]->Fill(dlen);
+           fSumd02Hist[index]->Fill(sumD02);
+           fPtMaxHist[index]->Fill(ptmax);
+           fDCAHist[index]->Fill(dca);
+           if(passTightCuts){
+            if(invMasspiKp>0. && invMasspKpi>0.){
+             if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp,0.5);
+             if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi,0.5);
+             }else{
+              if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
+              if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
+             }
+           }
+           
+         }else{
+           index=GetBackgroundHistoIndex(iPtBin);
+           if(invMasspiKp>0. && invMasspKpi>0.){
+            fMassHist[index]->Fill(invMasspiKp,0.5);
+            fMassHist[index]->Fill(invMasspKpi,0.5);
+           }else{
+            if(invMasspiKp>0.) fMassHist[index]->Fill(invMasspiKp);
+            if(invMasspKpi>0.) fMassHist[index]->Fill(invMasspKpi);
+          }
+           fCosPHist[index]->Fill(cosp);
+           fDLenHist[index]->Fill(dlen);
+           fSumd02Hist[index]->Fill(sumD02);
+           fPtMaxHist[index]->Fill(ptmax);
+           fDCAHist[index]->Fill(dca);
+           if(passTightCuts){
+           if(invMasspiKp>0. && invMasspKpi>0.){
+              fMassHistTC[index]->Fill(invMasspiKp,0.5);
+              fMassHistTC[index]->Fill(invMasspKpi,0.5);
+             }else{
+              if(invMasspiKp>0.) fMassHistTC[index]->Fill(invMasspiKp);
+              if(invMasspKpi>0.) fMassHistTC[index]->Fill(invMasspKpi);
+             }
+           }   
+         }
+       }
+      }
+      /*
+      //start OS analysis
+      if(labDp<0)fHistOSbkg->Fill(d->InvMassDplus());
+      fHistOS->Fill(d->InvMassDplus());
+      */
+      nOS++;
+    }
+    if(unsetvtx) d->UnsetOwnPrimaryVtx();
+  }
+  PostData(1,fOutput);    
+  return;
+}
+
+
+
+//________________________________________________________________________
+void AliAnalysisTaskSELambdac::Terminate(Option_t */*option*/)
+{
+  // Terminate analysis
+  //
+  if(fDebug > 1) printf("AnalysisTaskSELambdac: Terminate() \n");
+
+  fOutput = dynamic_cast<TList*> (GetOutputData(1));
+  if (!fOutput) {     
+    printf("ERROR: fOutput not available\n");
+    return;
+  }
+ fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
+
+ TString hisname;
+ Int_t index=0;
+
+ //Int_t indexLS=0;
+ for(Int_t i=0;i<fNPtBins;i++){
+    index=GetHistoIndex(i);
+    hisname.Form("hMassPt%d",i);
+    fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hCosPAllPt%d",i);
+    fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hDLenAllPt%d",i);
+    fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hSumd02AllPt%d",i);
+     fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hSigVertAllPt%d",i);
+     fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hPtMaxAllPt%d",i);
+     fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hDCAAllPt%d",i);
+     fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+     hisname.Form("hMassPt%dTC",i);
+    fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    
+    index=GetSignalHistoIndex(i);    
+    hisname.Form("hSigPt%d",i);
+    fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hCosPSigPt%d",i);
+    fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hDLenSigPt%d",i);
+    fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hSumd02SigPt%d",i);
+    fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hSigVertSigPt%d",i);
+    fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hPtMaxSigPt%d",i);
+    fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hDCASigPt%d",i);
+    fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    
+    hisname.Form("hSigPt%dTC",i);
+    fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    
+    index=GetBackgroundHistoIndex(i); 
+    hisname.Form("hBkgPt%d",i);
+    fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hCosPBkgPt%d",i);
+    fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hDLenBkgPt%d",i);
+    fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hSumd02BkgPt%d",i);
+    fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hSigVertBkgPt%d",i);
+    fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hPtMaxBkgPt%d",i);
+    fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hDCABkgPt%d",i);
+    fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+    hisname.Form("hBkgPt%dTC",i);
+    fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
+ }
+    fTPCSignal3Sigma=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3Sigma));
+    fTPCSignal3SigmaReK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaReK));
+    fTPCSignal3SigmaRep=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal3SigmaRep));
+    fTPCSignal2Sigma=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2Sigma));
+    fTPCSignal2SigmaReK=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaReK));
+    fTPCSignal2SigmaRep=dynamic_cast<TH2F*>(fOutput->FindObject(fTPCSignal2SigmaRep));
+
+  if(fFillNtuple){
+    fNtupleLambdac = dynamic_cast<TNtuple*>(GetOutputData(3));
+  }
+
+  //TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
+  //c1->cd();
+  //TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
+  //hMassPt->Draw();
+ return;
+}
+
+//________________________________________________________________________
+Int_t AliAnalysisTaskSELambdac::MatchToMCLambdac(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
+
+  Int_t lambdacLab[3]={0,0,0};
+  Int_t pdgs[3]={0,0,0};
+  for(Int_t i=0;i<3;i++){
+   AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
+   Int_t lab=daugh->GetLabel();
+   if(lab<0) return 0;
+   AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab);
+   if(!part) continue;
+   pdgs[i]=part->GetPdgCode();
+   Int_t partPdgcode = TMath::Abs(part->GetPdgCode());
+   if(partPdgcode==211 || partPdgcode==321 || partPdgcode==2212){
+        Int_t motherLabel=part->GetMother();
+            if(motherLabel<0) return 0;
+   AliAODMCParticle *motherPart = (AliAODMCParticle*)arrayMC->At(motherLabel);
+   if(!motherPart) continue;
+        Int_t motherPdg = TMath::Abs(motherPart->GetPdgCode());
+   if(motherPdg==4122) {
+        if(GetLambdacDaugh(motherPart,arrayMC)){lambdacLab[i]=motherLabel;continue;}
+             }
+   if(motherPdg==313 || motherPdg==2224 || motherPdg==3124){
+         Int_t GmotherLabel=motherPart->GetMother();
+              if(GmotherLabel<0) return 0;
+        AliAODMCParticle *GmotherPart = (AliAODMCParticle*)arrayMC->At(GmotherLabel);
+        if(!GmotherPart) continue;
+              Int_t GmotherPdg = TMath::Abs(GmotherPart->GetPdgCode());
+        if(GmotherPdg==4122) {
+               if(GetLambdacDaugh(GmotherPart,arrayMC)) {lambdacLab[i]=GmotherLabel;continue;}
+             }
+       }
+     }
+  }
+
+ if(lambdacLab[0]==lambdacLab[1] && lambdacLab[1]==lambdacLab[2]) {return lambdacLab[0];}
+ return 0;
+
+}
+//------------------------
+Bool_t AliAnalysisTaskSELambdac::GetLambdacDaugh(AliAODMCParticle *part,TClonesArray *arrayMC){
+
+ Int_t numberOfLambdac=0;
+ if(TMath::Abs(part->GetPdgCode())!=4122) return kFALSE;
+ Int_t daugh_tmp[2];
+ daugh_tmp[0]=part->GetDaughter(0);
+ daugh_tmp[1]=part->GetDaughter(1);
+ Int_t nDaugh = (Int_t)part->GetNDaughters();
+ if(nDaugh<2) return kFALSE;
+ if(nDaugh>3) return kFALSE;
+ AliAODMCParticle* pdaugh1 = (AliAODMCParticle*)arrayMC->At(part->GetDaughter(0));
+   if(!pdaugh1) {return kFALSE;}
+ Int_t number1 = TMath::Abs(pdaugh1->GetPdgCode());
+ AliAODMCParticle* pdaugh2 = (AliAODMCParticle*)arrayMC->At(part->GetDaughter(1));
+   if(!pdaugh2) {return kFALSE;}
+  Int_t number2 = TMath::Abs(pdaugh2->GetPdgCode());
+
+  if(nDaugh==3){
+   Int_t thirdDaugh=part->GetDaughter(1)-1;
+            AliAODMCParticle* pdaugh3 = (AliAODMCParticle*)arrayMC->At(thirdDaugh);
+  Int_t number3 = TMath::Abs(pdaugh3->GetPdgCode());
+   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)) {
+   numberOfLambdac++;
+   } 
+  }
+
+ if(nDaugh==2){
+
+         //Lambda in Lambda pi
+//  if((number1==211 && number2==3122)|| (number1==3122 && number2==211))return kFALSE;
+  
+  //Lambda resonant
+  
+  //Lambda -> p K*0
+  //
+  Int_t nfiglieK=0;
+
+   if((number1==2212 && number2==313)){
+     nfiglieK=pdaugh2->GetNDaughters();
+     if(nfiglieK!=2) return kFALSE;
+     AliAODMCParticle* pdaughK1 = (AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
+     AliAODMCParticle* pdaughK2 = (AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
+     if(!pdaughK1) return kFALSE;
+     if(!pdaughK2) return kFALSE;
+     if((TMath::Abs(pdaughK1->GetPdgCode())==211 && TMath::Abs(pdaughK2->GetPdgCode())==321) || (TMath::Abs(pdaughK1->GetPdgCode())==321 && TMath::Abs(pdaughK2->GetPdgCode())==211)) numberOfLambdac++;
+    }
+
+   if((number1==313 && number2==2212)){
+    nfiglieK=pdaugh1->GetNDaughters();
+    if(nfiglieK!=2) return kFALSE;
+    AliAODMCParticle* pdaughK1 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
+    AliAODMCParticle* pdaughK2 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
+    if(!pdaughK1) return kFALSE;
+    if(!pdaughK2) return kFALSE;
+     if((TMath::Abs(pdaughK1->GetPdgCode())==211 && TMath::Abs(pdaughK2->GetPdgCode())==321) || (TMath::Abs(pdaughK1->GetPdgCode())==321 && TMath::Abs(pdaughK2->GetPdgCode())==211)) numberOfLambdac++;
+   }
+
+   //Lambda -> Delta++ k
+   Int_t nfiglieDelta=0;
+   if(number1==321 && number2==2224){
+    nfiglieDelta=pdaugh2->GetNDaughters();
+    if(nfiglieDelta!=2) return kFALSE;
+    AliAODMCParticle *pdaughD1=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
+    AliAODMCParticle *pdaughD2=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
+    if(!pdaughD1) return kFALSE;
+    if(!pdaughD2) return kFALSE;
+    if((TMath::Abs(pdaughD1->GetPdgCode())==211 && TMath::Abs(pdaughD2->GetPdgCode())==2212) || (TMath::Abs(pdaughD1->GetPdgCode())==2212 && TMath::Abs(pdaughD2->GetPdgCode())==211)) numberOfLambdac++;
+   }
+   if(number1==2224 && number2==321){
+    nfiglieDelta=pdaugh1->GetNDaughters();
+    if(nfiglieDelta!=2) return kFALSE;
+    AliAODMCParticle* pdaughD1 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
+    AliAODMCParticle* pdaughD2 = (AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
+    if(!pdaughD1) return kFALSE;
+    if(!pdaughD2) return kFALSE;
+    if((TMath::Abs(pdaughD1->GetPdgCode())==211 && TMath::Abs(pdaughD2->GetPdgCode())==2212) || (TMath::Abs(pdaughD1->GetPdgCode())==2212 && TMath::Abs(pdaughD2->GetPdgCode())==211)) numberOfLambdac++;
+   }
+    
+
+   //Lambdac -> Lambda(1520) pi
+   Int_t nfiglieLa=0;
+   if(number1==3124 && number2==211){
+    nfiglieLa=pdaugh1->GetNDaughters();
+    if(nfiglieLa!=2) return kFALSE;
+    AliAODMCParticle *pdaughL1=(AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(0));
+    AliAODMCParticle *pdaughL2=(AliAODMCParticle*)arrayMC->At(pdaugh1->GetDaughter(1));
+    if(!pdaughL1) return kFALSE;
+    if(!pdaughL2) return kFALSE;
+    if((TMath::Abs(pdaughL1->GetPdgCode())==321 && TMath::Abs(pdaughL2->GetPdgCode())==2212) || (TMath::Abs(pdaughL1->GetPdgCode())==2212 && TMath::Abs(pdaughL2->GetPdgCode())==321)) numberOfLambdac++;
+   }
+   if(number1==211 && number2==3124){
+    nfiglieLa=pdaugh2->GetNDaughters();
+    if(nfiglieLa!=2) return kFALSE;
+    AliAODMCParticle *pdaughL1=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(0));
+    AliAODMCParticle *pdaughL2=(AliAODMCParticle*)arrayMC->At(pdaugh2->GetDaughter(1));
+    if(!pdaughL1) return kFALSE;
+    if(!pdaughL2) return kFALSE;
+    if((TMath::Abs(pdaughL1->GetPdgCode())==321 && TMath::Abs(pdaughL2->GetPdgCode())==2212) || (TMath::Abs(pdaughL1->GetPdgCode())==2212 && TMath::Abs(pdaughL2->GetPdgCode())==321)) numberOfLambdac++;
+   
+   }
+ }
+
+ if(numberOfLambdac>0) {return kTRUE;}
+         return kFALSE;
+}
+//-----------------------------
+Bool_t AliAnalysisTaskSELambdac::IspKpiMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
+
+   Int_t lab[3]={0,0,0},pdgs[3]={0,0,0};
+   for(Int_t i=0;i<3;i++){
+    AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
+    lab[i]=daugh->GetLabel();
+    if(lab[i]<0) return kFALSE;
+    AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab[i]);
+    if(!part) return kFALSE;
+    pdgs[i]=TMath::Abs(part->GetPdgCode());
+   }
+
+   if(pdgs[0]==2212 && pdgs[1]==321 && pdgs[2]==211) return kTRUE;
+
+   return kFALSE;
+}
+//-----------------------------
+Bool_t AliAnalysisTaskSELambdac::IspiKpMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC){
+
+   Int_t lab[3]={0,0,0},pdgs[3]={0,0,0};
+   for(Int_t i=0;i<3;i++){
+    AliAODTrack *daugh=(AliAODTrack*)d->GetDaughter(i);
+    lab[i]=daugh->GetLabel();
+    if(lab[i]<0) return kFALSE;
+    AliAODMCParticle *part= (AliAODMCParticle*)arrayMC->At(lab[i]);
+    if(!part) return kFALSE;
+    pdgs[i]=TMath::Abs(part->GetPdgCode());
+   }
+
+   if(pdgs[2]==2212 && pdgs[1]==321 && pdgs[0]==211) {return kTRUE;}
+
+   return kFALSE;
+}
+//--------------------------------------
+Bool_t AliAnalysisTaskSELambdac::IspiKpReal(AliAODRecoDecayHF3Prong *d){
+
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  Bool_t type[2]={kFALSE,kFALSE};
+  for(Int_t i=0;i<3;i++){
+//   Bool_t pid[3]={kFALSE,kFALSE,kFALSE};
+   AliAODTrack *track=(AliAODTrack*)d->GetDaughter(i);
+   AliAODPid *pidObjtrk=track->GetDetPid();
+//   pidObj->CombinedProbability(track,pid);
+
+//   Int_t i2=i-1;
+//   type[i2]=pid[i];
+   Double_t dedx=pidObjtrk->GetTPCsignal();
+   Double_t mom = pidObjtrk->GetTPCmomentum();
+   //printf("TPC mom= %f\n",mom);
+    pidObj->SetSigma(3.);
+   if(pidObj->IsProtonRaw(track,"TPC")){ 
+    fTPCSignal3Sigma->Fill(mom,dedx);
+ //   if(i==2) type[1]=kTRUE;
+    pidObj->SetSigma(1.);
+    //if(!pidObj->IsPionRaw(track,"TPC") && !pidObj->IsKaonRaw(track,"TPC"))
+    if(!pidObj->IsKaonRaw(track,"TPC")){
+     fTPCSignal3SigmaRep->Fill(mom,dedx);
+    }
+   }
+    pidObj->SetSigma(3.);
+   if(pidObj->IsKaonRaw(track,"TPC")) {
+    fTPCSignal3Sigma->Fill(mom,dedx);
+ //   if(i==1) type[0]=kTRUE;
+    pidObj->SetSigma(1.);
+    //if(!pidObj->IsPionRaw(track,"TPC") && !pidObj->IsProtonRaw(track,"TPC"))
+    if(!pidObj->IsPionRaw(track,"TPC")){
+     fTPCSignal3SigmaReK->Fill(mom,dedx);
+    }
+   }
+
+    pidObj->SetSigma(2.);
+   if(pidObj->IsProtonRaw(track,"TPC")){ 
+    if(i==2) type[1]=kTRUE;
+    fTPCSignal2Sigma->Fill(mom,dedx);
+    pidObj->SetSigma(1.);
+    //if(!pidObj->IsPionRaw(track,"TPC") && !pidObj->IsKaonRaw(track,"TPC"))
+    if(!pidObj->IsKaonRaw(track,"TPC")){
+     fTPCSignal2SigmaRep->Fill(mom,dedx);
+    }
+   }
+   if(pidObj->IsKaonRaw(track,"TPC")) {
+    if(i==1) type[0]=kTRUE;
+    fTPCSignal2Sigma->Fill(mom,dedx);
+    pidObj->SetSigma(1.);
+    //if(!pidObj->IsPionRaw(track,"TPC") && !pidObj->IsProtonRaw(track,"TPC"))
+    if(!pidObj->IsPionRaw(track,"TPC")){
+     fTPCSignal2SigmaReK->Fill(mom,dedx);
+    }
+   }
+  }
+
+  if(type[0] && type[1]) {delete pidObj;return kTRUE;}
+  delete pidObj;
+ return kFALSE;
+}
+//--------------------------------------
+Bool_t AliAnalysisTaskSELambdac::IspKpiReal(AliAODRecoDecayHF3Prong *d){
+
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  Bool_t type[2]={kFALSE,kFALSE};
+    pidObj->SetSigma(2.);
+  for(Int_t i=0;i<2;i++){
+   //Bool_t pid[3]={kFALSE,kFALSE,kFALSE};
+   AliAODTrack *track=(AliAODTrack*)d->GetDaughter(i);
+  // pidObj->CombinedProbability(track,pid);
+  // if(i==0) type[i]=pid[2];
+  // if(i==1) type[i]=pid[i];
+   if(pidObj->IsProtonRaw(track,"TPC") && i==0) type[1]=kTRUE;
+   if(pidObj->IsKaonRaw(track,"TPC") && i==1) type[0]=kTRUE;
+  }
+  if(type[0] && type[1]) {delete pidObj;return kTRUE;}
+  delete pidObj;
+ return kFALSE;
+}
+//------------------------------------
+Bool_t AliAnalysisTaskSELambdac::VertexingKF(AliAODRecoDecayHF3Prong *d,Int_t *pdgs,Double_t field){
+  
+   Int_t iprongs[3]={0,1,2};
+   Double_t mass[2]={0.,0.};
+ //topological constr
+   AliKFParticle *Lambdac=d->ApplyVertexingKF(iprongs,3,pdgs,kTRUE,field,mass);
+   if(!Lambdac) return kFALSE;
+  Double_t probTot=TMath::Prob(Lambdac->GetChi2(),Lambdac->GetNDF());
+  if(probTot<fCutsKF[0]) return kFALSE;
+
+  if(d->PtProng(0)<1.5 && d->Getd0Prong(0)>0.02) return kFALSE;
+    if(d->PtProng(1)<1.5 && d->Getd0Prong(1)>0.02) return kFALSE;
+      if(d->PtProng(2)<1.5 && d->Getd0Prong(2)>0.02) return kFALSE;
+ //mass constr for K*
+   Int_t ipRes[2];
+   Int_t pdgres[2];
+   mass[0]=0.8961;mass[1]=0.03;
+   if(TMath::Abs(pdgs[0])==211){
+    ipRes[0]=0;ipRes[1]=1;
+    pdgres[0]=pdgs[0];pdgres[1]=321;
+   }
+   if(TMath::Abs(pdgs[2])==211){
+    ipRes[0]=2;ipRes[1]=1;
+    pdgres[0]=pdgs[2];pdgres[1]=321;
+   }
+   AliKFParticle *KappaStar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+
+  Double_t probKstar=TMath::Prob(KappaStar->GetChi2(),KappaStar->GetNDF());
+  if(probKstar>fCutsKF[1]) {
+    AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
+    AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
+    AliKFParticle prong1(*esdProng1,pdgres[0]);
+    AliKFParticle prong2(*esdProng2,pdgres[1]);
+    if(KappaStar->GetPt()<fCutsKF[2] && prong1.GetAngle(prong2)>fCutsKF[3]) return kFALSE;
+  } 
+ //mass constr for Lambda
+   mass[0]=1.520;mass[1]=0.005;
+   if(TMath::Abs(pdgs[0])==2212){
+    ipRes[0]=0;ipRes[1]=1;
+    pdgres[0]=pdgs[0];pdgres[1]=pdgs[1];
+   }
+   if(TMath::Abs(pdgs[2])==2212){
+    ipRes[0]=2;ipRes[1]=1;
+    pdgres[0]=pdgs[2];pdgres[1]=pdgs[1];
+   }
+   AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+  Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
+  if(probLa>fCutsKF[4]) {
+    AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
+    AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
+    AliKFParticle prong1(*esdProng1,pdgres[0]);
+    AliKFParticle prong2(*esdProng2,pdgres[1]);
+    if(Lambda1520->GetPt()<fCutsKF[5] && prong1.GetAngle(prong2)>fCutsKF[6]) return kFALSE;
+  } 
+ //mass constr for Delta
+   mass[0]=1.2;mass[1]=0.15;
+   ipRes[0]=0;ipRes[1]=2;
+   pdgres[0]=pdgs[0];pdgres[2]=pdgs[2];
+   AliKFParticle *Delta=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+  Double_t probDelta=TMath::Prob(Delta->GetChi2(),Delta->GetNDF());
+  if(probDelta>fCutsKF[7]) {
+    AliAODTrack *esdProng1=(AliAODTrack*)d->GetDaughter(ipRes[0]);
+    AliAODTrack *esdProng2=(AliAODTrack*)d->GetDaughter(ipRes[1]);
+    AliKFParticle prong1(*esdProng1,pdgres[0]);
+    AliKFParticle prong2(*esdProng2,pdgres[1]);
+    if(Delta->GetPt()<fCutsKF[8] && prong1.GetAngle(prong2)>fCutsKF[9]) return kFALSE;
+  } 
+ return kTRUE;
+}
+//-------------------------------------
+Bool_t AliAnalysisTaskSELambdac::IspiKpResonant(AliAODRecoDecayHF3Prong *d,Double_t field){
+  
+ //if lambda* -> pk
+        Double_t mass[2]={1.520,0.005};
+        Int_t ipRes[2]={1,2};
+        Int_t pdgres[2]={321,2212};
+        AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+       Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
+       if(probLa>0.1) return kTRUE;
+ //K* -> kpi
+//        mass[0]=0.8961;mass[1]=0.03;
+//        ipRes[0]=0;ipRes[1]=1;
+//        pdgres[0]=211;pdgres[1]=321;
+//        AliKFParticle *Kstar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+//     Double_t probKa=TMath::Prob(Kstar->GetChi2(),Kstar->GetNDF());
+//     if(probKa>0.1) return kTRUE;
+
+ return kFALSE;
+
+}
+//-------------------------------------
+Bool_t AliAnalysisTaskSELambdac::IspKpiResonant(AliAODRecoDecayHF3Prong *d,Double_t field){
+  
+ //if lambda* -> pk
+        Double_t mass[2]={1.520,0.005};
+        Int_t ipRes[2]={0,1};
+        Int_t pdgres[2]={2212,321};
+        AliKFParticle *Lambda1520=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+       Double_t probLa=TMath::Prob(Lambda1520->GetChi2(),Lambda1520->GetNDF());
+       if(probLa>0.1) return kTRUE;
+ //K* -> kpi
+//        mass[0]=0.8961;mass[1]=0.03;
+//        ipRes[0]=1;ipRes[1]=2;
+//        pdgres[1]=211;pdgres[0]=321;
+//        AliKFParticle *Kstar=d->ApplyVertexingKF(ipRes,2,pdgres,kFALSE,field,mass);
+//     Double_t probKa=TMath::Prob(Kstar->GetChi2(),Kstar->GetNDF());
+//     if(probKa>0.1) return kTRUE;
+
+ return kFALSE;
+
+}
diff --git a/PWG3/vertexingHF/AliAnalysisTaskSELambdac.h b/PWG3/vertexingHF/AliAnalysisTaskSELambdac.h
new file mode 100644 (file)
index 0000000..021bbb5
--- /dev/null
@@ -0,0 +1,124 @@
+#ifndef ALIANALYSISTASKLAMBDAC_H
+#define ALIANALYSISTASKLAMBDAC_H
+
+/* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+//*************************************************************************
+// Class AliAnalysisTaskSELambdac
+// AliAnalysisTaskSE for the D+ candidates Invariant Mass Histogram and 
+//comparison of heavy-flavour decay candidates
+// to MC truth (kinematics stored in the AOD)
+// Renu Bala, bala@to.infn.it
+//*************************************************************************
+
+#include <TROOT.h>
+#include <TSystem.h>
+#include <TNtuple.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TArrayD.h>
+
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisVertexingHF.h"
+#include "AliAODMCParticle.h"
+#include "AliRDHFCutsLctopKpi.h"
+#include "AliRDHFCuts.h"
+#include "TClonesArray.h"
+
+class AliAnalysisTaskSELambdac : public AliAnalysisTaskSE
+{
+ public:
+
+  AliAnalysisTaskSELambdac();
+  AliAnalysisTaskSELambdac(const char *name, Bool_t fillNtuple,AliRDHFCutsLctopKpi *lccutsana, AliRDHFCutsLctopKpi *lccutsprod);
+  virtual ~AliAnalysisTaskSELambdac();
+
+  void SetReadMC(Bool_t readMC=kTRUE){fReadMC=readMC;}
+  void SetMCPid(){fMCPid=kTRUE;fReadMC=kTRUE;fRealPid=kFALSE;fResPid=kFALSE;return;}
+  void SetRealPid(){fRealPid=kTRUE;fMCPid=kFALSE;fResPid=kFALSE;return;}
+  void SetResonantPid(){fResPid=kTRUE;fRealPid=kFALSE;fMCPid=kFALSE;return;}
+  void SetCutsKF(Float_t cutsKF[10]){for(Int_t i=0;i<10;i++){fCutsKF[i]=cutsKF[i];}return;}
+  void SetUseKF(Bool_t useKF=kTRUE){fUseKF=useKF;}
+  void SetMassLimits(Float_t range);
+  void SetMassLimits(Float_t lowlimit, Float_t uplimit);
+  void SetPtBinLimit(Int_t n, Float_t *limitarray);
+  
+  Float_t GetUpperMassLimit(){return fUpmasslimit;}
+  Float_t GetLowerMassLimit(){return fLowmasslimit;}
+  Int_t GetNBinsPt(){return fNPtBins;}
+  Double_t GetPtBinLimit(Int_t ibin);
+  Bool_t IspiKpMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC);
+  Bool_t IspKpiMC(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC);
+  Bool_t IspiKpReal(AliAODRecoDecayHF3Prong *d);
+  Bool_t IspKpiReal(AliAODRecoDecayHF3Prong *d);
+  Bool_t IspiKpResonant(AliAODRecoDecayHF3Prong *d,Double_t field);
+  Bool_t IspKpiResonant(AliAODRecoDecayHF3Prong *d,Double_t field);
+  Bool_t VertexingKF(AliAODRecoDecayHF3Prong *d,Int_t *pdgs,Double_t field);
+  Int_t MatchToMCLambdac(AliAODRecoDecayHF3Prong *d,TClonesArray *arrayMC);
+  Bool_t GetLambdacDaugh(AliAODMCParticle *part, TClonesArray *arrayMC);
+
+  // Implementation of interface methods
+  virtual void UserCreateOutputObjects();
+  virtual void Init();
+  virtual void LocalInit() {Init();}
+  virtual void UserExec(Option_t *option);
+  virtual void Terminate(Option_t *option);
+    
+ private:
+
+  AliAnalysisTaskSELambdac(const AliAnalysisTaskSELambdac &source);
+  AliAnalysisTaskSELambdac& operator=(const AliAnalysisTaskSELambdac& source); 
+  Int_t GetHistoIndex(Int_t iPtBin) const { return iPtBin*3;}
+  Int_t GetSignalHistoIndex(Int_t iPtBin) const { return iPtBin*3+1;}
+  Int_t GetBackgroundHistoIndex(Int_t iPtBin) const { return iPtBin*3+2;}
+  Int_t GetLSHistoIndex(Int_t iPtBin)const { return iPtBin*5;}
+  enum {kMaxPtBins=10};
+
+  TList   *fOutput; //! list send on output slot 0
+  TH1F    *fHistNEvents; //!hist. for No. of events
+  TH1F *fMassHist[3*kMaxPtBins]; //!hist. for inv mass (LC)
+  TH1F*   fCosPHist[3*kMaxPtBins]; //!hist. for PointingAngle (LC)
+  TH1F*   fDLenHist[3*kMaxPtBins]; //!hist. for Dec Length (LC)
+  TH1F*   fSumd02Hist[3*kMaxPtBins]; //!hist. for sum d02 (LC)
+  TH1F*   fSigVertHist[3*kMaxPtBins]; //!hist. for sigVert (LC)
+  TH1F*   fPtMaxHist[3*kMaxPtBins]; //!hist. for Pt Max (LC)
+  TH1F*   fDCAHist[3*kMaxPtBins]; //!hist. for DCA (LC)
+  TH1F *fMassHistTC[3*kMaxPtBins]; //!hist. for inv mass (TC)
+  TH1F *fMassHistLS[5*kMaxPtBins];//!hist. for LS inv mass (LC)
+  TH1F *fCosPHistLS[3*kMaxPtBins];//!hist. for LS cuts variable 1 (LC)
+  TH1F *fDLenHistLS[3*kMaxPtBins];//!hist. for LS cuts variable 2 (LC)
+  TH1F *fSumd02HistLS[3*kMaxPtBins];//!hist. for LS cuts variable 3 (LC)
+  TH1F *fSigVertHistLS[3*kMaxPtBins];//!hist. for LS cuts variable 4 (LC)
+  TH1F *fPtMaxHistLS[3*kMaxPtBins];//!hist. for LS cuts variable 5 (LC)
+  TH1F *fDCAHistLS[3*kMaxPtBins];//!hist. for LS cuts variable 6 (LC)
+  TH1F *fMassHistLSTC[5*kMaxPtBins];//!hist. for LS inv mass (TC)
+  TH2F *fTPCSignal3Sigma;//!hist. for LS inv mass (TC)
+  TH2F *fTPCSignal3SigmaReK;//!hist. for LS inv mass (TC)
+  TH2F *fTPCSignal3SigmaRep;//!hist. for LS inv mass (TC)
+  TH2F *fTPCSignal2Sigma;//!hist. for LS inv mass (TC)
+  TH2F *fTPCSignal2SigmaReK;//!hist. for LS inv mass (TC)
+  TH2F *fTPCSignal2SigmaRep;//!hist. for LS inv mass (TC)
+  TNtuple *fNtupleLambdac; //! output ntuple
+  Float_t fUpmasslimit;  //upper inv mass limit for histos
+  Float_t fLowmasslimit; //lower inv mass limit for histos
+  Float_t fCutsKF[10];
+  Int_t fNPtBins; //number of bins in Pt for histograms
+  AliRDHFCutsLctopKpi *fRDCutsAnalysis; //Cuts for Analysis
+  AliRDHFCutsLctopKpi *fRDCutsProduction; //Production Cuts
+  TList *fListCuts; //list of cuts
+  Double_t fArrayBinLimits[kMaxPtBins+1]; //limits for the Pt bins
+  Bool_t fFillNtuple;   // flag for filling ntuple
+  Bool_t fReadMC;    //flag for access to MC
+  Bool_t fMCPid;    //flag for access to MC
+  Bool_t fRealPid;    //flag for access to MC
+  Bool_t fResPid;      //flag to do LS analysis
+  Bool_t fUseKF;      //flag to do LS analysis
+  AliAnalysisVertexingHF *fVHF;  // Vertexer heavy flavour (used to pass the cuts)
+  
+  ClassDef(AliAnalysisTaskSELambdac,1); // AliAnalysisTaskSE for the MC association of heavy-flavour decay candidates
+};
+
+#endif
+