Major update of the HFE package (comments inside the code
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEV0pid.cxx
index a2b1eed..658fe1f 100644 (file)
@@ -34,6 +34,9 @@
 #include "AliKFVertex.h"
 #include "AliVEvent.h"
 #include "AliVTrack.h"
+#include "AliMCParticle.h"
+#include "AliStack.h"
+#include "AliMCEvent.h"
 
 #include "AliHFEV0cuts.h"
 #include "AliHFEV0info.h"
@@ -46,6 +49,9 @@ ClassImp(AliHFEV0pid)
 AliHFEV0pid::AliHFEV0pid():
   TObject()
   , fInputEvent(NULL)
+  , fNtracks(0)
+  , fMCEvent(NULL)
+  , fMCon(kFALSE)
   , fPrimaryVertex(NULL)
   , fElectrons(NULL)
   , fPionsK0(NULL)
@@ -120,14 +126,15 @@ void AliHFEV0pid::InitQA(){
   fV0cuts->Init("V0cuts");
 
   if(!fQA){
-    fQA = new AliHFEcollection("v0pidQA", "QA histograms for V0 PID");
+    fQA = new AliHFEcollection("v0pidQA", "QA histograms for V0 selection");
 
     fQA->CreateTH1F("h_nV0s", "No. of found and accepted V0s", 5, -0.5, 4.5);
 
     // QA histograms for invariant mass
     fQA->CreateTH1F("h_InvMassGamma", "Gamma invariant mass; inv mass [GeV/c^{2}]; counts", 100, 0, 0.25);
     fQA->CreateTH1F("h_InvMassK0s", "K0s invariant mass; inv mass [GeV/c^{2}]; counts", 200, 0.4, 0.65);
-    fQA->CreateTH1F("h_InvMassLambda", "Lambda invariant mass; inv mass [GeV/c^{2}]; counts", 100, 1.05, 1.15);
+    fQA->CreateTH1F("h_InvMassL", "Lambda invariant mass; inv mass [GeV/c^{2}]; counts", 100, 1.05, 1.15);
+    fQA->CreateTH1F("h_InvMassAL", "Lambda invariant mass; inv mass [GeV/c^{2}]; counts", 100, 1.05, 1.15);
     
     // QA histograms for p distribution (of the daughters)
     fQA->CreateTH1F("h_P_electron", "P distribution of the gamma electrons; p (GeV/c); counts", 100, 0.1, 10);
@@ -138,8 +145,30 @@ void AliHFEV0pid::InitQA(){
     // QA pt of the V0
     fQA->CreateTH1F("h_Pt_Gamma", "Pt of the gamma conversion; p_{T} (GeV/c); counts", 100, 0, 10);
     fQA->CreateTH1F("h_Pt_K0", "Pt of the K0; p_{T} (GeV/c); counts", 100, 0, 10);
-    fQA->CreateTH1F("h_Pt_Lambda", "Pt of the Lambda; p_{T} (GeV/c); counts", 100, 0, 10);    
+    fQA->CreateTH1F("h_Pt_L", "Pt of the Lambda; p_{T} (GeV/c); counts", 100, 0, 10);    
+    fQA->CreateTH1F("h_Pt_AL", "Pt of the Lambda; p_{T} (GeV/c); counts", 100, 0, 10);    
     
+    // Armenteros plot V0 preselection
+    fQA->CreateTH2F("h_AP_all_V0s", "armenteros plot for all V0 candidates; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
+    fQA->CreateTH2F("h_AP_selected_V0s", "armenteros plot for all V0 candidates; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
+
+    //
+    // !!! MC plots !!!
+    // 
+    fQA->CreateTH2F("h_AP_MC_all_V0s", "armenteros plot for all MC tagged V0s; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
+    fQA->CreateTH2F("h_AP_MC_Gamma", "armenteros plot for MC tagged Gammas; #alpha; Q_{T}",  200, -1, 1, 200, 0, 0.25);
+    fQA->CreateTH2F("h_AP_MC_K0", "armenteros plot for MC tagged K0s; #alpha; Q_{T}",  200, -1, 1, 200, 0, 0.25);
+    fQA->CreateTH2F("h_AP_MC_Lambda", "armenteros plot for MC tagged Lambdas; #alpha; Q_{T}",  200, -1, 1, 200, 0, 0.25);
+    fQA->CreateTH2F("h_AP_MC_ALambda", "armenteros plot for MC tagged A-Lambdass; #alpha; Q_{T}",  200, -1, 1, 200, 0, 0.25);
+
+   
+    // armenteros plots for different V0 momenta - MC signal only
+    fQA->CreateTH2Fvector1(12, "h_AP_MC_Gamma_p", "armenteros plot for MC tagged Gammas; #alpha; Q_{T}",  200, -1., 1., 200, 0., 0.25);
+    fQA->CreateTH2Fvector1(12, "h_AP_MC_K0_p", "armenteros plot for MC tagged K0s; #alpha; Q_{T}",  200, -1., 1., 200, 0., 0.25);
+    fQA->CreateTH2Fvector1(12, "h_AP_MC_Lambda_p", "armenteros plot for MC tagged Lambdas; #alpha; Q_{T}",  200, -1., 1., 200, 0., 0.25);
+   //
+
     
   }
 }
@@ -154,11 +183,13 @@ void AliHFEV0pid::Process(AliVEvent * const inputEvent){
   
   Int_t nGamma = 0, nK0s = 0, nLambda = 0, nPhi = 0;
   fInputEvent = inputEvent;
+  fNtracks = fInputEvent->GetNumberOfTracks();
   fIndices->Init(fInputEvent->GetNumberOfV0s() * 2);
   fPrimaryVertex = new AliKFVertex(*(fInputEvent->GetPrimaryVertex()));
   if(!fPrimaryVertex) return;
   fV0cuts->SetInputEvent(fInputEvent);
   fV0cuts->SetPrimaryVertex(fPrimaryVertex);
+  if(fMCEvent) fV0cuts->SetMCEvent(fMCEvent);
   Int_t v0status = 0;
   for(Int_t iv0 = 0; iv0 < fInputEvent->GetNumberOfV0s(); iv0++){
     if(!TString(fInputEvent->IsA()->GetName()).CompareTo("AliESDEvent")){
@@ -174,14 +205,14 @@ void AliHFEV0pid::Process(AliVEvent * const inputEvent){
       AliAODv0 *aodV0 = (dynamic_cast<AliAODEvent *>(fInputEvent))->GetV0(iv0);
       if(aodV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
       v0status = ProcessV0(aodV0);
-      if(kUndef != v0status){
+      if(AliHFEV0cuts::kUndef != v0status){
       }
     }
     switch(v0status){
-    case kRecoGamma: nGamma++; break;
-    case kRecoK0s: nK0s++; break;
-    case kRecoPhi: nPhi++; break;  
-    case kRecoLambda: nLambda++; break;
+    case AliHFEV0cuts::kRecoGamma: nGamma++; break;
+    case AliHFEV0cuts::kRecoK0: nK0s++; break;
+    case AliHFEV0cuts::kRecoPhi: nPhi++; break;  
+    case AliHFEV0cuts::kRecoLambda: nLambda++; break;
     };
   }
 
@@ -210,33 +241,46 @@ Int_t AliHFEV0pid::ProcessV0(TObject *v0){
   AliVTrack* daughter[2];
   daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack((dynamic_cast<AliESDv0 *>(v0))->GetPindex()));
   daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack((dynamic_cast<AliESDv0 *>(v0))->GetNindex()));
-  if(!daughter[0] || !daughter[1]) return kUndef;
+  if(!daughter[0] || !daughter[1]) return AliHFEV0cuts::kUndef;
+
+  if(fMCEvent) fMCon = kTRUE;
+  //printf("-D: fMCEvent %x, fMCon: %i\n", fMCEvent, fMCon);
+
+  Int_t dMC[2] = {-1, -1};
+  Int_t idMC = AliHFEV0cuts::kUndef; 
 
-  
   if(IsESDanalysis()){
+    if(fMCon){
+      idMC = IdentifyV0(v0, dMC);
+      //printf("--D: V0 pid: %i, P: %i, N: %i\n", id, d[0], d[1]);
+      fV0cuts->SetCurrentV0id(idMC);
+      fV0cuts->SetDaughtersID(dMC);
+    }
+    // check common single track cuts
     for(Int_t i=0; i<2; ++i){
-      // check common single track cuts
-      if(!fV0cuts->TrackCutsCommon(dynamic_cast<AliESDtrack*>(daughter[i]))) return kUndef;
-    }    
+      if(!fV0cuts->TrackCutsCommon(dynamic_cast<AliESDtrack*>(daughter[i]))) return AliHFEV0cuts::kUndef;
+    }
     // check commom V0 cuts
-    if(!fV0cuts->V0CutsCommon(dynamic_cast<AliESDv0 *>(v0))) return kUndef;
+    if(!fV0cuts->V0CutsCommon(dynamic_cast<AliESDv0 *>(v0))) return AliHFEV0cuts::kUndef;
   }
 
+  // preselect the V0 candidates based on the Armenteros plot
+  Int_t id = PreselectV0((dynamic_cast<AliESDv0 *>(v0)), idMC);
 
   // store the resutls
-  if(IsGammaConv(v0)){
-    fQA->Fill("h_nV0s", kRecoGamma);
-    return kRecoGamma;
+  if(AliHFEV0cuts::kRecoGamma == id && IsGammaConv(v0)){
+    fQA->Fill("h_nV0s", AliHFEV0cuts::kRecoGamma);
+    return AliHFEV0cuts::kRecoGamma;
   }
-  else if(IsK0s(v0)){
-    fQA->Fill("h_nV0s", kRecoK0s);
-    return kRecoK0s;
+  else if(AliHFEV0cuts::kRecoK0 == id && IsK0s(v0)){
+    fQA->Fill("h_nV0s", AliHFEV0cuts::kRecoK0);
+    return AliHFEV0cuts::kRecoK0;
   }
-  else if(IsLambda(v0)){
-    fQA->Fill("h_nV0s", kRecoLambda);    
-    return kRecoLambda;
+  else if(AliHFEV0cuts::kRecoLambda == TMath::Abs(id) && IsLambda(v0)){
+    fQA->Fill("h_nV0s", AliHFEV0cuts::kRecoLambda);    
+    return AliHFEV0cuts::kRecoLambda;
   }
-  else return kUndef;
+  else return AliHFEV0cuts::kUndef;
     
 }
 //____________________________________________________________
@@ -252,6 +296,111 @@ void AliHFEV0pid::Flush(){
   fIndices->Flush();
 }
 //____________________________________________________________
+Int_t AliHFEV0pid::PreselectV0(AliESDv0 * const v0, Int_t idMC){
+  //
+  // Based on Armenteros plot preselet the possible identity of the V0 candidate
+  //
+
+  if(!v0) return -1;
+
+  // momentum dependent armenteros plots
+  ArmenterosPlotMC(v0, idMC);
+
+  // comupte the armenteros variables
+  Float_t ar[2];
+  fV0cuts->Armenteros(v0, ar);
+  // for clarity
+  const Float_t alpha = ar[0];
+  const Float_t qt = ar[1];
+  //printf(" -D: Alpha: %f, QT: %f \n", alpha, qt);
+
+  if(TMath::Abs(alpha) > 1) return AliHFEV0cuts::kUndef;
+
+  fQA->Fill("h_AP_all_V0s", alpha, qt);
+
+  // fill a few MC tagged histograms - AP plots
+  if(fMCEvent){
+    switch(idMC){
+    case AliHFEV0cuts::kRecoGamma :
+      fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
+      fQA->Fill("h_AP_MC_Gamma", alpha, qt);
+      break;
+    case AliHFEV0cuts::kRecoK0 :
+      fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
+    fQA->Fill("h_AP_MC_K0", alpha, qt);
+      break;
+    case AliHFEV0cuts::kRecoLambda :
+      fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
+      fQA->Fill("h_AP_MC_Lambda", alpha, qt);
+      break;
+    case AliHFEV0cuts::kRecoALambda :
+      fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
+      fQA->Fill("h_AP_MC_ALambda", alpha, qt);
+      break;
+    }
+  }
+
+  // Gamma cuts
+  const Double_t cutAlphaG = 0.35; 
+  const Double_t cutQTG = 0.05;
+  const Double_t cutAlphaG2[2] = {0.6, 0.8};
+  const Double_t cutQTG2 = 0.04;
+
+  // K0 cuts
+  const Float_t cutQTK0[2] = {0.1075, 0.215};
+  const Float_t cutAPK0[2] = {0.199, 0.8};   // parameters for curved QT cut
+  
+  // Lambda & A-Lambda cuts
+  const Float_t cutQTL = 0.03;
+  const Float_t cutAlphaL[2] = {0.35, 0.7};
+  const Float_t cutAlphaAL[2] = {-0.7,  -0.35};
+  const Float_t cutAPL[3] = {0.107, -0.69, 0.5};  // parameters fir curved QT cut
+
+
+  // Check for Gamma candidates
+  if(qt < cutQTG){
+    if( (TMath::Abs(alpha) < cutAlphaG) ){
+      fQA->Fill("h_AP_selected_V0s", alpha, qt);
+      return  AliHFEV0cuts::kRecoGamma;
+    }
+  }
+  // additional region - should help high pT gammas
+  if(qt < cutQTG2){
+    if( (TMath::Abs(alpha) > cutAlphaG2[0]) &&  (TMath::Abs(alpha) < cutAlphaG2[1]) ){
+      fQA->Fill("h_AP_selected_V0s", alpha, qt);
+      return  AliHFEV0cuts::kRecoGamma;
+    }
+  }
+
+
+  // Check for K0 candidates
+  Float_t q = cutAPK0[0] * TMath::Sqrt(TMath::Abs(1 - alpha*alpha/(cutAPK0[1]*cutAPK0[1])));
+  if( (qt > cutQTK0[0]) && (qt < cutQTK0[1]) && (qt > q) ){
+    fQA->Fill("h_AP_selected_V0s", alpha, qt);
+    return AliHFEV0cuts::kRecoK0;
+  }
+  
+  if( (alpha > 0) && (alpha > cutAlphaL[0])  && (alpha < cutAlphaL[1]) && (qt > cutQTL)){
+    q = cutAPL[0] * TMath::Sqrt(1 - ( (alpha + cutAPL[1]) * (alpha + cutAPL[1]))  / (cutAPL[2]*cutAPL[2]) );
+    if( qt < q  ){
+      fQA->Fill("h_AP_selected_V0s", alpha, qt);
+      return AliHFEV0cuts::kRecoLambda;
+    }
+  }
+
+  // Check for A-Lambda candidates
+  if( (alpha < 0) && (alpha > cutAlphaAL[0]) && (alpha < cutAlphaAL[1]) && (qt > cutQTL)){
+    q = cutAPL[0] * TMath::Sqrt(1 - ( (alpha - cutAPL[1]) * (alpha - cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) );
+    if( qt < q ){
+      fQA->Fill("h_AP_selected_V0s", alpha, qt);
+      return AliHFEV0cuts::kRecoLambda;
+    }
+  }
+  
+  return AliHFEV0cuts::kUndef;
+
+}
+//____________________________________________________________
 Bool_t AliHFEV0pid::IsGammaConv(TObject *v0){
   //
   // Identify Gamma
@@ -294,12 +443,12 @@ Bool_t AliHFEV0pid::IsGammaConv(TObject *v0){
   if(!fIndices->Find(daughter[0]->GetID())){
     AliDebug(1, Form("Gamma identified, daughter IDs: %d,%d", daughter[0]->GetID(), daughter[1]->GetID()));    
     fElectrons->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
-    fIndices->Add(daughter[0]->GetID(), AliHFEV0pid::kRecoElectron);
+    fIndices->Add(daughter[0]->GetID(), AliHFEV0cuts::kRecoElectron);
   }
   if(!fIndices->Find(daughter[1]->GetID())){
     AliDebug(1, Form("Gamma identified, daughter IDs: %d,%d", daughter[1]->GetID(), daughter[1]->GetID()));
     fElectrons->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
-    fIndices->Add(daughter[1]->GetID(), AliHFEV0pid::kRecoElectron);
+    fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoElectron);
   }
   fGammas->Add(v0);
   
@@ -346,12 +495,12 @@ Bool_t AliHFEV0pid::IsK0s(TObject *v0){
   if(!fIndices->Find(daughter[0]->GetID())){
     AliDebug(1, Form("Adding K0 Pion track with ID %d", daughter[0]->GetID()));
     fPionsK0->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
-    fIndices->Add(daughter[0]->GetID(), AliHFEV0pid::kRecoPionK0);
+    fIndices->Add(daughter[0]->GetID(), AliHFEV0cuts::kRecoPionK0);
   }
   if(!fIndices->Find(daughter[1]->GetID())){
     AliDebug(1, Form("Adding K0 Pion track with ID %d", daughter[1]->GetID()));
     fPionsK0->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
-    fIndices->Add(daughter[1]->GetID(), AliHFEV0pid::kRecoPionK0);
+    fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoPionK0);
   }
   fK0s->Add(v0);
   return kTRUE; 
@@ -390,12 +539,14 @@ Bool_t AliHFEV0pid::IsLambda(TObject *v0){
   Int_t pIndex = 0, nIndex = 0;
   Double_t invMass = 0.;
   Bool_t isLambda = kTRUE; // Lambda - kTRUE, Anti Lambda - kFALSE
+  Double_t mPt = 0.;
   Int_t v0id = -1;
   if(IsESDanalysis()){
     // ESD - cut V0
     AliESDv0 *esdV0 = dynamic_cast<AliESDv0 *>(v0);
     v0id = esdV0->GetLabel();
     if(!fV0cuts->LambdaCuts(esdV0,isLambda)) return kFALSE; 
+    mPt = esdV0->Pt();
     if(fV0cuts->CheckSigns(esdV0)){
       pIndex = esdV0->GetPindex();
       nIndex = esdV0->GetNindex();
@@ -422,24 +573,29 @@ Bool_t AliHFEV0pid::IsLambda(TObject *v0){
   
   // lambda
   if(isLambda){
+    fQA->Fill("h_Pt_L", mPt);
+    fQA->Fill("h_InvMassL", invMass);
+
     if(!fIndices->Find(daughter[0]->GetID())){
       fProtons->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
-      fIndices->Add(daughter[0]->GetID(), AliHFEV0pid::kRecoProton);
+      fIndices->Add(daughter[0]->GetID(), AliHFEV0cuts::kRecoProton);
     }
     if(!fIndices->Find(daughter[1]->GetID())){
       fPionsL->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
-      fIndices->Add(daughter[1]->GetID(), AliHFEV0pid::kRecoPionL);
+      fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoPionL);
     }
   }
   // antilambda
   else{
+    fQA->Fill("h_Pt_AL", mPt);
+    fQA->Fill("h_InvMassAL", invMass);
     if(!fIndices->Find(daughter [1]->GetID())){
       fProtons->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
-      fIndices->Add(daughter[1]->GetID(), AliHFEV0pid::kRecoProton);
+      fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoProton);
     }
     if(!fIndices->Find(daughter [0]->GetID())){
       fPionsL->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
-      fIndices->Add(daughter [0]->GetID(), AliHFEV0pid::kRecoPionL);
+      fIndices->Add(daughter [0]->GetID(), AliHFEV0cuts::kRecoPionL);
     }
   }
   if(isLambda) fLambdas->Add(v0);
@@ -449,6 +605,124 @@ Bool_t AliHFEV0pid::IsLambda(TObject *v0){
 }
 
 //____________________________________________________________
+Int_t AliHFEV0pid::IdentifyV0(TObject *esdV0, Int_t d[2]){
+  //
+  // for MC only, returns the V0 Id
+  //
+
+  //
+  // be carefull about changing the return values - they are used later selectively
+  // In particulra "-2" means that identity of either of the daughters could not be
+  // estimated
+  //
+
+  AliESDv0 *v0 = dynamic_cast<AliESDv0 *>(esdV0);
+  
+  if(!v0) return -1;
+  AliESDtrack* dN, *dP; 
+  Int_t iN, iP;
+  iN = iP = -1;
+  iP = v0->GetPindex();
+  iN = v0->GetNindex();
+  if(iN < 0 || iP < 0) return -1;
+  if(iN >= fNtracks || iP >= fNtracks) return -1;
+  dP = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iP));
+  dN = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iN));  
+  if(!dN || !dP) return -1;
+
+  // as of 26/10/2010
+  // there is still a problem with wrong assignment of positive and negative
+  // V0 daughter in V0 finder - a check is necessary
+  // if the V0 daughters are miss-assigned - swap their labels
+  Bool_t sign = fV0cuts->CheckSigns(v0);
+
+  // get the MC labels
+  Int_t lN, lP;
+  if(sign){
+    lN = dN->GetLabel();
+    lP = dP->GetLabel();
+  }
+  else{
+    lP = dN->GetLabel();
+    lN = dP->GetLabel();
+  }
+  if(lN < 0 || lP < 0) return -2;
+  // get the associated MC particles
+  AliMCParticle *mcP, *mcN;
+  mcP = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lP));
+  mcN = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lN));
+  if(!mcP || !mcN) return -2;
+  
+  // identify the daughter tracks and their mothers
+  Int_t pdgP, pdgN;
+  pdgP = TMath::Abs(mcP->PdgCode());
+  pdgN = TMath::Abs(mcN->PdgCode());
+  // store the daughter ID for later use
+  d[0] = pdgP;
+  d[1] = pdgN;
+  //printf(" -D: pdgP: %i, pdgN: %i\n", pdgP, pdgN);
+  // lablel of the mother particle
+  // -1 may mean it was a primary particle
+  Int_t lPm, lNm;
+  lPm = mcP->GetMother();
+  lNm = mcN->GetMother();
+  if(-1==lPm || -1==lNm) return -3;
+
+  // return if mothers are not the same particle
+  if(lPm != lNm) return -3;
+  // get MC mother particle - now we need only 1
+  AliMCParticle *m = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lPm));
+  if(!m) return -2;
+  // mother PDG
+  Int_t pdgM = m->PdgCode();
+
+  //   if(3122 == TMath::Abs(pdgM)){
+  //     printf("-D: v0 signs: %i\n", fV0cuts->CheckSigns(v0));
+  //     printf("-D: pdgM: %i, pdgN: %i, pdgP: %i \n", pdgM, pdgN, pdgP);
+  //   }
+  
+  // now check the mother and daughters identity
+  if(22 == TMath::Abs(pdgM) && 11 == pdgN && 11 == pdgP) return AliHFEV0cuts::kRecoGamma;
+  if(310 == TMath::Abs(pdgM) && 211 == pdgN && 211 == pdgP) return AliHFEV0cuts::kRecoK0;
+  if(-3122 == pdgM && 2212 == pdgN && 211 == pdgP) return AliHFEV0cuts::kRecoALambda;
+  if(3122 == pdgM && 211 == pdgN && 2212 == pdgP) return AliHFEV0cuts::kRecoLambda;
+    
+  return AliHFEV0cuts::kUndef;
+
+}
+//____________________________________________________________
+void   AliHFEV0pid::ArmenterosPlotMC(AliESDv0 * const v0, Int_t idMC){
+  //
+  // Armenteros plots as a function of Mohter Momentum
+  //
+  //const Float_t minP = 0.1;
+  //const Float_t maxP = 10.;
+  // approx log bins - over the 0.1 - 10 GeV/c
+  const Float_t bins[13] = {0.1, 0.1468, 0.2154, 0.3162, 0.4642, 0.6813, 1.0, 1.4678, 2.1544, 3.1623, 4.6416, 6.8129, 10.0};
+  
+  Float_t ar[2];
+  fV0cuts->Armenteros(v0, ar);
+  Float_t p = v0->P();
+  if( (p <=  bins[0]) || (p >= bins[12])) return;
+
+  Int_t pBin = 0;
+  Float_t tmp = bins[0];
+  while(tmp < p){
+    ++pBin;
+    tmp = bins[pBin];
+  }
+  pBin--;
+
+  if(AliHFEV0cuts::kRecoGamma == idMC) fQA->Fill("h_AP_MC_Gamma_p", pBin, ar[0], ar[1]);
+  if(AliHFEV0cuts::kRecoK0 == idMC) fQA->Fill("h_AP_MC_K0_p", pBin, ar[0], ar[1]);
+  if(AliHFEV0cuts::kRecoLambda == TMath::Abs(idMC)) fQA->Fill("h_AP_MC_Lambda_p", pBin, ar[0], ar[1]);
+  
+  
+
+}
+//____________________________________________________________
 AliHFEV0pid::AliHFEV0pidTrackIndex::AliHFEV0pidTrackIndex():
     fNElectrons(0)
   , fNPionsK0(0)
@@ -518,16 +792,16 @@ void AliHFEV0pid::AliHFEV0pidTrackIndex::Add(Int_t index, Int_t species){
   // Add new index to the list of identified particles
   //
   switch(species){
-    case AliHFEV0pid::kRecoElectron:
+    case AliHFEV0cuts::kRecoElectron:
       fIndexElectron[fNElectrons++] = index;
       break;
-    case AliHFEV0pid::kRecoPionK0:
+    case AliHFEV0cuts::kRecoPionK0:
       fIndexPionK0[fNPionsK0++] = index;
       break;
-    case AliHFEV0pid::kRecoPionL:
+    case AliHFEV0cuts::kRecoPionL:
       fIndexPionL[fNPionsL++] = index;
       break;
-    case AliHFEV0pid::kRecoProton:
+    case AliHFEV0cuts::kRecoProton:
       fIndexProton[fNProtons++] = index;
       break;
   };
@@ -541,19 +815,19 @@ Bool_t AliHFEV0pid::AliHFEV0pidTrackIndex::Find(Int_t index, Int_t species) cons
 
   Int_t *container = NULL; Int_t n = 0;
   switch(species){
-  case AliHFEV0pid::kRecoElectron:
+  case AliHFEV0cuts::kRecoElectron:
     container = fIndexElectron;
     n = fNElectrons;
     break;
-  case AliHFEV0pid::kRecoPionK0:
+  case AliHFEV0cuts::kRecoPionK0:
     container = fIndexPionK0;
     n = fNPionsK0;
     break;
-  case AliHFEV0pid::kRecoPionL:
+  case AliHFEV0cuts::kRecoPionL:
     container = fIndexPionL;
     n = fNPionsL;
     break;
-  case AliHFEV0pid::kRecoProton:
+  case AliHFEV0cuts::kRecoProton:
     container = fIndexProton;
     n = fNProtons;
     break;
@@ -575,10 +849,10 @@ Bool_t AliHFEV0pid::AliHFEV0pidTrackIndex::Find(Int_t index) const {
   // 
   // Find index in all samples
   //
-  if(Find(index, AliHFEV0pid::kRecoElectron)) return kTRUE;
-  else if(Find(index, AliHFEV0pid::kRecoPionK0)) return kTRUE;
-  else if(Find(index, AliHFEV0pid::kRecoPionL)) return kTRUE;
-  else return Find(index, AliHFEV0pid::kRecoProton);
+  if(Find(index, AliHFEV0cuts::kRecoElectron)) return kTRUE;
+  else if(Find(index, AliHFEV0cuts::kRecoPionK0)) return kTRUE;
+  else if(Find(index, AliHFEV0cuts::kRecoPionL)) return kTRUE;
+  else return Find(index, AliHFEV0cuts::kRecoProton);
 }
 
 //____________________________________________________________
@@ -586,7 +860,7 @@ TList *AliHFEV0pid::GetListOfQAhistograms(){
   //
   // Getter for V0 PID QA histograms
   //
-  
+
   TList *tmp = fV0cuts->GetList();
   tmp->SetName("V0cuts");
   fOutput->Add(tmp);
@@ -596,5 +870,11 @@ TList *AliHFEV0pid::GetListOfQAhistograms(){
     tmp->SetName("V0pid");
     fOutput->Add(tmp);
   } 
+  tmp = 0x0;
+  tmp = fV0cuts->GetListMC();
+  tmp->SetName("V0cutsMC");
+  //printf(" -D: adding MC V0 cuts stuff\n");
+  fOutput->Add(tmp);
+  
   return fOutput;
 }