]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Adding D meson efficiency maps and Coverity fixes (Sandro)
authorarossi <arossi@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 19 Jul 2013 13:08:41 +0000 (13:08 +0000)
committerarossi <arossi@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 19 Jul 2013 13:08:41 +0000 (13:08 +0000)
PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.cxx
PWGHF/correlationHF/AliAnalysisTaskDStarCorrelations.h
PWGHF/correlationHF/AliHFAssociatedTrackCuts.cxx
PWGHF/correlationHF/AliHFAssociatedTrackCuts.h
PWGHF/correlationHF/AliHFCorrelator.cxx
PWGHF/correlationHF/AliHFCorrelator.h
PWGHF/correlationHF/macros/AddTaskDStarCorrelations.C

index 8c9fa6f98091fbf3272f56bfd17aa107fa275b06..6e202c07a38283a5e0ee5bb059f31bcfd67775af 100644 (file)
@@ -75,12 +75,17 @@ fSystem(pp),
 fReco(kTRUE),\r
 fUseEfficiencyCorrection(kFALSE),\r
 fUseDmesonEfficiencyCorrection(kFALSE),\r
+fUseCentrality(kFALSE),\r
 fPhiBins(32),\r
 fEvents(0),\r
 fDebugLevel(0),\r
 fDisplacement(0),\r
+fDim(0),\r
+fNofPtBins(0),\r
 fDMesonSigmas(0),\r
 \r
+fD0Window(0),\r
+\r
 fOutput(0x0),\r
 fOutputMC(0x0),\r
 fCuts(0),\r
@@ -89,7 +94,7 @@ fUtils(0),
 fTracklets(0),\r
 fDeffMapvsPt(0),\r
 fDeffMapvsPtvsMult(0),\r
-fDeffMapvsPtvsMultvsEta(0)\r
+fDeffMapvsPtvsEta(0)\r
 {\r
     SetDim();\r
     // default constructor\r
@@ -115,7 +120,10 @@ fPhiBins(32),
 fEvents(0),\r
 fDebugLevel(0),\r
 fDisplacement(0),\r
+fDim(0),\r
+fNofPtBins(0),\r
 fDMesonSigmas(0),\r
+fD0Window(0),\r
 \r
 fOutput(0x0),\r
 fOutputMC(0x0),\r
@@ -125,11 +133,20 @@ fUtils(0),
 fTracklets(0),\r
 fDeffMapvsPt(0),\r
 fDeffMapvsPtvsMult(0),\r
-fDeffMapvsPtvsMultvsEta(0)\r
+fDeffMapvsPtvsEta(0)\r
 {\r
-  SetDim();\r
-  fCuts=cuts;\r
-  Info("AliAnalysisTaskDStarCorrelations","Calling Constructor");\r
+     Info("AliAnalysisTaskDStarCorrelations","Calling Constructor");\r
+  \r
+    SetDim();\r
+    if(fSystem == AA)  fUseCentrality = kTRUE; else fUseCentrality = kFALSE;\r
+  \r
+    fCuts=cuts;\r
+    fNofPtBins= fCuts->GetNPtBins();\r
+    //cout << "Enlarging the DZero window " << endl;\r
+    EnlargeDZeroMassWindow();\r
+   // cout << "Done" << endl;\r
+    \r
+   \r
   DefineInput(0, TChain::Class());\r
   \r
   DefineOutput(1,TList::Class()); // histos from data and MC\r
@@ -159,7 +176,7 @@ AliAnalysisTaskDStarCorrelations::~AliAnalysisTaskDStarCorrelations() {
        if(fCuts2) {delete fCuts2; fCuts2=0;}\r
     if(fDeffMapvsPt){delete fDeffMapvsPt; fDeffMapvsPt=0;}\r
     if(fDeffMapvsPtvsMult){delete fDeffMapvsPtvsMult; fDeffMapvsPtvsMult=0;}\r
-    if(fDeffMapvsPtvsMultvsEta){delete fDeffMapvsPtvsMultvsEta; fDeffMapvsPtvsMultvsEta=0;}\r
+    if(fDeffMapvsPtvsEta){delete fDeffMapvsPtvsEta; fDeffMapvsPtvsEta=0;}\r
 \r
 }\r
 \r
@@ -203,11 +220,12 @@ void AliAnalysisTaskDStarCorrelations::UserCreateOutputObjects(){
        DefineHistoForAnalysis();\r
     DefineThNSparseForAnalysis();\r
     \r
+\r
        fCounter = new AliNormalizationCounter(Form("%s",GetOutputSlot(5)->GetContainer()->GetName()));\r
        fCounter->Init();\r
        \r
     Double_t Pi = TMath::Pi();\r
-       fCorrelator = new AliHFCorrelator("Correlator",fCuts2,fSystem); // fCuts2 is the hadron cut object, fSystem to switch between pp or PbPb\r
+       fCorrelator = new AliHFCorrelator("Correlator",fCuts2,fUseCentrality); // fCuts2 is the hadron cut object, fSystem to switch between pp or PbPb\r
        fCorrelator->SetDeltaPhiInterval(  -0.5*Pi - Pi/fPhiBins, 1.5*Pi- Pi/fPhiBins); // set correct phi interval\r
        //fCorrelator->SetDeltaPhiInterval((-0.5)*Pi,(1.5)*Pi); // set correct phi interval\r
        fCorrelator->SetEventMixing(fmixing); //set kFALSE/kTRUE for mixing Off/On\r
@@ -231,118 +249,119 @@ void AliAnalysisTaskDStarCorrelations::UserCreateOutputObjects(){
 \r
 //_________________________________________________\r
 void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){\r
+  \r
+  \r
+  \r
+  \r
+  if(fDebugLevel){\r
+    \r
+    if(fReco) std::cout << "USING RECONSTRUCTION" << std::endl;\r
+    if(!fReco) std::cout << "USING MC TRUTH" << std::endl;\r
+    std::cout << " " << std::endl;\r
+    std::cout << "=================================================================================" << std::endl;\r
+    if(!fmixing){\r
+      if(fselect==1) std::cout << "TASK::Correlation with hadrons on SE "<< std::endl;\r
+      if(fselect==2) std::cout << "TASK::Correlation with kaons on SE "<< std::endl;\r
+      if(fselect==3) std::cout << "TASK::Correlation with kzeros on SE "<< std::endl;\r
+    }\r
+    if(fmixing){\r
+      if(fselect==1) std::cout << "TASK::Correlation with hadrons on ME "<< std::endl;\r
+      if(fselect==2) std::cout << "TASK::Correlation with kaons on ME "<< std::endl;\r
+      if(fselect==3) std::cout << "TASK::Correlation with kzeros on ME "<< std::endl;\r
+    }\r
     \r
-    \r
-    \r
-    \r
-       if(fDebugLevel){\r
-               \r
-               if(fReco) std::cout << "USING RECONSTRUCTION" << std::endl;\r
-               if(!fReco) std::cout << "USING MC TRUTH" << std::endl;\r
-        std::cout << " " << std::endl;\r
-        std::cout << "=================================================================================" << std::endl;\r
-        if(!fmixing){\r
-            if(fselect==1) std::cout << "TASK::Correlation with hadrons on SE "<< std::endl;\r
-            if(fselect==2) std::cout << "TASK::Correlation with kaons on SE "<< std::endl;\r
-            if(fselect==3) std::cout << "TASK::Correlation with kzeros on SE "<< std::endl;\r
-        }\r
-        if(fmixing){\r
-            if(fselect==1) std::cout << "TASK::Correlation with hadrons on ME "<< std::endl;\r
-            if(fselect==2) std::cout << "TASK::Correlation with kaons on ME "<< std::endl;\r
-            if(fselect==3) std::cout << "TASK::Correlation with kzeros on ME "<< std::endl;\r
-        }\r
-        \r
-    }// end if debug\r
-       \r
-    \r
-    \r
-       if (!fInputEvent) {\r
-               Error("UserExec","NO EVENT FOUND!");\r
+  }// end if debug\r
+  \r
+  \r
+  \r
+  if (!fInputEvent) {\r
+    Error("UserExec","NO EVENT FOUND!");\r
                return;\r
-       }\r
-       \r
-       AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);\r
-       if(!aodEvent){\r
-        AliError("AOD event not found!");\r
-        return;\r
-       }\r
-       \r
+  }\r
+  \r
+  AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);\r
+  if(!aodEvent){\r
+    AliError("AOD event not found!");\r
+    return;\r
+  }\r
+  \r
     fTracklets = aodEvent->GetTracklets();\r
-       \r
-       fEvents++; // event counter\r
-       ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(0);\r
+    \r
+    fEvents++; // event counter\r
+    ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(0);\r
     \r
     fCounter->StoreEvent(aodEvent,fCuts,fmontecarlo);\r
     \r
-       // load MC array\r
-       fmcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));\r
-       if(fmontecarlo && !fmcArray){\r
-        AliError("Array of MC particles not found");\r
-        return;\r
-       }\r
-       \r
+    // load MC array\r
+    fmcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));\r
+    if(fmontecarlo && !fmcArray){\r
+      AliError("Array of MC particles not found");\r
+      return;\r
+    }\r
+    \r
     \r
     \r
     \r
     // ********************************************** START EVENT SELECTION ****************************************************\r
     \r
