]> 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
modifications
[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);
a6df418c 1180 fCuts->SetCutITSpixel(AliHFEextraCuts::kAny);
1181 fCuts->SetCheckITSLayerStatus(kFALSE);
1182 fCuts->SetVertexRange(10.);
1183 fCuts->SetTOFPIDStep(kFALSE);
fbe402d7 1184 //fCuts->SetPtRange(2, 50);
1185 fCuts->SetPtRange(0.1, 50);
0f566acb 1186 //fCuts->SetMaxImpactParam(3.,3.);
1187 fCuts->SetMaxImpactParam(2.4,3.2); // standard in 2011
a6df418c 1188
bfc7c23b 1189 //--------Initialize correction Framework and Cuts
1190 fCFM = new AliCFManager;
1191 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
1192 fCFM->SetNStepParticle(kNcutSteps);
1193 for(Int_t istep = 0; istep < kNcutSteps; istep++)
1194 fCFM->SetParticleCutsList(istep, NULL);
1195
1196 if(!fCuts){
1197 AliWarning("Cuts not available. Default cuts will be used");
1198 fCuts = new AliHFEcuts;
1199 fCuts->CreateStandardCuts();
1200 }
1201 fCuts->Initialize(fCFM);
1202
1203 //---------Output Tlist
1204 fOutputList = new TList();
1205 fOutputList->SetOwner();
1206 fOutputList->Add(fPIDqa->MakeList("PIDQA"));
1207
1208 fNoEvents = new TH1F("fNoEvents","",4,-0.5,3.5) ;
1209 fOutputList->Add(fNoEvents);
1210
02cd82a0 1211 Int_t binsE[5] = {250, 100, 1000, 200, 10};
fb87d707 1212 Double_t xminE[5] = {1.0, -1, 0.0, 0, -0.5};
1213 Double_t xmaxE[5] = {3.5, 1, 100.0, 100, 9.5};
1214 fEMCAccE = new THnSparseD("fEMCAccE","EMC acceptance & E;#eta;#phi;Energy;Centrality;trugCondition;",5,binsE,xminE,xmaxE);
c82f1b9d 1215 //if(fqahist==1)fOutputList->Add(fEMCAccE);
bfc7c23b 1216
f4e0d2d5 1217 hEMCAccE = new TH2F("hEMCAccE","Cluster Energy",200,0,100,100,0,20);
1218 fOutputList->Add(hEMCAccE);
1219
bfc7c23b 1220 fTrkpt = new TH1F("fTrkpt","track pt",100,0,50);
1221 fOutputList->Add(fTrkpt);
1222
1223 fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50);
1224 fOutputList->Add(fTrackPtBefTrkCuts);
1225
1226 fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50);
1227 fOutputList->Add(fTrackPtAftTrkCuts);
1228
1229 fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10);
1230 fOutputList->Add(fTPCnsigma);
1231
1232 fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2);
1233 fOutputList->Add(fTrkEovPBef);
1234
1235 fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2);
1236 fOutputList->Add(fTrkEovPAft);
1237
42c75692 1238 fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,200,-10,10);
bfc7c23b 1239 fOutputList->Add(fdEdxBef);
1240
42c75692 1241 fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,200,-10,10);
bfc7c23b 1242 fOutputList->Add(fdEdxAft);
1243
02cd82a0 1244 fIncpT = new TH2F("fIncpT","HFE pid electro vs. centrality",200,0,100,100,0,50);
bfc7c23b 1245 fOutputList->Add(fIncpT);
1246
02cd82a0 1247 fIncpTM20 = new TH2F("fIncpTM20","HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
fb87d707 1248 fOutputList->Add(fIncpTM20);
4e01b68c 1249
5e38d973 1250 Int_t nBinspho[9] = { 10, 30, 600, 60, 50, 4, 40, 8, 30};
1251 Double_t minpho[9] = { 0., 0., -0.1, 40, 0, -0.5, 0,-1.5, 0};
1252 Double_t maxpho[9] = {100., 30., 0.5, 100, 1, 3.5, 2, 6.5, 30};
f09766dd 1253
b12bc641 1254 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 1255 if(fqahist==1)fOutputList->Add(fInvmassLS);
bfc7c23b 1256
b12bc641 1257 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 1258 if(fqahist==1)fOutputList->Add(fInvmassULS);
bfc7c23b 1259
8806ce6c 1260 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 1261 if(fqahist==1)fOutputList->Add(fInvmassLSmc);
8806ce6c 1262
1263 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 1264 if(fqahist==1)fOutputList->Add(fInvmassULSmc);
8806ce6c 1265
d89866fd 1266 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 1267 fInvmassLSreco->Sumw2();
1268 fOutputList->Add(fInvmassLSreco);
1269
d89866fd 1270 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 1271 fInvmassULSreco->Sumw2();
1272 fOutputList->Add(fInvmassULSreco);
1273
d89866fd 1274 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 1275 fInvmassLSmc0->Sumw2();
1276 fOutputList->Add(fInvmassLSmc0);
1277
d89866fd 1278 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 1279 fInvmassLSmc1->Sumw2();
1280 fOutputList->Add(fInvmassLSmc1);
1281
d89866fd 1282 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 1283 fInvmassLSmc2->Sumw2();
1284 fOutputList->Add(fInvmassLSmc2);
1285
d89866fd 1286 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 1287 fInvmassLSmc3->Sumw2();
1288 fOutputList->Add(fInvmassLSmc3);
1289
d89866fd 1290 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 1291 fInvmassULSmc0->Sumw2();
1292 fOutputList->Add(fInvmassULSmc0);
1293
d89866fd 1294 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 1295 fInvmassULSmc1->Sumw2();
1296 fOutputList->Add(fInvmassULSmc1);
1297
d89866fd 1298 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 1299 fInvmassULSmc2->Sumw2();
1300 fOutputList->Add(fInvmassULSmc2);
1301
d89866fd 1302 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 1303 fInvmassULSmc3->Sumw2();
93c189c5 1304 fOutputList->Add(fInvmassULSmc3);
1305
bfc7c23b 1306 fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
1307 fOutputList->Add(fOpeningAngleLS);
1308
1309 fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
1310 fOutputList->Add(fOpeningAngleULS);
1311
1312 fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50);
1313 fOutputList->Add(fPhotoElecPt);
1314
02cd82a0 1315 fPhoElecPt = new TH2F("fPhoElecPt", "Pho-inclusive electron pt",200,0,100,100,0,50);
bfc7c23b 1316 fOutputList->Add(fPhoElecPt);
1317
02cd82a0 1318 fPhoElecPtM20 = new TH2F("fPhoElecPtM20", "Pho-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 1319 fOutputList->Add(fPhoElecPtM20);
1320
0f566acb 1321 fPhoElecPtM20Mass = new TH2F("fPhoElecPtM20Mass", "Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1322 fPhoElecPtM20Mass->Sumw2();
1323 fOutputList->Add(fPhoElecPtM20Mass);
1324
02cd82a0 1325 fSameElecPt = new TH2F("fSameElecPt", "Same-inclusive electron pt",200,0,100,100,0,50);
f09766dd 1326 fOutputList->Add(fSameElecPt);
1327
02cd82a0 1328 fSameElecPtM20 = new TH2F("fSameElecPtM20", "Same-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 1329 fOutputList->Add(fSameElecPtM20);
1330
0f566acb 1331 fSameElecPtM20Mass = new TH2F("fSameElecPtM20Mass", "Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1332 fSameElecPtM20Mass->Sumw2();
1333 fOutputList->Add(fSameElecPtM20Mass);
1334
ebef1da9 1335 fSemiElecPtM20 = new TH2F("fSemiElecPtM20", "Semi-inclusive electron pt with M20",200,0,100,100,0,50);
1336 fOutputList->Add(fSemiElecPtM20);
1337
02cd82a0 1338 fCent = new TH1F("fCent","Centrality",200,0,100) ;
bfc7c23b 1339 fOutputList->Add(fCent);
1340
1341 // Make common binning
44be9e1d 1342 const Double_t kMinP = 0.;
d19af75d 1343 const Double_t kMaxP = 20.;
bfc7c23b 1344
05ac9660 1345 //+++ 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select)
1346 // 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)
1347 Int_t nBins[16] = { 100, 7, 60, 20, 90, 100, 25, 60, 60, 100, 40, 10, 250, 100, 100, 8};
1348 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};
1349 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};
1350 //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);
1351 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 1352 //if(fqahist==1)fOutputList->Add(fEleInfo);
bfc7c23b 1353
e1f0fb74 1354 // Make common binning
1355 Int_t nBinsEop[3] = { 10, 50, 100};
1356 Double_t minEop[3] = { 0, 0, -5};
1357 Double_t maxEop[3] = {100, 50, 5};
1358 fElenSigma= new THnSparseD("fElenSigma", "Electron nSigma; cent; pT [GeV/c]; nSigma", 3, nBinsEop, minEop, maxEop);
1359 fOutputList->Add(fElenSigma);
1360
1361
fb87d707 1362 //<--- Trigger info
feffe705 1363 /*
fb87d707 1364 fClsEBftTrigCut = new TH1F("fClsEBftTrigCut","cluster E before trigger selection",1000,0,100);
1365 fOutputList->Add(fClsEBftTrigCut);
1366
1367 fClsEAftTrigCut = new TH1F("fClsEAftTrigCut","cluster E if cls has 0 trigcut channel",1000,0,100);
1368 fOutputList->Add(fClsEAftTrigCut);
1369
1370 fClsEAftTrigCut1 = new TH1F("fClsEAftTrigCut1","cluster E if cls with trig channel",1000,0,100);
1371 fOutputList->Add(fClsEAftTrigCut1);
1372
1373 fClsEAftTrigCut2 = new TH1F("fClsEAftTrigCut2","cluster E if cls with trigcut channel",1000,0,100);
1374 fOutputList->Add(fClsEAftTrigCut2);
1375
1376 fClsEAftTrigCut3 = new TH1F("fClsEAftTrigCut3","cluster E if cls with trigcut channel + nCell>Ecorrect",1000,0,100);
1377 fOutputList->Add(fClsEAftTrigCut3);
1378
1379 fClsEAftTrigCut4 = new TH1F("fClsEAftTrigCut4","cluster E if cls with trigcut channel + nCell>Ecorrect + cls time cut",1000,0,100);
1380 fOutputList->Add(fClsEAftTrigCut4);
1381
1382 fClsETime = new TH2F("fClsETime", "Cls time vs E; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
1383 fOutputList->Add(fClsETime);
1384
1385 fClsETime1 = new TH2F("fClsETime1", "Cls time vs E if cls contains trigger channel; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
1386 fOutputList->Add(fClsETime1);
1387
1388 fTrigTimes = new TH1F("fTrigTimes", "Trigger time; time; N;",25,0,25);
1389 fOutputList->Add(fTrigTimes);
1390
42c75692 1391 fCellCheck = new TH2F("fCellCheck", "Cell vs E; E GeV; Cell ID",10,6,26,12000,0,12000);
1392 fOutputList->Add(fCellCheck);
feffe705 1393 */
e7d87aef 1394 //<---------- MC
1395
02cd82a0 1396 fInputHFEMC = new TH2F("fInputHFEMC","Input MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 1397 fOutputList->Add(fInputHFEMC);
1398
02cd82a0 1399 fInputAlle = new TH2F("fInputAlle","Input MC electro vs. centrality",200,0,100,100,0,50);
feffe705 1400 fOutputList->Add(fInputAlle);
1401
02cd82a0 1402 fIncpTMChfe = new TH2F("fIncpTMChfe","MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 1403 fOutputList->Add(fIncpTMChfe);
1404
02cd82a0 1405 fIncpTMChfeAll = new TH2F("fIncpTMChfeAll","MC Alle pid electro vs. centrality",200,0,100,100,0,50);
3db00c72 1406 fOutputList->Add(fIncpTMChfeAll);
1407
02cd82a0 1408 fIncpTMCM20hfe = new TH2F("fIncpTMCM20hfe","MC HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
e7d87aef 1409 fOutputList->Add(fIncpTMCM20hfe);
1410
02cd82a0 1411 fIncpTMCM20hfeAll = new TH2F("fIncpTMCM20hfeAll","MC Alle pid electro vs. centrality with M20",200,0,100,100,0,50);
3db00c72 1412 fOutputList->Add(fIncpTMCM20hfeAll);
1413
60544aea 1414 fIncpTMCM20hfeCheck = new TH2F("fIncpTMCM20hfeCheck","MC HFE pid electro vs. centrality with M20 Check",200,0,100,100,0,50);
1415 fOutputList->Add(fIncpTMCM20hfeCheck);
44be9e1d 1416
2b4460a6 1417 Int_t nBinspho2[5] = { 200, 100, 7, 3, 700};
987053ce 1418 Double_t minpho2[5] = { 0., 0., -2.5, -0.5, 0.};
2b4460a6 1419 Double_t maxpho2[5] = {100., 50., 4.5, 2.5, 70.};
bd6ee6dd 1420
1421 fInputHFEMC_weight = new THnSparseD("fInputHFEMC_weight", "MC HFE electron pt",5,nBinspho2,minpho2,maxpho2);
1422 fOutputList->Add(fInputHFEMC_weight);
1423
1424 fIncpTMCM20hfeCheck_weight = new THnSparseD("fIncpTMCM20hfeCheck_weight", "HFE electron pt with M20",5,nBinspho2,minpho2,maxpho2);
1425 fOutputList->Add(fIncpTMCM20hfeCheck_weight);
1426
987053ce 1427 fIncpTMCpho = new THnSparseD("fIncpTMCpho","MC Pho pid electro vs. centrality",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1428 fOutputList->Add(fIncpTMCpho);
1429
987053ce 1430 fIncpTMCM20pho = new THnSparseD("fIncpTMCM20pho","MC Pho pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1431 fOutputList->Add(fIncpTMCM20pho);
1432
987053ce 1433 fPhoElecPtMC = new THnSparseD("fPhoElecPtMC", "MC Pho-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1434 fOutputList->Add(fPhoElecPtMC);
1435
987053ce 1436 fPhoElecPtMCM20 = new THnSparseD("fPhoElecPtMCM20", "MC Pho-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
feffe705 1437 fOutputList->Add(fPhoElecPtMCM20);
e7d87aef 1438
0f566acb 1439 fPhoElecPtMCM20Mass = new TH2D("fPhoElecPtMCM20Mass", "MC Pho-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1440 fOutputList->Add(fPhoElecPtMCM20Mass);
1441
987053ce 1442 fSameElecPtMC = new THnSparseD("fSameElecPtMC", "MC Same-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1443 fOutputList->Add(fSameElecPtMC);
1444
987053ce 1445 fSameElecPtMCM20 = new THnSparseD("fSameElecPtMCM20", "MC Same-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1446 fOutputList->Add(fSameElecPtMCM20);
1447
0f566acb 1448 fSameElecPtMCM20Mass = new TH2D("fSameElecPtMCM20Mass", "MC Same-inclusive electron pt with M20 Mass",200,0,100,100,0,50);
1449 fOutputList->Add(fSameElecPtMCM20Mass);
1450
987053ce 1451 fIncpTMCM20pho_pi0e = new THnSparseD("fIncpTMCM20pho_pi0e","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1452 fIncpTMCM20pho_pi0e->Sumw2();
e4b0faf2 1453 fOutputList->Add(fIncpTMCM20pho_pi0e);
1454
987053ce 1455 fPhoElecPtMCM20_pi0e = new THnSparseD("fPhoElecPtMCM20_pi0e", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1456 fPhoElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1457 fOutputList->Add(fPhoElecPtMCM20_pi0e);
1458
987053ce 1459 fSameElecPtMCM20_pi0e = new THnSparseD("fSameElecPtMCM20_pi0e", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1460 fSameElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1461 fOutputList->Add(fSameElecPtMCM20_pi0e);
697ecf6b 1462 //
93c189c5 1463 fIncpTMCM20pho_eta = new THnSparseD("fIncpTMCM20pho_eta","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1464 fIncpTMCM20pho_eta->Sumw2();
1465 fOutputList->Add(fIncpTMCM20pho_eta);
1466
1467 fPhoElecPtMCM20_eta = new THnSparseD("fPhoElecPtMCM20_eta", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1468 fPhoElecPtMCM20_eta->Sumw2();
1469 fOutputList->Add(fPhoElecPtMCM20_eta);
1470
1471 fSameElecPtMCM20_eta = new THnSparseD("fSameElecPtMCM20_eta", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1472 fSameElecPtMCM20_eta->Sumw2();
1473 fOutputList->Add(fSameElecPtMCM20_eta);
697ecf6b 1474 // ------------
1475 fIncpTMCpho_pi0e_TPC = new THnSparseD("fIncpTMCpho_pi0e_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1476 fIncpTMCpho_pi0e_TPC->Sumw2();
1477 fOutputList->Add(fIncpTMCpho_pi0e_TPC);
1478
1479 fPhoElecPtMC_pi0e_TPC = new THnSparseD("fPhoElecPtMC_pi0e_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2);
1480 fPhoElecPtMC_pi0e_TPC->Sumw2();
1481 fOutputList->Add(fPhoElecPtMC_pi0e_TPC);
1482
1483 fSameElecPtMC_pi0e_TPC = new THnSparseD("fSameElecPtMC_pi0e_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1484 fSameElecPtMC_pi0e_TPC->Sumw2();
1485 fOutputList->Add(fSameElecPtMC_pi0e_TPC);
2f8ef315 1486
1487 fIncpTMCpho_eta_TPC = new THnSparseD("fIncpTMCpho_eta_TPC","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1488 fIncpTMCpho_eta_TPC->Sumw2();
1489 fOutputList->Add(fIncpTMCpho_eta_TPC);
1490
1491 fPhoElecPtMC_eta_TPC = new THnSparseD("fPhoElecPtMC_eta_TPC", "MC Pho-inclusive electron pt with pi0->e",5,nBinspho2,minpho2,maxpho2);
1492 fPhoElecPtMC_eta_TPC->Sumw2();
1493 fOutputList->Add(fPhoElecPtMC_eta_TPC);
1494
1495 fSameElecPtMC_eta_TPC = new THnSparseD("fSameElecPtMC_eta_TPC", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1496 fSameElecPtMC_eta_TPC->Sumw2();
1497 fOutputList->Add(fSameElecPtMC_eta_TPC);
1498
1499
697ecf6b 1500 //-------------
e4b0faf2 1501
a6df418c 1502 CheckNclust = new TH1D("CheckNclust","cluster check",200,0,200);
1503 fOutputList->Add(CheckNclust);
1504
1505 CheckNits = new TH1D("CheckNits","ITS cluster check",8,-0.5,7.5);
1506 fOutputList->Add(CheckNits);
ae88597c 1507
1508 CheckDCA = new TH2D("CheckDCA","DCA check",200,-5,5,200,-5,5);
1509 fOutputList->Add(CheckDCA);
68d718e7 1510 /*
d113d7cd 1511 Hpi0pTcheck = new TH2D("Hpi0pTcheck","Pi0 pT from Hijing",100,0,50,3,-0.5,2.5);
1512 fOutputList->Add(Hpi0pTcheck);
1513
19325033 1514 HETApTcheck = new TH2D("HETApTcheck","Eta pT from Hijing",100,0,50,3,-0.5,2.5);
1515 fOutputList->Add(HETApTcheck);
68d718e7 1516 */
1517
1518 Int_t nBinspho3[3] = { 200, 100, 3};
1519 Double_t minpho3[3] = { 0., 0., -0.5};
1520 Double_t maxpho3[3] = {100., 50., 2.5};
19325033 1521
68d718e7 1522 Hpi0pTcheck = new THnSparseD("Hpi0pTcheck","Pi0 pT from Hijing",3,nBinspho3,minpho3,maxpho3);
1523 fOutputList->Add(Hpi0pTcheck);
1524
1525 HETApTcheck = new THnSparseD("HETApTcheck","Eta pT from Hijing",3,nBinspho3,minpho3,maxpho3);
1526 fOutputList->Add(HETApTcheck);
1527 //--
e4b0faf2 1528 HphopTcheck = new TH2D("HphopTcheck","Pho pT from Hijing",100,0,50,3,-0.5,2.5);
1529 fOutputList->Add(HphopTcheck);
5e0e45b3 1530 //
1531 HDpTcheck = new TH2D("HDpTcheck","D pT from Hijing",100,0,50,3,-0.5,2.5);
1532 fOutputList->Add(HDpTcheck);
1533
1534 HBpTcheck = new TH2D("HBpTcheck","B pT from Hijing",100,0,50,3,-0.5,2.5);
1535 fOutputList->Add(HBpTcheck);
1536 //
e4b0faf2 1537
93c189c5 1538 fpTCheck = new TH1D("fpTCheck","pT check",500,0,50);
1539 fOutputList->Add(fpTCheck);
1540
50919258 1541 fMomDtoE = new TH2D("fMomDtoE","D->E pT correlations;e p_{T} GeV/c;D p_{T} GeV/c",400,0,40,400,0,40);
1542 fOutputList->Add(fMomDtoE);
d113d7cd 1543
70448e3b 1544 fLabelCheck = new TH2D("fLabelCheck","MC label",50,0,50,5,-1.5,3.5);
1545 fOutputList->Add(fLabelCheck);
1546
1547 fgeoFake = new TH2D("fgeoFake","Label==0 eta and phi",628,0,6.28,200,-1,1);
1548 fOutputList->Add(fgeoFake);
1549
d19af75d 1550 fFakeTrk0 = new TH2D("fFakeTrk0","fake trakcs",10,0,100,20,0,20);
1551 fOutputList->Add(fFakeTrk0);
1552
1553 fFakeTrk1 = new TH2D("fFakeTrk1","true all e a.f. eID",10,0,100,20,0,20);
1554 fOutputList->Add(fFakeTrk1);
1555
59d998de 1556 ftimingEle = new TH2D("ftimingEle","electron TOF",100,0,20,100,1e-7,1e-6);
1557 fOutputList->Add(ftimingEle);
1558
09022877 1559 // eta correction
1560 // note: parameters 01/31new.TPCnSigmaEtaDep
1561 // 70-90 delta_eta = 0.2
1562
1563 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};
1564 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)
1565 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)
1566 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)
1567 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)
1568 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)
1569 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)
1570 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)
1571
1572 fnSigEtaCorr[0] = new TGraphErrors(12,etaval,corr0); // 0-10
1573 fnSigEtaCorr[1] = new TGraphErrors(12,etaval,corr1); // 10-20
1574 fnSigEtaCorr[2] = new TGraphErrors(12,etaval,corr2); // 20-30
1575 fnSigEtaCorr[3] = new TGraphErrors(12,etaval,corr3); // 30-40
1576 fnSigEtaCorr[4] = new TGraphErrors(12,etaval,corr4); // 40-50
1577 fnSigEtaCorr[5] = new TGraphErrors(12,etaval,corr5); // 50-70
1578 fnSigEtaCorr[6] = new TGraphErrors(12,etaval,corr6); // 70-90
1579
6f4a9952 1580 fIncMaxE = new TH2D("fIncMaxE","Inc",10,0,100,10,0,100);
1581 fOutputList->Add(fIncMaxE);
09022877 1582
c82f1b9d 1583 fIncReco = new TH2D("fIncReco","Inc",10,0,100,100,0,500);
5e38d973 1584 fOutputList->Add(fIncReco);
1585
c82f1b9d 1586 fPhoReco = new TH2D("fPhoReco","Pho",10,0,100,100,0,500);
1587 fOutputList->Add(fPhoReco);
1588
1589 fSamReco = new TH2D("fSamReco","Same",10,0,100,100,0,500);
1590 fOutputList->Add(fSamReco);
1591
1592 fIncRecoMaxE = new TH2D("fIncRecoMaxE","Inc",10,0,100,100,0,500);
1593 fOutputList->Add(fIncRecoMaxE);
1594
1595 fPhoRecoMaxE = new TH2D("fPhoRecoMaxE","Pho",10,0,100,100,0,500);
eabc1527 1596 fOutputList->Add(fPhoRecoMaxE);
5e38d973 1597
c82f1b9d 1598 fSamRecoMaxE = new TH2D("fSamRecoMaxE","Same",10,0,100,100,0,500);
eabc1527 1599 fOutputList->Add(fSamRecoMaxE);
5e38d973 1600
466d4f7a 1601 fPhoVertexReco_TPC = new TH2D("fPhoVertexReco_TPC","photon production Vertex mass selection TPC",40,0,20,250,0,50);
1602 fPhoVertexReco_TPC->Sumw2();
1603 fOutputList->Add(fPhoVertexReco_TPC);
1604
1605 fPhoVertexReco_TPC_Invmass = new TH2D("fPhoVertexReco_TPC_Invmass","photon production Vertex mass selection TPC Invmass",40,0,20,250,0,50);
1606 fPhoVertexReco_TPC_Invmass->Sumw2();
1607 fOutputList->Add(fPhoVertexReco_TPC_Invmass);
05ac9660 1608
1609 fPhoVertexReco_EMCal = new TH2D("fPhoVertexReco_EMCal","photon production Vertex mass selection",40,0,20,250,0,50);
7d495dc5 1610 fPhoVertexReco_EMCal->Sumw2();
05ac9660 1611 fOutputList->Add(fPhoVertexReco_EMCal);
1612
1613 fPhoVertexReco_Invmass = new TH2D("fPhoVertexReco_Invmass","photon production Vertex mass selection",40,0,20,250,0,50);
7d495dc5 1614 fPhoVertexReco_Invmass->Sumw2();
05ac9660 1615 fOutputList->Add(fPhoVertexReco_Invmass);
1616
1617 fPhoVertexReco_step0= new TH2D("fPhoVertexReco_step0","photon production Vertex mass selection",40,0,20,250,0,50);
d26226a1 1618 fPhoVertexReco_step0->Sumw2();
05ac9660 1619 fOutputList->Add(fPhoVertexReco_step0);
1620
1621 fPhoVertexReco_step1= new TH2D("fPhoVertexReco_step1","photon production Vertex mass selection",40,0,20,250,0,50);
d26226a1 1622 fPhoVertexReco_step1->Sumw2();
05ac9660 1623 fOutputList->Add(fPhoVertexReco_step1);
1624
f2fa5a00 1625 fMatchV0_0 = new TH1D("fMatchV0_0","V0 match",100,0,20);
1626 fOutputList->Add(fMatchV0_0);
1627
1628 fMatchV0_1 = new TH1D("fMatchV0_1","V0 match",100,0,20);
1629 fOutputList->Add(fMatchV0_1);
1630
1631 fMatchMC_0 = new TH1D("fMatchMC_0","MC match",100,0,20);
1632 fOutputList->Add(fMatchMC_0);
1633
1634 fMatchMC_1 = new TH1D("fMatchMC_1","MC match",100,0,20);
1635 fOutputList->Add(fMatchMC_1);
1636
36e9774f 1637 fpair = new TH2D("fpair","pair of associate",100,0,20,21,-10.5,10.5);
1638 fOutputList->Add(fpair);
1639
c5a4ddbf 1640 fFakeRejection0 = new TH2D("fFakeRejection0","TPC PID",2,-0.5,1.5,100,0,20);
1c9a06c8 1641 fFakeRejection0->Sumw2();
c5a4ddbf 1642 fOutputList->Add(fFakeRejection0);
1643
1644 fFakeRejection1 = new TH2D("fFakeRejection1","TPC PID + Tr match",2,-0.5,1.5,100,0,20);
1c9a06c8 1645 fFakeRejection1->Sumw2();
c5a4ddbf 1646 fOutputList->Add(fFakeRejection1);
1647
1648 fFakeRejection2 = new TH2D("fFakeRejection2","TPC PID + Tr match + E/p",2,-0.5,1.5,100,0,20);
1c9a06c8 1649 fFakeRejection2->Sumw2();
c5a4ddbf 1650 fOutputList->Add(fFakeRejection2);
1651
32abf0b9 1652 EopFake = new TH2D("EopFake","negative track Eop",20,0,20,200,0,2);
9e88412a 1653 EopFake->Sumw2();
1654 fOutputList->Add(EopFake);
1655
32abf0b9 1656 EopTrue = new TH2D("EopTrue","true track Eop",20,0,20,200,0,2);
9e88412a 1657 EopTrue->Sumw2();
1658 fOutputList->Add(EopTrue);
1659
32abf0b9 1660 MatchFake = new TH2D("MatchFake","negative track Match",20,0,20,100,0,0.05);
9e88412a 1661 MatchFake->Sumw2();
1662 fOutputList->Add(MatchFake);
1663
93c64cd2 1664 MatchTrue = new TH2D("MatchTrue","true track Match",20,0,20,100,0,0.05);
9e88412a 1665 MatchTrue->Sumw2();
1666 fOutputList->Add(MatchTrue);
1667
ebf2f603 1668 MatchTrCheck = new TH2D("MatchTrCheck","Check Tr Match",2010,-10,2000,2010,-10,2000);
1669 fOutputList->Add(MatchTrCheck);
1670
1671 MatchTrEop = new TH2D("MatchTrEop","Check Tr Match E/p",200,0,2,200,0,2);
1672 fOutputList->Add(MatchTrEop);
1673
bfc7c23b 1674 PostData(1,fOutputList);
1675}
1676
1677//________________________________________________________________________
1678void AliAnalysisTaskHFECal::Terminate(Option_t *)
1679{
1680 // Info("Terminate");
1681 AliAnalysisTaskSE::Terminate();
1682}
1683
1684//________________________________________________________________________
1685Bool_t AliAnalysisTaskHFECal::ProcessCutStep(Int_t cutStep, AliVParticle *track)
1686{
1687 // Check single track cuts for a given cut step
1688 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
1689 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
1690 return kTRUE;
1691}
1692//_________________________________________
0f566acb 1693void 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 1694{
1695 //Identify non-heavy flavour electrons using Invariant mass method
1696
fb87d707 1697 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
feffe705 1698 fTrackCuts->SetRequireTPCRefit(kTRUE);
85985bb0 1699 fTrackCuts->SetRequireITSRefit(kTRUE);
feffe705 1700 fTrackCuts->SetEtaRange(-0.9,0.9);
bfc7c23b 1701 fTrackCuts->SetMaxChi2PerClusterTPC(3.5);
feffe705 1702 fTrackCuts->SetMinNClustersTPC(90);
bfc7c23b 1703
c6bca3f2 1704 //const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
6451f693 1705 Double_t bfield = fESD->GetMagneticField();
2f8ef315 1706 Double_t emass = 0.51*0.001; // (0.51 MeV)
6451f693 1707
93c189c5 1708 double ptEle = track->Pt(); //add
1709 if(ibgevent==0 && w > 0.0)
1710 {
1711 fpTCheck->Fill(ptEle,w);
1712 }
1713
bfc7c23b 1714 Bool_t flagPhotonicElec = kFALSE;
f09766dd 1715 Bool_t flagConvinatElec = kFALSE;
bfc7c23b 1716
8806ce6c 1717 int p1 = 0;
1718 if(mce==3)
1719 {
1720 Int_t label = TMath::Abs(track->GetLabel());
1721 TParticle* particle = stack->Particle(label);
1722 p1 = particle->GetFirstMother();
1723 }
1724
36e9774f 1725 int numULS = 0;
1726 int numLS = 0;
1727
5e65985c 1728 for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){
bfc7c23b 1729 AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
1730 if (!trackAsso) {
1731 printf("ERROR: Could not receive track %d\n", jTracks);
1732 continue;
1733 }
5e65985c 1734 if(itrack==jTracks)continue;
d86a815c 1735 int jbgevent = 0;
5e65985c 1736
8806ce6c 1737 int p2 = 0;
1738 if(mce==3)
1739 {
1740 Int_t label2 = TMath::Abs(trackAsso->GetLabel());
1741 TParticle* particle2 = stack->Particle(label2);
d86a815c 1742 Bool_t MChijing_ass = fMC->IsFromBGEvent(label2);
1743 if(MChijing_ass)jbgevent =1;
8806ce6c 1744 if(particle2->GetFirstMother()>-1)
1745 p2 = particle2->GetFirstMother();
1746 }
1747
f09766dd 1748 Double_t dEdxAsso = -999., ptPrim=-999., ptAsso=-999., openingAngle = -999.;
bfc7c23b 1749 Double_t mass=999., width = -999;
1750 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
1751
93c189c5 1752 //ptPrim = track->Pt();
1753 ptPrim = ptEle;
f09766dd 1754
bfc7c23b 1755 dEdxAsso = trackAsso->GetTPCsignal();
1756 ptAsso = trackAsso->Pt();
1757 Int_t chargeAsso = trackAsso->Charge();
1758 Int_t charge = track->Charge();
1759
b12bc641 1760
6451f693 1761 if(ptAsso <fMimpTassCut) continue;
bfc7c23b 1762 if(!fTrackCuts->AcceptTrack(trackAsso)) continue;
6b2cee73 1763 double fTPCnSigmaAss = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(trackAsso, AliPID::kElectron) : 1000;
1764 //cout << "fTPCnSigmaAss = " << fTPCnSigmaAss << endl;
1765 //cout << "fTPCnSigmaAss Cut = " << fMimNsigassCut << endl;
1766 if(fTPCnSigmaAss <fMimNsigassCut || fTPCnSigmaAss>5) continue;
1767 //cout << "fTPCnSigmaAss a.f. cut = " << fTPCnSigmaAss << endl;
bfc7c23b 1768
1769 Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
1770 if(charge>0) fPDGe1 = -11;
1771 if(chargeAsso>0) fPDGe2 = -11;
b12bc641 1772
1773 //printf("chargeAsso = %d\n",chargeAsso);
1774 //printf("charge = %d\n",charge);
bfc7c23b 1775 if(charge == chargeAsso) fFlagLS = kTRUE;
1776 if(charge != chargeAsso) fFlagULS = kTRUE;
1777
b12bc641 1778 //printf("fFlagLS = %d\n",fFlagLS);
1779 //printf("fFlagULS = %d\n",fFlagULS);
cccab38f 1780 //printf("\n");
b12bc641 1781
6451f693 1782 AliKFParticle::SetField(bfield);
bfc7c23b 1783 AliKFParticle ge1(*track, fPDGe1);
1784 AliKFParticle ge2(*trackAsso, fPDGe2);
2f8ef315 1785
1786 if(fSetMassNonlinear)
1787 {
1788 ge1.SetNonlinearMassConstraint(emass);
1789 ge2.SetNonlinearMassConstraint(emass);
1790 }
1791
bfc7c23b 1792 AliKFParticle recg(ge1, ge2);
1793
2f8ef315 1794 /*
d19af75d 1795 // vertex
bfc7c23b 1796 AliKFVertex primV(*pVtx);
1797 primV += recg;
8079a103 1798 primV -= ge1;
1799 primV -= ge2;
bfc7c23b 1800 recg.SetProductionVertex(primV);
2f8ef315 1801 */
1802
d19af75d 1803 // check chi2
1804 if(recg.GetNDF()<1) continue;
a009053a 1805 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
1806 Double_t chi2cut = 3.0;
f9447b3b 1807
a009053a 1808 // mass.
1809 if(fSetMassConstraint)
1810 {
1811 recg.SetMassConstraint(0,0.0001);
1812 chi2cut = 30.0;
1813 }
f9447b3b 1814 recg.GetMass(mass,width);
1815
3c475a9d 1816 if(fSetMassWidthCut && width>1e10)continue;
1817
a009053a 1818 // angle
bfc7c23b 1819 openingAngle = ge1.GetAngle(ge2);
1820 if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
1821 if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
0f566acb 1822
c2801a73 1823 double ishower = 0;
1824 if(shower>0.0 && shower<0.3)ishower = 1;
1825
0f566acb 1826 if(iCal==1)
1827 {
1828 double phoinfo[9];
1829 phoinfo[0] = cent;
1830 phoinfo[1] = ptPrim;
1831 phoinfo[2] = mass;
1832 phoinfo[3] = nSig;
1833 phoinfo[4] = openingAngle;
1834 phoinfo[5] = ishower;
1835 phoinfo[6] = ep;
1836 phoinfo[7] = mce;
1837 phoinfo[8] = ptAsso;
1838
1839 if(fFlagLS) fInvmassLS->Fill(phoinfo);
1840 if(fFlagULS) fInvmassULS->Fill(phoinfo);
1841 if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc->Fill(phoinfo,w);
1842 if(fFlagULS && ibgevent==0 && jbgevent==0)
1843 {
1844 fInvmassULSmc->Fill(phoinfo,w);
1845 }
1846 //printf("fInvmassCut %f\n",fInvmassCut);
1847 //printf("openingAngle %f\n",fOpeningAngleCut);
1848 }
b12bc641 1849
f9447b3b 1850 // angle cut
ffa14dda 1851 if(openingAngle > fOpeningAngleCut) continue;
a009053a 1852 // chi2 cut
6ec70f8c 1853 //if(TMath::Sqrt(TMath::Abs(chi2recg))>chi2cut) continue;
1854 if(chi2recg>chi2cut) continue;
1855
1856 if(fFlagLS ) fInvmassLSreco->Fill(ptPrim,mass);
1857 if(fFlagULS) fInvmassULSreco->Fill(ptPrim,mass);
1858
36e9774f 1859 // check double count
1860 if(mass<fInvmassCut && fFlagULS)numULS++;
1861 if(mass<fInvmassCut && fFlagLS)numLS++;
1862
d86a815c 1863 // for real data
93c189c5 1864 //printf("mce =%f\n",mce);
0f566acb 1865
d86a815c 1866 if(mce<-0.5) // mce==-1. is real
0f566acb 1867 {
1868 //printf("Real data\n");
1869
1870 // count
1871 if(mass<fInvmassCut && fFlagULS && shower>0.0 && shower<0.3 && iCal==1)fPhoElecPtM20Mass->Fill(cent,ptPrim);
1872 if(mass<fInvmassCut && fFlagLS && shower>0.0 && shower<0.3 && iCal==1)fSameElecPtM20Mass->Fill(cent,ptPrim);
1873
1874 // flag
1875 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
1876 flagPhotonicElec = kTRUE;
1877 }
1878 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec){
1879 flagConvinatElec = kTRUE;
1880 }
1881 }
d86a815c 1882 // for MC data
1883 else
0f566acb 1884 {
1885 // count
1886 if(mass<fInvmassCut && fFlagULS && shower>0.0 && shower<0.3 && mce>2.1 && iCal==1)fPhoElecPtMCM20Mass->Fill(cent,ptPrim,w);
1887 if(mass<fInvmassCut && fFlagLS && shower>0.0 && shower<0.3 && mce>2.1 && iCal==1)fSameElecPtMCM20Mass->Fill(cent,ptPrim,w);
1888
1889 //printf("MC data\n");
1890
1891 if(w>0.0)
1892 {
1893 //cout << "tagpi0 = " << tagpi0 << " ; tageta = " << tageta << endl;
1894 if(iCal==1)
1895 {
1896 if(fFlagLS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassLSmc0->Fill(ptPrim,mass);
1897 if(fFlagULS && ibgevent==0 && jbgevent==0 && tagpi0) fInvmassULSmc0->Fill(ptPrim,mass);
1898 if(fFlagLS && ibgevent==0 && jbgevent==0 && tageta) fInvmassLSmc1->Fill(ptPrim,mass);
1899 if(fFlagULS && ibgevent==0 && jbgevent==0 && tageta) fInvmassULSmc1->Fill(ptPrim,mass);
1900 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassLSmc2->Fill(ptPrim,mass);
1901 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tagpi0) fInvmassULSmc2->Fill(ptPrim,mass);
1902 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassLSmc3->Fill(ptPrim,mass);
1903 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2) && tageta) fInvmassULSmc3->Fill(ptPrim,mass);
1904 }
1905 }
1906
1907 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (ibgevent==jbgevent)){
1908 //if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (p1==p2)){ //<--- only MC train (55,56) v5-03-68-AN , 69 & v5-05-70-AN (till 74)
1909 flagPhotonicElec = kTRUE;
1910 }
1911 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec && (ibgevent==jbgevent)){
1912 flagConvinatElec = kTRUE;
1913 }
1914 }
d86a815c 1915
36e9774f 1916 } // end of associate loop
1917
6e1d7a8a 1918 if(numULS>0 || numLS>0)
1919 {
1920 int numPair = numULS-numLS;
0f566acb 1921 if(iCal==1)fpair->Fill(ptEle,numPair);
6e1d7a8a 1922 }
1923
0f566acb 1924 fFlagPhotonic = flagPhotonicElec;
1925 fFlagConvinat = flagConvinatElec;
6451f693 1926
1927}
1928
8806ce6c 1929//-------------------------------------------
46305ed6 1930
1931void AliAnalysisTaskHFECal::FindMother(TParticle* part, int &label, int &pid)
1932{
1933 //int label = 99999;
1934 //int pid = 99999;
1935
1936 if(part->GetFirstMother()>-1)
1937 {
1938 label = part->GetFirstMother();
1939 pid = stack->Particle(label)->GetPdgCode();
1940 }
70448e3b 1941 //cout << "Find Mother : label = " << label << " ; pid" << pid << endl;
46305ed6 1942}
1943
8806ce6c 1944double AliAnalysisTaskHFECal::GetMCweight(double mcPi0pT)
1945{
1946 double weight = 1.0;
1947
1948 if(mcPi0pT>0.0 && mcPi0pT<5.0)
1949 {
1950 weight = 0.323*mcPi0pT/(TMath::Exp(-1.6+0.767*mcPi0pT+0.0285*mcPi0pT*mcPi0pT));
1951 }
1952 else
1953 {
1954 weight = 115.0/(0.718*mcPi0pT*TMath::Power(mcPi0pT,3.65));
1955 }
1956 return weight;
1957}
fb87d707 1958
93c189c5 1959double AliAnalysisTaskHFECal::GetMCweightEta(double mcEtapT)
1960{
1961 double weight = 1.0;
1962
1963 weight = 223.3/TMath::Power((TMath::Exp(-0.17*mcEtapT-0.0322*mcEtapT*mcEtapT)+mcEtapT/1.69),5.65);
1964 return weight;
1965}
1966
1967
fb87d707 1968//_________________________________________
1969void AliAnalysisTaskHFECal::FindTriggerClusters()
1970{
bd6ee6dd 1971 //cout << "finding trigger patch" << endl;
fb87d707 1972 // constants
1973 const int nModuleCols = 2;
1974 const int nModuleRows = 5;
1975 const int nColsFeeModule = 48;
1976 const int nRowsFeeModule = 24;
1977 const int nColsFaltroModule = 24;
1978 const int nRowsFaltroModule = 12;
1979 //const int faltroWidthMax = 20;
1980
1981 // part 1, trigger extraction -------------------------------------
1982 Int_t globCol, globRow;
1983 //Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0, trigInCut=0;
1984 Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0;
1985
1986 //Int_t trigtimes[faltroWidthMax];
1987 Double_t cellTime[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1988 Double_t cellEnergy[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1989 //Double_t fTrigCutLow = 6;
1990 //Double_t fTrigCutHigh = 10;
1991 Double_t fTimeCutLow = 469e-09;
1992 Double_t fTimeCutHigh = 715e-09;
1993
1994 AliESDCaloTrigger * fCaloTrigger = fESD->GetCaloTrigger( "EMCAL" );
1995
1996 // erase trigger maps
1997 for(Int_t i = 0; i < nColsFaltroModule*nModuleCols; i++ )
1998 {
1999 for(Int_t j = 0; j < nRowsFaltroModule*nModuleRows; j++ )
2000 {
2001 ftriggersCut[i][j] = 0;
2002 ftriggers[i][j] = 0;
2003 ftriggersTime[i][j] = 0;
2004 }
2005 }
2006
2007 Int_t iglobCol=0, iglobRow=0;
2008 // go through triggers
2009 if( fCaloTrigger->GetEntries() > 0 )
2010 {
2011 // needs reset
2012 fCaloTrigger->Reset();
2013 while( fCaloTrigger->Next() )
2014 {
2015 fCaloTrigger->GetPosition( globCol, globRow );
2016 fCaloTrigger->GetNL0Times( ntimes );
2017 /*
2018 // no L0s
2019 if( ntimes < 1 ) continue;
2020 // get precise timings
2021 fCaloTrigger->GetL0Times( trigtimes );
2022 trigInCut = 0;
2023 for(Int_t i = 0; i < ntimes; i++ )
2024 {
2025 // save the first trigger time in channel
2026 if( i == 0 || triggersTime[globCol][globRow] > trigtimes[i] )
2027 triggersTime[globCol][globRow] = trigtimes[i];
2028 //printf("trigger times: %d\n",trigtimes[i]);
2029 // check if in cut
2030 if(trigtimes[i] > fTrigCutLow && trigtimes[i] < fTrigCutHigh )
2031 trigInCut = 1;
2032
2033 fTrigTimes->Fill(trigtimes[i]);
2034 }
2035 */
2036
2037 //L1 analysis from AliAnalysisTaskEMCALTriggerQA
2038 Int_t bit = 0;
2039 fCaloTrigger->GetTriggerBits(bit);
bd6ee6dd 2040 //cout << "bit = " << bit << endl;
fb87d707 2041
2042 Int_t ts = 0;
2043 fCaloTrigger->GetL1TimeSum(ts);
bd6ee6dd 2044 //cout << "ts = " << ts << endl;
fb87d707 2045 if (ts > 0)ftriggers[globCol][globRow] = 1;
2046 // number of triggered channels in event
2047 nTrigChannel++;
2048 // ... inside cut
2049 if(ts>0 && (bit >> 6 & 0x1))
2050 {
2051 iglobCol = globCol;
2052 iglobRow = globRow;
2053 nTrigChannelCut++;
bd6ee6dd 2054 //cout << "ts cut = " << ts << endl;
2055 //cout << "globCol = " << globCol << endl;
2056 //cout << "globRow = " << globRow << endl;
fb87d707 2057 ftriggersCut[globCol][globRow] = 1;
2058 }
2059
2060 } // calo trigger entries
2061 } // has calo trigger entries
2062
2063 // part 2 go through the clusters here -----------------------------------
bd6ee6dd 2064 //cout << " part 2 go through the clusters here ----------------------------------- " << endl;
fb87d707 2065 Int_t nCluster=0, nCell=0, iCell=0, gCell=0;
60d77596 2066 Short_t cellAddr, nSACell;
2067 Int_t mclabel;
fb87d707 2068 //Int_t nSACell, iSACell, mclabel;
2069 Int_t iSACell;
42c75692 2070 Double_t cellAmp=0, cellTimeT=0, clusterTime=0, efrac=0;
fb87d707 2071 Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta, gphi, geta, feta, fphi;
2072
2073 TRefArray *fCaloClusters = new TRefArray();
2074 fESD->GetEMCALClusters( fCaloClusters );
2075 nCluster = fCaloClusters->GetEntries();
2076
2077
2078 // save all cells times if there are clusters
2079 if( nCluster > 0 ){
2080 // erase time map
2081 for(Int_t i = 0; i < nColsFeeModule*nModuleCols; i++ ){
2082 for(Int_t j = 0; j < nRowsFeeModule*nModuleRows; j++ ){
2083 cellTime[i][j] = 0.;
2084 cellEnergy[i][j] = 0.;
2085 }
2086 }
2087
2088 // get cells
2089 AliESDCaloCells *fCaloCells = fESD->GetEMCALCells();
2090 //AliVCaloCells fCaloCells = fESD->GetEMCALCells();
2091 nSACell = fCaloCells->GetNumberOfCells();
2092 for(iSACell = 0; iSACell < nSACell; iSACell++ ){
2093 // get the cell info *fCal
2094 fCaloCells->GetCell( iSACell, cellAddr, cellAmp, cellTimeT , mclabel, efrac);
2095
2096 // get cell position
2097 fGeom->GetCellIndex( cellAddr, nSupMod, nModule, nIphi, nIeta );
2098 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
2099
2100 // convert co global phi eta
2101 gphi = iphi + nRowsFeeModule*(nSupMod/2);
2102 geta = ieta + nColsFeeModule*(nSupMod%2);
2103
2104 // save cell time and energy
2105 cellTime[geta][gphi] = cellTimeT;
2106 cellEnergy[geta][gphi] = cellAmp;
2107
2108 }
2109 }
2110
2111 Int_t nClusterTrig, nClusterTrigCut;
2112 UShort_t *cellAddrs;
2113 Double_t clsE=-999, clsEta=-999, clsPhi=-999;
2114 Float_t clsPos[3] = {0.,0.,0.};
2115
2116 for(Int_t icl=0; icl<fESD->GetNumberOfCaloClusters(); icl++)
2117 {
2118 AliESDCaloCluster *cluster = fESD->GetCaloCluster(icl);
2119 if(!cluster || !cluster->IsEMCAL()) continue;
2120
2121 // get cluster cells
2122 nCell = cluster->GetNCells();
2123
2124 // get cluster energy
2125 clsE = cluster->E();
2126
2127 // get cluster position
2128 cluster->GetPosition(clsPos);
2129 TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
2130 clsEta = clsPosVec.Eta();
2131 clsPhi = clsPosVec.Phi();
2132
2133 // get the cell addresses
2134 cellAddrs = cluster->GetCellsAbsId();
2135
2136 // check if the cluster contains cell, that was marked as triggered
2137 nClusterTrig = 0;
2138 nClusterTrigCut = 0;
2139
2140 // loop the cells to check, if cluser in acceptance
2141 // any cluster with a cell outside acceptance is not considered
2142 for( iCell = 0; iCell < nCell; iCell++ )
2143 {
42c75692 2144 // check hot cell
feffe705 2145 //if(clsE>6.0)fCellCheck->Fill(clsE,cellAddrs[iCell]);
42c75692 2146
fb87d707 2147 // get cell position
2148 fGeom->GetCellIndex( cellAddrs[iCell], nSupMod, nModule, nIphi, nIeta );
2149 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
2150
2151 // convert co global phi eta
2152 gphi = iphi + nRowsFeeModule*(nSupMod/2);
2153 geta = ieta + nColsFeeModule*(nSupMod%2);
2154
2155 if( cellTime[geta][gphi] > 0. ){
2156 clusterTime += cellTime[geta][gphi];
2157 gCell++;
2158 }
2159
2160 // get corresponding FALTRO
2161 fphi = gphi / 2;
2162 feta = geta / 2;
2163
bd6ee6dd 2164 //cout << "fphi = " << fphi << endl;
2165 //cout << "feta = " << feta << endl;
2166
fb87d707 2167 // try to match with a triggered
2168 if( ftriggers[feta][fphi]==1)
2169 { nClusterTrig++;
2170 }
2171 if( ftriggersCut[feta][fphi]==1)
2172 { nClusterTrigCut++;
2173 }
2174
bd6ee6dd 2175 //cout << "nClusterTrigCut : " << nClusterTrigCut << endl;
2176
fb87d707 2177 } // cells
2178
2179
2180 if( gCell > 0 )
2181 clusterTime = clusterTime / (Double_t)gCell;
2182 // fix the reconstriction code time 100ns jumps
2183 if( fESD->GetBunchCrossNumber() % 4 < 2 )
2184 clusterTime -= 0.0000001;
2185
feffe705 2186 //fClsETime->Fill(clsE,clusterTime);
2187 //fClsEBftTrigCut->Fill(clsE);
fb87d707 2188
2189 if(nClusterTrig>0){
feffe705 2190 //fClsETime1->Fill(clsE,clusterTime);
fb87d707 2191 }
2192
2193 if(nClusterTrig>0){
2194 cluster->SetChi2(1);
feffe705 2195 //fClsEAftTrigCut1->Fill(clsE);
fb87d707 2196 }
2197
2198 if(nClusterTrigCut>0){
2199 cluster->SetChi2(2);
feffe705 2200 //fClsEAftTrigCut2->Fill(clsE);
fb87d707 2201 }
2202
2203 if(nClusterTrigCut>0 && ( nCell > (1 + clsE / 3)))
2204 {
2205 cluster->SetChi2(3);
feffe705 2206 //fClsEAftTrigCut3->Fill(clsE);
fb87d707 2207 }
2208
2209 if(nClusterTrigCut>0 && (nCell > (1 + clsE / 3) )&&( clusterTime > fTimeCutLow && clusterTime < fTimeCutHigh ))
2210 {
2211 // cluster->SetChi2(4);
feffe705 2212 //fClsEAftTrigCut4->Fill(clsE);
fb87d707 2213 }
2214 if(nClusterTrigCut<1)
2215 {
2216 cluster->SetChi2(0);
2217
feffe705 2218 //fClsEAftTrigCut->Fill(clsE);
fb87d707 2219 }
2220
2221 } // clusters
2222}
2223
09022877 2224// <-------- only MC correction
8efacc22 2225double AliAnalysisTaskHFECal::MCEopMeanCorrection(double pTmc, float central)
8af58b0d 2226{
8efacc22 2227 TF1 *fcorr0 = new TF1("fcorr0","[0]*tanh([1]+[2]*x)");
2228 TF1 *fcorr1 = new TF1("fcorr1","[0]*tanh([1]+[2]*x)");
2229
8af58b0d 2230 double shift = 0.0;
8efacc22 2231
2232 if(central>0 && central<=10)
2233 {
2234 fcorr0->SetParameters(1.045,1.288,3.18e-01); //
2235 fcorr1->SetParameters(9.91e-01,3.466,2.344);
2236 }
09022877 2237 else if(central>10 && central<=20)
8efacc22 2238 {
2239 fcorr0->SetParameters(1.029,8.254e-01,4.07e-01);
2240 fcorr1->SetParameters(0.975,2.276,1.501e-01);
2241 }
2242 else if(central>20 && central<=30)
2243 {
2244 fcorr0->SetParameters(1.01,8.795e-01,3.904e-01);
2245 fcorr1->SetParameters(9.675e-01,1.654,2.583e-01);
2246 }
2247 else if(central>30 && central<=40)
2248 {
2249 fcorr0->SetParameters(1.00,1.466,2.305e-1);
2250 fcorr1->SetParameters(9.637e-01,1.473,2.754e-01);
2251 }
2252 else if(central>40 && central<=50)
2253 {
fd2126aa 2254 fcorr0->SetParameters(1.00,1.422,1.518e-01);
8efacc22 2255 fcorr1->SetParameters(9.59e-01,1.421,2.931e-01);
2256 }
2257
2258 else if(central>50 && central<=70)
2259 {
2260 fcorr0->SetParameters(0.989,2.495,2.167);
e1f0fb74 2261 fcorr1->SetParameters(0.961,1.734,1.438e-01);
8efacc22 2262 }
2263 else if(central>70 && central<=100)
2264 {
2265 fcorr0->SetParameters(0.981,-3.373,3.93327);
fd2126aa 2266 fcorr1->SetParameters(9.574e-01,1.698,1.58e-01);
8efacc22 2267 }
2268
2269
2270 shift = fcorr0->Eval(pTmc)-fcorr1->Eval(pTmc);
2271
89f41a30 2272 fcorr0->Delete();
2273 fcorr1->Delete();
2274
8af58b0d 2275 return shift;
2276}
fb87d707 2277
09022877 2278// <-------- only Data correction
2279double AliAnalysisTaskHFECal::NsigmaCorrection(double tmpeta, float central)
2280{
2281 int icent = 0;
2282
2283 if(central>=0 && central<10)
2284 {
2285 icent = 0;
2286 }
2287 else if(central>=10 && central<20)
2288 {
2289 icent = 1;
2290 }
2291 else if(central>=20 && central<30)
2292 {
2293 icent = 2;
2294 }
2295 else if(central>=30 && central<40)
2296 {
2297 icent = 3;
2298 }
2299 else if(central>=40 && central<50)
2300 {
2301 icent = 4;
2302 }
2303 else if(central>=50 && central<70)
2304 {
2305 icent = 5;
2306 }
2307 else
2308 {
2309 icent = 6;
2310 }
2311
2312 double shift = fnSigEtaCorr[icent]->Eval(tmpeta);
2313
2314 //cout << "eta correction"<< endl;
2315 //cout << "cent = "<< central<< endl;
2316 //cout << "icent = "<< icent << endl;
2317 //cout << "shift = "<< shift << endl;
2318
2319 return shift;
fb87d707 2320
09022877 2321}
fb87d707 2322