]> 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
Further TPC pid 2011
[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 ////////////////////
38be5083 1393
1394 UInt_t isEventSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
8b835997 1395 if(fTriggerUsed==0){
1396
38be5083 1397 // MB, semi-central and central
8b835997 1398
38be5083 1399 if ( !((isEventSelected & AliVEvent::kCentral) |
1400 (isEventSelected & AliVEvent::kSemiCentral) |
1401 (isEventSelected & AliVEvent::kMB)) ) return;
8b835997 1402
1403 }
38be5083 1404 else if(fTriggerUsed==1){
8b835997 1405
1406 // semi-central Ionut
1407
38be5083 1408 if ( !((isEventSelected & AliVEvent::kCentral) |
1409 (isEventSelected & AliVEvent::kSemiCentral) |
1410 (isEventSelected & AliVEvent::kMB)) ) return;
8b835997 1411
38be5083 1412 Bool_t isMB = (InputEvent()->GetTriggerMask() & (ULong64_t(1)<<1));
8b835997 1413 //Bool_t isCentral = (InputEvent()->GetTriggerMask() & (ULong64_t(1)<<4));
1414 Bool_t isSemiCentral = (InputEvent()->GetTriggerMask() & (ULong64_t(1)<<7));
1415
38be5083 1416 if(!(isSemiCentral | isMB)) return;
8b835997 1417
1418 }
38be5083 1419 else if(fTriggerUsed==2){
8b835997 1420
1421 // semi-central Andrea and Muons
1422
38be5083 1423 if ( !(isEventSelected & AliVEvent::kAny) ) return;
8b835997 1424
1425 //TString firedTriggerClasses = static_cast<const AliAODEvent*>(InputEvent())->GetFiredTriggerClasses();
1426 TString firedTriggerClasses = InputEvent()->GetFiredTriggerClasses();
1427
1428 if ( ! ( firedTriggerClasses.Contains("CVLN_B2-B-NOPF-ALLNOTRD") || firedTriggerClasses.Contains("CVLN_R1-B-NOPF-ALLNOTRD") || firedTriggerClasses.Contains("CSEMI_R1-B-NOPF-ALLNOTRD") ) ) return;
1429 }
c683985a 1430
1431 /////////////////
1432 // centrality
1433 /////////////////
1434
1435 //AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
1436 //if(!esd) return;
1437 AliCentrality *centrality = fInputEvent->GetCentrality();
1438 //AliDebug(2,"Got the centrality");
1439 if(!centrality) return;
1440 cntr = centrality->GetCentralityPercentile("V0M");
1441 if((0.0< cntr) && (cntr<5.0)) binct = 0.5;
1442 if((5.0< cntr) && (cntr<10.0)) binct = 1.5;
1443 if((10.0< cntr) && (cntr<20.0)) binct = 2.5;
1444 if((20.0< cntr) && (cntr<30.0)) binct = 3.5;
1445 if((30.0< cntr) && (cntr<40.0)) binct = 4.5;
1446 if((40.0< cntr) && (cntr<50.0)) binct = 5.5;
1447 if((50.0< cntr) && (cntr<60.0)) binct = 6.5;
1448 if((60.0< cntr) && (cntr<70.0)) binct = 7.5;
1449 if((70.0< cntr) && (cntr<80.0)) binct = 8.5;
1450 if((80.0< cntr) && (cntr<90.0)) binct = 9.5;
1451 if((90.0< cntr) && (cntr<100.0)) binct = 10.5;
1452
1453 if((0.< cntr) && (cntr < 20.)) binctt = 0.5;
1454 if((20.< cntr) && (cntr < 40.)) binctt = 1.5;
1455 if((40.< cntr) && (cntr < 80.)) binctt = 2.5;
1456
1457 if((0.0< cntr) && (cntr<5.0)) binctMore = 0.5;
1458 if((5.0< cntr) && (cntr<10.0)) binctMore = 1.5;
1459 if((10.0< cntr) && (cntr<15.0)) binctMore = 2.5;
1460 if((15.0< cntr) && (cntr<20.0)) binctMore = 3.5;
1461 if((20.0< cntr) && (cntr<25.0)) binctMore = 4.5;
1462 if((25.0< cntr) && (cntr<30.0)) binctMore = 5.5;
1463 if((30.0< cntr) && (cntr<35.0)) binctMore = 6.5;
1464 if((35.0< cntr) && (cntr<40.0)) binctMore = 7.5;
1465 if((40.0< cntr) && (cntr<45.0)) binctMore = 8.5;
1466 if((45.0< cntr) && (cntr<50.0)) binctMore = 9.5;
1467 if((50.0< cntr) && (cntr<55.0)) binctMore = 10.5;
1468 if((55.0< cntr) && (cntr<60.0)) binctMore = 11.5;
1469 if((60.0< cntr) && (cntr<65.0)) binctMore = 12.5;
1470 if((65.0< cntr) && (cntr<70.0)) binctMore = 13.5;
1471 if((70.0< cntr) && (cntr<75.0)) binctMore = 14.5;
1472 if((75.0< cntr) && (cntr<80.0)) binctMore = 15.5;
1473 if((80.0< cntr) && (cntr<85.0)) binctMore = 16.5;
1474 if((85.0< cntr) && (cntr<90.0)) binctMore = 17.5;
1475 if((90.0< cntr) && (cntr<95.0)) binctMore = 18.5;
1476 if((95.0< cntr) && (cntr<100.0)) binctMore = 19.5;
1477
1478 binctLess = cntr;
1479
1480
1481 if(binct > 11.0) return;
1482
1483 // centrality
1484 valuensparsea[3] = binct;
1485 valuensparseabis[1] = binct;
1486 valuensparsee[1] = binct;
1487 valuensparsef[1] = binctMore;
1488 valuensparsefsin[1] = binct;
1489 valuensparsefbis[3] = binctMore;
1490 valuensparsefbissin[3] = binct;
1491 valuensparseg[1] = binct;
1492 valuensparseh[1] = binct;
1493 valuefractioncont[1] = binct;
1494 valuensparsehprofile[1] = binct;
1495 valuecossinephiep[2] = binctMore;
1496 valuensparseMCSourceDeltaPhiMaps[0] = binct;
1497 valuedeltaphicontamination[1] = binct;
1498
1499 //////////////////////
1500 // run number
1501 //////////////////////
1502
1503 Int_t runnumber = fInputEvent->GetRunNumber();
1504 AliDebug(2,Form("Run number %d",runnumber));
1505
1506 if(!fPID->IsInitialized()){
1507 // Initialize PID with the given run number
1508 fPID->InitializePID(runnumber);
1509 }
1510 if(!fPIDTOFOnly->IsInitialized()){
1511 // Initialize PID with the given run number
1512 fPIDTOFOnly->InitializePID(runnumber);
1513 }
1514
1515 //
1516 if(!fPIDBackground->IsInitialized()){
1517 // Initialize PID with the given run number
1518 fPIDBackground->InitializePID(runnumber);
1519 }
1520
1521 fHFECuts->SetRecEvent(fInputEvent);
1522 if(mcEvent) fHFECuts->SetMCEvent(mcEvent);
1523
1524
1525 //////////
1526 // PID
1527 //////////
1528
1529 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
1530 if(!pidResponse){
1531 AliDebug(2,"No PID response set");
1532 return;
1533 }
1534 fPID->SetPIDResponse(pidResponse);
1535 fPIDTOFOnly->SetPIDResponse(pidResponse);
1536 fPIDBackground->SetPIDResponse(pidResponse);
1537 if(fMonitorPhotonic) fBackgroundSubtraction->InitRun(fInputEvent,pidResponse);
1538
1539 fHistEV->Fill(binctt,0.0);
1540
1541 //////////////////
1542 // Event cut
1543 //////////////////
1544 if(!fHFECuts->CheckEventCuts("fEvRecCuts", fInputEvent)) {
1545 AliDebug(2,"Does not pass the event cut");
1546 PostData(1, fListHist);
1547 return;
1548 }
1549
1550 fHistEV->Fill(binctt,1.0);
1551
1552
1553 ///////////////////////////////////////////////////////////
1554 // PileUpCut
1555 ///////////////////////////////////////////////////////////
1556
1557 Float_t multTPC(0.); // tpc mult estimate
1558 Float_t multGlob(0.); // global multiplicity
1559 const Int_t nGoodTracks = fInputEvent->GetNumberOfTracks();
1560 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill tpc mult
1561 AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*>(fInputEvent->GetTrack(iTracks));
1562 if (!trackAOD) continue;
1563 if (!(trackAOD->TestFilterBit(1))) continue;
1564 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;
1565 multTPC++;
1566 }
1567 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill global mult
1568 AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*>(fInputEvent->GetTrack(iTracks));
1569 if (!trackAOD) continue;
1570 if (!(trackAOD->TestFilterBit(16))) continue;
1571 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;
1572 Double_t b[2] = {-99., -99.};
1573 Double_t bCov[3] = {-99., -99., -99.};
1574 if (!(trackAOD->PropagateToDCA(fInputEvent->GetPrimaryVertex(), fInputEvent->GetMagneticField(), 100., b, bCov))) continue;
1575 if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
1576 multGlob++;
1577 } //track loop
1578
1579 Double_t pileup[4];
1580 pileup[0]=fInputEvent->GetCentrality()->GetCentralityPercentile("V0M");
1581 pileup[1]=fInputEvent->GetCentrality()->GetCentralityPercentile("TRK");
1582 pileup[2]=multTPC;
1583 pileup[3]=multGlob;
1584 fHistPileUp->Fill(pileup);
1585
1586 if(fPileUpCut){
1587 if (TMath::Abs(pileup[0]-pileup[1]) > 5) {
1588 AliDebug(2,"Does not pass the centrality correlation cut");
1589 return;
1590 }
1591 if(multTPC < (-36.81+1.48*multGlob) && multTPC > (63.03+1.78*multGlob)){
1592 AliDebug(2,"Does not pass the multiplicity correlation cut");
1593 return;
1594 }
1595 }
1596
1597 // AliVVZERO* vzeroData=fInputEvent->GetVZEROData();
1598 // Double_t mult[3],multV0A(0),multV0C(0);
1599 // for(Int_t i=0; i<32; ++i) {
1600 // multV0A += vzeroData->GetMultiplicityV0A(i);
1601 // multV0C += vzeroData->GetMultiplicityV0C(i);
1602 // }
1603
1604 // int ntrk=0;
1605 // for(Int_t k = 0; k < fInputEvent->GetNumberOfTracks(); k++){
1606 // AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);
1607 // if(!track) continue;
1608 // if(!(track->GetStatus()&AliVTrack::kITSrefit)) continue;
1609 // if(!(track->GetStatus()&AliVTrack::kTPCrefit)) continue;
1610 // ntrk++;
1611 // }
1612
1613 // mult[0]=fInputEvent->GetNumberOfTracks();
1614 // mult[1]=multV0A+multV0C;
1615 // mult[2]=binctMore;
1616 // fHistPileUp->Fill(mult);
1617
1618 // if(fUpperPileUpCut&&fLowerPileUpCut){
1619 // if((mult[0]<fLowerPileUpCut->Eval(mult[1])) ||
1620 // (mult[0]>fUpperPileUpCut->Eval(mult[1]))){
1621 // AliDebug(2,"Does not pass the pileup cut");
1622 // PostData(1, fListHist);
1623 // return;
1624 // }
1625 // }
1626
1627 ////////////////////////////////////
1628 // First method event plane
1629 ////////////////////////////////////
1630
1631 AliEventplane* vEPa = fInputEvent->GetEventplane();
1632 Float_t eventPlanea = 0.0;
1633 Float_t eventPlaneTPC = 0.0;
1634 Float_t eventPlaneV0A = 0.0;
1635 Float_t eventPlaneV0C = 0.0;
1636 Float_t eventPlaneV0 = 0.0;
1637 TVector2 *qTPC = 0x0;
1638 TVector2 *qsub1a = 0x0;
1639 TVector2 *qsub2a = 0x0;
1640 TVector2 qV0A,qV0C,qV0,*qAna;
1641
1642 // V0
1643
1644 if(fHFEVZEROEventPlane && (!fAODAnalysis)){
1645
1646 //AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
1647 //if(!esd) return;
1648
1649 fHFEVZEROEventPlane->ProcessEvent(fInputEvent);
1650
1651 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0A()+100) < 0.0000001) eventPlaneV0A = -100.0;
1652 else {
1653 eventPlaneV0A = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0A());
1654 if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();
1655 }
1656
1657 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0C()+100) < 0.0000001) eventPlaneV0C = -100.0;
1658 else {
1659 eventPlaneV0C = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0C());
1660 if(eventPlaneV0C > TMath::Pi()) eventPlaneV0C = eventPlaneV0C - TMath::Pi();
1661 }
1662
1663 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0()+100) < 0.0000001) eventPlaneV0 = -100.0;
1664 else {
1665 eventPlaneV0 = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0());
1666 if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();
1667 }
1668
1669 }
1670 else {
1671
1672 Double_t qVx, qVy; //TR: info
1673 eventPlaneV0 = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,10,2,qVx,qVy));
1674 if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();
1675 qV0.Set(qVx,qVy);
1676 eventPlaneV0A = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,8,2,qVx,qVy));
1677 if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();
1678 qV0A.Set(qVx,qVy);
1679 eventPlaneV0C = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,9,2,qVx,qVy));
1680 if(eventPlaneV0C > TMath::Pi()) eventPlaneV0C = eventPlaneV0C - TMath::Pi();
1681 qV0C.Set(qVx,qVy);
1682
1683 if(eventPlaneV0<-900) return;
1684 if(eventPlaneV0A<-900) return;
1685 if(eventPlaneV0C<-900) return;
1686
1687
1688 eventPlaneV0=TVector2::Phi_0_2pi(eventPlaneV0);
1689 eventPlaneV0A=TVector2::Phi_0_2pi(eventPlaneV0A);
1690 eventPlaneV0C=TVector2::Phi_0_2pi(eventPlaneV0C);
1691 }
1692
1693
1694 // TPC
1695
1696 qTPC = vEPa->GetQVector();
1697 Double_t qx = -1.0;
1698 Double_t qy = -1.0;
1699 if(qTPC) {
1700 qx = qTPC->X();
1701 qy = qTPC->Y();
1702 }
1703 TVector2 qVectorfortrack;
1704 qVectorfortrack.Set(qx,qy);
1705 eventPlaneTPC = TVector2::Phi_0_2pi(qVectorfortrack.Phi())/2.;
1706
1707 // Choose the one used for v2
1708
1709 if(fVZEROEventPlane){ //TR: info
1710 eventPlanea = eventPlaneV0;
1711 qAna = &qV0;
1712 }
1713 if(fVZEROEventPlaneA){
1714 eventPlanea = eventPlaneV0A;
1715 qAna = &qV0A;
1716 }
1717 if(fVZEROEventPlaneC){
1718 eventPlanea = eventPlaneV0C;
1719 qAna = &qV0C;
1720 }
1721 if(!fVZEROEventPlane){
1722 eventPlanea = eventPlaneTPC;
1723 qAna = &qV0C;
1724 }
1725
1726 valuecossinephiep[0] = TMath::Cos(2*eventPlanea);
1727 valuecossinephiep[1] = TMath::Sin(2*eventPlanea);
1728
1729 Float_t eventPlanesub1a = -100.0;
1730 Float_t eventPlanesub2a = -100.0;
1731 Double_t diffsub1sub2a = -100.0;
1732 Double_t diffsub1sub2asin = -100.0;
1733 Double_t diffsubasubb = -100.0;
1734 Double_t diffsubasubc = -100.0;
1735 Double_t diffsubbsubc = -100.0;
1736 Double_t diffsubasubbsin = -100.0;
1737 Double_t diffsubasubcsin = -100.0;
1738 Double_t diffsubbsubcsin = -100.0;
1739
1740 // two sub event TPC
1741 qsub1a = vEPa->GetQsub1();
1742 qsub2a = vEPa->GetQsub2();
1743
1744 /////////////////////////////////////////////////////////
1745 // Cut for event with event plane reconstructed by all
1746 ////////////////////////////////////////////////////////
1747
1748 if((!qTPC) || (!qsub1a) || (!qsub2a)) {
1749 AliDebug(2,"No event plane");
1750 return;
1751 }
1752
1753 eventPlanesub1a = TVector2::Phi_0_2pi(qsub1a->Phi())/2.;
1754 eventPlanesub2a = TVector2::Phi_0_2pi(qsub2a->Phi())/2.;
1755 diffsub1sub2a = TMath::Cos(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));
1756 diffsub1sub2asin = TMath::Sin(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));
1757
1758
1759 // if ( !fDebugStreamer ) {
1760 // //debug stream
1761 // TDirectory *backup = gDirectory;
1762 // fDebugStreamer = new TTreeSRedirector("TaskFlowTPCTOFEPSPdebug.root");
1763 // if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
1764 // }
1765
1766 // {
1767
1768 // double v0nrom = TMath::Sqrt(qV0.X()*qV0.X()+qV0.Y()*qV0.Y());
1769 // double v0Anrom = TMath::Sqrt(qV0A.X()*qV0A.X()+qV0A.Y()*qV0A.Y());
1770 // double v0Cnrom = TMath::Sqrt(qV0C.X()*qV0C.X()+qV0C.Y()*qV0C.Y());
1771 // double sub1nrom = TMath::Sqrt(qsub1a->X()*qsub1a->X()+qsub1a->Y()*qsub1a->Y());
1772 // double sub2nrom = TMath::Sqrt(qsub2a->X()*qsub2a->X()+qsub2a->Y()*qsub2a->Y());
1773
1774 // (* fDebugStreamer) << "UserExec" <<
1775 // "binct="<<binct<<
1776 // "qV0="<<v0nrom<<
1777 // "qV0A="<<v0Anrom<<
1778 // "qV0C="<<v0Cnrom<<
1779 // "qsub1a="<<sub1nrom<<
1780 // "qsub2a="<<sub2nrom<<
1781 // "\n";
1782 // }
1783
1784 // three sub events in case of VZEROA and VZEROC
1785 if(!fSP){
1786 diffsubasubb = TMath::Cos(2.*(eventPlaneV0A - eventPlaneV0C)); //TR:
1787 diffsubasubc = TMath::Cos(2.*(eventPlaneV0A - eventPlaneTPC)); //TR:
1788 diffsubbsubc = TMath::Cos(2.*(eventPlaneV0C - eventPlaneTPC)); //TR:
1789 }
1790 else{
1791 if(fVZEROEventPlaneA){
1792 diffsubasubb = qV0A.X()*qV0C.X()+qV0A.Y()*qV0C.Y();
1793 diffsubasubc = qV0A.X()*qTPC->X()+qV0A.Y()*qTPC->Y();
1794 diffsubbsubc = qV0C.X()*qTPC->X()+qV0C.Y()*qTPC->Y();
1795 }
1796 else if(fVZEROEventPlaneC){
1797 diffsubasubb = qV0C.X()*qV0A.X()+qV0C.Y()*qV0A.Y();
1798 diffsubasubc = qV0C.X()*qTPC->X()+qV0C.Y()*qTPC->Y();
1799 diffsubbsubc = qV0A.X()*qTPC->X()+qV0A.Y()*qTPC->Y();
1800 }
1801 }
1802
1803 diffsubasubbsin = TMath::Sin(2.*(eventPlaneV0A - eventPlaneV0C));
1804 diffsubasubcsin = TMath::Sin(2.*(eventPlaneV0A - eventPlaneTPC));
1805 diffsubbsubcsin = TMath::Sin(2.*(eventPlaneV0C - eventPlaneTPC));
1806 // three sub events in case of VZERO all
1807 if(fVZEROEventPlane && (!fVZEROEventPlaneA) && (!fVZEROEventPlaneC)) {
1808 if(!fSP){
1809 diffsubasubb = TMath::Cos(2.*(eventPlaneV0 - eventPlanesub1a)); //TR:
1810 diffsubasubc = TMath::Cos(2.*(eventPlaneV0 - eventPlanesub2a)); //TR:
1811 diffsubbsubc = TMath::Cos(2.*(eventPlanesub1a - eventPlanesub2a)); //TR:
1812 }
1813 else{
1814 diffsubasubb = qV0.X()*qsub1a->X()+qV0.Y()*qsub1a->Y();
1815 diffsubasubc = qV0.X()*qsub2a->X()+qV0.Y()*qsub2a->Y();
1816 diffsubbsubc = qsub1a->X()*qsub2a->X()+qsub1a->Y()*qsub2a->Y();
1817 }
1818
1819 diffsubasubbsin = TMath::Sin(2.*(eventPlaneV0 - eventPlanesub1a));
1820 diffsubasubcsin = TMath::Sin(2.*(eventPlaneV0 - eventPlanesub2a));
1821 diffsubbsubcsin = TMath::Sin(2.*(eventPlanesub1a - eventPlanesub2a));
1822 }
1823
1824 //////////////////////////////////////
1825 // AliFlowEvent and MC event plane
1826 /////////////////////////////////////
1827
1828 Int_t nbtracks = fInputEvent->GetNumberOfTracks();
1829 AliDebug(2,Form("Number of tracks %d",nbtracks));
1830
1831 if(fMonitorQCumulant) {
1832
1833 fcutsRP->SetEvent( InputEvent(), MCEvent());
1834 fcutsPOI->SetEvent( InputEvent(), MCEvent());
1835 if( fflowEvent ){
1836 fflowEvent->~AliFlowEvent();
1837 new(fflowEvent) AliFlowEvent(fcutsRP,fcutsPOI);
1838 }else fflowEvent = new AliFlowEvent(fcutsRP,fcutsPOI);
1839 if(mcEvent && mcEvent->GenEventHeader()) {
1840 fflowEvent->SetMCReactionPlaneAngle(mcEvent);
1841 //if reaction plane not set from elsewhere randomize it before adding flow
1842 //if (!fflowEvent->IsSetMCReactionPlaneAngle()) fflowEvent->SetMCReactionPlaneAngle(gRandom->Uniform(0.0,TMath::TwoPi()));
1843 mcReactionPlane = TVector2::Phi_0_2pi(fflowEvent->GetMCReactionPlaneAngle());
1844 if(mcReactionPlane > TMath::Pi()) mcReactionPlane = mcReactionPlane - TMath::Pi();
1845 AliDebug(2,Form("MC reaction plane %f",mcReactionPlane));
1846 }
1847 fflowEvent->SetReferenceMultiplicity( nbtracks );
1848 fflowEvent->DefineDeadZone(0,0,0,0);
1849 //fflowEvent.TagSubeventsInEta(-0.8,-0.1,0.1,0.8);
1850
1851 ////////////////
1852 // MC
1853 ///////////////
1854 if(fUseMCReactionPlane) {
1855 eventPlanea = mcReactionPlane;
1856 diffsub1sub2a = 0.0;
1857 }
1858 }
1859
1860
1861 //////////////////////
1862 // Fill Histos
1863 //////////////////////
1864
1865 fHistEV->Fill(binctt,2.0);
1866
1867 // Fill
1868 valuensparsea[0] = eventPlaneV0A;
1869 valuensparsea[1] = eventPlaneV0C;
1870 valuensparsea[2] = eventPlaneTPC;
1871 if(fVZEROEventPlane && (!fVZEROEventPlaneA) && (!fVZEROEventPlaneC)) {
1872 // case VZERO all
1873 valuensparsea[0] = eventPlaneV0;
1874 valuensparsea[1] = eventPlanesub1a;
1875 valuensparsea[2] = eventPlanesub2a;
1876 }
1877 fEventPlane->Fill(&valuensparsea[0]);
1878
1879 // Fill
1880 if(fMonitorEventPlane) fCosSin2phiep->Fill(&valuecossinephiep[0]);
1881
1882 if(!fVZEROEventPlane) {
1883 valuensparsef[0] = diffsub1sub2a;
1884 fCosRes->Fill(&valuensparsef[0]);
1885 valuensparsefsin[0] = diffsub1sub2asin;
1886 if(fMonitorEventPlane) fSinRes->Fill(&valuensparsefsin[0]);
1887 if(fMonitorEventPlane) {
1888 fProfileCosRes->Fill(valuensparsef[1],valuensparsef[0]);
1889 }
1890 }
1891 else {
1892 valuensparsefbis[0] = diffsubasubb;
1893 valuensparsefbis[1] = diffsubasubc;
1894 valuensparsefbis[2] = diffsubbsubc;
1895 fCosResabc->Fill(&valuensparsefbis[0]); //TR: info
1896 valuensparsefbissin[0] = diffsubasubbsin;
1897 valuensparsefbissin[1] = diffsubbsubcsin;
1898 valuensparsefbissin[2] = diffsubasubcsin;
1899 if(fMonitorEventPlane) fSinResabc->Fill(&valuensparsefbissin[0]);
1900 if(fMonitorEventPlane) {
1901 fProfileCosResab->Fill(valuensparsefbis[3],valuensparsefbis[0]);
1902 fProfileCosResac->Fill(valuensparsefbis[3],valuensparsefbis[1]);
1903 fProfileCosResbc->Fill(valuensparsefbis[3],valuensparsefbis[2]);
1904 }
1905 }
1906
1907 ////////////////////////////////////////
1908 // Loop to determine pool background
1909 /////////////////////////////////////////
1910 if(fMonitorPhotonic) {
1911
1912 fBackgroundSubtraction->FillPoolAssociatedTracks(fInputEvent,binct);
1913
1914 if( fArraytrack ){
1915 fArraytrack->~TArrayI();
1916 new(fArraytrack) TArrayI(nbtracks);
1917 }
1918 else {
1919 fArraytrack = new TArrayI(nbtracks);
1920 }
1921 fCounterPoolBackground = 0;
1922
1923 for(Int_t k = 0; k < nbtracks; k++){
1924
1925 AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);
1926 if(!track) continue;
1927
1928 // Track cuts
1929 Bool_t survivedbackground = kTRUE;
1930 if(fAODAnalysis) {
1931 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
1932 if(aodtrack) {
1933 AliESDtrack esdTrack(aodtrack);
1934 // set the TPC cluster info
1935 esdTrack.SetTPCClusterMap(aodtrack->GetTPCClusterMap());
1936 esdTrack.SetTPCSharedMap(aodtrack->GetTPCSharedMap());
1937 esdTrack.SetTPCPointsF(aodtrack->GetTPCNclsF());
1938 // needed to calculate the impact parameters
1939 AliAODEvent *aodeventu = dynamic_cast<AliAODEvent *>(fInputEvent);
1940 if(aodeventu) {
1941 AliAODVertex *vAOD = aodeventu->GetPrimaryVertex();
1942 Double_t bfield = aodeventu->GetMagneticField();
1943 Double_t pos[3],cov[6];
1944 vAOD->GetXYZ(pos);
1945 vAOD->GetCovarianceMatrix(cov);
1946 const AliESDVertex vESD(pos,cov,100.,100);
1947 esdTrack.RelateToVertex(&vESD,bfield,3.);
1948 }
1949 if(!fHFEBackgroundCuts->IsSelected(&esdTrack)) {
1950 survivedbackground = kFALSE;
1951 }
1952 }
1953 }
1954 else {
1955 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
1956 if(esdtrack) {
1957 if(!fHFEBackgroundCuts->IsSelected(esdtrack)) survivedbackground = kFALSE;
1958 }
1959 }
1960 // PID
1961 if(survivedbackground) {
1962 // PID track cuts
1963 AliHFEpidObject hfetrack2;
1964 if(!fAODAnalysis) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1965 else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);
1966 hfetrack2.SetRecTrack(track);
1967 hfetrack2.SetCentrality((Int_t)binct);
1968 AliDebug(2,Form("centrality %f and %d",binct,hfetrack2.GetCentrality()));
1969 hfetrack2.SetPbPb();
1970 if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundqa)) {
1971 fArraytrack->AddAt(k,fCounterPoolBackground);
1972 fCounterPoolBackground++;
1973 AliDebug(2,Form("fCounterPoolBackground %d, track %d",fCounterPoolBackground,k));
1974 }
1975 }
1976 }
1977 }
1978
1979
1980 //////////////////////////
1981 // Loop over track
1982 //////////////////////////
1983 for(Int_t k = 0; k < nbtracks; k++){
1984
1985 AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);
1986 if(!track) continue;
1987
1988 if(fAODAnalysis) {
1989 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
1990 if(!aodtrack){
1991 AliError("AOD track is not there");
1992 continue;
1993 }
1994 AliDebug(2,"Find AOD track on");
1995 if(!(aodtrack->TestFilterBit(fFilter))) continue; // Only process AOD tracks where the HFE is set
1996 }
1997
1998
1999 valuetrackingcuts[0] = track->Pt();
2000 valuetrackingcuts[1] = 0;
2001
2002 // RecKine: ITSTPC cuts
2003 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
2004 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
2005
2006 // RecPrim
2007 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
2008 valuetrackingcuts[1] = 1;
2009 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
2010
2011 // HFEcuts: ITS layers cuts
2012 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
2013 valuetrackingcuts[1] = 2;
2014 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
2015
2016 // HFE cuts: TOF and mismatch flag
2017 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTOF + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
2018 valuetrackingcuts[1] = 3;
2019 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
2020
2021 // HFE cuts: TPC PID cleanup
2022 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
2023 valuetrackingcuts[1] = 4;
2024 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
2025
2026 // HFEcuts: Nb of tracklets TRD0
2027 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
2028 valuetrackingcuts[1] = 5;
2029 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
2030
2031 AliDebug(2,"Survived");
2032
2033 /////////////////////////////////////////////////////////
2034 // Subtract candidate from TPC event plane
2035 ////////////////////////////////////////////////////////
2036 Float_t eventplanesubtracted = 0.0;
2037
2038 if(!fVZEROEventPlane) {
2039 // Subtract the tracks from the event plane
2040 Double_t qX = qTPC->X() - vEPa->GetQContributionX(track); //Modify the components: subtract the track you want to look at with your analysis
2041 Double_t qY = qTPC->Y() - vEPa->GetQContributionY(track); //Modify the components: subtract the track you want to look at with your analysis
2042 TVector2 newQVectorfortrack;
2043 newQVectorfortrack.Set(qX,qY);
2044 eventplanesubtracted = TVector2::Phi_0_2pi(newQVectorfortrack.Phi())/2;
2045 }
2046 else eventplanesubtracted = eventPlanea;
2047
2048 ///////////////////////////////////////////
2049 // Event plane
2050 //////////////////////////////////////////
2051 Bool_t fillEventPlane = kTRUE;
2052 if(!fVZEROEventPlane){
2053 if((!qsub1a) || (!qsub2a)) fillEventPlane = kFALSE;
2054 if(fSubEtaGapTPC) {
2055 if(track->Eta() < (- fEtaGap/2.)) eventplanesubtracted = eventPlanesub1a;
2056 else if(track->Eta() > (fEtaGap/2.)) eventplanesubtracted = eventPlanesub2a;
2057 else fillEventPlane = kFALSE;
2058 }
2059 }
2060
2061 ///////////////
2062 // MC
2063 //////////////
2064 if(fUseMCReactionPlane) {
2065 eventplanesubtracted = mcReactionPlane;
2066 fillEventPlane = kTRUE;
2067 }
2068
2069 //////////////////////////////////////////////////////////////////////////////
2070 ///////////////////////////AFTERBURNER
2071 Double_t phitrack = track->Phi();
2072 if (fAfterBurnerOn)
2073 {
2074 phitrack = GetPhiAfterAddV2(track->Phi(),mcReactionPlane);
2075 }
2076 //////////////////////////////////////////////////////////////////////////////
2077
2078
2079 ///////////////////////
2080 // Calculate deltaphi
2081 ///////////////////////
2082
2083 // Suppose phi track is between 0.0 and phi
2084 Double_t deltaphi = TVector2::Phi_0_2pi(phitrack - eventplanesubtracted);
2085 if(deltaphi > TMath::Pi()) deltaphi = deltaphi - TMath::Pi();
2086
2087 ////////////////////////////////
2088 // Determine the deltaphi bin
2089 ///////////////////////////////
2090
2091 // in-plane
2092 if(((deltaphi<(TMath::Pi()/4.)) && (deltaphi>0.0)) || ((deltaphi>(3*TMath::Pi()/4.)) && (deltaphi<TMath::Pi()))) valuedeltaphicontamination[0] = 0.5;
2093 // out-of-plane
2094 if((deltaphi>(TMath::Pi()/4.)) && (deltaphi<(3*TMath::Pi()/4.))) valuedeltaphicontamination[0] = 1.5;
2095
2096 ////////////////////////////////////////
2097 // Define variables
2098 ///////////////////////////////////////
2099
2100
2101 valuedeltaphicontamination[2] = track->Pt();
2102 valuensparsee[2] = track->Pt();
2103 valuensparsee[3] = track->Eta();
2104 valuensparseg[2] = track->Pt();
2105 valuensparseh[2] = track->Pt();
2106 valuefractioncont[0] = track->Pt();
2107 valuensparsehprofile[2] = track->Pt();
2108 valuensparseMCSourceDeltaPhiMaps[1] = track->Pt();
2109 if(track->Charge() > 0.0) {
2110 valuensparseg[3] = 0.2;
2111 valuensparseh[3] = 0.2;
2112 }
2113 else {
2114 valuensparseg[3] = -0.2;
2115 valuensparseh[3] = -0.2;
2116 }
2117 valuensparseh[4] = track->Eta();
2118 valuensparseg[4] = track->Eta();
2119
2120 AliDebug(2,Form("charge %d",(Int_t)track->Charge()));
2121
2122 ////////////////////////
2123 // Fill before PID
2124 ///////////////////////
2125
2126 if(fMonitorWithoutPID) {
2127
2128 valuensparseg[0] = deltaphi;
2129 if(fillEventPlane) fDeltaPhiMapsBeforePID->Fill(&valuensparseg[0]);
2130
2131 //
2132 valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));
2133 if(fillEventPlane) {
2134 fCosPhiMapsBeforePID->Fill(&valuensparseh[0]);
2135 }
2136 }
2137
2138 ////////////////////////
2139 // Apply PID
2140 ////////////////////////
2141 if(!fNoPID) {
2142 // Apply PID for Data
2143 if(!fMCPID) {
2144 // pid object
2145 AliHFEpidObject hfetrack;
2c0951cd 2146 if(!fAODAnalysis){
2147 hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1a319cb5 2148 if(fVariableMultiplicity==0) hfetrack.SetMulitplicity(cntr);
2149 //if(((AliESDEvent*)fInputEvent)->GetPrimaryVertexSPD()) hfetrack.SetMulitplicity(((AliESDEvent*)fInputEvent)->GetPrimaryVertexSPD()->GetNContributors());
38be5083 2150 if(fVariableMultiplicity==1)
8b835997 2151 hfetrack.SetMulitplicity(((AliESDEvent*)fInputEvent)->GetNumberOfESDTracks()/8.);
2c0951cd 2152 }else{
2153 hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);
1a319cb5 2154 if(fVariableMultiplicity==0) hfetrack.SetMulitplicity(cntr);
2155 //if(((AliAODEvent*)fInputEvent)->GetPrimaryVertexSPD()) hfetrack.SetMulitplicity(((AliAODEvent*)fInputEvent)->GetPrimaryVertexSPD()->GetNContributors());
38be5083 2156 if(fVariableMultiplicity==1)
8b835997 2157 hfetrack.SetMulitplicity(((AliAODEvent*)fInputEvent)->GetNumberOfESDTracks()/8.);
2c0951cd 2158 }
c683985a 2159 hfetrack.SetRecTrack(track);
2160 hfetrack.SetCentrality((Int_t)binct);
2161 AliDebug(2,Form("centrality %f and %d",binct,hfetrack.GetCentrality()));
2162 hfetrack.SetPbPb();
2163
2164 // Only TOF PID
2165 if(fMonitorContamination) {
2166 if(fPIDTOFOnly->IsSelected(&hfetrack,0x0,"recTrackCont",0x0)) {
2167 Float_t nsigma = pidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
2168 valuedeltaphicontamination[3] = nsigma;
2169 fDeltaPhiMapsContamination->Fill(&valuedeltaphicontamination[0]);
2170 }
2171 }
2172
2173 // Complete PID TOF+TPC
2174 if(!fPID->IsSelected(&hfetrack,0x0,"recTrackCont",fPIDqa)) {
2175 continue;
2176 }
2177 }
2178 else {
2179 if(!mcEvent) continue;
2180 if(!(mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;
2181 AliDebug(2,Form("PdgCode %d",TMath::Abs(mctrack->Particle()->GetPdgCode())));
2182 if(TMath::Abs(mctrack->Particle()->GetPdgCode())!=11) continue;
2183 }
2184 }
2185
2186
2187 /////////////////////////////////////////////////////////////////////////////
2188 // Add candidate to AliFlowEvent for POI and subtract from RP if needed
2189 ////////////////////////////////////////////////////////////////////////////
2190 if(fMonitorQCumulant) {
2191 Int_t idtrack = static_cast<AliVTrack*>(track)->GetID();
2192 Bool_t found = kFALSE;
2193 Int_t numberoffound = 0;
2194 AliDebug(2,Form("A: Number of tracks %d",fflowEvent->NumberOfTracks()));
2195 for(Int_t iRPs=0; iRPs< fflowEvent->NumberOfTracks(); iRPs++) {
2196 AliFlowTrack *iRP = (AliFlowTrack*) (fflowEvent->GetTrack(iRPs));
2197 if(!iRP) continue;
2198 //if(!iRP->InRPSelection()) continue;
2199 if( TMath::Abs(idtrack) == TMath::Abs(iRP->GetID()) ) {
2200 iRP->SetForPOISelection(kTRUE);
2201 found = kTRUE;
2202 numberoffound ++;
2203 }
2204 }
2205 AliDebug(2,Form("Found %d mal",numberoffound));
2206 if(!found) {
2207 AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*) MakeTrack(massElectron,track->Pt(),track->Phi(), track->Eta());
2208 sTrack->SetID(idtrack);
2209 fflowEvent->AddTrack(sTrack);
2210 AliDebug(2,"Add the track");
2211 }
2212 AliDebug(2,Form("B: Number of tracks %d",fflowEvent->NumberOfTracks()));
2213 }
2214
2215
2216 /////////////////////
2217 // Fill THnSparseF
2218 /////////////////////
2219
2220 //
2221 valuensparseabis[0] = eventplanesubtracted;
2222 if((fillEventPlane) && (fMonitorEventPlane)) fEventPlaneaftersubtraction->Fill(&valuensparseabis[0]);
2223
2224
2225 if(fMonitorEventPlane)
2226 {
2227 //
2228 valuensparsee[0] = TMath::Cos(2*phitrack);
2229 fCos2phie->Fill(&valuensparsee[0]);
2230 valuensparsee[0] = TMath::Sin(2*phitrack);
2231 fSin2phie->Fill(&valuensparsee[0]);
2232 //
2233 valuensparsee[0] = TMath::Cos(2*eventplanesubtracted);
2234 if(fillEventPlane) fCos2phiep->Fill(&valuensparsee[0]);
2235 valuensparsee[0] = TMath::Sin(2*eventplanesubtracted);
2236 if(fillEventPlane) fSin2phiep->Fill(&valuensparsee[0]);
2237 valuensparsee[0] = TMath::Sin(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));
2238 if(fillEventPlane) fSin2phiephiep->Fill(&valuensparsee[0]);
2239 //
2240 }
2241
2242 //
2243 valuensparseg[0] = deltaphi;
2244 if(fillEventPlane) fDeltaPhiMaps->Fill(&valuensparseg[0]);
2245
2246 //
2247 valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));
2248 if(fillEventPlane) {
2249 fCosPhiMaps->Fill(&valuensparseh[0]); //TR: fCosPhiQSum+=valuensparseh[0]*TMath:Sqrt(qAna->X()*qAna->X()+qAna->Y()*qAna->Y()); fCosPhiQN++;
2250 if((valuefractioncont[1] >=0) && (valuefractioncont[1] < 11)){
2251 if(fContamination[((Int_t)valuefractioncont[1])]){
58a496d1 2252 Double_t weight = 1.;
2253 if(fAsFunctionOfP) weight = fContamination[((Int_t)valuefractioncont[1])]->Eval(track->P());
2254 else weight = fContamination[((Int_t)valuefractioncont[1])]->Eval(track->Pt());
c683985a 2255 if(weight<0.0) weight=0.0;
2256 if(weight>1.0) weight=1.0;
2257 fFractionContamination->Fill(&valuefractioncont[0],weight);
2258 if(fv2contamination[((Int_t)valuefractioncont[1])]){
2259 Double_t v2 = fv2contamination[((Int_t)valuefractioncont[1])]->Eval(track->Pt());
2260 AliDebug(2,Form("value v2 %f, contamination %f and pt %f centrality %d\n",v2,weight,track->Pt(),(Int_t)valuefractioncont[1]));
2261 AliDebug(2,Form("Check for centrality 3: value v2 %f, contamination %f\n",fv2contamination[3]->Eval(track->Pt()),fContamination[3]->Eval(track->P())));
2262 AliDebug(2,Form("Check for centrality 4: value v2 %f, contamination %f\n",fv2contamination[4]->Eval(track->Pt()),fContamination[4]->Eval(track->P())));
2263 AliDebug(2,Form("Check for centrality 5: value v2 %f, contamination %f\n",fv2contamination[5]->Eval(track->Pt()),fContamination[5]->Eval(track->P())));
2264 fContaminationv2->Fill(valuefractioncont[1],valuefractioncont[0],v2,weight);
2265 }
58a496d1 2266 fContaminationmeanpt->Fill(valuefractioncont[1],valuefractioncont[0],TMath::Abs(track->Pt()));
c683985a 2267 }
2268 }
2269 if(fMonitorEventPlane) {
2270 if(fSP)
2271 valuensparseh[0] *= TMath::Sqrt(qAna->X()*qAna->X()+qAna->Y()*qAna->Y());
2272 fProfileCosPhiMaps->Fill(valuensparsehprofile[1],valuensparsehprofile[2],valuensparseh[0]); //TR: info
2273 }
2274 }
2275
2276 if(fMonitorPhotonic) {
2277 Int_t indexmother = -1;
2278 Int_t source = 1;
2279 if(mcthere) source = fBackgroundSubtraction->FindMother(mctrack->GetLabel(),indexmother);
2280 fBackgroundSubtraction->LookAtNonHFE(k, track, fInputEvent, 1, binct, deltaphi, source, indexmother);
2281
2282 if((!fAODAnalysis && mcthere) || !mcthere) {
2283 // background
2284 source = 0;
2285 indexmother = -1;
2286 source = FindMother(TMath::Abs(track->GetLabel()),mcEvent, indexmother);
2287 valuensparseMCSourceDeltaPhiMaps[2] = source;
2288 if(mcEvent) fMCSourceDeltaPhiMaps->Fill(&valuensparseMCSourceDeltaPhiMaps[0]);
2289 //LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);
2290 Int_t taggedvalue = LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);
2291 if(fMonitorPhotonic) {
2292 // No opposite charge partner found in the invariant mass choosen
2293 if((taggedvalue!=2) && (taggedvalue!=6)) {
2294 //fDeltaPhiMapsTaggedNonPhotonic->Fill(&valuensparseg[0]);
2295 //fCosPhiMapsTaggedNonPhotonic->Fill(&valuensparseh[0]);
2296 }
2297 // One opposite charge partner found in the invariant mass choosen
2298 if((taggedvalue==2) || (taggedvalue==6)) {
2299 fDeltaPhiMapsTaggedPhotonic->Fill(&valuensparseg[0]);
2300 //fCosPhiMapsTaggedPhotonic->Fill(&valuensparseh[0]);
2301 }
2302 // One same charge partner found in the invariant mass choosen
2303 if((taggedvalue==4) || (taggedvalue==6)) {
2304 fDeltaPhiMapsTaggedPhotonicLS->Fill(&valuensparseg[0]);
2305 //fCosPhiMapsTaggedPhotonicLS->Fill(&valuensparseh[0]);
2306 }
2307 }
2308 }
2309 }
2310
2311 }
2312
2313 //////////////////////////////////////////////////////////////////////////////
2314 ///////////////////////////AFTERBURNER
2315 if (fAfterBurnerOn & fMonitorQCumulant)
2316 {
2317 fflowEvent->AddFlow(fV1,fV2,fV3,fV4,fV5); //add flow
2318 fflowEvent->CloneTracks(fNonFlowNumberOfTrackClones); //add nonflow by cloning tracks
2319 }
2320 //////////////////////////////////////////////////////////////////////////////
2321
2322
2323
2324 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {
2325 // if((fBinCentralityLess[bincless]< cntr) && (cntr < fBinCentralityLess[bincless+1])) PostData(bincless+2,fflowEvent);
2326 //}
2327
2328 if(fMonitorPhotonic) {
2329 if(fArraytrack) {
2330 delete fArraytrack;
2331 fArraytrack = NULL;
2332 }
2333 }
2334
2335 if(fMonitorPhotonic) fBackgroundSubtraction->CountPoolAssociated(fInputEvent,binct);
2336
2337 PostData(1, fListHist);
2338
2339
2340
2341}
2342//______________________________________________________________________________
2343AliFlowCandidateTrack *AliAnalysisTaskFlowTPCTOFEPSP::MakeTrack( Double_t mass,
2344 Double_t pt, Double_t phi, Double_t eta) {
2345 //
2346 // Make Track (Not needed actually)
2347 //
2348
2349 AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
2350 sTrack->SetMass(mass);
2351 sTrack->SetPt(pt);
2352 sTrack->SetPhi(phi);
2353 sTrack->SetEta(eta);
2354 sTrack->SetForPOISelection(kTRUE);
2355 sTrack->SetForRPSelection(kFALSE);
2356 return sTrack;
2357}
2358//_________________________________________________________________________________
2359Double_t AliAnalysisTaskFlowTPCTOFEPSP::GetPhiAfterAddV2(Double_t phi,Double_t reactionPlaneAngle) const
2360{
2361 //
2362 // Adds v2, uses Newton-Raphson iteration
2363 //
2364 Double_t phiend=phi;
2365 Double_t phi0=phi;
2366 Double_t f=0.;
2367 Double_t fp=0.;
2368 Double_t phiprev=0.;
2369
2370 for (Int_t i=0; i<fMaxNumberOfIterations; i++)
2371 {
2372 phiprev=phiend; //store last value for comparison
2373 f = phiend-phi0+fV2*TMath::Sin(2.*(phiend-reactionPlaneAngle));
2374 fp = 1.0+2.0*fV2*TMath::Cos(2.*(phiend-reactionPlaneAngle)); //first derivative
2375 phiend -= f/fp;
2376 if (TMath::AreEqualAbs(phiprev,phiend,fPrecisionPhi)) break;
2377 }
2378 return phiend;
2379}
2380//_____________________________________________________________________________________________
2381Int_t AliAnalysisTaskFlowTPCTOFEPSP::LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *vEvent, AliMCEvent *mcEvent,Int_t binct,Double_t deltaphi,Int_t source,Int_t indexmother)
2382{
2383 //
2384 // Look At Non HFE
2385 //
2386
2387 // return -1 if nothing
2388 // return 2 if opposite charge within the mass range found
2389 // return 4 if like charge within the mass range found
2390 // return 6 if opposite charge and like charge within the mass range found
2391 //
2392
2393 Int_t taggedphotonic = -1;
2394
2395 Bool_t oppositetaggedphotonic = kFALSE;
2396 Bool_t sametaggedphotonic = kFALSE;
2397
2398 AliDebug(2,Form("fCounterPoolBackground %d in LookAtNonHFE!!!",fCounterPoolBackground));
2399 if(!fArraytrack) return taggedphotonic;
2400 AliDebug(2,Form("process track %d",iTrack1));
2401
2402 TVector3 v3D1;
2403 TVector3 v3D2;
2404
2405 Double_t valuensparseDeltaPhiMaps[5];
2406 Double_t valueangle[3];
2407
2408 valuensparseDeltaPhiMaps[1] = binct;
2409 valuensparseDeltaPhiMaps[2] = track1->Pt();
2410 valuensparseDeltaPhiMaps[0] = deltaphi;
2411 valuensparseDeltaPhiMaps[4] = source;
2412
2413 valueangle[2] = source;
2414 valueangle[1] = binct;
2415
2416 // Pdg code
2417 Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()),mcEvent);
2418 Int_t numberfound = 0;
2419
2420 //Magnetic Field
2421 Double_t bfield = vEvent->GetMagneticField();
2422
2423 // Get Primary vertex
2424 const AliVVertex *pVtx = vEvent->GetPrimaryVertex();
2425
2426 for(Int_t idex = 0; idex < fCounterPoolBackground; idex++)
2427 {
2428
2429 Int_t iTrack2 = fArraytrack->At(idex);
2430 AliDebug(2,Form("track %d",iTrack2));
2431 AliVTrack* track2 = (AliVTrack *) vEvent->GetTrack(iTrack2);
2432 if (!track2)
2433 {
2434 printf("ERROR: Could not receive track %d", iTrack2);
2435 continue;
2436 }
2437 if(iTrack2==iTrack1) continue;
2438 AliDebug(2,"Different");
2439
2440 // Reset the MC info
2441 valueangle[2] = source;
2442 valuensparseDeltaPhiMaps[4] = source;
2443
2444 // track cuts and PID already done
2445
2446 // if MC look
2447 Int_t pdg2 = -100;
2448 if(mcEvent) {
2449 Int_t source2 = 0;
2450 Int_t indexmother2 = -1;
2451 source2 = FindMother(TMath::Abs(track2->GetLabel()),mcEvent, indexmother2);
2452 pdg2 = CheckPdg(TMath::Abs(track2->GetLabel()),mcEvent);
2453 if(source2 >=0 ) {
2454 if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0)) {
2455 if(source == kElectronfromconversion) {
2456 valueangle[2] = kElectronfromconversionboth;
2457 valuensparseDeltaPhiMaps[4] = kElectronfromconversionboth;
2458 numberfound++;
2459 }
2460 if(source == kElectronfrompi0) {
2461 valueangle[2] = kElectronfrompi0both;
2462 valuensparseDeltaPhiMaps[4] = kElectronfrompi0both;
2463 }
2464 if(source == kElectronfrometa) {
2465 valueangle[2] = kElectronfrometaboth;
2466 valuensparseDeltaPhiMaps[4] = kElectronfrometaboth;
2467 }
2468 }
2469 }
2470 }
2471
2472 if(fAlgorithmMA && (!fAODAnalysis))
2473 {
2474 // tracks
2475 AliESDtrack *esdtrack2 = dynamic_cast<AliESDtrack *>(track2);
2476 AliESDtrack *esdtrack1 = dynamic_cast<AliESDtrack *>(track1);
2477 if((!esdtrack2) || (!esdtrack1)) continue;
2478
2479 //Variables
2480 Double_t p1[3];
2481 Double_t p2[3];
2482 Double_t xt1; //radial position track 1 at the DCA point
2483 Double_t xt2; //radial position track 2 at the DCA point
2484 //DCA track1-track2
2485 Double_t dca12 = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1);
2486
2487 // Cut dca
2488 if(dca12 > fMaxdca) continue;
2489
2490 //Momento of the track extrapolated to DCA track-track
2491 //Track1
2492 Bool_t hasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1);
2493 //Track2
2494 Bool_t hasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2);
2495
2496 if(!hasdcaT1 || !hasdcaT2) AliWarning("It could be a problem in the extrapolation");
2497
2498 //track1-track2 Invariant Mass
2499 Double_t eMass = 0.000510998910; //Electron mass in GeV
2500 Double_t pP1 = sqrt(p1[0]*p1[0]+p1[1]*p1[1]+p1[2]*p1[2]); //Track 1 momentum
2501 Double_t pP2 = sqrt(p2[0]*p2[0]+p2[1]*p2[1]+p2[2]*p2[2]); //Track 2 momentum
2502 Double_t eE1 = TMath::Sqrt(pP1*pP1+eMass*eMass);
2503 Double_t eE2 = TMath::Sqrt(pP2*pP2+eMass*eMass);
2504
2505 //TLorentzVector v1(p1[0],p1[1],p1[2],sqrt(eMass*eMass+pP1*pP1));
2506 //TLorentzVector v2(p2[0],p2[1],p2[2],sqrt(eMass*eMass+pP2*pP2));
2507 //Double_t imass = (v1+v2).M(); //Invariant Mass
2508 //Double_t angle3D = v1.Angle(v2.Vect()); //Opening Angle (Total Angle)
2509
2510 // daughter
2511 v3D1.SetXYZ(p1[0],p1[1],p1[2]);
2512 v3D2.SetXYZ(p2[0],p2[1],p2[2]);
2513 Double_t openingangle = TVector2::Phi_0_2pi(v3D2.Angle(v3D1));
2514
2515 // mother
2516 TVector3 motherrec = v3D1 + v3D2;
2517 Double_t invmass = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));
2518
2519 // xy
2520 //TVector3 vectordiff = v3D1 - v3D2;
2521 //Double_t diffphi = TVector2::Phi_0_2pi(vectordiff.Phi());
2522 //Double_t massxy = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(diffphi)));
2523
2524 // rz
2525 //Double_t difftheta = TVector2::Phi_0_2pi(vectordiff.Eta());
2526 //Double_t massrz = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(difftheta)));
2527
2528
2529 Float_t fCharge1 = track1->Charge();
2530 Float_t fCharge2 = track2->Charge();
2531
2532 // Fill Histo
2533 //valueangle[0] = diffphi;
2534 //valueangle[1] = difftheta;
2535 valueangle[0] = openingangle;
2536 if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);
2537 else fOppSignAngle->Fill(&valueangle[0]);
2538
2539 // Cut
2540 if(openingangle > fMaxopening3D) continue;
2541 //if(difftheta > fMaxopeningtheta) continue;
2542 //if(diffphi > fMaxopeningphi) continue;
2543
2544 // Invmass
2545 valuensparseDeltaPhiMaps[3] = invmass;
2546 if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2547 else fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2548
2549 // Cut
2550 if(invmass < fMaxInvmass) {
2551 if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;
2552 if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;
2553 }
2554
2555
2556 }
2557 else
2558 {
2559 Int_t fPDGtrack1 = 11;
2560 Int_t fPDGtrack2 = 11;
2561
2562 Float_t fCharge1 = track1->Charge();
2563 Float_t fCharge2 = track2->Charge();
2564
2565 if(fCharge1>0) fPDGtrack1 = -11;
2566 if(fCharge2>0) fPDGtrack2 = -11;
2567
2568 AliKFParticle ktrack1(*track1, fPDGtrack1);
2569 AliKFParticle ktrack2(*track2, fPDGtrack2);
2570 AliKFParticle recoGamma(ktrack1, ktrack2);
2571
2572 //Reconstruction Cuts
2573 if(recoGamma.GetNDF()<1) continue;
2574 Double_t chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();
2575 if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;
2576
2577 // DCA
2578 //Double_t dca12 = ktrack1.GetDistanceFromParticle(ktrack2);
2579 //if(dca12 > fMaxdca) continue;
2580
2581 // if set mass constraint
2582 if(fSetMassConstraint && pVtx) {
2583 AliKFVertex primV(*pVtx);
2584 primV += recoGamma;
2585 primV -= ktrack1;
2586 primV -= ktrack2;
2587 recoGamma.SetProductionVertex(primV);
2588 recoGamma.SetMassConstraint(0,0.0001);
2589 }
2590
2591 //Invariant Mass
2592 Double_t imass;
2593 Double_t width;
2594 recoGamma.GetMass(imass,width);
2595
2596 //Opening Angle (Total Angle)
2597 Double_t angle = ktrack1.GetAngle(ktrack2);
2598 valueangle[0] = angle;
2599 if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);
2600 else fOppSignAngle->Fill(&valueangle[0]);
2601
2602 // Cut
2603 if(angle > fMaxopening3D) continue;
2604
2605 // Invmass
2606 valuensparseDeltaPhiMaps[3] = imass;
2607 if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2608 else {
2609 fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2610 /*
2611 if(valueangle[2] == kElectronfromconversionboth) {
2612 printf("Reconstructed charge1 %f, charge 2 %f and invmass %f",fCharge1,fCharge2,imass);
2613 printf("MC charge1 %d, charge 2 %d",pdg1,pdg2);
2614 printf("DCA %f",dca12);
2615 printf("Number of found %d",numberfound);
2616 }
2617 */
2618 }
2619
2620 // Cut
2621 if(imass < fMaxInvmass) {
2622 if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;
2623 if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;
2624 }
2625 }
2626 }
2627
2628 if(oppositetaggedphotonic && sametaggedphotonic){
2629 taggedphotonic = 6;
2630 }
2631
2632 if(!oppositetaggedphotonic && sametaggedphotonic){
2633 taggedphotonic = 4;
2634 }
2635
2636 if(oppositetaggedphotonic && !sametaggedphotonic){
2637 taggedphotonic = 2;
2638 }
2639
2640
2641 return taggedphotonic;
2642}
2643//_________________________________________________________________________
2644Int_t AliAnalysisTaskFlowTPCTOFEPSP::FindMother(Int_t tr, AliMCEvent *mcEvent, Int_t &indexmother){
2645 //
2646 // Find the mother if MC
2647 //
2648
2649 if(!mcEvent) return 0;
2650
2651 Int_t pdg = CheckPdg(tr,mcEvent);
2652 if(TMath::Abs(pdg)!= 11) {
2653 indexmother = -1;
2654 return kNoElectron;
2655 }
2656
2657 indexmother = IsMotherGamma(tr,mcEvent);
2658 if(indexmother > 0) return kElectronfromconversion;
2659 indexmother = IsMotherPi0(tr,mcEvent);
2660 if(indexmother > 0) return kElectronfrompi0;
2661 indexmother = IsMotherC(tr,mcEvent);
2662 if(indexmother > 0) return kElectronfromC;
2663 indexmother = IsMotherB(tr,mcEvent);
2664 if(indexmother > 0) return kElectronfromB;
2665 indexmother = IsMotherEta(tr,mcEvent);
2666 if(indexmother > 0) return kElectronfrometa;
2667
2668 return kElectronfromother;
2669
2670
2671}
2672//____________________________________________________________________________________________________________
2673Int_t AliAnalysisTaskFlowTPCTOFEPSP::CheckPdg(Int_t tr, AliMCEvent* mcEvent) {
2674
2675 //
2676 // Return the pdg of the particle
2677 //
2678
2679
2680 Int_t pdgcode = -1;
2681 if(tr < 0) return pdgcode;
2682
2683 if(!mcEvent) return pdgcode;
2684
2685 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2686
2687
2688 if(mctrack->IsA() == AliMCParticle::Class()) {
2689 AliMCParticle *mctrackesd = NULL;
2690 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;
2691 pdgcode = mctrackesd->PdgCode();
2692 }
2693
2694 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2695 AliAODMCParticle *mctrackaod = NULL;
2696 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;
2697 pdgcode = mctrackaod->GetPdgCode();
2698 }
2699
2700 return pdgcode;
2701
2702
2703}
2704//____________________________________________________________________________________________________________
2705Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherGamma(Int_t tr, AliMCEvent* mcEvent) {
2706
2707 //
2708 // Return the lab of gamma mother or -1 if not gamma
2709 //
2710
2711 if(tr < 0) return -1;
2712 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2713
2714 if(mctrack->IsA() == AliMCParticle::Class()) {
2715 AliMCParticle *mctrackesd = NULL;
2716 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2717 TParticle *particle = 0x0;
2718 particle = mctrackesd->Particle();
2719 // Take mother
2720 if(!particle) return -1;
2721 Int_t imother = particle->GetFirstMother();
2722 if(imother < 0) return -1;
2723 AliMCParticle *mothertrack = NULL;
2724 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2725 TParticle * mother = mothertrack->Particle();
2726 if(!mother) return -1;
2727 // Check gamma
2728 Int_t pdg = mother->GetPdgCode();
2729 if(TMath::Abs(pdg) == 22) return imother;
2730 if(TMath::Abs(pdg) == 11) {
2731 return IsMotherGamma(imother,mcEvent);
2732 }
2733 return -1;
2734 }
2735
2736 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2737 AliAODMCParticle *mctrackaod = NULL;
2738 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2739 // Take mother
2740 Int_t imother = mctrackaod->GetMother();
2741 if(imother < 0) return -1;
2742 AliAODMCParticle *mothertrack = NULL;
2743 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2744 // Check gamma
2745 Int_t pdg = mothertrack->GetPdgCode();
2746 if(TMath::Abs(pdg) == 22) return imother;
2747 if(TMath::Abs(pdg) == 11) {
2748 return IsMotherGamma(imother,mcEvent);
2749 }
2750 return -1;
2751
2752 }
2753
2754 return -1;
2755
2756
2757}
2758//
2759//____________________________________________________________________________________________________________
2760Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherPi0(Int_t tr, AliMCEvent* mcEvent) {
2761
2762 //
2763 // Return the lab of pi0 mother or -1 if not pi0
2764 //
2765
2766 if(tr < 0) return -1;
2767 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2768
2769 if(mctrack->IsA() == AliMCParticle::Class()) {
2770 AliMCParticle *mctrackesd = NULL;
2771 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2772 TParticle *particle = 0x0;
2773 particle = mctrackesd->Particle();
2774 // Take mother
2775 if(!particle) return -1;
2776 Int_t imother = particle->GetFirstMother();
2777 if(imother < 0) return -1;
2778 AliMCParticle *mothertrack = NULL;
2779 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2780 TParticle * mother = mothertrack->Particle();
2781 if(!mother) return -1;
2782 // Check gamma
2783 Int_t pdg = mother->GetPdgCode();
2784 if(TMath::Abs(pdg) == 111) return imother;
2785 if(TMath::Abs(pdg) == 11) {
2786 return IsMotherPi0(imother,mcEvent);
2787 }
2788 return -1;
2789 }
2790
2791 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2792 AliAODMCParticle *mctrackaod = NULL;
2793 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2794 // Take mother
2795 Int_t imother = mctrackaod->GetMother();
2796 if(imother < 0) return -1;
2797 AliAODMCParticle *mothertrack = NULL;
2798 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2799 // Check gamma
2800 Int_t pdg = mothertrack->GetPdgCode();
2801 if(TMath::Abs(pdg) == 111) return imother;
2802 if(TMath::Abs(pdg) == 11) {
2803 return IsMotherPi0(imother,mcEvent);
2804 }
2805 return -1;
2806 }
2807
2808 return -1;
2809
2810}
2811//____________________________________________________________________________________________________________
2812Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherC(Int_t tr, AliMCEvent* mcEvent) {
2813
2814 //
2815 // Return the lab of signal mother or -1 if not signal
2816 //
2817
2818 if(tr < 0) return -1;
2819 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2820
2821 if(mctrack->IsA() == AliMCParticle::Class()) {
2822 AliMCParticle *mctrackesd = NULL;
2823 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2824 TParticle *particle = 0x0;
2825 particle = mctrackesd->Particle();
2826 // Take mother
2827 if(!particle) return -1;
2828 Int_t imother = particle->GetFirstMother();
2829 if(imother < 0) return -1;
2830 AliMCParticle *mothertrack = NULL;
2831 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2832 TParticle * mother = mothertrack->Particle();
2833 if(!mother) return -1;
2834 // Check gamma
2835 Int_t pdg = mother->GetPdgCode();
2836 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;
2837 if(TMath::Abs(pdg) == 11) {
2838 return IsMotherC(imother,mcEvent);
2839 }
2840 return -1;
2841 }
2842
2843 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2844 AliAODMCParticle *mctrackaod = NULL;
2845 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2846 // Take mother
2847 Int_t imother = mctrackaod->GetMother();
2848 if(imother < 0) return -1;
2849 AliAODMCParticle *mothertrack = NULL;
2850 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2851 // Check gamma
2852 Int_t pdg = mothertrack->GetPdgCode();
2853 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;
2854 if(TMath::Abs(pdg) == 11) {
2855 return IsMotherC(imother,mcEvent);
2856 }
2857 return -1;
2858 }
2859
2860 return -1;
2861
2862}
2863//____________________________________________________________________________________________________________
2864Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherB(Int_t tr, AliMCEvent* mcEvent) {
2865
2866 //
2867 // Return the lab of signal mother or -1 if not signal
2868 //
2869
2870 if(tr < 0) return -1;
2871 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2872
2873 if(mctrack->IsA() == AliMCParticle::Class()) {
2874 AliMCParticle *mctrackesd = NULL;
2875 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2876 TParticle *particle = 0x0;
2877 particle = mctrackesd->Particle();
2878 // Take mother
2879 if(!particle) return -1;
2880 Int_t imother = particle->GetFirstMother();
2881 if(imother < 0) return -1;
2882 AliMCParticle *mothertrack = NULL;
2883 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2884 TParticle * mother = mothertrack->Particle();
2885 if(!mother) return -1;
2886 // Check gamma
2887 Int_t pdg = mother->GetPdgCode();
2888 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;
2889 if(TMath::Abs(pdg) == 11) {
2890 return IsMotherB(imother,mcEvent);
2891 }
2892 return -1;
2893 }
2894
2895 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2896 AliAODMCParticle *mctrackaod = NULL;
2897 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2898 // Take mother
2899 Int_t imother = mctrackaod->GetMother();
2900 if(imother < 0) return -1;
2901 AliAODMCParticle *mothertrack = NULL;
2902 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2903 // Check gamma
2904 Int_t pdg = mothertrack->GetPdgCode();
2905 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;
2906 if(TMath::Abs(pdg) == 11) {
2907 return IsMotherB(imother,mcEvent);
2908 }
2909 return -1;
2910 }
2911
2912 return -1;
2913
2914}
2915//____________________________________________________________________________________________________________
2916Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherEta(Int_t tr, AliMCEvent* mcEvent) {
2917
2918 //
2919 // Return the lab of pi0 mother or -1 if not pi0
2920 //
2921
2922 if(tr < 0) return -1;
2923 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2924
2925 if(mctrack->IsA() == AliMCParticle::Class()) {
2926 AliMCParticle *mctrackesd = NULL;
2927 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2928 TParticle *particle = 0x0;
2929 particle = mctrackesd->Particle();
2930 // Take mother
2931 if(!particle) return -1;
2932 Int_t imother = particle->GetFirstMother();
2933 if(imother < 0) return -1;
2934 AliMCParticle *mothertrack = NULL;
2935 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2936 TParticle * mother = mothertrack->Particle();
2937 if(!mother) return -1;
2938 // Check gamma
2939 Int_t pdg = mother->GetPdgCode();
2940 if(TMath::Abs(pdg) == 221) return imother;
2941 if(TMath::Abs(pdg) == 11) {
2942 return IsMotherEta(imother,mcEvent);
2943 }
2944 return -1;
2945 }
2946
2947 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2948 AliAODMCParticle *mctrackaod = NULL;
2949 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2950 // Take mother
2951 Int_t imother = mctrackaod->GetMother();
2952 if(imother < 0) return -1;
2953 AliAODMCParticle *mothertrack = NULL;
2954 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2955 // Check gamma
2956 Int_t pdg = mothertrack->GetPdgCode();
2957 if(TMath::Abs(pdg) == 221) return imother;
2958 if(TMath::Abs(pdg) == 11) {
2959 return IsMotherEta(imother,mcEvent);
2960 }
2961 return -1;
2962 }
2963
2964 return -1;
2965
2966}