]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliAnalysisTaskHFECal.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
added EMCal QA task
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliAnalysisTaskHFECal.cxx
CommitLineData
bfc7c23b 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 with EMCal triggered events
17// Author: Shingo Sakai
18
19
20#include "TChain.h"
21#include "TTree.h"
22#include "TH2F.h"
8efacc22 23#include "TF1.h"
bfc7c23b 24#include "TMath.h"
25#include "TCanvas.h"
26#include "THnSparse.h"
27#include "TLorentzVector.h"
28#include "TString.h"
29#include "TFile.h"
09022877 30#include "TGraphErrors.h"
bfc7c23b 31
6451f693 32#include "TDatabasePDG.h"
33
bfc7c23b 34#include "AliAnalysisTask.h"
35#include "AliAnalysisManager.h"
36
37#include "AliESDEvent.h"
38#include "AliESDHandler.h"
39#include "AliAODEvent.h"
40#include "AliAODHandler.h"
41
e7d87aef 42#include "AliMCEventHandler.h"
43#include "AliMCEvent.h"
44#include "AliMCParticle.h"
45
bfc7c23b 46#include "AliAnalysisTaskHFECal.h"
47#include "TGeoGlobalMagField.h"
48#include "AliLog.h"
49#include "AliAnalysisTaskSE.h"
50#include "TRefArray.h"
51#include "TVector.h"
52#include "AliESDInputHandler.h"
53#include "AliESDpid.h"
54#include "AliESDtrackCuts.h"
55#include "AliPhysicsSelection.h"
56#include "AliESDCaloCluster.h"
57#include "AliAODCaloCluster.h"
58#include "AliEMCALRecoUtils.h"
59#include "AliEMCALGeometry.h"
60#include "AliGeomManager.h"
61#include "stdio.h"
62#include "TGeoManager.h"
63#include "iostream"
64#include "fstream"
65
66#include "AliEMCALTrack.h"
67#include "AliMagF.h"
68
69#include "AliKFParticle.h"
70#include "AliKFVertex.h"
71
72#include "AliPID.h"
73#include "AliPIDResponse.h"
74#include "AliHFEcontainer.h"
75#include "AliHFEcuts.h"
76#include "AliHFEpid.h"
77#include "AliHFEpidBase.h"
78#include "AliHFEpidQAmanager.h"
79#include "AliHFEtools.h"
80#include "AliCFContainer.h"
81#include "AliCFManager.h"
82
e7d87aef 83#include "AliStack.h"
84
bfc7c23b 85#include "AliCentrality.h"
86
c2aad3ae 87using namespace std;
88
bfc7c23b 89ClassImp(AliAnalysisTaskHFECal)
90//________________________________________________________________________
91AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name)
92 : AliAnalysisTaskSE(name)
93 ,fESD(0)
e7d87aef 94 ,fMC(0)
46305ed6 95 ,stack(0)
fb87d707 96 ,fGeom(0)
bfc7c23b 97 ,fOutputList(0)
c6bca3f2 98 ,fqahist(0)
bfc7c23b 99 ,fTrackCuts(0)
100 ,fCuts(0)
101 ,fIdentifiedAsOutInz(kFALSE)
102 ,fPassTheEventCut(kFALSE)
103 ,fRejectKinkMother(kFALSE)
e7d87aef 104 ,fmcData(kFALSE)
bfc7c23b 105 ,fVz(0.0)
106 ,fCFM(0)
107 ,fPID(0)
108 ,fPIDqa(0)
109 ,fOpeningAngleCut(0.1)
6451f693 110 ,fMimpTassCut(0.5)
6b2cee73 111 ,fMimNsigassCut(-3)
a009053a 112 ,fInvmassCut(0) // no mass
113 ,fSetMassConstraint(kTRUE)
3c475a9d 114 ,fSetMassWidthCut(kFALSE)
2f8ef315 115 ,fSetMassNonlinear(kFALSE)
6451f693 116 ,fSetKFpart(kTRUE)
bfc7c23b 117 ,fNoEvents(0)
118 ,fEMCAccE(0)
f4e0d2d5 119 ,hEMCAccE(0)
bfc7c23b 120 ,fTrkpt(0)
121 ,fTrkEovPBef(0)
122 ,fTrkEovPAft(0)
123 ,fdEdxBef(0)
124 ,fdEdxAft(0)
125 ,fIncpT(0)
fb87d707 126 ,fIncpTM20(0)
bfc7c23b 127 ,fInvmassLS(0)
128 ,fInvmassULS(0)
8806ce6c 129 ,fInvmassLSmc(0)
130 ,fInvmassULSmc(0)
6ec70f8c 131 ,fInvmassLSreco(0)
132 ,fInvmassULSreco(0)
93c189c5 133 ,fInvmassLSmc0(0)
134 ,fInvmassLSmc1(0)
135 ,fInvmassLSmc2(0)
136 ,fInvmassLSmc3(0)
137 ,fInvmassULSmc0(0)
138 ,fInvmassULSmc1(0)
139 ,fInvmassULSmc2(0)
140 ,fInvmassULSmc3(0)
bfc7c23b 141 ,fOpeningAngleLS(0)
142 ,fOpeningAngleULS(0)
143 ,fPhotoElecPt(0)
144 ,fPhoElecPt(0)
fb87d707 145 ,fPhoElecPtM20(0)
0f566acb 146 ,fPhoElecPtM20Mass(0)
f09766dd 147 ,fSameElecPt(0)
fb87d707 148 ,fSameElecPtM20(0)
0f566acb 149 ,fSameElecPtM20Mass(0)
ebef1da9 150 ,fSemiElecPtM20(0)
bfc7c23b 151 ,fTrackPtBefTrkCuts(0)
152 ,fTrackPtAftTrkCuts(0)
153 ,fTPCnsigma(0)
154 ,fCent(0)
155 ,fEleInfo(0)
e1f0fb74 156 ,fElenSigma(0)
feffe705 157 /*
fb87d707 158 ,fClsEBftTrigCut(0)
159 ,fClsEAftTrigCut(0)
160 ,fClsEAftTrigCut1(0)
161 ,fClsEAftTrigCut2(0)
162 ,fClsEAftTrigCut3(0)
163 ,fClsEAftTrigCut4(0)
164 ,fClsETime(0)
165 ,fClsETime1(0)
166 ,fTrigTimes(0)
42c75692 167 ,fCellCheck(0)
feffe705 168 */
e7d87aef 169 ,fInputHFEMC(0)
feffe705 170 ,fInputAlle(0)
e7d87aef 171 ,fIncpTMChfe(0)
3db00c72 172 ,fIncpTMChfeAll(0)
e7d87aef 173 ,fIncpTMCM20hfe(0)
3db00c72 174 ,fIncpTMCM20hfeAll(0)
60544aea 175 ,fIncpTMCM20hfeCheck(0)
bd6ee6dd 176 ,fInputHFEMC_weight(0)
177 ,fIncpTMCM20hfeCheck_weight(0)
e7d87aef 178 ,fIncpTMCpho(0)
179 ,fIncpTMCM20pho(0)
180 ,fPhoElecPtMC(0)
181 ,fPhoElecPtMCM20(0)
0f566acb 182 ,fPhoElecPtMCM20Mass(0)
e7d87aef 183 ,fSameElecPtMC(0)
184 ,fSameElecPtMCM20(0)
0f566acb 185 ,fSameElecPtMCM20Mass(0)
e4b0faf2 186 ,fIncpTMCM20pho_pi0e(0)
187 ,fPhoElecPtMCM20_pi0e(0)
188 ,fSameElecPtMCM20_pi0e(0)
93c189c5 189 ,fIncpTMCM20pho_eta(0)
190 ,fPhoElecPtMCM20_eta(0)
191 ,fSameElecPtMCM20_eta(0)
697ecf6b 192 ,fIncpTMCpho_pi0e_TPC(0)
193 ,fPhoElecPtMC_pi0e_TPC(0)
194 ,fSameElecPtMC_pi0e_TPC(0)
2f8ef315 195 ,fIncpTMCpho_eta_TPC(0)
196 ,fPhoElecPtMC_eta_TPC(0)
197 ,fSameElecPtMC_eta_TPC(0)
a6df418c 198 ,CheckNclust(0)
199 ,CheckNits(0)
ae88597c 200 ,CheckDCA(0)
d113d7cd 201 ,Hpi0pTcheck(0)
19325033 202 ,HETApTcheck(0)
bd6ee6dd 203 ,HphopTcheck(0)
5e0e45b3 204 ,HDpTcheck(0)
205 ,HBpTcheck(0)
93c189c5 206 ,fpTCheck(0)
50919258 207 ,fMomDtoE(0)
70448e3b 208 ,fLabelCheck(0)
209 ,fgeoFake(0)
d19af75d 210 ,fFakeTrk0(0)
211 ,fFakeTrk1(0)
5e38d973 212 ,ftimingEle(0)
6f4a9952 213 ,fIncMaxE(0)
5e38d973 214 ,fIncReco(0)
215 ,fPhoReco(0)
216 ,fSamReco(0)
c82f1b9d 217 ,fIncRecoMaxE(0)
218 ,fPhoRecoMaxE(0)
219 ,fSamRecoMaxE(0)
466d4f7a 220 ,fPhoVertexReco_TPC(0)
221 ,fPhoVertexReco_TPC_Invmass(0)
05ac9660 222 ,fPhoVertexReco_EMCal(0)
223 ,fPhoVertexReco_Invmass(0)
224 ,fPhoVertexReco_step0(0)
225 ,fPhoVertexReco_step1(0)
f2fa5a00 226 ,fMatchV0_0(0)
227 ,fMatchV0_1(0)
228 ,fMatchMC_0(0)
229 ,fMatchMC_1(0)
36e9774f 230 ,fpair(0)
c5a4ddbf 231 ,fFakeRejection0(0)
232 ,fFakeRejection1(0)
233 ,fFakeRejection2(0)
9e88412a 234 ,EopFake(0)
235 ,EopTrue(0)
236 ,MatchFake(0)
237 ,MatchTrue(0)
ebf2f603 238 ,MatchTrCheck(0)
239 ,MatchTrEop(0)
09022877 240 //,fnSigEtaCorr(NULL)
bfc7c23b 241{
242 //Named constructor
243
244 fPID = new AliHFEpid("hfePid");
245 fTrackCuts = new AliESDtrackCuts();
246
09022877 247 for(int i=0; i<7; i++)fnSigEtaCorr[i] = 0;
248
bfc7c23b 249 // Define input and output slots here
250 // Input slot #0 works with a TChain
251 DefineInput(0, TChain::Class());
252 // Output slot #0 id reserved by the base class for AOD
253 // Output slot #1 writes into a TH1 container
254 // DefineOutput(1, TH1I::Class());
255 DefineOutput(1, TList::Class());
256 // DefineOutput(3, TTree::Class());
257}
258
259//________________________________________________________________________
260AliAnalysisTaskHFECal::AliAnalysisTaskHFECal()
261 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskHFECal")
262 ,fESD(0)
e7d87aef 263 ,fMC(0)
46305ed6 264 ,stack(0)
fb87d707 265 ,fGeom(0)
bfc7c23b 266 ,fOutputList(0)
c6bca3f2 267 ,fqahist(0)
bfc7c23b 268 ,fTrackCuts(0)
269 ,fCuts(0)
270 ,fIdentifiedAsOutInz(kFALSE)
271 ,fPassTheEventCut(kFALSE)
272 ,fRejectKinkMother(kFALSE)
e7d87aef 273 ,fmcData(kFALSE)
bfc7c23b 274 ,fVz(0.0)
275 ,fCFM(0)
276 ,fPID(0)
277 ,fPIDqa(0)
278 ,fOpeningAngleCut(0.1)
6451f693 279 ,fMimpTassCut(0.5)
6b2cee73 280 ,fMimNsigassCut(-3)
a009053a 281 ,fInvmassCut(0) // no mass
282 ,fSetMassConstraint(kTRUE)
3c475a9d 283 ,fSetMassWidthCut(kFALSE)
2f8ef315 284 ,fSetMassNonlinear(kFALSE)
6451f693 285 ,fSetKFpart(kTRUE)
bfc7c23b 286 ,fNoEvents(0)
287 ,fEMCAccE(0)
f4e0d2d5 288 ,hEMCAccE(0)
bfc7c23b 289 ,fTrkpt(0)
290 ,fTrkEovPBef(0)
291 ,fTrkEovPAft(0)
292 ,fdEdxBef(0)
293 ,fdEdxAft(0)
294 ,fIncpT(0)
fb87d707 295 ,fIncpTM20(0)
bfc7c23b 296 ,fInvmassLS(0)
297 ,fInvmassULS(0)
8806ce6c 298 ,fInvmassLSmc(0)
299 ,fInvmassULSmc(0)
6ec70f8c 300 ,fInvmassLSreco(0)
301 ,fInvmassULSreco(0)
93c189c5 302 ,fInvmassLSmc0(0)
303 ,fInvmassLSmc1(0)
304 ,fInvmassLSmc2(0)
305 ,fInvmassLSmc3(0)
306 ,fInvmassULSmc0(0)
307 ,fInvmassULSmc1(0)
308 ,fInvmassULSmc2(0)
309 ,fInvmassULSmc3(0)
bfc7c23b 310 ,fOpeningAngleLS(0)
311 ,fOpeningAngleULS(0)
312 ,fPhotoElecPt(0)
313 ,fPhoElecPt(0)
fb87d707 314 ,fPhoElecPtM20(0)
0f566acb 315 ,fPhoElecPtM20Mass(0)
f09766dd 316 ,fSameElecPt(0)
fb87d707 317 ,fSameElecPtM20(0)
0f566acb 318 ,fSameElecPtM20Mass(0)
ebef1da9 319 ,fSemiElecPtM20(0)
bfc7c23b 320 ,fTrackPtBefTrkCuts(0)
321 ,fTrackPtAftTrkCuts(0)
322 ,fTPCnsigma(0)
323 ,fCent(0)
324 ,fEleInfo(0)
e1f0fb74 325 ,fElenSigma(0)
feffe705 326 /*
fb87d707 327 ,fClsEBftTrigCut(0)
328 ,fClsEAftTrigCut(0)
329 ,fClsEAftTrigCut1(0)
330 ,fClsEAftTrigCut2(0)
331 ,fClsEAftTrigCut3(0)
332 ,fClsEAftTrigCut4(0)
333 ,fClsETime(0)
334 ,fClsETime1(0)
335 ,fTrigTimes(0)
42c75692 336 ,fCellCheck(0)
feffe705 337 */
e7d87aef 338 ,fInputHFEMC(0)
feffe705 339 ,fInputAlle(0)
e7d87aef 340 ,fIncpTMChfe(0)
3db00c72 341 ,fIncpTMChfeAll(0)
e7d87aef 342 ,fIncpTMCM20hfe(0)
3db00c72 343 ,fIncpTMCM20hfeAll(0)
60544aea 344 ,fIncpTMCM20hfeCheck(0)
bd6ee6dd 345 ,fInputHFEMC_weight(0)
346 ,fIncpTMCM20hfeCheck_weight(0)
e7d87aef 347 ,fIncpTMCpho(0)
348 ,fIncpTMCM20pho(0)
349 ,fPhoElecPtMC(0)
350 ,fPhoElecPtMCM20(0)
0f566acb 351 ,fPhoElecPtMCM20Mass(0)
e7d87aef 352 ,fSameElecPtMC(0)
353 ,fSameElecPtMCM20(0)
0f566acb 354 ,fSameElecPtMCM20Mass(0)
e4b0faf2 355 ,fIncpTMCM20pho_pi0e(0)
356 ,fPhoElecPtMCM20_pi0e(0)
357 ,fSameElecPtMCM20_pi0e(0)
93c189c5 358 ,fIncpTMCM20pho_eta(0)
359 ,fPhoElecPtMCM20_eta(0)
360 ,fSameElecPtMCM20_eta(0)
697ecf6b 361 ,fIncpTMCpho_pi0e_TPC(0)
362 ,fPhoElecPtMC_pi0e_TPC(0)
363 ,fSameElecPtMC_pi0e_TPC(0)
2f8ef315 364 ,fIncpTMCpho_eta_TPC(0)
365 ,fPhoElecPtMC_eta_TPC(0)
366 ,fSameElecPtMC_eta_TPC(0)
a6df418c 367 ,CheckNclust(0)
368 ,CheckNits(0)
ae88597c 369 ,CheckDCA(0)
d113d7cd 370 ,Hpi0pTcheck(0)
19325033 371 ,HETApTcheck(0)
e4b0faf2 372 ,HphopTcheck(0)
5e0e45b3 373 ,HDpTcheck(0)
374 ,HBpTcheck(0)
93c189c5 375 ,fpTCheck(0)
50919258 376 ,fMomDtoE(0)
70448e3b 377 ,fLabelCheck(0)
378 ,fgeoFake(0)
d19af75d 379 ,fFakeTrk0(0)
380 ,fFakeTrk1(0)
59d998de 381 ,ftimingEle(0)
6f4a9952 382 ,fIncMaxE(0)
5e38d973 383 ,fIncReco(0)
384 ,fPhoReco(0)
385 ,fSamReco(0)
c82f1b9d 386 ,fIncRecoMaxE(0)
387 ,fPhoRecoMaxE(0)
388 ,fSamRecoMaxE(0)
466d4f7a 389 ,fPhoVertexReco_TPC(0)
390 ,fPhoVertexReco_TPC_Invmass(0)
05ac9660 391 ,fPhoVertexReco_EMCal(0)
392 ,fPhoVertexReco_Invmass(0)
393 ,fPhoVertexReco_step0(0)
394 ,fPhoVertexReco_step1(0)
f2fa5a00 395 ,fMatchV0_0(0)
396 ,fMatchV0_1(0)
397 ,fMatchMC_0(0)
398 ,fMatchMC_1(0)
36e9774f 399 ,fpair(0)
c5a4ddbf 400 ,fFakeRejection0(0)
401 ,fFakeRejection1(0)
402 ,fFakeRejection2(0)
9e88412a 403 ,EopFake(0)
404 ,EopTrue(0)
405 ,MatchFake(0)
406 ,MatchTrue(0)
ebf2f603 407 ,MatchTrCheck(0)
408 ,MatchTrEop(0)
09022877 409 //,fnSigEtaCorr(NULL)
bfc7c23b 410{
411 //Default constructor
412 fPID = new AliHFEpid("hfePid");
413
414 fTrackCuts = new AliESDtrackCuts();
415
09022877 416 for(int i=0; i<7; i++)fnSigEtaCorr[i] = 0;
417
bfc7c23b 418 // Constructor
419 // Define input and output slots here
420 // Input slot #0 works with a TChain
421 DefineInput(0, TChain::Class());
422 // Output slot #0 id reserved by the base class for AOD
423 // Output slot #1 writes into a TH1 container
424 // DefineOutput(1, TH1I::Class());
425 DefineOutput(1, TList::Class());
426 //DefineOutput(3, TTree::Class());
427}
428//_________________________________________
429
430AliAnalysisTaskHFECal::~AliAnalysisTaskHFECal()
431{
432 //Destructor
433
434 delete fOutputList;
fb87d707 435 delete fGeom;
bfc7c23b 436 delete fPID;
a6df418c 437 delete fCuts;
bfc7c23b 438 delete fCFM;
439 delete fPIDqa;
440 delete fTrackCuts;
441}
442//_________________________________________
443
444void AliAnalysisTaskHFECal::UserExec(Option_t*)
445{
446 //Main loop
447 //Called for each event
448
449 // create pointer to event
450 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
451 if (!fESD) {
452 printf("ERROR: fESD not available\n");
453 return;
454 }
455
456 if(!fCuts){
457 AliError("HFE cuts not available");
458 return;
459 }
460
461 if(!fPID->IsInitialized()){
462 // Initialize PID with the given run number
463 AliWarning("PID not initialised, get from Run no");
464 fPID->InitializePID(fESD->GetRunNumber());
465 }
466
e7d87aef 467 if(fmcData)fMC = MCEvent();
46305ed6 468 //AliStack* stack = NULL;
e7d87aef 469 if(fmcData && fMC)stack = fMC->Stack();
470
471 Float_t cent = -1.;
472 AliCentrality *centrality = fESD->GetCentrality();
473 cent = centrality->GetCentralityPercentile("V0M");
474
475 //---- fill MC track info
476 if(fmcData && fMC)
477 {
478 Int_t nParticles = stack->GetNtrack();
479 for (Int_t iParticle = 0; iParticle < nParticles; iParticle++) {
480 TParticle* particle = stack->Particle(iParticle);
481 int fPDG = particle->GetPdgCode();
482 double mcZvertex = fMC->GetPrimaryVertex()->GetZ();
483 double pTMC = particle->Pt();
feffe705 484 double proR = particle->R();
9b3495ff 485 double etaMC = particle->Eta();
486 if(fabs(etaMC)>0.6)continue;
e7d87aef 487 Bool_t mcInDtoE= kFALSE;
488 Bool_t mcInBtoE= kFALSE;
489
d113d7cd 490 Bool_t MChijing = fMC->IsFromBGEvent(iParticle);
fd426741 491 //if(!MChijing)printf("not MC hijing");
d113d7cd 492 int iHijing = 1;
493 if(!MChijing)iHijing = 0;
68d718e7 494 double mcphoinfo[5];
495 mcphoinfo[0] = cent;
496 mcphoinfo[1] = pTMC;
497 mcphoinfo[2] = iHijing;
68d718e7 498 if(fPDG==111)Hpi0pTcheck->Fill(mcphoinfo);
499 if(fPDG==221)HETApTcheck->Fill(mcphoinfo);
bd6ee6dd 500 if(fabs(fPDG)==411 || fabs(fPDG)==413 || fabs(fPDG)==421 || fabs(fPDG)==423 || fabs(fPDG)==431)HDpTcheck->Fill(pTMC,iHijing);
501 if(fabs(fPDG)==511 || fabs(fPDG)==513 || fabs(fPDG)==521 || fabs(fPDG)==523 || fabs(fPDG)==531)HBpTcheck->Fill(pTMC,iHijing);
d113d7cd 502
e4b0faf2 503 if(particle->GetFirstMother()>-1 && fPDG==22)
504 {
505 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
506 if(parentPID==111 || parentPID==221)HphopTcheck->Fill(pTMC,iHijing); // pi0->g & eta->g
507 }
508
e7d87aef 509 if(particle->GetFirstMother()>-1 && fabs(fPDG)==11)
510 {
511 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
bd6ee6dd 512 double pTMCparent = stack->Particle(particle->GetFirstMother())->Pt();
e7d87aef 513 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(fPDG)==11)mcInDtoE = kTRUE;
514 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(fPDG)==11)mcInBtoE = kTRUE;
bd6ee6dd 515 if((mcInBtoE || mcInDtoE) && fabs(mcZvertex)<10.0)
516 {
517 fInputHFEMC->Fill(cent,pTMC);
518 double mcinfo[5];
519 mcinfo[0] = cent;
520 mcinfo[1] = pTMC;
521 mcinfo[2] = 0.0;
522 mcinfo[3] = iHijing;
523 mcinfo[4] = pTMCparent;
524 fInputHFEMC_weight->Fill(mcinfo);
525 }
e7d87aef 526 }
527
5e0e45b3 528
feffe705 529 if(proR<7.0 && fabs(fPDG)==11)fInputAlle->Fill(cent,pTMC);
530
e7d87aef 531 }
532 }
533
bfc7c23b 534 fNoEvents->Fill(0);
535
536 Int_t fNOtrks = fESD->GetNumberOfTracks();
537 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
538
539 Double_t pVtxZ = -999;
540 pVtxZ = pVtx->GetZ();
541
542 if(TMath::Abs(pVtxZ)>10) return;
543 fNoEvents->Fill(1);
544
f4e0d2d5 545 if(fNOtrks<1) return;
bfc7c23b 546 fNoEvents->Fill(2);
547
548 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
549 if(!pidResponse){
550 AliDebug(1, "Using default PID Response");
551 pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
552 }
553
554 fPID->SetPIDResponse(pidResponse);
555
556 fCFM->SetRecEventInfo(fESD);
557
bfc7c23b 558 fCent->Fill(cent);
559
b12bc641 560 //if(cent>90.) return;
bfc7c23b 561
562 // Calorimeter info.
563
42c75692 564 FindTriggerClusters();
565
bfc7c23b 566 // make EMCAL array
6f4a9952 567 double maxE = 0.0;
bfc7c23b 568 for(Int_t iCluster=0; iCluster<fESD->GetNumberOfCaloClusters(); iCluster++)
569 {
570 AliESDCaloCluster *clust = fESD->GetCaloCluster(iCluster);
68d718e7 571 if(clust && clust->IsEMCAL())
bfc7c23b 572 {
573 double clustE = clust->E();
574 float emcx[3]; // cluster pos
575 clust->GetPosition(emcx);
576 TVector3 clustpos(emcx[0],emcx[1],emcx[2]);
577 double emcphi = clustpos.Phi();
578 double emceta = clustpos.Eta();
fb87d707 579 double calInfo[5];
580 calInfo[0] = emcphi; calInfo[1] = emceta; calInfo[2] = clustE; calInfo[3] = cent; calInfo[4] = clust->Chi2();
f4e0d2d5 581 //fEMCAccE->Fill(calInfo);
f4e0d2d5 582 hEMCAccE->Fill(cent,clustE);
6f4a9952 583 if(clustE>maxE)maxE = clustE;
bfc7c23b 584 }
585 }
586
587 // Track loop
588 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
589 AliESDtrack* track = fESD->GetTrack(iTracks);
590 if (!track) {
591 printf("ERROR: Could not receive track %d\n", iTracks);
592 continue;
593 }
e7d87aef 594
89f41a30 595 //--- Get MC informtion
596
46305ed6 597 int parentlabel = 99999;
598 int parentPID = 99999;
599 int grand_parentlabel = 99999;
600 int grand_parentPID = 99999;
e7d87aef 601 Bool_t mcPho = kFALSE;
602 Bool_t mcDtoE= kFALSE;
603 Bool_t mcBtoE= kFALSE;
93c189c5 604 Bool_t mcOrgPi0 = kFALSE;
605 Bool_t mcOrgEta = kFALSE;
e7d87aef 606 double mcele = -1.;
44be9e1d 607 double mcpT = 0.0;
b12bc641 608 double mcMompT = 0.0;
e4b0faf2 609 //double mcGrandMompT = 0.0;
610 double mcWeight = -10.0;
bd50bf24 611 double conv_proR = -1.0;
d113d7cd 612
613 int iHijing = 1;
70448e3b 614 int mcLabel = -1;
d113d7cd 615
e7d87aef 616 if(fmcData && fMC && stack)
617 {
feffe705 618 Int_t label = TMath::Abs(track->GetLabel());
7d495dc5 619 //mcLabel = track->GetLabel();
1dd85087 620 mcLabel = abs(track->GetLabel()); // check for conv. issue
bd6ee6dd 621
8efacc22 622 if(mcLabel>-1)
bd6ee6dd 623 {
624
625 Bool_t MChijing = fMC->IsFromBGEvent(label);
626 if(!MChijing)iHijing = 0;
627
628 TParticle* particle = stack->Particle(label);
629 int mcpid = particle->GetPdgCode();
630 mcpT = particle->Pt();
bd50bf24 631 conv_proR = particle->R();
bd6ee6dd 632 //printf("MCpid = %d",mcpid);
633 if(particle->GetFirstMother()>-1)
634 {
635 //int parentlabel = particle->GetFirstMother();
636 //int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
637 //mcMompT = stack->Particle(particle->GetFirstMother())->Pt();
638 FindMother(particle, parentlabel, parentPID);
639 mcMompT = stack->Particle(parentlabel)->Pt();
640 if((parentPID==22 || parentPID==111 || parentPID==221)&& fabs(mcpid)==11)mcPho = kTRUE;
641 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(mcpid)==11)mcDtoE = kTRUE;
642 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(mcpid)==11)mcBtoE = kTRUE;
643
644 // make D->e pT correlation
645 if(mcDtoE)fMomDtoE->Fill(mcpT,mcMompT);
646
647 //cout << "check PID = " << parentPID << endl;
648 //cout << "check pho = " << mcPho << endl;
649 //cout << "check D or B = " << mcDtoE << endl;
650 // pi->e (Dalitz)
651 if(parentPID==111 && fabs(mcpid)==11 && mcMompT>0.0)
652 {
653 //cout << "find pi0->e " << endl;
654 mcOrgPi0 = kTRUE;
655 mcWeight = GetMCweight(mcMompT);
656 }
657 // eta->e (Dalitz)
658 if(parentPID==221 && fabs(mcpid)==11 && mcMompT>0.0)
659 {
660 //cout << "find Eta->e " << endl;
661 mcOrgEta = kTRUE;
662 mcWeight = GetMCweightEta(mcMompT);
663 }
664
665 // access grand parent
666 TParticle* particle_parent = stack->Particle(parentlabel); // get parent pointer
667 //if(particle_parent->GetFirstMother()>-1 && parentPID==22 && fabs(mcpid)==11) // get grand parent g->e
668 if(particle_parent->GetFirstMother()>-1 && (parentPID==22 || parentPID==111) && fabs(mcpid)==11) // get grand parent g->e
669 {
670 //int grand_parentPID = stack->Particle(particle_parent->GetFirstMother())->GetPdgCode();
671 //double pTtmp = stack->Particle(particle_parent->GetFirstMother())->Pt();
672 FindMother(particle_parent, grand_parentlabel, grand_parentPID);
673 double mcGrandpT = stack->Particle(grand_parentlabel)->Pt();
674 if(grand_parentPID==111 && mcGrandpT>0.0)
675 {
676 // check eta->pi0 decay !
677 int grand2_parentlabel = 99999; int grand2_parentPID = 99999;
678 TParticle* particle_grand = stack->Particle(grand_parentlabel); // get parent pointer
679 FindMother(particle_grand, grand2_parentlabel, grand2_parentPID);
680 if(grand2_parentPID==221)
681 {
682 //cout << "find Eta->e " << endl;
683 double mcGrandpT2 = stack->Particle(grand2_parentlabel)->Pt();
684 mcOrgEta = kTRUE;
685 mcWeight = GetMCweight(mcGrandpT2);
686 mcMompT = mcGrandpT2;
687 }
688 else
689 {
690 //cout << "find pi0->e " << endl;
691 mcOrgPi0 = kTRUE;
692 mcWeight = GetMCweight(mcGrandpT);
693 mcMompT = mcGrandpT;
694 }
695 }
696
697 if(grand_parentPID==221 && mcGrandpT>0.0)
698 {
699 //cout << "find Eta->e " << endl;
700 mcOrgEta = kTRUE;
701 mcOrgPi0 = kFALSE;
702 mcWeight = GetMCweightEta(mcGrandpT);
703 mcMompT = mcGrandpT;
704 }
705 }
706 }
707
cc1605fe 708 //cout << "===================="<<endl;
709 //cout << "mcDtoE : " << mcDtoE << endl;
710 //cout << "mcBtoE : " << mcBtoE << endl;
711 //cout << "mcPho : " << mcPho << endl;
712
713 if(fabs(mcpid)==11)mcele= 0.;
714 //cout << "check e: " << mcele << endl;
bd6ee6dd 715 if(fabs(mcpid)==11 && mcDtoE)mcele= 1.;
cc1605fe 716 //cout << "check D->e: " << mcele << endl;
bd6ee6dd 717 if(fabs(mcpid)==11 && mcBtoE)mcele= 2.;
cc1605fe 718 //cout << "check B->e: " << mcele << endl;
bd6ee6dd 719 if(fabs(mcpid)==11 && mcPho)mcele= 3.;
cc1605fe 720 //cout << "check Pho->e: " << mcele << endl;
bd6ee6dd 721
09022877 722 //cout << "check PID " << endl;
fd2126aa 723 if(fabs(mcpid)!=11)
724 {
09022877 725 //cout << "!= 11" << endl;
726 //cout << mcpid << endl;
fd2126aa 727 }
728 if(mcele==-1)
729 {
09022877 730 //cout << "mcele==-1" << endl;
731 //cout << mcele << endl;
732 //cout << mcpid << endl;
fd2126aa 733 }
734
bd6ee6dd 735 } // end of mcLabel>-1
736
89f41a30 737 } // <------ end of MC info.
e7d87aef 738
2b4460a6 739 //cout << "Pi0 = " << mcOrgPi0 << " ; Eta = " << mcOrgEta << endl;
50919258 740 //printf("weight = %f\n",mcWeight);
741
9b3495ff 742 if(TMath::Abs(track->Eta())>0.6) continue;
fbe402d7 743 //if(TMath::Abs(track->Pt()<2.5)) continue;
744 if(TMath::Abs(track->Pt()<0.1)) continue;
bfc7c23b 745
05ac9660 746 int nITS = track->GetNcls(0);
747
bfc7c23b 748 fTrackPtBefTrkCuts->Fill(track->Pt());
f2111233 749
05ac9660 750 /*
f2fa5a00 751 UChar_t itsPixel = track->GetITSClusterMap();
05ac9660 752 cout << "nITS = " << nITS << endl;
753 if(itsPixel & BIT(0))cout << "1st layer hit" << endl;
754 if(itsPixel & BIT(1))cout << "2nd layer hit" << endl;
755 if(itsPixel & BIT(2))cout << "3rd layer hit" << endl;
756 if(itsPixel & BIT(3))cout << "4th layer hit" << endl;
757 if(itsPixel & BIT(4))cout << "5th layer hit" << endl;
758 */
759
bfc7c23b 760 // RecKine: ITSTPC cuts
761 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
762
763 //RecKink
764 if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
765 if(track->GetKinkIndex(0) != 0) continue;
766 }
767
768 // RecPrim
769 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
770
771 // HFEcuts: ITS layers cuts
6e1d7a8a 772 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
bfc7c23b 773
f09766dd 774 // HFE cuts: TPC PID cleanup
775 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
776
36e9774f 777 if(mcPho && iHijing==0)fPhoVertexReco_step0->Fill(track->Pt(),conv_proR); // check MC vertex
778 if(mcPho && iHijing==1)fPhoVertexReco_step1->Fill(track->Pt(),conv_proR); // check MC vertex
779
a6df418c 780 int nTPCcl = track->GetTPCNcls();
96167a04 781 //int nTPCclF = track->GetTPCNclsF(); // warnings
05ac9660 782 //int nITS = track->GetNcls(0);
783
05ac9660 784
bfc7c23b 785 fTrackPtAftTrkCuts->Fill(track->Pt());
786
787 Double_t mom = -999., eop=-999., pt = -999., dEdx=-999., fTPCnSigma=-10, phi=-999., eta=-999.;
feffe705 788 pt = track->Pt();
d26226a1 789 if(pt<2.5)continue;
fbe402d7 790 if(pt<0.1)continue;
bfc7c23b 791
6f4a9952 792 //Int_t charge = track->Charge();
bfc7c23b 793 fTrkpt->Fill(pt);
794 mom = track->P();
795 phi = track->Phi();
796 eta = track->Eta();
ae88597c 797 float dca_xy;
798 float dca_z;
799 track->GetImpactParameters(dca_xy,dca_z);
800
bfc7c23b 801 dEdx = track->GetTPCsignal();
802 fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
8af58b0d 803
09022877 804 //cout << "nSigma correctoon-----" << endl;
805 //cout << "org = " << fTPCnSigma << endl;
9eae9ea5 806 /*
09022877 807 if(!fmcData) // nsigma eta correction
808 {
809 double nSigexpCorr = NsigmaCorrection(eta,cent);
810 fTPCnSigma -= nSigexpCorr;
811 }
9eae9ea5 812 */
09022877 813 //cout << "correction = " << fTPCnSigma << endl;
814
89f41a30 815 //--- track cluster match
816
bfc7c23b 817 double ncells = -1.0;
818 double m20 = -1.0;
819 double m02 = -1.0;
820 double disp = -1.0;
821 double rmatch = -1.0;
822 double nmatch = -1.0;
89f41a30 823 //double oppstatus = 0.0;
59d998de 824 double emctof = 0.0;
6f4a9952 825 Bool_t MaxEmatch = kFALSE;
bfc7c23b 826
13055933 827 Int_t clsId = track->GetEMCALcluster();
ebf2f603 828
829 //cout << "match ID" << clsId << " ; " << clsId_Tender << endl;
830
831 if (clsId>=0){
bfc7c23b 832 AliESDCaloCluster *clust = fESD->GetCaloCluster(clsId);
ebf2f603 833
bfc7c23b 834 if(clust && clust->IsEMCAL()){
835
836 double clustE = clust->E();
6f4a9952 837 if(clustE==maxE)MaxEmatch = kTRUE;
bfc7c23b 838 eop = clustE/fabs(mom);
8efacc22 839
c5ec1e9b 840
bfc7c23b 841 //double clustT = clust->GetTOF();
842 ncells = clust->GetNCells();
843 m02 = clust->GetM02();
844 m20 = clust->GetM20();
845 disp = clust->GetDispersion();
846 double delphi = clust->GetTrackDx();
847 double deleta = clust->GetTrackDz();
848 rmatch = sqrt(pow(delphi,2)+pow(deleta,2));
849 nmatch = clust->GetNTracksMatched();
59d998de 850 emctof = clust->GetTOF();
851 //cout << "emctof = " << emctof << endl;
ebf2f603 852 //cout << "eop org = "<< eop << endl;
89f41a30 853 double eoporg = eop;
854 if(fmcData)
9eae9ea5 855 {
856 double mceopcorr = MCEopMeanCorrection(pt,cent);
857 eop += mceopcorr;
858 }
ebf2f603 859 //cout << "eop corr = " << eop << endl;
9eae9ea5 860
feffe705 861 double valdedx[16];
a6df418c 862 valdedx[0] = pt; valdedx[1] = nITS; valdedx[2] = phi; valdedx[3] = eta; valdedx[4] = fTPCnSigma;
05ac9660 863 valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = dca_xy, valdedx[8] = dca_z; valdedx[9] = m20; valdedx[10] = mcpT;
89f41a30 864 valdedx[11] = cent; valdedx[12] = dEdx; valdedx[13] = eoporg; valdedx[14] = nTPCcl;
feffe705 865 valdedx[15] = mcele;
89f41a30 866 fEleInfo->Fill(valdedx);
bfc7c23b 867
868 }
869 }
e1f0fb74 870
871 //Get Cal info PID response
9eae9ea5 872 /*
e1f0fb74 873 double eop2;
874 double ss[4];
875 Double_t nSigmaEop = fPID->GetPIDResponse()->NumberOfSigmasEMCAL(track,AliPID::kElectron,eop2,ss);
876 if(fTPCnSigma>-1.5 && fTPCnSigma<3.0 && nITS>2.5 && nTPCcl>100)
877 {
878 double valEop[3];
879 valEop[0] = cent;
880 valEop[1] = pt;
881 valEop[2] = nSigmaEop;
882 fElenSigma->Fill(valEop);
883 }
9eae9ea5 884 */
89f41a30 885 // --- tarck cut & e ID
e1f0fb74 886
a6df418c 887 if(nITS<2.5)continue;
888 if(nTPCcl<100)continue;
7d495dc5 889
7d495dc5 890
a6df418c 891 CheckNclust->Fill(nTPCcl);
892 CheckNits->Fill(nITS);
ae88597c 893 CheckDCA->Fill(dca_xy,dca_z);
f2111233 894 // check production vertex of photons
f2111233 895
a6df418c 896
42c75692 897 fdEdxBef->Fill(mom,fTPCnSigma);
bfc7c23b 898 fTPCnsigma->Fill(mom,fTPCnSigma);
f09766dd 899 if(fTPCnSigma >= -1.0 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,eop);
bfc7c23b 900
89f41a30 901 //--- track accepted by HFE
bfc7c23b 902
89f41a30 903 Int_t pidpassed = 1;
bfc7c23b 904 AliHFEpidObject hfetrack;
905 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
906 hfetrack.SetRecTrack(track);
3db00c72 907 double binct = 10.5;
908 if((0.0< cent) && (cent<5.0)) binct = 0.5;
909 if((5.0< cent) && (cent<10.0)) binct = 1.5;
910 if((10.0< cent) && (cent<20.0)) binct = 2.5;
911 if((20.0< cent) && (cent<30.0)) binct = 3.5;
912 if((30.0< cent) && (cent<40.0)) binct = 4.5;
913 if((40.0< cent) && (cent<50.0)) binct = 5.5;
914 if((50.0< cent) && (cent<60.0)) binct = 6.5;
915 if((60.0< cent) && (cent<70.0)) binct = 7.5;
916 if((70.0< cent) && (cent<80.0)) binct = 8.5;
917 if((80.0< cent) && (cent<90.0)) binct = 9.5;
918 if((90.0< cent) && (cent<100.0)) binct = 10.5;
919
920 hfetrack.SetCentrality((int)binct); //added
bfc7c23b 921 hfetrack.SetPbPb();
922 if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
923
f2fa5a00 924 if(pidpassed==0) continue; // nSigma rejection
89f41a30 925
926 //--- photonic ID
927
0f566acb 928 Bool_t fFlagPhotonicTPC = kFALSE;
929 Bool_t fFlagConvinatTPC = kFALSE;
930 SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicTPC,fFlagConvinatTPC,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta,0);
89f41a30 931
932 //--- check reco eff. with TPC
933 double phoval[5];
934 phoval[0] = cent;
935 phoval[1] = pt;
936 phoval[2] = fTPCnSigma;
937 phoval[3] = iHijing;
938 phoval[4] = mcMompT;
939
2f8ef315 940 if((fTPCnSigma >= -5.0 && fTPCnSigma <= 5) && (mcOrgPi0 || mcOrgEta))
89f41a30 941 {
942 if(iHijing==1)mcWeight = 1.0;
2f8ef315 943 if(mcOrgPi0)
944 {
945 fIncpTMCpho_pi0e_TPC->Fill(phoval,mcWeight);
0f566acb 946 if(fFlagPhotonicTPC) fPhoElecPtMC_pi0e_TPC->Fill(phoval,mcWeight);
947 if(fFlagConvinatTPC) fSameElecPtMC_pi0e_TPC->Fill(phoval,mcWeight);
466d4f7a 948 fPhoVertexReco_TPC->Fill(track->Pt(),conv_proR,mcWeight); // check MC vertex
949 if(fFlagPhotonicTPC)fPhoVertexReco_TPC_Invmass->Fill(track->Pt(),conv_proR,mcWeight); // check MC vertex
2f8ef315 950 }
951 if(mcOrgEta)
952 {
953 fIncpTMCpho_eta_TPC->Fill(phoval,mcWeight);
0f566acb 954 if(fFlagPhotonicTPC) fPhoElecPtMC_eta_TPC->Fill(phoval,mcWeight);
955 if(fFlagConvinatTPC) fSameElecPtMC_eta_TPC->Fill(phoval,mcWeight);
2f8ef315 956 }
89f41a30 957 }
958
0f566acb 959 //--- matching check
f2fa5a00 960 double emcphimim = 1.396;
961 double emcphimax = 3.14;
962 if(phi>emcphimim && phi<emcphimax)
963 {
0f566acb 964 if(fFlagPhotonicTPC)fMatchV0_0->Fill(pt); // data
d26226a1 965 if(mcele>2.1)fMatchMC_0->Fill(pt,mcWeight); // MC
966
967 if(eop>=0.0) // have a match
f2fa5a00 968 {
0f566acb 969 if(fFlagPhotonicTPC)fMatchV0_1->Fill(pt);
d26226a1 970 if(mcele>2.1)fMatchMC_1->Fill(pt,mcWeight);
f2fa5a00 971 }
972 }
973
c5a4ddbf 974 // check fake rejection
975 if(mcOrgPi0 || mcOrgEta)
976 {
1c9a06c8 977 double phiacc0 = 80.0/180.0*acos(-1);
978 double phiacc1 = 180.0/180.0*acos(-1);
c5a4ddbf 979 int TrStat = 0;
1c9a06c8 980 if(phi>phiacc0 && phi<phiacc1)
981 {
982 if(track->GetLabel()>0)TrStat = 1;
983 fFakeRejection0->Fill(TrStat,pt,mcWeight);
984 if(eop>-1.0)fFakeRejection1->Fill(TrStat,pt,mcWeight); // have match
985 if(eop>0.9 && eop<1.3)fFakeRejection2->Fill(TrStat,pt,mcWeight); // have PID
9e88412a 986
987 if(TrStat==0)
988 {
32abf0b9 989 EopFake->Fill(pt,eop,mcWeight);
990 MatchFake->Fill(pt,rmatch,mcWeight);
9e88412a 991 }
992 else
993 {
32abf0b9 994 EopTrue->Fill(pt,eop,mcWeight);
995 MatchTrue->Fill(pt,rmatch,mcWeight);
9e88412a 996 }
997
1c9a06c8 998 }
c5a4ddbf 999 }
1000
89f41a30 1001 //+++++++ E/p cut ++++++++++++++++
1002
1003 if(eop<0.9 || eop>1.3)continue;
0f566acb 1004
1005 Bool_t fFlagPhotonicElec = kFALSE;
1006 Bool_t fFlagConvinatElec = kFALSE;
1007 SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicElec,fFlagConvinatElec,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta,1);
1008
bfc7c23b 1009 fTrkEovPAft->Fill(pt,eop);
42c75692 1010 fdEdxAft->Fill(mom,fTPCnSigma);
bfc7c23b 1011
e7d87aef 1012 // Fill real data
fb87d707 1013 fIncpT->Fill(cent,pt);
bfc7c23b 1014 if(fFlagPhotonicElec) fPhoElecPt->Fill(cent,pt);
f09766dd 1015 if(fFlagConvinatElec) fSameElecPt->Fill(cent,pt);
fb87d707 1016
1017 if(m20>0.0 && m20<0.3)
1018 {
59d998de 1019 fIncpTM20->Fill(cent,pt);
1020 ftimingEle->Fill(pt,emctof);
fb87d707 1021 if(fFlagPhotonicElec) fPhoElecPtM20->Fill(cent,pt);
1022 if(fFlagConvinatElec) fSameElecPtM20->Fill(cent,pt);
ebef1da9 1023 if(!fFlagPhotonicElec) fSemiElecPtM20->Fill(cent,pt);
fb87d707 1024 }
feffe705 1025
5e38d973 1026
e7d87aef 1027 // MC
70448e3b 1028 // check label for electron candidiates
70448e3b 1029 int idlabel = 1;
1030 if(mcLabel==0)idlabel = 0;
1031 fLabelCheck->Fill(pt,idlabel);
1032 if(mcLabel==0)fgeoFake->Fill(phi,eta);
1033
d19af75d 1034 if(mcLabel<0 && m20>0.0 && m20<0.3 && fTPCnSigma>-1 && fTPCnSigma<3)
1035 {
1036 fFakeTrk0->Fill(cent,pt);
1037 }
1038
1039 if(mcele>-1) // select MC electrons
e7d87aef 1040 {
3db00c72 1041
1042 fIncpTMChfeAll->Fill(cent,pt);
1043 if(m20>0.0 && m20<0.3)fIncpTMCM20hfeAll->Fill(cent,pt);
d19af75d 1044 if(m20>0.0 && m20<0.3 && fTPCnSigma>-1 && fTPCnSigma<3)fFakeTrk1->Fill(cent,pt);
3db00c72 1045
e4b0faf2 1046 if(mcBtoE || mcDtoE) // select B->e & D->e
e7d87aef 1047 {
1048 fIncpTMChfe->Fill(cent,pt);
bd6ee6dd 1049 if(m20>0.0 && m20<0.3)
1050 {
09022877 1051 //cout << "MC label = " << mcLabel << endl;
bd6ee6dd 1052 fIncpTMCM20hfe->Fill(cent,pt);
1053 fIncpTMCM20hfeCheck->Fill(cent,mcpT);
1054 fIncpTMCM20hfeCheck_weight->Fill(phoval);
1055 }
e7d87aef 1056 }
e4b0faf2 1057
1058 if(mcPho) // select photonic electrons
e7d87aef 1059 {
44be9e1d 1060
1061 fIncpTMCpho->Fill(phoval);
1062 if(fFlagPhotonicElec) fPhoElecPtMC->Fill(phoval);
1063 if(fFlagConvinatElec) fSameElecPtMC->Fill(phoval);
e7d87aef 1064
1065 if(m20>0.0 && m20<0.3)
1066 {
44be9e1d 1067 fIncpTMCM20pho->Fill(phoval);
1068 if(fFlagPhotonicElec) fPhoElecPtMCM20->Fill(phoval);
1069 if(fFlagConvinatElec) fSameElecPtMCM20->Fill(phoval);
e4b0faf2 1070 // pi0->g->e
1071 if(mcWeight>-1)
1072 {
987053ce 1073 if(iHijing==1)mcWeight = 1.0;
93c189c5 1074 if(mcOrgPi0)
1075 {
697ecf6b 1076 fIncpTMCM20pho_pi0e->Fill(phoval,mcWeight);
1077 if(fFlagPhotonicElec) fPhoElecPtMCM20_pi0e->Fill(phoval,mcWeight);
1078 if(fFlagConvinatElec) fSameElecPtMCM20_pi0e->Fill(phoval,mcWeight);
05ac9660 1079
1080 // check production vertex
7d495dc5 1081 //fPhoVertexReco_EMCal->Fill(track->Pt(),conv_proR);
1082 //if(fFlagPhotonicElec)fPhoVertexReco_Invmass->Fill(track->Pt(),conv_proR);
1083 fPhoVertexReco_EMCal->Fill(track->Pt(),conv_proR,mcWeight);
1084 if(fFlagPhotonicElec)fPhoVertexReco_Invmass->Fill(track->Pt(),conv_proR,mcWeight);
05ac9660 1085
1086 }
bd50bf24 1087 // --- eta
93c189c5 1088 if(mcOrgEta)
1089 {
697ecf6b 1090 fIncpTMCM20pho_eta->Fill(phoval,mcWeight);
1091 if(fFlagPhotonicElec) fPhoElecPtMCM20_eta->Fill(phoval,mcWeight);
1092 if(fFlagConvinatElec) fSameElecPtMCM20_eta->Fill(phoval,mcWeight);
93c189c5 1093 }
bd50bf24 1094 // check production vertex
e4b0faf2 1095 }
e7d87aef 1096 }
1097 }
1098 }
bfc7c23b 1099 }
1100 PostData(1, fOutputList);
1101}
1102//_________________________________________
1103void AliAnalysisTaskHFECal::UserCreateOutputObjects()
1104{
42c75692 1105 //--- Check MC
1106
e7d87aef 1107 //Bool_t mcData = kFALSE;
42c75692 1108 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
1109 {
e7d87aef 1110 fmcData = kTRUE;
42c75692 1111 printf("+++++ MC Data available");
1112 }
e7d87aef 1113 if(fmcData)
42c75692 1114 {
1115 printf("++++++++= MC analysis \n");
1116 }
1117 else
1118 {
1119 printf("++++++++ real data analysis \n");
1120 }
1121
85985bb0 1122 printf("+++++++ QA hist %d",fqahist);
1123
fb87d707 1124 //---- Geometry
1125 fGeom = AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
1126
bfc7c23b 1127 //--------Initialize PID
42c75692 1128 //fPID->SetHasMCData(kFALSE);
e7d87aef 1129 fPID->SetHasMCData(fmcData);
bfc7c23b 1130 if(!fPID->GetNumberOfPIDdetectors())
1131 {
1132 fPID->AddDetector("TPC", 0);
89f41a30 1133 //fPID->AddDetector("EMCAL", 1); <--- apply PID selection in task
bfc7c23b 1134 }
1135
3db00c72 1136 Double_t params[4];
78ff1095 1137 const char *cutmodel;
3db00c72 1138 cutmodel = "pol0";
1139 params[0] = -1.0; //sigma min
9b3495ff 1140 double maxnSig = 3.0;
1141 if(fmcData)
1142 {
1143 params[0] = -5.0; //sigma min
1144 maxnSig = 5.0;
1145 }
1146
1147 for(Int_t a=0;a<11;a++)fPID->ConfigureTPCcentralityCut(a,cutmodel,params,maxnSig);
3db00c72 1148
bfc7c23b 1149 fPID->SortDetectors();
1150 fPIDqa = new AliHFEpidQAmanager();
1151 fPIDqa->Initialize(fPID);
a6df418c 1152
1153 //------- fcut --------------
1154 fCuts = new AliHFEcuts();
1155 fCuts->CreateStandardCuts();
1156 //fCuts->SetMinNClustersTPC(100);
1157 fCuts->SetMinNClustersTPC(90);
1158 fCuts->SetMinRatioTPCclusters(0.6);
1159 fCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
1160 //fCuts->SetMinNClustersITS(3);
1161 fCuts->SetMinNClustersITS(2);
05ac9660 1162 fCuts->SetProductionVertex(0,50,0,50);
a6df418c 1163 fCuts->SetCutITSpixel(AliHFEextraCuts::kAny);
1164 fCuts->SetCheckITSLayerStatus(kFALSE);
1165 fCuts->SetVertexRange(10.);
1166 fCuts->SetTOFPIDStep(kFALSE);
fbe402d7 1167 //fCuts->SetPtRange(2, 50);
1168 fCuts->SetPtRange(0.1, 50);
0f566acb 1169 //fCuts->SetMaxImpactParam(3.,3.);
1170 fCuts->SetMaxImpactParam(2.4,3.2); // standard in 2011
a6df418c 1171
bfc7c23b 1172 //--------Initialize correction Framework and Cuts
1173 fCFM = new AliCFManager;
1174 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
1175 fCFM->SetNStepParticle(kNcutSteps);
1176 for(Int_t istep = 0; istep < kNcutSteps; istep++)
1177 fCFM->SetParticleCutsList(istep, NULL);
1178
1179 if(!fCuts){
1180 AliWarning("Cuts not available. Default cuts will be used");
1181 fCuts = new AliHFEcuts;
1182 fCuts->CreateStandardCuts();
1183 }
1184 fCuts->Initialize(fCFM);
1185
1186 //---------Output Tlist
1187 fOutputList = new TList();
1188 fOutputList->SetOwner();
1189 fOutputList->Add(fPIDqa->MakeList("PIDQA"));
1190
1191 fNoEvents = new TH1F("fNoEvents","",4,-0.5,3.5) ;
1192 fOutputList->Add(fNoEvents);
1193
02cd82a0 1194 Int_t binsE[5] = {250, 100, 1000, 200, 10};
fb87d707 1195 Double_t xminE[5] = {1.0, -1, 0.0, 0, -0.5};
1196 Double_t xmaxE[5] = {3.5, 1, 100.0, 100, 9.5};
1197 fEMCAccE = new THnSparseD("fEMCAccE","EMC acceptance & E;#eta;#phi;Energy;Centrality;trugCondition;",5,binsE,xminE,xmaxE);
c82f1b9d 1198 //if(fqahist==1)fOutputList->Add(fEMCAccE);
bfc7c23b 1199
f4e0d2d5 1200 hEMCAccE = new TH2F("hEMCAccE","Cluster Energy",200,0,100,100,0,20);
1201 fOutputList->Add(hEMCAccE);
1202
bfc7c23b 1203 fTrkpt = new TH1F("fTrkpt","track pt",100,0,50);
1204 fOutputList->Add(fTrkpt);
1205
1206 fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50);
1207 fOutputList->Add(fTrackPtBefTrkCuts);
1208
1209 fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50);
1210 fOutputList->Add(fTrackPtAftTrkCuts);
1211
1212 fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10);
1213 fOutputList->Add(fTPCnsigma);
1214
1215 fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2);
1216 fOutputList->Add(fTrkEovPBef);
1217
1218 fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2);
1219 fOutputList->Add(fTrkEovPAft);
1220
42c75692 1221 fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,200,-10,10);
bfc7c23b 1222 fOutputList->Add(fdEdxBef);
1223
42c75692 1224 fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,200,-10,10);
bfc7c23b 1225 fOutputList->Add(fdEdxAft);
1226
02cd82a0 1227 fIncpT = new TH2F("fIncpT","HFE pid electro vs. centrality",200,0,100,100,0,50);
bfc7c23b 1228 fOutputList->Add(fIncpT);
1229
02cd82a0 1230 fIncpTM20 = new TH2F("fIncpTM20","HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
fb87d707 1231 fOutputList->Add(fIncpTM20);
4e01b68c 1232
5e38d973 1233 Int_t nBinspho[9] = { 10, 30, 600, 60, 50, 4, 40, 8, 30};
1234 Double_t minpho[9] = { 0., 0., -0.1, 40, 0, -0.5, 0,-1.5, 0};
1235 Double_t maxpho[9] = {100., 30., 0.5, 100, 1, 3.5, 2, 6.5, 30};
f09766dd 1236
b12bc641 1237 fInvmassLS = new THnSparseD("fInvmassLS", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2); nSigma; angle; m20cut; eop; Mcele;", 9, nBinspho,minpho, maxpho);
5cc5103f 1238 if(fqahist==1)fOutputList->Add(fInvmassLS);
bfc7c23b 1239
b12bc641 1240 fInvmassULS = new THnSparseD("fInvmassULS", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2); nSigma; angle; m20cut; eop; MCele", 9, nBinspho,minpho, maxpho);
5cc5103f 1241 if(fqahist==1)fOutputList->Add(fInvmassULS);
bfc7c23b 1242
8806ce6c 1243 fInvmassLSmc = new THnSparseD("fInvmassLSmc", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2); nSigma; angle; m20cut; eop; Mcele;", 9, nBinspho,minpho, maxpho);
5cc5103f 1244 if(fqahist==1)fOutputList->Add(fInvmassLSmc);
8806ce6c 1245
1246 fInvmassULSmc = new THnSparseD("fInvmassULSmc", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2); nSigma; angle; m20cut; eop; MCele", 9, nBinspho,minpho, maxpho);
5cc5103f 1247 if(fqahist==1)fOutputList->Add(fInvmassULSmc);
8806ce6c 1248
d89866fd 1249 fInvmassLSreco = new TH2D("fInvmassLSreco", "Inv mass of LS (e,e) reco; cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
6ec70f8c 1250 fInvmassLSreco->Sumw2();
1251 fOutputList->Add(fInvmassLSreco);
1252
d89866fd 1253 fInvmassULSreco = new TH2D("fInvmassULSreco", "Inv mass of ULS (e,e) reco; cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
6ec70f8c 1254 fInvmassULSreco->Sumw2();
1255 fOutputList->Add(fInvmassULSreco);
1256
d89866fd 1257 fInvmassLSmc0 = new TH2D("fInvmassLSmc0", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1258 fInvmassLSmc0->Sumw2();
1259 fOutputList->Add(fInvmassLSmc0);
1260
d89866fd 1261 fInvmassLSmc1 = new TH2D("fInvmassLSmc1", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1262 fInvmassLSmc1->Sumw2();
1263 fOutputList->Add(fInvmassLSmc1);
1264
d89866fd 1265 fInvmassLSmc2 = new TH2D("fInvmassLSmc2", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1266 fInvmassLSmc2->Sumw2();
1267 fOutputList->Add(fInvmassLSmc2);
1268
d89866fd 1269 fInvmassLSmc3 = new TH2D("fInvmassLSmc3", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1270 fInvmassLSmc3->Sumw2();
1271 fOutputList->Add(fInvmassLSmc3);
1272
d89866fd 1273 fInvmassULSmc0 = new TH2D("fInvmassULSmc0", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1274 fInvmassULSmc0->Sumw2();
1275 fOutputList->Add(fInvmassULSmc0);
1276
d89866fd 1277 fInvmassULSmc1 = new TH2D("fInvmassULSmc1", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1278 fInvmassULSmc1->Sumw2();
1279 fOutputList->Add(fInvmassULSmc1);
1280
d89866fd 1281 fInvmassULSmc2 = new TH2D("fInvmassULSmc2", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1282 fInvmassULSmc2->Sumw2();
1283 fOutputList->Add(fInvmassULSmc2);
1284
d89866fd 1285 fInvmassULSmc3 = new TH2D("fInvmassULSmc3", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,600,-0.1,0.5 );
93c189c5 1286 fInvmassULSmc3->Sumw2();
93c189c5 1287 fOutputList->Add(fInvmassULSmc3);
1288
bfc7c23b 1289 fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
1290 fOutputList->Add(fOpeningAngleLS);
1291
1292 fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
1293 fOutputList->Add(fOpeningAngleULS);
1294
1295 fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50);
1296 fOutputList->Add(fPhotoElecPt);
1297
02cd82a0 1298 fPhoElecPt = new TH2F("fPhoElecPt", "Pho-inclusive electron pt",200,0,100,100,0,50);
bfc7c23b 1299 fOutputList->Add(fPhoElecPt);
1300
02cd82a0 1301 fPhoElecPtM20 = new TH2F("fPhoElecPtM20", "Pho-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 1302 fOutputList->Add(fPhoElecPtM20);
1303
0f566acb 1304 fPhoElecPtM20Mass = new TH2F("fPhoElecPtM20Mass", "Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1305 fPhoElecPtM20Mass->Sumw2();
1306 fOutputList->Add(fPhoElecPtM20Mass);
1307
02cd82a0 1308 fSameElecPt = new TH2F("fSameElecPt", "Same-inclusive electron pt",200,0,100,100,0,50);
f09766dd 1309 fOutputList->Add(fSameElecPt);
1310
02cd82a0 1311 fSameElecPtM20 = new TH2F("fSameElecPtM20", "Same-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 1312 fOutputList->Add(fSameElecPtM20);
1313
0f566acb 1314 fSameElecPtM20Mass = new TH2F("fSameElecPtM20Mass", "Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1315 fSameElecPtM20Mass->Sumw2();
1316 fOutputList->Add(fSameElecPtM20Mass);
1317
ebef1da9 1318 fSemiElecPtM20 = new TH2F("fSemiElecPtM20", "Semi-inclusive electron pt with M20",200,0,100,100,0,50);
1319 fOutputList->Add(fSemiElecPtM20);
1320
02cd82a0 1321 fCent = new TH1F("fCent","Centrality",200,0,100) ;
bfc7c23b 1322 fOutputList->Add(fCent);
1323
1324 // Make common binning
44be9e1d 1325 const Double_t kMinP = 0.;
d19af75d 1326 const Double_t kMaxP = 20.;
bfc7c23b 1327
05ac9660 1328 //+++ 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select)
1329 // 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, dca_xy, dca_z, M20, Disp, Centrality, select)
1330 Int_t nBins[16] = { 100, 7, 60, 20, 90, 100, 25, 60, 60, 100, 40, 10, 250, 100, 100, 8};
1331 Double_t min[16] = {kMinP, -0.5, 1.0, -1.0, -5.0, 0, 0, -3.0, -3.0, 0.0, 0.0, 0, 0, 0, 80, -1.5};
1332 Double_t max[16] = {kMaxP, 6.5, 4.0, 1.0, 4.0, 2.0, 0.05, 3.0, 3.0, 1.0, 20.0, 100, 100, 2.0, 180, 6.5};
1333 //fEleInfo = new THnSparseD("fEleInfo", "Electron Info; pT [GeV/c]; TPC signal;phi;eta;nSig; E/p;Rmatch;Ncell;clsF;M20;mcpT;Centrality;charge;opp;same;trigCond;MCele", 16, nBins, min, max);
1334 fEleInfo = new THnSparseD("fEleInfo", "Electron Info; pT [GeV/c]; TPC signal;phi;eta;nSig; E/p;Rmatch;DCA_xy;DCA_z;M20;mcpT;Centrality;charge;opp;same;trigCond;MCele", 16, nBins, min, max);
f2fa5a00 1335 //if(fqahist==1)fOutputList->Add(fEleInfo);
bfc7c23b 1336
e1f0fb74 1337 // Make common binning
1338 Int_t nBinsEop[3] = { 10, 50, 100};
1339 Double_t minEop[3] = { 0, 0, -5};
1340 Double_t maxEop[3] = {100, 50, 5};
1341 fElenSigma= new THnSparseD("fElenSigma", "Electron nSigma; cent; pT [GeV/c]; nSigma", 3, nBinsEop, minEop, maxEop);
1342 fOutputList->Add(fElenSigma);
1343
1344
fb87d707 1345 //<--- Trigger info
feffe705 1346 /*
fb87d707 1347 fClsEBftTrigCut = new TH1F("fClsEBftTrigCut","cluster E before trigger selection",1000,0,100);
1348 fOutputList->Add(fClsEBftTrigCut);
1349
1350 fClsEAftTrigCut = new TH1F("fClsEAftTrigCut","cluster E if cls has 0 trigcut channel",1000,0,100);
1351 fOutputList->Add(fClsEAftTrigCut);
1352
1353 fClsEAftTrigCut1 = new TH1F("fClsEAftTrigCut1","cluster E if cls with trig channel",1000,0,100);
1354 fOutputList->Add(fClsEAftTrigCut1);
1355
1356 fClsEAftTrigCut2 = new TH1F("fClsEAftTrigCut2","cluster E if cls with trigcut channel",1000,0,100);
1357 fOutputList->Add(fClsEAftTrigCut2);
1358
1359 fClsEAftTrigCut3 = new TH1F("fClsEAftTrigCut3","cluster E if cls with trigcut channel + nCell>Ecorrect",1000,0,100);
1360 fOutputList->Add(fClsEAftTrigCut3);
1361
1362 fClsEAftTrigCut4 = new TH1F("fClsEAftTrigCut4","cluster E if cls with trigcut channel + nCell>Ecorrect + cls time cut",1000,0,100);
1363 fOutputList->Add(fClsEAftTrigCut4);
1364
1365 fClsETime = new TH2F("fClsETime", "Cls time vs E; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
1366 fOutputList->Add(fClsETime);
1367
1368 fClsETime1 = new TH2F("fClsETime1", "Cls time vs E if cls contains trigger channel; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
1369 fOutputList->Add(fClsETime1);
1370
1371 fTrigTimes = new TH1F("fTrigTimes", "Trigger time; time; N;",25,0,25);
1372 fOutputList->Add(fTrigTimes);
1373
42c75692 1374 fCellCheck = new TH2F("fCellCheck", "Cell vs E; E GeV; Cell ID",10,6,26,12000,0,12000);
1375 fOutputList->Add(fCellCheck);
feffe705 1376 */
e7d87aef 1377 //<---------- MC
1378
02cd82a0 1379 fInputHFEMC = new TH2F("fInputHFEMC","Input MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 1380 fOutputList->Add(fInputHFEMC);
1381
02cd82a0 1382 fInputAlle = new TH2F("fInputAlle","Input MC electro vs. centrality",200,0,100,100,0,50);
feffe705 1383 fOutputList->Add(fInputAlle);
1384
02cd82a0 1385 fIncpTMChfe = new TH2F("fIncpTMChfe","MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 1386 fOutputList->Add(fIncpTMChfe);
1387
02cd82a0 1388 fIncpTMChfeAll = new TH2F("fIncpTMChfeAll","MC Alle pid electro vs. centrality",200,0,100,100,0,50);
3db00c72 1389 fOutputList->Add(fIncpTMChfeAll);
1390
02cd82a0 1391 fIncpTMCM20hfe = new TH2F("fIncpTMCM20hfe","MC HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
e7d87aef 1392 fOutputList->Add(fIncpTMCM20hfe);
1393
02cd82a0 1394 fIncpTMCM20hfeAll = new TH2F("fIncpTMCM20hfeAll","MC Alle pid electro vs. centrality with M20",200,0,100,100,0,50);
3db00c72 1395 fOutputList->Add(fIncpTMCM20hfeAll);
1396
60544aea 1397 fIncpTMCM20hfeCheck = new TH2F("fIncpTMCM20hfeCheck","MC HFE pid electro vs. centrality with M20 Check",200,0,100,100,0,50);
1398 fOutputList->Add(fIncpTMCM20hfeCheck);
44be9e1d 1399
2b4460a6 1400 Int_t nBinspho2[5] = { 200, 100, 7, 3, 700};
987053ce 1401 Double_t minpho2[5] = { 0., 0., -2.5, -0.5, 0.};
2b4460a6 1402 Double_t maxpho2[5] = {100., 50., 4.5, 2.5, 70.};
bd6ee6dd 1403
1404 fInputHFEMC_weight = new THnSparseD("fInputHFEMC_weight", "MC HFE electron pt",5,nBinspho2,minpho2,maxpho2);
1405 fOutputList->Add(fInputHFEMC_weight);
1406
1407 fIncpTMCM20hfeCheck_weight = new THnSparseD("fIncpTMCM20hfeCheck_weight", "HFE electron pt with M20",5,nBinspho2,minpho2,maxpho2);
1408 fOutputList->Add(fIncpTMCM20hfeCheck_weight);
1409
987053ce 1410 fIncpTMCpho = new THnSparseD("fIncpTMCpho","MC Pho pid electro vs. centrality",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1411 fOutputList->Add(fIncpTMCpho);
1412
987053ce 1413 fIncpTMCM20pho = new THnSparseD("fIncpTMCM20pho","MC Pho pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1414 fOutputList->Add(fIncpTMCM20pho);
1415
987053ce 1416 fPhoElecPtMC = new THnSparseD("fPhoElecPtMC", "MC Pho-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1417 fOutputList->Add(fPhoElecPtMC);
1418
987053ce 1419 fPhoElecPtMCM20 = new THnSparseD("fPhoElecPtMCM20", "MC Pho-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
feffe705 1420 fOutputList->Add(fPhoElecPtMCM20);
e7d87aef 1421
0f566acb 1422 fPhoElecPtMCM20Mass = new TH2D("fPhoElecPtMCM20Mass", "MC Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1423 fOutputList->Add(fPhoElecPtMCM20Mass);
1424
987053ce 1425 fSameElecPtMC = new THnSparseD("fSameElecPtMC", "MC Same-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1426 fOutputList->Add(fSameElecPtMC);
1427
987053ce 1428 fSameElecPtMCM20 = new THnSparseD("fSameElecPtMCM20", "MC Same-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1429 fOutputList->Add(fSameElecPtMCM20);
1430
0f566acb 1431 fSameElecPtMCM20Mass = new TH2D("fSameElecPtMCM20Mass", "MC Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1432 fOutputList->Add(fSameElecPtMCM20Mass);
1433
987053ce 1434 fIncpTMCM20pho_pi0e = new THnSparseD("fIncpTMCM20pho_pi0e","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1435 fIncpTMCM20pho_pi0e->Sumw2();
e4b0faf2 1436 fOutputList->Add(fIncpTMCM20pho_pi0e);
1437
987053ce 1438 fPhoElecPtMCM20_pi0e = new THnSparseD("fPhoElecPtMCM20_pi0e", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1439 fPhoElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1440 fOutputList->Add(fPhoElecPtMCM20_pi0e);
1441
987053ce 1442 fSameElecPtMCM20_pi0e = new THnSparseD("fSameElecPtMCM20_pi0e", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1443 fSameElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1444 fOutputList->Add(fSameElecPtMCM20_pi0e);
697ecf6b 1445 //
93c189c5 1446 fIncpTMCM20pho_eta = new THnSparseD("fIncpTMCM20pho_eta","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1447 fIncpTMCM20pho_eta->Sumw2();
1448 fOutputList->Add(fIncpTMCM20pho_eta);
1449
1450 fPhoElecPtMCM20_eta = new THnSparseD("fPhoElecPtMCM20_eta", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1451 fPhoElecPtMCM20_eta->Sumw2();
1452 fOutputList->Add(fPhoElecPtMCM20_eta);
1453
1454 fSameElecPtMCM20_eta = new THnSparseD("fSameElecPtMCM20_eta", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1455 fSameElecPtMCM20_eta->Sumw2();
1456 fOutputList->Add(fSameElecPtMCM20_eta);
697ecf6b 1457 // ------------
1458 fIncpTMCpho_pi0e_TPC = new THnSparseD("fIncpTMCpho_pi0e_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1459 fIncpTMCpho_pi0e_TPC->Sumw2();
1460 fOutputList->Add(fIncpTMCpho_pi0e_TPC);
1461
1462 fPhoElecPtMC_pi0e_TPC = new THnSparseD("fPhoElecPtMC_pi0e_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2);
1463 fPhoElecPtMC_pi0e_TPC->Sumw2();
1464 fOutputList->Add(fPhoElecPtMC_pi0e_TPC);
1465
1466 fSameElecPtMC_pi0e_TPC = new THnSparseD("fSameElecPtMC_pi0e_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1467 fSameElecPtMC_pi0e_TPC->Sumw2();
1468 fOutputList->Add(fSameElecPtMC_pi0e_TPC);
2f8ef315 1469
1470 fIncpTMCpho_eta_TPC = new THnSparseD("fIncpTMCpho_eta_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1471 fIncpTMCpho_eta_TPC->Sumw2();
1472 fOutputList->Add(fIncpTMCpho_eta_TPC);
1473
1474 fPhoElecPtMC_eta_TPC = new THnSparseD("fPhoElecPtMC_eta_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2);
1475 fPhoElecPtMC_eta_TPC->Sumw2();
1476 fOutputList->Add(fPhoElecPtMC_eta_TPC);
1477
1478 fSameElecPtMC_eta_TPC = new THnSparseD("fSameElecPtMC_eta_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1479 fSameElecPtMC_eta_TPC->Sumw2();
1480 fOutputList->Add(fSameElecPtMC_eta_TPC);
1481
1482
697ecf6b 1483 //-------------
e4b0faf2 1484
a6df418c 1485 CheckNclust = new TH1D("CheckNclust","cluster check",200,0,200);
1486 fOutputList->Add(CheckNclust);
1487
1488 CheckNits = new TH1D("CheckNits","ITS cluster check",8,-0.5,7.5);
1489 fOutputList->Add(CheckNits);
ae88597c 1490
1491 CheckDCA = new TH2D("CheckDCA","DCA check",200,-5,5,200,-5,5);
1492 fOutputList->Add(CheckDCA);
68d718e7 1493 /*
d113d7cd 1494 Hpi0pTcheck = new TH2D("Hpi0pTcheck","Pi0 pT from Hijing",100,0,50,3,-0.5,2.5);
1495 fOutputList->Add(Hpi0pTcheck);
1496
19325033 1497 HETApTcheck = new TH2D("HETApTcheck","Eta pT from Hijing",100,0,50,3,-0.5,2.5);
1498 fOutputList->Add(HETApTcheck);
68d718e7 1499 */
1500
1501 Int_t nBinspho3[3] = { 200, 100, 3};
1502 Double_t minpho3[3] = { 0., 0., -0.5};
1503 Double_t maxpho3[3] = {100., 50., 2.5};
19325033 1504
68d718e7 1505 Hpi0pTcheck = new THnSparseD("Hpi0pTcheck","Pi0 pT from Hijing",3,nBinspho3,minpho3,maxpho3);
1506 fOutputList->Add(Hpi0pTcheck);
1507
1508 HETApTcheck = new THnSparseD("HETApTcheck","Eta pT from Hijing",3,nBinspho3,minpho3,maxpho3);
1509 fOutputList->Add(HETApTcheck);
1510 //--
e4b0faf2 1511 HphopTcheck = new TH2D("HphopTcheck","Pho pT from Hijing",100,0,50,3,-0.5,2.5);
1512 fOutputList->Add(HphopTcheck);
5e0e45b3 1513 //
1514 HDpTcheck = new TH2D("HDpTcheck","D pT from Hijing",100,0,50,3,-0.5,2.5);
1515 fOutputList->Add(HDpTcheck);
1516
1517 HBpTcheck = new TH2D("HBpTcheck","B pT from Hijing",100,0,50,3,-0.5,2.5);
1518 fOutputList->Add(HBpTcheck);
1519 //
e4b0faf2 1520
93c189c5 1521 fpTCheck = new TH1D("fpTCheck","pT check",500,0,50);
1522 fOutputList->Add(fpTCheck);
1523
50919258 1524 fMomDtoE = new TH2D("fMomDtoE","D->E pT correlations;e p_{T} GeV/c;D p_{T} GeV/c",400,0,40,400,0,40);
1525 fOutputList->Add(fMomDtoE);
d113d7cd 1526
70448e3b 1527 fLabelCheck = new TH2D("fLabelCheck","MC label",50,0,50,5,-1.5,3.5);
1528 fOutputList->Add(fLabelCheck);
1529
1530 fgeoFake = new TH2D("fgeoFake","Label==0 eta and phi",628,0,6.28,200,-1,1);
1531 fOutputList->Add(fgeoFake);
1532
d19af75d 1533 fFakeTrk0 = new TH2D("fFakeTrk0","fake trakcs",10,0,100,20,0,20);
1534 fOutputList->Add(fFakeTrk0);
1535
1536 fFakeTrk1 = new TH2D("fFakeTrk1","true all e a.f. eID",10,0,100,20,0,20);
1537 fOutputList->Add(fFakeTrk1);
1538
59d998de 1539 ftimingEle = new TH2D("ftimingEle","electron TOF",100,0,20,100,1e-7,1e-6);
1540 fOutputList->Add(ftimingEle);
1541
09022877 1542 // eta correction
1543 // note: parameters 01/31new.TPCnSigmaEtaDep
1544 // 70-90 delta_eta = 0.2
1545
1546 double etaval[12] = {-0.55,-0.45,-0.35,-0.25,-0.15,-0.05,0.05,0.15,0.25,0.35,0.45,0.55};
1547 double corr0[12]= {-0.569177,-0.528844,-0.391979,-0.165494,0.0283495,0.156171,0.266353,0.13103,-0.0250842,-0.274089,-0.45481,-0.536291}; // 0-10 (done)
1548 double corr1[12]= {-0.404742,-0.278953,-0.218069,0.00139927,0.191412,0.354403,0.524594,0.341778,0.244199,-0.112146,-0.160692,-0.352832}; // 10-20 (done)
1549 double corr2[12] = {-0.306007,-0.16821,-0.0248635,0.202233,0.447051,0.497197,0.712251,0.433482,0.337907,0.168426,-0.0693229,-0.0728351}; // 20-30 (done)
1550 double corr3[12] = {-0.13884,-0.0503553,0.104403,0.389773,0.50697,0.539048,0.751642,0.655636,0.518563,0.308156,0.0361159,-0.0491439}; // 30-40 (done)
1551 double corr4[12] = {-0.0319431,0.0808711,0.208774,0.443217,0.557762,0.61453,0.889519,0.808282,0.620394,0.267092,0.15241,-0.0458664}; // 40-50 (done)
1552 double corr5[12] = {-0.130625,0.0189124,0.190344,0.467431,0.546353,0.672251,0.731541,0.802101,0.437108,0.294081,0.193682,0.159074}; // 50-70(done)
1553 double corr6[12] = {0.0600197,0.0600197,0.358366,0.358366,0.973734,0.973734,0.759812,0.759812,0.667861,0.667861,0.415635,0.415635}; // 70-90(done)
1554
1555 fnSigEtaCorr[0] = new TGraphErrors(12,etaval,corr0); // 0-10
1556 fnSigEtaCorr[1] = new TGraphErrors(12,etaval,corr1); // 10-20
1557 fnSigEtaCorr[2] = new TGraphErrors(12,etaval,corr2); // 20-30
1558 fnSigEtaCorr[3] = new TGraphErrors(12,etaval,corr3); // 30-40
1559 fnSigEtaCorr[4] = new TGraphErrors(12,etaval,corr4); // 40-50
1560 fnSigEtaCorr[5] = new TGraphErrors(12,etaval,corr5); // 50-70
1561 fnSigEtaCorr[6] = new TGraphErrors(12,etaval,corr6); // 70-90
1562
6f4a9952 1563 fIncMaxE = new TH2D("fIncMaxE","Inc",10,0,100,10,0,100);
1564 fOutputList->Add(fIncMaxE);
09022877 1565
c82f1b9d 1566 fIncReco = new TH2D("fIncReco","Inc",10,0,100,100,0,500);
5e38d973 1567 fOutputList->Add(fIncReco);
1568
c82f1b9d 1569 fPhoReco = new TH2D("fPhoReco","Pho",10,0,100,100,0,500);
1570 fOutputList->Add(fPhoReco);
1571
1572 fSamReco = new TH2D("fSamReco","Same",10,0,100,100,0,500);
1573 fOutputList->Add(fSamReco);
1574
1575 fIncRecoMaxE = new TH2D("fIncRecoMaxE","Inc",10,0,100,100,0,500);
1576 fOutputList->Add(fIncRecoMaxE);
1577
1578 fPhoRecoMaxE = new TH2D("fPhoRecoMaxE","Pho",10,0,100,100,0,500);
eabc1527 1579 fOutputList->Add(fPhoRecoMaxE);
5e38d973 1580
c82f1b9d 1581 fSamRecoMaxE = new TH2D("fSamRecoMaxE","Same",10,0,100,100,0,500);
eabc1527 1582 fOutputList->Add(fSamRecoMaxE);
5e38d973 1583
466d4f7a 1584 fPhoVertexReco_TPC = new TH2D("fPhoVertexReco_TPC","photon production Vertex mass selection TPC",40,0,20,250,0,50);
1585 fPhoVertexReco_TPC->Sumw2();
1586 fOutputList->Add(fPhoVertexReco_TPC);
1587
1588 fPhoVertexReco_TPC_Invmass = new TH2D("fPhoVertexReco_TPC_Invmass","photon production Vertex mass selection TPC Invmass",40,0,20,250,0,50);
1589 fPhoVertexReco_TPC_Invmass->Sumw2();
1590 fOutputList->Add(fPhoVertexReco_TPC_Invmass);
05ac9660 1591
1592 fPhoVertexReco_EMCal = new TH2D("fPhoVertexReco_EMCal","photon production Vertex mass selection",40,0,20,250,0,50);
7d495dc5 1593 fPhoVertexReco_EMCal->Sumw2();
05ac9660 1594 fOutputList->Add(fPhoVertexReco_EMCal);
1595
1596 fPhoVertexReco_Invmass = new TH2D("fPhoVertexReco_Invmass","photon production Vertex mass selection",40,0,20,250,0,50);
7d495dc5 1597 fPhoVertexReco_Invmass->Sumw2();
05ac9660 1598 fOutputList->Add(fPhoVertexReco_Invmass);
1599
1600 fPhoVertexReco_step0= new TH2D("fPhoVertexReco_step0","photon production Vertex mass selection",40,0,20,250,0,50);
d26226a1 1601 fPhoVertexReco_step0->Sumw2();
05ac9660 1602 fOutputList->Add(fPhoVertexReco_step0);
1603
1604 fPhoVertexReco_step1= new TH2D("fPhoVertexReco_step1","photon production Vertex mass selection",40,0,20,250,0,50);
d26226a1 1605 fPhoVertexReco_step1->Sumw2();
05ac9660 1606 fOutputList->Add(fPhoVertexReco_step1);
1607
f2fa5a00 1608 fMatchV0_0 = new TH1D("fMatchV0_0","V0 match",100,0,20);
1609 fOutputList->Add(fMatchV0_0);
1610
1611 fMatchV0_1 = new TH1D("fMatchV0_1","V0 match",100,0,20);
1612 fOutputList->Add(fMatchV0_1);
1613
1614 fMatchMC_0 = new TH1D("fMatchMC_0","MC match",100,0,20);
1615 fOutputList->Add(fMatchMC_0);
1616
1617 fMatchMC_1 = new TH1D("fMatchMC_1","MC match",100,0,20);
1618 fOutputList->Add(fMatchMC_1);
1619
36e9774f 1620 fpair = new TH2D("fpair","pair of associate",100,0,20,21,-10.5,10.5);
1621 fOutputList->Add(fpair);
1622
c5a4ddbf 1623 fFakeRejection0 = new TH2D("fFakeRejection0","TPC PID",2,-0.5,1.5,100,0,20);
1c9a06c8 1624 fFakeRejection0->Sumw2();
c5a4ddbf 1625 fOutputList->Add(fFakeRejection0);
1626
1627 fFakeRejection1 = new TH2D("fFakeRejection1","TPC PID + Tr match",2,-0.5,1.5,100,0,20);
1c9a06c8 1628 fFakeRejection1->Sumw2();
c5a4ddbf 1629 fOutputList->Add(fFakeRejection1);
1630
1631 fFakeRejection2 = new TH2D("fFakeRejection2","TPC PID + Tr match + E/p",2,-0.5,1.5,100,0,20);
1c9a06c8 1632 fFakeRejection2->Sumw2();
c5a4ddbf 1633 fOutputList->Add(fFakeRejection2);
1634
32abf0b9 1635 EopFake = new TH2D("EopFake","negative track Eop",20,0,20,200,0,2);
9e88412a 1636 EopFake->Sumw2();
1637 fOutputList->Add(EopFake);
1638
32abf0b9 1639 EopTrue = new TH2D("EopTrue","true track Eop",20,0,20,200,0,2);
9e88412a 1640 EopTrue->Sumw2();
1641 fOutputList->Add(EopTrue);
1642
32abf0b9 1643 MatchFake = new TH2D("MatchFake","negative track Match",20,0,20,100,0,0.05);
9e88412a 1644 MatchFake->Sumw2();
1645 fOutputList->Add(MatchFake);
1646
93c64cd2 1647 MatchTrue = new TH2D("MatchTrue","true track Match",20,0,20,100,0,0.05);
9e88412a 1648 MatchTrue->Sumw2();
1649 fOutputList->Add(MatchTrue);
1650
ebf2f603 1651 MatchTrCheck = new TH2D("MatchTrCheck","Check Tr Match",2010,-10,2000,2010,-10,2000);
1652 fOutputList->Add(MatchTrCheck);
1653
1654 MatchTrEop = new TH2D("MatchTrEop","Check Tr Match E/p",200,0,2,200,0,2);
1655 fOutputList->Add(MatchTrEop);
1656
bfc7c23b 1657 PostData(1,fOutputList);
1658}
1659
1660//________________________________________________________________________
1661void AliAnalysisTaskHFECal::Terminate(Option_t *)
1662{
1663 // Info("Terminate");
1664 AliAnalysisTaskSE::Terminate();
1665}
1666
1667//________________________________________________________________________
1668Bool_t AliAnalysisTaskHFECal::ProcessCutStep(Int_t cutStep, AliVParticle *track)
1669{
1670 // Check single track cuts for a given cut step
1671 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
1672 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
1673 return kTRUE;
1674}
1675//_________________________________________
0f566acb 1676void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonic, Bool_t &fFlagConvinat, Double_t nSig, Double_t shower, Double_t ep, Double_t mce, Double_t w, Int_t ibgevent, Bool_t tagpi0, Bool_t tageta, Int_t iCal)
bfc7c23b 1677{
1678 //Identify non-heavy flavour electrons using Invariant mass method
1679
fb87d707 1680 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
feffe705 1681 fTrackCuts->SetRequireTPCRefit(kTRUE);
85985bb0 1682 fTrackCuts->SetRequireITSRefit(kTRUE);
feffe705 1683 fTrackCuts->SetEtaRange(-0.9,0.9);
bfc7c23b 1684 fTrackCuts->SetMaxChi2PerClusterTPC(3.5);
feffe705 1685 fTrackCuts->SetMinNClustersTPC(90);
bfc7c23b 1686
c6bca3f2 1687 //const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
6451f693 1688 Double_t bfield = fESD->GetMagneticField();
2f8ef315 1689 Double_t emass = 0.51*0.001; // (0.51 MeV)
6451f693 1690
93c189c5 1691 double ptEle = track->Pt(); //add
1692 if(ibgevent==0 && w > 0.0)
1693 {
1694 fpTCheck->Fill(ptEle,w);
1695 }
1696
bfc7c23b 1697 Bool_t flagPhotonicElec = kFALSE;
f09766dd 1698 Bool_t flagConvinatElec = kFALSE;
bfc7c23b 1699
8806ce6c 1700 int p1 = 0;
1701 if(mce==3)
1702 {
1703 Int_t label = TMath::Abs(track->GetLabel());
1704 TParticle* particle = stack->Particle(label);
1705 p1 = particle->GetFirstMother();
1706 }
1707
36e9774f 1708 int numULS = 0;
1709 int numLS = 0;
1710
5e65985c 1711 for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){
bfc7c23b 1712 AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
1713 if (!trackAsso) {
1714 printf("ERROR: Could not receive track %d\n", jTracks);
1715 continue;
1716 }
5e65985c 1717 if(itrack==jTracks)continue;
d86a815c 1718 int jbgevent = 0;
5e65985c 1719
8806ce6c 1720 int p2 = 0;
1721 if(mce==3)
1722 {
1723 Int_t label2 = TMath::Abs(trackAsso->GetLabel());
1724 TParticle* particle2 = stack->Particle(label2);
d86a815c 1725 Bool_t MChijing_ass = fMC->IsFromBGEvent(label2);
1726 if(MChijing_ass)jbgevent =1;
8806ce6c 1727 if(particle2->GetFirstMother()>-1)
1728 p2 = particle2->GetFirstMother();
1729 }
1730
f09766dd 1731 Double_t dEdxAsso = -999., ptPrim=-999., ptAsso=-999., openingAngle = -999.;
bfc7c23b 1732 Double_t mass=999., width = -999;
1733 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
1734
93c189c5 1735 //ptPrim = track->Pt();
1736 ptPrim = ptEle;
f09766dd 1737
bfc7c23b 1738 dEdxAsso = trackAsso->GetTPCsignal();
1739 ptAsso = trackAsso->Pt();
1740 Int_t chargeAsso = trackAsso->Charge();
1741 Int_t charge = track->Charge();
1742
b12bc641 1743
6451f693 1744 if(ptAsso <fMimpTassCut) continue;
bfc7c23b 1745 if(!fTrackCuts->AcceptTrack(trackAsso)) continue;
6b2cee73 1746 double fTPCnSigmaAss = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(trackAsso, AliPID::kElectron) : 1000;
1747 //cout << "fTPCnSigmaAss = " << fTPCnSigmaAss << endl;
1748 //cout << "fTPCnSigmaAss Cut = " << fMimNsigassCut << endl;
1749 if(fTPCnSigmaAss <fMimNsigassCut || fTPCnSigmaAss>5) continue;
1750 //cout << "fTPCnSigmaAss a.f. cut = " << fTPCnSigmaAss << endl;
bfc7c23b 1751
1752 Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
1753 if(charge>0) fPDGe1 = -11;
1754 if(chargeAsso>0) fPDGe2 = -11;
b12bc641 1755
1756 //printf("chargeAsso = %d\n",chargeAsso);
1757 //printf("charge = %d\n",charge);
bfc7c23b 1758 if(charge == chargeAsso) fFlagLS = kTRUE;
1759 if(charge != chargeAsso) fFlagULS = kTRUE;
1760
b12bc641 1761 //printf("fFlagLS = %d\n",fFlagLS);
1762 //printf("fFlagULS = %d\n",fFlagULS);
cccab38f 1763 //printf("\n");
b12bc641 1764
6451f693 1765 AliKFParticle::SetField(bfield);
bfc7c23b 1766 AliKFParticle ge1(*track, fPDGe1);
1767 AliKFParticle ge2(*trackAsso, fPDGe2);
2f8ef315 1768
1769 if(fSetMassNonlinear)
1770 {
1771 ge1.SetNonlinearMassConstraint(emass);
1772 ge2.SetNonlinearMassConstraint(emass);
1773 }
1774
bfc7c23b 1775 AliKFParticle recg(ge1, ge2);
1776
2f8ef315 1777 /*
d19af75d 1778 // vertex
bfc7c23b 1779 AliKFVertex primV(*pVtx);
1780 primV += recg;
8079a103 1781 primV -= ge1;
1782 primV -= ge2;
bfc7c23b 1783 recg.SetProductionVertex(primV);
2f8ef315 1784 */
1785
d19af75d 1786 // check chi2
1787 if(recg.GetNDF()<1) continue;
a009053a 1788 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
1789 Double_t chi2cut = 3.0;
f9447b3b 1790
a009053a 1791 // mass.
1792 if(fSetMassConstraint)
1793 {
1794 recg.SetMassConstraint(0,0.0001);
1795 chi2cut = 30.0;
1796 }
f9447b3b 1797 recg.GetMass(mass,width);
1798
3c475a9d 1799 if(fSetMassWidthCut && width>1e10)continue;
1800
a009053a 1801 // angle
bfc7c23b 1802 openingAngle = ge1.GetAngle(ge2);
1803 if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
1804 if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
0f566acb 1805
c2801a73 1806 double ishower = 0;
1807 if(shower>0.0 && shower<0.3)ishower = 1;
1808
0f566acb 1809 if(iCal==1)
1810 {
1811 double phoinfo[9];
1812 phoinfo[0] = cent;
1813 phoinfo[1] = ptPrim;
1814 phoinfo[2] = mass;
1815 phoinfo[3] = nSig;
1816 phoinfo[4] = openingAngle;
1817 phoinfo[5] = ishower;
1818 phoinfo[6] = ep;
1819 phoinfo[7] = mce;
1820 phoinfo[8] = ptAsso;
1821
1822 if(fFlagLS) fInvmassLS->Fill(phoinfo);
1823 if(fFlagULS) fInvmassULS->Fill(phoinfo);
1824 if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc->Fill(phoinfo,w);
1825 if(fFlagULS && ibgevent==0 && jbgevent==0)
1826 {
1827 fInvmassULSmc->Fill(phoinfo,w);
1828 }
1829 //printf("fInvmassCut %f\n",fInvmassCut);
1830 //printf("openingAngle %f\n",fOpeningAngleCut);
1831 }
b12bc641 1832
f9447b3b 1833 // angle cut
ffa14dda 1834 if(openingAngle > fOpeningAngleCut) continue;
a009053a 1835 // chi2 cut
6ec70f8c 1836 //if(TMath::Sqrt(TMath::Abs(chi2recg))>chi2cut) continue;
1837 if(chi2recg>chi2cut) continue;
1838
1839 if(fFlagLS ) fInvmassLSreco->Fill(ptPrim,mass);
1840 if(fFlagULS) fInvmassULSreco->Fill(ptPrim,mass);
1841
36e9774f 1842 // check double count
1843 if(mass<fInvmassCut && fFlagULS)numULS++;
1844 if(mass<fInvmassCut && fFlagLS)numLS++;
1845
d86a815c 1846 // for real data
93c189c5 1847 //printf("mce =%f\n",mce);
0f566acb 1848
d86a815c 1849 if(mce<-0.5) // mce==-1. is real
0f566acb 1850 {
1851 //printf("Real data\n");
1852
1853 // count
1854 if(mass<fInvmassCut && fFlagULS && shower>0.0 && shower<0.3 && iCal==1)fPhoElecPtM20Mass->Fill(cent,ptPrim);
1855 if(mass<fInvmassCut && fFlagLS && shower>0.0 && shower<0.3 && iCal==1)fSameElecPtM20Mass->Fill(cent,ptPrim);
1856
1857 // flag
1858 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
1859 flagPhotonicElec = kTRUE;
1860 }
1861 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec){
1862 flagConvinatElec = kTRUE;
1863 }
1864 }
d86a815c 1865 // for MC data
1866 else
0f566acb 1867 {
1868 // count
1869 if(mass<fInvmassCut && fFlagULS && shower>0.0 && shower<0.3 && mce>2.1 && iCal==1)fPhoElecPtMCM20Mass->Fill(cent,ptPrim,w);
1870 if(mass<fInvmassCut && fFlagLS && shower>0.0 && shower<0.3 && mce>2.1 && iCal==1)fSameElecPtMCM20Mass->Fill(cent,ptPrim,w);
1871
1872 //printf("MC data\n");
1873
1874 if(w>0.0)
1875 {
1876 //cout << "tagpi0 = " << tagpi0 << " ; tageta = " << tageta << endl;
1877 if(iCal==1)
1878 {
1879 if(fFlagLS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassLSmc0->Fill(ptPrim,mass);
1880 if(fFlagULS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassULSmc0->Fill(ptPrim,mass);
1881 if(fFlagLS && ibgevent==0 && jbgevent==0 && tageta) fInvmassLSmc1->Fill(ptPrim,mass);
1882 if(fFlagULS && ibgevent==0 && jbgevent==0 && tageta) fInvmassULSmc1->Fill(ptPrim,mass);
1883 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassLSmc2->Fill(ptPrim,mass);
1884 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassULSmc2->Fill(ptPrim,mass);
1885 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassLSmc3->Fill(ptPrim,mass);
1886 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassULSmc3->Fill(ptPrim,mass);
1887 }
1888 }
1889
1890 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (ibgevent==jbgevent)){
1891 //if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (p1==p2)){ //<--- only MC train (55,56) v5-03-68-AN , 69 & v5-05-70-AN (till 74)
1892 flagPhotonicElec = kTRUE;
1893 }
1894 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec && (ibgevent==jbgevent)){
1895 flagConvinatElec = kTRUE;
1896 }
1897 }
d86a815c 1898
36e9774f 1899 } // end of associate loop
1900
6e1d7a8a 1901 if(numULS>0 || numLS>0)
1902 {
1903 int numPair = numULS-numLS;
0f566acb 1904 if(iCal==1)fpair->Fill(ptEle,numPair);
6e1d7a8a 1905 }
1906
0f566acb 1907 fFlagPhotonic = flagPhotonicElec;
1908 fFlagConvinat = flagConvinatElec;
6451f693 1909
1910}
1911
8806ce6c 1912//-------------------------------------------
46305ed6 1913
1914void AliAnalysisTaskHFECal::FindMother(TParticle* part, int &label, int &pid)
1915{
1916 //int label = 99999;
1917 //int pid = 99999;
1918
1919 if(part->GetFirstMother()>-1)
1920 {
1921 label = part->GetFirstMother();
1922 pid = stack->Particle(label)->GetPdgCode();
1923 }
70448e3b 1924 //cout << "Find Mother : label = " << label << " ; pid" << pid << endl;
46305ed6 1925}
1926
8806ce6c 1927double AliAnalysisTaskHFECal::GetMCweight(double mcPi0pT)
1928{
1929 double weight = 1.0;
1930
1931 if(mcPi0pT>0.0 && mcPi0pT<5.0)
1932 {
1933 weight = 0.323*mcPi0pT/(TMath::Exp(-1.6+0.767*mcPi0pT+0.0285*mcPi0pT*mcPi0pT));
1934 }
1935 else
1936 {
1937 weight = 115.0/(0.718*mcPi0pT*TMath::Power(mcPi0pT,3.65));
1938 }
1939 return weight;
1940}
fb87d707 1941
93c189c5 1942double AliAnalysisTaskHFECal::GetMCweightEta(double mcEtapT)
1943{
1944 double weight = 1.0;
1945
1946 weight = 223.3/TMath::Power((TMath::Exp(-0.17*mcEtapT-0.0322*mcEtapT*mcEtapT)+mcEtapT/1.69),5.65);
1947 return weight;
1948}
1949
1950
fb87d707 1951//_________________________________________
1952void AliAnalysisTaskHFECal::FindTriggerClusters()
1953{
bd6ee6dd 1954 //cout << "finding trigger patch" << endl;
fb87d707 1955 // constants
1956 const int nModuleCols = 2;
1957 const int nModuleRows = 5;
1958 const int nColsFeeModule = 48;
1959 const int nRowsFeeModule = 24;
1960 const int nColsFaltroModule = 24;
1961 const int nRowsFaltroModule = 12;
1962 //const int faltroWidthMax = 20;
1963
1964 // part 1, trigger extraction -------------------------------------
1965 Int_t globCol, globRow;
1966 //Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0, trigInCut=0;
1967 Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0;
1968
1969 //Int_t trigtimes[faltroWidthMax];
1970 Double_t cellTime[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1971 Double_t cellEnergy[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1972 //Double_t fTrigCutLow = 6;
1973 //Double_t fTrigCutHigh = 10;
1974 Double_t fTimeCutLow = 469e-09;
1975 Double_t fTimeCutHigh = 715e-09;
1976
1977 AliESDCaloTrigger * fCaloTrigger = fESD->GetCaloTrigger( "EMCAL" );
1978
1979 // erase trigger maps
1980 for(Int_t i = 0; i < nColsFaltroModule*nModuleCols; i++ )
1981 {
1982 for(Int_t j = 0; j < nRowsFaltroModule*nModuleRows; j++ )
1983 {
1984 ftriggersCut[i][j] = 0;
1985 ftriggers[i][j] = 0;
1986 ftriggersTime[i][j] = 0;
1987 }
1988 }
1989
1990 Int_t iglobCol=0, iglobRow=0;
1991 // go through triggers
1992 if( fCaloTrigger->GetEntries() > 0 )
1993 {
1994 // needs reset
1995 fCaloTrigger->Reset();
1996 while( fCaloTrigger->Next() )
1997 {
1998 fCaloTrigger->GetPosition( globCol, globRow );
1999 fCaloTrigger->GetNL0Times( ntimes );
2000 /*
2001 // no L0s
2002 if( ntimes < 1 ) continue;
2003 // get precise timings
2004 fCaloTrigger->GetL0Times( trigtimes );
2005 trigInCut = 0;
2006 for(Int_t i = 0; i < ntimes; i++ )
2007 {
2008 // save the first trigger time in channel
2009 if( i == 0 || triggersTime[globCol][globRow] > trigtimes[i] )
2010 triggersTime[globCol][globRow] = trigtimes[i];
2011 //printf("trigger times: %d\n",trigtimes[i]);
2012 // check if in cut
2013 if(trigtimes[i] > fTrigCutLow && trigtimes[i] < fTrigCutHigh )
2014 trigInCut = 1;
2015
2016 fTrigTimes->Fill(trigtimes[i]);
2017 }
2018 */
2019
2020 //L1 analysis from AliAnalysisTaskEMCALTriggerQA
2021 Int_t bit = 0;
2022 fCaloTrigger->GetTriggerBits(bit);
bd6ee6dd 2023 //cout << "bit = " << bit << endl;
fb87d707 2024
2025 Int_t ts = 0;
2026 fCaloTrigger->GetL1TimeSum(ts);
bd6ee6dd 2027 //cout << "ts = " << ts << endl;
fb87d707 2028 if (ts > 0)ftriggers[globCol][globRow] = 1;
2029 // number of triggered channels in event
2030 nTrigChannel++;
2031 // ... inside cut
2032 if(ts>0 && (bit >> 6 & 0x1))
2033 {
2034 iglobCol = globCol;
2035 iglobRow = globRow;
2036 nTrigChannelCut++;
bd6ee6dd 2037 //cout << "ts cut = " << ts << endl;
2038 //cout << "globCol = " << globCol << endl;
2039 //cout << "globRow = " << globRow << endl;
fb87d707 2040 ftriggersCut[globCol][globRow] = 1;
2041 }
2042
2043 } // calo trigger entries
2044 } // has calo trigger entries
2045
2046 // part 2 go through the clusters here -----------------------------------
bd6ee6dd 2047 //cout << " part 2 go through the clusters here ----------------------------------- " << endl;
fb87d707 2048 Int_t nCluster=0, nCell=0, iCell=0, gCell=0;
60d77596 2049 Short_t cellAddr, nSACell;
2050 Int_t mclabel;
fb87d707 2051 //Int_t nSACell, iSACell, mclabel;
2052 Int_t iSACell;
42c75692 2053 Double_t cellAmp=0, cellTimeT=0, clusterTime=0, efrac=0;
fb87d707 2054 Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta, gphi, geta, feta, fphi;
2055
2056 TRefArray *fCaloClusters = new TRefArray();
2057 fESD->GetEMCALClusters( fCaloClusters );
2058 nCluster = fCaloClusters->GetEntries();
2059
2060
2061 // save all cells times if there are clusters
2062 if( nCluster > 0 ){
2063 // erase time map
2064 for(Int_t i = 0; i < nColsFeeModule*nModuleCols; i++ ){
2065 for(Int_t j = 0; j < nRowsFeeModule*nModuleRows; j++ ){
2066 cellTime[i][j] = 0.;
2067 cellEnergy[i][j] = 0.;
2068 }
2069 }
2070
2071 // get cells
2072 AliESDCaloCells *fCaloCells = fESD->GetEMCALCells();
2073 //AliVCaloCells fCaloCells = fESD->GetEMCALCells();
2074 nSACell = fCaloCells->GetNumberOfCells();
2075 for(iSACell = 0; iSACell < nSACell; iSACell++ ){
2076 // get the cell info *fCal
2077 fCaloCells->GetCell( iSACell, cellAddr, cellAmp, cellTimeT , mclabel, efrac);
2078
2079 // get cell position
2080 fGeom->GetCellIndex( cellAddr, nSupMod, nModule, nIphi, nIeta );
2081 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
2082
2083 // convert co global phi eta
2084 gphi = iphi + nRowsFeeModule*(nSupMod/2);
2085 geta = ieta + nColsFeeModule*(nSupMod%2);
2086
2087 // save cell time and energy
2088 cellTime[geta][gphi] = cellTimeT;
2089 cellEnergy[geta][gphi] = cellAmp;
2090
2091 }
2092 }
2093
2094 Int_t nClusterTrig, nClusterTrigCut;
2095 UShort_t *cellAddrs;
2096 Double_t clsE=-999, clsEta=-999, clsPhi=-999;
2097 Float_t clsPos[3] = {0.,0.,0.};
2098
2099 for(Int_t icl=0; icl<fESD->GetNumberOfCaloClusters(); icl++)
2100 {
2101 AliESDCaloCluster *cluster = fESD->GetCaloCluster(icl);
2102 if(!cluster || !cluster->IsEMCAL()) continue;
2103
2104 // get cluster cells
2105 nCell = cluster->GetNCells();
2106
2107 // get cluster energy
2108 clsE = cluster->E();
2109
2110 // get cluster position
2111 cluster->GetPosition(clsPos);
2112 TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
2113 clsEta = clsPosVec.Eta();
2114 clsPhi = clsPosVec.Phi();
2115
2116 // get the cell addresses
2117 cellAddrs = cluster->GetCellsAbsId();
2118
2119 // check if the cluster contains cell, that was marked as triggered
2120 nClusterTrig = 0;
2121 nClusterTrigCut = 0;
2122
2123 // loop the cells to check, if cluser in acceptance
2124 // any cluster with a cell outside acceptance is not considered
2125 for( iCell = 0; iCell < nCell; iCell++ )
2126 {
42c75692 2127 // check hot cell
feffe705 2128 //if(clsE>6.0)fCellCheck->Fill(clsE,cellAddrs[iCell]);
42c75692 2129
fb87d707 2130 // get cell position
2131 fGeom->GetCellIndex( cellAddrs[iCell], nSupMod, nModule, nIphi, nIeta );
2132 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
2133
2134 // convert co global phi eta
2135 gphi = iphi + nRowsFeeModule*(nSupMod/2);
2136 geta = ieta + nColsFeeModule*(nSupMod%2);
2137
2138 if( cellTime[geta][gphi] > 0. ){
2139 clusterTime += cellTime[geta][gphi];
2140 gCell++;
2141 }
2142
2143 // get corresponding FALTRO
2144 fphi = gphi / 2;
2145 feta = geta / 2;
2146
bd6ee6dd 2147 //cout << "fphi = " << fphi << endl;
2148 //cout << "feta = " << feta << endl;
2149
fb87d707 2150 // try to match with a triggered
2151 if( ftriggers[feta][fphi]==1)
2152 { nClusterTrig++;
2153 }
2154 if( ftriggersCut[feta][fphi]==1)
2155 { nClusterTrigCut++;
2156 }
2157
bd6ee6dd 2158 //cout << "nClusterTrigCut : " << nClusterTrigCut << endl;
2159
fb87d707 2160 } // cells
2161
2162
2163 if( gCell > 0 )
2164 clusterTime = clusterTime / (Double_t)gCell;
2165 // fix the reconstriction code time 100ns jumps
2166 if( fESD->GetBunchCrossNumber() % 4 < 2 )
2167 clusterTime -= 0.0000001;
2168
feffe705 2169 //fClsETime->Fill(clsE,clusterTime);
2170 //fClsEBftTrigCut->Fill(clsE);
fb87d707 2171
2172 if(nClusterTrig>0){
feffe705 2173 //fClsETime1->Fill(clsE,clusterTime);
fb87d707 2174 }
2175
2176 if(nClusterTrig>0){
2177 cluster->SetChi2(1);
feffe705 2178 //fClsEAftTrigCut1->Fill(clsE);
fb87d707 2179 }
2180
2181 if(nClusterTrigCut>0){
2182 cluster->SetChi2(2);
feffe705 2183 //fClsEAftTrigCut2->Fill(clsE);
fb87d707 2184 }
2185
2186 if(nClusterTrigCut>0 && ( nCell > (1 + clsE / 3)))
2187 {
2188 cluster->SetChi2(3);
feffe705 2189 //fClsEAftTrigCut3->Fill(clsE);
fb87d707 2190 }
2191
2192 if(nClusterTrigCut>0 && (nCell > (1 + clsE / 3) )&&( clusterTime > fTimeCutLow && clusterTime < fTimeCutHigh ))
2193 {
2194 // cluster->SetChi2(4);
feffe705 2195 //fClsEAftTrigCut4->Fill(clsE);
fb87d707 2196 }
2197 if(nClusterTrigCut<1)
2198 {
2199 cluster->SetChi2(0);
2200
feffe705 2201 //fClsEAftTrigCut->Fill(clsE);
fb87d707 2202 }
2203
2204 } // clusters
2205}
2206
09022877 2207// <-------- only MC correction
8efacc22 2208double AliAnalysisTaskHFECal::MCEopMeanCorrection(double pTmc, float central)
8af58b0d 2209{
8efacc22 2210 TF1 *fcorr0 = new TF1("fcorr0","[0]*tanh([1]+[2]*x)");
2211 TF1 *fcorr1 = new TF1("fcorr1","[0]*tanh([1]+[2]*x)");
2212
8af58b0d 2213 double shift = 0.0;
8efacc22 2214
2215 if(central>0 && central<=10)
2216 {
2217 fcorr0->SetParameters(1.045,1.288,3.18e-01); //
2218 fcorr1->SetParameters(9.91e-01,3.466,2.344);
2219 }
09022877 2220 else if(central>10 && central<=20)
8efacc22 2221 {
2222 fcorr0->SetParameters(1.029,8.254e-01,4.07e-01);
2223 fcorr1->SetParameters(0.975,2.276,1.501e-01);
2224 }
2225 else if(central>20 && central<=30)
2226 {
2227 fcorr0->SetParameters(1.01,8.795e-01,3.904e-01);
2228 fcorr1->SetParameters(9.675e-01,1.654,2.583e-01);
2229 }
2230 else if(central>30 && central<=40)
2231 {
2232 fcorr0->SetParameters(1.00,1.466,2.305e-1);
2233 fcorr1->SetParameters(9.637e-01,1.473,2.754e-01);
2234 }
2235 else if(central>40 && central<=50)
2236 {
fd2126aa 2237 fcorr0->SetParameters(1.00,1.422,1.518e-01);
8efacc22 2238 fcorr1->SetParameters(9.59e-01,1.421,2.931e-01);
2239 }
2240
2241 else if(central>50 && central<=70)
2242 {
2243 fcorr0->SetParameters(0.989,2.495,2.167);
e1f0fb74 2244 fcorr1->SetParameters(0.961,1.734,1.438e-01);
8efacc22 2245 }
2246 else if(central>70 && central<=100)
2247 {
2248 fcorr0->SetParameters(0.981,-3.373,3.93327);
fd2126aa 2249 fcorr1->SetParameters(9.574e-01,1.698,1.58e-01);
8efacc22 2250 }
2251
2252
2253 shift = fcorr0->Eval(pTmc)-fcorr1->Eval(pTmc);
2254
89f41a30 2255 fcorr0->Delete();
2256 fcorr1->Delete();
2257
8af58b0d 2258 return shift;
2259}
fb87d707 2260
09022877 2261// <-------- only Data correction
2262double AliAnalysisTaskHFECal::NsigmaCorrection(double tmpeta, float central)
2263{
2264 int icent = 0;
2265
2266 if(central>=0 && central<10)
2267 {
2268 icent = 0;
2269 }
2270 else if(central>=10 && central<20)
2271 {
2272 icent = 1;
2273 }
2274 else if(central>=20 && central<30)
2275 {
2276 icent = 2;
2277 }
2278 else if(central>=30 && central<40)
2279 {
2280 icent = 3;
2281 }
2282 else if(central>=40 && central<50)
2283 {
2284 icent = 4;
2285 }
2286 else if(central>=50 && central<70)
2287 {
2288 icent = 5;
2289 }
2290 else
2291 {
2292 icent = 6;
2293 }
2294
2295 double shift = fnSigEtaCorr[icent]->Eval(tmpeta);
2296
2297 //cout << "eta correction"<< endl;
2298 //cout << "cent = "<< central<< endl;
2299 //cout << "icent = "<< icent << endl;
2300 //cout << "shift = "<< shift << endl;
2301
2302 return shift;
fb87d707 2303
09022877 2304}
fb87d707 2305