-/**************************************************************************
- * 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&¢>=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&¢>=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
-#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