update of D0-HFE correlation
authorarossi <arossi@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 23 Jan 2013 09:05:55 +0000 (09:05 +0000)
committerarossi <arossi@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 23 Jan 2013 09:05:55 +0000 (09:05 +0000)
- using asymmetric cut on electron PID to reduce contamination of
electron selection by hadrons
- workaround for volatile external pointers in AliHFAssociatedTrackCuts
- correction of info and error messages
- removing carriage returns, accidentally added in r60367 on every line

PWGHF/correlationHF/AliAnalysisTaskDxHFECorrelation.cxx
PWGHF/correlationHF/AliDxHFECorrelation.cxx
PWGHF/correlationHF/macros/AddTaskDxHFECorrelation.C

index 2de1b3e..670551f 100644 (file)
-// $Id$\r
-\r
-//**************************************************************************\r
-//* This file is property of and copyright by the ALICE Project            * \r
-//* ALICE Experiment at CERN, All rights reserved.                         *\r
-//*                                                                        *\r
-//* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *\r
-//*                  Sedat Altinpinar <Sedat.Altinpinar@cern.ch>           *\r
-//*                  Hege Erdal       <hege.erdal@gmail.com>               *\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
-/// @file   AliAnalysisTaskDxHFECorrelation.cxx\r
-/// @author Sedat Altinpinar, Hege Erdal, Matthias Richter\r
-/// @date   2012-03-19\r
-/// @brief  AnalysisTask D0 - HFE correlation\r
-///\r
-\r
-#include "AliAnalysisTaskDxHFECorrelation.h"\r
-#include "AliDxHFECorrelation.h"\r
-#include "AliDxHFECorrelationMC.h"\r
-#include "AliDxHFEParticleSelectionD0.h"\r
-#include "AliDxHFEParticleSelectionMCD0.h"\r
-#include "AliDxHFEParticleSelectionEl.h"\r
-#include "AliDxHFEParticleSelectionMCEl.h"\r
-#include "AliDxHFEParticleSelection.h"\r
-#include "AliHFCorrelator.h"\r
-#include "AliAnalysisManager.h"\r
-#include "AliLog.h"\r
-#include "AliESDInputHandler.h"\r
-#include "AliAnalysisDataSlot.h"\r
-#include "AliAnalysisDataContainer.h"\r
-#include "AliAnalysisManager.h"\r
-#include "AliVertexerTracks.h"\r
-#include "AliAODHandler.h"\r
-#include "AliInputEventHandler.h"\r
-#include "AliAODEvent.h"\r
-#include "AliAODVertex.h"\r
-#include "AliAODTrack.h"\r
-#include "AliAODMCHeader.h"\r
-#include "AliAODMCParticle.h"\r
-#include "AliAODRecoDecayHF2Prong.h"\r
-#include "AliAODRecoCascadeHF.h"\r
-#include "AliRDHFCutsD0toKpi.h"\r
-#include "AliPID.h"\r
-#include "AliPIDResponse.h"\r
-#include "AliHFEcontainer.h"\r
-#include "AliHFEpid.h"\r
-#include "AliHFEpidBase.h"\r
-#include "AliHFEcuts.h"\r
-#include "AliHFEtools.h"\r
-#include "TObject.h"\r
-#include "TChain.h"\r
-#include "TSystem.h"\r
-#include "AliReducedParticle.h"\r
-#include "AliHFAssociatedTrackCuts.h" // initialization of event pool\r
-#include "TFile.h"\r
-#include <memory>\r
-\r
-using namespace std;\r
-\r
-/// ROOT macro for the implementation of ROOT specific class methods\r
-ClassImp(AliAnalysisTaskDxHFECorrelation)\r
-\r
-AliAnalysisTaskDxHFECorrelation::AliAnalysisTaskDxHFECorrelation(const char* opt)\r
-  : AliAnalysisTaskSE("AliAnalysisTaskDxHFECorrelation")\r
-  , fOutput(0)\r
-  , fOption(opt)\r
-  , fCorrelation(NULL)\r
-  , fD0s(NULL)\r
-  , fElectrons(NULL)\r
-  , fCutsD0(NULL)\r
-  , fCutsHFE(NULL)\r
-  , fCuts(NULL)\r
-  , fPID(NULL)\r
-  , fFillOnlyD0D0bar(0)\r
-  , fUseMC(kFALSE)\r
-  , fUseEventMixing(kFALSE)\r
-  , fSystem(0)\r
-  , fSelectedD0s(NULL)\r
-  , fSelectedElectrons(NULL)\r
-{\r
-  // constructor\r
-  //\r
-  //\r
-  DefineSlots();\r
-  fPID = new AliHFEpid("hfePid");\r
-\r
-}\r
-\r
-int AliAnalysisTaskDxHFECorrelation::DefineSlots()\r
-{\r
-  // define the data slots\r
-  DefineInput(0, TChain::Class());\r
-  DefineOutput(1, TList::Class());\r
-  DefineOutput(2,AliRDHFCutsD0toKpi::Class());\r
-  DefineOutput(3,AliHFEcuts::Class());\r
-  DefineOutput(4,AliHFAssociatedTrackCuts::Class());\r
-  return 0;\r
-}\r
-\r
-AliAnalysisTaskDxHFECorrelation::~AliAnalysisTaskDxHFECorrelation()\r
-{\r
-  // destructor\r
-  //\r
-  //\r
-\r
-  // histograms are in the output list and deleted when the output\r
-  // list is deleted by the TSelector dtor\r
-\r
-  if (fOutput && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {\r
-    delete fOutput;\r
-    fOutput = 0;\r
-  }\r
-  if (fD0s) delete fD0s;\r
-  fD0s=NULL;\r
-  if (fElectrons) delete fElectrons;\r
-  fElectrons=NULL;\r
-  if (fCorrelation) delete fCorrelation;\r
-  fCorrelation=NULL;\r
-  // external object, do not delete\r
-  fCutsD0=NULL;\r
-  // external object, do not delete\r
-  fCutsHFE=NULL;\r
-  if(fPID) delete fPID;\r
-  fPID=NULL;\r
-  if(fSelectedElectrons) delete fSelectedElectrons;\r
-  fSelectedElectrons=NULL;\r
-  if(fSelectedD0s) delete fSelectedD0s;\r
-  fSelectedD0s=NULL;\r
-\r
-\r
-}\r
-\r
-void AliAnalysisTaskDxHFECorrelation::UserCreateOutputObjects()\r
-{\r
-  // create result objects and add to output list\r
-  int iResult=0;\r
-\r
-  //Initialize PID for electron selection\r
-  if(!fPID->GetNumberOfPIDdetectors()) { \r
-    fPID->AddDetector("TOF",0);\r
-    fPID->AddDetector("TPC",1);\r
-  }\r
-  fPID->InitializePID();\r
-\r
-  fOutput = new TList;\r
-  fOutput->SetOwner();\r
-\r
-  // setting up for D0s\r
-  TString selectionD0Options;\r
-  switch (fFillOnlyD0D0bar) {\r
-  case 1: selectionD0Options+="FillOnlyD0 "; break;\r
-  case 2: selectionD0Options+="FillOnlyD0bar "; break;\r
-  default: selectionD0Options+="FillD0D0bar ";\r
-  }\r
-\r
-  if(fUseMC) fD0s=new AliDxHFEParticleSelectionMCD0(selectionD0Options);\r
-  else fD0s=new AliDxHFEParticleSelectionD0(selectionD0Options);\r
-  fD0s->SetCuts(fCutsD0);\r
-  iResult=fD0s->Init();\r
-  if (iResult<0) {\r
-    AliFatal(Form("initialization of worker class instance fD0s failed with error %d", iResult));\r
-  }\r
-\r
-  //Electrons\r
-  if(fUseMC) fElectrons=new AliDxHFEParticleSelectionMCEl;\r
-  else fElectrons=new AliDxHFEParticleSelectionEl;\r
-  fElectrons->SetCuts(fPID, AliDxHFEParticleSelectionEl::kCutPID);\r
-  fElectrons->SetCuts(fCutsHFE, AliDxHFEParticleSelectionEl::kCutHFE);\r
-  iResult=fElectrons->Init();\r
-  if (iResult<0) {\r
-    AliFatal(Form("initialization of worker class instance fElectrons failed with error %d", iResult));\r
-  }\r
-\r
-  //Correlation\r
-  if(fUseMC) fCorrelation=new AliDxHFECorrelationMC;\r
-  else fCorrelation=new AliDxHFECorrelation;\r
-  fCorrelation->SetCuts(fCuts);\r
-  // TODO: check if we can get rid of the mc flag in the correlation analysis class\r
-  // at the moment this is needed to pass on info to AliHFCorrelator\r
-  TString arguments;\r
-  if (fUseMC)          arguments+=" use-mc";\r
-  if (fUseEventMixing) arguments+=" event-mixing";\r
-  // TODO: fSystem is a boolean right now, needs to be changed to fit also p-Pb\r
-  if (!fSystem)         arguments+=" system=pp";\r
-  else                 arguments+=" system=Pb-Pb";\r
-  iResult=fCorrelation->Init(arguments);\r
-  if (iResult<0) {\r
-    AliFatal(Form("initialization of worker class instance fCorrelation failed with error %d", iResult));\r
-  }\r
-\r
-  // Fix for merging:\r
-  // Retrieving the individual objects created\r
-  // and storing them instead of fD0s, fElectrons etc.. \r
-  TList *list =(TList*)fCorrelation->GetControlObjects();\r
-  TObject *obj=NULL;\r
-\r
-  TIter next(list);\r
-  while((obj = next())){\r
-    fOutput->Add(obj);\r
-  }\r
-\r
-  list=(TList*)fD0s->GetControlObjects();\r
-  next=TIter(list);\r
-  while((obj= next())){\r
-    fOutput->Add(obj);\r
-  }\r
-\r
-  list=(TList*)fElectrons->GetControlObjects();\r
-  next=TIter(list);\r
-  while((obj = next()))\r
-    fOutput->Add(obj);\r
-\r
-  if (fCutsD0) {\r
-    // TODO: eliminate this copy\r
-    AliRDHFCutsD0toKpi* cuts=dynamic_cast<AliRDHFCutsD0toKpi*>(fCutsD0);\r
-    if (!cuts) {\r
-      AliFatal(Form("cut object %s is of incorrect type %s, expecting AliRDHFCutsD0toKpi", fCutsD0->GetName(), fCutsD0->ClassName()));\r
-      return;\r
-    }\r
-  }\r
-  // TODO: why copy? cleanup?\r
-  AliRDHFCutsD0toKpi* copyfCuts=new AliRDHFCutsD0toKpi(dynamic_cast<AliRDHFCutsD0toKpi&>(*fCutsD0));\r
-  const char* nameoutput=GetOutputSlot(2)->GetContainer()->GetName();\r
-  copyfCuts->SetName(nameoutput);\r
-\r
-  // all tasks must post data once for all outputs\r
-  PostData(1, fOutput);\r
-  PostData(2,copyfCuts);\r
-  PostData(3,fCutsHFE);\r
-  PostData(4,fCuts);\r
-\r
-}\r
-\r
-void AliAnalysisTaskDxHFECorrelation::UserExec(Option_t* /*option*/)\r
-{\r
-  // process the event\r
-  TObject* pInput=InputEvent();\r
-  if (!pInput) {\r
-    AliError("failed to get input");\r
-    return;\r
-  }\r
-  AliVEvent *pEvent = dynamic_cast<AliVEvent*>(pInput);\r
-  TClonesArray *inputArray=0;\r
-\r
-  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsAll);\r
-\r
-  if(!pEvent && AODEvent() && IsStandardAOD()) { //Not sure if this is needed.. Keep it for now. \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
-    pEvent = 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
-    \r
-    if(aodHandler->GetExtensions()) {\r
-      AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");\r
-      AliAODEvent* aodFromExt = ext->GetAOD();\r
-      inputArray=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");\r
-    }\r
-  } else if(pEvent) {\r
-    inputArray=(TClonesArray*)pEvent->GetList()->FindObject("D0toKpi");\r
-  }\r
-  if(!inputArray || !pEvent) {\r
-    AliError("Input branch not found!\n");\r
-    return;\r
-  }\r
-  // fix for temporary bug in ESDfilter\r
-  // the AODs with null vertex pointer didn't pass the PhysSel\r
-  if(!pEvent->GetPrimaryVertex() || TMath::Abs(pEvent->GetMagneticField())<0.001){\r
-    AliDebug(2,"Rejected at GetPrimaryvertex");\r
-    return;\r
-  }\r
-\r
-  AliRDHFCuts* cutsd0=dynamic_cast<AliRDHFCuts*>(fCutsD0);\r
-  if (!cutsd0) return; // Fatal thrown already in initialization\r
-\r
-  if(!cutsd0->IsEventSelected(pEvent)) {\r
-    AliDebug(2,"rejected at IsEventSelected");\r
-    return;\r
-  }\r
-\r
-  if(!fPID->IsInitialized()){ \r
-    // Initialize PID with the given run number\r
-    AliWarning("PID not initialised, get from Run no");\r
-    fPID->InitializePID(pEvent->GetRunNumber());\r
-  }\r
-\r
-  AliPIDResponse *pidResponse = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler())->GetPIDResponse();\r
-  if(!pidResponse){\r
-    AliDebug(1, "Using default PID Response");\r
-    pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class()); \r
-  }\r
\r
-  fPID->SetPIDResponse(pidResponse);\r
-\r
-  // Retrieving process from the AODMCHeader. \r
-  // TODO: Move it somewhere else? (keep it here for the moment since only need to read once pr event)\r
-  if(fUseMC){\r
-    AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(pEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));\r
-    \r
-    if (!mcHeader) {\r
-      AliError("Could not find MC Header in AOD");\r
-      return;\r
-    }\r
-    Int_t eventType = mcHeader->GetEventType();\r
-    fCorrelation->SetEventType(eventType);\r
-  }\r
-\r
-  Int_t nInD0toKpi = inputArray->GetEntriesFast();\r
-\r
-  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsSel);\r
-\r
-  if(fSelectedD0s) delete fSelectedD0s;\r
-  fSelectedD0s=(fD0s->Select(inputArray,pEvent));\r
-  \r
-  if(! fSelectedD0s) {\r
-    return;\r
-  }\r
-  Int_t nD0Selected = fSelectedD0s->GetEntriesFast();\r
-\r
-\r
-  /*std::auto_ptr<TObjArray> pSelectedD0s(fD0s->Select(inputArray,pEvent));\r
-  if(! pSelectedD0s.get()) {\r
-    return;\r
-  }\r
-  Int_t nD0Selected = pSelectedD0s->GetEntriesFast();*/\r
-\r
-  // When not using EventMixing, no need to go further if no D0s are found.\r
-  // For Event Mixing, need to store all found electrons in the pool\r
-  if(!fUseEventMixing && nD0Selected==0){\r
-    AliDebug(4,"No D0s found in this event");\r
-    return;\r
-  }\r
-\r
-  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsD0);\r
-\r
-  /* std::auto_ptr<TObjArray> pSelectedElectrons(fElectrons->Select(pEvent));\r
-  // note: the pointer is deleted automatically once the scope is left\r
-  // if the array should be published, the auto pointer must be released\r
-  // first, however some other cleanup will be necessary in that case\r
-  // probably a clone with a reduced AliVParticle implementation is\r
-  // appropriate.\r
-\r
-  if(! pSelectedElectrons.get()) {\r
-    return;\r
-  }\r
-\r
-  Int_t nElSelected =  pSelectedElectrons->GetEntriesFast();*/\r
-  if (fSelectedElectrons) delete fSelectedElectrons;\r
-  fSelectedElectrons=(fElectrons->Select(pEvent));\r
-  \r
-  if(! fSelectedElectrons) {\r
-    return;\r
-  }\r
-\r
-  Int_t nElSelected =  fSelectedElectrons->GetEntriesFast();\r
-\r
-\r
-  // No need to go further if no electrons are found, except for event mixing. Will here anyway correlate D0s with electrons from previous events\r
-  if(!fUseEventMixing && nElSelected==0){\r
-    AliDebug(4,"No electrons found in this event");\r
-    return;\r
-  }\r
-  if(nD0Selected==0 && nElSelected==0){\r
-    AliDebug(4,"Neither D0 nor electrons in this event");\r
-    return;\r
-  }\r
-  \r
-  AliDebug(4,Form("Number of D0->Kpi Start: %d , End: %d    Electrons Selected: %d\n", nInD0toKpi, nD0Selected, nElSelected));\r
-\r
-  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsD0e);\r
-\r
-  //int iResult=fCorrelation->Fill(pSelectedD0s.get(), pSelectedElectrons.get(), pEvent);\r
-  int iResult=fCorrelation->Fill(fSelectedD0s, fSelectedElectrons, pEvent);\r
-\r
-  if (iResult<0) {\r
-    AliFatal(Form("%s processing failed with error %d", fCorrelation->GetName(), iResult));\r
-  }\r
-\r
-  PostData(1, fOutput);\r
-  return;\r
-\r
-}\r
-\r
-void AliAnalysisTaskDxHFECorrelation::FinishTaskOutput()\r
-{\r
-  // end of the processing\r
-}\r
-\r
-void AliAnalysisTaskDxHFECorrelation::Terminate(Option_t *)\r
-{\r
-  // last action on the client\r
-  fOutput = dynamic_cast<TList*> (GetOutputData(1));\r
-  if (!fOutput) {\r
-    // looks like that is a valid condition if the task is run\r
-    // in mode "terminate"\r
-    return;\r
-  }\r
-}\r
+// $Id$
+
+//**************************************************************************
+//* This file is property of and copyright by the ALICE Project            * 
+//* ALICE Experiment at CERN, All rights reserved.                         *
+//*                                                                        *
+//* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
+//*                  Sedat Altinpinar <Sedat.Altinpinar@cern.ch>           *
+//*                  Hege Erdal       <hege.erdal@gmail.com>               *
+//*                                                                        *
+//* 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.                  *
+//**************************************************************************
+
+/// @file   AliAnalysisTaskDxHFECorrelation.cxx
+/// @author Sedat Altinpinar, Hege Erdal, Matthias Richter
+/// @date   2012-03-19
+/// @brief  AnalysisTask D0 - HFE correlation
+///
+
+#include "AliAnalysisTaskDxHFECorrelation.h"
+#include "AliDxHFECorrelation.h"
+#include "AliDxHFECorrelationMC.h"
+#include "AliDxHFEParticleSelectionD0.h"
+#include "AliDxHFEParticleSelectionMCD0.h"
+#include "AliDxHFEParticleSelectionEl.h"
+#include "AliDxHFEParticleSelectionMCEl.h"
+#include "AliDxHFEParticleSelection.h"
+#include "AliHFCorrelator.h"
+#include "AliAnalysisManager.h"
+#include "AliLog.h"
+#include "AliESDInputHandler.h"
+#include "AliAnalysisDataSlot.h"
+#include "AliAnalysisDataContainer.h"
+#include "AliAnalysisManager.h"
+#include "AliVertexerTracks.h"
+#include "AliAODHandler.h"
+#include "AliInputEventHandler.h"
+#include "AliAODEvent.h"
+#include "AliAODVertex.h"
+#include "AliAODTrack.h"
+#include "AliAODMCHeader.h"
+#include "AliAODMCParticle.h"
+#include "AliAODRecoDecayHF2Prong.h"
+#include "AliAODRecoCascadeHF.h"
+#include "AliRDHFCutsD0toKpi.h"
+#include "AliPID.h"
+#include "AliPIDResponse.h"
+#include "AliHFEcontainer.h"
+#include "AliHFEpid.h"
+#include "AliHFEpidBase.h"
+#include "AliHFEcuts.h"
+#include "AliHFEtools.h"
+#include "TObject.h"
+#include "TChain.h"
+#include "TSystem.h"
+#include "AliReducedParticle.h"
+#include "AliHFAssociatedTrackCuts.h" // initialization of event pool
+#include "TFile.h"
+#include <memory>
+
+using namespace std;
+
+/// ROOT macro for the implementation of ROOT specific class methods
+ClassImp(AliAnalysisTaskDxHFECorrelation)
+
+AliAnalysisTaskDxHFECorrelation::AliAnalysisTaskDxHFECorrelation(const char* opt)
+  : AliAnalysisTaskSE("AliAnalysisTaskDxHFECorrelation")
+  , fOutput(0)
+  , fOption(opt)
+  , fCorrelation(NULL)
+  , fD0s(NULL)
+  , fElectrons(NULL)
+  , fCutsD0(NULL)
+  , fCutsHFE(NULL)
+  , fCuts(NULL)
+  , fPID(NULL)
+  , fFillOnlyD0D0bar(0)
+  , fUseMC(kFALSE)
+  , fUseEventMixing(kFALSE)
+  , fSystem(0)
+  , fSelectedD0s(NULL)
+  , fSelectedElectrons(NULL)
+{
+  // constructor
+  //
+  //
+  DefineSlots();
+  fPID = new AliHFEpid("hfePid");
+
+}
+
+int AliAnalysisTaskDxHFECorrelation::DefineSlots()
+{
+  // define the data slots
+  DefineInput(0, TChain::Class());
+  DefineOutput(1, TList::Class());
+  DefineOutput(2,AliRDHFCutsD0toKpi::Class());
+  DefineOutput(3,AliHFEcuts::Class());
+  DefineOutput(4,AliHFAssociatedTrackCuts::Class());
+  return 0;
+}
+
+AliAnalysisTaskDxHFECorrelation::~AliAnalysisTaskDxHFECorrelation()
+{
+  // destructor
+  //
+  //
+
+  // histograms are in the output list and deleted when the output
+  // list is deleted by the TSelector dtor
+
+  if (fOutput && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
+    delete fOutput;
+    fOutput = 0;
+  }
+  if (fD0s) delete fD0s;
+  fD0s=NULL;
+  if (fElectrons) delete fElectrons;
+  fElectrons=NULL;
+  if (fCorrelation) delete fCorrelation;
+  fCorrelation=NULL;
+  // external object, do not delete
+  fCutsD0=NULL;
+  // external object, do not delete
+  fCutsHFE=NULL;
+  if(fPID) delete fPID;
+  fPID=NULL;
+  if(fSelectedElectrons) delete fSelectedElectrons;
+  fSelectedElectrons=NULL;
+  if(fSelectedD0s) delete fSelectedD0s;
+  fSelectedD0s=NULL;
+
+
+}
+
+void AliAnalysisTaskDxHFECorrelation::UserCreateOutputObjects()
+{
+  // create result objects and add to output list
+  int iResult=0;
+
+  //Initialize PID for electron selection
+  if(!fPID->GetNumberOfPIDdetectors()) { 
+    fPID->AddDetector("TOF",0);
+    fPID->AddDetector("TPC",1);
+  }
+  fPID->InitializePID();
+  const int paramSize=4;
+  Double_t params[paramSize];
+  memset(params, 0, sizeof(Double_t)*paramSize);
+  params[0]=1.;
+  fPID->ConfigureTPCdefaultCut(NULL, params, 3.);
+
+  fOutput = new TList;
+  fOutput->SetOwner();
+
+  // setting up for D0s
+  TString selectionD0Options;
+  switch (fFillOnlyD0D0bar) {
+  case 1: selectionD0Options+="FillOnlyD0 "; break;
+  case 2: selectionD0Options+="FillOnlyD0bar "; break;
+  default: selectionD0Options+="FillD0D0bar ";
+  }
+
+  if(fUseMC) fD0s=new AliDxHFEParticleSelectionMCD0(selectionD0Options);
+  else fD0s=new AliDxHFEParticleSelectionD0(selectionD0Options);
+  fD0s->SetCuts(fCutsD0);
+  iResult=fD0s->Init();
+  if (iResult<0) {
+    AliFatal(Form("initialization of worker class instance fD0s failed with error %d", iResult));
+  }
+
+  //Electrons
+  if(fUseMC) fElectrons=new AliDxHFEParticleSelectionMCEl;
+  else fElectrons=new AliDxHFEParticleSelectionEl;
+  fElectrons->SetCuts(fPID, AliDxHFEParticleSelectionEl::kCutPID);
+  fElectrons->SetCuts(fCutsHFE, AliDxHFEParticleSelectionEl::kCutHFE);
+  iResult=fElectrons->Init();
+  if (iResult<0) {
+    AliFatal(Form("initialization of worker class instance fElectrons failed with error %d", iResult));
+  }
+
+  //Correlation
+  if(fUseMC) fCorrelation=new AliDxHFECorrelationMC;
+  else fCorrelation=new AliDxHFECorrelation;
+  fCorrelation->SetCuts(fCuts);
+  // TODO: check if we can get rid of the mc flag in the correlation analysis class
+  // at the moment this is needed to pass on info to AliHFCorrelator
+  TString arguments;
+  if (fUseMC)          arguments+=" use-mc";
+  if (fUseEventMixing) arguments+=" event-mixing";
+  // TODO: fSystem is a boolean right now, needs to be changed to fit also p-Pb
+  if (!fSystem)         arguments+=" system=pp";
+  else                 arguments+=" system=Pb-Pb";
+  iResult=fCorrelation->Init(arguments);
+  if (iResult<0) {
+    AliFatal(Form("initialization of worker class instance fCorrelation failed with error %d", iResult));
+  }
+
+  // Fix for merging:
+  // Retrieving the individual objects created
+  // and storing them instead of fD0s, fElectrons etc.. 
+  TList *list =(TList*)fCorrelation->GetControlObjects();
+  TObject *obj=NULL;
+
+  TIter next(list);
+  while((obj = next())){
+    fOutput->Add(obj);
+  }
+
+  list=(TList*)fD0s->GetControlObjects();
+  next=TIter(list);
+  while((obj= next())){
+    fOutput->Add(obj);
+  }
+
+  list=(TList*)fElectrons->GetControlObjects();
+  next=TIter(list);
+  while((obj = next()))
+    fOutput->Add(obj);
+
+  if (fCutsD0) {
+    // TODO: eliminate this copy
+    AliRDHFCutsD0toKpi* cuts=dynamic_cast<AliRDHFCutsD0toKpi*>(fCutsD0);
+    if (!cuts) {
+      AliFatal(Form("cut object %s is of incorrect type %s, expecting AliRDHFCutsD0toKpi", fCutsD0->GetName(), fCutsD0->ClassName()));
+      return;
+    }
+  }
+  // TODO: why copy? cleanup?
+  AliRDHFCutsD0toKpi* copyfCuts=new AliRDHFCutsD0toKpi(dynamic_cast<AliRDHFCutsD0toKpi&>(*fCutsD0));
+  const char* nameoutput=GetOutputSlot(2)->GetContainer()->GetName();
+  copyfCuts->SetName(nameoutput);
+
+  // all tasks must post data once for all outputs
+  PostData(1, fOutput);
+  PostData(2,copyfCuts);
+  PostData(3,fCutsHFE);
+  PostData(4,fCuts);
+
+}
+
+void AliAnalysisTaskDxHFECorrelation::UserExec(Option_t* /*option*/)
+{
+  // process the event
+  TObject* pInput=InputEvent();
+  if (!pInput) {
+    AliError("failed to get input");
+    return;
+  }
+  AliVEvent *pEvent = dynamic_cast<AliVEvent*>(pInput);
+  TClonesArray *inputArray=0;
+
+  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsAll);
+
+  if(!pEvent && AODEvent() && IsStandardAOD()) { //Not sure if this is needed.. Keep it for now. 
+    // In case there is an AOD handler writing a standard AOD, use the AOD 
+    // event in memory rather than the input (ESD) event.    
+    pEvent = dynamic_cast<AliAODEvent*> (AODEvent());
+    // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
+    // have to taken from the AOD event hold by the AliAODExtension
+    AliAODHandler* aodHandler = (AliAODHandler*) 
+      ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
+    
+    if(aodHandler->GetExtensions()) {
+      AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
+      AliAODEvent* aodFromExt = ext->GetAOD();
+      inputArray=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
+    }
+  } else if(pEvent) {
+    inputArray=(TClonesArray*)pEvent->GetList()->FindObject("D0toKpi");
+  }
+  if(!inputArray || !pEvent) {
+    AliError("Input branch not found!\n");
+    return;
+  }
+  // fix for temporary bug in ESDfilter
+  // the AODs with null vertex pointer didn't pass the PhysSel
+  if(!pEvent->GetPrimaryVertex() || TMath::Abs(pEvent->GetMagneticField())<0.001){
+    AliDebug(2,"Rejected at GetPrimaryvertex");
+    return;
+  }
+
+  AliRDHFCuts* cutsd0=dynamic_cast<AliRDHFCuts*>(fCutsD0);
+  if (!cutsd0) return; // Fatal thrown already in initialization
+
+  if(!cutsd0->IsEventSelected(pEvent)) {
+    AliDebug(2,"rejected at IsEventSelected");
+    return;
+  }
+
+  if(!fPID->IsInitialized()){ 
+    // Initialize PID with the given run number
+    AliWarning("PID not initialised, get from Run no");
+    fPID->InitializePID(pEvent->GetRunNumber());
+  }
+
+  AliPIDResponse *pidResponse = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler())->GetPIDResponse();
+  if(!pidResponse){
+    AliDebug(1, "Using default PID Response");
+    pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class()); 
+  }
+  fPID->SetPIDResponse(pidResponse);
+
+  // Retrieving process from the AODMCHeader. 
+  // TODO: Move it somewhere else? (keep it here for the moment since only need to read once pr event)
+  if(fUseMC){
+    AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(pEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
+    
+    if (!mcHeader) {
+      AliError("Could not find MC Header in AOD");
+      return;
+    }
+    Int_t eventType = mcHeader->GetEventType();
+    fCorrelation->SetEventType(eventType);
+  }
+
+  Int_t nInD0toKpi = inputArray->GetEntriesFast();
+
+  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsSel);
+
+  if(fSelectedD0s) delete fSelectedD0s;
+  fSelectedD0s=(fD0s->Select(inputArray,pEvent));
+  
+  if(! fSelectedD0s) {
+    return;
+  }
+  Int_t nD0Selected = fSelectedD0s->GetEntriesFast();
+
+
+  /*std::auto_ptr<TObjArray> pSelectedD0s(fD0s->Select(inputArray,pEvent));
+  if(! pSelectedD0s.get()) {
+    return;
+  }
+  Int_t nD0Selected = pSelectedD0s->GetEntriesFast();*/
+
+  // When not using EventMixing, no need to go further if no D0s are found.
+  // For Event Mixing, need to store all found electrons in the pool
+  if(!fUseEventMixing && nD0Selected==0){
+    AliDebug(4,"No D0s found in this event");
+    return;
+  }
+
+  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsD0);
+
+  /* std::auto_ptr<TObjArray> pSelectedElectrons(fElectrons->Select(pEvent));
+  // note: the pointer is deleted automatically once the scope is left
+  // if the array should be published, the auto pointer must be released
+  // first, however some other cleanup will be necessary in that case
+  // probably a clone with a reduced AliVParticle implementation is
+  // appropriate.
+
+  if(! pSelectedElectrons.get()) {
+    return;
+  }
+
+  Int_t nElSelected =  pSelectedElectrons->GetEntriesFast();*/
+  if (fSelectedElectrons) delete fSelectedElectrons;
+  fSelectedElectrons=(fElectrons->Select(pEvent));
+  
+  if(! fSelectedElectrons) {
+    return;
+  }
+
+  Int_t nElSelected =  fSelectedElectrons->GetEntriesFast();
+
+
+  // No need to go further if no electrons are found, except for event mixing. Will here anyway correlate D0s with electrons from previous events
+  if(!fUseEventMixing && nElSelected==0){
+    AliDebug(4,"No electrons found in this event");
+    return;
+  }
+  if(nD0Selected==0 && nElSelected==0){
+    AliDebug(4,"Neither D0 nor electrons in this event");
+    return;
+  }
+  
+  AliDebug(4,Form("Number of D0->Kpi Start: %d , End: %d    Electrons Selected: %d\n", nInD0toKpi, nD0Selected, nElSelected));
+
+  fCorrelation->HistogramEventProperties(AliDxHFECorrelation::kEventsD0e);
+
+  //int iResult=fCorrelation->Fill(pSelectedD0s.get(), pSelectedElectrons.get(), pEvent);
+  int iResult=fCorrelation->Fill(fSelectedD0s, fSelectedElectrons, pEvent);
+
+  if (iResult<0) {
+    AliFatal(Form("%s processing failed with error %d", fCorrelation->GetName(), iResult));
+  }
+
+  PostData(1, fOutput);
+  return;
+
+}
+
+void AliAnalysisTaskDxHFECorrelation::FinishTaskOutput()
+{
+  // end of the processing
+}
+
+void AliAnalysisTaskDxHFECorrelation::Terminate(Option_t *)
+{
+  // last action on the client
+  fOutput = dynamic_cast<TList*> (GetOutputData(1));
+  if (!fOutput) {
+    // looks like that is a valid condition if the task is run
+    // in mode "terminate"
+    return;
+  }
+}
index b2fd08f..3129401 100644 (file)
-// $Id$\r
-\r
-//**************************************************************************\r
-//* This file is property of and copyright by the ALICE Project            * \r
-//* ALICE Experiment at CERN, All rights reserved.                         *\r
-//*                                                                        *\r
-//* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *\r
-//*                  Sedat Altinpinar <Sedat.Altinpinar@cern.ch>           *\r
-//*                  Hege Erdal       <hege.erdal@gmail.com>               *\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
-/// @file   AliDxHFECorrelation.cxx\r
-/// @author Sedat Altinpinar, Hege Erdal, Matthias Richter\r
-/// @date   2012-04-25\r
-/// @brief  Worker class for D0-HF electron correlation\r
-///\r
-\r
-#include "AliDxHFECorrelation.h"\r
-#include "AliVParticle.h"\r
-#include "AliLog.h"\r
-//#include "AliAnalysisCuts.h"         // required dependency libANALYSISalice.so\r
-//#include "AliFlowTrackSimple.h"      // required dependency libPWGflowBase.so\r
-//#include "AliFlowCandidateTrack.h"   // required dependency libPWGflowTasks.so\r
-//#include "AliCFContainer.h"          // required dependency libCORRFW.so\r
-#include "TObjArray.h"\r
-#include "AliHFCorrelator.h"\r
-#include "AliAODEvent.h"\r
-#include "AliAODVertex.h"\r
-#include "TH1D.h"\r
-#include "TH2D.h"\r
-#include "TH3D.h"\r
-#include "THnSparse.h"\r
-#include "TMath.h"\r
-#include "TFile.h"\r
-#include "TCanvas.h"\r
-#include "TDatabasePDG.h"\r
-#include "TLorentzVector.h"\r
-#include "AliReducedParticle.h"\r
-#include "AliDxHFEParticleSelection.h"\r
-#include <iostream>\r
-#include <cerrno>\r
-#include <memory>\r
-\r
-using namespace std;\r
-\r
-ClassImp(AliDxHFECorrelation)\r
-\r
-AliDxHFECorrelation::AliDxHFECorrelation(const char* name)\r
-  : TNamed(name?name:"AliDxHFECorrelation", "")\r
-  , fHistograms(NULL)  \r
-  , fControlObjects(NULL)\r
-  , fCorrProperties(NULL)\r
-  , fhEventControlCorr(NULL)\r
-  , fCuts(NULL)\r
-  , fUseMC(kFALSE)\r
-  , fCorrelator(NULL)\r
-  , fUseEventMixing(kFALSE)\r
-  , fSystem(0)\r
-  , fMinPhi(-TMath::Pi()/2)\r
-  , fMaxPhi(3*TMath::Pi()/2)\r
-  , fDeltaPhi(0)\r
-  , fDeltaEta(0)\r
-  , fDimThn(0)\r
-  , fCorrArray(NULL)\r
-  , fEventType(0)\r
-{\r
-  // default constructor\r
-  // \r
-  //\r
-\r
-}\r
-\r
-const char* AliDxHFECorrelation::fgkEventControlBinNames[]={\r
-  "nEventsAll",\r
-  "nEventsSelected",\r
-  "nEventsD0",\r
-  "nEventsD0e"\r
-};\r
-\r
-AliDxHFECorrelation::~AliDxHFECorrelation()\r
-{\r
-  // destructor\r
-  //\r
-  //\r
-  if (fHistograms) delete fHistograms;\r
-  fHistograms=NULL;\r
-\r
-  // NOTE: fControlObjects owns the object, and they are deleted in the\r
-  // destructor of TList\r
-  if (fControlObjects) delete fControlObjects;\r
-  fControlObjects=NULL;\r
-  fCorrProperties=NULL;\r
-  fhEventControlCorr=NULL;\r
-  if(fCorrelator) delete fCorrelator;\r
-  fCorrelator=NULL;\r
-  if(fCorrArray) delete fCorrArray;\r
-  fCorrArray=NULL;\r
-\r
-  // NOTE: the external object is deleted elsewhere\r
-  fCuts=NULL;\r
-}\r
-\r
-int AliDxHFECorrelation::Init(const char* arguments)\r
-{\r
-  //\r
-  // Will initialize thnsparse, histogram and AliHFCorrelator\r
-  //\r
-  AliInfo("Initializing correlation objects");\r
-  ParseArguments(arguments);\r
-\r
-  //----------------------------------------------\r
-  // Setting up THnSparse \r
-  fCorrProperties=DefineTHnSparse();\r
-  AddControlObject(fCorrProperties);\r
-\r
-  //----------------------------------------------\r
-  // Histogram for storing event information\r
-\r
-  std::auto_ptr<TH1D> hEventControl(new TH1D("hEventControlCorr", "hEventControlCorr", 10, 0, 10));\r
-  if (!hEventControl.get()) {\r
-    return -ENOMEM;\r
-  }\r
-  int iLabel=0;\r
-  for (iLabel=0; iLabel<kNEventControlLabels; iLabel++)\r
-    hEventControl->GetXaxis()->SetBinLabel(iLabel+1, fgkEventControlBinNames[iLabel]);\r
-\r
-  fhEventControlCorr=hEventControl.release();\r
-  AddControlObject(fhEventControlCorr);\r
-\r
-  //----------------------------------------------\r
-  // AliHFCorrelator for Event Mixing and correlation\r
-  // \r
-  // fCuts is the hadron cut object, fSystem to switch between pp or PbPb\r
-  AliHFAssociatedTrackCuts* cuts=dynamic_cast<AliHFAssociatedTrackCuts*>(fCuts);\r
-  if (!cuts) {\r
-    if (fCuts)\r
-      AliError(Form("cuts object of wrong type %s, required AliHFAssociatedTrackCuts", fCuts->ClassName()));\r
-    else\r
-      AliError("mandatory cuts object missing");\r
-    return -EINVAL;\r
-  }\r
-  if (cuts->GetNCentPoolBins()==0 || cuts->GetCentPoolBins()==NULL ||\r
-      cuts->GetNZvtxPoolBins()==0 || cuts->GetZvtxPoolBins()==NULL) {\r
-    // the bin information is used further downstream so it\r
-    // needs to be available in order to continue\r
-    AliError(Form("inavlid object %s: bin configuration is mandatory", cuts->GetName()));\r
-    cuts->Dump();\r
-    return -EINVAL;\r
-  }\r
-  fCorrelator = new AliHFCorrelator("Correlator", cuts, fSystem); \r
-  fCorrelator->SetDeltaPhiInterval(fMinPhi,fMaxPhi); //Correct Phi Interval\r
-  fCorrelator->SetEventMixing(fUseEventMixing);      // mixing Off/On \r
-  fCorrelator->SetAssociatedParticleType(AliHFCorrelator::kElectron);\r
-  // 0: don't calculate d0; 1: return d0; 2: return d0/d0err\r
-  fCorrelator->SetApplyDisplacementCut(kFALSE); \r
-  fCorrelator->SetUseMC(fUseMC);\r
-  fCorrelator->SetUseReco(kTRUE); // Reco/MCTruth\r
-  Bool_t pooldef = fCorrelator->DefineEventPool();\r
-       \r
-  if(!pooldef) AliInfo("Warning:: Event pool not defined properly");\r
-\r
-\r
-    // ============================= EVENT MIXING CHECKS ======================================\r
-  // TODO: Not sure if all 4 histos are needed. Keep for now.. \r
-  Int_t MaxNofEvents = cuts->GetMaxNEventsInPool();\r
-  Int_t MinNofTracks = cuts->GetMinNTracksInPool();\r
-  Int_t NofCentBins = cuts->GetNCentPoolBins();\r
-  const Double_t * CentBins = cuts->GetCentPoolBins();\r
-  const Double_t defaultCentBins[] = {0,100};\r
-  if (NofCentBins==0 || CentBins==NULL) {\r
-    NofCentBins=2;\r
-    CentBins=defaultCentBins;\r
-  }\r
-  Int_t NofZVrtxBins = cuts->GetNZvtxPoolBins();\r
-  const Double_t *ZVrtxBins = cuts->GetZvtxPoolBins();\r
-  const Double_t defaultZVrtxBins[] = {-10,10};\r
-  if (NofZVrtxBins==0 || ZVrtxBins==NULL) {\r
-    NofZVrtxBins=2;\r
-    ZVrtxBins=defaultZVrtxBins;\r
-  }\r
-\r
-  Int_t nofEventPropBins =0;\r
-\r
-  if(fSystem) nofEventPropBins = 100; // PbPb centrality\r
-  if(!fSystem) nofEventPropBins = NofCentBins; // pp multiplicity\r
-\r
-  Double_t minvalue = CentBins[0];\r
-  Double_t maxvalue = CentBins[NofCentBins];\r
-  Double_t Zminvalue = ZVrtxBins[0];\r
-  Double_t Zmaxvalue = ZVrtxBins[NofCentBins];\r
-\r
-  const Double_t Nevents[]={0,2*MaxNofEvents/10,4*MaxNofEvents/10,6*MaxNofEvents/10,8*MaxNofEvents/10,MaxNofEvents};\r
-  const Double_t * events = Nevents;\r
-\r
-  TH3D * EventsPerPoolBin = new TH3D("EventsPerPoolBin","Number of events in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins,5,events);\r
-  EventsPerPoolBin->GetXaxis()->SetTitle("Centrality/multiplicity ");\r
-  EventsPerPoolBin->GetYaxis()->SetTitle("Z vertex [cm]");\r
-  EventsPerPoolBin->GetZaxis()->SetTitle("Number of events in pool bin");\r
-  if(fUseEventMixing) AddControlObject(EventsPerPoolBin);\r
-\r
-  Int_t MaxNofTracks = (MaxNofEvents+1)*MinNofTracks;\r
-  Int_t Diff = MaxNofTracks-MinNofTracks;\r
-\r
-  Double_t Ntracks[]={MinNofTracks,MinNofTracks+Diff/5,MinNofTracks+2*Diff/5,MinNofTracks+3*Diff/5,MinNofTracks+4*Diff/5,MaxNofTracks};\r
-  Double_t  * trackN = Ntracks;\r
-\r
-  TH3D * NofTracksPerPoolBin = new TH3D("NofTracksPerPoolBin","Number of tracks in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins,5,trackN);\r
-  NofTracksPerPoolBin->GetXaxis()->SetTitle("Centrality/multiplicity ");\r
-  NofTracksPerPoolBin->GetYaxis()->SetTitle("Z vertex [cm]");\r
-  NofTracksPerPoolBin->GetZaxis()->SetTitle("Number of tracks per bin");\r
-\r
-  if(fUseEventMixing) AddControlObject(NofTracksPerPoolBin);\r
-\r
-  TH2D * NofPoolBinCalls = new TH2D("NofPoolBinCalls","Number of tracks in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins);\r
-  NofPoolBinCalls->GetXaxis()->SetTitle("Centrality/multiplicity ");\r
-  NofPoolBinCalls->GetYaxis()->SetTitle("Z vertex [cm]");\r
-  if(fUseEventMixing) AddControlObject(NofPoolBinCalls);\r
-\r
-  TH2D * EventProps = new TH2D("EventProps","Number of tracks in bin pool",nofEventPropBins,minvalue,maxvalue,100,Zminvalue,Zmaxvalue);\r
-  EventProps->GetXaxis()->SetTitle("Centrality/multiplicity ");\r
-  EventProps->GetYaxis()->SetTitle("Z vertex [cm]");\r
-  if(fUseEventMixing) AddControlObject(EventProps);\r
-\r
-  return 0;\r
-}\r
-\r
-int AliDxHFECorrelation::ParseArguments(const char* arguments)\r
-{\r
-  // parse arguments and set internal flags\r
-  TString strArguments(arguments);\r
-  auto_ptr<TObjArray> tokens(strArguments.Tokenize(" "));\r
-  if (!tokens.get()) return -ENOMEM;\r
-\r
-  TIter next(tokens.get());\r
-  TObject* token;\r
-  while ((token=next())) {\r
-    TString argument=token->GetName();\r
-   \r
-    if (argument.BeginsWith("event-mixing")) {\r
-      fUseEventMixing=true;\r
-      continue;\r
-    }\r
-      \r
-    if (argument.BeginsWith("use-mc")) {\r
-      fUseMC=true;\r
-      continue;\r
-    }\r
-    if (argument.BeginsWith("system=")) {\r
-      argument.ReplaceAll("system=", "");\r
-      if (argument.CompareTo("pp")==0) fSystem=0;\r
-      else if (argument.CompareTo("Pb-Pb")==0) fSystem=1;\r
-      else {\r
-       AliWarning(Form("can not set collision system, unknown parameter '%s'", argument.Data()));\r
-       // TODO: check what makes sense\r
-       fSystem=0;\r
-      }\r
-      continue;\r
-    }\r
-    AliWarning(Form("unknown argument '%s'", argument.Data()));\r
-      \r
-  }\r
-\r
-  return 0;\r
-}\r
-\r
-THnSparse* AliDxHFECorrelation::DefineTHnSparse()\r
-{\r
-  //\r
-  //Defines the THnSparse. For now, only calls CreateControlTHnSparse\r
-\r
-  // here is the only place to change the dimension\r
-  static const int sizeEventdphi = 7;  \r
-  InitTHnSparseArray(sizeEventdphi);\r
-  const double pi=TMath::Pi();\r
-\r
-  //TODO: add phi for electron??\r
-  //                                           0           1       2      3         4     5      6   \r
-  //                                         D0invmass   PtD0    PhiD0  PtbinD0    Pte   dphi   deta   \r
-  int         binsEventdphi[sizeEventdphi] = {   200,      1000,   100,    21,     1000,   100,    100};\r
-  double      minEventdphi [sizeEventdphi] = { 1.5648,      0,       0,     0,       0,  fMinPhi, -2};\r
-  double      maxEventdphi [sizeEventdphi] = { 2.1648,     100,   2*pi,    20,      100, fMaxPhi, 2};\r
-  const char* nameEventdphi[sizeEventdphi] = {\r
-    "D0InvMass",\r
-    "PtD0",\r
-    "PhiD0",\r
-    "PtBinD0",\r
-    "PtEl",\r
-    "#Delta#Phi", \r
-    "#Delta#eta"\r
-  };\r
-\r
-  TString name;\r
-  name.Form("%s info", GetName());\r
-\r
-\r
-  return CreateControlTHnSparse(name,sizeEventdphi,binsEventdphi,minEventdphi,maxEventdphi,nameEventdphi);\r
-\r
-}\r
-\r
-THnSparse* AliDxHFECorrelation::CreateControlTHnSparse(const char* name,\r
-                                                            int thnSize,\r
-                                                            int* thnBins,\r
-                                                            double* thnMin,\r
-                                                            double* thnMax,\r
-                                                            const char** binLabels) const\r
-{\r
-  //\r
-  // Creates THnSparse.\r
-  //\r
-\r
-  AliInfo("Setting up THnSparse");\r
-\r
-  std::auto_ptr<THnSparseD> th(new THnSparseD(name, name, thnSize, thnBins, thnMin, thnMax));\r
-  if (th.get()==NULL) {\r
-    return NULL;\r
-  }\r
-  for (int iLabel=0; iLabel<thnSize; iLabel++) {\r
-    th->GetAxis(iLabel)->SetTitle(binLabels[iLabel]);    \r
-   \r
-  }\r
-  return th.release();\r
-\r
-}\r
-\r
-int AliDxHFECorrelation::AddControlObject(TObject* pObj)\r
-{\r
-  AliInfo("Adding object");\r
-  /// add control object to list, the base class becomes owner of the object\r
-  if (!pObj) return -EINVAL;\r
-  if (!fControlObjects) {\r
-    fControlObjects=new TList;\r
-    if (!fControlObjects) return -ENOMEM;\r
-    fControlObjects->SetOwner();\r
-  }\r
-  if (fControlObjects->FindObject(pObj->GetName())) {\r
-    AliError(Form("ignoring duplicate object '%s' of type %s", pObj->GetName(), pObj->ClassName()));\r
-    return -EEXIST;\r
-  }\r
-  fControlObjects->Add(pObj);\r
-  return 0;\r
-}\r
-\r
-int AliDxHFECorrelation::HistogramEventProperties(int bin)\r
-{\r
-  /// histogram event properties\r
-  if (!fhEventControlCorr) return 0;\r
-  fhEventControlCorr->Fill(bin);\r
-\r
-  return 0;\r
-}\r
-\r
-int AliDxHFECorrelation::Fill(const TObjArray* triggerCandidates, const TObjArray* associatedTracks, const AliVEvent* pEvent)\r
-{\r
-  //\r
-  // will use AliHFCorrelator to process D0-electron pair and then fill THnSparse.\r
-  //\r
-  if (!triggerCandidates || !associatedTracks) return -EINVAL;\r
-  if (!fControlObjects) {\r
-    Init();\r
-  }\r
-  if (!fControlObjects) {\r
-    AliError("Initialisation failed, can not fill THnSparse");\r
-    return -ENODEV;\r
-  }\r
-  // set the event to be processed\r
-  // TODO: change the correlator class to take the const pointer\r
-  if (!fCorrelator) {\r
-    AliError("working class instance fCorrelator missing");\r
-    return -ENODEV;\r
-  }\r
-  fCorrelator->SetAODEvent(dynamic_cast<AliAODEvent*>(const_cast<AliVEvent*>(pEvent))); \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 1;\r
-  }\r
-\r
-  TIter itrigger(triggerCandidates);\r
-  TObject* otrigger=NULL;\r
-  int ctrigger=-1;\r
-\r
-  // For the moment this is very specific to D0-electron correlation. Should be \r
-  // changed to be less specific. \r
-  while ((otrigger=itrigger())!=NULL) {\r
-    // loop over trigger D0 particle\r
-    ctrigger++;\r
-    AliReducedParticle* ptrigger=dynamic_cast<AliReducedParticle*>(otrigger);\r
-    if (!ptrigger)  continue;\r
-\r
-    Double_t phiTrigger = ptrigger->Phi();\r
-    Double_t ptTrigger = ptrigger->Pt();\r
-    Double_t etaTrigger = ptrigger->Eta();\r
-\r
-    // set the phi of the D meson in the correct range\r
-    // TODO: Is this correct to do this??\r
-    phiTrigger = fCorrelator->SetCorrectPhiRange(phiTrigger);\r
-    // pass to the object the necessary trigger part parameters\r
-    fCorrelator->SetTriggerParticleProperties(ptTrigger,phiTrigger,etaTrigger); \r
-\r
-    Bool_t execPool = fCorrelator->ProcessEventPool();\r
-    if(fUseEventMixing && !execPool) {\r
-      AliInfo("Mixed event analysis: pool is not ready");\r
-      continue;\r
-    }\r
-    Int_t NofEventsinPool = 1;\r
-    if(fUseEventMixing) NofEventsinPool = fCorrelator->GetNofEventsInPool();\r
-               \r
-    // loop on events in the pool; if it is SE analysis, stops at one\r
-    for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){\r
-      Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix, associatedTracks);\r
-                       \r
-      if(!analyzetracks) {\r
-       AliInfo("AliHFCorrelator::Cannot process the track array");\r
-       continue;\r
-      }\r
-\r
-      Int_t NofTracks = fCorrelator->GetNofTracks();\r
-\r
-      // looping on track candidates\r
-      for(Int_t iTrack = 0; iTrack<NofTracks; iTrack++){ \r
-       Bool_t runcorrelation = fCorrelator->Correlate(iTrack);\r
-       if(!runcorrelation) continue;\r
-                       \r
-       fDeltaPhi = fCorrelator->GetDeltaPhi();\r
-       fDeltaEta = fCorrelator->GetDeltaEta();\r
-       \r
-       AliReducedParticle *assoc = fCorrelator->GetAssociatedParticle();\r
-\r
-       FillParticleProperties(ptrigger,assoc,ParticleProperties(),GetDimTHnSparse());\r
-       fCorrProperties->Fill(ParticleProperties());\r
-\r
-      } // loop over electron tracks in event\r
-    } // loop over events in pool\r
-  } // loop over trigger particle\r
-\r
-  Bool_t updated = fCorrelator->PoolUpdate(associatedTracks);\r
-  if(fUseEventMixing){\r
-    if(!updated) AliInfo("Pool was not updated");\r
-    else {\r
-      EventMixingChecks(pEvent);\r
-      AliInfo("Pool was updated");\r
-    }\r
-  }\r
-  return 0;\r
-}\r
-\r
-\r
-\r
-int AliDxHFECorrelation::FillParticleProperties(AliVParticle* tr, AliVParticle *as, Double_t* data, int dimension) const\r
-{\r
-  // fill the data array from the particle data\r
-  if (!data) return -EINVAL;\r
-  AliReducedParticle *ptrigger=(AliReducedParticle*)tr;\r
-  AliReducedParticle *assoc=(AliReducedParticle*)as;\r
-  if (!ptrigger || !assoc) return -ENODATA;\r
-  int i=0;\r
-  if (dimension!=GetDimTHnSparse()) {\r
-    // TODO: think about filling only the available data and throwing a warning\r
-    return -ENOSPC;\r
-  }\r
-  data[i++]=ptrigger->GetInvMass();\r
-  data[i++]=ptrigger->Pt();\r
-  data[i++]=ptrigger->Phi();\r
-  data[i++]=ptrigger->GetPtBin(); \r
-  data[i++]=assoc->Pt();\r
-  data[i++]=GetDeltaPhi();\r
-  data[i++]=GetDeltaEta();\r
-\r
-  return i;\r
-}\r
-\r
-void AliDxHFECorrelation::Clear(Option_t * /*option*/)\r
-{\r
-  /// overloaded from TObject: cleanup\r
-\r
-  // nothing to be done so far\r
-  return TObject::Clear();\r
-}\r
-\r
-void AliDxHFECorrelation::Print(Option_t */*option*/) const\r
-{\r
-  /// overloaded from TObject: print info\r
-  cout << "====================================================================" << endl;\r
-  TObject::Print();\r
-  if (fHistograms) {\r
-    fHistograms->Print();\r
-  }\r
-}\r
-\r
-void AliDxHFECorrelation::Draw(Option_t */*option*/)\r
-{\r
-  /// overloaded from TObject: draw histograms\r
-}\r
-\r
-TObject* AliDxHFECorrelation::FindObject(const char *name) const\r
-{\r
-  /// overloaded from TObject: find object by name\r
-  if (fControlObjects) {\r
-    return fControlObjects->FindObject(name);\r
-  }\r
-  return NULL;\r
-}\r
-\r
-TObject* AliDxHFECorrelation::FindObject(const TObject *obj) const\r
-{\r
-  /// overloaded from TObject: find object by pointer\r
-  if (fControlObjects) {\r
-    return fControlObjects->FindObject(obj);\r
-  }\r
-  return NULL;\r
-}\r
-\r
-void AliDxHFECorrelation::SaveAs(const char *filename, Option_t */*option*/) const\r
-{\r
-  /// overloaded from TObject: save to file\r
-  std::auto_ptr<TFile> output(TFile::Open(filename, "RECREATE"));\r
-  if (!output.get() || output->IsZombie()) {\r
-    AliError(Form("can not open file %s from writing", filename));\r
-    return;\r
-  }\r
-  output->cd();\r
-  if (fControlObjects) fControlObjects->Write();\r
-  output->Close();\r
-}\r
-\r
-AliDxHFECorrelation& AliDxHFECorrelation::operator+=(const AliDxHFECorrelation& other)\r
-{\r
-  /// add histograms from another instance\r
-  // TODO - need to change this to ThnSparse?\r
-  if (!fHistograms || !other.fHistograms) return *this;\r
-  \r
-  for (int i=0; i<kNofHistograms; i++) {\r
-    if (fHistograms->At(i)==NULL || other.fHistograms->At(i)==NULL) continue;\r
-    TH1* target=reinterpret_cast<TH1*>(fHistograms->At(i));\r
-    TH1* source=reinterpret_cast<TH1*>(other.fHistograms->At(i));\r
-    if (!target || !source) continue;\r
-    TString name(fHistograms->At(i)->GetName());\r
-    if (name.CompareTo(target->GetName())!=0) {\r
-      AliWarning(Form("skipping incompatible objects at position %d: %s vs %s", i, source->GetName(), target->GetName()));\r
-      continue;\r
-    }\r
-    if (source->IsA()!=target->IsA()) {\r
-      AliWarning(Form("skipping incompatible classes at position %d: %s vs %s", i, source->ClassName(), target->ClassName()));\r
-      continue;\r
-    }\r
-    target->Add(source);\r
-  }\r
-  return *this;\r
-}\r
-\r
-\r
-//____________________________  Run checks on event mixing ___________________________________________________\r
-void AliDxHFECorrelation::EventMixingChecks(const AliVEvent* pEvent){\r
-       \r
-  AliAODEvent *AOD= (AliAODEvent*)(pEvent);\r
-  AliCentrality *centralityObj = 0;\r
-  Int_t multiplicity = -1;\r
-  Double_t MultipOrCent = -1;\r
-       \r
-  // get the pool for event mixing\r
-  if(!fSystem){ // pp\r
-    multiplicity = AOD->GetNTracks();\r
-    MultipOrCent = multiplicity; // convert from Int_t to Double_t\r
-  }\r
-  if(fSystem){ // PbPb         \r
-    centralityObj = AOD->GetHeader()->GetCentralityP();\r
-    MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");\r
-    AliInfo(Form("Centrality is %f", MultipOrCent));\r
-  }\r
-       \r
-  AliAODVertex *vtx = AOD->GetPrimaryVertex();\r
-  Double_t zvertex = vtx->GetZ(); // zvertex\r
-\r
-  AliEventPool *pool = fCorrelator->GetPool();\r
-\r
-  ((TH2D*)fControlObjects->FindObject("NofPoolBinCalls"))->Fill(MultipOrCent,zvertex); // number of calls of pool\r
-  ((TH2D*)fControlObjects->FindObject("EventProps"))->Fill(MultipOrCent,zvertex); // event properties\r
-  ((TH3D*)fControlObjects->FindObject("EventsPerPoolBin"))->Fill(MultipOrCent,zvertex,pool->NTracksInPool()); // number of events in the pool\r
-  ((TH3D*)fControlObjects->FindObject("NofTracksPerPoolBin"))->Fill(MultipOrCent,zvertex,pool->GetCurrentNEvents()); // number of calls of pool\r
-}\r
-       \r
+// $Id$
+
+//**************************************************************************
+//* This file is property of and copyright by the ALICE Project            * 
+//* ALICE Experiment at CERN, All rights reserved.                         *
+//*                                                                        *
+//* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
+//*                  Sedat Altinpinar <Sedat.Altinpinar@cern.ch>           *
+//*                  Hege Erdal       <hege.erdal@gmail.com>               *
+//*                                                                        *
+//* 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.                  *
+//**************************************************************************
+
+/// @file   AliDxHFECorrelation.cxx
+/// @author Sedat Altinpinar, Hege Erdal, Matthias Richter
+/// @date   2012-04-25
+/// @brief  Worker class for D0-HF electron correlation
+///
+
+#include "AliDxHFECorrelation.h"
+#include "AliVParticle.h"
+#include "AliLog.h"
+//#include "AliAnalysisCuts.h"         // required dependency libANALYSISalice.so
+//#include "AliFlowTrackSimple.h"      // required dependency libPWGflowBase.so
+//#include "AliFlowCandidateTrack.h"   // required dependency libPWGflowTasks.so
+//#include "AliCFContainer.h"          // required dependency libCORRFW.so
+#include "TObjArray.h"
+#include "AliHFCorrelator.h"
+#include "AliAODEvent.h"
+#include "AliAODVertex.h"
+#include "TH1D.h"
+#include "TH2D.h"
+#include "TH3D.h"
+#include "THnSparse.h"
+#include "TMath.h"
+#include "TFile.h"
+#include "TCanvas.h"
+#include "TDatabasePDG.h"
+#include "TLorentzVector.h"
+#include "AliReducedParticle.h"
+#include "AliDxHFEParticleSelection.h"
+#include <iostream>
+#include <cerrno>
+#include <memory>
+
+using namespace std;
+
+ClassImp(AliDxHFECorrelation)
+
+AliDxHFECorrelation::AliDxHFECorrelation(const char* name)
+  : TNamed(name?name:"AliDxHFECorrelation", "")
+  , fHistograms(NULL)  
+  , fControlObjects(NULL)
+  , fCorrProperties(NULL)
+  , fhEventControlCorr(NULL)
+  , fCuts(NULL)
+  , fUseMC(kFALSE)
+  , fCorrelator(NULL)
+  , fUseEventMixing(kFALSE)
+  , fSystem(0)
+  , fMinPhi(-TMath::Pi()/2)
+  , fMaxPhi(3*TMath::Pi()/2)
+  , fDeltaPhi(0)
+  , fDeltaEta(0)
+  , fDimThn(0)
+  , fCorrArray(NULL)
+  , fEventType(0)
+{
+  // default constructor
+  // 
+  //
+
+}
+
+const char* AliDxHFECorrelation::fgkEventControlBinNames[]={
+  "nEventsAll",
+  "nEventsSelected",
+  "nEventsD0",
+  "nEventsD0e"
+};
+
+AliDxHFECorrelation::~AliDxHFECorrelation()
+{
+  // destructor
+  //
+  //
+  if (fHistograms) delete fHistograms;
+  fHistograms=NULL;
+
+  // NOTE: fControlObjects owns the object, and they are deleted in the
+  // destructor of TList
+  if (fControlObjects) delete fControlObjects;
+  fControlObjects=NULL;
+  fCorrProperties=NULL;
+  fhEventControlCorr=NULL;
+  if(fCorrelator) delete fCorrelator;
+  fCorrelator=NULL;
+  if(fCorrArray) delete fCorrArray;
+  fCorrArray=NULL;
+
+  // NOTE: the external object is deleted elsewhere
+  fCuts=NULL;
+}
+
+int AliDxHFECorrelation::Init(const char* arguments)
+{
+  //
+  // Will initialize thnsparse, histogram and AliHFCorrelator
+  //
+  AliInfo("Initializing correlation objects");
+  ParseArguments(arguments);
+
+  //----------------------------------------------
+  // Setting up THnSparse 
+  fCorrProperties=DefineTHnSparse();
+  AddControlObject(fCorrProperties);
+
+  //----------------------------------------------
+  // Histogram for storing event information
+
+  std::auto_ptr<TH1D> hEventControl(new TH1D("hEventControlCorr", "hEventControlCorr", 10, 0, 10));
+  if (!hEventControl.get()) {
+    return -ENOMEM;
+  }
+  int iLabel=0;
+  for (iLabel=0; iLabel<kNEventControlLabels; iLabel++)
+    hEventControl->GetXaxis()->SetBinLabel(iLabel+1, fgkEventControlBinNames[iLabel]);
+
+  fhEventControlCorr=hEventControl.release();
+  AddControlObject(fhEventControlCorr);
+
+  //----------------------------------------------
+  // AliHFCorrelator for Event Mixing and correlation
+  // 
+  // fCuts is the hadron cut object, fSystem to switch between pp or PbPb
+  AliHFAssociatedTrackCuts* cuts=dynamic_cast<AliHFAssociatedTrackCuts*>(fCuts);
+  if (!cuts) {
+    if (fCuts)
+      AliError(Form("cuts object of wrong type %s, required AliHFAssociatedTrackCuts", fCuts->ClassName()));
+    else
+      AliError("mandatory cuts object missing");
+    return -EINVAL;
+  }
+  if (cuts->GetNCentPoolBins()==0 || cuts->GetCentPoolBins()==NULL ||
+      cuts->GetNZvtxPoolBins()==0 || cuts->GetZvtxPoolBins()==NULL) {
+    // the bin information is used further downstream so it
+    // needs to be available in order to continue
+    AliError(Form("inavlid object %s: bin configuration is mandatory", cuts->GetName()));
+    cuts->Dump();
+    return -EINVAL;
+  }
+  cuts->Print("");
+  fCorrelator = new AliHFCorrelator("Correlator", cuts, fSystem); 
+  fCorrelator->SetDeltaPhiInterval(fMinPhi,fMaxPhi); //Correct Phi Interval
+  fCorrelator->SetEventMixing(fUseEventMixing);      // mixing Off/On 
+  fCorrelator->SetAssociatedParticleType(AliHFCorrelator::kElectron);
+  // 0: don't calculate d0; 1: return d0; 2: return d0/d0err
+  fCorrelator->SetApplyDisplacementCut(kFALSE); 
+  fCorrelator->SetUseMC(fUseMC);
+  fCorrelator->SetUseReco(kTRUE); // Reco/MCTruth
+  Bool_t pooldef = fCorrelator->DefineEventPool();
+       
+  if(!pooldef) AliInfo("Warning:: Event pool not defined properly");
+
+
+    // ============================= EVENT MIXING CHECKS ======================================
+  // TODO: Not sure if all 4 histos are needed. Keep for now.. 
+  Int_t MaxNofEvents = cuts->GetMaxNEventsInPool();
+  Int_t MinNofTracks = cuts->GetMinNTracksInPool();
+  Int_t NofCentBins = cuts->GetNCentPoolBins();
+  const Double_t * CentBins = cuts->GetCentPoolBins();
+  const Double_t defaultCentBins[] = {0,100};
+  if (NofCentBins==0 || CentBins==NULL) {
+    NofCentBins=1; // note: array dimension minus one, because of bin is bound by upper and lower
+    CentBins=defaultCentBins;
+  }
+  Int_t NofZVrtxBins = cuts->GetNZvtxPoolBins();
+  const Double_t *ZVrtxBins = cuts->GetZvtxPoolBins();
+  const Double_t defaultZVrtxBins[] = {-10,10};
+  if (NofZVrtxBins==0 || ZVrtxBins==NULL) {
+    NofZVrtxBins=1; // note: array dimension minus one, because of bin is bound by upper and lower
+    ZVrtxBins=defaultZVrtxBins;
+  }
+
+  Int_t nofEventPropBins =0;
+
+  if(fSystem) nofEventPropBins = 100; // PbPb centrality
+  if(!fSystem) nofEventPropBins = NofCentBins; // pp multiplicity
+
+  Double_t minvalue = CentBins[0];
+  Double_t maxvalue = CentBins[NofCentBins];
+  Double_t Zminvalue = ZVrtxBins[0];
+  Double_t Zmaxvalue = ZVrtxBins[NofCentBins];
+
+  const Double_t Nevents[]={0,2*MaxNofEvents/10,4*MaxNofEvents/10,6*MaxNofEvents/10,8*MaxNofEvents/10,MaxNofEvents};
+  const Double_t * events = Nevents;
+
+  TH3D * EventsPerPoolBin = new TH3D("EventsPerPoolBin","Number of events in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins,5,events);
+  EventsPerPoolBin->GetXaxis()->SetTitle("Centrality/multiplicity ");
+  EventsPerPoolBin->GetYaxis()->SetTitle("Z vertex [cm]");
+  EventsPerPoolBin->GetZaxis()->SetTitle("Number of events in pool bin");
+  if(fUseEventMixing) AddControlObject(EventsPerPoolBin);
+
+  Int_t MaxNofTracks = (MaxNofEvents+1)*MinNofTracks;
+  Int_t Diff = MaxNofTracks-MinNofTracks;
+
+  Double_t Ntracks[]={MinNofTracks,MinNofTracks+Diff/5,MinNofTracks+2*Diff/5,MinNofTracks+3*Diff/5,MinNofTracks+4*Diff/5,MaxNofTracks};
+  Double_t  * trackN = Ntracks;
+
+  TH3D * NofTracksPerPoolBin = new TH3D("NofTracksPerPoolBin","Number of tracks in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins,5,trackN);
+  NofTracksPerPoolBin->GetXaxis()->SetTitle("Centrality/multiplicity ");
+  NofTracksPerPoolBin->GetYaxis()->SetTitle("Z vertex [cm]");
+  NofTracksPerPoolBin->GetZaxis()->SetTitle("Number of tracks per bin");
+
+  if(fUseEventMixing) AddControlObject(NofTracksPerPoolBin);
+
+  TH2D * NofPoolBinCalls = new TH2D("NofPoolBinCalls","Number of tracks in bin pool",NofCentBins,CentBins,NofZVrtxBins,ZVrtxBins);
+  NofPoolBinCalls->GetXaxis()->SetTitle("Centrality/multiplicity ");
+  NofPoolBinCalls->GetYaxis()->SetTitle("Z vertex [cm]");
+  if(fUseEventMixing) AddControlObject(NofPoolBinCalls);
+
+  TH2D * EventProps = new TH2D("EventProps","Number of tracks in bin pool",nofEventPropBins,minvalue,maxvalue,100,Zminvalue,Zmaxvalue);
+  EventProps->GetXaxis()->SetTitle("Centrality/multiplicity ");
+  EventProps->GetYaxis()->SetTitle("Z vertex [cm]");
+  if(fUseEventMixing) AddControlObject(EventProps);
+
+  return 0;
+}
+
+int AliDxHFECorrelation::ParseArguments(const char* arguments)
+{
+  // parse arguments and set internal flags
+  TString strArguments(arguments);
+  auto_ptr<TObjArray> tokens(strArguments.Tokenize(" "));
+  if (!tokens.get()) return -ENOMEM;
+
+  TIter next(tokens.get());
+  TObject* token;
+  while ((token=next())) {
+    TString argument=token->GetName();
+   
+    if (argument.BeginsWith("event-mixing")) {
+      fUseEventMixing=true;
+      continue;
+    }
+      
+    if (argument.BeginsWith("use-mc")) {
+      fUseMC=true;
+      continue;
+    }
+    if (argument.BeginsWith("system=")) {
+      argument.ReplaceAll("system=", "");
+      if (argument.CompareTo("pp")==0) fSystem=0;
+      else if (argument.CompareTo("Pb-Pb")==0) fSystem=1;
+      else {
+       AliWarning(Form("can not set collision system, unknown parameter '%s'", argument.Data()));
+       // TODO: check what makes sense
+       fSystem=0;
+      }
+      continue;
+    }
+    AliWarning(Form("unknown argument '%s'", argument.Data()));
+      
+  }
+
+  return 0;
+}
+
+THnSparse* AliDxHFECorrelation::DefineTHnSparse()
+{
+  //
+  //Defines the THnSparse. For now, only calls CreateControlTHnSparse
+
+  // here is the only place to change the dimension
+  static const int sizeEventdphi = 7;  
+  InitTHnSparseArray(sizeEventdphi);
+  const double pi=TMath::Pi();
+
+  //TODO: add phi for electron??
+  //                                           0           1       2      3         4     5      6   
+  //                                         D0invmass   PtD0    PhiD0  PtbinD0    Pte   dphi   deta   
+  int         binsEventdphi[sizeEventdphi] = {   200,      1000,   100,    21,     1000,   100,    100};
+  double      minEventdphi [sizeEventdphi] = { 1.5648,      0,       0,     0,       0,  fMinPhi, -2};
+  double      maxEventdphi [sizeEventdphi] = { 2.1648,     100,   2*pi,    20,      100, fMaxPhi, 2};
+  const char* nameEventdphi[sizeEventdphi] = {
+    "D0InvMass",
+    "PtD0",
+    "PhiD0",
+    "PtBinD0",
+    "PtEl",
+    "#Delta#Phi", 
+    "#Delta#eta"
+  };
+
+  TString name;
+  name.Form("%s info", GetName());
+
+
+  return CreateControlTHnSparse(name,sizeEventdphi,binsEventdphi,minEventdphi,maxEventdphi,nameEventdphi);
+
+}
+
+THnSparse* AliDxHFECorrelation::CreateControlTHnSparse(const char* name,
+                                                            int thnSize,
+                                                            int* thnBins,
+                                                            double* thnMin,
+                                                            double* thnMax,
+                                                            const char** binLabels) const
+{
+  //
+  // Creates THnSparse.
+  //
+
+  AliInfo("Setting up THnSparse");
+
+  std::auto_ptr<THnSparseD> th(new THnSparseD(name, name, thnSize, thnBins, thnMin, thnMax));
+  if (th.get()==NULL) {
+    return NULL;
+  }
+  for (int iLabel=0; iLabel<thnSize; iLabel++) {
+    th->GetAxis(iLabel)->SetTitle(binLabels[iLabel]);    
+   
+  }
+  return th.release();
+
+}
+
+int AliDxHFECorrelation::AddControlObject(TObject* pObj)
+{
+  AliInfo("Adding object");
+  /// add control object to list, the base class becomes owner of the object
+  if (!pObj) return -EINVAL;
+  if (!fControlObjects) {
+    fControlObjects=new TList;
+    if (!fControlObjects) return -ENOMEM;
+    fControlObjects->SetOwner();
+  }
+  if (fControlObjects->FindObject(pObj->GetName())) {
+    AliError(Form("ignoring duplicate object '%s' of type %s", pObj->GetName(), pObj->ClassName()));
+    return -EEXIST;
+  }
+  fControlObjects->Add(pObj);
+  return 0;
+}
+
+int AliDxHFECorrelation::HistogramEventProperties(int bin)
+{
+  /// histogram event properties
+  if (!fhEventControlCorr) return 0;
+  fhEventControlCorr->Fill(bin);
+
+  return 0;
+}
+
+int AliDxHFECorrelation::Fill(const TObjArray* triggerCandidates, const TObjArray* associatedTracks, const AliVEvent* pEvent)
+{
+  //
+  // will use AliHFCorrelator to process D0-electron pair and then fill THnSparse.
+  //
+  if (!triggerCandidates || !associatedTracks) return -EINVAL;
+  if (!fControlObjects) {
+    Init();
+  }
+  if (!fControlObjects) {
+    AliError("Initialisation failed, can not fill THnSparse");
+    return -ENODEV;
+  }
+  // set the event to be processed
+  // TODO: change the correlator class to take the const pointer
+  if (!fCorrelator) {
+    AliError("working class instance fCorrelator missing");
+    return -ENODEV;
+  }
+  fCorrelator->SetAODEvent(dynamic_cast<AliAODEvent*>(const_cast<AliVEvent*>(pEvent))); 
+
+  Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing
+  if(!correlatorON) {
+    AliError("AliHFCorrelator didn't initialize the pool correctly or processed a bad event");
+    return 1;
+  }
+
+  TIter itrigger(triggerCandidates);
+  TObject* otrigger=NULL;
+  int ctrigger=-1;
+
+  // For the moment this is very specific to D0-electron correlation. Should be 
+  // changed to be less specific. 
+  while ((otrigger=itrigger())!=NULL) {
+    // loop over trigger D0 particle
+    ctrigger++;
+    AliReducedParticle* ptrigger=dynamic_cast<AliReducedParticle*>(otrigger);
+    if (!ptrigger)  continue;
+
+    Double_t phiTrigger = ptrigger->Phi();
+    Double_t ptTrigger = ptrigger->Pt();
+    Double_t etaTrigger = ptrigger->Eta();
+
+    // set the phi of the D meson in the correct range
+    // TODO: Is this correct to do this??
+    phiTrigger = fCorrelator->SetCorrectPhiRange(phiTrigger);
+    // pass to the object the necessary trigger part parameters
+    fCorrelator->SetTriggerParticleProperties(ptTrigger,phiTrigger,etaTrigger); 
+
+    Bool_t execPool = fCorrelator->ProcessEventPool();
+    if(fUseEventMixing && !execPool) {
+      AliInfo("Mixed event analysis: pool is not ready");
+      continue;
+    }
+    Int_t NofEventsinPool = 1;
+    if(fUseEventMixing) NofEventsinPool = fCorrelator->GetNofEventsInPool();
+               
+    // loop on events in the pool; if it is SE analysis, stops at one
+    for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){
+      Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix, associatedTracks);
+                       
+      if(!analyzetracks) {
+       AliError("AliHFCorrelator::Cannot process the track array");
+       continue;
+      }
+
+      Int_t NofTracks = fCorrelator->GetNofTracks();
+
+      // looping on track candidates
+      for(Int_t iTrack = 0; iTrack<NofTracks; iTrack++){ 
+       Bool_t runcorrelation = fCorrelator->Correlate(iTrack);
+       if(!runcorrelation) continue;
+                       
+       fDeltaPhi = fCorrelator->GetDeltaPhi();
+       fDeltaEta = fCorrelator->GetDeltaEta();
+       
+       AliReducedParticle *assoc = fCorrelator->GetAssociatedParticle();
+
+       FillParticleProperties(ptrigger,assoc,ParticleProperties(),GetDimTHnSparse());
+       fCorrProperties->Fill(ParticleProperties());
+
+      } // loop over electron tracks in event
+    } // loop over events in pool
+  } // loop over trigger particle
+
+  Bool_t updated = fCorrelator->PoolUpdate(associatedTracks);
+  if(fUseEventMixing){
+    if(!updated) AliInfo("Pool was not updated");
+    else {
+      EventMixingChecks(pEvent);
+      AliInfo("Pool was updated");
+    }
+  }
+  return 0;
+}
+
+
+
+int AliDxHFECorrelation::FillParticleProperties(AliVParticle* tr, AliVParticle *as, Double_t* data, int dimension) const
+{
+  // fill the data array from the particle data
+  if (!data) return -EINVAL;
+  AliReducedParticle *ptrigger=(AliReducedParticle*)tr;
+  AliReducedParticle *assoc=(AliReducedParticle*)as;
+  if (!ptrigger || !assoc) return -ENODATA;
+  int i=0;
+  if (dimension!=GetDimTHnSparse()) {
+    // TODO: think about filling only the available data and throwing a warning
+    return -ENOSPC;
+  }
+  data[i++]=ptrigger->GetInvMass();
+  data[i++]=ptrigger->Pt();
+  data[i++]=ptrigger->Phi();
+  data[i++]=ptrigger->GetPtBin(); 
+  data[i++]=assoc->Pt();
+  data[i++]=GetDeltaPhi();
+  data[i++]=GetDeltaEta();
+
+  return i;
+}
+
+void AliDxHFECorrelation::Clear(Option_t * /*option*/)
+{
+  /// overloaded from TObject: cleanup
+
+  // nothing to be done so far
+  return TObject::Clear();
+}
+
+void AliDxHFECorrelation::Print(Option_t */*option*/) const
+{
+  /// overloaded from TObject: print info
+  cout << "====================================================================" << endl;
+  TObject::Print();
+  if (fHistograms) {
+    fHistograms->Print();
+  }
+}
+
+void AliDxHFECorrelation::Draw(Option_t */*option*/)
+{
+  /// overloaded from TObject: draw histograms
+}
+
+TObject* AliDxHFECorrelation::FindObject(const char *name) const
+{
+  /// overloaded from TObject: find object by name
+  if (fControlObjects) {
+    return fControlObjects->FindObject(name);
+  }
+  return NULL;
+}
+
+TObject* AliDxHFECorrelation::FindObject(const TObject *obj) const
+{
+  /// overloaded from TObject: find object by pointer
+  if (fControlObjects) {
+    return fControlObjects->FindObject(obj);
+  }
+  return NULL;
+}
+
+void AliDxHFECorrelation::SaveAs(const char *filename, Option_t */*option*/) const
+{
+  /// overloaded from TObject: save to file
+  std::auto_ptr<TFile> output(TFile::Open(filename, "RECREATE"));
+  if (!output.get() || output->IsZombie()) {
+    AliError(Form("can not open file %s from writing", filename));
+    return;
+  }
+  output->cd();
+  if (fControlObjects) fControlObjects->Write();
+  output->Close();
+}
+
+AliDxHFECorrelation& AliDxHFECorrelation::operator+=(const AliDxHFECorrelation& other)
+{
+  /// add histograms from another instance
+  // TODO - need to change this to ThnSparse?
+  if (!fHistograms || !other.fHistograms) return *this;
+  
+  for (int i=0; i<kNofHistograms; i++) {
+    if (fHistograms->At(i)==NULL || other.fHistograms->At(i)==NULL) continue;
+    TH1* target=reinterpret_cast<TH1*>(fHistograms->At(i));
+    TH1* source=reinterpret_cast<TH1*>(other.fHistograms->At(i));
+    if (!target || !source) continue;
+    TString name(fHistograms->At(i)->GetName());
+    if (name.CompareTo(target->GetName())!=0) {
+      AliWarning(Form("skipping incompatible objects at position %d: %s vs %s", i, source->GetName(), target->GetName()));
+      continue;
+    }
+    if (source->IsA()!=target->IsA()) {
+      AliWarning(Form("skipping incompatible classes at position %d: %s vs %s", i, source->ClassName(), target->ClassName()));
+      continue;
+    }
+    target->Add(source);
+  }
+  return *this;
+}
+
+
+//____________________________  Run checks on event mixing ___________________________________________________
+void AliDxHFECorrelation::EventMixingChecks(const AliVEvent* pEvent){
+       
+  AliAODEvent *AOD= (AliAODEvent*)(pEvent);
+  AliCentrality *centralityObj = 0;
+  Int_t multiplicity = -1;
+  Double_t MultipOrCent = -1;
+       
+  // get the pool for event mixing
+  if(!fSystem){ // pp
+    multiplicity = AOD->GetNTracks();
+    MultipOrCent = multiplicity; // convert from Int_t to Double_t
+  }
+  if(fSystem){ // PbPb         
+    centralityObj = AOD->GetHeader()->GetCentralityP();
+    MultipOrCent = centralityObj->GetCentralityPercentileUnchecked("V0M");
+    AliInfo(Form("Centrality is %f", MultipOrCent));
+  }
+       
+  AliAODVertex *vtx = AOD->GetPrimaryVertex();
+  Double_t zvertex = vtx->GetZ(); // zvertex
+
+  AliEventPool *pool = fCorrelator->GetPool();
+
+  ((TH2D*)fControlObjects->FindObject("NofPoolBinCalls"))->Fill(MultipOrCent,zvertex); // number of calls of pool
+  ((TH2D*)fControlObjects->FindObject("EventProps"))->Fill(MultipOrCent,zvertex); // event properties
+  ((TH3D*)fControlObjects->FindObject("EventsPerPoolBin"))->Fill(MultipOrCent,zvertex,pool->NTracksInPool()); // number of events in the pool
+  ((TH3D*)fControlObjects->FindObject("NofTracksPerPoolBin"))->Fill(MultipOrCent,zvertex,pool->GetCurrentNEvents()); // number of calls of pool
+}
+       
index 7e49630..70f8528 100644 (file)
@@ -236,6 +236,14 @@ void AddTaskDxHFECorrelation(Bool_t bUseMC=kFALSE, TString analysisName="PWGHFco
   return;
 }
 
