Adding task for like-sign bkg for D0->Kpi (Carmelo)
authordainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 18 Jun 2009 11:40:34 +0000 (11:40 +0000)
committerdainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 18 Jun 2009 11:40:34 +0000 (11:40 +0000)
PWG3/CMake_libPWG3vertexingHF.txt
PWG3/PWG3vertexingHFLinkDef.h
PWG3/libPWG3vertexingHF.pkg
PWG3/vertexingHF/AddTaskBkgLikeSignD0.C [new file with mode: 0644]
PWG3/vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.cxx [new file with mode: 0644]
PWG3/vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.h [new file with mode: 0644]
PWG3/vertexingHF/RunAnalysisAODVertexingHF.C

index 9f900b3..3f76467 100644 (file)
@@ -11,6 +11,7 @@ set(SRCS
        vertexingHF/AliAnalysisTaskSECompareHF.cxx
        vertexingHF/AliAnalysisTaskSEDplus.cxx
        vertexingHF/AliAnalysisTaskSEBkgLikeSignJPSI.cxx
+       vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.cxx
        vertexingHF/AliCFHeavyFlavourTask.cxx
        vertexingHF/AliCFHeavyFlavourTaskMultiVar.cxx
        vertexingHF/AliCFHeavyFlavourTaskMultiVarMultiStep.cxx
index c031390..aefa12d 100644 (file)
@@ -22,6 +22,7 @@
 #pragma link C++ class AliSignificanceCalculator+;
 #pragma link C++ class AliHFMassFitter+;
 #pragma link C++ class AliAnalysisTaskSEBkgLikeSignJPSI+;
+#pragma link C++ class AliAnalysisTaskSEBkgLikeSignD0+;
 #pragma link C++ class AliAnalysisBtoJPSItoEle+;
 #pragma link C++ class AliAnalysisTaskSEBtoJPSItoEle+;
 #pragma link C++ class AliBtoJPSItoEleCDFfitFCN+;
index 5e9438a..07fc50c 100644 (file)
@@ -15,6 +15,7 @@ SRCS:=   vertexingHF/AliAODRecoDecayHF.cxx \
        vertexingHF/AliMultiDimVector.cxx vertexingHF/AliSignificanceCalculator.cxx \
        vertexingHF/AliHFMassFitter.cxx \
        vertexingHF/AliAnalysisTaskSEBkgLikeSignJPSI.cxx \
+       vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.cxx \
        vertexingHF/AliAnalysisBtoJPSItoEle.cxx \
        vertexingHF/AliAnalysisTaskSEBtoJPSItoEle.cxx \
        vertexingHF/AliBtoJPSItoEleCDFfitFCN.cxx \
diff --git a/PWG3/vertexingHF/AddTaskBkgLikeSignD0.C b/PWG3/vertexingHF/AddTaskBkgLikeSignD0.C
new file mode 100644 (file)
index 0000000..f374157
--- /dev/null
@@ -0,0 +1,36 @@
+AliAnalysisTaskSEBkgLikeSignD0 *AddTaskBkgLikeSignD0() 
+{
+  //
+  // Test macro for the AliAnalysisTaskSEBkgLikeSignD0 
+  // starting from AliAOD.root file with HF + Like Sign candidates.
+  // C.Di Giglio, carmelo.digiglio@ba.infn.it
+  //
+
+
+  // Get the pointer to the existing analysis manager via the static access method.
+  //==============================================================================
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+    ::Error("AddTaskBkgLikeSignD0", "No analysis manager to connect to.");
+    return NULL;
+  }   
+
+  // Like-sign background analysis task    
+  AliAnalysisTaskSEBkgLikeSignD0 *lsD0Task = new AliAnalysisTaskSEBkgLikeSignD0("CmpLikeSignD0Analysis");
+  lsD0Task->SetDebugLevel(2);
+
+  mgr->AddTask(lsD0Task);
+
+  //
+  // Create containers for input/output
+  AliAnalysisDataContainer *cinputLSD0 = mgr->CreateContainer("cinputLikeSignD0",TChain::Class(), 
+                                                         AliAnalysisManager::kInputContainer);
+  AliAnalysisDataContainer *coutputLSD0 = mgr->CreateContainer("coutputLikeSignD0",TList::Class(),
+                                                           AliAnalysisManager::kOutputContainer,
+                                                           "CmpLikeSignD0.root");
+
+  mgr->ConnectInput(lsD0Task,0,mgr->GetCommonInputContainer());
+  mgr->ConnectOutput(lsD0Task,1,coutputLSD0);
+
+  return lsD0Task;
+}
diff --git a/PWG3/vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.cxx b/PWG3/vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.cxx
new file mode 100644 (file)
index 0000000..c581cdf
--- /dev/null
@@ -0,0 +1,354 @@
+/**************************************************************************
+ * 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.                  *
+ **************************************************************************/
+
+///////////////////////////////////////////////////////////////////////////
+//
+// AliAnalysisTaskSE for reading both reconstructed D0->Kpi candidates
+// and like sign pairs and for drawing corresponding distributions
+//
+// Author: C.Di Giglio, carmelo.digiglio@ba.infn.it
+///////////////////////////////////////////////////////////////////////////
+
+#include <TSystem.h>
+#include <TROOT.h>
+#include <TClonesArray.h>
+#include <TList.h>
+#include <TH1F.h>
+
+#include "AliAODEvent.h"
+#include "AliAODVertex.h"
+#include "AliAODTrack.h"
+#include "AliAODRecoDecayHF2Prong.h"
+#include "AliAnalysisVertexingHF.h"
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisTaskSEBkgLikeSignD0.h"
+
+ClassImp(AliAnalysisTaskSEBkgLikeSignD0)
+
+//________________________________________________________________________
+AliAnalysisTaskSEBkgLikeSignD0::AliAnalysisTaskSEBkgLikeSignD0():
+AliAnalysisTaskSE(),
+fOutput(0), 
+fHistMassD0(0),
+fHistMassLS(0),
+fHistCtsD0(0),           
+fHistCtsLS(0),
+fHistCtsLSpos(0),
+fHistCtsLSneg(0),
+fHistCPtaD0(0),          
+fHistCPtaLS(0),
+fHistd0d0D0(0),          
+fHistd0d0LS(0),
+fHistDCAD0(0),           
+fHistDCALS(0),
+fVHF(0),
+fTotPosPairs(0),
+fTotNegPairs(0),
+fLsNormalization(1.)
+{
+  //
+  // Default constructor
+  //
+}
+
+//________________________________________________________________________
+AliAnalysisTaskSEBkgLikeSignD0::AliAnalysisTaskSEBkgLikeSignD0(const char *name):
+AliAnalysisTaskSE(name),
+fOutput(0),
+fHistMassD0(0),
+fHistMassLS(0),
+fHistCtsD0(0),
+fHistCtsLS(0),
+fHistCtsLSpos(0),
+fHistCtsLSneg(0),
+fHistCPtaD0(0),
+fHistCPtaLS(0),
+fHistd0d0D0(0),
+fHistd0d0LS(0),
+fHistDCAD0(0),
+fHistDCALS(0),
+fVHF(0),
+fTotPosPairs(0),
+fTotNegPairs(0),
+fLsNormalization(1.)
+{
+  //
+  // Standard constructor
+  //
+  // Output slot #1 writes into a TList container
+  DefineOutput(1,TList::Class());  //My private output
+}
+
+//________________________________________________________________________
+AliAnalysisTaskSEBkgLikeSignD0::~AliAnalysisTaskSEBkgLikeSignD0()
+{
+  // Destructor
+  if (fOutput) {
+    delete fOutput;
+    fOutput = 0;
+  }
+
+  if (fVHF) {
+    delete fVHF;
+    fVHF = 0;
+  }
+
+}
+//________________________________________________________________________
+void AliAnalysisTaskSEBkgLikeSignD0::Init()
+{
+  // Initialization
+
+  if(fDebug > 1) printf("AnalysisTaskSEBkgLikeSignD0::Init() \n");
+
+  gROOT->LoadMacro("ConfigVertexingHF.C");
+
+  fVHF = (AliAnalysisVertexingHF*)gROOT->ProcessLine("ConfigVertexingHF()");
+  fVHF->PrintStatus();
+
+  return;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSEBkgLikeSignD0::UserCreateOutputObjects()
+{
+  // Create the output container
+  //
+  if(fDebug > 1) printf("AnalysisTaskSEBkgLikeSignD0::UserCreateOutputObjects() \n");
+
+  // Several histograms are more conveniently managed in a TList
+  fOutput = new TList();
+  fOutput->SetOwner();
+
+  fHistMassD0 = new TH1F("fHistMassD0", "D0 invariant mass; M [GeV]; Entries",200,1.56,2.16);
+  fHistMassD0->Sumw2();
+  fHistMassD0->SetMinimum(0);
+  fOutput->Add(fHistMassD0);
+
+  fHistMassLS = new TH1F("fHistMassLS", "Like sign pairs invariant mass; M [GeV]; Entries",200,1.56,2.16);
+  fHistMassLS->Sumw2();
+  fHistMassLS->SetMinimum(0);
+  fOutput->Add(fHistMassLS);
+
+  fHistCtsD0 = new TH1F("fHistCtsD0", "D0 cosine of decay angle; Cos#Theta*; Entries",200,-1.,1.);
+  fHistCtsD0->Sumw2();
+  fHistCtsD0->SetMinimum(0);
+  fOutput->Add(fHistCtsD0);
+
+  fHistCtsLS = new TH1F("fHistCtsLS", "Like sign pairs cosine of decay angle; Cos#Theta*; Entries",200,-1.,1.);
+  fHistCtsLS->Sumw2();
+  fHistCtsLS->SetMinimum(0);
+  fOutput->Add(fHistCtsLS);
+
+  fHistCtsLSpos = new TH1F("fHistCtsLSpos", "Like sign ++ pairs cosine of decay angle; Cos#Theta*; Entries",200,-1.,1.);
+  fHistCtsLSpos->Sumw2();
+  fHistCtsLSpos->SetMinimum(0);
+  fOutput->Add(fHistCtsLSpos);
+
+  fHistCtsLSneg = new TH1F("fHistCtsLSneg", "Like sign -- pairs cosine of decay angle; Cos#Theta*; Entries",200,-1.,1.);
+  fHistCtsLSneg->Sumw2();
+  fHistCtsLSneg->SetMinimum(0);
+  fOutput->Add(fHistCtsLSneg);
+
+  fHistCPtaD0 = new TH1F("fHistCPtaD0", "D0 cosine of pointing angle; Cos#Theta_{point}; Entries",200,-1.,1.);
+  fHistCPtaD0->Sumw2();
+  fHistCPtaD0->SetMinimum(0);
+  fOutput->Add(fHistCPtaD0);
+
+  fHistCPtaLS = new TH1F("fHistCPtaLS", "Like sign pairs cosine of pointing angle; Cos#Theta_{point}; Entries",200,-1.,1.);
+  fHistCPtaLS->Sumw2();
+  fHistCPtaLS->SetMinimum(0);
+  fOutput->Add(fHistCPtaLS);
+
+  fHistd0d0D0 = new TH1F("fHistd0d0D0", "D0 product of impact parameters; d0xd0 [#mu m^{2}]; Entries",200,-100000.,100000.);
+  fHistd0d0D0->Sumw2(); 
+  fHistd0d0D0->SetMinimum(0);
+  fOutput->Add(fHistd0d0D0);
+
+  fHistd0d0LS = new TH1F("fHistd0d0LS", "Like sign pairs product of impact parameters; d0xd0 [#mu m^{2}]; Entries",200,-100000.,100000.);
+  fHistd0d0LS->Sumw2();
+  fHistd0d0LS->SetMinimum(0);
+  fOutput->Add(fHistd0d0LS);
+
+  fHistDCAD0 = new TH1F("fHistDCAD0", "D0 distance of closest approach; dca [10^{2}#mu m]; Entries",100,0.,5.);
+  fHistDCAD0->Sumw2(); 
+  fHistDCAD0->SetMinimum(0);
+  fOutput->Add(fHistDCAD0);
+
+  fHistDCALS = new TH1F("fHistDCALS", "Like sign pairs distance of closest approach; dca [10^{2}#mu m]; Entries",100,0.,5.);
+  fHistDCALS->Sumw2(); 
+  fHistDCALS->SetMinimum(0);
+  fOutput->Add(fHistDCALS);
+
+  return;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSEBkgLikeSignD0::UserExec(Option_t */*option*/)
+{
+  // Execute analysis for current event:
+  // heavy flavor candidates association to MC truth
+  
+  AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
+
+  // load heavy flavour vertices
+  TClonesArray *arrayVerticesHF =
+    (TClonesArray*)aod->GetList()->FindObject("VerticesHF");
+  if(!arrayVerticesHF) {
+    printf("AliAnalysisTaskSEBkgLikeSignD0::UserExec: VerticesHF branch not found!\n");
+    return;
+  }
+
+  // load D0->ee candidates                                                   
+  TClonesArray *arrayD0toKpi =
+    (TClonesArray*)aod->GetList()->FindObject("D0toKpi");
+  if(!arrayD0toKpi) {
+    printf("AliAnalysisTaskSEBkgLikeSignD0::UserExec: D0toKpi branch not found!\n");
+    return;
+  }
+
+  // load like sign candidates
+  TClonesArray *arrayLikeSign =
+    (TClonesArray*)aod->GetList()->FindObject("LikeSign2Prong");
+  if(!arrayLikeSign) {
+    printf("AliAnalysisTaskSEBkgLikeSignD0::UserExec: LikeSign2Prong branch not found!\n");
+    return;
+  }
+
+  // AOD primary vertex
+  AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
+
+  // make trkIDtoEntry register (temporary)
+  Int_t trkIDtoEntry[100000];
+  for(Int_t it=0;it<aod->GetNumberOfTracks();it++) {
+    AliAODTrack *track = aod->GetTrack(it);
+    trkIDtoEntry[track->GetID()]=it;
+  }
+
+  // loop over Like sign candidates
+  Int_t nPosPairs=0,nNegPairs=0;
+  Int_t nLikeSign = arrayLikeSign->GetEntriesFast();
+  printf("+++\n+++Number of like sign pairs ---> %d \n+++\n", nLikeSign);
+
+  for(Int_t iLikeSign = 0; iLikeSign < nLikeSign; iLikeSign++) {
+    AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)arrayLikeSign->UncheckedAt(iLikeSign);
+    Bool_t unsetvtx=kFALSE;
+    if(!d->GetOwnPrimaryVtx()) {
+        d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
+        unsetvtx=kTRUE;
+    }
+    Int_t okD0ls=0; Int_t okD0barls=0;
+    if(d->SelectD0(fVHF->GetD0toKpiCuts(),okD0ls,okD0barls)) {
+       if(okD0ls)fHistMassLS->Fill(d->InvMassD0(),0.5);
+       if(okD0barls)fHistMassLS->Fill(d->InvMassD0bar(),0.5);
+       fHistCPtaLS->Fill(d->CosPointingAngle());
+       fHistd0d0LS->Fill(1e8*d->Prodd0d0());
+       if(okD0ls)fHistCtsLS->Fill(d->CosThetaStarD0(),0.5);
+       if(okD0barls)fHistCtsLS->Fill(d->CosThetaStarD0bar(),0.5);
+       fHistDCALS->Fill(100*d->GetDCA());
+       PostData(1,fOutput);
+       AliAODTrack *trk0 = (AliAODTrack*)d->GetDaughter(0);
+       if(!trk0) {
+          trk0=aod->GetTrack(trkIDtoEntry[d->GetProngID(0)]);
+          printf("references to standard AOD not available \n");
+       }
+       if((trk0->Charge())==1) {
+          nPosPairs++;
+          fHistCtsLSpos->Fill(d->CosThetaStarD0());
+          PostData(1,fOutput);
+        } else {
+          nNegPairs++;
+          fHistCtsLSneg->Fill(d->CosThetaStarD0());
+          PostData(1,fOutput);
+        }
+      PostData(1,fOutput);
+      }
+    if(unsetvtx) d->UnsetOwnPrimaryVtx();
+  }
+
+  printf("------------ N. of positive pairs in Event ----- %d \n", nPosPairs);
+  printf("------------ N. of negative pairs in Event ----- %d \n", nNegPairs);
+
+  fTotPosPairs += nPosPairs;
+  fTotNegPairs += nNegPairs;
+
+  // loop over D0 candidates
+  Int_t nD0toKpi = arrayD0toKpi->GetEntriesFast();
+  printf("Number of like D0 -> Kpi candidates ---> %d \n", nD0toKpi);
+
+  for (Int_t iD0toKpi = 0; iD0toKpi < nD0toKpi; iD0toKpi++) {
+    AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)arrayD0toKpi->UncheckedAt(iD0toKpi);
+    Bool_t unsetvtx=kFALSE;
+    if(!d->GetOwnPrimaryVtx()) {
+      d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
+      unsetvtx=kTRUE;
+    }
+    Int_t okD0=0; Int_t okD0bar=0;
+    if(d->SelectD0(fVHF->GetD0toKpiCuts(),okD0,okD0bar)) {
+      if(okD0)fHistMassD0->Fill(d->InvMassD0(),0.5);
+      if(okD0bar)fHistMassD0->Fill(d->InvMassD0bar(),0.5);
+      if(okD0)fHistCtsD0->Fill(d->CosThetaStarD0(),0.5);
+      if(okD0bar)fHistCtsD0->Fill(d->CosThetaStarD0bar(),0.5);
+      fHistd0d0D0->Fill(1e8*d->Prodd0d0());
+      fHistCPtaD0->Fill(d->CosPointingAngle());
+      fHistDCAD0->Fill(100*d->GetDCA());
+      PostData(1,fOutput);
+    }
+   if(unsetvtx) d->UnsetOwnPrimaryVtx();
+  }
+
+  return;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskSEBkgLikeSignD0::Terminate(Option_t */*option*/)
+{
+  // Terminate analysis
+  //
+  if(fDebug > 1) printf("AnalysisTaskSEBkgLikeSignD0: Terminate() \n");
+
+  fOutput = dynamic_cast<TList*> (GetOutputData(1));
+  if (!fOutput) {     
+    printf("ERROR: fOutput not available\n");
+    return;
+  }
+
+  fLsNormalization = 2.*TMath::Sqrt(fTotPosPairs*fTotNegPairs); 
+
+  fHistMassD0 = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMassD0"));
+  fHistMassLS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMassLS"));
+  fHistCtsD0 = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCtsD0"));
+  fHistCtsLS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCtsLS"));
+  fHistCtsLSpos = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCtsLSpos"));
+  fHistCtsLSneg = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCtsLSneg"));
+  fHistCPtaD0 = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCPtaD0"));
+  fHistCPtaLS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistCPtaLS"));
+  fHistd0d0D0 = dynamic_cast<TH1F*>(fOutput->FindObject("fHistd0d0D0"));
+  fHistd0d0LS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistd0d0LS"));
+  fHistDCAD0 = dynamic_cast<TH1F*>(fOutput->FindObject("fHistDCAD0"));
+  fHistDCALS = dynamic_cast<TH1F*>(fOutput->FindObject("fHistDCALS"));
+
+  if(fLsNormalization>0) {
+    fHistMassLS->Scale((1/fLsNormalization)*fHistMassLS->GetEntries());
+    fHistCtsLS->Scale((1/fLsNormalization)*fHistCtsLS->GetEntries());
+    fHistCtsLSpos->Scale((1/fLsNormalization)*fHistCtsLSpos->GetEntries());
+    fHistCtsLSneg->Scale((1/fLsNormalization)*fHistCtsLSneg->GetEntries());
+    fHistCPtaLS->Scale((1/fLsNormalization)*fHistCPtaLS->GetEntries());
+    fHistd0d0LS->Scale((1/fLsNormalization)*fHistd0d0LS->GetEntries());
+    fHistDCALS->Scale((1/fLsNormalization)*fHistDCALS->GetEntries());
+  }
+
+  return;
+}
diff --git a/PWG3/vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.h b/PWG3/vertexingHF/AliAnalysisTaskSEBkgLikeSignD0.h
new file mode 100644 (file)
index 0000000..ecc8621
--- /dev/null
@@ -0,0 +1,66 @@
+#ifndef AliAnalysisTaskSEBkgLikeSignD0_H
+#define AliAnalysisTaskSEBkgLikeSignD0_H
+
+/* Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+//*************************************************************************
+// Class AliAnalysisTaskSEBkgLikeSignD0
+// AliAnalysisTaskSE for reading both reconstructed D0 -> Kpi candidates
+// and like sign pairs and for drawing corresponding distributions
+// Author: C.Di Giglio, carmelo.digiglio@ba.infn.it
+//*************************************************************************
+
+#include <TROOT.h>
+#include <TSystem.h>
+#include <TList.h>
+#include <TH1F.h>
+#include <TH2F.h>
+
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisVertexingHF.h"
+
+class AliAnalysisTaskSEBkgLikeSignD0 : public AliAnalysisTaskSE
+{
+ public:
+
+  AliAnalysisTaskSEBkgLikeSignD0();
+  AliAnalysisTaskSEBkgLikeSignD0(const char *name);
+  virtual ~AliAnalysisTaskSEBkgLikeSignD0();
+
+
+  // Implementation of interface methods
+  virtual void UserCreateOutputObjects();
+  virtual void Init();
+  virtual void LocalInit() {Init();}
+  virtual void UserExec(Option_t *option);
+  virtual void Terminate(Option_t *option);
+
+ private:
+
+  AliAnalysisTaskSEBkgLikeSignD0(const AliAnalysisTaskSEBkgLikeSignD0 &source);
+  AliAnalysisTaskSEBkgLikeSignD0& operator=(const AliAnalysisTaskSEBkgLikeSignD0& source); 
+
+  TList   *fOutput;                //! list send on output slot 0
+  TH1F    *fHistMassD0;            //! output histograms
+  TH1F    *fHistMassLS;            //!
+  TH1F    *fHistCtsD0;             //! Cosine of decay angle
+  TH1F    *fHistCtsLS;             //!
+  TH1F    *fHistCtsLSpos;          //!
+  TH1F    *fHistCtsLSneg;          //!
+  TH1F    *fHistCPtaD0;            //! Cosine of pointing angle
+  TH1F    *fHistCPtaLS;            //!
+  TH1F    *fHistd0d0D0;            //! Product of impact parameters
+  TH1F    *fHistd0d0LS;            //!
+  TH1F    *fHistDCAD0;             //! Distance of closest approach
+  TH1F    *fHistDCALS;             //! like-sign
+  AliAnalysisVertexingHF *fVHF;    // Vertexer heavy flavour (used to pass the cuts)
+
+  Int_t fTotPosPairs;              //
+  Int_t fTotNegPairs;              // normalization
+  Double_t fLsNormalization;       //
+  ClassDef(AliAnalysisTaskSEBkgLikeSignD0,1); // comparison of unlike-sign and like-sign background for D0->Kpi
+};
+
+#endif
index 1fdc9a5..2514747 100644 (file)
@@ -166,9 +166,13 @@ void RunAnalysisAODVertexingHF()
   //gROOT->LoadMacro(taskName.Data());\r
   //AliAnalysisTaskSESelectHF *seleTask = AddTaskSelectHF();\r
 \r
+  taskName="AddTaskBkgLikeSignD0.C"; taskName.Prepend(loadMacroPath.Data());\r
+  gROOT->LoadMacro(taskName.Data());\r
+  AliAnalysisTaskSEBkgLikeSignD0 *lsD0Task = AddTaskBkgLikeSignD0();\r
+\r
   taskName="AddTaskBkgLikeSignJPSI.C"; taskName.Prepend(loadMacroPath.Data());\r
   gROOT->LoadMacro(taskName.Data());\r
-  AliAnalysisTaskSEBkgLikeSignJPSI *lsTask = AddTaskBkgLikeSignJPSI();\r
+  AliAnalysisTaskSEBkgLikeSignJPSI *lsJPSITask = AddTaskBkgLikeSignJPSI();\r
 \r
   //taskName="AddTaskBtoJPSItoEle.C"; taskName.Prepend(loadMacroPath.Data());\r
   //gROOT->LoadMacro(taskName.Data());\r
@@ -228,6 +232,8 @@ AliAnalysisGrid* CreateAlienHandler(TString pluginmode="test",Bool_t useParFiles
    plugin->SetFriendChainName("AliAOD.VertexingHF.root");\r
    // ...then add run numbers to be considered\r
   plugin->AddRunNumber(529007);\r
+  //  or\r
+  //plugin->SetRunRange(529000,529007);\r
    // Method 2: Declare existing data files (raw collections, xml collections, root file)\r
    // If no path mentioned data is supposed to be in the work directory (see SetGridWorkingDir())\r
    // XML collections added via this method can be combined with the first method if\r