Update of Ds analysis (Sadhana, Francesco)

diff --git a/PWG3/vertexingHF/AliAnalysisTaskSEDs.cxx b/PWG3/vertexingHF/AliAnalysisTaskSEDs.cxx
index 22e7471f0e3a9198a08e23ed5f30a8120e108a5a..1576c65ca9e5eee34553c261874b12bc98dd0bba 100644 (file)
--- a/PWG3/vertexingHF/AliAnalysisTaskSEDs.cxx
+++ b/PWG3/vertexingHF/AliAnalysisTaskSEDs.cxx
@@ -51,6 +51,12 @@ AliAnalysisTaskSEDs::AliAnalysisTaskSEDs():
   AliAnalysisTaskSE(),
   fOutput(0), 
   fHistNEvents(0),
+  fPtVsMass(0),
+  fPtVsMassAC(0),
+  fYVsPt(0),
+  fYVsPtAC(0),
+  fYVsPtSig(0),
+  fYVsPtSigAC(0),
   fReadMC(kFALSE),
   fNPtBins(0),
   fListCuts(0),
@@ -67,6 +73,12 @@ AliAnalysisTaskSEDs::AliAnalysisTaskSEDs(const char *name, AliRDHFCutsDstoKKpi*
   AliAnalysisTaskSE(name),
   fOutput(0),
   fHistNEvents(0),
+  fPtVsMass(0),
+  fPtVsMassAC(0),
+  fYVsPt(0),
+  fYVsPtAC(0),
+  fYVsPtSig(0),
+  fYVsPtSigAC(0),
   fReadMC(kFALSE),
   fNPtBins(0),
   fListCuts(0),
@@ -116,11 +128,65 @@ AliAnalysisTaskSEDs::~AliAnalysisTaskSEDs()
     delete fOutput;
     fOutput = 0;
   }
+
+  if(fHistNEvents){
+    delete fHistNEvents;
+    fHistNEvents=0;
+  } 
+ 
   if (fListCuts) {
     delete fListCuts;
     fListCuts = 0;
   }
 