-       Bool_t isEvSel=fCuts->IsEventSelected(aodEvent);\r
+    Bool_t isEvSel=fCuts->IsEventSelected(aodEvent);\r
     \r
     if(!isEvSel) {\r
-        \r
-        if(fCuts->IsEventRejectedDueToPileupSPD()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(2);\r
-        if(fCuts->IsEventRejectedDueToCentrality()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(3);\r
-        if(fCuts->IsEventRejectedDueToNotRecoVertex()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(4);\r
-        if(fCuts->IsEventRejectedDueToVertexContributors()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(5);\r
-        if(fCuts->IsEventRejectedDueToZVertexOutsideFiducialRegion()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(6);\r
-        if(fCuts->IsEventRejectedDueToTrigger()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(7);\r
-        if(fCuts->IsEventRejectedDuePhysicsSelection()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(8);\r
-        \r
-        return;\r
+      \r
+      if(fCuts->IsEventRejectedDueToPileupSPD()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(2);\r
+      if(fCuts->IsEventRejectedDueToCentrality()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(3);\r
+      if(fCuts->IsEventRejectedDueToNotRecoVertex()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(4);\r
+      if(fCuts->IsEventRejectedDueToVertexContributors()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(5);\r
+      if(fCuts->IsEventRejectedDueToZVertexOutsideFiducialRegion()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(6);\r
+      if(fCuts->IsEventRejectedDueToTrigger()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(7);\r
+      if(fCuts->IsEventRejectedDuePhysicsSelection()) ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(8);\r
+      \r
+      return;\r
     }\r
     \r
     // added event selection for pA\r
     \r
     if(fSystem == pA){\r
-        \r
-        if(fUtils->IsFirstEventInChunk(aodEvent)) {\r
-            AliInfo("Rejecting the event - first in the chunk");\r
-            ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(9);\r
-            return;\r
-        }\r
-        if(!fUtils->IsVertexSelected2013pA(aodEvent)) {\r
-            ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(10);\r
-            AliInfo("Rejecting the event - bad vertex");\r
-            return;\r
+      \r
+      if(fUtils->IsFirstEventInChunk(aodEvent)) {\r
+       AliInfo("Rejecting the event - first in the chunk");\r
+       ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(9);\r
+       return;\r
         }\r
+      if(!fUtils->IsVertexSelected2013pA(aodEvent)) {\r
+       ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(10);\r
+       AliInfo("Rejecting the event - bad vertex");\r
+       return;\r
+      }\r
     }\r
     // ******************************** END event selections **************************************************\r
     \r
     AliCentrality *centralityObj = 0;\r
-       //Int_t multiplicity = -1;\r
-       Double_t MultipOrCent = -1;\r
-    \r
-    if(fSystem == AA ){                centralityObj = aodEvent->GetHeader()->GetCentralityP();\r
-               MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");\r
-               //AliInfo(Form("Centrality is %f", MultipOrCent));\r
-       }\r
+    //Int_t multiplicity = -1;\r
+    Double_t MultipOrCent = -1;\r
     \r
+    if(fUseCentrality){\r
+      /* if(fSystem == AA ){   */      centralityObj = aodEvent->GetHeader()->GetCentralityP();\r
+      MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");\r
+      //AliInfo(Form("Centrality is %f", MultipOrCent));\r
+    }\r
     \r
-    if(fSystem == pp || fSystem == pA){\r
+    if(!fUseCentrality){\r
+      /* if(fSystem == pp || fSystem == pA){*/\r
         //     MultipOrCent = fTracklets->GetNumberOfTracklets();\r
-        MultipOrCent = AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.,1.);\r
-               //AliInfo(Form("multiplicity is %f", MultipOrCent));\r
+      MultipOrCent = AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.,1.);\r
+      //AliInfo(Form("multiplicity is %f", MultipOrCent));\r
     }\r
     \r
     \r
     fCorrelator->SetAODEvent(aodEvent); // set the event to be processed\r
     \r
     ((TH1D*)fOutput->FindObject("NofEvents"))->Fill(1);\r
-       //\r
-       Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing\r
+    //\r
+    Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing\r
        if(!correlatorON) {\r
-               AliInfo("AliHFCorrelator didn't initialize the pool correctly or processed a bad event");\r
-               return;\r
+         AliInfo("AliHFCorrelator didn't initialize the pool correctly or processed a bad event");\r
+         return;\r
        }\r
        \r
        if(fmontecarlo) fCorrelator->SetMCArray(fmcArray);\r
@@ -351,59 +370,59 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){
        // check the event type\r
        // load MC header\r
        if(fmontecarlo){\r
-               AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));\r
-               if (fmontecarlo && !mcHeader) {\r
-                       AliError("Could not find MC Header in AOD");\r
-                       return;\r
-               }\r
-        \r
-               Bool_t isMCeventgood = kFALSE;\r
-        \r
-               \r
-               Int_t eventType = mcHeader->GetEventType();\r
-               Int_t NMCevents = fCuts2->GetNofMCEventType();\r
-               \r
-               for(Int_t k=0; k<NMCevents; k++){\r
-                       Int_t * MCEventType = fCuts2->GetMCEventType();\r
-                       \r
-                       if(eventType == MCEventType[k]) isMCeventgood= kTRUE;\r
-                       ((TH1D*)fOutputMC->FindObject("EventTypeMC"))->Fill(eventType);\r
-               }\r
-               \r
-               if(NMCevents && !isMCeventgood){\r
+         AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));\r
+         if (fmontecarlo && !mcHeader) {\r
+           AliError("Could not find MC Header in AOD");\r
+           return;\r
+         }\r
+         \r
+         Bool_t isMCeventgood = kFALSE;\r
+         \r
+         \r
+         Int_t eventType = mcHeader->GetEventType();\r
+         Int_t NMCevents = fCuts2->GetNofMCEventType();\r
+         \r
+         for(Int_t k=0; k<NMCevents; k++){\r
+           Int_t * MCEventType = fCuts2->GetMCEventType();\r
+           \r
+           if(eventType == MCEventType[k]) isMCeventgood= kTRUE;\r
+           ((TH1D*)fOutputMC->FindObject("EventTypeMC"))->Fill(eventType);\r
+         }\r
+         \r
+         if(NMCevents && !isMCeventgood){\r
             if(fDebugLevel)    std::cout << "The MC event " << eventType << " not interesting for this analysis: skipping" << std::endl;\r
-                       return;\r
-               }\r
-               \r
+           return;\r
+         }\r
+         \r
        } // end if montecarlo\r
        \r
-    \r
-    // checks on vertex and multiplicity of the event\r
-    AliAODVertex *vtx = aodEvent->GetPrimaryVertex();\r
+       \r
+       // checks on vertex and multiplicity of the event\r
+       AliAODVertex *vtx = aodEvent->GetPrimaryVertex();\r
        Double_t zVtxPosition = vtx->GetZ(); // zvertex\r
-    \r
-    \r
-    if(fFullmode) ((TH2F*)fOutput->FindObject("EventPropsCheck"))->Fill(MultipOrCent,zVtxPosition);\r
-    \r
-    \r
-    \r
+       \r
+       \r
+       if(fFullmode) ((TH2F*)fOutput->FindObject("EventPropsCheck"))->Fill(MultipOrCent,zVtxPosition);\r
+       \r
+       \r
+       \r
        // D* reconstruction\r
        TClonesArray *arrayDStartoD0pi=0;\r
        if(!aodEvent && AODEvent() && IsStandardAOD()) {\r
-               // In case there is an AOD handler writing a standard AOD, use the AOD\r
-               // event in memory rather than the input (ESD) event.\r
-               aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());\r
-               // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)\r
-               // have to taken from the AOD event hold by the AliAODExtension\r
-               AliAODHandler* aodHandler = (AliAODHandler*)\r
-               ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());\r
-               if(aodHandler->GetExtensions()) {\r
-                       AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");\r
-                       AliAODEvent *aodFromExt = ext->GetAOD();\r
-                       arrayDStartoD0pi=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");\r
-               }\r
+         // In case there is an AOD handler writing a standard AOD, use the AOD\r
+         // event in memory rather than the input (ESD) event.\r
+         aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());\r
+         // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)\r
+         // have to taken from the AOD event hold by the AliAODExtension\r
+         AliAODHandler* aodHandler = (AliAODHandler*)\r
+           ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());\r
+         if(aodHandler->GetExtensions()) {\r
+           AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");\r
+           AliAODEvent *aodFromExt = ext->GetAOD();\r
+           arrayDStartoD0pi=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");\r
+         }\r
        } else {\r
-               arrayDStartoD0pi=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");\r
+         arrayDStartoD0pi=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");\r
        }\r
        \r
        if(!aodEvent->GetPrimaryVertex() || TMath::Abs(aodEvent->GetMagneticField())<0.001) return;\r
@@ -420,80 +439,81 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){
        \r
        Double_t mPDGD0=1.8648;//TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
        Double_t mPDGDstar=2.01022;//TDatabasePDG::Instance()->GetParticle(413)->Mass();\r
-    \r
+       \r
        \r
        //MC tagging for DStar\r
        //D* and D0 prongs needed to MatchToMC method\r
        Int_t pdgDgDStartoD0pi[2]={421,211};\r
        Int_t pdgDgD0toKpi[2]={321,211};\r
-    \r
+       \r
        Bool_t isDStarCand = kFALSE;\r
-    Bool_t isFilled = kFALSE;\r
-    Bool_t isEventMixingFilledPeak = kFALSE;\r
-    Bool_t isEventMixingFilledSB = kFALSE;\r
+       Bool_t isFilled = kFALSE;\r
+       Bool_t isEventMixingFilledPeak = kFALSE;\r
+       Bool_t isEventMixingFilledSB = kFALSE;\r
        //loop on D* candidates\r
        \r
        Int_t looponDCands = 0;\r
        if(fReco) looponDCands = arrayDStartoD0pi->GetEntriesFast();\r
        if(!fReco) looponDCands = fmcArray->GetEntriesFast();\r
-    \r
-    Int_t nOfDStarCandidates = 0;\r
-    Int_t nOfSBCandidates = 0;\r
        \r
-    \r
-    \r
-    Double_t DmesonEfficiency = 1.;\r
-    Double_t DmesonWeight = 1.;\r
-    Bool_t isDfromB = kFALSE;\r
-    \r
-    \r
+       Int_t nOfDStarCandidates = 0;\r
+       Int_t nOfSBCandidates = 0;\r
+       \r
+       \r
+       \r
+       Double_t DmesonEfficiency = 1.;\r
+       Double_t DmesonWeight = 1.;\r
+       Bool_t isDfromB = kFALSE;\r
+       \r
+       \r
        for (Int_t iDStartoD0pi = 0; iDStartoD0pi<looponDCands; iDStartoD0pi++) {\r
-               isInPeak = kFALSE;\r
-               isInDStarSideBand = kFALSE;\r
-        isInDZeroSideBand = kFALSE;\r
-               isDStarMCtag = kFALSE;\r
-        isDfromB = kFALSE;\r
-               ptDStar = -123.4;\r
-               phiDStar = -999;\r
-               etaDStar = -56.;\r
-               invMassDZero = - 999;\r
-               deltainvMDStar = -998;\r
-               AliAODRecoCascadeHF* dstarD0pi;\r
-               AliAODRecoDecayHF2Prong* theD0particle;\r
-               AliAODMCParticle* DStarMC;\r
-               Short_t daughtercharge;\r
-               Int_t trackiddaugh0; // track id if it is reconstruction - label if it is montecarlo info\r
-               Int_t trackiddaugh1;\r
-               Int_t trackidsoftPi;\r
-               \r
-        // start the if reconstructed candidates from here ************************************************\r
-        \r
-               if(fReco){//// if reconstruction is applied\r
-                       dstarD0pi = (AliAODRecoCascadeHF*)arrayDStartoD0pi->At(iDStartoD0pi);\r
-                       if(!dstarD0pi->GetSecondaryVtx()) continue;\r
-                       theD0particle = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();\r
-                       if (!theD0particle) continue;\r
-            \r
-\r
-            \r
-                       // track quality cuts\r
-                       Int_t isTkSelected = fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kTracks); // quality cuts on tracks\r
-                       // region of interest + topological cuts + PID\r
-                       Int_t isSelected=fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kCandidate); //selected\r
+         isInPeak = kFALSE;\r
+         isInDStarSideBand = kFALSE;\r
+         isInDZeroSideBand = kFALSE;\r
+         isDStarMCtag = kFALSE;\r
+         isDfromB = kFALSE;\r
+         ptDStar = -123.4;\r
+         phiDStar = -999;\r
+         etaDStar = -56.;\r
+         invMassDZero = - 999;\r
+         deltainvMDStar = -998;\r
+         AliAODRecoCascadeHF* dstarD0pi;\r
+         AliAODRecoDecayHF2Prong* theD0particle;\r
+         AliAODMCParticle* DStarMC;\r
+         Short_t daughtercharge;\r
+         Int_t trackiddaugh0 = -1; // track id if it is reconstruction - label if it is montecarlo info\r
+         Int_t trackiddaugh1 = -1;\r
+         Int_t trackidsoftPi = -1;\r
+         \r
+         // start the if reconstructed candidates from here ************************************************\r
+         \r
+         if(fReco){//// if reconstruction is applied\r
+           dstarD0pi = (AliAODRecoCascadeHF*)arrayDStartoD0pi->At(iDStartoD0pi);\r
+           if(!dstarD0pi->GetSecondaryVtx()) continue;\r
+           theD0particle = (AliAODRecoDecayHF2Prong*)dstarD0pi->Get2Prong();\r
+           if (!theD0particle) continue;\r
+            \r
+           \r
+            \r
+           // track quality cuts\r
+           Int_t isTkSelected = fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kTracks); // quality cuts on tracks\r
+           // region of interest + topological cuts + PID\r
+           Int_t isSelected=fCuts->IsSelected(dstarD0pi,AliRDHFCuts::kCandidate); //selected\r
             \r
             //apply track selections\r
-                   if(!isTkSelected) continue;\r
-                       if(!isSelected) continue;\r
-                  if(!fCuts->IsInFiducialAcceptance(dstarD0pi->Pt(),dstarD0pi->YDstar())) continue;\r
+           if(!isTkSelected) continue;\r
+           if(!isSelected) continue;\r
+           //      if(!fCuts->IsInFiducialAcceptance(dstarD0pi->Pt(),dstarD0pi->YDstar())) continue;\r
             \r
             // new piece of code\r
             \r
             \r
-            Double_t rapidity = dstarD0pi->YDstar();\r
+           //            Double_t rapidity = dstarD0pi->YDstar();\r
             \r
+            DmesonEfficiency = fCuts2->GetTrigWeight(dstarD0pi->Pt(),MultipOrCent);\r
             \r
-            if(fDeffMapvsPt){\r
-                if(fDebugLevel)  std::cout << "Reading pT eff map from " << fDeffMapvsPt->GetName() <<  std::endl;\r
+           /* if(fDeffMapvsPt){\r
+              if(fDebugLevel)  std::cout << "Reading pT eff map from " << fDeffMapvsPt->GetName() <<  std::endl;\r
                 Int_t bin=fDeffMapvsPt->FindBin(dstarD0pi->Pt());\r
                 if(fDeffMapvsPt->IsBinUnderflow(bin)||fDeffMapvsPt->IsBinOverflow(bin)) DmesonEfficiency = 1.;\r
                 else DmesonEfficiency = fDeffMapvsPt->GetBinContent(bin);\r
@@ -501,174 +521,152 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){
                 \r
                 \r
                 \r
-            }\r
-            if(fDeffMapvsPtvsMult)\r
-            {\r
+               }\r
+               if(fDeffMapvsPtvsMult)\r
+               {\r
                 Int_t bin=fDeffMapvsPtvsMult->FindBin(dstarD0pi->Pt(),MultipOrCent);\r
                 if(fDeffMapvsPtvsMult->IsBinUnderflow(bin)||fDeffMapvsPtvsMult->IsBinOverflow(bin)) DmesonEfficiency = 1.;\r
                 else DmesonEfficiency = fDeffMapvsPtvsMult->GetBinContent(bin);\r
             }\r
-            if(fDeffMapvsPtvsMultvsEta){\r
-                Int_t bin=fDeffMapvsPtvsMultvsEta->FindBin(dstarD0pi->Pt(),rapidity,MultipOrCent);\r
-                if(fDeffMapvsPtvsMultvsEta->IsBinUnderflow(bin)||fDeffMapvsPtvsMultvsEta->IsBinOverflow(bin)) DmesonEfficiency = 1.;\r
-                else DmesonEfficiency = fDeffMapvsPtvsMultvsEta->GetBinContent(bin);\r
-                \r
+            if(fDeffMapvsPtvsEta){\r
+           Int_t bin=fDeffMapvsPtvsEta->FindBin(dstarD0pi->Pt(),rapidity);\r
+           if(fDeffMapvsPtvsEta->IsBinUnderflow(bin)||fDeffMapvsPtvsEta->IsBinOverflow(bin)) DmesonEfficiency = 1.;\r
+                else DmesonEfficiency = fDeffMapvsPtvsEta->GetBinContent(bin);\r
+               \r
                 \r
                 \r
                 \r
-            }\r
-            \r
+               }\r
+            */\r
             \r
             \r
             if(fUseDmesonEfficiencyCorrection)  DmesonWeight = 1./DmesonEfficiency;\r
             else DmesonWeight = 1.;\r
             \r
-            \r
             // continue;\r
             \r
-                       Int_t mcLabelDStar = -999;\r
-                       if(fmontecarlo){\r
-                               // find associated MC particle for D* ->D0toKpi\r
-                mcLabelDStar = dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,fmcArray,kFALSE);\r
-                if(mcLabelDStar>=0) isDStarMCtag = kTRUE;\r
-                       }\r
+           Int_t mcLabelDStar = -999;\r
+           if(fmontecarlo){\r
+             // find associated MC particle for D* ->D0toKpi\r
+             mcLabelDStar = dstarD0pi->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,fmcArray,kFALSE);\r
+             if(mcLabelDStar>=0) isDStarMCtag = kTRUE;\r
+           }\r
             \r
-                       ptDStar = dstarD0pi->Pt();\r
-                       phiDStar = dstarD0pi->Phi();\r
-                       etaDStar = dstarD0pi->Eta();\r
+           ptDStar = dstarD0pi->Pt();\r
+           phiDStar = dstarD0pi->Phi();\r
+           etaDStar = dstarD0pi->Eta();\r
             \r
             \r
-                       \r
-                       phiDStar = fCorrelator->SetCorrectPhiRange(phiDStar);\r
+           \r
+           phiDStar = fCorrelator->SetCorrectPhiRange(phiDStar);\r
             \r
             // set the phi of the D meson in the correct range\r
             \r
-                       Int_t ptbin=fCuts->PtBin(dstarD0pi->Pt());\r
+           Int_t ptbin=fCuts->PtBin(dstarD0pi->Pt());\r
             \r
-                       \r
-                       Double_t dmDStarWindow = 0.0019/3;// 0.0019 = 3 sigma\r
-                       Double_t mD0Window=0.074/3;\r
-            \r
-            \r
-            \r
-                       if (fSystem == pp || fSystem == pA){ // pp\r
-                               if (ptbin==1)  mD0Window = 0.026/3; //0.5-1\r
-                               if (ptbin==2)  mD0Window = 0.022/3; //1-2\r
-                               if (ptbin==3)  mD0Window = 0.024/3; //2-3\r
-                               if (ptbin==4)  mD0Window = 0.032/3;\r
-                               if (ptbin==5)  mD0Window = 0.032/3;\r
-                               if (ptbin==6)  mD0Window = 0.036/3;\r
-                               if (ptbin==7)  mD0Window = 0.036/3;\r
-                               if (ptbin==8)  mD0Window = 0.036/3;\r
-                               if (ptbin==9)  mD0Window = 0.058/3;\r
-                               if (ptbin==10) mD0Window = 0.058/3;\r
-                               if (ptbin>10)  mD0Window = 0.074/3;\r
-                       }\r
-                       if(fSystem == AA){// PbPb\r
-                               if (ptbin==0)  mD0Window = 0.032/3; //1-1\r
-                               if (ptbin==1)  mD0Window = 0.032/3; //2-3\r
-                               if (ptbin==2)  mD0Window = 0.032/3; //3-4\r
-                               if (ptbin==3)  mD0Window = 0.032/3; //4-5\r
-                               if (ptbin==4)  mD0Window = 0.036/3; //5-6\r
-                               if (ptbin==5)  mD0Window = 0.036/3; //6-8\r
-                               if (ptbin==6)  mD0Window = 0.055/3; //8-12\r
-                               if (ptbin==7)  mD0Window = 0.074/3; //12-16\r
-                               if (ptbin==8)  mD0Window = 0.074/3; //16-24\r
-                               if (ptbin==9)  mD0Window = 0.074/3; //24-35\r
-                       }\r
-            \r
-                       invMassDZero = dstarD0pi->InvMassD0();\r
-                       if(!fmixing && fFullmode) ((TH2F*)fOutput->FindObject("D0InvMass"))->Fill(ptDStar,invMassDZero);\r
-            \r
-                       deltainvMDStar = dstarD0pi->DeltaInvMass();\r
-            \r
-                       //good D0 candidates\r
-                       if (TMath::Abs(invMassDZero-mPDGD0)<fDMesonSigmas[1]*mD0Window){\r
-                \r
-                if(!fmixing)   ((TH3F*)fOutput->FindObject("DeltaInvMass"))->Fill(ptDStar,deltainvMDStar,MultipOrCent,DmesonWeight);\r
-                               \r
-                // good D*?\r
-                if(TMath::Abs(deltainvMDStar-(mPDGDstar-mPDGD0))<fDMesonSigmas[0]*dmDStarWindow){\r
-                    \r
-                    if(!fmixing)       ((TH1F*)fOutput->FindObject("RecoPtDStar"))->Fill(ptDStar,DmesonWeight);\r
-                    if(!fmixing)       ((TH2F*)fOutput->FindObject("PhiEtaTrigger"))->Fill(phiDStar,etaDStar);\r
-                    isInPeak = kTRUE;\r
-                    nOfDStarCandidates++;\r
-                               } // end Good D*\r
-                \r
+           \r
+           Double_t dmDStarWindow = 0.0019/3;// 0.0019 = 3 sigma\r
+           //  Double_t mD0Window=0.074/3;\r
+            \r
+            Double_t mD0Window= fD0Window[ptbin]/3;\r
+            //cout << "Check with new window " << fD0Window[ptbin]/3 << endl;\r
+            \r
+           \r
+            \r
+           invMassDZero = dstarD0pi->InvMassD0();\r
+           if(!fmixing && fFullmode) ((TH2F*)fOutput->FindObject("D0InvMass"))->Fill(ptDStar,invMassDZero);\r
+            \r
+           deltainvMDStar = dstarD0pi->DeltaInvMass();\r
+            \r
+           //good D0 candidates\r
+           if (TMath::Abs(invMassDZero-mPDGD0)<fDMesonSigmas[1]*mD0Window){\r
+             \r
+             if(!fmixing)      ((TH3F*)fOutput->FindObject("DeltaInvMass"))->Fill(ptDStar,deltainvMDStar,MultipOrCent,DmesonWeight);\r
+             \r
+             // good D*?\r
+             if(TMath::Abs(deltainvMDStar-(mPDGDstar-mPDGD0))<fDMesonSigmas[0]*dmDStarWindow){\r
+               \r
+               if(!fmixing)    ((TH1F*)fOutput->FindObject("RecoPtDStar"))->Fill(ptDStar,DmesonWeight);\r
+               if(!fmixing)    ((TH2F*)fOutput->FindObject("PhiEtaTrigger"))->Fill(phiDStar,etaDStar);\r
+               isInPeak = kTRUE;\r
+               nOfDStarCandidates++;\r
+             } // end Good D*\r
+              \r
                 //  D* sideband?\r
-                if((deltainvMDStar-(mPDGDstar-mPDGD0)>4*dmDStarWindow) && (deltainvMDStar-(mPDGDstar-mPDGD0)<8*dmDStarWindow)){\r
-                    isInDStarSideBand = kTRUE;\r
-                               } // end D* sideband\r
-                \r
+             if((deltainvMDStar-(mPDGDstar-mPDGD0)>4*dmDStarWindow) && (deltainvMDStar-(mPDGDstar-mPDGD0)<8*dmDStarWindow)){\r
+               isInDStarSideBand = kTRUE;\r
+             } // end D* sideband\r
+              \r
             }// end good D0 candidates\r
             \r
             //D0 sidebands\r
-                       if (TMath::Abs(invMassDZero-mPDGD0)>fDMesonSigmas[2]*mD0Window && TMath::Abs(invMassDZero-mPDGD0)<fDMesonSigmas[3]*mD0Window ){\r
-                               if(!fmixing)((TH3F*)fOutput->FindObject("bkgDeltaInvMass"))->Fill(ptDStar,deltainvMDStar,MultipOrCent,DmesonWeight);\r
-                               if(!fmixing && fFullmode)((TH2F*)fOutput->FindObject("D0InvMassinSB"))->Fill(ptDStar,invMassDZero,DmesonWeight);\r
-                \r
-                               if(TMath::Abs(deltainvMDStar-(mPDGDstar-mPDGD0))<fDMesonSigmas[0] *dmDStarWindow){ // is in DStar peak region?\r
-                    if(!fmixing)       ((TH1F*)fOutput->FindObject("RecoPtBkg"))->Fill(ptDStar,DmesonWeight);\r
-                                       isInDZeroSideBand = kTRUE;\r
-                    nOfSBCandidates++;\r
-                    if(!fmixing)       ((TH2F*)fOutput->FindObject("PhiEtaSideBand"))->Fill(phiDStar,etaDStar);\r
-                               }\r
-                \r
-                       }//end if sidebands\r
-            \r
-                       // getting the number of triggers in the MCtag D* case\r
+           if (TMath::Abs(invMassDZero-mPDGD0)>fDMesonSigmas[2]*mD0Window && TMath::Abs(invMassDZero-mPDGD0)<fDMesonSigmas[3]*mD0Window ){\r
+             if(!fmixing)((TH3F*)fOutput->FindObject("bkgDeltaInvMass"))->Fill(ptDStar,deltainvMDStar,MultipOrCent,DmesonWeight);\r
+             if(!fmixing && fFullmode)((TH2F*)fOutput->FindObject("D0InvMassinSB"))->Fill(ptDStar,invMassDZero,DmesonWeight);\r
+              \r
+             if(TMath::Abs(deltainvMDStar-(mPDGDstar-mPDGD0))<fDMesonSigmas[0] *dmDStarWindow){ // is in DStar peak region?\r
+               if(!fmixing)    ((TH1F*)fOutput->FindObject("RecoPtBkg"))->Fill(ptDStar,DmesonWeight);\r
+               isInDZeroSideBand = kTRUE;\r
+               nOfSBCandidates++;\r
+               if(!fmixing)    ((TH2F*)fOutput->FindObject("PhiEtaSideBand"))->Fill(phiDStar,etaDStar);\r
+             }\r
+              \r
+           }//end if sidebands\r
+            \r
+           // getting the number of triggers in the MCtag D* case\r
             if(fmontecarlo && isDStarMCtag) ((TH1F*)fOutput->FindObject("MCtagPtDStar"))->Fill(ptDStar);\r
             \r
             \r
-                       if(!isInPeak && !isInDStarSideBand && !isInDZeroSideBand) continue; // skip if it is not side band or peak event - SAVE CPU TIME\r
-                       \r
-                       \r
+           if(!isInPeak && !isInDStarSideBand && !isInDZeroSideBand) continue; // skip if it is not side band or peak event - SAVE CPU TIME\r
+           \r
+           \r
             // check properties of the events containing the D*\r
             \r
             \r
             \r
             if(isInPeak &&!isFilled) {\r
-                \r
-                if(fFullmode)  ((TH2F*)fOutput->FindObject("EventPropsCheckifDStar"))->Fill(MultipOrCent,zVtxPosition);\r
-                isFilled = kTRUE;\r
+             \r
+             if(fFullmode)  ((TH2F*)fOutput->FindObject("EventPropsCheckifDStar"))->Fill(MultipOrCent,zVtxPosition);\r
+             isFilled = kTRUE;\r
             }\r
             \r
             isDStarCand = kTRUE;\r
             \r
             // charge of the daughter od the\r
-                       daughtercharge = ((AliAODTrack*)dstarD0pi->GetBachelor())->Charge();\r
-            \r
-                       \r
-                       trackiddaugh0 = ((AliAODTrack*)theD0particle->GetDaughter(0))->GetID();\r
-                       trackiddaugh1 = ((AliAODTrack*)theD0particle->GetDaughter(1))->GetID();\r
-                       trackidsoftPi = ((AliAODTrack*)dstarD0pi->GetBachelor())->GetID();\r
-                       \r
-                       // end here the reco\r
-            \r
-                       \r
-               }// end of if for applied reconstruction to D*\r
-               \r
-               Int_t DStarLabel = -1;\r
-               \r
-               if(!fReco){ // use pure MC information\r
-                       \r
+           daughtercharge = ((AliAODTrack*)dstarD0pi->GetBachelor())->Charge();\r
+            \r
+           \r
+           trackiddaugh0 = ((AliAODTrack*)theD0particle->GetDaughter(0))->GetID();\r
+           trackiddaugh1 = ((AliAODTrack*)theD0particle->GetDaughter(1))->GetID();\r
+           trackidsoftPi = ((AliAODTrack*)dstarD0pi->GetBachelor())->GetID();\r
+           \r
+           // end here the reco\r
+            \r
+           \r
+         }// end of if for applied reconstruction to D*\r
+         \r
+         Int_t DStarLabel = -1;\r
+         \r
+         if(!fReco){ // use pure MC information\r
+           \r
             \r
             // check if DStar from B\r
             \r
             \r
             \r
-                       DStarMC = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iDStartoD0pi));\r
-                       if (!DStarMC) {\r
-                               AliWarning("Careful: DStar MC Particle not found in tree, skipping");\r
-                               continue;\r
-                       }\r
-                       DStarLabel = DStarMC->GetLabel();\r
-                       \r
-                       Int_t PDG =TMath::Abs(DStarMC->PdgCode());\r
-                       \r
-                       if(PDG !=413) continue; // skip if it is not a DStar\r
-            \r
+           DStarMC = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iDStartoD0pi));\r
+           if (!DStarMC) {\r
+             AliWarning("Careful: DStar MC Particle not found in tree, skipping");\r
+             continue;\r
+           }\r
+           DStarLabel = DStarMC->GetLabel();\r
+            if(DStarLabel>0)isDStarMCtag = kTRUE;\r
+           \r
+           Int_t PDG =TMath::Abs(DStarMC->PdgCode());\r
+           \r
+           if(PDG !=413) continue; // skip if it is not a DStar\r
+           // check fiducial acceptance\r
+            if(!fCuts->IsInFiducialAcceptance(DStarMC->Pt(),DStarMC->Y())) continue;\r
             \r
             //chech if DStar from B\r
             \r
@@ -676,248 +674,250 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){
             Int_t labelMother = DStarMC->GetMother();\r
             \r
             AliAODMCParticle * mother = dynamic_cast<AliAODMCParticle*>(fmcArray->At(labelMother));\r
-            \r
+           \r
+            if(!mother) continue;\r
             Int_t motherPDG =TMath::Abs(mother->PdgCode());\r
             \r
             \r
             \r
             if((motherPDG>=500 && motherPDG <600) || (motherPDG>=5000 && motherPDG<6000 ))\r
-            {isDfromB = kTRUE; }\r
+             {isDfromB = kTRUE; }\r
             \r
             \r
             //end check\r
-                       \r
-                       \r
-                       \r
-                       \r
-                       Bool_t isDZero = kFALSE;\r
-                       Bool_t isSoftPi = kFALSE;\r
-                       \r
-            \r
-            \r
-            \r
-            \r
-            \r
-            \r
-                       \r
-                       //check decay channel on MC ************************************************\r
-                       \r
-                               Int_t NDaugh = DStarMC->GetNDaughters();\r
-                               if(NDaugh != 2) continue; // skip decay channels w/0 2 prongs\r
-                               for(Int_t i=0; i<NDaugh;i++){ // loop on daughters\r
-                                       \r
-                                       Int_t daugh_label = DStarMC->GetDaughter(i);\r
-                                       AliAODMCParticle* mcDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daugh_label));\r
-                                       \r
-                                       Int_t daugh_pdg = TMath::Abs(mcDaughter->GetPdgCode());\r
-                                       if(fDebugLevel) std::cout << "Daughter " << i << " pdg code is " << daugh_pdg << std::endl;\r
-                                       \r
-                                       if(daugh_pdg == 421) {isDZero = kTRUE;\r
-                                               Int_t NDaughD0 = mcDaughter->GetNDaughters();\r
-                                               if(NDaughD0 != 2) continue; // skip decay channels w/0 2 prongs\r
-                                               trackiddaugh0 = mcDaughter->GetDaughter(0);\r
-                                               trackiddaugh1 = mcDaughter->GetDaughter(1);\r
-                                               \r
-                                               Bool_t isKaon = kFALSE;\r
-                                               Bool_t isPion = kFALSE;\r
-                                               \r
-                                               for(Int_t k=0;k<NDaughD0;k++){\r
-                                                       \r
-                                                       Int_t labelD0daugh = mcDaughter->GetDaughter(k);\r
-                                                       AliAODMCParticle* mcGrandDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(labelD0daugh));\r
-                                                       Int_t granddaugh_pdg = TMath::Abs(mcGrandDaughter->GetPdgCode());\r
-                                                       if(granddaugh_pdg==321) isKaon = kTRUE;\r
-                                                       if(granddaugh_pdg==211) isPion = kTRUE;\r
-                                               }\r
-                                               if(!isKaon || !isKaon) continue; // skip if not correct decay channel of D0\r
-                                       }\r
-                                       \r
-                                       if(daugh_pdg == 211) {\r
-                                               isSoftPi = kTRUE;\r
-                                               daughtercharge = mcDaughter->Charge();\r
-                                               trackidsoftPi = daugh_label;}\r
-                               }\r
-                               \r
-                               \r
-                               if(!isDZero || !isSoftPi) continue; // skip if not correct decay channel\r
-                               \r
-                       \r
-                       // end check decay channel\r
-                       \r
-                       ptDStar = DStarMC->Pt();\r
-                       phiDStar = DStarMC->Phi();\r
-                       etaDStar = DStarMC->Eta();\r
-            \r
-               } // end use pure MC information\r
-        \r
-        \r
-        // check if it is a DStar from B\r
-        \r
-        \r
-        \r
-        \r
-        \r
-               \r
-        \r
+           \r
+           \r
+           \r
+           \r
+           Bool_t isDZero = kFALSE;\r
+           Bool_t isSoftPi = kFALSE;\r
+           \r
+            \r
+            \r
+            \r
+            \r
+            \r
+            \r
+           \r
+           //check decay channel on MC ************************************************\r
+           \r
+           Int_t NDaugh = DStarMC->GetNDaughters();\r
+           if(NDaugh != 2) continue; // skip decay channels w/0 2 prongs\r
+           for(Int_t i=0; i<NDaugh;i++){ // loop on daughters\r
+             \r
+             Int_t daugh_label = DStarMC->GetDaughter(i);\r
+             AliAODMCParticle* mcDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(daugh_label));\r
+             if(!mcDaughter) continue;\r
+             Int_t daugh_pdg = TMath::Abs(mcDaughter->GetPdgCode());\r
+             if(fDebugLevel) std::cout << "Daughter " << i << " pdg code is " << daugh_pdg << std::endl;\r
+             \r
+             if(daugh_pdg == 421) {isDZero = kTRUE;\r
+               Int_t NDaughD0 = mcDaughter->GetNDaughters();\r
+               if(NDaughD0 != 2) continue; // skip decay channels w/0 2 prongs\r
+               trackiddaugh0 = mcDaughter->GetDaughter(0);\r
+               trackiddaugh1 = mcDaughter->GetDaughter(1);\r
                \r
-               fCorrelator->SetTriggerParticleProperties(ptDStar,phiDStar,etaDStar); // pass to the object the necessary trigger part parameters\r
-               \r
-               \r
-               fCorrelator->SetTriggerParticleDaughterCharge(daughtercharge);\r
-               \r
-               \r
-               // ************************************************ CORRELATION ANALYSIS STARTS HERE\r
-               \r
-        \r
\r
-        \r
-               Bool_t execPool = fCorrelator->ProcessEventPool();\r
+               Bool_t isKaon = kFALSE;\r
+               Bool_t isPion = kFALSE;\r
                \r
-        if(fmixing && !execPool) {\r
-                       AliInfo("Mixed event analysis: pool is not ready");\r
+               for(Int_t k=0;k<NDaughD0;k++){\r
+                 \r
+                 Int_t labelD0daugh = mcDaughter->GetDaughter(k);\r
+                 AliAODMCParticle* mcGrandDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(labelD0daugh));\r
+                 if(!mcGrandDaughter) continue;\r
+                 Int_t granddaugh_pdg = TMath::Abs(mcGrandDaughter->GetPdgCode());\r
+                 if(granddaugh_pdg==321) isKaon = kTRUE;\r
+                 if(granddaugh_pdg==211) isPion = kTRUE;\r
+               }\r
+               if(!isKaon || !isKaon) continue; // skip if not correct decay channel of D0\r
+             }\r
+             \r
+             if(daugh_pdg == 211) {\r
+               isSoftPi = kTRUE;\r
+               daughtercharge = mcDaughter->Charge();\r
+               trackidsoftPi = daugh_label;}\r
+           }\r
+           \r
+           \r
+           if(!isDZero || !isSoftPi) continue; // skip if not correct decay channel\r
+           \r
+           \r
+           // end check decay channel\r
+           \r
+           ptDStar = DStarMC->Pt();\r
+           phiDStar = DStarMC->Phi();\r
+           etaDStar = DStarMC->Eta();\r
+           \r
+         } // end use pure MC information\r
+         \r
+         \r
+         // check if it is a DStar from B\r
+         \r
+         \r
+         if(fmontecarlo && isDStarMCtag && !isDfromB) ((TH1D*)fOutputMC->FindObject("MCtagPtDStarfromCharm"))->Fill(ptDStar);\r
+         if(fmontecarlo && isDStarMCtag && isDfromB) ((TH1D*)fOutputMC->FindObject("MCtagPtDStarfromBeauty"))->Fill(ptDStar);\r
+         \r
+         \r
+         fCorrelator->SetTriggerParticleProperties(ptDStar,phiDStar,etaDStar); // pass to the object the necessary trigger part parameters\r
+         \r
+         \r
+         fCorrelator->SetTriggerParticleDaughterCharge(daughtercharge);\r
+         \r
+         \r
+         // ************************************************ CORRELATION ANALYSIS STARTS HERE\r
+         \r
+         \r
+         \r
+         \r
+         Bool_t execPool = fCorrelator->ProcessEventPool();\r
+         \r
+         \r
+         if(fmixing && !execPool) {\r
+           AliInfo("Mixed event analysis: pool is not ready");\r
             if(!isEventMixingFilledPeak && isInPeak)  {\r
-                ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(1);\r
-                isEventMixingFilledPeak = kTRUE;\r
+             ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(1);\r
+             isEventMixingFilledPeak = kTRUE;\r
             }\r
             if (!isEventMixingFilledSB && isInDZeroSideBand)  {\r
-                ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(3);\r
-                isEventMixingFilledSB=kTRUE;\r
+             ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(3);\r
+             isEventMixingFilledSB=kTRUE;\r
             }\r
             \r
-                       continue;\r
-               }\r
-        \r
-        // check event topology\r
-        if(fmixing&&execPool){\r
+           continue;\r
+         }\r
+         \r
+         // check event topology\r
+         if(fmixing&&execPool){\r
             // pool is ready - run checks on bins filling\r
             if(!isEventMixingFilledPeak && isInPeak)  {\r
-                ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(0);\r
-                if(fFullmode) EventMixingChecks(aodEvent);\r
-                isEventMixingFilledPeak = kTRUE;\r
+             ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(0);\r
+             if(fFullmode) EventMixingChecks(aodEvent);\r
+             isEventMixingFilledPeak = kTRUE;\r
             }\r
             \r
             if(!isEventMixingFilledSB && isInDZeroSideBand) {\r
-                ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(2);\r
-                isEventMixingFilledSB=kTRUE;\r
+             ((TH1D*)fOutput->FindObject("CheckPoolReadiness"))->Fill(2);\r
+             isEventMixingFilledSB=kTRUE;\r
             }\r
             \r
             \r
-        }\r
-               \r
-               Int_t NofEventsinPool = 1;\r
-               if(fmixing) NofEventsinPool = fCorrelator->GetNofEventsInPool();\r
-        \r
-        //  if(fmixing) cout << "Nof events in pool = " << NofEventsinPool << endl;\r
-        \r
-               for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, stops at one\r
-            \r
-                       Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix);\r
-                       \r
-                       if(!analyzetracks) {\r
-                               AliInfo("AliHFCorrelator::Cannot process the track array");\r
-                               continue;\r
-                       }\r
+         }\r
+         \r
+         Int_t NofEventsinPool = 1;\r
+         if(fmixing) NofEventsinPool = fCorrelator->GetNofEventsInPool();\r
+         \r
+         //  if(fmixing) cout << "Nof events in pool = " << NofEventsinPool << endl;\r
+         \r
+         for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, stops at one\r
+           \r
+           Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix);\r
+           \r
+           if(!analyzetracks) {\r
+             AliInfo("AliHFCorrelator::Cannot process the track array");\r
+             continue;\r
+           }\r
             \r
             //initialization of variables for correlations with leading particles\r
             Double_t DeltaPhiLeading = -999.;\r
-                       Double_t DeltaEtaLeading = -999.;\r
-                       //Double_t ptleading = -999.;\r
+           Double_t DeltaEtaLeading = -999.;\r
+           //Double_t ptleading = -999.;\r
             //     Int_t labelleading = -999;\r
             \r
             // Int_t crosscheck = (Int_t)((fCorrelator->GetTrackArray())->GetEntriesFast());\r
