]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWGHF/hfe/AliAnalysisTaskHFEFlow.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame_incremental
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (non-incremental) | history | HEAD
- Splines for LHC12i3 (also used for 12f1X, good statistics up to larger momenta...
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliAnalysisTaskHFEFlow.cxx
... / ...
CommitLineData
1/**************************************************************************\r
2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3* *\r
4* Author: The ALICE Off-line Project. *\r
5* Contributors are mentioned in the code where appropriate. *\r
6* *\r
7* Permission to use, copy, modify and distribute this software and its *\r
8* documentation strictly for non-commercial purposes is hereby granted *\r
9* without fee, provided that the above copyright notice appears in all *\r
10* copies and that both the copyright notice and this permission notice *\r
11* appear in the supporting documentation. The authors make no claims *\r
12* about the suitability of this software for any purpose. It is *\r
13* provided "as is" without express or implied warranty. *\r
14**************************************************************************/\r
15//\r
16// Flow task\r
17// \r
18// Authors:\r
19// Raphaelle Bailhache <R.Bailhache@gsi.de>\r
20// Theodor Rascanu <trascanu@stud.uni-frankfurt.de>\r
21//\r
22#include "TROOT.h"\r
23#include "TH1D.h"\r
24#include "TH2D.h"\r
25#include "TChain.h"\r
26#include "TVector2.h"\r
27#include "THnSparse.h"\r
28#include "TMath.h"\r
29#include "TRandom3.h"\r
30#include "TProfile.h"\r
31#include "TProfile2D.h"\r
32#include "TLorentzVector.h"\r
33#include "TParticle.h"\r
34#include "TF1.h"\r
35\r
36#include <TDirectory.h>\r
37#include <TTreeStream.h>\r
38\r
39#include "AliVEventHandler.h"\r
40#include "AliAnalysisTaskSE.h"\r
41#include "AliAnalysisManager.h"\r
42\r
43#include "AliVEvent.h"\r
44#include "AliESDInputHandler.h"\r
45#include "AliMCEvent.h"\r
46#include "AliESD.h"\r
47#include "AliESDEvent.h"\r
48#include "AliPID.h"\r
49#include "AliPIDResponse.h"\r
50#include "AliESDVZERO.h"\r
51#include "AliESDUtils.h"\r
52#include "AliMCParticle.h"\r
53#include "AliAODMCParticle.h"\r
54#include "AliAODEvent.h"\r
55#include "AliAODVertex.h"\r
56#include "AliAODTrack.h"\r
57#include "AliVTrack.h"\r
58#include "AliESDtrack.h"\r
59#include "AliESDtrackCuts.h"\r
60#include "AliAODTrack.h"\r
61#include "AliStack.h"\r
62#include "AliMCEvent.h"\r
63\r
64#include "AliFlowCandidateTrack.h"\r
65#include "AliFlowEvent.h"\r
66#include "AliFlowTrackCuts.h"\r
67#include "AliFlowVector.h"\r
68#include "AliFlowCommonConstants.h"\r
69#include "AliKFParticle.h"\r
70#include "AliKFVertex.h"\r
71\r
72#include "AliHFEcuts.h"\r
73#include "AliHFEpid.h"\r
74#include "AliHFEpidQAmanager.h"\r
75#include "AliHFEtools.h"\r
76#include "AliHFEVZEROEventPlane.h"\r
77\r
78#include "AliCentrality.h"\r
79#include "AliEventplane.h"\r
80#include "AliAnalysisTaskHFEFlow.h"\r
81#include "AliAODMCHeader.h"\r
82#include "TClonesArray.h"\r
83#include "AliHFENonPhotonicElectron.h"\r
84\r
85\r
86//____________________________________________________________________\r
87AliAnalysisTaskHFEFlow::AliAnalysisTaskHFEFlow() :\r
88 AliAnalysisTaskSE(),\r
89 fListHist(0x0), \r
90 fAODAnalysis(kFALSE),\r
91 fUseFilterAOD(kFALSE),\r
92 fApplyCut(kTRUE),\r
93 fFilter(1<<4),\r
94 fAODMCHeader(NULL),\r
95 fAODArrayMCInfo(NULL),\r
96 fBackgroundSubtraction(NULL),\r
97 fVZEROEventPlane(kFALSE),\r
98 fVZEROEventPlaneA(kFALSE),\r
99 fVZEROEventPlaneC(kFALSE),\r
100 fSubEtaGapTPC(kFALSE),\r
101 fEtaGap(0.0),\r
102 fNbBinsCentralityQCumulant(4),\r
103 fNbBinsPtQCumulant(12),\r
104 fMinPtQCumulant(0.2),\r
105 fMaxPtQCumulant(6.0),\r
106 fAfterBurnerOn(kFALSE),\r
107 fNonFlowNumberOfTrackClones(0),\r
108 fV1(0.),\r
109 fV2(0.),\r
110 fV3(0.),\r
111 fV4(0.),\r
112 fV5(0.),\r
113 fMaxNumberOfIterations(100),\r
114 fPrecisionPhi(0.001),\r
115 fUseMCReactionPlane(kFALSE),\r
116 fSP(kFALSE),\r
117 fMCPID(kFALSE),\r
118 fNoPID(kFALSE),\r
119 fChi2OverNDFCut(3.0),\r
120 fMaxdca(3.0),\r
121 fMaxopeningtheta(0.02),\r
122 fMaxopeningphi(0.1),\r
123 fMaxopening3D(0.1),\r
124 fMaxInvmass(0.1),\r
125 fSetMassConstraint(kFALSE),\r
126 fDebugLevel(0),\r
127 fMonitorEventPlane(kFALSE),\r
128 fMonitorContamination(kFALSE),\r
129 fMonitorPhotonic(kFALSE),\r
130 fMonitorWithoutPID(kFALSE),\r
131 fMonitorTrackCuts(kFALSE),\r
132 fMonitorQCumulant(kFALSE),\r
133 fcutsRP(0),\r
134 fcutsPOI(0),\r
135 fHFECuts(0),\r
136 fRejectKinkMother(kFALSE),\r
137 fPID(0),\r
138 fPIDTOFOnly(0),\r
139 fPIDqa(0),\r
140 fflowEvent(NULL),\r
141 fHFEBackgroundCuts(0),\r
142 fPIDBackground(0),\r
143 fPIDBackgroundqa(0),\r
144 fAlgorithmMA(kTRUE),\r
145 fArraytrack(NULL),\r
146 fCounterPoolBackground(0),\r
147 fHFEVZEROEventPlane(0x0),\r
148 fHistEV(0),\r
149 fEventPlane(0x0),\r
150 fEventPlaneaftersubtraction(0x0),\r
151 fFractionContamination(0x0),\r
152 fContaminationv2(0x0),\r
153 fCosSin2phiep(0x0),\r
154 fCos2phie(0x0),\r
155 fSin2phie(0x0),\r
156 fCos2phiep(0x0),\r
157 fSin2phiep(0x0),\r
158 fSin2phiephiep(0x0),\r
159 fCosResabc(0x0),\r
160 fSinResabc(0x0),\r
161 fProfileCosResab(0x0),\r
162 fProfileCosResac(0x0),\r
163 fProfileCosResbc(0x0),\r
164 fCosRes(0x0),\r
165 fSinRes(0x0),\r
166 fProfileCosRes(0x0),\r
167 fTrackingCuts(0x0),\r
168 fDeltaPhiMapsBeforePID(0x0),\r
169 fCosPhiMapsBeforePID(0x0),\r
170 fDeltaPhiMaps(0x0),\r
171 fDeltaPhiMapsContamination(0x0),\r
172 fCosPhiMaps(0x0),\r
173 fProfileCosPhiMaps(0x0),\r
174 fDeltaPhiMapsTaggedPhotonic(0x0),\r
175 //fCosPhiMapsTaggedPhotonic(0x0),\r
176 fDeltaPhiMapsTaggedNonPhotonic(0x0),\r
177 //fCosPhiMapsTaggedNonPhotonic(0x0),\r
178 fDeltaPhiMapsTaggedPhotonicLS(0x0),\r
179 //fCosPhiMapsTaggedPhotonicLS(0x0),\r
180 fMCSourceDeltaPhiMaps(0x0),\r
181 fOppSignDeltaPhiMaps(0x0),\r
182 fSameSignDeltaPhiMaps(0x0),\r
183 fOppSignAngle(0x0),\r
184 fSameSignAngle(0x0),\r
185 fDebugStreamer(0)\r
186{\r
187 // Constructor\r
188\r
189 for(Int_t k = 0; k < 10; k++) {\r
190 fBinCentralityLess[k] = 0.0;\r
191 }\r
192 for(Int_t k = 0; k < 11; k++) {\r
193 fContamination[k] = NULL;\r
194 fv2contamination[k] = NULL;\r
195 }\r
196 \r
197}\r
198//______________________________________________________________________________\r
199AliAnalysisTaskHFEFlow:: AliAnalysisTaskHFEFlow(const char *name) :\r
200 AliAnalysisTaskSE(name),\r
201 fListHist(0x0),\r
202 fAODAnalysis(kFALSE),\r
203 fUseFilterAOD(kFALSE),\r
204 fApplyCut(kTRUE),\r
205 fFilter(1<<4), \r
206 fAODMCHeader(NULL),\r
207 fAODArrayMCInfo(NULL),\r
208 fBackgroundSubtraction(NULL),\r
209 fVZEROEventPlane(kFALSE),\r
210 fVZEROEventPlaneA(kFALSE),\r
211 fVZEROEventPlaneC(kFALSE),\r
212 fSubEtaGapTPC(kFALSE),\r
213 fEtaGap(0.0),\r
214 fNbBinsCentralityQCumulant(4),\r
215 fNbBinsPtQCumulant(15),\r
216 fMinPtQCumulant(0.0),\r
217 fMaxPtQCumulant(6.0),\r
218 fAfterBurnerOn(kFALSE),\r
219 fNonFlowNumberOfTrackClones(0),\r
220 fV1(0.),\r
221 fV2(0.),\r
222 fV3(0.),\r
223 fV4(0.),\r
224 fV5(0.),\r
225 fMaxNumberOfIterations(100),\r
226 fPrecisionPhi(0.001),\r
227 fUseMCReactionPlane(kFALSE),\r
228 fSP(kFALSE),\r
229 fMCPID(kFALSE),\r
230 fNoPID(kFALSE),\r
231 fChi2OverNDFCut(3.0),\r
232 fMaxdca(3.0),\r
233 fMaxopeningtheta(0.02),\r
234 fMaxopeningphi(0.1),\r
235 fMaxopening3D(0.1),\r
236 fMaxInvmass(0.1),\r
237 fSetMassConstraint(kFALSE),\r
238 fDebugLevel(0),\r
239 fMonitorEventPlane(kFALSE),\r
240 fMonitorContamination(kFALSE),\r
241 fMonitorPhotonic(kFALSE),\r
242 fMonitorWithoutPID(kFALSE),\r
243 fMonitorTrackCuts(kFALSE),\r
244 fMonitorQCumulant(kFALSE),\r
245 fcutsRP(0),\r
246 fcutsPOI(0),\r
247 fHFECuts(0),\r
248 fRejectKinkMother(kFALSE),\r
249 fPID(0),\r
250 fPIDTOFOnly(0),\r
251 fPIDqa(0),\r
252 fflowEvent(NULL),\r
253 fHFEBackgroundCuts(0),\r
254 fPIDBackground(0),\r
255 fPIDBackgroundqa(0),\r
256 fAlgorithmMA(kTRUE), \r
257 fArraytrack(NULL),\r
258 fCounterPoolBackground(0),\r
259 fHFEVZEROEventPlane(0x0),\r
260 fHistEV(0),\r
261 fEventPlane(0x0),\r
262 fEventPlaneaftersubtraction(0x0),\r
263 fFractionContamination(0x0),\r
264 fContaminationv2(0x0),\r
265 fCosSin2phiep(0x0),\r
266 fCos2phie(0x0),\r
267 fSin2phie(0x0),\r
268 fCos2phiep(0x0),\r
269 fSin2phiep(0x0),\r
270 fSin2phiephiep(0x0),\r
271 fCosResabc(0x0),\r
272 fSinResabc(0x0),\r
273 fProfileCosResab(0x0),\r
274 fProfileCosResac(0x0),\r
275 fProfileCosResbc(0x0),\r
276 fCosRes(0x0),\r
277 fSinRes(0x0),\r
278 fProfileCosRes(0x0),\r
279 fTrackingCuts(0x0),\r
280 fDeltaPhiMapsBeforePID(0x0),\r
281 fCosPhiMapsBeforePID(0x0),\r
282 fDeltaPhiMaps(0x0),\r
283 fDeltaPhiMapsContamination(0x0),\r
284 fCosPhiMaps(0x0),\r
285 fProfileCosPhiMaps(0x0),\r
286 fDeltaPhiMapsTaggedPhotonic(0x0),\r
287 //fCosPhiMapsTaggedPhotonic(0x0),\r
288 fDeltaPhiMapsTaggedNonPhotonic(0x0),\r
289 //fCosPhiMapsTaggedNonPhotonic(0x0),\r
290 fDeltaPhiMapsTaggedPhotonicLS(0x0),\r
291 //fCosPhiMapsTaggedPhotonicLS(0x0),\r
292 fMCSourceDeltaPhiMaps(0x0),\r
293 fOppSignDeltaPhiMaps(0x0),\r
294 fSameSignDeltaPhiMaps(0x0),\r
295 fOppSignAngle(0x0),\r
296 fSameSignAngle(0x0),\r
297 fDebugStreamer(0)\r
298{\r
299 //\r
300 // named ctor\r
301 //\r
302 \r
303 for(Int_t k = 0; k < 10; k++) {\r
304 fBinCentralityLess[k] = 0.0;\r
305 }\r
306 fBinCentralityLess[0] = 0.0;\r
307 fBinCentralityLess[1] = 20.0;\r
308 fBinCentralityLess[2] = 40.0;\r
309 fBinCentralityLess[3] = 60.0;\r
310 fBinCentralityLess[4] = 80.0;\r
311\r
312 for(Int_t k = 0; k < 11; k++) {\r
313 fContamination[k] = NULL;\r
314 fv2contamination[k] = NULL;\r
315 }\r
316 \r
317 fPID = new AliHFEpid("hfePid");\r
318 fPIDqa = new AliHFEpidQAmanager;\r
319\r
320 fPIDBackground = new AliHFEpid("hfePidBackground");\r
321 fPIDBackgroundqa = new AliHFEpidQAmanager;\r
322\r
323 fPIDTOFOnly = new AliHFEpid("hfePidTOFOnly");\r
324\r
325 DefineInput(0,TChain::Class());\r
326 DefineOutput(1, TList::Class());\r
327 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {\r
328 // DefineOutput(bincless+2,AliFlowEventSimple::Class()); \r
329 //}\r
330 \r
331}\r
332//____________________________________________________________\r
333AliAnalysisTaskHFEFlow::AliAnalysisTaskHFEFlow(const AliAnalysisTaskHFEFlow &ref):\r
334 AliAnalysisTaskSE(ref),\r
335 fListHist(NULL),\r
336 fAODAnalysis(ref.fAODAnalysis), \r
337 fUseFilterAOD(ref.fUseFilterAOD),\r
338 fApplyCut(ref.fApplyCut),\r
339 fFilter(ref.fFilter),\r
340 fAODMCHeader(ref.fAODMCHeader),\r
341 fAODArrayMCInfo(ref.fAODArrayMCInfo),\r
342 fBackgroundSubtraction(ref.fBackgroundSubtraction),\r
343 fVZEROEventPlane(ref.fVZEROEventPlane),\r
344 fVZEROEventPlaneA(ref.fVZEROEventPlaneA),\r
345 fVZEROEventPlaneC(ref.fVZEROEventPlaneC),\r
346 fSubEtaGapTPC(ref.fSubEtaGapTPC),\r
347 fEtaGap(ref.fEtaGap),\r
348 fNbBinsCentralityQCumulant(ref.fNbBinsCentralityQCumulant),\r
349 fNbBinsPtQCumulant(ref.fNbBinsPtQCumulant),\r
350 fMinPtQCumulant(ref.fMinPtQCumulant),\r
351 fMaxPtQCumulant(ref.fMaxPtQCumulant),\r
352 fAfterBurnerOn(ref.fAfterBurnerOn),\r
353 fNonFlowNumberOfTrackClones(ref.fNonFlowNumberOfTrackClones),\r
354 fV1(ref.fV1),\r
355 fV2(ref.fV2),\r
356 fV3(ref.fV3),\r
357 fV4(ref.fV4),\r
358 fV5(ref.fV5),\r
359 fMaxNumberOfIterations(ref.fMaxNumberOfIterations),\r
360 fPrecisionPhi(ref.fPrecisionPhi),\r
361 fUseMCReactionPlane(ref.fUseMCReactionPlane),\r
362 fSP(ref.fSP),\r
363 fMCPID(ref.fMCPID),\r
364 fNoPID(ref.fNoPID),\r
365 fChi2OverNDFCut(ref.fChi2OverNDFCut),\r
366 fMaxdca(ref.fMaxdca),\r
367 fMaxopeningtheta(ref.fMaxopeningtheta),\r
368 fMaxopeningphi(ref.fMaxopeningphi),\r
369 fMaxopening3D(ref.fMaxopening3D),\r
370 fMaxInvmass(ref.fMaxInvmass),\r
371 fSetMassConstraint(ref.fSetMassConstraint),\r
372 fDebugLevel(ref.fDebugLevel),\r
373 fMonitorEventPlane(ref.fMonitorEventPlane),\r
374 fMonitorContamination(ref.fMonitorContamination),\r
375 fMonitorPhotonic(ref.fMonitorPhotonic),\r
376 fMonitorWithoutPID(ref.fMonitorWithoutPID),\r
377 fMonitorTrackCuts(ref.fMonitorTrackCuts),\r
378 fMonitorQCumulant(ref.fMonitorQCumulant),\r
379 fcutsRP(NULL),\r
380 fcutsPOI(NULL),\r
381 fHFECuts(NULL),\r
382 fRejectKinkMother(ref.fRejectKinkMother),\r
383 fPID(NULL),\r
384 fPIDTOFOnly(NULL),\r
385 fPIDqa(NULL),\r
386 fflowEvent(NULL),\r
387 fHFEBackgroundCuts(NULL),\r
388 fPIDBackground(NULL),\r
389 fPIDBackgroundqa(NULL),\r
390 fAlgorithmMA(ref.fAlgorithmMA),\r
391 fArraytrack(NULL),\r
392 fCounterPoolBackground(ref.fCounterPoolBackground),\r
393 fHFEVZEROEventPlane(NULL),\r
394 fHistEV(NULL),\r
395 fEventPlane(NULL),\r
396 fEventPlaneaftersubtraction(NULL),\r
397 fFractionContamination(NULL),\r
398 fContaminationv2(NULL),\r
399 fCosSin2phiep(NULL),\r
400 fCos2phie(NULL),\r
401 fSin2phie(NULL),\r
402 fCos2phiep(NULL),\r
403 fSin2phiep(NULL),\r
404 fSin2phiephiep(NULL),\r
405 fCosResabc(NULL),\r
406 fSinResabc(NULL),\r
407 fProfileCosResab(NULL),\r
408 fProfileCosResac(NULL),\r
409 fProfileCosResbc(NULL),\r
410 fCosRes(NULL),\r
411 fSinRes(NULL),\r
412 fProfileCosRes(NULL),\r
413 fTrackingCuts(NULL),\r
414 fDeltaPhiMapsBeforePID(NULL),\r
415 fCosPhiMapsBeforePID(NULL),\r
416 fDeltaPhiMaps(NULL),\r
417 fDeltaPhiMapsContamination(NULL),\r
418 fCosPhiMaps(NULL),\r
419 fProfileCosPhiMaps(NULL),\r
420 fDeltaPhiMapsTaggedPhotonic(NULL),\r
421 //fCosPhiMapsTaggedPhotonic(NULL),\r
422 fDeltaPhiMapsTaggedNonPhotonic(NULL),\r
423 //fCosPhiMapsTaggedNonPhotonic(NULL),\r
424 fDeltaPhiMapsTaggedPhotonicLS(NULL),\r
425 //fCosPhiMapsTaggedPhotonicLS(NULL),\r
426 fMCSourceDeltaPhiMaps(NULL),\r
427 fOppSignDeltaPhiMaps(NULL),\r
428 fSameSignDeltaPhiMaps(NULL),\r
429 fOppSignAngle(NULL),\r
430 fSameSignAngle(NULL),\r
431 fDebugStreamer(0)\r
432{\r
433 //\r
434 // Copy Constructor\r
435 //\r
436 ref.Copy(*this);\r
437}\r
438\r
439//____________________________________________________________\r
440AliAnalysisTaskHFEFlow &AliAnalysisTaskHFEFlow::operator=(const AliAnalysisTaskHFEFlow &ref){\r
441 //\r
442 // Assignment operator\r
443 //\r
444 if(this == &ref) \r
445 ref.Copy(*this);\r
446 return *this;\r
447}\r
448\r
449//____________________________________________________________\r
450void AliAnalysisTaskHFEFlow::Copy(TObject &o) const {\r
451 // \r
452 // Copy into object o\r
453 //\r
454 AliAnalysisTaskHFEFlow &target = dynamic_cast<AliAnalysisTaskHFEFlow &>(o);\r
455 target.fAODAnalysis = fAODAnalysis;\r
456 target.fUseFilterAOD = fUseFilterAOD;\r
457 target.fApplyCut = fApplyCut;\r
458 target.fFilter = fFilter;\r
459 target.fAODMCHeader = fAODMCHeader;\r
460 target.fAODArrayMCInfo = fAODArrayMCInfo;\r
461 target.fBackgroundSubtraction = fBackgroundSubtraction;\r
462 target.fVZEROEventPlane = fVZEROEventPlane;\r
463 target.fVZEROEventPlaneA = fVZEROEventPlaneA;\r
464 target.fVZEROEventPlaneC = fVZEROEventPlaneC;\r
465 target.fSubEtaGapTPC = fSubEtaGapTPC;\r
466 target.fEtaGap = fEtaGap;\r
467 target.fNbBinsCentralityQCumulant = fNbBinsCentralityQCumulant;\r
468 target.fNbBinsPtQCumulant = fNbBinsPtQCumulant;\r
469 target.fMinPtQCumulant = fMinPtQCumulant;\r
470 target.fMaxPtQCumulant = fMaxPtQCumulant;\r
471 target.fAfterBurnerOn = fAfterBurnerOn;\r
472 target.fNonFlowNumberOfTrackClones = fNonFlowNumberOfTrackClones;\r
473 target.fV1 = fV1;\r
474 target.fV2 = fV2;\r
475 target.fV3 = fV3;\r
476 target.fV4 = fV4;\r
477 target.fV5 = fV5;\r
478 target.fMaxNumberOfIterations = fMaxNumberOfIterations;\r
479 target.fPrecisionPhi = fPrecisionPhi;\r
480 target.fUseMCReactionPlane = fUseMCReactionPlane;\r
481 target.fSP = fSP;\r
482 target.fMCPID = fMCPID;\r
483 target.fNoPID = fNoPID;\r
484 target.fChi2OverNDFCut = fChi2OverNDFCut;\r
485 target.fMaxdca = fMaxdca;\r
486 target.fMaxopeningtheta = fMaxopeningtheta;\r
487 target.fMaxopeningphi = fMaxopeningphi;\r
488 target.fMaxopening3D = fMaxopening3D;\r
489 target.fMaxInvmass = fMaxInvmass;\r
490 target.fSetMassConstraint = fSetMassConstraint;\r
491 target.fAlgorithmMA = fAlgorithmMA;\r
492 target.fCounterPoolBackground = fCounterPoolBackground;\r
493 target.fDebugLevel = fDebugLevel;\r
494 target.fMonitorEventPlane = fMonitorEventPlane;\r
495 target.fMonitorContamination = fMonitorContamination;\r
496 target.fMonitorPhotonic = fMonitorPhotonic;\r
497 target.fMonitorWithoutPID = fMonitorWithoutPID;\r
498 target.fMonitorTrackCuts = fMonitorTrackCuts;\r
499 target.fMonitorQCumulant = fMonitorQCumulant;\r
500 target.fcutsRP = fcutsRP;\r
501 target.fcutsPOI = fcutsPOI;\r
502 target.fHFECuts = fHFECuts;\r
503 target.fRejectKinkMother = fRejectKinkMother;\r
504 target.fPID = fPID;\r
505 target.fPIDqa = fPIDqa;\r
506 target.fHFEVZEROEventPlane = fHFEVZEROEventPlane;\r
507 for(Int_t k = 0; k < 11; k++) {\r
508 target.fContamination[k] = fContamination[k];\r
509 target.fv2contamination[k] = fv2contamination[k];\r
510 }\r
511 target.fDebugStreamer=fDebugStreamer;\r
512}\r
513//____________________________________________________________\r
514AliAnalysisTaskHFEFlow::~AliAnalysisTaskHFEFlow(){\r
515 //\r
516 // Destructor\r
517 //\r
518 if(fArraytrack) delete fArraytrack;\r
519 if(fListHist) delete fListHist;\r
520 if(fcutsRP) delete fcutsRP;\r
521 if(fcutsPOI) delete fcutsPOI;\r
522 if(fHFECuts) delete fHFECuts;\r
523 if(fPID) delete fPID;\r
524 if(fPIDTOFOnly) delete fPIDTOFOnly;\r
525 //if(fPIDqa) delete fPIDqa;\r
526 if(fflowEvent) delete fflowEvent;\r
527 if(fHFEBackgroundCuts) delete fHFEBackgroundCuts;\r
528 if(fPIDBackground) delete fPIDBackground;\r
529 if(fBackgroundSubtraction) delete fBackgroundSubtraction;\r
530 //if(fPIDBackgroundqa) delete fPIDBackgroundqa;\r
531 //if(fHFEVZEROEventPlane) delete fHFEVZEROEventPlane;\r
532 if ( fDebugStreamer ) delete fDebugStreamer;\r
533\r
534}\r
535//________________________________________________________________________\r
536void AliAnalysisTaskHFEFlow::UserCreateOutputObjects()\r
537{\r
538\r
539 //********************\r
540 // Create histograms\r
541 //********************\r
542\r
543 //**************\r
544 // Cuts\r
545 //**************\r
546\r
547 //---------Data selection----------\r
548 //kMC, kGlobal, kTPCstandalone, kSPDtracklet, kPMD\r
549 //AliFlowTrackCuts::trackParameterType rptype = AliFlowTrackCuts::kGlobal;\r
550 //AliFlowTrackCuts::trackParameterType poitype = AliFlowTrackCuts::kGlobal;\r
551\r
552 //---------Parameter mixing--------\r
553 //kPure - no mixing, kTrackWithMCkine, kTrackWithMCPID, kTrackWithMCpt\r
554 //AliFlowTrackCuts::trackParameterMix rpmix = AliFlowTrackCuts::kPure;\r
555 //AliFlowTrackCuts::trackParameterMix poimix = AliFlowTrackCuts::kPure;\r
556\r
557 AliDebug(2,"test");\r
558\r
559 AliDebug(2,"AliAnalysisTaskHFEFlow: User create output objects");\r
560 \r
561 // AOD or ESD\r
562 AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
563 if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){\r
564 SetAODAnalysis(kTRUE);\r
565 AliDebug(2,"Put AOD analysis on");\r
566 } else {\r
567 SetAODAnalysis(kFALSE);\r
568 }\r
569\r
570 AliDebug(2,"AliAnalysisTaskHFEFlow: AOD ESD");\r
571\r
572 // RP TRACK CUTS:\r
573 fcutsRP = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();\r
574 fcutsRP->SetName("StandartTPC");\r
575 fcutsRP->SetEtaRange(-0.9,0.9);\r
576 fcutsRP->SetQA(kTRUE);\r
577 //TList *qaCutsRP = fcutsRP->GetQA();\r
578 //qaCutsRP->SetName("QA_StandartTPC_RP");\r
579\r
580 AliDebug(2,"AliAnalysisTaskHFEFlow: cutsRP");\r
581\r
582 //POI TRACK CUTS:\r
583 fcutsPOI = new AliFlowTrackCuts("dummy");\r
584 fcutsPOI->SetParamType(AliFlowTrackCuts::kGlobal);\r
585 fcutsPOI->SetPtRange(+1,-1); // select nothing QUICK\r
586 fcutsPOI->SetEtaRange(+1,-1); // select nothing VZERO\r
587\r
588 if( fflowEvent ){ \r
589 fflowEvent->~AliFlowEvent();\r
590 new(fflowEvent) AliFlowEvent(fcutsRP,fcutsPOI);\r
591 }\r
592 else fflowEvent = new AliFlowEvent(fcutsRP,fcutsPOI);\r
593 \r
594 AliDebug(2,"AliAnalysisTaskHFEFlow: cutsPOI");\r
595\r
596 // Flow\r
597 AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();\r
598 cc->SetNbinsMult(10000);\r
599 cc->SetMultMin(0);\r
600 cc->SetMultMax(10000.);\r
601 cc->SetNbinsPt(fNbBinsPtQCumulant);\r
602 cc->SetPtMin(fMinPtQCumulant);\r
603 cc->SetPtMax(fMaxPtQCumulant);\r
604 cc->SetNbinsPhi(180);\r
605 cc->SetPhiMin(0.0);\r
606 cc->SetPhiMax(TMath::TwoPi());\r
607 cc->SetNbinsEta(200);\r
608 cc->SetEtaMin(-0.9);\r
609 cc->SetEtaMax(+0.9);\r
610 cc->SetNbinsQ(500);\r
611 cc->SetQMin(0.0);\r
612 cc->SetQMax(3.0);\r
613\r
614 AliDebug(2,"AliAnalysisTaskHFEFlow: common constants");\r
615\r
616 \r
617 // HFE cuts\r
618\r
619 if(!fHFECuts){\r
620 fHFECuts = new AliHFEcuts;\r
621 fHFECuts->CreateStandardCuts();\r
622 }\r
623 fHFECuts->Initialize();\r
624 if(fAODAnalysis) fHFECuts->SetAOD(); \r
625\r
626 AliDebug(2,"AliAnalysisTaskHFEFlow: HFE cuts");\r
627\r
628\r
629 // PID HFE\r
630 //fPID->SetHasMCData(HasMCData());\r
631 if(!fPID) {\r
632 fPID =new AliHFEpid("hfePid");\r
633 AliDebug(2,"AliAnalysisTaskHFEFlow: pid init 0");\r
634 }\r
635 AliDebug(2,"AliAnalysisTaskHFEFlow: pid init 1");\r
636 if(!fPID->GetNumberOfPIDdetectors()) fPID->AddDetector("TPC", 0);\r
637 AliDebug(2,Form("AliAnalysisTaskHFEFlow: GetNumber of PID detectors %d",fPID->GetNumberOfPIDdetectors()));\r
638 fPID->InitializePID();\r
639 AliDebug(2,"Init ");\r
640 fPIDqa->Initialize(fPID);\r
641 AliDebug(2,"Init qa");\r
642 fPID->SortDetectors();\r
643 AliDebug(2,"Sort detectors");\r
644\r
645 AliDebug(2,"AliAnalysisTaskHFEFlow: pid and pidqa");\r
646\r
647 if(!fPIDTOFOnly->GetNumberOfPIDdetectors()) fPIDTOFOnly->AddDetector("TPC", 0);\r
648 fPIDTOFOnly->InitializePID();\r
649 fPIDTOFOnly->SortDetectors();\r
650\r
651 AliDebug(2,"AliAnalysisTaskHFEFlow: pidtof");\r
652\r
653 // HFE Background cuts\r
654\r
655 if(!fHFEBackgroundCuts){\r
656 fHFEBackgroundCuts = new AliESDtrackCuts();\r
657 fHFEBackgroundCuts->SetName("nackgroundcuts");\r
658 //Configure Default Track Cuts\r
659 fHFEBackgroundCuts->SetAcceptKinkDaughters(kFALSE);\r
660 fHFEBackgroundCuts->SetRequireTPCRefit(kTRUE);\r
661 fHFEBackgroundCuts->SetEtaRange(-0.9,0.9);\r
662 fHFEBackgroundCuts->SetRequireSigmaToVertex(kTRUE);\r
663 fHFEBackgroundCuts->SetMaxChi2PerClusterTPC(4.0);\r
664 fHFEBackgroundCuts->SetMinNClustersTPC(50);\r
665 fHFEBackgroundCuts->SetPtRange(0.3,1e10);\r
666 }\r
667 \r
668 AliDebug(2,"AliAnalysisTaskHFEFlow: hfe background");\r
669\r
670 // PID background HFE\r
671 if(!fPIDBackground->GetNumberOfPIDdetectors()) fPIDBackground->AddDetector("TPC", 0);\r
672 fPIDBackground->InitializePID();\r
673 fPIDBackgroundqa->Initialize(fPIDBackground);\r
674 fPIDBackground->SortDetectors();\r
675\r
676 AliDebug(2,"AliAnalysisTaskHFEFlow: pid background");\r
677\r
678 if (fMonitorPhotonic) {\r
679 if(!fBackgroundSubtraction) fBackgroundSubtraction = new AliHFENonPhotonicElectron();\r
680 fBackgroundSubtraction->Init();\r
681 }\r
682 \r
683\r
684\r
685 //**************************\r
686 // Bins for the THnSparse\r
687 //**************************\r
688\r
689 /*\r
690 Int_t nBinsPt = 44;\r
691 Double_t minPt = 0.1;\r
692 Double_t maxPt = 20.0;\r
693 Double_t binLimLogPt[nBinsPt+1];\r
694 Double_t binLimPt[nBinsPt+1];\r
695 for(Int_t i=0; i<=nBinsPt; i++) binLimLogPt[i]=(Double_t)TMath::Log10(minPt) + (TMath::Log10(maxPt)-TMath::Log10(minPt))/nBinsPt*(Double_t)i ;\r
696 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);\r
697 */\r
698\r
699 Int_t nBinsPt = 24;\r
700 Double_t binLimPt[25] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2,\r
701 1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 3., 3.5, 4., 5.,\r
702 6.};\r
703\r
704\r
705 Int_t nBinsPtPlus = fNbBinsPtQCumulant;\r
706 Double_t minPtPlus = fMinPtQCumulant;\r
707 Double_t maxPtPlus = fMaxPtQCumulant;\r
708 Double_t binLimPtPlus[nBinsPtPlus+1];\r
709 for(Int_t i=0; i<=nBinsPtPlus; i++) binLimPtPlus[i]=(Double_t)minPtPlus + (maxPtPlus-minPtPlus)/nBinsPtPlus*(Double_t)i ;\r
710\r
711 Int_t nBinsEta = 8;\r
712 Double_t minEta = -0.8;\r
713 Double_t maxEta = 0.8;\r
714 Double_t binLimEta[nBinsEta+1];\r
715 for(Int_t i=0; i<=nBinsEta; i++) binLimEta[i]=(Double_t)minEta + (maxEta-minEta)/nBinsEta*(Double_t)i ;\r
716\r
717 Int_t nBinsStep = 7;\r
718 Double_t minStep = 0.;\r
719 Double_t maxStep = 7.;\r
720 Double_t binLimStep[nBinsStep+1];\r
721 for(Int_t i=0; i<=nBinsStep; i++) binLimStep[i]=(Double_t)minStep + (maxStep-minStep)/nBinsStep*(Double_t)i ;\r
722\r
723 Int_t nBinsEtaLess = 2;\r
724 Double_t binLimEtaLess[nBinsEtaLess+1];\r
725 for(Int_t i=0; i<=nBinsEtaLess; i++) binLimEtaLess[i]=(Double_t)minEta + (maxEta-minEta)/nBinsEtaLess*(Double_t)i ;\r
726 \r
727 Int_t nBinsCos = 50;\r
728 Double_t minCos = -1.0;\r
729 Double_t maxCos = 1.0;\r
730 Double_t binLimCos[nBinsCos+1];\r
731 for(Int_t i=0; i<=nBinsCos; i++) binLimCos[i]=(Double_t)minCos + (maxCos-minCos)/nBinsCos*(Double_t)i ;\r
732\r
733 // Int_t nBinsCosSP = 50;\r
734 // Double_t minCosSP = -100.0;\r
735 // Double_t maxCosSP = 100.0;\r
736 // Double_t binLimCosSP[nBinsCosSP+1];\r
737 // for(Int_t i=0; i<=nBinsCosSP; i++) binLimCosSP[i]=(Double_t)minCosSP + (maxCosSP-minCosSP)/nBinsCosSP*(Double_t)i ;\r
738 \r
739 Int_t nBinsC = 11;\r
740 Double_t minC = 0.0;\r
741 Double_t maxC = 11.0;\r
742 Double_t binLimC[nBinsC+1];\r
743 for(Int_t i=0; i<=nBinsC; i++) binLimC[i]=(Double_t)minC + (maxC-minC)/nBinsC*(Double_t)i ;\r
744\r
745 Int_t nBinsCMore = 20;\r
746 Double_t minCMore = 0.0;\r
747 Double_t maxCMore = 20.0;\r
748 Double_t binLimCMore[nBinsCMore+1];\r
749 for(Int_t i=0; i<=nBinsCMore; i++) binLimCMore[i]=(Double_t)minCMore + (maxCMore-minCMore)/nBinsCMore*(Double_t)i ;\r
750\r
751 Int_t nBinsPhi = 8;\r
752 Double_t minPhi = 0.0;\r
753 Double_t maxPhi = TMath::Pi();\r
754 Double_t binLimPhi[nBinsPhi+1];\r
755 for(Int_t i=0; i<=nBinsPhi; i++) {\r
756 binLimPhi[i]=(Double_t)minPhi + (maxPhi-minPhi)/nBinsPhi*(Double_t)i ;\r
757 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));\r
758 }\r
759\r
760 Int_t nBinsPhiLess = 2.0;\r
761 Double_t minPhiLess = 0.0;\r
762 Double_t maxPhiLess = 2.0;\r
763 Double_t binLimPhiLess[nBinsPhiLess+1];\r
764 for(Int_t i=0; i<=nBinsPhiLess; i++) {\r
765 binLimPhiLess[i]=(Double_t)minPhiLess + (maxPhiLess-minPhiLess)/nBinsPhiLess*(Double_t)i ;\r
766 }\r
767\r
768 Int_t nBinsTPCdEdx = 140;\r
769 Double_t minTPCdEdx = -12.0;\r
770 Double_t maxTPCdEdx = 12.0;\r
771 Double_t binLimTPCdEdx[nBinsTPCdEdx+1];\r
772 for(Int_t i=0; i<=nBinsTPCdEdx; i++) {\r
773 binLimTPCdEdx[i]=(Double_t)minTPCdEdx + (maxTPCdEdx-minTPCdEdx)/nBinsTPCdEdx*(Double_t)i ;\r
774 }\r
775\r
776 Int_t nBinsAngle = 40;\r
777 Double_t minAngle = 0.0;\r
778 Double_t maxAngle = 1.0;\r
779 Double_t binLimAngle[nBinsAngle+1];\r
780 for(Int_t i=0; i<=nBinsAngle; i++) {\r
781 binLimAngle[i]=(Double_t)minAngle + (maxAngle-minAngle)/nBinsAngle*(Double_t)i ;\r
782 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));\r
783 }\r
784\r
785 Int_t nBinsCharge = 2;\r
786 Double_t minCharge = -1.0;\r
787 Double_t maxCharge = 1.0;\r
788 Double_t binLimCharge[nBinsCharge+1];\r
789 for(Int_t i=0; i<=nBinsCharge; i++) binLimCharge[i]=(Double_t)minCharge + (maxCharge-minCharge)/nBinsCharge*(Double_t)i ;\r
790\r
791 Int_t nBinsSource = 10;\r
792 Double_t minSource = 0.;\r
793 Double_t maxSource = 10.;\r
794 Double_t binLimSource[nBinsSource+1];\r
795 for(Int_t i=0; i<=nBinsSource; i++) binLimSource[i]=(Double_t)minSource + (maxSource-minSource)/nBinsSource*(Double_t)i ;\r
796\r
797 Int_t nBinsInvMass = 50;\r
798 Double_t minInvMass = 0.;\r
799 Double_t maxInvMass = 0.3;\r
800 Double_t binLimInvMass[nBinsInvMass+1];\r
801 for(Int_t i=0; i<=nBinsInvMass; i++) binLimInvMass[i]=(Double_t)minInvMass + (maxInvMass-minInvMass)/nBinsInvMass*(Double_t)i ;\r
802\r
803 AliDebug(2,"AliAnalysisTaskHFEFlow: variables");\r
804 \r
805 //******************\r
806 // Histograms\r
807 //******************\r
808 \r
809 fListHist = new TList();\r
810 fListHist->SetOwner();\r
811\r
812 AliDebug(2,"AliAnalysisTaskHFEFlow: list");\r
813\r
814 // Minimum histos\r
815\r
816 // Histos\r
817 fHistEV = new TH2D("fHistEV", "events", 3, 0, 3, 3, 0,3);\r
818 \r
819 AliDebug(2,"AliAnalysisTaskHFEFlow: histev");\r
820\r
821 // Event plane as function of phiep, centrality\r
822 const Int_t nDima=4;\r
823 Int_t nBina[nDima] = {nBinsPhi,nBinsPhi,nBinsPhi,nBinsC};\r
824 fEventPlane = new THnSparseF("EventPlane","EventPlane",nDima,nBina);\r
825 fEventPlane->SetBinEdges(0,binLimPhi);\r
826 fEventPlane->SetBinEdges(1,binLimPhi);\r
827 fEventPlane->SetBinEdges(2,binLimPhi);\r
828 fEventPlane->SetBinEdges(3,binLimC);\r
829 fEventPlane->Sumw2();\r
830\r
831 AliDebug(2,"AliAnalysisTaskHFEFlow: eventplane");\r
832\r
833 // Fraction of contamination, centrality\r
834 const Int_t nDimcont=2;\r
835 Int_t nBincont[nDimcont] = {nBinsPt,nBinsC};\r
836 fFractionContamination = new THnSparseF("Contamination","Contamination",nDimcont,nBincont);\r
837 fFractionContamination->SetBinEdges(0,binLimPt);\r
838 fFractionContamination->SetBinEdges(1,binLimC);\r
839 fFractionContamination->Sumw2();\r
840 // \r
841 fContaminationv2 = new TProfile2D("Contaminationv2","",nBinsC,binLimC,nBinsPt,binLimPt);\r
842 fContaminationv2->Sumw2();\r
843\r
844 AliDebug(2,"AliAnalysisTaskHFEFlow: fraction of contamination");\r
845\r
846 // Resolution cosres_abc centrality\r
847 const Int_t nDimfbis=4;\r
848 Int_t nBinfbis[nDimfbis] = {nBinsCos,nBinsCos,nBinsCos,nBinsCMore};\r
849 fCosResabc = new THnSparseF("CosRes_abc","CosRes_abc",nDimfbis,nBinfbis);\r
850 fCosResabc->SetBinEdges(0,binLimCos);\r
851 fCosResabc->SetBinEdges(1,binLimCos);\r
852 fCosResabc->SetBinEdges(2,binLimCos);\r
853 fCosResabc->SetBinEdges(3,binLimCMore);\r
854 fCosResabc->Sumw2();\r
855\r
856 AliDebug(2,"AliAnalysisTaskHFEFlow: cosresabc");\r
857\r
858 // Resolution cosres centrality\r
859 const Int_t nDimf=2;\r
860 Int_t nBinf[nDimf] = {nBinsCos, nBinsCMore};\r
861 fCosRes = new THnSparseF("CosRes","CosRes",nDimf,nBinf);\r
862 fCosRes->SetBinEdges(0,binLimCos);\r
863 fCosRes->SetBinEdges(1,binLimCMore);\r
864 fCosRes->Sumw2();\r
865\r
866 AliDebug(2,"AliAnalysisTaskHFEFlow: cosres");\r
867\r
868 // Maps delta phi\r
869 const Int_t nDimg=5;\r
870 Int_t nBing[nDimg] = {nBinsPhi,nBinsC,nBinsPt, nBinsCharge,nBinsEtaLess};\r
871 fDeltaPhiMaps = new THnSparseF("DeltaPhiMaps","DeltaPhiMaps",nDimg,nBing);\r
872 fDeltaPhiMaps->SetBinEdges(0,binLimPhi);\r
873 fDeltaPhiMaps->SetBinEdges(1,binLimC);\r
874 fDeltaPhiMaps->SetBinEdges(2,binLimPt);\r
875 fDeltaPhiMaps->SetBinEdges(3,binLimCharge);\r
876 fDeltaPhiMaps->SetBinEdges(4,binLimEtaLess);\r
877 fDeltaPhiMaps->Sumw2(); \r
878\r
879 AliDebug(2,"AliAnalysisTaskHFEFlow: deltaphimaps");\r
880\r
881 // Maps cos phi\r
882 const Int_t nDimh=5;\r
883 Int_t nBinh[nDimh] = {nBinsCos,nBinsC,nBinsPt,nBinsCharge,nBinsEtaLess};\r
884 fCosPhiMaps = new THnSparseF("CosPhiMaps","CosPhiMaps",nDimh,nBinh);\r
885 fCosPhiMaps->SetBinEdges(0,binLimCos);\r
886 fCosPhiMaps->SetBinEdges(1,binLimC);\r
887 fCosPhiMaps->SetBinEdges(2,binLimPt);\r
888 fCosPhiMaps->SetBinEdges(3,binLimCharge);\r
889 fCosPhiMaps->SetBinEdges(4,binLimEtaLess);\r
890 fCosPhiMaps->Sumw2();\r
891\r
892 AliDebug(2,"AliAnalysisTaskHFEFlow: cosphimaps");\r
893\r
894 //\r
895 // fMonitorEventPlane\r
896 //\r
897 //\r
898\r
899 if(fMonitorEventPlane) { \r
900 // Event Plane after subtraction as function of phiep, centrality, pt, eta\r
901 const Int_t nDimb=2;\r
902 Int_t nBinb[nDimb] = {nBinsPhi, nBinsC};\r
903 fEventPlaneaftersubtraction = new THnSparseF("EventPlane_aftersubtraction","EventPlane_aftersubtraction",nDimb,nBinb);\r
904 fEventPlaneaftersubtraction->SetBinEdges(0,binLimPhi);\r
905 fEventPlaneaftersubtraction->SetBinEdges(1,binLimC);\r
906 fEventPlaneaftersubtraction->Sumw2();\r
907\r
908 AliDebug(2,"AliAnalysisTaskHFEFlow: eventplane after sub");\r
909 \r
910 // Monitoring of the event Plane cos(2phi) sin(2phi) centrality\r
911 const Int_t nDimi=3;\r
912 Int_t nBini[nDimi] = {nBinsCos, nBinsCos, nBinsCMore};\r
913 fCosSin2phiep = new THnSparseF("CosSin2phiep","CosSin2phiep",nDimi,nBini);\r
914 fCosSin2phiep->SetBinEdges(0,binLimCos);\r
915 fCosSin2phiep->SetBinEdges(1,binLimCos);\r
916 fCosSin2phiep->SetBinEdges(2,binLimCMore);\r
917 fCosSin2phiep->Sumw2();\r
918\r
919 AliDebug(2,"AliAnalysisTaskHFEFlow: cossin2phiep");\r
920 \r
921 // Monitoring Event plane after subtraction of the track\r
922 const Int_t nDime=4;\r
923 Int_t nBine[nDime] = {nBinsCos, nBinsC, nBinsPt, nBinsEta};\r
924 fCos2phie = new THnSparseF("cos2phie","cos2phie",nDime,nBine);\r
925 fCos2phie->SetBinEdges(2,binLimPt);\r
926 fCos2phie->SetBinEdges(3,binLimEta);\r
927 fCos2phie->SetBinEdges(0,binLimCos);\r
928 fCos2phie->SetBinEdges(1,binLimC);\r
929 fCos2phie->Sumw2();\r
930 AliDebug(2,"AliAnalysisTaskHFEFlow: cos2phie");\r
931 fSin2phie = new THnSparseF("sin2phie","sin2phie",nDime,nBine);\r
932 fSin2phie->SetBinEdges(2,binLimPt);\r
933 fSin2phie->SetBinEdges(3,binLimEta);\r
934 fSin2phie->SetBinEdges(0,binLimCos);\r
935 fSin2phie->SetBinEdges(1,binLimC);\r
936 fSin2phie->Sumw2();\r
937 AliDebug(2,"AliAnalysisTaskHFEFlow: sin2phie");\r
938 fCos2phiep = new THnSparseF("cos2phiep","cos2phiep",nDime,nBine);\r
939 fCos2phiep->SetBinEdges(2,binLimPt);\r
940 fCos2phiep->SetBinEdges(3,binLimEta);\r
941 fCos2phiep->SetBinEdges(0,binLimCos);\r
942 fCos2phiep->SetBinEdges(1,binLimC);\r
943 fCos2phiep->Sumw2();\r
944 AliDebug(2,"AliAnalysisTaskHFEFlow: cos2phiep");\r
945 fSin2phiep = new THnSparseF("sin2phiep","sin2phiep",nDime,nBine);\r
946 fSin2phiep->SetBinEdges(2,binLimPt);\r
947 fSin2phiep->SetBinEdges(3,binLimEta);\r
948 fSin2phiep->SetBinEdges(0,binLimCos);\r
949 fSin2phiep->SetBinEdges(1,binLimC);\r
950 fSin2phiep->Sumw2();\r
951 AliDebug(2,"AliAnalysisTaskHFEFlow: sin2phiep");\r
952 fSin2phiephiep = new THnSparseF("sin2phie_phiep","sin2phie_phiep",nDime,nBine);\r
953 fSin2phiephiep->SetBinEdges(2,binLimPt);\r
954 fSin2phiephiep->SetBinEdges(3,binLimEta);\r
955 fSin2phiephiep->SetBinEdges(0,binLimCos);\r
956 fSin2phiephiep->SetBinEdges(1,binLimC);\r
957 fSin2phiephiep->Sumw2(); \r
958 AliDebug(2,"AliAnalysisTaskHFEFlow: sin2phiephiep");\r
959 \r
960 const Int_t nDimfbiss=4;\r
961 Int_t nBinfbiss[nDimfbiss] = {nBinsCos,nBinsCos,nBinsCos,nBinsC};\r
962 fSinResabc = new THnSparseF("SinRes_abc","SinRes_abc",nDimfbiss,nBinfbiss);\r
963 fSinResabc->SetBinEdges(0,binLimCos);\r
964 fSinResabc->SetBinEdges(1,binLimCos);\r
965 fSinResabc->SetBinEdges(2,binLimCos);\r
966 fSinResabc->SetBinEdges(3,binLimC);\r
967 fSinResabc->Sumw2();\r
968 AliDebug(2,"AliAnalysisTaskHFEFlow: sinresabc");\r
969 \r
970 // Profile cosres centrality with 3 subevents\r
971 fProfileCosResab = new TProfile("ProfileCosRes_a_b","ProfileCosRes_a_b",nBinsCMore,binLimCMore);\r
972 fProfileCosResab->Sumw2();\r
973 fProfileCosResac = new TProfile("ProfileCosRes_a_c","ProfileCosRes_a_c",nBinsCMore,binLimCMore);\r
974 fProfileCosResac->Sumw2();\r
975 fProfileCosResbc = new TProfile("ProfileCosRes_b_c","ProfileCosRes_b_c",nBinsCMore,binLimCMore);\r
976 fProfileCosResbc->Sumw2();\r
977 AliDebug(2,"AliAnalysisTaskHFEFlow: profilecosresbc");\r
978 \r
979 //\r
980 const Int_t nDimff=2;\r
981 Int_t nBinff[nDimff] = {nBinsCos, nBinsC};\r
982 fSinRes = new THnSparseF("SinRes","SinRes",nDimff,nBinff);\r
983 fSinRes->SetBinEdges(0,binLimCos);\r
984 fSinRes->SetBinEdges(1,binLimC);\r
985 fSinRes->Sumw2();\r
986 AliDebug(2,"AliAnalysisTaskHFEFlow: sinres");\r
987 \r
988 // Profile cosres centrality\r
989 fProfileCosRes = new TProfile("ProfileCosRes","ProfileCosRes",nBinsCMore,binLimCMore);\r
990 fProfileCosRes->Sumw2();\r
991 AliDebug(2,"AliAnalysisTaskHFEFlow: profilecosres");\r
992 \r
993 // Profile Maps cos phi\r
994 fProfileCosPhiMaps = new TProfile2D("ProfileCosPhiMaps","ProfileCosPhiMaps",nBinsC,binLimC,nBinsPt,binLimPt);\r
995 fProfileCosPhiMaps->Sumw2();\r
996 AliDebug(2,"AliAnalysisTaskHFEFlow: profilecosphimaps");\r
997\r
998 }\r
999 //\r
1000 // fMonitorTrackCuts\r
1001 //\r
1002\r
1003 if(fMonitorTrackCuts) {\r
1004 // Debugging tracking steps\r
1005 const Int_t nDimTrStep=2;\r
1006 Int_t nBinTrStep[nDimTrStep] = {nBinsPt,nBinsStep};\r
1007 fTrackingCuts = new THnSparseF("TrackingCuts","TrackingCuts",nDimTrStep,nBinTrStep);\r
1008 fTrackingCuts->SetBinEdges(0,binLimPt);\r
1009 fTrackingCuts->SetBinEdges(1,binLimStep);\r
1010 fTrackingCuts->Sumw2();\r
1011 AliDebug(2,"AliAnalysisTaskHFEFlow: trackingcuts");\r
1012 }\r
1013\r
1014 //\r
1015 // fMonitorContamination\r
1016 //\r
1017\r
1018 if(fMonitorContamination) { \r
1019 // Maps delta phi contamination\r
1020 const Int_t nDimgcont=4;\r
1021 Int_t nBingcont[nDimgcont] = {nBinsPhiLess,nBinsC,nBinsPt, nBinsTPCdEdx};\r
1022 fDeltaPhiMapsContamination = new THnSparseF("DeltaPhiMapsContamination","DeltaPhiMapsContamination",nDimgcont,nBingcont);\r
1023 fDeltaPhiMapsContamination->SetBinEdges(0,binLimPhiLess);\r
1024 fDeltaPhiMapsContamination->SetBinEdges(1,binLimC);\r
1025 fDeltaPhiMapsContamination->SetBinEdges(2,binLimPt);\r
1026 fDeltaPhiMapsContamination->SetBinEdges(3,binLimTPCdEdx);\r
1027 fDeltaPhiMapsContamination->Sumw2(); \r
1028 AliDebug(2,"AliAnalysisTaskHFEFlow: deltaphimapscontamination");\r
1029\r
1030 }\r
1031 //\r
1032 // fMonitorWithoutPID\r
1033 //\r
1034\r
1035 if(fMonitorWithoutPID) {\r
1036 //\r
1037 const Int_t nDimgb=3;\r
1038 Int_t nBingb[nDimgb] = {nBinsPhi,nBinsC,nBinsPt};\r
1039 \r
1040 fDeltaPhiMapsBeforePID = new THnSparseF("DeltaPhiMapsBeforePID","DeltaPhiMapsBeforePID",nDimgb,nBingb);\r
1041 fDeltaPhiMapsBeforePID->SetBinEdges(0,binLimPhi);\r
1042 fDeltaPhiMapsBeforePID->SetBinEdges(1,binLimC);\r
1043 fDeltaPhiMapsBeforePID->SetBinEdges(2,binLimPt);\r
1044 fDeltaPhiMapsBeforePID->Sumw2(); \r
1045 AliDebug(2,"AliAnalysisTaskHFEFlow: deltaphimapsbeforepid");\r
1046 \r
1047 const Int_t nDimhb=3;\r
1048 Int_t nBinhb[nDimhb] = {nBinsCos,nBinsC,nBinsPt};\r
1049 \r
1050 fCosPhiMapsBeforePID = new THnSparseF("CosPhiMapsBeforePID","CosPhiMapsBeforePID",nDimhb,nBinhb);\r
1051 fCosPhiMapsBeforePID->SetBinEdges(0,binLimCos);\r
1052 fCosPhiMapsBeforePID->SetBinEdges(1,binLimC);\r
1053 fCosPhiMapsBeforePID->SetBinEdges(2,binLimPt);\r
1054 fCosPhiMapsBeforePID->Sumw2();\r
1055 AliDebug(2,"AliAnalysisTaskHFEFlow: cosphimapsbeforepid");\r
1056 }\r
1057 //\r
1058 // fMonitorPhotonic\r
1059 //\r
1060\r
1061 if(fMonitorPhotonic) {\r
1062 \r
1063 const Int_t nDimgbp=3;\r
1064 Int_t nBingbp[nDimgbp] = {nBinsPhi,nBinsC,nBinsPt};\r
1065 \r
1066 fDeltaPhiMapsTaggedPhotonic = new THnSparseF("DeltaPhiMapsTaggedPhotonic","DeltaPhiMapsTaggedPhotonic",nDimgbp,nBingbp);\r
1067 fDeltaPhiMapsTaggedPhotonic->SetBinEdges(0,binLimPhi);\r
1068 fDeltaPhiMapsTaggedPhotonic->SetBinEdges(1,binLimC);\r
1069 fDeltaPhiMapsTaggedPhotonic->SetBinEdges(2,binLimPt);\r
1070 fDeltaPhiMapsTaggedPhotonic->Sumw2(); \r
1071 AliDebug(2,"AliAnalysisTaskHFEFlow: deltaphimapstaggedphotonic");\r
1072 \r
1073 fDeltaPhiMapsTaggedNonPhotonic = new THnSparseF("DeltaPhiMapsTaggedNonPhotonic","DeltaPhiMapsTaggedNonPhotonic",nDimgbp,nBingbp);\r
1074 fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(0,binLimPhi);\r
1075 fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(1,binLimC);\r
1076 fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(2,binLimPt);\r
1077 fDeltaPhiMapsTaggedNonPhotonic->Sumw2(); \r
1078 AliDebug(2,"AliAnalysisTaskHFEFlow: deltaphimapstaggednonphotonic");\r
1079 \r
1080 fDeltaPhiMapsTaggedPhotonicLS = new THnSparseF("DeltaPhiMapsTaggedPhotonicLS","DeltaPhiMapsTaggedPhotonicLS",nDimgbp,nBingbp);\r
1081 fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(0,binLimPhi);\r
1082 fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(1,binLimC);\r
1083 fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(2,binLimPt);\r
1084 fDeltaPhiMapsTaggedPhotonicLS->Sumw2(); \r
1085 AliDebug(2,"AliAnalysisTaskHFEFlow: deltaphimapstaggedphotonicls"); \r
1086\r
1087 /*\r
1088 const Int_t nDimhbp=3;\r
1089 Int_t nBinhbp[nDimhbp] = {nBinsCos,nBinsC,nBinsPt};\r
1090 \r
1091 fCosPhiMapsTaggedPhotonic = new THnSparseF("CosPhiMapsTaggedPhotonic","CosPhiMapsTaggedPhotonic",nDimhbp,nBinhbp);\r
1092 fCosPhiMapsTaggedPhotonic->SetBinEdges(0,binLimCos);\r
1093 fCosPhiMapsTaggedPhotonic->SetBinEdges(1,binLimC);\r
1094 fCosPhiMapsTaggedPhotonic->SetBinEdges(2,binLimPt);\r
1095 fCosPhiMapsTaggedPhotonic->Sumw2();\r
1096 AliDebug(2,"AliAnalysisTaskHFEFlow: cosphimapstaggedphotonic");\r
1097 \r
1098 fCosPhiMapsTaggedNonPhotonic = new THnSparseF("CosPhiMapsTaggedNonPhotonic","CosPhiMapsTaggedNonPhotonic",nDimhbp,nBinhbp);\r
1099 fCosPhiMapsTaggedNonPhotonic->SetBinEdges(0,binLimCos);\r
1100 fCosPhiMapsTaggedNonPhotonic->SetBinEdges(1,binLimC);\r
1101 fCosPhiMapsTaggedNonPhotonic->SetBinEdges(2,binLimPt);\r
1102 fCosPhiMapsTaggedNonPhotonic->Sumw2();\r
1103 AliDebug(2,"AliAnalysisTaskHFEFlow: cosphimapstaggednonphotonic");\r
1104 \r
1105 fCosPhiMapsTaggedPhotonicLS = new THnSparseF("CosPhiMapsTaggedPhotonicLS","CosPhiMapsTaggedPhotonicLS",nDimhbp,nBinhbp);\r
1106 fCosPhiMapsTaggedPhotonicLS->SetBinEdges(0,binLimCos);\r
1107 fCosPhiMapsTaggedPhotonicLS->SetBinEdges(1,binLimC);\r
1108 fCosPhiMapsTaggedPhotonicLS->SetBinEdges(2,binLimPt);\r
1109 fCosPhiMapsTaggedPhotonicLS->Sumw2();\r
1110 AliDebug(2,"AliAnalysisTaskHFEFlow: cosphimapstaggedphotonicls");\r
1111 */\r
1112 const Int_t nDimMCSource=3;\r
1113 Int_t nBinMCSource[nDimMCSource] = {nBinsC,nBinsPt,nBinsSource};\r
1114 fMCSourceDeltaPhiMaps = new THnSparseF("MCSourceDeltaPhiMaps","MCSourceDeltaPhiMaps",nDimMCSource,nBinMCSource);\r
1115 fMCSourceDeltaPhiMaps->SetBinEdges(0,binLimC);\r
1116 fMCSourceDeltaPhiMaps->SetBinEdges(1,binLimPt);\r
1117 fMCSourceDeltaPhiMaps->SetBinEdges(2,binLimSource);\r
1118 fMCSourceDeltaPhiMaps->Sumw2();\r
1119 AliDebug(2,"AliAnalysisTaskHFEFlow: mcsourcedeltaphimaps");\r
1120 \r
1121 // Maps invmass opposite\r
1122 const Int_t nDimOppSign=5;\r
1123 Int_t nBinOppSign[nDimOppSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource};\r
1124 fOppSignDeltaPhiMaps = new THnSparseF("OppSignDeltaPhiMaps","OppSignDeltaPhiMaps",nDimOppSign,nBinOppSign);\r
1125 fOppSignDeltaPhiMaps->SetBinEdges(0,binLimPhi);\r
1126 fOppSignDeltaPhiMaps->SetBinEdges(1,binLimC);\r
1127 fOppSignDeltaPhiMaps->SetBinEdges(2,binLimPt);\r
1128 fOppSignDeltaPhiMaps->SetBinEdges(3,binLimInvMass);\r
1129 fOppSignDeltaPhiMaps->SetBinEdges(4,binLimSource);\r
1130 fOppSignDeltaPhiMaps->Sumw2();\r
1131 AliDebug(2,"AliAnalysisTaskHFEFlow: oppsigndeltaphimaps");\r
1132 \r
1133 // Maps invmass same sign\r
1134 const Int_t nDimSameSign=5;\r
1135 Int_t nBinSameSign[nDimSameSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource};\r
1136 fSameSignDeltaPhiMaps = new THnSparseF("SameSignDeltaPhiMaps","SameSignDeltaPhiMaps",nDimSameSign,nBinSameSign);\r
1137 fSameSignDeltaPhiMaps->SetBinEdges(0,binLimPhi);\r
1138 fSameSignDeltaPhiMaps->SetBinEdges(1,binLimC);\r
1139 fSameSignDeltaPhiMaps->SetBinEdges(2,binLimPt);\r
1140 fSameSignDeltaPhiMaps->SetBinEdges(3,binLimInvMass);\r
1141 fSameSignDeltaPhiMaps->SetBinEdges(4,binLimSource);\r
1142 fSameSignDeltaPhiMaps->Sumw2();\r
1143 AliDebug(2,"AliAnalysisTaskHFEFlow: samesigndeltaphimaps");\r
1144 \r
1145 // Maps angle same sign\r
1146 const Int_t nDimAngleSameSign=3;\r
1147 Int_t nBinAngleSameSign[nDimAngleSameSign] = {nBinsAngle,nBinsC,nBinsSource};\r
1148 fSameSignAngle = new THnSparseF("SameSignAngleMaps","SameSignAngleMaps",nDimAngleSameSign,nBinAngleSameSign);\r
1149 fSameSignAngle->SetBinEdges(0,binLimAngle);\r
1150 fSameSignAngle->SetBinEdges(1,binLimC);\r
1151 fSameSignAngle->SetBinEdges(2,binLimSource);\r
1152 fSameSignAngle->Sumw2();\r
1153 AliDebug(2,"AliAnalysisTaskHFEFlow: samesignangle");\r
1154 \r
1155 // Maps angle opp sign\r
1156 const Int_t nDimAngleOppSign=3;\r
1157 Int_t nBinAngleOppSign[nDimAngleOppSign] = {nBinsAngle,nBinsC,nBinsSource};\r
1158 fOppSignAngle = new THnSparseF("OppSignAngleMaps","OppSignAngleMaps",nDimAngleOppSign,nBinAngleOppSign);\r
1159 fOppSignAngle->SetBinEdges(0,binLimAngle);\r
1160 fOppSignAngle->SetBinEdges(1,binLimC);\r
1161 fOppSignAngle->SetBinEdges(2,binLimSource);\r
1162 fOppSignAngle->Sumw2();\r
1163 AliDebug(2,"AliAnalysisTaskHFEFlow: oppsignangle");\r
1164\r
1165 }\r
1166\r
1167 //**************************\r
1168 // Add to the list\r
1169 //******************************\r
1170\r
1171 fListHist->Add(fHistEV);\r
1172 fListHist->Add(fEventPlane);\r
1173 fListHist->Add(fFractionContamination);\r
1174 fListHist->Add(fCosRes);\r
1175 fListHist->Add(fCosResabc);\r
1176 fListHist->Add(fCosPhiMaps);\r
1177 fListHist->Add(fDeltaPhiMaps);\r
1178 fListHist->Add(fPIDqa->MakeList("HFEpidQA"));\r
1179 fListHist->Add(fContaminationv2);\r
1180 AliDebug(2,"AliAnalysisTaskHFEFlow: add default");\r
1181\r
1182 if(fMonitorEventPlane) {\r
1183 fListHist->Add(fProfileCosRes);\r
1184 fListHist->Add(fProfileCosResab);\r
1185 fListHist->Add(fProfileCosResac);\r
1186 fListHist->Add(fProfileCosResbc);\r
1187 fListHist->Add(fCosSin2phiep);\r
1188 fListHist->Add(fCos2phie);\r
1189 fListHist->Add(fSin2phie);\r
1190 fListHist->Add(fCos2phiep);\r
1191 fListHist->Add(fSin2phiep);\r
1192 fListHist->Add(fSin2phiephiep);\r
1193 fListHist->Add(fSinRes);\r
1194 fListHist->Add(fSinResabc);\r
1195 fListHist->Add(fProfileCosPhiMaps);\r
1196 }\r
1197 AliDebug(2,"AliAnalysisTaskHFEFlow: add monitor");\r
1198\r
1199 if(fMonitorTrackCuts) fListHist->Add(fTrackingCuts);\r
1200\r
1201 AliDebug(2,"AliAnalysisTaskHFEFlow: add monitortrackcuts");\r
1202\r
1203 if(fMonitorContamination) {\r
1204 fListHist->Add(fDeltaPhiMapsContamination);\r
1205 }\r
1206 \r
1207 AliDebug(2,"AliAnalysisTaskHFEFlow: add deltaphimapscontamination");\r
1208\r
1209 if(fMonitorWithoutPID) {\r
1210 fListHist->Add(fDeltaPhiMapsBeforePID);\r
1211 fListHist->Add(fCosPhiMapsBeforePID);\r
1212 }\r
1213\r
1214 AliDebug(2,"AliAnalysisTaskHFEFlow: add without pid");\r
1215\r
1216 if(fMonitorPhotonic) {\r
1217 fListHist->Add(fPIDBackgroundqa->MakeList("HFEpidBackgroundQA"));\r
1218 fListHist->Add(fDeltaPhiMapsTaggedPhotonic);\r
1219 //fListHist->Add(fCosPhiMapsTaggedPhotonic);\r
1220 fListHist->Add(fDeltaPhiMapsTaggedNonPhotonic);\r
1221 //fListHist->Add(fCosPhiMapsTaggedNonPhotonic);\r
1222 fListHist->Add(fDeltaPhiMapsTaggedPhotonicLS);\r
1223 //fListHist->Add(fCosPhiMapsTaggedPhotonicLS);\r
1224 fListHist->Add(fMCSourceDeltaPhiMaps);\r
1225 fListHist->Add(fOppSignDeltaPhiMaps);\r
1226 fListHist->Add(fSameSignDeltaPhiMaps);\r
1227 fListHist->Add(fSameSignAngle);\r
1228 fListHist->Add(fOppSignAngle);\r
1229 fListHist->Add(fBackgroundSubtraction->GetListOutput());\r
1230 }\r
1231\r
1232 AliDebug(2,"AliAnalysisTaskHFEFlow: add photonic");\r
1233\r
1234 if(fHFEVZEROEventPlane && fMonitorEventPlane) fListHist->Add(fHFEVZEROEventPlane->GetOutputList());\r
1235 \r
1236 AliDebug(2,"AliAnalysisTaskHFEFlow: add event plane");\r
1237\r
1238 PostData(1, fListHist);\r
1239 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {\r
1240 // PostData(bincless+2,fflowEvent); \r
1241 //}\r
1242\r
1243 AliDebug(2,"AliAnalysisTaskHFEFlow: post");\r
1244\r
1245\r
1246}\r
1247 \r
1248//________________________________________________________________________\r
1249void AliAnalysisTaskHFEFlow::UserExec(Option_t */*option*/)\r
1250{\r
1251 //\r
1252 // Loop over event\r
1253 //\r
1254 \r
1255 Double_t massElectron = 0.000511;\r
1256 Double_t mcReactionPlane = 0.0;\r
1257 Bool_t eventplanedefined = kTRUE;\r
1258\r
1259 Float_t cntr = 0.0;\r
1260 Double_t binct = 11.5;\r
1261 Double_t binctMore = 20.5;\r
1262 Double_t binctLess = -0.5;\r
1263 Float_t binctt = -1.0;\r
1264 \r
1265 Double_t valuecossinephiep[3];\r
1266 Double_t valuensparsea[4];\r
1267 Double_t valuensparseabis[5];\r
1268 Double_t valuensparsee[4];\r
1269 Double_t valuensparsef[2];\r
1270 Double_t valuensparsefsin[2];\r
1271 Double_t valuensparsefbis[4];\r
1272 Double_t valuensparsefbissin[4];\r
1273 Double_t valuensparseg[5];\r
1274 Double_t valuensparseh[5];\r
1275 Double_t valuensparsehprofile[3];\r
1276 Double_t valuensparseMCSourceDeltaPhiMaps[3];\r
1277 Double_t valuetrackingcuts[2];\r
1278 Double_t valuedeltaphicontamination[4];\r
1279 Double_t valuefractioncont[2];\r
1280 \r
1281 AliMCEvent *mcEvent = MCEvent();\r
1282 AliMCParticle *mctrack = NULL;\r
1283\r
1284 // MC info\r
1285 Bool_t mcthere = kTRUE;\r
1286 if(fAODAnalysis) {\r
1287 AliAODEvent *aodE = dynamic_cast<AliAODEvent *>(fInputEvent);\r
1288 if(!aodE){\r
1289 // printf("testd\n");\r
1290 AliError("No AOD Event");\r
1291 return;\r
1292 }\r
1293 fAODMCHeader = dynamic_cast<AliAODMCHeader *>(fInputEvent->FindListObject(AliAODMCHeader::StdBranchName()));\r
1294 if(!fAODMCHeader){ \r
1295 mcthere = kFALSE;\r
1296 }\r
1297 fAODArrayMCInfo = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));\r
1298 if(!fAODArrayMCInfo){ \r
1299 mcthere = kFALSE;\r
1300 }\r
1301 else {\r
1302 fHFECuts->SetMCEvent(aodE);\r
1303 if(fMonitorPhotonic) fBackgroundSubtraction->SetAODArrayMCInfo(fAODArrayMCInfo);\r
1304 }\r
1305 }\r
1306 else {\r
1307 AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());\r
1308 if(!mcH) mcthere = kFALSE;\r
1309 else {\r
1310 if(fMonitorPhotonic) fBackgroundSubtraction->SetMCEvent(fMCEvent);\r
1311 }\r
1312 }\r
1313\r
1314 \r
1315 /////////////////\r
1316 // centrality\r
1317 /////////////////\r
1318 \r
1319 //AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
1320 //if(!esd) return;\r
1321 AliCentrality *centrality = fInputEvent->GetCentrality();\r
1322 //AliDebug(2,"Got the centrality");\r
1323 if(!centrality) return;\r
1324 cntr = centrality->GetCentralityPercentile("V0M");\r
1325 if((0.0< cntr) && (cntr<5.0)) binct = 0.5;\r
1326 if((5.0< cntr) && (cntr<10.0)) binct = 1.5;\r
1327 if((10.0< cntr) && (cntr<20.0)) binct = 2.5;\r
1328 if((20.0< cntr) && (cntr<30.0)) binct = 3.5;\r
1329 if((30.0< cntr) && (cntr<40.0)) binct = 4.5;\r
1330 if((40.0< cntr) && (cntr<50.0)) binct = 5.5;\r
1331 if((50.0< cntr) && (cntr<60.0)) binct = 6.5;\r
1332 if((60.0< cntr) && (cntr<70.0)) binct = 7.5;\r
1333 if((70.0< cntr) && (cntr<80.0)) binct = 8.5;\r
1334 if((80.0< cntr) && (cntr<90.0)) binct = 9.5;\r
1335 if((90.0< cntr) && (cntr<100.0)) binct = 10.5;\r
1336 \r
1337 if((0.< cntr) && (cntr < 20.)) binctt = 0.5;\r
1338 if((20.< cntr) && (cntr < 40.)) binctt = 1.5;\r
1339 if((40.< cntr) && (cntr < 80.)) binctt = 2.5;\r
1340\r
1341 if((0.0< cntr) && (cntr<5.0)) binctMore = 0.5;\r
1342 if((5.0< cntr) && (cntr<10.0)) binctMore = 1.5;\r
1343 if((10.0< cntr) && (cntr<15.0)) binctMore = 2.5;\r
1344 if((15.0< cntr) && (cntr<20.0)) binctMore = 3.5;\r
1345 if((20.0< cntr) && (cntr<25.0)) binctMore = 4.5;\r
1346 if((25.0< cntr) && (cntr<30.0)) binctMore = 5.5;\r
1347 if((30.0< cntr) && (cntr<35.0)) binctMore = 6.5;\r
1348 if((35.0< cntr) && (cntr<40.0)) binctMore = 7.5;\r
1349 if((40.0< cntr) && (cntr<45.0)) binctMore = 8.5;\r
1350 if((45.0< cntr) && (cntr<50.0)) binctMore = 9.5;\r
1351 if((50.0< cntr) && (cntr<55.0)) binctMore = 10.5;\r
1352 if((55.0< cntr) && (cntr<60.0)) binctMore = 11.5;\r
1353 if((60.0< cntr) && (cntr<65.0)) binctMore = 12.5;\r
1354 if((65.0< cntr) && (cntr<70.0)) binctMore = 13.5;\r
1355 if((70.0< cntr) && (cntr<75.0)) binctMore = 14.5;\r
1356 if((75.0< cntr) && (cntr<80.0)) binctMore = 15.5;\r
1357 if((80.0< cntr) && (cntr<85.0)) binctMore = 16.5;\r
1358 if((85.0< cntr) && (cntr<90.0)) binctMore = 17.5;\r
1359 if((90.0< cntr) && (cntr<95.0)) binctMore = 18.5;\r
1360 if((95.0< cntr) && (cntr<100.0)) binctMore = 19.5;\r
1361\r
1362 binctLess = cntr;\r
1363 \r
1364 \r
1365 if(binct > 11.0) return;\r
1366 \r
1367 // centrality\r
1368 valuensparsea[3] = binct; \r
1369 valuensparseabis[1] = binct; \r
1370 valuensparsee[1] = binct; \r
1371 valuensparsef[1] = binctMore; \r
1372 valuensparsefsin[1] = binct; \r
1373 valuensparsefbis[3] = binctMore; \r
1374 valuensparsefbissin[3] = binct; \r
1375 valuensparseg[1] = binct;\r
1376 valuensparseh[1] = binct; \r
1377 valuefractioncont[1] = binct;\r
1378 valuensparsehprofile[1] = binct; \r
1379 valuecossinephiep[2] = binctMore;\r
1380 valuensparseMCSourceDeltaPhiMaps[0] = binct;\r
1381 valuedeltaphicontamination[1] = binct;\r
1382 \r
1383 //////////////////////\r
1384 // run number\r
1385 //////////////////////\r
1386\r
1387 Int_t runnumber = fInputEvent->GetRunNumber();\r
1388 AliDebug(2,Form("Run number %d",runnumber));\r
1389 \r
1390 if(!fPID->IsInitialized()){\r
1391 // Initialize PID with the given run number\r
1392 fPID->InitializePID(runnumber);\r
1393 }\r
1394 if(!fPIDTOFOnly->IsInitialized()){\r
1395 // Initialize PID with the given run number\r
1396 fPIDTOFOnly->InitializePID(runnumber);\r
1397 }\r
1398\r
1399 //\r
1400 if(!fPIDBackground->IsInitialized()){\r
1401 // Initialize PID with the given run number\r
1402 fPIDBackground->InitializePID(runnumber);\r
1403 }\r
1404\r
1405 fHFECuts->SetRecEvent(fInputEvent);\r
1406 if(mcEvent) fHFECuts->SetMCEvent(mcEvent);\r
1407\r
1408\r
1409 //////////\r
1410 // PID\r
1411 //////////\r
1412 \r
1413 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();\r
1414 if(!pidResponse){\r
1415 AliDebug(2,"No PID response set");\r
1416 return;\r
1417 }\r
1418 fPID->SetPIDResponse(pidResponse);\r
1419 fPIDTOFOnly->SetPIDResponse(pidResponse);\r
1420 fPIDBackground->SetPIDResponse(pidResponse);\r
1421 if(fMonitorPhotonic) fBackgroundSubtraction->InitRun(fInputEvent,pidResponse);\r
1422\r
1423 fHistEV->Fill(binctt,0.0);\r
1424\r
1425 //////////////////\r
1426 // Event cut\r
1427 //////////////////\r
1428 if(!fHFECuts->CheckEventCuts("fEvRecCuts", fInputEvent)) {\r
1429 AliDebug(2,"Do not pass the event cut");\r
1430 PostData(1, fListHist);\r
1431 return;\r
1432 }\r
1433\r
1434 fHistEV->Fill(binctt,1.0);\r
1435\r
1436 //////////////////////////////////// \r
1437 // First method event plane\r
1438 ////////////////////////////////////\r
1439\r
1440 AliEventplane* vEPa = fInputEvent->GetEventplane();\r
1441 Float_t eventPlanea = 0.0;\r
1442 Float_t eventPlaneTPC = 0.0;\r
1443 Float_t eventPlaneV0A = 0.0;\r
1444 Float_t eventPlaneV0C = 0.0;\r
1445 Float_t eventPlaneV0 = 0.0;\r
1446 TVector2 *qTPC = 0x0;\r
1447 TVector2 *qsub1a = 0x0;\r
1448 TVector2 *qsub2a = 0x0;\r
1449 TVector2 qV0A,qV0C,qV0,*qAna;\r
1450 \r
1451 // V0\r
1452\r
1453 if(fHFEVZEROEventPlane && (!fAODAnalysis)){\r
1454\r
1455 //AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());\r
1456 //if(!esd) return;\r
1457 \r
1458 fHFEVZEROEventPlane->ProcessEvent(fInputEvent);\r
1459 \r
1460 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0A()+100) < 0.0000001) eventPlaneV0A = -100.0;\r
1461 else {\r
1462 eventPlaneV0A = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0A());\r
1463 if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();\r
1464 }\r
1465 \r
1466 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0C()+100) < 0.0000001) eventPlaneV0C = -100.0;\r
1467 else {\r
1468 eventPlaneV0C = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0C());\r
1469 if(eventPlaneV0C > TMath::Pi()) eventPlaneV0C = eventPlaneV0C - TMath::Pi();\r
1470 }\r
1471\r
1472 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0()+100) < 0.0000001) eventPlaneV0 = -100.0;\r
1473 else {\r
1474 eventPlaneV0 = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0());\r
1475 if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();\r
1476 }\r
1477 \r
1478 }\r
1479 else {\r
1480 \r
1481 Double_t qVx, qVy; //TR: info\r
1482 eventPlaneV0 = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,10,2,qVx,qVy));\r
1483 if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();\r
1484 qV0.Set(qVx,qVy);\r
1485 eventPlaneV0A = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,8,2,qVx,qVy));\r
1486 if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();\r
1487 qV0A.Set(qVx,qVy);\r
1488 eventPlaneV0C = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,9,2,qVx,qVy));\r
1489 if(eventPlaneV0C > TMath::Pi()) eventPlaneV0C = eventPlaneV0C - TMath::Pi();\r
1490 qV0C.Set(qVx,qVy);\r
1491 \r
1492 }\r
1493\r
1494 // TPC\r
1495\r
1496 qTPC = vEPa->GetQVector(); \r
1497 Double_t qx = -1.0;\r
1498 Double_t qy = -1.0;\r
1499 if(qTPC) {\r
1500 qx = qTPC->X();\r
1501 qy = qTPC->Y();\r
1502 } \r
1503 TVector2 qVectorfortrack;\r
1504 qVectorfortrack.Set(qx,qy);\r
1505 eventPlaneTPC = TVector2::Phi_0_2pi(qVectorfortrack.Phi())/2.; \r
1506\r
1507 // Choose the one used for v2\r
1508\r
1509 if(fVZEROEventPlane){ //TR: info\r
1510 eventPlanea = eventPlaneV0;\r
1511 qAna = &qV0;\r
1512 }\r
1513 if(fVZEROEventPlaneA){\r
1514 eventPlanea = eventPlaneV0A;\r
1515 qAna = &qV0A;\r
1516 }\r
1517 if(fVZEROEventPlaneC){\r
1518 eventPlanea = eventPlaneV0C;\r
1519 qAna = &qV0C;\r
1520 }\r
1521 if(!fVZEROEventPlane){\r
1522 eventPlanea = eventPlaneTPC;\r
1523 qAna = &qV0C;\r
1524 }\r
1525\r
1526 valuecossinephiep[0] = TMath::Cos(2*eventPlanea);\r
1527 valuecossinephiep[1] = TMath::Sin(2*eventPlanea);\r
1528\r
1529 Float_t eventPlanesub1a = -100.0;\r
1530 Float_t eventPlanesub2a = -100.0;\r
1531 Double_t diffsub1sub2a = -100.0;\r
1532 Double_t diffsub1sub2asin = -100.0;\r
1533 Double_t diffsubasubb = -100.0;\r
1534 Double_t diffsubasubc = -100.0;\r
1535 Double_t diffsubbsubc = -100.0;\r
1536 Double_t diffsubasubbsin = -100.0;\r
1537 Double_t diffsubasubcsin = -100.0;\r
1538 Double_t diffsubbsubcsin = -100.0;\r
1539\r
1540 // two sub event TPC\r
1541 qsub1a = vEPa->GetQsub1();\r
1542 qsub2a = vEPa->GetQsub2();\r
1543 if(qsub1a) eventPlanesub1a = TVector2::Phi_0_2pi(qsub1a->Phi())/2.;\r
1544 if(qsub2a) eventPlanesub2a = TVector2::Phi_0_2pi(qsub2a->Phi())/2.;\r
1545 if(qsub1a && qsub2a) {\r
1546 diffsub1sub2a = TMath::Cos(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));\r
1547 diffsub1sub2asin = TMath::Sin(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));\r
1548 }\r
1549\r
1550 // if ( !fDebugStreamer ) {\r
1551 // //debug stream\r
1552 // TDirectory *backup = gDirectory;\r
1553 // fDebugStreamer = new TTreeSRedirector("TaskHFEFlowdebug.root");\r
1554 // if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer\r
1555 // } \r
1556\r
1557 // {\r
1558\r
1559 // double v0nrom = TMath::Sqrt(qV0.X()*qV0.X()+qV0.Y()*qV0.Y());\r
1560 // double v0Anrom = TMath::Sqrt(qV0A.X()*qV0A.X()+qV0A.Y()*qV0A.Y());\r
1561 // double v0Cnrom = TMath::Sqrt(qV0C.X()*qV0C.X()+qV0C.Y()*qV0C.Y());\r
1562 // double sub1nrom = TMath::Sqrt(qsub1a->X()*qsub1a->X()+qsub1a->Y()*qsub1a->Y());\r
1563 // double sub2nrom = TMath::Sqrt(qsub2a->X()*qsub2a->X()+qsub2a->Y()*qsub2a->Y());\r
1564\r
1565 // (* fDebugStreamer) << "UserExec" <<\r
1566 // "binct="<<binct<<\r
1567 // "qV0="<<v0nrom<<\r
1568 // "qV0A="<<v0Anrom<<\r
1569 // "qV0C="<<v0Cnrom<<\r
1570 // "qsub1a="<<sub1nrom<<\r
1571 // "qsub2a="<<sub2nrom<<\r
1572 // "\n";\r
1573 // }\r
1574\r
1575 // three sub events in case of VZEROA and VZEROC\r
1576 if(!fSP){\r
1577 diffsubasubb = TMath::Cos(2.*(eventPlaneV0A - eventPlaneV0C)); //TR: \r
1578 diffsubasubc = TMath::Cos(2.*(eventPlaneV0A - eventPlaneTPC)); //TR: \r
1579 diffsubbsubc = TMath::Cos(2.*(eventPlaneV0C - eventPlaneTPC)); //TR: \r
1580 }\r
1581 else{\r
1582 if(fVZEROEventPlaneA){\r
1583 diffsubasubb = qV0A.X()*qV0C.X()+qV0A.Y()*qV0C.Y();\r
1584 diffsubasubc = qV0A.X()*qTPC->X()+qV0A.Y()*qTPC->Y();\r
1585 diffsubbsubc = qV0C.X()*qTPC->X()+qV0C.Y()*qTPC->Y();\r
1586 }\r
1587 else if(fVZEROEventPlaneC){\r
1588 diffsubasubb = qV0C.X()*qV0A.X()+qV0C.Y()*qV0A.Y();\r
1589 diffsubasubc = qV0C.X()*qTPC->X()+qV0C.Y()*qTPC->Y();\r
1590 diffsubbsubc = qV0A.X()*qTPC->X()+qV0A.Y()*qTPC->Y();\r
1591 }\r
1592 }\r
1593\r
1594 diffsubasubbsin = TMath::Sin(2.*(eventPlaneV0A - eventPlaneV0C));\r
1595 diffsubasubcsin = TMath::Sin(2.*(eventPlaneV0A - eventPlaneTPC));\r
1596 diffsubbsubcsin = TMath::Sin(2.*(eventPlaneV0C - eventPlaneTPC));\r
1597 // three sub events in case of VZERO all\r
1598 if(fVZEROEventPlane && (!fVZEROEventPlaneA) && (!fVZEROEventPlaneC)) {\r
1599 if(!fSP){\r
1600 diffsubasubb = TMath::Cos(2.*(eventPlaneV0 - eventPlanesub1a)); //TR: \r
1601 diffsubasubc = TMath::Cos(2.*(eventPlaneV0 - eventPlanesub2a)); //TR: \r
1602 diffsubbsubc = TMath::Cos(2.*(eventPlanesub1a - eventPlanesub2a)); //TR: \r
1603 }\r
1604 else{\r
1605 diffsubasubb = qV0.X()*qsub1a->X()+qV0.Y()*qsub1a->Y(); \r
1606 diffsubasubc = qV0.X()*qsub2a->X()+qV0.Y()*qsub2a->Y(); \r
1607 diffsubbsubc = qsub1a->X()*qsub2a->X()+qsub1a->Y()*qsub2a->Y();\r
1608 }\r
1609 \r
1610 diffsubasubbsin = TMath::Sin(2.*(eventPlaneV0 - eventPlanesub1a));\r
1611 diffsubasubcsin = TMath::Sin(2.*(eventPlaneV0 - eventPlanesub2a));\r
1612 diffsubbsubcsin = TMath::Sin(2.*(eventPlanesub1a - eventPlanesub2a));\r
1613 }\r
1614 \r
1615\r
1616 /////////////////////////////////////////////////////////\r
1617 // Cut for event with event plane reconstructed by all\r
1618 ////////////////////////////////////////////////////////\r
1619 \r
1620 //if(!fVZEROEventPlane) {\r
1621 // if(!qTPC) {\r
1622 // eventplanedefined = kFALSE;\r
1623 //PostData(1, fListHist);\r
1624 //return;\r
1625 //}\r
1626 //}\r
1627\r
1628 if((!qTPC) || (!qsub1a) || (!qsub2a)) {\r
1629 AliDebug(2,"No event plane");\r
1630 return;\r
1631 }\r
1632\r
1633 //////////////////////////////////////\r
1634 // AliFlowEvent and MC event plane\r
1635 /////////////////////////////////////\r
1636\r
1637 Int_t nbtracks = fInputEvent->GetNumberOfTracks();\r
1638 AliDebug(2,Form("Number of tracks %d",nbtracks));\r
1639\r
1640 if(fMonitorQCumulant) {\r
1641\r
1642 fcutsRP->SetEvent( InputEvent(), MCEvent());\r
1643 fcutsPOI->SetEvent( InputEvent(), MCEvent());\r
1644 if( fflowEvent ){ \r
1645 fflowEvent->~AliFlowEvent();\r
1646 new(fflowEvent) AliFlowEvent(fcutsRP,fcutsPOI);\r
1647 }else fflowEvent = new AliFlowEvent(fcutsRP,fcutsPOI);\r
1648 if(mcEvent && mcEvent->GenEventHeader()) {\r
1649 fflowEvent->SetMCReactionPlaneAngle(mcEvent);\r
1650 //if reaction plane not set from elsewhere randomize it before adding flow\r
1651 //if (!fflowEvent->IsSetMCReactionPlaneAngle()) fflowEvent->SetMCReactionPlaneAngle(gRandom->Uniform(0.0,TMath::TwoPi()));\r
1652 mcReactionPlane = TVector2::Phi_0_2pi(fflowEvent->GetMCReactionPlaneAngle());\r
1653 if(mcReactionPlane > TMath::Pi()) mcReactionPlane = mcReactionPlane - TMath::Pi();\r
1654 AliDebug(2,Form("MC reaction plane %f",mcReactionPlane));\r
1655 }\r
1656 fflowEvent->SetReferenceMultiplicity( nbtracks );\r
1657 fflowEvent->DefineDeadZone(0,0,0,0);\r
1658 //fflowEvent.TagSubeventsInEta(-0.8,-0.1,0.1,0.8);\r
1659 \r
1660 ////////////////\r
1661 // MC\r
1662 ///////////////\r
1663 if(fUseMCReactionPlane) {\r
1664 eventPlanea = mcReactionPlane;\r
1665 diffsub1sub2a = 0.0;\r
1666 }\r
1667 }\r
1668\r
1669 \r
1670 //////////////////////\r
1671 // Fill Histos\r
1672 //////////////////////\r
1673\r
1674 if(eventplanedefined) {\r
1675 \r
1676 fHistEV->Fill(binctt,2.0);\r
1677 \r
1678 // Fill\r
1679 valuensparsea[0] = eventPlaneV0A;\r
1680 valuensparsea[1] = eventPlaneV0C;\r
1681 valuensparsea[2] = eventPlaneTPC;\r
1682 if(fVZEROEventPlane && (!fVZEROEventPlaneA) && (!fVZEROEventPlaneC)) {\r
1683 // case VZERO all\r
1684 valuensparsea[0] = eventPlaneV0;\r
1685 valuensparsea[1] = eventPlanesub1a;\r
1686 valuensparsea[2] = eventPlanesub2a;\r
1687 } \r
1688 fEventPlane->Fill(&valuensparsea[0]);\r
1689\r
1690 // Fill\r
1691 if(fMonitorEventPlane) fCosSin2phiep->Fill(&valuecossinephiep[0]);\r
1692 \r
1693 if(!fVZEROEventPlane) {\r
1694 valuensparsef[0] = diffsub1sub2a;\r
1695 fCosRes->Fill(&valuensparsef[0]);\r
1696 valuensparsefsin[0] = diffsub1sub2asin;\r
1697 if(fMonitorEventPlane) fSinRes->Fill(&valuensparsefsin[0]);\r
1698 if(fMonitorEventPlane) {\r
1699 fProfileCosRes->Fill(valuensparsef[1],valuensparsef[0]);\r
1700 }\r
1701 }\r
1702 else {\r
1703 valuensparsefbis[0] = diffsubasubb;\r
1704 valuensparsefbis[1] = diffsubasubc; //TR: WTF!!!\r
1705 valuensparsefbis[2] = diffsubbsubc;\r
1706 fCosResabc->Fill(&valuensparsefbis[0]); //TR: info\r
1707 valuensparsefbissin[0] = diffsubasubbsin;\r
1708 valuensparsefbissin[1] = diffsubbsubcsin;\r
1709 valuensparsefbissin[2] = diffsubasubcsin;\r
1710 if(fMonitorEventPlane) fSinResabc->Fill(&valuensparsefbissin[0]);\r
1711 if(fMonitorEventPlane) {\r
1712 fProfileCosResab->Fill(valuensparsefbis[3],valuensparsefbis[0]);\r
1713 fProfileCosResac->Fill(valuensparsefbis[3],valuensparsefbis[1]); //TR: WTF!!!\r
1714 fProfileCosResbc->Fill(valuensparsefbis[3],valuensparsefbis[2]);\r
1715 }\r
1716 }\r
1717 \r
1718 }\r
1719 \r
1720 ////////////////////////////////////////\r
1721 // Loop to determine pool background\r
1722 /////////////////////////////////////////\r
1723 if(fMonitorPhotonic) {\r
1724\r
1725 fBackgroundSubtraction->FillPoolAssociatedTracks(fInputEvent,binct);\r
1726 \r
1727 if( fArraytrack ){ \r
1728 fArraytrack->~TArrayI();\r
1729 new(fArraytrack) TArrayI(nbtracks);\r
1730 }\r
1731 else { \r
1732 fArraytrack = new TArrayI(nbtracks);\r
1733 }\r
1734 fCounterPoolBackground = 0;\r
1735 \r
1736 for(Int_t k = 0; k < nbtracks; k++){\r
1737 \r
1738 AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);\r
1739 if(!track) continue;\r
1740 \r
1741 // Track cuts\r
1742 Bool_t survivedbackground = kTRUE;\r
1743 if(fAODAnalysis) {\r
1744 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);\r
1745 if(aodtrack) {\r
1746 AliESDtrack esdTrack(aodtrack);\r
1747 // set the TPC cluster info\r
1748 esdTrack.SetTPCClusterMap(aodtrack->GetTPCClusterMap());\r
1749 esdTrack.SetTPCSharedMap(aodtrack->GetTPCSharedMap());\r
1750 esdTrack.SetTPCPointsF(aodtrack->GetTPCNclsF());\r
1751 // needed to calculate the impact parameters\r
1752 AliAODEvent *aodeventu = dynamic_cast<AliAODEvent *>(fInputEvent);\r
1753 if(aodeventu) {\r
1754 AliAODVertex *vAOD = aodeventu->GetPrimaryVertex();\r
1755 Double_t bfield = aodeventu->GetMagneticField();\r
1756 Double_t pos[3],cov[6];\r
1757 vAOD->GetXYZ(pos);\r
1758 vAOD->GetCovarianceMatrix(cov);\r
1759 const AliESDVertex vESD(pos,cov,100.,100);\r
1760 esdTrack.RelateToVertex(&vESD,bfield,3.);\r
1761 } \r
1762 if(!fHFEBackgroundCuts->IsSelected(&esdTrack)) {\r
1763 survivedbackground = kFALSE;\r
1764 }\r
1765 }\r
1766 }\r
1767 else {\r
1768 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);\r
1769 if(esdtrack) {\r
1770 if(!fHFEBackgroundCuts->IsSelected(esdtrack)) survivedbackground = kFALSE;\r
1771 }\r
1772 }\r
1773 // PID\r
1774 if(survivedbackground) {\r
1775 // PID track cuts\r
1776 AliHFEpidObject hfetrack2;\r
1777 if(!fAODAnalysis) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);\r
1778 else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);\r
1779 hfetrack2.SetRecTrack(track);\r
1780 hfetrack2.SetCentrality((Int_t)binct);\r
1781 AliDebug(2,Form("centrality %f and %d",binct,hfetrack2.GetCentrality()));\r
1782 hfetrack2.SetPbPb();\r
1783 if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundqa)) {\r
1784 fArraytrack->AddAt(k,fCounterPoolBackground);\r
1785 fCounterPoolBackground++;\r
1786 AliDebug(2,Form("fCounterPoolBackground %d, track %d",fCounterPoolBackground,k));\r
1787 }\r
1788 }\r
1789 }\r
1790 }\r
1791\r
1792\r
1793 ///////////////////////////////////\r
1794 // Loop for kink mother AOD\r
1795 //////////////////////////////////\r
1796 // Look at kink mother in case of AOD\r
1797 Int_t numberofvertices = 1;\r
1798 AliAODEvent *aodevent = NULL;\r
1799 Int_t numberofmotherkink = 0; \r
1800 if(fAODAnalysis) {\r
1801 aodevent = dynamic_cast<AliAODEvent *>(fInputEvent);\r
1802 if(aodevent) {\r
1803 numberofvertices = aodevent->GetNumberOfVertices();\r
1804 }\r
1805 }\r
1806 Double_t listofmotherkink[numberofvertices];\r
1807 if(fAODAnalysis && aodevent) {\r
1808 for(Int_t ivertex=0; ivertex < numberofvertices; ivertex++) {\r
1809 AliAODVertex *aodvertex = aodevent->GetVertex(ivertex);\r
1810 if(!aodvertex) continue;\r
1811 if(aodvertex->GetType()==AliAODVertex::kKink) {\r
1812 AliAODTrack *mother = (AliAODTrack *) aodvertex->GetParent();\r
1813 if(!mother) continue;\r
1814 Int_t idmother = mother->GetID();\r
1815 listofmotherkink[numberofmotherkink] = idmother;\r
1816 AliDebug(2,Form("ID %d",idmother));\r
1817 numberofmotherkink++;\r
1818 }\r
1819 }\r
1820 }\r
1821 \r
1822 //////////////////////////\r
1823 // Loop over track\r
1824 //////////////////////////\r
1825 \r
1826 for(Int_t k = 0; k < nbtracks; k++){\r
1827 \r
1828 AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);\r
1829 if(!track) continue;\r
1830 \r
1831 if(fAODAnalysis) {\r
1832 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);\r
1833 if(!aodtrack){\r
1834 AliError("AOD track is not there");\r
1835 continue;\r
1836 } \r
1837 AliDebug(2,"Find AOD track on");\r
1838 if(fUseFilterAOD){\r
1839 if(!(aodtrack->TestFilterBit(fFilter))) continue; // Only process AOD tracks where the HFE is set\r
1840 }\r
1841 }\r
1842\r
1843 if(fApplyCut) {\r
1844\r
1845 valuetrackingcuts[0] = track->Pt(); \r
1846 valuetrackingcuts[1] = 0;\r
1847 \r
1848 // RecKine: ITSTPC cuts \r
1849 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
1850 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);\r
1851\r
1852 // Reject kink mother\r
1853 if(fRejectKinkMother) {\r
1854 if(fAODAnalysis) {\r
1855 Bool_t kinkmotherpass = kTRUE;\r
1856 for(Int_t kinkmother = 0; kinkmother < numberofmotherkink; kinkmother++) {\r
1857 if(track->GetID() == listofmotherkink[kinkmother]) {\r
1858 kinkmotherpass = kFALSE;\r
1859 continue;\r
1860 }\r
1861 }\r
1862 if(!kinkmotherpass) continue;\r
1863 }\r
1864 else {\r
1865 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);\r
1866 if(esdtrack){ \r
1867 if(esdtrack->GetKinkIndex(0) != 0) continue; \r
1868 } // Quick and dirty fix to reject both kink mothers and daughters\r
1869 }\r
1870 } \r
1871\r
1872 valuetrackingcuts[1] = 1; \r
1873 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
1874 // RecPrim\r
1875 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
1876 valuetrackingcuts[1] = 2; \r
1877 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
1878 \r
1879 // HFEcuts: ITS layers cuts\r
1880 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
1881 valuetrackingcuts[1] = 3; \r
1882 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
1883\r
1884 // HFE cuts: TOF PID and mismatch flag\r
1885 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTOF + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
1886 valuetrackingcuts[1] = 4; \r
1887 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
1888 \r
1889 // HFE cuts: TPC PID cleanup\r
1890 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
1891 valuetrackingcuts[1] = 5; \r
1892 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
1893 \r
1894 // HFEcuts: Nb of tracklets TRD0\r
1895 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;\r
1896 valuetrackingcuts[1] = 6; \r
1897 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]); \r
1898 \r
1899 }\r
1900 AliDebug(2,"Survived");\r
1901\r
1902 /////////////////////////////////////////////////////////\r
1903 // Subtract candidate from TPC event plane\r
1904 ////////////////////////////////////////////////////////\r
1905 Float_t eventplanesubtracted = 0.0; \r
1906\r
1907 //if(eventplanedefined && (!fVZEROEventPlane)) {\r
1908 if(!fVZEROEventPlane) {\r
1909 // Subtract the tracks from the event plane\r
1910 Double_t qX = qTPC->X() - vEPa->GetQContributionX(track); //Modify the components: subtract the track you want to look at with your analysis\r
1911 Double_t qY = qTPC->Y() - vEPa->GetQContributionY(track); //Modify the components: subtract the track you want to look at with your analysis\r
1912 TVector2 newQVectorfortrack;\r
1913 newQVectorfortrack.Set(qX,qY);\r
1914 eventplanesubtracted = TVector2::Phi_0_2pi(newQVectorfortrack.Phi())/2; \r
1915 }\r
1916 else eventplanesubtracted = eventPlanea;\r
1917\r
1918 ///////////////////////////////////////////\r
1919 // Event plane\r
1920 //////////////////////////////////////////\r
1921 Bool_t fillEventPlane = kTRUE;\r
1922 if(!fVZEROEventPlane){\r
1923 //if((!qsub1a) || (!qsub2a) || (!eventplanedefined)) fillEventPlane = kFALSE;\r
1924 if((!qsub1a) || (!qsub2a)) fillEventPlane = kFALSE;\r
1925 if(fSubEtaGapTPC) {\r
1926 if(track->Eta() < (- fEtaGap/2.)) eventplanesubtracted = eventPlanesub1a;\r
1927 else if(track->Eta() > (fEtaGap/2.)) eventplanesubtracted = eventPlanesub2a;\r
1928 else fillEventPlane = kFALSE;\r
1929 }\r
1930 }\r
1931\r
1932 ///////////////\r
1933 // MC\r
1934 //////////////\r
1935 if(fUseMCReactionPlane) {\r
1936 eventplanesubtracted = mcReactionPlane;\r
1937 fillEventPlane = kTRUE;\r
1938 }\r
1939 \r
1940 //////////////////////////////////////////////////////////////////////////////\r
1941 ///////////////////////////AFTERBURNER\r
1942 Double_t phitrack = track->Phi(); \r
1943 if (fAfterBurnerOn)\r
1944 {\r
1945 phitrack = GetPhiAfterAddV2(track->Phi(),mcReactionPlane);\r
1946 }\r
1947 //////////////////////////////////////////////////////////////////////////////\r
1948\r
1949\r
1950 ///////////////////////\r
1951 // Calculate deltaphi\r
1952 ///////////////////////\r
1953 \r
1954 // Suppose phi track is between 0.0 and phi\r
1955 Double_t deltaphi = TVector2::Phi_0_2pi(phitrack - eventplanesubtracted);\r
1956 if(deltaphi > TMath::Pi()) deltaphi = deltaphi - TMath::Pi();\r
1957\r
1958 ////////////////////////////////\r
1959 // Determine the deltaphi bin\r
1960 ///////////////////////////////\r
1961\r
1962 // in-plane\r
1963 if(((deltaphi<(TMath::Pi()/4.)) && (deltaphi>0.0)) || ((deltaphi>(3*TMath::Pi()/4.)) && (deltaphi<TMath::Pi()))) valuedeltaphicontamination[0] = 0.5;\r
1964 // out-of-plane\r
1965 if((deltaphi>(TMath::Pi()/4.)) && (deltaphi<(3*TMath::Pi()/4.))) valuedeltaphicontamination[0] = 1.5;\r
1966\r
1967 ////////////////////////////////////////\r
1968 // Define variables\r
1969 ///////////////////////////////////////\r
1970\r
1971\r
1972 valuedeltaphicontamination[2] = track->Pt();\r
1973 valuensparsee[2] = track->Pt();\r
1974 valuensparsee[3] = track->Eta(); \r
1975 valuensparseg[2] = track->Pt();\r
1976 valuensparseh[2] = track->Pt();\r
1977 valuefractioncont[0] = track->Pt();\r
1978 valuensparsehprofile[2] = track->Pt();\r
1979 valuensparseMCSourceDeltaPhiMaps[1] = track->Pt();\r
1980 if(track->Charge() > 0.0) {\r
1981 valuensparseg[3] = 0.2;\r
1982 valuensparseh[3] = 0.2;\r
1983 }\r
1984 else {\r
1985 valuensparseg[3] = -0.2;\r
1986 valuensparseh[3] = -0.2;\r
1987 }\r
1988 valuensparseh[4] = track->Eta();\r
1989 valuensparseg[4] = track->Eta();\r
1990\r
1991 AliDebug(2,Form("charge %d",(Int_t)track->Charge()));\r
1992\r
1993 ////////////////////////\r
1994 // Fill before PID\r
1995 ///////////////////////\r
1996 \r
1997 if(fMonitorWithoutPID) { \r
1998 \r
1999 valuensparseg[0] = deltaphi;\r
2000 if(fillEventPlane) fDeltaPhiMapsBeforePID->Fill(&valuensparseg[0]);\r
2001 \r
2002 //\r
2003 valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));\r
2004 if(fillEventPlane) {\r
2005 fCosPhiMapsBeforePID->Fill(&valuensparseh[0]);\r
2006 }\r
2007 }\r
2008 \r
2009 ////////////////////////\r
2010 // Apply PID\r
2011 ////////////////////////\r
2012 if(!fNoPID) {\r
2013 // Apply PID for Data\r
2014 if(!fMCPID) {\r
2015 // pid object\r
2016 AliHFEpidObject hfetrack;\r
2017 if(!fAODAnalysis) hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);\r
2018 else hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);\r
2019 hfetrack.SetRecTrack(track);\r
2020 hfetrack.SetCentrality((Int_t)binct);\r
2021 AliDebug(2,Form("centrality %f and %d",binct,hfetrack.GetCentrality()));\r
2022 hfetrack.SetPbPb();\r
2023\r
2024 // Only TOF PID\r
2025 if(fMonitorContamination) {\r
2026 if(fPIDTOFOnly->IsSelected(&hfetrack,0x0,"recTrackCont",0x0)) {\r
2027 Float_t nsigma = pidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);\r
2028 valuedeltaphicontamination[3] = nsigma;\r
2029 fDeltaPhiMapsContamination->Fill(&valuedeltaphicontamination[0]);\r
2030 }\r
2031 }\r
2032\r
2033 // Complete PID TOF+TPC\r
2034 if(!fPID->IsSelected(&hfetrack,0x0,"recTrackCont",fPIDqa)) {\r
2035 continue;\r
2036 }\r
2037 }\r
2038 else {\r
2039 if(!mcEvent) continue;\r
2040 if(!(mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;\r
2041 AliDebug(2,Form("PdgCode %d",TMath::Abs(mctrack->Particle()->GetPdgCode())));\r
2042 if(TMath::Abs(mctrack->Particle()->GetPdgCode())!=11) continue;\r
2043 }\r
2044 }\r
2045\r
2046\r
2047 /////////////////////////////////////////////////////////////////////////////\r
2048 // Add candidate to AliFlowEvent for POI and subtract from RP if needed\r
2049 ////////////////////////////////////////////////////////////////////////////\r
2050 if(fMonitorQCumulant) {\r
2051 Int_t idtrack = static_cast<AliVTrack*>(track)->GetID();\r
2052 Bool_t found = kFALSE;\r
2053 Int_t numberoffound = 0;\r
2054 AliDebug(2,Form("A: Number of tracks %d",fflowEvent->NumberOfTracks()));\r
2055 for(Int_t iRPs=0; iRPs< fflowEvent->NumberOfTracks(); iRPs++) {\r
2056 AliFlowTrack *iRP = (AliFlowTrack*) (fflowEvent->GetTrack(iRPs));\r
2057 //if(!iRP->InRPSelection()) continue;\r
2058 if( TMath::Abs(idtrack) == TMath::Abs(iRP->GetID()) ) {\r
2059 iRP->SetForPOISelection(kTRUE);\r
2060 found = kTRUE;\r
2061 numberoffound ++;\r
2062 }\r
2063 }\r
2064 AliDebug(2,Form("Found %d mal",numberoffound));\r
2065 if(!found) {\r
2066 AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*) MakeTrack(massElectron,track->Pt(),track->Phi(), track->Eta());\r
2067 sTrack->SetID(idtrack);\r
2068 fflowEvent->AddTrack(sTrack);\r
2069 AliDebug(2,"Add the track");\r
2070 }\r
2071 AliDebug(2,Form("B: Number of tracks %d",fflowEvent->NumberOfTracks()));\r
2072 }\r
2073 \r
2074 \r
2075 /////////////////////\r
2076 // Fill THnSparseF\r
2077 /////////////////////\r
2078\r
2079 //\r
2080 valuensparseabis[0] = eventplanesubtracted;\r
2081 if((fillEventPlane) && (fMonitorEventPlane)) fEventPlaneaftersubtraction->Fill(&valuensparseabis[0]);\r
2082 \r
2083\r
2084 if(fMonitorEventPlane) \r
2085 {\r
2086 //\r
2087 valuensparsee[0] = TMath::Cos(2*phitrack);\r
2088 fCos2phie->Fill(&valuensparsee[0]);\r
2089 valuensparsee[0] = TMath::Sin(2*phitrack);\r
2090 fSin2phie->Fill(&valuensparsee[0]);\r
2091 //\r
2092 valuensparsee[0] = TMath::Cos(2*eventplanesubtracted);\r
2093 if(fillEventPlane) fCos2phiep->Fill(&valuensparsee[0]);\r
2094 valuensparsee[0] = TMath::Sin(2*eventplanesubtracted);\r
2095 if(fillEventPlane) fSin2phiep->Fill(&valuensparsee[0]);\r
2096 valuensparsee[0] = TMath::Sin(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));\r
2097 if(fillEventPlane) fSin2phiephiep->Fill(&valuensparsee[0]);\r
2098 //\r
2099 }\r
2100\r
2101 // \r
2102 valuensparseg[0] = deltaphi;\r
2103 if(fillEventPlane) fDeltaPhiMaps->Fill(&valuensparseg[0]);\r
2104 \r
2105 //\r
2106 valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));\r
2107 if(fillEventPlane) {\r
2108 fCosPhiMaps->Fill(&valuensparseh[0]); //TR: fCosPhiQSum+=valuensparseh[0]*TMath:Sqrt(qAna->X()*qAna->X()+qAna->Y()*qAna->Y()); fCosPhiQN++;\r
2109 if((valuefractioncont[1] >=0) && (valuefractioncont[1] < 11)){\r
2110 if(fContamination[((Int_t)valuefractioncont[1])]){\r
2111 Double_t weight = fContamination[((Int_t)valuefractioncont[1])]->Eval(track->P());\r
2112 if(weight<0.0) weight=0.0;\r
2113 if(weight>1.0) weight=1.0;\r
2114 fFractionContamination->Fill(&valuefractioncont[0],weight);\r
2115 if(fv2contamination[((Int_t)valuefractioncont[1])]){\r
2116 Double_t v2 = fv2contamination[((Int_t)valuefractioncont[1])]->Eval(track->Pt());\r
2117 AliDebug(2,Form("value v2 %f, contamination %f and pt %f centrality %d\n",v2,weight,track->Pt(),(Int_t)valuefractioncont[1]));\r
2118 AliDebug(2,Form("Check for centrality 3: value v2 %f, contamination %f\n",fv2contamination[3]->Eval(track->Pt()),fContamination[3]->Eval(track->P())));\r
2119 AliDebug(2,Form("Check for centrality 4: value v2 %f, contamination %f\n",fv2contamination[4]->Eval(track->Pt()),fContamination[4]->Eval(track->P())));\r
2120 AliDebug(2,Form("Check for centrality 5: value v2 %f, contamination %f\n",fv2contamination[5]->Eval(track->Pt()),fContamination[5]->Eval(track->P())));\r
2121 fContaminationv2->Fill(valuefractioncont[1],valuefractioncont[0],v2,weight);\r
2122 }\r
2123 } \r
2124 }\r
2125 if(fMonitorEventPlane) {\r
2126 if(fSP)\r
2127 valuensparseh[0] *= TMath::Sqrt(qAna->X()*qAna->X()+qAna->Y()*qAna->Y());\r
2128 fProfileCosPhiMaps->Fill(valuensparsehprofile[1],valuensparsehprofile[2],valuensparseh[0]); //TR: info\r
2129 }\r
2130 }\r
2131 \r
2132 if(fMonitorPhotonic) {\r
2133 Int_t indexmother = -1;\r
2134 Int_t source = 1;\r
2135 if(mcthere) source = fBackgroundSubtraction->FindMother(mctrack->GetLabel(),indexmother);\r
2136 fBackgroundSubtraction->LookAtNonHFE(k, track, fInputEvent, 1, binct, deltaphi, source, indexmother);\r
2137 \r
2138 if((!fAODAnalysis && mcthere) || !mcthere) {\r
2139 // background\r
2140 source = 0;\r
2141 indexmother = -1;\r
2142 source = FindMother(TMath::Abs(track->GetLabel()),mcEvent, indexmother);\r
2143 valuensparseMCSourceDeltaPhiMaps[2] = source;\r
2144 if(mcEvent) fMCSourceDeltaPhiMaps->Fill(&valuensparseMCSourceDeltaPhiMaps[0]);\r
2145 //LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);\r
2146 Int_t taggedvalue = LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);\r
2147 if(fMonitorPhotonic) {\r
2148 // No opposite charge partner found in the invariant mass choosen\r
2149 if((taggedvalue!=2) && (taggedvalue!=6)) {\r
2150 //fDeltaPhiMapsTaggedNonPhotonic->Fill(&valuensparseg[0]);\r
2151 //fCosPhiMapsTaggedNonPhotonic->Fill(&valuensparseh[0]);\r
2152 }\r
2153 // One opposite charge partner found in the invariant mass choosen\r
2154 if((taggedvalue==2) || (taggedvalue==6)) {\r
2155 fDeltaPhiMapsTaggedPhotonic->Fill(&valuensparseg[0]);\r
2156 //fCosPhiMapsTaggedPhotonic->Fill(&valuensparseh[0]);\r
2157 }\r
2158 // One same charge partner found in the invariant mass choosen\r
2159 if((taggedvalue==4) || (taggedvalue==6)) {\r
2160 fDeltaPhiMapsTaggedPhotonicLS->Fill(&valuensparseg[0]);\r
2161 //fCosPhiMapsTaggedPhotonicLS->Fill(&valuensparseh[0]);\r
2162 }\r
2163 }\r
2164 }\r
2165 }\r
2166 \r
2167 }\r
2168\r
2169 //////////////////////////////////////////////////////////////////////////////\r
2170 ///////////////////////////AFTERBURNER\r
2171 if (fAfterBurnerOn & fMonitorQCumulant)\r
2172 {\r
2173 fflowEvent->AddFlow(fV1,fV2,fV3,fV4,fV5); //add flow\r
2174 fflowEvent->CloneTracks(fNonFlowNumberOfTrackClones); //add nonflow by cloning tracks\r
2175 }\r
2176 //////////////////////////////////////////////////////////////////////////////\r
2177\r
2178\r
2179\r
2180 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {\r
2181 // if((fBinCentralityLess[bincless]< cntr) && (cntr < fBinCentralityLess[bincless+1])) PostData(bincless+2,fflowEvent);\r
2182 //}\r
2183\r
2184 if(fMonitorPhotonic) {\r
2185 if(fArraytrack) {\r
2186 delete fArraytrack;\r
2187 fArraytrack = NULL;\r
2188 }\r
2189 }\r
2190 \r
2191 PostData(1, fListHist);\r
2192\r
2193 \r
2194}\r
2195//______________________________________________________________________________\r
2196AliFlowCandidateTrack *AliAnalysisTaskHFEFlow::MakeTrack( Double_t mass, \r
2197 Double_t pt, Double_t phi, Double_t eta) {\r
2198 //\r
2199 // Make Track (Not needed actually)\r
2200 //\r
2201\r
2202 AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();\r
2203 sTrack->SetMass(mass);\r
2204 sTrack->SetPt(pt);\r
2205 sTrack->SetPhi(phi);\r
2206 sTrack->SetEta(eta);\r
2207 sTrack->SetForPOISelection(kTRUE);\r
2208 sTrack->SetForRPSelection(kFALSE);\r
2209 return sTrack;\r
2210}\r
2211//_________________________________________________________________________________ \r
2212Double_t AliAnalysisTaskHFEFlow::GetPhiAfterAddV2(Double_t phi,Double_t reactionPlaneAngle) const\r
2213{\r
2214 //\r
2215 // Adds v2, uses Newton-Raphson iteration\r
2216 //\r
2217 Double_t phiend=phi;\r
2218 Double_t phi0=phi;\r
2219 Double_t f=0.;\r
2220 Double_t fp=0.;\r
2221 Double_t phiprev=0.;\r
2222\r
2223 for (Int_t i=0; i<fMaxNumberOfIterations; i++)\r
2224 {\r
2225 phiprev=phiend; //store last value for comparison\r
2226 f = phiend-phi0+fV2*TMath::Sin(2.*(phiend-reactionPlaneAngle));\r
2227 fp = 1.0+2.0*fV2*TMath::Cos(2.*(phiend-reactionPlaneAngle)); //first derivative\r
2228 phiend -= f/fp;\r
2229 if (TMath::AreEqualAbs(phiprev,phiend,fPrecisionPhi)) break;\r
2230 }\r
2231 return phiend;\r
2232}\r
2233//_____________________________________________________________________________________________\r
2234Int_t AliAnalysisTaskHFEFlow::LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *vEvent, AliMCEvent *mcEvent,Int_t binct,Double_t deltaphi,Int_t source,Int_t indexmother)\r
2235{ \r
2236 //\r
2237 // Look At Non HFE\r
2238 //\r
2239\r
2240 // return -1 if nothing\r
2241 // return 2 if opposite charge within the mass range found\r
2242 // return 4 if like charge within the mass range found\r
2243 // return 6 if opposite charge and like charge within the mass range found\r
2244 //\r
2245\r
2246 Int_t taggedphotonic = -1;\r
2247\r
2248 Bool_t oppositetaggedphotonic = kFALSE;\r
2249 Bool_t sametaggedphotonic = kFALSE;\r
2250\r
2251 AliDebug(2,Form("fCounterPoolBackground %d in LookAtNonHFE!!!",fCounterPoolBackground));\r
2252 if(!fArraytrack) return taggedphotonic;\r
2253 AliDebug(2,Form("process track %d",iTrack1));\r
2254 \r
2255 TVector3 v3D1;\r
2256 TVector3 v3D2;\r
2257\r
2258 Double_t valuensparseDeltaPhiMaps[5];\r
2259 Double_t valueangle[3];\r
2260\r
2261 valuensparseDeltaPhiMaps[1] = binct;\r
2262 valuensparseDeltaPhiMaps[2] = track1->Pt();\r
2263 valuensparseDeltaPhiMaps[0] = deltaphi;\r
2264 valuensparseDeltaPhiMaps[4] = source;\r
2265 \r
2266 valueangle[2] = source;\r
2267 valueangle[1] = binct;\r
2268\r
2269 // Pdg code\r
2270 Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()),mcEvent);\r
2271 Int_t numberfound = 0;\r
2272\r
2273 //Magnetic Field\r
2274 Double_t bfield = vEvent->GetMagneticField();\r
2275\r
2276 // Get Primary vertex\r
2277 const AliVVertex *pVtx = vEvent->GetPrimaryVertex();\r
2278 \r
2279 for(Int_t idex = 0; idex < fCounterPoolBackground; idex++) \r
2280 {\r
2281\r
2282 Int_t iTrack2 = fArraytrack->At(idex);\r
2283 AliDebug(2,Form("track %d",iTrack2));\r
2284 AliVTrack* track2 = (AliVTrack *) vEvent->GetTrack(iTrack2);\r
2285 if (!track2) \r
2286 {\r
2287 printf("ERROR: Could not receive track %d", iTrack2);\r
2288 continue;\r
2289 }\r
2290 if(iTrack2==iTrack1) continue;\r
2291 AliDebug(2,"Different");\r
2292\r
2293 // Reset the MC info\r
2294 valueangle[2] = source;\r
2295 valuensparseDeltaPhiMaps[4] = source;\r
2296\r
2297 // track cuts and PID already done\r
2298\r
2299 // if MC look\r
2300 Int_t pdg2 = -100;\r
2301 if(mcEvent) {\r
2302 Int_t source2 = 0;\r
2303 Int_t indexmother2 = -1;\r
2304 source2 = FindMother(TMath::Abs(track2->GetLabel()),mcEvent, indexmother2);\r
2305 pdg2 = CheckPdg(TMath::Abs(track2->GetLabel()),mcEvent);\r
2306 if(source2 >=0 ) {\r
2307 if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0)) {\r
2308 if(source == kElectronfromconversion) {\r
2309 valueangle[2] = kElectronfromconversionboth;\r
2310 valuensparseDeltaPhiMaps[4] = kElectronfromconversionboth;\r
2311 numberfound++;\r
2312 }\r
2313 if(source == kElectronfrompi0) {\r
2314 valueangle[2] = kElectronfrompi0both;\r
2315 valuensparseDeltaPhiMaps[4] = kElectronfrompi0both;\r
2316 }\r
2317 if(source == kElectronfrometa) {\r
2318 valueangle[2] = kElectronfrometaboth;\r
2319 valuensparseDeltaPhiMaps[4] = kElectronfrometaboth;\r
2320 }\r
2321 }\r
2322 }\r
2323 }\r
2324 \r
2325 if(fAlgorithmMA && (!fAODAnalysis))\r
2326 {\r
2327 // tracks\r
2328 AliESDtrack *esdtrack2 = dynamic_cast<AliESDtrack *>(track2); \r
2329 AliESDtrack *esdtrack1 = dynamic_cast<AliESDtrack *>(track1); \r
2330 if((!esdtrack2) || (!esdtrack1)) continue;\r
2331\r
2332 //Variables\r
2333 Double_t p1[3];\r
2334 Double_t p2[3];\r
2335 Double_t xt1; //radial position track 1 at the DCA point\r
2336 Double_t xt2; //radial position track 2 at the DCA point\r
2337 //DCA track1-track2\r
2338 Double_t dca12 = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1);\r
2339\r
2340 // Cut dca\r
2341 if(dca12 > fMaxdca) continue;\r
2342 \r
2343 //Momento of the track extrapolated to DCA track-track \r
2344 //Track1\r
2345 Bool_t hasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1);\r
2346 //Track2\r
2347 Bool_t hasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2);\r
2348 \r
2349 if(!hasdcaT1 || !hasdcaT2) AliWarning("It could be a problem in the extrapolation");\r
2350 \r
2351 //track1-track2 Invariant Mass\r
2352 Double_t eMass = 0.000510998910; //Electron mass in GeV\r
2353 Double_t pP1 = sqrt(p1[0]*p1[0]+p1[1]*p1[1]+p1[2]*p1[2]); //Track 1 momentum\r
2354 Double_t pP2 = sqrt(p2[0]*p2[0]+p2[1]*p2[1]+p2[2]*p2[2]); //Track 2 momentum\r
2355 Double_t eE1 = TMath::Sqrt(pP1*pP1+eMass*eMass);\r
2356 Double_t eE2 = TMath::Sqrt(pP2*pP2+eMass*eMass);\r
2357 \r
2358 //TLorentzVector v1(p1[0],p1[1],p1[2],sqrt(eMass*eMass+pP1*pP1));\r
2359 //TLorentzVector v2(p2[0],p2[1],p2[2],sqrt(eMass*eMass+pP2*pP2));\r
2360 //Double_t imass = (v1+v2).M(); //Invariant Mass\r
2361 //Double_t angle3D = v1.Angle(v2.Vect()); //Opening Angle (Total Angle)\r
2362 \r
2363 // daughter\r
2364 v3D1.SetXYZ(p1[0],p1[1],p1[2]);\r
2365 v3D2.SetXYZ(p2[0],p2[1],p2[2]);\r
2366 Double_t openingangle = TVector2::Phi_0_2pi(v3D2.Angle(v3D1));\r
2367 \r
2368 // mother\r
2369 TVector3 motherrec = v3D1 + v3D2;\r
2370 Double_t invmass = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));\r
2371 \r
2372 // xy\r
2373 //TVector3 vectordiff = v3D1 - v3D2;\r
2374 //Double_t diffphi = TVector2::Phi_0_2pi(vectordiff.Phi());\r
2375 //Double_t massxy = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(diffphi)));\r
2376\r
2377 // rz\r
2378 //Double_t difftheta = TVector2::Phi_0_2pi(vectordiff.Eta());\r
2379 //Double_t massrz = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(difftheta)));\r
2380 \r
2381\r
2382 Float_t fCharge1 = track1->Charge();\r
2383 Float_t fCharge2 = track2->Charge();\r
2384\r
2385 // Fill Histo\r
2386 //valueangle[0] = diffphi;\r
2387 //valueangle[1] = difftheta;\r
2388 valueangle[0] = openingangle;\r
2389 if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);\r
2390 else fOppSignAngle->Fill(&valueangle[0]);\r
2391\r
2392 // Cut\r
2393 if(openingangle > fMaxopening3D) continue;\r
2394 //if(difftheta > fMaxopeningtheta) continue;\r
2395 //if(diffphi > fMaxopeningphi) continue;\r
2396\r
2397 // Invmass\r
2398 valuensparseDeltaPhiMaps[3] = invmass;\r
2399 if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
2400 else fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
2401 \r
2402 // Cut\r
2403 if(invmass < fMaxInvmass) {\r
2404 if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;\r
2405 if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;\r
2406 }\r
2407\r
2408\r
2409 }\r
2410 else \r
2411 {\r
2412 Int_t fPDGtrack1 = 11; \r
2413 Int_t fPDGtrack2 = 11;\r
2414 \r
2415 Float_t fCharge1 = track1->Charge();\r
2416 Float_t fCharge2 = track2->Charge();\r
2417 \r
2418 if(fCharge1>0) fPDGtrack1 = -11;\r
2419 if(fCharge2>0) fPDGtrack2 = -11;\r
2420 \r
2421 AliKFParticle ktrack1(*track1, fPDGtrack1);\r
2422 AliKFParticle ktrack2(*track2, fPDGtrack2);\r
2423 AliKFParticle recoGamma(ktrack1, ktrack2);\r
2424 \r
2425 //Reconstruction Cuts\r
2426 if(recoGamma.GetNDF()<1) continue;\r
2427 Double_t chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();\r
2428 if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;\r
2429\r
2430 // DCA\r
2431 //Double_t dca12 = ktrack1.GetDistanceFromParticle(ktrack2);\r
2432 //if(dca12 > fMaxdca) continue; \r
2433\r
2434 // if set mass constraint\r
2435 if(fSetMassConstraint && pVtx) {\r
2436 AliKFVertex primV(*pVtx);\r
2437 primV += recoGamma;\r
2438 primV -= ktrack1;\r
2439 primV -= ktrack2;\r
2440 recoGamma.SetProductionVertex(primV);\r
2441 recoGamma.SetMassConstraint(0,0.0001);\r
2442 } \r
2443\r
2444 //Invariant Mass\r
2445 Double_t imass; \r
2446 Double_t width;\r
2447 recoGamma.GetMass(imass,width);\r
2448 \r
2449 //Opening Angle (Total Angle)\r
2450 Double_t angle = ktrack1.GetAngle(ktrack2);\r
2451 valueangle[0] = angle;\r
2452 if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);\r
2453 else fOppSignAngle->Fill(&valueangle[0]);\r
2454\r
2455 // Cut\r
2456 if(angle > fMaxopening3D) continue; \r
2457\r
2458 // Invmass\r
2459 valuensparseDeltaPhiMaps[3] = imass;\r
2460 if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
2461 else {\r
2462 fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);\r
2463 /*\r
2464 if(valueangle[2] == kElectronfromconversionboth) {\r
2465 printf("Reconstructed charge1 %f, charge 2 %f and invmass %f",fCharge1,fCharge2,imass);\r
2466 printf("MC charge1 %d, charge 2 %d",pdg1,pdg2);\r
2467 printf("DCA %f",dca12);\r
2468 printf("Number of found %d",numberfound);\r
2469 }\r
2470 */\r
2471 }\r
2472 \r
2473 // Cut\r
2474 if(imass < fMaxInvmass) {\r
2475 if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;\r
2476 if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;\r
2477 }\r
2478 }\r
2479 }\r
2480 \r
2481 if(oppositetaggedphotonic && sametaggedphotonic){\r
2482 taggedphotonic = 6;\r
2483 }\r
2484\r
2485 if(!oppositetaggedphotonic && sametaggedphotonic){\r
2486 taggedphotonic = 4;\r
2487 }\r
2488\r
2489 if(oppositetaggedphotonic && !sametaggedphotonic){\r
2490 taggedphotonic = 2;\r
2491 }\r
2492\r
2493 \r
2494 return taggedphotonic;\r
2495}\r
2496//_________________________________________________________________________\r
2497Int_t AliAnalysisTaskHFEFlow::FindMother(Int_t tr, AliMCEvent *mcEvent, Int_t &indexmother){\r
2498 //\r
2499 // Find the mother if MC\r
2500 //\r
2501\r
2502 if(!mcEvent) return 0;\r
2503\r
2504 Int_t pdg = CheckPdg(tr,mcEvent);\r
2505 if(TMath::Abs(pdg)!= 11) {\r
2506 indexmother = -1;\r
2507 return kNoElectron;\r
2508 }\r
2509 \r
2510 indexmother = IsMotherGamma(tr,mcEvent);\r
2511 if(indexmother > 0) return kElectronfromconversion;\r
2512 indexmother = IsMotherPi0(tr,mcEvent);\r
2513 if(indexmother > 0) return kElectronfrompi0;\r
2514 indexmother = IsMotherC(tr,mcEvent);\r
2515 if(indexmother > 0) return kElectronfromC;\r
2516 indexmother = IsMotherB(tr,mcEvent);\r
2517 if(indexmother > 0) return kElectronfromB;\r
2518 indexmother = IsMotherEta(tr,mcEvent);\r
2519 if(indexmother > 0) return kElectronfrometa;\r
2520 \r
2521 return kElectronfromother;\r
2522\r
2523\r
2524}\r
2525//____________________________________________________________________________________________________________\r
2526Int_t AliAnalysisTaskHFEFlow::CheckPdg(Int_t tr, AliMCEvent* mcEvent) {\r
2527\r
2528 //\r
2529 // Return the pdg of the particle\r
2530 //\r
2531\r
2532\r
2533 Int_t pdgcode = -1;\r
2534 if(tr < 0) return pdgcode;\r
2535\r
2536 if(!mcEvent) return pdgcode;\r
2537\r
2538 AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
2539 \r
2540 \r
2541 if(mctrack->IsA() == AliMCParticle::Class()) {\r
2542 AliMCParticle *mctrackesd = NULL;\r
2543 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;\r
2544 pdgcode = mctrackesd->PdgCode();\r
2545 }\r
2546\r
2547 if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
2548 AliAODMCParticle *mctrackaod = NULL;\r
2549 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;\r
2550 pdgcode = mctrackaod->GetPdgCode();\r
2551 }\r
2552 \r
2553 return pdgcode;\r
2554\r
2555 \r
2556}\r
2557//____________________________________________________________________________________________________________\r
2558Int_t AliAnalysisTaskHFEFlow::IsMotherGamma(Int_t tr, AliMCEvent* mcEvent) {\r
2559\r
2560 //\r
2561 // Return the lab of gamma mother or -1 if not gamma\r
2562 //\r
2563\r
2564 if(tr < 0) return -1;\r
2565 AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
2566 \r
2567 if(mctrack->IsA() == AliMCParticle::Class()) {\r
2568 AliMCParticle *mctrackesd = NULL;\r
2569 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2570 TParticle *particle = 0x0;\r
2571 particle = mctrackesd->Particle();\r
2572 // Take mother\r
2573 if(!particle) return -1;\r
2574 Int_t imother = particle->GetFirstMother(); \r
2575 if(imother < 0) return -1; \r
2576 AliMCParticle *mothertrack = NULL;\r
2577 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2578 TParticle * mother = mothertrack->Particle();\r
2579 if(!mother) return -1;\r
2580 // Check gamma \r
2581 Int_t pdg = mother->GetPdgCode();\r
2582 if(TMath::Abs(pdg) == 22) return imother;\r
2583 if(TMath::Abs(pdg) == 11) {\r
2584 return IsMotherGamma(imother,mcEvent);\r
2585 }\r
2586 return -1;\r
2587 }\r
2588\r
2589 if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
2590 AliAODMCParticle *mctrackaod = NULL;\r
2591 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2592 // Take mother\r
2593 Int_t imother = mctrackaod->GetMother();\r
2594 if(imother < 0) return -1; \r
2595 AliAODMCParticle *mothertrack = NULL;\r
2596 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2597 // Check gamma \r
2598 Int_t pdg = mothertrack->GetPdgCode();\r
2599 if(TMath::Abs(pdg) == 22) return imother;\r
2600 if(TMath::Abs(pdg) == 11) {\r
2601 return IsMotherGamma(imother,mcEvent);\r
2602 }\r
2603 return -1;\r
2604\r
2605 }\r
2606 \r
2607 return -1;\r
2608\r
2609 \r
2610}\r
2611//\r
2612//____________________________________________________________________________________________________________\r
2613Int_t AliAnalysisTaskHFEFlow::IsMotherPi0(Int_t tr, AliMCEvent* mcEvent) {\r
2614\r
2615 //\r
2616 // Return the lab of pi0 mother or -1 if not pi0\r
2617 //\r
2618\r
2619 if(tr < 0) return -1;\r
2620 AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
2621 \r
2622 if(mctrack->IsA() == AliMCParticle::Class()) {\r
2623 AliMCParticle *mctrackesd = NULL;\r
2624 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2625 TParticle *particle = 0x0;\r
2626 particle = mctrackesd->Particle();\r
2627 // Take mother\r
2628 if(!particle) return -1;\r
2629 Int_t imother = particle->GetFirstMother(); \r
2630 if(imother < 0) return -1; \r
2631 AliMCParticle *mothertrack = NULL;\r
2632 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2633 TParticle * mother = mothertrack->Particle();\r
2634 if(!mother) return -1;\r
2635 // Check gamma \r
2636 Int_t pdg = mother->GetPdgCode();\r
2637 if(TMath::Abs(pdg) == 111) return imother;\r
2638 if(TMath::Abs(pdg) == 11) {\r
2639 return IsMotherPi0(imother,mcEvent);\r
2640 }\r
2641 return -1;\r
2642 }\r
2643\r
2644 if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
2645 AliAODMCParticle *mctrackaod = NULL;\r
2646 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2647 // Take mother\r
2648 Int_t imother = mctrackaod->GetMother();\r
2649 if(imother < 0) return -1; \r
2650 AliAODMCParticle *mothertrack = NULL;\r
2651 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2652 // Check gamma \r
2653 Int_t pdg = mothertrack->GetPdgCode();\r
2654 if(TMath::Abs(pdg) == 111) return imother;\r
2655 if(TMath::Abs(pdg) == 11) {\r
2656 return IsMotherPi0(imother,mcEvent);\r
2657 }\r
2658 return -1;\r
2659 }\r
2660\r
2661 return -1;\r
2662 \r
2663}\r
2664//____________________________________________________________________________________________________________\r
2665Int_t AliAnalysisTaskHFEFlow::IsMotherC(Int_t tr, AliMCEvent* mcEvent) {\r
2666\r
2667 //\r
2668 // Return the lab of signal mother or -1 if not signal\r
2669 //\r
2670\r
2671 if(tr < 0) return -1;\r
2672 AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
2673 \r
2674 if(mctrack->IsA() == AliMCParticle::Class()) {\r
2675 AliMCParticle *mctrackesd = NULL;\r
2676 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2677 TParticle *particle = 0x0;\r
2678 particle = mctrackesd->Particle();\r
2679 // Take mother\r
2680 if(!particle) return -1;\r
2681 Int_t imother = particle->GetFirstMother(); \r
2682 if(imother < 0) return -1; \r
2683 AliMCParticle *mothertrack = NULL;\r
2684 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2685 TParticle * mother = mothertrack->Particle();\r
2686 if(!mother) return -1;\r
2687 // Check gamma \r
2688 Int_t pdg = mother->GetPdgCode();\r
2689 if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;\r
2690 if(TMath::Abs(pdg) == 11) {\r
2691 return IsMotherC(imother,mcEvent);\r
2692 }\r
2693 return -1;\r
2694 }\r
2695\r
2696 if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
2697 AliAODMCParticle *mctrackaod = NULL;\r
2698 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2699 // Take mother\r
2700 Int_t imother = mctrackaod->GetMother();\r
2701 if(imother < 0) return -1; \r
2702 AliAODMCParticle *mothertrack = NULL;\r
2703 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2704 // Check gamma \r
2705 Int_t pdg = mothertrack->GetPdgCode();\r
2706 if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;\r
2707 if(TMath::Abs(pdg) == 11) {\r
2708 return IsMotherC(imother,mcEvent);\r
2709 }\r
2710 return -1;\r
2711 }\r
2712\r
2713 return -1;\r
2714\r
2715}\r
2716//____________________________________________________________________________________________________________\r
2717Int_t AliAnalysisTaskHFEFlow::IsMotherB(Int_t tr, AliMCEvent* mcEvent) {\r
2718\r
2719 //\r
2720 // Return the lab of signal mother or -1 if not signal\r
2721 //\r
2722\r
2723 if(tr < 0) return -1;\r
2724 AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
2725 \r
2726 if(mctrack->IsA() == AliMCParticle::Class()) {\r
2727 AliMCParticle *mctrackesd = NULL;\r
2728 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2729 TParticle *particle = 0x0;\r
2730 particle = mctrackesd->Particle();\r
2731 // Take mother\r
2732 if(!particle) return -1;\r
2733 Int_t imother = particle->GetFirstMother(); \r
2734 if(imother < 0) return -1; \r
2735 AliMCParticle *mothertrack = NULL;\r
2736 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2737 TParticle * mother = mothertrack->Particle();\r
2738 if(!mother) return -1;\r
2739 // Check gamma \r
2740 Int_t pdg = mother->GetPdgCode();\r
2741 if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother; \r
2742 if(TMath::Abs(pdg) == 11) {\r
2743 return IsMotherB(imother,mcEvent);\r
2744 }\r
2745 return -1;\r
2746 }\r
2747\r
2748 if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
2749 AliAODMCParticle *mctrackaod = NULL;\r
2750 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2751 // Take mother\r
2752 Int_t imother = mctrackaod->GetMother();\r
2753 if(imother < 0) return -1; \r
2754 AliAODMCParticle *mothertrack = NULL;\r
2755 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2756 // Check gamma \r
2757 Int_t pdg = mothertrack->GetPdgCode();\r
2758 if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother;\r
2759 if(TMath::Abs(pdg) == 11) {\r
2760 return IsMotherB(imother,mcEvent);\r
2761 }\r
2762 return -1;\r
2763 }\r
2764\r
2765 return -1;\r
2766\r
2767}\r
2768//____________________________________________________________________________________________________________\r
2769Int_t AliAnalysisTaskHFEFlow::IsMotherEta(Int_t tr, AliMCEvent* mcEvent) {\r
2770\r
2771 //\r
2772 // Return the lab of pi0 mother or -1 if not pi0\r
2773 //\r
2774\r
2775 if(tr < 0) return -1;\r
2776 AliVParticle *mctrack = mcEvent->GetTrack(tr);\r
2777 \r
2778 if(mctrack->IsA() == AliMCParticle::Class()) {\r
2779 AliMCParticle *mctrackesd = NULL;\r
2780 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2781 TParticle *particle = 0x0;\r
2782 particle = mctrackesd->Particle();\r
2783 // Take mother\r
2784 if(!particle) return -1;\r
2785 Int_t imother = particle->GetFirstMother(); \r
2786 if(imother < 0) return -1; \r
2787 AliMCParticle *mothertrack = NULL;\r
2788 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2789 TParticle * mother = mothertrack->Particle();\r
2790 if(!mother) return -1;\r
2791 // Check gamma \r
2792 Int_t pdg = mother->GetPdgCode();\r
2793 if(TMath::Abs(pdg) == 221) return imother;\r
2794 if(TMath::Abs(pdg) == 11) {\r
2795 return IsMotherEta(imother,mcEvent);\r
2796 }\r
2797 return -1;\r
2798 }\r
2799\r
2800 if(mctrack->IsA() == AliAODMCParticle::Class()) {\r
2801 AliAODMCParticle *mctrackaod = NULL;\r
2802 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;\r
2803 // Take mother\r
2804 Int_t imother = mctrackaod->GetMother();\r
2805 if(imother < 0) return -1; \r
2806 AliAODMCParticle *mothertrack = NULL;\r
2807 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;\r
2808 // Check gamma \r
2809 Int_t pdg = mothertrack->GetPdgCode();\r
2810 if(TMath::Abs(pdg) == 221) return imother;\r
2811 if(TMath::Abs(pdg) == 11) {\r
2812 return IsMotherEta(imother,mcEvent);\r
2813 }\r
2814 return -1;\r
2815 }\r
2816\r
2817 return -1;\r
2818 \r
2819}\r