+  for(Int_t i=0;i<4*fNPtBins;i++){
+    
+    if(fMassHist[i]){ delete fMassHist[i]; fMassHist[i]=0;}
+    if(fMassHistCuts[i]){ delete fMassHistCuts[i]; fMassHistCuts[i]=0;}
+    if(fMassHistPhi[i]){ delete fMassHistPhi[i]; fMassHistPhi[i]=0;}
+    if(fMassHistCutsPhi[i]){ delete fMassHistCutsPhi[i]; fMassHistCutsPhi[i]=0;}
+    if(fMassHistK0st[i]){ delete fMassHistK0st[i]; fMassHistK0st[i]=0;}
+    if(fMassHistCutsK0st[i]){ delete fMassHistCutsK0st[i]; fMassHistCutsK0st[i]=0;}
+    if(fCosPHist[i]){ delete fCosPHist[i]; fCosPHist[i]=0;}
+    if(fDLenHist[i]){ delete fDLenHist[i]; fDLenHist[i]=0;}
+    if(fSumd02Hist[i]){ delete fSumd02Hist[i]; fSumd02Hist[i]=0;}
+    if(fSigVertHist[i]){ delete fSigVertHist[i]; fSigVertHist[i]=0;}
+    if(fPtMaxHist[i]){ delete fPtMaxHist[i]; fPtMaxHist[i]=0;}
+    if(fDCAHist[i]){ delete fDCAHist[i]; fDCAHist[i]=0;}
+    if(fPtProng0Hist[i]){ delete fPtProng0Hist[i]; fPtProng0Hist[i]=0;}
+    if(fPtProng1Hist[i]){ delete fPtProng1Hist[i]; fPtProng1Hist[i]=0;}
+    if(fPtProng2Hist[i]){ delete fPtProng2Hist[i]; fPtProng2Hist[i]=0;}
+    if(fDalitz[i]){ delete fDalitz[i]; fDalitz[i]=0;}
+    if(fDalitzPhi[i]){ delete fDalitzPhi[i]; fDalitzPhi[i]=0;}
+    if(fDalitzK0st[i]){ delete fDalitzK0st[i]; fDalitzK0st[i]=0;}
+
+  }
+
+  if(fPtVsMass){
+    delete fPtVsMass;
+    fPtVsMass=0;
+  }
+  if(fPtVsMassAC){
+    delete fPtVsMassAC;
+    fPtVsMassAC=0;
+  }
+  if(fYVsPt){
+    delete fYVsPt;
+    fYVsPt=0;
+  }
+  if(fYVsPtAC){
+    delete fYVsPtAC;
+    fYVsPtAC=0;
+  }
+  if(fYVsPtSig){
+    delete fYVsPtSig;
+    fYVsPtSig=0;
+  }
+  if(fYVsPtSigAC){
+    delete fYVsPtSigAC;
+    fYVsPtSigAC=0;
+  }
+
   if (fProdCuts) {
     delete fProdCuts;
     fProdCuts = 0;
@@ -165,7 +231,9 @@ void AliAnalysisTaskSEDs::UserCreateOutputObjects()
   Double_t massDs=TDatabasePDG::Instance()->GetParticle(431)->Mass();
   Int_t nInvMassBins=(Int_t)(fMassRange/fMassBinSize+0.5);
   if(nInvMassBins%2==1) nInvMassBins++;
+  // Double_t minMass=massDs-1.0*nInvMassBins*fMassBinSize;
   Double_t minMass=massDs-0.5*nInvMassBins*fMassBinSize;
+  //  Double_t maxMass=massDs+1.0*nInvMassBins*fMassBinSize;
   Double_t maxMass=massDs+0.5*nInvMassBins*fMassBinSize;
 
   TString hisname;
@@ -290,6 +358,19 @@ void AliAnalysisTaskSEDs::UserCreateOutputObjects()
   fHistNEvents->SetMinimum(0);
   fOutput->Add(fHistNEvents);
 
+  fPtVsMass=new TH2F("hPtVsMass","PtVsMass (prod. cuts)",nInvMassBins,minMass,maxMass,40,0.,20.);
+  fPtVsMassAC=new TH2F("hPtVsMassAC","PtVsMass (analysis cuts)",nInvMassBins,minMass,maxMass,40,0.,20.);  
+  fYVsPt=new TH2F("hYVsPt","YvsPt (prod. cuts)",40,0.,20.,80,-2.,2.);
+  fYVsPtAC=new TH2F("hYVsPtAC","YvsPt (analysis cuts)",40,0.,20.,80,-2.,2.);
+  fYVsPtSig=new TH2F("hYVsPtSig","YvsPt (MC, only sig., prod. cuts)",40,0.,20.,80,-2.,2.);
+  fYVsPtSigAC=new TH2F("hYVsPtSigAC","YvsPt (MC, only Sig, analysis cuts)",40,0.,20.,80,-2.,2.);
+
+  fOutput->Add(fPtVsMass);
+  fOutput->Add(fPtVsMassAC);
+  fOutput->Add(fYVsPt);
+  fOutput->Add(fYVsPtAC);
+  fOutput->Add(fYVsPtSig);
+  fOutput->Add(fYVsPtSigAC);
 
   return;
 }
@@ -406,6 +487,12 @@ void AliAnalysisTaskSEDs::UserExec(Option_t */*option*/)
       for(Int_t i=0;i<3;i++){
        if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
       }
+      
+      Double_t rapid=d->YDs(); 
+      fYVsPt->Fill(ptCand,rapid);
+      if(retCodeAnalysisCuts) fYVsPtAC->Fill(ptCand,rapid);
+ 
+      Bool_t isFidAcc=fAnalysisCuts->IsInFiducialAcceptance(ptCand,rapid);
 
       Int_t index=GetHistoIndex(iPtBin);
       Int_t type=0;
