]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG3/vertexingHF/AliRDHFCutsD0toKpi.cxx
Update (Alessandro)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliRDHFCutsD0toKpi.cxx
index 3953f712a6e47418e8f8d3413437e21230155f7f..67a01959123353966d3c4901e94afaa1fa840579 100644 (file)
@@ -13,6 +13,8 @@
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
+/* $Id$ */
+
 /////////////////////////////////////////////////////////////
 //
 // Class for cuts on AOD reconstructed D0->Kpi
@@ -38,44 +40,50 @@ ClassImp(AliRDHFCutsD0toKpi)
 //--------------------------------------------------------------------------
 AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const char* name) : 
 AliRDHFCuts(name),
-fUseSpecialCuts(kFALSE),
+fUseSpecialCuts(kTRUE),
 fLowPt(kTRUE),
-fDefaultPID(kTRUE),
+fDefaultPID(kFALSE),
 fUseKF(kFALSE)
 {
   //
   // Default Constructor
   //
-  Int_t nvars=9;
+  Int_t nvars=11;
   SetNVars(nvars);
-  TString varNames[9]={"inv. mass [GeV]",   
-                      "dca [cm]",
-                      "cosThetaStar", 
-                      "pTK [GeV/c]",
-                      "pTPi [GeV/c]",
-                      "d0K [cm]",
-                      "d0Pi [cm]",
-                      "d0d0 [cm^2]",
-                      "cosThetaPoint"};
-  Bool_t isUpperCut[9]={kTRUE,
-                       kTRUE,
-                       kTRUE,
-                       kFALSE,
-                       kFALSE,
-                       kTRUE,
-                       kTRUE,
-                       kTRUE,
-                       kFALSE};
+  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[9]={kFALSE,
-                   kTRUE,
-                   kTRUE,
-                   kFALSE,
-                   kFALSE,
-                   kFALSE,
-                   kFALSE,
-                   kTRUE,
-                   kTRUE};
+  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);
@@ -112,7 +120,7 @@ AliRDHFCutsD0toKpi &AliRDHFCutsD0toKpi::operator=(const AliRDHFCutsD0toKpi &sour
 
 
 //---------------------------------------------------------------------------
-void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters) {
+void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters,AliAODEvent *aod) {
   // 
   // Fills in vars the values of the variables 
   //
@@ -124,6 +132,18 @@ void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int
 
   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++;
@@ -190,6 +210,18 @@ void AliRDHFCutsD0toKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int
     vars[iter]=dd->CosPointingAngle();
   }
   
+  if(fVarsForOpt[9]){
+               iter++;
+               vars[iter]=TMath::Abs(dd->CosPointingAngleXY());
+       }
+  
+   if(fVarsForOpt[10]){
+               iter++;
+          vars[iter]=(dd->NormalizedDecayLengthXY()*(dd->P()/dd->Pt()));
+       }
+
+   if(cleanvtx)CleanOwnPrimaryVtx(dd,aod,origownvtx);
+
   return;
 }
 //---------------------------------------------------------------------------
@@ -197,7 +229,11 @@ Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEve
   //
   // Apply selection
   //
-  
+
+  fIsSelectedCuts=0;
+  fIsSelectedPID=0;
+
   if(!fCutsRD){
     cout<<"Cut matrice not inizialized. Exit..."<<endl;
     return 0;
@@ -210,108 +246,89 @@ Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEve
     return 0;
   }
 
