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