@@ -441,104 +528,124 @@ void AliAnalysisTaskSEDs::UserExec(Option_t */*option*/)
       if(isKKpi){
        index=GetHistoIndex(iPtBin);
        Double_t invMass=d->InvMassDsKKpi();
-       fMassHist[index]->Fill(invMass);
-       if(isPhi) fMassHistPhi[index]->Fill(invMass);
-       if(isK0star) fMassHistK0st[index]->Fill(invMass);
        Double_t massKK=d->InvMass2Prongs(0,1,321,321);
        Double_t massKp=d->InvMass2Prongs(1,2,321,211);
-       fDalitz[index]->Fill(massKK,massKp);
-       if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
-       if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
-       if(retCodeAnalysisCuts>0 && isKKpiAC){ 
-         fMassHistCuts[index]->Fill(invMass);
-         if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
-         if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+       if(isFidAcc){
+         fMassHist[index]->Fill(invMass);
+         fPtVsMass->Fill(invMass,ptCand);
+         if(isPhi) fMassHistPhi[index]->Fill(invMass);
+         if(isK0star) fMassHistK0st[index]->Fill(invMass);
+         fDalitz[index]->Fill(massKK,massKp);
+         if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
+         if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
+         if(retCodeAnalysisCuts>0 && isKKpiAC){ 
+           fMassHistCuts[index]->Fill(invMass);
+           fPtVsMassAC->Fill(invMass,ptCand);      
+           if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
+           if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+         }
        }
        if(fReadMC){
+         if(labDs>=0){
          Int_t labDau0=((AliAODTrack*)d->GetDaughter(0))->GetLabel();
          AliAODMCParticle* p=(AliAODMCParticle*)arrayMC->UncheckedAt(labDau0);
          Int_t pdgCode0=TMath::Abs(p->GetPdgCode());
-         if(labDs>=0){
+         //if(labDs>=0){
            if(pdgCode0==321) {   
              index=GetSignalHistoIndex(iPtBin);
+             fYVsPtSig->Fill(ptCand,rapid);
+             if(retCodeAnalysisCuts>0 && isKKpiAC) fYVsPtSigAC->Fill(ptCand,rapid);
            }else{
              index=GetReflSignalHistoIndex(iPtBin);
            }
          }else{
            index=GetBackgroundHistoIndex(iPtBin);
          }
-         fMassHist[index]->Fill(invMass);
-         if(isPhi) fMassHistPhi[index]->Fill(invMass);
-         if(isK0star) fMassHistK0st[index]->Fill(invMass);
-         fCosPHist[index]->Fill(cosp);
-         fDLenHist[index]->Fill(dlen);
-         fSigVertHist[index]->Fill(sigvert);
-         fSumd02Hist[index]->Fill(sumD02);
-         fPtMaxHist[index]->Fill(ptmax);
-         fPtCandHist[index]->Fill(ptCand);
-         fDCAHist[index]->Fill(dca);
-         fPtProng0Hist[index]->Fill(pt0);
-         fPtProng1Hist[index]->Fill(pt1);
-         fPtProng2Hist[index]->Fill(pt2);
-         fDalitz[index]->Fill(massKK,massKp);
-         if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
-         if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
-         if(retCodeAnalysisCuts>0 && isKKpiAC){
-           fMassHistCuts[index]->Fill(invMass);          
-           if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
-           if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+         if(isFidAcc){
+           fMassHist[index]->Fill(invMass);
+           if(isPhi) fMassHistPhi[index]->Fill(invMass);
+           if(isK0star) fMassHistK0st[index]->Fill(invMass);
+           fCosPHist[index]->Fill(cosp);
+           fDLenHist[index]->Fill(dlen);
+           fSigVertHist[index]->Fill(sigvert);
+           fSumd02Hist[index]->Fill(sumD02);
+           fPtMaxHist[index]->Fill(ptmax);
+           fPtCandHist[index]->Fill(ptCand);
+           fDCAHist[index]->Fill(dca);
+           fPtProng0Hist[index]->Fill(pt0);
+           fPtProng1Hist[index]->Fill(pt1);
+           fPtProng2Hist[index]->Fill(pt2);
+           fDalitz[index]->Fill(massKK,massKp);
+           if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
+           if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
+           if(retCodeAnalysisCuts>0 && isKKpiAC){
+             fMassHistCuts[index]->Fill(invMass);        
+             if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
+             if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+           }
          }
        }       
       }
       if(ispiKK){
        index=GetHistoIndex(iPtBin);
        Double_t invMass=d->InvMassDspiKK();
-       fMassHist[index]->Fill(invMass);
-       if(isPhi) fMassHistPhi[index]->Fill(invMass);
-       if(isK0star) fMassHistK0st[index]->Fill(invMass);       
        Double_t massKK=d->InvMass2Prongs(1,2,321,321);
        Double_t massKp=d->InvMass2Prongs(0,1,211,321);
-       fDalitz[index]->Fill(massKK,massKp);
-       if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
-       if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
-       if(retCodeAnalysisCuts>0 && ispiKKAC){
-         fMassHistCuts[index]->Fill(invMass);
-         if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
-         if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+
+
+       if(isFidAcc){
+         fMassHist[index]->Fill(invMass);
+         if(isPhi) fMassHistPhi[index]->Fill(invMass);
+         if(isK0star) fMassHistK0st[index]->Fill(invMass);     
+         fDalitz[index]->Fill(massKK,massKp);
+         if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
+         if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
+         if(retCodeAnalysisCuts>0 && ispiKKAC){
+           fMassHistCuts[index]->Fill(invMass);
+           fPtVsMassAC->Fill(invMass,ptCand);      
+           if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
+           if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+         }
        }
        if(fReadMC){
+          if(labDs>=0) {         
          Int_t labDau0=((AliAODTrack*)d->GetDaughter(0))->GetLabel();
          AliAODMCParticle* p=(AliAODMCParticle*)arrayMC->UncheckedAt(labDau0);
          Int_t pdgCode0=TMath::Abs(p->GetPdgCode());
-         if(labDs>=0) {          
+         // if(labDs>=0) {       
            if(pdgCode0==211) {   
              index=GetSignalHistoIndex(iPtBin);
+             fYVsPtSig->Fill(ptCand,rapid);
+             if(retCodeAnalysisCuts>0 && isKKpiAC) fYVsPtSigAC->Fill(ptCand,rapid);
            }else{
              index=GetReflSignalHistoIndex(iPtBin);
            }
          }else{
            index=GetBackgroundHistoIndex(iPtBin);
          }
-         fMassHist[index]->Fill(invMass);
-         if(isPhi) fMassHistPhi[index]->Fill(invMass);
-         if(isK0star) fMassHistK0st[index]->Fill(invMass);
-         fCosPHist[index]->Fill(cosp);
-         fDLenHist[index]->Fill(dlen);
-         fSigVertHist[index]->Fill(sigvert);
-         fSumd02Hist[index]->Fill(sumD02);
-         fPtMaxHist[index]->Fill(ptmax);
-         fPtCandHist[index]->Fill(ptCand);
-         fDCAHist[index]->Fill(dca);
-         fPtProng0Hist[index]->Fill(pt0);
-         fPtProng1Hist[index]->Fill(pt1);
-         fPtProng2Hist[index]->Fill(pt2);
-         fDalitz[index]->Fill(massKK,massKp);
-         if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
-         if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
-         if(retCodeAnalysisCuts>0 && ispiKKAC){ 
-           fMassHistCuts[index]->Fill(invMass);
-           if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
-           if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+         if(isFidAcc){
+           fMassHist[index]->Fill(invMass);
+           fPtVsMass->Fill(invMass,ptCand);
+           if(isPhi) fMassHistPhi[index]->Fill(invMass);
+           if(isK0star) fMassHistK0st[index]->Fill(invMass);
+           fCosPHist[index]->Fill(cosp);
+           fDLenHist[index]->Fill(dlen);
+           fSigVertHist[index]->Fill(sigvert);
+           fSumd02Hist[index]->Fill(sumD02);
+           fPtMaxHist[index]->Fill(ptmax);
+           fPtCandHist[index]->Fill(ptCand);
+           fDCAHist[index]->Fill(dca);
+           fPtProng0Hist[index]->Fill(pt0);
+           fPtProng1Hist[index]->Fill(pt1);
+           fPtProng2Hist[index]->Fill(pt2);
+           fDalitz[index]->Fill(massKK,massKp);
+           if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
+           if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
+           if(retCodeAnalysisCuts>0 && ispiKKAC){ 
+             fMassHistCuts[index]->Fill(invMass);
+             if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
+             if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
+           }
          }
        }
       }
@@ -564,6 +671,12 @@ void AliAnalysisTaskSEDs::Terminate(Option_t */*option*/)
     return;
   }
   fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("hNEvents"));
+  fYVsPt = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPt"));
+  fYVsPtAC = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPtAC"));
+  fYVsPtSig = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPtSig"));
+  fYVsPtSigAC = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPtSigAC"));
+  fPtVsMass = dynamic_cast<TH2F*>(fOutput->FindObject("hPtVsMass"));
+  fPtVsMassAC = dynamic_cast<TH2F*>(fOutput->FindObject("hPtVsMassAC"));
   fChanHist[0] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanAll"));
   fChanHist[1] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanSig"));
   fChanHist[2] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanBkg"));

diff --git a/PWG3/vertexingHF/AliAnalysisTaskSEDs.h b/PWG3/vertexingHF/AliAnalysisTaskSEDs.h
index 2ccfda274e4a1d45d9b777803b31b7e7f0aa3e7c..30818eb3149f38927899c68d485af1cbabb37437 100644 (file)
--- a/PWG3/vertexingHF/AliAnalysisTaskSEDs.h
+++ b/PWG3/vertexingHF/AliAnalysisTaskSEDs.h
@@ -76,6 +76,13 @@ class AliAnalysisTaskSEDs : public AliAnalysisTaskSE
   TH2F*   fDalitz[4*kMaxPtBins];      //! dalitz plot (sig,bkg,tot)
   TH2F*   fDalitzPhi[4*kMaxPtBins];   //! dalitz plot via phi (sig,bkg,tot)
   TH2F*   fDalitzK0st[4*kMaxPtBins];   //! dalitz plot via K0* (sig,bkg,tot)
+  TH2F *fPtVsMass;    //! hist. of pt vs. mass (prod. cuts)
+  TH2F *fPtVsMassAC;  //! hist. of pt vs. mass (analysis cuts)
+  TH2F *fYVsPt;       //! hist. of Y vs. Pt (prod. cuts)
+  TH2F *fYVsPtAC;     //! hist. of Y vs. Pt (analysis cuts)
+  TH2F *fYVsPtSig;    //! hist. of Y vs. Pt (MC, only sig, prod. cuts)
+  TH2F *fYVsPtSigAC;    //! hist. of Y vs. Pt (MC, only sig, analysis cuts)
+  
   Bool_t  fReadMC;                    //  flag for access to MC
   UChar_t fNPtBins;                   // number of Pt bins
   TList *fListCuts; //list of cuts
@@ -86,7 +93,7 @@ class AliAnalysisTaskSEDs : public AliAnalysisTaskSE
   AliRDHFCutsDstoKKpi *fProdCuts;     // Cuts for Analysis
   AliRDHFCutsDstoKKpi *fAnalysisCuts; // Cuts for Analysis
   
-  ClassDef(AliAnalysisTaskSEDs,5);    //  AliAnalysisTaskSE for Ds mass spectra
+  ClassDef(AliAnalysisTaskSEDs,6);    //  AliAnalysisTaskSE for Ds mass spectra
 };
 
 #endif

diff --git a/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.cxx b/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.cxx
index fecc85c0a2deda4c833680581629f83bf4590375..a2287a425d4a8f8fa5ca8af48d46f586347f5c37 100644 (file)
--- a/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.cxx
+++ b/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.cxx
@@ -85,6 +85,21 @@ AliRDHFCuts(name)
   SetVarsForOpt(7,forOpt);
   Float_t limits[2]={0,999999999.};
   SetPtBins(2,limits);
+  if(fPidHF)delete fPidHF;
+  fPidHF=new AliAODPidHF();
+  Double_t plim[2]={0.6,0.8};
+  Double_t nsigma[5]={2.,1.,2.,3.,0.};
+  
+  fPidHF->SetPLimit(plim);
+  fPidHF->SetAsym(kTRUE);
+  fPidHF->SetSigma(nsigma);
+  fPidHF->SetMatch(1);
+  fPidHF->SetTPC(1);
+  fPidHF->SetTOF(1);
+  fPidHF->SetITS(0);
+  fPidHF->SetTRD(0);
+  fPidHF->SetCompat(kTRUE);
+
 }
 //--------------------------------------------------------------------------
 AliRDHFCutsDstoKKpi::AliRDHFCutsDstoKKpi(const AliRDHFCutsDstoKKpi &source) :
@@ -234,6 +249,29 @@ void AliRDHFCutsDstoKKpi::GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,In
   
   return;
 }
+//---------------------------------------------------------------------------
+Bool_t AliRDHFCutsDstoKKpi::IsInFiducialAcceptance(Double_t pt, Double_t y) const
+{
+  //
+  // Checking if Ds is in fiducial acceptance region 
+  //
+
+  if(pt > 5.) {
+    // applying cut for pt > 5 GeV
+    AliDebug(2,Form("pt of Ds = %f (> 5), cutting at |y| < 0.8",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 Ds = %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 AliRDHFCutsDstoKKpi::IsSelectedPID(AliAODRecoDecayHF *rd) {
   // PID selection
@@ -281,7 +319,7 @@ Int_t AliRDHFCutsDstoKKpi::IsSelectedPID(AliAODRecoDecayHF *rd) {
 }
 
 //---------------------------------------------------------------------------
-Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel) {
+Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel, AliAODEvent* aod) {
   //
   // Apply selection
   //
@@ -318,11 +356,33 @@ Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel) {
   Bool_t okPidDsKKpi=returnvaluePID&1;
   Bool_t okPidDspiKK=returnvaluePID&2;
  
-  Int_t returnvalue=1;
+  Int_t returnvalueCuts=15;
   // selection on candidate
   if(selectionLevel==AliRDHFCuts::kAll || 
      selectionLevel==AliRDHFCuts::kCandidate) {
-    
+    //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");
+       return 0;
+      }
+      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;
+       }
+       return 0;
+      }
+      //set recalculed primary vertex
+      d->SetOwnPrimaryVtx(recvtx);
+      delete recvtx; recvtx=NULL;
+    }
 
     Int_t okDsKKpi=1;
     Int_t okDspiKK=1;
@@ -337,13 +397,13 @@ Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel) {
     Double_t mDspiKK=d->InvMassDspiKK();
     if(TMath::Abs(mDsKKpi-mDsPDG)>fCutsRD[GetGlobalIndex(0,ptbin)]) okDsKKpi = 0;
     if(TMath::Abs(mDspiKK-mDsPDG)>fCutsRD[GetGlobalIndex(0,ptbin)]) okDspiKK = 0;
-    if(!okDsKKpi && !okDspiKK) return 0;
-    if(okPidDsKKpi && !okDsKKpi)  return 0;
-    if(okPidDspiKK && !okDspiKK) return 0;
+    if(!okDsKKpi && !okDspiKK) returnvalueCuts=0;
+    if(okPidDsKKpi && !okDsKKpi)  returnvalueCuts=0;
+    if(okPidDspiKK && !okDspiKK) returnvalueCuts=0;
 
     //single track
     if(TMath::Abs(d->PtProng(1)) < fCutsRD[GetGlobalIndex(1,ptbin)] || 
-       TMath::Abs(d->Getd0Prong(1))<fCutsRD[GetGlobalIndex(3,ptbin)]) return 0;
+       TMath::Abs(d->Getd0Prong(1))<fCutsRD[GetGlobalIndex(3,ptbin)]) returnvalueCuts=0;
     if(okDsKKpi){
       if(TMath::Abs(d->PtProng(0)) < fCutsRD[GetGlobalIndex(1,ptbin)] || 
         TMath::Abs(d->Getd0Prong(0))<fCutsRD[GetGlobalIndex(3,ptbin)]) okDsKKpi=0;
@@ -356,7 +416,7 @@ Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel) {
       if(TMath::Abs(d->PtProng(2)) < fCutsRD[GetGlobalIndex(1,ptbin)] || 
         TMath::Abs(d->Getd0Prong(2))<fCutsRD[GetGlobalIndex(3,ptbin)]) okDspiKK=0;
     }
-    if(!okDsKKpi && !okDspiKK) return 0;
+    if(!okDsKKpi && !okDspiKK) returnvalueCuts=0;
 
     // cuts on resonant decays (via Phi or K0*)
     Double_t mPhiPDG = TDatabasePDG::Instance()->GetParticle(333)->Mass();
@@ -375,32 +435,44 @@ Int_t AliRDHFCutsDstoKKpi::IsSelected(TObject* obj,Int_t selectionLevel) {
       if(TMath::Abs(mass12phi-mPhiPDG)<fCutsRD[GetGlobalIndex(12,ptbin)]) okMassPhi=1;
       if(!okMassPhi && !okMassK0star) okDspiKK=0;
     }
-    if(!okDsKKpi && !okDspiKK) return 0;
+    if(!okDsKKpi && !okDspiKK) returnvalueCuts=0;
 
     // Cuts on track pairs
-    for(Int_t i=0;i<3;i++) if(d->GetDCA(i)>fCutsRD[GetGlobalIndex(11,ptbin)])  return 0;
+    for(Int_t i=0;i<3;i++) if(d->GetDCA(i)>fCutsRD[GetGlobalIndex(11,ptbin)])  returnvalueCuts=0;
     if(d->GetDist12toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)] || 
-       d->GetDist23toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)]) return 0;
+       d->GetDist23toPrim()<fCutsRD[GetGlobalIndex(5,ptbin)]) returnvalueCuts=0;
 
 
     // Cuts on candidate triplet
