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