]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
update Dplus correlation task
authorebruna <ebruna@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 5 Sep 2013 12:49:58 +0000 (12:49 +0000)
committerebruna <ebruna@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 5 Sep 2013 12:49:58 +0000 (12:49 +0000)
PWGHF/correlationHF/AliAnalysisTaskSEDplusCorrelations.cxx
PWGHF/correlationHF/AliAnalysisTaskSEDplusCorrelations.h
PWGHF/correlationHF/macros/AddTaskDplusCorr.C

index d07ec29d97a7beed73f93d719c1f35cd1c01b329..2c39583a97203c2e702541d4ad77a73d4c11caf1 100644 (file)
@@ -1,3 +1,4 @@
+
 /**************************************************************************
  * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
  *                                                                        *
@@ -42,6 +43,7 @@
 #include "AliAODRecoDecayHF3Prong.h"
 #include "AliAnalysisVertexingHF.h"
 #include "AliAnalysisTaskSE.h"
+#include "AliVertexingHFUtils.h"
 #include "AliAnalysisTaskSEDplusCorrelations.h"
 #include "AliNormalizationCounter.h"
 #include "AliVParticle.h"
@@ -64,50 +66,36 @@ AliAnalysisTaskSE(),
   fSelect(0),
   fDisplacement(0),
   fHistNEvents(0),
-  fMCSources(0),
-  fK0Origin(0),
-  fKaonOrigin(0),
-  fInvMassK0S(0),
+  fTrig(kTRUE),
   fEventMix(0),
-  fPtVsMass(0),
-  fPtVsMassTC(0),
-  fYVsPt(0),
-  fYVsPtTC(0),
-  fYVsPtSig(0),
-  fYVsPtSigTC(0),
   fUpmasslimit(1.965),
   fLowmasslimit(1.765),
   fNPtBins(0),
   fBinWidth(0.002),
   fListCuts(0),
-  fListCutsAsso(0), 
+//fListCutsAsso(0), 
   fRDCutsAnalysis(0),
   fCuts(0),
   fCounter(0),
   fReadMC(kFALSE),
   fUseBit(kTRUE),
   fMixing(kFALSE),
-  fSystem(kFALSE)
+  fSystem(kFALSE),
+  fReco(kTRUE)
 {
   // Default constructor
-   
+  
   for(Int_t i=0;i<3*kMaxPtBins;i++){
-    fMassVsdPhiHistHad[i]=0;
-    fMassVsdEtaHistHad[i]=0;
-    fMassVsdPhiHistKaon[i]=0;
-    fMassVsdEtaHistKaon[i]=0;
-    fMassVsdPhiHistKshort[i]=0;
-    fMassVsdEtaHistKshort[i]=0;
-    fMassVsdPhiHistLeadHad[i]=0;
-    fMassVsdEtaHistLeadHad[i]=0;
-    fMassHistK0S[i]=0;
-    fLeadPt[i]=0;
-    fMassHist[i]=0;
-    fMassHistTC[i]=0;
-    fMassHistTCPlus[i]=0;
-    fMassHistTCMinus[i]=0;
+   
+    if(fMassHistK0S[i])fMassHistK0S[i]=0;
+    if(fLeadPt[i])fLeadPt[i]=0;
+    if(fMassHist[i])fMassHist[i]=0;
+    if(fMassHistTC[i])fMassHistTC[i]=0;
+    if(fMassHistOrigC[i])fMassHistOrigC[i]=0;
+    if(fMassHistOrigB[i])fMassHistOrigB[i]=0;
+    if(fMassHistMC[i])fMassHistMC[i]=0;
   }
-
+  
   for(Int_t i=0;i<kMaxPtBins+1;i++){
     fArrayBinLimits[i]=0;
   }
@@ -122,53 +110,42 @@ AliAnalysisTaskSEDplusCorrelations::AliAnalysisTaskSEDplusCorrelations(const cha
   fSelect(0),
   fDisplacement(0),
   fHistNEvents(0),
-  fMCSources(0),
-  fK0Origin(0),
-  fKaonOrigin(0),
-  fInvMassK0S(0),
+  fTrig(kTRUE),
   fEventMix(0),
-  fPtVsMass(0),
-  fPtVsMassTC(0),
-  fYVsPt(0),
-  fYVsPtTC(0),
-  fYVsPtSig(0),
-  fYVsPtSigTC(0),
   fUpmasslimit(1.965),
   fLowmasslimit(1.765),
   fNPtBins(0),
   fBinWidth(0.002),
   fListCuts(0),
-  fListCutsAsso(0),
+  //fListCutsAsso(0),
   fRDCutsAnalysis(dpluscutsana),
   fCuts(asstrkcuts),
   fCounter(0),
   fReadMC(kFALSE),
   fUseBit(kTRUE),
   fMixing(kFALSE),
-  fSystem(kFALSE)
+  fSystem(kFALSE),
+  fReco(kTRUE)
 {
   // 
   // Standrd constructor
   //
   fNPtBins=fRDCutsAnalysis->GetNPtBins();
     
+    
   for(Int_t i=0;i<3*kMaxPtBins;i++){
-    fMassVsdPhiHistHad[i]=0;
-    fMassVsdEtaHistHad[i]=0;
-    fMassVsdPhiHistKaon[i]=0;
-    fMassVsdEtaHistKaon[i]=0;
-    fMassVsdPhiHistKshort[i]=0;
-    fMassVsdEtaHistKshort[i]=0;
-    fMassVsdPhiHistLeadHad[i]=0;
-    fMassVsdEtaHistLeadHad[i]=0;
-    fMassHistK0S[i]=0;
-    fLeadPt[i]=0;
-    fMassHist[i]=0;
-    fMassHistTC[i]=0;
-    fMassHistTCPlus[i]=0;
-    fMassHistTCMinus[i]=0;
-  }
    
+    if(fMassHistK0S[i])fMassHistK0S[i]=0;
+    if(fLeadPt[i])fLeadPt[i]=0;
+    if(fMassHist[i])fMassHist[i]=0;
+    if(fMassHistTC[i])fMassHistTC[i]=0;
+    if(fMassHistOrigC[i])fMassHistOrigC[i]=0;
+    if(fMassHistOrigB[i])fMassHistOrigB[i]=0;
+    if(fMassHistMC[i])fMassHistMC[i]=0;
+   
+    
+    
+  }
   
   for(Int_t i=0;i<kMaxPtBins+1;i++){
     fArrayBinLimits[i]=0;
@@ -179,11 +156,10 @@ AliAnalysisTaskSEDplusCorrelations::AliAnalysisTaskSEDplusCorrelations(const cha
   // Output slot #1 writes into a TList container
   DefineOutput(1,TList::Class());  //My private output
   // Output slot #2 writes cut to private output
-  //  DefineOutput(2,AliRDHFCutsDplusCorrelationstoKpipi::Class());
   DefineOutput(2,TList::Class());
   // Output slot #3 writes cut to private output
   DefineOutput(3,AliNormalizationCounter::Class());
-  DefineOutput(4,AliHFAssociatedTrackCuts::Class());
+  //DefineOutput(4,AliHFAssociatedTrackCuts::Class());
   
  
 }
@@ -194,12 +170,13 @@ AliAnalysisTaskSEDplusCorrelations::~AliAnalysisTaskSEDplusCorrelations()
   //
   // Destructor
   //
-  delete fOutput;
-  delete fListCuts;
-  delete fRDCutsAnalysis;
-  delete fCuts;
-  delete fCounter;
-  delete fCorrelator; 
+  if(fOutput) {delete fOutput; fOutput = 0;}
+  if(fListCuts) {delete fListCuts; fListCuts = 0;}
+  if(fRDCutsAnalysis) {delete fRDCutsAnalysis; fRDCutsAnalysis = 0;}
+  if(fCuts) {delete fCuts; fCuts = 0;}
+  if(fCounter) {delete fCounter; fCounter = 0;}
+  if(fCorrelator) {delete fCorrelator; fCorrelator=0;}
+
 
 
 }  
@@ -252,8 +229,9 @@ void AliAnalysisTaskSEDplusCorrelations::Init(){
   //
   if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations::Init() \n");
   
-  //PostData(2,fRDCutsloose);//we should then put those cuts in a tlist if we have more than 1
+
   fListCuts=new TList();
+
   // fListCutsAsso=new TList();
   
   AliRDHFCutsDplustoKpipi *analysis = new AliRDHFCutsDplustoKpipi(*fRDCutsAnalysis);
@@ -263,12 +241,12 @@ void AliAnalysisTaskSEDplusCorrelations::Init(){
   //trkcuts->SetName("Assotrkcuts");
   
   fListCuts->Add(analysis);
-  //fListCuts->Add(trkcuts);
+  //fListCutsAsso->Add(trkcuts);
 
    
 
   PostData(2,fListCuts);
-  PostData(4,fCuts);
+  //≈ß PostData(4,fListCutsAsso);
   
   return;
 }
@@ -278,19 +256,31 @@ void AliAnalysisTaskSEDplusCorrelations::UserCreateOutputObjects()
 {
   // Create the output container
   //
+
   if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations::UserCreateOutputObjects() \n");
-  // correlator creation and definition
+  // Correlator creation and definition
 
   Double_t Pi = TMath::Pi();
+
+  // Set up the correlator
+
+
   fCorrelator = new AliHFCorrelator("Correlator",fCuts,fSystem); // fCuts is the hadron cut object, fSystem to switch between pp or PbPb
+
   fCorrelator->SetDeltaPhiInterval((-0.5-1./32)*Pi,(1.5-1./32)*Pi); // set correct phi interval
-  //fCorrelator->SetDeltaPhiInterval(-1.57,4.71);
+
   fCorrelator->SetEventMixing(fMixing); //set kFALSE/kTRUE for mixing Off/On
+
   fCorrelator->SetAssociatedParticleType(fSelect); // set 1/2/3 for hadron/kaons/kzeros
+
   fCorrelator->SetApplyDisplacementCut(fDisplacement); //set kFALSE/kTRUE for using the displacement cut
+  
   fCorrelator->SetUseMC(fReadMC);
   fCorrelator->SetPIDmode(2);
-
+  fCorrelator->SetUseReco(fReco); // to analyse reco/kine
+  // For Event Mixing
   Bool_t pooldef = fCorrelator->DefineEventPool();
   
   if(!pooldef) AliInfo("Warning:: Event pool not defined properly");
@@ -306,53 +296,11 @@ void AliAnalysisTaskSEDplusCorrelations::UserCreateOutputObjects()
   Int_t nbins=GetNBinsHistos();
 
 
-   Int_t nbinsphi = 32;
-   Double_t philow = -0.5*Pi - Pi/32; // shift the bin by half the width so that at 0 is it the bin center
-   Double_t phiup = 1.5*Pi - Pi/32;    
-
-   Int_t nbinseta = 16;
-   Double_t etalow = -1.6; // shift the bin by half the width so that at 0 is it the bin center
-   Double_t etaup = +1.6;      
-     
-
   
   for(Int_t i=0;i<fNPtBins;i++){
 
     index=GetHistoIndex(i);
-
-
-    hisname.Form("hMassVsdPhiHad%d",i);
-    fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaHad%d",i);
-    fMassVsdEtaHistHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiKaon%d",i);
-    fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistKaon[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaKaon%d",i);
-    fMassVsdEtaHistKaon[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistKaon[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiK0%d",i);
-    fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistKshort[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaK0%d",i);
-    fMassVsdEtaHistKshort[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistKshort[index]->Sumw2();
     
-    hisname.Form("hMassVsdPhiLeadHad%d",i);
-    fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistLeadHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaLeadHad%d",i);
-    fMassVsdEtaHistLeadHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistLeadHad[index]->Sumw2();
-
     hisname.Form("hMassPtK0S%d",i);
     fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
     fMassHistK0S[index]->Sumw2();
@@ -361,61 +309,19 @@ void AliAnalysisTaskSEDplusCorrelations::UserCreateOutputObjects()
     fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
     fLeadPt[index]->Sumw2();
 
-       
+          
     hisname.Form("hMassPt%d",i);
     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
     fMassHist[index]->Sumw2();
-    
-   
+
+      
     hisname.Form("hMassPt%dTC",i);
     fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
     fMassHistTC[index]->Sumw2();
 
-    hisname.Form("hMassPt%dTCPlus",i);
-    fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTCPlus[index]->Sumw2();
-
-    hisname.Form("hMassPt%dTCMinus",i);
-    fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTCMinus[index]->Sumw2();
-
-
-    
+       
     index=GetSignalHistoIndex(i); 
-
-    hisname.Form("hMassVsdPhiHadSig%d",i);
-    fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaHadSig%d",i);
-    fMassVsdEtaHistHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiKaonSig%d",i);
-    fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistKaon[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaKaonSig%d",i);
-    fMassVsdEtaHistKaon[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistKaon[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiK0Sig%d",i);
-    fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistKshort[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaK0Sig%d",i);
-    fMassVsdEtaHistKshort[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistKshort[index]->Sumw2();
     
-    hisname.Form("hMassVsdPhiLeadHadSig%d",i);
-    fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistLeadHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaLeadHadSig%d",i);
-    fMassVsdEtaHistLeadHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistLeadHad[index]->Sumw2();
-
     hisname.Form("hSigPtK0S%d",i);
     fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
     fMassHistK0S[index]->Sumw2();
@@ -428,168 +334,106 @@ void AliAnalysisTaskSEDplusCorrelations::UserCreateOutputObjects()
     fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
     fMassHist[index]->Sumw2();
 
-    hisname.Form("hSigPt%dTC",i);
-    fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTC[index]->Sumw2();
-    hisname.Form("hSigPt%dTCPlus",i);
-    fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTCPlus[index]->Sumw2();
-    hisname.Form("hSigPt%dTCMinus",i);
-    fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTCMinus[index]->Sumw2();
+    hisname.Form("hSigOrigCPt%d",i);
+    fMassHistOrigC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
+    fMassHistOrigC[index]->Sumw2();
 
+    hisname.Form("hSigOrigBPt%d",i);
+    fMassHistOrigB[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
+    fMassHistOrigB[index]->Sumw2();
 
+      
+    hisname.Form("hSigMCPt%d",i);
+    fMassHistMC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
+    fMassHistMC[index]->Sumw2();
    
-    index=GetBackgroundHistoIndex(i);
-
-    hisname.Form("hMassVsdPhiBkgHad%d",i);
-    fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaBkgHad%d",i);
-    fMassVsdEtaHistHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiBkgKaon%d",i);
-    fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistKaon[index]->Sumw2();
-
-    hisname.Form("hMassVsdEtaBkgKaon%d",i);
-    fMassVsdEtaHistKaon[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistKaon[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiBkgKshort%d",i);
-    fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistKshort[index]->Sumw2();
-
-
-    hisname.Form("hMassVsdPhiBkgKshort%d",i);
-    fMassVsdEtaHistKshort[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistKshort[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiBkgLeadHad%d",i);
-    fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
-    fMassVsdPhiHistLeadHad[index]->Sumw2();
-
-    hisname.Form("hMassVsdPhiBkgKshort%d",i);
-    fMassVsdEtaHistLeadHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
-    fMassVsdEtaHistLeadHad[index]->Sumw2();
-    
-    hisname.Form("hBkgPtK0S%d",i);
-    fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
-    fMassHistK0S[index]->Sumw2();
-
-    hisname.Form("hLeadBkgPt%d",i);
-    fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
-    fLeadPt[index]->Sumw2();
-
-    hisname.Form("hBkgPt%dTC",i);
-    fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTC[index]->Sumw2();
-
-    hisname.Form("hBkgPt%dTCPlus",i);
-    fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTCPlus[index]->Sumw2();
-
-    hisname.Form("hBkgPt%dTCMinus",i);
-    fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
-    fMassHistTCMinus[index]->Sumw2();
-  }
+     }
     
 
   for(Int_t i=0; i<3*fNPtBins; i++){
-    fOutput->Add(fMassVsdPhiHistHad[i]);
-    fOutput->Add(fMassVsdEtaHistHad[i]);
-    fOutput->Add(fMassVsdPhiHistKaon[i]);
-    fOutput->Add(fMassVsdEtaHistKaon[i]);
-    fOutput->Add(fMassVsdPhiHistKshort[i]);
-    fOutput->Add(fMassVsdEtaHistKshort[i]);
-    fOutput->Add(fMassVsdPhiHistLeadHad[i]);
-    fOutput->Add(fMassVsdEtaHistLeadHad[i]);
+   
     fOutput->Add(fMassHistK0S[i]);
     fOutput->Add(fLeadPt[i]);
     fOutput->Add(fMassHist[i]);
     fOutput->Add(fMassHistTC[i]);
-    fOutput->Add(fMassHistTCPlus[i]);
-    fOutput->Add(fMassHistTCMinus[i]);
-    
-
-
+    fOutput->Add(fMassHistOrigC[i]);
+    fOutput->Add(fMassHistOrigB[i]);
+    fOutput->Add(fMassHistMC[i]);
+   
   }
-  fInvMassK0S = new TH2F("K0S","K0S", 500,0.3,0.8,500,0,50);
-  fInvMassK0S->GetXaxis()->SetTitle("Invariant Mass (#pi #pi) (GeV/c^{2})");
-  fInvMassK0S->GetYaxis()->SetTitle("K0S pt (GeV/c)");
-  fOutput->Add(fInvMassK0S);
-
-  fMCSources = new TH1F("MCSources","Origin of associated particles in MC", 10, -0.5, 9.5);
-  fMCSources->GetXaxis()->SetBinLabel(1,"All ");
-  fMCSources->GetXaxis()->SetBinLabel(2," from Heavy flavour");
-  fMCSources->GetXaxis()->SetBinLabel(3," from c->D");
-  fMCSources->GetXaxis()->SetBinLabel(4," from b->D");
-  fMCSources->GetXaxis()->SetBinLabel(5," from b->B");
-  if(fReadMC) fOutput->Add(fMCSources);
-
-  fK0Origin = new TH1F("K0Origin","Origin of K0 ", 10, -0.5, 5.5);
-  fK0Origin->GetXaxis()->SetBinLabel(1,"All K0s");
-  fK0Origin->GetXaxis()->SetBinLabel(2,"K0s from Heavy flavour");
-  fK0Origin->GetXaxis()->SetBinLabel(3,"K0s from Charm");
-  fK0Origin->GetXaxis()->SetBinLabel(4,"K0s from Beauty");
-  if(fReadMC) fOutput->Add(fK0Origin);
-
-  fKaonOrigin = new TH1F("K0Origin","Origin of Kaon ", 10, -0.5, 5.5);
-  fKaonOrigin->GetXaxis()->SetBinLabel(1,"All Kaons");
-  fKaonOrigin->GetXaxis()->SetBinLabel(2,"Kaons from Heavy flavour");
-  fKaonOrigin->GetXaxis()->SetBinLabel(3,"Kaons from Charm");
-  fKaonOrigin->GetXaxis()->SetBinLabel(4,"Kaons from Beauty");
-  if(fReadMC) fOutput->Add(fKaonOrigin);
   
-       
   
-  fEventMix = new TH2F("EventMixingCheck","EventMixingCheck",5,-0.5,4.5,7,-0.5,6.5);
+  TH1D * EventTypeMC = new TH1D("EventTypeMC","EventTypeMC",100,-0.5,99.5);
+  if(fReadMC) fOutput->Add(EventTypeMC);
+
+       
+  //Event Mixing histos (for control plots)
+
+  Int_t NumberEvents = fCuts->GetMaxNEventsInPool();
+  Int_t NofTracks = fCuts->GetMinNTracksInPool();
+  Int_t NofCentBins = fCuts->GetNCentPoolBins();
+  Double_t * CentBins = fCuts->GetCentPoolBins();
+  Int_t NumberZVtxBins = fCuts->GetNZvtxPoolBins();
+  Double_t *ZVtxBins = fCuts->GetZvtxPoolBins();
+  
+  
+  
+  fEventMix = new TH2D("EventMixingCheck","EventMixingCheck",100,0,600,100,-15,15);
+  
+  fEventMix->GetXaxis()->SetTitle("Multiplicity ");
+  fEventMix->GetYaxis()->SetTitle("Z vertex [cm]");
+  
   if(fMixing)fOutput->Add(fEventMix);
-       
   
-  fHistNEvents = new TH1F("fHistNEvents", "number of events ",11,-0.5,10.5);
+  TH3D * EventsInPool = new TH3D("EventsInPool","Number of events in  pool",NofCentBins,0,600,NumberZVtxBins,-15,15,1000000,0,1000000);
+  
+  EventsInPool->GetXaxis()->SetTitle("Multiplicity ");
+  EventsInPool->GetYaxis()->SetTitle("Z vertex [cm]");
+  EventsInPool->GetZaxis()->SetTitle("Number of events in pool");
+  if(fMixing) fOutput->Add(EventsInPool);
+  
+  Int_t MaxNofTracks = (NumberEvents+1)*NofTracks;
+  
+  
+  TH3D * NofTracksInPool = new TH3D("NofTracksInPool","Number of tracks in  pool",NofCentBins,0,500,NumberZVtxBins,-15,15,MaxNofTracks,0,MaxNofTracks);
+  NofTracksInPool->GetXaxis()->SetTitle("Multiplicity ");
+  NofTracksInPool->GetYaxis()->SetTitle("Z vertex [cm]");
+  NofTracksInPool->GetZaxis()->SetTitle("Number of tracks");
+  
+  if(fMixing) fOutput->Add(NofTracksInPool);
+  
+  TH2D * NofPoolBinCalls = new TH2D("NofPoolBinCalls","Calls per pool bin",NofCentBins,CentBins,NumberZVtxBins,ZVtxBins);
+  NofPoolBinCalls->GetXaxis()->SetTitle("Multiplicity ");
+  NofPoolBinCalls->GetYaxis()->SetTitle("Z vertex [cm]");
+  if(fMixing) fOutput->Add(NofPoolBinCalls);
+  
+       
+  
+  fHistNEvents = new TH1F("fHistNEvents", "number of events ",10,-0.5,10.5);
   fHistNEvents->GetXaxis()->SetBinLabel(1,"nEventsAnal");
   fHistNEvents->GetXaxis()->SetBinLabel(2,"nEvents accepted");
-  fHistNEvents->GetXaxis()->SetBinLabel(3,"Rejected due to trigger");
-  fHistNEvents->GetXaxis()->SetBinLabel(4,"Rejected pileup events");
-  fHistNEvents->GetXaxis()->SetBinLabel(5,"Rejected due to centrality"); 
-  fHistNEvents->GetXaxis()->SetBinLabel(6,"Rejected due to vtxz");
-  fHistNEvents->GetXaxis()->SetBinLabel(7,"Rejected due to Physics Sel");
-  fHistNEvents->GetXaxis()->SetBinLabel(8,"Total no. of candidate");
-  fHistNEvents->GetXaxis()->SetBinLabel(9,"no. of cand wo bitmask");
-  fHistNEvents->GetXaxis()->SetBinLabel(10,"D+ after loose cuts");
-  fHistNEvents->GetXaxis()->SetBinLabel(11,"D+ after tight cuts");
+  fHistNEvents->GetXaxis()->SetBinLabel(3,"nEvent with good vertex");
+  fHistNEvents->GetXaxis()->SetBinLabel(4,"Total number of candidates");
+  fHistNEvents->GetXaxis()->SetBinLabel(5,"Number without bitmask"); 
+  fHistNEvents->GetXaxis()->SetBinLabel(6,"Number after Physics Selection");
+  fHistNEvents->GetXaxis()->SetBinLabel(7,"Total number in Fiducial Accpt");
+  fHistNEvents->GetXaxis()->SetBinLabel(8,"Total no. of good candidates");
+
  
   fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);  
   fHistNEvents->Sumw2();
   fHistNEvents->SetMinimum(0);
   fOutput->Add(fHistNEvents);
-
-  fPtVsMass=new TH2F("hPtVsMass","PtVsMass (prod. cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
-  fPtVsMassTC=new TH2F("hPtVsMassTC","PtVsMass (analysis cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);  
-  fYVsPt=new TH2F("hYVsPt","YvsPt (prod. cuts)",40,0.,20.,80,-2.,2.);
-  fYVsPtTC=new TH2F("hYVsPtTC","YvsPt (analysis cuts)",40,0.,20.,80,-2.,2.);
-  fYVsPtSig=new TH2F("hYVsPtSig","YvsPt (MC, only sig., prod. cuts)",40,0.,20.,80,-2.,2.);
-  fYVsPtSigTC=new TH2F("hYVsPtSigTC","YvsPt (MC, only Sig, analysis cuts)",40,0.,20.,80,-2.,2.);
-
-  fOutput->Add(fPtVsMass);
-  fOutput->Add(fPtVsMassTC);
-  fOutput->Add(fYVsPt);
-  fOutput->Add(fYVsPtTC);
-  fOutput->Add(fYVsPtSig);
-  fOutput->Add(fYVsPtSigTC);
-
-
+  
   // Counter for Normalization
   TString normName="NormalizationCounter";
   AliAnalysisDataContainer *cont = GetOutputSlot(3)->GetContainer();
+  
   if(cont)normName=(TString)cont->GetName();
   fCounter = new AliNormalizationCounter(normName.Data());
   fCounter->Init();
 
-  
+  CreateCorrelationObjs(); 
 
   PostData(1,fOutput);      
   PostData(3,fCounter);      
@@ -603,11 +447,14 @@ void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
   TClonesArray *array3Prong = 0;
 
+    if(fReco) std::cout << "USING RECONSTRUCTION" << std::endl;
+   if(!fReco) std::cout << "USING MC TRUTH" << std::endl;
+  
   if(!fMixing){
     if(fSelect==1) cout << "TASK::Correlation with hadrons on SE "<< endl;
     if(fSelect==2) cout << "TASK::Correlation with kaons on SE "<< endl;
     if(fSelect==3) cout << "TASK::Correlation with kzeros on SE "<< endl;
-  }
+  }  
   if(fMixing){
     if(fSelect==1) cout << "TASK::Correlation with hadrons on ME "<< endl;
     if(fSelect==2) cout << "TASK::Correlation with kaons on ME "<< endl;
@@ -636,8 +483,10 @@ void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
 
   
   // the AODs with null vertex pointer didn't pass the PhysSel
+  
   if(!aod->GetPrimaryVertex()||TMath::Abs(aod->GetMagneticField())<0.001) return;
   fCounter->StoreEvent(aod,fRDCutsAnalysis,fReadMC);
+  
   fHistNEvents->Fill(0); // count event
   
 
@@ -646,20 +495,31 @@ void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
 
   
   PostData(1,fOutput);
+  
   if(!isEvSel)return;
+  
   fHistNEvents->Fill(1);
 
   // set PIDResponse for associated tracks
+  
   fCorrelator->SetPidAssociated(); 
 
   TClonesArray *arrayMC=0;
   AliAODMCHeader *mcHeader=0;
+  
   // AOD primary vertex
   AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
-  //    vtx1->Print();
+  
+  
   TString primTitle = vtx1->GetTitle();
-  //if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0)fHistNEvents->Fill(2);
-   
+  
+  if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0){
+  fHistNEvents->Fill(2);
+  }
+  
+  
+  
   // load MC particles
   if(fReadMC){
     
@@ -675,222 +535,493 @@ void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
       printf("AliAnalysisTaskSEDplusCorrelations::UserExec: MC header branch not found!\n");
       return;
     }
+
+       Bool_t isMCeventgood = kFALSE;
+       
+               
+       Int_t eventType = mcHeader->GetEventType();
+       Int_t NMCevents = fCuts->GetNofMCEventType();
+              
+       
+               for(Int_t k=0; k<NMCevents; k++){
+               Int_t * MCEventType = fCuts->GetMCEventType();
+
+               if(eventType == MCEventType[k]) isMCeventgood= kTRUE;
+
+       }
+       ((TH1D*)fOutput->FindObject("EventTypeMC"))->Fill(eventType);
+               
+               if(NMCevents && !isMCeventgood){
+                 
+               std::cout << "The MC event " << eventType << " not interesting for this analysis: skipping" << std::endl;
+               return; 
+               }
+    
   }
 
   //HFCorrelators initialization (for this event)
 
-  fCorrelator->SetAODEvent(aod); // set the event to be processedfCorrelator->
+  fCorrelator->SetAODEvent(aod); // set the event to be processed
+
   Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing
+
   if(!correlatorON) {
     AliInfo("AliHFCorrelator didn't initialize the pool correctly or processed a bad event");
     return;
   }
   if(fReadMC) fCorrelator->SetMCArray(arrayMC);
 
-
+  Int_t n3Prong = 0;
  
-  Int_t n3Prong = array3Prong->GetEntriesFast();
+  if(fReco) n3Prong = array3Prong->GetEntriesFast();
+
+  
+  if(!fReco) n3Prong = arrayMC->GetEntriesFast();  //for MC kine
+  
+  
    printf("Number of D+->Kpipi: %d and of tracks: %d\n",n3Prong,aod->GetNumberOfTracks());  
-  Int_t nOS=0;
+  
   Int_t index;
   Int_t pdgDgDplustoKpipi[3]={321,211,211};
+  
   Int_t nSelectedloose=0,nSelectedtight=0;
 
+  Double_t ptCand;
+  Double_t phiCand;
+  Double_t etaCand;
+  Double_t invMass = -1;
+  
+  Int_t nIDs[3] = {-9999999};
      
-   
+  Bool_t isDplus = kFALSE;
+  
+  Int_t labDp = -1;
+  
   for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
+
+    AliAODMCParticle* DplusMC;
+
+    if(fReco){
     AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
-    fHistNEvents->Fill(7);
+  
+    fHistNEvents->Fill(3);
+
+
+    
+    if(d->Pt()<2.0) continue;  // Dplus below 2.0 not considered.
+
+    
+
+    //      cout<<multipli<<"    multi"<<endl;
+    
     
+
     if(fUseBit && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts)){
-      fHistNEvents->Fill(8);
+      fHistNEvents->Fill(4);
       continue;
     }
 
     Int_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kAll,aod);
-     
-     if(!fRDCutsAnalysis->GetIsSelectedCuts()) continue;
-      
-    Bool_t unsetvtx=kFALSE;
-    if(!d->GetOwnPrimaryVtx()){
-      d->SetOwnPrimaryVtx(vtx1);
-      unsetvtx=kTRUE;
-    }
-
-      
-    Double_t ptCand = d->Pt();
-    Int_t iPtBin = fRDCutsAnalysis->PtBin(ptCand);
-    Bool_t recVtx=kFALSE;
-    AliAODVertex *origownvtx=0x0; 
-    if(fRDCutsAnalysis->GetIsPrimaryWithoutDaughters()){
-      if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());   
-      if(fRDCutsAnalysis->RecalcOwnPrimaryVtx(d,aod))recVtx=kTRUE;
-      else fRDCutsAnalysis->CleanOwnPrimaryVtx(d,aod,origownvtx);
-    }
-      
-     
-      
-                
-    Int_t labDp=-1;
-    Bool_t isDplus = kFALSE;
-      
+    
+      if(!fRDCutsAnalysis->GetIsSelectedCuts()) continue;
+    
+      fHistNEvents->Fill(5);
+    
+    
     if(fReadMC){
       labDp = d->MatchToMC(411,arrayMC,3,pdgDgDplustoKpipi);
+      
+      //  cout<<labDp<<"***"<<endl;
       if(labDp>=0){
        isDplus = kTRUE;
-      }
-    }
 
-          
-    Double_t invMass=d->InvMassDplus();
+      }
+    }    
+    
+    invMass=d->InvMassDplus();
     Double_t rapid=d->YDplus();
-    fYVsPt->Fill(ptCand,rapid);
-    if(passTightCuts>0.) {fYVsPtTC->Fill(ptCand,rapid);nOS++;}
-    //printf("****************: %d and of tracks: %d\n",n3Prong,nOS);
+    etaCand = d->Eta();
+    phiCand = d->Phi();
+    ptCand = d->Pt();
+    
     Bool_t isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,rapid);
+
     if(isFidAcc){
-      fPtVsMass->Fill(invMass,ptCand);
-      if(passTightCuts>0) fPtVsMassTC->Fill(invMass,ptCand);
+      fHistNEvents->Fill(6);
+      nSelectedloose++;}
+      if(!isFidAcc) continue;         
+      if(!passTightCuts)continue; 
+     
+    fHistNEvents->Fill(7);
+    
+    nSelectedtight++;
+    
+    
+    nIDs[0] = ((AliAODTrack*)d->GetDaughter(0))->GetID();
+    nIDs[1] = ((AliAODTrack*)d->GetDaughter(1))->GetID();
+    nIDs[2] = ((AliAODTrack*)d->GetDaughter(2))->GetID();
+    
     }
-            
+
+    Int_t labDplus=-1;                         
+    
+    if(!fReco){
+           
+      DplusMC = dynamic_cast<AliAODMCParticle*>(arrayMC->At(i3Prong));
+    
+      if (!DplusMC) {
+       AliWarning("Dplus MC Particle not found "); 
+       
+       continue;
+      }
+      labDplus = DplusMC->GetLabel();
+      
+      Int_t PDG =TMath::Abs(DplusMC->PdgCode()); 
+      
+      if(PDG !=411) continue; // not a Dplus
+      ptCand = DplusMC->Pt();
+      phiCand = DplusMC->Phi();
+      etaCand = DplusMC->Eta();
+
+      Bool_t isFidAccMC =fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,DplusMC->Y());
+      
+      if(!isFidAccMC)continue;
            
+    }                          
+    
+    //cout << "PHI D+ = " << phiCand << endl;  
+    // cout << "ETA D+ = " << etaCand << endl;
+    // cout << "PT D+ = " << ptCand << endl;
+    
+   
+    
+    
+        Int_t iPtBin = fRDCutsAnalysis->PtBin(ptCand);                             
+    
     if(iPtBin>=0){
+      
       index=GetHistoIndex(iPtBin);
-      if(isFidAcc){
-       fHistNEvents->Fill(9);
-       nSelectedloose++;
-       fMassHist[index]->Fill(invMass);
+      
+     cout<<"*****"<<iPtBin<<endl;
+      Double_t effTrig;
+      if (fTrig){
        
-       // loop for tight cuts
-               if(passTightCuts>0){   
-         fHistNEvents->Fill(10);
-         nSelectedtight++;
-         fMassHistTC[index]->Fill(invMass);
+       Double_t multipli  = AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aod,-1,1.); 
+        
+       effTrig = fCuts->GetTrigWeight(ptCand,multipli);
 
-         //Dplus info
+       //      cout<<"*****"<<effTrig<<"  "<<multipli<<"   "<<iPtBin<<endl;
+       
+      }
+      else
+       {
+         
+        effTrig = 1.0;
+       
+       }
+     
+      
+      Double_t trigweight = 1/effTrig ;
+
+      // cout<<"****"<<trigweight<<"  ***"<<endl;
+      
+      if(fReco && !fMixing){
+       
+       fMassHist[index]->Fill(invMass); //without trig weight
+       
+       fMassHistTC[index]->Fill(invMass,trigweight); //with trig wt
+      }
+     
+      if( fReco && fReadMC && isDplus && !fMixing) {
+       
+       index=GetSignalHistoIndex(iPtBin);
+       
+        fMassHistMC[index]->Fill(invMass,trigweight);
 
-         Double_t phiDplus = fCorrelator->SetCorrectPhiRange(d->Phi());
-         fCorrelator->SetTriggerParticleProperties(d->Pt(),phiDplus,d->Eta());
+         if(labDp>=0){
+     
+        AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
+        Int_t DplusSource = CheckOrigin(arrayMC,partDp);
+        
+         if(DplusSource == 4){ // is from charm 
            
-         Int_t nIDs[3] = {-9999999};
-         nIDs[0] = ((AliAODTrack*)d->GetDaughter(0))->GetID();
-         nIDs[1] = ((AliAODTrack*)d->GetDaughter(1))->GetID();
-         nIDs[2] = ((AliAODTrack*)d->GetDaughter(2))->GetID();
-
-         Double_t ptlead = 0;
-         Double_t philead = 0;
-         Double_t etalead = 0;
-         Double_t refpt = 0;
+            fMassHistOrigC[index]->Fill(invMass,trigweight);
+         }
+
+         if(DplusSource == 5){ // is from beauty
            
+           fMassHistOrigB[index]->Fill(invMass,trigweight);
 
 
-         Bool_t execPool = fCorrelator->ProcessEventPool();
+         }     
+        }
+      }
+    
+      
+      if(!fReco && fReadMC && !fMixing){
 
-         //     printf("*************: %d\n",execPool);
-         if(fMixing && !execPool) {
-           AliInfo("Mixed event analysis: pool is not ready");
-           continue;
-         }
-         Int_t NofEventsinPool = 1;
-         if(fMixing) {
-           NofEventsinPool = fCorrelator->GetNofEventsInPool();
+       
+       if(labDplus>=0){
+
+         fMassHist[index]->Fill(1.869);
+       
+         AliAODMCParticle *kineDp = (AliAODMCParticle*)arrayMC->At(labDplus);
+         
+         Int_t DplusSource = CheckOrigin(arrayMC,kineDp);
+         
+         
+         if(DplusSource==4){ // is from charm
+           
+            ((TH1F*)fOutput->FindObject(Form("histOrig_Dplus_Bin%d",iPtBin)))->Fill(0.);
          }
-               
-               
-               
-         for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, stops at one
 
-           Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix);
-           if(!analyzetracks) {
-             AliInfo("AliHFCorrelator::Cannot process the track array");
-             continue;
-           }
-                
-           //start the track loop
-                 
-           // Int_t NofTracks = fCorrelator->GetNofTracks();
-
-           //cout<<"*******"<<NofTracks<<endl;
-
-           for (Int_t iTrack = 0;iTrack<fCorrelator->GetNofTracks();iTrack++) {                       
-             Bool_t runcorrelation = fCorrelator->Correlate(iTrack);
-                      
-             if(!runcorrelation) continue;
-             Double_t DeltaPhi = fCorrelator->GetDeltaPhi();
-             Double_t DeltaEta = fCorrelator->GetDeltaEta();
-                      
-             AliReducedParticle* redpart = fCorrelator->GetAssociatedParticle();
-             Double_t phiHad=redpart->Phi();
-             Double_t ptHad=redpart->Pt();
-             Double_t etaHad=redpart->Eta();
-             Int_t label = redpart->GetLabel();
-             Int_t trackid = redpart->GetID();
-             phiHad = fCorrelator->SetCorrectPhiRange(phiHad);
-
-                      
-             //  discard the dplus daughters
-             if (!fMixing){
-               if( trackid == nIDs[0] || trackid == nIDs[1] || trackid == nIDs[2]) continue;
-             }
-             // discard the negative id tracks 
-             if(trackid < 0) continue;
+         if(DplusSource==5){ // is from beauty
 
-                   
-             FillCorrelations(d,DeltaPhi,DeltaEta,index,fSelect);
-                   
-             // For leading particle
-                   
-             if (ptHad > refpt) {
-               refpt = ptHad; ptlead = ptHad;
-               philead = d->Phi() - phiHad;
-               etalead = d->Eta() - etaHad;
-               if (philead < (-1)*TMath::Pi()/2) philead += 2*TMath::Pi();
-               if (philead > 3*TMath::Pi()/2) philead -= 2*TMath::Pi();
-                     
-             }
+         
+           ((TH1F*)fOutput->FindObject(Form("histOrig_Dplus_Bin%d",iPtBin)))->Fill(1.); 
+         }     
+       }
+      }
 
-             // montecarlo
-              
-             if(fReadMC && isDplus) {
+          
+         //Dplus info
 
-               index=GetSignalHistoIndex(iPtBin);
+      Double_t phiDplus = fCorrelator->SetCorrectPhiRange(phiCand);
+      fCorrelator->SetTriggerParticleProperties(ptCand,phiDplus,etaCand);
+           
+         
+      
+      Double_t ptlead = 0;
+      Double_t philead = 0;
+      Double_t etalead = 0;
+      Double_t refpt = 0;
+     
+      Int_t LeadSource = 0;
+      
+      
+      
+      Bool_t execPool = fCorrelator->ProcessEventPool();
+      
+      //printf("*************: %d\n",execPool);
+      
+      if(fMixing && !execPool) {
+       AliInfo("Mixed event analysis: pool is not ready");
+       continue;
+      }
+
+      Int_t NofEventsinPool = 1;
+      if(fMixing) {
+       NofEventsinPool = fCorrelator->GetNofEventsInPool();
 
-               Bool_t* partSource = fCuts->IsMCpartFromHF(label,arrayMC); // check source of associated particle (hadron/kaon/K0)
-               FillMCCorrelations(d,DeltaPhi,DeltaEta,index,partSource,fSelect);   
-               delete partSource;
+       //cout<<"*********"<<NofEventsinPool<<endl;
+      }
+      
+      
+      
+      for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, it stops at 1
+       
+       Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix);
+       if(!analyzetracks) {
+         AliInfo("AliHFCorrelator::Cannot process the track array");
+         continue;
+       }
+       
                
-                        
-             } // readMC 
+       //Int_t NofTracks = fCorrelator->GetNofTracks();
+
+       //cout<<"***number of tracks****"<<fCorrelator->GetNofTracks()<<endl;
+       
+       for (Int_t iTrack = 0;iTrack<fCorrelator->GetNofTracks();iTrack++) {                   
+         Bool_t runcorrelation = fCorrelator->Correlate(iTrack);
+         
+         if(!runcorrelation) continue;
+         
+         Double_t DeltaPhi = fCorrelator->GetDeltaPhi();
+        
+         
+         Double_t DeltaEta = fCorrelator->GetDeltaEta();
+         
+         AliReducedParticle* redpart = fCorrelator->GetAssociatedParticle();
+         Double_t phiHad=redpart->Phi();
+         Double_t ptHad=redpart->Pt();
+         Double_t etaHad=redpart->Eta();
+         Int_t label = redpart->GetLabel();
+         Int_t trackid = redpart->GetID();
+
+         Double_t effi = redpart->GetWeight();
+         Double_t weight = (1/effi)*(trigweight);
+         
+         phiHad = fCorrelator->SetCorrectPhiRange(phiHad);
+         //  cout<<effi<<"******"<<endl;
+
 
-           }//count good tracks
 
-           // For leading particle                               
-           fMassVsdPhiHistLeadHad[index]->Fill(invMass,philead);       
-           fMassVsdEtaHistLeadHad[index]->Fill(invMass,philead,etalead);
+         if (!fMixing && fReco){
+           if( trackid == nIDs[0] || trackid == nIDs[1] || trackid == nIDs[2]) continue;  //discard the Dplus daughters
+         }
+
+         if(fReco && !fReadMC)
+           {
+             FillCorrelations(ptHad,invMass,DeltaPhi,DeltaEta,iPtBin,fSelect,weight);
+           }
+
+        
+           if(fReco && fReadMC && isDplus){
+                   
+            if(label<0) continue;
+            
+            AliAODMCParticle* mchad = (AliAODMCParticle*)arrayMC->At(label);
+            
+            if (!mchad) continue;
+            
+            if (!mchad->IsPhysicalPrimary()) continue; //reject the Reco track if correspondent Kine track is not primary
+            
+            if(labDp < 0) continue;
+            
+            AliAODMCParticle *reDp = (AliAODMCParticle*)arrayMC->At(labDp);
+            
+            Int_t RDSource = CheckOrigin(arrayMC,reDp);
+            
+            Int_t MCSource = CheckTrackOrigin(arrayMC,(AliAODMCParticle*)arrayMC->At(label)) ;// check source of associated particle (hadron/kaon/K0)
 
-           fLeadPt[index]->Fill(ptlead);
             
-           if(fReadMC && isDplus) {
-             index=GetSignalHistoIndex(iPtBin);
-             fMassVsdPhiHistLeadHad[index]->Fill(invMass,philead);
-             fMassVsdEtaHistLeadHad[index]->Fill(invMass,philead,etalead);     
-             fLeadPt[index]->Fill(ptlead);
-                    
+           FillMCRCorrelations(ptHad,invMass,DeltaPhi,DeltaEta,iPtBin,MCSource,RDSource,weight);   
+            
            }
-             
-         }//jmix
+               
+             // montecarlo kine
+             if( !fReco && fReadMC){
                
-         }// tc 
-      }//fid
+               if(label<0) continue;
                
+               if(TMath::Abs(etaHad) > 0.8) continue;
                
+               AliAODMCParticle *hadMC = (AliAODMCParticle*)arrayMC->At(label);
+               if(!hadMC) continue;
+               if (!hadMC->IsPhysicalPrimary()) continue;
+               if(labDplus<0) continue;
+                 
+               AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDplus);
+               if(IsDDaughter(partDp,hadMC,arrayMC)) continue;
+                 
                
+               AliAODMCParticle *trueDp = (AliAODMCParticle*)arrayMC->At(labDplus);
+               Int_t DSource = CheckOrigin(arrayMC,trueDp);
+               
+               Int_t PartSource = CheckTrackOrigin(arrayMC,(AliAODMCParticle*)arrayMC->At(label)) ;// check source of associated particle (hadron/kaon/K0)
                
                
-    }
-    if(recVtx)fRDCutsAnalysis->CleanOwnPrimaryVtx(d,aod,origownvtx);
+               FillMCTruthCorrelations(ptHad,DeltaPhi,DeltaEta,iPtBin,PartSource,DSource,fSelect);   
+             }
+             
+             
+             
+             // For leading particle
+             
+                 if (ptHad > refpt) {
+               //refpt = ptHad; ptlead = ptHad;
+
+               if(fReadMC){
+                if(label<0) continue;
+                if(!(AliAODMCParticle*)arrayMC->At(label))continue;
+                
+                LeadSource = CheckTrackOrigin(arrayMC,(AliAODMCParticle*)arrayMC->At(label)) ;// check source of associated particle (hadron/kaon/K0)
+               }
+
+               philead = phiCand - phiHad;
+               etalead = etaCand - etaHad;
+               if (philead < (-1)*TMath::Pi()/2) philead += 2*TMath::Pi();
+               
+               if (philead > 3*TMath::Pi()/2) philead -= 2*TMath::Pi();
+               
+               refpt = ptHad; ptlead = ptHad;
+             }
+       }
+
+       /*if (fReco){
+        Double_t fillSparseLeadDplus[3] = {philead,invMass,etalead};
+       }
+        else{
+         Double_t fillSparseLeadDplus[3] = {philead,1.87,etalead};
+         }
+
+
+               if(fReco && !fReadMC && !fMixing){
+
+       ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+       }
+
+       if(fReco && fReadMC && isDplus && !fMixing){
+
+       ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+       
+       
+
+       if(LeadSource>=1&&LeadSource<=4){ // is from charm 
+       
+         ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_From_c_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+         
+       }       
+       if(LeadSource>=4&&LeadSource<=8){ // is from beauty
+       
+
+         ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_From_b_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+
+         
+       }       
+       
+    // from NHF
+       if(LeadSource==0){ // is from charm ->D
+         
+
+         ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_From_nhf_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+       }       
+       
+       fLeadPt[index]->Fill(ptlead);
+       }
+       
+       if(!fReco && fReadMC && !fMixing){
+
+         index=GetSignalHistoIndex(iPtBin);
+         
+         
+               ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+         
+         if(LeadSource>=1&&LeadSource<=4){ // is from charm
+ ((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_From_c_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+
+           
+         }     
+         if(LeadSource>=4&&LeadSource<=8){ // is from beauty
+
+((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_From_b_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+
+         }     
+         
+         // from NHF
+         if(LeadSource==0){ // is from charm ->D
+           
+((THnSparseF*)fOutput->FindObject(Form("hPhi_Lead_From_nhf_Bin%d",iPtBin)))->Fill(fillSparseLeadDplus,eweight);
+           
+
+               }
+         fLeadPt[index]->Fill(ptlead);
+         
+         }*/
+       
+         } //jmix
       