-    if(d->GetSigmaVert()>fCutsRD[GetGlobalIndex(6,ptbin)]) return 0;
-    if(d->DecayLength()<fCutsRD[GetGlobalIndex(7,ptbin)]) return 0;
+    if(d->GetSigmaVert()>fCutsRD[GetGlobalIndex(6,ptbin)]) returnvalueCuts=0;
+    if(d->DecayLength()<fCutsRD[GetGlobalIndex(7,ptbin)]) returnvalueCuts=0;
     if(TMath::Abs(d->PtProng(0))<fCutsRD[GetGlobalIndex(8,ptbin)] && 
        TMath::Abs(d->PtProng(1))<fCutsRD[GetGlobalIndex(8,ptbin)] && 
-       TMath::Abs(d->PtProng(2))<fCutsRD[GetGlobalIndex(8,ptbin)]) return 0;
-    if(d->CosPointingAngle()< fCutsRD[GetGlobalIndex(9,ptbin)])return 0;
+       TMath::Abs(d->PtProng(2))<fCutsRD[GetGlobalIndex(8,ptbin)]) returnvalueCuts=0;
+    if(d->CosPointingAngle()< fCutsRD[GetGlobalIndex(9,ptbin)])returnvalueCuts=0;
     Double_t sum2=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
-    if(sum2<fCutsRD[GetGlobalIndex(10,ptbin)])return 0;
-
-    returnvalue=0;
+    if(sum2<fCutsRD[GetGlobalIndex(10,ptbin)])returnvalueCuts=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");
+    }
+   
+    if(returnvalueCuts == 0) return 0;
+    Int_t returnvalue=0;
     if(okDsKKpi) returnvalue+=1;
     if(okDspiKK) returnvalue+=2;
     if(okMassPhi) returnvalue+=4;
     if(okMassK0star) returnvalue+=8;
 
