]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliAnalysisTaskFlowTPCEMCalEP.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Yvonne for the TPC-TOF MB pPb analysis
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliAnalysisTaskFlowTPCEMCalEP.cxx
CommitLineData
1d08b6e8 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// Class for heavy-flavour electron v2 with EMCal triggered events
17// Author: Denise Godoy
18
19
20#include "TChain.h"
21#include "TTree.h"
22#include "TH2F.h"
23#include "TMath.h"
24#include "TCanvas.h"
25#include "THnSparse.h"
26#include "TLorentzVector.h"
27#include "TString.h"
28#include "TFile.h"
29
30#include "AliAnalysisTask.h"
31#include "AliAnalysisManager.h"
32
33#include "AliESDEvent.h"
34#include "AliESDHandler.h"
35#include "AliAODEvent.h"
36#include "AliAODHandler.h"
37
914042c2 38#include "AliAnalysisTaskFlowTPCEMCalEP.h"
1d08b6e8 39#include "TGeoGlobalMagField.h"
40#include "AliLog.h"
41#include "AliAnalysisTaskSE.h"
42#include "TRefArray.h"
43#include "TVector.h"
44#include "AliESDInputHandler.h"
45#include "AliESDpid.h"
46#include "AliESDtrackCuts.h"
47#include "AliPhysicsSelection.h"
48#include "AliESDCaloCluster.h"
49#include "AliAODCaloCluster.h"
50#include "AliEMCALRecoUtils.h"
51#include "AliEMCALGeometry.h"
52#include "AliGeomManager.h"
53#include "stdio.h"
54#include "TGeoManager.h"
55#include "iostream"
56#include "fstream"
57
58#include "AliEMCALTrack.h"
59#include "AliMagF.h"
60
61#include "AliKFParticle.h"
62#include "AliKFVertex.h"
63
3a94d6fe 64#include "AliMCEventHandler.h"
65#include "AliMCEvent.h"
66#include "AliMCParticle.h"
67#include "AliStack.h"
68
1d08b6e8 69#include "AliPID.h"
70#include "AliPIDResponse.h"
71#include "AliHFEcontainer.h"
72#include "AliHFEcuts.h"
73#include "AliHFEpid.h"
74#include "AliHFEpidBase.h"
75#include "AliHFEpidQAmanager.h"
76#include "AliHFEtools.h"
77#include "AliCFContainer.h"
78#include "AliCFManager.h"
79
80#include "AliEventplane.h"
81#include "AliCentrality.h"
82
914042c2 83ClassImp(AliAnalysisTaskFlowTPCEMCalEP)
1d08b6e8 84//________________________________________________________________________
914042c2 85AliAnalysisTaskFlowTPCEMCalEP::AliAnalysisTaskFlowTPCEMCalEP(const char *name)
1d08b6e8 86 : AliAnalysisTaskSE(name)
87 ,fESD(0)
3a94d6fe 88 ,fMC(0)
1d08b6e8 89 ,fOutputList(0)
90 ,fTrackCuts(0)
91 ,fCuts(0)
92 ,fIdentifiedAsOutInz(kFALSE)
93 ,fPassTheEventCut(kFALSE)
94 ,fRejectKinkMother(kFALSE)
3a94d6fe 95 ,fIsMC(kFALSE)
1d08b6e8 96 ,fVz(0.0)
97 ,fCFM(0)
98 ,fPID(0)
99 ,fPIDqa(0)
100 ,fOpeningAngleCut(0.1)
101 ,fInvmassCut(0.01)
102 ,fNoEvents(0)
103 ,fTrkpt(0)
104 ,fTrkEovPBef(0)
105 ,fTrkEovPAft(0)
106 ,fdEdxBef(0)
107 ,fdEdxAft(0)
108 ,fInvmassLS(0)
109 ,fInvmassULS(0)
110 ,fOpeningAngleLS(0)
111 ,fOpeningAngleULS(0)
112 ,fPhotoElecPt(0)
113 ,fSemiInclElecPt(0)
3a94d6fe 114 ,fMCphotoElecPt(0)
1d08b6e8 115 ,fTrackPtBefTrkCuts(0)
116 ,fTrackPtAftTrkCuts(0)
117 ,fTPCnsigma(0)
118 ,fCent(0)
44fdf7ac 119 ,fevPlaneV0A(0)
120 ,fevPlaneV0C(0)
121 ,fevPlaneTPC(0)
1d08b6e8 122 ,fTPCsubEPres(0)
123 ,fEPres(0)
124 ,fCorr(0)
125 ,feTPCV2(0)
126 ,feV2(0)
127 ,fphoteV2(0)
128 ,fChargPartV2(0)
6522b3b5 129 ,fGammaWeight(0)
130 ,fPi0Weight(0)
131 ,fEtaWeight(0)
def484dc 132 ,fD0Weight(0)
133 ,fDplusWeight(0)
134 ,fDminusWeight(0)
937b2ff2 135 ,fD0_e(0)
1c424962 136 ,fTot_pi0e(0)
137 ,fPhot_pi0e(0)
138 ,fPhotBCG_pi0e(0)
139 ,fTot_etae(0)
140 ,fPhot_etae(0)
141 ,fPhotBCG_etae(0)
1d08b6e8 142{
143 //Named constructor
937b2ff2 144
145 for(Int_t k = 0; k < 6; k++) {
146 fDe[k]= NULL;
147 fD0e[k]= NULL;
148 fDpluse[k]= NULL;
149 fDminuse[k]= NULL;
150 }
151
1d08b6e8 152 fPID = new AliHFEpid("hfePid");
153 fTrackCuts = new AliESDtrackCuts();
154
155 // Define input and output slots here
156 // Input slot #0 works with a TChain
157 DefineInput(0, TChain::Class());
158 // Output slot #0 id reserved by the base class for AOD
159 // Output slot #1 writes into a TH1 container
160 // DefineOutput(1, TH1I::Class());
161 DefineOutput(1, TList::Class());
162 // DefineOutput(3, TTree::Class());
163}
164
165//________________________________________________________________________
914042c2 166AliAnalysisTaskFlowTPCEMCalEP::AliAnalysisTaskFlowTPCEMCalEP()
1d08b6e8 167 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisElecHadCorrel")
168 ,fESD(0)
3a94d6fe 169 ,fMC(0)
1d08b6e8 170 ,fOutputList(0)
171 ,fTrackCuts(0)
172 ,fCuts(0)
173 ,fIdentifiedAsOutInz(kFALSE)
174 ,fPassTheEventCut(kFALSE)
175 ,fRejectKinkMother(kFALSE)
3a94d6fe 176 ,fIsMC(kFALSE)
1d08b6e8 177 ,fVz(0.0)
178 ,fCFM(0)
179 ,fPID(0)
180 ,fPIDqa(0)
181 ,fOpeningAngleCut(0.1)
182 ,fInvmassCut(0.01)
183 ,fNoEvents(0)
184 ,fTrkpt(0)
185 ,fTrkEovPBef(0)
186 ,fTrkEovPAft(0)
187 ,fdEdxBef(0)
188 ,fdEdxAft(0)
189 ,fInvmassLS(0)
190 ,fInvmassULS(0)
191 ,fOpeningAngleLS(0)
192 ,fOpeningAngleULS(0)
193 ,fPhotoElecPt(0)
194 ,fSemiInclElecPt(0)
3a94d6fe 195 ,fMCphotoElecPt(0)
1d08b6e8 196 ,fTrackPtBefTrkCuts(0)
197 ,fTrackPtAftTrkCuts(0)
198 ,fTPCnsigma(0)
199 ,fCent(0)
44fdf7ac 200 ,fevPlaneV0A(0)
201 ,fevPlaneV0C(0)
202 ,fevPlaneTPC(0)
1d08b6e8 203 ,fTPCsubEPres(0)
204 ,fEPres(0)
205 ,fCorr(0)
206 ,feTPCV2(0)
207 ,feV2(0)
208 ,fphoteV2(0)
209 ,fChargPartV2(0)
6522b3b5 210 ,fGammaWeight(0)
211 ,fPi0Weight(0)
212 ,fEtaWeight(0)
def484dc 213 ,fD0Weight(0)
214 ,fDplusWeight(0)
215 ,fDminusWeight(0)
937b2ff2 216 ,fD0_e(0)
1c424962 217 ,fTot_pi0e(0)
218 ,fPhot_pi0e(0)
219 ,fPhotBCG_pi0e(0)
220 ,fTot_etae(0)
221 ,fPhot_etae(0)
222 ,fPhotBCG_etae(0)
1d08b6e8 223{
224 //Default constructor
937b2ff2 225
226 for(Int_t k = 0; k < 6; k++) {
227 fDe[k]= NULL;
228 fD0e[k]= NULL;
229 fDpluse[k]= NULL;
230 fDminuse[k]= NULL;
231 }
232
1d08b6e8 233 fPID = new AliHFEpid("hfePid");
234
235 fTrackCuts = new AliESDtrackCuts();
236
237 // Constructor
238 // Define input and output slots here
239 // Input slot #0 works with a TChain
240 DefineInput(0, TChain::Class());
241 // Output slot #0 id reserved by the base class for AOD
242 // Output slot #1 writes into a TH1 container
243 // DefineOutput(1, TH1I::Class());
244 DefineOutput(1, TList::Class());
245 //DefineOutput(3, TTree::Class());
246}
247//_________________________________________
248
914042c2 249AliAnalysisTaskFlowTPCEMCalEP::~AliAnalysisTaskFlowTPCEMCalEP()
1d08b6e8 250{
251 //Destructor
252
253 delete fOutputList;
254 delete fPID;
255 delete fCFM;
256 delete fPIDqa;
257 delete fTrackCuts;
258}
259//_________________________________________
260
914042c2 261void AliAnalysisTaskFlowTPCEMCalEP::UserExec(Option_t*)
1d08b6e8 262{
263 //Main loop
264 //Called for each event
265
266 // create pointer to event
267 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
268 if (!fESD) {
269 printf("ERROR: fESD not available\n");
270 return;
271 }
272
273 if(!fCuts){
274 AliError("HFE cuts not available");
275 return;
276 }
277
278 if(!fPID->IsInitialized()){
279 // Initialize PID with the given run number
280 AliWarning("PID not initialised, get from Run no");
281 fPID->InitializePID(fESD->GetRunNumber());
282 }
283
3a94d6fe 284 if(fIsMC)fMC = MCEvent();
285 AliStack* stack = NULL;
286 if(fIsMC && fMC) stack = fMC->Stack();
1d08b6e8 287
288 Int_t fNOtrks = fESD->GetNumberOfTracks();
289 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
290
291 Double_t pVtxZ = -999;
292 pVtxZ = pVtx->GetZ();
293
294 if(TMath::Abs(pVtxZ)>10) return;
295 fNoEvents->Fill(0);
296
297 if(fNOtrks<2) return;
298
299 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
300 if(!pidResponse){
301 AliDebug(1, "Using default PID Response");
302 pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
303 }
304
305 fPID->SetPIDResponse(pidResponse);
306
307 fCFM->SetRecEventInfo(fESD);
308
309 Float_t cent = -1.;
310 AliCentrality *centrality = fESD->GetCentrality();
311 cent = centrality->GetCentralityPercentile("V0M");
312 fCent->Fill(cent);
6522b3b5 313
1d08b6e8 314 //Event planes
315
316 Double_t evPlaneV0A = TVector2::Phi_0_2pi(fESD->GetEventplane()->GetEventplane("V0A",fESD,2));
317 if(evPlaneV0A > TMath::Pi()) evPlaneV0A = evPlaneV0A - TMath::Pi();
5945f829 318 fevPlaneV0A->Fill(evPlaneV0A,cent);
1d08b6e8 319
320 Double_t evPlaneV0C = TVector2::Phi_0_2pi(fESD->GetEventplane()->GetEventplane("V0C",fESD,2));
321 if(evPlaneV0C > TMath::Pi()) evPlaneV0C = evPlaneV0C - TMath::Pi();
5945f829 322 fevPlaneV0C->Fill(evPlaneV0C,cent);
1d08b6e8 323
324 AliEventplane* esdTPCep = fESD->GetEventplane();
44fdf7ac 325 TVector2 *standardQ = 0x0;
1d08b6e8 326 Double_t qx = -999., qy = -999.;
44fdf7ac 327 standardQ = esdTPCep->GetQVector();
ee89e744 328 if(!standardQ)return;
329
330 qx = standardQ->X();
331 qy = standardQ->Y();
332
1d08b6e8 333 TVector2 qVectorfortrack;
334 qVectorfortrack.Set(qx,qy);
335 Float_t evPlaneTPC = TVector2::Phi_0_2pi(qVectorfortrack.Phi())/2.;
5945f829 336 fevPlaneTPC->Fill(evPlaneTPC,cent);
44fdf7ac 337
1d08b6e8 338 TVector2 *qsub1a = esdTPCep->GetQsub1();
339 TVector2 *qsub2a = esdTPCep->GetQsub2();
340 Double_t evPlaneResTPC = -999.;
341 if(qsub1a && qsub2a)
342 {
343 evPlaneResTPC = TMath::Cos(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));
344 }
345
346 fTPCsubEPres->Fill(evPlaneResTPC,cent);
347
348 Double_t evPlaneRes[4]={GetCos2DeltaPhi(evPlaneV0A,evPlaneV0C),GetCos2DeltaPhi(evPlaneV0A,evPlaneTPC),GetCos2DeltaPhi(evPlaneV0C,evPlaneTPC),cent};
349 fEPres->Fill(evPlaneRes);
350
1c424962 351 if(fIsMC && fMC && stack && cent>30 && cent<50){
6522b3b5 352 Int_t nParticles = stack->GetNtrack();
353 for (Int_t iParticle = 0; iParticle < nParticles; iParticle++) {
354 TParticle* particle = stack->Particle(iParticle);
355 int fPDG = particle->GetPdgCode();
356 double pTMC = particle->Pt();
357 double etaMC = particle->Eta();
358 if(fabs(etaMC)>0.7)continue;
359
360 Bool_t MChijing = fMC->IsFromBGEvent(iParticle);
361 int iHijing = 1;
362 if(!MChijing)iHijing = 0;
363
364 if(fPDG==111)fPi0Weight->Fill(pTMC,iHijing);//pi0
365 if(fPDG==221)fEtaWeight->Fill(pTMC,iHijing);//eta
def484dc 366 if(fPDG==421)fD0Weight->Fill(pTMC,iHijing);//D0
367 if(fPDG==411)fDplusWeight->Fill(pTMC,iHijing);//D+
368 if(fPDG==-411)fDminusWeight->Fill(pTMC,iHijing);//D-
369
370
6522b3b5 371 Int_t idMother = particle->GetFirstMother();
372 if (idMother>0){
373 TParticle *mother = stack->Particle(idMother);
374 int motherPDG = mother->GetPdgCode();
375 if(fPDG==22 && motherPDG!=111 && motherPDG!=221)fGammaWeight->Fill(pTMC,iHijing);//gamma
376 }
377
378 }
379 }
380
381
1d08b6e8 382 // Track loop
383 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
384 AliESDtrack* track = fESD->GetTrack(iTracks);
385 if (!track) {
386 printf("ERROR: Could not receive track %d\n", iTracks);
387 continue;
388 }
389
390 if(TMath::Abs(track->Eta())>0.7) continue;
391
392 fTrackPtBefTrkCuts->Fill(track->Pt());
b1a1875f 393
1d08b6e8 394 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
395
1d08b6e8 396 if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
397 if(track->GetKinkIndex(0) != 0) continue;
398 }
399
1d08b6e8 400 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
401
1d08b6e8 402 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
403
1d08b6e8 404 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
405
406 fTrackPtAftTrkCuts->Fill(track->Pt());
407
1c424962 408 Double_t clsE = -999., p = -999., EovP=-999., pt = -999., dEdx=-999., fTPCnSigma=0, phi=-999., wclsE = -999., wEovP = -999., m02= -999., m20= -999.;
32d22164 409
1d08b6e8 410 pt = track->Pt();
44fdf7ac 411 if(pt<2) continue;
1d08b6e8 412 fTrkpt->Fill(pt);
3a94d6fe 413
1d08b6e8 414 Int_t clsId = track->GetEMCALcluster();
415 if (clsId>0){
416 AliESDCaloCluster *cluster = fESD->GetCaloCluster(clsId);
417 if(cluster && cluster->IsEMCAL()){
418 clsE = cluster->E();
1c424962 419 m20 = cluster->GetM20();
420 m02 = cluster->GetM02();
1d08b6e8 421 }
422 }
b1a1875f 423
1d08b6e8 424 p = track->P();
425 phi = track->Phi();
426 dEdx = track->GetTPCsignal();
427 EovP = clsE/p;
b1a1875f 428 wEovP = wclsE/p;
1d08b6e8 429 fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
430 fdEdxBef->Fill(p,dEdx);
431 fTPCnsigma->Fill(p,fTPCnSigma);
32d22164 432
433 //Remove electron candidate from the event plane
6522b3b5 434 Float_t evPlaneCorrTPC = evPlaneTPC;
32d22164 435 if(dEdx>70 && dEdx<90){
436 Double_t qX = standardQ->X() - esdTPCep->GetQContributionX(track);
437 Double_t qY = standardQ->Y() - esdTPCep->GetQContributionY(track);
438 TVector2 newQVectorfortrack;
439 newQVectorfortrack.Set(qX,qY);
440 evPlaneCorrTPC = TVector2::Phi_0_2pi(newQVectorfortrack.Phi())/2;
441 }
442
06d74a7b 443 Bool_t fFlagPhotonicElec = kFALSE;
444 Bool_t fFlagPhotonicElecBCG = kFALSE;
5945f829 445
06d74a7b 446 SelectPhotonicElectron(iTracks,track, fFlagPhotonicElec, fFlagPhotonicElecBCG);
447
1c424962 448 Double_t corr[12]={phi,fTPCnSigma,cent,pt,EovP,GetDeltaPhi(phi,evPlaneCorrTPC),GetDeltaPhi(phi,evPlaneV0A),GetDeltaPhi(phi,evPlaneV0C),fFlagPhotonicElec,fFlagPhotonicElecBCG,m02,m20};
1d08b6e8 449 fCorr->Fill(corr);
06d74a7b 450
32d22164 451
452 Int_t whichFirstMother = 0, whichSecondMother = 0, whichThirdMother = 0;
6522b3b5 453 Int_t whichPart = -99;
32d22164 454 Int_t partPDG = -99, motherPDG = -99, secondMotherPDG = -99, thirdMotherPDG = -99;
1c424962 455 Double_t partPt = -99. , motherPt = -99., secondMotherPt = -99.,thirdMotherPt = -99.;
456 Double_t weight = 1.;
def484dc 457 Double_t Dweight = 1.;
32d22164 458 Bool_t MChijing;
459
1c424962 460 Bool_t pi0Decay= kFALSE;
461 Bool_t etaDecay= kFALSE;
462
32d22164 463 if(fIsMC && fMC && stack){
06d74a7b 464 Int_t label = track->GetLabel();
465 if(label>0){
466 TParticle *particle = stack->Particle(label);
467 if(particle){
32d22164 468 partPDG = particle->GetPdgCode();
469 partPt = particle->Pt();
470
6522b3b5 471 if (TMath::Abs(partPDG)==11) whichPart = 0; //electron
472 if (partPDG==22) whichPart = 3; //gamma
473 if (partPDG==111) whichPart = 2; //pi0
474 if (partPDG==221) whichPart = 1; //eta
475
32d22164 476 MChijing = fMC->IsFromBGEvent(label);
06d74a7b 477
32d22164 478 int iHijing = 1;
479 if(!MChijing) iHijing = 0; // 0 if enhanced sample
480
06d74a7b 481 Int_t idMother = particle->GetFirstMother();
482 if (idMother>0){
32d22164 483 TParticle *mother = stack->Particle(idMother);
484 motherPt = mother->Pt();
485 motherPDG = mother->GetPdgCode();
486
32d22164 487 if (motherPDG==22) whichFirstMother = 3; //gamma
488 if (motherPDG==111) whichFirstMother = 2; //pi0
489 if (motherPDG==221) whichFirstMother = 1; //eta
490
491 Int_t idSecondMother = particle->GetSecondMother();
492 if (idSecondMother>0){
493 TParticle *secondMother = stack->Particle(idSecondMother);
494 secondMotherPt = secondMother->Pt();
495 secondMotherPDG = secondMother->GetPdgCode();
496
497 if (secondMotherPDG==111) whichSecondMother = 2; //pi0
498 if (secondMotherPDG==221) whichSecondMother = 1; //eta
499
500 Int_t idThirdMother = secondMother->GetFirstMother();
501 if (idThirdMother>0){
502 TParticle *thirdMother = stack->Particle(idThirdMother);
503 thirdMotherPt = thirdMother->Pt();
504 thirdMotherPDG = thirdMother->GetPdgCode();
505
506 if (thirdMotherPDG==221) whichThirdMother = 1; //eta
507 }
508 }
def484dc 509
510
511 if (cent>30 && cent<50 && TMath::Abs(partPDG)==11){
937b2ff2 512
513 if (motherPDG==421 || TMath::Abs(motherPDG)==411){ // D
908103bf 514 Double_t phi_D = -99., Deltaphi_De=-99.;
515 phi_D = mother->Phi();
516 Deltaphi_De = phi_D - phi;
517
518 fD0_e->Fill(pt,Deltaphi_De);
519
937b2ff2 520 Dweight = GetDweight(0,motherPt);
521 if(iHijing==1) Dweight = 1.0;
522
523 if (motherPt>=2 && motherPt<3) fDe[0]->Fill(pt,Dweight);
524 if (motherPt>=3 && motherPt<4) fDe[1]->Fill(pt,Dweight);
525 if (motherPt>=4 && motherPt<6) fDe[2]->Fill(pt,Dweight);
526 if (motherPt>=6 && motherPt<8) fDe[3]->Fill(pt,Dweight);
527 if (motherPt>=8 && motherPt<12) fDe[4]->Fill(pt,Dweight);
528 if (motherPt>=12 && motherPt<16) fDe[5]->Fill(pt,Dweight);
529 }
530 if (motherPDG==421){ // D0
908103bf 531
937b2ff2 532 Dweight = GetDweight(1,motherPt);
533 if(iHijing==1) Dweight = 1.0;
534
535 if (motherPt>=2 && motherPt<3) fD0e[0]->Fill(pt,Dweight);
536 if (motherPt>=3 && motherPt<4) fD0e[1]->Fill(pt,Dweight);
537 if (motherPt>=4 && motherPt<6) fD0e[2]->Fill(pt,Dweight);
538 if (motherPt>=6 && motherPt<8) fD0e[3]->Fill(pt,Dweight);
539 if (motherPt>=8 && motherPt<12) fD0e[4]->Fill(pt,Dweight);
540 if (motherPt>=12 && motherPt<16) fD0e[5]->Fill(pt,Dweight);
541 }
542 if (motherPDG==411){ // D+
543 Dweight = GetDweight(2,motherPt);
544 if(iHijing==1) Dweight = 1.0;
545
546 if (motherPt>=2 && motherPt<3) fDpluse[0]->Fill(pt,Dweight);
547 if (motherPt>=3 && motherPt<4) fDpluse[1]->Fill(pt,Dweight);
548 if (motherPt>=4 && motherPt<6) fDpluse[2]->Fill(pt,Dweight);
549 if (motherPt>=6 && motherPt<8) fDpluse[3]->Fill(pt,Dweight);
550 if (motherPt>=8 && motherPt<12) fDpluse[4]->Fill(pt,Dweight);
551 if (motherPt>=12 && motherPt<16) fDpluse[5]->Fill(pt,Dweight);
552 }
553 if (motherPDG==-411){ //D-
554 Dweight = GetDweight(3,motherPt);
555 if(iHijing==1) Dweight = 1.0;
556
557 if (motherPt>=2 && motherPt<3) fDminuse[0]->Fill(pt,Dweight);
558 if (motherPt>=3 && motherPt<4) fDminuse[1]->Fill(pt,Dweight);
559 if (motherPt>=4 && motherPt<6) fDminuse[2]->Fill(pt,Dweight);
560 if (motherPt>=6 && motherPt<8) fDminuse[3]->Fill(pt,Dweight);
561 if (motherPt>=8 && motherPt<12) fDminuse[4]->Fill(pt,Dweight);
562 if (motherPt>=12 && motherPt<16) fDminuse[5]->Fill(pt,Dweight);
563 }
def484dc 564 }
565
1c424962 566 //pi0 decay
567 if (TMath::Abs(partPDG)==11 && motherPDG==111 && secondMotherPDG!=221){ //not eta -> pi0 -> e
568 weight = GetPi0weight(motherPt);
569 pi0Decay = kTRUE;
570 }
571 if (TMath::Abs(partPDG)==11 && motherPDG==22 && secondMotherPDG==111 && thirdMotherPDG!=221){ //not eta -> pi0 -> gamma -> e
572 weight = GetPi0weight(secondMotherPt);
573 pi0Decay = kTRUE;
574 }
575
576 //eta decay
577 if (TMath::Abs(partPDG)==11 && motherPDG==221){ //eta -> e
578 weight = GetEtaweight(motherPt);
579 etaDecay = kTRUE;
580 }
581 if (TMath::Abs(partPDG)==11 && motherPDG==111 && secondMotherPDG==221){ //eta -> pi0 -> e
582 weight = GetEtaweight(secondMotherPt);
583 etaDecay = kTRUE;
584 }
585 if (TMath::Abs(partPDG)==11 && motherPDG==22 && secondMotherPDG==221){ //eta -> gamma -> e
586 weight = GetEtaweight(secondMotherPt);
587 etaDecay = kTRUE;
588 }
589 if (TMath::Abs(partPDG)==11 && motherPDG==22 && secondMotherPDG==111 && thirdMotherPDG==221){ //eta -> pi0 -> gamma -> e
590 weight = GetEtaweight(thirdMotherPt);
591 etaDecay = kTRUE;
592 }
593
def484dc 594 if (cent>30 && cent<50 && fTPCnSigma>-1 && fTPCnSigma<3 && EovP>1 && EovP<1.3 && (motherPDG==22 || motherPDG==111 || motherPDG==221)){
1c424962 595 if(iHijing==1) weight = 1.0;
596
597 if (pi0Decay){
598 fTot_pi0e->Fill(partPt,weight);
599 if(fFlagPhotonicElec) fPhot_pi0e->Fill(partPt,weight);
600 if(fFlagPhotonicElecBCG) fPhotBCG_pi0e->Fill(partPt,weight);
601 }
602
603 if (etaDecay){
604 fTot_etae->Fill(partPt,weight);
605 if(fFlagPhotonicElec) fPhot_etae->Fill(partPt,weight);
606 if(fFlagPhotonicElecBCG) fPhotBCG_etae->Fill(partPt,weight);
607 }
608
609 }
610
611
6522b3b5 612 Double_t mc[15]={EovP,fTPCnSigma,partPt,fFlagPhotonicElec,fFlagPhotonicElecBCG,whichPart,cent,pt,whichFirstMother,whichSecondMother,whichThirdMother,iHijing,motherPt,secondMotherPt,thirdMotherPt};
1c424962 613
614// fMCphotoElecPt->Fill(mc);
615 if (motherPDG==22 || motherPDG==111 || motherPDG==221) fMCphotoElecPt->Fill(mc);// mother = gamma, pi0, eta
06d74a7b 616 }
617 }
618 }
619 }
06d74a7b 620
1d08b6e8 621
622 if(fTPCnSigma >= 1.5 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,EovP);
623 Int_t pidpassed = 1;
624
625 //--- track accepted
626 AliHFEpidObject hfetrack;
627 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
628 hfetrack.SetRecTrack(track);
629 hfetrack.SetPbPb();
630 if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
631
44fdf7ac 632 Double_t corrV2[7]={phi,cent,pt,EovP,GetCos2DeltaPhi(phi,evPlaneTPC),GetCos2DeltaPhi(phi,evPlaneV0A),GetCos2DeltaPhi(phi,evPlaneV0C)};
1d08b6e8 633 fChargPartV2->Fill(corrV2);
b3bd01b3 634
635 if(fTPCnSigma >= -0.5){
32d22164 636 Double_t correctedV2[5]={cent,pt,GetCos2DeltaPhi(phi,evPlaneCorrTPC),GetCos2DeltaPhi(phi,evPlaneV0A),GetCos2DeltaPhi(phi,evPlaneV0C)};
b3bd01b3 637 feTPCV2->Fill(correctedV2);
638 }
639
640 if(pidpassed==0) continue;
641
32d22164 642 Double_t correctedV2[5]={cent,pt,GetCos2DeltaPhi(phi,evPlaneCorrTPC),GetCos2DeltaPhi(phi,evPlaneV0A),GetCos2DeltaPhi(phi,evPlaneV0C)};
1d08b6e8 643
1d08b6e8 644 feV2->Fill(correctedV2);
645
646 fTrkEovPAft->Fill(pt,EovP);
647 fdEdxAft->Fill(p,dEdx);
648
3a94d6fe 649 if(fFlagPhotonicElec){
650 fphoteV2->Fill(correctedV2);
651 fPhotoElecPt->Fill(pt);
652 }
653
654 if(!fFlagPhotonicElec) fSemiInclElecPt->Fill(pt);
655
1d08b6e8 656 }
657 PostData(1, fOutputList);
658}
659//_________________________________________
914042c2 660void AliAnalysisTaskFlowTPCEMCalEP::UserCreateOutputObjects()
1d08b6e8 661{
3a94d6fe 662 //--- Check MC
663 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()){
664 fIsMC = kTRUE;
665 printf("+++++ MC Data available");
666 }
1d08b6e8 667 //--------Initialize PID
3a94d6fe 668 fPID->SetHasMCData(fIsMC);
669
1d08b6e8 670 if(!fPID->GetNumberOfPIDdetectors())
671 {
672 fPID->AddDetector("TPC", 0);
673 fPID->AddDetector("EMCAL", 1);
674 }
675
676 fPID->SortDetectors();
677 fPIDqa = new AliHFEpidQAmanager();
678 fPIDqa->Initialize(fPID);
679
680 //--------Initialize correction Framework and Cuts
681 fCFM = new AliCFManager;
682 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
683 fCFM->SetNStepParticle(kNcutSteps);
684 for(Int_t istep = 0; istep < kNcutSteps; istep++)
685 fCFM->SetParticleCutsList(istep, NULL);
686
687 if(!fCuts){
688 AliWarning("Cuts not available. Default cuts will be used");
689 fCuts = new AliHFEcuts;
690 fCuts->CreateStandardCuts();
691 }
692 fCuts->Initialize(fCFM);
693
694 //---------Output Tlist
695 fOutputList = new TList();
696 fOutputList->SetOwner();
697 fOutputList->Add(fPIDqa->MakeList("PIDQA"));
698
699 fNoEvents = new TH1F("fNoEvents","",1,0,1) ;
700 fOutputList->Add(fNoEvents);
701
702 fTrkpt = new TH1F("fTrkpt","track pt",100,0,50);
703 fOutputList->Add(fTrkpt);
704
705 fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50);
706 fOutputList->Add(fTrackPtBefTrkCuts);
707
708 fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50);
709 fOutputList->Add(fTrackPtAftTrkCuts);
710
711 fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10);
712 fOutputList->Add(fTPCnsigma);
713
714 fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2);
715 fOutputList->Add(fTrkEovPBef);
716
717 fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2);
718 fOutputList->Add(fTrkEovPAft);
719
720 fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,150,0,150);
721 fOutputList->Add(fdEdxBef);
722
723 fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,150,0,150);
724 fOutputList->Add(fdEdxAft);
725
726 fInvmassLS = new TH1F("fInvmassLS", "Inv mass of LS (e,e); mass(GeV/c^2); counts;", 500,0,0.5);
727 fOutputList->Add(fInvmassLS);
728
729 fInvmassULS = new TH1F("fInvmassULS", "Inv mass of ULS (e,e); mass(GeV/c^2); counts;", 500,0,0.5);
730 fOutputList->Add(fInvmassULS);
731
732 fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
733 fOutputList->Add(fOpeningAngleLS);
734
735 fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
736 fOutputList->Add(fOpeningAngleULS);
737
738 fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50);
739 fOutputList->Add(fPhotoElecPt);
740
741 fSemiInclElecPt = new TH1F("fSemiInclElecPt", "Semi-inclusive electron pt",100,0,50);
742 fOutputList->Add(fSemiInclElecPt);
743
744 fCent = new TH1F("fCent","Centrality",100,0,100) ;
745 fOutputList->Add(fCent);
746
5945f829 747 fevPlaneV0A = new TH2F("fevPlaneV0A","V0A EP",100,0,TMath::Pi(),90,0,90);
44fdf7ac 748 fOutputList->Add(fevPlaneV0A);
749
5945f829 750 fevPlaneV0C = new TH2F("fevPlaneV0C","V0C EP",100,0,TMath::Pi(),90,0,90);
44fdf7ac 751 fOutputList->Add(fevPlaneV0C);
752
5945f829 753 fevPlaneTPC = new TH2F("fevPlaneTPC","TPC EP",100,0,TMath::Pi(),90,0,90);
44fdf7ac 754 fOutputList->Add(fevPlaneTPC);
755
1d08b6e8 756 fTPCsubEPres = new TH2F("fTPCsubEPres","TPC subevent plane resolution",100,-1,1,90,0,90);
757 fOutputList->Add(fTPCsubEPres);
758
759 Int_t binsv1[4]={100,100,100,90}; // V0A-V0C, V0A-TPC, V0C-TPC, cent
760 Double_t xminv1[4]={-1,-1,-1,0};
761 Double_t xmaxv1[4]={1,1,1,90};
762 fEPres = new THnSparseD ("fEPres","EP resolution",4,binsv1,xminv1,xmaxv1);
763 fOutputList->Add(fEPres);
764
06d74a7b 765 //phi,fTPCnSigma,cent,pt,EovP,GetDeltaPhi(phi,evPlaneTPC),GetDeltaPhi(phi,evPlaneV0A),GetDeltaPhi(phi,evPlaneV0C),fFlagPhotonicElec,fFlagPhotonicElecBCG,m20,m02
1c424962 766 Int_t binsv2[12]={100,100,90,100,100,100,100,100,3,3,100,100};
767 Double_t xminv2[12]={0,-3.5,0,0,0,0,0,0,-1,-1,0,0};
768 Double_t xmaxv2[12]={2*TMath::Pi(),3.5,90,50,3,TMath::Pi(),TMath::Pi(),TMath::Pi(),2,2,1,1};
769 fCorr = new THnSparseD ("fCorr","Correlations",12,binsv2,xminv2,xmaxv2);
1d08b6e8 770 fOutputList->Add(fCorr);
771
772 Int_t binsv3[5]={90,100,100,100,100}; // cent, pt, TPCcos2DeltaPhi, V0Acos2DeltaPhi, V0Ccos2DeltaPhi
773 Double_t xminv3[5]={0,0,-1,-1,-1};
774 Double_t xmaxv3[5]={90,50,1,1,1};
775 feV2 = new THnSparseD ("feV2","inclusive electron v2",5,binsv3,xminv3,xmaxv3);
776 fOutputList->Add(feV2);
777
778 Int_t binsv4[5]={90,100,100,100,100}; // cent, pt, TPCdeltaPhi, V0AdeltaPhi, V0CdeltaPhi
779 Double_t xminv4[5]={0,0,-1,-1,-1};
780 Double_t xmaxv4[5]={90,50,1,1,1};
781 fphoteV2 = new THnSparseD ("fphoteV2","photonic electron v2",5,binsv4,xminv4,xmaxv4);
782 fOutputList->Add(fphoteV2);
783
44fdf7ac 784 Int_t binsv5[7]={100,90,100,100,100,100,100}; // phi, cent, pt, EovP, TPCdeltaPhi, V0AdeltaPhi, V0CdeltaPhi
785 Double_t xminv5[7]={0,0,0,0,-1,-1,-1};
786 Double_t xmaxv5[7]={2*TMath::Pi(),90,50,3,1,1,1};
787 fChargPartV2 = new THnSparseD ("fChargPartV2","Charged particle v2",7,binsv5,xminv5,xmaxv5);
1d08b6e8 788 fOutputList->Add(fChargPartV2);
789
790 Int_t binsv6[5]={90,100,100,100,100}; // cent, pt, TPCdeltaPhi, V0AdeltaPhi, V0CdeltaPhi
791 Double_t xminv6[5]={0,0,-1,-1,-1};
792 Double_t xmaxv6[5]={90,50,1,1,1};
793 feTPCV2 = new THnSparseD ("feTPCV2","inclusive electron v2 (TPC)",5,binsv6,xminv6,xmaxv6);
794 fOutputList->Add(feTPCV2);
795
6522b3b5 796 //EovP,fTPCnSigma,partPt,fFlagPhotonicElec,fFlagPhotonicElecBCG,whichPart,cent,pt,firstMother,secondMother,thirdMother,iHijing,motherPt,secondMotherPt,thirdMotherPt
797 Int_t binsv7[15]={100,100,100,3,3,5,90,100,5,5,5,3,100,100,100};
32d22164 798 Double_t xminv7[15]={0,-3.5,0,-1,-1,-1,0,0,-1,-1,-1,-1,0,0,0};
6522b3b5 799 Double_t xmaxv7[15]={3,3.5,50,2,2,4,90,50,4,4,4,2,50,50,50};
32d22164 800 fMCphotoElecPt = new THnSparseD ("fMCphotoElecPt", "pt distribution (MC)",15,binsv7,xminv7,xmaxv7);
3a94d6fe 801 fOutputList->Add(fMCphotoElecPt);
6522b3b5 802
803 fGammaWeight = new TH2F("fGammaWeight", "Gamma weight",100,0,50,3,-1,2);
804 fOutputList->Add(fGammaWeight);
805
806 fPi0Weight = new TH2F("fPi0Weight", "Pi0 weight",100,0,50,3,-1,2);
807 fOutputList->Add(fPi0Weight);
808
809 fEtaWeight = new TH2F("fEtaWeight", "Eta weight",100,0,50,3,-1,2);
810 fOutputList->Add(fEtaWeight);
1c424962 811
def484dc 812 fD0Weight = new TH2F("fD0Weight", "D0 weight",100,0,50,3,-1,2);
813 fOutputList->Add(fD0Weight);
814
815 fDplusWeight = new TH2F("fDplusWeight", "D+ weight",100,0,50,3,-1,2);
816 fOutputList->Add(fDplusWeight);
817
818 fDminusWeight = new TH2F("fDminusWeight", "D- weight",100,0,50,3,-1,2);
819 fOutputList->Add(fDminusWeight);
820
908103bf 821 fD0_e = new TH2F("fD0_e", "D0 vs e",100,0,50,200,-6.3,6.3);
937b2ff2 822 fOutputList->Add(fD0_e);
def484dc 823
937b2ff2 824 for(Int_t k = 0; k < 6; k++) {
825
826 TString De_name = Form("fDe%d",k);
827 TString D0e_name = Form("fD0e%d",k);
828 TString Dpluse_name = Form("fDpluse%d",k);
829 TString Dminuse_name = Form("fDminuse%d",k);
830
831 fDe[k] = new TH1F((const char*)De_name,"",100,0,50);
832 fD0e[k] = new TH1F((const char*)D0e_name,"",100,0,50);
833 fDpluse[k] = new TH1F((const char*)Dpluse_name,"",100,0,50);
834 fDminuse[k] = new TH1F((const char*)Dminuse_name,"",100,0,50);
835
836 fOutputList->Add(fDe[k]);
837 fOutputList->Add(fD0e[k]);
838 fOutputList->Add(fDpluse[k]);
839 fOutputList->Add(fDminuse[k]);
840
841 }
def484dc 842
1c424962 843 int nbin_v2 = 8;
844 double bin_v2[9] = {2,2.5,3,4,6,8,10,13,18};
845
846 fTot_pi0e = new TH1F("fTot_pi0e","fTot_pi0e",nbin_v2,bin_v2);
847 fOutputList->Add(fTot_pi0e);
6522b3b5 848
1c424962 849 fPhot_pi0e = new TH1F("fPhot_pi0e","fPhot_pi0e",nbin_v2,bin_v2);
850 fOutputList->Add(fPhot_pi0e);
851
852 fPhotBCG_pi0e = new TH1F("fPhotBCG_pi0e","fPhotBCG_pi0e",nbin_v2,bin_v2);
853 fOutputList->Add(fPhotBCG_pi0e);
854
855 fTot_etae = new TH1F("fTot_etae","fTot_etae",nbin_v2,bin_v2);
856 fOutputList->Add(fTot_etae);
857
858 fPhot_etae = new TH1F("fPhot_etae","fPhot_etae",nbin_v2,bin_v2);
859 fOutputList->Add(fPhot_etae);
860
861 fPhotBCG_etae = new TH1F("fPhotBCG_etae","fPhotBCG_etae",nbin_v2,bin_v2);
862 fOutputList->Add(fPhotBCG_etae);
863
1d08b6e8 864 PostData(1,fOutputList);
865}
866
867//________________________________________________________________________
914042c2 868void AliAnalysisTaskFlowTPCEMCalEP::Terminate(Option_t *)
1d08b6e8 869{
870 // Info("Terminate");
871 AliAnalysisTaskSE::Terminate();
872}
873
874//________________________________________________________________________
914042c2 875Bool_t AliAnalysisTaskFlowTPCEMCalEP::ProcessCutStep(Int_t cutStep, AliVParticle *track)
1d08b6e8 876{
877 // Check single track cuts for a given cut step
878 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
879 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
880 return kTRUE;
881}
882//_________________________________________
914042c2 883void AliAnalysisTaskFlowTPCEMCalEP::SelectPhotonicElectron(Int_t iTracks,AliESDtrack *track,Bool_t &fFlagPhotonicElec, Bool_t &fFlagPhotonicElecBCG)
1d08b6e8 884{
885 //Identify non-heavy flavour electrons using Invariant mass method
886
887 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
888 fTrackCuts->SetRequireTPCRefit(kTRUE);
5945f829 889 fTrackCuts->SetRequireITSRefit(kTRUE);
1d08b6e8 890 fTrackCuts->SetEtaRange(-0.7,0.7);
891 fTrackCuts->SetRequireSigmaToVertex(kTRUE);
892 fTrackCuts->SetMaxChi2PerClusterTPC(3.5);
893 fTrackCuts->SetMinNClustersTPC(100);
894
895 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
896
897 Bool_t flagPhotonicElec = kFALSE;
06d74a7b 898 Bool_t flagPhotonicElecBCG = kFALSE;
1d08b6e8 899
b1a1875f 900 for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){
901
5945f829 902 if(jTracks==iTracks) continue;
b1a1875f 903
1d08b6e8 904 AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
905 if (!trackAsso) {
906 printf("ERROR: Could not receive track %d\n", jTracks);
907 continue;
908 }
909
5945f829 910 Double_t dEdxAsso = -999., ptAsso=-999.;
06d74a7b 911 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
912 Double_t openingAngle = -999., mass=999., width = -999;
1d08b6e8 913
914 dEdxAsso = trackAsso->GetTPCsignal();
915 ptAsso = trackAsso->Pt();
916 Int_t chargeAsso = trackAsso->Charge();
917 Int_t charge = track->Charge();
918
06d74a7b 919 if(ptAsso <0.5) continue;
1d08b6e8 920 if(!fTrackCuts->AcceptTrack(trackAsso)) continue;
06d74a7b 921 if(dEdxAsso <65 || dEdxAsso>100) continue;
1d08b6e8 922
923 Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
924 if(charge>0) fPDGe1 = -11;
925 if(chargeAsso>0) fPDGe2 = -11;
926
927 if(charge == chargeAsso) fFlagLS = kTRUE;
928 if(charge != chargeAsso) fFlagULS = kTRUE;
929
930 AliKFParticle ge1(*track, fPDGe1);
931 AliKFParticle ge2(*trackAsso, fPDGe2);
932 AliKFParticle recg(ge1, ge2);
933
934 if(recg.GetNDF()<1) continue;
935 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
936 if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
937
938 AliKFVertex primV(*pVtx);
939 primV += recg;
940 recg.SetProductionVertex(primV);
941
942 recg.SetMassConstraint(0,0.0001);
943
944 openingAngle = ge1.GetAngle(ge2);
945 if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
946 if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
947
5945f829 948 if(openingAngle > fOpeningAngleCut) continue;
06d74a7b 949
950 recg.GetMass(mass,width);
5945f829 951
1d08b6e8 952 if(fFlagLS) fInvmassLS->Fill(mass);
953 if(fFlagULS) fInvmassULS->Fill(mass);
954
5945f829 955 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec) flagPhotonicElec = kTRUE;
06d74a7b 956 if(mass<fInvmassCut && fFlagLS && !flagPhotonicElecBCG) flagPhotonicElecBCG = kTRUE;
1d08b6e8 957
958 }
959 fFlagPhotonicElec = flagPhotonicElec;
06d74a7b 960 fFlagPhotonicElecBCG = flagPhotonicElecBCG;
1d08b6e8 961
962}
963//_________________________________________
914042c2 964Double_t AliAnalysisTaskFlowTPCEMCalEP::GetCos2DeltaPhi(Double_t phiA,Double_t phiB) const
1d08b6e8 965{
966 //Get cos[2(phi-psi_EP)] or cos[2(psi_subEP1 - psi_subEP2)]
967 Double_t dPhi = TVector2::Phi_0_2pi(phiA - phiB);
968 if(dPhi > TMath::Pi()) dPhi = dPhi - TMath::Pi();
969 Double_t cos2DeltaPhi = TMath::Cos(2*dPhi);
970
971 return cos2DeltaPhi;
972}
1d08b6e8 973//_________________________________________
914042c2 974Double_t AliAnalysisTaskFlowTPCEMCalEP::GetDeltaPhi(Double_t phiA,Double_t phiB) const
1d08b6e8 975{
976 //Get phi-psi_EP
977 Double_t dPhi = TVector2::Phi_0_2pi(phiA - phiB);
978 if(dPhi > TMath::Pi()) dPhi = dPhi - TMath::Pi();
979
980 return dPhi;
981}
1c424962 982//_________________________________________
983Double_t AliAnalysisTaskFlowTPCEMCalEP::GetPi0weight(Double_t mcPi0pT) const
984{
937b2ff2 985 //Get Pi0 weight
1c424962 986 double weight = 1.0;
987
988 if(mcPi0pT>0.0 && mcPi0pT<5.0)
989 {
990
def484dc 991// weight = (2.877091*mcPi0pT)/(TMath::Power(0.706963+mcPi0pT/3.179309,17.336628)*exp(-mcPi0pT));
992 weight = (2.392024*mcPi0pT)/(TMath::Power(0.688810+mcPi0pT/3.005145,16.811845)*exp(-mcPi0pT));
1c424962 993 }
994 else
995 {
def484dc 996// weight = (0.0004*mcPi0pT)/TMath::Power(-0.176181+mcPi0pT/3.989747,5.629235);
997 weight = (0.000186*mcPi0pT)/TMath::Power(-0.606279+mcPi0pT/3.158642,4.365540);
1c424962 998 }
999 return weight;
1000}
1001//_________________________________________
1002Double_t AliAnalysisTaskFlowTPCEMCalEP::GetEtaweight(Double_t mcEtapT) const
1003{
937b2ff2 1004 //Get eta weight
1c424962 1005 double weight = 1.0;
1006
def484dc 1007// weight = (0.818052*mcEtapT)/(TMath::Power(0.358651+mcEtapT/2.878631,9.494043));
1008 weight = (0.622703*mcEtapT)/(TMath::Power(0.323045+mcEtapT/2.736407,9.180356));
1c424962 1009
1010 return weight;
1011}
def484dc 1012//_________________________________________
937b2ff2 1013Double_t AliAnalysisTaskFlowTPCEMCalEP::GetDweight(Int_t whichD, Double_t mcDpT) const
def484dc 1014{
937b2ff2 1015 //get D weights
def484dc 1016 double weight = 1.0;
937b2ff2 1017
1018 if (whichD == 0) weight = 0.271583*TMath::Landau(mcDpT,3.807103,1.536753,0); // D
1019 if (whichD == 1) weight = 0.300771*TMath::Landau(mcDpT,3.725771,1.496980,0); // D0
1020 if (whichD == 2) weight = 0.247280*TMath::Landau(mcDpT,3.746811,1.607551,0); // D+
1021 if (whichD == 3) weight = 0.249410*TMath::Landau(mcDpT,3.611508,1.632196,0); //D-
1022
def484dc 1023 return weight;
1024}