]> 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
Merge branch 'feature-movesplit'
[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
3b68b3fb 534
535 Bool_t SelColl = kTRUE;
536 //cout <<"check trigger : " << GetCollisionCandidates() << endl;
537 //cout <<"check kAny : " << AliVEvent::kAny<< endl;
538 if(GetCollisionCandidates()==AliVEvent::kAny)
539 {
540 //cout <<"kAny selection"<< endl;
541 SelColl = kFALSE;
542 TString firedTrigger;
543 firedTrigger = fESD->GetFiredTriggerClasses();
544 if(firedTrigger.Contains("CVLN_B2-B-NOPF-ALLNOTRD") || firedTrigger.Contains("CVLN_R1-B-NOPF-ALLNOTRD") || firedTrigger.Contains("CSEMI_R1-B-NOPF-ALLNOTRD"))SelColl=kTRUE;
545
546 //cout << "SemiCentral ? " << SelColl << endl;
547 }
548
549 if(!SelColl)return;
550
bfc7c23b 551 fNoEvents->Fill(0);
552
553 Int_t fNOtrks = fESD->GetNumberOfTracks();
554 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
555
556 Double_t pVtxZ = -999;
557 pVtxZ = pVtx->GetZ();
558
559 if(TMath::Abs(pVtxZ)>10) return;
560 fNoEvents->Fill(1);
561
f4e0d2d5 562 if(fNOtrks<1) return;
bfc7c23b 563 fNoEvents->Fill(2);
564
565 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
566 if(!pidResponse){
567 AliDebug(1, "Using default PID Response");
568 pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
569 }
570
571 fPID->SetPIDResponse(pidResponse);
572
573 fCFM->SetRecEventInfo(fESD);
574
bfc7c23b 575 fCent->Fill(cent);
576
b12bc641 577 //if(cent>90.) return;
bfc7c23b 578
579 // Calorimeter info.
580
42c75692 581 FindTriggerClusters();
582
bfc7c23b 583 // make EMCAL array
6f4a9952 584 double maxE = 0.0;
bfc7c23b 585 for(Int_t iCluster=0; iCluster<fESD->GetNumberOfCaloClusters(); iCluster++)
586 {
587 AliESDCaloCluster *clust = fESD->GetCaloCluster(iCluster);
68d718e7 588 if(clust && clust->IsEMCAL())
bfc7c23b 589 {
590 double clustE = clust->E();
591 float emcx[3]; // cluster pos
592 clust->GetPosition(emcx);
593 TVector3 clustpos(emcx[0],emcx[1],emcx[2]);
594 double emcphi = clustpos.Phi();
595 double emceta = clustpos.Eta();
fb87d707 596 double calInfo[5];
597 calInfo[0] = emcphi; calInfo[1] = emceta; calInfo[2] = clustE; calInfo[3] = cent; calInfo[4] = clust->Chi2();
f4e0d2d5 598 //fEMCAccE->Fill(calInfo);
f4e0d2d5 599 hEMCAccE->Fill(cent,clustE);
6f4a9952 600 if(clustE>maxE)maxE = clustE;
bfc7c23b 601 }
602 }
603
604 // Track loop
605 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
606 AliESDtrack* track = fESD->GetTrack(iTracks);
607 if (!track) {
608 printf("ERROR: Could not receive track %d\n", iTracks);
609 continue;
610 }
e7d87aef 611
89f41a30 612 //--- Get MC informtion
613
46305ed6 614 int parentlabel = 99999;
615 int parentPID = 99999;
616 int grand_parentlabel = 99999;
617 int grand_parentPID = 99999;
e7d87aef 618 Bool_t mcPho = kFALSE;
619 Bool_t mcDtoE= kFALSE;
620 Bool_t mcBtoE= kFALSE;
93c189c5 621 Bool_t mcOrgPi0 = kFALSE;
622 Bool_t mcOrgEta = kFALSE;
e7d87aef 623 double mcele = -1.;
44be9e1d 624 double mcpT = 0.0;
b12bc641 625 double mcMompT = 0.0;
e4b0faf2 626 //double mcGrandMompT = 0.0;
627 double mcWeight = -10.0;
bd50bf24 628 double conv_proR = -1.0;
d113d7cd 629
630 int iHijing = 1;
70448e3b 631 int mcLabel = -1;
d113d7cd 632
e7d87aef 633 if(fmcData && fMC && stack)
634 {
feffe705 635 Int_t label = TMath::Abs(track->GetLabel());
7d495dc5 636 //mcLabel = track->GetLabel();
1dd85087 637 mcLabel = abs(track->GetLabel()); // check for conv. issue
bd6ee6dd 638
8efacc22 639 if(mcLabel>-1)
bd6ee6dd 640 {
641
642 Bool_t MChijing = fMC->IsFromBGEvent(label);
643 if(!MChijing)iHijing = 0;
644
645 TParticle* particle = stack->Particle(label);
646 int mcpid = particle->GetPdgCode();
647 mcpT = particle->Pt();
bd50bf24 648 conv_proR = particle->R();
bd6ee6dd 649 //printf("MCpid = %d",mcpid);
650 if(particle->GetFirstMother()>-1)
651 {
652 //int parentlabel = particle->GetFirstMother();
653 //int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
654 //mcMompT = stack->Particle(particle->GetFirstMother())->Pt();
655 FindMother(particle, parentlabel, parentPID);
656 mcMompT = stack->Particle(parentlabel)->Pt();
657 if((parentPID==22 || parentPID==111 || parentPID==221)&& fabs(mcpid)==11)mcPho = kTRUE;
658 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(mcpid)==11)mcDtoE = kTRUE;
659 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(mcpid)==11)mcBtoE = kTRUE;
660
661 // make D->e pT correlation
662 if(mcDtoE)fMomDtoE->Fill(mcpT,mcMompT);
663
664 //cout << "check PID = " << parentPID << endl;
665 //cout << "check pho = " << mcPho << endl;
666 //cout << "check D or B = " << mcDtoE << endl;
667 // pi->e (Dalitz)
668 if(parentPID==111 && fabs(mcpid)==11 && mcMompT>0.0)
669 {
670 //cout << "find pi0->e " << endl;
671 mcOrgPi0 = kTRUE;
672 mcWeight = GetMCweight(mcMompT);
673 }
674 // eta->e (Dalitz)
675 if(parentPID==221 && fabs(mcpid)==11 && mcMompT>0.0)
676 {
677 //cout << "find Eta->e " << endl;
678 mcOrgEta = kTRUE;
679 mcWeight = GetMCweightEta(mcMompT);
680 }
681
682 // access grand parent
683 TParticle* particle_parent = stack->Particle(parentlabel); // get parent pointer
684 //if(particle_parent->GetFirstMother()>-1 && parentPID==22 && fabs(mcpid)==11) // get grand parent g->e
685 if(particle_parent->GetFirstMother()>-1 && (parentPID==22 || parentPID==111) && fabs(mcpid)==11) // get grand parent g->e
686 {
687 //int grand_parentPID = stack->Particle(particle_parent->GetFirstMother())->GetPdgCode();
688 //double pTtmp = stack->Particle(particle_parent->GetFirstMother())->Pt();
689 FindMother(particle_parent, grand_parentlabel, grand_parentPID);
690 double mcGrandpT = stack->Particle(grand_parentlabel)->Pt();
691 if(grand_parentPID==111 && mcGrandpT>0.0)
692 {
693 // check eta->pi0 decay !
694 int grand2_parentlabel = 99999; int grand2_parentPID = 99999;
695 TParticle* particle_grand = stack->Particle(grand_parentlabel); // get parent pointer
696 FindMother(particle_grand, grand2_parentlabel, grand2_parentPID);
697 if(grand2_parentPID==221)
698 {
699 //cout << "find Eta->e " << endl;
700 double mcGrandpT2 = stack->Particle(grand2_parentlabel)->Pt();
701 mcOrgEta = kTRUE;
702 mcWeight = GetMCweight(mcGrandpT2);
703 mcMompT = mcGrandpT2;
704 }
705 else
706 {
707 //cout << "find pi0->e " << endl;
708 mcOrgPi0 = kTRUE;
709 mcWeight = GetMCweight(mcGrandpT);
710 mcMompT = mcGrandpT;
711 }
712 }
713
714 if(grand_parentPID==221 && mcGrandpT>0.0)
715 {
716 //cout << "find Eta->e " << endl;
717 mcOrgEta = kTRUE;
718 mcOrgPi0 = kFALSE;
719 mcWeight = GetMCweightEta(mcGrandpT);
720 mcMompT = mcGrandpT;
721 }
722 }
723 }
724
cc1605fe 725 //cout << "===================="<<endl;
726 //cout << "mcDtoE : " << mcDtoE << endl;
727 //cout << "mcBtoE : " << mcBtoE << endl;
728 //cout << "mcPho : " << mcPho << endl;
729
730 if(fabs(mcpid)==11)mcele= 0.;
731 //cout << "check e: " << mcele << endl;
bd6ee6dd 732 if(fabs(mcpid)==11 && mcDtoE)mcele= 1.;
cc1605fe 733 //cout << "check D->e: " << mcele << endl;
bd6ee6dd 734 if(fabs(mcpid)==11 && mcBtoE)mcele= 2.;
cc1605fe 735 //cout << "check B->e: " << mcele << endl;
bd6ee6dd 736 if(fabs(mcpid)==11 && mcPho)mcele= 3.;
cc1605fe 737 //cout << "check Pho->e: " << mcele << endl;
bd6ee6dd 738
09022877 739 //cout << "check PID " << endl;
fd2126aa 740 if(fabs(mcpid)!=11)
741 {
09022877 742 //cout << "!= 11" << endl;
743 //cout << mcpid << endl;
fd2126aa 744 }
745 if(mcele==-1)
746 {
09022877 747 //cout << "mcele==-1" << endl;
748 //cout << mcele << endl;
749 //cout << mcpid << endl;
fd2126aa 750 }
751
bd6ee6dd 752 } // end of mcLabel>-1
753
89f41a30 754 } // <------ end of MC info.
e7d87aef 755
2b4460a6 756 //cout << "Pi0 = " << mcOrgPi0 << " ; Eta = " << mcOrgEta << endl;
50919258 757 //printf("weight = %f\n",mcWeight);
758
9b3495ff 759 if(TMath::Abs(track->Eta())>0.6) continue;
fbe402d7 760 //if(TMath::Abs(track->Pt()<2.5)) continue;
761 if(TMath::Abs(track->Pt()<0.1)) continue;
bfc7c23b 762
05ac9660 763 int nITS = track->GetNcls(0);
764
bfc7c23b 765 fTrackPtBefTrkCuts->Fill(track->Pt());
f2111233 766
05ac9660 767 /*
f2fa5a00 768 UChar_t itsPixel = track->GetITSClusterMap();
05ac9660 769 cout << "nITS = " << nITS << endl;
770 if(itsPixel & BIT(0))cout << "1st layer hit" << endl;
771 if(itsPixel & BIT(1))cout << "2nd layer hit" << endl;
772 if(itsPixel & BIT(2))cout << "3rd layer hit" << endl;
773 if(itsPixel & BIT(3))cout << "4th layer hit" << endl;
774 if(itsPixel & BIT(4))cout << "5th layer hit" << endl;
775 */
776
bfc7c23b 777 // RecKine: ITSTPC cuts
778 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
779
780 //RecKink
781 if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
782 if(track->GetKinkIndex(0) != 0) continue;
783 }
784
785 // RecPrim
786 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
787
788 // HFEcuts: ITS layers cuts
6e1d7a8a 789 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
bfc7c23b 790
f09766dd 791 // HFE cuts: TPC PID cleanup
792 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
793
36e9774f 794 if(mcPho && iHijing==0)fPhoVertexReco_step0->Fill(track->Pt(),conv_proR); // check MC vertex
795 if(mcPho && iHijing==1)fPhoVertexReco_step1->Fill(track->Pt(),conv_proR); // check MC vertex
796
a6df418c 797 int nTPCcl = track->GetTPCNcls();
96167a04 798 //int nTPCclF = track->GetTPCNclsF(); // warnings
05ac9660 799 //int nITS = track->GetNcls(0);
800
05ac9660 801
bfc7c23b 802 fTrackPtAftTrkCuts->Fill(track->Pt());
803
804 Double_t mom = -999., eop=-999., pt = -999., dEdx=-999., fTPCnSigma=-10, phi=-999., eta=-999.;
feffe705 805 pt = track->Pt();
d26226a1 806 if(pt<2.5)continue;
fbe402d7 807 if(pt<0.1)continue;
bfc7c23b 808
6f4a9952 809 //Int_t charge = track->Charge();
bfc7c23b 810 fTrkpt->Fill(pt);
811 mom = track->P();
812 phi = track->Phi();
813 eta = track->Eta();
ae88597c 814 float dca_xy;
815 float dca_z;
816 track->GetImpactParameters(dca_xy,dca_z);
817
bfc7c23b 818 dEdx = track->GetTPCsignal();
819 fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
8af58b0d 820
09022877 821 //cout << "nSigma correctoon-----" << endl;
822 //cout << "org = " << fTPCnSigma << endl;
9eae9ea5 823 /*
09022877 824 if(!fmcData) // nsigma eta correction
825 {
826 double nSigexpCorr = NsigmaCorrection(eta,cent);
827 fTPCnSigma -= nSigexpCorr;
828 }
9eae9ea5 829 */
09022877 830 //cout << "correction = " << fTPCnSigma << endl;
831
89f41a30 832 //--- track cluster match
833
bfc7c23b 834 double ncells = -1.0;
835 double m20 = -1.0;
836 double m02 = -1.0;
837 double disp = -1.0;
838 double rmatch = -1.0;
839 double nmatch = -1.0;
89f41a30 840 //double oppstatus = 0.0;
59d998de 841 double emctof = 0.0;
6f4a9952 842 Bool_t MaxEmatch = kFALSE;
bfc7c23b 843
13055933 844 Int_t clsId = track->GetEMCALcluster();
ebf2f603 845
846 //cout << "match ID" << clsId << " ; " << clsId_Tender << endl;
847
848 if (clsId>=0){
bfc7c23b 849 AliESDCaloCluster *clust = fESD->GetCaloCluster(clsId);
ebf2f603 850
bfc7c23b 851 if(clust && clust->IsEMCAL()){
852
853 double clustE = clust->E();
6f4a9952 854 if(clustE==maxE)MaxEmatch = kTRUE;
bfc7c23b 855 eop = clustE/fabs(mom);
8efacc22 856
c5ec1e9b 857
bfc7c23b 858 //double clustT = clust->GetTOF();
859 ncells = clust->GetNCells();
860 m02 = clust->GetM02();
861 m20 = clust->GetM20();
862 disp = clust->GetDispersion();
863 double delphi = clust->GetTrackDx();
864 double deleta = clust->GetTrackDz();
865 rmatch = sqrt(pow(delphi,2)+pow(deleta,2));
866 nmatch = clust->GetNTracksMatched();
59d998de 867 emctof = clust->GetTOF();
868 //cout << "emctof = " << emctof << endl;
ebf2f603 869 //cout << "eop org = "<< eop << endl;
89f41a30 870 double eoporg = eop;
871 if(fmcData)
9eae9ea5 872 {
873 double mceopcorr = MCEopMeanCorrection(pt,cent);
874 eop += mceopcorr;
875 }
ebf2f603 876 //cout << "eop corr = " << eop << endl;
9eae9ea5 877
feffe705 878 double valdedx[16];
a6df418c 879 valdedx[0] = pt; valdedx[1] = nITS; valdedx[2] = phi; valdedx[3] = eta; valdedx[4] = fTPCnSigma;
05ac9660 880 valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = dca_xy, valdedx[8] = dca_z; valdedx[9] = m20; valdedx[10] = mcpT;
89f41a30 881 valdedx[11] = cent; valdedx[12] = dEdx; valdedx[13] = eoporg; valdedx[14] = nTPCcl;
feffe705 882 valdedx[15] = mcele;
89f41a30 883 fEleInfo->Fill(valdedx);
bfc7c23b 884
885 }
886 }
e1f0fb74 887
888 //Get Cal info PID response
9eae9ea5 889 /*
e1f0fb74 890 double eop2;
891 double ss[4];
892 Double_t nSigmaEop = fPID->GetPIDResponse()->NumberOfSigmasEMCAL(track,AliPID::kElectron,eop2,ss);
893 if(fTPCnSigma>-1.5 && fTPCnSigma<3.0 && nITS>2.5 && nTPCcl>100)
894 {
895 double valEop[3];
896 valEop[0] = cent;
897 valEop[1] = pt;
898 valEop[2] = nSigmaEop;
899 fElenSigma->Fill(valEop);
900 }
9eae9ea5 901 */
89f41a30 902 // --- tarck cut & e ID
e1f0fb74 903
a6df418c 904 if(nITS<2.5)continue;
905 if(nTPCcl<100)continue;
7d495dc5 906
7d495dc5 907
a6df418c 908 CheckNclust->Fill(nTPCcl);
909 CheckNits->Fill(nITS);
ae88597c 910 CheckDCA->Fill(dca_xy,dca_z);
f2111233 911 // check production vertex of photons
f2111233 912
a6df418c 913
42c75692 914 fdEdxBef->Fill(mom,fTPCnSigma);
bfc7c23b 915 fTPCnsigma->Fill(mom,fTPCnSigma);
f09766dd 916 if(fTPCnSigma >= -1.0 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,eop);
bfc7c23b 917
89f41a30 918 //--- track accepted by HFE
bfc7c23b 919
89f41a30 920 Int_t pidpassed = 1;
bfc7c23b 921 AliHFEpidObject hfetrack;
922 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
923 hfetrack.SetRecTrack(track);
3db00c72 924 double binct = 10.5;
925 if((0.0< cent) && (cent<5.0)) binct = 0.5;
926 if((5.0< cent) && (cent<10.0)) binct = 1.5;
927 if((10.0< cent) && (cent<20.0)) binct = 2.5;
928 if((20.0< cent) && (cent<30.0)) binct = 3.5;
929 if((30.0< cent) && (cent<40.0)) binct = 4.5;
930 if((40.0< cent) && (cent<50.0)) binct = 5.5;
931 if((50.0< cent) && (cent<60.0)) binct = 6.5;
932 if((60.0< cent) && (cent<70.0)) binct = 7.5;
933 if((70.0< cent) && (cent<80.0)) binct = 8.5;
934 if((80.0< cent) && (cent<90.0)) binct = 9.5;
935 if((90.0< cent) && (cent<100.0)) binct = 10.5;
936
937 hfetrack.SetCentrality((int)binct); //added
bfc7c23b 938 hfetrack.SetPbPb();
939 if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
940
f2fa5a00 941 if(pidpassed==0) continue; // nSigma rejection
89f41a30 942
943 //--- photonic ID
944
0f566acb 945 Bool_t fFlagPhotonicTPC = kFALSE;
946 Bool_t fFlagConvinatTPC = kFALSE;
947 SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicTPC,fFlagConvinatTPC,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta,0);
89f41a30 948
949 //--- check reco eff. with TPC
950 double phoval[5];
951 phoval[0] = cent;
952 phoval[1] = pt;
953 phoval[2] = fTPCnSigma;
954 phoval[3] = iHijing;
955 phoval[4] = mcMompT;
956
2f8ef315 957 if((fTPCnSigma >= -5.0 && fTPCnSigma <= 5) && (mcOrgPi0 || mcOrgEta))
89f41a30 958 {
959 if(iHijing==1)mcWeight = 1.0;
2f8ef315 960 if(mcOrgPi0)
961 {
962 fIncpTMCpho_pi0e_TPC->Fill(phoval,mcWeight);
0f566acb 963 if(fFlagPhotonicTPC) fPhoElecPtMC_pi0e_TPC->Fill(phoval,mcWeight);
964 if(fFlagConvinatTPC) fSameElecPtMC_pi0e_TPC->Fill(phoval,mcWeight);
466d4f7a 965 fPhoVertexReco_TPC->Fill(track->Pt(),conv_proR,mcWeight); // check MC vertex
966 if(fFlagPhotonicTPC)fPhoVertexReco_TPC_Invmass->Fill(track->Pt(),conv_proR,mcWeight); // check MC vertex
2f8ef315 967 }
968 if(mcOrgEta)
969 {
970 fIncpTMCpho_eta_TPC->Fill(phoval,mcWeight);
0f566acb 971 if(fFlagPhotonicTPC) fPhoElecPtMC_eta_TPC->Fill(phoval,mcWeight);
972 if(fFlagConvinatTPC) fSameElecPtMC_eta_TPC->Fill(phoval,mcWeight);
2f8ef315 973 }
89f41a30 974 }
975
0f566acb 976 //--- matching check
f2fa5a00 977 double emcphimim = 1.396;
978 double emcphimax = 3.14;
979 if(phi>emcphimim && phi<emcphimax)
980 {
0f566acb 981 if(fFlagPhotonicTPC)fMatchV0_0->Fill(pt); // data
d26226a1 982 if(mcele>2.1)fMatchMC_0->Fill(pt,mcWeight); // MC
983
984 if(eop>=0.0) // have a match
f2fa5a00 985 {
0f566acb 986 if(fFlagPhotonicTPC)fMatchV0_1->Fill(pt);
d26226a1 987 if(mcele>2.1)fMatchMC_1->Fill(pt,mcWeight);
f2fa5a00 988 }
989 }
990
c5a4ddbf 991 // check fake rejection
992 if(mcOrgPi0 || mcOrgEta)
993 {
1c9a06c8 994 double phiacc0 = 80.0/180.0*acos(-1);
995 double phiacc1 = 180.0/180.0*acos(-1);
c5a4ddbf 996 int TrStat = 0;
1c9a06c8 997 if(phi>phiacc0 && phi<phiacc1)
998 {
999 if(track->GetLabel()>0)TrStat = 1;
1000 fFakeRejection0->Fill(TrStat,pt,mcWeight);
1001 if(eop>-1.0)fFakeRejection1->Fill(TrStat,pt,mcWeight); // have match
1002 if(eop>0.9 && eop<1.3)fFakeRejection2->Fill(TrStat,pt,mcWeight); // have PID
9e88412a 1003
1004 if(TrStat==0)
1005 {
32abf0b9 1006 EopFake->Fill(pt,eop,mcWeight);
1007 MatchFake->Fill(pt,rmatch,mcWeight);
9e88412a 1008 }
1009 else
1010 {
32abf0b9 1011 EopTrue->Fill(pt,eop,mcWeight);
1012 MatchTrue->Fill(pt,rmatch,mcWeight);
9e88412a 1013 }
1014
1c9a06c8 1015 }
c5a4ddbf 1016 }
1017
89f41a30 1018 //+++++++ E/p cut ++++++++++++++++
1019
1020 if(eop<0.9 || eop>1.3)continue;
0f566acb 1021
1022 Bool_t fFlagPhotonicElec = kFALSE;
1023 Bool_t fFlagConvinatElec = kFALSE;
1024 SelectPhotonicElectron(iTracks,cent,track,fFlagPhotonicElec,fFlagConvinatElec,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing,mcOrgPi0,mcOrgEta,1);
1025
bfc7c23b 1026 fTrkEovPAft->Fill(pt,eop);
42c75692 1027 fdEdxAft->Fill(mom,fTPCnSigma);
bfc7c23b 1028
e7d87aef 1029 // Fill real data
fb87d707 1030 fIncpT->Fill(cent,pt);
bfc7c23b 1031 if(fFlagPhotonicElec) fPhoElecPt->Fill(cent,pt);
f09766dd 1032 if(fFlagConvinatElec) fSameElecPt->Fill(cent,pt);
fb87d707 1033
1034 if(m20>0.0 && m20<0.3)
1035 {
59d998de 1036 fIncpTM20->Fill(cent,pt);
1037 ftimingEle->Fill(pt,emctof);
fb87d707 1038 if(fFlagPhotonicElec) fPhoElecPtM20->Fill(cent,pt);
1039 if(fFlagConvinatElec) fSameElecPtM20->Fill(cent,pt);
ebef1da9 1040 if(!fFlagPhotonicElec) fSemiElecPtM20->Fill(cent,pt);
fb87d707 1041 }
feffe705 1042
5e38d973 1043
e7d87aef 1044 // MC
70448e3b 1045 // check label for electron candidiates
70448e3b 1046 int idlabel = 1;
1047 if(mcLabel==0)idlabel = 0;
1048 fLabelCheck->Fill(pt,idlabel);
1049 if(mcLabel==0)fgeoFake->Fill(phi,eta);
1050
d19af75d 1051 if(mcLabel<0 && m20>0.0 && m20<0.3 && fTPCnSigma>-1 && fTPCnSigma<3)
1052 {
1053 fFakeTrk0->Fill(cent,pt);
1054 }
1055
1056 if(mcele>-1) // select MC electrons
e7d87aef 1057 {
3db00c72 1058
1059 fIncpTMChfeAll->Fill(cent,pt);
1060 if(m20>0.0 && m20<0.3)fIncpTMCM20hfeAll->Fill(cent,pt);
d19af75d 1061 if(m20>0.0 && m20<0.3 && fTPCnSigma>-1 && fTPCnSigma<3)fFakeTrk1->Fill(cent,pt);
3db00c72 1062
e4b0faf2 1063 if(mcBtoE || mcDtoE) // select B->e & D->e
e7d87aef 1064 {
1065 fIncpTMChfe->Fill(cent,pt);
bd6ee6dd 1066 if(m20>0.0 && m20<0.3)
1067 {
09022877 1068 //cout << "MC label = " << mcLabel << endl;
bd6ee6dd 1069 fIncpTMCM20hfe->Fill(cent,pt);
1070 fIncpTMCM20hfeCheck->Fill(cent,mcpT);
1071 fIncpTMCM20hfeCheck_weight->Fill(phoval);
1072 }
e7d87aef 1073 }
e4b0faf2 1074
1075 if(mcPho) // select photonic electrons
e7d87aef 1076 {
44be9e1d 1077
1078 fIncpTMCpho->Fill(phoval);
1079 if(fFlagPhotonicElec) fPhoElecPtMC->Fill(phoval);
1080 if(fFlagConvinatElec) fSameElecPtMC->Fill(phoval);
e7d87aef 1081
1082 if(m20>0.0 && m20<0.3)
1083 {
44be9e1d 1084 fIncpTMCM20pho->Fill(phoval);
1085 if(fFlagPhotonicElec) fPhoElecPtMCM20->Fill(phoval);
1086 if(fFlagConvinatElec) fSameElecPtMCM20->Fill(phoval);
e4b0faf2 1087 // pi0->g->e
1088 if(mcWeight>-1)
1089 {
987053ce 1090 if(iHijing==1)mcWeight = 1.0;
93c189c5 1091 if(mcOrgPi0)
1092 {
697ecf6b 1093 fIncpTMCM20pho_pi0e->Fill(phoval,mcWeight);
1094 if(fFlagPhotonicElec) fPhoElecPtMCM20_pi0e->Fill(phoval,mcWeight);
1095 if(fFlagConvinatElec) fSameElecPtMCM20_pi0e->Fill(phoval,mcWeight);
05ac9660 1096
1097 // check production vertex
7d495dc5 1098 //fPhoVertexReco_EMCal->Fill(track->Pt(),conv_proR);
1099 //if(fFlagPhotonicElec)fPhoVertexReco_Invmass->Fill(track->Pt(),conv_proR);
1100 fPhoVertexReco_EMCal->Fill(track->Pt(),conv_proR,mcWeight);
1101 if(fFlagPhotonicElec)fPhoVertexReco_Invmass->Fill(track->Pt(),conv_proR,mcWeight);
05ac9660 1102
1103 }
bd50bf24 1104 // --- eta
93c189c5 1105 if(mcOrgEta)
1106 {
697ecf6b 1107 fIncpTMCM20pho_eta->Fill(phoval,mcWeight);
1108 if(fFlagPhotonicElec) fPhoElecPtMCM20_eta->Fill(phoval,mcWeight);
1109 if(fFlagConvinatElec) fSameElecPtMCM20_eta->Fill(phoval,mcWeight);
93c189c5 1110 }
bd50bf24 1111 // check production vertex
e4b0faf2 1112 }
e7d87aef 1113 }
1114 }
1115 }
bfc7c23b 1116 }
1117 PostData(1, fOutputList);
1118}
1119//_________________________________________
1120void AliAnalysisTaskHFECal::UserCreateOutputObjects()
1121{
42c75692 1122 //--- Check MC
1123
e7d87aef 1124 //Bool_t mcData = kFALSE;
42c75692 1125 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
1126 {
e7d87aef 1127 fmcData = kTRUE;
42c75692 1128 printf("+++++ MC Data available");
1129 }
e7d87aef 1130 if(fmcData)
42c75692 1131 {
1132 printf("++++++++= MC analysis \n");
1133 }
1134 else
1135 {
1136 printf("++++++++ real data analysis \n");
1137 }
1138
85985bb0 1139 printf("+++++++ QA hist %d",fqahist);
1140
fb87d707 1141 //---- Geometry
1142 fGeom = AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
1143
bfc7c23b 1144 //--------Initialize PID
42c75692 1145 //fPID->SetHasMCData(kFALSE);
e7d87aef 1146 fPID->SetHasMCData(fmcData);
bfc7c23b 1147 if(!fPID->GetNumberOfPIDdetectors())
1148 {
1149 fPID->AddDetector("TPC", 0);
89f41a30 1150 //fPID->AddDetector("EMCAL", 1); <--- apply PID selection in task
bfc7c23b 1151 }
1152
3db00c72 1153 Double_t params[4];
78ff1095 1154 const char *cutmodel;
3db00c72 1155 cutmodel = "pol0";
1156 params[0] = -1.0; //sigma min
9b3495ff 1157 double maxnSig = 3.0;
1158 if(fmcData)
1159 {
1160 params[0] = -5.0; //sigma min
1161 maxnSig = 5.0;
1162 }
1163
1164 for(Int_t a=0;a<11;a++)fPID->ConfigureTPCcentralityCut(a,cutmodel,params,maxnSig);
3db00c72 1165
bfc7c23b 1166 fPID->SortDetectors();
1167 fPIDqa = new AliHFEpidQAmanager();
1168 fPIDqa->Initialize(fPID);
a6df418c 1169
1170 //------- fcut --------------
1171 fCuts = new AliHFEcuts();
1172 fCuts->CreateStandardCuts();
1173 //fCuts->SetMinNClustersTPC(100);
1174 fCuts->SetMinNClustersTPC(90);
1175 fCuts->SetMinRatioTPCclusters(0.6);
1176 fCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
1177 //fCuts->SetMinNClustersITS(3);
1178 fCuts->SetMinNClustersITS(2);
05ac9660 1179 fCuts->SetProductionVertex(0,50,0,50);
4d84ab70 1180 //fCuts->SetCutITSpixel(AliHFEextraCuts::kAny); // default
1181 fCuts->SetCutITSpixel(AliHFEextraCuts::kBoth);
a6df418c 1182 fCuts->SetCheckITSLayerStatus(kFALSE);
1183 fCuts->SetVertexRange(10.);
1184 fCuts->SetTOFPIDStep(kFALSE);
fbe402d7 1185 //fCuts->SetPtRange(2, 50);
1186 fCuts->SetPtRange(0.1, 50);
0f566acb 1187 //fCuts->SetMaxImpactParam(3.,3.);
1188 fCuts->SetMaxImpactParam(2.4,3.2); // standard in 2011
a6df418c 1189
bfc7c23b 1190 //--------Initialize correction Framework and Cuts
1191 fCFM = new AliCFManager;
1192 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
1193 fCFM->SetNStepParticle(kNcutSteps);
1194 for(Int_t istep = 0; istep < kNcutSteps; istep++)
1195 fCFM->SetParticleCutsList(istep, NULL);
1196
1197 if(!fCuts){
1198 AliWarning("Cuts not available. Default cuts will be used");
1199 fCuts = new AliHFEcuts;
1200 fCuts->CreateStandardCuts();
1201 }
1202 fCuts->Initialize(fCFM);
1203
1204 //---------Output Tlist
1205 fOutputList = new TList();
1206 fOutputList->SetOwner();
1207 fOutputList->Add(fPIDqa->MakeList("PIDQA"));
1208
1209 fNoEvents = new TH1F("fNoEvents","",4,-0.5,3.5) ;
1210 fOutputList->Add(fNoEvents);
1211
02cd82a0 1212 Int_t binsE[5] = {250, 100, 1000, 200, 10};
fb87d707 1213 Double_t xminE[5] = {1.0, -1, 0.0, 0, -0.5};
1214 Double_t xmaxE[5] = {3.5, 1, 100.0, 100, 9.5};
1215 fEMCAccE = new THnSparseD("fEMCAccE","EMC acceptance & E;#eta;#phi;Energy;Centrality;trugCondition;",5,binsE,xminE,xmaxE);
c82f1b9d 1216 //if(fqahist==1)fOutputList->Add(fEMCAccE);
bfc7c23b 1217
f4e0d2d5 1218 hEMCAccE = new TH2F("hEMCAccE","Cluster Energy",200,0,100,100,0,20);
1219 fOutputList->Add(hEMCAccE);
1220
bfc7c23b 1221 fTrkpt = new TH1F("fTrkpt","track pt",100,0,50);
1222 fOutputList->Add(fTrkpt);
1223
1224 fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50);
1225 fOutputList->Add(fTrackPtBefTrkCuts);
1226
1227 fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50);
1228 fOutputList->Add(fTrackPtAftTrkCuts);
1229
1230 fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10);
1231 fOutputList->Add(fTPCnsigma);
1232
1233 fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2);
1234 fOutputList->Add(fTrkEovPBef);
1235
1236 fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2);
1237 fOutputList->Add(fTrkEovPAft);
1238
42c75692 1239 fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,200,-10,10);
bfc7c23b 1240 fOutputList->Add(fdEdxBef);
1241
42c75692 1242 fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,200,-10,10);
bfc7c23b 1243 fOutputList->Add(fdEdxAft);
1244
02cd82a0 1245 fIncpT = new TH2F("fIncpT","HFE pid electro vs. centrality",200,0,100,100,0,50);
bfc7c23b 1246 fOutputList->Add(fIncpT);
1247
02cd82a0 1248 fIncpTM20 = new TH2F("fIncpTM20","HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
fb87d707 1249 fOutputList->Add(fIncpTM20);
4e01b68c 1250
5e38d973 1251 Int_t nBinspho[9] = { 10, 30, 600, 60, 50, 4, 40, 8, 30};
1252 Double_t minpho[9] = { 0., 0., -0.1, 40, 0, -0.5, 0,-1.5, 0};
1253 Double_t maxpho[9] = {100., 30., 0.5, 100, 1, 3.5, 2, 6.5, 30};
f09766dd 1254
b12bc641 1255 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 1256 if(fqahist==1)fOutputList->Add(fInvmassLS);
bfc7c23b 1257
b12bc641 1258 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 1259 if(fqahist==1)fOutputList->Add(fInvmassULS);
bfc7c23b 1260
8806ce6c 1261 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 1262 if(fqahist==1)fOutputList->Add(fInvmassLSmc);
8806ce6c 1263
1264 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 1265 if(fqahist==1)fOutputList->Add(fInvmassULSmc);
8806ce6c 1266
d89866fd 1267 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 1268 fInvmassLSreco->Sumw2();
1269 fOutputList->Add(fInvmassLSreco);
1270
d89866fd 1271 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 1272 fInvmassULSreco->Sumw2();
1273 fOutputList->Add(fInvmassULSreco);
1274
d89866fd 1275 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 1276 fInvmassLSmc0->Sumw2();
1277 fOutputList->Add(fInvmassLSmc0);
1278
d89866fd 1279 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 1280 fInvmassLSmc1->Sumw2();
1281 fOutputList->Add(fInvmassLSmc1);
1282
d89866fd 1283 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 1284 fInvmassLSmc2->Sumw2();
1285 fOutputList->Add(fInvmassLSmc2);
1286
d89866fd 1287 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 1288 fInvmassLSmc3->Sumw2();
1289 fOutputList->Add(fInvmassLSmc3);
1290
d89866fd 1291 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 1292 fInvmassULSmc0->Sumw2();
1293 fOutputList->Add(fInvmassULSmc0);
1294
d89866fd 1295 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 1296 fInvmassULSmc1->Sumw2();
1297 fOutputList->Add(fInvmassULSmc1);
1298
d89866fd 1299 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 1300 fInvmassULSmc2->Sumw2();
1301 fOutputList->Add(fInvmassULSmc2);
1302
d89866fd 1303 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 1304 fInvmassULSmc3->Sumw2();
93c189c5 1305 fOutputList->Add(fInvmassULSmc3);
1306
bfc7c23b 1307 fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
1308 fOutputList->Add(fOpeningAngleLS);
1309
1310 fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
1311 fOutputList->Add(fOpeningAngleULS);
1312
1313 fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50);
1314 fOutputList->Add(fPhotoElecPt);
1315
02cd82a0 1316 fPhoElecPt = new TH2F("fPhoElecPt", "Pho-inclusive electron pt",200,0,100,100,0,50);
bfc7c23b 1317 fOutputList->Add(fPhoElecPt);
1318
02cd82a0 1319 fPhoElecPtM20 = new TH2F("fPhoElecPtM20", "Pho-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 1320 fOutputList->Add(fPhoElecPtM20);
1321
0f566acb 1322 fPhoElecPtM20Mass = new TH2F("fPhoElecPtM20Mass", "Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1323 fPhoElecPtM20Mass->Sumw2();
1324 fOutputList->Add(fPhoElecPtM20Mass);
1325
02cd82a0 1326 fSameElecPt = new TH2F("fSameElecPt", "Same-inclusive electron pt",200,0,100,100,0,50);
f09766dd 1327 fOutputList->Add(fSameElecPt);
1328
02cd82a0 1329 fSameElecPtM20 = new TH2F("fSameElecPtM20", "Same-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 1330 fOutputList->Add(fSameElecPtM20);
1331
0f566acb 1332 fSameElecPtM20Mass = new TH2F("fSameElecPtM20Mass", "Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1333 fSameElecPtM20Mass->Sumw2();
1334 fOutputList->Add(fSameElecPtM20Mass);
1335
ebef1da9 1336 fSemiElecPtM20 = new TH2F("fSemiElecPtM20", "Semi-inclusive electron pt with M20",200,0,100,100,0,50);
1337 fOutputList->Add(fSemiElecPtM20);
1338
02cd82a0 1339 fCent = new TH1F("fCent","Centrality",200,0,100) ;
bfc7c23b 1340 fOutputList->Add(fCent);
1341
1342 // Make common binning
44be9e1d 1343 const Double_t kMinP = 0.;
d19af75d 1344 const Double_t kMaxP = 20.;
bfc7c23b 1345
05ac9660 1346 //+++ 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select)
1347 // 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)
1348 Int_t nBins[16] = { 100, 7, 60, 20, 90, 100, 25, 60, 60, 100, 40, 10, 250, 100, 100, 8};
1349 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};
1350 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};
1351 //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);
1352 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 1353 //if(fqahist==1)fOutputList->Add(fEleInfo);
bfc7c23b 1354
e1f0fb74 1355 // Make common binning
1356 Int_t nBinsEop[3] = { 10, 50, 100};
1357 Double_t minEop[3] = { 0, 0, -5};
1358 Double_t maxEop[3] = {100, 50, 5};
1359 fElenSigma= new THnSparseD("fElenSigma", "Electron nSigma; cent; pT [GeV/c]; nSigma", 3, nBinsEop, minEop, maxEop);
1360 fOutputList->Add(fElenSigma);
1361
1362
fb87d707 1363 //<--- Trigger info
feffe705 1364 /*
fb87d707 1365 fClsEBftTrigCut = new TH1F("fClsEBftTrigCut","cluster E before trigger selection",1000,0,100);
1366 fOutputList->Add(fClsEBftTrigCut);
1367
1368 fClsEAftTrigCut = new TH1F("fClsEAftTrigCut","cluster E if cls has 0 trigcut channel",1000,0,100);
1369 fOutputList->Add(fClsEAftTrigCut);
1370
1371 fClsEAftTrigCut1 = new TH1F("fClsEAftTrigCut1","cluster E if cls with trig channel",1000,0,100);
1372 fOutputList->Add(fClsEAftTrigCut1);
1373
1374 fClsEAftTrigCut2 = new TH1F("fClsEAftTrigCut2","cluster E if cls with trigcut channel",1000,0,100);
1375 fOutputList->Add(fClsEAftTrigCut2);
1376
1377 fClsEAftTrigCut3 = new TH1F("fClsEAftTrigCut3","cluster E if cls with trigcut channel + nCell>Ecorrect",1000,0,100);
1378 fOutputList->Add(fClsEAftTrigCut3);
1379
1380 fClsEAftTrigCut4 = new TH1F("fClsEAftTrigCut4","cluster E if cls with trigcut channel + nCell>Ecorrect + cls time cut",1000,0,100);
1381 fOutputList->Add(fClsEAftTrigCut4);
1382
1383 fClsETime = new TH2F("fClsETime", "Cls time vs E; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
1384 fOutputList->Add(fClsETime);
1385
1386 fClsETime1 = new TH2F("fClsETime1", "Cls time vs E if cls contains trigger channel; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
1387 fOutputList->Add(fClsETime1);
1388
1389 fTrigTimes = new TH1F("fTrigTimes", "Trigger time; time; N;",25,0,25);
1390 fOutputList->Add(fTrigTimes);
1391
42c75692 1392 fCellCheck = new TH2F("fCellCheck", "Cell vs E; E GeV; Cell ID",10,6,26,12000,0,12000);
1393 fOutputList->Add(fCellCheck);
feffe705 1394 */
e7d87aef 1395 //<---------- MC
1396
02cd82a0 1397 fInputHFEMC = new TH2F("fInputHFEMC","Input MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 1398 fOutputList->Add(fInputHFEMC);
1399
02cd82a0 1400 fInputAlle = new TH2F("fInputAlle","Input MC electro vs. centrality",200,0,100,100,0,50);
feffe705 1401 fOutputList->Add(fInputAlle);
1402
02cd82a0 1403 fIncpTMChfe = new TH2F("fIncpTMChfe","MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 1404 fOutputList->Add(fIncpTMChfe);
1405
02cd82a0 1406 fIncpTMChfeAll = new TH2F("fIncpTMChfeAll","MC Alle pid electro vs. centrality",200,0,100,100,0,50);
3db00c72 1407 fOutputList->Add(fIncpTMChfeAll);
1408
02cd82a0 1409 fIncpTMCM20hfe = new TH2F("fIncpTMCM20hfe","MC HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
e7d87aef 1410 fOutputList->Add(fIncpTMCM20hfe);
1411
02cd82a0 1412 fIncpTMCM20hfeAll = new TH2F("fIncpTMCM20hfeAll","MC Alle pid electro vs. centrality with M20",200,0,100,100,0,50);
3db00c72 1413 fOutputList->Add(fIncpTMCM20hfeAll);
1414
60544aea 1415 fIncpTMCM20hfeCheck = new TH2F("fIncpTMCM20hfeCheck","MC HFE pid electro vs. centrality with M20 Check",200,0,100,100,0,50);
1416 fOutputList->Add(fIncpTMCM20hfeCheck);
44be9e1d 1417
2b4460a6 1418 Int_t nBinspho2[5] = { 200, 100, 7, 3, 700};
987053ce 1419 Double_t minpho2[5] = { 0., 0., -2.5, -0.5, 0.};
2b4460a6 1420 Double_t maxpho2[5] = {100., 50., 4.5, 2.5, 70.};
bd6ee6dd 1421
1422 fInputHFEMC_weight = new THnSparseD("fInputHFEMC_weight", "MC HFE electron pt",5,nBinspho2,minpho2,maxpho2);
1423 fOutputList->Add(fInputHFEMC_weight);
1424
1425 fIncpTMCM20hfeCheck_weight = new THnSparseD("fIncpTMCM20hfeCheck_weight", "HFE electron pt with M20",5,nBinspho2,minpho2,maxpho2);
1426 fOutputList->Add(fIncpTMCM20hfeCheck_weight);
1427
987053ce 1428 fIncpTMCpho = new THnSparseD("fIncpTMCpho","MC Pho pid electro vs. centrality",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1429 fOutputList->Add(fIncpTMCpho);
1430
987053ce 1431 fIncpTMCM20pho = new THnSparseD("fIncpTMCM20pho","MC Pho pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1432 fOutputList->Add(fIncpTMCM20pho);
1433
987053ce 1434 fPhoElecPtMC = new THnSparseD("fPhoElecPtMC", "MC Pho-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1435 fOutputList->Add(fPhoElecPtMC);
1436
987053ce 1437 fPhoElecPtMCM20 = new THnSparseD("fPhoElecPtMCM20", "MC Pho-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
feffe705 1438 fOutputList->Add(fPhoElecPtMCM20);
e7d87aef 1439
0f566acb 1440 fPhoElecPtMCM20Mass = new TH2D("fPhoElecPtMCM20Mass", "MC Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1441 fOutputList->Add(fPhoElecPtMCM20Mass);
1442
987053ce 1443 fSameElecPtMC = new THnSparseD("fSameElecPtMC", "MC Same-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1444 fOutputList->Add(fSameElecPtMC);
1445
987053ce 1446 fSameElecPtMCM20 = new THnSparseD("fSameElecPtMCM20", "MC Same-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1447 fOutputList->Add(fSameElecPtMCM20);
1448
0f566acb 1449 fSameElecPtMCM20Mass = new TH2D("fSameElecPtMCM20Mass", "MC Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1450 fOutputList->Add(fSameElecPtMCM20Mass);
1451
987053ce 1452 fIncpTMCM20pho_pi0e = new THnSparseD("fIncpTMCM20pho_pi0e","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1453 fIncpTMCM20pho_pi0e->Sumw2();
e4b0faf2 1454 fOutputList->Add(fIncpTMCM20pho_pi0e);
1455
987053ce 1456 fPhoElecPtMCM20_pi0e = new THnSparseD("fPhoElecPtMCM20_pi0e", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1457 fPhoElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1458 fOutputList->Add(fPhoElecPtMCM20_pi0e);
1459
987053ce 1460 fSameElecPtMCM20_pi0e = new THnSparseD("fSameElecPtMCM20_pi0e", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1461 fSameElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1462 fOutputList->Add(fSameElecPtMCM20_pi0e);
697ecf6b 1463 //
93c189c5 1464 fIncpTMCM20pho_eta = new THnSparseD("fIncpTMCM20pho_eta","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1465 fIncpTMCM20pho_eta->Sumw2();
1466 fOutputList->Add(fIncpTMCM20pho_eta);
1467
1468 fPhoElecPtMCM20_eta = new THnSparseD("fPhoElecPtMCM20_eta", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1469 fPhoElecPtMCM20_eta->Sumw2();
1470 fOutputList->Add(fPhoElecPtMCM20_eta);
1471
1472 fSameElecPtMCM20_eta = new THnSparseD("fSameElecPtMCM20_eta", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1473 fSameElecPtMCM20_eta->Sumw2();
1474 fOutputList->Add(fSameElecPtMCM20_eta);
697ecf6b 1475 // ------------
1476 fIncpTMCpho_pi0e_TPC = new THnSparseD("fIncpTMCpho_pi0e_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1477 fIncpTMCpho_pi0e_TPC->Sumw2();
1478 fOutputList->Add(fIncpTMCpho_pi0e_TPC);
1479
1480 fPhoElecPtMC_pi0e_TPC = new THnSparseD("fPhoElecPtMC_pi0e_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2);
1481 fPhoElecPtMC_pi0e_TPC->Sumw2();
1482 fOutputList->Add(fPhoElecPtMC_pi0e_TPC);
1483
1484 fSameElecPtMC_pi0e_TPC = new THnSparseD("fSameElecPtMC_pi0e_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1485 fSameElecPtMC_pi0e_TPC->Sumw2();
1486 fOutputList->Add(fSameElecPtMC_pi0e_TPC);
2f8ef315 1487
1488 fIncpTMCpho_eta_TPC = new THnSparseD("fIncpTMCpho_eta_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1489 fIncpTMCpho_eta_TPC->Sumw2();
1490 fOutputList->Add(fIncpTMCpho_eta_TPC);
1491
1492 fPhoElecPtMC_eta_TPC = new THnSparseD("fPhoElecPtMC_eta_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2);
1493 fPhoElecPtMC_eta_TPC->Sumw2();
1494 fOutputList->Add(fPhoElecPtMC_eta_TPC);
1495
1496 fSameElecPtMC_eta_TPC = new THnSparseD("fSameElecPtMC_eta_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1497 fSameElecPtMC_eta_TPC->Sumw2();
1498 fOutputList->Add(fSameElecPtMC_eta_TPC);
1499
1500
697ecf6b 1501 //-------------
e4b0faf2 1502
a6df418c 1503 CheckNclust = new TH1D("CheckNclust","cluster check",200,0,200);
1504 fOutputList->Add(CheckNclust);
1505
1506 CheckNits = new TH1D("CheckNits","ITS cluster check",8,-0.5,7.5);
1507 fOutputList->Add(CheckNits);
ae88597c 1508
1509 CheckDCA = new TH2D("CheckDCA","DCA check",200,-5,5,200,-5,5);
1510 fOutputList->Add(CheckDCA);
68d718e7 1511 /*
d113d7cd 1512 Hpi0pTcheck = new TH2D("Hpi0pTcheck","Pi0 pT from Hijing",100,0,50,3,-0.5,2.5);
1513 fOutputList->Add(Hpi0pTcheck);
1514
19325033 1515 HETApTcheck = new TH2D("HETApTcheck","Eta pT from Hijing",100,0,50,3,-0.5,2.5);
1516 fOutputList->Add(HETApTcheck);
68d718e7 1517 */
1518
1519 Int_t nBinspho3[3] = { 200, 100, 3};
1520 Double_t minpho3[3] = { 0., 0., -0.5};
1521 Double_t maxpho3[3] = {100., 50., 2.5};
19325033 1522
68d718e7 1523 Hpi0pTcheck = new THnSparseD("Hpi0pTcheck","Pi0 pT from Hijing",3,nBinspho3,minpho3,maxpho3);
1524 fOutputList->Add(Hpi0pTcheck);
1525
1526 HETApTcheck = new THnSparseD("HETApTcheck","Eta pT from Hijing",3,nBinspho3,minpho3,maxpho3);
1527 fOutputList->Add(HETApTcheck);
1528 //--
e4b0faf2 1529 HphopTcheck = new TH2D("HphopTcheck","Pho pT from Hijing",100,0,50,3,-0.5,2.5);
1530 fOutputList->Add(HphopTcheck);
5e0e45b3 1531 //
1532 HDpTcheck = new TH2D("HDpTcheck","D pT from Hijing",100,0,50,3,-0.5,2.5);
1533 fOutputList->Add(HDpTcheck);
1534
1535 HBpTcheck = new TH2D("HBpTcheck","B pT from Hijing",100,0,50,3,-0.5,2.5);
1536 fOutputList->Add(HBpTcheck);
1537 //
e4b0faf2 1538
93c189c5 1539 fpTCheck = new TH1D("fpTCheck","pT check",500,0,50);
1540 fOutputList->Add(fpTCheck);
1541
50919258 1542 fMomDtoE = new TH2D("fMomDtoE","D->E pT correlations;e p_{T} GeV/c;D p_{T} GeV/c",400,0,40,400,0,40);
1543 fOutputList->Add(fMomDtoE);
d113d7cd 1544
70448e3b 1545 fLabelCheck = new TH2D("fLabelCheck","MC label",50,0,50,5,-1.5,3.5);
1546 fOutputList->Add(fLabelCheck);
1547
1548 fgeoFake = new TH2D("fgeoFake","Label==0 eta and phi",628,0,6.28,200,-1,1);
1549 fOutputList->Add(fgeoFake);
1550
d19af75d 1551 fFakeTrk0 = new TH2D("fFakeTrk0","fake trakcs",10,0,100,20,0,20);
1552 fOutputList->Add(fFakeTrk0);
1553
1554 fFakeTrk1 = new TH2D("fFakeTrk1","true all e a.f. eID",10,0,100,20,0,20);
1555 fOutputList->Add(fFakeTrk1);
1556
59d998de 1557 ftimingEle = new TH2D("ftimingEle","electron TOF",100,0,20,100,1e-7,1e-6);
1558 fOutputList->Add(ftimingEle);
1559
09022877 1560 // eta correction
1561 // note: parameters 01/31new.TPCnSigmaEtaDep
1562 // 70-90 delta_eta = 0.2
1563
1564 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};
1565 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)
1566 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)
1567 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)
1568 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)
1569 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)
1570 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)
1571 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)
1572
1573 fnSigEtaCorr[0] = new TGraphErrors(12,etaval,corr0); // 0-10
1574 fnSigEtaCorr[1] = new TGraphErrors(12,etaval,corr1); // 10-20
1575 fnSigEtaCorr[2] = new TGraphErrors(12,etaval,corr2); // 20-30
1576 fnSigEtaCorr[3] = new TGraphErrors(12,etaval,corr3); // 30-40
1577 fnSigEtaCorr[4] = new TGraphErrors(12,etaval,corr4); // 40-50
1578 fnSigEtaCorr[5] = new TGraphErrors(12,etaval,corr5); // 50-70
1579 fnSigEtaCorr[6] = new TGraphErrors(12,etaval,corr6); // 70-90
1580
6f4a9952 1581 fIncMaxE = new TH2D("fIncMaxE","Inc",10,0,100,10,0,100);
1582 fOutputList->Add(fIncMaxE);
09022877 1583
c82f1b9d 1584 fIncReco = new TH2D("fIncReco","Inc",10,0,100,100,0,500);
5e38d973 1585 fOutputList->Add(fIncReco);
1586
c82f1b9d 1587 fPhoReco = new TH2D("fPhoReco","Pho",10,0,100,100,0,500);
1588 fOutputList->Add(fPhoReco);
1589
1590 fSamReco = new TH2D("fSamReco","Same",10,0,100,100,0,500);
1591 fOutputList->Add(fSamReco);
1592
1593 fIncRecoMaxE = new TH2D("fIncRecoMaxE","Inc",10,0,100,100,0,500);
1594 fOutputList->Add(fIncRecoMaxE);
1595
1596 fPhoRecoMaxE = new TH2D("fPhoRecoMaxE","Pho",10,0,100,100,0,500);
eabc1527 1597 fOutputList->Add(fPhoRecoMaxE);
5e38d973 1598
c82f1b9d 1599 fSamRecoMaxE = new TH2D("fSamRecoMaxE","Same",10,0,100,100,0,500);
eabc1527 1600 fOutputList->Add(fSamRecoMaxE);
5e38d973 1601
466d4f7a 1602 fPhoVertexReco_TPC = new TH2D("fPhoVertexReco_TPC","photon production Vertex mass selection TPC",40,0,20,250,0,50);
1603 fPhoVertexReco_TPC->Sumw2();
1604 fOutputList->Add(fPhoVertexReco_TPC);
1605
1606 fPhoVertexReco_TPC_Invmass = new TH2D("fPhoVertexReco_TPC_Invmass","photon production Vertex mass selection TPC Invmass",40,0,20,250,0,50);
1607 fPhoVertexReco_TPC_Invmass->Sumw2();
1608 fOutputList->Add(fPhoVertexReco_TPC_Invmass);
05ac9660 1609
1610 fPhoVertexReco_EMCal = new TH2D("fPhoVertexReco_EMCal","photon production Vertex mass selection",40,0,20,250,0,50);
7d495dc5 1611 fPhoVertexReco_EMCal->Sumw2();
05ac9660 1612 fOutputList->Add(fPhoVertexReco_EMCal);
1613
1614 fPhoVertexReco_Invmass = new TH2D("fPhoVertexReco_Invmass","photon production Vertex mass selection",40,0,20,250,0,50);
7d495dc5 1615 fPhoVertexReco_Invmass->Sumw2();
05ac9660 1616 fOutputList->Add(fPhoVertexReco_Invmass);
1617
1618 fPhoVertexReco_step0= new TH2D("fPhoVertexReco_step0","photon production Vertex mass selection",40,0,20,250,0,50);
d26226a1 1619 fPhoVertexReco_step0->Sumw2();
05ac9660 1620 fOutputList->Add(fPhoVertexReco_step0);
1621
1622 fPhoVertexReco_step1= new TH2D("fPhoVertexReco_step1","photon production Vertex mass selection",40,0,20,250,0,50);
d26226a1 1623 fPhoVertexReco_step1->Sumw2();
05ac9660 1624 fOutputList->Add(fPhoVertexReco_step1);
1625
f2fa5a00 1626 fMatchV0_0 = new TH1D("fMatchV0_0","V0 match",100,0,20);
1627 fOutputList->Add(fMatchV0_0);
1628
1629 fMatchV0_1 = new TH1D("fMatchV0_1","V0 match",100,0,20);
1630 fOutputList->Add(fMatchV0_1);
1631
1632 fMatchMC_0 = new TH1D("fMatchMC_0","MC match",100,0,20);
1633 fOutputList->Add(fMatchMC_0);
1634
1635 fMatchMC_1 = new TH1D("fMatchMC_1","MC match",100,0,20);
1636 fOutputList->Add(fMatchMC_1);
1637
36e9774f 1638 fpair = new TH2D("fpair","pair of associate",100,0,20,21,-10.5,10.5);
1639 fOutputList->Add(fpair);
1640
c5a4ddbf 1641 fFakeRejection0 = new TH2D("fFakeRejection0","TPC PID",2,-0.5,1.5,100,0,20);
1c9a06c8 1642 fFakeRejection0->Sumw2();
c5a4ddbf 1643 fOutputList->Add(fFakeRejection0);
1644
1645 fFakeRejection1 = new TH2D("fFakeRejection1","TPC PID + Tr match",2,-0.5,1.5,100,0,20);
1c9a06c8 1646 fFakeRejection1->Sumw2();
c5a4ddbf 1647 fOutputList->Add(fFakeRejection1);
1648
1649 fFakeRejection2 = new TH2D("fFakeRejection2","TPC PID + Tr match + E/p",2,-0.5,1.5,100,0,20);
1c9a06c8 1650 fFakeRejection2->Sumw2();
c5a4ddbf 1651 fOutputList->Add(fFakeRejection2);
1652
32abf0b9 1653 EopFake = new TH2D("EopFake","negative track Eop",20,0,20,200,0,2);
9e88412a 1654 EopFake->Sumw2();
1655 fOutputList->Add(EopFake);
1656
32abf0b9 1657 EopTrue = new TH2D("EopTrue","true track Eop",20,0,20,200,0,2);
9e88412a 1658 EopTrue->Sumw2();
1659 fOutputList->Add(EopTrue);
1660
32abf0b9 1661 MatchFake = new TH2D("MatchFake","negative track Match",20,0,20,100,0,0.05);
9e88412a 1662 MatchFake->Sumw2();
1663 fOutputList->Add(MatchFake);
1664
93c64cd2 1665 MatchTrue = new TH2D("MatchTrue","true track Match",20,0,20,100,0,0.05);
9e88412a 1666 MatchTrue->Sumw2();
1667 fOutputList->Add(MatchTrue);
1668
ebf2f603 1669 MatchTrCheck = new TH2D("MatchTrCheck","Check Tr Match",2010,-10,2000,2010,-10,2000);
1670 fOutputList->Add(MatchTrCheck);
1671
1672 MatchTrEop = new TH2D("MatchTrEop","Check Tr Match E/p",200,0,2,200,0,2);
1673 fOutputList->Add(MatchTrEop);
1674
bfc7c23b 1675 PostData(1,fOutputList);
1676}
1677
1678//________________________________________________________________________
1679void AliAnalysisTaskHFECal::Terminate(Option_t *)
1680{
1681 // Info("Terminate");
1682 AliAnalysisTaskSE::Terminate();
1683}
1684
1685//________________________________________________________________________
1686Bool_t AliAnalysisTaskHFECal::ProcessCutStep(Int_t cutStep, AliVParticle *track)
1687{
1688 // Check single track cuts for a given cut step
1689 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
1690 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
1691 return kTRUE;
1692}
1693//_________________________________________
0f566acb 1694void 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 1695{
1696 //Identify non-heavy flavour electrons using Invariant mass method
1697
fb87d707 1698 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
feffe705 1699 fTrackCuts->SetRequireTPCRefit(kTRUE);
85985bb0 1700 fTrackCuts->SetRequireITSRefit(kTRUE);
feffe705 1701 fTrackCuts->SetEtaRange(-0.9,0.9);
bfc7c23b 1702 fTrackCuts->SetMaxChi2PerClusterTPC(3.5);
feffe705 1703 fTrackCuts->SetMinNClustersTPC(90);
bfc7c23b 1704
c6bca3f2 1705 //const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
6451f693 1706 Double_t bfield = fESD->GetMagneticField();
2f8ef315 1707 Double_t emass = 0.51*0.001; // (0.51 MeV)
6451f693 1708
93c189c5 1709 double ptEle = track->Pt(); //add
1710 if(ibgevent==0 && w > 0.0)
1711 {
1712 fpTCheck->Fill(ptEle,w);
1713 }
1714
bfc7c23b 1715 Bool_t flagPhotonicElec = kFALSE;
f09766dd 1716 Bool_t flagConvinatElec = kFALSE;
bfc7c23b 1717
8806ce6c 1718 int p1 = 0;
1719 if(mce==3)
1720 {
1721 Int_t label = TMath::Abs(track->GetLabel());
1722 TParticle* particle = stack->Particle(label);
1723 p1 = particle->GetFirstMother();
1724 }
1725
36e9774f 1726 int numULS = 0;
1727 int numLS = 0;
1728
5e65985c 1729 for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){
bfc7c23b 1730 AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
1731 if (!trackAsso) {
1732 printf("ERROR: Could not receive track %d\n", jTracks);
1733 continue;
1734 }
5e65985c 1735 if(itrack==jTracks)continue;
d86a815c 1736 int jbgevent = 0;
5e65985c 1737
8806ce6c 1738 int p2 = 0;
1739 if(mce==3)
1740 {
1741 Int_t label2 = TMath::Abs(trackAsso->GetLabel());
1742 TParticle* particle2 = stack->Particle(label2);
d86a815c 1743 Bool_t MChijing_ass = fMC->IsFromBGEvent(label2);
1744 if(MChijing_ass)jbgevent =1;
8806ce6c 1745 if(particle2->GetFirstMother()>-1)
1746 p2 = particle2->GetFirstMother();
1747 }
1748
f09766dd 1749 Double_t dEdxAsso = -999., ptPrim=-999., ptAsso=-999., openingAngle = -999.;
bfc7c23b 1750 Double_t mass=999., width = -999;
1751 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
1752
93c189c5 1753 //ptPrim = track->Pt();
1754 ptPrim = ptEle;
f09766dd 1755
bfc7c23b 1756 dEdxAsso = trackAsso->GetTPCsignal();
1757 ptAsso = trackAsso->Pt();
1758 Int_t chargeAsso = trackAsso->Charge();
1759 Int_t charge = track->Charge();
1760
b12bc641 1761
6451f693 1762 if(ptAsso <fMimpTassCut) continue;
bfc7c23b 1763 if(!fTrackCuts->AcceptTrack(trackAsso)) continue;
6b2cee73 1764 double fTPCnSigmaAss = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(trackAsso, AliPID::kElectron) : 1000;
1765 //cout << "fTPCnSigmaAss = " << fTPCnSigmaAss << endl;
1766 //cout << "fTPCnSigmaAss Cut = " << fMimNsigassCut << endl;
1767 if(fTPCnSigmaAss <fMimNsigassCut || fTPCnSigmaAss>5) continue;
1768 //cout << "fTPCnSigmaAss a.f. cut = " << fTPCnSigmaAss << endl;
bfc7c23b 1769
1770 Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
1771 if(charge>0) fPDGe1 = -11;
1772 if(chargeAsso>0) fPDGe2 = -11;
b12bc641 1773
1774 //printf("chargeAsso = %d\n",chargeAsso);
1775 //printf("charge = %d\n",charge);
bfc7c23b 1776 if(charge == chargeAsso) fFlagLS = kTRUE;
1777 if(charge != chargeAsso) fFlagULS = kTRUE;
1778
b12bc641 1779 //printf("fFlagLS = %d\n",fFlagLS);
1780 //printf("fFlagULS = %d\n",fFlagULS);
cccab38f 1781 //printf("\n");
b12bc641 1782
6451f693 1783 AliKFParticle::SetField(bfield);
bfc7c23b 1784 AliKFParticle ge1(*track, fPDGe1);
1785 AliKFParticle ge2(*trackAsso, fPDGe2);
2f8ef315 1786
1787 if(fSetMassNonlinear)
1788 {
1789 ge1.SetNonlinearMassConstraint(emass);
1790 ge2.SetNonlinearMassConstraint(emass);
1791 }
1792
bfc7c23b 1793 AliKFParticle recg(ge1, ge2);
1794
2f8ef315 1795 /*
d19af75d 1796 // vertex
bfc7c23b 1797 AliKFVertex primV(*pVtx);
1798 primV += recg;
8079a103 1799 primV -= ge1;
1800 primV -= ge2;
bfc7c23b 1801 recg.SetProductionVertex(primV);
2f8ef315 1802 */
1803
d19af75d 1804 // check chi2
1805 if(recg.GetNDF()<1) continue;
a009053a 1806 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
1807 Double_t chi2cut = 3.0;
f9447b3b 1808
a009053a 1809 // mass.
1810 if(fSetMassConstraint)
1811 {
1812 recg.SetMassConstraint(0,0.0001);
1813 chi2cut = 30.0;
1814 }
f9447b3b 1815 recg.GetMass(mass,width);
1816
3c475a9d 1817 if(fSetMassWidthCut && width>1e10)continue;
1818
a009053a 1819 // angle
bfc7c23b 1820 openingAngle = ge1.GetAngle(ge2);
1821 if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
1822 if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
0f566acb 1823
c2801a73 1824 double ishower = 0;
1825 if(shower>0.0 && shower<0.3)ishower = 1;
1826
0f566acb 1827 if(iCal==1)
1828 {
1829 double phoinfo[9];
1830 phoinfo[0] = cent;
1831 phoinfo[1] = ptPrim;
1832 phoinfo[2] = mass;
1833 phoinfo[3] = nSig;
1834 phoinfo[4] = openingAngle;
1835 phoinfo[5] = ishower;
1836 phoinfo[6] = ep;
1837 phoinfo[7] = mce;
1838 phoinfo[8] = ptAsso;
1839
1840 if(fFlagLS) fInvmassLS->Fill(phoinfo);
1841 if(fFlagULS) fInvmassULS->Fill(phoinfo);
1842 if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc->Fill(phoinfo,w);
1843 if(fFlagULS && ibgevent==0 && jbgevent==0)
1844 {
1845 fInvmassULSmc->Fill(phoinfo,w);
1846 }
1847 //printf("fInvmassCut %f\n",fInvmassCut);
1848 //printf("openingAngle %f\n",fOpeningAngleCut);
1849 }
b12bc641 1850
f9447b3b 1851 // angle cut
ffa14dda 1852 if(openingAngle > fOpeningAngleCut) continue;
a009053a 1853 // chi2 cut
6ec70f8c 1854 //if(TMath::Sqrt(TMath::Abs(chi2recg))>chi2cut) continue;
1855 if(chi2recg>chi2cut) continue;
1856
1857 if(fFlagLS ) fInvmassLSreco->Fill(ptPrim,mass);
1858 if(fFlagULS) fInvmassULSreco->Fill(ptPrim,mass);
1859
36e9774f 1860 // check double count
1861 if(mass<fInvmassCut && fFlagULS)numULS++;
1862 if(mass<fInvmassCut && fFlagLS)numLS++;
1863
d86a815c 1864 // for real data
93c189c5 1865 //printf("mce =%f\n",mce);
0f566acb 1866
d86a815c 1867 if(mce<-0.5) // mce==-1. is real
0f566acb 1868 {
1869 //printf("Real data\n");
1870
1871 // count
1872 if(mass<fInvmassCut && fFlagULS && shower>0.0 && shower<0.3 && iCal==1)fPhoElecPtM20Mass->Fill(cent,ptPrim);
1873 if(mass<fInvmassCut && fFlagLS && shower>0.0 && shower<0.3 && iCal==1)fSameElecPtM20Mass->Fill(cent,ptPrim);
1874
1875 // flag
1876 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
1877 flagPhotonicElec = kTRUE;
1878 }
1879 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec){
1880 flagConvinatElec = kTRUE;
1881 }
1882 }
d86a815c 1883 // for MC data
1884 else
0f566acb 1885 {
1886 // count
1887 if(mass<fInvmassCut && fFlagULS && shower>0.0 && shower<0.3 && mce>2.1 && iCal==1)fPhoElecPtMCM20Mass->Fill(cent,ptPrim,w);
1888 if(mass<fInvmassCut && fFlagLS && shower>0.0 && shower<0.3 && mce>2.1 && iCal==1)fSameElecPtMCM20Mass->Fill(cent,ptPrim,w);
1889
1890 //printf("MC data\n");
1891
1892 if(w>0.0)
1893 {
1894 //cout << "tagpi0 = " << tagpi0 << " ; tageta = " << tageta << endl;
1895 if(iCal==1)
1896 {
1897 if(fFlagLS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassLSmc0->Fill(ptPrim,mass);
1898 if(fFlagULS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassULSmc0->Fill(ptPrim,mass);
1899 if(fFlagLS && ibgevent==0 && jbgevent==0 && tageta) fInvmassLSmc1->Fill(ptPrim,mass);
1900 if(fFlagULS && ibgevent==0 && jbgevent==0 && tageta) fInvmassULSmc1->Fill(ptPrim,mass);
1901 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassLSmc2->Fill(ptPrim,mass);
1902 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassULSmc2->Fill(ptPrim,mass);
1903 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassLSmc3->Fill(ptPrim,mass);
1904 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassULSmc3->Fill(ptPrim,mass);
1905 }
1906 }
1907
1908 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (ibgevent==jbgevent)){
1909 //if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (p1==p2)){ //<--- only MC train (55,56) v5-03-68-AN , 69 & v5-05-70-AN (till 74)
1910 flagPhotonicElec = kTRUE;
1911 }
1912 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec && (ibgevent==jbgevent)){
1913 flagConvinatElec = kTRUE;
1914 }
1915 }
d86a815c 1916
36e9774f 1917 } // end of associate loop
1918
6e1d7a8a 1919 if(numULS>0 || numLS>0)
1920 {
1921 int numPair = numULS-numLS;
0f566acb 1922 if(iCal==1)fpair->Fill(ptEle,numPair);
6e1d7a8a 1923 }
1924
0f566acb 1925 fFlagPhotonic = flagPhotonicElec;
1926 fFlagConvinat = flagConvinatElec;
6451f693 1927
1928}
1929
8806ce6c 1930//-------------------------------------------
46305ed6 1931
1932void AliAnalysisTaskHFECal::FindMother(TParticle* part, int &label, int &pid)
1933{
1934 //int label = 99999;
1935 //int pid = 99999;
1936
1937 if(part->GetFirstMother()>-1)
1938 {
1939 label = part->GetFirstMother();
1940 pid = stack->Particle(label)->GetPdgCode();
1941 }
70448e3b 1942 //cout << "Find Mother : label = " << label << " ; pid" << pid << endl;
46305ed6 1943}
1944
8806ce6c 1945double AliAnalysisTaskHFECal::GetMCweight(double mcPi0pT)
1946{
1947 double weight = 1.0;
1948
1949 if(mcPi0pT>0.0 && mcPi0pT<5.0)
1950 {
1951 weight = 0.323*mcPi0pT/(TMath::Exp(-1.6+0.767*mcPi0pT+0.0285*mcPi0pT*mcPi0pT));
1952 }
1953 else
1954 {
1955 weight = 115.0/(0.718*mcPi0pT*TMath::Power(mcPi0pT,3.65));
1956 }
1957 return weight;
1958}
fb87d707 1959
93c189c5 1960double AliAnalysisTaskHFECal::GetMCweightEta(double mcEtapT)
1961{
1962 double weight = 1.0;
1963
1964 weight = 223.3/TMath::Power((TMath::Exp(-0.17*mcEtapT-0.0322*mcEtapT*mcEtapT)+mcEtapT/1.69),5.65);
1965 return weight;
1966}
1967
1968
fb87d707 1969//_________________________________________
1970void AliAnalysisTaskHFECal::FindTriggerClusters()
1971{
bd6ee6dd 1972 //cout << "finding trigger patch" << endl;
fb87d707 1973 // constants
1974 const int nModuleCols = 2;
1975 const int nModuleRows = 5;
1976 const int nColsFeeModule = 48;
1977 const int nRowsFeeModule = 24;
1978 const int nColsFaltroModule = 24;
1979 const int nRowsFaltroModule = 12;
1980 //const int faltroWidthMax = 20;
1981
1982 // part 1, trigger extraction -------------------------------------
1983 Int_t globCol, globRow;
1984 //Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0, trigInCut=0;
1985 Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0;
1986
1987 //Int_t trigtimes[faltroWidthMax];
1988 Double_t cellTime[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1989 Double_t cellEnergy[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1990 //Double_t fTrigCutLow = 6;
1991 //Double_t fTrigCutHigh = 10;
1992 Double_t fTimeCutLow = 469e-09;
1993 Double_t fTimeCutHigh = 715e-09;
1994
1995 AliESDCaloTrigger * fCaloTrigger = fESD->GetCaloTrigger( "EMCAL" );
1996
1997 // erase trigger maps
1998 for(Int_t i = 0; i < nColsFaltroModule*nModuleCols; i++ )
1999 {
2000 for(Int_t j = 0; j < nRowsFaltroModule*nModuleRows; j++ )
2001 {
2002 ftriggersCut[i][j] = 0;
2003 ftriggers[i][j] = 0;
2004 ftriggersTime[i][j] = 0;
2005 }
2006 }
2007
2008 Int_t iglobCol=0, iglobRow=0;
2009 // go through triggers
2010 if( fCaloTrigger->GetEntries() > 0 )
2011 {
2012 // needs reset
2013 fCaloTrigger->Reset();
2014 while( fCaloTrigger->Next() )
2015 {
2016 fCaloTrigger->GetPosition( globCol, globRow );
2017 fCaloTrigger->GetNL0Times( ntimes );
2018 /*
2019 // no L0s
2020 if( ntimes < 1 ) continue;
2021 // get precise timings
2022 fCaloTrigger->GetL0Times( trigtimes );
2023 trigInCut = 0;
2024 for(Int_t i = 0; i < ntimes; i++ )
2025 {
2026 // save the first trigger time in channel
2027 if( i == 0 || triggersTime[globCol][globRow] > trigtimes[i] )
2028 triggersTime[globCol][globRow] = trigtimes[i];
2029 //printf("trigger times: %d\n",trigtimes[i]);
2030 // check if in cut
2031 if(trigtimes[i] > fTrigCutLow && trigtimes[i] < fTrigCutHigh )
2032 trigInCut = 1;
2033
2034 fTrigTimes->Fill(trigtimes[i]);
2035 }
2036 */
2037
2038 //L1 analysis from AliAnalysisTaskEMCALTriggerQA
2039 Int_t bit = 0;
2040 fCaloTrigger->GetTriggerBits(bit);
bd6ee6dd 2041 //cout << "bit = " << bit << endl;
fb87d707 2042
2043 Int_t ts = 0;
2044 fCaloTrigger->GetL1TimeSum(ts);
bd6ee6dd 2045 //cout << "ts = " << ts << endl;
fb87d707 2046 if (ts > 0)ftriggers[globCol][globRow] = 1;
2047 // number of triggered channels in event
2048 nTrigChannel++;
2049 // ... inside cut
2050 if(ts>0 && (bit >> 6 & 0x1))
2051 {
2052 iglobCol = globCol;
2053 iglobRow = globRow;
2054 nTrigChannelCut++;
bd6ee6dd 2055 //cout << "ts cut = " << ts << endl;
2056 //cout << "globCol = " << globCol << endl;
2057 //cout << "globRow = " << globRow << endl;
fb87d707 2058 ftriggersCut[globCol][globRow] = 1;
2059 }
2060
2061 } // calo trigger entries
2062 } // has calo trigger entries
2063
2064 // part 2 go through the clusters here -----------------------------------
bd6ee6dd 2065 //cout << " part 2 go through the clusters here ----------------------------------- " << endl;
fb87d707 2066 Int_t nCluster=0, nCell=0, iCell=0, gCell=0;
60d77596 2067 Short_t cellAddr, nSACell;
2068 Int_t mclabel;
fb87d707 2069 //Int_t nSACell, iSACell, mclabel;
2070 Int_t iSACell;
42c75692 2071 Double_t cellAmp=0, cellTimeT=0, clusterTime=0, efrac=0;
fb87d707 2072 Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta, gphi, geta, feta, fphi;
2073
2074 TRefArray *fCaloClusters = new TRefArray();
2075 fESD->GetEMCALClusters( fCaloClusters );
2076 nCluster = fCaloClusters->GetEntries();
2077
2078
2079 // save all cells times if there are clusters
2080 if( nCluster > 0 ){
2081 // erase time map
2082 for(Int_t i = 0; i < nColsFeeModule*nModuleCols; i++ ){
2083 for(Int_t j = 0; j < nRowsFeeModule*nModuleRows; j++ ){
2084 cellTime[i][j] = 0.;
2085 cellEnergy[i][j] = 0.;
2086 }
2087 }
2088
2089 // get cells
2090 AliESDCaloCells *fCaloCells = fESD->GetEMCALCells();
2091 //AliVCaloCells fCaloCells = fESD->GetEMCALCells();
2092 nSACell = fCaloCells->GetNumberOfCells();
2093 for(iSACell = 0; iSACell < nSACell; iSACell++ ){
2094 // get the cell info *fCal
2095 fCaloCells->GetCell( iSACell, cellAddr, cellAmp, cellTimeT , mclabel, efrac);
2096
2097 // get cell position
2098 fGeom->GetCellIndex( cellAddr, nSupMod, nModule, nIphi, nIeta );
2099 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
2100
2101 // convert co global phi eta
2102 gphi = iphi + nRowsFeeModule*(nSupMod/2);
2103 geta = ieta + nColsFeeModule*(nSupMod%2);
2104
2105 // save cell time and energy
2106 cellTime[geta][gphi] = cellTimeT;
2107 cellEnergy[geta][gphi] = cellAmp;
2108
2109 }
2110 }
2111
2112 Int_t nClusterTrig, nClusterTrigCut;
2113 UShort_t *cellAddrs;
2114 Double_t clsE=-999, clsEta=-999, clsPhi=-999;
2115 Float_t clsPos[3] = {0.,0.,0.};
2116
2117 for(Int_t icl=0; icl<fESD->GetNumberOfCaloClusters(); icl++)
2118 {
2119 AliESDCaloCluster *cluster = fESD->GetCaloCluster(icl);
2120 if(!cluster || !cluster->IsEMCAL()) continue;
2121
2122 // get cluster cells
2123 nCell = cluster->GetNCells();
2124
2125 // get cluster energy
2126 clsE = cluster->E();
2127
2128 // get cluster position
2129 cluster->GetPosition(clsPos);
2130 TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
2131 clsEta = clsPosVec.Eta();
2132 clsPhi = clsPosVec.Phi();
2133
2134 // get the cell addresses
2135 cellAddrs = cluster->GetCellsAbsId();
2136
2137 // check if the cluster contains cell, that was marked as triggered
2138 nClusterTrig = 0;
2139 nClusterTrigCut = 0;
2140
2141 // loop the cells to check, if cluser in acceptance
2142 // any cluster with a cell outside acceptance is not considered
2143 for( iCell = 0; iCell < nCell; iCell++ )
2144 {
42c75692 2145 // check hot cell
feffe705 2146 //if(clsE>6.0)fCellCheck->Fill(clsE,cellAddrs[iCell]);
42c75692 2147
fb87d707 2148 // get cell position
2149 fGeom->GetCellIndex( cellAddrs[iCell], nSupMod, nModule, nIphi, nIeta );
2150 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
2151
2152 // convert co global phi eta
2153 gphi = iphi + nRowsFeeModule*(nSupMod/2);
2154 geta = ieta + nColsFeeModule*(nSupMod%2);
2155
2156 if( cellTime[geta][gphi] > 0. ){
2157 clusterTime += cellTime[geta][gphi];
2158 gCell++;
2159 }
2160
2161 // get corresponding FALTRO
2162 fphi = gphi / 2;
2163 feta = geta / 2;
2164
bd6ee6dd 2165 //cout << "fphi = " << fphi << endl;
2166 //cout << "feta = " << feta << endl;
2167
fb87d707 2168 // try to match with a triggered
2169 if( ftriggers[feta][fphi]==1)
2170 { nClusterTrig++;
2171 }
2172 if( ftriggersCut[feta][fphi]==1)
2173 { nClusterTrigCut++;
2174 }
2175
bd6ee6dd 2176 //cout << "nClusterTrigCut : " << nClusterTrigCut << endl;
2177
fb87d707 2178 } // cells
2179
2180
2181 if( gCell > 0 )
2182 clusterTime = clusterTime / (Double_t)gCell;
2183 // fix the reconstriction code time 100ns jumps
2184 if( fESD->GetBunchCrossNumber() % 4 < 2 )
2185 clusterTime -= 0.0000001;
2186
feffe705 2187 //fClsETime->Fill(clsE,clusterTime);
2188 //fClsEBftTrigCut->Fill(clsE);
fb87d707 2189
2190 if(nClusterTrig>0){
feffe705 2191 //fClsETime1->Fill(clsE,clusterTime);
fb87d707 2192 }
2193
2194 if(nClusterTrig>0){
2195 cluster->SetChi2(1);
feffe705 2196 //fClsEAftTrigCut1->Fill(clsE);
fb87d707 2197 }
2198
2199 if(nClusterTrigCut>0){
2200 cluster->SetChi2(2);
feffe705 2201 //fClsEAftTrigCut2->Fill(clsE);
fb87d707 2202 }
2203
2204 if(nClusterTrigCut>0 && ( nCell > (1 + clsE / 3)))
2205 {
2206 cluster->SetChi2(3);
feffe705 2207 //fClsEAftTrigCut3->Fill(clsE);
fb87d707 2208 }
2209
2210 if(nClusterTrigCut>0 && (nCell > (1 + clsE / 3) )&&( clusterTime > fTimeCutLow && clusterTime < fTimeCutHigh ))
2211 {
2212 // cluster->SetChi2(4);
feffe705 2213 //fClsEAftTrigCut4->Fill(clsE);
fb87d707 2214 }
2215 if(nClusterTrigCut<1)
2216 {
2217 cluster->SetChi2(0);
2218
feffe705 2219 //fClsEAftTrigCut->Fill(clsE);
fb87d707 2220 }
2221
2222 } // clusters
2223}
2224
09022877 2225// <-------- only MC correction
8efacc22 2226double AliAnalysisTaskHFECal::MCEopMeanCorrection(double pTmc, float central)
8af58b0d 2227{
8efacc22 2228 TF1 *fcorr0 = new TF1("fcorr0","[0]*tanh([1]+[2]*x)");
2229 TF1 *fcorr1 = new TF1("fcorr1","[0]*tanh([1]+[2]*x)");
2230
8af58b0d 2231 double shift = 0.0;
8efacc22 2232
2233 if(central>0 && central<=10)
2234 {
2235 fcorr0->SetParameters(1.045,1.288,3.18e-01); //
2236 fcorr1->SetParameters(9.91e-01,3.466,2.344);
2237 }
09022877 2238 else if(central>10 && central<=20)
8efacc22 2239 {
2240 fcorr0->SetParameters(1.029,8.254e-01,4.07e-01);
2241 fcorr1->SetParameters(0.975,2.276,1.501e-01);
2242 }
2243 else if(central>20 && central<=30)
2244 {
2245 fcorr0->SetParameters(1.01,8.795e-01,3.904e-01);
2246 fcorr1->SetParameters(9.675e-01,1.654,2.583e-01);
2247 }
2248 else if(central>30 && central<=40)
2249 {
2250 fcorr0->SetParameters(1.00,1.466,2.305e-1);
2251 fcorr1->SetParameters(9.637e-01,1.473,2.754e-01);
2252 }
2253 else if(central>40 && central<=50)
2254 {
fd2126aa 2255 fcorr0->SetParameters(1.00,1.422,1.518e-01);
8efacc22 2256 fcorr1->SetParameters(9.59e-01,1.421,2.931e-01);
2257 }
2258
2259 else if(central>50 && central<=70)
2260 {
2261 fcorr0->SetParameters(0.989,2.495,2.167);
e1f0fb74 2262 fcorr1->SetParameters(0.961,1.734,1.438e-01);
8efacc22 2263 }
2264 else if(central>70 && central<=100)
2265 {
2266 fcorr0->SetParameters(0.981,-3.373,3.93327);
fd2126aa 2267 fcorr1->SetParameters(9.574e-01,1.698,1.58e-01);
8efacc22 2268 }
2269
2270
2271 shift = fcorr0->Eval(pTmc)-fcorr1->Eval(pTmc);
2272
89f41a30 2273 fcorr0->Delete();
2274 fcorr1->Delete();
2275
8af58b0d 2276 return shift;
2277}
fb87d707 2278
09022877 2279// <-------- only Data correction
2280double AliAnalysisTaskHFECal::NsigmaCorrection(double tmpeta, float central)
2281{
2282 int icent = 0;
2283
2284 if(central>=0 && central<10)
2285 {
2286 icent = 0;
2287 }
2288 else if(central>=10 && central<20)
2289 {
2290 icent = 1;
2291 }
2292 else if(central>=20 && central<30)
2293 {
2294 icent = 2;
2295 }
2296 else if(central>=30 && central<40)
2297 {
2298 icent = 3;
2299 }
2300 else if(central>=40 && central<50)
2301 {
2302 icent = 4;
2303 }
2304 else if(central>=50 && central<70)
2305 {
2306 icent = 5;
2307 }
2308 else
2309 {
2310 icent = 6;
2311 }
2312
2313 double shift = fnSigEtaCorr[icent]->Eval(tmpeta);
2314
2315 //cout << "eta correction"<< endl;
2316 //cout << "cent = "<< central<< endl;
2317 //cout << "icent = "<< icent << endl;
2318 //cout << "shift = "<< shift << endl;
2319
2320 return shift;
fb87d707 2321
09022877 2322}
fb87d707 2323