+    return returnvalue;
   }
-  return returnvalue;
+  return returnvalueCuts;
 
 }
 //---------------------------------------------------------------------------

diff --git a/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.h b/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.h
index 7b78bb52eee995d33a652c4bd8c82e29296c6522..b73e5c55eee3c0f0e7b2d57dee1c11669168f40d 100644 (file)
--- a/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.h
+++ b/PWG3/vertexingHF/AliRDHFCutsDstoKKpi.h
@@ -25,9 +25,15 @@ class AliRDHFCutsDstoKKpi : public AliRDHFCuts
   virtual void GetCutVarsForOpt(AliAODRecoDecayHF *d,Float_t *vars,Int_t nvars,Int_t *pdgdaughters);
 
   using AliRDHFCuts::IsSelected;
-  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel);
+  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel){
+    return IsSelected(obj,selectionLevel,0x0);
+  }
+  virtual Int_t IsSelected(TObject* obj,Int_t selectionLevel,AliAODEvent* aod);
+
+
   virtual Int_t IsSelectedPID(AliAODRecoDecayHF *rd);
-  
+   
+  virtual Bool_t IsInFiducialAcceptance(Double_t pt,Double_t y) const;
   Float_t GetMassCut(Int_t iPtBin=0) const { return (GetCuts() ? fCutsRD[GetGlobalIndex(0,iPtBin)] : 1.e6);}
   Float_t GetDCACut(Int_t iPtBin=0) const { return (GetCuts() ? fCutsRD[GetGlobalIndex(11,iPtBin)] : 1.e6);}