From bfd68286113d591e2f0ea4c27985808a57a9c52f Mon Sep 17 00:00:00 2001
From: morsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Date: Wed, 11 Mar 2009 10:58:57 +0000
Subject: [PATCH] A. Dainese added a method to match the reconstructed
 particles to the MC (returns the label of the corresponding
 AliAODMCParticle). He prefers to have this as a method of AliAODRecoDecay,
 because it is not completely trivial and the same piece of code was getting
 duplicated in several different places.

---
 STEER/AliAODRecoDecay.cxx | 97 +++++++++++++++++++++++++++++++++++++++
 STEER/AliAODRecoDecay.h   |  9 ++++
 2 files changed, 106 insertions(+)

diff --git a/STEER/AliAODRecoDecay.cxx b/STEER/AliAODRecoDecay.cxx
index fb93e84c389..372c3468e91 100644
--- a/STEER/AliAODRecoDecay.cxx
+++ b/STEER/AliAODRecoDecay.cxx
@@ -22,8 +22,10 @@
 
 #include <TDatabasePDG.h>
 #include <TVector3.h>
+#include <TClonesArray.h>
 #include "AliLog.h"
 #include "AliVTrack.h"
+#include "AliAODMCParticle.h"
 #include "AliAODRecoDecay.h"
 
 ClassImp(AliAODRecoDecay)
@@ -472,6 +474,101 @@ Double_t AliAODRecoDecay::InvMass2Prongs(Int_t ip1,Int_t ip2,
 
   return mass;
 }
+//----------------------------------------------------------------------------
+Int_t AliAODRecoDecay::MatchToMC(Int_t pdgabs,TClonesArray *mcArray) const
+{
+  //
+  // Check if this candidate is matched to a MC signal
+  // If no, return -1
+  // If yes, return label (>=0) of the AliAODMCParticle
+  // 
+
+  if(!GetNDaughters()) {
+    AliError("No daughters available");
+    return -1;
+  }
+  
+  Int_t *dgLabels = new Int_t[GetNDaughters()];
+
+  // loop on daughters and write the labels
+  for(Int_t i=0; i<GetNDaughters(); i++) {
+    AliVTrack *trk = (AliVTrack*)GetDaughter(i);
+    dgLabels[i] = trk->GetLabel();
+  }
+
+  Int_t labMother = MatchToMC(pdgabs,mcArray,dgLabels);
+
+  delete [] dgLabels; dgLabels=NULL;
+
+  return labMother;
+}
+//----------------------------------------------------------------------------
+Int_t AliAODRecoDecay::MatchToMC(Int_t pdgabs,TClonesArray *mcArray,
+				 Int_t *dgLabels) const
+{
+  //
+  // Check if this candidate is matched to a MC signal
+  // If no, return -1
+  // If yes, return label (>=0) of the AliAODMCParticle
+  // 
+
+  if(!GetNDaughters()) {
+    AliError("No daughters available");
+    return -1;
+  }
+  
+  Int_t *labMom = new Int_t[GetNDaughters()];
+  Int_t i,lab,labMother,pdgMother;
+  AliAODMCParticle *part=0;
+
+  // loop on daughter labels
+  for(i=0; i<GetNDaughters(); i++) {
+    labMom[i]=-1;
+    lab = dgLabels[i];
+    if(lab<0) {
+      printf("daughter with negative label\n");
+      continue;
+    }
+    part = (AliAODMCParticle*)mcArray->At(lab);
+    if(!part) { 
+      printf("no MC particle\n");
+      continue;
+    }
+    while(part->GetMother()>=0) {
+      labMother=part->GetMother();
+      part = (AliAODMCParticle*)mcArray->At(labMother);
+      if(!part) {
+	printf("no MC mother particle\n");
+	break;
+      }
+      pdgMother = TMath::Abs(part->GetPdgCode());
+      if(pdgMother==pdgabs) {
+	labMom[i]=labMother;
+	break;
+      }
+    }
+  }
+
+  // check if the candidate is signal
+  Bool_t isSignal=kTRUE;
+  labMother=labMom[0];
+  // all labels have to be the same and !=-1
+  for(i=0; i<GetNDaughters(); i++) {
+    if(labMom[i]==-1 || labMom[i]!=labMother) isSignal=kFALSE;
+  }
+
+  delete [] labMom; labMom=NULL;
+
+  if(!isSignal) return -1;
+
+  // check that the mother decayed in <GetNDaughters()> prongs
+  part = (AliAODMCParticle*)mcArray->At(labMother);
+  Int_t ndg = TMath::Abs(part->GetDaughter(1)-part->GetDaughter(0))+1;
+
+  if(ndg!=GetNDaughters()) return -1;
+
+  return labMother;
+}
 //---------------------------------------------------------------------------
 void AliAODRecoDecay::Print(Option_t* /*option*/) const 
 {
diff --git a/STEER/AliAODRecoDecay.h b/STEER/AliAODRecoDecay.h
index 1554b95f328..a711d0ff34b 100644
--- a/STEER/AliAODRecoDecay.h
+++ b/STEER/AliAODRecoDecay.h
@@ -11,6 +11,7 @@
 
 #include <TMath.h>
 #include <TRef.h>
+#include <TClonesArray.h>
 #include "AliAODVertex.h"
 #include "AliAODTrack.h"
 #include "AliVTrack.h"
@@ -46,6 +47,13 @@ class AliAODRecoDecay : public AliVTrack {
   Short_t  GetCharge() const {return fCharge;}
   void     SetCharge(Short_t charge=0) {fCharge=charge;}
 
+  // Match to MC signal:
+  // check if this candidate is matched to a MC signal
+  // If no, return -1
+  // If yes, return label (>=0) of the AliAODMCParticle
+  Int_t    MatchToMC(Int_t pdgabs,TClonesArray *mcArray) const;
+  Int_t    MatchToMC(Int_t pdgabs,TClonesArray *mcArray,Int_t *dgLabels) const;
+
   // PID
   void      SetPID(Int_t nprongs,Double_t *pid);
   Double_t *GetPID() const { return fPID; }
@@ -348,3 +356,4 @@ inline void AliAODRecoDecay::GetPIDProng(Int_t ip,Double_t *pid) const
 
 
 #endif
+
-- 
2.39.3