-    if(unsetvtx) d->UnsetOwnPrimaryVtx();
+       
+    } //pt bin 
+    
+    
     
   }
 
@@ -906,10 +1037,12 @@ void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
   PostData(1,fOutput); 
   PostData(2,fListCuts);
   PostData(3,fCounter);    
+  //  PostData(4,fListCutsAsso);    
   return;
 }
 
 //________________________________________________________________________
+
 void AliAnalysisTaskSEDplusCorrelations::Terminate(Option_t */*option*/)
 {
   // Terminate analysis
@@ -921,121 +1054,610 @@ void AliAnalysisTaskSEDplusCorrelations::Terminate(Option_t */*option*/)
     printf("ERROR: fOutput not available\n");
     return;
   }
-  fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
-
-  TString hisname;
-  Int_t index=0;
-
-  for(Int_t i=0;i<fNPtBins;i++){
-    index=GetHistoIndex(i);
-    
-    hisname.Form("hMassPt%dTC",i);
-    fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
-  } 
-    
-  TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
-  c1->cd();
-  TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
-  hMassPt->SetLineColor(kBlue);
-  hMassPt->SetXTitle("M[GeV/c^{2}]"); 
-  hMassPt->Draw();
  
   return;
 }
 
 
 //________________________________________________________________________
