]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliHFENonPhotonicElectron.cxx
Update of npe Pb-Pb code
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliHFENonPhotonicElectron.cxx
CommitLineData
76d0b522 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 *
afb48e1d 9
76d0b522 10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 * *
16 *************************************************************************************/
17
18/*************************************************************************************
19 * *
20 * Class for the Selection of Non-Heavy-Flavour-Electrons trought *
21 * the invariant mass method. The selection can be done from two *
22 * different algorithms, which can be choosed calling the function *
23 * "SetAlgorithm(TString Algorithm)". *
24 * *
25 * Authors: R.Bailhache, C.A.Schmidt *
26 * *
27 *************************************************************************************/
28
29#include "TVector2.h"
30#include "THnSparse.h"
31#include "TMath.h"
32#include "TLorentzVector.h"
33#include "TParticle.h"
34#include "TList.h"
35#include "TDatabasePDG.h"
36
37#include "AliVEvent.h"
38#include "AliMCEvent.h"
39#include "AliESDEvent.h"
40#include "AliMCParticle.h"
41#include "AliAODMCParticle.h"
42#include "AliAODEvent.h"
43#include "AliAODVertex.h"
44#include "AliAODTrack.h"
45#include "AliVTrack.h"
46#include "AliESDtrack.h"
47#include "AliESDtrackCuts.h"
48#include "AliPIDResponse.h"
8a9b2231 49#include "AliPID.h"
465ff082 50#include "AliStack.h"
76d0b522 51
52#include "AliKFParticle.h"
53#include "AliKFVertex.h"
54
55#include "AliHFEcuts.h"
56#include "AliHFEpid.h"
57#include "AliHFEpidQAmanager.h"
58#include "AliHFEtools.h"
afb48e1d 59#include "AliHFEmcQA.h"
76d0b522 60
61#include "AliHFENonPhotonicElectron.h"
62
63ClassImp(AliHFENonPhotonicElectron)
465ff082 64 //________________________________________________________________________
65 AliHFENonPhotonicElectron::AliHFENonPhotonicElectron(const char *name, const Char_t *title)
66 :TNamed (name, title)
67 ,fIsAOD (kFALSE)
68 ,fMCEvent (NULL)
69 ,fAODArrayMCInfo (NULL)
70 ,fLevelBack(-1)
71 ,fHFEBackgroundCuts (NULL)
72 ,fPIDBackground (0x0)
73 ,fPIDBackgroundQA (0)
74 ,fkPIDRespons (NULL)
75 ,fPtBinning()
76 ,fEtaBinning()
77 ,fInvMassBinning()
78 ,fStudyRadius (kFALSE)
79 ,fAlgorithmMA (kTRUE)
80 ,fChi2OverNDFCut (3.0)
81 ,fMaxDCA (3.0)
82 // ,fMaxOpeningTheta (0.02)
83 // ,fMaxOpeningPhi (0.1)
84 ,fMaxOpening3D (TMath::Pi())
85 ,fMaxInvMass (1000)
86 ,fSetMassConstraint (kFALSE)
87 ,fSelectCategory1tracks(kTRUE)
88 ,fSelectCategory2tracks(kFALSE)
89 ,fITSmeanShift(0.)
90 ,fITSnSigmaHigh(3.)
91 ,fITSnSigmaLow(-3.)
92 ,fminPt(0.1)
93 ,fEtaDalitzWeightFactor(1.0)
94 ,fArraytrack (NULL)
95 ,fCounterPoolBackground (0)
96 ,fnumberfound (0)
97 ,fListOutput (NULL)
98 ,fAssElectron (NULL)
99 ,fIncElectron (NULL)
100 ,fUSign (NULL)
101 ,fLSign (NULL)
102 ,fUSmatches(NULL)
103 ,fLSmatches(NULL)
104 ,fHnsigmaITS(NULL)
105 ,fWeightsSource(NULL)
106 ,fIncElectronRadius(NULL)
107 ,fRecElectronRadius(NULL)
108 // ,fUSignAngle (NULL)
109 // ,fLSignAngle (NULL)
110 ,fAnaPairGen(kFALSE)
111 ,fNumberofGenerations(1)
112 ,fDisplayMCStack(kFALSE)
76d0b522 113{
465ff082 114 //
115 // Constructor
116 //
117 fPIDBackground = new AliHFEpid("hfePidBackground");
118 fPIDBackgroundQA = new AliHFEpidQAmanager;
76d0b522 119}
120
121//________________________________________________________________________
465ff082 122 AliHFENonPhotonicElectron::AliHFENonPhotonicElectron()
123 :TNamed ()
124 ,fIsAOD (kFALSE)
125 ,fMCEvent (NULL)
126 ,fAODArrayMCInfo (NULL)
127 ,fLevelBack(-1)
128 ,fHFEBackgroundCuts (NULL)
129 ,fPIDBackground (0x0)
130 ,fPIDBackgroundQA (0)
131 ,fkPIDRespons (NULL)
132 ,fPtBinning()
133 ,fEtaBinning()
134 ,fInvMassBinning()
135 ,fStudyRadius (kFALSE)
136 ,fAlgorithmMA (kTRUE)
137 ,fChi2OverNDFCut (3.0)
138 ,fMaxDCA (3.0)
139 // ,fMaxOpeningTheta (0.02)
140 // ,fMaxOpeningPhi (0.1)
141 ,fMaxOpening3D (TMath::TwoPi())
142 ,fMaxInvMass (1000)
143 ,fSetMassConstraint (kFALSE)
144 ,fSelectCategory1tracks(kTRUE)
145 ,fSelectCategory2tracks(kFALSE)
146 ,fITSmeanShift(0.)
147 ,fITSnSigmaHigh(3.)
148 ,fITSnSigmaLow(-3.)
149 ,fminPt(0.1)
150 ,fEtaDalitzWeightFactor(1.0)
151 ,fArraytrack (NULL)
152 ,fCounterPoolBackground (0)
153 ,fnumberfound (0)
154 ,fListOutput (NULL)
155 ,fAssElectron (NULL)
156 ,fIncElectron (NULL)
157 ,fUSign (NULL)
158 ,fLSign (NULL)
159 ,fUSmatches(NULL)
160 ,fLSmatches(NULL)
161 ,fHnsigmaITS(NULL)
162 ,fWeightsSource(NULL)
163 ,fIncElectronRadius(NULL)
164 ,fRecElectronRadius(NULL)
165 // ,fUSignAngle (NULL)
166 // ,fLSignAngle (NULL)
167 ,fAnaPairGen(kFALSE)
168 ,fNumberofGenerations(1)
169 ,fDisplayMCStack(kFALSE)
76d0b522 170{
465ff082 171 //
172 // Constructor
173 //
174 fPIDBackground = new AliHFEpid("hfePidBackground");
175 fPIDBackgroundQA = new AliHFEpidQAmanager;
76d0b522 176}
177
178//________________________________________________________________________
465ff082 179 AliHFENonPhotonicElectron::AliHFENonPhotonicElectron(const AliHFENonPhotonicElectron &ref)
180 :TNamed(ref)
181 ,fIsAOD (ref.fIsAOD)
182 ,fMCEvent (NULL)
183 ,fAODArrayMCInfo (NULL)
184 ,fLevelBack (ref.fLevelBack)
185 ,fHFEBackgroundCuts (ref.fHFEBackgroundCuts)
186 ,fPIDBackground (ref.fPIDBackground)
187 ,fPIDBackgroundQA (ref.fPIDBackgroundQA)
188 ,fkPIDRespons (ref.fkPIDRespons)
189 ,fPtBinning (ref.fPtBinning)
190 ,fEtaBinning (ref.fEtaBinning)
191 ,fInvMassBinning (ref.fInvMassBinning)
192 ,fStudyRadius (ref.fStudyRadius)
193 ,fAlgorithmMA (ref.fAlgorithmMA)
194 ,fChi2OverNDFCut (ref.fChi2OverNDFCut)
195 ,fMaxDCA (ref.fMaxDCA)
196 // ,fMaxOpeningTheta (ref.fMaxOpeningTheta)
197 // ,fMaxOpeningPhi (ref.fMaxOpeningPhi)
198 ,fMaxOpening3D (ref.fMaxOpening3D)
199 ,fMaxInvMass (ref.fMaxInvMass)
200 ,fSetMassConstraint (ref.fSetMassConstraint)
201 ,fSelectCategory1tracks(ref.fSelectCategory1tracks)
202 ,fSelectCategory2tracks(ref.fSelectCategory2tracks)
203 ,fITSmeanShift(ref.fITSmeanShift)
204 ,fITSnSigmaHigh(ref.fITSnSigmaHigh)
205 ,fITSnSigmaLow(ref.fITSnSigmaLow)
206 ,fminPt(ref.fminPt)
207 ,fEtaDalitzWeightFactor(ref.fEtaDalitzWeightFactor)
208 ,fArraytrack (NULL)
209 ,fCounterPoolBackground (0)
210 ,fnumberfound (0)
211 ,fListOutput (ref.fListOutput)
212 ,fAssElectron (ref.fAssElectron)
213 ,fIncElectron (ref.fIncElectron)
214 ,fUSign (ref.fUSign)
215 ,fLSign (ref.fLSign)
216 ,fUSmatches(ref.fUSmatches)
217 ,fLSmatches(ref.fLSmatches)
218 ,fHnsigmaITS(ref.fHnsigmaITS)
219 ,fWeightsSource(ref.fWeightsSource)
220 ,fIncElectronRadius(ref.fIncElectronRadius)
221 ,fRecElectronRadius(ref.fRecElectronRadius)
222 // ,fUSignAngle (ref.fUSignAngle)
223 // ,fLSignAngle (ref.fLSignAngle)
224 ,fAnaPairGen(kFALSE)
225 ,fNumberofGenerations(1)
226 ,fDisplayMCStack(kFALSE)
76d0b522 227{
465ff082 228 //
229 // Copy Constructor
230 //
231 ref.Copy(*this);
76d0b522 232}
233
234//____________________________________________________________
235AliHFENonPhotonicElectron &AliHFENonPhotonicElectron::operator=(const AliHFENonPhotonicElectron &ref){
465ff082 236 //
237 // Assignment operator
238 //
239 if(this == &ref) ref.Copy(*this);
240 return *this;
76d0b522 241}
242
243//_________________________________________
244AliHFENonPhotonicElectron::~AliHFENonPhotonicElectron()
245{
465ff082 246 //
247 // Destructor
248 //
249 if(fArraytrack) delete fArraytrack;
250 //if(fHFEBackgroundCuts) delete fHFEBackgroundCuts;
251 if(fPIDBackground) delete fPIDBackground;
252 if(fPIDBackgroundQA) delete fPIDBackgroundQA;
76d0b522 253}
254
255//_____________________________________________________________________________________________
256void AliHFENonPhotonicElectron::Init()
257{
465ff082 258 //
259 // Init
260 //
261
262 //printf("Analysis Mode for AliHFENonPhotonicElectron: %s Analysis\n", fIsAOD ? "AOD" : "ESD");
263
264 if(!fListOutput) fListOutput = new TList;
265 fListOutput->SetName("HFENonPhotonicElectron");
266 fListOutput->SetOwner();
267
268 if(!fHFEBackgroundCuts) fHFEBackgroundCuts = new AliHFEcuts();
269 if(fIsAOD) fHFEBackgroundCuts->SetAOD();
270 fHFEBackgroundCuts->Initialize();
271 if(fHFEBackgroundCuts->IsQAOn()) {
272 fListOutput->Add(fHFEBackgroundCuts->GetQAhistograms());
273 }
274
275 // Initialize PID
276 if(!fPIDBackground) fPIDBackground = new AliHFEpid("default pid");
277 if(fMCEvent || fAODArrayMCInfo) fPIDBackground->SetHasMCData(kTRUE); // does nothing since the fMCEvent are set afterwards at the moment
278 if(!fPIDBackground->GetNumberOfPIDdetectors())
279 {
280 //fPIDBackground->AddDetector("TOF", 0);
281 fPIDBackground->AddDetector("TPC", 0);
282 }
283 AliInfo("PID Background QA switched on");
284 fPIDBackgroundQA->Initialize(fPIDBackground);
285 fListOutput->Add(fPIDBackgroundQA->MakeList("HFENP_PID_Background"));
286 fPIDBackground->SortDetectors();
287
288 const Int_t kBinsPtDefault = 35;
289 Double_t binLimPtDefault[kBinsPtDefault+1] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2, 1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 4.5, 5., 5.5, 6., 7., 8., 10., 12., 14., 16., 18., 20.};
290 const Int_t kBinsEtaInclusiveDefault = 8;
291 Double_t binLimEtaInclusiveDefault[kBinsEtaInclusiveDefault+1] = {-0.8, -0.6, -0.4, -0.2, 0., 0.2, 0.4, 0.6, 0.8};
292 const Int_t kBinsEtaAssociated = 15;
293 Double_t binLimEtaAssociat[kBinsEtaAssociated+1] = {-1.5,-1.3,-1.1,-0.9,-0.7,-0.5,-0.3,-0.1,0.1,0.3,0.5,0.7,0.9,1.1,1.3,1.5};
294
295 if(!fPtBinning.GetSize()) fPtBinning.Set(kBinsPtDefault+1, binLimPtDefault);
296 if(!fEtaBinning.GetSize()) fEtaBinning.Set(kBinsEtaInclusiveDefault+1, binLimEtaInclusiveDefault);
297
298 //Int_t nBinsP = 400;
299 //Double_t minP = 0.0;
300 //Double_t maxP = 20.0;
301 //Double_t binLimP[nBinsP+1];
302 //for(Int_t i=0; i<=nBinsP; i++) binLimP[i]=(Double_t)minP + (maxP-minP)/nBinsP*(Double_t)i ;
303
304 Int_t nBinsC = 11;
305 Double_t minC = 0.0;
306 Double_t maxC = 11.0;
307 Double_t binLimC[nBinsC+1];
308 for(Int_t i=0; i<=nBinsC; i++) binLimC[i]=(Double_t)minC + (maxC-minC)/nBinsC*(Double_t)i ;
309
310 Int_t nBinsSource = 12;
311 Double_t minSource = 0.;
312 Double_t maxSource = 12.;
313 Double_t binLimSource[nBinsSource+1];
314 for(Int_t i=0; i<=nBinsSource; i++) binLimSource[i]=(Double_t)minSource + (maxSource-minSource)/nBinsSource*(Double_t)i ;
315
316 Int_t nBinsInvMass = 30;
317 Double_t minInvMass = 0.;
318 Double_t maxInvMass = 0.3;
319 Double_t binLimInvMass[nBinsInvMass+1];
320 for(Int_t i=0; i<=nBinsInvMass; i++) binLimInvMass[i]=(Double_t)minInvMass + (maxInvMass-minInvMass)/nBinsInvMass*(Double_t)i ;
321 if(!fInvMassBinning.GetSize()) fInvMassBinning.Set(nBinsInvMass+1, binLimInvMass);
322
323 Int_t nBinsPhi = 8;
324 Double_t minPhi = 0.0;
325 Double_t maxPhi = TMath::Pi();
326 Double_t binLimPhi[nBinsPhi+1];
327 for(Int_t i=0; i<=nBinsPhi; i++)
328 {
329 binLimPhi[i]=(Double_t)minPhi + (maxPhi-minPhi)/nBinsPhi*(Double_t)i ;
330 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
331 }
332
333 //now used for pair generation
334 Int_t nBinsGen = fNumberofGenerations;
335 Double_t minGen = 0;
336 Double_t maxGen = fNumberofGenerations;
337 Double_t binLimGen[nBinsGen+1];
338 for(Int_t i=0; i<=nBinsGen; i++)
339 {
340 binLimGen[i]=(Double_t)minGen + (maxGen-minGen)/nBinsGen*(Double_t)i ;
341 AliDebug(2,Form("bin Generation is %f for %d",binLimGen[i],i));
342 }
343
344 // Constrain histograms
345 const Int_t nDimSingle=4;
346 const Int_t nDimPair=9;
347 Int_t nBinPair[nDimPair] = {nBinsPhi,nBinsC,fPtBinning.GetSize()-1,fInvMassBinning.GetSize()-1,nBinsSource,nBinsGen,fPtBinning.GetSize()-1,fEtaBinning.GetSize()-1,kBinsEtaAssociated};
348
349 // Associated Electron
350 Int_t nBinAssElectron[nDimSingle] = {nBinsC,fPtBinning.GetSize()-1,nBinsSource,kBinsEtaAssociated};
351 fAssElectron = new THnSparseF("fAssElectron","fAssElectron",nDimSingle,nBinAssElectron);
352 fAssElectron->SetBinEdges(0,binLimC);
353 fAssElectron->SetBinEdges(1,fPtBinning.GetArray());
354 fAssElectron->SetBinEdges(2,binLimSource);
355 fAssElectron->SetBinEdges(3,binLimEtaAssociat);
356 fAssElectron->Sumw2();
357 AliDebug(2,"AliHFENonPhotonicElectron: fAssElectron");
358
359 // Inclusive Electron
360 Int_t nBinIncElectron[nDimSingle] = {nBinsC,fPtBinning.GetSize()-1,nBinsSource,fEtaBinning.GetSize()-1};
361 fIncElectron = new THnSparseF("fIncElectron","fIncElectron",nDimSingle,nBinIncElectron);
362 fIncElectron->SetBinEdges(0,binLimC);
363 fIncElectron->SetBinEdges(1,fPtBinning.GetArray());
364 fIncElectron->SetBinEdges(2,binLimSource);
365 fIncElectron->SetBinEdges(3,fEtaBinning.GetArray());
366 fIncElectron->Sumw2();
367 AliDebug(2,"AliHFENonPhotonicElectron: fIncElectron");
368
369 // ee invariant mass Unlike Sign
370 fUSign = new THnSparseF("fUSign","fUSign",nDimPair,nBinPair);
371 fUSign->SetBinEdges(0,binLimPhi);
372 fUSign->SetBinEdges(1,binLimC);
373 fUSign->SetBinEdges(2,fPtBinning.GetArray());
374 fUSign->SetBinEdges(3,fInvMassBinning.GetArray());
375 fUSign->SetBinEdges(4,binLimSource);
376 fUSign->SetBinEdges(5,binLimGen);
377 fUSign->SetBinEdges(6,fPtBinning.GetArray());
378 fUSign->SetBinEdges(7,fEtaBinning.GetArray());
379 fUSign->SetBinEdges(8,binLimEtaAssociat);
380 fUSign->Sumw2();
381 AliDebug(2,"AliHFENonPhotonicElectron: fUSign");
382
383 // ee invariant mass Like Sign
384 fLSign = new THnSparseF("fLSign","fLSign",nDimPair,nBinPair);
385 fLSign->SetBinEdges(0,binLimPhi);
386 fLSign->SetBinEdges(1,binLimC);
387 fLSign->SetBinEdges(2,fPtBinning.GetArray());
388 fLSign->SetBinEdges(3,fInvMassBinning.GetArray());
389 fLSign->SetBinEdges(4,binLimSource);
390 fLSign->SetBinEdges(5,binLimGen);
391 fLSign->SetBinEdges(6,fPtBinning.GetArray());
392 fLSign->SetBinEdges(7,fEtaBinning.GetArray());
393 fLSign->SetBinEdges(8,binLimEtaAssociat);
394 fLSign->Sumw2();
395 AliDebug(2,"AliHFENonPhotonicElectron: fLSign");
396
397 // Histograms counting the number of like sign / unlike sign matches per inclusive track
398 const Int_t nBinsMatches = 50;
399 Double_t binLimMatches[nBinsMatches+1];
400 for(Int_t ib = 0; ib <= nBinsMatches; ib++) binLimMatches[ib] = ib;
401 const Int_t nDimMatches = 3; // centrality, pt_inc, number of matches
402 const Int_t nBinsMatchHist[nDimMatches] = {nBinsC, fPtBinning.GetSize()-1, nBinsMatches};
403 fUSmatches = new THnSparseF("fUSmatches", "fUSmatches", nDimMatches, nBinsMatchHist);
404 fUSmatches->SetBinEdges(0,binLimC);
405 fUSmatches->SetBinEdges(1,fPtBinning.GetArray());
406 fUSmatches->SetBinEdges(2,binLimMatches);
407
408 fLSmatches = new THnSparseF("fLSmatches", "fLSmatches", nDimMatches, nBinsMatchHist);
409 fLSmatches->SetBinEdges(0,binLimC);
410 fLSmatches->SetBinEdges(1,fPtBinning.GetArray());
411 fLSmatches->SetBinEdges(2,binLimMatches);
412
413 // Histograms for radius studies
414 Int_t nBinsradius = 13;
415 Double_t minradius = 0.0;
416 Double_t maxradius = 25.0;
417 Double_t binLimradius[nBinsradius+1];
418 for(Int_t i=0; i<=nBinsradius; i++) binLimradius[i]=(Double_t)minradius + (maxradius-minradius)/nBinsradius*(Double_t)i ;
419 const Int_t nDimIncElectronRadius = 4; // centrality, pt_inc, radius
420 const Int_t nBinsIncElectronRadius[nDimIncElectronRadius] = {nBinsC, fPtBinning.GetSize()-1, nBinsradius, nBinsSource};
421 fIncElectronRadius = new THnSparseF("fIncElectronRadius", "fIncElectronRadius", nDimIncElectronRadius, nBinsIncElectronRadius);
422 fIncElectronRadius->SetBinEdges(0,binLimC);
423 fIncElectronRadius->SetBinEdges(1,fPtBinning.GetArray());
424 fIncElectronRadius->SetBinEdges(2,binLimradius);
425 fIncElectronRadius->SetBinEdges(3,binLimSource);
426 fIncElectronRadius->Sumw2();
427
428 fRecElectronRadius = new THnSparseF("fRecElectronRadius", "fRecElectronRadius", nDimIncElectronRadius, nBinsIncElectronRadius);
429 fRecElectronRadius->SetBinEdges(0,binLimC);
430 fRecElectronRadius->SetBinEdges(1,fPtBinning.GetArray());
431 fRecElectronRadius->SetBinEdges(2,binLimradius);
432 fRecElectronRadius->SetBinEdges(3,binLimSource);
433 fRecElectronRadius->Sumw2();
434
435 /*
436 // ee angle Unlike Sign
437 const Int_t nDimUSignAngle=3;
438 Int_t nBinUSignAngle[nDimUSignAngle] = {nBinsAngle,nBinsC,nBinsSource};
439 fUSignAngle = new THnSparseF("fUSignAngle","fUSignAngle",nDimUSignAngle,nBinUSignAngle);
440 fUSignAngle->SetBinEdges(0,binLimAngle);
441 fUSignAngle->SetBinEdges(1,binLimC);
442 fUSignAngle->SetBinEdges(2,binLimSource);
443 fUSignAngle->Sumw2();
444 AliDebug(2,"AliHFENonPhotonicElectron: fUSignAngle");
445
446 // ee angle Like Sign
447 const Int_t nDimLSignAngle=3;
448 Int_t nBinLSignAngle[nDimLSignAngle] = {nBinsAngle,nBinsC,nBinsSource};
449 fLSignAngle = new THnSparseF("fLSignAngle","fLSignAngle",nDimLSignAngle,nBinLSignAngle);
450 fLSignAngle->SetBinEdges(0,binLimAngle);
451 fLSignAngle->SetBinEdges(1,binLimC);
452 fLSignAngle->SetBinEdges(2,binLimSource);
453 fLSignAngle->Sumw2();
454 AliDebug(2,"AliHFENonPhotonicElectron: fLSignAngle");
455 */
456
457 // control histogram for ITS PID
458 fHnsigmaITS = new TH2F("fHnsigmaITS", "Number of sigmas in the ITS", 30, 0., 0.3, 1200, -10., 10.);
459
460 // control histogram for weights sources
461 fWeightsSource = new TH2F("fWeightsSource", "Source code for weights", 13, -1.5, 11.5, 29, -1.5, 27.5);
462
463 fListOutput->Add(fAssElectron);
464 fListOutput->Add(fIncElectron);
465 fListOutput->Add(fUSign);
466 fListOutput->Add(fLSign);
467 fListOutput->Add(fUSmatches);
468 fListOutput->Add(fLSmatches);
469 fListOutput->Add(fHnsigmaITS);
470 fListOutput->Add(fWeightsSource);
471 fListOutput->Add(fIncElectronRadius);
472 fListOutput->Add(fRecElectronRadius);
473 // fListOutput->Add(fUSignAngle);
474 // fListOutput->Add(fLSignAngle);
76d0b522 475
476}
477
afb48e1d 478//_____________________________________________________________________________________________
479void AliHFENonPhotonicElectron::SetWithWeights(Int_t levelBack)
480{
465ff082 481 //
482 // Init the HFE level
483 //
484 if(levelBack >= 0) fLevelBack = levelBack;
afb48e1d 485
486}
487
488//_____________________________________________________________________________________________
489void AliHFENonPhotonicElectron::SetMCEvent(AliMCEvent *mcEvent)
490{
465ff082 491 //
492 // Pass the mcEvent
493 //
494
495 fMCEvent = mcEvent;
496
afb48e1d 497}
498
499//_____________________________________________________________________________________________
500void AliHFENonPhotonicElectron::SetAODArrayMCInfo(TClonesArray *aodArrayMCInfo)
501{
465ff082 502 //
503 // Pass the mcEvent info
504 //
505
506 fAODArrayMCInfo = aodArrayMCInfo;
afb48e1d 507
afb48e1d 508}
509
76d0b522 510//_____________________________________________________________________________________________
663a2523 511void AliHFENonPhotonicElectron::InitRun(const AliVEvent *inputEvent,const AliPIDResponse *pidResponse)
76d0b522 512{
465ff082 513 //
514 // Init run
515 //
516
517 if(!pidResponse)
518 {
519 AliDebug(1, "Using default PID Response");
520 Bool_t hasmc = kFALSE;
521 if(fMCEvent || fAODArrayMCInfo) hasmc=kTRUE;
522 pidResponse = AliHFEtools::GetDefaultPID(hasmc, inputEvent->IsA() == AliESDEvent::Class());
523 }
524
525 if(!fPIDBackground) return;
526 fPIDBackground->SetPIDResponse(pidResponse);
527
528 if(!fPIDBackground->IsInitialized())
529 {
530 // Initialize PID with the given run number
531 fPIDBackground->InitializePID(inputEvent->GetRunNumber());
532 }
76d0b522 533
534}
535
536//_____________________________________________________________________________________________
8a9b2231 537Int_t AliHFENonPhotonicElectron::FillPoolAssociatedTracks(AliVEvent *inputEvent, Int_t binct)
76d0b522 538{
4437a0d2 539 //
465ff082 540 // Fill the pool of associated tracks
541 // Return the number of associated tracks
542 //
543
544 fnumberfound = 0;
545
546 fHFEBackgroundCuts->SetRecEvent(inputEvent);
547 Int_t nbtracks = inputEvent->GetNumberOfTracks();
548
549 if( fArraytrack ){
550 fArraytrack->~TArrayI();
551 new(fArraytrack) TArrayI(nbtracks);
552 } else {
553 fArraytrack = new TArrayI(nbtracks);
76d0b522 554 }
76d0b522 555
465ff082 556 fCounterPoolBackground = 0;
557
558 Bool_t isSelected(kFALSE);
559 Bool_t isAOD = (dynamic_cast<AliAODEvent *>(inputEvent) != NULL);
560 AliDebug(2, Form("isAOD: %s", isAOD ? "yes" : "no"));
561 for(Int_t k = 0; k < nbtracks; k++) {
562 AliVTrack *track = (AliVTrack *) inputEvent->GetTrack(k);
563 if(!track) continue;
564
565 //
566 isSelected = kFALSE;
567 if(fSelectCategory1tracks && FilterCategory1Track(track, isAOD, binct)) isSelected = kTRUE;
568 else if(fSelectCategory2tracks && FilterCategory2Track(track, isAOD)) isSelected = kTRUE;
569
570 if(isSelected){
571 AliDebug(2,Form("fCounterPoolBackground %d, track %d",fCounterPoolBackground,k));
572 fArraytrack->AddAt(k,fCounterPoolBackground);
573 fCounterPoolBackground++;
574 }
575 } // loop tracks
576
577 //printf(Form("Associated Pool: Tracks %d, fCounterPoolBackground %d \n", nbtracks, fCounterPoolBackground));
76d0b522 578
465ff082 579 return fCounterPoolBackground;
8a9b2231 580
581}
582
583//_____________________________________________________________________________________________
584Int_t AliHFENonPhotonicElectron::CountPoolAssociated(AliVEvent *inputEvent, Int_t binct)
585{
465ff082 586 //
587 // Count the pool of assiocated tracks
588 //
8a9b2231 589
590
465ff082 591 if(fnumberfound > 0) //!count only events with an inclusive electron
592 {
593 Double_t valueAssElectron[4] = {(Double_t) binct, -1, -1}; //Centrality Pt Source
594 Int_t iTrack2 = 0;
595 Int_t indexmother2 = -1;
596 AliVTrack *track2 = 0x0;
8a9b2231 597
465ff082 598 for(Int_t ii = 0; ii < fCounterPoolBackground; ii++){
599 iTrack2 = fArraytrack->At(ii);
600 AliDebug(2,Form("track %d",iTrack2));
601 track2 = (AliVTrack *)inputEvent->GetTrack(iTrack2);
8a9b2231 602
465ff082 603 if(!track2){
604 //printf("ERROR: Could not receive track %d", iTrack2);
605 continue;
606 }
8a9b2231 607
465ff082 608 // if MC look
609 if(fMCEvent || fAODArrayMCInfo) valueAssElectron[2] = FindMother(TMath::Abs(track2->GetLabel()), indexmother2) ;
8a9b2231 610
465ff082 611 fkPIDRespons = fPIDBackground->GetPIDResponse();
8a9b2231 612
465ff082 613 valueAssElectron[1] = track2->Pt() ;
614 valueAssElectron[3] = track2->Eta() ;
8a9b2231 615
465ff082 616 fAssElectron->Fill( valueAssElectron) ;
617 }
618 //printf(Form("Associated Pool: fCounterPoolBackground %d \n", fCounterPoolBackground));
8a9b2231 619 }
465ff082 620 return fnumberfound;
76d0b522 621}
622
623//_____________________________________________________________________________________________
afb48e1d 624Int_t AliHFENonPhotonicElectron::LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *vEvent, Double_t weight, Int_t binct, Double_t deltaphi, Int_t source, Int_t indexmother,Int_t mcQAsource)
76d0b522 625{
465ff082 626 //
627 // Look At Non HFE
628 //
629
630 /***********************************************************************************
631 * *
632 * iTrack1: index of the tagged electrons in AliVEvent *
633 * track1: tagged electron *
634 * vEvent: event *
635 * weight: weight in pt if not realistic *
636 * binct: centrality bin *
637 * deltaphi: phi-phi event plane for v2 *
638 * source: MC sources *
639 * indexmother: MC index mother *
640 * *
641 * *
642 * return -1 if nothing *
643 * return 2 if opposite charge within the mass range *
644 * return 4 if like charge within the mass range *
645 * return 6 if opposite & like charge within the mass range *
646 * *
647 ***********************************************************************************/
648
649 //printf("weight %f and source %d\n",weight,source);
650
651
652 AliAODEvent *aodeventu = dynamic_cast<AliAODEvent*>(vEvent);
653 Int_t taggedphotonic = -1;
654
655 AliDebug(2,Form("fCounterPoolBackground %d in LookAtNonHFE!!!",fCounterPoolBackground));
656 if(!fArraytrack) return taggedphotonic;
657 AliDebug(2,Form("process track %d",iTrack1));
658 AliDebug(1,Form("Inclusive source is %d\n", source));
659
660 fkPIDRespons = fPIDBackground->GetPIDResponse();
661
662 //Set Fill-Arrays for THnSparse
663 Double_t valueIncElectron[4] = { (Double_t) binct, track1->Pt(), (Double_t)source, track1->Eta()}; //Centrality Pt Source P
664 Double_t valueSign[9] = { deltaphi, (Double_t)binct, track1->Pt(), -1, (Double_t)source, -1, -1, track1->Eta(), -1}; //DeltaPhi Centrality Pt InvariantMass Source Angle Pt
665 //Double_t valueAngle[3] = { -1, binct, source};
666 Double_t valueradius[4] = { (Double_t)binct, track1->Pt(), 0.,(Double_t)source};
667
668
669 Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()));
670 Double_t radius = Radius(TMath::Abs(track1->GetLabel()));
671 AliKFParticle::SetField(vEvent->GetMagneticField());
672 AliKFVertex primV(*(vEvent->GetPrimaryVertex()));
673 valueradius[2] = radius;
674
675 AliVTrack *track2(NULL);
676 Int_t iTrack2 = 0;
677 Int_t indexmother2 = -1;
678 Int_t pdg2 = -100;
679 Int_t source2 = -1;
680 Float_t fCharge2 = 0;
681
682 // count number of matches with opposite/same sign track in the given mass range
683 Int_t countsMatchLikesign(0),
684 countsMatchUnlikesign(0);
685
686 // Values to fill
687 Double_t angle(-1.);
688 Double_t invmass(-1);
76d0b522 689
465ff082 690 Float_t fCharge1 = track1->Charge(); //Charge from track1
691
692 Bool_t kUSignPhotonic = kFALSE;
693 Bool_t kLSignPhotonic = kFALSE;
694
695 //
696 //TEMPORARY SWITCH
697 //changes the weighting factor for electrons from Eta Dalitz decays
698 //
699 if(mcQAsource == AliHFEmcQA::kEta) weight *= fEtaDalitzWeightFactor;
700
701 // --------------------------
702 // | Switches for Generations |
703 // --------------------------
704 // Number of Generations
705 // Display MC stack for likesign pairs with common mother in generation 1.
706
707 //! FILL Inclusive Electron
708 fWeightsSource->Fill(source,mcQAsource);
709 fIncElectron->Fill(valueIncElectron,weight);
710 fnumberfound++;
711 if(fStudyRadius) fIncElectronRadius->Fill(valueradius,weight);
712
713 //printf(Form("Inclusive Pool: TrackNr. %d, fnumberfound %d \n", iTrack1, fnumberfound));
714
715 for(Int_t idex = 0; idex < fCounterPoolBackground; idex++){
716 iTrack2 = fArraytrack->At(idex);
717 AliDebug(2,Form("track %d",iTrack2));
718 track2 = (AliVTrack *)vEvent->GetTrack(iTrack2);
719
720 if(!track2){
721 //printf("ERROR: Could not receive track %d", iTrack2);
722 continue;
c3e32eae 723 }
76d0b522 724
465ff082 725 fCharge2 = track2->Charge(); //Charge from track2
726
727 // Reset the MC info
728 //valueAngle[2] = source;
729 valueradius[3] = source;
730 valueSign[4] = source;
731 valueSign[6] = track2->Pt();
732 valueSign[8] = track2->Eta();
733
734 // track cuts and PID already done
735
736 // Checking if it is the same Track!
737 if(iTrack2==iTrack1) continue;
738 AliDebug(2,"Different");
739
740 // if MC look
741 if(fMCEvent || fAODArrayMCInfo){
742 AliDebug(2, "Checking for source");
743 source2 = FindMother(TMath::Abs(track2->GetLabel()), indexmother2);
744 AliDebug(2, Form("source is %d", source2));
745 AliDebug(2, Form("sourceindex is %i", indexmother2));
746 AliDebug(2, Form("getlabel: %i", track2->GetLabel()));
747 pdg2 = CheckPdg(TMath::Abs(track2->GetLabel()));
748
749 if(source == kElectronfromconversion){
750 AliDebug(2, Form("Electron from conversion (source %d), paired with source %d", source, source2));
751 AliDebug(2, Form("Index of the mothers: incl %d, associated %d", indexmother, indexmother2));
752 AliDebug(2, Form("PDGs: incl %d, associated %d", pdg1, pdg2));
753 }
754
755 if(source2 >=0 ){
756 AliDebug(1,"------------------------------------ \n");
757 if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0)){
758 AliDebug(1, "Real pair");
759 switch(source){
760 case kElectronfromconversion:
761 valueSign[4] = kElectronfromconversionboth;
762 valueradius[3] = kElectronfromconversionboth;
763 break;
764 case kElectronfrompi0:
765 valueSign[4] = kElectronfrompi0both;
766 valueradius[3] = kElectronfrompi0both;
767 break;
768 case kElectronfrometa:
769 valueSign[4] = kElectronfrometaboth;
770 valueradius[3] = kElectronfrometaboth;
771 break;
772 case kElectronfromomega:
773 valueSign[4] = kElectronfromomegaboth;
774 valueradius[3] = kElectronfromomegaboth;
775 break;
776 };
777 }
778
779 if(fAnaPairGen){
780
781 Int_t MotherArray1[fNumberofGenerations];
782 Int_t MotherArray2[fNumberofGenerations];
783 for(Int_t i = 0; i < fNumberofGenerations;++i){
784 MotherArray1[i]=-1;
785 MotherArray2[i]=-1;
786 }
787 FillMotherArray(TMath::Abs(track1->GetLabel()),0,MotherArray1,fNumberofGenerations);
788 FillMotherArray(TMath::Abs(track2->GetLabel()),0,MotherArray2,fNumberofGenerations);
789 AliDebug(2,Form(" indextrack inclusive: %i pdg: %i || indextrack assoc: %i pdg: %i \n",track1->GetLabel(),pdg1 , track2->GetLabel(),pdg2));
790 AliDebug(2,Form(" Mother Gen 1: %i || Mother Gen 1: %i \n", MotherArray1[0],MotherArray2[0]));
791 AliDebug(2,Form(" Mother Gen 2: %i || Mother Gen 2: %i \n", MotherArray1[1],MotherArray2[1]));
792 AliDebug(2,Form(" Mother Gen 3: %i || Mother Gen 3: %i \n", MotherArray1[2],MotherArray2[2]));
793 AliDebug(2,Form(" Mother Gen 4: %i || Mother Gen 4: %i \n", MotherArray1[3],MotherArray2[3]));
794 valueSign[5] = FindGeneration(MotherArray1,MotherArray2,fNumberofGenerations);
795 }
796 }
797 }
798
799 if(fAlgorithmMA){
800 // Use TLorentzVector
801 if(!MakePairDCA(track1, track2, vEvent, (aodeventu != NULL), invmass, angle)) continue;
802 } else {
803 // Use AliKF package
804 if(!MakePairKF(track1, track2, primV, invmass, angle)) continue;
805 }
76d0b522 806
465ff082 807 valueSign[3] = invmass;
808 // valueSign[5] = angle;
76d0b522 809
465ff082 810 //if((fCharge1*fCharge2)>0.0) fLSignAngle->Fill(&valueAngle[0],weight);
811 //else fUSignAngle->Fill(&valueAngle[0],weight);
76d0b522 812
465ff082 813 if(angle > fMaxOpening3D) continue; //! Cut on Opening Angle
814 if(invmass > fMaxInvMass) continue; //! Cut on Invariant Mass
8a9b2231 815
465ff082 816 if((fCharge1*fCharge2)>0.0){
817 if(invmass < 1.0){
818 fLSign->Fill( valueSign, weight);
819 //if switched on produces mcstackdump for likesign with commonmother in gen 1
820 if(fDisplayMCStack && fMCEvent && valueSign[5] == 1) {
821 AliStack* stack = fMCEvent->Stack();
822 if(!stack) AliError("No Stack");
823 else stack->DumpPStack();
824 }
825
826 }
827 // count like-sign background matched pairs per inclusive based on mass cut
828 if(invmass < 0.14) countsMatchLikesign++;
829 AliDebug(1, "-> Selected Like sign");
830 } else {
831 if(invmass < 1.0){
832 fUSign->Fill( valueSign, weight);
833 }
834 // count unlike-sign matched pairs per inclusive based on mass cut
835 if(invmass < 0.14) {
836 countsMatchUnlikesign++;
837 if(fStudyRadius) fRecElectronRadius->Fill(valueradius,weight);
838 }
839 AliDebug(1, "-> Selected Unlike sign");
840
841 }
4437a0d2 842
76d0b522 843
465ff082 844 if((fCharge1*fCharge2)>0.0) kLSignPhotonic=kTRUE;
845 else kUSignPhotonic=kTRUE;
846 }
847
848 // Fill counted
849 Double_t valCountsLS[3] = {(Double_t)binct, track1->Pt(),(Double_t)countsMatchLikesign},
850 valCountsUS[3] = {(Double_t)binct, track1->Pt(),(Double_t)countsMatchUnlikesign};
851 fUSmatches->Fill(valCountsUS);
852 fLSmatches->Fill(valCountsLS);
c3e32eae 853
465ff082 854 if( kUSignPhotonic && kLSignPhotonic) taggedphotonic = 6;
855 if(!kUSignPhotonic && kLSignPhotonic) taggedphotonic = 4;
856 if( kUSignPhotonic && !kLSignPhotonic) taggedphotonic = 2;
76d0b522 857
465ff082 858 AliDebug(1,"------------------------------------ \n");
859 return taggedphotonic;
76d0b522 860}
861
862//_________________________________________________________________________
39427568 863Int_t AliHFENonPhotonicElectron::FindMother(Int_t tr, Int_t &indexmother) const {
465ff082 864 //
865 // Find the mother if MC
866 //
867
868 if(!fMCEvent && !fAODArrayMCInfo) return -1;
869
870 Int_t pdg = CheckPdg(tr);
871 if(TMath::Abs(pdg)!= 11)
872 {
873 indexmother = -1;
874 return kNoElectron;
875 }
876
877 indexmother = IsMotherGamma(tr);
878 if(indexmother > 0) return kElectronfromconversion;
879 indexmother = IsMotherPi0(tr);
880 if(indexmother > 0) return kElectronfrompi0;
881 indexmother = IsMotherC(tr);
882 if(indexmother > 0) return kElectronfromC;
883 indexmother = IsMotherB(tr);
884 if(indexmother > 0) return kElectronfromB;
885 indexmother = IsMotherEta(tr);
886 if(indexmother > 0) return kElectronfrometa;
887 indexmother = IsMotherOmega(tr);
888 if(indexmother > 0) return kElectronfromomega;
889
890 return kElectronfromother;
76d0b522 891}
892
893//________________________________________________________________________________________________
39427568 894Int_t AliHFENonPhotonicElectron::CheckPdg(Int_t tr) const {
76d0b522 895
465ff082 896 //
897 // Return the pdg of the particle
898 //
76d0b522 899
465ff082 900 Int_t pdgcode = -1;
901 if(tr < 0) return pdgcode;
76d0b522 902
465ff082 903 AliMCParticle *mctrackesd = NULL; AliAODMCParticle *mctrackaod = NULL;
904 if(fMCEvent){
905 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
906 if(mctrack){
907 if((mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) pdgcode = mctrackesd->PdgCode();
908 else if((mctrackaod = dynamic_cast<AliAODMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) pdgcode = mctrackaod->GetPdgCode();
909 }
910 } else if(fAODArrayMCInfo) {
911 if(tr < fAODArrayMCInfo->GetEntriesFast()){
912 mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
913 if(mctrackaod) return pdgcode = mctrackaod->GetPdgCode();
914 }
39427568 915 }
76d0b522 916
465ff082 917 return pdgcode;
76d0b522 918}
919
afb48e1d 920//________________________________________________________________________________________________
921Double_t AliHFENonPhotonicElectron::Radius(Int_t tr) const {
922
465ff082 923 //
924 // Return the production vertex radius
925 //
afb48e1d 926
465ff082 927 Double_t radius = 0.;
928 if(tr < 0) return radius;
afb48e1d 929
465ff082 930 AliMCParticle *mctrackesd = NULL; AliAODMCParticle *mctrackaod = NULL;
931 if(fMCEvent){
932 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
933 if(mctrack){
934 if((mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) radius = TMath::Sqrt(mctrackesd->Xv()*mctrackesd->Xv()+mctrackesd->Yv()*mctrackesd->Yv());
935 else if((mctrackaod = dynamic_cast<AliAODMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) radius = TMath::Sqrt(mctrackaod->Xv()*mctrackaod->Xv()+mctrackaod->Yv()*mctrackaod->Yv());
936 }
937 } else if(fAODArrayMCInfo) {
938 if(tr < fAODArrayMCInfo->GetEntriesFast()){
939 mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
940 if(mctrackaod) return radius = TMath::Sqrt(mctrackaod->Xv()*mctrackaod->Xv()+mctrackaod->Yv()*mctrackaod->Yv());
941 }
afb48e1d 942 }
afb48e1d 943
465ff082 944 return radius;
afb48e1d 945}
946
76d0b522 947//_______________________________________________________________________________________________
39427568 948Int_t AliHFENonPhotonicElectron::GetMotherPDG(Int_t tr, Int_t &motherIndex) const {
465ff082 949 //
950 // Returns the mother PDG of the track (return value) and the index of the
951 // mother track
952 //
953 if(tr < 0) return -1;
954
955 Int_t pdg(-1);
956 AliMCParticle *mctrackesd(NULL); AliAODMCParticle *mctrackaod(NULL);
957
958 motherIndex = -1;
959 if(fMCEvent) {
960 AliDebug(2, "Using MC Event");
961 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
962 if(mctrack){
963 if((mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))){
964 // Case ESD
965 TParticle *particle = mctrackesd->Particle();
966
967 // Take mother
968 if(particle){
969 motherIndex = particle->GetFirstMother();
970 if(motherIndex >= 0){
971 AliMCParticle *mothertrack = NULL;
972 if((mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(motherIndex))))){
973 TParticle * mother = mothertrack->Particle();
974 pdg = mother->GetPdgCode();
975 }
976 }
977 }
978 } else if((mctrackaod = dynamic_cast<AliAODMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))){
979 // Case AOD
980 // Take mother
981 motherIndex = mctrackaod->GetMother();
982 if(motherIndex >= 0){
983 AliAODMCParticle *mothertrack = dynamic_cast<AliAODMCParticle *>(fMCEvent->GetTrack(motherIndex));
984 if(mothertrack) pdg = mothertrack->GetPdgCode();
985 }
39427568 986 }
39427568 987 }
465ff082 988 } else if(fAODArrayMCInfo) {
989 AliDebug(2, "Using AOD list");
990 if(tr < fAODArrayMCInfo->GetEntriesFast()){
991 mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
992
993 // Take mother
994 if(mctrackaod){
995 motherIndex = mctrackaod->GetMother();
996 if(motherIndex >= 0 && motherIndex < fAODArrayMCInfo->GetEntriesFast()){
997 AliAODMCParticle *mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(motherIndex)));
998 if(mothertrack) pdg = mothertrack->GetPdgCode();
999 }
1000 }
39427568 1001 }
76d0b522 1002 }
465ff082 1003 return pdg;
39427568 1004}
76d0b522 1005
39427568 1006//_______________________________________________________________________________________________
1007Int_t AliHFENonPhotonicElectron::IsMotherGamma(Int_t tr) const {
1008
465ff082 1009 //
1010 // Return the lab of gamma mother or -1 if not gamma
1011 //
1012
1013 Int_t imother(-1), pdg(-1);
1014 pdg = GetMotherPDG(tr, imother);
1015
1016 // Check gamma
1017 if(imother >= 0){
1018 if(TMath::Abs(pdg) == 22){
1019 AliDebug(2, "Gamma Mother selected");
1020 return imother;
1021 }
1022 if(TMath::Abs(pdg) == 11){
1023 AliDebug(2, "Mother is electron - look further in hierarchy");
1024 return IsMotherGamma(imother);
1025 }
1026 AliDebug(2, "Nothing selected");
1027 return -1;
76d0b522 1028 }
465ff082 1029 AliDebug(2, "Not mother");
76d0b522 1030 return -1;
76d0b522 1031}
1032
1033//________________________________________________________________________________________________
39427568 1034Int_t AliHFENonPhotonicElectron::IsMotherPi0(Int_t tr) const {
76d0b522 1035
465ff082 1036 //
1037 // Return the lab of pi0 mother or -1 if not pi0
1038 //
76d0b522 1039
465ff082 1040 Int_t imother(-1), pdg(-1);
1041 pdg = GetMotherPDG(tr, imother);
76d0b522 1042
465ff082 1043 // Check pi0
1044 if(imother >= 0){
1045 if(TMath::Abs(pdg) == 111){
1046 AliDebug(2, "Pi0 Mother selected");
1047 return imother;
1048 }
1049 if(TMath::Abs(pdg) == 11){
1050 AliDebug(2, "Mother is electron - look further in hierarchy");
1051 return IsMotherPi0(imother);
1052 }
1053 AliDebug(2, "Nothing selected");
1054 return -1;
76d0b522 1055 }
465ff082 1056 AliDebug(2, "Not mother");
76d0b522 1057 return -1;
76d0b522 1058}
1059//________________________________________________________________________________________________
39427568 1060Int_t AliHFENonPhotonicElectron::IsMotherC(Int_t tr) const {
76d0b522 1061
465ff082 1062 //
1063 // Return the lab of signal mother or -1 if not from C
1064 //
76d0b522 1065
465ff082 1066 Int_t imother(-1), pdg(-1);
1067 pdg = GetMotherPDG(tr, imother);
76d0b522 1068
1069 // Check C
465ff082 1070 if(imother >= 0){
1071 if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)){
1072 AliDebug(2, "Charm Mother selected");
1073 return imother;
1074 }
1075 if(TMath::Abs(pdg) == 11){
1076 AliDebug(2, "Mother is electron - look further in hierarchy");
1077 return IsMotherC(imother);
1078 }
1079 AliDebug(2, "Nothing selected");
1080 return -1;
76d0b522 1081 }
465ff082 1082 AliDebug(2, "Not mother");
76d0b522 1083 return -1;
76d0b522 1084}
39427568 1085
76d0b522 1086//_______________________________________________________________________________________________
39427568 1087Int_t AliHFENonPhotonicElectron::IsMotherB(Int_t tr) const {
76d0b522 1088
465ff082 1089 //
1090 // Return the lab of signal mother or -1 if not B
1091 //
76d0b522 1092
465ff082 1093 Int_t imother(-1), pdg(-1);
1094 pdg = GetMotherPDG(tr, imother);
76d0b522 1095
465ff082 1096 // Check B
1097 if(imother >= 0){
1098 if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)){
1099 AliDebug(2, "Bottom Mother selected");
1100 return imother;
1101 }
1102 if(TMath::Abs(pdg) == 11){
1103 return IsMotherB(imother);
1104 AliDebug(2, "Mother is electron - look further in hierarchy");
1105 }
1106 AliDebug(2, "Nothing selected");
1107 return -1;
76d0b522 1108 }
465ff082 1109 AliDebug(2, "Not mother");
76d0b522 1110 return -1;
76d0b522 1111}
1112
1113//_______________________________________________________________________________________________
39427568 1114Int_t AliHFENonPhotonicElectron::IsMotherEta(Int_t tr) const {
76d0b522 1115
465ff082 1116 //
1117 // Return the lab of eta mother or -1 if not eta
1118 //
76d0b522 1119
465ff082 1120 Int_t imother(-1), pdg(-1);
1121 pdg = GetMotherPDG(tr, imother);
76d0b522 1122
465ff082 1123 // Check eta
1124 if(imother >= 0){
1125 if(TMath::Abs(pdg) == 221){
1126 AliDebug(2, "Eta mother selected");
1127 return imother;
1128 }
1129 if(TMath::Abs(pdg) == 11){
1130 AliDebug(2, "Mother is electron - look further in hierarchy");
1131 return IsMotherEta(imother);
1132 }
1133 AliDebug(2, "Nothing selected");
1134 return -1;
76d0b522 1135 }
465ff082 1136 AliDebug(2, "Not mother");
1137 return -1;
1138}
1139
1140//_______________________________________________________________________________________________
1141Int_t AliHFENonPhotonicElectron::IsMotherOmega(Int_t tr) const {
1142
1143 //
1144 // Return the lab of omega mother or -1 if not omega
1145 //
1146
1147 Int_t imother(-1), pdg(-1);
1148 pdg = GetMotherPDG(tr, imother);
1149
1150 // Check omega
1151 if(imother >= 0){
1152 if(TMath::Abs(pdg) == 223){
1153 AliDebug(2, "Omega mother selected");
1154 return imother;
1155 }
1156 if(TMath::Abs(pdg) == 11){
1157 AliDebug(2, "Mother is electron - look further in hierarchy");
1158 return IsMotherOmega(imother);
1159 }
1160 AliDebug(2, "Nothing selected");
1161 return -1;
76d0b522 1162 }
465ff082 1163 AliDebug(2, "Not mother");
76d0b522 1164 return -1;
76d0b522 1165}
4437a0d2 1166
1167//_______________________________________________________________________________________________
1168Bool_t AliHFENonPhotonicElectron::MakePairDCA(const AliVTrack *inclusive, const AliVTrack *associated, AliVEvent *vEvent, Bool_t isAOD, Double_t &invMass, Double_t &angle) const {
465ff082 1169 //
1170 // Make Pairs of electrons using TLorentzVector
1171 //
1172 Double_t eMass = TDatabasePDG::Instance()->GetParticle(11)->Mass(); //Electron mass in GeV
1173 Double_t bfield = vEvent->GetMagneticField();
1174
1175 AliESDtrack *esdtrack1, *esdtrack2;
1176 if(isAOD){
1177 // call copy constructor for AODs
1178 esdtrack1 = new AliESDtrack(inclusive);
1179 esdtrack2 = new AliESDtrack(associated);
1180 } else {
1181 // call copy constructor for ESDs
1182 esdtrack1 = new AliESDtrack(*(static_cast<const AliESDtrack *>(inclusive)));
1183 esdtrack2 = new AliESDtrack(*(static_cast<const AliESDtrack *>(associated)));
1184 }
1185 if((!esdtrack1) || (!esdtrack2)){
1186 delete esdtrack1;
1187 delete esdtrack2;
1188 return kFALSE;
1189 }
1190
1191 Double_t xt1 = 0; //radial position track 1 at the DCA point
1192 Double_t xt2 = 0; //radial position track 2 at the DCA point
1193 Double_t dca = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1); //DCA track1-track2
1194 if(dca > fMaxDCA){
1195 // Apply DCA cut already in the function
1196 delete esdtrack1;
1197 delete esdtrack2;
1198 return kFALSE;
1199 }
1200
1201 //Momenta of the track extrapolated to DCA track-track
1202 Double_t p1[3] = {0,0,0};
1203 Double_t p2[3] = {0,0,0};
1204 Bool_t kHasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1); //Track1
1205 Bool_t kHasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2); //Track2
1206 if(!kHasdcaT1 || !kHasdcaT2) AliWarning("It could be a problem in the extrapolation");
1207
1208 TLorentzVector electron1, electron2, mother;
1209 electron1.SetXYZM(p1[0], p1[1], p1[2], eMass);
1210 electron2.SetXYZM(p2[0], p2[1], p2[2], eMass);
1211
1212 mother = electron1 + electron2;
1213 invMass = mother.M();
1214 angle = TVector2::Phi_0_2pi(electron1.Angle(electron2.Vect()));
1215
4437a0d2 1216 delete esdtrack1;
1217 delete esdtrack2;
465ff082 1218 return kTRUE;
4437a0d2 1219}
1220
1221//_______________________________________________________________________________________________
1222Bool_t AliHFENonPhotonicElectron::MakePairKF(const AliVTrack *inclusive, const AliVTrack *associated, AliKFVertex &primV, Double_t &invMass, Double_t &angle) const {
465ff082 1223 //
1224 // Make pairs of electrons using the AliKF package
1225 //
1226
1227 //printf("AOD HFE non photonic\n");
1228
1229 Int_t fPDGtrack1 = 11;
1230 if(inclusive->Charge()>0) fPDGtrack1 = -11;
1231 Int_t fPDGtrack2 = 11;
1232 if(associated->Charge()>0) fPDGtrack2 = -11;
1233
1234 AliKFParticle ktrack1(*inclusive, fPDGtrack1);
1235 AliKFParticle ktrack2(*associated, fPDGtrack2);
1236 AliKFParticle recoGamma(ktrack1,ktrack2);
1237
1238 if(recoGamma.GetNDF()<1) return kFALSE; //! Cut on Reconstruction
1239
1240 Double_t chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();
1241 if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) return kFALSE;
1242
1243 if(fSetMassConstraint){
1244 primV += recoGamma;
1245 primV -= ktrack1;
1246 primV -= ktrack2;
1247 recoGamma.SetProductionVertex(primV);
1248 recoGamma.SetMassConstraint(0,0.0001);
1249 }
1250
1251 Double_t width(0.);
1252 recoGamma.GetMass(invMass,width);
1253 angle = ktrack1.GetAngle(ktrack2);
1254 return kTRUE;
4437a0d2 1255}
1256
1257//_______________________________________________________________________________________________
1258Bool_t AliHFENonPhotonicElectron::FilterCategory1Track(const AliVTrack * const track, Bool_t isAOD, Int_t binct){
465ff082 1259 //
1260 // Selection of good associated tracks for the pool
1261 // selection is done using strong cuts
1262 // Tracking in the TPC and the ITS is a minimal requirement
1263 //
1264 Bool_t survivedbackground = kTRUE;
1265
1266 if(!fHFEBackgroundCuts->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *) track)) survivedbackground = kFALSE;
1267 AliDebug(3, Form("First cut: %s\n", survivedbackground == kTRUE ? "yes" : "no"));
1268 if(!fHFEBackgroundCuts->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, (TObject *) track)) survivedbackground = kFALSE;
1269 AliDebug(3, Form("Second cut: %s\n", survivedbackground == kTRUE ? "yes" : "no"));
1270
1271 if(survivedbackground){
1272 // PID track cuts
1273 AliHFEpidObject hfetrack2;
1274
1275 if(!isAOD) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1276 else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);
1277
1278 hfetrack2.SetRecTrack(track);
1279 if(binct>-1){
1280 hfetrack2.SetCentrality((Int_t)binct);
1281 AliDebug(3,Form("centrality %d and %d",binct,hfetrack2.GetCentrality()));
1282 hfetrack2.SetPbPb();
1283 }
4437a0d2 1284
465ff082 1285 if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundQA)){
1286 survivedbackground = kTRUE;
1287 } else survivedbackground = kFALSE;
1288 }
1289 AliDebug(3, Form("PID: %s\n", survivedbackground == kTRUE ? "yes" : "no"));
1290 return survivedbackground;
4437a0d2 1291}
1292
1293//_______________________________________________________________________________________________
1294Bool_t AliHFENonPhotonicElectron::FilterCategory2Track(const AliVTrack * const track, Bool_t isAOD){
465ff082 1295 //
1296 // Selection of category 2 tracks: These tracks are exclusively low pt tracks below
1297 // 300 MeV/c which have at least 2 hits in the 4 outer ITS layers and are identified as
1298 // electron candidates by the ITS
1299 //
1300 if(TMath::Abs(track->Pt()) > 0.3) return kFALSE;
1301 if(TMath::Abs(track->Pt()) < fminPt) return kFALSE;
1302 Int_t nclustersITS(0), nclustersOuter(0);
1303 if(isAOD){
1304 const AliAODTrack *aodtrack = static_cast<const AliAODTrack *>(track);
1305 if(!(aodtrack->TestFilterBit(AliAODTrack::kTrkGlobalNoDCA) || aodtrack->TestFilterBit(AliAODTrack::kTrkITSsa))) return kFALSE;
1306 if(!aodtrack->IsOn(AliAODTrack::kITSrefit)) return kFALSE;
1307 nclustersITS = aodtrack->GetITSNcls();
1308 for(Int_t ily = 2; ily < 5; ily++)
1309 if(aodtrack->HasPointOnITSLayer(ily)) nclustersOuter++;
1310 } else {
1311 const AliESDtrack *esdtrack = static_cast<const AliESDtrack *>(track);
1312 if(esdtrack->GetStatus() & AliESDtrack::kITSpureSA) return kFALSE;
1313 if(esdtrack->GetStatus() & AliESDtrack::kTPCin) return kFALSE;
1314 if(!(esdtrack->GetStatus() & AliESDtrack::kITSrefit)) return kFALSE;
1315 nclustersITS = esdtrack->GetITSclusters(NULL);
1316 for(Int_t ily = 2; ily < 5; ily++)
1317 if(esdtrack->HasPointOnITSLayer(ily)) nclustersOuter++;
1318 }
1319 if(nclustersITS < 4) return kFALSE;
1320 if(nclustersOuter < 3) return kFALSE;
1321
1322 // Do ITS PID
1323 Double_t nsigmaITS = fPIDBackground->GetPIDResponse()->NumberOfSigmasITS(track, AliPID::kElectron);
1324 fHnsigmaITS->Fill(track->Pt(), nsigmaITS);
1325 if((nsigmaITS - fITSmeanShift) > fITSnSigmaHigh ) return kFALSE;
1326 if((nsigmaITS - fITSmeanShift) < fITSnSigmaLow ) return kFALSE;
1327 // if global track, we apply also TPC PID
1328 return kTRUE;
1329}
1330
1331//_______________________________________________________________________________________________
1332void AliHFENonPhotonicElectron::FillMotherArray(Int_t tr, Int_t index, Int_t a[], Int_t NumberofGenerations){
1333 //
1334 // Fill an array of mothers for a particle of a given trackid
1335 //
1336 Int_t imother(-1),pdg(-1);
1337 pdg = GetMotherPDG(tr,imother);
1338 if(index < NumberofGenerations){
1339 if(TMath::Abs(pdg) == 11){
1340 AliDebug(2, "Mother is electron, look further to fill motherarray");
1341 FillMotherArray(imother,index,a,NumberofGenerations);
1342 } else {
1343 a[index] = imother;
1344 index++;
1345 FillMotherArray(imother,index,a,NumberofGenerations);
1346 }
1347 }
1348}
1349
1350
1351//_______________________________________________________________________________________________
1352Int_t AliHFENonPhotonicElectron::FindGeneration(Int_t a[], Int_t b[], Int_t NumberofGenerations){
1353 //
1354 // Find the common mother in 2 given arrays and return its generation.
1355 //
1356 Int_t Generation = 0;
1357 for(Int_t iGeneration = 0; iGeneration < NumberofGenerations; iGeneration++){
1358 for(Int_t jGeneration = 0; jGeneration < NumberofGenerations; jGeneration++){
1359 if(a[iGeneration] == b[jGeneration] && a[iGeneration]!= -1 && b[jGeneration]!= -1){
1360 AliDebug(3,Form("* compared a and b %i == %i \n", a[iGeneration], b[jGeneration]));
1361 if(iGeneration > jGeneration){
1362 Generation = iGeneration + 1;
1363 } else {
1364 Generation = jGeneration + 1;
1365 }
1366 AliDebug(1,Form(" Found common mother gen: %i \n",Generation));
1367 return Generation;
1368 }
1369 }
1370 }
1371 AliDebug(1,Form(" Found common mother gen: %i \n",Generation));
1372 return Generation;
4437a0d2 1373}
465ff082 1374
1375