]> 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
Add three particle correlation loop for correction of auto-correlaton bias
[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"
23#include "TMath.h"
24#include "TCanvas.h"
25#include "THnSparse.h"
26#include "TLorentzVector.h"
27#include "TString.h"
28#include "TFile.h"
29
30#include "AliAnalysisTask.h"
31#include "AliAnalysisManager.h"
32
33#include "AliESDEvent.h"
34#include "AliESDHandler.h"
35#include "AliAODEvent.h"
36#include "AliAODHandler.h"
37
e7d87aef 38#include "AliMCEventHandler.h"
39#include "AliMCEvent.h"
40#include "AliMCParticle.h"
41
bfc7c23b 42#include "AliAnalysisTaskHFECal.h"
43#include "TGeoGlobalMagField.h"
44#include "AliLog.h"
45#include "AliAnalysisTaskSE.h"
46#include "TRefArray.h"
47#include "TVector.h"
48#include "AliESDInputHandler.h"
49#include "AliESDpid.h"
50#include "AliESDtrackCuts.h"
51#include "AliPhysicsSelection.h"
52#include "AliESDCaloCluster.h"
53#include "AliAODCaloCluster.h"
54#include "AliEMCALRecoUtils.h"
55#include "AliEMCALGeometry.h"
56#include "AliGeomManager.h"
57#include "stdio.h"
58#include "TGeoManager.h"
59#include "iostream"
60#include "fstream"
61
62#include "AliEMCALTrack.h"
63#include "AliMagF.h"
64
65#include "AliKFParticle.h"
66#include "AliKFVertex.h"
67
68#include "AliPID.h"
69#include "AliPIDResponse.h"
70#include "AliHFEcontainer.h"
71#include "AliHFEcuts.h"
72#include "AliHFEpid.h"
73#include "AliHFEpidBase.h"
74#include "AliHFEpidQAmanager.h"
75#include "AliHFEtools.h"
76#include "AliCFContainer.h"
77#include "AliCFManager.h"
78
e7d87aef 79#include "AliStack.h"
80
bfc7c23b 81#include "AliCentrality.h"
82
83ClassImp(AliAnalysisTaskHFECal)
84//________________________________________________________________________
85AliAnalysisTaskHFECal::AliAnalysisTaskHFECal(const char *name)
86 : AliAnalysisTaskSE(name)
87 ,fESD(0)
e7d87aef 88 ,fMC(0)
fb87d707 89 ,fGeom(0)
bfc7c23b 90 ,fOutputList(0)
85985bb0 91 ,fqahist(1)
bfc7c23b 92 ,fTrackCuts(0)
93 ,fCuts(0)
94 ,fIdentifiedAsOutInz(kFALSE)
95 ,fPassTheEventCut(kFALSE)
96 ,fRejectKinkMother(kFALSE)
e7d87aef 97 ,fmcData(kFALSE)
bfc7c23b 98 ,fVz(0.0)
99 ,fCFM(0)
100 ,fPID(0)
101 ,fPIDqa(0)
102 ,fOpeningAngleCut(0.1)
103 ,fInvmassCut(0.01)
104 ,fNoEvents(0)
105 ,fEMCAccE(0)
106 ,fTrkpt(0)
107 ,fTrkEovPBef(0)
108 ,fTrkEovPAft(0)
109 ,fdEdxBef(0)
110 ,fdEdxAft(0)
111 ,fIncpT(0)
fb87d707 112 ,fIncpTM20(0)
bfc7c23b 113 ,fInvmassLS(0)
114 ,fInvmassULS(0)
8806ce6c 115 ,fInvmassLSmc(0)
116 ,fInvmassULSmc(0)
93c189c5 117 ,fInvmassLSmc0(0)
118 ,fInvmassLSmc1(0)
119 ,fInvmassLSmc2(0)
120 ,fInvmassLSmc3(0)
121 ,fInvmassULSmc0(0)
122 ,fInvmassULSmc1(0)
123 ,fInvmassULSmc2(0)
124 ,fInvmassULSmc3(0)
bfc7c23b 125 ,fOpeningAngleLS(0)
126 ,fOpeningAngleULS(0)
127 ,fPhotoElecPt(0)
128 ,fPhoElecPt(0)
fb87d707 129 ,fPhoElecPtM20(0)
f09766dd 130 ,fSameElecPt(0)
fb87d707 131 ,fSameElecPtM20(0)
bfc7c23b 132 ,fTrackPtBefTrkCuts(0)
133 ,fTrackPtAftTrkCuts(0)
134 ,fTPCnsigma(0)
135 ,fCent(0)
136 ,fEleInfo(0)
feffe705 137 /*
fb87d707 138 ,fClsEBftTrigCut(0)
139 ,fClsEAftTrigCut(0)
140 ,fClsEAftTrigCut1(0)
141 ,fClsEAftTrigCut2(0)
142 ,fClsEAftTrigCut3(0)
143 ,fClsEAftTrigCut4(0)
144 ,fClsETime(0)
145 ,fClsETime1(0)
146 ,fTrigTimes(0)
42c75692 147 ,fCellCheck(0)
feffe705 148 */
e7d87aef 149 ,fInputHFEMC(0)
feffe705 150 ,fInputAlle(0)
e7d87aef 151 ,fIncpTMChfe(0)
3db00c72 152 ,fIncpTMChfeAll(0)
e7d87aef 153 ,fIncpTMCM20hfe(0)
3db00c72 154 ,fIncpTMCM20hfeAll(0)
60544aea 155 ,fIncpTMCM20hfeCheck(0)
e7d87aef 156 ,fIncpTMCpho(0)
157 ,fIncpTMCM20pho(0)
158 ,fPhoElecPtMC(0)
159 ,fPhoElecPtMCM20(0)
160 ,fSameElecPtMC(0)
161 ,fSameElecPtMCM20(0)
e4b0faf2 162 ,fIncpTMCM20pho_pi0e(0)
163 ,fPhoElecPtMCM20_pi0e(0)
164 ,fSameElecPtMCM20_pi0e(0)
93c189c5 165 ,fIncpTMCM20pho_eta(0)
166 ,fPhoElecPtMCM20_eta(0)
167 ,fSameElecPtMCM20_eta(0)
a6df418c 168 ,CheckNclust(0)
169 ,CheckNits(0)
d113d7cd 170 ,Hpi0pTcheck(0)
19325033 171 ,HETApTcheck(0)
e4b0faf2 172 ,HphopTcheck(0)
93c189c5 173 ,fpTCheck(0)
50919258 174 ,fMomDtoE(0)
bfc7c23b 175{
176 //Named constructor
177
178 fPID = new AliHFEpid("hfePid");
179 fTrackCuts = new AliESDtrackCuts();
180
181 // Define input and output slots here
182 // Input slot #0 works with a TChain
183 DefineInput(0, TChain::Class());
184 // Output slot #0 id reserved by the base class for AOD
185 // Output slot #1 writes into a TH1 container
186 // DefineOutput(1, TH1I::Class());
187 DefineOutput(1, TList::Class());
188 // DefineOutput(3, TTree::Class());
189}
190
191//________________________________________________________________________
192AliAnalysisTaskHFECal::AliAnalysisTaskHFECal()
193 : AliAnalysisTaskSE("DefaultAnalysis_AliAnalysisTaskHFECal")
194 ,fESD(0)
e7d87aef 195 ,fMC(0)
fb87d707 196 ,fGeom(0)
bfc7c23b 197 ,fOutputList(0)
85985bb0 198 ,fqahist(1)
bfc7c23b 199 ,fTrackCuts(0)
200 ,fCuts(0)
201 ,fIdentifiedAsOutInz(kFALSE)
202 ,fPassTheEventCut(kFALSE)
203 ,fRejectKinkMother(kFALSE)
e7d87aef 204 ,fmcData(kFALSE)
bfc7c23b 205 ,fVz(0.0)
206 ,fCFM(0)
207 ,fPID(0)
208 ,fPIDqa(0)
209 ,fOpeningAngleCut(0.1)
210 ,fInvmassCut(0.01)
211 ,fNoEvents(0)
212 ,fEMCAccE(0)
213 ,fTrkpt(0)
214 ,fTrkEovPBef(0)
215 ,fTrkEovPAft(0)
216 ,fdEdxBef(0)
217 ,fdEdxAft(0)
218 ,fIncpT(0)
fb87d707 219 ,fIncpTM20(0)
bfc7c23b 220 ,fInvmassLS(0)
221 ,fInvmassULS(0)
8806ce6c 222 ,fInvmassLSmc(0)
223 ,fInvmassULSmc(0)
93c189c5 224 ,fInvmassLSmc0(0)
225 ,fInvmassLSmc1(0)
226 ,fInvmassLSmc2(0)
227 ,fInvmassLSmc3(0)
228 ,fInvmassULSmc0(0)
229 ,fInvmassULSmc1(0)
230 ,fInvmassULSmc2(0)
231 ,fInvmassULSmc3(0)
bfc7c23b 232 ,fOpeningAngleLS(0)
233 ,fOpeningAngleULS(0)
234 ,fPhotoElecPt(0)
235 ,fPhoElecPt(0)
fb87d707 236 ,fPhoElecPtM20(0)
f09766dd 237 ,fSameElecPt(0)
fb87d707 238 ,fSameElecPtM20(0)
bfc7c23b 239 ,fTrackPtBefTrkCuts(0)
240 ,fTrackPtAftTrkCuts(0)
241 ,fTPCnsigma(0)
242 ,fCent(0)
243 ,fEleInfo(0)
feffe705 244 /*
fb87d707 245 ,fClsEBftTrigCut(0)
246 ,fClsEAftTrigCut(0)
247 ,fClsEAftTrigCut1(0)
248 ,fClsEAftTrigCut2(0)
249 ,fClsEAftTrigCut3(0)
250 ,fClsEAftTrigCut4(0)
251 ,fClsETime(0)
252 ,fClsETime1(0)
253 ,fTrigTimes(0)
42c75692 254 ,fCellCheck(0)
feffe705 255 */
e7d87aef 256 ,fInputHFEMC(0)
feffe705 257 ,fInputAlle(0)
e7d87aef 258 ,fIncpTMChfe(0)
3db00c72 259 ,fIncpTMChfeAll(0)
e7d87aef 260 ,fIncpTMCM20hfe(0)
3db00c72 261 ,fIncpTMCM20hfeAll(0)
60544aea 262 ,fIncpTMCM20hfeCheck(0)
e7d87aef 263 ,fIncpTMCpho(0)
264 ,fIncpTMCM20pho(0)
265 ,fPhoElecPtMC(0)
266 ,fPhoElecPtMCM20(0)
267 ,fSameElecPtMC(0)
268 ,fSameElecPtMCM20(0)
e4b0faf2 269 ,fIncpTMCM20pho_pi0e(0)
270 ,fPhoElecPtMCM20_pi0e(0)
271 ,fSameElecPtMCM20_pi0e(0)
93c189c5 272 ,fIncpTMCM20pho_eta(0)
273 ,fPhoElecPtMCM20_eta(0)
274 ,fSameElecPtMCM20_eta(0)
a6df418c 275 ,CheckNclust(0)
276 ,CheckNits(0)
d113d7cd 277 ,Hpi0pTcheck(0)
19325033 278 ,HETApTcheck(0)
e4b0faf2 279 ,HphopTcheck(0)
93c189c5 280 ,fpTCheck(0)
50919258 281 ,fMomDtoE(0)
bfc7c23b 282{
283 //Default constructor
284 fPID = new AliHFEpid("hfePid");
285
286 fTrackCuts = new AliESDtrackCuts();
287
288 // Constructor
289 // Define input and output slots here
290 // Input slot #0 works with a TChain
291 DefineInput(0, TChain::Class());
292 // Output slot #0 id reserved by the base class for AOD
293 // Output slot #1 writes into a TH1 container
294 // DefineOutput(1, TH1I::Class());
295 DefineOutput(1, TList::Class());
296 //DefineOutput(3, TTree::Class());
297}
298//_________________________________________
299
300AliAnalysisTaskHFECal::~AliAnalysisTaskHFECal()
301{
302 //Destructor
303
304 delete fOutputList;
fb87d707 305 delete fGeom;
bfc7c23b 306 delete fPID;
a6df418c 307 delete fCuts;
bfc7c23b 308 delete fCFM;
309 delete fPIDqa;
310 delete fTrackCuts;
311}
312//_________________________________________
313
314void AliAnalysisTaskHFECal::UserExec(Option_t*)
315{
316 //Main loop
317 //Called for each event
318
319 // create pointer to event
320 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
321 if (!fESD) {
322 printf("ERROR: fESD not available\n");
323 return;
324 }
325
326 if(!fCuts){
327 AliError("HFE cuts not available");
328 return;
329 }
330
331 if(!fPID->IsInitialized()){
332 // Initialize PID with the given run number
333 AliWarning("PID not initialised, get from Run no");
334 fPID->InitializePID(fESD->GetRunNumber());
335 }
336
e7d87aef 337 if(fmcData)fMC = MCEvent();
338 AliStack* stack = NULL;
339 if(fmcData && fMC)stack = fMC->Stack();
340
341 Float_t cent = -1.;
342 AliCentrality *centrality = fESD->GetCentrality();
343 cent = centrality->GetCentralityPercentile("V0M");
344
345 //---- fill MC track info
346 if(fmcData && fMC)
347 {
348 Int_t nParticles = stack->GetNtrack();
349 for (Int_t iParticle = 0; iParticle < nParticles; iParticle++) {
350 TParticle* particle = stack->Particle(iParticle);
351 int fPDG = particle->GetPdgCode();
352 double mcZvertex = fMC->GetPrimaryVertex()->GetZ();
353 double pTMC = particle->Pt();
feffe705 354 double proR = particle->R();
9b3495ff 355 double etaMC = particle->Eta();
356 if(fabs(etaMC)>0.6)continue;
e7d87aef 357 Bool_t mcInDtoE= kFALSE;
358 Bool_t mcInBtoE= kFALSE;
359
d113d7cd 360 Bool_t MChijing = fMC->IsFromBGEvent(iParticle);
fd426741 361 //if(!MChijing)printf("not MC hijing");
d113d7cd 362 int iHijing = 1;
363 if(!MChijing)iHijing = 0;
364 if(fPDG==111)Hpi0pTcheck->Fill(pTMC,iHijing);
19325033 365 if(fPDG==221)HETApTcheck->Fill(pTMC,iHijing);
d113d7cd 366
e4b0faf2 367 if(particle->GetFirstMother()>-1 && fPDG==22)
368 {
369 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
370 if(parentPID==111 || parentPID==221)HphopTcheck->Fill(pTMC,iHijing); // pi0->g & eta->g
371 }
372
e7d87aef 373 if(particle->GetFirstMother()>-1 && fabs(fPDG)==11)
374 {
375 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
376 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(fPDG)==11)mcInDtoE = kTRUE;
377 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(fPDG)==11)mcInBtoE = kTRUE;
378 if((mcInBtoE || mcInDtoE) && fabs(mcZvertex)<10.0)fInputHFEMC->Fill(cent,pTMC);
379 }
380
feffe705 381 if(proR<7.0 && fabs(fPDG)==11)fInputAlle->Fill(cent,pTMC);
382
e7d87aef 383 }
384 }
385
bfc7c23b 386 fNoEvents->Fill(0);
387
388 Int_t fNOtrks = fESD->GetNumberOfTracks();
389 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
390
391 Double_t pVtxZ = -999;
392 pVtxZ = pVtx->GetZ();
393
394 if(TMath::Abs(pVtxZ)>10) return;
395 fNoEvents->Fill(1);
396
397 if(fNOtrks<2) return;
398 fNoEvents->Fill(2);
399
400 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
401 if(!pidResponse){
402 AliDebug(1, "Using default PID Response");
403 pidResponse = AliHFEtools::GetDefaultPID(kFALSE, fInputEvent->IsA() == AliAODEvent::Class());
404 }
405
406 fPID->SetPIDResponse(pidResponse);
407
408 fCFM->SetRecEventInfo(fESD);
409
e7d87aef 410 //Float_t cent = -1.;
411 //AliCentrality *centrality = fESD->GetCentrality();
412 //cent = centrality->GetCentralityPercentile("V0M");
bfc7c23b 413 fCent->Fill(cent);
414
b12bc641 415 //if(cent>90.) return;
bfc7c23b 416
417 // Calorimeter info.
418
42c75692 419 FindTriggerClusters();
420
bfc7c23b 421 // make EMCAL array
422 for(Int_t iCluster=0; iCluster<fESD->GetNumberOfCaloClusters(); iCluster++)
423 {
424 AliESDCaloCluster *clust = fESD->GetCaloCluster(iCluster);
425 if(clust->IsEMCAL())
426 {
427 double clustE = clust->E();
428 float emcx[3]; // cluster pos
429 clust->GetPosition(emcx);
430 TVector3 clustpos(emcx[0],emcx[1],emcx[2]);
431 double emcphi = clustpos.Phi();
432 double emceta = clustpos.Eta();
fb87d707 433 double calInfo[5];
434 calInfo[0] = emcphi; calInfo[1] = emceta; calInfo[2] = clustE; calInfo[3] = cent; calInfo[4] = clust->Chi2();
85985bb0 435 //if(clustE>1.5)fEMCAccE->Fill(calInfo);
44be9e1d 436 //if(fqahist==1 && clustE>1.5)fEMCAccE->Fill(calInfo);
bfc7c23b 437 }
438 }
439
440 // Track loop
441 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
442 AliESDtrack* track = fESD->GetTrack(iTracks);
443 if (!track) {
444 printf("ERROR: Could not receive track %d\n", iTracks);
445 continue;
446 }
e7d87aef 447
448 Bool_t mcPho = kFALSE;
449 Bool_t mcDtoE= kFALSE;
450 Bool_t mcBtoE= kFALSE;
93c189c5 451 Bool_t mcOrgPi0 = kFALSE;
452 Bool_t mcOrgEta = kFALSE;
e7d87aef 453 double mcele = -1.;
44be9e1d 454 double mcpT = 0.0;
b12bc641 455 double mcMompT = 0.0;
e4b0faf2 456 //double mcGrandMompT = 0.0;
457 double mcWeight = -10.0;
d113d7cd 458
459 int iHijing = 1;
460
e7d87aef 461 if(fmcData && fMC && stack)
462 {
feffe705 463 Int_t label = TMath::Abs(track->GetLabel());
d113d7cd 464
465 Bool_t MChijing = fMC->IsFromBGEvent(label);
466 if(!MChijing)iHijing = 0;
467
e7d87aef 468 TParticle* particle = stack->Particle(label);
469 int mcpid = particle->GetPdgCode();
44be9e1d 470 mcpT = particle->Pt();
9b3495ff 471 //printf("MCpid = %d",mcpid);
e7d87aef 472 if(particle->GetFirstMother()>-1)
473 {
e4b0faf2 474 int parentlabel = particle->GetFirstMother();
e7d87aef 475 int parentPID = stack->Particle(particle->GetFirstMother())->GetPdgCode();
b12bc641 476 mcMompT = stack->Particle(particle->GetFirstMother())->Pt();
e7d87aef 477 if((parentPID==22 || parentPID==111 || parentPID==221)&& fabs(mcpid)==11)mcPho = kTRUE;
478 if((fabs(parentPID)==411 || fabs(parentPID)==413 || fabs(parentPID)==421 || fabs(parentPID)==423 || fabs(parentPID)==431)&& fabs(mcpid)==11)mcDtoE = kTRUE;
479 if((fabs(parentPID)==511 || fabs(parentPID)==513 || fabs(parentPID)==521 || fabs(parentPID)==523 || fabs(parentPID)==531)&& fabs(mcpid)==11)mcBtoE = kTRUE;
e4b0faf2 480
50919258 481 // make D->e pT correlation
482 if(mcDtoE)fMomDtoE->Fill(mcpT,mcMompT);
483
484 // pi->e (Dalitz)
485 if(parentPID==111 && fabs(mcpid)==11 && mcMompT>0.0)
e4b0faf2 486 {
93c189c5 487 mcOrgPi0 = kTRUE;
8806ce6c 488 mcWeight = GetMCweight(mcMompT);
e4b0faf2 489 }
93c189c5 490 // eta->e (Dalitz)
491 if(parentPID==221 && fabs(mcpid)==11 && mcMompT>0.0)
492 {
493 mcOrgEta = kTRUE;
494 mcWeight = GetMCweightEta(mcMompT);
495 }
496
e4b0faf2 497
498 // access grand parent pi0->g->e
499 TParticle* particle_parent = stack->Particle(parentlabel); // get parent pointer
500 if(particle_parent->GetFirstMother()>-1 && parentPID==22 && fabs(mcpid)==11) // get grand parent g->e
501 {
502 int grand_parentPID = stack->Particle(particle_parent->GetFirstMother())->GetPdgCode();
93c189c5 503 double pTtmp = stack->Particle(particle_parent->GetFirstMother())->Pt();
504 if(grand_parentPID==111 && pTtmp>0.0)
e4b0faf2 505 {
93c189c5 506 mcOrgPi0 = kTRUE;
11704995 507 mcWeight = GetMCweight(pTtmp);
93c189c5 508 }
509 if(grand_parentPID==221 && pTtmp>0.0)
510 {
511 mcOrgEta = kTRUE;
512 mcWeight = GetMCweightEta(pTtmp);
e4b0faf2 513 }
514 }
e7d87aef 515 }
9b3495ff 516 if(fabs(mcpid)==11 && mcDtoE)mcele= 1.;
517 if(fabs(mcpid)==11 && mcBtoE)mcele= 2.;
518 if(fabs(mcpid)==11 && mcPho)mcele= 3.;
e7d87aef 519 }
520
50919258 521 //printf("weight = %f\n",mcWeight);
522
9b3495ff 523 if(TMath::Abs(track->Eta())>0.6) continue;
f09766dd 524 if(TMath::Abs(track->Pt()<2.0)) continue;
bfc7c23b 525
526 fTrackPtBefTrkCuts->Fill(track->Pt());
527 // RecKine: ITSTPC cuts
528 if(!ProcessCutStep(AliHFEcuts::kStepRecKineITSTPC, track)) continue;
529
530 //RecKink
531 if(fRejectKinkMother) { // Quick and dirty fix to reject both kink mothers and daughters
532 if(track->GetKinkIndex(0) != 0) continue;
533 }
534
535 // RecPrim
536 if(!ProcessCutStep(AliHFEcuts::kStepRecPrim, track)) continue;
537
538 // HFEcuts: ITS layers cuts
539 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsITS, track)) continue;
540
f09766dd 541 // HFE cuts: TPC PID cleanup
542 if(!ProcessCutStep(AliHFEcuts::kStepHFEcutsTPC, track)) continue;
543
a6df418c 544 int nTPCcl = track->GetTPCNcls();
6e92c428 545 int nTPCclF = track->GetTPCNclsF();
a6df418c 546 int nITS = track->GetNcls(0);
bfc7c23b 547
548 fTrackPtAftTrkCuts->Fill(track->Pt());
549
550 Double_t mom = -999., eop=-999., pt = -999., dEdx=-999., fTPCnSigma=-10, phi=-999., eta=-999.;
feffe705 551 pt = track->Pt();
552 if(pt<2.0)continue;
bfc7c23b 553
554 // Track extrapolation
555
bfc7c23b 556 Int_t charge = track->Charge();
557 fTrkpt->Fill(pt);
558 mom = track->P();
559 phi = track->Phi();
560 eta = track->Eta();
561 dEdx = track->GetTPCsignal();
562 fTPCnSigma = fPID->GetPIDResponse() ? fPID->GetPIDResponse()->NumberOfSigmasTPC(track, AliPID::kElectron) : 1000;
563
564 double ncells = -1.0;
565 double m20 = -1.0;
566 double m02 = -1.0;
567 double disp = -1.0;
568 double rmatch = -1.0;
569 double nmatch = -1.0;
f09766dd 570 double oppstatus = 0.0;
bfc7c23b 571
f09766dd 572 Bool_t fFlagPhotonicElec = kFALSE;
573 Bool_t fFlagConvinatElec = kFALSE;
bfc7c23b 574
575 Int_t clsId = track->GetEMCALcluster();
576 if (clsId>0){
577 AliESDCaloCluster *clust = fESD->GetCaloCluster(clsId);
578 if(clust && clust->IsEMCAL()){
579
580 double clustE = clust->E();
581 eop = clustE/fabs(mom);
582 //double clustT = clust->GetTOF();
583 ncells = clust->GetNCells();
584 m02 = clust->GetM02();
585 m20 = clust->GetM20();
586 disp = clust->GetDispersion();
587 double delphi = clust->GetTrackDx();
588 double deleta = clust->GetTrackDz();
589 rmatch = sqrt(pow(delphi,2)+pow(deleta,2));
590 nmatch = clust->GetNTracksMatched();
591
c2801a73 592 if(fTPCnSigma>-1.5 && fTPCnSigma<3.0)
593 {
8806ce6c 594 SelectPhotonicElectron(iTracks,cent,track,stack,fFlagPhotonicElec,fFlagConvinatElec,fTPCnSigma,m20,eop,mcele,mcWeight,iHijing);
c2801a73 595 }
596 if(fFlagPhotonicElec)oppstatus = 1.0;
597 if(fFlagConvinatElec)oppstatus = 2.0;
598 if(fFlagPhotonicElec && fFlagConvinatElec)oppstatus = 3.0;
599
feffe705 600 double valdedx[16];
a6df418c 601 valdedx[0] = pt; valdedx[1] = nITS; valdedx[2] = phi; valdedx[3] = eta; valdedx[4] = fTPCnSigma;
6e92c428 602 valdedx[5] = eop; valdedx[6] = rmatch; valdedx[7] = ncells, valdedx[8] = nTPCclF; valdedx[9] = m20; valdedx[10] = mcpT;
a6df418c 603 valdedx[11] = cent; valdedx[12] = charge; valdedx[13] = oppstatus; valdedx[14] = nTPCcl;
feffe705 604 valdedx[15] = mcele;
85985bb0 605 if(fqahist==1)fEleInfo->Fill(valdedx);
bfc7c23b 606
607
608 }
609 }
a6df418c 610
611 if(nITS<2.5)continue;
612 if(nTPCcl<100)continue;
613
614 CheckNclust->Fill(nTPCcl);
615 CheckNits->Fill(nITS);
616
42c75692 617 fdEdxBef->Fill(mom,fTPCnSigma);
bfc7c23b 618 fTPCnsigma->Fill(mom,fTPCnSigma);
f09766dd 619 if(fTPCnSigma >= -1.0 && fTPCnSigma <= 3)fTrkEovPBef->Fill(pt,eop);
bfc7c23b 620
bfc7c23b 621 Int_t pidpassed = 1;
622
623
624 //--- track accepted
625 AliHFEpidObject hfetrack;
626 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
627 hfetrack.SetRecTrack(track);
3db00c72 628 double binct = 10.5;
629 if((0.0< cent) && (cent<5.0)) binct = 0.5;
630 if((5.0< cent) && (cent<10.0)) binct = 1.5;
631 if((10.0< cent) && (cent<20.0)) binct = 2.5;
632 if((20.0< cent) && (cent<30.0)) binct = 3.5;
633 if((30.0< cent) && (cent<40.0)) binct = 4.5;
634 if((40.0< cent) && (cent<50.0)) binct = 5.5;
635 if((50.0< cent) && (cent<60.0)) binct = 6.5;
636 if((60.0< cent) && (cent<70.0)) binct = 7.5;
637 if((70.0< cent) && (cent<80.0)) binct = 8.5;
638 if((80.0< cent) && (cent<90.0)) binct = 9.5;
639 if((90.0< cent) && (cent<100.0)) binct = 10.5;
640
641 hfetrack.SetCentrality((int)binct); //added
bfc7c23b 642 hfetrack.SetPbPb();
643 if(!fPID->IsSelected(&hfetrack, NULL, "", fPIDqa)) pidpassed = 0;
644
645 if(pidpassed==0) continue;
646
647 fTrkEovPAft->Fill(pt,eop);
42c75692 648 fdEdxAft->Fill(mom,fTPCnSigma);
bfc7c23b 649
e7d87aef 650 // Fill real data
fb87d707 651 fIncpT->Fill(cent,pt);
bfc7c23b 652 if(fFlagPhotonicElec) fPhoElecPt->Fill(cent,pt);
f09766dd 653 if(fFlagConvinatElec) fSameElecPt->Fill(cent,pt);
fb87d707 654
655 if(m20>0.0 && m20<0.3)
656 {
657 fIncpTM20->Fill(cent,pt);
658 if(fFlagPhotonicElec) fPhoElecPtM20->Fill(cent,pt);
659 if(fFlagConvinatElec) fSameElecPtM20->Fill(cent,pt);
660 }
feffe705 661
e7d87aef 662 // MC
e4b0faf2 663 if(mcele>0) // select MC electrons
e7d87aef 664 {
3db00c72 665
666 fIncpTMChfeAll->Fill(cent,pt);
667 if(m20>0.0 && m20<0.3)fIncpTMCM20hfeAll->Fill(cent,pt);
668
e4b0faf2 669 if(mcBtoE || mcDtoE) // select B->e & D->e
e7d87aef 670 {
671 fIncpTMChfe->Fill(cent,pt);
672 if(m20>0.0 && m20<0.3)fIncpTMCM20hfe->Fill(cent,pt);
60544aea 673 if(m20>0.0 && m20<0.3)fIncpTMCM20hfeCheck->Fill(cent,mcpT);
e7d87aef 674 }
e4b0faf2 675
676 if(mcPho) // select photonic electrons
e7d87aef 677 {
987053ce 678 double phoval[5];
44be9e1d 679 phoval[0] = cent;
680 phoval[1] = pt;
681 phoval[2] = fTPCnSigma;
d113d7cd 682 phoval[3] = iHijing;
987053ce 683 phoval[4] = mcMompT;
44be9e1d 684
685 fIncpTMCpho->Fill(phoval);
686 if(fFlagPhotonicElec) fPhoElecPtMC->Fill(phoval);
687 if(fFlagConvinatElec) fSameElecPtMC->Fill(phoval);
e7d87aef 688
689 if(m20>0.0 && m20<0.3)
690 {
44be9e1d 691 fIncpTMCM20pho->Fill(phoval);
692 if(fFlagPhotonicElec) fPhoElecPtMCM20->Fill(phoval);
693 if(fFlagConvinatElec) fSameElecPtMCM20->Fill(phoval);
e4b0faf2 694 // pi0->g->e
695 if(mcWeight>-1)
696 {
987053ce 697 if(iHijing==1)mcWeight = 1.0;
93c189c5 698 if(mcOrgPi0)
699 {
700 fIncpTMCM20pho_pi0e->Fill(phoval,mcWeight);
701 if(fFlagPhotonicElec) fPhoElecPtMCM20_pi0e->Fill(phoval,mcWeight);
702 if(fFlagConvinatElec) fSameElecPtMCM20_pi0e->Fill(phoval,mcWeight);
703 }
704 if(mcOrgEta)
705 {
706 fIncpTMCM20pho_eta->Fill(phoval,mcWeight);
707 if(fFlagPhotonicElec) fPhoElecPtMCM20_eta->Fill(phoval,mcWeight);
708 if(fFlagConvinatElec) fSameElecPtMCM20_eta->Fill(phoval,mcWeight);
709 }
710 // --- eta
e4b0faf2 711 }
e7d87aef 712 }
713 }
714 }
bfc7c23b 715 }
716 PostData(1, fOutputList);
717}
718//_________________________________________
719void AliAnalysisTaskHFECal::UserCreateOutputObjects()
720{
42c75692 721 //--- Check MC
722
e7d87aef 723 //Bool_t mcData = kFALSE;
42c75692 724 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())
725 {
e7d87aef 726 fmcData = kTRUE;
42c75692 727 printf("+++++ MC Data available");
728 }
e7d87aef 729 if(fmcData)
42c75692 730 {
731 printf("++++++++= MC analysis \n");
732 }
733 else
734 {
735 printf("++++++++ real data analysis \n");
736 }
737
85985bb0 738 printf("+++++++ QA hist %d",fqahist);
739
fb87d707 740 //---- Geometry
741 fGeom = AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1");
742
bfc7c23b 743 //--------Initialize PID
42c75692 744 //fPID->SetHasMCData(kFALSE);
e7d87aef 745 fPID->SetHasMCData(fmcData);
bfc7c23b 746 if(!fPID->GetNumberOfPIDdetectors())
747 {
748 fPID->AddDetector("TPC", 0);
749 fPID->AddDetector("EMCAL", 1);
750 }
751
3db00c72 752 Double_t params[4];
78ff1095 753 const char *cutmodel;
3db00c72 754 cutmodel = "pol0";
755 params[0] = -1.0; //sigma min
9b3495ff 756 double maxnSig = 3.0;
757 if(fmcData)
758 {
759 params[0] = -5.0; //sigma min
760 maxnSig = 5.0;
761 }
762
763 for(Int_t a=0;a<11;a++)fPID->ConfigureTPCcentralityCut(a,cutmodel,params,maxnSig);
3db00c72 764
bfc7c23b 765 fPID->SortDetectors();
766 fPIDqa = new AliHFEpidQAmanager();
767 fPIDqa->Initialize(fPID);
a6df418c 768
769 //------- fcut --------------
770 fCuts = new AliHFEcuts();
771 fCuts->CreateStandardCuts();
772 //fCuts->SetMinNClustersTPC(100);
773 fCuts->SetMinNClustersTPC(90);
774 fCuts->SetMinRatioTPCclusters(0.6);
775 fCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
776 //fCuts->SetMinNClustersITS(3);
777 fCuts->SetMinNClustersITS(2);
778 fCuts->SetCutITSpixel(AliHFEextraCuts::kAny);
779 fCuts->SetCheckITSLayerStatus(kFALSE);
780 fCuts->SetVertexRange(10.);
781 fCuts->SetTOFPIDStep(kFALSE);
782 fCuts->SetPtRange(2, 50);
783 fCuts->SetMaxImpactParam(3.,3.);
784
bfc7c23b 785 //--------Initialize correction Framework and Cuts
786 fCFM = new AliCFManager;
787 const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
788 fCFM->SetNStepParticle(kNcutSteps);
789 for(Int_t istep = 0; istep < kNcutSteps; istep++)
790 fCFM->SetParticleCutsList(istep, NULL);
791
792 if(!fCuts){
793 AliWarning("Cuts not available. Default cuts will be used");
794 fCuts = new AliHFEcuts;
795 fCuts->CreateStandardCuts();
796 }
797 fCuts->Initialize(fCFM);
798
799 //---------Output Tlist
800 fOutputList = new TList();
801 fOutputList->SetOwner();
802 fOutputList->Add(fPIDqa->MakeList("PIDQA"));
803
804 fNoEvents = new TH1F("fNoEvents","",4,-0.5,3.5) ;
805 fOutputList->Add(fNoEvents);
806
02cd82a0 807 Int_t binsE[5] = {250, 100, 1000, 200, 10};
fb87d707 808 Double_t xminE[5] = {1.0, -1, 0.0, 0, -0.5};
809 Double_t xmaxE[5] = {3.5, 1, 100.0, 100, 9.5};
810 fEMCAccE = new THnSparseD("fEMCAccE","EMC acceptance & E;#eta;#phi;Energy;Centrality;trugCondition;",5,binsE,xminE,xmaxE);
85985bb0 811 if(fqahist==1)fOutputList->Add(fEMCAccE);
bfc7c23b 812
813 fTrkpt = new TH1F("fTrkpt","track pt",100,0,50);
814 fOutputList->Add(fTrkpt);
815
816 fTrackPtBefTrkCuts = new TH1F("fTrackPtBefTrkCuts","track pt before track cuts",100,0,50);
817 fOutputList->Add(fTrackPtBefTrkCuts);
818
819 fTrackPtAftTrkCuts = new TH1F("fTrackPtAftTrkCuts","track pt after track cuts",100,0,50);
820 fOutputList->Add(fTrackPtAftTrkCuts);
821
822 fTPCnsigma = new TH2F("fTPCnsigma", "TPC - n sigma",100,0,50,200,-10,10);
823 fOutputList->Add(fTPCnsigma);
824
825 fTrkEovPBef = new TH2F("fTrkEovPBef","track E/p before HFE pid",100,0,50,100,0,2);
826 fOutputList->Add(fTrkEovPBef);
827
828 fTrkEovPAft = new TH2F("fTrkEovPAft","track E/p after HFE pid",100,0,50,100,0,2);
829 fOutputList->Add(fTrkEovPAft);
830
42c75692 831 fdEdxBef = new TH2F("fdEdxBef","track dEdx vs p before HFE pid",100,0,50,200,-10,10);
bfc7c23b 832 fOutputList->Add(fdEdxBef);
833
42c75692 834 fdEdxAft = new TH2F("fdEdxAft","track dEdx vs p after HFE pid",100,0,50,200,-10,10);
bfc7c23b 835 fOutputList->Add(fdEdxAft);
836
02cd82a0 837 fIncpT = new TH2F("fIncpT","HFE pid electro vs. centrality",200,0,100,100,0,50);
bfc7c23b 838 fOutputList->Add(fIncpT);
839
02cd82a0 840 fIncpTM20 = new TH2F("fIncpTM20","HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
fb87d707 841 fOutputList->Add(fIncpTM20);
4e01b68c 842
02cd82a0 843 Int_t nBinspho[9] = { 200, 100, 500, 12, 50, 4, 200, 8, 100};
b12bc641 844 Double_t minpho[9] = { 0., 0., 0., -2.5, 0, -0.5, 0,-1.5, 0};
845 Double_t maxpho[9] = {100., 50., 0.5, 3.5, 1, 3.5, 2, 6.5, 50};
f09766dd 846
b12bc641 847 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);
bfc7c23b 848 fOutputList->Add(fInvmassLS);
849
b12bc641 850 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);
bfc7c23b 851 fOutputList->Add(fInvmassULS);
852
8806ce6c 853 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);
854 fOutputList->Add(fInvmassLSmc);
855
856 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);
857 fOutputList->Add(fInvmassULSmc);
858
93c189c5 859 fInvmassLSmc0 = new TH2D("fInvmassLSmc0", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
860 fInvmassLSmc0->Sumw2();
861 fOutputList->Add(fInvmassLSmc0);
862
863 fInvmassLSmc1 = new TH2D("fInvmassLSmc1", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
864 fInvmassLSmc1->Sumw2();
865 fOutputList->Add(fInvmassLSmc1);
866
867 fInvmassLSmc2 = new TH2D("fInvmassLSmc2", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
868 fInvmassLSmc2->Sumw2();
869 fOutputList->Add(fInvmassLSmc2);
870
871 fInvmassLSmc3 = new TH2D("fInvmassLSmc3", "Inv mass of LS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
872 fInvmassLSmc3->Sumw2();
873 fOutputList->Add(fInvmassLSmc3);
874
875 fInvmassULSmc0 = new TH2D("fInvmassULSmc0", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
876 fInvmassULSmc0->Sumw2();
877 fOutputList->Add(fInvmassULSmc0);
878
879 fInvmassULSmc1 = new TH2D("fInvmassULSmc1", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
880 fInvmassULSmc1->Sumw2();
881 fOutputList->Add(fInvmassULSmc1);
882
883 fInvmassULSmc2 = new TH2D("fInvmassULSmc2", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
884 fInvmassULSmc2->Sumw2();
885 fOutputList->Add(fInvmassULSmc2);
886
887 fInvmassULSmc3 = new TH2D("fInvmassULSmc3", "Inv mass of ULS (e,e); cent; p_{T} (GeV/c); mass(GeV/c^2)",20,0,20,500,0,0.5 );
888 fInvmassULSmc3->Sumw2();
889 fOutputList->Add(fInvmassULSmc2);
890 fOutputList->Add(fInvmassULSmc3);
891
bfc7c23b 892 fOpeningAngleLS = new TH1F("fOpeningAngleLS","Opening angle for LS pairs",100,0,1);
893 fOutputList->Add(fOpeningAngleLS);
894
895 fOpeningAngleULS = new TH1F("fOpeningAngleULS","Opening angle for ULS pairs",100,0,1);
896 fOutputList->Add(fOpeningAngleULS);
897
898 fPhotoElecPt = new TH1F("fPhotoElecPt", "photonic electron pt",100,0,50);
899 fOutputList->Add(fPhotoElecPt);
900
02cd82a0 901 fPhoElecPt = new TH2F("fPhoElecPt", "Pho-inclusive electron pt",200,0,100,100,0,50);
bfc7c23b 902 fOutputList->Add(fPhoElecPt);
903
02cd82a0 904 fPhoElecPtM20 = new TH2F("fPhoElecPtM20", "Pho-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 905 fOutputList->Add(fPhoElecPtM20);
906
02cd82a0 907 fSameElecPt = new TH2F("fSameElecPt", "Same-inclusive electron pt",200,0,100,100,0,50);
f09766dd 908 fOutputList->Add(fSameElecPt);
909
02cd82a0 910 fSameElecPtM20 = new TH2F("fSameElecPtM20", "Same-inclusive electron pt with M20",200,0,100,100,0,50);
fb87d707 911 fOutputList->Add(fSameElecPtM20);
912
02cd82a0 913 fCent = new TH1F("fCent","Centrality",200,0,100) ;
bfc7c23b 914 fOutputList->Add(fCent);
915
916 // Make common binning
44be9e1d 917 const Double_t kMinP = 0.;
bfc7c23b 918 const Double_t kMaxP = 50.;
a6df418c 919 //const Double_t kTPCSigMim = 40.;
920 //const Double_t kTPCSigMax = 140.;
bfc7c23b 921
922 // 1st histogram: TPC dEdx with/without EMCAL (p, pT, TPC Signal, phi, eta, Sig, e/p, ,match, cell, M02, M20, Disp, Centrality, select)
a6df418c 923 Int_t nBins[16] = { 250, 10, 60, 20, 100, 300, 50, 40, 200, 200, 250, 200, 3, 5, 100, 8};
924 Double_t min[16] = {kMinP, -0.5, 1.0, -1.0, -6.0, 0, 0, 0, 0.0, 0.0, 0.0, 0, -1.5, -0.5, 80, -1.5};
6e92c428 925 Double_t max[16] = {kMaxP, 9.5, 4.0, 1.0, 4.0, 3.0, 0.1, 40, 200, 2.0, 50.0, 100, 1.5, 4.5, 180, 6.5};
926 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);
85985bb0 927 if(fqahist==1)fOutputList->Add(fEleInfo);
bfc7c23b 928
fb87d707 929 //<--- Trigger info
feffe705 930 /*
fb87d707 931 fClsEBftTrigCut = new TH1F("fClsEBftTrigCut","cluster E before trigger selection",1000,0,100);
932 fOutputList->Add(fClsEBftTrigCut);
933
934 fClsEAftTrigCut = new TH1F("fClsEAftTrigCut","cluster E if cls has 0 trigcut channel",1000,0,100);
935 fOutputList->Add(fClsEAftTrigCut);
936
937 fClsEAftTrigCut1 = new TH1F("fClsEAftTrigCut1","cluster E if cls with trig channel",1000,0,100);
938 fOutputList->Add(fClsEAftTrigCut1);
939
940 fClsEAftTrigCut2 = new TH1F("fClsEAftTrigCut2","cluster E if cls with trigcut channel",1000,0,100);
941 fOutputList->Add(fClsEAftTrigCut2);
942
943 fClsEAftTrigCut3 = new TH1F("fClsEAftTrigCut3","cluster E if cls with trigcut channel + nCell>Ecorrect",1000,0,100);
944 fOutputList->Add(fClsEAftTrigCut3);
945
946 fClsEAftTrigCut4 = new TH1F("fClsEAftTrigCut4","cluster E if cls with trigcut channel + nCell>Ecorrect + cls time cut",1000,0,100);
947 fOutputList->Add(fClsEAftTrigCut4);
948
949 fClsETime = new TH2F("fClsETime", "Cls time vs E; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
950 fOutputList->Add(fClsETime);
951
952 fClsETime1 = new TH2F("fClsETime1", "Cls time vs E if cls contains trigger channel; E; time;",1000,0,100,1000,-0.0000002,0.0000009);
953 fOutputList->Add(fClsETime1);
954
955 fTrigTimes = new TH1F("fTrigTimes", "Trigger time; time; N;",25,0,25);
956 fOutputList->Add(fTrigTimes);
957
42c75692 958 fCellCheck = new TH2F("fCellCheck", "Cell vs E; E GeV; Cell ID",10,6,26,12000,0,12000);
959 fOutputList->Add(fCellCheck);
feffe705 960 */
e7d87aef 961 //<---------- MC
962
02cd82a0 963 fInputHFEMC = new TH2F("fInputHFEMC","Input MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 964 fOutputList->Add(fInputHFEMC);
965
02cd82a0 966 fInputAlle = new TH2F("fInputAlle","Input MC electro vs. centrality",200,0,100,100,0,50);
feffe705 967 fOutputList->Add(fInputAlle);
968
02cd82a0 969 fIncpTMChfe = new TH2F("fIncpTMChfe","MC HFE pid electro vs. centrality",200,0,100,100,0,50);
e7d87aef 970 fOutputList->Add(fIncpTMChfe);
971
02cd82a0 972 fIncpTMChfeAll = new TH2F("fIncpTMChfeAll","MC Alle pid electro vs. centrality",200,0,100,100,0,50);
3db00c72 973 fOutputList->Add(fIncpTMChfeAll);
974
02cd82a0 975 fIncpTMCM20hfe = new TH2F("fIncpTMCM20hfe","MC HFE pid electro vs. centrality with M20",200,0,100,100,0,50);
e7d87aef 976 fOutputList->Add(fIncpTMCM20hfe);
977
02cd82a0 978 fIncpTMCM20hfeAll = new TH2F("fIncpTMCM20hfeAll","MC Alle pid electro vs. centrality with M20",200,0,100,100,0,50);
3db00c72 979 fOutputList->Add(fIncpTMCM20hfeAll);
980
60544aea 981 fIncpTMCM20hfeCheck = new TH2F("fIncpTMCM20hfeCheck","MC HFE pid electro vs. centrality with M20 Check",200,0,100,100,0,50);
982 fOutputList->Add(fIncpTMCM20hfeCheck);
44be9e1d 983
987053ce 984 Int_t nBinspho2[5] = { 200, 100, 7, 3, 100};
985 Double_t minpho2[5] = { 0., 0., -2.5, -0.5, 0.};
986 Double_t maxpho2[5] = {100., 50., 4.5, 2.5, 50.};
44be9e1d 987
987053ce 988 fIncpTMCpho = new THnSparseD("fIncpTMCpho","MC Pho pid electro vs. centrality",5,nBinspho2,minpho2,maxpho2);
e7d87aef 989 fOutputList->Add(fIncpTMCpho);
990
987053ce 991 fIncpTMCM20pho = new THnSparseD("fIncpTMCM20pho","MC Pho pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 992 fOutputList->Add(fIncpTMCM20pho);
993
987053ce 994 fPhoElecPtMC = new THnSparseD("fPhoElecPtMC", "MC Pho-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 995 fOutputList->Add(fPhoElecPtMC);
996
987053ce 997 fPhoElecPtMCM20 = new THnSparseD("fPhoElecPtMCM20", "MC Pho-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
feffe705 998 fOutputList->Add(fPhoElecPtMCM20);
e7d87aef 999
987053ce 1000 fSameElecPtMC = new THnSparseD("fSameElecPtMC", "MC Same-inclusive electron pt",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1001 fOutputList->Add(fSameElecPtMC);
1002
987053ce 1003 fSameElecPtMCM20 = new THnSparseD("fSameElecPtMCM20", "MC Same-inclusive electron pt with M20",5,nBinspho2,minpho2,maxpho2);
e7d87aef 1004 fOutputList->Add(fSameElecPtMCM20);
1005
987053ce 1006 fIncpTMCM20pho_pi0e = new THnSparseD("fIncpTMCM20pho_pi0e","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1007 fIncpTMCM20pho_pi0e->Sumw2();
e4b0faf2 1008 fOutputList->Add(fIncpTMCM20pho_pi0e);
1009
987053ce 1010 fPhoElecPtMCM20_pi0e = new THnSparseD("fPhoElecPtMCM20_pi0e", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1011 fPhoElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1012 fOutputList->Add(fPhoElecPtMCM20_pi0e);
1013
987053ce 1014 fSameElecPtMCM20_pi0e = new THnSparseD("fSameElecPtMCM20_pi0e", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1015 fSameElecPtMCM20_pi0e->Sumw2();
e4b0faf2 1016 fOutputList->Add(fSameElecPtMCM20_pi0e);
93c189c5 1017 //
1018 fIncpTMCM20pho_eta = new THnSparseD("fIncpTMCM20pho_eta","MC Pho pi0->e pid electro vs. centrality with M20",5,nBinspho2,minpho2,maxpho2);
1019 fIncpTMCM20pho_eta->Sumw2();
1020 fOutputList->Add(fIncpTMCM20pho_eta);
1021
1022 fPhoElecPtMCM20_eta = new THnSparseD("fPhoElecPtMCM20_eta", "MC Pho-inclusive electron pt with M20 pi0->e",5,nBinspho2,minpho2,maxpho2);
1023 fPhoElecPtMCM20_eta->Sumw2();
1024 fOutputList->Add(fPhoElecPtMCM20_eta);
1025
1026 fSameElecPtMCM20_eta = new THnSparseD("fSameElecPtMCM20_eta", "MC Same-inclusive electron pt pi0->e",5,nBinspho2,minpho2,maxpho2);
1027 fSameElecPtMCM20_eta->Sumw2();
1028 fOutputList->Add(fSameElecPtMCM20_eta);
1029
e4b0faf2 1030
a6df418c 1031 CheckNclust = new TH1D("CheckNclust","cluster check",200,0,200);
1032 fOutputList->Add(CheckNclust);
1033
1034 CheckNits = new TH1D("CheckNits","ITS cluster check",8,-0.5,7.5);
1035 fOutputList->Add(CheckNits);
e7d87aef 1036
d113d7cd 1037 Hpi0pTcheck = new TH2D("Hpi0pTcheck","Pi0 pT from Hijing",100,0,50,3,-0.5,2.5);
1038 fOutputList->Add(Hpi0pTcheck);
1039
19325033 1040 HETApTcheck = new TH2D("HETApTcheck","Eta pT from Hijing",100,0,50,3,-0.5,2.5);
1041 fOutputList->Add(HETApTcheck);
1042
e4b0faf2 1043 HphopTcheck = new TH2D("HphopTcheck","Pho pT from Hijing",100,0,50,3,-0.5,2.5);
1044 fOutputList->Add(HphopTcheck);
1045
93c189c5 1046 fpTCheck = new TH1D("fpTCheck","pT check",500,0,50);
1047 fOutputList->Add(fpTCheck);
1048
50919258 1049 fMomDtoE = new TH2D("fMomDtoE","D->E pT correlations;e p_{T} GeV/c;D p_{T} GeV/c",400,0,40,400,0,40);
1050 fOutputList->Add(fMomDtoE);
d113d7cd 1051
bfc7c23b 1052 PostData(1,fOutputList);
1053}
1054
1055//________________________________________________________________________
1056void AliAnalysisTaskHFECal::Terminate(Option_t *)
1057{
1058 // Info("Terminate");
1059 AliAnalysisTaskSE::Terminate();
1060}
1061
1062//________________________________________________________________________
1063Bool_t AliAnalysisTaskHFECal::ProcessCutStep(Int_t cutStep, AliVParticle *track)
1064{
1065 // Check single track cuts for a given cut step
1066 const Int_t kMCOffset = AliHFEcuts::kNcutStepsMCTrack;
1067 if(!fCFM->CheckParticleCuts(cutStep + kMCOffset, track)) return kFALSE;
1068 return kTRUE;
1069}
1070//_________________________________________
f09766dd 1071//void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec)
c2801a73 1072//void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, Bool_t &fFlagPhotonicElec, Bool_t &fFlagConvinatElec, Double_t nSig)
8806ce6c 1073void AliAnalysisTaskHFECal::SelectPhotonicElectron(Int_t itrack, Double_t cent, AliESDtrack *track, AliStack *stack ,Bool_t &fFlagPhotonicElec, Bool_t &fFlagConvinatElec, Double_t nSig, Double_t shower, Double_t ep, Double_t mce, Double_t w, Int_t ibgevent)
bfc7c23b 1074{
1075 //Identify non-heavy flavour electrons using Invariant mass method
1076
fb87d707 1077 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
feffe705 1078 fTrackCuts->SetRequireTPCRefit(kTRUE);
85985bb0 1079 fTrackCuts->SetRequireITSRefit(kTRUE);
feffe705 1080 fTrackCuts->SetEtaRange(-0.9,0.9);
f09766dd 1081 //fTrackCuts->SetRequireSigmaToVertex(kTRUE);
bfc7c23b 1082 fTrackCuts->SetMaxChi2PerClusterTPC(3.5);
feffe705 1083 fTrackCuts->SetMinNClustersTPC(90);
bfc7c23b 1084
1085 const AliESDVertex *pVtx = fESD->GetPrimaryVertex();
1086
93c189c5 1087 double ptEle = track->Pt(); //add
1088 if(ibgevent==0 && w > 0.0)
1089 {
1090 fpTCheck->Fill(ptEle,w);
1091 }
1092
bfc7c23b 1093 Bool_t flagPhotonicElec = kFALSE;
f09766dd 1094 Bool_t flagConvinatElec = kFALSE;
bfc7c23b 1095
8806ce6c 1096 int p1 = 0;
1097 if(mce==3)
1098 {
1099 Int_t label = TMath::Abs(track->GetLabel());
1100 TParticle* particle = stack->Particle(label);
1101 p1 = particle->GetFirstMother();
1102 }
1103
5e65985c 1104 //for(Int_t jTracks = itrack+1; jTracks<fESD->GetNumberOfTracks(); jTracks++){
1105 for(Int_t jTracks = 0; jTracks<fESD->GetNumberOfTracks(); jTracks++){
bfc7c23b 1106 AliESDtrack* trackAsso = fESD->GetTrack(jTracks);
1107 if (!trackAsso) {
1108 printf("ERROR: Could not receive track %d\n", jTracks);
1109 continue;
1110 }
5e65985c 1111 if(itrack==jTracks)continue;
d86a815c 1112 int jbgevent = 0;
5e65985c 1113
8806ce6c 1114 int p2 = 0;
1115 if(mce==3)
1116 {
1117 Int_t label2 = TMath::Abs(trackAsso->GetLabel());
1118 TParticle* particle2 = stack->Particle(label2);
d86a815c 1119 Bool_t MChijing_ass = fMC->IsFromBGEvent(label2);
1120 if(MChijing_ass)jbgevent =1;
8806ce6c 1121 if(particle2->GetFirstMother()>-1)
1122 p2 = particle2->GetFirstMother();
1123 }
1124
f09766dd 1125 Double_t dEdxAsso = -999., ptPrim=-999., ptAsso=-999., openingAngle = -999.;
bfc7c23b 1126 Double_t mass=999., width = -999;
1127 Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
1128
93c189c5 1129 //ptPrim = track->Pt();
1130 ptPrim = ptEle;
f09766dd 1131
bfc7c23b 1132 dEdxAsso = trackAsso->GetTPCsignal();
1133 ptAsso = trackAsso->Pt();
1134 Int_t chargeAsso = trackAsso->Charge();
1135 Int_t charge = track->Charge();
1136
b12bc641 1137
3db00c72 1138 if(ptAsso <0.5) continue;
bfc7c23b 1139 if(!fTrackCuts->AcceptTrack(trackAsso)) continue;
f09766dd 1140 if(dEdxAsso <65 || dEdxAsso>100) continue; //11a pass1
bfc7c23b 1141
1142 Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
1143 if(charge>0) fPDGe1 = -11;
1144 if(chargeAsso>0) fPDGe2 = -11;
b12bc641 1145
1146 //printf("chargeAsso = %d\n",chargeAsso);
1147 //printf("charge = %d\n",charge);
bfc7c23b 1148 if(charge == chargeAsso) fFlagLS = kTRUE;
1149 if(charge != chargeAsso) fFlagULS = kTRUE;
1150
b12bc641 1151 //printf("fFlagLS = %d\n",fFlagLS);
1152 //printf("fFlagULS = %d\n",fFlagULS);
1153 //printf("\n");
1154
bfc7c23b 1155 AliKFParticle ge1(*track, fPDGe1);
1156 AliKFParticle ge2(*trackAsso, fPDGe2);
1157 AliKFParticle recg(ge1, ge2);
1158
1159 if(recg.GetNDF()<1) continue;
1160 Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
1161 if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
1162
d86a815c 1163 // for v5-03-70-AN comment
bfc7c23b 1164 AliKFVertex primV(*pVtx);
1165 primV += recg;
1166 recg.SetProductionVertex(primV);
1167
d113d7cd 1168 recg.SetMassConstraint(0,0.0001);
d86a815c 1169
bfc7c23b 1170 openingAngle = ge1.GetAngle(ge2);
1171 if(fFlagLS) fOpeningAngleLS->Fill(openingAngle);
1172 if(fFlagULS) fOpeningAngleULS->Fill(openingAngle);
1173
bfc7c23b 1174
1175 recg.GetMass(mass,width);
1176
c2801a73 1177 double ishower = 0;
1178 if(shower>0.0 && shower<0.3)ishower = 1;
1179
b12bc641 1180 double phoinfo[9];
f09766dd 1181 phoinfo[0] = cent;
1182 phoinfo[1] = ptPrim;
1183 phoinfo[2] = mass;
3db00c72 1184 phoinfo[3] = nSig;
4e01b68c 1185 phoinfo[4] = openingAngle;
c2801a73 1186 phoinfo[5] = ishower;
1187 phoinfo[6] = ep;
1188 phoinfo[7] = mce;
b12bc641 1189 phoinfo[8] = ptAsso;
f09766dd 1190
1191 if(fFlagLS) fInvmassLS->Fill(phoinfo);
1192 if(fFlagULS) fInvmassULS->Fill(phoinfo);
d86a815c 1193 if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc->Fill(phoinfo,w);
93c189c5 1194 if(fFlagULS && ibgevent==0 && jbgevent==0)
1195 {
1196 fInvmassULSmc->Fill(phoinfo,w);
1197 }
b12bc641 1198 //printf("fInvmassCut %f\n",fInvmassCut);
1199 //printf("openingAngle %f\n",fOpeningAngleCut);
1200
ffa14dda 1201 if(openingAngle > fOpeningAngleCut) continue;
bfc7c23b 1202
d86a815c 1203 // for real data
93c189c5 1204 //printf("mce =%f\n",mce);
d86a815c 1205 if(mce<-0.5) // mce==-1. is real
1206 {
93c189c5 1207 //printf("Real data\n");
d86a815c 1208 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec){
1209 //if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (p1==p2)){ <--- only MC train (55,56) v5-03-68-AN & 69 for check
1210 flagPhotonicElec = kTRUE;
1211 }
1212 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec){
1213 flagConvinatElec = kTRUE;
1214 }
1215 }
1216 // for MC data
1217 else
1218 {
93c189c5 1219 //printf("MC data\n");
1220
1221 if(w>0.0)
1222 {
1223 if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc0->Fill(ptPrim,mass,w);
1224 if(fFlagULS && ibgevent==0 && jbgevent==0) fInvmassULSmc0->Fill(ptPrim,mass,w);
1225 if(fFlagLS && ibgevent==0 && jbgevent==0) fInvmassLSmc1->Fill(ptPrim,mass);
1226 if(fFlagULS && ibgevent==0 && jbgevent==0) fInvmassULSmc1->Fill(ptPrim,mass);
1227 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2)) fInvmassLSmc2->Fill(ptPrim,mass,w);
1228 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2)) fInvmassULSmc2->Fill(ptPrim,mass,w);
1229 if(fFlagLS && ibgevent==0 && jbgevent==0 && (p1==p2)) fInvmassLSmc3->Fill(ptPrim,mass);
1230 if(fFlagULS && ibgevent==0 && jbgevent==0 && (p1==p2)) fInvmassULSmc3->Fill(ptPrim,mass);
1231 }
d86a815c 1232 if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (ibgevent==jbgevent)){
1233 //if(mass<fInvmassCut && fFlagULS && !flagPhotonicElec && (p1==p2)){ <--- only MC train (55,56) v5-03-68-AN & 69 for check
1234 flagPhotonicElec = kTRUE;
1235 }
1236 if(mass<fInvmassCut && fFlagLS && !flagConvinatElec && (ibgevent==jbgevent)){
1237 flagConvinatElec = kTRUE;
1238 }
1239 }
1240
bfc7c23b 1241 }
1242 fFlagPhotonicElec = flagPhotonicElec;
f09766dd 1243 fFlagConvinatElec = flagConvinatElec;
bfc7c23b 1244
1245}
8806ce6c 1246//-------------------------------------------
1247//void AliAnalysisTaskHFECal::GetMCweight(double mcPi0pT)
1248double AliAnalysisTaskHFECal::GetMCweight(double mcPi0pT)
1249{
1250 double weight = 1.0;
1251
1252 if(mcPi0pT>0.0 && mcPi0pT<5.0)
1253 {
1254 weight = 0.323*mcPi0pT/(TMath::Exp(-1.6+0.767*mcPi0pT+0.0285*mcPi0pT*mcPi0pT));
1255 }
1256 else
1257 {
1258 weight = 115.0/(0.718*mcPi0pT*TMath::Power(mcPi0pT,3.65));
1259 }
1260 return weight;
1261}
fb87d707 1262
93c189c5 1263double AliAnalysisTaskHFECal::GetMCweightEta(double mcEtapT)
1264{
1265 double weight = 1.0;
1266
1267 weight = 223.3/TMath::Power((TMath::Exp(-0.17*mcEtapT-0.0322*mcEtapT*mcEtapT)+mcEtapT/1.69),5.65);
1268 return weight;
1269}
1270
1271
fb87d707 1272//_________________________________________
1273void AliAnalysisTaskHFECal::FindTriggerClusters()
1274{
1275 // constants
1276 const int nModuleCols = 2;
1277 const int nModuleRows = 5;
1278 const int nColsFeeModule = 48;
1279 const int nRowsFeeModule = 24;
1280 const int nColsFaltroModule = 24;
1281 const int nRowsFaltroModule = 12;
1282 //const int faltroWidthMax = 20;
1283
1284 // part 1, trigger extraction -------------------------------------
1285 Int_t globCol, globRow;
1286 //Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0, trigInCut=0;
1287 Int_t ntimes=0, nTrigChannel=0, nTrigChannelCut=0;
1288
1289 //Int_t trigtimes[faltroWidthMax];
1290 Double_t cellTime[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1291 Double_t cellEnergy[nColsFeeModule*nModuleCols][nRowsFeeModule*nModuleRows];
1292 //Double_t fTrigCutLow = 6;
1293 //Double_t fTrigCutHigh = 10;
1294 Double_t fTimeCutLow = 469e-09;
1295 Double_t fTimeCutHigh = 715e-09;
1296
1297 AliESDCaloTrigger * fCaloTrigger = fESD->GetCaloTrigger( "EMCAL" );
1298
1299 // erase trigger maps
1300 for(Int_t i = 0; i < nColsFaltroModule*nModuleCols; i++ )
1301 {
1302 for(Int_t j = 0; j < nRowsFaltroModule*nModuleRows; j++ )
1303 {
1304 ftriggersCut[i][j] = 0;
1305 ftriggers[i][j] = 0;
1306 ftriggersTime[i][j] = 0;
1307 }
1308 }
1309
1310 Int_t iglobCol=0, iglobRow=0;
1311 // go through triggers
1312 if( fCaloTrigger->GetEntries() > 0 )
1313 {
1314 // needs reset
1315 fCaloTrigger->Reset();
1316 while( fCaloTrigger->Next() )
1317 {
1318 fCaloTrigger->GetPosition( globCol, globRow );
1319 fCaloTrigger->GetNL0Times( ntimes );
1320 /*
1321 // no L0s
1322 if( ntimes < 1 ) continue;
1323 // get precise timings
1324 fCaloTrigger->GetL0Times( trigtimes );
1325 trigInCut = 0;
1326 for(Int_t i = 0; i < ntimes; i++ )
1327 {
1328 // save the first trigger time in channel
1329 if( i == 0 || triggersTime[globCol][globRow] > trigtimes[i] )
1330 triggersTime[globCol][globRow] = trigtimes[i];
1331 //printf("trigger times: %d\n",trigtimes[i]);
1332 // check if in cut
1333 if(trigtimes[i] > fTrigCutLow && trigtimes[i] < fTrigCutHigh )
1334 trigInCut = 1;
1335
1336 fTrigTimes->Fill(trigtimes[i]);
1337 }
1338 */
1339
1340 //L1 analysis from AliAnalysisTaskEMCALTriggerQA
1341 Int_t bit = 0;
1342 fCaloTrigger->GetTriggerBits(bit);
1343
1344 Int_t ts = 0;
1345 fCaloTrigger->GetL1TimeSum(ts);
1346 if (ts > 0)ftriggers[globCol][globRow] = 1;
1347 // number of triggered channels in event
1348 nTrigChannel++;
1349 // ... inside cut
1350 if(ts>0 && (bit >> 6 & 0x1))
1351 {
1352 iglobCol = globCol;
1353 iglobRow = globRow;
1354 nTrigChannelCut++;
1355 ftriggersCut[globCol][globRow] = 1;
1356 }
1357
1358 } // calo trigger entries
1359 } // has calo trigger entries
1360
1361 // part 2 go through the clusters here -----------------------------------
1362 Int_t nCluster=0, nCell=0, iCell=0, gCell=0;
1363 Short_t cellAddr, nSACell, mclabel;
1364 //Int_t nSACell, iSACell, mclabel;
1365 Int_t iSACell;
42c75692 1366 Double_t cellAmp=0, cellTimeT=0, clusterTime=0, efrac=0;
fb87d707 1367 Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta, gphi, geta, feta, fphi;
1368
1369 TRefArray *fCaloClusters = new TRefArray();
1370 fESD->GetEMCALClusters( fCaloClusters );
1371 nCluster = fCaloClusters->GetEntries();
1372
1373
1374 // save all cells times if there are clusters
1375 if( nCluster > 0 ){
1376 // erase time map
1377 for(Int_t i = 0; i < nColsFeeModule*nModuleCols; i++ ){
1378 for(Int_t j = 0; j < nRowsFeeModule*nModuleRows; j++ ){
1379 cellTime[i][j] = 0.;
1380 cellEnergy[i][j] = 0.;
1381 }
1382 }
1383
1384 // get cells
1385 AliESDCaloCells *fCaloCells = fESD->GetEMCALCells();
1386 //AliVCaloCells fCaloCells = fESD->GetEMCALCells();
1387 nSACell = fCaloCells->GetNumberOfCells();
1388 for(iSACell = 0; iSACell < nSACell; iSACell++ ){
1389 // get the cell info *fCal
1390 fCaloCells->GetCell( iSACell, cellAddr, cellAmp, cellTimeT , mclabel, efrac);
1391
1392 // get cell position
1393 fGeom->GetCellIndex( cellAddr, nSupMod, nModule, nIphi, nIeta );
1394 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
1395
1396 // convert co global phi eta
1397 gphi = iphi + nRowsFeeModule*(nSupMod/2);
1398 geta = ieta + nColsFeeModule*(nSupMod%2);
1399
1400 // save cell time and energy
1401 cellTime[geta][gphi] = cellTimeT;
1402 cellEnergy[geta][gphi] = cellAmp;
1403
1404 }
1405 }
1406
1407 Int_t nClusterTrig, nClusterTrigCut;
1408 UShort_t *cellAddrs;
1409 Double_t clsE=-999, clsEta=-999, clsPhi=-999;
1410 Float_t clsPos[3] = {0.,0.,0.};
1411
1412 for(Int_t icl=0; icl<fESD->GetNumberOfCaloClusters(); icl++)
1413 {
1414 AliESDCaloCluster *cluster = fESD->GetCaloCluster(icl);
1415 if(!cluster || !cluster->IsEMCAL()) continue;
1416
1417 // get cluster cells
1418 nCell = cluster->GetNCells();
1419
1420 // get cluster energy
1421 clsE = cluster->E();
1422
1423 // get cluster position
1424 cluster->GetPosition(clsPos);
1425 TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
1426 clsEta = clsPosVec.Eta();
1427 clsPhi = clsPosVec.Phi();
1428
1429 // get the cell addresses
1430 cellAddrs = cluster->GetCellsAbsId();
1431
1432 // check if the cluster contains cell, that was marked as triggered
1433 nClusterTrig = 0;
1434 nClusterTrigCut = 0;
1435
1436 // loop the cells to check, if cluser in acceptance
1437 // any cluster with a cell outside acceptance is not considered
1438 for( iCell = 0; iCell < nCell; iCell++ )
1439 {
42c75692 1440 // check hot cell
feffe705 1441 //if(clsE>6.0)fCellCheck->Fill(clsE,cellAddrs[iCell]);
42c75692 1442
fb87d707 1443 // get cell position
1444 fGeom->GetCellIndex( cellAddrs[iCell], nSupMod, nModule, nIphi, nIeta );
1445 fGeom->GetCellPhiEtaIndexInSModule( nSupMod,nModule, nIphi, nIeta, iphi, ieta);
1446
1447 // convert co global phi eta
1448 gphi = iphi + nRowsFeeModule*(nSupMod/2);
1449 geta = ieta + nColsFeeModule*(nSupMod%2);
1450
1451 if( cellTime[geta][gphi] > 0. ){
1452 clusterTime += cellTime[geta][gphi];
1453 gCell++;
1454 }
1455
1456 // get corresponding FALTRO
1457 fphi = gphi / 2;
1458 feta = geta / 2;
1459
1460 // try to match with a triggered
1461 if( ftriggers[feta][fphi]==1)
1462 { nClusterTrig++;
1463 }
1464 if( ftriggersCut[feta][fphi]==1)
1465 { nClusterTrigCut++;
1466 }
1467
1468 } // cells
1469
1470
1471 if( gCell > 0 )
1472 clusterTime = clusterTime / (Double_t)gCell;
1473 // fix the reconstriction code time 100ns jumps
1474 if( fESD->GetBunchCrossNumber() % 4 < 2 )
1475 clusterTime -= 0.0000001;
1476
feffe705 1477 //fClsETime->Fill(clsE,clusterTime);
1478 //fClsEBftTrigCut->Fill(clsE);
fb87d707 1479
1480 if(nClusterTrig>0){
feffe705 1481 //fClsETime1->Fill(clsE,clusterTime);
fb87d707 1482 }
1483
1484 if(nClusterTrig>0){
1485 cluster->SetChi2(1);
feffe705 1486 //fClsEAftTrigCut1->Fill(clsE);
fb87d707 1487 }
1488
1489 if(nClusterTrigCut>0){
1490 cluster->SetChi2(2);
feffe705 1491 //fClsEAftTrigCut2->Fill(clsE);
fb87d707 1492 }
1493
1494 if(nClusterTrigCut>0 && ( nCell > (1 + clsE / 3)))
1495 {
1496 cluster->SetChi2(3);
feffe705 1497 //fClsEAftTrigCut3->Fill(clsE);
fb87d707 1498 }
1499
1500 if(nClusterTrigCut>0 && (nCell > (1 + clsE / 3) )&&( clusterTime > fTimeCutLow && clusterTime < fTimeCutHigh ))
1501 {
1502 // cluster->SetChi2(4);
feffe705 1503 //fClsEAftTrigCut4->Fill(clsE);
fb87d707 1504 }
1505 if(nClusterTrigCut<1)
1506 {
1507 cluster->SetChi2(0);
1508
feffe705 1509 //fClsEAftTrigCut->Fill(clsE);
fb87d707 1510 }
1511
1512 } // clusters
1513}
1514
1515
1516
1517
1518
1519
1520