]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliHFENonPhotonicElectron.cxx
Sorry, I will manage....change back
[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 *
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
17/*************************************************************************************
18 * *
19 * Class for the Selection of Non-Heavy-Flavour-Electrons trought *
20 * the invariant mass method. The selection can be done from two *
21 * different algorithms, which can be choosed calling the function *
22 * "SetAlgorithm(TString Algorithm)". *
23 * *
24 * Authors: R.Bailhache, C.A.Schmidt *
25 * *
26 *************************************************************************************/
27
28#include "TVector2.h"
29#include "THnSparse.h"
30#include "TMath.h"
31#include "TLorentzVector.h"
32#include "TParticle.h"
33#include "TList.h"
34#include "TDatabasePDG.h"
35
36#include "AliVEvent.h"
37#include "AliMCEvent.h"
38#include "AliESDEvent.h"
39#include "AliMCParticle.h"
40#include "AliAODMCParticle.h"
41#include "AliAODEvent.h"
42#include "AliAODVertex.h"
43#include "AliAODTrack.h"
44#include "AliVTrack.h"
45#include "AliESDtrack.h"
46#include "AliESDtrackCuts.h"
47#include "AliPIDResponse.h"
48
49#include "AliKFParticle.h"
50#include "AliKFVertex.h"
51
52#include "AliHFEcuts.h"
53#include "AliHFEpid.h"
54#include "AliHFEpidQAmanager.h"
55#include "AliHFEtools.h"
56
57#include "AliHFENonPhotonicElectron.h"
58
59ClassImp(AliHFENonPhotonicElectron)
60//________________________________________________________________________
61AliHFENonPhotonicElectron::AliHFENonPhotonicElectron(const char *name, const Char_t *title)
62 :TNamed(name, title)
63 ,fMCEvent(NULL)
64 ,fAODArrayMCInfo(NULL)
65 ,fHFEBackgroundCuts(0)
66 ,fPIDBackground(0x0)
67 ,fPIDBackgroundQA(0)
68 ,fAlgorithmMA(kTRUE)
69 ,fUseFilterAOD(kTRUE)
70 ,fFilter(-1)
71 ,fChi2OverNDFCut(3.0)
72 ,fMaxDCA(3.0)
73 ,fMaxOpeningTheta(0.02)
74 ,fMaxOpeningPhi(0.1)
75 ,fMaxOpening3D(TMath::Pi())
76 ,fMaxInvMass(0.6)
77 ,fSetMassConstraint(kFALSE)
78 ,fArraytrack(NULL)
79 ,fCounterPoolBackground(0)
80 ,fListOutput(NULL)
81 ,fMCSourceee(NULL)
82 ,fUSignee(NULL)
83 ,fLSignee(NULL)
84 ,fUSignAngle(NULL)
85 ,fLSignAngle(NULL)
86{
87 //
88 // Constructor
89 //
90 fPIDBackground = new AliHFEpid("hfePidBackground");
91 fPIDBackgroundQA = new AliHFEpidQAmanager;
92}
93
94//________________________________________________________________________
95AliHFENonPhotonicElectron::AliHFENonPhotonicElectron()
96 :TNamed()
97 ,fMCEvent(0)
98 ,fAODArrayMCInfo(NULL)
99 ,fHFEBackgroundCuts(NULL)
100 ,fPIDBackground(0x0)
101 ,fPIDBackgroundQA(0)
102 ,fAlgorithmMA(kTRUE)
103 ,fUseFilterAOD(kTRUE)
104 ,fFilter(-1)
105 ,fChi2OverNDFCut(3.0)
106 ,fMaxDCA(3.0)
107 ,fMaxOpeningTheta(0.02)
108 ,fMaxOpeningPhi(0.1)
109 ,fMaxOpening3D(TMath::TwoPi())
110 ,fMaxInvMass(0.6)
111 ,fSetMassConstraint(kFALSE)
112 ,fArraytrack(NULL)
113 ,fCounterPoolBackground(0)
114 ,fListOutput(NULL)
115 ,fMCSourceee(NULL)
116 ,fUSignee(NULL)
117 ,fLSignee(NULL)
118 ,fUSignAngle(NULL)
119 ,fLSignAngle(NULL)
120{
121 //
122 // Constructor
123 //
124 fPIDBackground = new AliHFEpid("hfePidBackground");
125 fPIDBackgroundQA = new AliHFEpidQAmanager;
126}
127
128//________________________________________________________________________
129AliHFENonPhotonicElectron::AliHFENonPhotonicElectron(const AliHFENonPhotonicElectron &ref)
130 :TNamed(ref)
131 ,fMCEvent(NULL)
132 ,fAODArrayMCInfo(NULL)
133 ,fHFEBackgroundCuts(ref.fHFEBackgroundCuts)
134 ,fPIDBackground(ref.fPIDBackground)
135 ,fPIDBackgroundQA(ref.fPIDBackgroundQA)
136 ,fAlgorithmMA(ref.fAlgorithmMA)
137 ,fUseFilterAOD(ref.fUseFilterAOD)
138 ,fFilter(ref.fFilter)
139 ,fChi2OverNDFCut(ref.fChi2OverNDFCut)
140 ,fMaxDCA(ref.fMaxDCA)
141 ,fMaxOpeningTheta(ref.fMaxOpeningTheta)
142 ,fMaxOpeningPhi(ref.fMaxOpeningPhi)
143 ,fMaxOpening3D(ref.fMaxOpening3D)
144 ,fMaxInvMass(ref.fMaxInvMass)
145 ,fSetMassConstraint(ref.fSetMassConstraint)
146 ,fArraytrack(NULL)
147 ,fCounterPoolBackground(0)
148 ,fListOutput(ref.fListOutput)
149 ,fMCSourceee(ref.fMCSourceee)
150 ,fUSignee(ref.fUSignee)
151 ,fLSignee(ref.fLSignee)
152 ,fUSignAngle(ref.fUSignAngle)
153 ,fLSignAngle(ref.fLSignAngle)
154{
155 //
156 // Copy Constructor
157 //
158 ref.Copy(*this);
159}
160
161//____________________________________________________________
162AliHFENonPhotonicElectron &AliHFENonPhotonicElectron::operator=(const AliHFENonPhotonicElectron &ref){
163 //
164 // Assignment operator
165 //
166 if(this == &ref) ref.Copy(*this);
167 return *this;
168}
169
170//_________________________________________
171AliHFENonPhotonicElectron::~AliHFENonPhotonicElectron()
172{
173 //
174 // Destructor
175 //
176 if(fArraytrack) delete fArraytrack;
177 if(fHFEBackgroundCuts) delete fHFEBackgroundCuts;
178 if(fPIDBackground) delete fPIDBackground;
179}
180
181//_____________________________________________________________________________________________
182void AliHFENonPhotonicElectron::Init()
183{
184 //
185 // Init
186 //
187
188 if(!fListOutput) fListOutput = new TList;
189 fListOutput->SetName("HFENonPhotonicElectron");
190 fListOutput->SetOwner();
191
192 if(!fHFEBackgroundCuts) fHFEBackgroundCuts = new AliESDtrackCuts();
193
194 // Initialize PID
195 if(!fPIDBackground) fPIDBackground = new AliHFEpid("default pid");
196 if(fMCEvent || fAODArrayMCInfo) fPIDBackground->SetHasMCData(kTRUE);
197 if(!fPIDBackground->GetNumberOfPIDdetectors()) fPIDBackground->AddDetector("TPC", 0);
198 AliInfo("PID Background QA switched on");
199 fPIDBackgroundQA->Initialize(fPIDBackground);
200 fListOutput->Add(fPIDBackgroundQA->MakeList("HFENP_PID_Background"));
201 fPIDBackground->SortDetectors();
202
203 Int_t nBinsPt = 24;
204 Double_t binLimPt[25] = {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, 3., 3.5, 4., 5., 6.};
205
206 Int_t nBinsC = 11;
207 Double_t minC = 0.0;
208 Double_t maxC = 11.0;
209 Double_t binLimC[nBinsC+1];
210 for(Int_t i=0; i<=nBinsC; i++) binLimC[i]=(Double_t)minC + (maxC-minC)/nBinsC*(Double_t)i ;
211
212 Int_t nBinsSource = 10;
213 Double_t minSource = 0.;
214 Double_t maxSource = 10.;
215 Double_t binLimSource[nBinsSource+1];
216 for(Int_t i=0; i<=nBinsSource; i++) binLimSource[i]=(Double_t)minSource + (maxSource-minSource)/nBinsSource*(Double_t)i ;
217
218 Int_t nBinsInvMass = 60;
219 Double_t minInvMass = 0.;
220 Double_t maxInvMass = 0.6;
221 Double_t binLimInvMass[nBinsInvMass+1];
222 for(Int_t i=0; i<=nBinsInvMass; i++) binLimInvMass[i]=(Double_t)minInvMass + (maxInvMass-minInvMass)/nBinsInvMass*(Double_t)i ;
223
224 Int_t nBinsPhi = 180;
225 Double_t minPhi = 0.0;
226 Double_t maxPhi = TMath::Pi();
227 Double_t binLimPhi[nBinsPhi+1];
228 for(Int_t i=0; i<=nBinsPhi; i++)
229 {
230 binLimPhi[i]=(Double_t)minPhi + (maxPhi-minPhi)/nBinsPhi*(Double_t)i ;
231 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
232 }
233
234 Int_t nBinsAngle = 180;
235 Double_t minAngle = 0.0;
236 Double_t maxAngle = TMath::Pi();
237 Double_t binLimAngle[nBinsAngle+1];
238 for(Int_t i=0; i<=nBinsAngle; i++)
239 {
240 binLimAngle[i]=(Double_t)minAngle + (maxAngle-minAngle)/nBinsAngle*(Double_t)i ;
241 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
242 }
243
244 const Int_t nDimMCSource=3;
245 Int_t nBinMCSource[nDimMCSource] = {nBinsC,nBinsPt,nBinsSource};
246 fMCSourceee = new THnSparseF("fMCSourceee","fMCSourceee",nDimMCSource,nBinMCSource);
247 fMCSourceee->SetBinEdges(0,binLimC);
248 fMCSourceee->SetBinEdges(1,binLimPt);
249 fMCSourceee->SetBinEdges(2,binLimSource);
250/* fMCSourceee->GetAxis(0)->SetTitle("Multiplicity (arb.unit)");
251 fMCSourceee->GetAxis(1)->SetTitle("#it{p}_{mc,T} (GeV/c)");
252 fMCSourceee->GetAxis(2)->SetTitle("MC");*/
253 fMCSourceee->Sumw2();
254 AliDebug(2,"AliHFENonPhotonicElectron: fMCSourceee");
255
256 // ee invariant mass Unlike Sign
257 const Int_t nDimUSign=6;
258 Int_t nBinUSign[nDimUSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource,nBinsAngle};
259 fUSignee = new THnSparseF("fUSignee","fUSignee",nDimUSign,nBinUSign);
260 fUSignee->SetBinEdges(0,binLimPhi);
261 fUSignee->SetBinEdges(1,binLimC);
262 fUSignee->SetBinEdges(2,binLimPt);
263 fUSignee->SetBinEdges(3,binLimInvMass);
264 fUSignee->SetBinEdges(4,binLimSource);
265 fUSignee->SetBinEdges(5,binLimAngle);
266/* fUSignee->GetAxis(0)->SetTitle("#Delta #phi");
267 fUSignee->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
268 fUSignee->GetAxis(2)->SetTitle("#it{p}_{T} (GeV/c)");
269 fUSignee->GetAxis(3)->SetTitle("m_{e^{#pm}e^{#mp}}} (GeV/#it{c}^{2}");
270 fUSignee->GetAxis(4)->SetTitle("MC");
271 fUSignee->GetAxis(5)->SetTitle("#alpha");*/
272 fUSignee->Sumw2();
273 AliDebug(2,"AliHFENonPhotonicElectron: fUSignee");
274
275 // ee invariant mass Like Sign
276 const Int_t nDimLSign=6;
277 Int_t nBinLSign[nDimLSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource,nBinsAngle};
278 fLSignee = new THnSparseF("fLSignee","fLSignee",nDimLSign,nBinLSign);
279 fLSignee->SetBinEdges(0,binLimPhi);
280 fLSignee->SetBinEdges(1,binLimC);
281 fLSignee->SetBinEdges(2,binLimPt);
282 fLSignee->SetBinEdges(3,binLimInvMass);
283 fLSignee->SetBinEdges(4,binLimSource);
284 fLSignee->SetBinEdges(5,binLimAngle);
285/* fLSignee->GetAxis(0)->SetTitle("#Delta #phi");
286 fLSignee->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
287 fLSignee->GetAxis(2)->SetTitle("#it{p}_{T} (GeV/c)");
288 fLSignee->GetAxis(3)->SetTitle("m_{e^{#pm}e^{#pm}}} (GeV/#it{c}^{2}");
289 fLSignee->GetAxis(4)->SetTitle("MC");
290 fLSignee->GetAxis(5)->SetTitle("#alpha");*/
291 fLSignee->Sumw2();
292 AliDebug(2,"AliHFENonPhotonicElectron: fLSignee");
293
294 // ee angle Unlike Sign
295 const Int_t nDimUSignAngle=3;
296 Int_t nBinUSignAngle[nDimUSignAngle] = {nBinsAngle,nBinsC,nBinsSource};
297 fUSignAngle = new THnSparseF("fUSignAngle","fUSignAngle",nDimUSignAngle,nBinUSignAngle);
298 fUSignAngle->SetBinEdges(0,binLimAngle);
299 fUSignAngle->SetBinEdges(1,binLimC);
300 fUSignAngle->SetBinEdges(2,binLimSource);
301/* fUSignAngle->GetAxis(0)->SetTitle("#alpha");
302 fUSignAngle->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
303 fUSignAngle->GetAxis(2)->SetTitle("MC");*/
304 fUSignAngle->Sumw2();
305 AliDebug(2,"AliHFENonPhotonicElectron: fUSignAngle");
306
307 // ee angle Like Sign
308 const Int_t nDimLSignAngle=3;
309 Int_t nBinLSignAngle[nDimLSignAngle] = {nBinsAngle,nBinsC,nBinsSource};
310 fLSignAngle = new THnSparseF("fLSignAngle","fLSignAngle",nDimLSignAngle,nBinLSignAngle);
311 fLSignAngle->SetBinEdges(0,binLimAngle);
312 fLSignAngle->SetBinEdges(1,binLimC);
313 fLSignAngle->SetBinEdges(2,binLimSource);
314/* fLSignAngle->GetAxis(0)->SetTitle("#alpha");
315 fLSignAngle->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
316 fLSignAngle->GetAxis(2)->SetTitle("MC");*/
317 fLSignAngle->Sumw2();
318 AliDebug(2,"AliHFENonPhotonicElectron: fLSignAngle");
319
320
321 fListOutput->Add(fMCSourceee);
322 fListOutput->Add(fUSignee);
323 fListOutput->Add(fLSignee);
324 fListOutput->Add(fUSignAngle);
325 fListOutput->Add(fLSignAngle);
326
327}
328
329//_____________________________________________________________________________________________
663a2523 330void AliHFENonPhotonicElectron::InitRun(const AliVEvent *inputEvent,const AliPIDResponse *pidResponse)
76d0b522 331{
332 //
333 // Init run
334 //
335
336 if(!pidResponse)
337 {
338 AliDebug(1, "Using default PID Response");
339 Bool_t hasmc = kFALSE;
340 if(fMCEvent || fAODArrayMCInfo) hasmc=kTRUE;
341 pidResponse = AliHFEtools::GetDefaultPID(hasmc, inputEvent->IsA() == AliESDEvent::Class());
342 }
343
344 if(!fPIDBackground) return;
345 fPIDBackground->SetPIDResponse(pidResponse);
346
347 if(!fPIDBackground->IsInitialized())
348 {
349 // Initialize PID with the given run number
350 fPIDBackground->InitializePID(inputEvent->GetRunNumber());
351 }
352
353}
354
355//_____________________________________________________________________________________________
356Int_t AliHFENonPhotonicElectron::FillPoolAssociatedTracks(AliVEvent *inputEvent,Int_t binct)
357{
358 //
359 // Fill the pool of associated tracks
360 // Return the number of associated tracks
361 //
362
363 Int_t nbtracks = inputEvent->GetNumberOfTracks();
364
365 if( fArraytrack )
366 {
367 fArraytrack->~TArrayI();
368 new(fArraytrack) TArrayI(nbtracks);
369 }
370 else
371 {
372 fArraytrack = new TArrayI(nbtracks);
373 }
374
375 fCounterPoolBackground = 0;
376
377 for(Int_t k = 0; k < nbtracks; k++)
378 {
379 AliVTrack *track = (AliVTrack *) inputEvent->GetTrack(k);
380 if(!track) continue;
381
382 // Track cuts
383 Bool_t survivedbackground = kTRUE;
384 AliAODEvent *aodeventu = dynamic_cast<AliAODEvent *>(inputEvent);
385
386 if(aodeventu)
387 {
388 // AOD
389 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
390 if(aodtrack)
391 {
392 // filter
393 if(fUseFilterAOD)
394 {
395 if(!(aodtrack->TestFilterBit(fFilter))) survivedbackground = kFALSE;
396 }
397
398 // additional cuts
399 if(survivedbackground)
400 {
401 AliESDtrack esdTrack(aodtrack);
402
403 // set the TPC cluster info
404 esdTrack.SetTPCClusterMap(aodtrack->GetTPCClusterMap());
405 esdTrack.SetTPCSharedMap(aodtrack->GetTPCSharedMap());
406 esdTrack.SetTPCPointsF(aodtrack->GetTPCNclsF());
407 AliAODVertex *vAOD = aodeventu->GetPrimaryVertex();
408 Double_t bfield = aodeventu->GetMagneticField();
409 Double_t pos[3],cov[6];
410 vAOD->GetXYZ(pos);
411 vAOD->GetCovarianceMatrix(cov);
412 const AliESDVertex vESD(pos,cov,100.,100);
413 esdTrack.RelateToVertex(&vESD,bfield,3.);
414
415 if(!fHFEBackgroundCuts->IsSelected(&esdTrack))
416 {
417 survivedbackground = kFALSE;
418 }
419 }
420 }
421 }
422 else
423 {
424 // ESD
425 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
426 if(esdtrack)
427 {
428 if(!fHFEBackgroundCuts->IsSelected(esdtrack)) survivedbackground = kFALSE;
429 }
430 }
431
432 // PID
433 if(survivedbackground)
434 {
435 // PID track cuts
436 AliHFEpidObject hfetrack2;
437
438 if(!aodeventu) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);
439 else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);
440
441 hfetrack2.SetRecTrack(track);
442 if(binct>-1)
443 {
444 hfetrack2.SetCentrality((Int_t)binct);
445 AliDebug(2,Form("centrality %d and %d",binct,hfetrack2.GetCentrality()));
446 hfetrack2.SetPbPb();
447 }
448
449 if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundQA))
450 {
451 fArraytrack->AddAt(k,fCounterPoolBackground);
452 fCounterPoolBackground++;
453 AliDebug(2,Form("fCounterPoolBackground %d, track %d",fCounterPoolBackground,k));
454 }
455 }
456 } // loop tracks
457
458 //printf(Form("Test pool: nbtrack %d, binct %d \n fCounterPoolBackground %d \n",nbtracks,binct,fCounterPoolBackground));
459
460 return fCounterPoolBackground;
461}
462
463//_____________________________________________________________________________________________
464Int_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)
465{
466 //
467 // Look At Non HFE
468 //
469
470 /***********************************************************************************
471 * *
472 * iTrack1: index of the tagged electrons in AliVEvent *
473 * track1: tagged electron *
474 * vEvent: event *
475 * weight: weight in pt if not realistic *
476 * binct: centrality bin *
477 * deltaphi: phi-phi event plane for v2 *
478 * source: MC sources *
479 * indexmother: MC index mother *
480 * *
481 * *
482 * return -1 if nothing *
483 * return 2 if opposite charge within the mass range *
484 * return 4 if like charge within the mass range *
485 * return 6 if opposite & like charge within the mass range *
486 * *
487 ***********************************************************************************/
488
489 AliAODEvent *aodeventu = dynamic_cast<AliAODEvent*>(vEvent);
490 Int_t taggedphotonic = -1;
491
492 AliDebug(2,Form("fCounterPoolBackground %d in LookAtNonHFE!!!",fCounterPoolBackground));
493 if(!fArraytrack) return taggedphotonic;
494 AliDebug(2,Form("process track %d",iTrack1));
495
496 //Set Fill-Arrays for THnSparse
497 Double_t value_MCSource[3] = { binct, track1->Pt(), source}; //Centrality Pt Source
498 Double_t valueee[6] = { deltaphi, binct, track1->Pt(), -1, source,-1}; //DeltaPhi Centrality Pt InvariantMass Source
499 Double_t valueAngle[3] = { -1, binct, source}; //Angle Centrality Source
500
501 Bool_t USignPhotonic = kFALSE;
502 Bool_t LSignPhotonic = kFALSE;
503 Bool_t hasdcaT1 = kFALSE;
504 Bool_t hasdcaT2 = kFALSE;
505
506 Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()));
507 Int_t pdg2 = -100;
508 Int_t numberfound = 0;
509 Int_t iTrack2 = 0;
510 Int_t source2 = 0;
511 Int_t indexmother2 = -1;
512 Int_t fPDGtrack1 = 0;
513 Int_t fPDGtrack2 = 0;
514
515 Double_t eMass = TDatabasePDG::Instance()->GetParticle(11)->Mass(); //Electron mass in GeV
516 Double_t bfield = vEvent->GetMagneticField();
517 Double_t xt1 = 0; //radial position track 1 at the DCA point
518 Double_t xt2 = 0; //radial position track 2 at the DCA point
519 Double_t p1[3] = {0,0,0};
520 Double_t p2[3] = {0,0,0};
521 Double_t dca12 = 0;
522
523 Double_t angleESD = 0;
524 Double_t invmassESD = 0;
525
526 Double_t chi2OverNDF = 0;
527 Double_t width = 0;
528 Double_t angleAOD = 0;
529 Double_t invmassAOD = 0;
530
531 Float_t fCharge1 = 0;
532 Float_t fCharge2 = 0;
533
534 TLorentzVector electron1;
535 TLorentzVector electron2;
536 TLorentzVector mother;
537
538 AliVTrack *track2;
539 AliESDtrack *esdtrack1;
540 AliESDtrack *esdtrack2;
541 AliKFParticle *ktrack1;
542 AliKFParticle *ktrack2;
543 AliKFVertex primV(*(vEvent->GetPrimaryVertex()));
544
545 // FILL MCsource TODO: correct?
546 fMCSourceee->Fill(&value_MCSource[0],weight);
547
548
549 for(Int_t idex = 0; idex < fCounterPoolBackground; idex++)
550 {
551 iTrack2 = fArraytrack->At(idex);
552 AliDebug(2,Form("track %d",iTrack2));
553 track2 = (AliVTrack *)vEvent->GetTrack(iTrack2);
554
555 if(!track2)
556 {
557 printf("ERROR: Could not receive track %d", iTrack2);
558 continue;
559 }
560
561 if(iTrack2==iTrack1) continue;
562 AliDebug(2,"Different");
563
564 fCharge1 = track1->Charge(); //Charge from track1
565 fCharge2 = track2->Charge(); //Charge from track2
566
567 // Reset the MC info
568 valueAngle[2] = source;
569 valueee[4] = source;
570
571 // track cuts and PID already done
572
573 // if MC look
574 if(fMCEvent || fAODArrayMCInfo)
575 {
576 source2 = FindMother(TMath::Abs(track2->GetLabel()), indexmother2);
577 pdg2 = CheckPdg(TMath::Abs(track2->GetLabel()));
578
579 if(source2 >=0 )
580 {
581 if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0))
582 {
583 if(source == kElectronfromconversion)
584 {
585 valueAngle[2] = kElectronfromconversionboth;
586 valueee[4] = kElectronfromconversionboth;
587 numberfound++;
588 }
589
590 if(source == kElectronfrompi0)
591 {
592 valueAngle[2] = kElectronfrompi0both;
593 valueee[4] = kElectronfrompi0both;
594 }
595
596 if(source == kElectronfrometa)
597 {
598 valueAngle[2] = kElectronfrometaboth;
599 valueee[4] = kElectronfrometaboth;
600 }
601 }
602 }
603 }
604
605 if(fAlgorithmMA && (!aodeventu))
606 {
607 /** *
608 * ESD-Analysis *
609 ** */
610
611 esdtrack1 = dynamic_cast<AliESDtrack *>(track1); //ESD-track1
612 esdtrack2 = dynamic_cast<AliESDtrack *>(track2); //ESD-track2
613 if((!esdtrack1) || (!esdtrack2)) continue;
614
615 dca12 = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1); //DCA track1-track2
616
617 if(dca12 > fMaxDCA) continue; //! Cut on DCA
618
619 //Momento of the track extrapolated to DCA track-track
620 hasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1); //Track1
621 hasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2); //Track2
622 if(!hasdcaT1 || !hasdcaT2) AliWarning("It could be a problem in the extrapolation");
623
624 electron1.SetXYZM(esdtrack1->Px(), esdtrack1->Py(), esdtrack1->Pz(), eMass);
625 electron2.SetXYZM(esdtrack2->Px(), esdtrack2->Py(), esdtrack2->Pz(), eMass);
626
627 mother = electron1 + electron2;
628 invmassESD = mother.M();
629 angleESD = TVector2::Phi_0_2pi(electron1.Angle(electron2.Vect()));
630
631 valueAngle[0] = angleESD;
632 valueee[5] = angleESD;
633
634 if((fCharge1*fCharge2)>0.0) fLSignAngle->Fill(&valueAngle[0],weight);
635 else fUSignAngle->Fill(&valueAngle[0],weight);
636
637 if(angleESD > fMaxOpening3D) continue; //! Cut on Opening Angle
638 valueee[3] = invmassESD;
639
640 if((fCharge1*fCharge2)>0.0) fLSignee->Fill(&valueee[0],weight);
641 else fUSignee->Fill(&valueee[0],weight);
642
643 if(invmassESD > fMaxInvMass) continue; //! Cut on Invariant Mass
644
645 if((fCharge1*fCharge2)>0.0) LSignPhotonic=kTRUE;
646 else USignPhotonic=kTRUE;
647 }
648 else
649 {
650 /** *
651 * AliKF-Analysis *
652 ** */
653
654 fPDGtrack1 = 11;
655 fPDGtrack2 = 11;
656
657 if(fCharge1>0) fPDGtrack1 = -11;
658 if(fCharge2>0) fPDGtrack2 = -11;
659
660 ktrack1 = new AliKFParticle(*track1, fPDGtrack1);
661 ktrack2 = new AliKFParticle(*track2, fPDGtrack2);
662 AliKFParticle recoGamma(*ktrack1,*ktrack2);
663
664 if(recoGamma.GetNDF()<1) continue; //! Cut on Reconstruction
665
666 chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();
667 if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;
668
669 // DCA
670 //Double_t dca12 = ktrack1.GetDistanceFromParticle(ktrack2);
671 //if(dca12 > fMaxDCA) continue;
672
673 // if set mass constraint
674 if(fSetMassConstraint) //&& pVtx)
675 {
676 primV += recoGamma;
677 primV -= *ktrack1;
678 primV -= *ktrack2;
679 recoGamma.SetProductionVertex(primV);
680 recoGamma.SetMassConstraint(0,0.0001);
681 }
682
683 recoGamma.GetMass(invmassAOD,width);
684 angleAOD = ktrack1->GetAngle(*ktrack2);
685
686 valueAngle[0] = angleAOD;
687 valueee[5] = angleAOD;
688
689 if((fCharge1*fCharge2)>0.0) fLSignAngle->Fill(&valueAngle[0],weight);
690 else fUSignAngle->Fill(&valueAngle[0],weight);
691
692 if(angleAOD > fMaxOpening3D) continue; //! Cut on Opening Angle
693
694 valueee[3] = invmassAOD;
695
696 if((fCharge1*fCharge2)>0.0) fLSignee->Fill(&valueee[0],weight);
697 else fUSignee->Fill(&valueee[0],weight);
698
699 if(invmassAOD > fMaxInvMass) continue; //! Cut on Invariant Mass
700
701 if((fCharge1*fCharge2)>0.0) LSignPhotonic=kTRUE;
702 else USignPhotonic=kTRUE;
703 }
704 }
705
706 if( USignPhotonic && LSignPhotonic) taggedphotonic = 6;
707 if(!USignPhotonic && LSignPhotonic) taggedphotonic = 4;
708 if( USignPhotonic && !LSignPhotonic) taggedphotonic = 2;
709
710 return taggedphotonic;
711}
712
713//_________________________________________________________________________
714Int_t AliHFENonPhotonicElectron::FindMother(Int_t tr, Int_t &indexmother){
715 //
716 // Find the mother if MC
717 //
718
719 if(!fMCEvent && !fAODArrayMCInfo) return 0;
720
721 Int_t pdg = CheckPdg(tr);
722 if(TMath::Abs(pdg)!= 11)
723 {
724 indexmother = -1;
725 return kNoElectron;
726 }
727
728 indexmother = IsMotherGamma(tr);
729 if(indexmother > 0) return kElectronfromconversion;
730 indexmother = IsMotherPi0(tr);
731 if(indexmother > 0) return kElectronfrompi0;
732 indexmother = IsMotherC(tr);
733 if(indexmother > 0) return kElectronfromC;
734 indexmother = IsMotherB(tr);
735 if(indexmother > 0) return kElectronfromB;
736 indexmother = IsMotherEta(tr);
737 if(indexmother > 0) return kElectronfrometa;
738
739 return kElectronfromother;
740}
741
742//________________________________________________________________________________________________
743Int_t AliHFENonPhotonicElectron::CheckPdg(Int_t tr) {
744
745 //
746 // Return the pdg of the particle
747 //
748
749 Int_t pdgcode = -1;
750 if(tr < 0) return pdgcode;
751
752 if(fMCEvent)
753 {
754 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
755 if(!mctrack) return -1;
756 AliMCParticle *mctrackesd = NULL;
757 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;
758 pdgcode = mctrackesd->PdgCode();
759 }
760
761 if(fAODArrayMCInfo)
762 {
763 if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
764 AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
765 if(!mctrackaod) return pdgcode;
766 pdgcode = mctrackaod->GetPdgCode();
767 }
768
769 return pdgcode;
770}
771
772//_______________________________________________________________________________________________
773Int_t AliHFENonPhotonicElectron::IsMotherGamma(Int_t tr) {
774
775 //
776 // Return the lab of gamma mother or -1 if not gamma
777 //
778
779 if(tr < 0) return -1;
780
781 if(fMCEvent)
782 {
783 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
784 if(!mctrack) return -1;
785 AliMCParticle *mctrackesd = NULL;
786 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
787 TParticle *particle = 0x0;
788 particle = mctrackesd->Particle();
789
790 // Take mother
791 if(!particle) return -1;
792 Int_t imother = particle->GetFirstMother();
793 if(imother < 0) return -1;
794 AliMCParticle *mothertrack = NULL;
795 if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
796 TParticle * mother = mothertrack->Particle();
797 if(!mother) return -1;
798
799 // Check gamma
800 Int_t pdg = mother->GetPdgCode();
801 if(TMath::Abs(pdg) == 22) return imother;
802 if(TMath::Abs(pdg) == 11)
803 {
804 return IsMotherGamma(imother);
805 }
806
807 return -1;
808 }
809
810 if(fAODArrayMCInfo)
811 {
812 if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
813 AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
814 if(!mctrackaod) return -1;
815
816 // Take mother
817 Int_t imother = mctrackaod->GetMother();
818 if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
819 AliAODMCParticle *mothertrack = NULL;
820 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
821
822 // Check gamma
823 Int_t pdg = mothertrack->GetPdgCode();
824 if(TMath::Abs(pdg) == 22) return imother;
825 if(TMath::Abs(pdg) == 11)
826 {
827 return IsMotherGamma(imother);
828 }
829
830 return -1;
831 }
832
833 return -1;
834}
835
836//________________________________________________________________________________________________
837Int_t AliHFENonPhotonicElectron::IsMotherPi0(Int_t tr) {
838
839 //
840 // Return the lab of pi0 mother or -1 if not pi0
841 //
842
843 if(tr < 0) return -1;
844
845 if(fMCEvent)
846 {
847 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
848 if(!mctrack) return -1;
849 AliMCParticle *mctrackesd = NULL;
850 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
851 TParticle *particle = 0x0;
852 particle = mctrackesd->Particle();
853 // Take mother
854 if(!particle) return -1;
855 Int_t imother = particle->GetFirstMother();
856 if(imother < 0) return -1;
857 AliMCParticle *mothertrack = NULL;
858 if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
859 TParticle * mother = mothertrack->Particle();
860 if(!mother) return -1;
861 // Check pi0
862 Int_t pdg = mother->GetPdgCode();
863 if(TMath::Abs(pdg) == 111) return imother;
864 if(TMath::Abs(pdg) == 11)
865 {
866 return IsMotherPi0(imother);
867 }
868
869 return -1;
870 }
871
872 if(fAODArrayMCInfo) {
873
874 if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
875 AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
876 if(!mctrackaod) return -1;
877
878 // Take mother
879 Int_t imother = mctrackaod->GetMother();
880 if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
881 AliAODMCParticle *mothertrack = NULL;
882 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
883 // Check pi0
884 Int_t pdg = mothertrack->GetPdgCode();
885 if(TMath::Abs(pdg) == 111) return imother;
886 if(TMath::Abs(pdg) == 11)
887 {
888 return IsMotherPi0(imother);
889 }
890
891 return -1;
892 }
893
894 return -1;
895}
896//________________________________________________________________________________________________
897Int_t AliHFENonPhotonicElectron::IsMotherC(Int_t tr) {
898
899 //
900 // Return the lab of signal mother or -1 if not from C
901 //
902
903 if(tr < 0) return -1;
904
905 if(fMCEvent)
906 {
907 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
908 if(!mctrack) return -1;
909 AliMCParticle *mctrackesd = NULL;
910 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
911 TParticle *particle = 0x0;
912 particle = mctrackesd->Particle();
913
914 // Take mother
915 if(!particle) return -1;
916 Int_t imother = particle->GetFirstMother();
917 if(imother < 0) return -1;
918 AliMCParticle *mothertrack = NULL;
919 if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
920 TParticle * mother = mothertrack->Particle();
921 if(!mother) return -1;
922
923 // Check C
924 Int_t pdg = mother->GetPdgCode();
925 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)) return imother;
926 if(TMath::Abs(pdg) == 11)
927 {
928 return IsMotherC(imother);
929 }
930
931 return -1;
932 }
933
934 if(fAODArrayMCInfo)
935 {
936 if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
937 AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
938 if(!mctrackaod) return -1;
939
940 // Take mother
941 Int_t imother = mctrackaod->GetMother();
942 if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
943 AliAODMCParticle *mothertrack = NULL;
944 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
945
946 // Check C
947 Int_t pdg = mothertrack->GetPdgCode();
948 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)) return imother;
949 if(TMath::Abs(pdg) == 11)
950 {
951 return IsMotherC(imother);
952 }
953
954 return -1;
955 }
956
957 return -1;
958}
959//_______________________________________________________________________________________________
960Int_t AliHFENonPhotonicElectron::IsMotherB(Int_t tr) {
961
962 //
963 // Return the lab of signal mother or -1 if not B
964 //
965
966 if(tr < 0) return -1;
967
968 if(fMCEvent)
969 {
970 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
971 if(!mctrack) return -1;
972 AliMCParticle *mctrackesd = NULL;
973 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
974 TParticle *particle = 0x0;
975 particle = mctrackesd->Particle();
976
977 // Take mother
978 if(!particle) return -1;
979 Int_t imother = particle->GetFirstMother();
980 if(imother < 0) return -1;
981 AliMCParticle *mothertrack = NULL;
982 if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
983 TParticle * mother = mothertrack->Particle();
984 if(!mother) return -1;
985
986 // Check B
987 Int_t pdg = mother->GetPdgCode();
988 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)) return imother;
989 if(TMath::Abs(pdg) == 11)
990 {
991 return IsMotherB(imother);
992 }
993
994 return -1;
995 }
996
997 if(fAODArrayMCInfo)
998 {
999 if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
1000 AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
1001 if(!mctrackaod) return -1;
1002
1003 // Take mother
1004 Int_t imother = mctrackaod->GetMother();
1005 if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
1006 AliAODMCParticle *mothertrack = NULL;
1007 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
1008 // Check B
1009 Int_t pdg = mothertrack->GetPdgCode();
1010 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)) return imother;
1011 if(TMath::Abs(pdg) == 11)
1012 {
1013 return IsMotherB(imother);
1014 }
1015
1016 return -1;
1017 }
1018
1019 return -1;
1020}
1021
1022//_______________________________________________________________________________________________
1023Int_t AliHFENonPhotonicElectron::IsMotherEta(Int_t tr) {
1024
1025 //
1026 // Return the lab of eta mother or -1 if not eta
1027 //
1028
1029 if(tr < 0) return -1;
1030
1031 if(fMCEvent)
1032 {
1033 AliVParticle *mctrack = fMCEvent->GetTrack(tr);
1034 if(!mctrack) return -1;
1035 AliMCParticle *mctrackesd = NULL;
1036 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
1037 TParticle *particle = 0x0;
1038 particle = mctrackesd->Particle();
1039
1040 // Take mother
1041 if(!particle) return -1;
1042 Int_t imother = particle->GetFirstMother();
1043 if(imother < 0) return -1;
1044 AliMCParticle *mothertrack = NULL;
1045 if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
1046 TParticle * mother = mothertrack->Particle();
1047 if(!mother) return -1;
1048
1049 // Check eta
1050 Int_t pdg = mother->GetPdgCode();
1051 if(TMath::Abs(pdg) == 221) return imother;
1052 if(TMath::Abs(pdg) == 11)
1053 {
1054 return IsMotherEta(imother);
1055 }
1056
1057 return -1;
1058 }
1059
1060 if(fAODArrayMCInfo)
1061 {
1062 if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
1063 AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
1064 if(!mctrackaod) return -1;
1065
1066 // Take mother
1067 Int_t imother = mctrackaod->GetMother();
1068 if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
1069 AliAODMCParticle *mothertrack = NULL;
1070 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
1071
1072 // Check eta
1073 Int_t pdg = mothertrack->GetPdgCode();
1074 if(TMath::Abs(pdg) == 221) return imother;
1075 if(TMath::Abs(pdg) == 11)
1076 {
1077 return IsMotherEta(imother);
1078 }
1079
1080 return -1;
1081 }
1082
1083 return -1;
1084}