AliPhysics  master (3d17d9d)
AliAnalysisTaskEmcalHfeTagging.cxx
Go to the documentation of this file.
1 //
2 // HFE tagged jet shape analysis task.
3 //
4 // Authors: D. Caffarri, L. Cunqueiro (jet); D. Godoy (HFE)
5 //
6 #include <TClonesArray.h>
7 #include <TH1F.h>
8 #include <TH2F.h>
9 #include <TH3F.h>
10 #include <THnSparse.h>
11 #include <TTree.h>
12 #include <TList.h>
13 #include <TLorentzVector.h>
14 #include <TProfile.h>
15 #include <TChain.h>
16 #include <TSystem.h>
17 #include <TFile.h>
18 #include <TKey.h>
19 #include <AliAnalysisDataSlot.h>
20 #include <AliAnalysisDataContainer.h>
21 #include "TMatrixD.h"
22 #include "TMatrixDSym.h"
23 #include "TMatrixDSymEigen.h"
24 #include "TVector3.h"
25 #include "TVector2.h"
26 #include "AliVCluster.h"
27 #include "AliVTrack.h"
28 #include "AliEmcalJet.h"
29 #include "AliRhoParameter.h"
30 #include "AliLog.h"
31 #include "AliEmcalParticle.h"
32 #include "AliMCEvent.h"
33 #include "AliGenPythiaEventHeader.h"
34 #include "AliAODMCHeader.h"
35 #include "AliAODHandler.h"
36 #include "AliAODEvent.h"
37 #include "AliMCEvent.h"
38 #include "AliAnalysisManager.h"
39 #include "AliJetContainer.h"
40 #include "AliParticleContainer.h"
41 #include "AliEmcalPythiaInfo.h"
42 #include "TRandom3.h"
43 #include "AliAODInputHandler.h"
44 #include "AliPIDResponse.h"
45 #include "AliCFContainer.h"
46 #include "AliCFManager.h"
47 #include "AliMultiEventInputHandler.h"
48 #include "AliKFParticle.h"
49 #include "AliKFVertex.h"
50 #include "AliMultSelection.h"
51 #include "AliAnalysisUtils.h"
52 
53 #include "AliAODEvent.h"
55 
56 using std::cout;
57 using std::endl;
58 
60 
61 //________________________________________________________________________
64 fAOD(0),
65 fVevent(0),
66 fpidResponse(0),
67 fTracksTender(0),
68 pVtx(0),
69 spdVtx(0),
70 fMC(0),
71 fStack(0),
72 fMCparticle(0),
73 fMCarray(0),
74 fMCheader(0),
75 fContainer(0),
76 fMinFractionShared(0),
77 fJetShapeType(kData),
78 fJetShapeSub(kNoSub),
79 fJetSelection(kInclusive),
80 fPtThreshold(-9999.),
81 fRMatching(0.2),
82 fSelectedShapes(0),
83 fminpTTrig(20.),
84 fmaxpTTrig(50.),
85 fangWindowRecoil(0.6),
86 fSemigoodCorrect(0),
87 fHolePos(0),
88 fHoleWidth(0),
89 fCentSelectOn(kTRUE),
90 fCentMin(0),
91 fCentMax(10),
92 fOneConstSelectOn(kFALSE),
93 fDerivSubtrOrder(0),
94 fMCweight(0),
95 fRunNumber(0),
96 fAssPtCut(0.1),
97 fITSncut(3),
98 fAssTPCnCut(60),
99 fTPCnCut(100),
100 fAssITSrefitCut(kTRUE),
101 fUseTender(kFALSE),
102 fSigmaTOFcut(3.),
103 fSigmaTPCcut(-1.),
104 fDcaXYcut(1.),
105 fDcaZcut(2.),
106 fIMcut(0.1),
107 fNeventV0(0),
108 fNeventT0(0),
109 fh2ResponseUW(0x0),
110 fh2ResponseW(0x0),
111 fPhiJetCorr6(0x0),
112 fPhiJetCorr7(0x0),
113 fEtaJetCorr6(0x0),
114 fEtaJetCorr7(0x0),
115 fPtJetCorr(0x0),
116 fPtJet(0x0),
117 fPtGenJet(0),
118 fPhiJet(0x0),
119 fEtaJet(0x0),
120 fEtaPhiJet(0x0),
121 fAreaJet(0x0),
122 fJetProbDensityDetPart(0x0),
123 fJetProbDensityPartDet(0x0),
124 fNbOfConstvspT(0x0),
125 fnTPCnTOFnocut(0),
126 fnTPCnocutP(0),
127 fnTOFnocutP(0),
128 fnTPCcutP(0),
129 fnTPCcutPt(0),
130 fnULSmLSpairsPerElectron(0),
131 fnPartPerJet(0),
132 fnElecOverPartPerJet(0),
133 fnInclElecPerJet(0),
134 fnPhotElecPerJet(0),
135 fnIncSubPhotElecPerJet(0),
136 fnTrueElecPerJet(0),
137 fnTrueHFElecPerJet(0),
138 fnTruePElecPerJet(0),
139 fPi0PtGen(0),
140 fPi0PtEnh(0),
141 fEtaPtGen(0),
142 fEtaPtEnh(0),
143 fGenHfePt(0),
144 fGenPePt(0),
145 fPtP(0),
146 fptJetIE(0),
147 fptJetPE(0),
148 fptJetHFE(0),
149 fptRecPE(0),
150 fptTruePE(0),
151 fptWrongPE(0),
152 fPtTrack(0),
153 fPhiTrack(0x0),
154 fEtaTrack(0x0),
155 fEtaPhiTrack(0x0),
156 fPhiRecElecTPC(0x0),
157 fEtaRecElecTPC(0x0),
158 fEtaPhiRecElecTPC(0x0),
159 fPhiRecElecEMCal(0x0),
160 fEtaRecElecEMCal(0x0),
161 fEtaPhiRecElecEMCal(0x0),
162 fPhiTrueElec(0x0),
163 fEtaTrueElec(0x0),
164 fEtaPhiTrueElec(0x0),
165 fnEovPelecNoTPCcut(0x0),
166 fnEovPelecTPCcut(0x0),
167 fnEovPelecTPCEMCalcut(0x0),
168 fnEovPbackg(0x0),
169 fnClsE(0x0),
170 fnM20(0x0),
171 fnM02(0x0),
172 fnClsTime(0x0),
173 fAngULS(0x0),
174 fAngLS(0x0),
175 fAngChargPart(0x0),
176 fAngHadron(0x0),
177 fAngIncElec(0x0),
178 fAngPhotElec(0x0),
179 fAngElecFromD(0x0),
180 fAngElecFromB(0x0),
181 fAngElecFromDFromB(0x0),
182 fAngD(0x0),
183 fAngB(0x0),
184 fAngCharm(0x0),
185 fAngBeauty(0x0),
186 fAngQuark(0x0),
187 fAngGluon(0x0),
188 fDispULS(0x0),
189 fDispLS(0x0),
190 fDispChargPart(0x0),
191 fDispHadron(0x0),
192 fDispIncElec(0x0),
193 fDispPhotElec(0x0),
194 fDispElecFromD(0x0),
195 fDispElecFromB(0x0),
196 fDispElecFromDFromB(0x0),
197 fDispD(0x0),
198 fDispB(0x0),
199 fDispCharm(0x0),
200 fDispBeauty(0x0),
201 fDispQuark(0x0),
202 fDispGluon(0x0),
203 fTreeObservableTagging(0)
204 
205 {
206  for(Int_t i=0;i<26;i++){
207  fShapesVar[i]=0;
208  }
209 
210  for(Int_t i=0;i<5;i++){
211  fptTrueHFEeffTPCTOF[i] = NULL;
212  fptTrueHFEeffEMCal[i] = NULL;
213  fnEovPelecTPCsscut[i] = NULL;
214  }
215 
216  for(Int_t i=0;i<2;i++){
217  for(Int_t j=0;j<5;j++){
218  fptTrueHFEeffTPCTOFwJetPt[i][j] = NULL;
219  fptTrueHFEeffEMCalwJetPt[i][j] = NULL;
220  }
221  }
222 
223 
224  for(Int_t i=0;i<5;i++){
225  fInvmassLS[i] = NULL;
226  fInvmassULS[i] = NULL;
227  fRecPEJetPt[i] = NULL;
228  fTotPEJetPt[i] = NULL;
229  fRecPEAng[i] = NULL;
230  fTotPEAng[i] = NULL;
231  for(Int_t j=0;j<18;j++){
232  fnTPCSigma[i][j] = NULL;
233  }
234  }
235 
236  for(Int_t i=0;i<6;i++){
237  fRecPEDisp[i] = NULL;
238  fTotPEDisp[i] = NULL;
239  }
240 
241 
242  SetMakeGeneralHistograms(kTRUE);
243  DefineOutput(1, TList::Class());
244  DefineOutput(2, TTree::Class());
245 }
246 
247 //________________________________________________________________________
249 AliAnalysisTaskEmcalJet(name, kTRUE),
250 fAOD(0),
251 fVevent(0),
252 fpidResponse(0),
253 fTracksTender(0),
254 pVtx(0),
255 spdVtx(0),
256 fMC(0),
257 fStack(0),
258 fMCparticle(0),
259 fMCarray(0),
260 fMCheader(0),
261 fContainer(0),
262 fMinFractionShared(0),
263 fJetShapeType(kData),
264 fJetShapeSub(kNoSub),
265 fJetSelection(kInclusive),
266 fPtThreshold(-9999.),
267 fRMatching(0.2),
268 fSelectedShapes(0),
269 fminpTTrig(20.),
270 fmaxpTTrig(50.),
271 fangWindowRecoil(0.6),
272 fSemigoodCorrect(0),
273 fHolePos(0),
274 fHoleWidth(0),
275 fCentSelectOn(kTRUE),
276 fCentMin(0),
277 fCentMax(10),
278 fOneConstSelectOn(kFALSE),
279 fDerivSubtrOrder(0),
280 fMCweight(0),
281 fRunNumber(0),
282 fAssPtCut(0.1),
283 fITSncut(3),
284 fAssTPCnCut(60),
285 fTPCnCut(100),
286 fAssITSrefitCut(kTRUE),
287 fUseTender(kFALSE),
288 fSigmaTOFcut(3.),
289 fSigmaTPCcut(-1.),
290 fDcaXYcut(1.),
291 fDcaZcut(2.),
292 fIMcut(0.1),
293 fNeventV0(0),
294 fNeventT0(0),
295 fh2ResponseUW(0x0),
296 fh2ResponseW(0x0),
297 fPhiJetCorr6(0x0),
298 fPhiJetCorr7(0x0),
299 fEtaJetCorr6(0x0),
300 fEtaJetCorr7(0x0),
301 fPtJetCorr(0x0),
302 fPtJet(0x0),
303 fPtGenJet(0),
304 fPhiJet(0x0),
305 fEtaJet(0x0),
306 fEtaPhiJet(0x0),
307 fAreaJet(0x0),
308 fJetProbDensityDetPart(0x0),
309 fJetProbDensityPartDet(0x0),
310 fNbOfConstvspT(0x0),
311 fnTPCnTOFnocut(0),
312 fnTPCnocutP(0),
313 fnTOFnocutP(0),
314 fnTPCcutP(0),
315 fnTPCcutPt(0),
316 fnULSmLSpairsPerElectron(0),
317 fnPartPerJet(0),
318 fnElecOverPartPerJet(0),
319 fnInclElecPerJet(0),
320 fnPhotElecPerJet(0),
321 fnIncSubPhotElecPerJet(0),
322 fnTrueElecPerJet(0),
323 fnTrueHFElecPerJet(0),
324 fnTruePElecPerJet(0),
325 fPi0PtGen(0),
326 fPi0PtEnh(0),
327 fEtaPtGen(0),
328 fEtaPtEnh(0),
329 fGenHfePt(0),
330 fGenPePt(0),
331 fPtP(0),
332 fptJetIE(0),
333 fptJetPE(0),
334 fptJetHFE(0),
335 fptRecPE(0),
336 fptTruePE(0),
337 fptWrongPE(0),
338 fPtTrack(0),
339 fPhiTrack(0x0),
340 fEtaTrack(0x0),
341 fEtaPhiTrack(0x0),
342 fPhiRecElecTPC(0x0),
343 fEtaRecElecTPC(0x0),
344 fEtaPhiRecElecTPC(0x0),
345 fPhiRecElecEMCal(0x0),
346 fEtaRecElecEMCal(0x0),
347 fEtaPhiRecElecEMCal(0x0),
348 fPhiTrueElec(0x0),
349 fEtaTrueElec(0x0),
350 fEtaPhiTrueElec(0x0),
351 fnEovPelecNoTPCcut(0x0),
352 fnEovPelecTPCcut(0x0),
353 fnEovPelecTPCEMCalcut(0x0),
354 fnEovPbackg(0x0),
355 fnClsE(0x0),
356 fnM20(0x0),
357 fnM02(0x0),
358 fnClsTime(0x0),
359 fAngULS(0x0),
360 fAngLS(0x0),
361 fAngChargPart(0x0),
362 fAngHadron(0x0),
363 fAngIncElec(0x0),
364 fAngPhotElec(0x0),
365 fAngElecFromD(0x0),
366 fAngElecFromB(0x0),
367 fAngElecFromDFromB(0x0),
368 fAngD(0x0),
369 fAngB(0x0),
370 fAngCharm(0x0),
371 fAngBeauty(0x0),
372 fAngQuark(0x0),
373 fAngGluon(0x0),
374 fDispULS(0x0),
375 fDispLS(0x0),
376 fDispChargPart(0x0),
377 fDispHadron(0x0),
378 fDispIncElec(0x0),
379 fDispPhotElec(0x0),
380 fDispElecFromD(0x0),
381 fDispElecFromB(0x0),
382 fDispElecFromDFromB(0x0),
383 fDispD(0x0),
384 fDispB(0x0),
385 fDispCharm(0x0),
386 fDispBeauty(0x0),
387 fDispQuark(0x0),
388 fDispGluon(0x0),
389 fTreeObservableTagging(0)
390 
391 {
392  // Standard constructor.
393  for(Int_t i=0;i<26;i++){
394  fShapesVar[i]=0;
395  }
396 
397  for(Int_t i=0;i<5;i++){
398  fptTrueHFEeffTPCTOF[i] = NULL;
399  fptTrueHFEeffEMCal[i] = NULL;
400  fnEovPelecTPCsscut[i] = NULL;
401  }
402 
403  for(Int_t i=0;i<2;i++){
404  for(Int_t j=0;j<5;j++){
405  fptTrueHFEeffTPCTOFwJetPt[i][j] = NULL;
406  fptTrueHFEeffEMCalwJetPt[i][j] = NULL;
407  }
408  }
409 
410  for(Int_t i=0;i<5;i++){
411  fInvmassLS[i] = NULL;
412  fInvmassULS[i] = NULL;
413  fRecPEJetPt[i] = NULL;
414  fTotPEJetPt[i] = NULL;
415  fRecPEAng[i] = NULL;
416  fTotPEAng[i] = NULL;
417  for(Int_t j=0;j<18;j++){
418  fnTPCSigma[i][j] = NULL;
419  }
420  }
421 
422  for(Int_t i=0;i<6;i++){
423  fRecPEDisp[i] = NULL;
424  fTotPEDisp[i] = NULL;
425  }
426 
427 
429  DefineOutput(1, TList::Class());
430  DefineOutput(2, TTree::Class());
431 }
432 
433 //________________________________________________________________________
435 {
436  // Destructor.
437 }
438 
439 //________________________________________________________________________
441 {
442  // Create user output.
443 
445 
446  Bool_t oldStatus = TH1::AddDirectoryStatus();
447  TH1::AddDirectory(kFALSE);
448 
449  Double_t ptRange[34] = {0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2, 1.3, 1.4,
450  1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 5.,
451  6., 8., 10., 12., 14., 16., 19., 22., 26., 30.,
452  35., 40., 45., 50.};
453 
454  int nbin = 59;
455  double xbins[60] = {0.01,0.1,0.12,0.14,0.16,0.18,0.2,0.25,0.3,0.35,
456  0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,
457  0.9,0.95,1,1.1,1.2,1.3,1.4,1.5,1.6,1.7,
458  1.8,1.9,2,2.2,2.4,2.6,2.8,3,3.2,3.4,
459  3.6,3.8,4,4.5,5,5.5,6,6.5,7,8,
460  9,10,11,12,13,14,15,16,18,20};
461 
462  Double_t bin_JetPt[6] = {5.,20.,40.,60.,80.,120.};
463  Double_t bin_g[9] = {0.,0.02,0.03,0.04,0.05,0.06,0.07,0.08,0.12};//g
464  Double_t bin_ptd[7] = {0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 1.0};//pTD
465 
466  fNeventV0 = new TH1F("fNeventV0","Number of Events (V0)",5,-0.5,4.5);
467  fOutput->Add(fNeventV0);
468 
469  fNeventT0 = new TH1F("fNeventT0","Number of Events (T0)",5,-0.5,4.5);
470  fOutput->Add(fNeventT0);
471 
472  fh2ResponseUW= new TH2F("fh2ResponseUW", "fh2ResponseUW", 100, 0, 200, 100, 0, 200);
473  fOutput->Add(fh2ResponseUW);
474 
475  fh2ResponseW= new TH2F("fh2ResponseW", "fh2ResponseW", 100, 0, 200, 100, 0, 200);
476  fOutput->Add(fh2ResponseW);
477 
478  fPhiJetCorr6= new TH2F("fPhiJetCorr6", "fPhiJetCorr6", 50, 0, 2*TMath::Pi(), 50, 0, 2*TMath::Pi());
479  fOutput->Add(fPhiJetCorr6);
480 
481  fEtaJetCorr6= new TH2F("fEtaJetCorr6", "fEtaJetCorr6", 50, -1.5, 1.5, 50, -1.5, 1.5);
482  fOutput->Add(fEtaJetCorr6);
483 
484  fPhiJetCorr7= new TH2F("fPhiJetCorr7", "fPhiJetCorr7", 50, 0, 2*TMath::Pi(), 50, 0, 2*TMath::Pi());
485  fOutput->Add(fPhiJetCorr7);
486 
487  fEtaJetCorr7= new TH2F("fEtaJetCorr7", "fEtaJetCorr7", 50, -1.5, 1.5, 50, -1.5, 1.5);
488  fOutput->Add(fEtaJetCorr7);
489 
490  fPtJetCorr= new TH2F("fPtJetCorr", "fPtJetCorr", 100, 0, 200, 100, 0, 200);
491  fOutput->Add(fPtJetCorr);
492 
493  fPtJet= new TH1F("fPtJet", "fPtJet", 100, 0, 200);
494  fOutput->Add(fPtJet);
495 
496  fPtGenJet= new TH1F("fPtGenJet", "fPtGenJet", 100, 0, 200);
497  fOutput->Add(fPtGenJet);
498 
499  fPhiJet= new TH2F("fPhiJet", "fPhiJet", 100, 0, 200, 100, 0, TMath::TwoPi());
500  fOutput->Add(fPhiJet);
501 
502  fEtaJet= new TH2F("fEtaJet", "fEtaJet", 100, 0, 200, 100, -1,1);
503  fOutput->Add(fEtaJet);
504 
505  fEtaPhiJet= new TH2F("fEtaPhiJet", "fEtaPhiJet", 100, 0, TMath::TwoPi(), 100, -1,1);
506  fOutput->Add(fEtaPhiJet);
507 
508  fAreaJet= new TH2F("fAreaJet", "fAreaJet", 100, 0, 200, 100, 0,1.5);
509  fOutput->Add(fAreaJet);
510 
511  fJetProbDensityDetPart=new TH2F("fJetProbDensityDetPart", "fJetProbDensityDetPart",100,-2,2,100, 0,200);
513 
514  fJetProbDensityPartDet=new TH2F("fJetProbDensityPartDet", "fJetProbDensityPartDet",100,-2,2,100, 0,200);
516 
517  fNbOfConstvspT=new TH2F("fNbOfConstvspT", "fNbOfConstvspT", 100, 0, 100, 200, 0, 200);
518  fOutput->Add(fNbOfConstvspT);
519 
520  fnTPCnTOFnocut=new TH2F("fnTPCnTOFnocut", "fnTPCnTOFnocut", 100, -10, 10, 100, -10, 10);
521  fOutput->Add(fnTPCnTOFnocut);
522 
523  fnTPCnocutP=new TH2F("fnTPCnocutP", "fnTPCnocutP", 50, 0, 5, 200, -10, 10);
524  fOutput->Add(fnTPCnocutP);
525 
526  fnTOFnocutP=new TH2F("fnTOFnocutP", "fnTOFnocutP", 50, 0, 5, 200, -10, 10);
527  fOutput->Add(fnTOFnocutP);
528 
529  fnTPCcutP=new TH2F("fnTPCcutP", "fnTPCcutP", 50, 0, 5, 200, -10, 10);
530  fOutput->Add(fnTPCcutP);
531 
532  fnTPCcutPt=new TH2F("fnTPCcutPt", "fnTPCcutPt", 50, 0, 5, 200, -10, 10);
533  fOutput->Add(fnTPCcutPt);
534 
535  for(Int_t i=0;i<5;i++){
536  for(Int_t j=0;j<18;j++){
537  fnTPCSigma[i][j] = new TH1F(Form("fnTPCSigma%d%d",i,j),Form("fnTPCSigma%d%d",i,j), 100,-15,15);
538  fOutput->Add(fnTPCSigma[i][j]);
539  }
540  }
541 
542  fnULSmLSpairsPerElectron =new TH2F("fnULSmLSpairsPerElectron", "fnULSmLSpairsPerElectron", 50, 0, 5, 100, 0, 100);
544 
545  for(Int_t i=0;i<5;i++){
546  fInvmassLS[i] = new TH2F(Form("fInvmassLS%d",i), Form("fInvmassLS%d",i), 33,ptRange, 100, 0, 0.5);
547  fOutput->Add(fInvmassLS[i]);
548  }
549 
550  for(Int_t i=0;i<5;i++){
551  fInvmassULS[i] = new TH2F(Form("fInvmassULS%d",i), Form("fInvmassULS%d",i), 33,ptRange, 100, 0, 0.5);
552  fOutput->Add(fInvmassULS[i]);
553  }
554 
555  fnPartPerJet=new TH1F("fnPartPerJet", "fnPartPerJet", 50,0,50);
556  fOutput->Add(fnPartPerJet);
557 
558  fnElecOverPartPerJet=new TH1F("fnElecOverPartPerJet", "fnElecOverPartPerJet", 50,0,0.1);
560 
561  fnInclElecPerJet=new TH1F("fnInclElecPerJet", "fnInclElecPerJet", 50,0,50);
563 
564  fnPhotElecPerJet=new TH1F("fnPhotElecPerJet", "fnPhotElecPerJet", 50,0,50);
566 
567  fnIncSubPhotElecPerJet=new TH1F("fnIncSubPhotElecPerJet", "fnIncSubPhotElecPerJet", 51,-1,50);
569 
570  fnTrueElecPerJet=new TH1F("fnTrueElecPerJet", "fnTrueElecPerJet", 50,0,50);
572 
573  fnTrueHFElecPerJet=new TH1F("fnTrueHFElecPerJet", "fnTrueHFElecPerJet", 50,0,50);
575 
576  fnTruePElecPerJet=new TH1F("fnTruePElecPerJet", "fnTruePElecPerJet", 50,0,50);
578 
579  fPi0PtGen = new TH1F("fPi0PtGen","fPi0PtGen",nbin,xbins);
580  fOutput->Add(fPi0PtGen);
581 
582  fPi0PtEnh = new TH1F("fPi0PtEnh","fPi0PtEnh",nbin,xbins);
583  fOutput->Add(fPi0PtEnh);
584 
585  fEtaPtGen = new TH1F("fEtaPtGen", "fEtaPtGen",nbin,xbins);
586  fOutput->Add(fEtaPtGen);
587 
588  fEtaPtEnh = new TH1F("fEtaPtEnh", "fEtaPtEnh",nbin,xbins);
589  fOutput->Add(fEtaPtEnh);
590 
591  fGenHfePt = new TH1F("fGenHfePt","fGenHfePt",nbin,xbins);
592  fOutput->Add(fGenHfePt);
593 
594  fGenPePt = new TH1F("fGenPePt", "fGenPePt",nbin,xbins);
595  fOutput->Add(fGenPePt);
596 
597  for(Int_t i=0;i<5;i++){
598  fRecPEJetPt[i] = new TH1F(Form("fRecPEJetPt%d",i), Form("fRecPEJetPt%d",i), 33,ptRange);
599  fOutput->Add(fRecPEJetPt[i]);
600  }
601 
602  for(Int_t i=0;i<5;i++){
603  fTotPEJetPt[i] = new TH1F(Form("fTotPEJetPt%d",i), Form("fTotPEJetPt%d",i), 33,ptRange);
604  fOutput->Add(fTotPEJetPt[i]);
605  }
606 
607  for(Int_t i=0;i<5;i++){
608  fRecPEAng[i] = new TH1F(Form("fRecPEAng%d",i), Form("fRecPEAng%d",i), 33,ptRange);
609  fOutput->Add(fRecPEAng[i]);
610  }
611 
612  for(Int_t i=0;i<5;i++){
613  fTotPEAng[i] = new TH1F(Form("fTotPEAng%d",i), Form("fTotPEAng%d",i), 33,ptRange);
614  fOutput->Add(fTotPEAng[i]);
615  }
616 
617  for(Int_t i=0;i<6;i++){
618  fRecPEDisp[i] = new TH1F(Form("fRecPEDisp%d",i), Form("fRecPEDisp%d",i), 33,ptRange);
619  fOutput->Add(fRecPEDisp[i]);
620  }
621 
622  for(Int_t i=0;i<6;i++){
623  fTotPEDisp[i] = new TH1F(Form("fTotPEDisp%d",i), Form("fTotPEDisp%d",i), 33,ptRange);
624  fOutput->Add(fTotPEDisp[i]);
625  }
626 
627  fPtP=new TH2F("fPtP", "fPtP", 33,ptRange,33,ptRange);
628  fOutput->Add(fPtP);
629 
630  fptJetIE= new TH1F("fptJetIE", "fptJetIE", 100, 0, 200);
631  fOutput->Add(fptJetIE);
632 
633  fptJetPE= new TH1F("fptJetPE", "fptJetPE", 100, 0, 200);
634  fOutput->Add(fptJetPE);
635 
636  fptJetHFE= new TH1F("fptJetHFE", "fptJetHFE", 100, 0, 200);
637  fOutput->Add(fptJetHFE);
638 
639  fptRecPE= new TH1F("fptRecPE", "fptRecPE", 33,ptRange);
640  fOutput->Add(fptRecPE);
641 
642  fptTruePE= new TH1F("fptTruePE", "fptTruePE", 33,ptRange);
643  fOutput->Add(fptTruePE);
644 
645  for(Int_t i=0;i<5;i++){
646  fptTrueHFEeffTPCTOF[i] = new TH1F(Form("fptTrueHFEeffTPCTOF%d",i), Form("fptTrueHFEeffTPCTOF%d",i), 33,ptRange);
647  fOutput->Add(fptTrueHFEeffTPCTOF[i]);
648  }
649 
650  for(Int_t i=0;i<2;i++){
651  for(Int_t j=0;j<5;j++){
652  fptTrueHFEeffTPCTOFwJetPt[i][j] = new TH1F(Form("fptTrueHFEeffTPCTOFwJetPt%d%d",i,j), Form("fptTrueHFEeffTPCTOFwJetPt%d%d",i,j), 33,ptRange);
654  }
655  }
656 
657 
658  for(Int_t i=0;i<5;i++){
659  fptTrueHFEeffEMCal[i] = new TH1F(Form("fptTrueHFEeffEMCal%d",i), Form("fptTrueHFEeffEMCal%d",i), 33,ptRange);
660  fOutput->Add(fptTrueHFEeffEMCal[i]);
661  }
662 
663  for(Int_t i=0;i<2;i++){
664  for(Int_t j=0;j<5;j++){
665  fptTrueHFEeffEMCalwJetPt[i][j] = new TH1F(Form("fptTrueHFEeffEMCalwJetPt%d%d",i,j), Form("fptTrueHFEeffEMCalwJetPt%d%d",i,j), 33,ptRange);
666  fOutput->Add(fptTrueHFEeffEMCalwJetPt[i][j]);
667  }
668  }
669 
670  fptWrongPE= new TH1F("fptWrongPE", "fptWrongPE", 33,ptRange);
671  fOutput->Add(fptWrongPE);
672 
673  fPtTrack= new TH1F("fPtTrack", "fPtTrack", 100, 0, 200);
674  fOutput->Add(fPtTrack);
675 
676  fPhiTrack= new TH2F("fPhiTrack", "fPhiTrack", 100, 0, 200, 100, 0, TMath::TwoPi());
677  fOutput->Add(fPhiTrack);
678 
679  fEtaTrack= new TH2F("fEtaTrack", "fEtaTrack", 100, 0, 200, 100, -1,1);
680  fOutput->Add(fEtaTrack);
681 
682  fEtaPhiTrack= new TH2F("fEtaPhiTrack", "fEtaPhiTrack", 100, 0, TMath::TwoPi(), 100, -1,1);
683  fOutput->Add(fEtaPhiTrack);
684 
685  fPhiRecElecTPC= new TH2F("fPhiRecElecTPC", "fPhiRecElecTPC", 100, 0, 50, 100, 0, TMath::TwoPi());
686  fOutput->Add(fPhiRecElecTPC);
687 
688  fEtaRecElecTPC= new TH2F("fEtaRecElecTPC", "fEtaRecElecTPC", 100, 0, 50, 100, -1,1);
689  fOutput->Add(fEtaRecElecTPC);
690 
691  fEtaPhiRecElecTPC= new TH2F("fEtaPhiRecElecTPC", "fEtaPhiRecElecTPC", 100, 0, TMath::TwoPi(), 100, -1,1);
693 
694  fPhiRecElecEMCal= new TH2F("fPhiRecElecEMCal", "fPhiRecElecEMCal", 100, 0, 50, 100, 0, TMath::TwoPi());
696 
697  fEtaRecElecEMCal= new TH2F("fEtaRecElecEMCal", "fEtaRecElecEMCal", 100, 0, 50, 100, -1,1);
699 
700  fEtaPhiRecElecEMCal= new TH2F("fEtaPhiRecElecEMCal", "fEtaPhiRecElecEMCal", 100, 0, TMath::TwoPi(), 100, -1,1);
702 
703  fPhiTrueElec= new TH2F("fPhiTrueElec", "fPhiTrueElec", 100, 0, 50, 100, 0, TMath::TwoPi());
704  fOutput->Add(fPhiTrueElec);
705 
706  fEtaTrueElec= new TH2F("fEtaTrueElec", "fEtaTrueElec", 100, 0, 50, 100, -1,1);
707  fOutput->Add(fEtaTrueElec);
708 
709  fEtaPhiTrueElec= new TH2F("fEtaPhiTrueElec", "fEtaPhiTrueElec", 100, 0, TMath::TwoPi(), 100, -1,1);
710  fOutput->Add(fEtaPhiTrueElec);
711 
712  fnEovPelecNoTPCcut = new TH2F("fnEovPelecNoTPCcut", "fnEovPelecNoTPCcut", 100, 0, 100, 40,0,2);
714 
715  fnEovPelecTPCcut = new TH2F("fnEovPelecTPCcut", "fnEovPelecTPCcut", 100, 0, 100, 40,0,2);
717 
718  fnEovPelecTPCEMCalcut = new TH2F("fnEovPelecTPCEMCalcut", "fnEovPelecTPCEMCalcut", 100, 0, 100,40,0,2);
720 
721  for(Int_t i=0;i<5;i++){
722  fnEovPelecTPCsscut[i] = new TH2F(Form("fnEovPelecTPCsscut%d",i), Form("fnEovPelecTPCsscut%d",i), 33,ptRange,40,0,2);
723  fOutput->Add(fnEovPelecTPCsscut[i]);
724  }
725 
726  fnEovPbackg = new TH2F("fnEovPbackg", "fnEovPbackg", 100, 0, 100,40,0,2);
727  fOutput->Add(fnEovPbackg);
728 
729  fnClsE = new TH2F("fnClsE", "fnClsE", 100, 0, 100, 100, 0,100);
730  fOutput->Add(fnClsE);
731 
732  fnM20 = new TH2F("fnM20", "fnM20", 100, 0, 100, 100, 0,2);
733  fOutput->Add(fnM20);
734 
735  fnM02 = new TH2F("fnM02", "fnM02", 100, 0, 100, 100, 0,2);
736  fOutput->Add(fnM02);
737 
738  fnClsTime = new TH2F("fnClsTime", "fnClsTime", 100, 0, 100, 100, -200,200);
739  fOutput->Add(fnClsTime);
740 
741  fAngULS = new TH2F("fAngULS", "fAngULS", 5, bin_JetPt, 8, bin_g);
742  fOutput->Add(fAngULS);
743 
744  fAngLS = new TH2F("fAngLS", "fAngLS", 5, bin_JetPt, 8, bin_g);
745  fOutput->Add(fAngLS);
746 
747  fAngChargPart = new TH2F("fAngChargPart", "fAngChargPart", 5, bin_JetPt, 8, bin_g);
748  fOutput->Add(fAngChargPart);
749 
750  fAngHadron = new TH2F("fAngHadron", "fAngHadron", 5, bin_JetPt, 8, bin_g);
751  fOutput->Add(fAngHadron);
752 
753  fAngIncElec = new TH2F("fAngIncElec", "fAngIncElec", 5, bin_JetPt, 8, bin_g);
754  fOutput->Add(fAngIncElec);
755 
756  fAngPhotElec = new TH2F("fAngPhotElec", "fAngPhotElec", 5, bin_JetPt, 8, bin_g);
757  fOutput->Add(fAngPhotElec);
758 
759  fAngElecFromD = new TH2F("fAngElecFromD", "fAngElecFromD", 5, bin_JetPt, 8, bin_g);
760  fOutput->Add(fAngElecFromD);
761 
762  fAngElecFromB = new TH2F("fAngElecFromB", "fAngElecFromB", 5, bin_JetPt, 8, bin_g);
763  fOutput->Add(fAngElecFromB);
764 
765  fAngElecFromDFromB = new TH2F("fAngElecFromDFromB", "fAngElecFromDFromB", 5, bin_JetPt, 8, bin_g);
767 
768  fAngD = new TH2F("fAngD", "fAngD", 5, bin_JetPt, 8, bin_g);
769  fOutput->Add(fAngD);
770 
771  fAngB = new TH2F("fAngB", "fAngB", 5, bin_JetPt, 8, bin_g);
772  fOutput->Add(fAngB);
773 
774  fAngCharm = new TH2F("fAngCharm", "fAngCharm", 5, bin_JetPt, 8, bin_g);
775  fOutput->Add(fAngCharm);
776 
777  fAngBeauty = new TH2F("fAngBeauty", "fAngBeauty", 5, bin_JetPt, 8, bin_g);
778  fOutput->Add(fAngBeauty);
779 
780  fAngQuark = new TH2F("fAngQuark", "fAngQuark", 5, bin_JetPt, 8, bin_g);
781  fOutput->Add(fAngQuark);
782 
783  fAngGluon = new TH2F("fAngGluon", "fAngGluon", 5, bin_JetPt, 8, bin_g);
784  fOutput->Add(fAngGluon);
785 
786  fDispULS = new TH2F("fDispULS", "fDispULS", 5, bin_JetPt, 6, bin_ptd);
787  fOutput->Add(fDispULS);
788 
789  fDispLS = new TH2F("fDispLS", "fDispLS", 5, bin_JetPt, 6, bin_ptd);
790  fOutput->Add(fDispLS);
791 
792  fDispChargPart = new TH2F("fDispChargPart", "fDispChargPart", 5, bin_JetPt, 6, bin_ptd);
793  fOutput->Add(fDispChargPart);
794 
795  fDispHadron = new TH2F("fDispHadron", "fDispHadron", 5, bin_JetPt, 6, bin_ptd);
796  fOutput->Add(fDispHadron);
797 
798  fDispIncElec = new TH2F("fDispIncElec", "fDispIncElec", 5, bin_JetPt, 6, bin_ptd);
799  fOutput->Add(fDispIncElec);
800 
801  fDispPhotElec = new TH2F("fDispPhotElec", "fDispPhotElec", 5, bin_JetPt, 6, bin_ptd);
802  fOutput->Add(fDispPhotElec);
803 
804  fDispElecFromD = new TH2F("fDispElecFromD", "fDispElecFromD", 5, bin_JetPt, 6, bin_ptd);
805  fOutput->Add(fDispElecFromD);
806 
807  fDispElecFromB = new TH2F("fDispElecFromB", "fDispElecFromB", 5, bin_JetPt, 6, bin_ptd);
808  fOutput->Add(fDispElecFromB);
809 
810  fDispElecFromDFromB = new TH2F("fDispElecFromDFromB", "fDispElecFromDFromB", 5, bin_JetPt, 6, bin_ptd);
812 
813  fDispD = new TH2F("fDispD", "fDispD", 5, bin_JetPt, 6, bin_ptd);
814  fOutput->Add(fDispD);
815 
816  fDispB = new TH2F("fDispB", "fDispB", 5, bin_JetPt, 6, bin_ptd);
817  fOutput->Add(fDispB);
818 
819  fDispCharm = new TH2F("fDispCharm", "fDispCharm", 5, bin_JetPt, 6, bin_ptd);
820  fOutput->Add(fDispCharm);
821 
822  fDispBeauty = new TH2F("fDispBeauty", "fDispBeauty", 5, bin_JetPt, 6, bin_ptd);
823  fOutput->Add(fDispBeauty);
824 
825  fDispQuark = new TH2F("fDispQuark", "fDispQuark", 5, bin_JetPt, 6, bin_ptd);
826  fOutput->Add(fDispQuark);
827 
828  fDispGluon = new TH2F("fDispGluon", "fDispGluon", 5, bin_JetPt, 6, bin_ptd);
829  fOutput->Add(fDispGluon);
830 
831  // =========== Switch on Sumw2 for all histos ===========
832  for (Int_t i=0; i<fOutput->GetEntries(); ++i) {
833  TH1 *h = dynamic_cast<TH1*>(fOutput->At(i));
834  if (h){
835  h->Sumw2();
836  continue;
837  }
838  THnSparse *hn = dynamic_cast<THnSparse*>(fOutput->At(i));
839  if(hn)hn->Sumw2();
840  }
841 
842 
843  TH1::AddDirectory(oldStatus);
844  const Int_t nVar = 26;
845 
846  fTreeObservableTagging = new TTree(GetOutputSlot(2)->GetContainer()->GetName(), GetOutputSlot(2)->GetContainer()->GetName());
847 
848  TString *fShapesVarNames = new TString [nVar];
849 
850  fShapesVarNames[0] = "partonCode";
851  fShapesVarNames[1] = "ptJet";
852  fShapesVarNames[2] = "ptDJet";
853  fShapesVarNames[3] = "mJet";
854  // fShapesVarNames[4] = "nbOfConst";
855  fShapesVarNames[4] = "angularity";
856  fShapesVarNames[5] = "circularity";
857  fShapesVarNames[6] = "lesub";
858  fShapesVarNames[7] = "coronna";
859 
860  fShapesVarNames[8] = "ptJetMatch";
861  fShapesVarNames[9] = "ptDJetMatch";
862  fShapesVarNames[10] = "mJetMatch";
863  // fShapesVarNames[12] = "nbOfConstMatch";
864  fShapesVarNames[11] = "angularityMatch";
865  fShapesVarNames[12] = "circularityMatch";
866  fShapesVarNames[13] = "lesubMatch";
867  fShapesVarNames[14] = "coronnaMatch";
868  fShapesVarNames[15] = "weightPythia";
869  //fShapesVarNames[14]="ntrksEvt";
870  fShapesVarNames[16] = "rhoVal";
871  fShapesVarNames[17] = "rhoMassVal";
872  fShapesVarNames[18] = "ptUnsub";
873  fShapesVarNames[19] = "ptTrueHFE";
874  fShapesVarNames[20] = "ptElec";
875  fShapesVarNames[21] = "nInclElec";
876  fShapesVarNames[22] = "nPhotElec";
877  fShapesVarNames[23] = "hasElec";
878  fShapesVarNames[24] = "nTrueElectronsMC";
879  fShapesVarNames[25] = "nTrueHFElecMC";
880 
881 
882  for(Int_t ivar=0; ivar < nVar; ivar++){
883  //cout<<"looping over variables"<<endl;
884  fTreeObservableTagging->Branch(fShapesVarNames[ivar].Data(), &fShapesVar[ivar], Form("%s/F", fShapesVarNames[ivar].Data()));
885  }
886 
887  PostData(1,fOutput);
888  PostData(2,fTreeObservableTagging);
889 
890  delete [] fShapesVarNames;
891 }
892 
893 //________________________________________________________________________
895 {
896  // // Run analysis code here, if needed. It will be executed before FillHistograms().
897 
898  fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
899  if (!fAOD){
900  printf("ERROR: fAOD not available\n");
901  return kFALSE;
902  }
903 
904  fVevent = dynamic_cast<AliVEvent*>(InputEvent());
905  if(!fVevent){
906  printf("ERROR: fVEvent not available\n");
907  return kFALSE;
908  }
909 
910  fRunNumber = fVevent->GetRunNumber();
911 
912  //remaining event selection
913  pVtx = fVevent->GetPrimaryVertex();
914  spdVtx = fVevent->GetPrimaryVertexSPD();
915 
916  fMCarray = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
917  fMCheader = dynamic_cast<AliAODMCHeader*>(fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
918 
919  fNeventV0->Fill(0);
920  if (!fMCheader){
921  TString firedTriggerClasses = static_cast<const AliAODEvent*>(InputEvent())->GetFiredTriggerClasses();
922  if ((firedTriggerClasses.Contains("C0TVX-B-NOPF-CENT"))){
923  fNeventT0->Fill(0);
924  }
925  }
926 
927  Int_t NpureMC = -1, NpureMCproc = -1;
928 
929  if (fMCheader){
930  TList *lh=fMCheader->GetCocktailHeaders();
931  if(lh){
932  for(int igene=0; igene<lh->GetEntries(); igene++){
933  AliGenEventHeader* gh=(AliGenEventHeader*)lh->At(igene);
934  if(gh){
935  // cout << "<------- imc = "<< gh->GetName() << endl;
936  if(igene==0)NpureMC = gh->NProduced(); // generate by PYTHIA or HIJING
937  NpureMCproc += gh->NProduced();
938  }
939  }
940  }
941  }
942  // pileup rejection with SPD multiple vertices and track multivertexer (using default configuration)
943 
944  Bool_t isPileupfromSPDmulbins = kFALSE, isPileupFromMV = kFALSE;
945  isPileupfromSPDmulbins = fAOD->IsPileupFromSPDInMultBins();
946  if (isPileupfromSPDmulbins) return kFALSE;
947  fNeventV0->Fill(1);
948 
949  AliAnalysisUtils utils;
950  utils.SetMinPlpContribMV(5);
951  utils.SetMaxPlpChi2MV(5.);
952  utils.SetMinWDistMV(15);
953  utils.SetCheckPlpFromDifferentBCMV(kFALSE);
954  isPileupFromMV = utils.IsPileUpMV(fAOD);
955  if (isPileupFromMV) return kFALSE;
956  fNeventV0->Fill(2);
957 
958  // Selection of pi0 and eta in MC to compute the weight
959 
960  Bool_t isPrimary = kFALSE, isFromLMdecay = kTRUE, isFromHFdecay=kTRUE, hasMother = kTRUE;
961 
962  Double_t MCweight = 1.;
963  Int_t iDecay = 0;
964 
965  if(fMCarray){
966  //Int_t nParticles = fMCarray->GetEntries();
967  for (Int_t iParticle = 0; iParticle < NpureMCproc; iParticle++) {
968  AliAODMCParticle* particle = (AliAODMCParticle*) fMCarray->At(iParticle);
969  int fPDG = particle->GetPdgCode();
970  Double_t pTMC = particle->Pt();
971  Double_t phiMC = particle->Phi();
972  Double_t etaMC = particle->Eta();
973 
974  Bool_t iEnhance = kFALSE;
975  if(iParticle>=NpureMC)iEnhance = kTRUE;
976 
977  //if (TMath::Abs(etaMC)>1.2)continue;
978 
979  isPrimary = IsPrimary(particle);
980  isFromLMdecay = IsFromLMdecay(particle);
981  isFromHFdecay = IsFromHFdecay(particle);
982  hasMother = HasMother(particle);
983 
984  MCweight = 1.;
985  iDecay = 0;
986 
987  GetWeightAndDecay(particle,iDecay,MCweight);
988 
989 
990  if (TMath::Abs(etaMC)<0.8 && TMath::Abs(fPDG)==11){
991  fPhiTrueElec->Fill(pTMC,phiMC);
992  fEtaTrueElec->Fill(pTMC,etaMC);
993  if (pTMC > 0.5) fEtaPhiTrueElec->Fill(phiMC,etaMC);
994 
995  if (isFromHFdecay) fGenHfePt->Fill(pTMC);
996  if (iDecay>0 && iDecay<7) fGenPePt->Fill(pTMC);
997  }
998 
999 
1000  Double_t yMC = particle->Y();
1001  if (TMath::Abs(yMC)>1.0)continue;
1002 
1003  if (isPrimary){
1004  if (!hasMother || (!isFromLMdecay && !isFromHFdecay)){
1005  if(fPDG==111 && iEnhance==kFALSE) fPi0PtGen->Fill(pTMC);
1006  if(fPDG==111 && iEnhance==kTRUE) fPi0PtEnh->Fill(pTMC);
1007 
1008  if(fPDG==221 && iEnhance==kFALSE) fEtaPtGen->Fill(pTMC);
1009  if(fPDG==221 && iEnhance==kTRUE) fEtaPtEnh->Fill(pTMC);
1010  }
1011  }
1012  }
1013  }//MC
1014  return kTRUE;
1015 }
1016 //________________________________________________________________________
1018 {
1019  // Fill histograms.
1020  //cout<<"base container"<<endl;
1021  AliEmcalJet* jet1 = NULL;
1022  AliJetContainer *jetCont = GetJetContainer(0);
1023 
1024  Float_t kWeight=1;
1025  if (fCentSelectOn)
1026  if ((fCent>fCentMax) || (fCent<fCentMin)) return 0;
1027 
1028  // list of kink mothers
1029  Int_t nVertices = 1;
1030  nVertices = fAOD->GetNumberOfVertices();
1031  Double_t listofmotherkink[nVertices];
1032  Int_t nMotherKink = 0;
1033  for(Int_t ivertex=0; ivertex < nVertices; ivertex++) {
1034  AliAODVertex *vertex = fAOD->GetVertex(ivertex);
1035  if(!vertex) continue;
1036  if(vertex->GetType()==AliAODVertex::kKink) {
1037  AliAODTrack *mother = (AliAODTrack *) vertex->GetParent();
1038  if(!mother) continue;
1039  Int_t idmother = mother->GetID();
1040  listofmotherkink[nMotherKink] = idmother;
1041  nMotherKink++;
1042  }
1043  }
1044 
1045  Float_t rhoVal=0, rhoMassVal = 0.;
1046 
1047  if(jetCont) {
1048  jetCont->ResetCurrentID();
1050  //rho
1051  AliRhoParameter* rhoParam = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject("RhoSparseR020"));
1052  if (!rhoParam) {
1053  Printf("%s: Could not retrieve rho %s (some histograms will be filled with zero)!", GetName(), jetCont->GetRhoName().Data());
1054  } else rhoVal = rhoParam->GetVal();
1055  //rhom
1056  AliRhoParameter* rhomParam = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject("RhoMassSparseR020"));
1057  if (!rhomParam) {
1058  Printf("%s: Could not retrieve rho_m %s (some histograms will be filled with zero)!", GetName(), jetCont->GetRhoMassName().Data());
1059  } else rhoMassVal = rhomParam->GetVal();
1060  }
1061 
1062  while((jet1 = jetCont->GetNextAcceptJet())) {
1063  if (!jet1) continue;
1064  AliEmcalJet* jet2 = 0x0;
1065  AliEmcalJet* jet3 = 0x0;
1066  fPtJet->Fill(jet1->Pt());
1067  fPhiJet->Fill(jet1->Pt(),jet1->Phi());
1068  fEtaJet->Fill(jet1->Pt(),jet1->Eta());
1069  if (jet1->Pt() > 5.) fEtaPhiJet->Fill(jet1->Phi(),jet1->Eta());
1070  fAreaJet->Fill(jet1->Pt(),jet1->Area());
1071  AliEmcalJet *jetUS = NULL;
1072  Int_t ifound=0, ilab=-1;
1073 
1074  fShapesVar[0] = 0.;
1076  //AliJetContainer *jetContTrue = GetJetContainer(1);
1077  AliJetContainer *jetContUS = GetJetContainer(2);
1078 
1080  for(Int_t i = 0; i<jetContUS->GetNJets(); i++) {
1081  jetUS = jetContUS->GetJet(i);
1082  if(jetUS->GetLabel()==jet1->GetLabel()) {
1083  ifound++;
1084  if(ifound==1) ilab = i;
1085  }
1086  }
1087  if(ilab==-1) continue;
1088  jetUS=jetContUS->GetJet(ilab);
1089  jet2=jetUS->ClosestJet();
1090  }
1091 
1093  if (!jet2) {
1094  Printf("jet2 does not exist, returning");
1095  continue;
1096  }
1097 
1098  //AliJetContainer *jetContPart=GetJetContainer(3);
1099  jet3=jet2->ClosestJet();
1100 
1101  if(!jet3){
1102  Printf("jet3 does not exist, returning");
1103  continue;
1104  }
1105  //cout<<"jet 3 exists"<<jet3->Pt()<<endl;
1106 
1107 
1108  fh2ResponseUW->Fill(jet1->Pt(),jet2->Pt());
1109 
1110  Double_t fraction=0;
1111  if(!(fJetShapeSub==kConstSub)) fraction = jetCont->GetFractionSharedPt(jet1);
1112  if(fJetShapeSub==kConstSub) fraction = jetContUS->GetFractionSharedPt(jetUS);
1113  //if (fraction > 0.1) cout<<"***** hey a jet matched with fraction"<<fraction<<" "<<jet1->Pt()<<" "<<jet2->Pt()<<" "<<fCent<<endl;
1114 
1115  if(fraction<fMinFractionShared) continue;
1116  //InputEvent()->Print();
1117 
1118  }
1119 
1121  jet3 = jet1->ClosestJet();
1122  if (!jet3) continue;
1123 
1124  fh2ResponseUW->Fill(jet1->Pt(),jet3->Pt());
1125 
1126  Double_t probDensDetPart = -999., probDensPartDet = -999.;
1127 
1128  if (jet1->Pt()>0) probDensPartDet = (jet3->Pt()-jet1->Pt())/jet1->Pt();
1129  if (jet3->Pt()>0) probDensDetPart = (jet1->Pt()-jet3->Pt())/jet3->Pt();
1130 
1131  fJetProbDensityDetPart->Fill(probDensDetPart,jet3->Pt());
1132  fJetProbDensityPartDet->Fill(probDensPartDet,jet1->Pt());
1133 
1134  }
1135 
1136  if (fJetShapeType == kGenOnTheFly){
1137  const AliEmcalPythiaInfo *partonsInfo = 0x0;
1138  partonsInfo = GetPythiaInfo();
1139  Double_t jp1=RelativePhi(jet1->Phi(),partonsInfo->GetPartonPhi6());
1140  Double_t detap1=(jet1->Eta())-(partonsInfo->GetPartonEta6());
1141  kWeight=partonsInfo->GetPythiaEventWeight();
1142  fh2ResponseW->Fill(jet1->Pt(),jet1->Pt(),kWeight);
1143 
1144  Float_t dRp1 = TMath::Sqrt(jp1 * jp1 + detap1 * detap1);
1145  fEtaJetCorr6->Fill(jet1->Eta(), partonsInfo->GetPartonEta6());
1146  fPhiJetCorr6->Fill(jet1->Phi(), partonsInfo->GetPartonPhi6());
1147  if(dRp1 < fRMatching) {
1148  fShapesVar[0] = partonsInfo->GetPartonFlag6();
1149  fPtJetCorr ->Fill(partonsInfo->GetPartonPt6(), jet1->Pt());
1150  }
1151  else {
1152  jp1=RelativePhi(jet1->Phi(),partonsInfo->GetPartonPhi7());
1153  detap1=(jet1->Eta())-(partonsInfo->GetPartonEta7());
1154  dRp1 = TMath::Sqrt(jp1 * jp1 + detap1 * detap1);
1155  fEtaJetCorr7->Fill(jet1->Eta(), partonsInfo->GetPartonEta7());
1156  fPhiJetCorr7->Fill(jet1->Phi(), partonsInfo->GetPartonPhi7());
1157  if(dRp1 < fRMatching) {
1158  fShapesVar[0] = partonsInfo->GetPartonFlag7();
1159  fPtJetCorr->Fill(partonsInfo->GetPartonPt7(), jet1->Pt());
1160  }
1161  else fShapesVar[0]=0;
1162  }
1163  }
1164 
1165  Double_t ptSubtracted = 0;
1166  if (fJetShapeSub==kConstSub) ptSubtracted= jet1->Pt();
1167 
1168  else if (fJetShapeSub==kDerivSub) {
1169  ptSubtracted=jet1->Pt()-GetRhoVal(0)*jet1->Area();
1170  }
1171 
1172  else if (fJetShapeSub==kNoSub) {
1173  if ((fJetShapeType==kData) || (fJetShapeType==kDetEmbPartPythia)) ptSubtracted=jet1->Pt()-GetRhoVal(0)*jet1->Area();
1174  else if ((fJetShapeType==kPythiaDef) || (fJetShapeType==kMCTrue) || (fJetShapeType==kGenOnTheFly)) ptSubtracted= jet1->Pt();
1175  }
1176 
1177  //Printf("ptSubtracted=%f,fPtThreshold =%f ", ptSubtracted, fPtThreshold);
1178  if (ptSubtracted < fPtThreshold) continue;
1179 
1180  if (fOneConstSelectOn == kTRUE) fNbOfConstvspT->Fill(GetJetNumberOfConstituents(jet1,0), ptSubtracted);
1181  if ((fCentSelectOn == kFALSE) && (jet1->GetNumberOfTracks() <= 1)) continue;
1182 
1183  fShapesVar[1] = ptSubtracted;
1184  fShapesVar[2] = GetJetpTD(jet1,0);
1185  fShapesVar[3] = GetJetMass(jet1,0);
1186  fShapesVar[4] = GetJetAngularity(jet1,0);
1187  fShapesVar[5] = GetJetCircularity(jet1,0);
1188  fShapesVar[6] = GetJetLeSub(jet1,0);
1189  fShapesVar[7] = GetJetCoronna(jet1,0);
1190 
1191  Float_t ptMatch=0., ptDMatch=0., massMatch=0., angulMatch=0.,circMatch=0., lesubMatch=0., coronnaMatch=0;
1192  //Float constMatch=0., sigma2Match=0.;
1193  Int_t kMatched = 0;
1194 
1195  if (fJetShapeType==kPythiaDef) {
1196  kMatched =1;
1197  if(fJetShapeSub==kConstSub) kMatched = 3;
1198 
1199  if (!jet3) {
1200  Printf("jet3 does not exist, returning");
1201  continue;
1202  }
1203 
1204  ptMatch=jet3->Pt();
1205  ptDMatch=GetJetpTD(jet3, kMatched);
1206  massMatch=GetJetMass(jet3,kMatched);
1207  //constMatch=1.*GetJetNumberOfConstituents(jet2,kMatched);
1208  angulMatch=GetJetAngularity(jet3, kMatched);
1209  circMatch=GetJetCircularity(jet3, kMatched);
1210  lesubMatch=GetJetLeSub(jet3, kMatched);
1211  coronnaMatch=GetJetCoronna(jet3,kMatched);
1212  //sigma2Match = GetSigma2(jet2, kMatched);
1213  }
1214 
1216  if(fJetShapeSub==kConstSub) kMatched = 3;
1217  if(fJetShapeSub==kDerivSub) kMatched = 2;
1218  ptMatch=jet3->Pt();
1219  ptDMatch=GetJetpTD(jet3, kMatched);
1220  massMatch=GetJetMass(jet3,kMatched);
1221  // constMatch=1.*GetJetNumberOfConstituents(jet3,kMatched);
1222  angulMatch=GetJetAngularity(jet3, kMatched);
1223  circMatch=GetJetCircularity(jet3, kMatched);
1224  lesubMatch=GetJetLeSub(jet3, kMatched);
1225  coronnaMatch = GetJetCoronna(jet3, kMatched);
1226  }
1227 
1228 
1229 
1231  kMatched = 0;
1232  ptMatch=0.;
1233  ptDMatch=0.;
1234  massMatch=0.;
1235  //constMatch=0.;
1236  angulMatch=0.;
1237  circMatch=0.;
1238  lesubMatch=0.;
1239  coronnaMatch =0.;
1240 
1241  }
1242 
1243  fShapesVar[8] = ptMatch;
1244  fShapesVar[9] = ptDMatch;
1245  fShapesVar[10] = massMatch;
1246  fShapesVar[11] = angulMatch;
1247  fShapesVar[12] = circMatch;
1248  fShapesVar[13] = lesubMatch;
1249  fShapesVar[14] = coronnaMatch;
1250  fShapesVar[15] = kWeight;
1251  fShapesVar[16] = rhoVal;
1252  fShapesVar[17] = rhoMassVal;
1253  fShapesVar[18] = jet1->Pt();
1254 
1255  Int_t nInclusiveElectrons = 0, nPhotonicElectrons = 0, nTrueElectronsMC= 0, nTrueHFElecMC= 0;
1256  Double_t pElec = 0., ptElec = 0.;
1257  Bool_t hasElectrons = kFALSE;
1258 
1259  GetNumberOfElectrons(jet1, 0,nMotherKink,listofmotherkink,nInclusiveElectrons,nPhotonicElectrons,pElec,ptElec,hasElectrons);
1260 
1261  fAngChargPart->Fill(jet1->Pt(),GetJetAngularity(jet1,0));
1262  fDispChargPart->Fill(jet1->Pt(),GetJetpTD(jet1,0));
1263 
1264  if(nInclusiveElectrons==1){
1265  fAngIncElec->Fill(jet1->Pt(),GetJetAngularity(jet1,0));
1266  fDispIncElec->Fill(jet1->Pt(),GetJetpTD(jet1,0));
1267  }
1268 
1269  if(!hasElectrons){
1270  fAngHadron->Fill(jet1->Pt(),GetJetAngularity(jet1,0));
1271  fDispHadron->Fill(jet1->Pt(),GetJetpTD(jet1,0));
1272  }
1273 
1274 
1275 
1276  // generated HFE jets
1277 
1278  AliVParticle *vp1 = 0x0;
1279  Int_t pdgCode = 0;
1280  Int_t elecCounter = 0;
1281  Double_t ptTrueHFE = -1.;
1282 
1283  if(fMCarray){
1284 
1285  GetNumberOfTrueElectrons(jet1, 0,nMotherKink,listofmotherkink,nTrueElectronsMC,nTrueHFElecMC, ptTrueHFE);
1286 
1287  for(UInt_t i = 0; i < jet1->GetNumberOfTracks(); i++) {
1288  vp1 = static_cast<AliVParticle*>(jet1->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
1289 
1290  if (!vp1){
1291  Printf("AliVParticle associated to constituent not found");
1292  continue;
1293  }
1294 
1295  pdgCode = vp1->PdgCode();
1296  if (TMath::Abs(pdgCode)==11) elecCounter++;
1297  }
1298  if (elecCounter==1) fPtGenJet->Fill(jet1->Pt());
1299  }
1300 
1301  fShapesVar[19] = ptTrueHFE;
1302  fShapesVar[20] = ptElec;
1303  fShapesVar[21] = nInclusiveElectrons;
1304  fShapesVar[22] = nPhotonicElectrons;
1305  fShapesVar[23] = hasElectrons;
1306  fShapesVar[24] = nTrueElectronsMC;
1307  fShapesVar[25] = nTrueHFElecMC;
1308 
1309  fTreeObservableTagging->Fill();
1310 
1311 
1312  }//jet loop
1313  }
1314  return kTRUE;
1315 }
1316 
1317 //________________________________________________________________________
1319  //calc subtracted jet mass
1320  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
1321  if (fDerivSubtrOrder == 1) return jet->GetShapeProperties()->GetFirstOrderSubtracted();
1322  else return jet->GetShapeProperties()->GetSecondOrderSubtracted();
1323  else
1324  return jet->M();
1325 }
1326 
1327 //________________________________________________________________________
1328 void AliAnalysisTaskEmcalHfeTagging::GetNumberOfElectrons(AliEmcalJet *jet, Int_t jetContNb , Int_t nMother, Double_t listMother[] , Int_t &nIncElec, Int_t &nPhotElec, Double_t &pElec, Double_t &ptElec, Bool_t &hasElec){
1329  // count the number of inclusive electrons per jet and per event
1330 
1331  AliVParticle *vp1 = 0x0;
1332  Int_t nIE=0, nPairs=0, nPE=0, sub = -1;
1333  Float_t ratioElec = -1.;
1334  Double_t pte=0., pe=0.;
1335  Bool_t hasElecCand = kFALSE;
1336 
1337  Double_t jetPt = jet->Pt();
1338 
1339  Double_t ptRange[34] = {0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2, 1.3, 1.4,
1340  1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.5, 4., 5.,
1341  6., 8., 10., 12., 14., 16., 19., 22., 26., 30.,
1342  35., 40., 45., 50.};
1343  Double_t ptJetRange[6] = {5,20,40,60,80,120};
1344 
1345  AliJetContainer *jetCont = GetJetContainer(jetContNb);
1346  if (jet->GetNumberOfTracks()){
1347 
1348  fnPartPerJet->Fill(jet->GetNumberOfTracks());
1349  fpidResponse = fInputHandler->GetPIDResponse();
1350  if(!fpidResponse){
1351  printf("ERROR: pid response not available\n");
1352  }
1353 
1354  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
1355  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
1356 
1357  if (!vp1){
1358  Printf("AliVParticle associated to constituent not found");
1359  continue;
1360  }
1361 
1362  AliVTrack *vtrack = dynamic_cast<AliVTrack*>(vp1);
1363  if (!vtrack) {
1364  printf("ERROR: Could not receive track%d\n", i);
1365  continue;
1366  }
1367 
1368  AliAODTrack *track = dynamic_cast<AliAODTrack*>(vtrack);
1369  if (!track) continue;
1370 
1371  // track cuts for electron identification
1372  Bool_t passTrackCut = kFALSE;
1373  passTrackCut = InclElecTrackCuts(pVtx,track,nMother,listMother);
1374  if (!passTrackCut) continue;
1375 
1376  Double_t p=-9., pt=-9., fTPCnSigma=-99., fTOFnSigma=-99., phi = -9., eta = -99.;
1377  p = track->P();
1378  pt = track->Pt();
1379  phi = track->Phi();
1380  eta = track->Eta();
1381 
1382  fPtTrack->Fill(pt);
1383  fPhiTrack->Fill(pt,phi);
1384  fEtaTrack->Fill(pt,eta);
1385  if (pt > 0.5) fEtaPhiTrack->Fill(phi,eta);
1386 
1387  nPairs=0;
1388 
1389  fTPCnSigma = fpidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
1390  fTOFnSigma = fpidResponse->NumberOfSigmasTOF(track, AliPID::kElectron);
1391 
1392  fnTPCnTOFnocut->Fill(fTPCnSigma,fTOFnSigma);
1393  fnTPCnocutP->Fill(p,fTPCnSigma);
1394  fnTOFnocutP->Fill(p,fTOFnSigma);
1395  if (TMath::Abs(fTOFnSigma)<fSigmaTOFcut){
1396  fnTPCcutP->Fill(p,fTPCnSigma);
1397  fnTPCcutPt->Fill(pt,fTPCnSigma);
1398 
1399  for (Int_t l=0;l<5;l++){// pt jet range
1400  for(Int_t k=0;k<18;k++){// pt electron range (up to 4 GeV/c)
1401  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1] && p>=ptRange[k] && p<ptRange[k+1]) fnTPCSigma[l][k]->Fill(fTPCnSigma);
1402  }
1403  }
1404  }
1405 
1406 
1407  if(TMath::Abs(fTPCnSigma)<3.5) hasElecCand = kTRUE;
1408 
1409  if (TMath::Abs(fTOFnSigma)<fSigmaTOFcut) fPtP->Fill(pt,p);
1410 
1411  if (fTPCnSigma>fSigmaTPCcut && fTPCnSigma<3 && TMath::Abs(fTOFnSigma)<fSigmaTOFcut && pt >= 0.5 && pt < 4){
1412  fPhiRecElecTPC->Fill(pt,phi);
1413  fEtaRecElecTPC->Fill(pt,eta);
1414  fEtaPhiRecElecTPC->Fill(phi,eta);
1415 
1416  nIE++;
1417  pte=pt;
1418  pe=p;
1419  nPairs = GetNumberOfPairs(jet,track,pVtx,nMother,listMother);
1420  if (nPairs>0) nPE++;
1421  }
1422 
1423  // Electron ID with EMCal
1424 
1425  Double_t clsE = -9., m02 = -9., m20 = -9., clsTime = -9, EovP = -9.;
1426 
1427  Double_t emcphimim = 1.39;
1428  Double_t emcphimax = 3.265;
1429 
1430  Int_t clsId = track->GetEMCALcluster();
1431 
1432  AliVCluster *cluster=0x0;
1433 
1434  if (clsId>=0){
1435 
1436  cluster = (AliVCluster*)fVevent->GetCaloCluster(clsId);
1437 
1438  if(cluster && cluster->IsEMCAL()){
1439 
1440  Float_t emcx[3]; // cluster pos
1441  cluster->GetPosition(emcx);
1442  TVector3 clustpos(emcx[0],emcx[1],emcx[2]);
1443  Double_t emcphi = clustpos.Phi();
1444 
1445  if(emcphi < 0) emcphi = emcphi+(2*TMath::Pi());
1446 
1447  if (emcphi>emcphimim && emcphi<emcphimax){
1448  clsE = cluster->E();
1449  m20 = cluster->GetM20();
1450  m02 = cluster->GetM02();
1451  clsTime = cluster->GetTOF()*1e+9; // ns
1452 
1453  EovP = clsE/p;
1454  }
1455  }
1456  }
1457 
1458  fnEovPelecNoTPCcut->Fill(pt,EovP);
1459 
1460 
1461  if (fTPCnSigma>-1.5 && fTPCnSigma<3) fnEovPelecTPCcut->Fill(pt,EovP);
1462 
1463  if (fTPCnSigma>-1.5 && fTPCnSigma<3 && m20 > 0.01 && m20 < 0.35){
1464  fnEovPelecTPCEMCalcut->Fill(pt,EovP);
1465 
1466  for (Int_t l=0;l<5;l++){// pt jet range
1467  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1]) fnEovPelecTPCsscut[l]->Fill(pt,EovP);
1468  }
1469  }
1470 
1471  if (fTPCnSigma < -3.5) fnEovPbackg->Fill(pt,EovP);
1472  fnClsE->Fill(pt,clsE);
1473  fnM20->Fill(pt,m20);
1474  fnM02->Fill(pt,m02);
1475  fnClsTime->Fill(pt,clsTime);
1476 
1477  if (fTPCnSigma>-1.5 && fTPCnSigma<3 && EovP>0.9 && EovP<1.3 && m20 > 0.01 && m20 < 0.35 && pt >= 4 && pt < 50){
1478  fPhiRecElecEMCal->Fill(pt,phi);
1479  fEtaRecElecEMCal->Fill(pt,eta);
1480  fEtaPhiRecElecEMCal->Fill(phi,eta);
1481 
1482  nIE++;
1483  pte=pt;
1484  pe=p;
1485  nPairs = GetNumberOfPairs(jet,track,pVtx,nMother,listMother);
1486  if (nPairs>0) nPE++;
1487  }
1488  }//tagged track
1489  }
1490 
1491  sub = nIE - nPE;
1492  if (jet->GetNumberOfTracks()>0) ratioElec = nIE/jet->GetNumberOfTracks();
1493  fnInclElecPerJet->Fill(nIE);
1494  fnPhotElecPerJet->Fill(nPE);
1495  fnIncSubPhotElecPerJet->Fill(sub);
1496  fnElecOverPartPerJet->Fill(ratioElec);
1497 
1498  nIncElec = nIE;
1499  nPhotElec = nPE;
1500  pElec = pe;//to be used for jets with only one IE
1501  ptElec = pte;//to be used for jets with only one IE
1502  hasElec = hasElecCand;
1503 }
1504 //________________________________________________________________________
1505 void AliAnalysisTaskEmcalHfeTagging::GetNumberOfTrueElectrons(AliEmcalJet *jet, Int_t jetContNb , Int_t nMother, Double_t listMother[] , Int_t &nTrueElec, Int_t &nTrueHFElec, Double_t &ptTrueHFElec){
1506  // count the number of inclusive and HF electrons per jet and per event (true MC)
1507 
1508  AliVParticle *vp1 = 0x0;
1509  Int_t nIE=0, nHFE=0, nPE=0, nPairs=0, iDecay = 0, nDmeson = 0, nBmeson = 0, nElecFromB = 0, nElecFromD = 0, nElecFromDfromB = 0, nQuark = 0, nGluon = 0, nBeauty = 0, nCharm = 0;
1510  Double_t p=-9., pt=-9., fTPCnSigma=-99., fTOFnSigma=-99., MCweight = 1., eta = -99., phi = -99., pte=0.;
1511 
1512  Double_t ptJetRange[6] = {5.,20.,40.,60.,80.,120.};
1513  Double_t angRange[6] = {0.,0.02,0.04,0.06,0.08,0.12};
1514  Double_t dispRange[7] = {0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 1.0};
1515 
1516  Bool_t isFromHFdecay=kFALSE;
1517  Bool_t isFromLMdecay=kFALSE;
1518  Bool_t passTrackCut = kFALSE;
1519 
1520  Double_t jetPt = jet->Pt();
1521  Double_t ang = GetJetAngularity(jet,0);
1522  Double_t disp = GetJetpTD(jet,0);
1523 
1524  AliJetContainer *jetCont = GetJetContainer(jetContNb);
1525  if (jet->GetNumberOfTracks()){
1526 
1527  fpidResponse = fInputHandler->GetPIDResponse();
1528  if(!fpidResponse){
1529  printf("ERROR: pid response not available\n");
1530  }
1531 
1532  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
1533  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
1534 
1535  if (!vp1){
1536  Printf("AliVParticle associated to constituent not found");
1537  continue;
1538  }
1539 
1540  AliVTrack *vtrack = dynamic_cast<AliVTrack*>(vp1);
1541  if (!vtrack) {
1542  printf("ERROR: Could not receive track%d\n", i);
1543  continue;
1544  }
1545 
1546  AliAODTrack *track = dynamic_cast<AliAODTrack*>(vtrack);
1547  if (!track) continue;
1548 
1549  passTrackCut = kFALSE;
1550  isFromHFdecay=kFALSE;
1551  isFromLMdecay=kFALSE;
1552 
1553  p = track->P();
1554  pt = track->Pt();
1555  eta = track->Eta();
1556  phi = track->Phi();
1557 
1558  fTPCnSigma = fpidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
1559  fTOFnSigma = fpidResponse->NumberOfSigmasTOF(track, AliPID::kElectron);
1560 
1561  nPairs=0;
1562  MCweight = 1.;
1563  iDecay = 0;
1564 
1565  // Electron ID with EMCal
1566 
1567  Double_t EovP = -9., clsE = -9., m20 = -9.;
1568 
1569  Double_t emcphimim = 1.39;
1570  Double_t emcphimax = 3.265;
1571 
1572  Int_t clsId = track->GetEMCALcluster();
1573  if (clsId>0){
1574  AliVCluster *cluster=0x0;
1575  cluster = (AliVCluster*)fVevent->GetCaloCluster(clsId);
1576 
1577  if(cluster && cluster->IsEMCAL() && phi > emcphimim && phi < emcphimax){
1578  clsE = cluster->E();
1579  m20 = cluster->GetM20();
1580  }
1581  }
1582 
1583  EovP = clsE/p;
1584 
1585  if(fMCarray){
1586  Int_t label = track->GetLabel();
1587  if(label!=0){
1588  fMCparticle = (AliAODMCParticle*) fMCarray->At(TMath::Abs(track->GetLabel()));
1589  if(fMCparticle){
1590  Int_t partPDG = TMath::Abs(fMCparticle->GetPdgCode());
1591 
1592  if (((partPDG/100)%10) == 4) nDmeson++;
1593  if (((partPDG/100)%10) == 5) nBmeson++;
1594 
1595 
1596  Int_t idMother = fMCparticle->GetMother();
1597 
1598  if (idMother>0){
1599  AliAODMCParticle* mother = (AliAODMCParticle*) fMCarray->At(idMother);
1600  Int_t motherPDG = TMath::Abs(mother->GetPdgCode());
1601 
1602  if ((partPDG==11) && ((motherPDG/100)%10) == 5) nElecFromB++;
1603 
1604  if (((motherPDG/100)%10) == 4) nDmeson++;
1605  if (((motherPDG/100)%10) == 5) nBmeson++;
1606 
1607 
1608  Int_t idSecondMother = mother->GetMother();
1609  if (idSecondMother>0){
1610  AliAODMCParticle* secondMother = (AliAODMCParticle*) fMCarray->At(idSecondMother);
1611  Int_t secondMotherPDG = TMath::Abs(secondMother->GetPdgCode());
1612 
1613  if ((partPDG==11) && (((secondMotherPDG/100)%10) != 5) && (((motherPDG/100)%10) == 4))
1614  nElecFromD++;
1615  if ((partPDG==11) && (((secondMotherPDG/100)%10) == 5) && (((motherPDG/100)%10) == 4))
1616  nElecFromDfromB++;
1617 
1618  if (((secondMotherPDG/100)%10) == 4) nDmeson++;
1619  if (((secondMotherPDG/100)%10) == 5) nBmeson++;
1620 
1621 
1622  Int_t idThirdMother = secondMother->GetMother();
1623  if (idThirdMother>0){
1624  AliAODMCParticle* thirdMother = (AliAODMCParticle*) fMCarray->At(idThirdMother);
1625  Int_t thirdMotherPDG = TMath::Abs(thirdMother->GetPdgCode());
1626 
1627  if ((partPDG==11) && (((motherPDG/100)%10) == 5) && secondMotherPDG == 5){
1628  if (thirdMotherPDG < 9) nQuark++;
1629  if (thirdMotherPDG == 21) nGluon++;
1630  }
1631 
1632  if ((partPDG==11) && (((motherPDG/100)%10) == 4) && secondMotherPDG == 4){
1633  if (thirdMotherPDG < 9) nQuark++;
1634  if (thirdMotherPDG == 21) nGluon++;
1635  }
1636 
1637  if (((thirdMotherPDG/100)%10) == 4) nDmeson++;
1638  if (((thirdMotherPDG/100)%10) == 5) nBmeson++;
1639 
1640 
1641  Int_t idForthMother = thirdMother->GetMother();
1642  if (idForthMother>0){
1643  AliAODMCParticle* forthMother = (AliAODMCParticle*) fMCarray->At(idForthMother);
1644  Int_t forthMotherPDG = TMath::Abs(forthMother->GetPdgCode());
1645 
1646  if ((partPDG==11) && (((secondMotherPDG/100)%10) == 5) && (((motherPDG/100)%10) == 4) && thirdMotherPDG == 5) {
1647  if (forthMotherPDG < 9) nQuark++;
1648  if (forthMotherPDG == 21) nGluon++;
1649  }
1650 
1651  if (((forthMotherPDG/100)%10) == 4) nDmeson++;
1652  if (((forthMotherPDG/100)%10) == 5) nBmeson++;
1653 
1654  if (thirdMotherPDG == 4 || forthMotherPDG == 4) nCharm++;
1655  if (thirdMotherPDG == 5 || forthMotherPDG == 5) nBeauty++;
1656 
1657  }
1658  }//3rd mother
1659  }//2nd mother
1660  }//1st mother
1661 
1662  GetWeightAndDecay(fMCparticle,iDecay,MCweight);
1663  isFromHFdecay = IsFromHFdecay(fMCparticle);
1664  isFromLMdecay = IsFromLMdecay(fMCparticle);
1665 
1666  if ((partPDG==11) && isFromHFdecay){
1667  fptTrueHFEeffTPCTOF[0]->Fill(pt);
1668  fptTrueHFEeffEMCal[0]->Fill(pt);
1669 
1670  for (Int_t l=0;l<5;l++){// pt jet range
1671  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1]){
1672  fptTrueHFEeffTPCTOFwJetPt[0][l]->Fill(pt);
1673  fptTrueHFEeffEMCalwJetPt[0][l]->Fill(pt);
1674  }
1675  }
1676  }
1677 
1678  // track cuts
1679  passTrackCut = InclElecTrackCuts(pVtx,track,nMother,listMother);
1680  if (!passTrackCut) continue;
1681 
1682  if ((partPDG==11) && isFromHFdecay){
1683  fptTrueHFEeffTPCTOF[1]->Fill(pt);
1684  fptTrueHFEeffEMCal[1]->Fill(pt);
1685 
1686  if (fTPCnSigma>fSigmaTPCcut && fTPCnSigma<3)
1687  fptTrueHFEeffTPCTOF[2]->Fill(pt);
1688 
1689  if (TMath::Abs(fTOFnSigma)<fSigmaTOFcut)
1690  fptTrueHFEeffTPCTOF[3]->Fill(pt);
1691 
1692  if (fTPCnSigma>fSigmaTPCcut && fTPCnSigma<3 && TMath::Abs(fTOFnSigma)<fSigmaTOFcut)
1693  fptTrueHFEeffTPCTOF[4]->Fill(pt);
1694 
1695  if (fTPCnSigma>-1.5 && fTPCnSigma<3 && m20 > 0.01 && m20 < 0.35)
1696  fptTrueHFEeffEMCal[2]->Fill(pt);
1697 
1698  if (fTPCnSigma>-1.5 && fTPCnSigma<3 && EovP>0.9 && EovP<1.3)
1699  fptTrueHFEeffEMCal[3]->Fill(pt);
1700 
1701  if (fTPCnSigma>-1.5 && fTPCnSigma<3 && EovP>0.9 && EovP<1.3 && m20 > 0.01 && m20 < 0.35){
1702 
1703  fptTrueHFEeffEMCal[4]->Fill(pt);
1704 
1705  for (Int_t l=0;l<5;l++){// pt jet range
1706  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1]){
1707  fptTrueHFEeffTPCTOFwJetPt[1][l]->Fill(pt);
1708  fptTrueHFEeffEMCalwJetPt[1][l]->Fill(pt);
1709  }
1710  }
1711  }
1712  }
1713 
1714  if (partPDG == 11){
1715  nIE++;
1716  pte=pt;
1717 
1718  if (isFromHFdecay) nHFE++;
1719  if (iDecay>0 && iDecay<7) nPE++;
1720  }
1721 
1722 
1723  // TPC-TOF
1724  if (fTPCnSigma>fSigmaTPCcut && fTPCnSigma<3 && TMath::Abs(fTOFnSigma)<fSigmaTOFcut && pt >= 0.5 && pt < 4){
1725  nPairs = GetNumberOfPairs(jet,track,pVtx,nMother,listMother);
1726  if (nPairs>0 && iDecay>0 && iDecay<7) fptRecPE->Fill(pt);
1727  if (nPairs>0 && iDecay==0) fptWrongPE->Fill(pt);
1728  if (iDecay>0 && iDecay<7) fptTruePE->Fill(pt);
1729 
1730  for (Int_t l=0;l<5;l++){// pt jet range
1731  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1]){
1732  if (iDecay>0 && iDecay<7) fTotPEJetPt[l]->Fill(pt,MCweight);
1733  if (nPairs>0 && iDecay>0 && iDecay<7) fRecPEJetPt[l]->Fill(pt,MCweight);
1734  }
1735  }
1736 
1737  for (Int_t l=0;l<5;l++){// angularity
1738  if (ang>=angRange[l] && ang<angRange[l+1]){
1739  if (iDecay>0 && iDecay<7) fTotPEAng[l]->Fill(pt,MCweight);
1740  if (nPairs>0 && iDecay>0 && iDecay<7) fRecPEAng[l]->Fill(pt,MCweight);
1741  }
1742  }
1743 
1744  for (Int_t l=0;l<6;l++){// dispersion
1745  if (disp>=dispRange[l] && disp<dispRange[l+1]){
1746  if (iDecay>0 && iDecay<7) fTotPEDisp[l]->Fill(pt,MCweight);
1747  if (nPairs>0 && iDecay>0 && iDecay<7) fRecPEDisp[l]->Fill(pt,MCweight);
1748  }
1749  }
1750  }// PID cuts
1751 
1752  //EMCal
1753  if (fTPCnSigma>-1.5 && fTPCnSigma<3 && EovP>0.9 && EovP<1.3 && m20 > 0.01 && m20 < 0.35 && pt >= 4 && pt < 50 ){
1754  nPairs = GetNumberOfPairs(jet,track,pVtx,nMother,listMother);
1755  if (nPairs>0 && iDecay>0 && iDecay<7) fptRecPE->Fill(pt);
1756  if (nPairs>0 && iDecay==0) fptWrongPE->Fill(pt);
1757  if (iDecay>0 && iDecay<7) fptTruePE->Fill(pt);
1758 
1759  for (Int_t l=0;l<5;l++){// pt jet range
1760  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1]){
1761  if (iDecay>0 && iDecay<7) fTotPEJetPt[l]->Fill(pt,MCweight);
1762  if (nPairs>0 && iDecay>0 && iDecay<7) fRecPEJetPt[l]->Fill(pt,MCweight);
1763  }
1764  }
1765 
1766  for (Int_t l=0;l<5;l++){// angularity
1767  if (ang>=angRange[l] && ang<angRange[l+1]){
1768  if (iDecay>0 && iDecay<7) fTotPEAng[l]->Fill(pt,MCweight);
1769  if (nPairs>0 && iDecay>0 && iDecay<7) fRecPEAng[l]->Fill(pt,MCweight);
1770  }
1771  }
1772 
1773  for (Int_t l=0;l<6;l++){// dispersion
1774  if (disp>=dispRange[l] && disp<dispRange[l+1]){
1775  if (iDecay>0 && iDecay<7) fTotPEDisp[l]->Fill(pt,MCweight);
1776  if (nPairs>0 && iDecay>0 && iDecay<7) fRecPEDisp[l]->Fill(pt,MCweight);
1777  }
1778  }
1779  }// PID cuts
1780  }
1781  }
1782  }
1783  }//track loop
1784  }
1785 
1786  if (nIE==1) fptJetIE->Fill(jet->Pt());
1787  if (nIE==1 && nPE==1) fptJetPE->Fill(jet->Pt());
1788  if (nIE==1 && nHFE==1) fptJetHFE->Fill(jet->Pt());
1789 
1790  fnTrueElecPerJet->Fill(nIE);
1791  fnTrueHFElecPerJet->Fill(nHFE);
1792  fnTruePElecPerJet->Fill(nPE);
1793 
1794  nTrueElec = nIE;
1795  nTrueHFElec = nHFE;
1796 
1797  if (nIE==1 && nHFE==1) ptTrueHFElec = pte;
1798 
1799  if (nCharm > 1){
1800  fAngCharm->Fill(jet->Pt(),GetJetAngularity(jet,0));
1801  fDispCharm->Fill(jet->Pt(),GetJetpTD(jet,0));
1802  }
1803 
1804  if (nBeauty > 1){
1805  fAngBeauty->Fill(jet->Pt(),GetJetAngularity(jet,0));
1806  fDispBeauty->Fill(jet->Pt(),GetJetpTD(jet,0));
1807  }
1808 
1809  if (nDmeson > 1){
1810  fAngD->Fill(jet->Pt(),GetJetAngularity(jet,0));
1811  fDispD->Fill(jet->Pt(),GetJetpTD(jet,0));
1812  }
1813 
1814  if (nBmeson > 1){
1815  fAngB->Fill(jet->Pt(),GetJetAngularity(jet,0));
1816  fDispB->Fill(jet->Pt(),GetJetpTD(jet,0));
1817  }
1818 
1819  if (nElecFromB == 1){
1820  fAngElecFromB->Fill(jet->Pt(),GetJetAngularity(jet,0));
1821  fDispElecFromB->Fill(jet->Pt(),GetJetpTD(jet,0));
1822  }
1823 
1824  if (nElecFromD == 1){
1825  fAngElecFromD->Fill(jet->Pt(),GetJetAngularity(jet,0));
1826  fDispElecFromD->Fill(jet->Pt(),GetJetpTD(jet,0));
1827  }
1828 
1829  if (nElecFromDfromB == 1){
1830  fAngElecFromDFromB->Fill(jet->Pt(),GetJetAngularity(jet,0));
1831  fDispElecFromDFromB->Fill(jet->Pt(),GetJetpTD(jet,0));
1832  }
1833 
1834  if (nQuark == 1){
1835  fAngQuark->Fill(jet->Pt(),GetJetAngularity(jet,0));
1836  fDispQuark->Fill(jet->Pt(),GetJetpTD(jet,0));
1837  }
1838  if (nGluon == 1){
1839  fAngGluon->Fill(jet->Pt(),GetJetAngularity(jet,0));
1840  fDispGluon->Fill(jet->Pt(),GetJetpTD(jet,0));
1841  }
1842 }
1843 
1844 //_________________________________________
1845 Int_t AliAnalysisTaskEmcalHfeTagging::GetNumberOfPairs(AliEmcalJet *jet, AliAODTrack *track,const AliVVertex *pVtx, Int_t nMother, Double_t listMother[])
1846 {
1847  //Get number of ULS and LS pairs per event
1848 
1849  Int_t ntracks = 0;
1850  ntracks = fVevent->GetNumberOfTracks();
1851 
1852  Int_t nULSpairs = 0;
1853  Int_t nLSpairs = 0;
1854  Int_t sub = -1;
1855 
1856  Double_t ptJetRange[6] = {5,20,40,60,80,120};
1857  Double_t jetPt = jet->Pt();
1858 
1859  for(Int_t jTracks = 0; jTracks < ntracks; jTracks++){
1860 
1861  AliVParticle* Vassotrack = fVevent->GetTrack(jTracks);
1862 
1863  if (!Vassotrack) {
1864  printf("ERROR: Could not receive associated track %d\n", jTracks);
1865  continue;
1866  }
1867  AliAODTrack *trackAsso = dynamic_cast<AliAODTrack*>(Vassotrack);
1868  if(!trackAsso) continue;
1869 
1870  if((track->GetID())==(trackAsso->GetID())) continue;
1871 
1872  // track cuts
1873  Bool_t passAssoTrackCut = kFALSE;
1874  passAssoTrackCut = PhotElecTrackCuts(pVtx,trackAsso,nMother,listMother);
1875  if (!passAssoTrackCut) continue;
1876 
1877  // tagged particle
1878  Double_t p=-9.,pt=-9.,eta =-9.,phi=-9.;
1879  Int_t charge = 0;
1880 
1881  pt = track->Pt();
1882  p = track->P();
1883  phi = track->Phi();
1884  eta = track->Eta();
1885  charge = track->Charge();
1886 
1887 
1888  // associated particle variables
1889  Double_t pAsso=-9.,ptAsso=-9.,etaAsso =-9.,fITSnSigmaAsso=-9.,fTPCnSigmaAsso=-9.,phiAsso=-9.;
1890  Int_t chargeAsso = 0;
1891 
1892  ptAsso = trackAsso->Pt();
1893  pAsso = trackAsso->P();
1894  phiAsso = trackAsso->Phi();
1895  etaAsso = trackAsso->Eta();
1896  chargeAsso = trackAsso->Charge();
1897 
1898 
1899  // looser PID cuts
1900  fITSnSigmaAsso = fpidResponse->NumberOfSigmasITS(trackAsso, AliPID::kElectron);
1901  fTPCnSigmaAsso = fpidResponse->NumberOfSigmasTPC(trackAsso, AliPID::kElectron);
1902 
1903  if(TMath::Abs(fTPCnSigmaAsso)>3.5) continue;
1904 
1905  // invariant mass
1906  Bool_t fFlagLS=kFALSE, fFlagULS=kFALSE;
1907  Double_t openingAngle = -999., mass=999., width = -999;
1908 
1909  Int_t fPDGe1 = 11; Int_t fPDGe2 = 11;
1910  if(charge>0) fPDGe1 = -11;
1911  if(chargeAsso>0) fPDGe2 = -11;
1912 
1913  if(charge == chargeAsso) fFlagLS = kTRUE;
1914  if(charge != chargeAsso) fFlagULS = kTRUE;
1915 
1916  AliKFParticle::SetField(fVevent->GetMagneticField());
1917  AliKFParticle ge1(*track, fPDGe1);
1918  AliKFParticle ge2(*trackAsso, fPDGe2);
1919  AliKFParticle recg(ge1, ge2);
1920 
1921  if(recg.GetNDF()<1) continue;
1922  Double_t chi2recg = recg.GetChi2()/recg.GetNDF();
1923  if(TMath::Sqrt(TMath::Abs(chi2recg))>3.) continue;
1924 
1925  openingAngle = ge1.GetAngle(ge2);
1926  //if(openingAngle > fOpeningAngleCut) continue;
1927 
1928  Int_t MassCorrect=-9;
1929  MassCorrect = recg.GetMass(mass,width);
1930 
1931  for (Int_t l=0;l<5;l++){// pt jet range
1932  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1] && fFlagULS) fInvmassULS[l]->Fill(pt,mass);
1933  if (jetPt>=ptJetRange[l] && jetPt<ptJetRange[l+1] && fFlagLS) fInvmassLS[l]->Fill(pt,mass);
1934  }
1935 
1936  if(mass<fIMcut && fFlagULS){
1937  nULSpairs++;
1938  fAngULS->Fill(jet->Pt(),GetJetAngularity(jet,0));
1939  fDispULS->Fill(jet->Pt(),GetJetpTD(jet,0));
1940  }
1941 
1942  if(mass<fIMcut && fFlagLS){
1943  nLSpairs++;
1944  fAngLS->Fill(jet->Pt(),GetJetAngularity(jet,0));
1945  fDispLS->Fill(jet->Pt(),GetJetpTD(jet,0));
1946  }
1947 
1948  }//track loop
1949 
1950  sub = nULSpairs-nLSpairs;
1951  fnULSmLSpairsPerElectron->Fill(track->Pt(),sub);
1952 
1953  if (sub>0){
1954  fAngPhotElec->Fill(jet->Pt(),GetJetAngularity(jet,0));
1955  fDispPhotElec->Fill(jet->Pt(),GetJetpTD(jet,0));
1956  }
1957 
1958  return sub;
1959 }
1960 //_________________________________________
1962 {
1963  // Check if the mother comes from heavy-flavour decays
1964  Bool_t isHFdecay = kFALSE;
1965 
1966  Int_t idMother = particle->GetMother();
1967  if (idMother>0){
1968  AliAODMCParticle* mother = (AliAODMCParticle*) fMCarray->At(idMother);
1969  Int_t motherPDG = TMath::Abs(mother->GetPdgCode());
1970 
1971  // c decays
1972  if( (((motherPDG/1000)%10) == 4) || (((motherPDG/100)%10) == 4) ) isHFdecay = kTRUE;
1973 
1974  // b decays
1975  if( (((motherPDG/1000)%10) == 5) || (((motherPDG/100)%10) == 5) ) isHFdecay = kTRUE;
1976  }
1977 
1978  return isHFdecay;
1979 }
1980 //_________________________________________
1982 {
1983  // Check if mother comes from light-meson decays
1984  Bool_t isLMdecay = kFALSE;
1985 
1986  Int_t idMother = particle->GetMother();
1987  if (idMother>0){
1988  AliAODMCParticle* mother = (AliAODMCParticle*) fMCarray->At(idMother);
1989  Int_t motherPDG = TMath::Abs(mother->GetPdgCode());
1990 
1991  if(motherPDG == 111 || motherPDG == 221 || motherPDG==223 || motherPDG==333 || motherPDG==331 ||
1992  motherPDG==113 || motherPDG==213 || motherPDG==313 || motherPDG==323) isLMdecay = kTRUE;
1993  }
1994 
1995  return isLMdecay;
1996 }
1997 //_________________________________________
1999 {
2000  // Check if particle is primary
2001  Bool_t isprimary = kFALSE;
2002  if (particle->IsPrimary()) isprimary = kTRUE;
2003 
2004  return isprimary;
2005 }
2006 //_________________________________________
2008 {
2009  // Check if the particle has mother
2010  Bool_t hasMother = kTRUE;
2011 
2012  Int_t idMother = particle->GetMother();
2013  if (idMother < 0) hasMother = kFALSE;
2014 
2015  return hasMother;
2016 }
2017 //_________________________________________
2019 {
2020  //Get Pi0 weight
2021  double weight = 1.;
2022  double parPi0_enh[5] = {0,0,0,0,0};
2023 
2024  if (fRunNumber >= 252235 && fRunNumber <= 264347){//LHC17h8b (o/p)
2025  parPi0_enh[0] = 11.4919;
2026  parPi0_enh[1] = 0.021972;
2027  parPi0_enh[2] = 0.133756;
2028  parPi0_enh[3] = 1.67114;
2029  parPi0_enh[4] = 5.80327;
2030  }
2031 
2032  if (fRunNumber >= 256504 && fRunNumber <= 259888){//LHC18f4b (k/l)
2033  parPi0_enh[0] = 17.2349;
2034  parPi0_enh[1] = 0.101117;
2035  parPi0_enh[2] = 0.165525;
2036  parPi0_enh[3] = 1.52289;
2037  parPi0_enh[4] = 5.72523;
2038  }
2039 
2040  if (fMCweight==1) weight = parPi0_enh[0] / TMath::Power((exp(parPi0_enh[1]*mcPi0pT - parPi0_enh[2]*mcPi0pT*mcPi0pT) + (mcPi0pT/parPi0_enh[3])) , parPi0_enh[4]);
2041 
2042  return weight;
2043 }
2044 //_________________________________________
2046 {
2047  //Get Pi0 weight
2048  double weight = 1.;
2049  double parEta_enh[5] = {0,0,0,0,0};
2050 
2051  if (fRunNumber >= 252235 && fRunNumber <= 264347){//LHC17h8b (o/p)
2052  parEta_enh[0] = 3.89483;
2053  parEta_enh[1] = 2.8197;
2054  parEta_enh[2] = 0.0662309;
2055  parEta_enh[3] = 0.00563749;
2056  parEta_enh[4] = 0.435908;
2057  }
2058 
2059  if (fRunNumber >= 256504 && fRunNumber <= 259888){//LHC18f4b (k/l)
2060  parEta_enh[0] = 0.328435;
2061  parEta_enh[1] = -0.450847;
2062  parEta_enh[2] = 0.00972497;
2063  parEta_enh[3] = 3.13337;
2064  parEta_enh[4] = 5.91296;
2065  }
2066 
2067  if (fMCweight==1) weight = parEta_enh[0] / TMath::Power((exp(parEta_enh[1]*mcEtapT - parEta_enh[2]*mcEtapT*mcEtapT) + (mcEtapT/parEta_enh[3])), parEta_enh[4]);
2068 
2069  return weight;
2070 }
2071 //________________________________________________________________________
2072 void AliAnalysisTaskEmcalHfeTagging::GetWeightAndDecay(AliAODMCParticle *particle, Int_t &decay, Double_t &weight)
2073 {
2074  //Get decay channel and weight for pi0 and eta (MC with enhanced signal)
2075  Double_t w = 1.;
2076  Int_t d = 0;
2077  Int_t partPDG = particle->GetPdgCode();
2078 
2079  if (TMath::Abs(partPDG)==11){
2080  Int_t idMother = particle->GetMother();
2081 
2082  if (idMother>0){
2083  AliAODMCParticle* mother = (AliAODMCParticle*) fMCarray->At(idMother);
2084  Int_t motherPDG = mother->GetPdgCode();
2085  Double_t motherPt = mother->Pt();
2086 
2087  Bool_t isMotherPrimary = IsPrimary(mother);
2088  Bool_t isMotherFromHF = IsFromHFdecay(mother);
2089  Bool_t isMotherFromLM = IsFromLMdecay(mother);
2090  Bool_t hasSecondMother = HasMother(mother);
2091 
2092  if (motherPDG==111 && isMotherPrimary){// pi0 -> e
2093  if (!hasSecondMother || (!isMotherFromHF && !isMotherFromLM)){
2094  d = 1;
2095  w = GetPi0weight(motherPt);
2096  }
2097  }
2098 
2099  if (motherPDG==221 && isMotherPrimary){// eta -> e
2100  if (!hasSecondMother || (!isMotherFromHF && !isMotherFromLM)){
2101  d = 2;
2102  w = GetEtaweight(motherPt);
2103  }
2104  }
2105 
2106  Int_t idSecondMother = mother->GetMother();
2107 
2108  if (idSecondMother>0){
2109  AliAODMCParticle* secondMother = (AliAODMCParticle*) fMCarray->At(idSecondMother);
2110  Int_t secondMotherPDG = secondMother->GetPdgCode();
2111  Double_t secondMotherPt = secondMother->Pt();
2112 
2113  Bool_t isSecondMotherPrimary = IsPrimary(secondMother);
2114  Bool_t isSecondMotherFromHF = IsFromHFdecay(secondMother);
2115  Bool_t isSecondMotherFromLM = IsFromLMdecay(secondMother);
2116  Bool_t hasThirdMother = HasMother(secondMother);
2117 
2118  if (motherPDG==22 && secondMotherPDG==111 && isSecondMotherPrimary){ //pi0 -> g -> e
2119  if (!hasThirdMother || (!isSecondMotherFromHF && !isSecondMotherFromLM)){
2120  d = 3;
2121  w = GetPi0weight(secondMotherPt);
2122  }
2123  }
2124 
2125  if (motherPDG==22 && secondMotherPDG==221 && isSecondMotherPrimary){ //eta -> g -> e
2126  if (!hasThirdMother || (!isSecondMotherFromHF && !isSecondMotherFromLM)){
2127  d = 4;
2128  w = GetEtaweight(secondMotherPt);
2129  }
2130  }
2131 
2132  if (motherPDG==111 && secondMotherPDG==221 && isSecondMotherPrimary){ //eta -> pi0 -> e
2133  if (!hasThirdMother || (!isSecondMotherFromHF && !isSecondMotherFromLM)){
2134  d = 5;
2135  w = GetEtaweight(secondMotherPt);
2136  }
2137  }
2138 
2139  Int_t idThirdMother = secondMother->GetMother();
2140 
2141  if (idThirdMother>0){
2142  AliAODMCParticle* thirdMother = (AliAODMCParticle*) fMCarray->At(idThirdMother);
2143  Int_t thirdMotherPDG = thirdMother->GetPdgCode();
2144  Double_t thirdMotherPt = thirdMother->Pt();
2145 
2146  Bool_t isThirdMotherPrimary = IsPrimary(thirdMother);
2147  Bool_t isThirdMotherFromHF = IsFromHFdecay(thirdMother);
2148  Bool_t isThirdMotherFromLM = IsFromLMdecay(thirdMother);
2149  Bool_t hasFourthMother = HasMother(thirdMother);
2150 
2151  if (motherPDG==22 && secondMotherPDG==111 && thirdMotherPDG==221 && isThirdMotherPrimary){//eta->pi0->g-> e
2152  if (!hasFourthMother || (!isThirdMotherFromHF && !isThirdMotherFromLM)){
2153  d = 6;
2154  w = GetEtaweight(thirdMotherPt);
2155  }
2156  }
2157  }//third mother
2158  }//second mother
2159  }//mother
2160  }// if electron
2161  decay = d;
2162  weight = w;
2163 }
2164 //________________________________________________________________________
2166 
2167  AliJetContainer *jetCont = GetJetContainer(jetContNb);
2168  if (!jet->GetNumberOfTracks())
2169  return 0;
2170  Double_t den=0.;
2171  Double_t num = 0.;
2172  AliVParticle *vp1 = 0x0;
2173  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
2174  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
2175 
2176  if (!vp1){
2177  Printf("AliVParticle associated to constituent not found");
2178  continue;
2179  }
2180 
2181  Double_t dphi = RelativePhi(vp1->Phi(),jet->Phi());
2182  Double_t dr2 = (vp1->Eta()-jet->Eta())*(vp1->Eta()-jet->Eta()) + dphi*dphi;
2183  Double_t dr = TMath::Sqrt(dr2);
2184  num=num+vp1->Pt()*dr;
2185  den=den+vp1->Pt();
2186  }
2187  return num/den;
2188 }
2189 
2190 //________________________________________________________________________
2192 
2193  if((fJetShapeSub==kDerivSub) && (jetContNb==0))
2196  else
2197  return Angularity(jet, jetContNb);
2198 
2199 }
2200 
2201 //____________________________________________________________________________
2202 
2204 
2205  AliTrackContainer *PartCont = NULL;
2206  AliParticleContainer *PartContMC = NULL;
2207 
2208 
2209  if (fJetShapeSub==kConstSub){
2211  else PartCont = GetTrackContainer(1);
2212  }
2213  else{
2215  else PartCont = GetTrackContainer(0);
2216  }
2217 
2218  TClonesArray *TracksArray = NULL;
2219  TClonesArray *TracksArrayMC = NULL;
2220 
2221  if (fJetShapeType == AliAnalysisTaskEmcalHfeTagging::kGenOnTheFly) TracksArrayMC = PartContMC->GetArray();
2222  else TracksArray = PartCont->GetArray();
2223 
2225  if(!PartContMC || !TracksArrayMC) return -2;
2226  }
2227  else {
2228  if(!PartCont || !TracksArray) return -2;
2229  }
2230 
2231 
2232  AliAODTrack *Track = 0x0;
2233  Float_t sumpt=0;
2234  Int_t NTracks=0;
2235  if (fJetShapeType == AliAnalysisTaskEmcalHfeTagging::kGenOnTheFly) NTracks = TracksArrayMC->GetEntriesFast();
2236  else NTracks = TracksArray->GetEntriesFast();
2237 
2238  for(Int_t i=0; i < NTracks; i++){
2240  if((Track = static_cast<AliAODTrack*>(PartContMC->GetAcceptParticle(i)))){
2241  if (!Track) continue;
2242  if(TMath::Abs(Track->Eta())>0.9) continue;
2243  Double_t dphi = RelativePhi(Track->Phi(),jet->Phi());
2244  Double_t dr2 = (Track->Eta()-jet->Eta())*(Track->Eta()-jet->Eta()) + dphi*dphi;
2245  Double_t dr = TMath::Sqrt(dr2);
2246  if((dr>=0.8) && (dr<1)) sumpt=sumpt+Track->Pt();
2247  }
2248  }
2249  else{
2250  if((Track = static_cast<AliAODTrack*>(PartCont->GetAcceptTrack(i)))){
2251  if (!Track) continue;
2252  if(TMath::Abs(Track->Eta())>0.9) continue;
2253  if (Track->Pt()<0.15) continue;
2254  Double_t dphi = RelativePhi(Track->Phi(),jet->Phi());
2255  Double_t dr2 = (Track->Eta()-jet->Eta())*(Track->Eta()-jet->Eta()) + dphi*dphi;
2256  Double_t dr = TMath::Sqrt(dr2);
2257  if((dr>=0.8) && (dr<1)) sumpt=sumpt+Track->Pt();
2258 
2259  }
2260  }
2261  }
2262  return sumpt;
2263 }
2264 
2265 //________________________________________________________________________
2267 
2268  if((fJetShapeSub==kDerivSub) && (jetContNb==0)) return -2;
2269  else
2270  return Coronna(jet, jetContNb);
2271 
2272 }
2273 
2274 //________________________________________________________________________
2276 
2277  AliJetContainer *jetCont = GetJetContainer(jetContNb);
2278  if (!jet->GetNumberOfTracks())
2279  return 0;
2280  Double_t den=0.;
2281  Double_t num = 0.;
2282  AliVParticle *vp1 = 0x0;
2283  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
2284  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
2285 
2286  if (!vp1){
2287  Printf("AliVParticle associated to constituent not found");
2288  continue;
2289  }
2290 
2291  num=num+vp1->Pt()*vp1->Pt();
2292  den=den+vp1->Pt();
2293  }
2294  return TMath::Sqrt(num)/den;
2295 }
2296 
2297 //________________________________________________________________________
2299  //calc subtracted jet mass
2300  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
2302  else return jet->GetShapeProperties()->GetSecondOrderSubtractedpTD();
2303  else
2304  return PTD(jet, jetContNb);
2305 
2306 }
2307 
2308 //_____________________________________________________________________________
2310 
2311  AliJetContainer *jetCont = GetJetContainer(jetContNb);
2312  if (!jet->GetNumberOfTracks())
2313  return 0;
2314  Double_t mxx = 0.;
2315  Double_t myy = 0.;
2316  Double_t mxy = 0.;
2317  int nc = 0;
2318  Double_t sump2 = 0.;
2319  Double_t pxjet=jet->Px();
2320  Double_t pyjet=jet->Py();
2321  Double_t pzjet=jet->Pz();
2322 
2323 
2324  //2 general normalized vectors perpendicular to the jet
2325  TVector3 ppJ1(pxjet, pyjet, pzjet);
2326  TVector3 ppJ3(- pxjet* pzjet, - pyjet * pzjet, pxjet * pxjet + pyjet * pyjet);
2327  ppJ3.SetMag(1.);
2328  TVector3 ppJ2(-pyjet, pxjet, 0);
2329  ppJ2.SetMag(1.);
2330  AliVParticle *vp1 = 0x0;
2331  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
2332  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
2333 
2334  if (!vp1){
2335  Printf("AliVParticle associated to constituent not found");
2336  continue;
2337  }
2338 
2339  TVector3 pp(vp1->Px(), vp1->Py(), vp1->Pz());
2340 
2341  //local frame
2342  TVector3 pLong = pp.Dot(ppJ1) / ppJ1.Mag2() * ppJ1;
2343  TVector3 pPerp = pp - pLong;
2344  //projection onto the two perpendicular vectors defined above
2345 
2346  Float_t ppjX = pPerp.Dot(ppJ2);
2347  Float_t ppjY = pPerp.Dot(ppJ3);
2348  Float_t ppjT = TMath::Sqrt(ppjX * ppjX + ppjY * ppjY);
2349  if(ppjT<=0) return 0;
2350 
2351  mxx += (ppjX * ppjX / ppjT);
2352  myy += (ppjY * ppjY / ppjT);
2353  mxy += (ppjX * ppjY / ppjT);
2354  nc++;
2355  sump2 += ppjT;}
2356 
2357  if(nc<2) return 0;
2358  if(sump2==0) return 0;
2359  // Sphericity Matrix
2360  Double_t ele[4] = {mxx / sump2, mxy / sump2, mxy / sump2, myy / sump2};
2361  TMatrixDSym m0(2,ele);
2362 
2363  // Find eigenvectors
2364  TMatrixDSymEigen m(m0);
2365  TVectorD eval(2);
2366  TMatrixD evecm = m.GetEigenVectors();
2367  eval = m.GetEigenValues();
2368  // Largest eigenvector
2369  int jev = 0;
2370  // cout<<eval[0]<<" "<<eval[1]<<endl;
2371  if (eval[0] < eval[1]) jev = 1;
2372  TVectorD evec0(2);
2373  // Principle axis
2374  evec0 = TMatrixDColumn(evecm, jev);
2375  Double_t compx=evec0[0];
2376  Double_t compy=evec0[1];
2377  TVector2 evec(compx, compy);
2378  Double_t circ=0;
2379  if(jev==1) circ=2*eval[0];
2380  if(jev==0) circ=2*eval[1];
2381 
2382  return circ;
2383 
2384 
2385 
2386 }
2387 
2388 //________________________________________________________________________
2390  //calc subtracted jet mass
2391 
2392  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
2395  else
2396  return Circularity(jet, jetContNb);
2397 
2398 }
2399 
2400 //________________________________________________________________________
2402 
2403  AliJetContainer *jetCont = GetJetContainer(jetContNb);
2404  if (!jet->GetNumberOfTracks())
2405  return 0;
2406  Double_t den=0.;
2407  Double_t num = 0.;
2408  AliVParticle *vp1 = 0x0;
2409  AliVParticle *vp2 = 0x0;
2410  std::vector<int> ordindex;
2412  //Printf("Nbof const = %d", jet->GetNumberOfTracks());
2413  //Printf("ordindex[0] = %d, ordindex[1] = %d", ordindex[0], ordindex[1]);
2414 
2415  if(ordindex.size()<2) return -1;
2416 
2417  vp1 = static_cast<AliVParticle*>(jet->TrackAt(ordindex[0], jetCont->GetParticleContainer()->GetArray()));
2418  if (!vp1){
2419  Printf("AliVParticle associated to Leading constituent not found");
2420  return -1;
2421  }
2422 
2423  vp2 = static_cast<AliVParticle*>(jet->TrackAt(ordindex[1], jetCont->GetParticleContainer()->GetArray()));
2424  if (!vp2){
2425  Printf("AliVParticle associated to Subleading constituent not found");
2426  return -1;
2427  }
2428 
2429 
2430  num=vp1->Pt();
2431  den=vp2->Pt();
2432  //Printf("vp1->Pt() =%f, vp2->Pt() =%f", vp1->Pt(), vp2->Pt());
2433 
2434  return num-den;
2435 }
2436 
2437 //________________________________________________________________________
2439  //calc subtracted jet mass
2440 
2441  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
2443  else return jet->GetShapeProperties()->GetSecondOrderSubtractedLeSub();
2444  else
2445  return LeSub(jet, jetContNb);
2446 
2447 }
2448 
2449 //________________________________________________________________________
2451  //calc subtracted jet mass
2452 
2453  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
2456  else
2457  return jet->GetNumberOfTracks();
2458 
2459 }
2460 
2461 
2462 //______________________________________________________________________________
2464 
2465  AliJetContainer *jetCont = GetJetContainer(jetContNb);
2466  if (!jet->GetNumberOfTracks())
2467  return 0;
2468  Double_t mxx = 0.;
2469  Double_t myy = 0.;
2470  Double_t mxy = 0.;
2471  int nc = 0;
2472  Double_t sump2 = 0.;
2473 
2474  AliVParticle *vp1 = 0x0;
2475  for(UInt_t i = 0; i < jet->GetNumberOfTracks(); i++) {
2476  vp1 = static_cast<AliVParticle*>(jet->TrackAt(i, jetCont->GetParticleContainer()->GetArray()));
2477 
2478  if (!vp1){
2479  Printf("AliVParticle associated to constituent not found");
2480  continue;
2481  }
2482 
2483  Double_t ppt=vp1->Pt();
2484  Double_t dphi = RelativePhi(vp1->Phi(),jet->Phi());
2485 
2486  Double_t deta = vp1->Eta()-jet->Eta();
2487  mxx += ppt*ppt*deta*deta;
2488  myy += ppt*ppt*dphi*dphi;
2489  mxy -= ppt*ppt*deta*TMath::Abs(dphi);
2490  nc++;
2491  sump2 += ppt*ppt;
2492 
2493  }
2494  if(nc<2) return 0;
2495  if(sump2==0) return 0;
2496  // Sphericity Matrix
2497  Double_t ele[4] = {mxx , mxy , mxy , myy };
2498  TMatrixDSym m0(2,ele);
2499 
2500  // Find eigenvectors
2501  TMatrixDSymEigen m(m0);
2502  TVectorD eval(2);
2503  TMatrixD evecm = m.GetEigenVectors();
2504  eval = m.GetEigenValues();
2505  // Largest eigenvector
2506  int jev = 0;
2507  // cout<<eval[0]<<" "<<eval[1]<<endl;
2508  if (eval[0] < eval[1]) jev = 1;
2509  TVectorD evec0(2);
2510  // Principle axis
2511  evec0 = TMatrixDColumn(evecm, jev);
2512  Double_t compx=evec0[0];
2513  Double_t compy=evec0[1];
2514  TVector2 evec(compx, compy);
2515  Double_t sig=0;
2516  if(jev==1) sig=TMath::Sqrt(TMath::Abs(eval[0])/sump2);
2517  if(jev==0) sig=TMath::Sqrt(TMath::Abs(eval[1])/sump2);
2518 
2519  return sig;
2520 
2521 }
2522 
2523 //________________________________________________________________________
2525  //calc subtracted jet mass
2526 
2527  if((fJetShapeSub==kDerivSub)&&(jetContNb==0))
2529  else return jet->GetShapeProperties()->GetSecondOrderSubtractedSigma2();
2530  else
2531  return Sigma2(jet, jetContNb);
2532 
2533 }
2534 
2535 //________________________________________________________________________
2537 
2538  AliTrackContainer *PartCont = NULL;
2539  AliParticleContainer *PartContMC = NULL;
2540 
2541 
2542  if (fJetShapeSub==kConstSub){
2544  else PartCont = GetTrackContainer(1);
2545  }
2546  else{
2548  else PartCont = GetTrackContainer(0);
2549  }
2550 
2551  TClonesArray *TracksArray = NULL;
2552  TClonesArray *TracksArrayMC = NULL;
2553 
2554  if (fJetShapeType == AliAnalysisTaskEmcalHfeTagging::kGenOnTheFly) TracksArrayMC = PartContMC->GetArray();
2555  else TracksArray = PartCont->GetArray();
2556 
2558  if(!PartContMC || !TracksArrayMC) return -99999;
2559  }
2560  else {
2561  if(!PartCont || !TracksArray) return -99999;
2562  }
2563 
2564 
2565  AliAODTrack *Track = 0x0;
2566 
2567 
2568 
2569  //TList *trackList = new TList();
2570  Int_t triggers[100];
2571  for (Int_t iTrigger=0; iTrigger<100; iTrigger++) triggers[iTrigger] = 0;
2572  Int_t iTT = 0;
2573  Int_t NTracks=0;
2574  if (fJetShapeType == AliAnalysisTaskEmcalHfeTagging::kGenOnTheFly) NTracks = TracksArrayMC->GetEntriesFast();
2575  else NTracks = TracksArray->GetEntriesFast();
2576 
2577  for(Int_t i=0; i < NTracks; i++){
2579  if((Track = static_cast<AliAODTrack*>(PartContMC->GetAcceptParticle(i)))){
2580  if (!Track) continue;
2581  if(TMath::Abs(Track->Eta())>0.9) continue;
2582  if (Track->Pt()<0.15) continue;
2583  if ((Track->Pt() >= minpT) && (Track->Pt()< maxpT)) {
2584  triggers[iTT] = i;
2585  iTT++;
2586  }
2587  }
2588  }
2589  else{
2590  if((Track = static_cast<AliAODTrack*>(PartCont->GetAcceptTrack(i)))){
2591  if (!Track) continue;
2592  if(TMath::Abs(Track->Eta())>0.9) continue;
2593  if (Track->Pt()<0.15) continue;
2594  if ((Track->Pt() >= minpT) && (Track->Pt()< maxpT)) {
2595  triggers[iTT] = i;
2596  iTT++;
2597  }
2598  }
2599  }
2600  }
2601 
2602 
2603  if (iTT == 0) return -99999;
2604  Int_t nbRn = 0, index = 0 ;
2605  TRandom3* random = new TRandom3(0);
2606  nbRn = random->Integer(iTT);
2607  index = triggers[nbRn];
2608  //Printf("iTT Total= %d, nbRn = %d, Index = %d",iTT, nbRn, index );
2609  return index;
2610 
2611 }
2612 
2613 //__________________________________________________________________________________
2615 
2616  if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
2617  else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
2618  if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
2619  else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
2620  double dphi = mphi-vphi;
2621  if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
2622  else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
2623  return dphi;//dphi in [-Pi, Pi]
2624 }
2625 
2626 
2627 //________________________________________________________________________
2629  //
2630  // retrieve event objects
2631  //
2633  return kFALSE;
2634 
2635  return kTRUE;
2636 }
2637 
2638 //_________________________________________
2639 Bool_t AliAnalysisTaskEmcalHfeTagging::InclElecTrackCuts(const AliVVertex *pVietx,AliAODTrack *ietrack,Int_t nMother, Double_t listMother[])
2640 {
2641  // track cuts for inclusive electrons
2642 
2643  if(!ietrack->TestFilterMask(AliAODTrack::kTrkGlobalNoDCA)) return kFALSE;
2644  if(TMath::Abs(ietrack->Eta())>0.7) return kFALSE;
2645  if(ietrack->GetTPCNcls() < fTPCnCut) return kFALSE;
2646  if (ietrack->GetITSNcls() < fITSncut) return kFALSE;
2647  if(!ietrack->IsOn(AliAODTrack::kITSrefit)) return kFALSE;
2648  if(!ietrack->IsOn(AliAODTrack::kTPCrefit)) return kFALSE;
2649  if(!(ietrack->HasPointOnITSLayer(0) && ietrack->HasPointOnITSLayer(1))) return kFALSE;
2650  if(ietrack->GetTPCFoundFraction() < 0.6) return kFALSE;
2651  if(ietrack->Chi2perNDF() > 4) return kFALSE;
2652 
2653  Bool_t kinkmotherpass = kTRUE;
2654  for(Int_t kinkmother = 0; kinkmother < nMother; kinkmother++) {
2655  if(ietrack->GetID() == listMother[kinkmother]) {
2656  kinkmotherpass = kFALSE;
2657  continue;
2658  }
2659  }
2660  if(!kinkmotherpass) return kFALSE;
2661 
2662  Double_t d0z0[2]={-999,-999}, cov[3];
2663 
2664  if(ietrack->PropagateToDCA(pVietx, fVevent->GetMagneticField(), 20., d0z0, cov))
2665  if(TMath::Abs(d0z0[0]) > fDcaXYcut || TMath::Abs(d0z0[1]) > fDcaZcut) return kFALSE;
2666 
2667  return kTRUE;
2668 
2669 }
2670 //_________________________________________
2671 Bool_t AliAnalysisTaskEmcalHfeTagging::PhotElecTrackCuts(const AliVVertex *pVaetx,AliAODTrack *aetrack,Int_t nMother, Double_t listMother[])
2672 {
2673  // track cuts for associate tracks of photonic electrons
2674 
2675  if(!aetrack->TestFilterMask(AliAODTrack::kTrkTPCOnly)) return kFALSE;
2676  if(aetrack->Pt() < fAssPtCut) return kFALSE;
2677  if(TMath::Abs(aetrack->Eta())>0.9) return kFALSE;
2678  if(aetrack->GetTPCNcls() < fAssTPCnCut) return kFALSE;
2679  if (fAssITSrefitCut && !(aetrack->IsOn(AliAODTrack::kITSrefit))) return kFALSE;
2680  if(!aetrack->IsOn(AliAODTrack::kTPCrefit)) return kFALSE;
2681 
2682  return kTRUE;
2683 
2684 }
2685 
2686 //_______________________________________________________________________
2688 {
2689  // Called once at the end of the analysis.
2690 
2691  // fTreeObservableTagging = dynamic_cast<TTree*>(GetOutputData(1));
2692  // if (!fTreeObservableTagging){
2693  // Printf("ERROR: fTreeObservableTagging not available");
2694  // return;
2695  // }
2696 
2697 }
2698 
Int_t charge
Float_t Angularity(AliEmcalJet *jet, Int_t jetContNb)
Double_t GetFirstOrderSubtractedAngularity() const
Float_t GetJetLeSub(AliEmcalJet *jet, Int_t jetContNb)
Double_t Area() const
Definition: AliEmcalJet.h:130
Bool_t IsFromHFdecay(AliAODMCParticle *particle)
void GetNumberOfTrueElectrons(AliEmcalJet *jet, Int_t jetContNb, Int_t nMother, Double_t listMother[], Int_t &nTrueElec, Int_t &nTrueHFElec, Double_t &ptTrueHFElec)
const TString & GetRhoName() const
double Double_t
Definition: External.C:58
Double_t GetSecondOrderSubtractedSigma2() const
Definition: External.C:236
Bool_t Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
Bool_t FillHistograms()
Function filling histograms.
Float_t GetPartonEta6() const
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:341
AliJetContainer * GetJetContainer(Int_t i=0) const
Float_t GetPartonEta7() const
Double_t GetSecondOrderSubtractedConstituent() const
Double_t Eta() const
Definition: AliEmcalJet.h:121
Double_t Py() const
Definition: AliEmcalJet.h:107
Double_t Phi() const
Definition: AliEmcalJet.h:117
Double_t mass
Float_t GetPythiaEventWeight() const
Container with name, TClonesArray and cuts for particles.
Double_t GetPi0weight(Double_t mcPi0pT) const
Int_t GetLabel() const
Definition: AliEmcalJet.h:124
Float_t GetPartonPhi7() const
Bool_t RetrieveEventObjects()
Retrieve common objects from event.
Bool_t PhotElecTrackCuts(const AliVVertex *pVtx, AliAODTrack *aetrack, Int_t nMother, Double_t listMother[])
Float_t Circularity(AliEmcalJet *jet, Int_t jetContNb)
Bool_t HasMother(AliAODMCParticle *particle)
Int_t GetNumberOfPairs(AliEmcalJet *jet, AliAODTrack *track, const AliVVertex *pVtx, Int_t nMother, Double_t listMother[])
Float_t GetJetCircularity(AliEmcalJet *jet, Int_t jetContNb)
Int_t GetPartonFlag7() const
Container for particles within the EMCAL framework.
Float_t Sigma2(AliEmcalJet *jet, Int_t jetContNb)
Float_t GetPartonPhi6() const
Int_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:160
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:139
TString kData
Declare data MC or deltaAOD.
Double_t Px() const
Definition: AliEmcalJet.h:106
Float_t PTD(AliEmcalJet *jet, Int_t jetContNb)
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
const TString & GetRhoMassName() const
Int_t GetPartonFlag6() const
Double_t GetEtaweight(Double_t mcEtapT) const
AliParticleContainer * GetParticleContainer() const
Double_t GetFirstOrderSubtractedConstituent() const
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
Bool_t IsPrimary(AliAODMCParticle *particle)
Int_t GetNJets() const
Float_t GetJetMass(AliEmcalJet *jet, Int_t jetContNb)
Double_t GetSecondOrderSubtractedAngularity() const
Bool_t InclElecTrackCuts(const AliVVertex *pVtx, AliAODTrack *ietrack, Int_t nMother, Double_t listMother[])
Double_t RelativePhi(Double_t mphi, Double_t vphi)
virtual AliVParticle * GetAcceptParticle(Int_t i=-1) const
Double_t fCent
!event centrality
TClonesArray * GetArray() const
AliEmcalJet * GetNextAcceptJet()
void GetNumberOfElectrons(AliEmcalJet *jet, Int_t jetContNb, Int_t nMother, Double_t listMother[], Int_t &nIncElec, Int_t &nPhotElec, Double_t &pElec, Double_t &ptElec, Bool_t &hasElec)
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
Float_t LeSub(AliEmcalJet *jet, Int_t jetContNb)
Double_t Pt() const
Definition: AliEmcalJet.h:109
Double_t GetRhoVal(Int_t i=0) const
void GetWeightAndDecay(AliAODMCParticle *particle, Int_t &decay, Double_t &weight)
Double_t GetFirstOrderSubtractedCircularity() const
void Track()
decay
Definition: HFPtSpectrum.C:42
AliEmcalList * fOutput
!output list
Float_t GetJetCoronna(AliEmcalJet *jet, Int_t jetContNb)
Double_t GetSecondOrderSubtractedCircularity() const
AliTrackContainer * GetTrackContainer(Int_t i=0) const
Store some informaion about a Pythia eventThis class is used to store some information about a Pythia...
void SetMakeGeneralHistograms(Bool_t g)
Enable general histograms.
Float_t GetJetAngularity(AliEmcalJet *jet, Int_t jetContNb)
Float_t GetSigma2(AliEmcalJet *jet, Int_t jetContNb)
virtual AliVTrack * GetAcceptTrack(Int_t i=-1) const
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
const AliEmcalPythiaInfo * GetPythiaInfo() const
Double_t Pz() const
Definition: AliEmcalJet.h:108
Double_t GetFractionSharedPt(const AliEmcalJet *jet, AliParticleContainer *cont2=0x0) const
Declaration of class AliEmcalPythiaInfo.
const char Option_t
Definition: External.C:48
void UserCreateOutputObjects()
Main initialization function on the worker.
Double_t GetFirstOrderSubtractedSigma2() const
bool Bool_t
Definition: External.C:53
Bool_t IsFromLMdecay(AliAODMCParticle *particle)
AliEmcalJetShapeProperties * GetShapeProperties() const
Definition: AliEmcalJet.h:375
Int_t SelectTrigger(Float_t minpT, Float_t maxpT)
Float_t GetJetpTD(AliEmcalJet *jet, Int_t jetContNb)
Float_t Coronna(AliEmcalJet *jet, Int_t jetContNb)
Float_t GetPartonPt7() const
Double_t M() const
Definition: AliEmcalJet.h:120
void ResetCurrentID(Int_t i=-1)
Reset the iterator to a given index.
Container for jet within the EMCAL jet framework.
Definition: External.C:196
std::vector< int > GetPtSortedTrackConstituentIndexes(TClonesArray *tracks) const
Float_t GetPartonPt6() const
Float_t GetJetNumberOfConstituents(AliEmcalJet *jet, Int_t jetContNb)
AliEmcalJet * GetJet(Int_t i) const