-  // selection on daughter tracks 
-  if(selectionLevel==AliRDHFCuts::kAll || 
-     selectionLevel==AliRDHFCuts::kTracks) {
-    if(!AreDaughtersSelected(d)) 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;
  
   // returnvalue: 0 not sel, 1 only D0, 2 only D0bar, 3 both
   Int_t returnvaluePID=3;
   Int_t returnvalueCuts=3;
 
-  // selection on PID 
-  if(selectionLevel==AliRDHFCuts::kAll || 
-     selectionLevel==AliRDHFCuts::kCandidate ||
-     selectionLevel==AliRDHFCuts::kPID) {
-    returnvaluePID = IsSelectedPID(d);
-  }
-
-
-
-
   // selection on candidate
   if(selectionLevel==AliRDHFCuts::kAll || 
      selectionLevel==AliRDHFCuts::kCandidate) {
 
-    // go to selection with Kalman vertexing, if requested
-    if(fUseKF) {
-      returnvalueCuts = IsSelectedKF(d,aod);
-    } else {
+    if(!fUseKF) {
 
       //recalculate vertex w/o daughters
       AliAODVertex *origownvtx=0x0;
-      AliAODVertex *recvtx=0x0;
-      
-      if(fRemoveDaughtersFromPrimary) {
-       if(!aod) {
-         AliError("Can not remove daughters from vertex without AOD event");
+      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());
-       recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
-       if(!recvtx){
-         AliDebug(2,"Removal of daughter tracks failed");
-         //recvtx=d->GetPrimaryVtx();
-         if(origownvtx){
-           delete origownvtx;
-           origownvtx=NULL;
-         }
+       if(!SetMCPrimaryVtx(d,aod)) {
+         CleanOwnPrimaryVtx(d,aod,origownvtx);
          return 0;
        }
-       //set recalculed primary vertex
-       d->SetOwnPrimaryVtx(recvtx);
-       delete recvtx; recvtx=NULL;
       }
       
-      
       Double_t pt=d->Pt();
       
       Int_t okD0=0,okD0bar=0;
       
       Int_t ptbin=PtBin(pt);
       if (ptbin==-1) {
-       if(origownvtx){
-         d->SetOwnPrimaryVtx(origownvtx);
-         delete origownvtx;
-         origownvtx=NULL;
-       }
-       else d->UnsetOwnPrimaryVtx();
+       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->PtProng(1) < fCutsRD[GetGlobalIndex(3,ptbin)] || d->PtProng(0) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0 = 0;
-      if(d->PtProng(0) < fCutsRD[GetGlobalIndex(3,ptbin)] || d->PtProng(1) < fCutsRD[GetGlobalIndex(4,ptbin)]) okD0bar = 0;
-      if(!okD0 && !okD0bar) returnvalueCuts=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)  returnvalueCuts=0;
+      if(!okD0 && !okD0bar)  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
       
-      if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)])  returnvalueCuts=0;
+      if(d->GetDCA() > fCutsRD[GetGlobalIndex(1,ptbin)])  {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
       
-      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)  returnvalueCuts=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)   returnvalueCuts=0;
-    
-      if(d->Prodd0d0() > fCutsRD[GetGlobalIndex(7,ptbin)])  returnvalueCuts=0;
+      if(!okD0 && !okD0bar)   {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
     
-      if(d->CosPointingAngle() < fCutsRD[GetGlobalIndex(8,ptbin)])  returnvalueCuts=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()*(d->P()/d->Pt()));
+      if (normalDecayLengXY < fCutsRD[GetGlobalIndex(10, ptbin)]) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
       
       if (returnvalueCuts!=0) {
        if (okD0) returnvalueCuts=1; //cuts passed as D0
@@ -322,25 +339,42 @@ Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEve
       // call special cuts
       Int_t special=1;
       if(fUseSpecialCuts) special=IsSelectedSpecialCuts(d);
-      if(!special) returnvalueCuts=0;
+      if(!special) {CleanOwnPrimaryVtx(d,aod,origownvtx); return 0;}
 
       // unset recalculated primary vertex when not needed any more
-      if(origownvtx) {
-       d->SetOwnPrimaryVtx(origownvtx);
-       delete origownvtx;
-       origownvtx=NULL;
-      } else if(fRemoveDaughtersFromPrimary) {
-       d->UnsetOwnPrimaryVtx();
-       AliDebug(3,"delete new vertex\n");
-      }
+      CleanOwnPrimaryVtx(d,aod,origownvtx);
+
+    } else {
+      // go to selection with Kalman vertexing, if requested
+      returnvalueCuts = IsSelectedKF(d,aod);
+    }
 
-    } // if(fUseKF)
+    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 CombineSelectionLevels(3,returnvalueCuts,returnvaluePID);
+  return returnvalueComb;
 }
 
 //------------------------------------------------------------------------------------------
@@ -649,10 +683,11 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d)
   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;
+      checkPIDInfo[daught]=kFALSE; 
+      continue;
     }
 
     // identify kaon
@@ -670,26 +705,27 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d)
      }
 
 
-   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(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()<2.){
      Double_t sigmaTPC[3]={3.,2.,0.};
@@ -740,7 +776,7 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d)
     return 0;
   }
 
-  if(d->Pt()<2.){
+  if(fLowPt && d->Pt()<2.){
     if(fLowPt) fPidHF->SetSigmaForTPC(sigma_tmp);
     if(combinedPID[0][0]<=0&&combinedPID[1][0]<=0){
       fWhyRejection=32;// reject cases where the Kaon is not identified
@@ -884,6 +920,7 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedPIDdefault(AliAODRecoDecayHF* d)
      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;
        }
@@ -1036,15 +1073,13 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedSpecialCuts(AliAODRecoDecayHF *d) const
   for(Int_t i=0;i<2/*prongs*/;i++){
     if(TMath::Abs(d->Normalizedd0Prong(i))<normd0Cut) return 0; //normd0Cut not passed
   }
-  if(d->DecayLength()<decLengthCut)  return 0; //decLengthCut not passed
-  if(d->NormalizedDecayLength()<normDecLengthCut)  return 0; //decLengthCut 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
@@ -1086,9 +1121,259 @@ void AliRDHFCutsD0toKpi::SetUseKF(Bool_t useKF)
                    kFALSE,
                    kTRUE,
                    kTRUE,
-                    kFALSE,
-                    kFALSE};
+                   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);
+
+  
+  // 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);
+
+
+  PrintAll();
+
+  delete pidObj;
+  pidObj=NULL;
+
+  return;
+
+}
+
+
+void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010() {
+  //
+  //PRELIMINARY CUTS USED FOR 2010 PbPb analysis
+  //... EVOLVING SOON 
+  // 
+  
+  SetName("D0toKpiCutsStandard");
+  SetTitle("Standard Cuts for D0 analysis in PbPb2010 run");
+  
+  // PILE UP REJECTION
+  //SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
+  // CENTRALITY SELECTION
+  SetMinCentrality(0.);
+  SetMaxCentrality(80.);
+  SetUseCentrality(AliRDHFCuts::kCentV0M);
+
+
+  // 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.8,1.e10);
+
+  esdTrackCuts->SetMaxDCAToVertexXY(1.);  
+  esdTrackCuts->SetMaxDCAToVertexZ(1.);
+  esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0100*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");  
+
+
+  AddTrackCuts(esdTrackCuts);
+
+  // 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,-50000.*1E-8,0.85,0.,5.},/* pt<0.5*/
+                                                 {0.400,400.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-50000.*1E-8,0.85,0.,5.},/* 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.,5.},/* 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.998,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.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,5.},/* 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.998,5.},/* 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.998,5.},/* 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.998,5.}};/* 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);
+  
+  SetPidHF(pidObj);
+  SetUsePID(kTRUE);
+  SetUseDefaultPID(kFALSE);
+
+
+  PrintAll();
+
+
+  delete pidObj;
+  pidObj=NULL;
+
+  return;
+
+}