Added possibility to request SPD kFirst for candidates below a certain pt (Zaida...
authordainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 12 Nov 2011 00:15:29 +0000 (00:15 +0000)
committerdainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 12 Nov 2011 00:15:29 +0000 (00:15 +0000)
PWG3/vertexingHF/AliRDHFCuts.cxx
PWG3/vertexingHF/AliRDHFCuts.h

index 0b7d4ad..befa66e 100644 (file)
-/**************************************************************************
- * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
- *                                                                        *
- * Author: The ALICE Off-line Project.                                    *
- * Contributors are mentioned in the code where appropriate.              *
- *                                                                        *
- * Permission to use, copy, modify and distribute this software and its   *
- * documentation strictly for non-commercial purposes is hereby granted   *
- * without fee, provided that the above copyright notice appears in all   *
- * copies and that both the copyright notice and this permission notice   *
- * appear in the supporting documentation. The authors make no claims     *
- * about the suitability of this software for any purpose. It is          *
- * provided "as is" without express or implied warranty.                  *
- **************************************************************************/
-
-/* $Id$ */
-
-/////////////////////////////////////////////////////////////
-//
-// Base class for cuts on AOD reconstructed heavy-flavour decay
-//
-// Author: A.Dainese, andrea.dainese@pd.infn.it
-/////////////////////////////////////////////////////////////
-#include <Riostream.h>
-
-#include "AliVEvent.h"
-#include "AliESDEvent.h"
-#include "AliAODEvent.h"
-#include "AliVVertex.h"
-#include "AliESDVertex.h"
-#include "AliLog.h"
-#include "AliAODVertex.h"
-#include "AliESDtrack.h"
-#include "AliAODTrack.h"
-#include "AliESDtrackCuts.h"
-#include "AliCentrality.h"
-#include "AliAODRecoDecayHF.h"
-#include "AliAnalysisVertexingHF.h"
-#include "AliAODMCHeader.h"
-#include "AliAODMCParticle.h"
-#include "AliRDHFCuts.h"
-#include "AliAnalysisManager.h"
-#include "AliInputEventHandler.h"
-#include "AliPIDResponse.h"
-
-ClassImp(AliRDHFCuts)
-
-//--------------------------------------------------------------------------
-AliRDHFCuts::AliRDHFCuts(const Char_t* name, const Char_t* title) : 
-AliAnalysisCuts(name,title),
-fMinVtxType(3),
-fMinVtxContr(1),
-fMaxVtxRedChi2(1e6),
-fMaxVtxZ(10.),
-fMinSPDMultiplicity(0),
-fTriggerMask(0),
-fTriggerClass("CINT1"),
-fTrackCuts(0),
-fnPtBins(1),
-fnPtBinLimits(1),
-fPtBinLimits(0),
-fnVars(1),
-fVarNames(0),
-fnVarsForOpt(0),
-fVarsForOpt(0),
-fGlobalIndex(1),
-fCutsRD(0),
-fIsUpperCut(0),
-fUsePID(kFALSE),
-fUseAOD049(kFALSE),
-fPidHF(0),
-fWhyRejection(0),
-fEvRejectionBits(0),
-fRemoveDaughtersFromPrimary(kFALSE),
-fUseMCVertex(kFALSE),
-fUsePhysicsSelection(kTRUE),
-fOptPileup(0),
-fMinContrPileup(3),
-fMinDzPileup(0.6),
-fUseCentrality(0),
-fMinCentrality(0.),
-fMaxCentrality(100.),
-fFixRefs(kFALSE),
-fIsSelectedCuts(0),
-fIsSelectedPID(0),
-fMinPtCand(-1.),
-fMaxPtCand(100000.),
-fKeepSignalMC(kFALSE)
-{
-  //
-  // Default Constructor
-  //
-}
-//--------------------------------------------------------------------------
-AliRDHFCuts::AliRDHFCuts(const AliRDHFCuts &source) :
-  AliAnalysisCuts(source),
-  fMinVtxType(source.fMinVtxType),
-  fMinVtxContr(source.fMinVtxContr),
-  fMaxVtxRedChi2(source.fMaxVtxRedChi2),
-  fMaxVtxZ(source.fMaxVtxZ),
-  fMinSPDMultiplicity(source.fMinSPDMultiplicity),
-  fTriggerMask(source.fTriggerMask),
-  fTriggerClass(source.fTriggerClass),
-  fTrackCuts(0),
-  fnPtBins(source.fnPtBins),
-  fnPtBinLimits(source.fnPtBinLimits),
-  fPtBinLimits(0),
-  fnVars(source.fnVars),
-  fVarNames(0),
-  fnVarsForOpt(source.fnVarsForOpt),
-  fVarsForOpt(0),
-  fGlobalIndex(source.fGlobalIndex),
-  fCutsRD(0),
-  fIsUpperCut(0),
-  fUsePID(source.fUsePID),
-  fUseAOD049(source.fUseAOD049),
-  fPidHF(0),
-  fWhyRejection(source.fWhyRejection),
-  fEvRejectionBits(source.fEvRejectionBits),
-  fRemoveDaughtersFromPrimary(source.fRemoveDaughtersFromPrimary),
-  fUseMCVertex(source.fUseMCVertex),
-  fUsePhysicsSelection(source.fUsePhysicsSelection),
-  fOptPileup(source.fOptPileup),
-  fMinContrPileup(source.fMinContrPileup),
-  fMinDzPileup(source.fMinDzPileup),
-  fUseCentrality(source.fUseCentrality),
-  fMinCentrality(source.fMinCentrality),
-  fMaxCentrality(source.fMaxCentrality),
-  fFixRefs(source.fFixRefs),
-  fIsSelectedCuts(source.fIsSelectedCuts),
-  fIsSelectedPID(source.fIsSelectedPID),
-  fMinPtCand(source.fMinPtCand),
-  fMaxPtCand(source.fMaxPtCand),
-  fKeepSignalMC(source.fKeepSignalMC)
-{
-  //
-  // Copy constructor
-  //
-  cout<<"Copy constructor"<<endl;
-  if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
-  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
-  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
-  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
-  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
-  if(source.fPidHF) SetPidHF(source.fPidHF);
-  PrintAll();
-
-}
-//--------------------------------------------------------------------------
-AliRDHFCuts &AliRDHFCuts::operator=(const AliRDHFCuts &source)
-{
-  //
-  // assignment operator
-  //
-  if(&source == this) return *this;
-
-  AliAnalysisCuts::operator=(source);
-
-  fMinVtxType=source.fMinVtxType;
-  fMinVtxContr=source.fMinVtxContr;
-  fMaxVtxRedChi2=source.fMaxVtxRedChi2;
-  fMaxVtxZ=source.fMaxVtxZ;
-  fMinSPDMultiplicity=source.fMinSPDMultiplicity;
-  fTriggerMask=source.fTriggerMask;
-  fTriggerClass=source.fTriggerClass;
-  fnPtBins=source.fnPtBins;
-  fnVars=source.fnVars;
-  fGlobalIndex=source.fGlobalIndex;
-  fnVarsForOpt=source.fnVarsForOpt;
-  fUsePID=source.fUsePID;
-  fUseAOD049=source.fUseAOD049;
-  SetPidHF(source.GetPidHF());
-  fWhyRejection=source.fWhyRejection;
-  fEvRejectionBits=source.fEvRejectionBits;
-  fRemoveDaughtersFromPrimary=source.fRemoveDaughtersFromPrimary;
-  fUseMCVertex=source.fUseMCVertex;
-  fUsePhysicsSelection=source.fUsePhysicsSelection;
-  fOptPileup=source.fOptPileup;
-  fMinContrPileup=source.fMinContrPileup;
-  fMinDzPileup=source.fMinDzPileup;
-  fUseCentrality=source.fUseCentrality;
-  fMinCentrality=source.fMinCentrality;
-  fMaxCentrality=source.fMaxCentrality;
-  fFixRefs=source.fFixRefs;
-  fIsSelectedCuts=source.fIsSelectedCuts;
-  fIsSelectedPID=source.fIsSelectedPID;
-  fMinPtCand=source.fMinPtCand;
-  fMaxPtCand=source.fMaxPtCand;
-  fKeepSignalMC=source.fKeepSignalMC;
-
-  if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
-  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
-  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
-  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
-  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
-  PrintAll();
-
-  return *this;
-}
-//--------------------------------------------------------------------------
-AliRDHFCuts::~AliRDHFCuts() {
-  //  
-  // Default Destructor
-  //
-  if(fPtBinLimits) {delete [] fPtBinLimits; fPtBinLimits=0;}
-  if(fVarNames) {delete [] fVarNames; fVarNames=0;}
-  if(fVarsForOpt) {delete [] fVarsForOpt; fVarsForOpt=0;}
-  if(fCutsRD) {
-    delete [] fCutsRD;
-    fCutsRD=0;
-  }
-  if(fIsUpperCut) {delete [] fIsUpperCut; fIsUpperCut=0;}
-  if(fPidHF){ 
-    delete fPidHF;
-    fPidHF=0;
-  }
-}
-//---------------------------------------------------------------------------
-Int_t AliRDHFCuts::IsEventSelectedInCentrality(AliVEvent *event) {
-  //
-  // Centrality selection
-  //
-  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){    
-    AliWarning("Centrality estimator not valid");    
-    return 3;  
-  }else{    
-    Float_t centvalue=GetCentrality((AliAODEvent*)event);          
-    if (centvalue<-998.){//-999 if no centralityP
-      return 0;    
-    }else{      
-      if (centvalue<fMinCentrality || centvalue>fMaxCentrality){
-       return 2;      
-      }    
-    } 
-  }  
-  return 0;
-}
-//---------------------------------------------------------------------------
-Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) {
-  //
-  // Event selection
-  // 
-  //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;
-
-  fWhyRejection=0;
-  fEvRejectionBits=0;
-  Bool_t accept=kTRUE;
-
-  // check if it's MC
-  Bool_t isMC=kFALSE;
-  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
-  if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}
-
-  // settings for the TPC dE/dx BB parameterization
-  if(fPidHF && fPidHF->GetOldPid()) {
-    // pp, from LHC10d onwards
-    if((event->GetRunNumber()>121693 && event->GetRunNumber()<136851) ||
-       event->GetRunNumber()>139517) fPidHF->SetOnePad(kTRUE);
-    // pp, 2011 low energy run
-    if((event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860)){
-      fPidHF->SetppLowEn2011(kTRUE);
-      fPidHF->SetOnePad(kFALSE);
-    }
-    // PbPb LHC10h
-    if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517) fPidHF->SetPbPb(kTRUE);
-    // MC
-    if(isMC) fPidHF->SetMC(kTRUE);
-    if(isMC && (event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860))
-      fPidHF->SetMClowenpp2011(kTRUE);
-  } 
-  else if(fPidHF && !fPidHF->GetOldPid()) {
-   if(fPidHF->GetPidResponse()==0x0){
-    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
-    AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
-    AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
-    fPidHF->SetPidResponse(pidResp);
-   }
-  }
-
-
-  // trigger class
-  TString firedTriggerClasses=((AliAODEvent*)event)->GetFiredTriggerClasses();
-  // don't do for MC and for PbPb 2010 data
-  if(!isMC && (event->GetRunNumber()<136851 || event->GetRunNumber()>139517)) {
-    if(!firedTriggerClasses.Contains(fTriggerClass.Data())) {
-      fWhyRejection=5;
-      fEvRejectionBits+=1<<kNotSelTrigger;
-      accept=kFALSE;
-    }
-  }
-
-  // TEMPORARY FIX FOR REFERENCES
-  // Fix references to daughter tracks
-  //  if(fFixRefs) {
-  //    AliAnalysisVertexingHF *fixer = new AliAnalysisVertexingHF();
-  //    fixer->FixReferences((AliAODEvent*)event);
-  //    delete fixer;
-  //  }
-  //
-
-
-  // physics selection requirements
-  if(fUsePhysicsSelection){
-    Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
-    if(!isSelected) {
-      if(accept) fWhyRejection=7;
-      fEvRejectionBits+=1<<kPhysicsSelection;
-      accept=kFALSE;
-    }
-  }
-
-  // vertex requirements
-   
-  const AliVVertex *vertex = event->GetPrimaryVertex();
-
-  if(!vertex){
-    accept=kFALSE;
-    fEvRejectionBits+=1<<kNoVertex;
-  }else{
-    TString title=vertex->GetTitle();
-    if(title.Contains("Z") && fMinVtxType>1){
-      accept=kFALSE;
-      fEvRejectionBits+=1<<kNoVertex;
-    }
-    else if(title.Contains("3D") && fMinVtxType>2){
-      accept=kFALSE;
-      fEvRejectionBits+=1<<kNoVertex;
-    }
-    if(vertex->GetNContributors()<fMinVtxContr){
-      accept=kFALSE;
-      fEvRejectionBits+=1<<kTooFewVtxContrib;
-    }
-    if(TMath::Abs(vertex->GetZ())>fMaxVtxZ) {
-      fEvRejectionBits+=1<<kZVtxOutFid;
-      if(accept) fWhyRejection=6;
-      accept=kFALSE;
-    } 
-  }
-
-
-  // pile-up rejection
-  if(fOptPileup==kRejectPileupEvent){
-    Int_t cutc=(Int_t)fMinContrPileup;
-    Double_t cutz=(Double_t)fMinDzPileup;
-    if(event->IsPileupFromSPD(cutc,cutz,3.,2.,10.)) {
-      if(accept) fWhyRejection=1;
-      fEvRejectionBits+=1<<kPileupSPD;
-      accept=kFALSE;
-    }
-  }
-
-  // centrality selection
-  if (fUseCentrality!=kCentOff) {  
-    Int_t rejection=IsEventSelectedInCentrality(event);    
-    if(rejection>1){      
-      if(accept) fWhyRejection=rejection;      
-      fEvRejectionBits+=1<<kOutsideCentrality;
-      accept=kFALSE;
-    }
-  }
-
-  return accept;
-}
-//---------------------------------------------------------------------------
-Bool_t AliRDHFCuts::AreDaughtersSelected(AliAODRecoDecayHF *d) const {
-  //
-  // Daughter tracks selection
-  // 
-  if(!fTrackCuts) return kTRUE;
-  Int_t ndaughters = d->GetNDaughters();
-  AliAODVertex *vAOD = d->GetPrimaryVtx();
-  Double_t pos[3],cov[6];
-  vAOD->GetXYZ(pos);
-  vAOD->GetCovarianceMatrix(cov);
-  const AliESDVertex vESD(pos,cov,100.,100);
-  
-  Bool_t retval=kTRUE;
-  
-  for(Int_t idg=0; idg<ndaughters; idg++) {
-    AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(idg);
-    if(!dgTrack) {retval = kFALSE; continue;}
-    //printf("charge %d\n",dgTrack->Charge());
-    if(dgTrack->Charge()==0) continue; // it's not a track, but a V0
-
-    if(!IsDaughterSelected(dgTrack,&vESD,fTrackCuts)) retval = kFALSE;
-  }
-  
-  return retval;
-}
-//---------------------------------------------------------------------------
-Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts) const {
-  //
-  // Convert to ESDtrack, relate to vertex and check cuts
-  //
-  if(!cuts) return kTRUE;
-
-  Bool_t retval=kTRUE;
-
-  // convert to ESD track here
-  AliESDtrack esdTrack(track);
-  // set the TPC cluster info
-  esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());
-  esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());
-  esdTrack.SetTPCPointsF(track->GetTPCNclsF());
-  // needed to calculate the impact parameters
-  esdTrack.RelateToVertex(primary,0.,3.); 
-
-  if(!cuts->IsSelected(&esdTrack)) retval = kFALSE;
-  if(fOptPileup==kRejectTracksFromPileupVertex){
-    // to be implemented
-    // we need either to have here the AOD Event, 
-    // or to have the pileup vertex object
-  }
-  return retval; 
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits) {
-  // Set the pt bins
-
-  if(fPtBinLimits) {
-    delete [] fPtBinLimits;
-    fPtBinLimits = NULL;
-    printf("Changing the pt bins\n");
-  }
-
-  if(nPtBinLimits != fnPtBins+1){
-    cout<<"Warning: ptBinLimits dimention "<<nPtBinLimits<<" != nPtBins+1 ("<<fnPtBins+1<<")\nSetting nPtBins to "<<nPtBinLimits-1<<endl;
-    SetNPtBins(nPtBinLimits-1);
-  }
-
-  fnPtBinLimits = nPtBinLimits;
-  SetGlobalIndex();
-  //cout<<"Changing also Global Index -> "<<fGlobalIndex<<endl;
-  fPtBinLimits = new Float_t[fnPtBinLimits];
-  for(Int_t ib=0; ib<nPtBinLimits; ib++) fPtBinLimits[ib]=ptBinLimits[ib];
-
-  return;
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut){
-  // Set the variable names
-
-  if(fVarNames) {
-    delete [] fVarNames;
-    fVarNames = NULL;
-    //printf("Changing the variable names\n");
-  }
-  if(nVars!=fnVars){
-    printf("Wrong number of variables: it has to be %d\n",fnVars);
-    return;
-  }
-  //fnVars=nVars;
-  fVarNames = new TString[nVars];
-  fIsUpperCut = new Bool_t[nVars];
-  for(Int_t iv=0; iv<nVars; iv++) {
-    fVarNames[iv] = varNames[iv];
-    fIsUpperCut[iv] = isUpperCut[iv];
-  }
-
-  return;
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::SetVarsForOpt(Int_t nVars,Bool_t *forOpt) {
-  // Set the variables to be used for cuts optimization
-
-  if(fVarsForOpt) {
-    delete [] fVarsForOpt;
-    fVarsForOpt = NULL;
-    //printf("Changing the variables for cut optimization\n");
-  }
-  
-  if(nVars==0){//!=fnVars) {
-    printf("%d not accepted as number of variables: it has to be %d\n",nVars,fnVars);
-    return;
-  } 
-  
-  fnVarsForOpt = 0;
-  fVarsForOpt = new Bool_t[fnVars];
-  for(Int_t iv=0; iv<fnVars; iv++) {
-    fVarsForOpt[iv]=forOpt[iv];
-    if(fVarsForOpt[iv]) fnVarsForOpt++;
-  }
-
-  return;
-}
-
-//---------------------------------------------------------------------------
-void AliRDHFCuts::SetUseCentrality(Int_t flag) {
-  //
-  // set centrality estimator  
-  //
-  fUseCentrality=flag;
-  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");
-  return;
-}
-
-
-//---------------------------------------------------------------------------
-void AliRDHFCuts::SetCuts(Int_t nVars,Int_t nPtBins,Float_t **cutsRD) {
-  //
-  // store the cuts
-  //
-  if(nVars!=fnVars) {
-    printf("Wrong number of variables: it has to be %d\n",fnVars);
-    AliFatal("exiting");
-  } 
-  if(nPtBins!=fnPtBins) {
-    printf("Wrong number of pt bins: it has to be %d\n",fnPtBins);
-    AliFatal("exiting");
-  } 
-
-  if(!fCutsRD)  fCutsRD = new Float_t[fGlobalIndex];
-  
-
-  for(Int_t iv=0; iv<fnVars; iv++) {
-
-    for(Int_t ib=0; ib<fnPtBins; ib++) {
-
-      //check
-      if(GetGlobalIndex(iv,ib)>=fGlobalIndex) {
-       cout<<"Overflow, exit..."<<endl;
-       return;
-      }
-
-      fCutsRD[GetGlobalIndex(iv,ib)] = cutsRD[iv][ib];
-
-    }
-  }
-  return;
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::SetCuts(Int_t glIndex,Float_t* cutsRDGlob){
-  //
-  // store the cuts
-  //
-  if(glIndex != fGlobalIndex){
-    cout<<"Wrong array size: it has to be "<<fGlobalIndex<<endl;
-    AliFatal("exiting");
-  }
-  if(!fCutsRD)  fCutsRD = new Float_t[fGlobalIndex];
-
-  for(Int_t iGl=0;iGl<fGlobalIndex;iGl++){
-    fCutsRD[iGl] = cutsRDGlob[iGl];
-  }
-  return;
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::PrintAll() const {
-  //
-  // print all cuts values
-  // 
-
-  printf("Minimum vtx type %d\n",fMinVtxType);
-  printf("Minimum vtx contr %d\n",fMinVtxContr);
-  printf("Max vtx red chi2 %f\n",fMaxVtxRedChi2);
-  printf("Min SPD mult %d\n",fMinSPDMultiplicity);
-  printf("Use PID %d\n",(Int_t)fUsePID);
-  printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);
-  printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");
-  if(fOptPileup==1) printf(" -- Reject pileup event");
-  if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");
-  if(fUseCentrality>0) {
-    TString estimator="";
-    if(fUseCentrality==1) estimator = "V0";
-    if(fUseCentrality==2) estimator = "Tracks";
-    if(fUseCentrality==3) estimator = "Tracklets";
-    if(fUseCentrality==4) estimator = "SPD clusters outer"; 
-    printf("Centrality class considered: %.1f-%.1f, estimated with %s",fMinCentrality,fMaxCentrality,estimator.Data());
-  }
-
-  if(fVarNames){
-    cout<<"Array of variables"<<endl;
-    for(Int_t iv=0;iv<fnVars;iv++){
-      cout<<fVarNames[iv]<<"\t";
-    }
-    cout<<endl;
-  }
-  if(fVarsForOpt){
-    cout<<"Array of optimization"<<endl;
-    for(Int_t iv=0;iv<fnVars;iv++){
-      cout<<fVarsForOpt[iv]<<"\t";
-    }
-    cout<<endl;
-  }
-  if(fIsUpperCut){
-    cout<<"Array of upper/lower cut"<<endl;
-   for(Int_t iv=0;iv<fnVars;iv++){
-     cout<<fIsUpperCut[iv]<<"\t";
-   }
-   cout<<endl;
-  }
-  if(fPtBinLimits){
-    cout<<"Array of ptbin limits"<<endl;
-    for(Int_t ib=0;ib<fnPtBinLimits;ib++){
-      cout<<fPtBinLimits[ib]<<"\t";
-    }
-    cout<<endl;
-  }
-  if(fCutsRD){
-    cout<<"Matrix of cuts"<<endl;
-   for(Int_t iv=0;iv<fnVars;iv++){
-     for(Int_t ib=0;ib<fnPtBins;ib++){
-       cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"\t";
-     } 
-     cout<<endl;
-   }
-   cout<<endl;
-  }
-  return;
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::GetCuts(Float_t**& cutsRD) const{
-  //
-  // get the cuts
-  //
-
-  //cout<<"Give back a "<<fnVars<<"x"<<fnPtBins<<" matrix."<<endl;
-
-
-  Int_t iv,ib;
-  if(!cutsRD) {
-    //cout<<"Initialization..."<<endl;
-    cutsRD=new Float_t*[fnVars];
-    for(iv=0; iv<fnVars; iv++) {
-      cutsRD[iv] = new Float_t[fnPtBins];
-    }
-  }
-  
-  for(Int_t iGlobal=0; iGlobal<fGlobalIndex; iGlobal++) {
-    GetVarPtIndex(iGlobal,iv,ib);
-    cutsRD[iv][ib] = fCutsRD[iGlobal];
-  }
-
-  return;
-}
-
-//---------------------------------------------------------------------------
-Int_t AliRDHFCuts::GetGlobalIndex(Int_t iVar,Int_t iPtBin) const{
-  //
-  // give the global index from variable and pt bin
-  //
-  return iPtBin*fnVars+iVar;
-}
-
-//---------------------------------------------------------------------------
-void AliRDHFCuts::GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const {
-  //
-  //give the index of the variable and of the pt bin from the global index
-  //
-  iPtBin=(Int_t)iGlob/fnVars;
-  iVar=iGlob%fnVars;
-
-  return;
-}
-
-//---------------------------------------------------------------------------
-Int_t AliRDHFCuts::PtBin(Double_t pt) const {
-  //
-  //give the pt bin where the pt lies.
-  //
-  Int_t ptbin=-1;
-  if(pt<fPtBinLimits[0])return ptbin;
-  for (Int_t i=0;i<fnPtBins;i++){
-    if(pt<fPtBinLimits[i+1]) {
-      ptbin=i;
-      break;
-    }
-  }
-  return ptbin;
-}
-//-------------------------------------------------------------------
-Float_t AliRDHFCuts::GetCutValue(Int_t iVar,Int_t iPtBin) const {
-  // 
-  // Give the value of cut set for the variable iVar and the pt bin iPtBin
-  //
-  if(!fCutsRD){
-    cout<<"Cuts not iniziaisez yet"<<endl;
-    return 0;
-  }
-  return fCutsRD[GetGlobalIndex(iVar,iPtBin)];
-}
-//-------------------------------------------------------------------
-Float_t AliRDHFCuts::GetCentrality(AliAODEvent* aodEvent,AliRDHFCuts::ECentrality estimator) {
-  //
-  // Get centrality percentile
-  //
-
-  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aodEvent)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
-  if(mcArray) {fUseAOD049=kFALSE;}
-
-  AliAODHeader *header=aodEvent->GetHeader();
-  AliCentrality *centrality=header->GetCentralityP();
-  Float_t cent=-999.;
-  Bool_t isSelRun=kFALSE;
-  Int_t selRun[5]={138364, 138826, 138828, 138836, 138871};
-  if(!centrality) return cent;
-  else{
-    if (estimator==kCentV0M){
-      cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));
-      if(cent<0){
-       Int_t quality = centrality->GetQuality();
-       if(quality<=1){
-         cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
-       }else{
-         Int_t runnum=aodEvent->GetRunNumber();
-         for(Int_t ir=0;ir<5;ir++){
-           if(runnum==selRun[ir]){
-             isSelRun=kTRUE;
-             break;
-           }
-         }
-         if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
-       }
-      }
-
-      //temporary fix for AOD049 outliers
-      if(fUseAOD049&&cent>=0){
-       Float_t v0=0;
-       AliAODVZERO* aodV0 = aodEvent->GetVZEROData();
-       v0+=aodV0->GetMTotV0A();
-       v0+=aodV0->GetMTotV0C();
-       if(cent==0&&v0<19500)return -1;//filtering issue
-       Float_t tkl = (Float_t)(aodEvent->GetTracklets()->GetNumberOfTracklets());
-       Float_t val= 1.30552 +  0.147931 * v0;
-       Float_t tklSigma[101]={176.644, 156.401, 153.789, 153.015, 142.476, 137.951, 136.127, 129.852, 127.436, 124.86, 120.788, 115.611, 113.172, 110.496, 109.127, 104.421, 102.479, 99.9766, 97.5152, 94.0654, 92.4602, 89.3364, 87.1342, 83.3497, 82.6216, 81.1084, 78.0793, 76.1234, 72.9434, 72.1334, 68.0056, 68.2755, 66.0376, 62.9666, 62.4274, 59.65, 58.3776, 56.6361, 54.5184, 53.4224, 51.932, 50.8922, 48.2848, 47.912, 46.5717, 43.4114, 43.2083, 41.3065, 40.1863, 38.5255, 37.2851, 37.5396, 34.4949, 33.8366, 31.8043, 31.7412, 30.8392, 30.0274, 28.8793, 27.6398, 26.6488, 25.0183, 25.1489, 24.4185, 22.9107, 21.2002, 21.6977, 20.1242, 20.4963, 19.0235, 19.298, 17.4103, 16.868, 15.2939, 15.2939, 16.0295, 14.186, 14.186, 15.2173, 12.9504, 12.9504, 12.9504, 15.264, 12.3674, 12.3674, 12.3674, 12.3674, 12.3674, 18.3811, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544};
-       if ( TMath::Abs(tkl-val) > 6.*tklSigma[(Int_t)cent] )return -1;//outlier
-      }
-    }
-    else {
-      if (estimator==kCentTRK) {
-       cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));
-       if(cent<0){
-         Int_t quality = centrality->GetQuality();
-         if(quality<=1){
-           cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
-         }else{
-           Int_t runnum=aodEvent->GetRunNumber();
-           for(Int_t ir=0;ir<5;ir++){
-             if(runnum==selRun[ir]){
-               isSelRun=kTRUE;
-               break;
-             }
-           }
-           if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
-         }
-       }
-      }
-      else{
-       if (estimator==kCentTKL){
-         cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));
-         if(cent<0){
-           Int_t quality = centrality->GetQuality();
-           if(quality<=1){
-             cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
-           }else{
-             Int_t runnum=aodEvent->GetRunNumber();
-             for(Int_t ir=0;ir<5;ir++){
-               if(runnum==selRun[ir]){
-                 isSelRun=kTRUE;
-                 break;
-           }
-             }
-             if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
-           }   
-         }
-       }
-       else{
-         if (estimator==kCentCL1){
-           cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));
-           if(cent<0){
-             Int_t quality = centrality->GetQuality();
-             if(quality<=1){
-               cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
-             }else{
-               Int_t runnum=aodEvent->GetRunNumber();
-               for(Int_t ir=0;ir<5;ir++){
-                 if(runnum==selRun[ir]){
-                   isSelRun=kTRUE;
-                   break;
-                 }
-               }
-               if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
-             }
-           }
-         }
-         else {
-           AliWarning("Centrality estimator not valid");
-           
-         }
-       }
-      }
-    } 
-  }
-  return cent;
-}
-//-------------------------------------------------------------------
-Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {
-  //
-  // Compare two cuts objects
-  //
-
-  Bool_t areEqual=kTRUE;
-
-  if(fMinVtxType!=obj->fMinVtxType) { printf("Minimum vtx type %d  %d\n",fMinVtxType,obj->fMinVtxType); areEqual=kFALSE;}
-
-  if(fMinVtxContr!=obj->fMinVtxContr) { printf("Minimum vtx contr %d  %d\n",fMinVtxContr,obj->fMinVtxContr); areEqual=kFALSE;}
-
-  if(TMath::Abs(fMaxVtxRedChi2-obj->fMaxVtxRedChi2)>1.e-10) {   printf("Max vtx red chi2 %f  %f\n",fMaxVtxRedChi2,obj->fMaxVtxRedChi2);areEqual=kFALSE;}
-
-  if(fMinSPDMultiplicity!=obj->fMinSPDMultiplicity) {  printf("Min SPD mult %d\n  %d",fMinSPDMultiplicity,obj->fMinSPDMultiplicity);areEqual=kFALSE;}
-
-  if(fUsePID!=obj->fUsePID) { printf("Use PID %d  %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}
-
-  if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d  %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}
-  if(fTrackCuts){
-    if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d  %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
-
-    if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d  %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
-
-    if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f  %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
-
-    if(fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)!=obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)) {printf("ClusterReq SPD %d  %d\n",fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD),obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)); areEqual=kFALSE;}
-  }
-
-  if(fCutsRD) {
-   for(Int_t iv=0;iv<fnVars;iv++) {
-     for(Int_t ib=0;ib<fnPtBins;ib++) {
-       if(TMath::Abs(fCutsRD[GetGlobalIndex(iv,ib)]-obj->fCutsRD[GetGlobalIndex(iv,ib)])>1.e-10) {
-        cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"   "<<obj->fCutsRD[GetGlobalIndex(iv,ib)]<<"\n";
-        areEqual=kFALSE;
-       }
-     }
-   }
-  }
-
-  return areEqual;
-}
-//---------------------------------------------------------------------------
-void AliRDHFCuts::MakeTable() const {
-  //
-  // print cuts values in table format
-  // 
-
-       TString ptString = "pT range";
-       if(fVarNames && fPtBinLimits && fCutsRD){
-               TString firstLine(Form("*       %-15s",ptString.Data()));
-               for (Int_t ivar=0; ivar<fnVars; ivar++){
-                       firstLine+=Form("*    %-15s  ",fVarNames[ivar].Data());
-                       if (ivar == fnVars){
-                               firstLine+="*\n";
-                       }
-               }
-               Printf("%s",firstLine.Data());
-               
-               for (Int_t ipt=0; ipt<fnPtBins; ipt++){
-                       TString line;
-                       if (ipt==fnPtBins-1){
-                               line=Form("*  %5.1f < pt < inf    ",fPtBinLimits[ipt]);
-                       }
-                       else{
-                               line=Form("*  %5.1f < pt < %4.1f   ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);
-                       }
-                       for (Int_t ivar=0; ivar<fnVars; ivar++){
-                               line+=Form("*     %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);
-                       }
-                       Printf("%s",line.Data());
-               }
-
-       }
-
-
-  return;
-}
-//--------------------------------------------------------------------------
-Bool_t AliRDHFCuts::RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d,
-                                       AliAODEvent *aod) const
-{
-  //
-  // Recalculate primary vertex without daughters
-  //
-
-  if(!aod) {
-    AliError("Can not remove daughters from vertex without AOD event");
-    return 0;
-  }   
-
-  AliAODVertex *recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
-  if(!recvtx){
-    AliDebug(2,"Removal of daughter tracks failed");
-    return kFALSE;
-  }
-
-
-  //set recalculed primary vertex
-  d->SetOwnPrimaryVtx(recvtx);
-  delete recvtx;
-
-  return kTRUE;
-}
-//--------------------------------------------------------------------------
-Bool_t AliRDHFCuts::SetMCPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const
-{
-  //
-  // Recalculate primary vertex without daughters
-  //
-
-  if(!aod) {
-    AliError("Can not get MC vertex without AOD event");
-    return kFALSE;
-  }   
-
-  // load MC header
-  AliAODMCHeader *mcHeader = 
-    (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
-  if(!mcHeader) {
-    AliError("Can not get MC vertex without AODMCHeader event");
-    return kFALSE;
-  }
-  Double_t pos[3];
-  Double_t covmatrix[6]={0.,0.,0.,0.,0.,0.};
-  mcHeader->GetVertex(pos);
-  AliAODVertex *recvtx=new AliAODVertex(pos,covmatrix);
-
-  if(!recvtx){
-    AliDebug(2,"Removal of daughter tracks failed");
-    return kFALSE;
-  }
-
-  //set recalculed primary vertex
-  d->SetOwnPrimaryVtx(recvtx);
-
-  d->RecalculateImpPars(recvtx,aod);
-
-  delete recvtx;
-
-  return kTRUE;
-}
-//--------------------------------------------------------------------------
-void AliRDHFCuts::CleanOwnPrimaryVtx(AliAODRecoDecayHF *d,
-                                    AliAODEvent *aod,
-                                    AliAODVertex *origownvtx) const
-{
-  //
-  // Clean-up own primary vertex if needed
-  //
-
-  if(fRemoveDaughtersFromPrimary || fUseMCVertex) {
-    d->UnsetOwnPrimaryVtx();
-    if(origownvtx) {
-      d->SetOwnPrimaryVtx(origownvtx);
-      delete origownvtx; origownvtx=NULL;
-    }
-    d->RecalculateImpPars(d->GetPrimaryVtx(),aod);
-  } else {
-    if(origownvtx) {
-      delete origownvtx; origownvtx=NULL;
-    }
-  }
-  return;
-}
-//--------------------------------------------------------------------------
-Bool_t AliRDHFCuts::IsSignalMC(AliAODRecoDecay *d,AliAODEvent *aod,Int_t pdg) const 
-{
-  //
-  // Checks if this candidate is matched to MC signal
-  // 
-
-  if(!aod) return kFALSE;
-
-  // get the MC array
-  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aod)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
-
-  if(!mcArray) return kFALSE;
-
-  // try to match  
-  Int_t label = d->MatchToMC(pdg,mcArray);
-  
-  if(label>=0) {
-    //printf("MATCH!\n");
-    return kTRUE;
-  }
-
-  return kFALSE;
-}
-
-
+/**************************************************************************\r
+ * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+\r
+/* $Id$ */\r
+\r
+/////////////////////////////////////////////////////////////\r
+//\r
+// Base class for cuts on AOD reconstructed heavy-flavour decay\r
+//\r
+// Author: A.Dainese, andrea.dainese@pd.infn.it\r
+/////////////////////////////////////////////////////////////\r
+#include <Riostream.h>\r
+\r
+#include "AliVEvent.h"\r
+#include "AliESDEvent.h"\r
+#include "AliAODEvent.h"\r
+#include "AliVVertex.h"\r
+#include "AliESDVertex.h"\r
+#include "AliLog.h"\r
+#include "AliAODVertex.h"\r
+#include "AliESDtrack.h"\r
+#include "AliAODTrack.h"\r
+#include "AliESDtrackCuts.h"\r
+#include "AliCentrality.h"\r
+#include "AliAODRecoDecayHF.h"\r
+#include "AliAnalysisVertexingHF.h"\r
+#include "AliAODMCHeader.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliRDHFCuts.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliInputEventHandler.h"\r
+#include "AliPIDResponse.h"\r
+\r
+ClassImp(AliRDHFCuts)\r
+\r
+//--------------------------------------------------------------------------\r
+AliRDHFCuts::AliRDHFCuts(const Char_t* name, const Char_t* title) : \r
+AliAnalysisCuts(name,title),\r
+fMinVtxType(3),\r
+fMinVtxContr(1),\r
+fMaxVtxRedChi2(1e6),\r
+fMaxVtxZ(10.),\r
+fMinSPDMultiplicity(0),\r
+fTriggerMask(0),\r
+fTriggerClass("CINT1"),\r
+fTrackCuts(0),\r
+fnPtBins(1),\r
+fnPtBinLimits(1),\r
+fPtBinLimits(0),\r
+fnVars(1),\r
+fVarNames(0),\r
+fnVarsForOpt(0),\r
+fVarsForOpt(0),\r
+fGlobalIndex(1),\r
+fCutsRD(0),\r
+fIsUpperCut(0),\r
+fUsePID(kFALSE),\r
+fUseAOD049(kFALSE),\r
+fPidHF(0),\r
+fWhyRejection(0),\r
+fEvRejectionBits(0),\r
+fRemoveDaughtersFromPrimary(kFALSE),\r
+fUseMCVertex(kFALSE),\r
+fUsePhysicsSelection(kTRUE),\r
+fOptPileup(0),\r
+fMinContrPileup(3),\r
+fMinDzPileup(0.6),\r
+fUseCentrality(0),\r
+fMinCentrality(0.),\r
+fMaxCentrality(100.),\r
+fFixRefs(kFALSE),\r
+fIsSelectedCuts(0),\r
+fIsSelectedPID(0),\r
+fMinPtCand(-1.),\r
+fMaxPtCand(100000.),\r
+fKeepSignalMC(kFALSE),\r
+fIsCandTrackSPDFirst(kFALSE),\r
+fMaxPtCandTrackSPDFirst(0.)\r
+{\r
+  //\r
+  // Default Constructor\r
+  //\r
+}\r
+//--------------------------------------------------------------------------\r
+AliRDHFCuts::AliRDHFCuts(const AliRDHFCuts &source) :\r
+  AliAnalysisCuts(source),\r
+  fMinVtxType(source.fMinVtxType),\r
+  fMinVtxContr(source.fMinVtxContr),\r
+  fMaxVtxRedChi2(source.fMaxVtxRedChi2),\r
+  fMaxVtxZ(source.fMaxVtxZ),\r
+  fMinSPDMultiplicity(source.fMinSPDMultiplicity),\r
+  fTriggerMask(source.fTriggerMask),\r
+  fTriggerClass(source.fTriggerClass),\r
+  fTrackCuts(0),\r
+  fnPtBins(source.fnPtBins),\r
+  fnPtBinLimits(source.fnPtBinLimits),\r
+  fPtBinLimits(0),\r
+  fnVars(source.fnVars),\r
+  fVarNames(0),\r
+  fnVarsForOpt(source.fnVarsForOpt),\r
+  fVarsForOpt(0),\r
+  fGlobalIndex(source.fGlobalIndex),\r
+  fCutsRD(0),\r
+  fIsUpperCut(0),\r
+  fUsePID(source.fUsePID),\r
+  fUseAOD049(source.fUseAOD049),\r
+  fPidHF(0),\r
+  fWhyRejection(source.fWhyRejection),\r
+  fEvRejectionBits(source.fEvRejectionBits),\r
+  fRemoveDaughtersFromPrimary(source.fRemoveDaughtersFromPrimary),\r
+  fUseMCVertex(source.fUseMCVertex),\r
+  fUsePhysicsSelection(source.fUsePhysicsSelection),\r
+  fOptPileup(source.fOptPileup),\r
+  fMinContrPileup(source.fMinContrPileup),\r
+  fMinDzPileup(source.fMinDzPileup),\r
+  fUseCentrality(source.fUseCentrality),\r
+  fMinCentrality(source.fMinCentrality),\r
+  fMaxCentrality(source.fMaxCentrality),\r
+  fFixRefs(source.fFixRefs),\r
+  fIsSelectedCuts(source.fIsSelectedCuts),\r
+  fIsSelectedPID(source.fIsSelectedPID),\r
+  fMinPtCand(source.fMinPtCand),\r
+  fMaxPtCand(source.fMaxPtCand),\r
+  fKeepSignalMC(source.fKeepSignalMC),\r
+  fIsCandTrackSPDFirst(source.fIsCandTrackSPDFirst),\r
+  fMaxPtCandTrackSPDFirst(source.fMaxPtCandTrackSPDFirst)\r
+{\r
+  //\r
+  // Copy constructor\r
+  //\r
+  cout<<"Copy constructor"<<endl;\r
+  if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());\r
+  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);\r
+  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);\r
+  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);\r
+  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);\r
+  if(source.fPidHF) SetPidHF(source.fPidHF);\r
+  PrintAll();\r
+\r
+}\r
+//--------------------------------------------------------------------------\r
+AliRDHFCuts &AliRDHFCuts::operator=(const AliRDHFCuts &source)\r
+{\r
+  //\r
+  // assignment operator\r
+  //\r
+  if(&source == this) return *this;\r
+\r
+  AliAnalysisCuts::operator=(source);\r
+\r
+  fMinVtxType=source.fMinVtxType;\r
+  fMinVtxContr=source.fMinVtxContr;\r
+  fMaxVtxRedChi2=source.fMaxVtxRedChi2;\r
+  fMaxVtxZ=source.fMaxVtxZ;\r
+  fMinSPDMultiplicity=source.fMinSPDMultiplicity;\r
+  fTriggerMask=source.fTriggerMask;\r
+  fTriggerClass=source.fTriggerClass;\r
+  fnPtBins=source.fnPtBins;\r
+  fnVars=source.fnVars;\r
+  fGlobalIndex=source.fGlobalIndex;\r
+  fnVarsForOpt=source.fnVarsForOpt;\r
+  fUsePID=source.fUsePID;\r
+  fUseAOD049=source.fUseAOD049;\r
+  SetPidHF(source.GetPidHF());\r
+  fWhyRejection=source.fWhyRejection;\r
+  fEvRejectionBits=source.fEvRejectionBits;\r
+  fRemoveDaughtersFromPrimary=source.fRemoveDaughtersFromPrimary;\r
+  fUseMCVertex=source.fUseMCVertex;\r
+  fUsePhysicsSelection=source.fUsePhysicsSelection;\r
+  fOptPileup=source.fOptPileup;\r
+  fMinContrPileup=source.fMinContrPileup;\r
+  fMinDzPileup=source.fMinDzPileup;\r
+  fUseCentrality=source.fUseCentrality;\r
+  fMinCentrality=source.fMinCentrality;\r
+  fMaxCentrality=source.fMaxCentrality;\r
+  fFixRefs=source.fFixRefs;\r
+  fIsSelectedCuts=source.fIsSelectedCuts;\r
+  fIsSelectedPID=source.fIsSelectedPID;\r
+  fMinPtCand=source.fMinPtCand;\r
+  fMaxPtCand=source.fMaxPtCand;\r
+  fKeepSignalMC=source.fKeepSignalMC;\r
+  fIsCandTrackSPDFirst=source.fIsCandTrackSPDFirst;\r
+  fMaxPtCandTrackSPDFirst=source.fMaxPtCandTrackSPDFirst;\r
+\r
+  if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());\r
+  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);\r
+  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);\r
+  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);\r
+  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);\r
+  PrintAll();\r
+\r
+  return *this;\r
+}\r
+//--------------------------------------------------------------------------\r
+AliRDHFCuts::~AliRDHFCuts() {\r
+  //  \r
+  // Default Destructor\r
+  //\r
+  if(fPtBinLimits) {delete [] fPtBinLimits; fPtBinLimits=0;}\r
+  if(fVarNames) {delete [] fVarNames; fVarNames=0;}\r
+  if(fVarsForOpt) {delete [] fVarsForOpt; fVarsForOpt=0;}\r
+  if(fCutsRD) {\r
+    delete [] fCutsRD;\r
+    fCutsRD=0;\r
+  }\r
+  if(fIsUpperCut) {delete [] fIsUpperCut; fIsUpperCut=0;}\r
+  if(fPidHF){ \r
+    delete fPidHF;\r
+    fPidHF=0;\r
+  }\r
+}\r
+//---------------------------------------------------------------------------\r
+Int_t AliRDHFCuts::IsEventSelectedInCentrality(AliVEvent *event) {\r
+  //\r
+  // Centrality selection\r
+  //\r
+  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){    \r
+    AliWarning("Centrality estimator not valid");    \r
+    return 3;  \r
+  }else{    \r
+    Float_t centvalue=GetCentrality((AliAODEvent*)event);          \r
+    if (centvalue<-998.){//-999 if no centralityP\r
+      return 0;    \r
+    }else{      \r
+      if (centvalue<fMinCentrality || centvalue>fMaxCentrality){\r
+       return 2;      \r
+      }    \r
+    } \r
+  }  \r
+  return 0;\r
+}\r
+//---------------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) {\r
+  //\r
+  // Event selection\r
+  // \r
+  //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;\r
+\r
+  fWhyRejection=0;\r
+  fEvRejectionBits=0;\r
+  Bool_t accept=kTRUE;\r
+\r
+  // check if it's MC\r
+  Bool_t isMC=kFALSE;\r
+  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());\r
+  if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}\r
+\r
+  // settings for the TPC dE/dx BB parameterization\r
+  if(fPidHF && fPidHF->GetOldPid()) {\r
+    // pp, from LHC10d onwards\r
+    if((event->GetRunNumber()>121693 && event->GetRunNumber()<136851) ||\r
+       event->GetRunNumber()>139517) fPidHF->SetOnePad(kTRUE);\r
+    // pp, 2011 low energy run\r
+    if((event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860)){\r
+      fPidHF->SetppLowEn2011(kTRUE);\r
+      fPidHF->SetOnePad(kFALSE);\r
+    }\r
+    // PbPb LHC10h\r
+    if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517) fPidHF->SetPbPb(kTRUE);\r
+    // MC\r
+    if(isMC) fPidHF->SetMC(kTRUE);\r
+    if(isMC && (event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860))\r
+      fPidHF->SetMClowenpp2011(kTRUE);\r
+  } \r
+  else if(fPidHF && !fPidHF->GetOldPid()) {\r
+   if(fPidHF->GetPidResponse()==0x0){\r
+    AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
+    AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();\r
+    AliPIDResponse *pidResp=inputHandler->GetPIDResponse();\r
+    fPidHF->SetPidResponse(pidResp);\r
+   }\r
+  }\r
+\r
+\r
+  // trigger class\r
+  TString firedTriggerClasses=((AliAODEvent*)event)->GetFiredTriggerClasses();\r
+  // don't do for MC and for PbPb 2010 data\r
+  if(!isMC && (event->GetRunNumber()<136851 || event->GetRunNumber()>139517)) {\r
+    if(!firedTriggerClasses.Contains(fTriggerClass.Data())) {\r
+      fWhyRejection=5;\r
+      fEvRejectionBits+=1<<kNotSelTrigger;\r
+      accept=kFALSE;\r
+    }\r
+  }\r
+\r
+  // TEMPORARY FIX FOR REFERENCES\r
+  // Fix references to daughter tracks\r
+  //  if(fFixRefs) {\r
+  //    AliAnalysisVertexingHF *fixer = new AliAnalysisVertexingHF();\r
+  //    fixer->FixReferences((AliAODEvent*)event);\r
+  //    delete fixer;\r
+  //  }\r
+  //\r
+\r
+\r
+  // physics selection requirements\r
+  if(fUsePhysicsSelection){\r
+    Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);\r
+    if(!isSelected) {\r
+      if(accept) fWhyRejection=7;\r
+      fEvRejectionBits+=1<<kPhysicsSelection;\r
+      accept=kFALSE;\r
+    }\r
+  }\r
+\r
+  // vertex requirements\r
+   \r
+  const AliVVertex *vertex = event->GetPrimaryVertex();\r
+\r
+  if(!vertex){\r
+    accept=kFALSE;\r
+    fEvRejectionBits+=1<<kNoVertex;\r
+  }else{\r
+    TString title=vertex->GetTitle();\r
+    if(title.Contains("Z") && fMinVtxType>1){\r
+      accept=kFALSE;\r
+      fEvRejectionBits+=1<<kNoVertex;\r
+    }\r
+    else if(title.Contains("3D") && fMinVtxType>2){\r
+      accept=kFALSE;\r
+      fEvRejectionBits+=1<<kNoVertex;\r
+    }\r
+    if(vertex->GetNContributors()<fMinVtxContr){\r
+      accept=kFALSE;\r
+      fEvRejectionBits+=1<<kTooFewVtxContrib;\r
+    }\r
+    if(TMath::Abs(vertex->GetZ())>fMaxVtxZ) {\r
+      fEvRejectionBits+=1<<kZVtxOutFid;\r
+      if(accept) fWhyRejection=6;\r
+      accept=kFALSE;\r
+    } \r
+  }\r
+\r
+\r
+  // pile-up rejection\r
+  if(fOptPileup==kRejectPileupEvent){\r
+    Int_t cutc=(Int_t)fMinContrPileup;\r
+    Double_t cutz=(Double_t)fMinDzPileup;\r
+    if(event->IsPileupFromSPD(cutc,cutz,3.,2.,10.)) {\r
+      if(accept) fWhyRejection=1;\r
+      fEvRejectionBits+=1<<kPileupSPD;\r
+      accept=kFALSE;\r
+    }\r
+  }\r
+\r
+  // centrality selection\r
+  if (fUseCentrality!=kCentOff) {  \r
+    Int_t rejection=IsEventSelectedInCentrality(event);    \r
+    if(rejection>1){      \r
+      if(accept) fWhyRejection=rejection;      \r
+      fEvRejectionBits+=1<<kOutsideCentrality;\r
+      accept=kFALSE;\r
+    }\r
+  }\r
+\r
+  return accept;\r
+}\r
+//---------------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::AreDaughtersSelected(AliAODRecoDecayHF *d) const {\r
+  //\r
+  // Daughter tracks selection\r
+  // \r
+  if(!fTrackCuts) return kTRUE;\r
\r
+  Int_t ndaughters = d->GetNDaughters();\r
+  AliAODVertex *vAOD = d->GetPrimaryVtx();\r
+  Double_t pos[3],cov[6];\r
+  vAOD->GetXYZ(pos);\r
+  vAOD->GetCovarianceMatrix(cov);\r
+  const AliESDVertex vESD(pos,cov,100.,100);\r
+  \r
+  Bool_t retval=kTRUE;\r
+  \r
+  for(Int_t idg=0; idg<ndaughters; idg++) {\r
+    AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(idg);\r
+    if(!dgTrack) {retval = kFALSE; continue;}\r
+    //printf("charge %d\n",dgTrack->Charge());\r
+    if(dgTrack->Charge()==0) continue; // it's not a track, but a V0\r
+\r
+    if(fIsCandTrackSPDFirst && d->Pt()<fMaxPtCandTrackSPDFirst)\r
+      { if(!dgTrack->HasPointOnITSLayer(0)) { retval = kFALSE; continue; } }\r
+\r
+    if(!IsDaughterSelected(dgTrack,&vESD,fTrackCuts)) retval = kFALSE;\r
+  }\r
+  \r
+  return retval;\r
+}\r
+//---------------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts) const {\r
+  //\r
+  // Convert to ESDtrack, relate to vertex and check cuts\r
+  //\r
+  if(!cuts) return kTRUE;\r
+\r
+  Bool_t retval=kTRUE;\r
+\r
+  // convert to ESD track here\r
+  AliESDtrack esdTrack(track);\r
+  // set the TPC cluster info\r
+  esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());\r
+  esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());\r
+  esdTrack.SetTPCPointsF(track->GetTPCNclsF());\r
+  // needed to calculate the impact parameters\r
+  esdTrack.RelateToVertex(primary,0.,3.); \r
+\r
+  if(!cuts->IsSelected(&esdTrack)) retval = kFALSE;\r
\r
+  if(fOptPileup==kRejectTracksFromPileupVertex){\r
+    // to be implemented\r
+    // we need either to have here the AOD Event, \r
+    // or to have the pileup vertex object\r
+  }\r
+  return retval; \r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits) {\r
+  // Set the pt bins\r
+\r
+  if(fPtBinLimits) {\r
+    delete [] fPtBinLimits;\r
+    fPtBinLimits = NULL;\r
+    printf("Changing the pt bins\n");\r
+  }\r
+\r
+  if(nPtBinLimits != fnPtBins+1){\r
+    cout<<"Warning: ptBinLimits dimention "<<nPtBinLimits<<" != nPtBins+1 ("<<fnPtBins+1<<")\nSetting nPtBins to "<<nPtBinLimits-1<<endl;\r
+    SetNPtBins(nPtBinLimits-1);\r
+  }\r
+\r
+  fnPtBinLimits = nPtBinLimits;\r
+  SetGlobalIndex();\r
+  //cout<<"Changing also Global Index -> "<<fGlobalIndex<<endl;\r
+  fPtBinLimits = new Float_t[fnPtBinLimits];\r
+  for(Int_t ib=0; ib<nPtBinLimits; ib++) fPtBinLimits[ib]=ptBinLimits[ib];\r
+\r
+  return;\r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut){\r
+  // Set the variable names\r
+\r
+  if(fVarNames) {\r
+    delete [] fVarNames;\r
+    fVarNames = NULL;\r
+    //printf("Changing the variable names\n");\r
+  }\r
+  if(nVars!=fnVars){\r
+    printf("Wrong number of variables: it has to be %d\n",fnVars);\r
+    return;\r
+  }\r
+  //fnVars=nVars;\r
+  fVarNames = new TString[nVars];\r
+  fIsUpperCut = new Bool_t[nVars];\r
+  for(Int_t iv=0; iv<nVars; iv++) {\r
+    fVarNames[iv] = varNames[iv];\r
+    fIsUpperCut[iv] = isUpperCut[iv];\r
+  }\r
+\r
+  return;\r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::SetVarsForOpt(Int_t nVars,Bool_t *forOpt) {\r
+  // Set the variables to be used for cuts optimization\r
+\r
+  if(fVarsForOpt) {\r
+    delete [] fVarsForOpt;\r
+    fVarsForOpt = NULL;\r
+    //printf("Changing the variables for cut optimization\n");\r
+  }\r
+  \r
+  if(nVars==0){//!=fnVars) {\r
+    printf("%d not accepted as number of variables: it has to be %d\n",nVars,fnVars);\r
+    return;\r
+  } \r
+  \r
+  fnVarsForOpt = 0;\r
+  fVarsForOpt = new Bool_t[fnVars];\r
+  for(Int_t iv=0; iv<fnVars; iv++) {\r
+    fVarsForOpt[iv]=forOpt[iv];\r
+    if(fVarsForOpt[iv]) fnVarsForOpt++;\r
+  }\r
+\r
+  return;\r
+}\r
+\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::SetUseCentrality(Int_t flag) {\r
+  //\r
+  // set centrality estimator  \r
+  //\r
+  fUseCentrality=flag;\r
+  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");\r
\r
+  return;\r
+}\r
+\r
+\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::SetCuts(Int_t nVars,Int_t nPtBins,Float_t **cutsRD) {\r
+  //\r
+  // store the cuts\r
+  //\r
+  if(nVars!=fnVars) {\r
+    printf("Wrong number of variables: it has to be %d\n",fnVars);\r
+    AliFatal("exiting");\r
+  } \r
+  if(nPtBins!=fnPtBins) {\r
+    printf("Wrong number of pt bins: it has to be %d\n",fnPtBins);\r
+    AliFatal("exiting");\r
+  } \r
+\r
+  if(!fCutsRD)  fCutsRD = new Float_t[fGlobalIndex];\r
+  \r
+\r
+  for(Int_t iv=0; iv<fnVars; iv++) {\r
+\r
+    for(Int_t ib=0; ib<fnPtBins; ib++) {\r
+\r
+      //check\r
+      if(GetGlobalIndex(iv,ib)>=fGlobalIndex) {\r
+       cout<<"Overflow, exit..."<<endl;\r
+       return;\r
+      }\r
+\r
+      fCutsRD[GetGlobalIndex(iv,ib)] = cutsRD[iv][ib];\r
+\r
+    }\r
+  }\r
+  return;\r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::SetCuts(Int_t glIndex,Float_t* cutsRDGlob){\r
+  //\r
+  // store the cuts\r
+  //\r
+  if(glIndex != fGlobalIndex){\r
+    cout<<"Wrong array size: it has to be "<<fGlobalIndex<<endl;\r
+    AliFatal("exiting");\r
+  }\r
+  if(!fCutsRD)  fCutsRD = new Float_t[fGlobalIndex];\r
+\r
+  for(Int_t iGl=0;iGl<fGlobalIndex;iGl++){\r
+    fCutsRD[iGl] = cutsRDGlob[iGl];\r
+  }\r
+  return;\r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::PrintAll() const {\r
+  //\r
+  // print all cuts values\r
+  // \r
+\r
+  printf("Minimum vtx type %d\n",fMinVtxType);\r
+  printf("Minimum vtx contr %d\n",fMinVtxContr);\r
+  printf("Max vtx red chi2 %f\n",fMaxVtxRedChi2);\r
+  printf("Min SPD mult %d\n",fMinSPDMultiplicity);\r
+  printf("Use PID %d\n",(Int_t)fUsePID);\r
+  printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);\r
+  printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");\r
+  if(fOptPileup==1) printf(" -- Reject pileup event");\r
+  if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");\r
+  if(fUseCentrality>0) {\r
+    TString estimator="";\r
+    if(fUseCentrality==1) estimator = "V0";\r
+    if(fUseCentrality==2) estimator = "Tracks";\r
+    if(fUseCentrality==3) estimator = "Tracklets";\r
+    if(fUseCentrality==4) estimator = "SPD clusters outer"; \r
+    printf("Centrality class considered: %.1f-%.1f, estimated with %s",fMinCentrality,fMaxCentrality,estimator.Data());\r
+  }\r
+  if(fIsCandTrackSPDFirst) printf("Check for candidates with pt < %2.2f, that daughters fullfill kFirst criteria\n",fMaxPtCandTrackSPDFirst);\r
+\r
+  if(fVarNames){\r
+    cout<<"Array of variables"<<endl;\r
+    for(Int_t iv=0;iv<fnVars;iv++){\r
+      cout<<fVarNames[iv]<<"\t";\r
+    }\r
+    cout<<endl;\r
+  }\r
+  if(fVarsForOpt){\r
+    cout<<"Array of optimization"<<endl;\r
+    for(Int_t iv=0;iv<fnVars;iv++){\r
+      cout<<fVarsForOpt[iv]<<"\t";\r
+    }\r
+    cout<<endl;\r
+  }\r
+  if(fIsUpperCut){\r
+    cout<<"Array of upper/lower cut"<<endl;\r
+   for(Int_t iv=0;iv<fnVars;iv++){\r
+     cout<<fIsUpperCut[iv]<<"\t";\r
+   }\r
+   cout<<endl;\r
+  }\r
+  if(fPtBinLimits){\r
+    cout<<"Array of ptbin limits"<<endl;\r
+    for(Int_t ib=0;ib<fnPtBinLimits;ib++){\r
+      cout<<fPtBinLimits[ib]<<"\t";\r
+    }\r
+    cout<<endl;\r
+  }\r
+  if(fCutsRD){\r
+    cout<<"Matrix of cuts"<<endl;\r
+   for(Int_t iv=0;iv<fnVars;iv++){\r
+     for(Int_t ib=0;ib<fnPtBins;ib++){\r
+       cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"\t";\r
+     } \r
+     cout<<endl;\r
+   }\r
+   cout<<endl;\r
+  }\r
+  return;\r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::GetCuts(Float_t**& cutsRD) const{\r
+  //\r
+  // get the cuts\r
+  //\r
+\r
+  //cout<<"Give back a "<<fnVars<<"x"<<fnPtBins<<" matrix."<<endl;\r
+\r
+\r
+  Int_t iv,ib;\r
+  if(!cutsRD) {\r
+    //cout<<"Initialization..."<<endl;\r
+    cutsRD=new Float_t*[fnVars];\r
+    for(iv=0; iv<fnVars; iv++) {\r
+      cutsRD[iv] = new Float_t[fnPtBins];\r
+    }\r
+  }\r
+  \r
+  for(Int_t iGlobal=0; iGlobal<fGlobalIndex; iGlobal++) {\r
+    GetVarPtIndex(iGlobal,iv,ib);\r
+    cutsRD[iv][ib] = fCutsRD[iGlobal];\r
+  }\r
+\r
+  return;\r
+}\r
+\r
+//---------------------------------------------------------------------------\r
+Int_t AliRDHFCuts::GetGlobalIndex(Int_t iVar,Int_t iPtBin) const{\r
+  //\r
+  // give the global index from variable and pt bin\r
+  //\r
+  return iPtBin*fnVars+iVar;\r
+}\r
+\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const {\r
+  //\r
+  //give the index of the variable and of the pt bin from the global index\r
+  //\r
+  iPtBin=(Int_t)iGlob/fnVars;\r
+  iVar=iGlob%fnVars;\r
+\r
+  return;\r
+}\r
+\r
+//---------------------------------------------------------------------------\r
+Int_t AliRDHFCuts::PtBin(Double_t pt) const {\r
+  //\r
+  //give the pt bin where the pt lies.\r
+  //\r
+  Int_t ptbin=-1;\r
+  if(pt<fPtBinLimits[0])return ptbin;\r
+  for (Int_t i=0;i<fnPtBins;i++){\r
+    if(pt<fPtBinLimits[i+1]) {\r
+      ptbin=i;\r
+      break;\r
+    }\r
+  }\r
+  return ptbin;\r
+}\r
+//-------------------------------------------------------------------\r
+Float_t AliRDHFCuts::GetCutValue(Int_t iVar,Int_t iPtBin) const {\r
+  // \r
+  // Give the value of cut set for the variable iVar and the pt bin iPtBin\r
+  //\r
+  if(!fCutsRD){\r
+    cout<<"Cuts not iniziaisez yet"<<endl;\r
+    return 0;\r
+  }\r
+  return fCutsRD[GetGlobalIndex(iVar,iPtBin)];\r
+}\r
+//-------------------------------------------------------------------\r
+Float_t AliRDHFCuts::GetCentrality(AliAODEvent* aodEvent,AliRDHFCuts::ECentrality estimator) {\r
+  //\r
+  // Get centrality percentile\r
+  //\r
+\r
+  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aodEvent)->GetList()->FindObject(AliAODMCParticle::StdBranchName());\r
+  if(mcArray) {fUseAOD049=kFALSE;}\r
+\r
+  AliAODHeader *header=aodEvent->GetHeader();\r
+  AliCentrality *centrality=header->GetCentralityP();\r
+  Float_t cent=-999.;\r
+  Bool_t isSelRun=kFALSE;\r
+  Int_t selRun[5]={138364, 138826, 138828, 138836, 138871};\r
+  if(!centrality) return cent;\r
+  else{\r
+    if (estimator==kCentV0M){\r
+      cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));\r
+      if(cent<0){\r
+       Int_t quality = centrality->GetQuality();\r
+       if(quality<=1){\r
+         cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");\r
+       }else{\r
+         Int_t runnum=aodEvent->GetRunNumber();\r
+         for(Int_t ir=0;ir<5;ir++){\r
+           if(runnum==selRun[ir]){\r
+             isSelRun=kTRUE;\r
+             break;\r
+           }\r
+         }\r
+         if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");\r
+       }\r
+      }\r
+\r
+      //temporary fix for AOD049 outliers\r
+      if(fUseAOD049&&cent>=0){\r
+       Float_t v0=0;\r
+       AliAODVZERO* aodV0 = aodEvent->GetVZEROData();\r
+       v0+=aodV0->GetMTotV0A();\r
+       v0+=aodV0->GetMTotV0C();\r
+       if(cent==0&&v0<19500)return -1;//filtering issue\r
+       Float_t tkl = (Float_t)(aodEvent->GetTracklets()->GetNumberOfTracklets());\r
+       Float_t val= 1.30552 +  0.147931 * v0;\r
+       Float_t tklSigma[101]={176.644, 156.401, 153.789, 153.015, 142.476, 137.951, 136.127, 129.852, 127.436, 124.86, 120.788, 115.611, 113.172, 110.496, 109.127, 104.421, 102.479, 99.9766, 97.5152, 94.0654, 92.4602, 89.3364, 87.1342, 83.3497, 82.6216, 81.1084, 78.0793, 76.1234, 72.9434, 72.1334, 68.0056, 68.2755, 66.0376, 62.9666, 62.4274, 59.65, 58.3776, 56.6361, 54.5184, 53.4224, 51.932, 50.8922, 48.2848, 47.912, 46.5717, 43.4114, 43.2083, 41.3065, 40.1863, 38.5255, 37.2851, 37.5396, 34.4949, 33.8366, 31.8043, 31.7412, 30.8392, 30.0274, 28.8793, 27.6398, 26.6488, 25.0183, 25.1489, 24.4185, 22.9107, 21.2002, 21.6977, 20.1242, 20.4963, 19.0235, 19.298, 17.4103, 16.868, 15.2939, 15.2939, 16.0295, 14.186, 14.186, 15.2173, 12.9504, 12.9504, 12.9504, 15.264, 12.3674, 12.3674, 12.3674, 12.3674, 12.3674, 18.3811, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544};\r
+       if ( TMath::Abs(tkl-val) > 6.*tklSigma[(Int_t)cent] )return -1;//outlier\r
+      }\r
+    }\r
+    else {\r
+      if (estimator==kCentTRK) {\r
+       cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));\r
+       if(cent<0){\r
+         Int_t quality = centrality->GetQuality();\r
+         if(quality<=1){\r
+           cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");\r
+         }else{\r
+           Int_t runnum=aodEvent->GetRunNumber();\r
+           for(Int_t ir=0;ir<5;ir++){\r
+             if(runnum==selRun[ir]){\r
+               isSelRun=kTRUE;\r
+               break;\r
+             }\r
+           }\r
+           if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");\r
+         }\r
+       }\r
+      }\r
+      else{\r
+       if (estimator==kCentTKL){\r
+         cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));\r
+         if(cent<0){\r
+           Int_t quality = centrality->GetQuality();\r
+           if(quality<=1){\r
+             cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");\r
+           }else{\r
+             Int_t runnum=aodEvent->GetRunNumber();\r
+             for(Int_t ir=0;ir<5;ir++){\r
+               if(runnum==selRun[ir]){\r
+                 isSelRun=kTRUE;\r
+                 break;\r
+           }\r
+             }\r
+             if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");\r
+           }   \r
+         }\r
+       }\r
+       else{\r
+         if (estimator==kCentCL1){\r
+           cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));\r
+           if(cent<0){\r
+             Int_t quality = centrality->GetQuality();\r
+             if(quality<=1){\r
+               cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");\r
+             }else{\r
+               Int_t runnum=aodEvent->GetRunNumber();\r
+               for(Int_t ir=0;ir<5;ir++){\r
+                 if(runnum==selRun[ir]){\r
+                   isSelRun=kTRUE;\r
+                   break;\r
+                 }\r
+               }\r
+               if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");\r
+             }\r
+           }\r
+         }\r
+         else {\r
+           AliWarning("Centrality estimator not valid");\r
+           \r
+         }\r
+       }\r
+      }\r
+    } \r
+  }\r
+  return cent;\r
+}\r
+//-------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {\r
+  //\r
+  // Compare two cuts objects\r
+  //\r
+\r
+  Bool_t areEqual=kTRUE;\r
+\r
+  if(fMinVtxType!=obj->fMinVtxType) { printf("Minimum vtx type %d  %d\n",fMinVtxType,obj->fMinVtxType); areEqual=kFALSE;}\r
+\r
+  if(fMinVtxContr!=obj->fMinVtxContr) { printf("Minimum vtx contr %d  %d\n",fMinVtxContr,obj->fMinVtxContr); areEqual=kFALSE;}\r
+\r
+  if(TMath::Abs(fMaxVtxRedChi2-obj->fMaxVtxRedChi2)>1.e-10) {   printf("Max vtx red chi2 %f  %f\n",fMaxVtxRedChi2,obj->fMaxVtxRedChi2);areEqual=kFALSE;}\r
+\r
+  if(fMinSPDMultiplicity!=obj->fMinSPDMultiplicity) {  printf("Min SPD mult %d\n  %d",fMinSPDMultiplicity,obj->fMinSPDMultiplicity);areEqual=kFALSE;}\r
+\r
+  if(fUsePID!=obj->fUsePID) { printf("Use PID %d  %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}\r
+\r
+  if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d  %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}\r
+  if(fTrackCuts){\r
+    if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d  %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}\r
+\r
+    if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d  %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}\r
+\r
+    if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f  %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}\r
+\r
+    if(fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)!=obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)) {printf("ClusterReq SPD %d  %d\n",fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD),obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)); areEqual=kFALSE;}\r
+  }\r
+\r
+  if(fCutsRD) {\r
+   for(Int_t iv=0;iv<fnVars;iv++) {\r
+     for(Int_t ib=0;ib<fnPtBins;ib++) {\r
+       if(TMath::Abs(fCutsRD[GetGlobalIndex(iv,ib)]-obj->fCutsRD[GetGlobalIndex(iv,ib)])>1.e-10) {\r
+        cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"   "<<obj->fCutsRD[GetGlobalIndex(iv,ib)]<<"\n";\r
+        areEqual=kFALSE;\r
+       }\r
+     }\r
+   }\r
+  }\r
+\r
+  return areEqual;\r
+}\r
+//---------------------------------------------------------------------------\r
+void AliRDHFCuts::MakeTable() const {\r
+  //\r
+  // print cuts values in table format\r
+  // \r
+\r
+       TString ptString = "pT range";\r
+       if(fVarNames && fPtBinLimits && fCutsRD){\r
+               TString firstLine(Form("*       %-15s",ptString.Data()));\r
+               for (Int_t ivar=0; ivar<fnVars; ivar++){\r
+                       firstLine+=Form("*    %-15s  ",fVarNames[ivar].Data());\r
+                       if (ivar == fnVars){\r
+                               firstLine+="*\n";\r
+                       }\r
+               }\r
+               Printf("%s",firstLine.Data());\r
+               \r
+               for (Int_t ipt=0; ipt<fnPtBins; ipt++){\r
+                       TString line;\r
+                       if (ipt==fnPtBins-1){\r
+                               line=Form("*  %5.1f < pt < inf    ",fPtBinLimits[ipt]);\r
+                       }\r
+                       else{\r
+                               line=Form("*  %5.1f < pt < %4.1f   ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);\r
+                       }\r
+                       for (Int_t ivar=0; ivar<fnVars; ivar++){\r
+                               line+=Form("*     %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);\r
+                       }\r
+                       Printf("%s",line.Data());\r
+               }\r
+\r
+       }\r
+\r
+\r
+  return;\r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d,\r
+                                       AliAODEvent *aod) const\r
+{\r
+  //\r
+  // Recalculate primary vertex without daughters\r
+  //\r
+\r
+  if(!aod) {\r
+    AliError("Can not remove daughters from vertex without AOD event");\r
+    return 0;\r
+  }   \r
+\r
+  AliAODVertex *recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);\r
+  if(!recvtx){\r
+    AliDebug(2,"Removal of daughter tracks failed");\r
+    return kFALSE;\r
+  }\r
+\r
+\r
+  //set recalculed primary vertex\r
+  d->SetOwnPrimaryVtx(recvtx);\r
+  delete recvtx;\r
+\r
+  return kTRUE;\r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::SetMCPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const\r
+{\r
+  //\r
+  // Recalculate primary vertex without daughters\r
+  //\r
+\r
+  if(!aod) {\r
+    AliError("Can not get MC vertex without AOD event");\r
+    return kFALSE;\r
+  }   \r
+\r
+  // load MC header\r
+  AliAODMCHeader *mcHeader = \r
+    (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());\r
+  if(!mcHeader) {\r
+    AliError("Can not get MC vertex without AODMCHeader event");\r
+    return kFALSE;\r
+  }\r
+  Double_t pos[3];\r
+  Double_t covmatrix[6]={0.,0.,0.,0.,0.,0.};\r
+  mcHeader->GetVertex(pos);\r
+  AliAODVertex *recvtx=new AliAODVertex(pos,covmatrix);\r
+\r
+  if(!recvtx){\r
+    AliDebug(2,"Removal of daughter tracks failed");\r
+    return kFALSE;\r
+  }\r
+\r
+  //set recalculed primary vertex\r
+  d->SetOwnPrimaryVtx(recvtx);\r
+\r
+  d->RecalculateImpPars(recvtx,aod);\r
+\r
+  delete recvtx;\r
+\r
+  return kTRUE;\r
+}\r
+//--------------------------------------------------------------------------\r
+void AliRDHFCuts::CleanOwnPrimaryVtx(AliAODRecoDecayHF *d,\r
+                                    AliAODEvent *aod,\r
+                                    AliAODVertex *origownvtx) const\r
+{\r
+  //\r
+  // Clean-up own primary vertex if needed\r
+  //\r
+\r
+  if(fRemoveDaughtersFromPrimary || fUseMCVertex) {\r
+    d->UnsetOwnPrimaryVtx();\r
+    if(origownvtx) {\r
+      d->SetOwnPrimaryVtx(origownvtx);\r
+      delete origownvtx; origownvtx=NULL;\r
+    }\r
+    d->RecalculateImpPars(d->GetPrimaryVtx(),aod);\r
+  } else {\r
+    if(origownvtx) {\r
+      delete origownvtx; origownvtx=NULL;\r
+    }\r
+  }\r
+  return;\r
+}\r
+//--------------------------------------------------------------------------\r
+Bool_t AliRDHFCuts::IsSignalMC(AliAODRecoDecay *d,AliAODEvent *aod,Int_t pdg) const \r
+{\r
+  //\r
+  // Checks if this candidate is matched to MC signal\r
+  // \r
+\r
+  if(!aod) return kFALSE;\r
+\r
+  // get the MC array\r
+  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aod)->GetList()->FindObject(AliAODMCParticle::StdBranchName());\r
+\r
+  if(!mcArray) return kFALSE;\r
+\r
+  // try to match  \r
+  Int_t label = d->MatchToMC(pdg,mcArray);\r
+  \r
+  if(label>=0) {\r
+    //printf("MATCH!\n");\r
+    return kTRUE;\r
+  }\r
+\r
+  return kFALSE;\r
+}\r
+\r
+\r
index 1a993f7..4e16f0e 100644 (file)
-#ifndef ALIRDHFCUTS_H
-#define ALIRDHFCUTS_H
-/* Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice                               */
-
-//***********************************************************
-// Class AliRDHFCuts
-// base class for cuts on AOD reconstructed heavy-flavour decays
-// Author: A.Dainese, andrea.dainese@pd.infn.it
-//***********************************************************
-
-#include <TString.h>
-
-#include "AliAnalysisCuts.h"
-#include "AliESDtrackCuts.h"
-#include "AliAODPidHF.h"
-#include "AliAODEvent.h"
-#include "AliVEvent.h"
-
-class AliAODTrack;
-class AliAODRecoDecayHF;
-class AliESDVertex;
-
-class AliRDHFCuts : public AliAnalysisCuts 
-{
- public:
-
-  enum ECentrality {kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid};
-  enum ESelLevel {kAll,kTracks,kPID,kCandidate};
-  enum EPileup {kNoPileupSelection,kRejectPileupEvent,kRejectTracksFromPileupVertex};
-  enum ESele {kD0toKpiCuts,kD0toKpiPID,kD0fromDstarCuts,kD0fromDstarPID,kDplusCuts,kDplusPID,kDsCuts,kDsPID,kLcCuts,kLcPID,kDstarCuts,kDstarPID};
-  enum ERejBits {kNotSelTrigger,kNoVertex,kTooFewVtxContrib,kZVtxOutFid,kPileupSPD,kOutsideCentrality,kPhysicsSelection};
-  AliRDHFCuts(const Char_t* name="RDHFCuts", const Char_t* title="");
-  
-  virtual ~AliRDHFCuts();
-  
-  AliRDHFCuts(const AliRDHFCuts& source);
-  AliRDHFCuts& operator=(const AliRDHFCuts& source); 
-
-  virtual void SetStandardCutsPP2010() {return;}  
-  virtual void SetStandardCutsPbPb2010() {return;}  
-
-
-  void SetMinCentrality(Float_t minCentrality=0.) {fMinCentrality=minCentrality;} 
-  void SetMaxCentrality(Float_t maxCentrality=100.) {fMaxCentrality=maxCentrality;} 
-  void SetMinVtxType(Int_t type=3) {fMinVtxType=type;}  
-  void SetMinVtxContr(Int_t contr=1) {fMinVtxContr=contr;}  
-  void SetMaxVtxRdChi2(Float_t chi2=1e6) {fMaxVtxRedChi2=chi2;}  
-  void SetMaxVtxZ(Float_t z=1e6) {fMaxVtxZ=z;}  
-  void SetMinSPDMultiplicity(Int_t mult=0) {fMinSPDMultiplicity=mult;}  
-  void SetTriggerMask(ULong64_t mask=0) {fTriggerMask=mask;} 
-  void SetTriggerClass(TString trclass) {fTriggerClass=trclass;} 
-  void SetVarsForOpt(Int_t nVars,Bool_t *forOpt);
-  void SetGlobalIndex(){fGlobalIndex=fnVars*fnPtBins;}
-  void SetGlobalIndex(Int_t nVars,Int_t nptBins){fnVars=nVars; fnPtBins=nptBins; SetGlobalIndex();}
-  void SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut);  
-  void SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits);
-  void SetCuts(Int_t nVars,Int_t nPtBins,Float_t** cutsRD);
-  void SetCuts(Int_t glIndex, Float_t* cutsRDGlob);
-  void AddTrackCuts(const AliESDtrackCuts *cuts) 
-          {delete fTrackCuts; fTrackCuts=new AliESDtrackCuts(*cuts); return;}
-  void SetUsePID(Bool_t flag=kTRUE) {fUsePID=flag; return;}
-  void SetUseAOD049(Bool_t flag=kTRUE) {fUseAOD049=flag; return;}
-  void SetUseCentrality(Int_t flag=1);    // see enum below
-  void SetPidHF(AliAODPidHF* pidObj) {
-    if(fPidHF) delete fPidHF;
-    fPidHF=new AliAODPidHF(*pidObj);
-  }
-  void SetRemoveDaughtersFromPrim(Bool_t removeDaughtersPrim) {fRemoveDaughtersFromPrimary=removeDaughtersPrim;}
-  void SetMinPtCandidate(Double_t ptCand=-1.) {fMinPtCand=ptCand; return;}
-  void SetMaxPtCandidate(Double_t ptCand=1000.) {fMaxPtCand=ptCand; return;}
-  void SetOptPileup(Int_t opt=0){
-    // see enum below
-    fOptPileup=opt;
-  }
-  void ConfigurePileupCuts(Int_t minContrib=3, Float_t minDz=0.6){
-    fMinContrPileup=minContrib;
-    fMinDzPileup=minDz;
-  }
-
-
-  AliAODPidHF* GetPidHF() const {return fPidHF;}
-  Float_t *GetPtBinLimits() const {return fPtBinLimits;}
-  Int_t   GetNPtBins() const {return fnPtBins;}
-  Int_t   GetNVars() const {return fnVars;} 
-  TString *GetVarNames() const {return fVarNames;} 
-  Bool_t  *GetVarsForOpt() const {return fVarsForOpt;} 
-  Int_t   GetNVarsForOpt() const {return fnVarsForOpt;}
-  const Float_t *GetCuts() const {return fCutsRD;} 
-  void    GetCuts(Float_t**& cutsRD) const;
-  Float_t GetCutValue(Int_t iVar,Int_t iPtBin) const;
-  Double_t GetMaxVtxZ() const {return fMaxVtxZ;}  
-  Float_t GetCentrality(AliAODEvent* aodEvent){return GetCentrality(aodEvent,(AliRDHFCuts::ECentrality)fUseCentrality);}
-  Float_t GetCentrality(AliAODEvent* aodEvent, AliRDHFCuts::ECentrality estimator);
-  Bool_t  *GetIsUpperCut() const {return fIsUpperCut;}
-  AliESDtrackCuts *GetTrackCuts() const {return fTrackCuts;}
-  virtual AliESDtrackCuts *GetTrackCutsSoftPi() const {return 0;}
-  virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) = 0;
-  virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent * /*aod*/)
-            {return GetCutVarsForOpt(d,vars,nvars,pdgdaughters);}
-  Int_t   GetGlobalIndex(Int_t iVar,Int_t iPtBin) const;
-  void    GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const;
-  Bool_t  GetIsUsePID() const {return fUsePID;}
-  Bool_t  GetUseAOD049() const {return fUseAOD049;}
-  Bool_t  GetIsPrimaryWithoutDaughters() const {return fRemoveDaughtersFromPrimary;}
-  Bool_t GetOptPileUp() const {return fOptPileup;}
-  Int_t GetUseCentrality() const {return fUseCentrality;}
-  Float_t GetMinCentrality() const {return fMinCentrality;}
-  Float_t GetMaxCentrality() const {return fMaxCentrality;}
-  Double_t GetMinPtCandidate() const {return fMinPtCand;}
-  Double_t GetMaxPtCandidate() const {return fMaxPtCand;}
-  Bool_t IsSelected(TObject *obj) {return IsSelected(obj,AliRDHFCuts::kAll);}
-  Bool_t IsSelected(TList *list) {if(!list) return kTRUE; return kFALSE;}
-  Int_t  IsEventSelectedInCentrality(AliVEvent *event);
-  Bool_t IsEventSelected(AliVEvent *event);
-  Bool_t AreDaughtersSelected(AliAODRecoDecayHF *rd) const;
-  Bool_t IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts) const;
-  virtual Int_t IsSelectedPID(AliAODRecoDecayHF * /*rd*/) {return 1;}
-
-  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel) = 0;
-  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* /*aod*/)
-                {return IsSelected(obj,selectionLevel);}
-  Int_t PtBin(Double_t pt) const;
-  void PrintAll()const;
-
-  virtual Bool_t IsInFiducialAcceptance(Double_t /*pt*/,Double_t /*y*/) const {return kTRUE;}
-
-  void SetWhyRejection(Int_t why) {fWhyRejection=why; return;}
-  Int_t GetWhyRejection() const {return fWhyRejection;}
-  UInt_t GetEventRejectionBitMap() const {return fEvRejectionBits;}
-  Bool_t IsEventRejectedDueToTrigger() const {
-    return fEvRejectionBits&(1<<kNotSelTrigger);
-  }
-  Bool_t IsEventRejectedDueToNotRecoVertex() const {
-    return fEvRejectionBits&(1<<kNoVertex);
-  }
-  Bool_t IsEventRejectedDueToVertexContributors() const {
-    return fEvRejectionBits&(1<<kTooFewVtxContrib);
-  }
-  Bool_t IsEventRejectedDueToZVertexOutsideFiducialRegion() const {
-    return fEvRejectionBits&(1<<kZVtxOutFid);
-  }
-  Bool_t IsEventRejectedDueToPileupSPD() const {
-    return fEvRejectionBits&(1<<kPileupSPD);
-  }
-  Bool_t IsEventRejectedDueToCentrality() const {
-    return fEvRejectionBits&(1<<kOutsideCentrality);
-  }
-  Bool_t IsEventRejectedDuePhysicsSelection() const {
-    return fEvRejectionBits&(1<<kPhysicsSelection);
-  }
-
-
-  void SetFixRefs(Bool_t fix=kTRUE) {fFixRefs=fix; return;}
-  void SetUsePhysicsSelection(Bool_t use=kTRUE){fUsePhysicsSelection=use; return;}
-
-
-
-  Bool_t CompareCuts(const AliRDHFCuts *obj) const;
-  void MakeTable()const;
-
-  Int_t GetIsSelectedCuts() const {return fIsSelectedCuts;}
-  Int_t GetIsSelectedPID() const  {return fIsSelectedPID;}
-
-  void SetUseMCVertex() { fUseMCVertex=kTRUE; }
-  Bool_t GetUseMCVertex() const { return fUseMCVertex; }
-
-  Bool_t RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const;
-  Bool_t SetMCPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const;
-  void   CleanOwnPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod,AliAODVertex *origownvtx) const;
-
-  Bool_t CountEventForNormalization() const 
-  { if(fWhyRejection==0) {return kTRUE;} else {return kFALSE;} }
-
-  void SetKeepSignalMC() {fKeepSignalMC=kTRUE; return;}
-
-
- protected:
-
-  void SetNPtBins(Int_t nptBins){fnPtBins=nptBins;}
-  void SetNVars(Int_t nVars){fnVars=nVars;}
-
-  Bool_t IsSignalMC(AliAODRecoDecay *d,AliAODEvent *aod,Int_t pdg) const;
-
-  // cuts on the event
-  Int_t fMinVtxType; // 0: not cut; 1: SPDZ; 2: SPD3D; 3: Tracks
-  Int_t fMinVtxContr;   // minimum vertex contributors
-  Float_t fMaxVtxRedChi2; // maximum chi2/ndf
-  Float_t fMaxVtxZ; // maximum |z| of primary vertex
-  Int_t fMinSPDMultiplicity; // SPD multiplicity
-  ULong64_t fTriggerMask; // trigger mask
-  TString  fTriggerClass; // trigger class
-  // quality cuts on the daughter tracks
-  AliESDtrackCuts *fTrackCuts; // tracks for daughter tracks (AOD converted to ESD on the flight!)
-  // cuts on the candidate
-  Int_t fnPtBins;  // number of pt bins for cuts
-  Int_t fnPtBinLimits; // "number of limits", that is fnPtBins+1
-  Float_t* fPtBinLimits; //[fnPtBinLimits]  pt bins
-  Int_t fnVars;    // number of cut vars for candidates
-  TString *fVarNames; //[fnVars] names of the variables
-  Int_t fnVarsForOpt;    // number of cut vars to be optimized for candidates
-  Bool_t *fVarsForOpt; //[fnVars] kTRUE for vars to be used in optimization
-  Int_t fGlobalIndex; // fnVars*fnPtBins
-  Float_t *fCutsRD; //[fGlobalIndex] the cuts values
-  Bool_t  *fIsUpperCut; //[fnVars] use > or < to select
-  Bool_t fUsePID; // enable PID usage (off by default)
-  Bool_t fUseAOD049; // enable AOD049 centrality cleanup
-  AliAODPidHF *fPidHF; // PID for heavy flavours manager
-  Int_t fWhyRejection; // used to code the step at which candidate was rejected
-  UInt_t fEvRejectionBits; //bit map storing the full info about event rejection
-  Bool_t fRemoveDaughtersFromPrimary; // flag to switch on the removal of duaghters from the primary vertex computation
-  Bool_t fUseMCVertex; // use MC primary vertex 
-  Bool_t fUsePhysicsSelection; // use Physics selection criteria
-  Int_t  fOptPileup;      // option for pielup selection
-  Int_t  fMinContrPileup; // min. n. of tracklets in pileup vertex
-  Float_t fMinDzPileup;   // min deltaz between main and pileup vertices
-  Int_t   fUseCentrality; // off =0 (default)
-                          // 1 = V0 
-                          // 2 = Tracks
-                          // 3 = Tracklets
-                          // 4 = SPD clusters outer 
-  Float_t fMinCentrality; // minimum centrality for selected events
-  Float_t fMaxCentrality; // maximum centrality for selected events
-  Bool_t  fFixRefs;       // fix the daughter track references 
-  Int_t  fIsSelectedCuts; // outcome of cuts selection
-  Int_t  fIsSelectedPID;  // outcome of PID selection
-  Double_t fMinPtCand; // minimum pt of the candidate
-  Double_t fMaxPtCand; // minimum pt of the candidate
-  Bool_t  fKeepSignalMC; // IsSelected returns always kTRUE for MC signal
-
-  ClassDef(AliRDHFCuts,19);  // base class for cuts on AOD reconstructed heavy-flavour decays
-};
-
-#endif
+#ifndef ALIRDHFCUTS_H\r
+#define ALIRDHFCUTS_H\r
+/* Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *\r
+ * See cxx source for full Copyright notice                               */\r
+\r
+//***********************************************************\r
+// Class AliRDHFCuts\r
+// base class for cuts on AOD reconstructed heavy-flavour decays\r
+// Author: A.Dainese, andrea.dainese@pd.infn.it\r
+//***********************************************************\r
+\r
+#include <TString.h>\r
+\r
+#include "AliAnalysisCuts.h"\r
+#include "AliESDtrackCuts.h"\r
+#include "AliAODPidHF.h"\r
+#include "AliAODEvent.h"\r
+#include "AliVEvent.h"\r
+\r
+class AliAODTrack;\r
+class AliAODRecoDecayHF;\r
+class AliESDVertex;\r
+\r
+class AliRDHFCuts : public AliAnalysisCuts \r
+{\r
+ public:\r
+\r
+  enum ECentrality {kCentOff,kCentV0M,kCentTRK,kCentTKL,kCentCL1,kCentInvalid};\r
+  enum ESelLevel {kAll,kTracks,kPID,kCandidate};\r
+  enum EPileup {kNoPileupSelection,kRejectPileupEvent,kRejectTracksFromPileupVertex};\r
+  enum ESele {kD0toKpiCuts,kD0toKpiPID,kD0fromDstarCuts,kD0fromDstarPID,kDplusCuts,kDplusPID,kDsCuts,kDsPID,kLcCuts,kLcPID,kDstarCuts,kDstarPID};\r
+  enum ERejBits {kNotSelTrigger,kNoVertex,kTooFewVtxContrib,kZVtxOutFid,kPileupSPD,kOutsideCentrality,kPhysicsSelection};\r
+  AliRDHFCuts(const Char_t* name="RDHFCuts", const Char_t* title="");\r
+  \r
+  virtual ~AliRDHFCuts();\r
+  \r
+  AliRDHFCuts(const AliRDHFCuts& source);\r
+  AliRDHFCuts& operator=(const AliRDHFCuts& source); \r
+\r
+  virtual void SetStandardCutsPP2010() {return;}  \r
+  virtual void SetStandardCutsPbPb2010() {return;}  \r
+\r
+\r
+  void SetMinCentrality(Float_t minCentrality=0.) {fMinCentrality=minCentrality;} \r
+  void SetMaxCentrality(Float_t maxCentrality=100.) {fMaxCentrality=maxCentrality;} \r
+  void SetMinVtxType(Int_t type=3) {fMinVtxType=type;}  \r
+  void SetMinVtxContr(Int_t contr=1) {fMinVtxContr=contr;}  \r
+  void SetMaxVtxRdChi2(Float_t chi2=1e6) {fMaxVtxRedChi2=chi2;}  \r
+  void SetMaxVtxZ(Float_t z=1e6) {fMaxVtxZ=z;}  \r
+  void SetMinSPDMultiplicity(Int_t mult=0) {fMinSPDMultiplicity=mult;}  \r
+  void SetTriggerMask(ULong64_t mask=0) {fTriggerMask=mask;} \r
+  void SetTriggerClass(TString trclass) {fTriggerClass=trclass;} \r
+  void SetVarsForOpt(Int_t nVars,Bool_t *forOpt);\r
+  void SetGlobalIndex(){fGlobalIndex=fnVars*fnPtBins;}\r
+  void SetGlobalIndex(Int_t nVars,Int_t nptBins){fnVars=nVars; fnPtBins=nptBins; SetGlobalIndex();}\r
+  void SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut);  \r
+  void SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits);\r
+  void SetCuts(Int_t nVars,Int_t nPtBins,Float_t** cutsRD);\r
+  void SetCuts(Int_t glIndex, Float_t* cutsRDGlob);\r
+  void AddTrackCuts(const AliESDtrackCuts *cuts) \r
+          {delete fTrackCuts; fTrackCuts=new AliESDtrackCuts(*cuts); return;}\r
+  void SetUsePID(Bool_t flag=kTRUE) {fUsePID=flag; return;}\r
+  void SetUseAOD049(Bool_t flag=kTRUE) {fUseAOD049=flag; return;}\r
+  void SetUseCentrality(Int_t flag=1);    // see enum below\r
+  void SetPidHF(AliAODPidHF* pidObj) {\r
+    if(fPidHF) delete fPidHF;\r
+    fPidHF=new AliAODPidHF(*pidObj);\r
+  }\r
+  void SetRemoveDaughtersFromPrim(Bool_t removeDaughtersPrim) {fRemoveDaughtersFromPrimary=removeDaughtersPrim;}\r
+  void SetMinPtCandidate(Double_t ptCand=-1.) {fMinPtCand=ptCand; return;}\r
+  void SetMaxPtCandidate(Double_t ptCand=1000.) {fMaxPtCand=ptCand; return;}\r
+  void SetOptPileup(Int_t opt=0){\r
+    // see enum below\r
+    fOptPileup=opt;\r
+  }\r
+  void ConfigurePileupCuts(Int_t minContrib=3, Float_t minDz=0.6){\r
+    fMinContrPileup=minContrib;\r
+    fMinDzPileup=minDz;\r
+  }\r
+\r
+\r
+  AliAODPidHF* GetPidHF() const {return fPidHF;}\r
+  Float_t *GetPtBinLimits() const {return fPtBinLimits;}\r
+  Int_t   GetNPtBins() const {return fnPtBins;}\r
+  Int_t   GetNVars() const {return fnVars;} \r
+  TString *GetVarNames() const {return fVarNames;} \r
+  Bool_t  *GetVarsForOpt() const {return fVarsForOpt;} \r
+  Int_t   GetNVarsForOpt() const {return fnVarsForOpt;}\r
+  const Float_t *GetCuts() const {return fCutsRD;} \r
+  void    GetCuts(Float_t**& cutsRD) const;\r
+  Float_t GetCutValue(Int_t iVar,Int_t iPtBin) const;\r
+  Double_t GetMaxVtxZ() const {return fMaxVtxZ;}  \r
+  Float_t GetCentrality(AliAODEvent* aodEvent){return GetCentrality(aodEvent,(AliRDHFCuts::ECentrality)fUseCentrality);}\r
+  Float_t GetCentrality(AliAODEvent* aodEvent, AliRDHFCuts::ECentrality estimator);\r
+  Bool_t  *GetIsUpperCut() const {return fIsUpperCut;}\r
+  AliESDtrackCuts *GetTrackCuts() const {return fTrackCuts;}\r
+  virtual AliESDtrackCuts *GetTrackCutsSoftPi() const {return 0;}\r
+  virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) = 0;\r
+  virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent * /*aod*/)\r
+            {return GetCutVarsForOpt(d,vars,nvars,pdgdaughters);}\r
+  Int_t   GetGlobalIndex(Int_t iVar,Int_t iPtBin) const;\r
+  void    GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const;\r
+  Bool_t  GetIsUsePID() const {return fUsePID;}\r
+  Bool_t  GetUseAOD049() const {return fUseAOD049;}\r
+  Bool_t  GetIsPrimaryWithoutDaughters() const {return fRemoveDaughtersFromPrimary;}\r
+  Bool_t GetOptPileUp() const {return fOptPileup;}\r
+  Int_t GetUseCentrality() const {return fUseCentrality;}\r
+  Float_t GetMinCentrality() const {return fMinCentrality;}\r
+  Float_t GetMaxCentrality() const {return fMaxCentrality;}\r
+  Double_t GetMinPtCandidate() const {return fMinPtCand;}\r
+  Double_t GetMaxPtCandidate() const {return fMaxPtCand;}\r
+  Bool_t IsSelected(TObject *obj) {return IsSelected(obj,AliRDHFCuts::kAll);}\r
+  Bool_t IsSelected(TList *list) {if(!list) return kTRUE; return kFALSE;}\r
+  Int_t  IsEventSelectedInCentrality(AliVEvent *event);\r
+  Bool_t IsEventSelected(AliVEvent *event);\r
+  Bool_t AreDaughtersSelected(AliAODRecoDecayHF *rd) const;\r
+  Bool_t IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts) const;\r
+  virtual Int_t IsSelectedPID(AliAODRecoDecayHF * /*rd*/) {return 1;}\r
+\r
+  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel) = 0;\r
+  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* /*aod*/)\r
+                {return IsSelected(obj,selectionLevel);}\r
+  Int_t PtBin(Double_t pt) const;\r
+  void PrintAll()const;\r
+\r
+  virtual Bool_t IsInFiducialAcceptance(Double_t /*pt*/,Double_t /*y*/) const {return kTRUE;}\r
+\r
+  void SetWhyRejection(Int_t why) {fWhyRejection=why; return;}\r
+  Int_t GetWhyRejection() const {return fWhyRejection;}\r
+  UInt_t GetEventRejectionBitMap() const {return fEvRejectionBits;}\r
+  Bool_t IsEventRejectedDueToTrigger() const {\r
+    return fEvRejectionBits&(1<<kNotSelTrigger);\r
+  }\r
+  Bool_t IsEventRejectedDueToNotRecoVertex() const {\r
+    return fEvRejectionBits&(1<<kNoVertex);\r
+  }\r
+  Bool_t IsEventRejectedDueToVertexContributors() const {\r
+    return fEvRejectionBits&(1<<kTooFewVtxContrib);\r
+  }\r
+  Bool_t IsEventRejectedDueToZVertexOutsideFiducialRegion() const {\r
+    return fEvRejectionBits&(1<<kZVtxOutFid);\r
+  }\r
+  Bool_t IsEventRejectedDueToPileupSPD() const {\r
+    return fEvRejectionBits&(1<<kPileupSPD);\r
+  }\r
+  Bool_t IsEventRejectedDueToCentrality() const {\r
+    return fEvRejectionBits&(1<<kOutsideCentrality);\r
+  }\r
+  Bool_t IsEventRejectedDuePhysicsSelection() const {\r
+    return fEvRejectionBits&(1<<kPhysicsSelection);\r
+  }\r
+\r
+\r
+  void SetFixRefs(Bool_t fix=kTRUE) {fFixRefs=fix; return;}\r
+  void SetUsePhysicsSelection(Bool_t use=kTRUE){fUsePhysicsSelection=use; return;}\r
+\r
+\r
+\r
+  Bool_t CompareCuts(const AliRDHFCuts *obj) const;\r
+  void MakeTable()const;\r
+\r
+  Int_t GetIsSelectedCuts() const {return fIsSelectedCuts;}\r
+  Int_t GetIsSelectedPID() const  {return fIsSelectedPID;}\r
+\r
+  void SetUseMCVertex() { fUseMCVertex=kTRUE; }\r
+  Bool_t GetUseMCVertex() const { return fUseMCVertex; }\r
+\r
+  Bool_t RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const;\r
+  Bool_t SetMCPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const;\r
+  void   CleanOwnPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod,AliAODVertex *origownvtx) const;\r
+\r
+  Bool_t CountEventForNormalization() const \r
+  { if(fWhyRejection==0) {return kTRUE;} else {return kFALSE;} }\r
+\r
+  void SetKeepSignalMC() {fKeepSignalMC=kTRUE; return;}\r
+\r
+  // Flag and pt-maximum to check if the candidate daughters fulfill the kFirst criteria\r
+  void SetSelectCandTrackSPDFirst( Bool_t flag, Double_t ptmax )\r
+  { fIsCandTrackSPDFirst=flag; fMaxPtCandTrackSPDFirst=ptmax; }\r
+  Bool_t IsSelectCandTrackSPDFirst() const { return fIsCandTrackSPDFirst; }\r
+  Double_t IsMaxCandTrackSPDFirst() const { return fMaxPtCandTrackSPDFirst; }\r
+\r
+\r
+ protected:\r
+\r
+  void SetNPtBins(Int_t nptBins){fnPtBins=nptBins;}\r
+  void SetNVars(Int_t nVars){fnVars=nVars;}\r
+\r
+  Bool_t IsSignalMC(AliAODRecoDecay *d,AliAODEvent *aod,Int_t pdg) const;\r
+\r
+  // cuts on the event\r
+  Int_t fMinVtxType; // 0: not cut; 1: SPDZ; 2: SPD3D; 3: Tracks\r
+  Int_t fMinVtxContr;   // minimum vertex contributors\r
+  Float_t fMaxVtxRedChi2; // maximum chi2/ndf\r
+  Float_t fMaxVtxZ; // maximum |z| of primary vertex\r
+  Int_t fMinSPDMultiplicity; // SPD multiplicity\r
+  ULong64_t fTriggerMask; // trigger mask\r
+  TString  fTriggerClass; // trigger class\r
+  // quality cuts on the daughter tracks\r
+  AliESDtrackCuts *fTrackCuts; // tracks for daughter tracks (AOD converted to ESD on the flight!)\r
+  // cuts on the candidate\r
+  Int_t fnPtBins;  // number of pt bins for cuts\r
+  Int_t fnPtBinLimits; // "number of limits", that is fnPtBins+1\r
+  Float_t* fPtBinLimits; //[fnPtBinLimits]  pt bins\r
+  Int_t fnVars;    // number of cut vars for candidates\r
+  TString *fVarNames; //[fnVars] names of the variables\r
+  Int_t fnVarsForOpt;    // number of cut vars to be optimized for candidates\r
+  Bool_t *fVarsForOpt; //[fnVars] kTRUE for vars to be used in optimization\r
+  Int_t fGlobalIndex; // fnVars*fnPtBins\r
+  Float_t *fCutsRD; //[fGlobalIndex] the cuts values\r
+  Bool_t  *fIsUpperCut; //[fnVars] use > or < to select\r
+  Bool_t fUsePID; // enable PID usage (off by default)\r
+  Bool_t fUseAOD049; // enable AOD049 centrality cleanup\r
+  AliAODPidHF *fPidHF; // PID for heavy flavours manager\r
+  Int_t fWhyRejection; // used to code the step at which candidate was rejected\r
+  UInt_t fEvRejectionBits; //bit map storing the full info about event rejection\r
+  Bool_t fRemoveDaughtersFromPrimary; // flag to switch on the removal of duaghters from the primary vertex computation\r
+  Bool_t fUseMCVertex; // use MC primary vertex \r
+  Bool_t fUsePhysicsSelection; // use Physics selection criteria\r
+  Int_t  fOptPileup;      // option for pielup selection\r
+  Int_t  fMinContrPileup; // min. n. of tracklets in pileup vertex\r
+  Float_t fMinDzPileup;   // min deltaz between main and pileup vertices\r
+  Int_t   fUseCentrality; // off =0 (default)\r
+                          // 1 = V0 \r
+                          // 2 = Tracks\r
+                          // 3 = Tracklets\r
+                          // 4 = SPD clusters outer \r
+  Float_t fMinCentrality; // minimum centrality for selected events\r
+  Float_t fMaxCentrality; // maximum centrality for selected events\r
+  Bool_t  fFixRefs;       // fix the daughter track references \r
+  Int_t  fIsSelectedCuts; // outcome of cuts selection\r
+  Int_t  fIsSelectedPID;  // outcome of PID selection\r
+  Double_t fMinPtCand; // minimum pt of the candidate\r
+  Double_t fMaxPtCand; // minimum pt of the candidate\r
+  Bool_t  fKeepSignalMC; // IsSelected returns always kTRUE for MC signal\r
+  Bool_t fIsCandTrackSPDFirst; // flag to select the track kFirst criteria for pt < ptlimit\r
+  Double_t fMaxPtCandTrackSPDFirst; // maximum pt of the candidate for which to check if the daughters fulfill kFirst criteria\r
+\r
+\r
+  ClassDef(AliRDHFCuts,20);  // base class for cuts on AOD reconstructed heavy-flavour decays\r
+};\r
+\r
+#endif\r