default D0 cuts for the low energy run (Zaida)
authorprino <prino@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 1 May 2012 14:39:20 +0000 (14:39 +0000)
committerprino <prino@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 1 May 2012 14:39:20 +0000 (14:39 +0000)
PWGHF/vertexingHF/AliRDHFCutsD0toKpi.cxx
PWGHF/vertexingHF/AliRDHFCutsD0toKpi.h

index 20d2c7f..143f45f 100644 (file)
-/**************************************************************************\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
-// Class for cuts on AOD reconstructed D0->Kpi\r
-//\r
-// Author: A.Dainese, andrea.dainese@pd.infn.it\r
-/////////////////////////////////////////////////////////////\r
-\r
-#include <TDatabasePDG.h>\r
-#include <Riostream.h>\r
-\r
-#include "AliRDHFCutsD0toKpi.h"\r
-#include "AliAODRecoDecayHF2Prong.h"\r
-#include "AliAODTrack.h"\r
-#include "AliESDtrack.h"\r
-#include "AliAODPid.h"\r
-#include "AliTPCPIDResponse.h"\r
-#include "AliAODVertex.h"\r
-#include "AliKFVertex.h"\r
-#include "AliKFParticle.h"\r
-\r
-ClassImp(AliRDHFCutsD0toKpi)\r
-\r
-//--------------------------------------------------------------------------\r
-AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const char* name) : \r
-AliRDHFCuts(name),\r
-fUseSpecialCuts(kFALSE),\r
-fLowPt(kTRUE),\r
-fDefaultPID(kFALSE),\r
-fUseKF(kFALSE),\r
-fPtLowPID(2.),\r
-fPtMaxSpecialCuts(9999.)\r
-{\r
-  //\r
-  // Default Constructor\r
-  //\r
-  Int_t nvars=11;\r
-  SetNVars(nvars);\r
-  TString varNames[11]={"inv. mass [GeV]",   \r
-                       "dca [cm]",\r
-                       "cosThetaStar", \r
-                       "pTK [GeV/c]",\r
-                       "pTPi [GeV/c]",\r
-                       "d0K [cm]",\r
-                       "d0Pi [cm]",\r
-                       "d0d0 [cm^2]",\r
-                       "cosThetaPoint",\r
-                       "|cosThetaPointXY|", \r
-                       "NormDecayLenghtXY"};\r
-  Bool_t isUpperCut[11]={kTRUE,\r
-                        kTRUE,\r
-                        kTRUE,\r
-                        kFALSE,\r
-                        kFALSE,\r
-                        kTRUE,\r
-                        kTRUE,\r
-                        kTRUE,\r
-                        kFALSE,\r
-                        kFALSE, \r
-                        kFALSE};\r
-  SetVarNames(nvars,varNames,isUpperCut);\r
-  Bool_t forOpt[11]={kFALSE,\r
-                    kTRUE,\r
-                    kTRUE,\r
-                    kFALSE,\r
-                    kFALSE,\r
-                    kFALSE,\r
-                    kFALSE,\r
-                    kTRUE,\r
-                    kTRUE,\r
-                    kFALSE,\r
-                    kFALSE};\r
-  SetVarsForOpt(4,forOpt);\r
-  Float_t limits[2]={0,999999999.};\r
-  SetPtBins(2,limits);\r
-\r
-}\r
-//--------------------------------------------------------------------------\r
-AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const AliRDHFCutsD0toKpi &source) :\r
-  AliRDHFCuts(source),   \r
-  fUseSpecialCuts(source.fUseSpecialCuts),\r
-  fLowPt(source.fLowPt),\r
-  fDefaultPID(source.fDefaultPID),\r
-  fUseKF(source.fUseKF),\r
-  fPtLowPID(source.fPtLowPID),\r
-  fPtMaxSpecialCuts(source.fPtMaxSpecialCuts)\r
-{\r
-  //\r
-  // Copy constructor\r
-  //\r
-\r
-}\r
-//--------------------------------------------------------------------------\r
-AliRDHFCutsD0toKpi &AliRDHFCutsD0toKpi::operator=(const AliRDHFCutsD0toKpi &source)\r
-{\r
-  //\r
-  // assignment operator\r
-  //\r
-  if(&source == this) return *this;\r
-\r
-  AliRDHFCuts::operator=(source); \r
-  fUseSpecialCuts=source.fUseSpecialCuts;\r
-  fLowPt=source.fLowPt;\r
-  fDefaultPID=source.fDefaultPID;\r
-  fUseKF=source.fUseKF;\r
-  fPtLowPID=source.fPtLowPID;\r
-  fPtMaxSpecialCuts=source.fPtMaxSpecialCuts;\r
-\r
-  return *this;\r
-}\r
-\r
-\r
-//---------------------------------------------------------------------------\r
-void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent *aod) {\r
-  // \r
-  // Fills in vars the values of the variables \r
-  //\r
-\r
-  if(nvars!=fnVarsForOpt) {\r
-    printf("AliRDHFCutsD0toKpi::GetCutsVarsForOpt: wrong number of variables\n");\r
-    return;\r
-  }\r
-\r
-  AliAODRecoDecayHF2Prong *dd = (AliAODRecoDecayHF2Prong*)d;\r
\r
-  //recalculate vertex w/o daughters\r
-  Bool_t cleanvtx=kFALSE;\r
-  AliAODVertex *origownvtx=0x0;\r
-  if(fRemoveDaughtersFromPrimary) {\r
-    if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());\r
-    cleanvtx=kTRUE;\r
-    if(!RecalcOwnPrimaryVtx(dd,aod)) {\r
-      CleanOwnPrimaryVtx(dd,aod,origownvtx);\r
-      cleanvtx=kFALSE;\r
-    }\r
-  }\r
-\r
-  Int_t iter=-1;\r
-  if(fVarsForOpt[0]){\r
-    iter++;\r
-    if(TMath::Abs(pdgdaughters[0])==211) {\r
-      vars[iter]=dd->InvMassD0();\r
-    } else {\r
-      vars[iter]=dd->InvMassD0bar();\r
-    }\r
-  }\r
-  if(fVarsForOpt[1]){\r
-    iter++;\r
-    vars[iter]=dd->GetDCA();\r
-  }\r
-  if(fVarsForOpt[2]){\r
-    iter++;\r
-    if(TMath::Abs(pdgdaughters[0])==211) {\r
-      vars[iter] = dd->CosThetaStarD0();\r
-    } else {\r
-      vars[iter] = dd->CosThetaStarD0bar();\r
-    }\r
-  }\r
-  if(fVarsForOpt[3]){\r
-    iter++;\r
-   if(TMath::Abs(pdgdaughters[0])==321) {\r
-     vars[iter]=dd->PtProng(0);\r
-   }\r
-   else{\r
-     vars[iter]=dd->PtProng(1);\r
-   }\r
-  }\r
-  if(fVarsForOpt[4]){\r
-    iter++;\r
-   if(TMath::Abs(pdgdaughters[0])==211) {\r
-     vars[iter]=dd->PtProng(0);\r
-   }\r
-   else{\r
-     vars[iter]=dd->PtProng(1);\r
-   }\r
-  }\r
-  if(fVarsForOpt[5]){\r
-    iter++;\r
-    if(TMath::Abs(pdgdaughters[0])==321) {\r
-     vars[iter]=dd->Getd0Prong(0);\r
-   }\r
-   else{\r
-     vars[iter]=dd->Getd0Prong(1);\r
-   }\r
-  }\r
-  if(fVarsForOpt[6]){\r
-    iter++;\r
-     if(TMath::Abs(pdgdaughters[0])==211) {\r
-     vars[iter]=dd->Getd0Prong(0);\r
-   }\r
-   else{\r
-     vars[iter]=dd->Getd0Prong(1);\r
-   }\r
-  }\r
-  if(fVarsForOpt[7]){\r
-    iter++;\r
-    vars[iter]= dd->Prodd0d0();\r
-  }\r
-  if(fVarsForOpt[8]){\r
-    iter++;\r
-    vars[iter]=dd->CosPointingAngle();\r
-  }\r
-  \r
-  if(fVarsForOpt[9]){\r
-               iter++;\r
-               vars[iter]=TMath::Abs(dd->CosPointingAngleXY());\r
-       }\r
-  \r
-   if(fVarsForOpt[10]){\r
-               iter++;\r
-          vars[iter]=dd->NormalizedDecayLengthXY();\r
-       }\r
-\r
-   if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);\r
-\r
-  return;\r
-}\r
-//---------------------------------------------------------------------------\r
-Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* aod) {\r
-  //\r
-  // Apply selection\r
-  //\r
\r
-\r
-  fIsSelectedCuts=0;\r
-  fIsSelectedPID=0;\r
-\r
-  if(!fCutsRD){\r
-    cout<<"Cut matrice not inizialized. Exit..."<<endl;\r
-    return 0;\r
-  }\r
-  //PrintAll();\r
-  AliAODRecoDecayHF2Prong* d=(AliAODRecoDecayHF2Prong*)obj;\r
-\r
-  if(!d){\r
-    cout<<"AliAODRecoDecayHF2Prong null"<<endl;\r
-    return 0;\r
-  }\r
-\r
-  if(fKeepSignalMC) if(IsSignalMC(d,aod,421)) return 3;\r
-\r
-  Double_t ptD=d->Pt();\r
-  if(ptD<fMinPtCand) return 0;\r
-  if(ptD>fMaxPtCand) return 0;\r
-\r
-  if(d->HasBadDaughters()) return 0;\r
-\r
-  // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both\r
-  Int_t returnvaluePID=3;\r
-  Int_t returnvalueCuts=3;\r
-\r
-  // selection on candidate\r
-  if(selectionLevel==AliRDHFCuts::kAll || \r
-     selectionLevel==AliRDHFCuts::kCandidate) {\r
-\r
-    if(!fUseKF) {\r
-\r
-      //recalculate vertex w/o daughters\r
-      AliAODVertex *origownvtx=0x0;\r
-      if(fRemoveDaughtersFromPrimary && !fUseMCVertex) {\r
-       if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());\r
-       if(!RecalcOwnPrimaryVtx(d,aod)) { \r
-         CleanOwnPrimaryVtx(d,aod,origownvtx);\r
-         return 0;\r
-       }\r
-      }\r
-\r
-      if(fUseMCVertex) {\r
-       if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());\r
-       if(!SetMCPrimaryVtx(d,aod)) {\r
-         CleanOwnPrimaryVtx(d,aod,origownvtx);\r
-         return 0;\r
-       }\r
-      }\r
-      \r
-      Double_t pt=d->Pt();\r
-      \r
-      Int_t okD0=0,okD0bar=0;\r
-      \r
-      Int_t ptbin=PtBin(pt);\r
-      if (ptbin==-1) {\r
-       CleanOwnPrimaryVtx(d,aod,origownvtx);\r
-       return 0;\r
-      }\r
-\r
-      Double_t mD0,mD0bar,ctsD0,ctsD0bar;\r
-      okD0=1; okD0bar=1;\r
-      \r
-      Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
-\r
-      d->InvMassD0(mD0,mD0bar);\r
-      if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;\r
-      if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)])  okD0bar = 0;\r
-      if(!okD0 && !okD0bar)  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-\r
-      if(d->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-    \r
-      \r
-      if(d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(3,ptbin)]*fCutsRD[GetGlobalIndex(3,ptbin)] || d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(4,ptbin)]*fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;\r
-      if(d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(3,ptbin)]*fCutsRD[GetGlobalIndex(3,ptbin)] || d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(4,ptbin)]*fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;\r
-      if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-      \r
-      \r
-      if(TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] || \r
-        TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;\r
-      if(TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||\r
-        TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;\r
-      if(!okD0 && !okD0bar)  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-      \r
-      if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-      \r
-    \r
-      d->CosThetaStarD0(ctsD0,ctsD0bar);\r
-      if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0; \r
-      if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;\r
-      if(!okD0 && !okD0bar)   {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-    \r
-      if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-\r
-      \r
-      if(TMath::Abs(d->CosPointingAngleXY()) < fCutsRD[GetGlobalIndex(9,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-       \r
-      Double_t normalDecayLengXY=d->NormalizedDecayLengthXY();\r
-      if (normalDecayLengXY < fCutsRD[GetGlobalIndex(10, ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-      \r
-      if (returnvalueCuts!=0) {\r
-       if (okD0) returnvalueCuts=1; //cuts passed as D0\r
-       if (okD0bar) returnvalueCuts=2; //cuts passed as D0bar\r
-       if (okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar\r
-      }\r
-\r
-      // call special cuts\r
-      Int_t special=1;\r
-      if(fUseSpecialCuts && (pt<fPtMaxSpecialCuts)) special=IsSelectedSpecialCuts(d);\r
-      if(!special) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}\r
-\r
-      // unset recalculated primary vertex when not needed any more\r
-      CleanOwnPrimaryVtx(d,aod,origownvtx);\r
-\r
-    } else {\r
-      // go to selection with Kalman vertexing, if requested\r
-      returnvalueCuts = IsSelectedKF(d,aod);\r
-    }\r
-\r
-    fIsSelectedCuts=returnvalueCuts;\r
-    if(!returnvalueCuts) return 0;\r
-  }\r
-\r
-\r
-  // selection on PID \r
-  if(selectionLevel==AliRDHFCuts::kAll || \r
-     selectionLevel==AliRDHFCuts::kCandidate ||\r
-     selectionLevel==AliRDHFCuts::kPID) {\r
-    returnvaluePID = IsSelectedPID(d);\r
-    fIsSelectedPID=returnvaluePID;\r
-    if(!returnvaluePID) return 0;\r
-  }\r
-\r
-  Int_t returnvalueComb=CombineSelectionLevels(3,returnvalueCuts,returnvaluePID);\r
-\r
-  if(!returnvalueComb) return 0;\r
-\r
-  // selection on daughter tracks \r
-  if(selectionLevel==AliRDHFCuts::kAll || \r
-     selectionLevel==AliRDHFCuts::kTracks) {\r
-    if(!AreDaughtersSelected(d)) return 0;\r
-  }\r
\r
-  //  cout<<"Pid = "<<returnvaluePID<<endl;\r
-  return returnvalueComb;\r
-}\r
-\r
-//------------------------------------------------------------------------------------------\r
-Int_t AliRDHFCutsD0toKpi::IsSelectedKF(AliAODRecoDecayHF2Prong *d,\r
-                                      AliAODEvent* aod) const {\r
-  //\r
-  // Apply selection using KF-vertexing\r
-  //\r
-       \r
-  AliAODTrack *track0 = (AliAODTrack*)d->GetDaughter(0);\r
-  AliAODTrack *track1 = (AliAODTrack*)d->GetDaughter(1); \r
-   \r
-  if(!track0 || !track1) {\r
-    cout<<"one or two D0 daughters missing!"<<endl;\r
-    return 0;\r
-  }\r
-\r
-  // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both\r
-  Int_t returnvalueCuts=3;\r
-        \r
-  // candidate selection with AliKF\r
-  AliKFParticle::SetField(aod->GetMagneticField()); // set the magnetic field\r
-  \r
-  Int_t okD0=0,okD0bar=0;\r
-  okD0=1; okD0bar=1;\r
-  \r
-  // convert tracks into AliKFParticles\r
-  \r
-  AliKFParticle negPiKF(*track1,-211); // neg pion kandidate\r
-  AliKFParticle negKKF(*track1,-321); // neg kaon kandidate\r
-  AliKFParticle posPiKF(*track0,211); // pos pion kandidate\r
-  AliKFParticle posKKF(*track0,321); // pos kaon kandidate\r
-  \r
-  // build D0 candidates\r
-  \r
-  AliKFParticle d0c(negKKF,posPiKF); // D0 candidate\r
-  AliKFParticle ad0c(posKKF,negPiKF); // D0bar candidate\r
-  \r
-  // create kf primary vertices\r
-  \r
-  AliAODVertex *vtx1 = aod->GetPrimaryVertex();\r
-  AliKFVertex primVtx(*vtx1); \r
-  AliKFVertex aprimVtx(*vtx1);\r
-  \r
-  if(primVtx.GetNContributors()<=0) okD0 = 0;\r
-  if(aprimVtx.GetNContributors()<=0) okD0bar = 0;\r
-  if(!okD0 && !okD0bar) returnvalueCuts=0;\r
-       \r
-  // calculate mass\r
-       \r
-  Double_t d0mass = d0c.GetMass();\r
-  Double_t ad0mass = ad0c.GetMass();\r
-       \r
-  // calculate P of D0 and D0bar\r
-  Double_t d0P = d0c.GetP();\r
-  Double_t d0Px = d0c.GetPx();\r
-  Double_t d0Py = d0c.GetPy();\r
-  Double_t d0Pz = d0c.GetPz();\r
-  Double_t ad0P = ad0c.GetP(); \r
-  Double_t ad0Px = ad0c.GetPx();\r
-  Double_t ad0Py = ad0c.GetPy();\r
-  Double_t ad0Pz = ad0c.GetPz();\r
-  \r
-  //calculate Pt of D0 and D0bar\r
-       \r
-  Double_t pt=d0c.GetPt(); \r
-  Double_t apt=ad0c.GetPt();\r
-       \r
-  // remove D0 daughters from primary vertices (if used in vertex fit) and add D0-candidates\r
-  \r
-  if(track0->GetUsedForPrimVtxFit()) {\r
-    primVtx -= posPiKF; \r
-    aprimVtx -= posKKF;\r
-  }\r
-  \r
-  if(track1->GetUsedForPrimVtxFit()) { \r
-    primVtx -= negKKF; \r
-    aprimVtx -= negPiKF;\r
-  }\r
-  \r
-  primVtx += d0c;\r
-  aprimVtx += ad0c;\r
-  \r
-  if(primVtx.GetNContributors()<=0) okD0 = 0;\r
-  if(aprimVtx.GetNContributors()<=0) okD0bar = 0;\r
-  if(!okD0 && !okD0bar) returnvalueCuts=0;\r
-  \r
-  //calculate cut variables\r
-  \r
-  // calculate impact params of daughters w.r.t recalculated vertices\r
-  \r
-  Double_t impactPi = posPiKF.GetDistanceFromVertexXY(primVtx);\r
-  Double_t aimpactPi = negPiKF.GetDistanceFromVertexXY(aprimVtx);\r
-  Double_t impactKa = negKKF.GetDistanceFromVertexXY(primVtx);\r
-  Double_t aimpactKa = posKKF.GetDistanceFromVertexXY(aprimVtx);\r
-       \r
-  // calculate Product of Impact Params\r
-       \r
-  Double_t prodParam = impactPi*impactKa;\r
-  Double_t aprodParam = aimpactPi*aimpactKa;\r
-       \r
-  // calculate cosine of pointing angles\r
-       \r
-  TVector3 mom(d0c.GetPx(),d0c.GetPy(),d0c.GetPz());\r
-  TVector3 fline(d0c.GetX()-primVtx.GetX(),\r
-                d0c.GetY()-primVtx.GetY(),\r
-                d0c.GetZ()-primVtx.GetZ());\r
-  Double_t pta = mom.Angle(fline);\r
-  Double_t cosP = TMath::Cos(pta); // cosine of pta for D0 candidate\r
-  \r
-  TVector3 amom(ad0c.GetPx(),ad0c.GetPy(),ad0c.GetPz());\r
-  TVector3 afline(ad0c.GetX()-aprimVtx.GetX(),\r
-                 ad0c.GetY()-aprimVtx.GetY(),\r
-                 ad0c.GetZ()-aprimVtx.GetZ());\r
-  Double_t apta = amom.Angle(afline);\r
-  Double_t acosP = TMath::Cos(apta); // cosine of pta for D0bar candidate\r
-  \r
-  // calculate P of Pions at Decay Position of D0 and D0bar candidates\r
-  negKKF.TransportToParticle(d0c);\r
-  posPiKF.TransportToParticle(d0c);\r
-  posKKF.TransportToParticle(ad0c);\r
-  negPiKF.TransportToParticle(ad0c);\r
-  \r
-  Double_t pxPi =  posPiKF.GetPx();\r
-  Double_t pyPi =  posPiKF.GetPy();\r
-  Double_t pzPi =  posPiKF.GetPz();\r
-  Double_t ptPi =  posPiKF.GetPt();\r
-  \r
-  Double_t apxPi =  negPiKF.GetPx();\r
-  Double_t apyPi =  negPiKF.GetPy();\r
-  Double_t apzPi =  negPiKF.GetPz();\r
-  Double_t aptPi =  negPiKF.GetPt();\r
-  \r
-  // calculate Pt of Kaons at Decay Position of D0 and D0bar candidates\r
-  \r
-  Double_t ptK =  negKKF.GetPt();\r
-  Double_t aptK =  posKKF.GetPt();\r
-       \r
-  //calculate cos(thetastar)\r
-  Double_t massvtx = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
-  Double_t massp[2];\r
-  massp[0] = TDatabasePDG::Instance()->GetParticle(321)->Mass();\r
-  massp[1] = TDatabasePDG::Instance()->GetParticle(211)->Mass();\r
-  Double_t pStar = TMath::Sqrt(TMath::Power(massvtx*massvtx-massp[0]*massp[0]-massp[1]*massp[1],2.)\r
-                              -4.*massp[0]*massp[0]*massp[1]*massp[1])/(2.*massvtx);\r
-  \r
-  // cos(thetastar) for D0 and Pion\r
-  \r
-  Double_t d0E = TMath::Sqrt(massvtx*massvtx + d0P*d0P);\r
-  Double_t beta = d0P/d0E;\r
-  Double_t gamma = d0E/massvtx;\r
-  TVector3 momPi(pxPi,pyPi,pzPi);\r
-  TVector3 momTot(d0Px,d0Py,d0Pz);\r
-  Double_t q1 = momPi.Dot(momTot)/momTot.Mag();\r
-  Double_t cts = (q1/gamma-beta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar;     \r
-       \r
-  // cos(thetastar) for D0bar and Pion \r
-  \r
-  Double_t ad0E = TMath::Sqrt(massvtx*massvtx + ad0P*ad0P);\r
-  Double_t abeta = ad0P/ad0E;\r
-  Double_t agamma = ad0E/massvtx;\r
-  TVector3 amomPi(apxPi,apyPi,apzPi);\r
-  TVector3 amomTot(ad0Px,ad0Py,ad0Pz);\r
-  Double_t aq1 = amomPi.Dot(amomTot)/amomTot.Mag();\r
-  Double_t acts = (aq1/agamma-abeta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar; \r
-  \r
-  // calculate reduced Chi2 for the full D0 fit\r
-  d0c.SetProductionVertex(primVtx);\r
-  ad0c.SetProductionVertex(aprimVtx);\r
-  negKKF.SetProductionVertex(d0c);\r
-  posPiKF.SetProductionVertex(d0c);\r
-  posKKF.SetProductionVertex(ad0c);\r
-  negPiKF.SetProductionVertex(ad0c);\r
-  d0c.TransportToProductionVertex();\r
-  ad0c.TransportToProductionVertex();\r
-       \r
-  // calculate the decay length\r
-  Double_t decayLengthD0 = d0c.GetDecayLength();\r
-  Double_t adecayLengthD0 = ad0c.GetDecayLength();\r
-  \r
-  Double_t chi2D0 = 50.;\r
-  if(d0c.GetNDF() > 0 && d0c.GetChi2() >= 0) {\r
-    chi2D0 = d0c.GetChi2()/d0c.GetNDF();\r
-  }\r
-  \r
-  Double_t achi2D0 = 50.;\r
-  if(ad0c.GetNDF() > 0 && ad0c.GetChi2() >= 0) {\r
-    achi2D0 = ad0c.GetChi2()/ad0c.GetNDF();\r
-  }\r
-       \r
-  // Get the Pt-bins\r
-  Int_t ptbin=PtBin(pt);\r
-  Int_t aptbin=PtBin(apt);\r
-\r
-  if(ptbin < 0) okD0 = 0;\r
-  if(aptbin < 0) okD0bar = 0;\r
-  if(!okD0 && !okD0bar) returnvalueCuts=0;\r
-  \r
-  if(ptK < fCutsRD[GetGlobalIndex(3,ptbin)] || ptPi < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;\r
-  if(aptK < fCutsRD[GetGlobalIndex(3,aptbin)] || aptPi < fCutsRD[GetGlobalIndex(4,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar) returnvalueCuts=0;\r
-  \r
-  \r
-  if(TMath::Abs(impactKa) > fCutsRD[GetGlobalIndex(5,ptbin)] || \r
-     TMath::Abs(impactPi) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;\r
-  \r
-  if(TMath::Abs(aimpactKa) > fCutsRD[GetGlobalIndex(5,aptbin)] ||\r
-     TMath::Abs(aimpactPi) > fCutsRD[GetGlobalIndex(6,aptbin)]) okD0bar = 0;\r
-  \r
-  if(!okD0 && !okD0bar)  returnvalueCuts=0;\r
-  \r
-  // for the moment via the standard method due to bug in AliKF\r
-  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) okD0 = 0;\r
-  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar) returnvalueCuts=0;\r
-    \r
-    \r
-  if(TMath::Abs(d0mass-massvtx) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;\r
-  if(TMath::Abs(ad0mass-massvtx) > fCutsRD[GetGlobalIndex(0,aptbin)])  okD0bar = 0;\r
-  if(!okD0 && !okD0bar)  returnvalueCuts=0;\r
-  \r
-  \r
-  if(TMath::Abs(cts) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0; \r
-  if(TMath::Abs(acts) > fCutsRD[GetGlobalIndex(2,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar)   returnvalueCuts=0;\r
-  \r
-  if(prodParam  > fCutsRD[GetGlobalIndex(7,ptbin)]) okD0 = 0;\r
-  if(aprodParam > fCutsRD[GetGlobalIndex(7,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar)  returnvalueCuts=0;\r
-    \r
-  if(cosP  < fCutsRD[GetGlobalIndex(8,ptbin)]) okD0 = 0; \r
-  if(acosP < fCutsRD[GetGlobalIndex(8,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar)  returnvalueCuts=0;\r
-       \r
-  if(chi2D0  > fCutsRD[GetGlobalIndex(10,ptbin)]) okD0 = 0; \r
-  if(achi2D0 > fCutsRD[GetGlobalIndex(10,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar)  returnvalueCuts=0;\r
-       \r
-  if(decayLengthD0  < fCutsRD[GetGlobalIndex(9,ptbin)]) okD0 = 0; \r
-  if(adecayLengthD0 < fCutsRD[GetGlobalIndex(9,aptbin)]) okD0bar = 0;\r
-  if(!okD0 && !okD0bar)  returnvalueCuts=0;\r
-    \r
-  if(returnvalueCuts!=0) {\r
-    if(okD0) returnvalueCuts=1; //cuts passed as D0\r
-    if(okD0bar) returnvalueCuts=2; //cuts passed as D0bar\r
-    if(okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar\r
-  }\r
-\r
-  return returnvalueCuts;  \r
-}\r
-\r
-//---------------------------------------------------------------------------\r
-\r
-Bool_t AliRDHFCutsD0toKpi::IsInFiducialAcceptance(Double_t pt, Double_t y) const\r
-{\r
-  //\r
-  // Checking if D0 is in fiducial acceptance region \r
-  //\r
-\r
-  if(pt > 5.) {\r
-    // applying cut for pt > 5 GeV\r
-    AliDebug(2,Form("pt of D0 = %f (> 5), cutting at |y| < 0.8\n",pt)); \r
-    if (TMath::Abs(y) > 0.8){\r
-      return kFALSE;\r
-    }\r
-  } else {\r
-    // appliying smooth cut for pt < 5 GeV\r
-    Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5; \r
-    Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;                \r
-    AliDebug(2,Form("pt of D0 = %f (< 5), cutting  according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY)); \r
-    if (y < minFiducialY || y > maxFiducialY){\r
-      return kFALSE;\r
-    }\r
-  }\r
-\r
-  return kTRUE;\r
-}\r
-//---------------------------------------------------------------------------\r
-Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d) \r
-{\r
-  // ############################################################\r
-  //\r
-  // Apply PID selection\r
-  //\r
-  //\r
-  // ############################################################\r
-\r
-  if(!fUsePID) return 3;\r
-  if(fDefaultPID) return IsSelectedPIDdefault(d);\r
-  fWhyRejection=0;\r
-  Int_t isD0D0barPID[2]={1,2};\r
-  Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value; \r
-  //                                                                                                 same for prong 2\r
-  //                                               values convention -1 = discarded \r
-  //                                                                  0 = not identified (but compatible) || No PID (->hasPID flag)\r
-  //                                                                  1 = identified\r
-  // PID search:   pion (TPC) or not K (TOF), Kaon hypothesis for both \r
-  // Initial hypothesis: unknwon (but compatible) \r
-  combinedPID[0][0]=0;  // prima figlia, Kaon\r
-  combinedPID[0][1]=0;  // prima figlia, pione\r
-  combinedPID[1][0]=0;  // seconda figlia, Kaon\r
-  combinedPID[1][1]=0;  // seconda figlia, pion\r
-  \r
-  Bool_t checkPIDInfo[2]={kTRUE,kTRUE};\r
-  Double_t sigma_tmp[3]={fPidHF->GetSigma(0),fPidHF->GetSigma(1),fPidHF->GetSigma(2)};\r
-  for(Int_t daught=0;daught<2;daught++){\r
-    //Loop con prongs\r
-    AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);\r
-    if(fPidHF->IsTOFPiKexcluded(aodtrack,5.)) return 0; \r
-    \r
-    if(!(fPidHF->CheckStatus(aodtrack,"TPC")) && !(fPidHF->CheckStatus(aodtrack,"TOF"))) {\r
-      checkPIDInfo[daught]=kFALSE; \r
-      continue;\r
-    }\r
-\r
-    // identify kaon\r
-    combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);\r
-\r
-    // identify pion\r
-\r
-    if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {\r
-     combinedPID[daught][1]=0;\r
-    }else{\r
-      fPidHF->SetTOF(kFALSE);\r
-      combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);\r
-      fPidHF->SetTOF(kTRUE);\r
-      fPidHF->SetCompat(kTRUE);\r
-     }\r
-\r
-\r
-    if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){ // if not a K- and not a pi- both D0 and D0bar excluded\r
-      isD0D0barPID[0]=0;\r
-      isD0D0barPID[1]=0;\r
-    }\r
-    else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){\r
-      if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded\r
-      else isD0D0barPID[0]=0;// if K+ D0 excluded\r
-    }\r
-    /*    else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){\r
-         isD0D0barPID[0]=0;\r
-         isD0D0barPID[1]=0;\r
-         }\r
-    */\r
-    else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){ \r
-      if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded\r
-      else isD0D0barPID[0]=0;//  not a D0 if K+ or if pi+ excluded\r
-    }\r
-    else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){\r
-      if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded\r
-      else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded\r
-    }\r
-\r
-    if(fLowPt && d->Pt()<fPtLowPID){\r
-     Double_t sigmaTPC[3]={3.,2.,0.};\r
-     fPidHF->SetSigmaForTPC(sigmaTPC);\r
-    // identify kaon\r
-    combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);\r
-\r
-    Double_t ptProng=aodtrack->P();\r
-\r
-    if(ptProng<0.6){\r
-     fPidHF->SetCompat(kFALSE);\r
-     combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);\r
-     fPidHF->SetCompat(kTRUE);\r
-    }\r
-\r
-    if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {\r
-     combinedPID[daught][1]=0;\r
-    }else{\r
-      fPidHF->SetTOF(kFALSE);\r
-      Double_t sigmaTPCpi[3]={3.,3.,0.};\r
-      fPidHF->SetSigmaForTPC(sigmaTPCpi);\r
-      combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);\r
-      fPidHF->SetTOF(kTRUE);\r
-       if(ptProng<0.8){\r
-        Bool_t isTPCpion=fPidHF->IsPionRaw(aodtrack,"TPC");\r
-        if(isTPCpion){\r
-         combinedPID[daught][1]=1;\r
-        }else{\r
-         combinedPID[daught][1]=-1;\r
-        }\r
-      }\r
-    }\r
-\r
-   }\r
-   fPidHF->SetSigmaForTPC(sigma_tmp);\r
-  }// END OF LOOP ON DAUGHTERS\r
-\r
-   if(!checkPIDInfo[0] && !checkPIDInfo[1]) {\r
-    if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);\r
-    return 0;\r
-   }\r
-\r
-\r
-  // FURTHER PID REQUEST (both daughter info is needed)\r
-  if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){\r
-    fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found\r
-    if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);\r
-    return 0;\r
-  }\r
-\r
-  if(fLowPt && d->Pt()<fPtLowPID){    \r
-    if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){\r
-      fWhyRejection=32;// reject cases where the Kaon is not identified\r
-      fPidHF->SetSigmaForTPC(sigma_tmp);\r
-      return 0;\r
-    }\r
-  }\r
-    if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);\r
-\r
-  //  cout<<"Why? "<<fWhyRejection<<endl;  \r
-  return isD0D0barPID[0]+isD0D0barPID[1];\r
-}\r
-//---------------------------------------------------------------------------\r
-Int_t AliRDHFCutsD0toKpi::IsSelectedPIDdefault(AliAODRecoDecayHF* d) \r
-{\r
-  // ############################################################\r
-  //\r
-  // Apply PID selection\r
-  //\r
-  //\r
-  // temporary selection: PID AS USED FOR D0 by Andrea Rossi (up to 28/06/2010)\r
-  //\r
-  // d must be a AliAODRecoDecayHF2Prong object\r
-  // returns 0 if both D0 and D0bar are rejectecd\r
-  //         1 if D0 is accepted while D0bar is rejected\r
-  //         2 if D0bar is accepted while D0 is rejected\r
-  //         3 if both are accepted\r
-  // fWhyRejection variable (not returned for the moment, print it if needed)\r
-  //               keeps some information on why a candidate has been \r
-  //               rejected according to the following (unfriendly?) scheme \r
-  //             if more rejection cases are considered interesting, just add numbers\r
-  //\r
-  //      TO BE CONSIDERED WITH A GRAIN OF SALT (the order in which cut are applied is relevant) \r
-  //              from 20 to 30: "detector" selection (PID acceptance)                                             \r
-  //                                                  26: TPC refit\r
-  //                                                  27: ITS refit\r
-  //                                                  28: no (TOF||TPC) pid information (no kTOFpid,kTOFout,kTIME,kTPCpid,...)\r
-  //\r
-  //              from 30 to 40: PID selection\r
-  //                                                  31: no Kaon compatible tracks found between daughters\r
-  //                                                  32: no Kaon identified tracks found (strong sel. at low momenta)\r
-  //                                                  33: both mass hypotheses are rejected \r
-  //                  \r
-  // ############################################################\r
-\r
-  if(!fUsePID) return 3;\r
-  fWhyRejection=0;\r
-  Int_t isD0D0barPID[2]={1,2};\r
-  Double_t nsigmaTPCpi=-1., nsigmaTPCK=-1.; //used for TPC pid\r
-  Double_t tofSig,times[5];// used fot TOF pid\r
-  Int_t hasPID[2]={2,2};// flag to count how many detectors give PID info for the daughters\r
-  Int_t isKaonPionTOF[2][2],isKaonPionTPC[2][2];\r
-  Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value; \r
-  //                                                                                                 same for prong 2\r
-  //                                               values convention -1 = discarded \r
-  //                                                                  0 = not identified (but compatible) || No PID (->hasPID flag)\r
-  //                                                                  1 = identified\r
-  // PID search:   pion (TPC) or not K (TOF), Kaon hypothesis for both \r
-  // Initial hypothesis: unknwon (but compatible) \r
-  isKaonPionTOF[0][0]=0;\r
-  isKaonPionTOF[0][1]=0;\r
-  isKaonPionTOF[1][0]=0;\r
-  isKaonPionTOF[1][1]=0;\r
-  \r
-  isKaonPionTPC[0][0]=0;\r
-  isKaonPionTPC[0][1]=0;\r
-  isKaonPionTPC[1][0]=0;\r
-  isKaonPionTPC[1][1]=0;\r
-  \r
-  combinedPID[0][0]=0;\r
-  combinedPID[0][1]=0;\r
-  combinedPID[1][0]=0;\r
-  combinedPID[1][1]=0;\r
-  \r
-  \r
\r
-  for(Int_t daught=0;daught<2;daught++){\r
-    //Loop con prongs\r
-    \r
-    // ########### Step 0- CHECKING minimal PID "ACCEPTANCE" ####################\r
-\r
-    AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught); \r
-   \r
-    if(!(aodtrack->GetStatus()&AliESDtrack::kTPCrefit)){\r
-      fWhyRejection=26;\r
-      return 0;\r
-    } \r
-    if(!(aodtrack->GetStatus()&AliESDtrack::kITSrefit)){\r
-      fWhyRejection=27;\r
-      return 0;\r
-    } \r
-    \r
-    AliAODPid *pid=aodtrack->GetDetPid();\r
-    if(!pid) {\r
-      //delete esdtrack;\r
-      hasPID[daught]--;\r
-      continue;\r
-    }\r
-  \r
-    // ########### Step 1- Check of TPC and TOF response ####################\r
-\r
-    Double_t ptrack=aodtrack->P();\r
-    //#################### TPC PID #######################\r
-     if (!(aodtrack->GetStatus()&AliESDtrack::kTPCpid )){\r
-       // NO TPC PID INFO FOR THIS TRACK \r
-       hasPID[daught]--;\r
-     }\r
-     else {\r
-       static AliTPCPIDResponse theTPCpid;\r
-       AliAODPid *pidObj = aodtrack->GetDetPid();\r
-       Double_t ptProng=pidObj->GetTPCmomentum();\r
-       nsigmaTPCpi = theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kPion);\r
-       nsigmaTPCK =  theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kKaon);\r
-       //if(ptrack<0.6){\r
-       if(ptProng<0.6){\r
-        if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;\r
-        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;\r
-        if(TMath::Abs(nsigmaTPCpi)<2.)isKaonPionTPC[daught][1]=1;\r
-        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;\r
-       }\r
-       //else if(ptrack<.8){\r
-       else if(ptProng<.8){\r
-        if(TMath::Abs(nsigmaTPCK)<1.)isKaonPionTPC[daught][0]=1;\r
-        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;\r
-        if(TMath::Abs(nsigmaTPCpi)<1.)isKaonPionTPC[daught][1]=1;\r
-        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;\r
-       }     \r
-       else {\r
-        //     if(nsigmaTPCK>-2.&&nsigmaTPCK<1.)isKaonPionTPC[daught][0]=1;\r
-        if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;\r
-        //if(nsigmaTPCpi>-1.&&nsigmaTPCpi<2.)isKaonPionTPC[daught][1]=1;\r
-        if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;\r
-       }\r
-     }\r
-    \r
-    \r
-    // ##### TOF PID: do not ask nothing for pion/protons ############\r
-     if(!((aodtrack->GetStatus()&AliESDtrack::kTOFpid)&&(aodtrack->GetStatus()&AliESDtrack::kTOFout)&&(aodtrack->GetStatus()&AliESDtrack::kTIME))){\r
-       // NO TOF PID INFO FOR THIS TRACK      \r
-       hasPID[daught]--;\r
-     }\r
-     else{\r
-       tofSig=pid->GetTOFsignal(); \r
-       pid->GetIntegratedTimes(times);\r
-       if((tofSig-times[3])>5.*160.)return 0;// PROTON REJECTION\r
-       if(TMath::Abs(tofSig-times[3])>3.*160.){\r
-        isKaonPionTOF[daught][0]=-1;\r
-       }\r
-       else {   \r
-        if(ptrack<1.5){\r
-          isKaonPionTOF[daught][0]=1;\r
-        }\r
-       }\r
-     }\r
-     \r
-     //######### Step 2: COMBINE TOF and TPC PID ###############\r
-     // we apply the following convention: if TPC and TOF disagree (discarded Vs identified) -> unknown\r
-     combinedPID[daught][0]=isKaonPionTOF[daught][0]+isKaonPionTPC[daught][0];\r
-     combinedPID[daught][1]=isKaonPionTOF[daught][1]+isKaonPionTPC[daught][1];\r
-     \r
-     \r
-     //######### Step 3:   USE PID INFO     \r
-     \r
-     if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){// if not a K- and not a pi- both D0 and D0bar excluded\r
-       isD0D0barPID[0]=0;\r
-       isD0D0barPID[1]=0;\r
-     }\r
-     else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){// if in conflict (both pi- and K-), if k for both TPC and TOF -> is K\r
-       if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded\r
-       else isD0D0barPID[0]=0;// if K+ D0 excluded\r
-     }\r
-     else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){// if in conflict (both pi- and K-) and k- only for TPC or TOF -> reject\r
-       isD0D0barPID[0]=0;\r
-       isD0D0barPID[1]=0;\r
-     }\r
-     else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){\r
-       if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded\r
-       else isD0D0barPID[0]=0;//  not a D0 if K+ or if pi+ excluded\r
-     }\r
-     else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){\r
-       if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded\r
-      else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded\r
-     }\r
-     \r
-     // ##########  ALSO DIFFERENT TPC PID REQUEST FOR LOW pt D0: request of K identification      ###############################\r
-     // ########## more tolerant criteria for single particle ID-> more selective criteria for D0   ##############################\r
-     // ###############                     NOT OPTIMIZED YET                                  ###################################\r
-     if(d->Pt()<2.){\r
-       isKaonPionTPC[daught][0]=0;\r
-       isKaonPionTPC[daught][1]=0;\r
-       AliAODPid *pidObj = aodtrack->GetDetPid();\r
-       Double_t ptProng=pidObj->GetTPCmomentum();\r
-       //if(ptrack<0.6){\r
-       if(ptProng<0.6){\r
-        if(TMath::Abs(nsigmaTPCK)<3.)isKaonPionTPC[daught][0]=1;\r
-        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;\r
-        if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;\r
-        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;\r
-     }\r
-       //else if(ptrack<.8){\r
-       else if(ptProng<.8){\r
-        if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;\r
-        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;\r
-        if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;\r
-        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;\r
-       }     \r
-       else {\r
-        if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;\r
-        if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;\r
-       }\r
-     }\r
-     \r
-  }// END OF LOOP ON DAUGHTERS\r
-  \r
-  // FURTHER PID REQUEST (both daughter info is needed)\r
-  if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){\r
-    fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found\r
-    return 0;\r
-  }\r
-  else if(hasPID[0]==0&&hasPID[1]==0){\r
-    fWhyRejection=28;// reject cases in which no PID info is available  \r
-    return 0;\r
-  }\r
-  if(d->Pt()<2.){\r
-    // request of K identification at low D0 pt\r
-    combinedPID[0][0]=0;\r
-    combinedPID[0][1]=0;\r
-    combinedPID[1][0]=0;\r
-    combinedPID[1][1]=0;\r
-    \r
-    combinedPID[0][0]=isKaonPionTOF[0][0]+isKaonPionTPC[0][0];\r
-    combinedPID[0][1]=isKaonPionTOF[0][1]+isKaonPionTPC[0][1];\r
-    combinedPID[1][0]=isKaonPionTOF[1][0]+isKaonPionTPC[1][0];\r
-    combinedPID[1][1]=isKaonPionTOF[1][1]+isKaonPionTPC[1][1];\r
-    \r
-    if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){\r
-      fWhyRejection=32;// reject cases where the Kaon is not identified\r
-      return 0;\r
-    }\r
-  }\r
-\r
-  //  cout<<"Why? "<<fWhyRejection<<endl;  \r
-  return isD0D0barPID[0]+isD0D0barPID[1];\r
-}\r
-\r
-\r
-\r
-//---------------------------------------------------------------------------\r
-Int_t AliRDHFCutsD0toKpi::CombineSelectionLevels(Int_t selectionvalTrack,\r
-                                                Int_t selectionvalCand,\r
-                                                Int_t selectionvalPID) const\r
-{\r
-  //\r
-  // This method combines the tracks, PID and cuts selection results\r
-  //\r
-  if(selectionvalTrack==0) return 0;\r
-\r
-  Int_t returnvalue;\r
-\r
-  switch(selectionvalPID) {\r
-  case 0:\r
-    returnvalue=0;\r
-    break;\r
-  case 1:\r
-    returnvalue=((selectionvalCand==1 || selectionvalCand==3) ? 1 : 0);\r
-    break;\r
-  case 2:\r
-    returnvalue=((selectionvalCand==2 || selectionvalCand==3) ? 2 : 0);\r
-    break;\r
-  case 3:\r
-    returnvalue=selectionvalCand;\r
-    break;\r
-  default:\r
-    returnvalue=0;\r
-    break;\r
-  }\r
-\r
-  return returnvalue;\r
-}\r
-//----------------------------------------------------------------------------\r
-Int_t AliRDHFCutsD0toKpi::IsSelectedSpecialCuts(AliAODRecoDecayHF *d) const\r
-{\r
-  //\r
-  // Note: this method is temporary\r
-  // Additional cuts on decay lenght and lower cut for d0 norm are applied using vertex without candidate's daughters\r
-  //\r
-\r
-  //apply cuts\r
-\r
-  Float_t normDecLengthCut=1.,decLengthCut=TMath::Min(d->P()*0.0066+0.01,0.06/*cm*/), normd0Cut=0.5;\r
-  // "decay length" expo law with tau' = beta*gamma*ctau= p/m*ctau =p*0.0123/1.864~p*0.0066\r
-  // decay lenght > ctau' implies to retain (1-1/e) (for signal without considering detector resolution), \r
-\r
-  Int_t returnvalue=3; //cut passed\r
-  for(Int_t i=0;i<2/*prongs*/;i++){\r
-    if(TMath::Abs(d->Normalizedd0Prong(i))<normd0Cut) return 0; //normd0Cut not passed\r
-  }\r
-  if(d->DecayLength2()<decLengthCut*decLengthCut)  return 0; //decLengthCut not passed\r
-  if(d->NormalizedDecayLength2()<normDecLengthCut*normDecLengthCut)  return 0; //decLengthCut not passed\r
-       \r
-  return returnvalue;\r
-}\r
-\r
-//----------------------------------------------\r
-void AliRDHFCutsD0toKpi::SetUseKF(Bool_t useKF)\r
-{\r
-  //switch on candidate selection via AliKFparticle\r
-  if(!useKF) return;\r
-  if(useKF){\r
-    fUseKF=useKF;\r
-    Int_t nvarsKF=11;\r
-    SetNVars(nvarsKF);\r
-    TString varNamesKF[11]={"inv. mass [GeV]",   \r
-                           "dca [cm]",\r
-                           "cosThetaStar", \r
-                           "pTK [GeV/c]",\r
-                           "pTPi [GeV/c]",\r
-                           "d0K [cm]",\r
-                           "d0Pi [cm]",\r
-                           "d0d0 [cm^2]",\r
-                           "cosThetaPoint"\r
-                           "DecayLength[cm]",\r
-                           "RedChi2"};\r
-    Bool_t isUpperCutKF[11]={kTRUE,\r
-                            kTRUE,\r
-                            kTRUE,\r
-                            kFALSE,\r
-                            kFALSE,\r
-                            kTRUE,\r
-                            kTRUE,\r
-                            kTRUE,\r
-                            kFALSE,\r
-                            kFALSE,\r
-                            kTRUE};\r
-    SetVarNames(nvarsKF,varNamesKF,isUpperCutKF);\r
-    SetGlobalIndex();\r
-    Bool_t forOpt[11]={kFALSE,\r
-                   kTRUE,\r
-                   kTRUE,\r
-                   kFALSE,\r
-                   kFALSE,\r
-                   kFALSE,\r
-                   kFALSE,\r
-                   kTRUE,\r
-                   kTRUE,\r
-                   kFALSE,\r
-                   kFALSE};\r
-    SetVarsForOpt(4,forOpt);\r
-  }\r
-  return;\r
-}\r
-\r
-\r
-void AliRDHFCutsD0toKpi::SetStandardCutsPP2010() {\r
-  //\r
-  //STANDARD CUTS USED FOR 2010 pp analysis \r
-  //dca cut will be enlarged soon to 400 micron\r
-  //\r
-  \r
-  SetName("D0toKpiCutsStandard");\r
-  SetTitle("Standard Cuts for D0 analysis");\r
-  \r
-  // PILE UP REJECTION\r
-  SetOptPileup(AliRDHFCuts::kRejectPileupEvent);\r
-\r
-  // EVENT CUTS\r
-  SetMinVtxContr(1);\r
- // MAX Z-VERTEX CUT\r
-  SetMaxVtxZ(10.);\r
-  \r
-  // TRACKS ON SINGLE TRACKS\r
-  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");\r
-  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);\r
-  esdTrackCuts->SetRequireTPCRefit(kTRUE);\r
-  esdTrackCuts->SetRequireITSRefit(kTRUE);\r
-  //  esdTrackCuts->SetMinNClustersITS(4);\r
-  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);\r
-  esdTrackCuts->SetMinDCAToVertexXY(0.);\r
-  esdTrackCuts->SetEtaRange(-0.8,0.8);\r
-  esdTrackCuts->SetPtRange(0.3,1.e10);\r
-  \r
-  AddTrackCuts(esdTrackCuts);\r
-\r
-  \r
-  const Int_t nptbins =14;\r
-  const Double_t ptmax = 9999.;\r
-  const Int_t nvars=11;\r
-  Float_t ptbins[nptbins+1];\r
-  ptbins[0]=0.;\r
-  ptbins[1]=0.5;       \r
-  ptbins[2]=1.;\r
-  ptbins[3]=2.;\r
-  ptbins[4]=3.;\r
-  ptbins[5]=4.;\r
-  ptbins[6]=5.;\r
-  ptbins[7]=6.;\r
-  ptbins[8]=7.;\r
-  ptbins[9]=8.;\r
-  ptbins[10]=12.;\r
-  ptbins[11]=16.;\r
-  ptbins[12]=20.;\r
-  ptbins[13]=24.;\r
-  ptbins[14]=ptmax;\r
-\r
-  SetGlobalIndex(nvars,nptbins);\r
-  SetPtBins(nptbins+1,ptbins);\r
-  \r
-  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* pt<0.5*/\r
-                                                 {0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* 0.5<pt<1*/\r
-                                                 {0.400,300.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.80,0.,0.},/* 1<pt<2 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.,0.},/* 2<pt<3 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 3<pt<4 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 4<pt<5 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 5<pt<6 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 6<pt<7 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-7000.*1E-8,0.85,0.,0.},/* 7<pt<8 */\r
-                                                 {0.400,300.*1E-4,0.9,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.,0.},/* 8<pt<12 */\r
-                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,10000.*1E-8,0.85,0.,0.},/* 12<pt<16 */\r
-                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 16<pt<20 */\r
-                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 20<pt<24 */\r
-                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.}};/* pt>24 */\r
-  \r
-  \r
-  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts\r
-  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];\r
-  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];\r
-  \r
-  for (Int_t ibin=0;ibin<nptbins;ibin++){\r
-    for (Int_t ivar = 0; ivar<nvars; ivar++){\r
-      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      \r
-    }\r
-  }\r
-  \r
-  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);\r
-  SetUseSpecialCuts(kTRUE);\r
-  SetRemoveDaughtersFromPrim(kTRUE);\r
-  \r
-  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];\r
-  delete [] cutsMatrixTransposeStand;\r
-  cutsMatrixTransposeStand=NULL;\r
-\r
-  // PID SETTINGS\r
-  AliAODPidHF* pidObj=new AliAODPidHF();\r
-  //pidObj->SetName("pid4D0");\r
-  Int_t mode=1;\r
-  const Int_t nlims=2;\r
-  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]\r
-  Bool_t compat=kTRUE; //effective only for this mode\r
-  Bool_t asym=kTRUE;\r
-  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella\r
-  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC\r
-  pidObj->SetMatch(mode);\r
-  pidObj->SetPLimit(plims,nlims);\r
-  pidObj->SetSigma(sigmas);\r
-  pidObj->SetCompat(compat);\r
-  pidObj->SetTPC(kTRUE);\r
-  pidObj->SetTOF(kTRUE);\r
-  pidObj->SetPCompatTOF(1.5);\r
-  pidObj->SetSigmaForTPCCompat(3.);\r
-  pidObj->SetSigmaForTOFCompat(3.);\r
-  \r
-\r
-  SetPidHF(pidObj);\r
-  SetUsePID(kTRUE);\r
-  SetUseDefaultPID(kFALSE);\r
-\r
-  SetLowPt(kFALSE);\r
-\r
-  PrintAll();\r
-\r
-  delete pidObj;\r
-  pidObj=NULL;\r
-\r
-  return;\r
-\r
-}\r
-\r
-\r
-void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010() {\r
-  //\r
-  // Final CUTS USED FOR 2010 PbPb analysis of central events\r
-  // REMEMBER TO EVENTUALLY SWITCH ON \r
-  //          SetUseAOD049(kTRUE);\r
-  // \r
-  \r
-  SetName("D0toKpiCutsStandard");\r
-  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run");\r
-  \r
-  \r
-  // CENTRALITY SELECTION\r
-  SetMinCentrality(0.);\r
-  SetMaxCentrality(80.);\r
-  SetUseCentrality(AliRDHFCuts::kCentV0M);\r
-\r
-\r
-  \r
-  // EVENT CUTS\r
-  SetMinVtxContr(1);\r
-  // MAX Z-VERTEX CUT\r
-  SetMaxVtxZ(10.);\r
-  // PILE UP\r
-  SetOptPileup(AliRDHFCuts::kNoPileupSelection);\r
-  // PILE UP REJECTION\r
-  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);\r
-\r
-  // TRACKS ON SINGLE TRACKS\r
-  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");\r
-  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);\r
-  esdTrackCuts->SetRequireTPCRefit(kTRUE);\r
-  esdTrackCuts->SetRequireITSRefit(kTRUE);\r
-  //  esdTrackCuts->SetMinNClustersITS(4);\r
-  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);\r
-  esdTrackCuts->SetMinDCAToVertexXY(0.);\r
-  esdTrackCuts->SetEtaRange(-0.8,0.8);\r
-  esdTrackCuts->SetPtRange(0.7,1.e10);\r
-\r
-  esdTrackCuts->SetMaxDCAToVertexXY(1.);  \r
-  esdTrackCuts->SetMaxDCAToVertexZ(1.);\r
-  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0075*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");  \r
-\r
-\r
-  AddTrackCuts(esdTrackCuts);\r
-  // SPD k FIRST for Pt<3 GeV/c\r
-  SetSelectCandTrackSPDFirst(kTRUE, 3); \r
-\r
-  // CANDIDATE CUTS  \r
-  const Int_t nptbins =13;\r
-  const Double_t ptmax = 9999.;\r
-  const Int_t nvars=11;\r
-  Float_t ptbins[nptbins+1];\r
-  ptbins[0]=0.;\r
-  ptbins[1]=0.5;       \r
-  ptbins[2]=1.;\r
-  ptbins[3]=2.;\r
-  ptbins[4]=3.;\r
-  ptbins[5]=4.;\r
-  ptbins[6]=5.;\r
-  ptbins[7]=6.;\r
-  ptbins[8]=8.;\r
-  ptbins[9]=12.;\r
-  ptbins[10]=16.;\r
-  ptbins[11]=20.;\r
-  ptbins[12]=24.;\r
-  ptbins[13]=ptmax;\r
-\r
-  SetGlobalIndex(nvars,nptbins);\r
-  SetPtBins(nptbins+1,ptbins);\r
-  SetMinPtCandidate(2.);\r
-\r
-  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.85,0.99,2.},/* pt<0.5*/\r
-                                                 {0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-35000.*1E-8,0.9,0.99,2.},/* 0.5<pt<1*/\r
-                                                 {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-43000.*1E-8,0.85,0.995,2.},/* 1<pt<2 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-45000.*1E-8,0.95,0.998,7.},/* 2<pt<3 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-36000.*1E-8,0.95,0.998,5.},/* 3<pt<4 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-27000.*1E-8,0.95,0.998,5.},/* 4<pt<5 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-21000.*1E-8,0.92,0.998,5.},/* 5<pt<6 */\r
-                                                 {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-14000.*1E-8,0.88,0.998,5.},/* 6<pt<8 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.998,5.},/* 8<pt<12 */\r
-                                                 {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.83,0.998,8.},/* 12<pt<16 */\r
-                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.82,0.995,6.},/* 16<pt<20 */\r
-                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.81,0.995,6.},/* 20<pt<24 */\r
-                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.8,0.99,2.}};/* pt>24 */\r
-  \r
-  \r
-  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts\r
-  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];\r
-  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];\r
-  \r
-  for (Int_t ibin=0;ibin<nptbins;ibin++){\r
-    for (Int_t ivar = 0; ivar<nvars; ivar++){\r
-      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      \r
-    }\r
-  }\r
-  \r
-  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);\r
-  SetUseSpecialCuts(kTRUE);\r
-  SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb\r
-  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];\r
-  delete [] cutsMatrixTransposeStand;\r
-  cutsMatrixTransposeStand=NULL;\r
-  \r
-  // PID SETTINGS\r
-  AliAODPidHF* pidObj=new AliAODPidHF();\r
-  //pidObj->SetName("pid4D0");\r
-  Int_t mode=1;\r
-  const Int_t nlims=2;\r
-  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]\r
-  Bool_t compat=kTRUE; //effective only for this mode\r
-  Bool_t asym=kTRUE;\r
-  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella\r
-  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC\r
-  pidObj->SetMatch(mode);\r
-  pidObj->SetPLimit(plims,nlims);\r
-  pidObj->SetSigma(sigmas);\r
-  pidObj->SetCompat(compat);\r
-  pidObj->SetTPC(kTRUE);\r
-  pidObj->SetTOF(kTRUE);\r
-  pidObj->SetPCompatTOF(2.);\r
-  pidObj->SetSigmaForTPCCompat(3.);\r
-  pidObj->SetSigmaForTOFCompat(3.);  \r
-\r
-\r
-  SetPidHF(pidObj);\r
-  SetUsePID(kTRUE);\r
-  SetUseDefaultPID(kFALSE);\r
-\r
-\r
-  PrintAll();\r
-\r
-\r
-  delete pidObj;\r
-  pidObj=NULL;\r
-\r
-  return;\r
-\r
-}\r
-\r
-void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010Peripherals() {\r
-  // CUTS USED FOR D0 RAA for the CENTRALITY RANGE 20-80%\r
-\r
-  \r
-  SetName("D0toKpiCutsStandard");\r
-  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run, for peripheral events");\r
-  \r
-  \r
-  // CENTRALITY SELECTION\r
-  SetMinCentrality(40.);\r
-  SetMaxCentrality(80.);\r
-  SetUseCentrality(AliRDHFCuts::kCentV0M);\r
-  \r
-  // EVENT CUTS\r
-  SetMinVtxContr(1);\r
-  // MAX Z-VERTEX CUT\r
-  SetMaxVtxZ(10.);\r
-  // PILE UP\r
-  SetOptPileup(AliRDHFCuts::kNoPileupSelection);\r
-  // PILE UP REJECTION\r
-  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);\r
-\r
-  // TRACKS ON SINGLE TRACKS\r
-  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");\r
-  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);\r
-  esdTrackCuts->SetRequireTPCRefit(kTRUE);\r
-  esdTrackCuts->SetRequireITSRefit(kTRUE);\r
-  //  esdTrackCuts->SetMinNClustersITS(4);\r
-  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);\r
-  esdTrackCuts->SetMinDCAToVertexXY(0.);\r
-  esdTrackCuts->SetEtaRange(-0.8,0.8);\r
-  esdTrackCuts->SetPtRange(0.5,1.e10);\r
-\r
-  esdTrackCuts->SetMaxDCAToVertexXY(1.);  \r
-  esdTrackCuts->SetMaxDCAToVertexZ(1.);\r
-  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0025*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");  \r
-\r
-\r
-  AddTrackCuts(esdTrackCuts);\r
-  // SPD k FIRST for Pt<3 GeV/c\r
-  SetSelectCandTrackSPDFirst(kTRUE, 3); \r
-\r
-  // CANDIDATE CUTS  \r
-  const Int_t nptbins =13;\r
-  const Double_t ptmax = 9999.;\r
-  const Int_t nvars=11;\r
-  Float_t ptbins[nptbins+1];\r
-  ptbins[0]=0.;\r
-  ptbins[1]=0.5;       \r
-  ptbins[2]=1.;\r
-  ptbins[3]=2.;\r
-  ptbins[4]=3.;\r
-  ptbins[5]=4.;\r
-  ptbins[6]=5.;\r
-  ptbins[7]=6.;\r
-  ptbins[8]=8.;\r
-  ptbins[9]=12.;\r
-  ptbins[10]=16.;\r
-  ptbins[11]=20.;\r
-  ptbins[12]=24.;\r
-  ptbins[13]=ptmax;\r
-\r
-  SetGlobalIndex(nvars,nptbins);\r
-  SetPtBins(nptbins+1,ptbins);\r
-  SetMinPtCandidate(0.);\r
-\r
-  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-20000.*1E-8,0.7,0.993,2.},/* pt<0.5*/\r
-                                                 {0.400,400.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.85,0.993,2.},/* 0.5<pt<1*/\r
-                                                 {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.85,0.995,6.},/* 1<pt<2 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.95,0.991,5.},/* 2<pt<3 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-36000.*1E-8,0.95,0.993,5.},/* 3<pt<4 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-27000.*1E-8,0.95,0.995,5.},/* 4<pt<5 */\r
-                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-21000.*1E-8,0.92,0.998,5.},/* 5<pt<6 */\r
-                                                 {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-14000.*1E-8,0.88,0.998,5.},/* 6<pt<8 */\r
-                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.995,5.},/* 8<pt<12 */\r
-                                                 {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.83,0.995,4.},/* 12<pt<16 */\r
-                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.82,0.995,4.},/* 16<pt<20 */\r
-                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.81,0.995,4.},/* 20<pt<24 */\r
-                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.8,0.995,4.}};/* pt>24 */\r
-  \r
-  \r
-  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts\r
-  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];\r
-  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];\r
-  \r
-  for (Int_t ibin=0;ibin<nptbins;ibin++){\r
-    for (Int_t ivar = 0; ivar<nvars; ivar++){\r
-      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      \r
-    }\r
-  }\r
-  \r
-  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);\r
-  SetUseSpecialCuts(kTRUE);\r
-  SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb\r
-  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];\r
-  delete [] cutsMatrixTransposeStand;\r
-  cutsMatrixTransposeStand=NULL;\r
-  \r
-  // PID SETTINGS\r
-  AliAODPidHF* pidObj=new AliAODPidHF();\r
-  //pidObj->SetName("pid4D0");\r
-  Int_t mode=1;\r
-  const Int_t nlims=2;\r
-  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]\r
-  Bool_t compat=kTRUE; //effective only for this mode\r
-  Bool_t asym=kTRUE;\r
-  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella\r
-  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC\r
-  pidObj->SetMatch(mode);\r
-  pidObj->SetPLimit(plims,nlims);\r
-  pidObj->SetSigma(sigmas);\r
-  pidObj->SetCompat(compat);\r
-  pidObj->SetTPC(kTRUE);\r
-  pidObj->SetTOF(kTRUE);\r
-  pidObj->SetPCompatTOF(2.);\r
-  pidObj->SetSigmaForTPCCompat(3.);\r
-  pidObj->SetSigmaForTOFCompat(3.);  \r
-\r
-  SetPidHF(pidObj);\r
-  SetUsePID(kTRUE);\r
-  SetUseDefaultPID(kFALSE);\r
-\r
-  SetLowPt(kTRUE,2.);\r
-\r
-  PrintAll();\r
-\r
-\r
-  delete pidObj;\r
-  pidObj=NULL;\r
-\r
-  return;\r
-  \r
-}\r
-\r
-\r
-void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2011() {\r
-  \r
-  // Default 2010 PbPb cut object\r
-  SetStandardCutsPbPb2010();\r
-\r
-  //\r
-  // Enable all 2011 PbPb run triggers\r
-  //  \r
-  SetTriggerClass("");\r
-  ResetMaskAndEnableMBTrigger();\r
-  EnableCentralTrigger();\r
-  EnableSemiCentralTrigger();\r
-\r
-}\r
+/**************************************************************************
+ * 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$ */
+
+/////////////////////////////////////////////////////////////
+//
+// Class for cuts on AOD reconstructed D0->Kpi
+//
+// Author: A.Dainese, andrea.dainese@pd.infn.it
+/////////////////////////////////////////////////////////////
+
+#include <TDatabasePDG.h>
+#include <Riostream.h>
+
+#include "AliRDHFCutsD0toKpi.h"
+#include "AliAODRecoDecayHF2Prong.h"
+#include "AliAODTrack.h"
+#include "AliESDtrack.h"
+#include "AliAODPid.h"
+#include "AliTPCPIDResponse.h"
+#include "AliAODVertex.h"
+#include "AliKFVertex.h"
+#include "AliKFParticle.h"
+
+ClassImp(AliRDHFCutsD0toKpi)
+
+//--------------------------------------------------------------------------
+AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const char* name) : 
+AliRDHFCuts(name),
+fUseSpecialCuts(kFALSE),
+fLowPt(kTRUE),
+fDefaultPID(kFALSE),
+fUseKF(kFALSE),
+fPtLowPID(2.),
+fPtMaxSpecialCuts(9999.)
+{
+  //
+  // Default Constructor
+  //
+  Int_t nvars=11;
+  SetNVars(nvars);
+  TString varNames[11]={"inv. mass [GeV]",   
+                       "dca [cm]",
+                       "cosThetaStar", 
+                       "pTK [GeV/c]",
+                       "pTPi [GeV/c]",
+                       "d0K [cm]",
+                       "d0Pi [cm]",
+                       "d0d0 [cm^2]",
+                       "cosThetaPoint",
+                       "|cosThetaPointXY|", 
+                       "NormDecayLenghtXY"};
+  Bool_t isUpperCut[11]={kTRUE,
+                        kTRUE,
+                        kTRUE,
+                        kFALSE,
+                        kFALSE,
+                        kTRUE,
+                        kTRUE,
+                        kTRUE,
+                        kFALSE,
+                        kFALSE, 
+                        kFALSE};
+  SetVarNames(nvars,varNames,isUpperCut);
+  Bool_t forOpt[11]={kFALSE,
+                    kTRUE,
+                    kTRUE,
+                    kFALSE,
+                    kFALSE,
+                    kFALSE,
+                    kFALSE,
+                    kTRUE,
+                    kTRUE,
+                    kFALSE,
+                    kFALSE};
+  SetVarsForOpt(4,forOpt);
+  Float_t limits[2]={0,999999999.};
+  SetPtBins(2,limits);
+
+}
+//--------------------------------------------------------------------------
+AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const AliRDHFCutsD0toKpi &source) :
+  AliRDHFCuts(source),   
+  fUseSpecialCuts(source.fUseSpecialCuts),
+  fLowPt(source.fLowPt),
+  fDefaultPID(source.fDefaultPID),
+  fUseKF(source.fUseKF),
+  fPtLowPID(source.fPtLowPID),
+  fPtMaxSpecialCuts(source.fPtMaxSpecialCuts)
+{
+  //
+  // Copy constructor
+  //
+
+}
+//--------------------------------------------------------------------------
+AliRDHFCutsD0toKpi &AliRDHFCutsD0toKpi::operator=(const AliRDHFCutsD0toKpi &source)
+{
+  //
+  // assignment operator
+  //
+  if(&source == this) return *this;
+
+  AliRDHFCuts::operator=(source); 
+  fUseSpecialCuts=source.fUseSpecialCuts;
+  fLowPt=source.fLowPt;
+  fDefaultPID=source.fDefaultPID;
+  fUseKF=source.fUseKF;
+  fPtLowPID=source.fPtLowPID;
+  fPtMaxSpecialCuts=source.fPtMaxSpecialCuts;
+
+  return *this;
+}
+
+
+//---------------------------------------------------------------------------
+void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent *aod) {
+  // 
+  // Fills in vars the values of the variables 
+  //
+
+  if(nvars!=fnVarsForOpt) {
+    printf("AliRDHFCutsD0toKpi::GetCutsVarsForOpt: wrong number of variables\n");
+    return;
+  }
+
+  AliAODRecoDecayHF2Prong *dd = (AliAODRecoDecayHF2Prong*)d;
+  //recalculate vertex w/o daughters
+  Bool_t cleanvtx=kFALSE;
+  AliAODVertex *origownvtx=0x0;
+  if(fRemoveDaughtersFromPrimary) {
+    if(dd->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*dd->GetOwnPrimaryVtx());
+    cleanvtx=kTRUE;
+    if(!RecalcOwnPrimaryVtx(dd,aod)) {
+      CleanOwnPrimaryVtx(dd,aod,origownvtx);
+      cleanvtx=kFALSE;
+    }
+  }
+
+  Int_t iter=-1;
+  if(fVarsForOpt[0]){
+    iter++;
+    if(TMath::Abs(pdgdaughters[0])==211) {
+      vars[iter]=dd->InvMassD0();
+    } else {
+      vars[iter]=dd->InvMassD0bar();
+    }
+  }
+  if(fVarsForOpt[1]){
+    iter++;
+    vars[iter]=dd->GetDCA();
+  }
+  if(fVarsForOpt[2]){
+    iter++;
+    if(TMath::Abs(pdgdaughters[0])==211) {
+      vars[iter] = dd->CosThetaStarD0();
+    } else {
+      vars[iter] = dd->CosThetaStarD0bar();
+    }
+  }
+  if(fVarsForOpt[3]){
+    iter++;
+   if(TMath::Abs(pdgdaughters[0])==321) {
+     vars[iter]=dd->PtProng(0);
+   }
+   else{
+     vars[iter]=dd->PtProng(1);
+   }
+  }
+  if(fVarsForOpt[4]){
+    iter++;
+   if(TMath::Abs(pdgdaughters[0])==211) {
+     vars[iter]=dd->PtProng(0);
+   }
+   else{
+     vars[iter]=dd->PtProng(1);
+   }
+  }
+  if(fVarsForOpt[5]){
+    iter++;
+    if(TMath::Abs(pdgdaughters[0])==321) {
+     vars[iter]=dd->Getd0Prong(0);
+   }
+   else{
+     vars[iter]=dd->Getd0Prong(1);
+   }
+  }
+  if(fVarsForOpt[6]){
+    iter++;
+     if(TMath::Abs(pdgdaughters[0])==211) {
+     vars[iter]=dd->Getd0Prong(0);
+   }
+   else{
+     vars[iter]=dd->Getd0Prong(1);
+   }
+  }
+  if(fVarsForOpt[7]){
+    iter++;
+    vars[iter]= dd->Prodd0d0();
+  }
+  if(fVarsForOpt[8]){
+    iter++;
+    vars[iter]=dd->CosPointingAngle();
+  }
+  
+  if(fVarsForOpt[9]){
+               iter++;
+               vars[iter]=TMath::Abs(dd->CosPointingAngleXY());
+       }
+  
+   if(fVarsForOpt[10]){
+               iter++;
+          vars[iter]=dd->NormalizedDecayLengthXY();
+       }
+
+   if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
+
+  return;
+}
+//---------------------------------------------------------------------------
+Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* aod) {
+  //
+  // Apply selection
+  //
+
+  fIsSelectedCuts=0;
+  fIsSelectedPID=0;
+
+  if(!fCutsRD){
+    cout<<"Cut matrice not inizialized. Exit..."<<endl;
+    return 0;
+  }
+  //PrintAll();
+  AliAODRecoDecayHF2Prong* d=(AliAODRecoDecayHF2Prong*)obj;
+
+  if(!d){
+    cout<<"AliAODRecoDecayHF2Prong null"<<endl;
+    return 0;
+  }
+
+  if(fKeepSignalMC) if(IsSignalMC(d,aod,421)) return 3;
+
+  Double_t ptD=d->Pt();
+  if(ptD<fMinPtCand) return 0;
+  if(ptD>fMaxPtCand) return 0;
+
+  if(d->HasBadDaughters()) return 0;
+
+  // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
+  Int_t returnvaluePID=3;
+  Int_t returnvalueCuts=3;
+
+  // selection on candidate
+  if(selectionLevel==AliRDHFCuts::kAll || 
+     selectionLevel==AliRDHFCuts::kCandidate) {
+
+    if(!fUseKF) {
+
+      //recalculate vertex w/o daughters
+      AliAODVertex *origownvtx=0x0;
+      if(fRemoveDaughtersFromPrimary && !fUseMCVertex) {
+       if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
+       if(!RecalcOwnPrimaryVtx(d,aod)) { 
+         CleanOwnPrimaryVtx(d,aod,origownvtx);
+         return 0;
+       }
+      }
+
+      if(fUseMCVertex) {
+       if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
+       if(!SetMCPrimaryVtx(d,aod)) {
+         CleanOwnPrimaryVtx(d,aod,origownvtx);
+         return 0;
+       }
+      }
+      
+      Double_t pt=d->Pt();
+      
+      Int_t okD0=0,okD0bar=0;
+      
+      Int_t ptbin=PtBin(pt);
+      if (ptbin==-1) {
+       CleanOwnPrimaryVtx(d,aod,origownvtx);
+       return 0;
+      }
+
+      Double_t mD0,mD0bar,ctsD0,ctsD0bar;
+      okD0=1; okD0bar=1;
+      
+      Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
+
+      d->InvMassD0(mD0,mD0bar);
+      if(TMath::Abs(mD0-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
+      if(TMath::Abs(mD0bar-mD0PDG) > fCutsRD[GetGlobalIndex(0,ptbin)])  okD0bar = 0;
+      if(!okD0 && !okD0bar)  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+
+      if(d->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+    
+      
+      if(d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(3,ptbin)]*fCutsRD[GetGlobalIndex(3,ptbin)] || d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(4,ptbin)]*fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
+      if(d->Pt2Prong(0) < fCutsRD[GetGlobalIndex(3,ptbin)]*fCutsRD[GetGlobalIndex(3,ptbin)] || d->Pt2Prong(1) < fCutsRD[GetGlobalIndex(4,ptbin)]*fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;
+      if(!okD0 && !okD0bar) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+      
+      
+      if(TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)] || 
+        TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
+      if(TMath::Abs(d->Getd0Prong(0)) > fCutsRD[GetGlobalIndex(6,ptbin)] ||
+        TMath::Abs(d->Getd0Prong(1)) > fCutsRD[GetGlobalIndex(5,ptbin)]) okD0bar = 0;
+      if(!okD0 && !okD0bar)  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+      
+      if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+      
+    
+      d->CosThetaStarD0(ctsD0,ctsD0bar);
+      if(TMath::Abs(ctsD0) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0; 
+      if(TMath::Abs(ctsD0bar) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0bar = 0;
+      if(!okD0 && !okD0bar)   {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+    
+      if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+
+      
+      if(TMath::Abs(d->CosPointingAngleXY()) < fCutsRD[GetGlobalIndex(9,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+       
+      Double_t normalDecayLengXY=d->NormalizedDecayLengthXY();
+      if (normalDecayLengXY < fCutsRD[GetGlobalIndex(10, ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+      
+      if (returnvalueCuts!=0) {
+       if (okD0) returnvalueCuts=1; //cuts passed as D0
+       if (okD0bar) returnvalueCuts=2; //cuts passed as D0bar
+       if (okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar
+      }
+
+      // call special cuts
+      Int_t special=1;
+      if(fUseSpecialCuts && (pt<fPtMaxSpecialCuts)) special=IsSelectedSpecialCuts(d);
+      if(!special) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
+
+      // unset recalculated primary vertex when not needed any more
+      CleanOwnPrimaryVtx(d,aod,origownvtx);
+
+    } else {
+      // go to selection with Kalman vertexing, if requested
+      returnvalueCuts = IsSelectedKF(d,aod);
+    }
+
+    fIsSelectedCuts=returnvalueCuts;
+    if(!returnvalueCuts) return 0;
+  }
+
+
+  // selection on PID 
+  if(selectionLevel==AliRDHFCuts::kAll || 
+     selectionLevel==AliRDHFCuts::kCandidate ||
+     selectionLevel==AliRDHFCuts::kPID) {
+    returnvaluePID = IsSelectedPID(d);
+    fIsSelectedPID=returnvaluePID;
+    if(!returnvaluePID) return 0;
+  }
+
+  Int_t returnvalueComb=CombineSelectionLevels(3,returnvalueCuts,returnvaluePID);
+
+  if(!returnvalueComb) return 0;
+
+  // selection on daughter tracks 
+  if(selectionLevel==AliRDHFCuts::kAll || 
+     selectionLevel==AliRDHFCuts::kTracks) {
+    if(!AreDaughtersSelected(d)) return 0;
+  }
+  //  cout<<"Pid = "<<returnvaluePID<<endl;
+  return returnvalueComb;
+}
+
+//------------------------------------------------------------------------------------------
+Int_t AliRDHFCutsD0toKpi::IsSelectedKF(AliAODRecoDecayHF2Prong *d,
+                                      AliAODEvent* aod) const {
+  //
+  // Apply selection using KF-vertexing
+  //
+       
+  AliAODTrack *track0 = (AliAODTrack*)d->GetDaughter(0);
+  AliAODTrack *track1 = (AliAODTrack*)d->GetDaughter(1); 
+   
+  if(!track0 || !track1) {
+    cout<<"one or two D0 daughters missing!"<<endl;
+    return 0;
+  }
+
+  // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
+  Int_t returnvalueCuts=3;
+        
+  // candidate selection with AliKF
+  AliKFParticle::SetField(aod->GetMagneticField()); // set the magnetic field
+  
+  Int_t okD0=0,okD0bar=0;
+  okD0=1; okD0bar=1;
+  
+  // convert tracks into AliKFParticles
+  
+  AliKFParticle negPiKF(*track1,-211); // neg pion kandidate
+  AliKFParticle negKKF(*track1,-321); // neg kaon kandidate
+  AliKFParticle posPiKF(*track0,211); // pos pion kandidate
+  AliKFParticle posKKF(*track0,321); // pos kaon kandidate
+  
+  // build D0 candidates
+  
+  AliKFParticle d0c(negKKF,posPiKF); // D0 candidate
+  AliKFParticle ad0c(posKKF,negPiKF); // D0bar candidate
+  
+  // create kf primary vertices
+  
+  AliAODVertex *vtx1 = aod->GetPrimaryVertex();
+  AliKFVertex primVtx(*vtx1); 
+  AliKFVertex aprimVtx(*vtx1);
+  
+  if(primVtx.GetNContributors()<=0) okD0 = 0;
+  if(aprimVtx.GetNContributors()<=0) okD0bar = 0;
+  if(!okD0 && !okD0bar) returnvalueCuts=0;
+       
+  // calculate mass
+       
+  Double_t d0mass = d0c.GetMass();
+  Double_t ad0mass = ad0c.GetMass();
+       
+  // calculate P of D0 and D0bar
+  Double_t d0P = d0c.GetP();
+  Double_t d0Px = d0c.GetPx();
+  Double_t d0Py = d0c.GetPy();
+  Double_t d0Pz = d0c.GetPz();
+  Double_t ad0P = ad0c.GetP(); 
+  Double_t ad0Px = ad0c.GetPx();
+  Double_t ad0Py = ad0c.GetPy();
+  Double_t ad0Pz = ad0c.GetPz();
+  
+  //calculate Pt of D0 and D0bar
+       
+  Double_t pt=d0c.GetPt(); 
+  Double_t apt=ad0c.GetPt();
+       
+  // remove D0 daughters from primary vertices (if used in vertex fit) and add D0-candidates
+  
+  if(track0->GetUsedForPrimVtxFit()) {
+    primVtx -= posPiKF; 
+    aprimVtx -= posKKF;
+  }
+  
+  if(track1->GetUsedForPrimVtxFit()) { 
+    primVtx -= negKKF; 
+    aprimVtx -= negPiKF;
+  }
+  
+  primVtx += d0c;
+  aprimVtx += ad0c;
+  
+  if(primVtx.GetNContributors()<=0) okD0 = 0;
+  if(aprimVtx.GetNContributors()<=0) okD0bar = 0;
+  if(!okD0 && !okD0bar) returnvalueCuts=0;
+  
+  //calculate cut variables
+  
+  // calculate impact params of daughters w.r.t recalculated vertices
+  
+  Double_t impactPi = posPiKF.GetDistanceFromVertexXY(primVtx);
+  Double_t aimpactPi = negPiKF.GetDistanceFromVertexXY(aprimVtx);
+  Double_t impactKa = negKKF.GetDistanceFromVertexXY(primVtx);
+  Double_t aimpactKa = posKKF.GetDistanceFromVertexXY(aprimVtx);
+       
+  // calculate Product of Impact Params
+       
+  Double_t prodParam = impactPi*impactKa;
+  Double_t aprodParam = aimpactPi*aimpactKa;
+       
+  // calculate cosine of pointing angles
+       
+  TVector3 mom(d0c.GetPx(),d0c.GetPy(),d0c.GetPz());
+  TVector3 fline(d0c.GetX()-primVtx.GetX(),
+                d0c.GetY()-primVtx.GetY(),
+                d0c.GetZ()-primVtx.GetZ());
+  Double_t pta = mom.Angle(fline);
+  Double_t cosP = TMath::Cos(pta); // cosine of pta for D0 candidate
+  
+  TVector3 amom(ad0c.GetPx(),ad0c.GetPy(),ad0c.GetPz());
+  TVector3 afline(ad0c.GetX()-aprimVtx.GetX(),
+                 ad0c.GetY()-aprimVtx.GetY(),
+                 ad0c.GetZ()-aprimVtx.GetZ());
+  Double_t apta = amom.Angle(afline);
+  Double_t acosP = TMath::Cos(apta); // cosine of pta for D0bar candidate
+  
+  // calculate P of Pions at Decay Position of D0 and D0bar candidates
+  negKKF.TransportToParticle(d0c);
+  posPiKF.TransportToParticle(d0c);
+  posKKF.TransportToParticle(ad0c);
+  negPiKF.TransportToParticle(ad0c);
+  
+  Double_t pxPi =  posPiKF.GetPx();
+  Double_t pyPi =  posPiKF.GetPy();
+  Double_t pzPi =  posPiKF.GetPz();
+  Double_t ptPi =  posPiKF.GetPt();
+  
+  Double_t apxPi =  negPiKF.GetPx();
+  Double_t apyPi =  negPiKF.GetPy();
+  Double_t apzPi =  negPiKF.GetPz();
+  Double_t aptPi =  negPiKF.GetPt();
+  
+  // calculate Pt of Kaons at Decay Position of D0 and D0bar candidates
+  
+  Double_t ptK =  negKKF.GetPt();
+  Double_t aptK =  posKKF.GetPt();
+       
+  //calculate cos(thetastar)
+  Double_t massvtx = TDatabasePDG::Instance()->GetParticle(421)->Mass();
+  Double_t massp[2];
+  massp[0] = TDatabasePDG::Instance()->GetParticle(321)->Mass();
+  massp[1] = TDatabasePDG::Instance()->GetParticle(211)->Mass();
+  Double_t pStar = TMath::Sqrt(TMath::Power(massvtx*massvtx-massp[0]*massp[0]-massp[1]*massp[1],2.)
+                              -4.*massp[0]*massp[0]*massp[1]*massp[1])/(2.*massvtx);
+  
+  // cos(thetastar) for D0 and Pion
+  
+  Double_t d0E = TMath::Sqrt(massvtx*massvtx + d0P*d0P);
+  Double_t beta = d0P/d0E;
+  Double_t gamma = d0E/massvtx;
+  TVector3 momPi(pxPi,pyPi,pzPi);
+  TVector3 momTot(d0Px,d0Py,d0Pz);
+  Double_t q1 = momPi.Dot(momTot)/momTot.Mag();
+  Double_t cts = (q1/gamma-beta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar;     
+       
+  // cos(thetastar) for D0bar and Pion 
+  
+  Double_t ad0E = TMath::Sqrt(massvtx*massvtx + ad0P*ad0P);
+  Double_t abeta = ad0P/ad0E;
+  Double_t agamma = ad0E/massvtx;
+  TVector3 amomPi(apxPi,apyPi,apzPi);
+  TVector3 amomTot(ad0Px,ad0Py,ad0Pz);
+  Double_t aq1 = amomPi.Dot(amomTot)/amomTot.Mag();
+  Double_t acts = (aq1/agamma-abeta*TMath::Sqrt(pStar*pStar+massp[1]*massp[1]))/pStar; 
+  
+  // calculate reduced Chi2 for the full D0 fit
+  d0c.SetProductionVertex(primVtx);
+  ad0c.SetProductionVertex(aprimVtx);
+  negKKF.SetProductionVertex(d0c);
+  posPiKF.SetProductionVertex(d0c);
+  posKKF.SetProductionVertex(ad0c);
+  negPiKF.SetProductionVertex(ad0c);
+  d0c.TransportToProductionVertex();
+  ad0c.TransportToProductionVertex();
+       
+  // calculate the decay length
+  Double_t decayLengthD0 = d0c.GetDecayLength();
+  Double_t adecayLengthD0 = ad0c.GetDecayLength();
+  
+  Double_t chi2D0 = 50.;
+  if(d0c.GetNDF() > 0 && d0c.GetChi2() >= 0) {
+    chi2D0 = d0c.GetChi2()/d0c.GetNDF();
+  }
+  
+  Double_t achi2D0 = 50.;
+  if(ad0c.GetNDF() > 0 && ad0c.GetChi2() >= 0) {
+    achi2D0 = ad0c.GetChi2()/ad0c.GetNDF();
+  }
+       
+  // Get the Pt-bins
+  Int_t ptbin=PtBin(pt);
+  Int_t aptbin=PtBin(apt);
+
+  if(ptbin < 0) okD0 = 0;
+  if(aptbin < 0) okD0bar = 0;
+  if(!okD0 && !okD0bar) returnvalueCuts=0;
+  
+  if(ptK < fCutsRD[GetGlobalIndex(3,ptbin)] || ptPi < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
+  if(aptK < fCutsRD[GetGlobalIndex(3,aptbin)] || aptPi < fCutsRD[GetGlobalIndex(4,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar) returnvalueCuts=0;
+  
+  
+  if(TMath::Abs(impactKa) > fCutsRD[GetGlobalIndex(5,ptbin)] || 
+     TMath::Abs(impactPi) > fCutsRD[GetGlobalIndex(6,ptbin)]) okD0 = 0;
+  
+  if(TMath::Abs(aimpactKa) > fCutsRD[GetGlobalIndex(5,aptbin)] ||
+     TMath::Abs(aimpactPi) > fCutsRD[GetGlobalIndex(6,aptbin)]) okD0bar = 0;
+  
+  if(!okD0 && !okD0bar)  returnvalueCuts=0;
+  
+  // for the moment via the standard method due to bug in AliKF
+  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)]) okD0 = 0;
+  if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar) returnvalueCuts=0;
+    
+    
+  if(TMath::Abs(d0mass-massvtx) > fCutsRD[GetGlobalIndex(0,ptbin)]) okD0 = 0;
+  if(TMath::Abs(ad0mass-massvtx) > fCutsRD[GetGlobalIndex(0,aptbin)])  okD0bar = 0;
+  if(!okD0 && !okD0bar)  returnvalueCuts=0;
+  
+  
+  if(TMath::Abs(cts) > fCutsRD[GetGlobalIndex(2,ptbin)]) okD0 = 0; 
+  if(TMath::Abs(acts) > fCutsRD[GetGlobalIndex(2,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar)   returnvalueCuts=0;
+  
+  if(prodParam  > fCutsRD[GetGlobalIndex(7,ptbin)]) okD0 = 0;
+  if(aprodParam > fCutsRD[GetGlobalIndex(7,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar)  returnvalueCuts=0;
+    
+  if(cosP  < fCutsRD[GetGlobalIndex(8,ptbin)]) okD0 = 0; 
+  if(acosP < fCutsRD[GetGlobalIndex(8,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar)  returnvalueCuts=0;
+       
+  if(chi2D0  > fCutsRD[GetGlobalIndex(10,ptbin)]) okD0 = 0; 
+  if(achi2D0 > fCutsRD[GetGlobalIndex(10,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar)  returnvalueCuts=0;
+       
+  if(decayLengthD0  < fCutsRD[GetGlobalIndex(9,ptbin)]) okD0 = 0; 
+  if(adecayLengthD0 < fCutsRD[GetGlobalIndex(9,aptbin)]) okD0bar = 0;
+  if(!okD0 && !okD0bar)  returnvalueCuts=0;
+    
+  if(returnvalueCuts!=0) {
+    if(okD0) returnvalueCuts=1; //cuts passed as D0
+    if(okD0bar) returnvalueCuts=2; //cuts passed as D0bar
+    if(okD0 && okD0bar) returnvalueCuts=3; //cuts passed as D0 and D0bar
+  }
+
+  return returnvalueCuts;  
+}
+
+//---------------------------------------------------------------------------
+
+Bool_t AliRDHFCutsD0toKpi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
+{
+  //
+  // Checking if D0 is in fiducial acceptance region 
+  //
+
+  if(pt > 5.) {
+    // applying cut for pt > 5 GeV
+    AliDebug(2,Form("pt of D0 = %f (> 5), cutting at |y| < 0.8\n",pt)); 
+    if (TMath::Abs(y) > 0.8){
+      return kFALSE;
+    }
+  } else {
+    // appliying smooth cut for pt < 5 GeV
+    Double_t maxFiducialY = -0.2/15*pt*pt+1.9/15*pt+0.5; 
+    Double_t minFiducialY = 0.2/15*pt*pt-1.9/15*pt-0.5;                
+    AliDebug(2,Form("pt of D0 = %f (< 5), cutting  according to the fiducial zone [%f, %f]\n",pt,minFiducialY,maxFiducialY)); 
+    if (y < minFiducialY || y > maxFiducialY){
+      return kFALSE;
+    }
+  }
+
+  return kTRUE;
+}
+//---------------------------------------------------------------------------
+Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d) 
+{
+  // ############################################################
+  //
+  // Apply PID selection
+  //
+  //
+  // ############################################################
+
+  if(!fUsePID) return 3;
+  if(fDefaultPID) return IsSelectedPIDdefault(d);
+  fWhyRejection=0;
+  Int_t isD0D0barPID[2]={1,2};
+  Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value; 
+  //                                                                                                 same for prong 2
+  //                                               values convention -1 = discarded 
+  //                                                                  0 = not identified (but compatible) || No PID (->hasPID flag)
+  //                                                                  1 = identified
+  // PID search:   pion (TPC) or not K (TOF), Kaon hypothesis for both 
+  // Initial hypothesis: unknwon (but compatible) 
+  combinedPID[0][0]=0;  // prima figlia, Kaon
+  combinedPID[0][1]=0;  // prima figlia, pione
+  combinedPID[1][0]=0;  // seconda figlia, Kaon
+  combinedPID[1][1]=0;  // seconda figlia, pion
+  
+  Bool_t checkPIDInfo[2]={kTRUE,kTRUE};
+  Double_t sigma_tmp[3]={fPidHF->GetSigma(0),fPidHF->GetSigma(1),fPidHF->GetSigma(2)};
+  for(Int_t daught=0;daught<2;daught++){
+    //Loop con prongs
+    AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);
+    if(fPidHF->IsTOFPiKexcluded(aodtrack,5.)) return 0; 
+    
+    if(!(fPidHF->CheckStatus(aodtrack,"TPC")) && !(fPidHF->CheckStatus(aodtrack,"TOF"))) {
+      checkPIDInfo[daught]=kFALSE; 
+      continue;
+    }
+
+    // identify kaon
+    combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
+
+    // identify pion
+
+    if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {
+     combinedPID[daught][1]=0;
+    }else{
+      fPidHF->SetTOF(kFALSE);
+      combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
+      fPidHF->SetTOF(kTRUE);
+      fPidHF->SetCompat(kTRUE);
+     }
+
+
+    if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){ // if not a K- and not a pi- both D0 and D0bar excluded
+      isD0D0barPID[0]=0;
+      isD0D0barPID[1]=0;
+    }
+    else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){
+      if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded
+      else isD0D0barPID[0]=0;// if K+ D0 excluded
+    }
+    /*    else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){
+         isD0D0barPID[0]=0;
+         isD0D0barPID[1]=0;
+         }
+    */
+    else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){ 
+      if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded
+      else isD0D0barPID[0]=0;//  not a D0 if K+ or if pi+ excluded
+    }
+    else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){
+      if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded
+      else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded
+    }
+
+    if(fLowPt && d->Pt()<fPtLowPID){
+     Double_t sigmaTPC[3]={3.,2.,0.};
+     fPidHF->SetSigmaForTPC(sigmaTPC);
+    // identify kaon
+    combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
+
+    Double_t ptProng=aodtrack->P();
+
+    if(ptProng<0.6){
+     fPidHF->SetCompat(kFALSE);
+     combinedPID[daught][0]=fPidHF->MakeRawPid(aodtrack,3);
+     fPidHF->SetCompat(kTRUE);
+    }
+
+    if(!(fPidHF->CheckStatus(aodtrack,"TPC"))) {
+     combinedPID[daught][1]=0;
+    }else{
+      fPidHF->SetTOF(kFALSE);
+      Double_t sigmaTPCpi[3]={3.,3.,0.};
+      fPidHF->SetSigmaForTPC(sigmaTPCpi);
+      combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
+      fPidHF->SetTOF(kTRUE);
+       if(ptProng<0.8){
+        Bool_t isTPCpion=fPidHF->IsPionRaw(aodtrack,"TPC");
+        if(isTPCpion){
+         combinedPID[daught][1]=1;
+        }else{
+         combinedPID[daught][1]=-1;
+        }
+      }
+    }
+
+   }
+   fPidHF->SetSigmaForTPC(sigma_tmp);
+  }// END OF LOOP ON DAUGHTERS
+
+   if(!checkPIDInfo[0] && !checkPIDInfo[1]) {
+    if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
+    return 0;
+   }
+
+
+  // FURTHER PID REQUEST (both daughter info is needed)
+  if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){
+    fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found
+    if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
+    return 0;
+  }
+
+  if(fLowPt && d->Pt()<fPtLowPID){    
+    if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
+      fWhyRejection=32;// reject cases where the Kaon is not identified
+      fPidHF->SetSigmaForTPC(sigma_tmp);
+      return 0;
+    }
+  }
+    if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
+
+  //  cout<<"Why? "<<fWhyRejection<<endl;  
+  return isD0D0barPID[0]+isD0D0barPID[1];
+}
+//---------------------------------------------------------------------------
+Int_t AliRDHFCutsD0toKpi::IsSelectedPIDdefault(AliAODRecoDecayHF* d) 
+{
+  // ############################################################
+  //
+  // Apply PID selection
+  //
+  //
+  // temporary selection: PID AS USED FOR D0 by Andrea Rossi (up to 28/06/2010)
+  //
+  // d must be a AliAODRecoDecayHF2Prong object
+  // returns 0 if both D0 and D0bar are rejectecd
+  //         1 if D0 is accepted while D0bar is rejected
+  //         2 if D0bar is accepted while D0 is rejected
+  //         3 if both are accepted
+  // fWhyRejection variable (not returned for the moment, print it if needed)
+  //               keeps some information on why a candidate has been 
+  //               rejected according to the following (unfriendly?) scheme 
+  //             if more rejection cases are considered interesting, just add numbers
+  //
+  //      TO BE CONSIDERED WITH A GRAIN OF SALT (the order in which cut are applied is relevant) 
+  //              from 20 to 30: "detector" selection (PID acceptance)                                             
+  //                                                  26: TPC refit
+  //                                                  27: ITS refit
+  //                                                  28: no (TOF||TPC) pid information (no kTOFpid,kTOFout,kTIME,kTPCpid,...)
+  //
+  //              from 30 to 40: PID selection
+  //                                                  31: no Kaon compatible tracks found between daughters
+  //                                                  32: no Kaon identified tracks found (strong sel. at low momenta)
+  //                                                  33: both mass hypotheses are rejected 
+  //                  
+  // ############################################################
+
+  if(!fUsePID) return 3;
+  fWhyRejection=0;
+  Int_t isD0D0barPID[2]={1,2};
+  Double_t nsigmaTPCpi=-1., nsigmaTPCK=-1.; //used for TPC pid
+  Double_t tofSig,times[5];// used fot TOF pid
+  Int_t hasPID[2]={2,2};// flag to count how many detectors give PID info for the daughters
+  Int_t isKaonPionTOF[2][2],isKaonPionTPC[2][2];
+  Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value; 
+  //                                                                                                 same for prong 2
+  //                                               values convention -1 = discarded 
+  //                                                                  0 = not identified (but compatible) || No PID (->hasPID flag)
+  //                                                                  1 = identified
+  // PID search:   pion (TPC) or not K (TOF), Kaon hypothesis for both 
+  // Initial hypothesis: unknwon (but compatible) 
+  isKaonPionTOF[0][0]=0;
+  isKaonPionTOF[0][1]=0;
+  isKaonPionTOF[1][0]=0;
+  isKaonPionTOF[1][1]=0;
+  
+  isKaonPionTPC[0][0]=0;
+  isKaonPionTPC[0][1]=0;
+  isKaonPionTPC[1][0]=0;
+  isKaonPionTPC[1][1]=0;
+  
+  combinedPID[0][0]=0;
+  combinedPID[0][1]=0;
+  combinedPID[1][0]=0;
+  combinedPID[1][1]=0;
+  
+  
+  for(Int_t daught=0;daught<2;daught++){
+    //Loop con prongs
+    
+    // ########### Step 0- CHECKING minimal PID "ACCEPTANCE" ####################
+
+    AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught); 
+   
+    if(!(aodtrack->GetStatus()&AliESDtrack::kTPCrefit)){
+      fWhyRejection=26;
+      return 0;
+    } 
+    if(!(aodtrack->GetStatus()&AliESDtrack::kITSrefit)){
+      fWhyRejection=27;
+      return 0;
+    } 
+    
+    AliAODPid *pid=aodtrack->GetDetPid();
+    if(!pid) {
+      //delete esdtrack;
+      hasPID[daught]--;
+      continue;
+    }
+  
+    // ########### Step 1- Check of TPC and TOF response ####################
+
+    Double_t ptrack=aodtrack->P();
+    //#################### TPC PID #######################
+     if (!(aodtrack->GetStatus()&AliESDtrack::kTPCpid )){
+       // NO TPC PID INFO FOR THIS TRACK 
+       hasPID[daught]--;
+     }
+     else {
+       static AliTPCPIDResponse theTPCpid;
+       AliAODPid *pidObj = aodtrack->GetDetPid();
+       Double_t ptProng=pidObj->GetTPCmomentum();
+       nsigmaTPCpi = theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kPion);
+       nsigmaTPCK =  theTPCpid.GetNumberOfSigmas(ptProng,(Float_t)pid->GetTPCsignal(),(Int_t)aodtrack->GetTPCClusterMap().CountBits(),AliPID::kKaon);
+       //if(ptrack<0.6){
+       if(ptProng<0.6){
+        if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;
+        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
+        if(TMath::Abs(nsigmaTPCpi)<2.)isKaonPionTPC[daught][1]=1;
+        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
+       }
+       //else if(ptrack<.8){
+       else if(ptProng<.8){
+        if(TMath::Abs(nsigmaTPCK)<1.)isKaonPionTPC[daught][0]=1;
+        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
+        if(TMath::Abs(nsigmaTPCpi)<1.)isKaonPionTPC[daught][1]=1;
+        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
+       }     
+       else {
+        //     if(nsigmaTPCK>-2.&&nsigmaTPCK<1.)isKaonPionTPC[daught][0]=1;
+        if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
+        //if(nsigmaTPCpi>-1.&&nsigmaTPCpi<2.)isKaonPionTPC[daught][1]=1;
+        if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
+       }
+     }
+    
+    
+    // ##### TOF PID: do not ask nothing for pion/protons ############
+     if(!((aodtrack->GetStatus()&AliESDtrack::kTOFpid)&&(aodtrack->GetStatus()&AliESDtrack::kTOFout)&&(aodtrack->GetStatus()&AliESDtrack::kTIME))){
+       // NO TOF PID INFO FOR THIS TRACK      
+       hasPID[daught]--;
+     }
+     else{
+       tofSig=pid->GetTOFsignal(); 
+       pid->GetIntegratedTimes(times);
+       if((tofSig-times[3])>5.*160.)return 0;// PROTON REJECTION
+       if(TMath::Abs(tofSig-times[3])>3.*160.){
+        isKaonPionTOF[daught][0]=-1;
+       }
+       else {   
+        if(ptrack<1.5){
+          isKaonPionTOF[daught][0]=1;
+        }
+       }
+     }
+     
+     //######### Step 2: COMBINE TOF and TPC PID ###############
+     // we apply the following convention: if TPC and TOF disagree (discarded Vs identified) -> unknown
+     combinedPID[daught][0]=isKaonPionTOF[daught][0]+isKaonPionTPC[daught][0];
+     combinedPID[daught][1]=isKaonPionTOF[daught][1]+isKaonPionTPC[daught][1];
+     
+     
+     //######### Step 3:   USE PID INFO     
+     
+     if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){// if not a K- and not a pi- both D0 and D0bar excluded
+       isD0D0barPID[0]=0;
+       isD0D0barPID[1]=0;
+     }
+     else if(combinedPID[daught][0]==2&&combinedPID[daught][1]>=1){// if in conflict (both pi- and K-), if k for both TPC and TOF -> is K
+       if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;//if K- D0bar excluded
+       else isD0D0barPID[0]=0;// if K+ D0 excluded
+     }
+     else if(combinedPID[daught][0]==1&&combinedPID[daught][1]>=1){// if in conflict (both pi- and K-) and k- only for TPC or TOF -> reject
+       isD0D0barPID[0]=0;
+       isD0D0barPID[1]=0;
+     }
+     else if(combinedPID[daught][0]>=1||combinedPID[daught][1]<=-1){
+       if(aodtrack->Charge()==-1)isD0D0barPID[1]=0;// not a D0bar if K- or if pi- excluded
+       else isD0D0barPID[0]=0;//  not a D0 if K+ or if pi+ excluded
+     }
+     else if(combinedPID[daught][0]<=-1||combinedPID[daught][1]>=1){
+       if(aodtrack->Charge()==-1)isD0D0barPID[0]=0;// not a D0 if pi- or if K- excluded
+      else isD0D0barPID[1]=0;// not a D0bar if pi+ or if K+ excluded
+     }
+     
+     // ##########  ALSO DIFFERENT TPC PID REQUEST FOR LOW pt D0: request of K identification      ###############################
+     // ########## more tolerant criteria for single particle ID-> more selective criteria for D0   ##############################
+     // ###############                     NOT OPTIMIZED YET                                  ###################################
+     if(d->Pt()<2.){
+       isKaonPionTPC[daught][0]=0;
+       isKaonPionTPC[daught][1]=0;
+       AliAODPid *pidObj = aodtrack->GetDetPid();
+       Double_t ptProng=pidObj->GetTPCmomentum();
+       //if(ptrack<0.6){
+       if(ptProng<0.6){
+        if(TMath::Abs(nsigmaTPCK)<3.)isKaonPionTPC[daught][0]=1;
+        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
+        if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;
+        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
+     }
+       //else if(ptrack<.8){
+       else if(ptProng<.8){
+        if(TMath::Abs(nsigmaTPCK)<2.)isKaonPionTPC[daught][0]=1;
+        else if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
+        if(TMath::Abs(nsigmaTPCpi)<3.)isKaonPionTPC[daught][1]=1;
+        else if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
+       }     
+       else {
+        if(TMath::Abs(nsigmaTPCK)>3.)isKaonPionTPC[daught][0]=-1;
+        if(TMath::Abs(nsigmaTPCpi)>3.)isKaonPionTPC[daught][1]=-1;
+       }
+     }
+     
+  }// END OF LOOP ON DAUGHTERS
+  
+  // FURTHER PID REQUEST (both daughter info is needed)
+  if(combinedPID[0][0]<=-1&&combinedPID[1][0]<=-1){
+    fWhyRejection=31;// reject cases in which no kaon-compatible tracks are found
+    return 0;
+  }
+  else if(hasPID[0]==0&&hasPID[1]==0){
+    fWhyRejection=28;// reject cases in which no PID info is available  
+    return 0;
+  }
+  if(d->Pt()<2.){
+    // request of K identification at low D0 pt
+    combinedPID[0][0]=0;
+    combinedPID[0][1]=0;
+    combinedPID[1][0]=0;
+    combinedPID[1][1]=0;
+    
+    combinedPID[0][0]=isKaonPionTOF[0][0]+isKaonPionTPC[0][0];
+    combinedPID[0][1]=isKaonPionTOF[0][1]+isKaonPionTPC[0][1];
+    combinedPID[1][0]=isKaonPionTOF[1][0]+isKaonPionTPC[1][0];
+    combinedPID[1][1]=isKaonPionTOF[1][1]+isKaonPionTPC[1][1];
+    
+    if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
+      fWhyRejection=32;// reject cases where the Kaon is not identified
+      return 0;
+    }
+  }
+
+  //  cout<<"Why? "<<fWhyRejection<<endl;  
+  return isD0D0barPID[0]+isD0D0barPID[1];
+}
+
+
+
+//---------------------------------------------------------------------------
+Int_t AliRDHFCutsD0toKpi::CombineSelectionLevels(Int_t selectionvalTrack,
+                                                Int_t selectionvalCand,
+                                                Int_t selectionvalPID) const
+{
+  //
+  // This method combines the tracks, PID and cuts selection results
+  //
+  if(selectionvalTrack==0) return 0;
+
+  Int_t returnvalue;
+
+  switch(selectionvalPID) {
+  case 0:
+    returnvalue=0;
+    break;
+  case 1:
+    returnvalue=((selectionvalCand==1 || selectionvalCand==3) ? 1 : 0);
+    break;
+  case 2:
+    returnvalue=((selectionvalCand==2 || selectionvalCand==3) ? 2 : 0);
+    break;
+  case 3:
+    returnvalue=selectionvalCand;
+    break;
+  default:
+    returnvalue=0;
+    break;
+  }
+
+  return returnvalue;
+}
+//----------------------------------------------------------------------------
+Int_t AliRDHFCutsD0toKpi::IsSelectedSpecialCuts(AliAODRecoDecayHF *d) const
+{
+  //
+  // Note: this method is temporary
+  // Additional cuts on decay lenght and lower cut for d0 norm are applied using vertex without candidate's daughters
+  //
+
+  //apply cuts
+
+  Float_t normDecLengthCut=1.,decLengthCut=TMath::Min(d->P()*0.0066+0.01,0.06/*cm*/), normd0Cut=0.5;
+  // "decay length" expo law with tau' = beta*gamma*ctau= p/m*ctau =p*0.0123/1.864~p*0.0066
+  // decay lenght > ctau' implies to retain (1-1/e) (for signal without considering detector resolution), 
+
+  Int_t returnvalue=3; //cut passed
+  for(Int_t i=0;i<2/*prongs*/;i++){
+    if(TMath::Abs(d->Normalizedd0Prong(i))<normd0Cut) return 0; //normd0Cut not passed
+  }
+  if(d->DecayLength2()<decLengthCut*decLengthCut)  return 0; //decLengthCut not passed
+  if(d->NormalizedDecayLength2()<normDecLengthCut*normDecLengthCut)  return 0; //decLengthCut not passed
+       
+  return returnvalue;
+}
+
+//----------------------------------------------
+void AliRDHFCutsD0toKpi::SetUseKF(Bool_t useKF)
+{
+  //switch on candidate selection via AliKFparticle
+  if(!useKF) return;
+  if(useKF){
+    fUseKF=useKF;
+    Int_t nvarsKF=11;
+    SetNVars(nvarsKF);
+    TString varNamesKF[11]={"inv. mass [GeV]",   
+                           "dca [cm]",
+                           "cosThetaStar", 
+                           "pTK [GeV/c]",
+                           "pTPi [GeV/c]",
+                           "d0K [cm]",
+                           "d0Pi [cm]",
+                           "d0d0 [cm^2]",
+                           "cosThetaPoint"
+                           "DecayLength[cm]",
+                           "RedChi2"};
+    Bool_t isUpperCutKF[11]={kTRUE,
+                            kTRUE,
+                            kTRUE,
+                            kFALSE,
+                            kFALSE,
+                            kTRUE,
+                            kTRUE,
+                            kTRUE,
+                            kFALSE,
+                            kFALSE,
+                            kTRUE};
+    SetVarNames(nvarsKF,varNamesKF,isUpperCutKF);
+    SetGlobalIndex();
+    Bool_t forOpt[11]={kFALSE,
+                   kTRUE,
+                   kTRUE,
+                   kFALSE,
+                   kFALSE,
+                   kFALSE,
+                   kFALSE,
+                   kTRUE,
+                   kTRUE,
+                   kFALSE,
+                   kFALSE};
+    SetVarsForOpt(4,forOpt);
+  }
+  return;
+}
+
+
+void AliRDHFCutsD0toKpi::SetStandardCutsPP2010() {
+  //
+  //STANDARD CUTS USED FOR 2010 pp analysis 
+  //dca cut will be enlarged soon to 400 micron
+  //
+  
+  SetName("D0toKpiCutsStandard");
+  SetTitle("Standard Cuts for D0 analysis");
+  
+  // PILE UP REJECTION
+  SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
+
+  // EVENT CUTS
+  SetMinVtxContr(1);
+ // MAX Z-VERTEX CUT
+  SetMaxVtxZ(10.);
+  
+  // TRACKS ON SINGLE TRACKS
+  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  //  esdTrackCuts->SetMinNClustersITS(4);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
+  esdTrackCuts->SetMinDCAToVertexXY(0.);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+  esdTrackCuts->SetPtRange(0.3,1.e10);
+  
+  AddTrackCuts(esdTrackCuts);
+
+  
+  const Int_t nptbins =14;
+  const Double_t ptmax = 9999.;
+  const Int_t nvars=11;
+  Float_t ptbins[nptbins+1];
+  ptbins[0]=0.;
+  ptbins[1]=0.5;       
+  ptbins[2]=1.;
+  ptbins[3]=2.;
+  ptbins[4]=3.;
+  ptbins[5]=4.;
+  ptbins[6]=5.;
+  ptbins[7]=6.;
+  ptbins[8]=7.;
+  ptbins[9]=8.;
+  ptbins[10]=12.;
+  ptbins[11]=16.;
+  ptbins[12]=20.;
+  ptbins[13]=24.;
+  ptbins[14]=ptmax;
+
+  SetGlobalIndex(nvars,nptbins);
+  SetPtBins(nptbins+1,ptbins);
+  
+  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* pt<0.5*/
+                                                 {0.400,350.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,0.},/* 0.5<pt<1*/
+                                                 {0.400,300.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.80,0.,0.},/* 1<pt<2 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.,0.},/* 2<pt<3 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 3<pt<4 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 4<pt<5 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 5<pt<6 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 6<pt<7 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-7000.*1E-8,0.85,0.,0.},/* 7<pt<8 */
+                                                 {0.400,300.*1E-4,0.9,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.,0.},/* 8<pt<12 */
+                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,10000.*1E-8,0.85,0.,0.},/* 12<pt<16 */
+                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 16<pt<20 */
+                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.},/* 20<pt<24 */
+                                                 {0.400,300.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,999999.*1E-8,0.85,0.,0.}};/* pt>24 */
+  
+  
+  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
+  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
+  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
+  
+  for (Int_t ibin=0;ibin<nptbins;ibin++){
+    for (Int_t ivar = 0; ivar<nvars; ivar++){
+      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      
+    }
+  }
+  
+  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
+  SetUseSpecialCuts(kTRUE);
+  SetRemoveDaughtersFromPrim(kTRUE);
+  
+  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
+  delete [] cutsMatrixTransposeStand;
+  cutsMatrixTransposeStand=NULL;
+
+  // PID SETTINGS
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  //pidObj->SetName("pid4D0");
+  Int_t mode=1;
+  const Int_t nlims=2;
+  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
+  Bool_t compat=kTRUE; //effective only for this mode
+  Bool_t asym=kTRUE;
+  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
+  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
+  pidObj->SetMatch(mode);
+  pidObj->SetPLimit(plims,nlims);
+  pidObj->SetSigma(sigmas);
+  pidObj->SetCompat(compat);
+  pidObj->SetTPC(kTRUE);
+  pidObj->SetTOF(kTRUE);
+  pidObj->SetPCompatTOF(1.5);
+  pidObj->SetSigmaForTPCCompat(3.);
+  pidObj->SetSigmaForTOFCompat(3.);
+  
+
+  SetPidHF(pidObj);
+  SetUsePID(kTRUE);
+  SetUseDefaultPID(kFALSE);
+
+  SetLowPt(kFALSE);
+
+  PrintAll();
+
+  delete pidObj;
+  pidObj=NULL;
+
+  return;
+
+}
+
+
+void AliRDHFCutsD0toKpi::SetStandardCutsPP2011_276TeV() {
+  //
+  // STANDARD CUTS USED FOR 2011 pp analysis at 2.76TeV
+  //
+  
+  SetName("D0toKpiCutsStandard");
+  SetTitle("Standard Cuts for D0 analysis in pp2011 at 2.76TeV run");
+
+  //
+  // Track cuts
+  //
+  AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  //default
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
+                                        AliESDtrackCuts::kAny); 
+ // default is kBoth, otherwise kAny
+  esdTrackCuts->SetMinDCAToVertexXY(0.);
+  esdTrackCuts->SetPtRange(0.3,1.e10);
+
+  esdTrackCuts->SetMaxDCAToVertexXY(1.);
+  esdTrackCuts->SetMaxDCAToVertexZ(1.);
+  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0075*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");
+
+  AddTrackCuts(esdTrackCuts);
+
+
+  const Int_t nvars=11;
+  const Int_t nptbins=13;
+  Float_t* ptbins;
+  ptbins=new Float_t[nptbins+1];
+  ptbins[0]=0.;
+  ptbins[1]=0.5;
+  ptbins[2]=1.;
+  ptbins[3]=2.;
+  ptbins[4]=3.;
+  ptbins[5]=4.;
+  ptbins[6]=5.;
+  ptbins[7]=6.;
+  ptbins[8]=8.;
+  ptbins[9]=12.;
+  ptbins[10]=16.;
+  ptbins[11]=20.;
+  ptbins[12]=24.;
+  ptbins[13]=9999.;
+
+  SetPtBins(nptbins+1,ptbins);
+
+  Float_t** rdcutsvalmine;
+  rdcutsvalmine=new Float_t*[nvars];
+  for(Int_t iv=0;iv<nvars;iv++){
+    rdcutsvalmine[iv]=new Float_t[nptbins];
+  }
+       
+  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,0.04,0.75,0.3,0.3,1000.*1E-4,1000.*1E-4,0.,0.85,0.,0.},/* pt<0.5*/
+                                                 {0.400,0.04,0.75,0.3,0.3,1000.*1E-4,1000.*1E-4,0.,0.85,0.,0.},/* 0.5<pt<1*/
+                                                 {0.400,0.03,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.8,0.,0.},/* 1<pt<2 */
+                                                 {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.0003,0.9,0.,0.},/* 2<pt<3 */
+                                                 {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.0002,0.9,0.,0.},/* 3<pt<4 */
+                                                 {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.00015,0.9,0.,0.},/* 4<pt<5 */
+                                                 {0.400,0.03,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.0001,0.9,0.,0.},/* 5<pt<6 */
+                                                 {0.400,0.09,0.85,0.7,0.7,1000.*1E-4,1000.*1E-4,0.,0.85,0.,0.},/* 6<pt<8 */
+                                                 {0.400,0.06,0.85,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.00001,0.85,0.,0.},/* 8<pt<12 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 12<pt<16 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 16<pt<20 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 20<pt<24 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.}};/* pt>24 */
+    
+  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
+  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
+  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
+  for (Int_t ibin=0;ibin<nptbins;ibin++){
+    for (Int_t ivar = 0; ivar<nvars; ivar++){
+      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];
+    }
+  }
+  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
+
+
+  //pid settings
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  Int_t mode=1;
+  const Int_t nlims=2;
+  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
+  Bool_t compat=kTRUE; //effective only for this mode
+  Bool_t asym=kTRUE;
+  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
+  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
+  pidObj->SetMatch(mode);
+  pidObj->SetPLimit(plims,nlims);
+  pidObj->SetSigma(sigmas);
+  pidObj->SetCompat(compat);
+  pidObj->SetTPC(kTRUE);
+  pidObj->SetTOF(kTRUE);
+  SetPidHF(pidObj);
+
+  SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
+
+  SetUsePID(kTRUE);
+
+  //activate pileup rejection (for pp)
+  SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
+
+  //Do recalculate the vertex
+  SetRemoveDaughtersFromPrim(kTRUE);
+
+  // Cut on the zvtx
+  SetMaxVtxZ(10.);
+  
+  // Use the kFirst selection for tracks with candidate D with pt<2
+  SetSelectCandTrackSPDFirst(kTRUE,2.);
+
+  // Use special cuts for D candidates with pt<2
+  SetUseSpecialCuts(kTRUE);
+  SetMaximumPtSpecialCuts(2.);
+
+  PrintAll();
+
+  delete pidObj;
+  pidObj=NULL;
+
+  return;
+}
+
+
+void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010() {
+  //
+  // Final CUTS USED FOR 2010 PbPb analysis of central events
+  // REMEMBER TO EVENTUALLY SWITCH ON 
+  //          SetUseAOD049(kTRUE);
+  // 
+  
+  SetName("D0toKpiCutsStandard");
+  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run");
+  
+  
+  // CENTRALITY SELECTION
+  SetMinCentrality(0.);
+  SetMaxCentrality(80.);
+  SetUseCentrality(AliRDHFCuts::kCentV0M);
+
+
+  
+  // EVENT CUTS
+  SetMinVtxContr(1);
+  // MAX Z-VERTEX CUT
+  SetMaxVtxZ(10.);
+  // PILE UP
+  SetOptPileup(AliRDHFCuts::kNoPileupSelection);
+  // PILE UP REJECTION
+  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
+
+  // TRACKS ON SINGLE TRACKS
+  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  //  esdTrackCuts->SetMinNClustersITS(4);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
+  esdTrackCuts->SetMinDCAToVertexXY(0.);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+  esdTrackCuts->SetPtRange(0.7,1.e10);
+
+  esdTrackCuts->SetMaxDCAToVertexXY(1.);  
+  esdTrackCuts->SetMaxDCAToVertexZ(1.);
+  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0075*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");  
+
+
+  AddTrackCuts(esdTrackCuts);
+  // SPD k FIRST for Pt<3 GeV/c
+  SetSelectCandTrackSPDFirst(kTRUE, 3); 
+
+  // CANDIDATE CUTS  
+  const Int_t nptbins =13;
+  const Double_t ptmax = 9999.;
+  const Int_t nvars=11;
+  Float_t ptbins[nptbins+1];
+  ptbins[0]=0.;
+  ptbins[1]=0.5;       
+  ptbins[2]=1.;
+  ptbins[3]=2.;
+  ptbins[4]=3.;
+  ptbins[5]=4.;
+  ptbins[6]=5.;
+  ptbins[7]=6.;
+  ptbins[8]=8.;
+  ptbins[9]=12.;
+  ptbins[10]=16.;
+  ptbins[11]=20.;
+  ptbins[12]=24.;
+  ptbins[13]=ptmax;
+
+  SetGlobalIndex(nvars,nptbins);
+  SetPtBins(nptbins+1,ptbins);
+  SetMinPtCandidate(2.);
+
+  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.85,0.99,2.},/* pt<0.5*/
+                                                 {0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-35000.*1E-8,0.9,0.99,2.},/* 0.5<pt<1*/
+                                                 {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-43000.*1E-8,0.85,0.995,2.},/* 1<pt<2 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-45000.*1E-8,0.95,0.998,7.},/* 2<pt<3 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-36000.*1E-8,0.95,0.998,5.},/* 3<pt<4 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-27000.*1E-8,0.95,0.998,5.},/* 4<pt<5 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-21000.*1E-8,0.92,0.998,5.},/* 5<pt<6 */
+                                                 {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-14000.*1E-8,0.88,0.998,5.},/* 6<pt<8 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.998,5.},/* 8<pt<12 */
+                                                 {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.83,0.998,8.},/* 12<pt<16 */
+                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.82,0.995,6.},/* 16<pt<20 */
+                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.81,0.995,6.},/* 20<pt<24 */
+                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,-1000.*1E-8,0.8,0.99,2.}};/* pt>24 */
+  
+  
+  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
+  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
+  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
+  
+  for (Int_t ibin=0;ibin<nptbins;ibin++){
+    for (Int_t ivar = 0; ivar<nvars; ivar++){
+      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      
+    }
+  }
+  
+  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
+  SetUseSpecialCuts(kTRUE);
+  SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb
+  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
+  delete [] cutsMatrixTransposeStand;
+  cutsMatrixTransposeStand=NULL;
+  
+  // PID SETTINGS
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  //pidObj->SetName("pid4D0");
+  Int_t mode=1;
+  const Int_t nlims=2;
+  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
+  Bool_t compat=kTRUE; //effective only for this mode
+  Bool_t asym=kTRUE;
+  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
+  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
+  pidObj->SetMatch(mode);
+  pidObj->SetPLimit(plims,nlims);
+  pidObj->SetSigma(sigmas);
+  pidObj->SetCompat(compat);
+  pidObj->SetTPC(kTRUE);
+  pidObj->SetTOF(kTRUE);
+  pidObj->SetPCompatTOF(2.);
+  pidObj->SetSigmaForTPCCompat(3.);
+  pidObj->SetSigmaForTOFCompat(3.);  
+
+
+  SetPidHF(pidObj);
+  SetUsePID(kTRUE);
+  SetUseDefaultPID(kFALSE);
+
+
+  PrintAll();
+
+
+  delete pidObj;
+  pidObj=NULL;
+
+  return;
+
+}
+
+void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010Peripherals() {
+  // CUTS USED FOR D0 RAA for the CENTRALITY RANGE 20-80%
+
+  
+  SetName("D0toKpiCutsStandard");
+  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run, for peripheral events");
+  
+  
+  // CENTRALITY SELECTION
+  SetMinCentrality(40.);
+  SetMaxCentrality(80.);
+  SetUseCentrality(AliRDHFCuts::kCentV0M);
+  
+  // EVENT CUTS
+  SetMinVtxContr(1);
+  // MAX Z-VERTEX CUT
+  SetMaxVtxZ(10.);
+  // PILE UP
+  SetOptPileup(AliRDHFCuts::kNoPileupSelection);
+  // PILE UP REJECTION
+  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
+
+  // TRACKS ON SINGLE TRACKS
+  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","default");
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  //  esdTrackCuts->SetMinNClustersITS(4);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
+  esdTrackCuts->SetMinDCAToVertexXY(0.);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+  esdTrackCuts->SetPtRange(0.5,1.e10);
+
+  esdTrackCuts->SetMaxDCAToVertexXY(1.);  
+  esdTrackCuts->SetMaxDCAToVertexZ(1.);
+  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0025*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");  
+
+
+  AddTrackCuts(esdTrackCuts);
+  // SPD k FIRST for Pt<3 GeV/c
+  SetSelectCandTrackSPDFirst(kTRUE, 3); 
+
+  // CANDIDATE CUTS  
+  const Int_t nptbins =13;
+  const Double_t ptmax = 9999.;
+  const Int_t nvars=11;
+  Float_t ptbins[nptbins+1];
+  ptbins[0]=0.;
+  ptbins[1]=0.5;       
+  ptbins[2]=1.;
+  ptbins[3]=2.;
+  ptbins[4]=3.;
+  ptbins[5]=4.;
+  ptbins[6]=5.;
+  ptbins[7]=6.;
+  ptbins[8]=8.;
+  ptbins[9]=12.;
+  ptbins[10]=16.;
+  ptbins[11]=20.;
+  ptbins[12]=24.;
+  ptbins[13]=ptmax;
+
+  SetGlobalIndex(nvars,nptbins);
+  SetPtBins(nptbins+1,ptbins);
+  SetMinPtCandidate(0.);
+
+  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,400.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-20000.*1E-8,0.7,0.993,2.},/* pt<0.5*/
+                                                 {0.400,400.*1E-4,0.8,0.5,0.5,1000.*1E-4,1000.*1E-4,-25000.*1E-8,0.85,0.993,2.},/* 0.5<pt<1*/
+                                                 {0.400,400.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.85,0.995,6.},/* 1<pt<2 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-40000.*1E-8,0.95,0.991,5.},/* 2<pt<3 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-36000.*1E-8,0.95,0.993,5.},/* 3<pt<4 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-27000.*1E-8,0.95,0.995,5.},/* 4<pt<5 */
+                                                 {0.400,250.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-21000.*1E-8,0.92,0.998,5.},/* 5<pt<6 */
+                                                 {0.400,270.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-14000.*1E-8,0.88,0.998,5.},/* 6<pt<8 */
+                                                 {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.85,0.995,5.},/* 8<pt<12 */
+                                                 {0.400,350.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.83,0.995,4.},/* 12<pt<16 */
+                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.82,0.995,4.},/* 16<pt<20 */
+                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.81,0.995,4.},/* 20<pt<24 */
+                                                 {0.400,400.*1E-4,1.0,0.7,0.7,1000.*1E-4,1000.*1E-4,0.*1E-8,0.8,0.995,4.}};/* pt>24 */
+  
+  
+  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
+  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
+  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
+  
+  for (Int_t ibin=0;ibin<nptbins;ibin++){
+    for (Int_t ivar = 0; ivar<nvars; ivar++){
+      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      
+    }
+  }
+  
+  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
+  SetUseSpecialCuts(kTRUE);
+  SetRemoveDaughtersFromPrim(kFALSE);// THIS IS VERY IMPORTANT! TOO SLOW IN PbPb
+  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
+  delete [] cutsMatrixTransposeStand;
+  cutsMatrixTransposeStand=NULL;
+  
+  // PID SETTINGS
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  //pidObj->SetName("pid4D0");
+  Int_t mode=1;
+  const Int_t nlims=2;
+  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
+  Bool_t compat=kTRUE; //effective only for this mode
+  Bool_t asym=kTRUE;
+  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
+  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
+  pidObj->SetMatch(mode);
+  pidObj->SetPLimit(plims,nlims);
+  pidObj->SetSigma(sigmas);
+  pidObj->SetCompat(compat);
+  pidObj->SetTPC(kTRUE);
+  pidObj->SetTOF(kTRUE);
+  pidObj->SetPCompatTOF(2.);
+  pidObj->SetSigmaForTPCCompat(3.);
+  pidObj->SetSigmaForTOFCompat(3.);  
+
+  SetPidHF(pidObj);
+  SetUsePID(kTRUE);
+  SetUseDefaultPID(kFALSE);
+
+  SetLowPt(kTRUE,2.);
+
+  PrintAll();
+
+
+  delete pidObj;
+  pidObj=NULL;
+
+  return;
+  
+}
+
+
+void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2011() {
+  
+  // Default 2010 PbPb cut object
+  SetStandardCutsPbPb2010();
+
+  //
+  // Enable all 2011 PbPb run triggers
+  //  
+  SetTriggerClass("");
+  ResetMaskAndEnableMBTrigger();
+  EnableCentralTrigger();
+  EnableSemiCentralTrigger();
+
+}
index 98e4321..152a16e 100644 (file)
@@ -44,6 +44,7 @@ class AliRDHFCutsD0toKpi : public AliRDHFCuts
   Int_t CombineSelectionLevels(Int_t selectionvalTrack,Int_t selectionvalCand,Int_t selectionvalPID)const;
   virtual Bool_t IsInFiducialAcceptance(Double_t pt,Double_t y) const;
   virtual void SetStandardCutsPP2010();
+  virtual void SetStandardCutsPP2011_276TeV();
   virtual void SetStandardCutsPbPb2010();
   virtual void SetStandardCutsPbPb2011();
   void SetStandardCutsPbPb2010Peripherals();