-                       \r
-                       \r
-                       Int_t NofTracks = fCorrelator->GetNofTracks();\r
-                       \r
-                       if(isInPeak && fFullmode) ((TH1D*)fOutput->FindObject("NofTracksInPeak"))->Fill(NofTracks);\r
-                       if(isInDZeroSideBand && fFullmode) ((TH1D*)fOutput->FindObject("NofTracksInSB"))->Fill(NofTracks);\r
-                       \r
-                       \r
-                       \r
+           \r
+           \r
+           Int_t NofTracks = fCorrelator->GetNofTracks();\r
+            \r
+           \r
+           if(isInPeak && fFullmode) ((TH1D*)fOutput->FindObject("NofTracksInPeak"))->Fill(NofTracks);\r
+           if(isInDZeroSideBand && fFullmode) ((TH1D*)fOutput->FindObject("NofTracksInSB"))->Fill(NofTracks);\r
+           \r
+           \r
+           \r
             Double_t arraytofill[6];\r
             Double_t MCarraytofill[7];\r
             \r
             \r
             Double_t weight;\r
             \r
-                       for(Int_t iTrack = 0; iTrack<NofTracks; iTrack++){ // looping on track candidates\r
-                               Bool_t runcorrelation = fCorrelator->Correlate(iTrack);\r
-                               if(!runcorrelation) continue;\r
-                \r
-                               Double_t DeltaPhi = fCorrelator->GetDeltaPhi();\r
-                               Double_t DeltaEta = fCorrelator->GetDeltaEta();\r
-                               \r
-                AliReducedParticle * hadron = fCorrelator->GetAssociatedParticle();\r
-                               if(!hadron) {/*cout << "No Hadron" << endl;*/ continue;}\r
-                \r
-                               Double_t ptHad = hadron->Pt();\r
-                               Double_t phiHad = hadron->Phi();\r
-                               Double_t etaHad = hadron->Eta();\r
-                               Int_t label = hadron->GetLabel();\r
-                               Int_t trackid = hadron->GetID();\r
-                Double_t efficiency = hadron->GetWeight();\r
-                \r
-                weight = 1;\r
-                if(fUseEfficiencyCorrection && efficiency){\r
-                    //weight = 1./efficiency;\r
-                    weight = DmesonWeight * (1./efficiency);\r
-                }\r
-                // weight = DmesonWeight * (1./efficiency);\r
-                \r
-                phiHad = fCorrelator->SetCorrectPhiRange(phiHad);\r
-                \r
-                \r
-                if(fFullmode) ((TH2F*)fOutput->FindObject("WeightChecks"))->Fill(ptHad,efficiency);\r
-                \r
-                arraytofill[0] = DeltaPhi;\r
-                arraytofill[1] = DeltaEta;\r
-                arraytofill[2] = ptDStar;\r
-                arraytofill[3] = MultipOrCent;\r
-                arraytofill[4] = ptHad;\r
-                arraytofill[5] = zVtxPosition;\r
-                \r
-                \r
-                \r
-                \r
-                MCarraytofill[0] = DeltaPhi;\r
-                MCarraytofill[1] = DeltaEta;\r
-                MCarraytofill[2] = ptDStar;\r
-                MCarraytofill[3] = etaDStar;\r
-                MCarraytofill[4] = ptHad;\r
-                               \r
-                Bool_t isDdaughter = kFALSE;\r
-                               if(fmontecarlo){\r
-                                       if(label<0 && fFullmode) ((TH2D*)fOutputMC->FindObject("TrackLabels"))->Fill(0.,NofTracks);\r
-                                       if(label>=0 && fFullmode) ((TH2D*)fOutputMC->FindObject("TrackLabels"))->Fill(1.,NofTracks);\r
-                                       if(label<0) continue; // skip track with wrong label\r
-                               }\r
-                               if(!fmixing && fReco){ // skip D* Daughetrs if it is reconstruced DStar\r
-                                       if(trackid == trackiddaugh0) continue;\r
-                                       if(trackid == trackiddaugh1) continue;\r
-                                       if(trackid == trackidsoftPi) continue;\r
-                               }\r
-                               if(!fmixing && !fReco && fmontecarlo){  // skip D* Daughetrs if it is Pure MCDStar\r
-                    // if(label == trackiddaugh0) continue;\r
-                    // if(label == trackiddaugh1) continue;\r
+           for(Int_t iTrack = 0; iTrack<NofTracks; iTrack++){ // looping on track candidates\r
+             Bool_t runcorrelation = fCorrelator->Correlate(iTrack);\r
+             if(!runcorrelation) continue;\r
+              \r
+             Double_t DeltaPhi = fCorrelator->GetDeltaPhi();\r
+             Double_t DeltaEta = fCorrelator->GetDeltaEta();\r
+             \r
+             AliReducedParticle * hadron = fCorrelator->GetAssociatedParticle();\r
+             if(!hadron) {/*cout << "No Hadron" << endl;*/ continue;}\r
+              \r
+             Double_t ptHad = hadron->Pt();\r
+             Double_t phiHad = hadron->Phi();\r
+             Double_t etaHad = hadron->Eta();\r
+             Int_t label = hadron->GetLabel();\r
+             Int_t trackid = hadron->GetID();\r
+             Double_t efficiency = hadron->GetWeight();\r
+              \r
+             weight = 1;\r
+             if(fUseEfficiencyCorrection && efficiency){\r
+               //weight = 1./efficiency;\r
+               weight = DmesonWeight * (1./efficiency);\r
+             }\r
+             // weight = DmesonWeight * (1./efficiency);\r
+              \r
+             phiHad = fCorrelator->SetCorrectPhiRange(phiHad);\r
+              \r
+              \r
+             if(fFullmode) ((TH2F*)fOutput->FindObject("WeightChecks"))->Fill(ptHad,efficiency);\r
+              \r
+             arraytofill[0] = DeltaPhi;\r
+             arraytofill[1] = DeltaEta;\r
+             arraytofill[2] = ptDStar;\r
+             arraytofill[3] = MultipOrCent;\r
+             arraytofill[4] = ptHad;\r
+             arraytofill[5] = zVtxPosition;\r
+              \r
+              \r
+              \r
+              \r
+             MCarraytofill[0] = DeltaPhi;\r
+             MCarraytofill[1] = DeltaEta;\r
+             MCarraytofill[2] = ptDStar;\r
+             MCarraytofill[3] = etaDStar;\r
+             MCarraytofill[4] = ptHad;\r
+             \r
+             Bool_t isDdaughter = kFALSE;\r
+             if(fmontecarlo){\r
+               if(label<0 && fFullmode) ((TH2D*)fOutputMC->FindObject("TrackLabels"))->Fill(0.,NofTracks);\r
+               if(label>=0 && fFullmode) ((TH2D*)fOutputMC->FindObject("TrackLabels"))->Fill(1.,NofTracks);\r
+               if(label<0) continue; // skip track with wrong label\r
+             }\r
+             if(!fmixing && fReco){ // skip D* Daughetrs if it is reconstruced DStar\r
+               if(trackid == trackiddaugh0) continue;\r
+               if(trackid == trackiddaugh1) continue;\r
+               if(trackid == trackidsoftPi) continue;\r
+             }\r
+             if(!fmixing && !fReco && fmontecarlo){  // skip D* Daughetrs if it is Pure MCDStar\r
+               //      if(label == trackiddaugh0) continue;\r
+               //      if(label == trackiddaugh1) continue;\r
                     // if(label == trackidsoftPi) continue;\r
-                                       Int_t hadronlabel = label;\r
-                                       for(Int_t k=0; k<4;k++){ // go back 4 generations and check the mothers\r
-                                               if(DStarLabel<0){ break;}\r
-                                               if(hadronlabel<0) { break;}\r
-                                               AliAODMCParticle* mcParticle = dynamic_cast<AliAODMCParticle*>(fmcArray->At(hadronlabel));\r
-                                               if(!mcParticle) {AliInfo("NO MC PARTICLE"); break;}\r
-                                               \r
-                                               hadronlabel = mcParticle->GetMother();\r
-                                               if(hadronlabel == DStarLabel) isDdaughter = kTRUE;\r
-                                       }\r
-                                       \r
-                                       \r
-                                       if(isDdaughter && fDebugLevel){\r
-                                               std::cout << "It is the D* daughter with label " << label << std::endl;\r
-                                               std::cout << "Daughter 0 label = " << trackiddaugh0 << std::endl;\r
-                                               std::cout << "Daughter 1 label = " << trackiddaugh1 << std::endl;\r
-                                               std::cout << "Soft pi label = " << trackidsoftPi << std::endl;\r
-                                       }\r
+               Int_t hadronlabel = label;\r
+               for(Int_t k=0; k<4;k++){ // go back 4 generations and check the mothers\r
+                 if(DStarLabel<0){ break;}\r
+                 if(hadronlabel<0) { break;}\r
+                 AliAODMCParticle* mcParticle = dynamic_cast<AliAODMCParticle*>(fmcArray->At(hadronlabel));\r
+                 if(!mcParticle) {AliInfo("NO MC PARTICLE"); break;}\r
+                 \r
+                 hadronlabel = mcParticle->GetMother();\r
+                 if(hadronlabel == DStarLabel) isDdaughter = kTRUE;\r
+               }\r
+               \r
+               \r
+               if(isDdaughter && fDebugLevel){\r
+                 std::cout << "It is the D* daughter with label " << label << std::endl;\r
+                 std::cout << "Daughter 0 label = " << trackiddaugh0 << std::endl;\r
+                 std::cout << "Daughter 1 label = " << trackiddaugh1 << std::endl;\r
+                 std::cout << "Soft pi label = " << trackidsoftPi << std::endl;\r
+               }\r
                     \r
                                        if(isDdaughter) continue; // skip if track is from DStar\r
                                }\r
@@ -939,6 +939,8 @@ void AliAnalysisTaskDStarCorrelations::UserExec(Option_t *){
                                        if(!fReco && TMath::Abs(etaHad)>0.9) continue; // makes sure you study the correlation on MC  truth only if particles are in acceptance\r
                     ((THnSparseF*)fOutputMC->FindObject("MCDStarCorrelationsDStarHadron"))->Fill(MCarraytofill);\r
                     \r
+                    delete[] PartSource;\r
+                    \r
                                }\r
                                if(isInPeak)  {\r
                     \r
@@ -1110,7 +1112,7 @@ void AliAnalysisTaskDStarCorrelations::DefineThNSparseForAnalysis(){
     \r
     \r
 }\r
-\r
+//__________________________________________________________________________________________________\r
 void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){\r
        \r
        Double_t Pi = TMath::Pi();\r
@@ -1160,6 +1162,12 @@ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){
        \r
        TH1F *RecoPtBkg = new TH1F("RecoPtBkg","RECO pt distribution side bands",50,0,50);\r
        if(!fmixing) fOutput->Add(RecoPtBkg);\r
+    \r
+    TH1D *MCtagPtDStarfromCharm = new TH1D("MCtagPtDStarfromCharm","RECO pt of MCtagged DStars from charm",50,0,50);\r
+    if(fmontecarlo) fOutputMC->Add(MCtagPtDStarfromCharm);\r
+    \r
+    TH1D *MCtagPtDStarfromBeauty = new TH1D("MCtagPtDStarfromBeauty","RECO pt of MCtagged DStars from beauty",50,0,50);\r
+    if(fmontecarlo) fOutputMC->Add(MCtagPtDStarfromBeauty);\r
        \r
        TH1F *MCtagPtDStar = new TH1F("MCtagPtDStar","RECO pt of MCtagged DStars side bands",50,0,50);\r
        if(!fmixing) fOutput->Add(MCtagPtDStar);\r
@@ -1420,9 +1428,11 @@ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){
        //Double_t Nevents[]={0,2*MaxNofEvents/10,4*MaxNofEvents/10,6*MaxNofEvents/10,8*MaxNofEvents/10,MaxNofEvents};\r
     // Double_t Nevents[]={0,2*MaxNofEvents/10,4*MaxNofEvents/10,6*MaxNofEvents/10,8*MaxNofEvents/10,MaxNofEvents};\r
     // Double_t * events = Nevents;\r
-    Double_t * events = new Double_t[2];\r
-    events[0] = 0;\r
-       events[1] = 1000000;\r
+    Double_t eventsv[] ={0,1000000};\r
+    //Double_t * events = new Double_t[2];\r
+   // events[0] = 0;\r
+//     events[1] = 1000000;\r
+    Double_t *events = eventsv;\r
     Int_t Nevents = 1000000;\r
     //  TH3D * EventsPerPoolBin = new TH3D("EventsPerPoolBin","Number of events in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins,Nevents,events);\r
     \r
@@ -1469,7 +1479,46 @@ void AliAnalysisTaskDStarCorrelations::DefineHistoForAnalysis(){
        \r
 }\r
 \r
+//__________________________________________________________________________________________________\r
+void AliAnalysisTaskDStarCorrelations::EnlargeDZeroMassWindow(){\r
+    \r
 \r
+       //Float_t* ptbins = fCuts->GetPtBinLimits();\r
+    \r
+    fD0Window = new Float_t[fNofPtBins];\r
+    \r
+    AliInfo("Enlarging the D0 mass windows from cut object\n"); \r
+    Int_t nvars = fCuts->GetNVars();\r
+    \r
+    Float_t** rdcutsvalmine;\r
+       rdcutsvalmine=new Float_t*[nvars];\r
+       for(Int_t iv=0;iv<nvars;iv++){\r
+               rdcutsvalmine[iv]=new Float_t[fNofPtBins];\r
+       }\r
+    \r
+    \r
+    for (Int_t k=0;k<nvars;k++){\r
+               for (Int_t j=0;j<fNofPtBins;j++){\r
+            \r
+            // enlarge D0 window\r
+            if(k==0)   {fD0Window[j] =fCuts->GetCutValue(0,j);\r
+                rdcutsvalmine[k][j] = 5.* fCuts->GetCutValue(0,j);\r
+                cout << "the set window = " << fD0Window[j] << " for ptbin " << j << endl;}\r
+                   else rdcutsvalmine[k][j] =fCuts->GetCutValue(k,j);\r
+                       \r
+            // set same windows\r
+            //rdcutsvalmine[k][j] =oldCuts->GetCutValue(k,j);\r
+        }\r
+       }\r
+       \r
+       fCuts->SetCuts(nvars,fNofPtBins,rdcutsvalmine);\r
+    \r
+   AliInfo("\n New windows set\n");     \r
+    fCuts->PrintAll();\r
+    \r
+    \r
+    \r
+}\r
 \r
 \r
 //____________________________  Run checks on event mixing ___________________________________________________\r
index 7cefe51b9aac5f7b39d8f4363300c5ef3e85b021..ade088c547c83dc7569d294c08edd5e6be710095 100644 (file)
@@ -40,7 +40,7 @@
 #include "AliHFCorrelator.h"\r
 #include <THnSparse.h>\r
 #include "AliAnalysisUtils.h"\r
-\r
+#include "AliVertexingHFUtils.h"\r
 class TParticle ;\r
 class TClonesArray ;\r
 class AliAODMCParticle;\r
@@ -73,6 +73,7 @@ class AliAnalysisTaskDStarCorrelations : public AliAnalysisTaskSE
   \r
   void DefineThNSparseForAnalysis();\r
   void DefineHistoForAnalysis();\r
+  void EnlargeDZeroMassWindow();\r
   \r
   \r
   // checker for event mixing\r
@@ -98,8 +99,8 @@ class AliAnalysisTaskDStarCorrelations : public AliAnalysisTaskSE
   void SetDim(){fDim = 4;\r
     fDMesonSigmas = new Float_t[4];}\r
   void SetDeffMapvsPt(TH1D * map){fDeffMapvsPt = map;}\r
-  void SetDeffMapvsPtvsMult(TH2D * map){fDeffMapvsPtvsMult = map;}\r
-  void SetDeffMapvsPtvsMultvsEta(TH3D * map){fDeffMapvsPtvsMultvsEta = map;}\r
+  void SetDeffMapvsPtvsMult(TH2D * map){fDeffMapvsPtvsMult = (TH2D*)map;}\r
+  void SetDeffMapvsPtvsMultvsEta(TH2D * map){fDeffMapvsPtvsEta = map;}\r
   void SetNofPhiBins(Int_t nbins){fPhiBins = nbins;}\r
   \r
   \r
@@ -124,12 +125,15 @@ private:
   Bool_t fReco; // use reconstruction or MC truth\r
   Bool_t fUseEfficiencyCorrection; // boolean variable to use or not the efficiency correction\r
   Bool_t fUseDmesonEfficiencyCorrection; // boolean flag for the use of Dmeson efficiency correction\r
+    Bool_t fUseCentrality;// boolean to switch in between centrality or multiplicity\r
   Int_t fPhiBins;\r
   Int_t fEvents; //! number of event\r
   Int_t fDebugLevel; //! debug level\r
   Int_t fDisplacement; // set 0 for no displacement cut, 1 for absolute d0, 2 for d0/sigma_d0\r
   Int_t fDim;//\r
+    Int_t fNofPtBins;\r
   Float_t *fDMesonSigmas;//[fDim]\r
+ Float_t * fD0Window;  //[fNofPtBins]\r
   \r
   \r
   \r
@@ -142,9 +146,9 @@ private:
   \r
   TH1D * fDeffMapvsPt; // histo for Deff mappin\r
   TH2D * fDeffMapvsPtvsMult; // histo for Deff mappin\r
-  TH3D * fDeffMapvsPtvsMultvsEta; // histo for Deff mappin\r
+  TH2D * fDeffMapvsPtvsEta; // histo for Deff mappin\r
   \r
-  ClassDef(AliAnalysisTaskDStarCorrelations,4); // class for D meson correlations\r
+  ClassDef(AliAnalysisTaskDStarCorrelations,5); // class for D meson correlations\r
   \r
 };\r
 \r
index 018df5b5e979b908f874a9e8b25eea7caf7148cb..9217091586961f80db2146c0641cf20304703358 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-2010, 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.                  *
- **************************************************************************/
-
-/* $Id$ */
-
-////////////////////////////////////////////////////////////////////////
-//
-// Base class for cuts on Associated tracks for HF Correlation analysis
-//
-// Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch
-////////////////////////////////////////////////////////////////////////
-#include <Riostream.h>
-#include "AliHFAssociatedTrackCuts.h"
-#include "AliAODPidHF.h"
-#include "AliESDtrackCuts.h"
-#include "AliESDtrack.h"
-#include "AliESDVertex.h"
-#include "AliAODv0.h"
-#include "AliAODVertex.h"
-#include "AliAODMCParticle.h"
-#include "AliAnalysisManager.h"
-#include "AliInputEventHandler.h"
-#include "TString.h"
-
-using std::cout;
-using std::endl;
-
-ClassImp(AliHFAssociatedTrackCuts)
-
-//--------------------------------------------------------------------------
-AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts():
-AliAnalysisCuts(),
-fESDTrackCuts(0),
-fPidObj(0),
-  fEffWeights(0),
-fPoolMaxNEvents(0), 
-fPoolMinNTracks(0), 
-fMinEventsToMix(0),
-fNzVtxBins(0), 
-fNzVtxBinsDim(0), 
-fZvtxBins(0), 
-fNCentBins(0), 
-fNCentBinsDim(0), 
-fCentBins(0),
-
-fNofMCEventType(0),
-fMCEventType(0),
-
-fNTrackCuts(0),
-fAODTrackCuts(0),
-fTrackCutsNames(0),
-fNvZeroCuts(0),
-fAODvZeroCuts(0),
-fvZeroCutsNames(0),
-fBit(-1),
-fCharge(0),
-fDescription("")
-
-{
-       //
-       //default constructor
-       //
-       //
-       //default constructor
-       //
-       
-}
-
-//--------------------------------------------------------------------------
-AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const char* name, const char* title):
-AliAnalysisCuts(name,title),
-fESDTrackCuts(0),
-fPidObj(0),
-fEffWeights(0),
-fPoolMaxNEvents(0), 
-fPoolMinNTracks(0), 
-fMinEventsToMix(0),
-fNzVtxBins(0), 
-fNzVtxBinsDim(0), 
-fZvtxBins(0), 
-fNCentBins(0), 
-fNCentBinsDim(0), 
-fCentBins(0),
-
-fNofMCEventType(0),
-fMCEventType(0),
-
-fNTrackCuts(0),
-fAODTrackCuts(0),
-fTrackCutsNames(0),
-fNvZeroCuts(0),
-fAODvZeroCuts(0),
-fvZeroCutsNames(0),
-fBit(-1),
-fCharge(0),
-fDescription("")
-
-{
-       //
-       //default constructor
-       //
-       
-}
-//--------------------------------------------------------------------------
-AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts &source) :
-AliAnalysisCuts(source),
-fESDTrackCuts(source.fESDTrackCuts),
-fPidObj(source.fPidObj),
-fEffWeights(source.fEffWeights),
-fPoolMaxNEvents(source.fPoolMaxNEvents), 
-fPoolMinNTracks(source.fPoolMinNTracks), 
-fMinEventsToMix(source.fMinEventsToMix),
-fNzVtxBins(source.fNzVtxBins), 
-fNzVtxBinsDim(source.fNzVtxBinsDim), 
-fZvtxBins(source.fZvtxBins), 
-fNCentBins(source.fNCentBins), 
-fNCentBinsDim(source.fNCentBinsDim), 
-fCentBins(source.fCentBins),
-
-fNofMCEventType(source.fNofMCEventType),
-fMCEventType(source.fMCEventType),
-
-fNTrackCuts(source.fNTrackCuts),
-fAODTrackCuts(source.fAODTrackCuts),
-fTrackCutsNames(source.fTrackCutsNames),
-fNvZeroCuts(source.fNvZeroCuts),
-fAODvZeroCuts(source.fAODvZeroCuts),
-fvZeroCutsNames(source.fvZeroCutsNames),
-fBit(source.fBit),
-fCharge(source.fCharge),
-fDescription(source.fDescription)
-{
-       //
-       // copy constructor
-       //
-       
-
-        AliInfo("AliHFAssociatedTrackCuts::Copy constructor ");
-       if(source.fESDTrackCuts) AddTrackCuts(source.fESDTrackCuts);
-       if(source.fAODTrackCuts) SetAODTrackCuts(source.fAODTrackCuts);
-       if(source.fAODvZeroCuts) SetAODvZeroCuts(source.fAODvZeroCuts);
-       if(source.fPidObj) SetPidHF(source.fPidObj);
-       if(source.fEffWeights) SetEfficiencyWeightMap(source.fEffWeights);
-}
-//--------------------------------------------------------------------------
-AliHFAssociatedTrackCuts &AliHFAssociatedTrackCuts::operator=(const AliHFAssociatedTrackCuts &source)
-{
-       //
-       // assignment operator
-       //      
-       if(&source == this) return *this;
-       
-       AliAnalysisCuts::operator=(source);
-       fESDTrackCuts=source.fESDTrackCuts;
-       fPidObj=source.fPidObj;
-       fEffWeights=source.fEffWeights;
-       fNTrackCuts=source.fNTrackCuts;
-       fAODTrackCuts=source.fAODTrackCuts;
-       fTrackCutsNames=source.fTrackCutsNames;
-       fNvZeroCuts=source.fNvZeroCuts;
-       fAODvZeroCuts=source.fAODvZeroCuts;
-       fvZeroCutsNames=source.fvZeroCutsNames;
-       fBit=source.fBit;
-       fCharge=source.fCharge;
-       
-       return *this;
-
-}
-
-
-//--------------------------------------------------------------------------
-AliHFAssociatedTrackCuts::~AliHFAssociatedTrackCuts()
-{
-       if(fESDTrackCuts) {delete fESDTrackCuts; fESDTrackCuts = 0;}
-       if(fPidObj) {delete fPidObj; fPidObj = 0;}
-       if(fEffWeights){delete fEffWeights;fEffWeights=0;}
-       if(fZvtxBins) {delete[] fZvtxBins; fZvtxBins=0;} 
-       if(fCentBins) {delete[] fCentBins; fCentBins=0;}
-       if(fAODTrackCuts) {delete[] fAODTrackCuts; fAODTrackCuts=0;}
-       if(fTrackCutsNames) {delete[] fTrackCutsNames; fTrackCutsNames=0;}
-       if(fAODvZeroCuts){delete[] fAODvZeroCuts; fAODvZeroCuts=0;}
-       if(fvZeroCutsNames) {delete[] fvZeroCutsNames; fvZeroCutsNames=0;}
-
-       
-}
-//--------------------------------------------------------------------------
-Bool_t AliHFAssociatedTrackCuts::IsInAcceptance()
-{
-       printf("Careful: method AliHFAssociatedTrackCuts::IsInAcceptance is not implemented yet \n");
-       return kFALSE;
-}
-//--------------------------------------------------------------------------
-Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track,const AliESDVertex *primary,const Double_t magfield)
-{
-  
-  AliESDtrack esdtrack(track);
-  if(primary){// needed to calculate impact parameters
-    // needed to calculate the impact parameters
-    esdtrack.RelateToVertex(primary,magfield,3.); 
-  }
-  // set the TPC cluster info
-  esdtrack.SetTPCClusterMap(track->GetTPCClusterMap());
-  esdtrack.SetTPCSharedMap(track->GetTPCSharedMap());
-  esdtrack.SetTPCPointsF(track->GetTPCNclsF());
-  
-  if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE;
-  
-  if(fBit>-1 && !track->TestFilterBit(fBit)) return kFALSE; // check the filter bit
-  
-  return kTRUE;
-       
-}
-
-//--------------------------------------------------------------------------
-Bool_t AliHFAssociatedTrackCuts::CheckHadronKinematic(Double_t pt, Double_t d0) 
-{
-       
-       
-       
-       if(pt < fAODTrackCuts[0]) return kFALSE;
-       if(pt > fAODTrackCuts[1]) return kFALSE;
-       if(d0 < fAODTrackCuts[2]) return kFALSE;
-       if(d0 > fAODTrackCuts[3]) return kFALSE;
-       
-       return kTRUE;
-
-       
-}
-//--------------------------------------------------------------------------
-
-Bool_t AliHFAssociatedTrackCuts::Charge(Short_t charge, AliAODTrack* track) 
-{// charge is the charge to compare to (for example, a daughter of a D meson)
-       
-       if(!fCharge) return kTRUE; // if fCharge is set to 0 (no selection on the charge), returns always true
-       if(track->Charge()!= fCharge*charge) return kFALSE;
-       return kTRUE;
-}
-
-//--------------------------------------------------------------------------
-Bool_t AliHFAssociatedTrackCuts::CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method)
-{
-       Bool_t isKaon = kFALSE;
-       
-       if(useMc) { // on MC
-               Int_t hadLabel = track->GetLabel();
-               if(hadLabel < 0) return kFALSE;
-               AliAODMCParticle* hadron = dynamic_cast<AliAODMCParticle*>(mcArray->At(hadLabel)); 
-               if(hadron){
-                 Int_t pdg = TMath::Abs(hadron->GetPdgCode()); 
-                 if (pdg == 321) isKaon = kTRUE;
-               }
-       }
-       
-       if(!useMc) { // on DATA
-          switch(method) {
-         case(1): {
-               Bool_t isKTPC=kFALSE;
-               Bool_t isPiTPC=kFALSE;
-               Bool_t isPTPC=kFALSE;
-               Bool_t isKTOF=kFALSE;
-               Bool_t isPiTOF=kFALSE;
-               Bool_t isPTOF=kFALSE;
-               
-               Bool_t KaonHyp = kFALSE;
-               Bool_t PionHyp = kFALSE;
-               Bool_t ProtonHyp = kFALSE;
-               
-               if(fPidObj->CheckStatus(track,"TOF")) {
-                       isKTOF=fPidObj->IsKaonRaw(track,"TOF");
-                       isPiTOF=fPidObj->IsPionRaw(track,"TOF");
-                       isPTOF=fPidObj->IsProtonRaw(track,"TOF");
-               }
-               if(fPidObj->CheckStatus(track,"TPC")){
-                       isKTPC=fPidObj->IsKaonRaw(track,"TPC");
-                       isPiTPC=fPidObj->IsPionRaw(track,"TPC");
-                       isPTPC=fPidObj->IsProtonRaw(track,"TPC");               
-               }
-               
-               if (isKTOF && isKTPC) KaonHyp = kTRUE;
-               if (isPiTOF && isPiTPC) PionHyp = kTRUE;
-               if (isPTOF && isPTPC) ProtonHyp = kTRUE;
-               
-               if(KaonHyp && !PionHyp && !ProtonHyp) isKaon = kTRUE; 
-               break;
-             }
-      case(2): {
-               if(fPidObj->MakeRawPid(track,3)>=1) isKaon = kTRUE;
-               break;
-             }
-      }
-       }
-       
-       return isKaon;
-       
-}
-//--------------------------------------------------------------------------
-Bool_t AliHFAssociatedTrackCuts::IsKZeroSelected(AliAODv0 *vzero, AliAODVertex *vtx1)
-{
-       
-       if(vzero->DcaV0Daughters()>fAODvZeroCuts[0]) return kFALSE;
-       if(vzero->Chi2V0()>fAODvZeroCuts[1]) return kFALSE;
-       if(vzero->DecayLength(vtx1) < fAODvZeroCuts[2]) return kFALSE;
-       if(vzero->DecayLength(vtx1) > fAODvZeroCuts[3]) return kFALSE;
-       if(vzero->OpenAngleV0() > fAODvZeroCuts[4]) return kFALSE;
-       if(vzero->Pt() < fAODvZeroCuts[5]) return kFALSE;
-       if(TMath::Abs(vzero->Eta()) > fAODvZeroCuts[6]) return kFALSE;
-
-       
-       return kTRUE;
-}
-//--------------------------------------------------------------------------
-Bool_t *AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray){
-  // Check origin in MC
-
-  AliAODMCParticle* mcParticle;
-  Int_t pdgCode = -1;
-       
-  Bool_t isCharmy = kFALSE;
-  Bool_t isBeauty = kFALSE;
-  Bool_t isD = kFALSE;
-  Bool_t isB = kFALSE;
-    
-     Bool_t *originvect = new Bool_t[4];
-    
-    originvect[0] = kFALSE;
-       originvect[1] = kFALSE;
-       originvect[2] = kFALSE;
-       originvect[3] = kFALSE;
-
-       if (label<0) return originvect;
-  
-       while(pdgCode!=2212){ // loops back to the collision to check the particle source
-
-    mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));
-    if(!mcParticle) {AliError("NO MC PARTICLE"); break;}
-    pdgCode =  TMath::Abs(mcParticle->GetPdgCode());
-
-    label = mcParticle->GetMother();
-
-
-    if((pdgCode>=400 && pdgCode <500) || (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;
-    if((pdgCode>=500 && pdgCode <600) || (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}
-
-
-    if(pdgCode == 4) isCharmy = kTRUE;
-    if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}
-    if(label<0) break;
-
-  }
-
-       
-       originvect[0] = isCharmy;
-       originvect[1] = isBeauty;
-       originvect[2] = isD;
-       originvect[3] = isB;
-    
-  return originvect;
-}
-
-//--------------------------------------------------------------------------
-Bool_t AliHFAssociatedTrackCuts::InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const {
-       //
-       // Calculates invmass of track+D0 and rejects if compatible with D*
-       // (to remove pions from D*)
-       // 
-       Double_t nsigma = 3.;
-       
-       Double_t mD0, mD0bar;
-       d->InvMassD0(mD0,mD0bar);
-       
-       Double_t invmassDstar1 = 0, invmassDstar2 = 0; 
-       Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0
-       Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar
-       Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px())
-       +(d->Py()+track->Py())*(d->Py()+track->Py())
-       +(d->Pz()+track->Pz())*(d->Pz()+track->Pz());
-       
-       switch(hypD0) {
-               case 1:
-                       invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
-                       if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
-                       break;
-               case 2:
-                       invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
-                       if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
-                       break;
-               case 3:
-                       invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);
-                       invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);
-                       if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
-                       if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;
-                       break;
-       }
-       
-       return kTRUE;
-}
-//________________________________________________________
-void AliHFAssociatedTrackCuts::SetMCEventTypes(Int_t *MCEventTypeArray)
-// set the array of event types you want to process in MonteCarlo (gluon splitting, pair production etc.)
-{
-       if(!fMCEventType) fMCEventType = new Int_t[fNofMCEventType];
-       
-       for(Int_t k=0; k<fNofMCEventType; k++){
-               fMCEventType[k] = MCEventTypeArray[k];
-       }
-       return; 
-}
-
-//________________________________________________________
-void AliHFAssociatedTrackCuts::SetAODTrackCuts(Float_t *cutsarray)
-{
-       if(!fAODTrackCuts) fAODTrackCuts = new Float_t[fNTrackCuts];
-       for(Int_t i =0; i<fNTrackCuts; i++){
-               fAODTrackCuts[i] = cutsarray[i];
-       }
-       SetTrackCutsNames();
-       return;
-}
-//________________________________________________________
-void AliHFAssociatedTrackCuts::SetTrackCutsNames(/*TString *namearray*/){
-       
-       fTrackCutsNames = new TString[4];
-       fTrackCutsNames[0]= "associated track:: pt min [GeV/c]................: ";
-       fTrackCutsNames[1]= "associated track:: pt max [GeV/c]................: ";
-       fTrackCutsNames[2]= "associated track:: d0 min [cm]...................: ";
-       fTrackCutsNames[3]= "associated track:: d0 max [cm]...................: ";
-       
-
-       
-       return;
-}
-//--------------------------------------------------------------------------
-void AliHFAssociatedTrackCuts::SetAODvZeroCuts(Float_t *cutsarray)
-{
-       
-
-       if(!fAODvZeroCuts) fAODvZeroCuts = new Float_t[fNvZeroCuts];
-       for(Int_t i =0; i<fNvZeroCuts; i++){
-               fAODvZeroCuts[i] = cutsarray[i];
-       }
-       SetvZeroCutsNames();
-       return;
-}
-//--------------------------------------------------------------------------
-void AliHFAssociatedTrackCuts::SetvZeroCutsNames(/*TString *namearray*/){
-       
-       fvZeroCutsNames = new TString[7];
-       fvZeroCutsNames[0] = "vZero:: max DCA between two daughters [cm].......: ";
-       fvZeroCutsNames[1] = "vZero:: max fit Chi Square between two daughters.: ";
-       fvZeroCutsNames[2] = "vZero:: min decay length [cm]....................: ";
-       fvZeroCutsNames[3] = "vZero:: max decay length [cm]....................: ";
-       fvZeroCutsNames[4] = "vZero:: max opening angle between daughters [rad]: ";
-       fvZeroCutsNames[5] = "vZero:: pt min [Gev/c]...........................: ";
-       fvZeroCutsNames[6] = "vZero:: |Eta| range <............................: ";
-       
-       
-       return;
-}
-
-//--------------------------------------------------------------------------
-void AliHFAssociatedTrackCuts::SetPidAssociated()
-{
-  //setting PidResponse
-  if(fPidObj->GetOldPid()==kFALSE && fPidObj->GetPidResponse()==0x0){
-    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
-    AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
-    AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
-    fPidObj->SetPidResponse(pidResp);
-  }
-}
-
-void AliHFAssociatedTrackCuts::Print(Option_t *option) const
-{
-  /// overloaded from TObject: print info
-  if (strcmp(option, "parameters")==0) {
-    PrintPoolParameters();
-    return;
-  } else if (strcmp(option, "selectedMC")==0) {
-    PrintSelectedMCevents();
-    return;
-  }
-  PrintAll();
-}
-
-//--------------------------------------------------------------------------
-void AliHFAssociatedTrackCuts::PrintAll() const
-{
-       
-       if(fDescription){
-         printf("=================================================");
-         printf("\nAdditional description\n");
-         std::cout << fDescription << std::endl;
-         printf("\n");
-       }
-       printf("\n=================================================");
-       if(fESDTrackCuts){
-         printf("\nCuts for the associated track: \n \n");
-         
-         printf("ITS Refit........................................: %s\n",fESDTrackCuts->GetRequireITSRefit() ? "Yes" : "No");
-         printf("TPC Refit........................................: %s\n",fESDTrackCuts->GetRequireTPCRefit() ? "Yes" : "No");
-         printf("ITS SA...........................................: %s\n",fESDTrackCuts->GetRequireITSStandAlone() ? "Yes" : "No");
-         printf("TPC SA...........................................: %s\n",fESDTrackCuts->GetRequireTPCStandAlone() ? "Yes" : "No");
-         printf("Min number of ITS clusters.......................: %d\n",fESDTrackCuts->GetMinNClustersITS());
-         printf("Min number of TPC clusters.......................: %d\n",fESDTrackCuts->GetMinNClusterTPC());
-         Int_t spd = fESDTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD);
-         if(spd==0) std::cout <<  "SPD..............................................: kOff"  << std::endl;
-         if(spd==1) std::cout <<  "SPD..............................................: kNone"  << std::endl;
-         if(spd==2) std::cout <<  "SPD..............................................: kAny"  << std::endl;
-         if(spd==3) std::cout <<  "SPD..............................................: kFirst"  << std::endl;
-         if(spd==4) std::cout <<  "SPD..............................................: kOnlyFirst"  << std::endl;
-         if(spd==5) std::cout <<  "SPD..............................................: kSecond"  << std::endl;
-         if(spd==6) std::cout <<  "SPD..............................................: kOnlySecond"  << std::endl;
-         if(spd==7) std::cout <<  "SPD..............................................: kBoth"  << std::endl;
-       }
-       else printf("\nNo Cuts for Associated Tracks\n");
-       std::cout <<  "Filter Bit.......................................: " << fBit  << std::endl;
-       std::cout <<  "Charge...........................................: " << fCharge  << std::endl;
-       
-       if(fAODTrackCuts){
-         for(Int_t j=0;j<fNTrackCuts;j++){
-           std::cout << fTrackCutsNames[j] << fAODTrackCuts[j] << std::endl;
-         }
-       }
-
-       if(fAODvZeroCuts){
-         printf("\n");
-         printf("=================================================");
-         printf("\nCuts for the K0 candidates: \n \n");
-         for(Int_t k=0;k<fNvZeroCuts;k++){
-           std::cout << fvZeroCutsNames[k] <<  fAODvZeroCuts[k] << std::endl;
-         }
-       }
-       else printf("\nNo Cuts for the K0 candidates\n");
-       std::cout << " " << std::endl;
-       PrintPoolParameters();
-       PrintSelectedMCevents();
-
-}
-
-//--------------------------------------------------------------------------
-void AliHFAssociatedTrackCuts::PrintPoolParameters() const
-{   
-       printf("=================================================");
-       printf("\nEvent Pool settings: \n \n");
-       
-       printf("Number of zVtx Bins: %d\n", fNzVtxBins);
-       printf("\nzVtx Bins:\n");
-       //Double_t zVtxbinLims[fNzVtxBins+1] = fNzVtxBins;
-       for(Int_t k=0; k<fNzVtxBins; k++){
-               printf("Bin %d..............................................: %.1f - %.1f cm\n", k, fZvtxBins[k], fZvtxBins[k+1]);      
-       }
-       printf("\n");
-       printf("\nNumber of Centrality(multiplicity) Bins: %d\n", fNCentBins);
-       printf("\nCentrality(multiplicity) Bins:\n");
-       for(Int_t k=0; k<fNCentBins; k++){
-               printf("Bin %d..............................................: %.1f - %.1f\n", k, fCentBins[k], fCentBins[k+1]);
-       }
-
-       
-       
-}
-Double_t AliHFAssociatedTrackCuts::GetTrackWeight(Double_t pt, Double_t eta,Double_t zvtx){
-  if(!fEffWeights)return 1.;
-  
-  Int_t bin=fEffWeights->FindBin(pt,eta,zvtx);
-  if(fEffWeights->IsBinUnderflow(bin)||fEffWeights->IsBinOverflow(bin))return 1.;
-  return fEffWeights->GetBinContent(bin);
-
-}
-//--------------------------------------------------------------------------
-void AliHFAssociatedTrackCuts::PrintSelectedMCevents() const
-{
-       printf("\n=================================================");
-       
-       printf("\nSelected MC events: \n \n");
-       printf("Number of selected events: %d\n",fNofMCEventType);
-       
-       for(Int_t k=0; k<fNofMCEventType; k++){
-       if(fMCEventType[k]==28) printf("=> Flavour excitation \n");     
-       if(fMCEventType[k]==53) printf("=> Pair creation \n");  
-       if(fMCEventType[k]==68) printf("=> Gluon splitting \n");        
-       }
-       
-       printf("\n");
-       for(Int_t k=0; k<fNofMCEventType; k++){
-               printf("MC process code %d \n",fMCEventType[k]);                
-       }
-       
-       printf("\n");
-       
-       
-               
-       
-}
-
-
+/**************************************************************************\r
+ * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id$ */\r
+\r
+////////////////////////////////////////////////////////////////////////\r
+//\r
+// Base class for cuts on Associated tracks for HF Correlation analysis\r
+//\r
+// Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch\r
+////////////////////////////////////////////////////////////////////////\r
+#include <Riostream.h>\r
+#include "AliHFAssociatedTrackCuts.h"\r
+#include "AliAODPidHF.h"\r
+#include "AliESDtrackCuts.h"\r
+#include "AliESDtrack.h"\r
+#include "AliESDVertex.h"\r
+#include "AliAODv0.h"\r
+#include "AliAODVertex.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliInputEventHandler.h"\r
+#include "TString.h"\r
+\r
+using std::cout;\r
+using std::endl;\r
+\r
+ClassImp(AliHFAssociatedTrackCuts)\r
+\r
+//--------------------------------------------------------------------------\r
+AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts():\r
+AliAnalysisCuts(),\r
+fESDTrackCuts(0),\r
+fPidObj(0),\r
+  fEffWeights(0),\r
+\r
+fTrigEffWeightsvspt(0),\r
+fTrigEffWeightsvsptB(0),\r
+fTrigEffWeights(0),\r
+fTrigEffWeightsB(0),\r
+fPoolMaxNEvents(0), \r
+fPoolMinNTracks(0), \r
+fMinEventsToMix(0),\r
+fNzVtxBins(0), \r
+fNzVtxBinsDim(0), \r
+fZvtxBins(0), \r
+fNCentBins(0), \r
+fNCentBinsDim(0), \r
+fCentBins(0),\r
+\r
+fNofMCEventType(0),\r
+fMCEventType(0),\r
+\r
+fNTrackCuts(0),\r
+fAODTrackCuts(0),\r
+fTrackCutsNames(0),\r
+fNvZeroCuts(0),\r
+fAODvZeroCuts(0),\r
+fvZeroCutsNames(0),\r
+fBit(-1),\r
+fCharge(0),\r
+fDescription("")\r
+\r
+{\r
+       //\r
+       //default constructor\r
+       //\r
+       //\r
+       //default constructor\r
+       //\r
+       \r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const char* name, const char* title):\r
+AliAnalysisCuts(name,title),\r
+fESDTrackCuts(0),\r
+fPidObj(0),\r
+fEffWeights(0),\r
+fTrigEffWeightsvspt(0),\r
+fTrigEffWeightsvsptB(0),\r
+fTrigEffWeights(0),\r
+fTrigEffWeightsB(0),\r
+fPoolMaxNEvents(0), \r
+fPoolMinNTracks(0), \r
+fMinEventsToMix(0),\r
+fNzVtxBins(0), \r
+fNzVtxBinsDim(0), \r
+fZvtxBins(0), \r
+fNCentBins(0), \r
+fNCentBinsDim(0), \r
+fCentBins(0),\r
+\r
+fNofMCEventType(0),\r
+fMCEventType(0),\r
+\r
+fNTrackCuts(0),\r
+fAODTrackCuts(0),\r
+fTrackCutsNames(0),\r
+fNvZeroCuts(0),\r
+fAODvZeroCuts(0),\r
+fvZeroCutsNames(0),\r
+fBit(-1),\r
+fCharge(0),\r
+fDescription("")\r
+\r
+{\r
+       //\r
+       //default constructor\r
+       //\r
+       \r
+}\r
+//--------------------------------------------------------------------------\r
+AliHFAssociatedTrackCuts::AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts &source) :\r
+AliAnalysisCuts(source),\r
+fESDTrackCuts(source.fESDTrackCuts),\r
+fPidObj(source.fPidObj),\r
+fEffWeights(source.fEffWeights),\r
+fTrigEffWeightsvspt(source.fTrigEffWeightsvspt),\r
+fTrigEffWeightsvsptB(source.fTrigEffWeightsvsptB),\r
+fTrigEffWeights(source.fTrigEffWeights),\r
+fTrigEffWeightsB(source.fTrigEffWeightsB),\r
+\r
+fPoolMaxNEvents(source.fPoolMaxNEvents), \r
+fPoolMinNTracks(source.fPoolMinNTracks), \r
+fMinEventsToMix(source.fMinEventsToMix),\r
+fNzVtxBins(source.fNzVtxBins), \r
+fNzVtxBinsDim(source.fNzVtxBinsDim), \r
+fZvtxBins(source.fZvtxBins), \r
+fNCentBins(source.fNCentBins), \r
+fNCentBinsDim(source.fNCentBinsDim), \r
+fCentBins(source.fCentBins),\r
+\r
+fNofMCEventType(source.fNofMCEventType),\r
+fMCEventType(source.fMCEventType),\r
+\r
+fNTrackCuts(source.fNTrackCuts),\r
+fAODTrackCuts(source.fAODTrackCuts),\r
+fTrackCutsNames(source.fTrackCutsNames),\r
+fNvZeroCuts(source.fNvZeroCuts),\r
+fAODvZeroCuts(source.fAODvZeroCuts),\r
+fvZeroCutsNames(source.fvZeroCutsNames),\r
+fBit(source.fBit),\r
+fCharge(source.fCharge),\r
+fDescription(source.fDescription)\r
+{\r
+       //\r
+       // copy constructor\r
+       //\r
+       \r
+\r
+        AliInfo("AliHFAssociatedTrackCuts::Copy constructor ");\r
+       if(source.fESDTrackCuts) AddTrackCuts(source.fESDTrackCuts);\r
+       if(source.fAODTrackCuts) SetAODTrackCuts(source.fAODTrackCuts);\r
+       if(source.fAODvZeroCuts) SetAODvZeroCuts(source.fAODvZeroCuts);\r
+       if(source.fPidObj) SetPidHF(source.fPidObj);\r
+       if(source.fEffWeights) SetEfficiencyWeightMap(source.fEffWeights);\r
+    if(source.fTrigEffWeightsvspt) SetTriggerEffWeightMapvspt(source.fTrigEffWeightsvspt);\r
+    if(source.fTrigEffWeightsvsptB) SetTriggerEffWeightMapvsptB(source.fTrigEffWeightsvspt);\r
+    if(source.fTrigEffWeights) SetTriggerEffWeightMap(source.fTrigEffWeights);\r
+    if(source.fTrigEffWeightsB)SetTriggerEffWeightMapB(source.fTrigEffWeightsB);\r
+   \r
+    \r
+}\r
+//--------------------------------------------------------------------------\r
+AliHFAssociatedTrackCuts &AliHFAssociatedTrackCuts::operator=(const AliHFAssociatedTrackCuts &source)\r
+{\r
+       //\r
+       // assignment operator\r
+       //      \r
+       if(&source == this) return *this;\r
+       \r
+       AliAnalysisCuts::operator=(source);\r
+       fESDTrackCuts=source.fESDTrackCuts;\r
+       fPidObj=source.fPidObj;\r
+       fEffWeights=source.fEffWeights;\r
+    fTrigEffWeightsvspt=source.fTrigEffWeightsvspt;\r
+       fTrigEffWeightsvsptB=source.fTrigEffWeightsvsptB;\r
+    fTrigEffWeights=source.fTrigEffWeights;\r
+       fTrigEffWeightsB=source.fTrigEffWeightsB;\r
+       fNTrackCuts=source.fNTrackCuts;\r
+       fAODTrackCuts=source.fAODTrackCuts;\r
+       fTrackCutsNames=source.fTrackCutsNames;\r
+       fNvZeroCuts=source.fNvZeroCuts;\r
+       fAODvZeroCuts=source.fAODvZeroCuts;\r
+       fvZeroCutsNames=source.fvZeroCutsNames;\r
+       fBit=source.fBit;\r
+       fCharge=source.fCharge;\r
+       \r
+       return *this;\r
+\r
+}\r
+\r
+\r
+//--------------------------------------------------------------------------\r
+AliHFAssociatedTrackCuts::~AliHFAssociatedTrackCuts()\r
+{\r
+       if(fESDTrackCuts) {delete fESDTrackCuts; fESDTrackCuts = 0;}\r
+       if(fPidObj) {delete fPidObj; fPidObj = 0;}\r
+       if(fEffWeights){delete fEffWeights;fEffWeights=0;}\r
+    if(fTrigEffWeightsvspt){delete fTrigEffWeightsvspt;fTrigEffWeightsvspt=0;}\r
+       if(fTrigEffWeightsvsptB){delete fTrigEffWeightsvsptB;fTrigEffWeightsvsptB=0;}\r
+    if(fTrigEffWeights){delete fTrigEffWeights;fTrigEffWeights=0;}\r
+       if(fTrigEffWeightsB){delete fTrigEffWeightsB;fTrigEffWeightsB=0;}\r
+       if(fZvtxBins) {delete[] fZvtxBins; fZvtxBins=0;} \r
+       if(fCentBins) {delete[] fCentBins; fCentBins=0;}\r
+       if(fAODTrackCuts) {delete[] fAODTrackCuts; fAODTrackCuts=0;}\r
+       if(fTrackCutsNames) {delete[] fTrackCutsNames; fTrackCutsNames=0;}\r
+       if(fAODvZeroCuts){delete[] fAODvZeroCuts; fAODvZeroCuts=0;}\r
+       if(fvZeroCutsNames) {delete[] fvZeroCutsNames; fvZeroCutsNames=0;}\r
+\r
+       \r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t AliHFAssociatedTrackCuts::IsInAcceptance()\r
+{\r
+       printf("Careful: method AliHFAssociatedTrackCuts::IsInAcceptance is not implemented yet \n");\r
+       return kFALSE;\r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t AliHFAssociatedTrackCuts::IsHadronSelected(AliAODTrack * track,const AliESDVertex *primary,const Double_t magfield)\r
+{\r
+  \r
+  AliESDtrack esdtrack(track);\r
+  if(primary){// needed to calculate impact parameters\r
+    // needed to calculate the impact parameters\r
+    esdtrack.RelateToVertex(primary,magfield,3.); \r
+  }\r
+  // set the TPC cluster info\r
+  esdtrack.SetTPCClusterMap(track->GetTPCClusterMap());\r
+  esdtrack.SetTPCSharedMap(track->GetTPCSharedMap());\r
+  esdtrack.SetTPCPointsF(track->GetTPCNclsF());\r
+  \r
+  if(!fESDTrackCuts->IsSelected(&esdtrack)) return kFALSE;\r
+  \r
+  if(fBit>-1 && !track->TestFilterBit(fBit)) return kFALSE; // check the filter bit\r
+  \r
+  return kTRUE;\r
+       \r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+Bool_t AliHFAssociatedTrackCuts::CheckHadronKinematic(Double_t pt, Double_t d0) \r
+{\r
+       \r
+       \r
+       \r
+       if(pt < fAODTrackCuts[0]) return kFALSE;\r
+       if(pt > fAODTrackCuts[1]) return kFALSE;\r
+       if(d0 < fAODTrackCuts[2]) return kFALSE;\r
+       if(d0 > fAODTrackCuts[3]) return kFALSE;\r
+       \r
+       return kTRUE;\r
+\r
+       \r
+}\r
+//--------------------------------------------------------------------------\r
+\r
+Bool_t AliHFAssociatedTrackCuts::Charge(Short_t charge, AliAODTrack* track) \r
+{// charge is the charge to compare to (for example, a daughter of a D meson)\r
+       \r
+       if(!fCharge) return kTRUE; // if fCharge is set to 0 (no selection on the charge), returns always true\r
+       if(track->Charge()!= fCharge*charge) return kFALSE;\r
+       return kTRUE;\r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+Bool_t AliHFAssociatedTrackCuts::CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method)\r
+{\r
+       Bool_t isKaon = kFALSE;\r
+       \r
+       if(useMc) { // on MC\r
+               Int_t hadLabel = track->GetLabel();\r
+               if(hadLabel < 0) return kFALSE;\r
+               AliAODMCParticle* hadron = dynamic_cast<AliAODMCParticle*>(mcArray->At(hadLabel)); \r
+               if(hadron){\r
+                 Int_t pdg = TMath::Abs(hadron->GetPdgCode()); \r
+                 if (pdg == 321) isKaon = kTRUE;\r
+               }\r
+       }\r
+       \r
+       if(!useMc) { // on DATA\r
+          switch(method) {\r
+         case(1): {\r
+               Bool_t isKTPC=kFALSE;\r
+               Bool_t isPiTPC=kFALSE;\r
+               Bool_t isPTPC=kFALSE;\r
+               Bool_t isKTOF=kFALSE;\r
+               Bool_t isPiTOF=kFALSE;\r
+               Bool_t isPTOF=kFALSE;\r
+               \r
+               Bool_t KaonHyp = kFALSE;\r
+               Bool_t PionHyp = kFALSE;\r
+               Bool_t ProtonHyp = kFALSE;\r
+               \r
+               if(fPidObj->CheckStatus(track,"TOF")) {\r
+                       isKTOF=fPidObj->IsKaonRaw(track,"TOF");\r
+                       isPiTOF=fPidObj->IsPionRaw(track,"TOF");\r
+                       isPTOF=fPidObj->IsProtonRaw(track,"TOF");\r
+               }\r
+               if(fPidObj->CheckStatus(track,"TPC")){\r
+                       isKTPC=fPidObj->IsKaonRaw(track,"TPC");\r
+                       isPiTPC=fPidObj->IsPionRaw(track,"TPC");\r
+                       isPTPC=fPidObj->IsProtonRaw(track,"TPC");               \r
+               }\r
+               \r
+               if (isKTOF && isKTPC) KaonHyp = kTRUE;\r
+               if (isPiTOF && isPiTPC) PionHyp = kTRUE;\r
+               if (isPTOF && isPTPC) ProtonHyp = kTRUE;\r
+               \r
+               if(KaonHyp && !PionHyp && !ProtonHyp) isKaon = kTRUE; \r
+               break;\r
+             }\r
+      case(2): {\r
+               if(fPidObj->MakeRawPid(track,3)>=1) isKaon = kTRUE;\r
+               break;\r
+             }\r
+      }\r
+       }\r
+       \r
+       return isKaon;\r
+       \r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t AliHFAssociatedTrackCuts::IsKZeroSelected(AliAODv0 *vzero, AliAODVertex *vtx1)\r
+{\r
+       \r
+       if(vzero->DcaV0Daughters()>fAODvZeroCuts[0]) return kFALSE;\r
+       if(vzero->Chi2V0()>fAODvZeroCuts[1]) return kFALSE;\r
+       if(vzero->DecayLength(vtx1) < fAODvZeroCuts[2]) return kFALSE;\r
+       if(vzero->DecayLength(vtx1) > fAODvZeroCuts[3]) return kFALSE;\r
+       if(vzero->OpenAngleV0() > fAODvZeroCuts[4]) return kFALSE;\r
+       if(vzero->Pt() < fAODvZeroCuts[5]) return kFALSE;\r
+       if(TMath::Abs(vzero->Eta()) > fAODvZeroCuts[6]) return kFALSE;\r
+\r
+       \r
+       return kTRUE;\r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t *AliHFAssociatedTrackCuts::IsMCpartFromHF(Int_t label, TClonesArray*mcArray){\r
+  // Check origin in MC\r
+\r
+  AliAODMCParticle* mcParticle;\r
+  Int_t pdgCode = -1;\r
+       \r
+  Bool_t isCharmy = kFALSE;\r
+  Bool_t isBeauty = kFALSE;\r
+  Bool_t isD = kFALSE;\r
+  Bool_t isB = kFALSE;\r
+    \r
+     Bool_t *originvect = new Bool_t[4];\r
+    \r
+    originvect[0] = kFALSE;\r
+       originvect[1] = kFALSE;\r
+       originvect[2] = kFALSE;\r
+       originvect[3] = kFALSE;\r
+\r
+       if (label<0) return originvect;\r
+  \r
+       while(pdgCode!=2212){ // loops back to the collision to check the particle source\r
+\r
+    mcParticle = dynamic_cast<AliAODMCParticle*>(mcArray->At(label));\r
+    if(!mcParticle) {AliError("NO MC PARTICLE"); break;}\r
+    pdgCode =  TMath::Abs(mcParticle->GetPdgCode());\r
+\r
+    label = mcParticle->GetMother();\r
+\r
+\r
+    if((pdgCode>=400 && pdgCode <500) || (pdgCode>=4000 && pdgCode<5000 )) isD = kTRUE;\r
+    if((pdgCode>=500 && pdgCode <600) || (pdgCode>=5000 && pdgCode<6000 )) {isD = kFALSE; isB = kTRUE;}\r
+\r
+\r
+    if(pdgCode == 4) isCharmy = kTRUE;\r
+    if(pdgCode == 5) {isBeauty = kTRUE; isCharmy = kFALSE;}\r
+    if(label<0) break;\r
+\r
+  }\r
+\r
+       \r
+       originvect[0] = isCharmy;\r
+       originvect[1] = isBeauty;\r
+       originvect[2] = isD;\r
+       originvect[3] = isB;\r
\r
+    \r
+  return originvect;\r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+Bool_t AliHFAssociatedTrackCuts::InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const {\r
+       //\r
+       // Calculates invmass of track+D0 and rejects if compatible with D*\r
+       // (to remove pions from D*)\r
+       // \r
+       Double_t nsigma = 3.;\r
+       \r
+       Double_t mD0, mD0bar;\r
+       d->InvMassD0(mD0,mD0bar);\r
+       \r
+       Double_t invmassDstar1 = 0, invmassDstar2 = 0; \r
+       Double_t e1Pi = d->EProng(0,211), e2K = d->EProng(1,321); //hyp 1 (pi,K) - D0\r
+       Double_t e1K = d->EProng(0,321), e2Pi = d->EProng(1,211); //hyp 2 (K,pi) - D0bar\r
+       Double_t psum2 = (d->Px()+track->Px())*(d->Px()+track->Px())\r
+       +(d->Py()+track->Py())*(d->Py()+track->Py())\r
+       +(d->Pz()+track->Pz())*(d->Pz()+track->Pz());\r
+       \r
+       switch(hypD0) {\r
+               case 1:\r
+                       invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);\r
+                       if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;\r
+                       break;\r
+               case 2:\r
+                       invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);\r
+                       if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;\r
+                       break;\r
+               case 3:\r
+                       invmassDstar1 = TMath::Sqrt(pow(e1Pi+e2K+track->E(0.1396),2.)-psum2);\r
+                       invmassDstar2 = TMath::Sqrt(pow(e2Pi+e1K+track->E(0.1396),2.)-psum2);\r
+                       if ((TMath::Abs(invmassDstar1-mD0)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;\r
+                       if ((TMath::Abs(invmassDstar2-mD0bar)-0.14543) < nsigma*800.*pow(10.,-6.)) return kFALSE;\r
+                       break;\r
+       }\r
+       \r
+       return kTRUE;\r
+}\r
+//________________________________________________________\r
+void AliHFAssociatedTrackCuts::SetMCEventTypes(Int_t *MCEventTypeArray)\r
+// set the array of event types you want to process in MonteCarlo (gluon splitting, pair production etc.)\r
+{\r
+       if(!fMCEventType) fMCEventType = new Int_t[fNofMCEventType];\r
+       \r
+       for(Int_t k=0; k<fNofMCEventType; k++){\r
+               fMCEventType[k] = MCEventTypeArray[k];\r
+       }\r
+       return; \r
+}\r
+\r
+//________________________________________________________\r
+void AliHFAssociatedTrackCuts::SetAODTrackCuts(Float_t *cutsarray)\r
+{\r
+       if(!fAODTrackCuts) fAODTrackCuts = new Float_t[fNTrackCuts];\r
+       for(Int_t i =0; i<fNTrackCuts; i++){\r
+               fAODTrackCuts[i] = cutsarray[i];\r
+       }\r
+       SetTrackCutsNames();\r
+       return;\r
+}\r
+//________________________________________________________\r
+void AliHFAssociatedTrackCuts::SetTrackCutsNames(/*TString *namearray*/){\r
+       \r
+       fTrackCutsNames = new TString[4];\r
+       fTrackCutsNames[0]= "associated track:: pt min [GeV/c]................: ";\r
+       fTrackCutsNames[1]= "associated track:: pt max [GeV/c]................: ";\r
+       fTrackCutsNames[2]= "associated track:: d0 min [cm]...................: ";\r
+       fTrackCutsNames[3]= "associated track:: d0 max [cm]...................: ";\r
+       \r
+\r
+       \r
+       return;\r
+}\r
+//--------------------------------------------------------------------------\r
+void AliHFAssociatedTrackCuts::SetAODvZeroCuts(Float_t *cutsarray)\r
+{\r
+       \r
+\r
+       if(!fAODvZeroCuts) fAODvZeroCuts = new Float_t[fNvZeroCuts];\r
+       for(Int_t i =0; i<fNvZeroCuts; i++){\r
+               fAODvZeroCuts[i] = cutsarray[i];\r
+       }\r
+       SetvZeroCutsNames();\r
+       return;\r
+}\r
+//--------------------------------------------------------------------------\r
+void AliHFAssociatedTrackCuts::SetvZeroCutsNames(/*TString *namearray*/){\r
+       \r
+       fvZeroCutsNames = new TString[7];\r
+       fvZeroCutsNames[0] = "vZero:: max DCA between two daughters [cm].......: ";\r
+       fvZeroCutsNames[1] = "vZero:: max fit Chi Square between two daughters.: ";\r
+       fvZeroCutsNames[2] = "vZero:: min decay length [cm]....................: ";\r
+       fvZeroCutsNames[3] = "vZero:: max decay length [cm]....................: ";\r
+       fvZeroCutsNames[4] = "vZero:: max opening angle between daughters [rad]: ";\r
+       fvZeroCutsNames[5] = "vZero:: pt min [Gev/c]...........................: ";\r
+       fvZeroCutsNames[6] = "vZero:: |Eta| range <............................: ";\r
+       \r
+       \r
+       return;\r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+void AliHFAssociatedTrackCuts::SetPidAssociated()\r
+{\r
+  //setting PidResponse\r
+  if(fPidObj->GetOldPid()==kFALSE && fPidObj->GetPidResponse()==0x0){\r
+    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
+    AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();\r
+    AliPIDResponse *pidResp=inputHandler->GetPIDResponse();\r
+    fPidObj->SetPidResponse(pidResp);\r
+  }\r
+}\r
+//--------------------------------------------------------------------------\r
+\r
+void AliHFAssociatedTrackCuts::Print(Option_t *option) const\r
+{\r
+  /// overloaded from TObject: print info\r
+  if (strcmp(option, "parameters")==0) {\r
+    PrintPoolParameters();\r
+    return;\r
+  } else if (strcmp(option, "selectedMC")==0) {\r
+    PrintSelectedMCevents();\r
+    return;\r
+  }\r
+  PrintAll();\r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+void AliHFAssociatedTrackCuts::PrintAll() const\r
+{\r
+       \r
+       if(fDescription){\r
+         printf("=================================================");\r
+         printf("\nAdditional description\n");\r
+         std::cout << fDescription << std::endl;\r
+         printf("\n");\r
+       }\r
+       printf("\n=================================================");\r
+       if(fESDTrackCuts){\r
+         printf("\nCuts for the associated track: \n \n");\r
+         \r
+         printf("ITS Refit........................................: %s\n",fESDTrackCuts->GetRequireITSRefit() ? "Yes" : "No");\r
+         printf("TPC Refit........................................: %s\n",fESDTrackCuts->GetRequireTPCRefit() ? "Yes" : "No");\r
+         printf("ITS SA...........................................: %s\n",fESDTrackCuts->GetRequireITSStandAlone() ? "Yes" : "No");\r
+         printf("TPC SA...........................................: %s\n",fESDTrackCuts->GetRequireTPCStandAlone() ? "Yes" : "No");\r
+         printf("Min number of ITS clusters.......................: %d\n",fESDTrackCuts->GetMinNClustersITS());\r
+         printf("Min number of TPC clusters.......................: %d\n",fESDTrackCuts->GetMinNClusterTPC());\r
+         Int_t spd = fESDTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD);\r
+         if(spd==0) std::cout <<  "SPD..............................................: kOff"  << std::endl;\r
+         if(spd==1) std::cout <<  "SPD..............................................: kNone"  << std::endl;\r
+         if(spd==2) std::cout <<  "SPD..............................................: kAny"  << std::endl;\r
+         if(spd==3) std::cout <<  "SPD..............................................: kFirst"  << std::endl;\r
+         if(spd==4) std::cout <<  "SPD..............................................: kOnlyFirst"  << std::endl;\r
+         if(spd==5) std::cout <<  "SPD..............................................: kSecond"  << std::endl;\r
+         if(spd==6) std::cout <<  "SPD..............................................: kOnlySecond"  << std::endl;\r
+         if(spd==7) std::cout <<  "SPD..............................................: kBoth"  << std::endl;\r
+       }\r
+       else printf("\nNo Cuts for Associated Tracks\n");\r
+       std::cout <<  "Filter Bit.......................................: " << fBit  << std::endl;\r
+       std::cout <<  "Charge...........................................: " << fCharge  << std::endl;\r
+       \r
+       if(fAODTrackCuts){\r
+         for(Int_t j=0;j<fNTrackCuts;j++){\r
+           std::cout << fTrackCutsNames[j] << fAODTrackCuts[j] << std::endl;\r
+         }\r
+       }\r
+\r
+       if(fAODvZeroCuts){\r
+         printf("\n");\r
+         printf("=================================================");\r
+         printf("\nCuts for the K0 candidates: \n \n");\r
+         for(Int_t k=0;k<fNvZeroCuts;k++){\r
+           std::cout << fvZeroCutsNames[k] <<  fAODvZeroCuts[k] << std::endl;\r
+         }\r
+       }\r
+       else printf("\nNo Cuts for the K0 candidates\n");\r
+       std::cout << " " << std::endl;\r
+       PrintPoolParameters();\r
+       PrintSelectedMCevents();\r
+\r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+void AliHFAssociatedTrackCuts::PrintPoolParameters() const\r
+{   \r
+       printf("=================================================");\r
+       printf("\nEvent Pool settings: \n \n");\r
+       \r
+       printf("Number of zVtx Bins: %d\n", fNzVtxBins);\r
+       printf("\nzVtx Bins:\n");\r
+       //Double_t zVtxbinLims[fNzVtxBins+1] = fNzVtxBins;\r
+       for(Int_t k=0; k<fNzVtxBins; k++){\r
+               printf("Bin %d..............................................: %.1f - %.1f cm\n", k, fZvtxBins[k], fZvtxBins[k+1]);      \r
+       }\r
+       printf("\n");\r
+       printf("\nNumber of Centrality(multiplicity) Bins: %d\n", fNCentBins);\r
+       printf("\nCentrality(multiplicity) Bins:\n");\r
+       for(Int_t k=0; k<fNCentBins; k++){\r
+               printf("Bin %d..............................................: %.1f - %.1f\n", k, fCentBins[k], fCentBins[k+1]);\r
+       }\r
+\r
+       \r
+       \r
+}\r
+//--------------------------------------------------------------------------\r
+\r
+Double_t AliHFAssociatedTrackCuts::GetTrackWeight(Double_t pt, Double_t eta,Double_t zvtx){\r
+  if(!fEffWeights)return 1.;\r
+  \r
+  Int_t bin=fEffWeights->FindBin(pt,eta,zvtx);\r
+  if(fEffWeights->IsBinUnderflow(bin)||fEffWeights->IsBinOverflow(bin))return 1.;\r
+  return fEffWeights->GetBinContent(bin);\r
+\r
+}\r
+\r
+\r
+//--------------------------------------------------------------------------\r
+Double_t AliHFAssociatedTrackCuts::GetTrigWeight(Double_t pt, Double_t mult){\r
+    \r
+    \r
+    \r
+    if(fTrigEffWeightsvspt){\r
+       Int_t bin=fTrigEffWeightsvspt->FindBin(pt);\r
+        if(fTrigEffWeightsvspt->IsBinUnderflow(bin)||fTrigEffWeightsvspt->IsBinOverflow(bin))return 1.;\r
+        return fTrigEffWeightsvspt->GetBinContent(bin);\r
+        \r
+    }\r
+    \r
+    if(fTrigEffWeights){\r
+        Int_t bin=fTrigEffWeights->FindBin(pt,mult);\r
+        if(fTrigEffWeights->IsBinUnderflow(bin)||fTrigEffWeights->IsBinOverflow(bin))return 1.;\r
+        return fTrigEffWeights->GetBinContent(bin);\r
+        \r
+    }\r
+    \r
+    //if(!fTrigEffWeights && !fTrigEffWeightsvspt)return 1.;\r
+    \r
+    return 1.;\r
+    \r
+}\r
+\r
+//--------------------------------------------------------------------------\r
+Double_t AliHFAssociatedTrackCuts::GetTrigWeightB(Double_t pt, Double_t mult){\r
+    \r
+    if(fTrigEffWeightsvsptB){\r
+        Int_t bin=fTrigEffWeightsvsptB->FindBin(pt);\r
+        if(fTrigEffWeightsvsptB->IsBinUnderflow(bin)||fTrigEffWeightsvsptB->IsBinOverflow(bin))return 1.;\r
+        return fTrigEffWeightsvsptB->GetBinContent(bin);\r
+        \r
+    }\r
+    \r
+    if(fTrigEffWeightsB){\r
+        Int_t bin=fTrigEffWeightsB->FindBin(pt,mult);\r
+        if(fTrigEffWeightsB->IsBinUnderflow(bin)||fTrigEffWeightsB->IsBinOverflow(bin))return 1.;\r
+        return fTrigEffWeightsB->GetBinContent(bin);\r
+        \r
+    }\r
+    \r
+ //   if(!fTrigEffWeightsB && !fTrigEffWeightsvsptB)return 1.;\r
+    return 1;\r
+}\r
+//--------------------------------------------------------------------------\r
+void AliHFAssociatedTrackCuts::PrintSelectedMCevents() const\r
+{\r
+       printf("\n=================================================");\r
+       \r
+       printf("\nSelected MC events: \n \n");\r
+       printf("Number of selected events: %d\n",fNofMCEventType);\r
+       \r
+       for(Int_t k=0; k<fNofMCEventType; k++){\r
+       if(fMCEventType[k]==28) printf("=> Flavour excitation \n");     \r
+       if(fMCEventType[k]==53) printf("=> Pair creation \n");  \r
+       if(fMCEventType[k]==68) printf("=> Gluon splitting \n");        \r
+       }\r
+       \r
+       printf("\n");\r
+       for(Int_t k=0; k<fNofMCEventType; k++){\r
+               printf("MC process code %d \n",fMCEventType[k]);                \r
+       }\r
+       \r
+       printf("\n");\r
+       \r
+       \r
+               \r
+       \r
+}\r
+\r
+\r
index c93db09d58bc8bca77f3b6f28d488c1dde6d5d17..27614e94e9001a63616b0d6fefa4d4b879479247 100644 (file)
-#ifndef AliHFAssociatedTrackCuts_H
-#define AliHFAssociatedTrackCuts_H
-/**************************************************************************
- * Copyright(c) 1998-2010, 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.                  *
- **************************************************************************/
-
-/* $Id$ */
-
-////////////////////////////////////////////////////////////////////////
-//
-// Base class for cuts on Associated tracks for HF Correlation analysis
-//
-// Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch
-////////////////////////////////////////////////////////////////////////
-
-#include <TString.h>
-#include "AliAnalysisCuts.h"
-#include "AliESDtrackCuts.h"
-#include "AliESDVertex.h"
-#include "AliAODPidHF.h"
-#include "AliAODEvent.h"
-#include "AliAODRecoDecayHF2Prong.h"
-#include <TClonesArray.h>
-#include <TH3D.h>
-
-
-class AliAODTrack;
-class AliAODEvent;
-
-
-//
-class AliHFAssociatedTrackCuts : public AliAnalysisCuts
-{
-       public: 
-       AliHFAssociatedTrackCuts();
-       AliHFAssociatedTrackCuts(const char* name, const char* title);
-       
-       
-       AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts& source);
-       AliHFAssociatedTrackCuts& operator=(const AliHFAssociatedTrackCuts& source);
-       
-       virtual ~AliHFAssociatedTrackCuts(); // destructor
-       Bool_t IsSelected(TList*  list) {if(list) return kTRUE; return kFALSE;};
-       Bool_t IsSelected(TObject*  obj) {if(obj) return kTRUE; return kFALSE;};
-       Bool_t IsInAcceptance();
-       Bool_t IsHadronSelected(AliAODTrack * track,const AliESDVertex *primary=0x0,const Double_t magfield=0);
-       Bool_t CheckHadronKinematic(Double_t pt, Double_t d0); 
-       Bool_t Charge(Short_t charge, AliAODTrack* track);
-       Bool_t CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method=1);
-       Bool_t IsKZeroSelected(AliAODv0 *vzero, AliAODVertex *vtx1);
-       Bool_t *IsMCpartFromHF(Int_t label, TClonesArray*mcArray);
-       Bool_t InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const;
-       void SetPidAssociated();        
-       
-       // getters
-       AliESDtrackCuts * GetESDTrackCuts() const {return fESDTrackCuts;}
-       AliAODPidHF * GetPIDObject() const {return fPidObj;}
-       TH3D * GetEfficiencyWeight() const {return fEffWeights;}
-       
-       Int_t GetMaxNEventsInPool() const {return fPoolMaxNEvents;}
-       Int_t GetMinNTracksInPool() const {return fPoolMinNTracks;}
-       Int_t GetMinEventsToMix() const {return fMinEventsToMix;}
-       Int_t GetNZvtxPoolBins() const {return fNzVtxBins;}
-       Double_t *GetZvtxPoolBins() const {return fZvtxBins;}
-       Int_t GetNCentPoolBins() const {return fNCentBins;}
-       Double_t *GetCentPoolBins() const {return fCentBins;}
-       
-       Int_t GetNofMCEventType() const {return fNofMCEventType;}
-       Int_t *GetMCEventType() const {return fMCEventType;}
-       
-       Int_t GetNTrackCuts() const {return fNTrackCuts;}
-       Float_t* GetAODTrackCuts() const {return fAODTrackCuts;}
-       TString * GetTrackCutNames() const {return fTrackCutsNames;}
-       Int_t GetNvZeroCuts() const {return fNvZeroCuts;}
-       Float_t * GetAODvZeroCuts() const {return fAODvZeroCuts;}
-       TString * GetvZeroCutNames() const {return fvZeroCutsNames;}
-       Int_t GetFilterBit() const {return fBit;}
-       Short_t GetCharge() const {return fCharge;}
-       TString GetDescription() const {return fDescription;}
-       
-       
-       
-       void AddTrackCuts(const AliESDtrackCuts *cuts) {
-               delete fESDTrackCuts; 
-               fESDTrackCuts=new AliESDtrackCuts(*cuts); 
-               return;
-       }
-       
-       void AddDescription(TString description){fDescription=description;}
-       
-       //setters
-       //event pool settings
-       void SetMaxNEventsInPool(Int_t events){fPoolMaxNEvents=events;}
-       void SetMinNTracksInPool(Int_t tracks){fPoolMinNTracks=tracks;}
-       void SetMinEventsToMix(Int_t events){fMinEventsToMix=events;}
-       
-       void SetNofPoolBins(Int_t Nzvtxbins, Int_t Ncentbins){
-               fNzVtxBins=Nzvtxbins;
-               fNzVtxBinsDim=Nzvtxbins+1;
-               
-           fNCentBins=Ncentbins;
-               fNCentBinsDim=Ncentbins+1;
-       }
-       
-       void SetPoolBins(Double_t *ZvtxBins, Double_t* CentBins){
-               fZvtxBins=ZvtxBins; 
-               fCentBins=CentBins;
-       }
-       
-       // set MC events to process
-       
-       void SetNofMCEventTypes(Int_t k) {fNofMCEventType=k;}
-       void SetMCEventTypes(Int_t *MCEventTypeArray);
-       
-       //cut settings
-       void SetAODTrackCuts(Float_t *cutsarray);
-       void SetTrackCutsNames(/*TString *namearray*/);
-       void SetAODvZeroCuts(Float_t *cutsarray);
-       void SetvZeroCutsNames(/*TString *namearray*/);
-       void SetPidHF(AliAODPidHF* pid) {fPidObj = pid; return;}
-       void SetCharge(Short_t charge) {fCharge = charge;}
-       void SetFilterBit(Int_t bit) {fBit = bit;}
-       void SetEfficiencyWeightMap(TH3D *hMap){if(fEffWeights)delete fEffWeights;fEffWeights=(TH3D*)hMap->Clone();}
-       Double_t GetTrackWeight(Double_t pt, Double_t eta,Double_t zvtx);
-       void Print(Option_t *option) const;
-       virtual void PrintAll() const;
-       virtual void PrintPoolParameters() const;
-       virtual void PrintSelectedMCevents() const;
-
-       
-       
-       
-       void SetNVarsTrack(Int_t nVars){fNTrackCuts=nVars;}
-       void SetNVarsVzero(Int_t nVars){fNvZeroCuts=nVars;}
-       
-       
-       
-private:
-       AliESDtrackCuts *fESDTrackCuts; // track cut object
-       AliAODPidHF * fPidObj;     /// PID object
-       TH3D *fEffWeights;     // weight map (pt,eta,zvtx) to account for single track efficiency  
-       Int_t fPoolMaxNEvents; // set maximum number of events in the pool
-       Int_t fPoolMinNTracks; // se minimum number of tracks in the pool
-       Int_t fMinEventsToMix; // set the minimum number of events you wanna mix
-       
-       Int_t fNzVtxBins; // number of z vrtx bins
-       Int_t fNzVtxBinsDim; // number of z vrtx bins +1 : necessary to initialize correctly the array
-       Double_t* fZvtxBins; // [fNzVtxBinsDim]
-       
-       
-       Int_t fNCentBins; //number of centrality bins
-       Int_t fNCentBinsDim; //number of centrality bins bins +1 : necessary to initialize correctly the array
-       Double_t* fCentBins; // [fNCentBinsDim]
-       
-       Int_t fNofMCEventType;// number of event types to be selected in MC simultaneously;
-       Int_t *fMCEventType;//[fNofMCEventType]
-       
-       Int_t fNTrackCuts;     // array dimension
-       Float_t* fAODTrackCuts;//[fNTrackCuts]
-       TString * fTrackCutsNames;//[fNTrackCuts]
-       Int_t fNvZeroCuts;// array dimension
-       Float_t *fAODvZeroCuts;//[fNvZeroCuts]
-       TString * fvZeroCutsNames;//[fNvZeroCuts]
-       Int_t fBit; // filterBit
-       Short_t fCharge; // charge (+1 or -1)
-       TString fDescription; // additional description to the cuts
-       
-       
-       ClassDef(AliHFAssociatedTrackCuts,5);
-};
-
-
-#endif
+#ifndef AliHFAssociatedTrackCuts_H\r
+#define AliHFAssociatedTrackCuts_H\r
+/**************************************************************************\r
+ * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id$ */\r
+\r
+////////////////////////////////////////////////////////////////////////\r
+//\r
+// Base class for cuts on Associated tracks for HF Correlation analysis\r
+//\r
+// Author: S.Bjelogrlic (Utrecht) sandro.bjelogrlic@cern.ch\r
+////////////////////////////////////////////////////////////////////////\r
+\r
+#include <TString.h>\r
+#include "AliAnalysisCuts.h"\r
+#include "AliESDtrackCuts.h"\r
+#include "AliESDVertex.h"\r
+#include "AliAODPidHF.h"\r
+#include "AliAODEvent.h"\r
+#include "AliAODRecoDecayHF2Prong.h"\r
+#include <TClonesArray.h>\r
+#include <TH3D.h>\r
+\r
+\r
+class AliAODTrack;\r
+class AliAODEvent;\r
+\r
+\r
+//\r
+class AliHFAssociatedTrackCuts : public AliAnalysisCuts\r
+{\r
+       public: \r
+       AliHFAssociatedTrackCuts();\r
+       AliHFAssociatedTrackCuts(const char* name, const char* title);\r
+       \r
+       \r
+       AliHFAssociatedTrackCuts(const AliHFAssociatedTrackCuts& source);\r
+       AliHFAssociatedTrackCuts& operator=(const AliHFAssociatedTrackCuts& source);\r
+       \r
+       virtual ~AliHFAssociatedTrackCuts(); // destructor\r
+       Bool_t IsSelected(TList*  list) {if(list) return kTRUE; return kFALSE;};\r
+       Bool_t IsSelected(TObject*  obj) {if(obj) return kTRUE; return kFALSE;};\r
+       Bool_t IsInAcceptance();\r
+       Bool_t IsHadronSelected(AliAODTrack * track,const AliESDVertex *primary=0x0,const Double_t magfield=0);\r
+       Bool_t CheckHadronKinematic(Double_t pt, Double_t d0); \r
+       Bool_t Charge(Short_t charge, AliAODTrack* track);\r
+       Bool_t CheckKaonCompatibility(AliAODTrack * track, Bool_t useMc, TClonesArray* mcArray, Int_t method=1);\r
+       Bool_t IsKZeroSelected(AliAODv0 *vzero, AliAODVertex *vtx1);\r
+       Bool_t *IsMCpartFromHF(Int_t label, TClonesArray*mcArray);\r
+       Bool_t InvMassDstarRejection(AliAODRecoDecayHF2Prong* d, AliAODTrack *track, Int_t hypD0) const;\r
+       void SetPidAssociated();        \r
+       \r
+       // getters\r
+       AliESDtrackCuts * GetESDTrackCuts() const {return fESDTrackCuts;}\r
+       AliAODPidHF * GetPIDObject() const {return fPidObj;}\r
+       TH3D * GetEfficiencyWeight() const {return fEffWeights;}\r
+       \r
+       Int_t GetMaxNEventsInPool() const {return fPoolMaxNEvents;}\r
+       Int_t GetMinNTracksInPool() const {return fPoolMinNTracks;}\r
+       Int_t GetMinEventsToMix() const {return fMinEventsToMix;}\r
+       Int_t GetNZvtxPoolBins() const {return fNzVtxBins;}\r
+       Double_t *GetZvtxPoolBins() const {return fZvtxBins;}\r
+       Int_t GetNCentPoolBins() const {return fNCentBins;}\r
+       Double_t *GetCentPoolBins() const {return fCentBins;}\r
+       \r
+       Int_t GetNofMCEventType() const {return fNofMCEventType;}\r
+       Int_t *GetMCEventType() const {return fMCEventType;}\r
+       \r
+       Int_t GetNTrackCuts() const {return fNTrackCuts;}\r
+       Float_t* GetAODTrackCuts() const {return fAODTrackCuts;}\r
+       TString * GetTrackCutNames() const {return fTrackCutsNames;}\r
+       Int_t GetNvZeroCuts() const {return fNvZeroCuts;}\r
+       Float_t * GetAODvZeroCuts() const {return fAODvZeroCuts;}\r
+       TString * GetvZeroCutNames() const {return fvZeroCutsNames;}\r
+       Int_t GetFilterBit() const {return fBit;}\r
+       Short_t GetCharge() const {return fCharge;}\r
+       TString GetDescription() const {return fDescription;}\r
+    \r
+\r
+       \r
+       \r
+       \r
+       void AddTrackCuts(const AliESDtrackCuts *cuts) {\r
+               delete fESDTrackCuts; \r
+               fESDTrackCuts=new AliESDtrackCuts(*cuts); \r
+               return;\r
+       }\r
+       \r
+       void AddDescription(TString description){fDescription=description;}\r
+       \r
+       //setters\r
+       //event pool settings\r
+       void SetMaxNEventsInPool(Int_t events){fPoolMaxNEvents=events;}\r
+       void SetMinNTracksInPool(Int_t tracks){fPoolMinNTracks=tracks;}\r
+       void SetMinEventsToMix(Int_t events){fMinEventsToMix=events;}\r
+       \r
+       void SetNofPoolBins(Int_t Nzvtxbins, Int_t Ncentbins){\r
+               fNzVtxBins=Nzvtxbins;\r
+               fNzVtxBinsDim=Nzvtxbins+1;\r
+               \r
+           fNCentBins=Ncentbins;\r
+               fNCentBinsDim=Ncentbins+1;\r
+       }\r
+       \r
+       void SetPoolBins(Double_t *ZvtxBins, Double_t* CentBins){\r
+               fZvtxBins=ZvtxBins; \r
+               fCentBins=CentBins;\r
+       }\r
+       \r
+       // set MC events to process\r
+       \r
+       void SetNofMCEventTypes(Int_t k) {fNofMCEventType=k;}\r
+       void SetMCEventTypes(Int_t *MCEventTypeArray);\r
+       \r
+       //cut settings\r
+       void SetAODTrackCuts(Float_t *cutsarray);\r
+       void SetTrackCutsNames(/*TString *namearray*/);\r
+       void SetAODvZeroCuts(Float_t *cutsarray);\r
+       void SetvZeroCutsNames(/*TString *namearray*/);\r
+       void SetPidHF(AliAODPidHF* pid) {fPidObj = pid; return;}\r
+       void SetCharge(Short_t charge) {fCharge = charge;}\r
+       void SetFilterBit(Int_t bit) {fBit = bit;}\r
+       void SetEfficiencyWeightMap(TH3D *hMap){if(fEffWeights)delete fEffWeights;fEffWeights=(TH3D*)hMap->Clone();}\r
+    \r
+    void SetTriggerEffWeightMapvspt(TH1D* hTrigMap) {if(fTrigEffWeightsvspt) delete fTrigEffWeightsvspt; fTrigEffWeightsvspt=(TH1D*)hTrigMap->Clone();}\r
+       void SetTriggerEffWeightMapvsptB(TH1D* hTrigMapB) {if(fTrigEffWeightsvsptB) delete fTrigEffWeightsvsptB; fTrigEffWeightsvsptB=(TH1D*)hTrigMapB->Clone();}\r
+    \r
+    void SetTriggerEffWeightMap(TH2D* hTrigMap) {if(fTrigEffWeights) delete fTrigEffWeights; fTrigEffWeights=(TH2D*)hTrigMap->Clone();}\r
+       void SetTriggerEffWeightMapB(TH2D* hTrigMapB) {if(fTrigEffWeightsB) delete fTrigEffWeightsB; fTrigEffWeightsB=(TH2D*)hTrigMapB->Clone();}\r
+    \r
+    \r
+       Double_t GetTrackWeight(Double_t pt, Double_t eta,Double_t zvtx);\r
+    Double_t GetTrigWeight(Double_t pt, Double_t mult=0);\r
+       Double_t GetTrigWeightB(Double_t pt, Double_t mult=0);\r
+    \r
+    Bool_t IsTrackEffMap(){if(fEffWeights) return kTRUE; else return kFALSE;}\r
+    Bool_t IsTrigEffMap1D(){ if(fTrigEffWeightsvspt) return kTRUE; else return kFALSE;}\r
+    Bool_t IsTrigEffMap1DB(){ if(fTrigEffWeightsvsptB) return kTRUE; else return kFALSE;}\r
+    Bool_t IsTrigEffMap2D(){ if(fTrigEffWeights) return kTRUE; else return kFALSE;}\r
+    Bool_t IsTrigEffMap2DB(){ if(fTrigEffWeightsB) return kTRUE; else return kFALSE;}\r
+    \r
+    \r
+       void Print(Option_t *option) const;\r
+       virtual void PrintAll() const;\r
+       virtual void PrintPoolParameters() const;\r
+       virtual void PrintSelectedMCevents() const;\r
+\r
+       \r
+       \r
+       \r
+       void SetNVarsTrack(Int_t nVars){fNTrackCuts=nVars;}\r
+       void SetNVarsVzero(Int_t nVars){fNvZeroCuts=nVars;}\r
+       \r
+       \r
+       \r
+private:\r
+       AliESDtrackCuts *fESDTrackCuts; // track cut object\r
+       AliAODPidHF * fPidObj;     /// PID object\r
+       TH3D *fEffWeights;     // weight map (pt,eta,zvtx) to account for single track efficiency\r
+    TH1D *fTrigEffWeightsvspt;     // weight map (pt,mult) to account for trigger efficiency (on data, from c)\r
+       TH1D *fTrigEffWeightsvsptB;     // weight map (pt,mult) to account for trigger efficiency (from b)\r
+    TH2D *fTrigEffWeights;     // weight map (pt,mult) to account for trigger efficiency (on data, from c)\r
+       TH2D *fTrigEffWeightsB;     // weight map (pt,mult) to account for trigger efficiency (from b)\r
+       Int_t fPoolMaxNEvents; // set maximum number of events in the pool\r
+       Int_t fPoolMinNTracks; // se minimum number of tracks in the pool\r
+       Int_t fMinEventsToMix; // set the minimum number of events you wanna mix\r
+       \r
+       Int_t fNzVtxBins; // number of z vrtx bins\r
+       Int_t fNzVtxBinsDim; // number of z vrtx bins +1 : necessary to initialize correctly the array\r
+       Double_t* fZvtxBins; // [fNzVtxBinsDim]\r
+       \r
+       \r
+       Int_t fNCentBins; //number of centrality bins\r
+       Int_t fNCentBinsDim; //number of centrality bins bins +1 : necessary to initialize correctly the array\r
+       Double_t* fCentBins; // [fNCentBinsDim]\r
+       \r
+       Int_t fNofMCEventType;// number of event types to be selected in MC simultaneously;\r
+       Int_t *fMCEventType;//[fNofMCEventType]\r
+       \r
+       Int_t fNTrackCuts;     // array dimension\r
+       Float_t* fAODTrackCuts;//[fNTrackCuts]\r
+       TString * fTrackCutsNames;//[fNTrackCuts]\r
+       Int_t fNvZeroCuts;// array dimension\r
+       Float_t *fAODvZeroCuts;//[fNvZeroCuts]\r
+       TString * fvZeroCutsNames;//[fNvZeroCuts]\r
+       Int_t fBit; // filterBit\r
+       Short_t fCharge; // charge (+1 or -1)\r
+       TString fDescription; // additional description to the cuts\r
+       \r
+       \r
+       ClassDef(AliHFAssociatedTrackCuts,6);\r
+};\r
+\r
+\r
+#endif\r
index 607d96d21e447188189906691a3cb12621ca2737..4bbb7b76f5e7cc6095c8dee03583ee5903faf459 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-2009, 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.                  *
- **************************************************************************/
-//
-//
-//             Base class for Heavy Flavour Correlations Analysis
-//             Single Event and Mixed Event Analysis are implemented
-//
-//-----------------------------------------------------------------------
-//          
-//
-//                                                Author S.Bjelogrlic
-//                         Utrecht University 
-//                      sandro.bjelogrlic@cern.ch
-//
-//-----------------------------------------------------------------------
-
-/* $Id$ */
-
-#include <TParticle.h>
-#include <TVector3.h>
-#include <TChain.h>
-#include "TROOT.h"
-#include "AliHFCorrelator.h"
-#include "AliRDHFCutsDStartoKpipi.h"
-#include "AliHFAssociatedTrackCuts.h"
-#include "AliEventPoolManager.h"
-#include "AliReducedParticle.h"
-#include "AliCentrality.h"
-#include "AliAODMCParticle.h"
-
-using std::cout;
-using std::endl;
-
-//_____________________________________________________
-AliHFCorrelator::AliHFCorrelator() :
-//
-// default constructor
-//
-TNamed(),
-fPoolMgr(0x0),         
-fPool(0x0),
-fhadcuts(0x0),
-fAODEvent(0x0),
-fAssociatedTracks(0x0),
-fmcArray(0x0),
-fReducedPart(0x0),
-fD0cand(0x0), 
-fhypD0(0), 
-fDCharge(0),
-
-fmixing(kFALSE),
-fmontecarlo(kFALSE),
-fsystem(kFALSE),
-fUseReco(kTRUE),
-fselect(kUndefined),
-
-fUseImpactParameter(0),
-fPIDmode(0),
-
-fNofTracks(0),
-fPoolContent(0),
-
-fPhiMin(0),
-fPhiMax(0),
-
-fPtTrigger(0),
-fPhiTrigger(0),
-fEtaTrigger(0),
-
-
-fDeltaPhi(0),
-fDeltaEta(0),
-fk0InvMass(0)
-
-{
-       // default constructor  
-}
-
-
-
-//_____________________________________________________
-AliHFCorrelator::AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb) :
-TNamed(name,"title"),
-fPoolMgr(0x0),         
-fPool(0x0),
-fhadcuts(0x0),
-fAODEvent(0x0),
-fAssociatedTracks(0x0),
-fmcArray(0x0),
-fReducedPart(0x0),
-fD0cand(0x0), 
-fhypD0(0),
-fDCharge(0),
-
-fmixing(kFALSE),
-fmontecarlo(kFALSE),
-fsystem(ppOrPbPb),
-fUseReco(kTRUE),
-fselect(kUndefined),
-fUseImpactParameter(0),
-fPIDmode(0),
-
-fNofTracks(0),
-fPoolContent(0),
-
-fPhiMin(0),
-fPhiMax(0),
-
-fPtTrigger(0),
-fPhiTrigger(0),
-fEtaTrigger(0),
-
-
-fDeltaPhi(0),
-fDeltaEta(0),
-fk0InvMass(0)
-{
-       fhadcuts = cuts; 
-}
-
-//_____________________________________________________
-AliHFCorrelator::~AliHFCorrelator() 
-{
-//
-// destructor
-//     
-       
-       if(fPoolMgr)  {delete fPoolMgr; fPoolMgr=0;}       
-       if(fPool) {delete fPool; fPool=0;}
-       if(fhadcuts) {delete fhadcuts; fhadcuts=0;}
-       if(fAODEvent) {delete fAODEvent; fAODEvent=0;}
-       if(fAssociatedTracks) {delete fAssociatedTracks; fAssociatedTracks=0;}
-       if(fmcArray) {delete fmcArray; fmcArray=0;}
-       if(fReducedPart) {delete fReducedPart; fReducedPart=0;}
-       if(fD0cand) {delete fD0cand; fD0cand=0;}
-       
-       
-       if(fNofTracks) fNofTracks = 0;
-       
-       if(fPhiMin) fPhiMin = 0;
-       if(fPhiMax) fPhiMax = 0;
-       
-       if(fPtTrigger) fPtTrigger=0;
-       if(fPhiTrigger) fPhiTrigger=0;
-       if(fEtaTrigger) fEtaTrigger=0;
-       
-       if(fDeltaPhi) fDeltaPhi=0;
-       if(fDeltaEta) fDeltaEta=0;
-       
-       if(fk0InvMass) fk0InvMass=0;
-       
-}
-//_____________________________________________________
-Bool_t AliHFCorrelator::DefineEventPool(){
-       // definition of the Pool Manager for Event Mixing
-       
-
-       Int_t MaxNofEvents = fhadcuts->GetMaxNEventsInPool();
-       Int_t MinNofTracks = fhadcuts->GetMinNTracksInPool();
-       Int_t NofCentBins = fhadcuts->GetNCentPoolBins();
-       Double_t * CentBins = fhadcuts->GetCentPoolBins();
-       Int_t NofZVrtxBins = fhadcuts->GetNZvtxPoolBins();
-       Double_t *ZVrtxBins = fhadcuts->GetZvtxPoolBins();
-               
-                       
-       fPoolMgr = new AliEventPoolManager(MaxNofEvents, MinNofTracks, NofCentBins, CentBins, NofZVrtxBins, ZVrtxBins);
-       if(!fPoolMgr) return kFALSE;
-       return kTRUE;
-}
-//_____________________________________________________
-Bool_t AliHFCorrelator::Initialize(){
-       
-    //  std::cout << "AliHFCorrelator::Initialize"<< std::endl;
-  AliInfo("AliHFCorrelator::Initialize") ;
-  if(!fAODEvent){
-    AliInfo("No AOD event") ;
-    return kFALSE;
-  }
-    //std::cout << "No AOD event" << std::endl;
-       
-       AliCentrality *centralityObj = 0;
-       Int_t multiplicity = -1;
-       Double_t MultipOrCent = -1;
-       
-       // initialize the pool for event mixing
-       if(!fsystem){ // pp
-       multiplicity = fAODEvent->GetNTracks();
-               MultipOrCent = multiplicity; // convert from Int_t to Double_t
-       }
-       if(fsystem){ // PbPb
-               
-               centralityObj = fAODEvent->GetHeader()->GetCentralityP();
-               MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");
-               AliInfo(Form("Centrality is %f", MultipOrCent));
-       }
-       
-       AliAODVertex *vtx = fAODEvent->GetPrimaryVertex();
-       Double_t zvertex = vtx->GetZ(); // zvertex
-       Double_t * CentBins = fhadcuts->GetCentPoolBins();
-       Double_t poolmin=CentBins[0];
-       Double_t poolmax=CentBins[fhadcuts->GetNCentPoolBins()];
-
-       
-               if(TMath::Abs(zvertex)>=10 || MultipOrCent>poolmax || MultipOrCent < poolmin) {
-               if(!fsystem)AliInfo(Form("pp Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,MultipOrCent));
-               if(fsystem) AliInfo(Form("PbPb Event with Zvertex = %.2f cm and centrality = %.1f  out of pool bounds, SKIPPING",zvertex,MultipOrCent));
-
-                       return kFALSE;
-               }
-       
-       fPool = fPoolMgr->GetEventPool(MultipOrCent, zvertex);
-       
-       if (!fPool){
-               AliInfo(Form("No pool found for multiplicity = %f, zVtx = %f cm", MultipOrCent, zvertex));
-           return kFALSE;
-       }
-       //fPool->PrintInfo();
-       return kTRUE;
-}
-
-//_____________________________________________________
-Bool_t AliHFCorrelator::ProcessEventPool(){
-        // analysis on Mixed Events
-       //cout << "AliHFCorrelator::ProcessEventPool"<< endl;
-               if(!fmixing) return kFALSE;
-               if(!fPool->IsReady()) return kFALSE;
-               if(fPool->GetCurrentNEvents()<fhadcuts->GetMinEventsToMix()) return kFALSE;
-       //      fPool->PrintInfo();
-               fPoolContent = fPool->GetCurrentNEvents();
-               
-               return kTRUE;
-       
-}
-
-//_____________________________________________________
-Bool_t AliHFCorrelator::ProcessAssociatedTracks(Int_t EventLoopIndex, const TObjArray* associatedTracks){
-  // TODO: memory leak needs to be fixed, for every call, a new array
-  // is allocated, but the pointer immediately lost. The cleanup is
-  // not straightforward as in the case of event mixing the pointer
-  // will be an external array which must not be deleted.
-       fAssociatedTracks = new TObjArray();
-
-       if(!fmixing){ // analysis on Single Event
-               
-               
-               
-               if(fselect==kHadron || fselect ==kKaon) fAssociatedTracks = AcceptAndReduceTracks(fAODEvent);
-               if(fselect==kKZero) {fAssociatedTracks = AcceptAndReduceKZero(fAODEvent);}      
-               if(fselect==kElectron && associatedTracks) fAssociatedTracks=new TObjArray(*associatedTracks);
-               
-       }
-       
-       if(fmixing) { // analysis on Mixed Events
-               
-                       
-               fAssociatedTracks = fPool->GetEvent(EventLoopIndex);
-                               
-                               
-                       
-
-       } // end if mixing
-       
-       if(!fAssociatedTracks) return kFALSE;
-       
-       fNofTracks = fAssociatedTracks->GetEntriesFast(); 
-               
-       return kTRUE;
-       
-}
-//_____________________________________________________
-Bool_t AliHFCorrelator::Correlate(Int_t loopindex){
-
-       if(loopindex >= fNofTracks) return kFALSE;
-       if(!fAssociatedTracks) return kFALSE;
-       
-       fReducedPart = (AliReducedParticle*)fAssociatedTracks->At(loopindex);
-       
-
-       fDeltaPhi = SetCorrectPhiRange(fPhiTrigger - fReducedPart->Phi());
-       
-       fDeltaEta = fEtaTrigger - fReducedPart->Eta();
-
-       return kTRUE;
-       
-}
-               
-//_____________________________________________________
-Bool_t AliHFCorrelator::PoolUpdate(const TObjArray* associatedTracks){
-
-       if(!fmixing) return kFALSE;
-       if(!fPool) return kFALSE;
-       if(fmixing) { // update the pool for Event Mixing
-               TObjArray* objArr = NULL;
-               if(fselect==kHadron || fselect==kKaon) objArr = (TObjArray*)AcceptAndReduceTracks(fAODEvent);
-               else if(fselect==kKZero) objArr = (TObjArray*)AcceptAndReduceKZero(fAODEvent);
-               else if(fselect==kElectron && associatedTracks){
-                 objArr = new TObjArray(*associatedTracks);
-               }
-               else return kFALSE;
-               if(objArr->GetEntriesFast()>0) fPool->UpdatePool(objArr); // updating the pool only if there are entries in the array
-       }
-               
-       return kTRUE;
-       
-}
-               
-//_____________________________________________________
-Double_t AliHFCorrelator::SetCorrectPhiRange(Double_t phi){
-       Double_t pi = TMath::Pi();
-       
-       if(phi<fPhiMin) phi = phi + 2*pi;
-       if(phi>fPhiMax) phi = phi - 2*pi;
-       
-       return phi;
-}
-
-//_____________________________________________________
-TObjArray*  AliHFCorrelator::AcceptAndReduceTracks(AliAODEvent* inputEvent){
-
-  Double_t weight=1.;
-  Int_t nTracks = inputEvent->GetNTracks();
-  AliAODVertex * vtx = inputEvent->GetPrimaryVertex();
-  Double_t pos[3],cov[6];
-  vtx->GetXYZ(pos);
-  vtx->GetCovarianceMatrix(cov);
-  const AliESDVertex vESD(pos,cov,100.,100);
-  
-  Double_t Bz = inputEvent->GetMagneticField();
-       
-  
-  TObjArray* tracksClone = new TObjArray;
-  tracksClone->SetOwner(kTRUE);
-  
-  //*******************************************************
-  // use reconstruction
-  if(fUseReco){
-    for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {
-      AliAODTrack* track = inputEvent->GetTrack(iTrack);
-      if (!track) continue;
-      if(!fhadcuts->IsHadronSelected(track,&vESD,Bz)) continue; // apply ESD level selections
-      if(!fhadcuts->Charge(fDCharge,track)) continue; // apply selection on charge, if required
-
-      Double_t pT = track->Pt();
-      
-      //compute impact parameter
-      Double_t d0z0[2],covd0z0[3];
-      Double_t d0=-999999.;
-      if(fUseImpactParameter) track->PropagateToDCA(vtx,Bz,100,d0z0,covd0z0);
-      else d0z0[0] = 1. ; // random number - be careful with the cuts you applied
-      
-      if(fUseImpactParameter==1) d0 = TMath::Abs(d0z0[0]); // use impact parameter
-      if(fUseImpactParameter==2) { // use impact parameter over resolution
-       if(TMath::Abs(covd0z0[0])>0.00000001) d0 = TMath::Abs(d0z0[0])/TMath::Sqrt(covd0z0[0]); 
-       else d0 = -1.; // if the resoultion is Zero, rejects the track - to be on the safe side
-       
-      }
-      
-      if(fmontecarlo) {// THIS TO BE CHECKED
-        Int_t hadLabel = track->GetLabel();
-       if(hadLabel < 0) continue;      
-      }
-      
-      if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts
-      Bool_t rejectsoftpi = kTRUE;// TO BE CHECKED: DO WE WANT IT TO kTRUE AS A DEFAULT?
-      if(fD0cand && !fmixing) rejectsoftpi = fhadcuts->InvMassDstarRejection(fD0cand,track,fhypD0); // TO BE CHECKED: WHY NOT FOR EM?
-      
-      
-      if(fselect ==kKaon){     
-       if(!fhadcuts->CheckKaonCompatibility(track,fmontecarlo,fmcArray,fPIDmode)) continue; // check if it is a Kaon - data and MC
-      }
-      weight=fhadcuts->GetTrackWeight(pT,track->Eta(),pos[2]);
-      tracksClone->Add(new AliReducedParticle(track->Eta(), track->Phi(), pT,track->GetLabel(),track->GetID(),d0,rejectsoftpi,track->Charge(),weight));
-    } // end loop on tracks
-  } // end if use reconstruction kTRUE
-  
-  //*******************************************************
-  
-  //use MC truth
-  if(fmontecarlo && !fUseReco){
-    
-    for (Int_t iPart=0; iPart<fmcArray->GetEntriesFast(); iPart++) { 
-      AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iPart));
-      if (!mcPart) {
-       AliWarning("MC Particle not found in tree, skipping"); 
-       continue;
-      }
-      if(!mcPart->Charge()) continue; // consider only charged tracks
-      
-      Int_t PDG =TMath::Abs(mcPart->PdgCode()); 
-if(fselect ==kHadron) {if(!((PDG==321)||(PDG==211)||(PDG==2212)||(PDG==13)||(PDG==11))) continue;} // select only if kaon, pion, proton, muon or electron
-      else if(fselect ==kKaon) {if(!(PDG==321)) continue;} // select only if kaon
-      else if(fselect ==kElectron) {if(!(PDG==11)) continue;} // select only if electron
-
-      Double_t pT = mcPart->Pt();
-      Double_t d0 =1; // set 1 fot the moment - no displacement calculation implemented yet
-      if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts
-      
-      tracksClone->Add(new AliReducedParticle(mcPart->Eta(), mcPart->Phi(), pT,iPart,-1,d0,kFALSE,mcPart->Charge()));
-    }
-    
-  } // end if use  MC truth
-  
-  
-  return tracksClone;
-}
-
-//_____________________________________________________
-TObjArray*  AliHFCorrelator::AcceptAndReduceKZero(AliAODEvent* inputEvent){
-       
-       Int_t nOfVZeros = inputEvent->GetNumberOfV0s();
-       TObjArray* KZeroClone = new TObjArray;
-       AliAODVertex *vertex1 = (AliAODVertex*)inputEvent->GetPrimaryVertex();
-
- // use reconstruction         
-     if(fUseReco){
-       Int_t v0label = -1;
-       Int_t pdgDgK0toPipi[2] = {211,211};
-       Double_t mPDGK0=0.497614;//TDatabasePDG::Instance()->GetParticle(310)->Mass();
-       const Int_t size = inputEvent->GetNumberOfV0s();
-       Int_t prevnegID[size];
-       Int_t prevposID[size];
-       for(Int_t iv0 =0; iv0< nOfVZeros; iv0++){// loop on all v0 candidates
-               AliAODv0 *v0 = (static_cast<AliAODEvent*> (inputEvent))->GetV0(iv0);
-               if(!v0) {
-                 AliInfo(Form("is not a v0 at step, %d", iv0)) ;
-                 //cout << "is not a v0 at step " << iv0 << endl;
-                 continue;
-               }
-               
-               if(!fhadcuts->IsKZeroSelected(v0,vertex1)) continue; // check if it is a k0
-           
-               // checks to avoid double counting
-               Int_t negID = -999;
-               Int_t posID = -998;
-               //Int_t a = 0;
-               prevnegID[iv0] = -997;
-               prevposID[iv0]= -996;
-               negID = v0->GetNegID();
-               posID = v0->GetPosID();
-               Bool_t DoubleCounting = kFALSE;
-               
-               for(Int_t k=0; k<iv0; k++){
-                       if(((negID==prevnegID[k])&&(posID==prevposID[k]))||((negID==prevposID[k])&&(posID==prevnegID[k]))){
-                               DoubleCounting = kTRUE;
-                               //a=k;
-                               break;
-                       }//end if
-               } // end for
-               
-               if(DoubleCounting) continue;
-               else{ 
-                       prevposID[iv0] = posID; 
-                       prevnegID[iv0] = negID;
-               }
-               
-               if(fmontecarlo) v0label = v0->MatchToMC(310,fmcArray, 0, pdgDgK0toPipi); //match a K0 short
-               Double_t k0pt = v0->Pt();
-               Double_t k0eta = v0->Eta();
-               Double_t k0Phi = v0->Phi();
-           fk0InvMass = v0->MassK0Short();     
-               
-               //if (loopindex == 0) {
-               //      if(!plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectra"))->Fill(k0InvMass,k0pt); // spectra for all k0
-               //      if(plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectraifHF"))->Fill(k0InvMass,k0pt); // spectra for k0 in association with a D*
-               //}
-               // if there are more D* candidates per event, loopindex == 0 makes sure you fill the mass spectra only once!
-               
-               if(TMath::Abs(fk0InvMass-mPDGK0)>3*0.004) continue; // select candidates within 3 sigma
-               KZeroClone->Add(new AliReducedParticle(k0eta,k0Phi,k0pt,v0label));
-               
-       }
-     } // end if use reconstruction kTRUE
-       
-
-
-//*********************************************************************//
-     //use MC truth
-     if(fmontecarlo && !fUseReco){
-               
-               for (Int_t iPart=0; iPart<fmcArray->GetEntriesFast(); iPart++) { 
-                       AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iPart));
-                       if (!mcPart) {
-                               AliWarning("MC Particle not found in tree, skipping"); 
-                               continue;
-                       }
-                       
-                       Int_t PDG =TMath::Abs(mcPart->PdgCode()); 
-                       if(!(PDG==310)) continue; // select only if k0 short
-               
-                       Double_t pT = mcPart->Pt();
-            Double_t d0 =1; // set 1 fot the moment - no displacement calculation implemented yet
-                       if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts
-                       
-                       KZeroClone->Add(new AliReducedParticle(mcPart->Eta(), mcPart->Phi(), pT,iPart,-1,d0,kTRUE,mcPart->Charge()));
-               }
-               
-       } // end if use  MC truth
-
-
-
-       return KZeroClone;
-}
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+//\r
+//\r
+//             Base class for Heavy Flavour Correlations Analysis\r
+//             Single Event and Mixed Event Analysis are implemented\r
+//\r
+//-----------------------------------------------------------------------\r
+//          \r
+//\r
+//                                                Author S.Bjelogrlic\r
+//                         Utrecht University \r
+//                      sandro.bjelogrlic@cern.ch\r
+//\r
+//-----------------------------------------------------------------------\r
+\r
+/* $Id$ */\r
+\r
+#include <TParticle.h>\r
+#include <TVector3.h>\r
+#include <TChain.h>\r
+#include "TROOT.h"\r
+#include "AliHFCorrelator.h"\r
+#include "AliRDHFCutsDStartoKpipi.h"\r
+#include "AliHFAssociatedTrackCuts.h"\r
+#include "AliEventPoolManager.h"\r
+#include "AliReducedParticle.h"\r
+#include "AliCentrality.h"\r
+#include "AliAODMCParticle.h"\r
+\r
+using std::cout;\r
+using std::endl;\r
+\r
+//_____________________________________________________\r
+AliHFCorrelator::AliHFCorrelator() :\r
+//\r
+// default constructor\r
+//\r
+TNamed(),\r
+fPoolMgr(0x0),         \r
+fPool(0x0),\r
+fhadcuts(0x0),\r
+fAODEvent(0x0),\r
+fAssociatedTracks(0x0),\r
+fmcArray(0x0),\r
+fReducedPart(0x0),\r
+fD0cand(0x0), \r
+fhypD0(0), \r
+fDCharge(0),\r
+\r
+fmixing(kFALSE),\r
+fmontecarlo(kFALSE),\r
+fsystem(kFALSE),\r
+fUseReco(kTRUE),\r
+fselect(kUndefined),\r
+\r
+fUseImpactParameter(0),\r
+fPIDmode(0),\r
+\r
+fNofTracks(0),\r
+fPoolContent(0),\r
+\r
+fPhiMin(0),\r
+fPhiMax(0),\r
+\r
+fPtTrigger(0),\r
+fPhiTrigger(0),\r
+fEtaTrigger(0),\r
+\r
+\r
+fDeltaPhi(0),\r
+fDeltaEta(0),\r
+fk0InvMass(0)\r
+\r
+{\r
+       // default constructor  \r
+}\r
+\r
+\r
+\r
+//_____________________________________________________\r
+AliHFCorrelator::AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb) :\r
+TNamed(name,"title"),\r
+fPoolMgr(0x0),         \r
+fPool(0x0),\r
+fhadcuts(0x0),\r
+fAODEvent(0x0),\r
+fAssociatedTracks(0x0),\r
+fmcArray(0x0),\r
+fReducedPart(0x0),\r
+fD0cand(0x0), \r
+fhypD0(0),\r
+fDCharge(0),\r
+\r
+fmixing(kFALSE),\r
+fmontecarlo(kFALSE),\r
+fsystem(ppOrPbPb),\r
+fUseReco(kTRUE),\r
+fselect(kUndefined),\r
+fUseImpactParameter(0),\r
+fPIDmode(0),\r
+\r
+fNofTracks(0),\r
+fPoolContent(0),\r
+\r
+fPhiMin(0),\r
+fPhiMax(0),\r
+\r
+fPtTrigger(0),\r
+fPhiTrigger(0),\r
+fEtaTrigger(0),\r
+\r
+\r
+fDeltaPhi(0),\r
+fDeltaEta(0),\r
+fk0InvMass(0)\r
+{\r
+       fhadcuts = cuts; \r
+}\r
+\r
+//_____________________________________________________\r
+AliHFCorrelator::~AliHFCorrelator() \r
+{\r
+//\r
+// destructor\r
+//     \r
+       \r
+       if(fPoolMgr)  {delete fPoolMgr; fPoolMgr=0;}       \r
+       if(fPool) {delete fPool; fPool=0;}\r
+       if(fhadcuts) {delete fhadcuts; fhadcuts=0;}\r
+       if(fAODEvent) {delete fAODEvent; fAODEvent=0;}\r
+       if(fAssociatedTracks) {delete fAssociatedTracks; fAssociatedTracks=0;}\r
+       if(fmcArray) {delete fmcArray; fmcArray=0;}\r
+       if(fReducedPart) {delete fReducedPart; fReducedPart=0;}\r
+       if(fD0cand) {delete fD0cand; fD0cand=0;}\r
+       \r
+       \r
+       if(fNofTracks) fNofTracks = 0;\r
+       \r
+       if(fPhiMin) fPhiMin = 0;\r
+       if(fPhiMax) fPhiMax = 0;\r
+       \r
+       if(fPtTrigger) fPtTrigger=0;\r
+       if(fPhiTrigger) fPhiTrigger=0;\r
+       if(fEtaTrigger) fEtaTrigger=0;\r
+       \r
+       if(fDeltaPhi) fDeltaPhi=0;\r
+       if(fDeltaEta) fDeltaEta=0;\r
+       \r
+       if(fk0InvMass) fk0InvMass=0;\r
+       \r
+}\r
+//_____________________________________________________\r
+Bool_t AliHFCorrelator::DefineEventPool(){\r
+       // definition of the Pool Manager for Event Mixing\r
+       \r
+\r
+       Int_t MaxNofEvents = fhadcuts->GetMaxNEventsInPool();\r
+       Int_t MinNofTracks = fhadcuts->GetMinNTracksInPool();\r
+       Int_t NofCentBins = fhadcuts->GetNCentPoolBins();\r
+       Double_t * CentBins = fhadcuts->GetCentPoolBins();\r
+       Int_t NofZVrtxBins = fhadcuts->GetNZvtxPoolBins();\r
+       Double_t *ZVrtxBins = fhadcuts->GetZvtxPoolBins();\r
+               \r
+                       \r
+       fPoolMgr = new AliEventPoolManager(MaxNofEvents, MinNofTracks, NofCentBins, CentBins, NofZVrtxBins, ZVrtxBins);\r
+       if(!fPoolMgr) return kFALSE;\r
+       return kTRUE;\r
+}\r
+//_____________________________________________________\r
+Bool_t AliHFCorrelator::Initialize(){\r
+       \r
+    //  std::cout << "AliHFCorrelator::Initialize"<< std::endl;\r
+//  AliInfo("AliHFCorrelator::Initialize") ;\r
+  if(!fAODEvent){\r
+    AliInfo("No AOD event") ;\r
+    return kFALSE;\r
+  }\r
+    //std::cout << "No AOD event" << std::endl;\r
+       \r
+       AliCentrality *centralityObj = 0;\r
+       //Int_t multiplicity = -1;\r
+       Double_t MultipOrCent = -1;\r
+       \r
+       // initialize the pool for event mixing\r
+       if(!fsystem){ // pp\r
+       //multiplicity = fAODEvent->GetNTracks();\r
+        MultipOrCent = AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(fAODEvent,-1.,1.);\r
+       //      MultipOrCent = multiplicity; // convert from Int_t to Double_t\r
+     //   AliInfo(Form("Multiplicity is %f", MultipOrCent));\r
+       }\r
+       if(fsystem){ // PbPb\r
+               \r
+               centralityObj = fAODEvent->GetHeader()->GetCentralityP();\r
+               MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");\r
+//             AliInfo(Form("Centrality is %f", MultipOrCent));\r
+       }\r
+       \r
+       AliAODVertex *vtx = fAODEvent->GetPrimaryVertex();\r
+       Double_t zvertex = vtx->GetZ(); // zvertex\r
+       Double_t * CentBins = fhadcuts->GetCentPoolBins();\r
+       Double_t poolmin=CentBins[0];\r
+       Double_t poolmax=CentBins[fhadcuts->GetNCentPoolBins()];\r
+\r
+       \r
+               if(TMath::Abs(zvertex)>=10 || MultipOrCent>poolmax || MultipOrCent < poolmin) {\r
+               if(!fsystem)AliInfo(Form("pp or pA Event with Zvertex = %.2f cm and multiplicity = %.0f out of pool bounds, SKIPPING",zvertex,MultipOrCent));\r
+               if(fsystem) AliInfo(Form("PbPb Event with Zvertex = %.2f cm and centrality = %.1f  out of pool bounds, SKIPPING",zvertex,MultipOrCent));\r
+\r
+                       return kFALSE;\r
+               }\r
+       \r
+       fPool = fPoolMgr->GetEventPool(MultipOrCent, zvertex);\r
+       \r
+       if (!fPool){\r
+               AliInfo(Form("No pool found for multiplicity = %f, zVtx = %f cm", MultipOrCent, zvertex));\r
+           return kFALSE;\r
+       }\r
+       //fPool->PrintInfo();\r
+       return kTRUE;\r
+}\r
+\r
+//_____________________________________________________\r
+Bool_t AliHFCorrelator::ProcessEventPool(){\r
+        // analysis on Mixed Events\r
+       //cout << "AliHFCorrelator::ProcessEventPool"<< endl;\r
+               if(!fmixing) return kFALSE;\r
+               if(!fPool->IsReady()) return kFALSE;\r
+               if(fPool->GetCurrentNEvents()<fhadcuts->GetMinEventsToMix()) return kFALSE;\r
+       //      fPool->PrintInfo();\r
+               fPoolContent = fPool->GetCurrentNEvents();\r
+               \r
+               return kTRUE;\r
+       \r
+}\r
+\r
+//_____________________________________________________\r
+Bool_t AliHFCorrelator::ProcessAssociatedTracks(Int_t EventLoopIndex, const TObjArray* associatedTracks){\r
+  // TODO: memory leak needs to be fixed, for every call, a new array\r
+  // is allocated, but the pointer immediately lost. The cleanup is\r
+  // not straightforward as in the case of event mixing the pointer\r
+  // will be an external array which must not be deleted.\r
+       fAssociatedTracks = new TObjArray();\r
+\r
+       if(!fmixing){ // analysis on Single Event\r
+               \r
+               \r
+               \r
+               if(fselect==kHadron || fselect ==kKaon) fAssociatedTracks = AcceptAndReduceTracks(fAODEvent);\r
+               if(fselect==kKZero) {fAssociatedTracks = AcceptAndReduceKZero(fAODEvent);}      \r
+               if(fselect==kElectron && associatedTracks) fAssociatedTracks=new TObjArray(*associatedTracks);\r
+               \r
+       }\r
+       \r
+       if(fmixing) { // analysis on Mixed Events\r
+               \r
+                       \r
+               fAssociatedTracks = fPool->GetEvent(EventLoopIndex);\r
+                               \r
+                               \r
+                       \r
+\r
+       } // end if mixing\r
+       \r
+       if(!fAssociatedTracks) return kFALSE;\r
+       \r
+       fNofTracks = fAssociatedTracks->GetEntriesFast(); \r
+               \r
+       return kTRUE;\r
+       \r
+}\r
+//_____________________________________________________\r
+Bool_t AliHFCorrelator::Correlate(Int_t loopindex){\r
+\r
+       if(loopindex >= fNofTracks) return kFALSE;\r
+       if(!fAssociatedTracks) return kFALSE;\r
+       \r
+       fReducedPart = (AliReducedParticle*)fAssociatedTracks->At(loopindex);\r
+       \r
+\r
+       fDeltaPhi = SetCorrectPhiRange(fPhiTrigger - fReducedPart->Phi());\r
+       \r
+       fDeltaEta = fEtaTrigger - fReducedPart->Eta();\r
+\r
+       return kTRUE;\r
+       \r
+}\r
+               \r
+//_____________________________________________________\r
+Bool_t AliHFCorrelator::PoolUpdate(const TObjArray* associatedTracks){\r
+\r
+       if(!fmixing) return kFALSE;\r
+       if(!fPool) return kFALSE;\r
+       if(fmixing) { // update the pool for Event Mixing\r
+               TObjArray* objArr = NULL;\r
+               if(fselect==kHadron || fselect==kKaon) objArr = (TObjArray*)AcceptAndReduceTracks(fAODEvent);\r
+               else if(fselect==kKZero) objArr = (TObjArray*)AcceptAndReduceKZero(fAODEvent);\r
+               else if(fselect==kElectron && associatedTracks){\r
+                 objArr = new TObjArray(*associatedTracks);\r
+               }\r
+               else return kFALSE;\r
+               if(objArr->GetEntriesFast()>0) fPool->UpdatePool(objArr); // updating the pool only if there are entries in the array\r
+       }\r
+               \r
+       return kTRUE;\r
+       \r
+}\r
+               \r
+//_____________________________________________________\r
+Double_t AliHFCorrelator::SetCorrectPhiRange(Double_t phi){\r
+       Double_t pi = TMath::Pi();\r
+       \r
+       if(phi<fPhiMin) phi = phi + 2*pi;\r
+       if(phi>fPhiMax) phi = phi - 2*pi;\r
+       \r
+       return phi;\r
+}\r
+\r
+//_____________________________________________________\r
+TObjArray*  AliHFCorrelator::AcceptAndReduceTracks(AliAODEvent* inputEvent){\r
+\r
+  Double_t weight=1.;\r
+  Int_t nTracks = inputEvent->GetNTracks();\r
+  AliAODVertex * vtx = inputEvent->GetPrimaryVertex();\r
+  Double_t pos[3],cov[6];\r
+  vtx->GetXYZ(pos);\r
+  vtx->GetCovarianceMatrix(cov);\r
+  const AliESDVertex vESD(pos,cov,100.,100);\r
+  \r
+  Double_t Bz = inputEvent->GetMagneticField();\r
+       \r
+  \r
+  TObjArray* tracksClone = new TObjArray;\r
+  tracksClone->SetOwner(kTRUE);\r
+  \r
+  //*******************************************************\r
+  // use reconstruction\r
+  if(fUseReco){\r
+    for (Int_t iTrack=0; iTrack<nTracks; ++iTrack) {\r
+      AliAODTrack* track = inputEvent->GetTrack(iTrack);\r
+      if (!track) continue;\r
+      if(!fhadcuts->IsHadronSelected(track,&vESD,Bz)) continue; // apply ESD level selections\r
+      if(!fhadcuts->Charge(fDCharge,track)) continue; // apply selection on charge, if required\r
+\r
+      Double_t pT = track->Pt();\r
+      \r
+      //compute impact parameter\r
+      Double_t d0z0[2],covd0z0[3];\r
+      Double_t d0=-999999.;\r
+      if(fUseImpactParameter) track->PropagateToDCA(vtx,Bz,100,d0z0,covd0z0);\r
+      else d0z0[0] = 1. ; // random number - be careful with the cuts you applied\r
+      \r
+      if(fUseImpactParameter==1) d0 = TMath::Abs(d0z0[0]); // use impact parameter\r
+      if(fUseImpactParameter==2) { // use impact parameter over resolution\r
+       if(TMath::Abs(covd0z0[0])>0.00000001) d0 = TMath::Abs(d0z0[0])/TMath::Sqrt(covd0z0[0]); \r
+       else d0 = -1.; // if the resoultion is Zero, rejects the track - to be on the safe side\r
+       \r
+      }\r
+      \r
+      if(fmontecarlo) {// THIS TO BE CHECKED\r
+        Int_t hadLabel = track->GetLabel();\r
+       if(hadLabel < 0) continue;      \r
+      }\r
+      \r
+      if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts\r
+      Bool_t rejectsoftpi = kTRUE;// TO BE CHECKED: DO WE WANT IT TO kTRUE AS A DEFAULT?\r
+      if(fD0cand && !fmixing) rejectsoftpi = fhadcuts->InvMassDstarRejection(fD0cand,track,fhypD0); // TO BE CHECKED: WHY NOT FOR EM?\r
+      \r
+      \r
+      if(fselect ==kKaon){     \r
+       if(!fhadcuts->CheckKaonCompatibility(track,fmontecarlo,fmcArray,fPIDmode)) continue; // check if it is a Kaon - data and MC\r
+      }\r
+      weight=fhadcuts->GetTrackWeight(pT,track->Eta(),pos[2]);\r
+      tracksClone->Add(new AliReducedParticle(track->Eta(), track->Phi(), pT,track->GetLabel(),track->GetID(),d0,rejectsoftpi,track->Charge(),weight));\r
+    } // end loop on tracks\r
+  } // end if use reconstruction kTRUE\r
+  \r
+  //*******************************************************\r
+  \r
+  //use MC truth\r
+  if(fmontecarlo && !fUseReco){\r
+    \r
+    for (Int_t iPart=0; iPart<fmcArray->GetEntriesFast(); iPart++) { \r
+      AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iPart));\r
+      if (!mcPart) {\r
+       AliWarning("MC Particle not found in tree, skipping"); \r
+       continue;\r
+      }\r
+      if(!mcPart->Charge()) continue; // consider only charged tracks\r
+      \r
+      Int_t PDG =TMath::Abs(mcPart->PdgCode()); \r
+if(fselect ==kHadron) {if(!((PDG==321)||(PDG==211)||(PDG==2212)||(PDG==13)||(PDG==11))) continue;} // select only if kaon, pion, proton, muon or electron\r
+      else if(fselect ==kKaon) {if(!(PDG==321)) continue;} // select only if kaon\r
+      else if(fselect ==kElectron) {if(!(PDG==11)) continue;} // select only if electron\r
+\r
+      Double_t pT = mcPart->Pt();\r
+      Double_t d0 =1; // set 1 fot the moment - no displacement calculation implemented yet\r
+      if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts\r
+      \r
+      tracksClone->Add(new AliReducedParticle(mcPart->Eta(), mcPart->Phi(), pT,iPart,-1,d0,kFALSE,mcPart->Charge()));\r
+    }\r
+    \r
+  } // end if use  MC truth\r
+  \r
+  \r
+  return tracksClone;\r
+}\r
+\r
+//_____________________________________________________\r
+TObjArray*  AliHFCorrelator::AcceptAndReduceKZero(AliAODEvent* inputEvent){\r
+       \r
+       Int_t nOfVZeros = inputEvent->GetNumberOfV0s();\r
+       TObjArray* KZeroClone = new TObjArray;\r
+       AliAODVertex *vertex1 = (AliAODVertex*)inputEvent->GetPrimaryVertex();\r
+\r
+ // use reconstruction         \r
+     if(fUseReco){\r
+       Int_t v0label = -1;\r
+       Int_t pdgDgK0toPipi[2] = {211,211};\r
+       Double_t mPDGK0=0.497614;//TDatabasePDG::Instance()->GetParticle(310)->Mass();\r
+       const Int_t size = inputEvent->GetNumberOfV0s();\r
+       Int_t prevnegID[size];\r
+       Int_t prevposID[size];\r
+       for(Int_t iv0 =0; iv0< nOfVZeros; iv0++){// loop on all v0 candidates\r
+               AliAODv0 *v0 = (static_cast<AliAODEvent*> (inputEvent))->GetV0(iv0);\r
+               if(!v0) {\r
+                 AliInfo(Form("is not a v0 at step, %d", iv0)) ;\r
+                 //cout << "is not a v0 at step " << iv0 << endl;\r
+                 continue;\r
+               }\r
+               \r
+               if(!fhadcuts->IsKZeroSelected(v0,vertex1)) continue; // check if it is a k0\r
+           \r
+               // checks to avoid double counting\r
+               Int_t negID = -999;\r
+               Int_t posID = -998;\r
+               //Int_t a = 0;\r
+               prevnegID[iv0] = -997;\r
+               prevposID[iv0]= -996;\r
+               negID = v0->GetNegID();\r
+               posID = v0->GetPosID();\r
+               Bool_t DoubleCounting = kFALSE;\r
+               \r
+               for(Int_t k=0; k<iv0; k++){\r
+                       if(((negID==prevnegID[k])&&(posID==prevposID[k]))||((negID==prevposID[k])&&(posID==prevnegID[k]))){\r
+                               DoubleCounting = kTRUE;\r
+                               //a=k;\r
+                               break;\r
+                       }//end if\r
+               } // end for\r
+               \r
+               if(DoubleCounting) continue;\r
+               else{ \r
+                       prevposID[iv0] = posID; \r
+                       prevnegID[iv0] = negID;\r
+               }\r
+               \r
+               if(fmontecarlo) v0label = v0->MatchToMC(310,fmcArray, 0, pdgDgK0toPipi); //match a K0 short\r
+               Double_t k0pt = v0->Pt();\r
+               Double_t k0eta = v0->Eta();\r
+               Double_t k0Phi = v0->Phi();\r
+           fk0InvMass = v0->MassK0Short();     \r
+               \r
+               //if (loopindex == 0) {\r
+               //      if(!plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectra"))->Fill(k0InvMass,k0pt); // spectra for all k0\r
+               //      if(plotassociation) ((TH2F*)fOutput->FindObject("KZeroSpectraifHF"))->Fill(k0InvMass,k0pt); // spectra for k0 in association with a D*\r
+               //}\r
+               // if there are more D* candidates per event, loopindex == 0 makes sure you fill the mass spectra only once!\r
+               \r
+               if(TMath::Abs(fk0InvMass-mPDGK0)>3*0.004) continue; // select candidates within 3 sigma\r
+               KZeroClone->Add(new AliReducedParticle(k0eta,k0Phi,k0pt,v0label));\r
+               \r
+       }\r
+     } // end if use reconstruction kTRUE\r
+       \r
+\r
+\r
+//*********************************************************************//\r
+     //use MC truth\r
+     if(fmontecarlo && !fUseReco){\r
+               \r
+               for (Int_t iPart=0; iPart<fmcArray->GetEntriesFast(); iPart++) { \r
+                       AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iPart));\r
+                       if (!mcPart) {\r
+                               AliWarning("MC Particle not found in tree, skipping"); \r
+                               continue;\r
+                       }\r
+                       \r
+                       Int_t PDG =TMath::Abs(mcPart->PdgCode()); \r
+                       if(!(PDG==310)) continue; // select only if k0 short\r
+               \r
+                       Double_t pT = mcPart->Pt();\r
+            Double_t d0 =1; // set 1 fot the moment - no displacement calculation implemented yet\r
+                       if(!fhadcuts->CheckHadronKinematic(pT,d0)) continue; // apply kinematic cuts\r
+                       \r
+                       KZeroClone->Add(new AliReducedParticle(mcPart->Eta(), mcPart->Phi(), pT,iPart,-1,d0,kTRUE,mcPart->Charge()));\r
+               }\r
+               \r
+       } // end if use  MC truth\r
+\r
+\r
+\r
+       return KZeroClone;\r
+}\r
index 3c2b64df7b6b217b3ff4a8d263d92bdf087b79d3..4ca6f5c239422b58880522e2362134b889ca0898 100644 (file)
-#ifndef AliHFCorrelator_H
-#define AliHFCorrelator_H
-
-/**************************************************************************
- * Copyright(c) 1998-2009, 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.                  *
- **************************************************************************/
-
-//
-//             Base class for Heavy Flavour Correlations Analysis
-//             Single Event and Mixed Event Analysis are implemented
-//-----------------------------------------------------------------------
-//          
-//
-//                                                Author S.Bjelogrlic
-//                         Utrecht University 
-//                      sandro.bjelogrlic@cern.ch
-//
-//-----------------------------------------------------------------------
-
-/* $Id$ */
-
-#include "AliHFAssociatedTrackCuts.h"
-#include "AliEventPoolManager.h"
-#include "AliVParticle.h"
-#include "AliReducedParticle.h"
-
-
-class AliHFCorrelator : public TNamed
-{
-       
- public:
-       
-       AliHFCorrelator();
-       AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb);
-       virtual ~AliHFCorrelator();
-       
-       // enum for setting which associated particle type to work with
-       enum{
-         kUndefined=0,
-         kHadron=1,
-         kKaon,
-         kKZero,
-         kElectron
-       };
-
-       //setters
-       void SetDeltaPhiInterval (Double_t min, Double_t max){
-               fPhiMin = min; fPhiMax = max;
-               if(TMath::Abs(fPhiMin-fPhiMax) != 2*TMath::Pi()) AliInfo("AliHFCorrelator::Warning: the delta phi interval is not set to 2 Pi");
-       }
-       void SetEventMixing(Bool_t mixON){fmixing=mixON;}
-       void SetTriggerParticleProperties(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig)
-       {fPtTrigger = ptTrig; fPhiTrigger = phiTrig; fEtaTrigger = etaTrig;}
-       void SetTriggerParticleDaughterCharge(Short_t charge) {fDCharge=charge;}
-       
-       
-       void SetAssociatedParticleType(Int_t type){fselect = type;}
-       void SetAODEvent(AliAODEvent* inputevent){fAODEvent = inputevent;}
-       void SetMCArray(TClonesArray* mcArray){fmcArray = mcArray;}
-       void SetUseMC(Bool_t useMC){fmontecarlo = useMC;}
-       void SetApplyDisplacementCut(Int_t applycut){fUseImpactParameter = applycut;}
-       void SetPIDmode(Int_t mode){fPIDmode = mode;}
-       
-       void SetD0Properties(AliAODRecoDecayHF2Prong* d, Int_t D0hyp)
-       {fD0cand = d; fhypD0 = D0hyp;}
-       
-       void SetUseReco(Bool_t useReco) {fUseReco = useReco;}
-       
-       
-       
-        Bool_t DefineEventPool(); // Definition of the Event pool parameters
-       Bool_t Initialize(); // function that initlize everything for the analysis      
-       Bool_t ProcessEventPool(); // processes the event pool
-       Bool_t ProcessAssociatedTracks(Int_t EventLoopIndex, const TObjArray* associatedTracks=NULL); //
-       Bool_t Correlate(Int_t loopindex); // function that computes the correlations between the trigger particle and the track n. loopindex
-       Bool_t PoolUpdate(const TObjArray* associatedTracks=NULL);// updates the event pool
-       Double_t SetCorrectPhiRange(Double_t phi); // sets all the angles in the correct range
-       void SetPidAssociated() {fhadcuts->SetPidAssociated();}
-
-       //getters
-       AliEventPool* GetPool() {return fPool;}
-       TObjArray * GetTrackArray(){return fAssociatedTracks;}
-       AliHFAssociatedTrackCuts* GetSelectionCuts() {return fhadcuts;}
-       AliReducedParticle* GetAssociatedParticle() {return fReducedPart;}
-       
-       Int_t GetNofTracks(){return fNofTracks;}
-       Int_t GetNofEventsInPool(){return fPoolContent;}
-
-       Double_t GetDeltaPhi(){return fDeltaPhi;} // Delta Phi, needs to be called after the method correlate 
-       Double_t GetDeltaEta(){return fDeltaEta;} // Delta Eta
-       
-       Double_t GetAssociatedKZeroInvariantmass(){return fk0InvMass;}
-       
-       
-       
-       // methods to reduce the tracks to correlate with track selection cuts applied here
-       TObjArray*  AcceptAndReduceTracks(AliAODEvent* inputEvent); // selecting hadrons and kaons
-       TObjArray*  AcceptAndReduceKZero(AliAODEvent* inputEvent); // selecting kzeros
-       
-       
- private:
-
-       AliHFCorrelator(const AliHFCorrelator& vtxr);
-       AliHFCorrelator& operator=(const AliHFCorrelator& vtxr );
-
-       AliEventPoolManager* fPoolMgr;         //! event pool manager
-       AliEventPool * fPool; //! Pool for event mixing
-       AliHFAssociatedTrackCuts* fhadcuts;//! hadron cuts
-       AliAODEvent * fAODEvent;//! AOD Event
-       TObjArray* fAssociatedTracks; // Array of associated tracks
-       TClonesArray* fmcArray; //mcarray
-       AliReducedParticle * fReducedPart; // reduced AOD particle;
-       AliAODRecoDecayHF2Prong* fD0cand; //D0 candidate
-       Int_t fhypD0; //hypothesis necessary for
-       Int_t fDCharge; // charge of a daughter of the D meson
-       
-       Bool_t fmixing;// switch for event mixing
-       Bool_t fmontecarlo; // switch for MonteCarlo
-       Bool_t fsystem; // select pp (kFALSE) or PbPb (kTRUE)
-       Bool_t fUseReco; // switch to use reconstruction (kTRUE) or MC truth (kFALSE)
-       
-       Int_t fselect; // 1 for hadrons, 2 for kaons, 3 for KZeros
-       Int_t fUseImpactParameter; // switch to use the impact parameter cut
-       Int_t fPIDmode; // set the PID mode for Kaon identification
-       
-       Int_t fNofTracks; // number of tracks in track array
-       Int_t fPoolContent; //  n of events in pool
-       
-       Double_t fPhiMin; // min for phi
-       Double_t fPhiMax; // max for phi
-       
-       Double_t fPtTrigger; // pt of the trigger D meson
-       Double_t fPhiTrigger; // phi of the trigger D meson
-       Double_t fEtaTrigger; // Eta of the trigger D meson
-
-       
-       Double_t fDeltaPhi; // delta phi between D meson and associated track
-       Double_t fDeltaEta; // delta eta between D meson and associated track
-       
-       Double_t fk0InvMass; // KZero invariant mass
-       
-       
-       ClassDef(AliHFCorrelator,3); // class for HF correlations       
-};
-
-
-
-
-
-#endif
+#ifndef AliHFCorrelator_H\r
+#define AliHFCorrelator_H\r
+\r
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+//\r
+//             Base class for Heavy Flavour Correlations Analysis\r
+//             Single Event and Mixed Event Analysis are implemented\r
+//-----------------------------------------------------------------------\r
+//          \r
+//\r
+//                                                Author S.Bjelogrlic\r
+//                         Utrecht University \r
+//                      sandro.bjelogrlic@cern.ch\r
+//\r
+//-----------------------------------------------------------------------\r
+\r
+/* $Id$ */\r
+\r
+#include "AliHFAssociatedTrackCuts.h"\r
+#include "AliEventPoolManager.h"\r
+#include "AliVParticle.h"\r
+#include "AliReducedParticle.h"\r
+#include "AliVertexingHFUtils.h"\r
+\r
+\r
+class AliHFCorrelator : public TNamed\r
+{\r
+       \r
+ public:\r
+       \r
+       AliHFCorrelator();\r
+       AliHFCorrelator(const Char_t* name, AliHFAssociatedTrackCuts *cuts, Bool_t ppOrPbPb);\r
+       virtual ~AliHFCorrelator();\r
+       \r
+       // enum for setting which associated particle type to work with\r
+       enum{\r
+         kUndefined=0,\r
+         kHadron=1,\r
+         kKaon,\r
+         kKZero,\r
+         kElectron\r
+       };\r
+\r
+       //setters\r
+       void SetDeltaPhiInterval (Double_t min, Double_t max){\r
+               fPhiMin = min; fPhiMax = max;\r
+               if(TMath::Abs(fPhiMin-fPhiMax) != 2*TMath::Pi()) AliInfo("AliHFCorrelator::Warning: the delta phi interval is not set to 2 Pi");\r
+       }\r
+       void SetEventMixing(Bool_t mixON){fmixing=mixON;}\r
+       void SetTriggerParticleProperties(Double_t ptTrig, Double_t phiTrig, Double_t etaTrig)\r
+       {fPtTrigger = ptTrig; fPhiTrigger = phiTrig; fEtaTrigger = etaTrig;}\r
+       void SetTriggerParticleDaughterCharge(Short_t charge) {fDCharge=charge;}\r
+       \r
+       \r
+       void SetAssociatedParticleType(Int_t type){fselect = type;}\r
+       void SetAODEvent(AliAODEvent* inputevent){fAODEvent = inputevent;}\r
+       void SetMCArray(TClonesArray* mcArray){fmcArray = mcArray;}\r
+       void SetUseMC(Bool_t useMC){fmontecarlo = useMC;}\r
+       void SetApplyDisplacementCut(Int_t applycut){fUseImpactParameter = applycut;}\r
+       void SetPIDmode(Int_t mode){fPIDmode = mode;}\r
+       \r
+       void SetD0Properties(AliAODRecoDecayHF2Prong* d, Int_t D0hyp)\r
+       {fD0cand = d; fhypD0 = D0hyp;}\r
+       \r
+       void SetUseReco(Bool_t useReco) {fUseReco = useReco;}\r
+       \r
+       \r
+       \r
+        Bool_t DefineEventPool(); // Definition of the Event pool parameters\r
+       Bool_t Initialize(); // function that initlize everything for the analysis      \r
+       Bool_t ProcessEventPool(); // processes the event pool\r
+       Bool_t ProcessAssociatedTracks(Int_t EventLoopIndex, const TObjArray* associatedTracks=NULL); //\r
+       Bool_t Correlate(Int_t loopindex); // function that computes the correlations between the trigger particle and the track n. loopindex\r
+       Bool_t PoolUpdate(const TObjArray* associatedTracks=NULL);// updates the event pool\r
+       Double_t SetCorrectPhiRange(Double_t phi); // sets all the angles in the correct range\r
+       void SetPidAssociated() {fhadcuts->SetPidAssociated();}\r
+\r
+       //getters\r
+       AliEventPool* GetPool() {return fPool;}\r
+       TObjArray * GetTrackArray(){return fAssociatedTracks;}\r
+       AliHFAssociatedTrackCuts* GetSelectionCuts() {return fhadcuts;}\r
+       AliReducedParticle* GetAssociatedParticle() {return fReducedPart;}\r
+       \r
+       Int_t GetNofTracks(){return fNofTracks;}\r
+       Int_t GetNofEventsInPool(){return fPoolContent;}\r
+\r
+       Double_t GetDeltaPhi(){return fDeltaPhi;} // Delta Phi, needs to be called after the method correlate \r
+       Double_t GetDeltaEta(){return fDeltaEta;} // Delta Eta\r
+       \r
+       Double_t GetAssociatedKZeroInvariantmass(){return fk0InvMass;}\r
+       \r
+       \r
+       \r
+       // methods to reduce the tracks to correlate with track selection cuts applied here\r
+       TObjArray*  AcceptAndReduceTracks(AliAODEvent* inputEvent); // selecting hadrons and kaons\r
+       TObjArray*  AcceptAndReduceKZero(AliAODEvent* inputEvent); // selecting kzeros\r
+       \r
+       \r
+ private:\r
+\r
+       AliHFCorrelator(const AliHFCorrelator& vtxr);\r
+       AliHFCorrelator& operator=(const AliHFCorrelator& vtxr );\r
+\r
+       AliEventPoolManager* fPoolMgr;         //! event pool manager\r
+       AliEventPool * fPool; //! Pool for event mixing\r
+       AliHFAssociatedTrackCuts* fhadcuts;//! hadron cuts\r
+       AliAODEvent * fAODEvent;//! AOD Event\r
+       TObjArray* fAssociatedTracks; // Array of associated tracks\r
+       TClonesArray* fmcArray; //mcarray\r
+       AliReducedParticle * fReducedPart; // reduced AOD particle;\r
+       AliAODRecoDecayHF2Prong* fD0cand; //D0 candidate\r
+       Int_t fhypD0; //hypothesis necessary for\r
+       Int_t fDCharge; // charge of a daughter of the D meson\r
+       \r
+       Bool_t fmixing;// switch for event mixing\r
+       Bool_t fmontecarlo; // switch for MonteCarlo\r
+       Bool_t fsystem; // select pp (kFALSE) or PbPb (kTRUE)\r
+       Bool_t fUseReco; // switch to use reconstruction (kTRUE) or MC truth (kFALSE)\r
+       \r
+       Int_t fselect; // 1 for hadrons, 2 for kaons, 3 for KZeros\r
+       Int_t fUseImpactParameter; // switch to use the impact parameter cut\r
+       Int_t fPIDmode; // set the PID mode for Kaon identification\r
+       \r
+       Int_t fNofTracks; // number of tracks in track array\r
+       Int_t fPoolContent; //  n of events in pool\r
+       \r
+       Double_t fPhiMin; // min for phi\r
+       Double_t fPhiMax; // max for phi\r
+       \r
+       Double_t fPtTrigger; // pt of the trigger D meson\r
+       Double_t fPhiTrigger; // phi of the trigger D meson\r
+       Double_t fEtaTrigger; // Eta of the trigger D meson\r
+\r
+       \r
+       Double_t fDeltaPhi; // delta phi between D meson and associated track\r
+       Double_t fDeltaEta; // delta eta between D meson and associated track\r
+       \r
+       Double_t fk0InvMass; // KZero invariant mass\r
+       \r
+       \r
+       ClassDef(AliHFCorrelator,3); // class for HF correlations       \r
+};\r
+\r
+\r
+\r
+\r
+\r
+#endif\r
index 6cf9fcb6dfa169da7eaf865a2d829c960d3548e2..08f97f2ff9456b3e203eede0d93a70811e719df0 100644 (file)
@@ -6,7 +6,7 @@
 AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(AliAnalysisTaskDStarCorrelations::CollSyst syst,\r
                                                                      Bool_t theMCon, Bool_t mixing, Bool_t UseReco = kTRUE, Bool_t fullmode = kFALSE ,Bool_t UseEffic = kFALSE, Bool_t UseDEffic = kFALSE ,\r
                                                                      Int_t trackselect =1, Int_t usedispl =0, Int_t nbins, Float_t DStarSigma, Float_t D0Sigma, Float_t D0SBSigmaLow, Float_t D0SBSigmaHigh, \r
-                                                                     TString DStarCutsFile, TString TrackCutsFile, TString DStarEffMap, TString TrackEffMap,\r
+                                                                     TString DStarCutsFile, TString TrackCutsFile,\r
                                                                      Int_t tasknumber = 0)\r
 {\r
 \r
@@ -68,8 +68,7 @@ AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(AliAnalysisTaskDStarC
     \r
     cout << "DStar cut object:     " << DStarCutsFile << endl;\r
     cout << "Tracks cut object:    " << TrackCutsFile << endl;\r
-    cout << "Track Efficiency map: " << TrackEffMap << endl;\r
-    cout << "Dmeson Efficiency map:" << DStarEffMap << endl;\r
+  \r
     \r
     cout << "==========================================================" << endl;\r
     //gSystem->Sleep(2000);\r
@@ -111,86 +110,13 @@ AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(AliAnalysisTaskDStarC
        AliHFAssociatedTrackCuts* corrCuts=new AliHFAssociatedTrackCuts();\r
        corrCuts = (AliHFAssociatedTrackCuts*)filecuts2->Get("AssociatedCuts");\r
        corrCuts->SetName("AssociatedCuts");\r
-//     corrCuts->PrintAll();\r
+       corrCuts->PrintAll();\r
        if(!corrCuts){\r
                cout<<"Specific associated track cuts not found"<<endl;\r
                return;\r
        }\r
        \r
-       \r
-// ******************************** OPENING THE EFFICIENCY MAPS  ************************************\r
-   \r
-    cout << "Getting Efficiency map object from file \n" << TrackEffMap.Data() << "\n "<<  endl;\r
-    TFile* effFile=new TFile(TrackEffMap.Data());\r
-    if(UseEffic && !effFile->IsOpen()){\r
-        cout<<"Input file not found for efficiency! Exiting..."<<endl;\r
-        return;\r
-    }\r
-    \r
-    \r
-    \r
-    TCanvas *c = (TCanvas*)effFile->Get("c");\r
-    if(!c) {cout << "No canvas !" << endl;}\r
-    TH3D *h3D = (TH3D*)c->FindObject("heff_rebin");\r
\r
-    \r
 \r
-    if(!h3D) {cout << "No Single track efficiency histo " << endl; ;}\r
-       \r
-    corrCuts->SetEfficiencyWeightMap(h3D);\r
-    \r
-    \r
-    \r
-    \r
-    \r
-    // ******************************** OPENING THE DMeson EFFICIENCY MAPS  ************************************\r
-    \r
-    cout << "Getting Efficiency map object from file \n" << DStarEffMap.Data() << "\n "<<  endl;\r
-    TFile* DeffFile=new TFile(DStarEffMap.Data());\r
-    \r
-    if(!DeffFile){\r
-        cout << "No D efficiencymapinput" << endl;\r
-        gSystem->Sleep(1000);\r
-    }\r
-    \r
-    if(UseDEffic && !DeffFile->IsOpen()){\r
-        cout<<"Input file not found for Dmeson efficiency! Exiting..."<<endl;\r
-        gSystem->Sleep(1000);\r
-        return;\r
-    }\r
-    \r
-    \r
- //    gSystem->Sleep(7000);\r
-    TString canvasname = "DMesonEfficiency";\r
-  //  if(DStarEffMap == "~/Analysis/HFCorrelations/DStarEfficiencies/DStarEfficiency_reco_genacc.root") canvasname = "c1";\r
-  //  else canvasname ="correctefficiency";\r
-    \r
-    TString histoname = "DMesonEffMap";\r
-    \r
-    \r
-   // if(DStarEffMap == "~/Analysis/HFCorrelations/DStarEfficiencies/DStarEfficiency_reco_genacc.root") histoname = "Reco_GenAcc";\r
-   // else histoname ="effacc";\r
-      \r
-    \r
-    \r
-    TCanvas *c1 = (TCanvas*)DeffFile->Get(canvasname.Data());\r
-    if (!c1) cout << "No canvas with name " << canvasname <<  endl;\r
-    TH1D *effMap = (TH1D*)c1->FindObject(histoname.Data());\r
-     if (!effMap) cout << "No map " << endl;\r
-  //  TH1D *effMap = (TH1D*)c1->FindObject("EffMap");\r
-    TH2D *effMap2d = (TH2D*)c1->FindObject("EffMap");\r
-  \r
-  //  effMap->Draw("ep");\r
-  //  return;\r
-    \r
-    \r
-    if(!c1) {cout << "No canvas for D eff" << endl; ;}\r
-    if(!effMap) {cout << "No histo for D eff " << endl;}\r
-    if(!effMap2d) {cout << "No 2D histo for D eff " << endl;}\r
-       \r
-    \r
-   // gSystem->Sleep(6000);\r
-       \r
 // ******************************** SELECTING THE MC PROCESS  ************************************\r
 \r
        TString selectMCproc = "";\r
@@ -220,8 +146,8 @@ AliAnalysisTaskDStarCorrelations *AddTaskDStarCorrelations(AliAnalysisTaskDStarC
        \r
     \r
    // cout << "Adding efficiency map to Assoc track cut object \n" << endl;\r
-    if(effMap) task->SetDeffMapvsPt(effMap);\r
-    if(effMap2d) task->SetDeffMapvsPtvsMult(effMap2d);\r
+  //  if(effMap) task->SetDeffMapvsPt(effMap);\r
+  //  if(effMap2d) task->SetDeffMapvsPtvsMult(effMap2d);\r
     \r
     task->SetNofPhiBins(nbins);\r
        task->SetMonteCarlo(theMCon);\r