-void AliAnalysisTaskSEDplusCorrelations::FillCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t sel) const{
-  //Filling histogams
+void AliAnalysisTaskSEDplusCorrelations::FillCorrelations(Double_t ptTrack, Double_t mass, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t sel, Double_t eweight) const{
   
-  Double_t invMass=d->InvMassDplus();
-       
 
+  //Filling THnSparse for correlations , only for charged tracks
+  
   if(sel==1){          
-    fMassVsdPhiHistHad[ind]->Fill(invMass,deltaPhi);
-    fMassVsdEtaHistHad[ind]->Fill(invMass,deltaPhi,deltaEta);
+      
+    if(!fMixing){
+      
+    Double_t ptLim_Sparse = ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->GetAxis(2)->GetXmax(); //all plots have same axes...
+    Double_t EtaLim_Sparse = ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->GetAxis(3)->GetXmax();
+    
+    if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;   if(deltaEta > EtaLim_Sparse) deltaEta = EtaLim_Sparse-0.01;
+    if(deltaEta < -EtaLim_Sparse) deltaEta = -EtaLim_Sparse+0.01;
+  
+    //variables for filling histos
+  
+    Double_t fillSparseDplus[4] = {deltaPhi,mass,ptTrack,deltaEta};
+   
+   
+
+   
+      //sparse fill for data (tracks) + weighted
+    
+    ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->Fill(fillSparseDplus,eweight);
+          
+    }
+    
+    if(fMixing) {
+    
+
+    //variables for filling histos
+    
+      Double_t fillSparseDplusMix[3] = {deltaPhi,mass,deltaEta};
+
+
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ind)))->Fill(fillSparseDplusMix,eweight);
+    }
+
   }
