Update of the HFE package
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliHFEsignalCuts.cxx
1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 *                                                                        *
4 * Author: The ALICE Off-line Project.                                    *
5 * Contributors are mentioned in the code where appropriate.              *
6 *                                                                        *
7 * Permission to use, copy, modify and distribute this software and its   *
8 * documentation strictly for non-commercial purposes is hereby granted   *
9 * without fee, provided that the above copyright notice appears in all   *
10 * copies and that both the copyright notice and this permission notice   *
11 * appear in the supporting documentation. The authors make no claims     *
12 * about the suitability of this software for any purpose. It is          *
13 * provided "as is" without express or implied warranty.                  *
14 **************************************************************************/
15 //
16 // Signal cuts
17 // Checks whether a particle (reconstructed or MC) is coming from MC Signal
18 // For more information see implementation file
19 //
20 // Autor:
21 //   Markus Fasel <M.Fasel@gsi.de>
22 //
23 #include <TClass.h>
24 #include <TMath.h>
25 #include <TParticle.h>
26 #include <TString.h>
27
28 #include "AliAODTrack.h"
29 #include "AliAODMCParticle.h"
30 #include "AliESDtrack.h"
31 #include "AliLog.h"
32 #include "AliMCEvent.h"
33 #include "AliMCParticle.h"
34 #include "AliVParticle.h"
35
36 #include "AliHFEsignalCuts.h"
37 #include "AliHFEmcQA.h"
38
39 ClassImp(AliHFEsignalCuts)
40
41 //____________________________________________________________
42 AliHFEsignalCuts::AliHFEsignalCuts():
43   AliAnalysisCuts(),
44   fMC(NULL),
45   fMCQA(NULL)
46 {
47   //
48   // Dummy constructor
49   //
50 }
51
52 //____________________________________________________________
53 AliHFEsignalCuts::AliHFEsignalCuts(const Char_t *name, const Char_t *title):
54   AliAnalysisCuts(name, title),
55   fMC(NULL),
56   fMCQA(NULL)
57 {
58   //
59   // Default constructor
60   //
61   fMCQA = new AliHFEmcQA;
62   if(fMCQA) fMCQA->Init();
63 }
64
65 //____________________________________________________________
66 AliHFEsignalCuts::AliHFEsignalCuts(const AliHFEsignalCuts &ref):
67   AliAnalysisCuts(ref),
68   fMC(ref.fMC),
69   fMCQA(ref.fMCQA)
70 {
71   //
72   // Copy constructor
73   //
74 }
75
76 //____________________________________________________________
77 AliHFEsignalCuts &AliHFEsignalCuts::operator=(const AliHFEsignalCuts &ref){
78   //
79   // Assignment operator
80   //
81   if(this != &ref){
82     fMC = ref.fMC; 
83     fMCQA = ref.fMCQA; 
84   }
85   return *this;
86 }
87
88 //____________________________________________________________
89 AliHFEsignalCuts::~AliHFEsignalCuts(){
90   //
91   // Destructor
92   //
93   if(fMCQA) delete fMCQA;
94 }
95
96 //____________________________________________________________
97 void AliHFEsignalCuts::SetMCEvent(AliMCEvent *mc){ 
98   //
99   // Set mc event
100   //
101   fMC = mc; 
102   if(fMCQA) fMCQA->SetMCEvent(mc);
103 }
104
105 //____________________________________________________________
106 Bool_t AliHFEsignalCuts::IsSelected(TObject *o){
107   //
108   // Define signal as electron coming from charm or beauty
109   // @TODO: Implement setter so that also either of them can be defined
110   // as signal alone
111   
112
113   return IsCharmElectron(o) || IsBeautyElectron(o);
114 /*  
115   //saving time?
116   Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
117   if(esources>0)printf("esources= %d\n",esources);
118   if(esources == AliHFEmcQA::kDirectCharm || esources == AliHFEmcQA::kDirectBeauty || esources == AliHFEmcQA::kBeautyCharm)  // 1: direct D->e, 2: B->e 3: B->D->e
119     return kTRUE;
120   else
121     return kFALSE;
122 */
123
124 }
125
126 //____________________________________________________________
127 Bool_t AliHFEsignalCuts::IsCharmElectron(const TObject * const o) const {
128   //
129   // Check if mother is coming from Charm
130   //
131   if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
132   Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
133   if(esources == AliHFEmcQA::kDirectCharm)  // 1: direct D->e
134     return kTRUE;
135   else
136     return kFALSE;
137 }
138
139 //____________________________________________________________
140 Bool_t AliHFEsignalCuts::IsBeautyElectron(const TObject * const o) const {
141   //
142   // Check if mother is coming from Beauty
143   //
144   if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
145   Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
146   if(esources == AliHFEmcQA::kDirectBeauty || esources == AliHFEmcQA::kBeautyCharm)  // 2: B->e 3: B->D->e
147     return kTRUE;
148   else
149     return kFALSE;
150 }
151
152 //____________________________________________________________
153 Bool_t AliHFEsignalCuts::IsGammaElectron(const TObject * const o) const {
154   //
155   // Check for MC if the electron is coming from Gamma  
156   //
157   if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
158   Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
159   if(esources >= AliHFEmcQA::kGammaPi0 && esources <= AliHFEmcQA::kGammaRho0 )  // 4: conversion electrons
160   //if(esources == AliHFEmcQA::kGammaPi0 || esources == AliHFEmcQA::kGammaEta || esources == AliHFEmcQA::kGammaOmega || esources == AliHFEmcQA::kGammaPhi || esources == AliHFEmcQA::kGammaEtaPrime || esources == AliHFEmcQA::kGammaRho0 )  // 4: conversion electrons
161     return kTRUE;
162   else
163     return kFALSE;
164 }
165
166 //____________________________________________________________
167 Bool_t AliHFEsignalCuts::IsNonHFElectron(const TObject * const o) const {
168   //
169   // Check for MC if the electron is coming from NonHFE except for conversion
170   //
171   if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
172   Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));  
173   if(esources == AliHFEmcQA:: kPi0 || esources == AliHFEmcQA::kEta || esources == AliHFEmcQA::kOmega || esources == AliHFEmcQA::kPhi || esources == AliHFEmcQA::kEtaPrime || esources == AliHFEmcQA::kRho0)  // 4: conversion electrons
174     return kTRUE;
175   else
176     return kFALSE;
177 }
178 //____________________________________________________________
179
180 /*
181 //____________________________________________________________
182 Bool_t AliHFEsignalCuts::IsCharmElectron(const TObject * const o) const {
183   //
184   // Check if mother is coming from Charm
185   //
186   if(TMath::Abs(GetTrackPDG(dynamic_cast<const AliVParticle *>(o))) != 11) return kFALSE;
187   Int_t motherpdg = TMath::Abs(GetMotherPDG(dynamic_cast<const AliVParticle *>(o)));
188   AliDebug(1, Form("Mother PDG %d\n", motherpdg));
189
190   if((motherpdg % 1000) / 100 == 4) return kTRUE;    // charmed meson, 3rd position in pdg code == 4
191   if(motherpdg / 1000 == 4) return kTRUE;            // charmed baryon, 4th position in pdg code == 4
192   AliDebug(1, "No Charm\n");
193   return kFALSE;
194 }
195
196 //____________________________________________________________
197 Bool_t AliHFEsignalCuts::IsBeautyElectron(const TObject * const o) const {
198   //
199   // Check if mother is coming from Beauty
200   //
201   if(TMath::Abs(GetTrackPDG(dynamic_cast<const AliVParticle *>(o))) != 11) return kFALSE;
202   Int_t motherpdg = TMath::Abs(GetMotherPDG(dynamic_cast<const AliVParticle *>(o)));
203   AliDebug(1, Form("Mother PDG %d\n", motherpdg));
204
205   if((motherpdg % 1000) / 100 == 5) return kTRUE;   // beauty meson, 3rd position in pdg code == 5
206   if(motherpdg / 1000 == 5) return kTRUE;           // beauty baryon, 4th position in pdg code == 5   
207   AliDebug(1, "No Beauty\n");
208   return kFALSE;
209 }
210
211 //____________________________________________________________
212 Bool_t AliHFEsignalCuts::IsGammaElectron(const TObject * const o) const {
213   //
214   // Check for MC if the electron is coming from Gamma
215   //
216   if(TMath::Abs(GetTrackPDG(dynamic_cast<const AliVParticle *>(o))) != 11) return kFALSE;
217   Int_t motherpdg = TMath::Abs(GetMotherPDG(dynamic_cast<const AliVParticle *>(o)));
218   AliDebug(1, Form("Mother PDG %d\n", motherpdg));
219
220   if(motherpdg!=22){
221     AliDebug(1, "No Gamma");
222     return kFALSE;
223   } else { 
224     AliDebug(1, "Gamma");
225     return kTRUE;
226   }
227 }
228 */
229
230 //____________________________________________________________
231 Int_t AliHFEsignalCuts::GetMotherPDG(const AliVParticle * const track) const {
232   //
233   // Get Mother Pdg code for reconstructed respectively MC tracks
234   // 
235   if(!fMC){
236     AliDebug(1, "No MC Event Available\n");
237     return 0;
238   }
239   const AliVParticle *motherParticle = NULL, *mctrack = NULL;
240   TString objectType = track->IsA()->GetName();
241   if(objectType.CompareTo("AliESDtrack") == 0 || objectType.CompareTo("AliAODTrack") == 0){
242     // Reconstructed track
243     if(track->GetLabel())
244       mctrack = fMC->GetTrack(TMath::Abs(track->GetLabel()));
245   } else {
246     // MCParticle
247     mctrack = track;
248   }
249
250   if(!mctrack) return 0;
251   
252   Int_t motherPDG = 0;
253   if(TString(mctrack->IsA()->GetName()).CompareTo("AliMCParticle") == 0){
254     // case MC Particle
255     const AliMCParticle *esdmctrack = dynamic_cast<const AliMCParticle *>(mctrack);
256     if(esdmctrack) motherParticle = fMC->GetTrack(esdmctrack->Particle()->GetFirstMother());
257     if(motherParticle){
258       const AliMCParticle *esdmcmother = dynamic_cast<const AliMCParticle *>(motherParticle);
259       if(esdmcmother) motherPDG = TMath::Abs(esdmcmother->Particle()->GetPdgCode());
260     }
261   } else {
262     // case AODMCParticle
263     const AliAODMCParticle *aodmctrack = dynamic_cast<const AliAODMCParticle *>(mctrack);
264     if(aodmctrack) motherParticle = fMC->GetTrack(aodmctrack->GetMother());
265     if(motherParticle){
266       const AliAODMCParticle *aodmcmother = dynamic_cast<const AliAODMCParticle *>(motherParticle);
267       if(aodmcmother) motherPDG = TMath::Abs(aodmcmother->GetPdgCode());
268     }
269   }
270   return motherPDG;
271 }
272
273 //____________________________________________________________
274 Int_t AliHFEsignalCuts::GetTrackPDG(const AliVParticle * const track) const {
275         //
276         // Return PDG code of a particle itself
277         //
278   if(!fMC){
279     AliDebug(1, "No MC Event Available\n");
280     return 0;
281   }
282         TString sourcetype = track->IsA()->GetName();
283         const AliVParticle *mctrack = NULL;
284         if(!sourcetype.CompareTo("AliESDtrack") || !sourcetype.CompareTo("AliAODTrack")){
285                 mctrack = fMC->GetTrack(TMath::Abs(track->GetLabel()));
286         } else  mctrack = track;
287         if(!mctrack) return 0;
288
289         TString mctype = mctrack->IsA()->GetName();
290         Int_t trackPdg = 0;
291         if(!mctype.CompareTo("AliMCParticle")){
292                 const AliMCParticle *esdmc = dynamic_cast<const AliMCParticle *>(mctrack);
293                 if(esdmc) trackPdg = esdmc->Particle()->GetPdgCode();
294         } else {
295                 const AliAODMCParticle *aodmc = dynamic_cast< const AliAODMCParticle *>(mctrack);
296                 if(aodmc) trackPdg = aodmc->GetPdgCode();
297         }
298         return trackPdg;
299 }
300
301 //____________________________________________________________
302 Int_t AliHFEsignalCuts::GetElecSource(const AliVParticle * const track) const {
303         //
304         // Return PDG code of a particle itself
305         //
306         
307   if(!fMC){
308     AliDebug(1, "No MC Event Available\n");
309     return 0;
310   }
311   if(!fMCQA){
312     AliDebug(1, "No MCQA Available\n");
313     return 0;
314   }
315   if(!track){
316     AliDebug(1, "Track not Available\n");
317     return 0;
318   }
319
320   TClass *tracktype;
321   const AliVParticle *mctrack = NULL;
322   TParticle *mcpart = NULL;
323   if((tracktype = track->IsA()) == AliESDtrack::Class() || tracktype == AliAODTrack::Class()){
324     mctrack = fMC->GetTrack(TMath::Abs(track->GetLabel()));
325   } else  mctrack = track;
326   if(!mctrack) return 0;
327
328   Int_t eSource = 0;
329   if(mctrack->IsA() == AliMCParticle::Class()){
330     const AliMCParticle *esdmc = dynamic_cast<const AliMCParticle *>(mctrack);
331     if(esdmc){
332       mcpart = esdmc->Particle();
333       eSource=fMCQA->GetElecSource(mcpart);
334     }
335   } else {
336     return -1;
337   }
338   return eSource;
339 }