+// Note: AliHFAssociatedTrackCuts keeps an instance of the external
+// pointer, the arrays thus need to be global
+// TODO: try a proper implementation of AliHFAssociatedTrackCuts later
+Double_t MBins[]={0,20,40,60,80,500};
+Double_t * MultiplicityBins = MBins;
+Double_t ZBins[]={-10,-5,-2.5,2.5,5,10};
+Double_t *ZVrtxBins = ZBins;
+
 AliAnalysisCuts* createDefaultPoolConfig()
 {
   AliHFAssociatedTrackCuts* HFCorrelationCuts=new AliHFAssociatedTrackCuts();
@@ -246,11 +254,7 @@ AliAnalysisCuts* createDefaultPoolConfig()
   HFCorrelationCuts->SetMaxNEventsInPool(200);
   HFCorrelationCuts->SetMinNTracksInPool(100);
   HFCorrelationCuts->SetMinEventsToMix(8);
-  HFCorrelationCuts->SetNofPoolBins(5,5);
-  Double_t MBins[]={0,20,40,60,80,500};
-  Double_t * MultiplicityBins = MBins;
-  Double_t ZBins[]={-10,-5,-2.5,2.5,5,10};
-  Double_t *ZVrtxBins = ZBins;
+  HFCorrelationCuts->SetNofPoolBins(5,5); // Note: the arrays have dimension x+1
   HFCorrelationCuts->SetPoolBins(ZVrtxBins,MultiplicityBins);
 
   TString description = "Info on Pool for EventMixing";