+
   if(sel==2){
-    fMassVsdPhiHistKaon[ind]->Fill(invMass,deltaPhi);
-    fMassVsdEtaHistKaon[ind]->Fill(invMass,deltaPhi,deltaEta);
+
   }
   if(sel==3){
-    fMassVsdPhiHistKshort[ind]->Fill(invMass,deltaPhi);
-    fMassVsdEtaHistKshort[ind]->Fill(invMass,deltaPhi,deltaEta);
+
   }
-        
+  
   return;
 }
 
 
 
 //________________________________________________________________________
-void AliAnalysisTaskSEDplusCorrelations::FillMCCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind,Bool_t* mcSource, Int_t sel) const{
-  // Filling histos with MC information
 
-  Double_t invMass=d->InvMassDplus();
+
+
+void AliAnalysisTaskSEDplusCorrelations::FillMCTruthCorrelations( Double_t ptTrack,Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t mcSource, Int_t origDplus, Int_t sel) const{
+
+   // Filling histos with MC Kine information
+
+     Double_t invm = 1.876;
+     
+   Double_t ptLim_Sparse = ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->GetAxis(2)->GetXmax(); //all plots have same axes...
+   Double_t EtaLim_Sparse = ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->GetAxis(3)->GetXmax();
+    if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
+    if(deltaEta > EtaLim_Sparse) deltaEta = EtaLim_Sparse-0.01;
+    if(deltaEta < -EtaLim_Sparse) deltaEta = -EtaLim_Sparse+0.01;
+
+    
   
+    //variables for filling histos
+     
+    Double_t fillSparseDplus[4] = {deltaPhi,invm,ptTrack,deltaEta};
+    
 
   if(sel==1){
-    fMassVsdPhiHistHad[ind]->Fill(invMass,deltaPhi);
-    fMassVsdEtaHistHad[ind]->Fill(invMass,deltaPhi,deltaEta);
+    
+    if(!fMixing){
+    
+    ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->Fill(fillSparseDplus);
+    
+    if(origDplus==4&&mcSource>=1&&mcSource<=3){ // is from charm 
+      
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_HFC_Bin%d",ind)))->Fill(fillSparseDplus);
+    }
+    if(origDplus==5&&mcSource>=4&&mcSource<=8){ // is from beauty 
+      
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_HFB_Bin%d",ind)))->Fill(fillSparseDplus);
+    }
+    
+    if(mcSource==0){ // is from  NHF
+      
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_NHF_Bin%d",ind)))->Fill(fillSparseDplus);
+    }
+    
+    }
+    
+    if(fMixing) {
 
-       
-    fMCSources->Fill(0);
+
+    //variables for filling histos
+      Double_t fillSparseDplusMix[3] = {deltaPhi,invm,deltaEta};
+      
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ind)))->Fill(fillSparseDplusMix);
+    }  
+    
+  }
+  
+  if(sel==2){
+    
+  }
+  
+  if(sel==3){
+    
+  }
+  
+  
+  return;
+}
+
+
+void AliAnalysisTaskSEDplusCorrelations::FillMCRCorrelations(Double_t ptTrack, Double_t mass, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t mcSource, Int_t origDplus, Double_t eweight) const{
+  
+  // Filling histos with MC information
+
+    
+   Double_t ptLim_Sparse = ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->GetAxis(2)->GetXmax(); //all plots have same axes...
+  Double_t EtaLim_Sparse = ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->GetAxis(3)->GetXmax();
+    if(ptTrack > ptLim_Sparse) ptTrack = ptLim_Sparse-0.01;
+     if(deltaEta > EtaLim_Sparse) deltaEta = EtaLim_Sparse-0.01;
+    if(deltaEta < -EtaLim_Sparse) deltaEta = -EtaLim_Sparse+0.01;
+
+
+  if(!fMixing){
+    Double_t fillSparseDplus[4] = {deltaPhi,mass,ptTrack,deltaEta};
+    
+
+    
+    ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d",ind)))->Fill(fillSparseDplus,eweight);
+  
+     if(origDplus==4&&mcSource>=1&&mcSource<=3){ // is from charm 
+
+    ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_HFC_Bin%d",ind)))->Fill(fillSparseDplus,eweight);
+     }
+     if(origDplus==5&&mcSource>=4&&mcSource<=8){ // is from beauty 
+
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_HFB_Bin%d",ind)))->Fill(fillSparseDplus,eweight);
+     }
+
+     if(mcSource==0){ // is from  NHF
+
+      ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_NHF_Bin%d",ind)))->Fill(fillSparseDplus,eweight);
+     }
+  }
+         
+     if(fMixing) {
+    
+    //variables for filling histos
+    Double_t fillSparseDplusMix[3] = {deltaPhi,mass,deltaEta};
+
+    ((THnSparseF*)fOutput->FindObject(Form("hPhi_Charg_Bin%d_EvMix",ind)))->Fill(fillSparseDplusMix,eweight);
+     }
+
+
+  return;
+}
+
+Bool_t AliAnalysisTaskSEDplusCorrelations::IsDDaughter(AliAODMCParticle* d, AliAODMCParticle* track, TClonesArray* mcArray) const {
+
+  //Daughter removal in MCKine case
+  Bool_t isDaughter = kFALSE;
+  Int_t labelDplus = d->GetLabel();
+  
+  Int_t mother = track->GetMother();
+
+  
+  while (mother > 0){
+    AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(mcArray->At(mother)); 
+    if (mcMoth){
+      if (mcMoth->GetLabel() == labelDplus) isDaughter = kTRUE;
+      mother = mcMoth->GetMother(); //goes back to last generation
+    } else{
+      AliError(" mother particle not found!");
+      break;
+    }
+  }
+
+  return isDaughter;
+}
+
+
+Int_t AliAnalysisTaskSEDplusCorrelations::CheckTrackOrigin(TClonesArray* arrayMC, AliAODMCParticle *mcPartCandidate) const {           
+  //
+  // checks on particle (#) origin:
+  // 0) Not HF
+  // 1) D->#
+  // 2) D->X->#
+  // 3) c hadronization
+  // 4) B->#
+  // 5) B->X-># (X!=D)
+  // 6) B->D->#
+  // 7) B->D->X->#
+  // 8) b hadronization
+  //
+  if(fDebug>2) printf("AliAnalysisTaskSEDplusCorrelations::CheckTrkOrigin() \n");
        
-    if(mcSource[2]&&mcSource[0]){ // is from charm ->D
-      fMCSources->Fill(1);
-      fMCSources->Fill(2);
+  Int_t pdgGranma = 0;
+  Int_t mother = 0;
+  mother = mcPartCandidate->GetMother();
+  Int_t istep = 0;
+  Int_t abspdgGranma =0;
+  Bool_t isFromB=kFALSE;
+  Bool_t isDdaugh=kFALSE;
+  Bool_t isDchaindaugh=kFALSE;
+  Bool_t isBdaugh=kFALSE;
+  Bool_t isBchaindaugh=kFALSE;
+  Bool_t isQuarkFound=kFALSE;
+
+  while (mother > 0){
+    istep++;
+    AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
+    if (mcMoth){
+      pdgGranma = mcMoth->GetPdgCode();
+      abspdgGranma = TMath::Abs(pdgGranma);
+      if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
+       isBchaindaugh=kTRUE;
+        if(istep==1) isBdaugh=kTRUE;
+      }
+      if ((abspdgGranma > 400 && abspdgGranma < 500) || (abspdgGranma > 4000 && abspdgGranma < 5000)){
+       isDchaindaugh=kTRUE;
+        if(istep==1) isDdaugh=kTRUE;
+      }
+      if(abspdgGranma==4 || abspdgGranma==5) {isQuarkFound=kTRUE; if(abspdgGranma==5) isFromB = kTRUE;}
+      mother = mcMoth->GetMother();
+    }else{
+      AliError("Failed casting the mother particle!");
+      break;
     }
-       if(mcSource[2]&&mcSource[1]){ // is from beauty -> D
-      fMCSources->Fill(1);
-      fMCSources->Fill(3);
-    }  
-    if(mcSource[3]&&mcSource[1]){ // is from beauty ->B
-      fMCSources->Fill(1);
-      fMCSources->Fill(4);
+  }
+
+  //decides what to return based on the flag status
+  if(isQuarkFound) {
+    if(!isFromB) {  //charm
+      if(isDdaugh) return 1; //charm immediate
+      else if(isDchaindaugh) return 2; //charm chain
+      else return 3; //charm hadronization
+    }
+    else { //beauty
+      if(isBdaugh) return 4; //b immediate
+      else if(isBchaindaugh) { //b chain
+        if(isDchaindaugh) {
+          if(isDdaugh) return 6; //d immediate
+          return 7; //d chain
+          }
+        else return 5; //b, not d
+      }
+      else return 8; //b hadronization
     }
   }
+  else return 0; //no HF quark
+}
 
-  if(sel==2){
-    fMassVsdPhiHistKaon[ind]->Fill(invMass,deltaPhi);
-    fMassVsdEtaHistKaon[ind]->Fill(invMass,deltaPhi,deltaEta);
-    fKaonOrigin->Fill(0);
-    if(mcSource[2]&&mcSource[0]){ // is from charm ->D
-      fKaonOrigin->Fill(1);
-      fKaonOrigin->Fill(2);
-    }  
-    if(mcSource[2]&&mcSource[1]){ // is from beauty -> D
-      fKaonOrigin->Fill(1);
-      fKaonOrigin->Fill(3);
-    }  
-    if(mcSource[3]&&mcSource[1]){ // is from beauty ->B
-      fKaonOrigin->Fill(1);
-      fKaonOrigin->Fill(4);
+Int_t AliAnalysisTaskSEDplusCorrelations::CheckOrigin(TClonesArray* arrayMC,  AliAODMCParticle *mcPartCandidate) const {               
+  //
+  // checking whether the mother of the particles come from a charm or a bottom quark
+  //
+       
+  Int_t pdgGranma = 0;
+  Int_t mother = 0;
+  mother = mcPartCandidate->GetMother();
+  Int_t istep = 0;
+  Int_t abspdgGranma =0;
+  Bool_t isFromB=kFALSE;
+  Bool_t isQuarkFound=kFALSE;
+  while (mother >0 ){
+    istep++;
+    AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
+    if (mcGranma){
+      pdgGranma = mcGranma->GetPdgCode();
+      abspdgGranma = TMath::Abs(pdgGranma);
+      if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
+       isFromB=kTRUE;
+      }
+      if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
+      mother = mcGranma->GetMother();
+    }else{
+      AliError("Failed casting the mother particle!");
+      break;
     }
   }
-  if(sel==3){
-    fMassVsdPhiHistKshort[ind]->Fill(invMass,deltaPhi);
-    fMassVsdEtaHistKshort[ind]->Fill(invMass,deltaPhi,deltaEta);
-    fK0Origin->Fill(0);
-    if(mcSource[2]&&mcSource[0]){ // is from charm ->D
-      fK0Origin->Fill(1);
-      fK0Origin->Fill(2);
-    }  
-    if(mcSource[2]&&mcSource[1]){ // is from beauty -> D
-      fK0Origin->Fill(1);
-      fK0Origin->Fill(3);
+  
+  if(isFromB) return 5;
+  else return 4;
+}
+
+
+void AliAnalysisTaskSEDplusCorrelations::CreateCorrelationObjs() {
+//
+
+  TString namePlot = "";
+  Int_t nbinsmass=GetNBinsHistos();
+  //These for limits in THnSparse (one per bin, same limits). 
+  //Vars: DeltaPhi, InvMass, PtTrack, Displacement, DeltaEta
+  Int_t nBinsPhi[4] = {32,nbinsmass,30,16};
+  Double_t binMinPhi[4] = {-TMath::Pi()/2.+TMath::Pi()/32.,fLowmasslimit,0.,-1.6};  //is the minimum for all the bins
+  Double_t binMaxPhi[4] = {3*TMath::Pi()/2.+TMath::Pi()/32.,fUpmasslimit,30.0,1.6};  //is the maximum for all the bins
+
+  //Vars: DeltaPhi, InvMass, DeltaEta
+  Int_t nBinsMix[3] = {32,nbinsmass,16};
+  Double_t binMinMix[3] = {-TMath::Pi()/2.+TMath::Pi()/32.,fLowmasslimit,-1.6};  //is the minimum for all the bins
+  Double_t binMaxMix[3] = {3*TMath::Pi()/2.+TMath::Pi()/32.,fUpmasslimit,1.6};  //is the maximum for all the bins
+
+  for(Int_t i=0;i<fNPtBins;i++){
+
+    if(!fMixing) {
+      //THnSparse plots: correlations for various invariant mass (MC and data)
+      namePlot="hPhi_Charg_Bin";
+      namePlot+=i;
+
+      THnSparseF *hPhi = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+      hPhi->Sumw2();
+      fOutput->Add(hPhi);
+
+      namePlot="hPhi_Kaon_Bin";
+      namePlot+=i;
+
+      THnSparseF *hPhiK = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+      hPhiK->Sumw2();
+      fOutput->Add(hPhiK);
+
+      namePlot="hPhi_K0_Bin";
+      namePlot+=i;
+
+      THnSparseF *hPhiK0 = new THnSparseF(namePlot.Data(), "Azimuthal correlation; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+      hPhiK0->Sumw2();
+      fOutput->Add(hPhiK0);
+  
+      //histos for c/b origin for D+ (MC only)
+      if (fReadMC) {
+
+        //generic origin for tracks
+        namePlot="hPhi_Charg_HFC_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiH_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiH_c->Sumw2();
+        fOutput->Add(hPhiH_c);
+
+        namePlot="hPhi_Kaon_HFC_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiK_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiK_c->Sumw2();
+        fOutput->Add(hPhiK_c);
+
+        namePlot="hPhi_K0_HFC_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiK0_c = new THnSparseF(namePlot.Data(), "Azimuthal correlation - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiK0_c->Sumw2();
+        fOutput->Add(hPhiK0_c);
+  
+        namePlot="hPhi_Charg_HFB_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiH_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiH_b->Sumw2();
+        fOutput->Add(hPhiH_b);
+
+        namePlot="hPhi_Kaon_HFB_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiK_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiK_b->Sumw2();
+        fOutput->Add(hPhiK_b);
+
+        namePlot="hPhi_K0_HFB_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiK0_b = new THnSparseF(namePlot.Data(), "Azimuthal correlation - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiK0_b->Sumw2();
+        fOutput->Add(hPhiK0_b);
+
+       namePlot="hPhi_Charg_NHF_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiH_NHF = new THnSparseF(namePlot.Data(), "Azimuthal correlation - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiH_NHF->Sumw2();
+        fOutput->Add(hPhiH_NHF);
+
+        namePlot="hPhi_Kaon_NHF_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiK_Non = new THnSparseF(namePlot.Data(), "Azimuthal correlation - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiK_Non->Sumw2();
+        fOutput->Add(hPhiK_Non);
+
+        namePlot="hPhi_K0_NHF_Bin";
+        namePlot+=i;
+
+        THnSparseF *hPhiK0_Non = new THnSparseF(namePlot.Data(), "Azimuthal correlation - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",4,nBinsPhi,binMinPhi,binMaxPhi);
+        hPhiK0_Non->Sumw2();
+        fOutput->Add(hPhiK0_Non);
+      }
+
+      //leading hadron correlations
+      namePlot="hPhi_Lead_Bin";
+      namePlot+=i;
+
+      THnSparseF *hCorrLead = new THnSparseF(namePlot.Data(), "Leading particle correlations; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+      hCorrLead->Sumw2();
+      fOutput->Add(hCorrLead);
+
+      if (fReadMC) {
+        namePlot="hPhi_Lead_From_c_Bin";
+        namePlot+=i;
+
+        THnSparseF *hCorrLead_c = new THnSparseF(namePlot.Data(), "Leading particle correlations - c origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+        hCorrLead_c->Sumw2();
+        fOutput->Add(hCorrLead_c);
+  
+        namePlot="hPhi_Lead_From_b_Bin";
+        namePlot+=i;
+  
+        THnSparseF *hCorrLead_b = new THnSparseF(namePlot.Data(), "Leading particle correlations - b origin; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+        hCorrLead_b->Sumw2();
+        fOutput->Add(hCorrLead_b);
+  
+
+
+        namePlot="hPhi_Lead_From_nhf_Bin";
+        namePlot+=i;
+  
+        THnSparseF *hCorrLead_Non = new THnSparseF(namePlot.Data(), "Leading particle correlations - Non HF; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+        hCorrLead_Non->Sumw2();
+        fOutput->Add(hCorrLead_Non);
+      }
+      
+           
+
+    //pT distribution histos
+    namePlot = "hist_Pt_Charg_Bin"; namePlot+=i;
+    TH1F *hPtHad = new TH1F(namePlot.Data(), "Charged track pT (in D+ events); p_{T} (GeV/c)",240,0.,20.);
+    hPtHad->SetMinimum(0);
+    fOutput->Add(hPtHad);
+
+    namePlot = "hist_Pt_Kcharg_Bin"; namePlot+=i;
+    TH1F *hPtH = new TH1F(namePlot.Data(), "Hadrons pT (in D+ events); p_{T} (GeV/c)",240,0.,12.);
+    hPtH->SetMinimum(0);
+    fOutput->Add(hPtH);
+
+    namePlot = "hist_Pt_K0_Bin"; namePlot+=i;
+    TH1F *hPtK = new TH1F(namePlot.Data(), "Kaons pT (in D+ events); p_{T} (GeV/c)",240,0.,12.);
+    hPtK->SetMinimum(0);
+    fOutput->Add(hPtK);
+
+     
     }
-    if(mcSource[3]&&mcSource[1]){ // is from beauty ->B
-      fK0Origin->Fill(1);
-      fK0Origin->Fill(4);
+
+    if(fMixing) {
+      //THnSparse plots for event mixing!
+
+
+      namePlot="hPhi_K0_Bin";
+      namePlot+=i;namePlot+="_EvMix";
+
+      THnSparseF *hPhiK_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+      hPhiK_EvMix->Sumw2();
+      fOutput->Add(hPhiK_EvMix);
+
+      namePlot="hPhi_Kcharg_Bin";
+      namePlot+=i;namePlot+="_EvMix";
+  
+      THnSparseF *hPhiH_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+      hPhiK_EvMix->Sumw2();
+      fOutput->Add(hPhiH_EvMix);
+
+      namePlot="hPhi_Charg_Bin";
+      namePlot+=i;namePlot+="_EvMix";
+
+      THnSparseF *hPhiC_EvMix = new THnSparseF(namePlot.Data(), "Az. corr. EvMix; #Delta#phi; Inv. Mass (GeV/c^{2}); p_{t} (GeV/c)",3,nBinsMix,binMinMix,binMaxMix);
+      hPhiC_EvMix->Sumw2();
+      fOutput->Add(hPhiC_EvMix);
     }
-         
   }
 
-  return;
-}
+  //out of bin loop
+  if(!fMixing) {
+    TH1F *hCountC = new TH1F("hist_Count_Charg", "Charged track counter; # Tracks",100,0.,100.);
+    hCountC->SetMinimum(0);
+    fOutput->Add(hCountC);
+
+    TH1F *hCountH = new TH1F("hist_Count_Kcharg", "Hadrons counter; # Tracks",20,0.,20.);
+    hCountH->SetMinimum(0);
+    fOutput->Add(hCountH);
+
+    TH1F *hCountK = new TH1F("hist_Count_K0", "Kaons counter; # Tracks",20,0.,20.);
+    hCountK->SetMinimum(0);
+    fOutput->Add(hCountK);
+  }
+
+  if (fReadMC) {
+    TH1D *hEventTypeMC = new TH1D("EventTypeMC","EventTypeMC",100,-0.5,99.5);
+    fOutput->Add(hEventTypeMC); 
+  }
+
+  // if (fFillGlobal) { //all-events plots
+    //pt distributions
 
+   if(!fMixing) {
+    TH1F *hPtHAll = new TH1F("hist_Pt_Charg_AllEv", "Charged track pT (All); p_{T} (GeV/c)",240,0.,12.);
+    hPtHAll->SetMinimum(0);
+    fOutput->Add(hPtHAll);
+
+    TH1F *hPtKAll = new TH1F("hist_Pt_Kaons_AllEv", "Kaons pT (All); p_{T} (GeV/c)",240,0.,12.);
+    hPtKAll->SetMinimum(0);
+    fOutput->Add(hPtKAll);
+
+    TH1F *hPtK0All = new TH1F("hist_Pt_K0_AllEv", "K0 pT (All); p_{T} (GeV/c)",240,0.,12.);
+    hPtK0All->SetMinimum(0);
+    fOutput->Add(hPtK0All);
+
+   
+      //phi distributions
+      TH1F *hPhiDistHAll = new TH1F("hist_PhiDistr_Charg", "Charged track phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
+      hPhiDistHAll->SetMinimum(0);
+      fOutput->Add(hPhiDistHAll);
+
+      TH1F *hPhiDistKAll = new TH1F("hist_PhiDistr_Kcharg", "Kaons phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
+      hPhiDistKAll->SetMinimum(0);
+      fOutput->Add(hPhiDistKAll);
+
+      TH1F *hPhiDistK0All = new TH1F("hist_PhiDistr_K0", "K0 phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
+      hPhiDistK0All->SetMinimum(0);
+      fOutput->Add(hPhiDistK0All);
+
+      TH1F *hPhiDistDAll = new TH1F("hist_PhiDistr_Dplus", "D^{+} phi distr. (All); p_{T} (GeV/c)",64,0,6.283);
+      hPhiDistDAll->SetMinimum(0);
+      fOutput->Add(hPhiDistDAll);
+   }
+    //}
+
+  //for MC analysis only
+  for(Int_t i=0;i<fNPtBins;i++) {
+
+    if (fReadMC && !fMixing) {
+
+      
+      //origin of tracks histos
+      namePlot="histOrig_Had_Bin";  namePlot+=i;
+      TH1F *hOrigin_had = new TH1F(namePlot.Data(), "Origin of associated tracksin MC",10,-0.5,9.5);
+      hOrigin_had->SetMinimum(0);
+      
+      hOrigin_had->GetXaxis()->SetBinLabel(1,"All ");
+      hOrigin_had->GetXaxis()->SetBinLabel(2," from Heavy flavour");
+      hOrigin_had->GetXaxis()->SetBinLabel(3," from c->D");
+      hOrigin_had->GetXaxis()->SetBinLabel(4," from b->D");
+      hOrigin_had->GetXaxis()->SetBinLabel(5," from b->B");
+      hOrigin_had->GetXaxis()->SetBinLabel(6," from NHF");
+      if(fReadMC) fOutput->Add(hOrigin_had);
+
+      namePlot="histOrig_Kaon_Bin";  namePlot+=i;
+      TH1F *hOrigin_kaon = new TH1F(namePlot.Data(), "Origin of kaons in MC",10,-0.5,9.5);
+      hOrigin_kaon->SetMinimum(0);
+      hOrigin_kaon->GetXaxis()->SetBinLabel(1,"All Kaons");
+      hOrigin_kaon->GetXaxis()->SetBinLabel(2,"Kaons from Heavy flavour");
+      hOrigin_kaon->GetXaxis()->SetBinLabel(3,"Kaons from Charm");
+      hOrigin_kaon->GetXaxis()->SetBinLabel(4,"Kaons from Beauty");
+      hOrigin_kaon->GetXaxis()->SetBinLabel(5,"Kaons from NHF");
+      if(fReadMC) fOutput->Add(hOrigin_kaon);
+      
+
+      namePlot="histOrig_K0_Bin";  namePlot+=i;
+      TH1F *hOrigin_K0 = new TH1F(namePlot.Data(), "Origin of Kshorts in MC",10,-0.5,9.5);
+      hOrigin_K0->SetMinimum(0);
+
+      hOrigin_K0->GetXaxis()->SetBinLabel(1,"All K0s");
+      hOrigin_K0->GetXaxis()->SetBinLabel(2,"K0s from Heavy flavour");
+      hOrigin_K0->GetXaxis()->SetBinLabel(3,"K0s from Charm");
+      hOrigin_K0->GetXaxis()->SetBinLabel(4,"K0s from Beauty");
+      hOrigin_K0->GetXaxis()->SetBinLabel(5,"K0s from NHF");
+      if(fReadMC) fOutput->Add(hOrigin_K0);
+    }
+
+    if (fReadMC) {
+      //origin of Dplus histos
+      namePlot="histOrig_Dplus_Bin";  namePlot+=i;
+      TH1F *hOrigin_Dplus = new TH1F(namePlot.Data(),"Origin of D+ in MC",2,0.,2.);
+      hOrigin_Dplus->SetMinimum(0);
+      hOrigin_Dplus->GetXaxis()->SetBinLabel(1,"From c");
+      hOrigin_Dplus->GetXaxis()->SetBinLabel(2,"From b");
+      fOutput->Add(hOrigin_Dplus);
+    }
+  }
+}
 
 
 
index 0a718b93af3431bf28052e788c1c3e64926e7b7f..cf8f165aaa69db41610c4021d0ca5f6cfca7ff3e 100644 (file)
@@ -54,21 +54,37 @@ class AliAnalysisTaskSEDplusCorrelations : public AliAnalysisTaskSE
   void SetReadMC(Bool_t readMC=kTRUE){fReadMC=readMC;}
   void SetEventMix(Bool_t mixing){fMixing=mixing;}
   
-  // void SetUseStrangeness(Bool_t uses=kTRUE){fUseStrangeness=uses;}
+  void CreateCorrelationObjs();
+
+  //  standard way used in Dplus task  
   void SetMassLimits(Float_t range);
   void SetMassLimits(Float_t lowlimit, Float_t uplimit);
   void SetBinWidth(Float_t w);
+
   void SetUseBit(Bool_t dols=kTRUE){fUseBit=dols;}
+
   void SetCorrelator(Int_t l) {fSelect = l;} // select 1 for hadrons, 2 for Kaons, 3 for Kzeros
   void SetUseDisplacement(Int_t m) {fDisplacement=m;} // select 0 for no displ, 1 for abs displ, 2 for d0/sigma_d0
   void SetSystem(Bool_t system){fSystem=system;} // select between pp (kFALSE) or PbPb (kTRUE)
 
+  void SetUseReconstruction(Bool_t reco){fReco = reco;}
+
+  void SetTrigEfficiency(Bool_t trigeff = kTRUE) {fTrig = trigeff;}
+
+  void FillCorrelations(Double_t ptTrack,Double_t mass, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t sel, Double_t eweight) const;
+
+  void FillMCTruthCorrelations(Double_t ptTrack, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t mcSource, Int_t origDplus, Int_t sel) const;
+  void FillMCRCorrelations(Double_t ptTrack,Double_t mass, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t mcSource, Int_t origDplus,Double_t eweight) const;
+  Bool_t IsDDaughter(AliAODMCParticle* d, AliAODMCParticle* track, TClonesArray* mcArray) const;
   
-  void FillCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t sel) const;
-  void FillMCCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Bool_t* mcSource,Int_t sel) const;
+  Int_t CheckTrackOrigin(TClonesArray* arrayMC, AliAODMCParticle* mcPartCandidate) const ;             
   
+  Int_t CheckOrigin(TClonesArray* arrayMC, AliAODMCParticle* mcPartCandidate) const ;          
   
 
+
+
   
   
   Float_t GetUpperMassLimit(){return fUpmasslimit;}
@@ -78,6 +94,7 @@ class AliAnalysisTaskSEDplusCorrelations : public AliAnalysisTaskSE
   Int_t GetNBinsHistos();
   
     // Implementation of interface methods
+  
   virtual void UserCreateOutputObjects();
   virtual void Init();
   virtual void LocalInit() {Init();}
@@ -89,70 +106,71 @@ class AliAnalysisTaskSEDplusCorrelations : public AliAnalysisTaskSE
   AliAnalysisTaskSEDplusCorrelations(const AliAnalysisTaskSEDplusCorrelations &source);
   AliAnalysisTaskSEDplusCorrelations& operator=(const AliAnalysisTaskSEDplusCorrelations& source); 
 
-  // AliEventPoolManager*     fPoolMgr;         //! event pool manager
-  // TObject* fHandler; //! Analysis Handler
-  //AliHFCorrelator * fCorrelator; // object for correlations
+  
 
   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;}
-  Float_t GetStrangenessWeights(const AliAODRecoDecayHF3Prong* d, TClonesArray* arrayMC, Float_t factor[3]) const;
-  Float_t GetTrueImpactParameter(AliAODMCHeader *mcHeader, TClonesArray* arrayMC, AliAODMCParticle *partD0) const ;
+   
   enum {kMaxPtBins=20};
+  
+  //For dplus efficiency
 
   TList  *fOutput; //! list send on output slot 0
-  //  TObject* fHandler; //! Analysis Handler
-  //AliEventPoolManager*     fPoolMgr;         //! event pool manager
   AliHFCorrelator* fCorrelator; // object for correlations
-  Int_t fSelect; // select what to correlate with a D* 1-chargedtracks,2-chargedkaons,3-k0s
+  Int_t fSelect; // select what to correlate with a Dplus 1-chargedtracks,2-chargedkaons,3-k0s
   Int_t fDisplacement; // set 0 for no displacement cut, 1 for absolute d0, 2 for d0/sigma_d0
   TH1F *fHistNEvents; //!hist. for No. of events
-  TH1F *fMCSources; //!hist. for No. of events
-  TH1F *fK0Origin; //!hist. for No. of events
-  TH1F *fKaonOrigin; //!hist. for No. of events
-  TH2F *fInvMassK0S; //!hist. for D- inv mass (TC)
-  TH2F *fEventMix; //!hist. for event mixing
-  TH2F *fMassVsdPhiHistHad[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH3D *fMassVsdEtaHistHad[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH2F *fMassVsdPhiHistKaon[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH3D *fMassVsdEtaHistKaon[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH2F *fMassVsdPhiHistKshort[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH3D *fMassVsdEtaHistKshort[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH2F *fMassVsdPhiHistLeadHad[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH3D *fMassVsdEtaHistLeadHad[3*kMaxPtBins]; //!hist. for inv mass (LC)
-  TH1F *fMassHistK0S[3*kMaxPtBins]; //!hist. for inv mass (LC)
   
+  Bool_t fTrig;  // flag for using trig eff         
 
+  TH2D *fEventMix; //!hist. for event mixing
+   
+  TH1F *fMassHistK0S[3*kMaxPtBins]; //!hist. for inv mass (LC)
+  
   TH1F *fLeadPt[3*kMaxPtBins]; //!hist. for D- inv mass (TC)
+
+  TH1F *fPtSig[3*kMaxPtBins]; //!hist. for D- inv mass (TC)
+  
   TH1F *fMassHist[3*kMaxPtBins]; //!hist. for inv mass (LC)
   TH1F *fMassHistTC[3*kMaxPtBins]; //!hist. for inv mass (TC)
-  TH1F *fMassHistTCPlus[3*kMaxPtBins]; //!hist. for D+ inv mass (TC)
-  TH1F *fMassHistTCMinus[3*kMaxPtBins]; //!hist. for D- inv mass (TC)
-  // TH2F *fInvMassK0S; //!hist. for D- inv mass (TC)
-  TH2F *fPtVsMass;    //! hist. of pt vs. mass (prod. cuts)
-  TH2F *fPtVsMassTC;  //! hist. of pt vs. mass (analysis cuts)
-  TH2F *fYVsPt;       //! hist. of Y vs. Pt (prod. cuts)
-  TH2F *fYVsPtTC;     //! hist. of Y vs. Pt (analysis cuts)
-  TH2F *fYVsPtSig;    //! hist. of Y vs. Pt (MC, only sig, prod. cuts)
-  TH2F *fYVsPtSigTC;    //! hist. of Y vs. Pt (MC, only sig, analysis cuts)
+  TH1F *fMassHistOrigC[3*kMaxPtBins]; //!hist. for inv mass (TC)
+  TH1F *fMassHistOrigB[3*kMaxPtBins]; //!hist. for inv mass (TC)
+  TH1F *fMassHistMC[3*kMaxPtBins]; //!hist. for inv mass (TC)
+
   Float_t fUpmasslimit;  //upper inv mass limit for histos
   Float_t fLowmasslimit; //lower inv mass limit for histos
+  
   Int_t fNPtBins; //Number of Pt Bins
+  
   Float_t fBinWidth;//width of one bin in output histos
+
   TList *fListCuts; //list of cuts
-  TList *fListCutsAsso; //list of cuts
+  
+  //TList *fListCutsAsso; //list of cuts
+
   AliRDHFCutsDplustoKpipi *fRDCutsAnalysis; //Cuts for Analysis
+  
   AliHFAssociatedTrackCuts *fCuts; 
+  
   AliNormalizationCounter *fCounter;//!Counter for normalization
+  
   Double_t fArrayBinLimits[kMaxPtBins+1]; //limits for the Pt bins
+  
   Bool_t fReadMC;    //flag for access to MC
   //  Bool_t fUseStrangeness;//flag to enhance strangeness in MC to fit to data
   Bool_t fUseBit;      // flag to use bitmask
   Bool_t fMixing;      // flag to use bitmask
   
   Bool_t fSystem; //
-  
-  ClassDef(AliAnalysisTaskSEDplusCorrelations,2); // AliAnalysisTaskSE for the MC association of heavy-flavour decay candidates
+  Bool_t fReco; // use reconstruction or MC truth
+  ClassDef(AliAnalysisTaskSEDplusCorrelations,3); // AliAnalysisTaskSE for the MC association of heavy-flavour decay candidates
 };
 
 #endif
index 0598479ccf135d66168b896c3a13274b4a0050f1..f59554f23bf2d092de83e4b02884180eb36b721c 100644 (file)
@@ -1,11 +1,11 @@
-AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorr(Bool_t system=kFALSE,
-                                                    Bool_t readMC=kFALSE,
-                                                    Bool_t mixing=kFALSE,
-                                                    Int_t select=2,
-                                                    Int_t usedisp=0,
-                                                    TString finDirname="Loose",
-                                                    TString filename="",
-                                                    TString finAnObjname="AnalysisCuts")
+AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorrWT(Bool_t system=kFALSE,                                                                                                                                Bool_t readMC=kFALSE,
+                                            Bool_t ifTrig=kTRUE,            
+                                            Bool_t mixing=kFALSE,
+                                            Int_t select=1,
+                                            Bool_t reco=kTRUE,                                                              Int_t usedisp=0,
+                                            TString finDirname="",
+                                            TString filename="",
+                                            TString finAnObjname="AnalysisCuts")
 {
 
     //                                                                                                                                    
@@ -32,7 +32,7 @@ AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorr(Bool_t system=kFALSE,
       }
   }
   
-  TFile* filecuts1=new TFile("AssocPartCutsDplus.root");
+  TFile* filecuts1=new TFile("AssocPartCuts.root");
          if(!filecuts1->IsOpen()){
                  cout<<"Input file1 not found: exit"<<endl;
                  return;
@@ -44,6 +44,7 @@ AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorr(Bool_t system=kFALSE,
   
   AliRDHFCutsDplustoKpipi* analysiscuts=new AliRDHFCutsDplustoKpipi();
   if(stdcuts) {
+    cout<<analysiscuts->GetNPtBins()<<endl;
     if(system==0) analysiscuts->SetStandardCutsPP2010();
     else if(system==1){
       analysiscuts->SetStandardCutsPbPb2011();
@@ -64,15 +65,18 @@ AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorr(Bool_t system=kFALSE,
   
     
   AliAnalysisTaskSEDplusCorrelations *dplusTask = new AliAnalysisTaskSEDplusCorrelations("DplusAnalysis",analysiscuts,corrCuts);
+
+
   dplusTask->SetReadMC(readMC);
   dplusTask->SetEventMix(mixing);
   dplusTask->SetCorrelator(select);
   dplusTask->SetUseDisplacement(usedisp);
+  dplusTask->SetTrigEfficiency(ifTrig);
   dplusTask->SetDebugLevel(0);
   dplusTask->SetMassLimits(0.2);
   dplusTask->SetUseBit(kTRUE);
   dplusTask->SetSystem(kFALSE);
-   
+dplusTask->SetUseReconstruction(kTRUE);   
 
   //  if (system==0) dplusTask->SetDoImpactParameterHistos(kTRUE);
 
@@ -84,25 +88,24 @@ AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorr(Bool_t system=kFALSE,
   TString outname = "coutputDplus";
   TString cutsname = "coutputDplusCuts";
   TString normname = "coutputDplusNorm";
-  TString trackcutsname = "coutputTrackCuts";
+  //TString assocutsname = "coutputAssoCuts";
   inname += finDirname.Data();
   outname += finDirname.Data();
   cutsname += finDirname.Data();
   normname += finDirname.Data();
-  trackcutsname += finDirname.Data();
+  //assocutsname += finDirname.Data();
   TString centr=Form("%.0f%.0f",analysiscuts->GetMinCentrality(),analysiscuts->GetMaxCentrality());
   inname += centr;
   outname += centr;
   cutsname += centr;
   normname += centr;
-  trackcutsname += centr;
-  
+  // assocutsname += centr;
 
 
   AliAnalysisDataContainer *cinputDplus = mgr->CreateContainer(inname,TChain::Class(),
                                                               AliAnalysisManager::kInputContainer);
   TString outputfile = AliAnalysisManager::GetCommonFileName();
-  outputfile += ":PWG3_D2H_InvMassDplus";
+  outputfile += ":PWGHF_CorrHF_Dplus";
   
   AliAnalysisDataContainer *coutputDplusCuts = mgr->CreateContainer(cutsname,TList::Class(),
                                                                    AliAnalysisManager::kOutputContainer,
@@ -115,17 +118,21 @@ AliAnalysisTaskSEDplusCorrelations *AddTaskDplusCorr(Bool_t system=kFALSE,
                                                                AliAnalysisManager::kOutputContainer,
                                                                outputfile.Data());
   
-AliAnalysisDataContainer *coutputTrackCuts = mgr->CreateContainer(trackcutsname,AliHFAssociatedTrackCuts::Class(),
-                                                                 AliAnalysisManager::kOutputContainer,                                                                                                          outputfile.Data());
+
+  //  AliAnalysisDataContainer *coutputAssoCuts = mgr->CreateContainer(assocutsname,TList::Class(),
+  //                                                              AliAnalysisManager::kOutputContainer,
+  //                                                              outputfile.Data());
+    
+
+  mgr->ConnectInput(dplusTask,0,mgr->GetCommonInputContainer());
   
- mgr->ConnectInput(dplusTask,0,mgr->GetCommonInputContainer());
   mgr->ConnectOutput(dplusTask,1,coutputDplus);
   
   mgr->ConnectOutput(dplusTask,2,coutputDplusCuts);
-  
+
   mgr->ConnectOutput(dplusTask,3,coutputDplusNorm);  
-  mgr->ConnectOutput(dplusTask,4,coutputTrackCuts); 
+
+  // mgr->ConnectOutput(dplusTask,4,coutputAssoCuts);
   
   return dplusTask;
 }