1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 // Raphaelle Bailhache <R.Bailhache@gsi.de>
20 // Theodor Rascanu <trascanu@stud.uni-frankfurt.de>
27 #include "THnSparse.h"
31 #include "TProfile2D.h"
32 #include "TLorentzVector.h"
33 #include "TParticle.h"
37 #include <TDirectory.h>
38 #include <TTreeStream.h>
40 #include "AliVEventHandler.h"
41 #include "AliAnalysisTaskSE.h"
42 #include "AliAnalysisManager.h"
44 #include "AliVEvent.h"
45 #include "AliESDInputHandler.h"
46 #include "AliMCEvent.h"
48 #include "AliESDEvent.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"
63 #include "AliMCEvent.h"
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"
73 #include "AliHFEcuts.h"
74 #include "AliHFEpid.h"
75 #include "AliHFEpidQAmanager.h"
76 #include "AliHFEtools.h"
77 #include "AliHFEVZEROEventPlane.h"
79 #include "AliCentrality.h"
80 #include "AliEventplane.h"
81 #include "AliAnalysisTaskFlowTPCTOFEPSP.h"
82 #include "AliAODMCHeader.h"
83 #include "TClonesArray.h"
84 #include "AliHFENonPhotonicElectron.h"
86 ClassImp(AliAnalysisTaskFlowTPCTOFEPSP)
88 //____________________________________________________________________
89 AliAnalysisTaskFlowTPCTOFEPSP::AliAnalysisTaskFlowTPCTOFEPSP() :
95 fAODArrayMCInfo(NULL),
96 fBackgroundSubtraction(NULL),
97 fVZEROEventPlane(kFALSE),
98 fVZEROEventPlaneA(kFALSE),
99 fVZEROEventPlaneC(kFALSE),
100 fSubEtaGapTPC(kFALSE),
103 fNbBinsCentralityQCumulant(4),
104 fNbBinsPtQCumulant(12),
105 fMinPtQCumulant(0.2),
106 fMaxPtQCumulant(6.0),
107 fAfterBurnerOn(kFALSE),
108 fNonFlowNumberOfTrackClones(0),
114 fMaxNumberOfIterations(100),
115 fPrecisionPhi(0.001),
116 fUseMCReactionPlane(kFALSE),
120 fChi2OverNDFCut(3.0),
122 fMaxopeningtheta(0.02),
126 fSetMassConstraint(kFALSE),
127 fMonitorEventPlane(kFALSE),
128 fMonitorContamination(kFALSE),
129 fMonitorPhotonic(kFALSE),
130 fMonitorWithoutPID(kFALSE),
131 fMonitorTrackCuts(kFALSE),
132 fMonitorQCumulant(kFALSE),
140 fAsFunctionOfP(kTRUE),
141 fHFEBackgroundCuts(0),
146 fCounterPoolBackground(0),
147 fHFEVZEROEventPlane(0x0),
152 fEventPlaneaftersubtraction(0x0),
153 fFractionContamination(0x0),
154 fContaminationv2(0x0),
155 fContaminationmeanpt(0x0),
164 fProfileCosResab(0x0),
165 fProfileCosResac(0x0),
166 fProfileCosResbc(0x0),
171 fDeltaPhiMapsBeforePID(0x0),
172 fCosPhiMapsBeforePID(0x0),
174 fDeltaPhiMapsContamination(0x0),
176 fProfileCosPhiMaps(0x0),
177 fDeltaPhiMapsTaggedPhotonic(0x0),
178 //fCosPhiMapsTaggedPhotonic(0x0),
179 fDeltaPhiMapsTaggedNonPhotonic(0x0),
180 //fCosPhiMapsTaggedNonPhotonic(0x0),
181 fDeltaPhiMapsTaggedPhotonicLS(0x0),
182 //fCosPhiMapsTaggedPhotonicLS(0x0),
183 fMCSourceDeltaPhiMaps(0x0),
184 fOppSignDeltaPhiMaps(0x0),
185 fSameSignDeltaPhiMaps(0x0),
192 for(Int_t k = 0; k < 10; k++) {
193 fBinCentralityLess[k] = 0.0;
195 for(Int_t k = 0; k < 11; k++) {
196 fContamination[k] = NULL;
197 fv2contamination[k] = NULL;
201 //______________________________________________________________________________
202 AliAnalysisTaskFlowTPCTOFEPSP:: AliAnalysisTaskFlowTPCTOFEPSP(const char *name) :
203 AliAnalysisTaskSE(name),
205 fAODAnalysis(kFALSE),
208 fAODArrayMCInfo(NULL),
209 fBackgroundSubtraction(NULL),
210 fVZEROEventPlane(kFALSE),
211 fVZEROEventPlaneA(kFALSE),
212 fVZEROEventPlaneC(kFALSE),
213 fSubEtaGapTPC(kFALSE),
216 fNbBinsCentralityQCumulant(4),
217 fNbBinsPtQCumulant(15),
218 fMinPtQCumulant(0.0),
219 fMaxPtQCumulant(6.0),
220 fAfterBurnerOn(kFALSE),
221 fNonFlowNumberOfTrackClones(0),
227 fMaxNumberOfIterations(100),
228 fPrecisionPhi(0.001),
229 fUseMCReactionPlane(kFALSE),
233 fChi2OverNDFCut(3.0),
235 fMaxopeningtheta(0.02),
239 fSetMassConstraint(kFALSE),
240 fMonitorEventPlane(kFALSE),
241 fMonitorContamination(kFALSE),
242 fMonitorPhotonic(kFALSE),
243 fMonitorWithoutPID(kFALSE),
244 fMonitorTrackCuts(kFALSE),
245 fMonitorQCumulant(kFALSE),
253 fAsFunctionOfP(kTRUE),
254 fHFEBackgroundCuts(0),
259 fCounterPoolBackground(0),
260 fHFEVZEROEventPlane(0x0),
265 fEventPlaneaftersubtraction(0x0),
266 fFractionContamination(0x0),
267 fContaminationv2(0x0),
268 fContaminationmeanpt(0x0),
277 fProfileCosResab(0x0),
278 fProfileCosResac(0x0),
279 fProfileCosResbc(0x0),
284 fDeltaPhiMapsBeforePID(0x0),
285 fCosPhiMapsBeforePID(0x0),
287 fDeltaPhiMapsContamination(0x0),
289 fProfileCosPhiMaps(0x0),
290 fDeltaPhiMapsTaggedPhotonic(0x0),
291 //fCosPhiMapsTaggedPhotonic(0x0),
292 fDeltaPhiMapsTaggedNonPhotonic(0x0),
293 //fCosPhiMapsTaggedNonPhotonic(0x0),
294 fDeltaPhiMapsTaggedPhotonicLS(0x0),
295 //fCosPhiMapsTaggedPhotonicLS(0x0),
296 fMCSourceDeltaPhiMaps(0x0),
297 fOppSignDeltaPhiMaps(0x0),
298 fSameSignDeltaPhiMaps(0x0),
307 for(Int_t k = 0; k < 10; k++) {
308 fBinCentralityLess[k] = 0.0;
310 fBinCentralityLess[0] = 0.0;
311 fBinCentralityLess[1] = 20.0;
312 fBinCentralityLess[2] = 40.0;
313 fBinCentralityLess[3] = 60.0;
314 fBinCentralityLess[4] = 80.0;
316 for(Int_t k = 0; k < 11; k++) {
317 fContamination[k] = NULL;
318 fv2contamination[k] = NULL;
321 fPID = new AliHFEpid("hfePid");
322 fPIDqa = new AliHFEpidQAmanager;
324 fPIDBackground = new AliHFEpid("hfePidBackground");
325 fPIDBackgroundqa = new AliHFEpidQAmanager;
327 fPIDTOFOnly = new AliHFEpid("hfePidTOFOnly");
329 DefineInput(0,TChain::Class());
330 DefineOutput(1, TList::Class());
331 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {
332 // DefineOutput(bincless+2,AliFlowEventSimple::Class());
336 //____________________________________________________________
337 AliAnalysisTaskFlowTPCTOFEPSP::AliAnalysisTaskFlowTPCTOFEPSP(const AliAnalysisTaskFlowTPCTOFEPSP &ref):
338 AliAnalysisTaskSE(ref),
340 fAODAnalysis(ref.fAODAnalysis),
341 fFilter(ref.fFilter),
342 fAODMCHeader(ref.fAODMCHeader),
343 fAODArrayMCInfo(ref.fAODArrayMCInfo),
344 fBackgroundSubtraction(ref.fBackgroundSubtraction),
345 fVZEROEventPlane(ref.fVZEROEventPlane),
346 fVZEROEventPlaneA(ref.fVZEROEventPlaneA),
347 fVZEROEventPlaneC(ref.fVZEROEventPlaneC),
348 fSubEtaGapTPC(ref.fSubEtaGapTPC),
349 fEtaGap(ref.fEtaGap),
350 fPtBinning(ref.fPtBinning),
351 fNbBinsCentralityQCumulant(ref.fNbBinsCentralityQCumulant),
352 fNbBinsPtQCumulant(ref.fNbBinsPtQCumulant),
353 fMinPtQCumulant(ref.fMinPtQCumulant),
354 fMaxPtQCumulant(ref.fMaxPtQCumulant),
355 fAfterBurnerOn(ref.fAfterBurnerOn),
356 fNonFlowNumberOfTrackClones(ref.fNonFlowNumberOfTrackClones),
362 fMaxNumberOfIterations(ref.fMaxNumberOfIterations),
363 fPrecisionPhi(ref.fPrecisionPhi),
364 fUseMCReactionPlane(ref.fUseMCReactionPlane),
368 fChi2OverNDFCut(ref.fChi2OverNDFCut),
369 fMaxdca(ref.fMaxdca),
370 fMaxopeningtheta(ref.fMaxopeningtheta),
371 fMaxopeningphi(ref.fMaxopeningphi),
372 fMaxopening3D(ref.fMaxopening3D),
373 fMaxInvmass(ref.fMaxInvmass),
374 fSetMassConstraint(ref.fSetMassConstraint),
375 fMonitorEventPlane(ref.fMonitorEventPlane),
376 fMonitorContamination(ref.fMonitorContamination),
377 fMonitorPhotonic(ref.fMonitorPhotonic),
378 fMonitorWithoutPID(ref.fMonitorWithoutPID),
379 fMonitorTrackCuts(ref.fMonitorTrackCuts),
380 fMonitorQCumulant(ref.fMonitorQCumulant),
388 fAsFunctionOfP(ref.fAsFunctionOfP),
389 fHFEBackgroundCuts(NULL),
390 fPIDBackground(NULL),
391 fPIDBackgroundqa(NULL),
392 fAlgorithmMA(ref.fAlgorithmMA),
394 fCounterPoolBackground(ref.fCounterPoolBackground),
395 fHFEVZEROEventPlane(NULL),
400 fEventPlaneaftersubtraction(NULL),
401 fFractionContamination(NULL),
402 fContaminationv2(NULL),
403 fContaminationmeanpt(0x0),
409 fSin2phiephiep(NULL),
412 fProfileCosResab(NULL),
413 fProfileCosResac(NULL),
414 fProfileCosResbc(NULL),
417 fProfileCosRes(NULL),
419 fDeltaPhiMapsBeforePID(NULL),
420 fCosPhiMapsBeforePID(NULL),
422 fDeltaPhiMapsContamination(NULL),
424 fProfileCosPhiMaps(NULL),
425 fDeltaPhiMapsTaggedPhotonic(NULL),
426 //fCosPhiMapsTaggedPhotonic(NULL),
427 fDeltaPhiMapsTaggedNonPhotonic(NULL),
428 //fCosPhiMapsTaggedNonPhotonic(NULL),
429 fDeltaPhiMapsTaggedPhotonicLS(NULL),
430 //fCosPhiMapsTaggedPhotonicLS(NULL),
431 fMCSourceDeltaPhiMaps(NULL),
432 fOppSignDeltaPhiMaps(NULL),
433 fSameSignDeltaPhiMaps(NULL),
435 fSameSignAngle(NULL),
442 for(Int_t k = 0; k < 10; k++) {
443 fBinCentralityLess[k] = 0.0;
445 for(Int_t k = 0; k < 11; k++) {
446 fContamination[k] = NULL;
447 fv2contamination[k] = NULL;
454 //____________________________________________________________
455 AliAnalysisTaskFlowTPCTOFEPSP &AliAnalysisTaskFlowTPCTOFEPSP::operator=(const AliAnalysisTaskFlowTPCTOFEPSP &ref){
457 // Assignment operator
464 //____________________________________________________________
465 void AliAnalysisTaskFlowTPCTOFEPSP::Copy(TObject &o) const {
467 // Copy into object o
469 AliAnalysisTaskFlowTPCTOFEPSP &target = dynamic_cast<AliAnalysisTaskFlowTPCTOFEPSP &>(o);
470 target.fListHist = fListHist;
471 target.fAODAnalysis = fAODAnalysis;
472 target.fFilter = fFilter;
473 target.fAODMCHeader = fAODMCHeader;
474 target.fAODArrayMCInfo = fAODArrayMCInfo;
475 target.fBackgroundSubtraction = fBackgroundSubtraction;
476 target.fVZEROEventPlane = fVZEROEventPlane;
477 target.fVZEROEventPlaneA = fVZEROEventPlaneA;
478 target.fVZEROEventPlaneC = fVZEROEventPlaneC;
479 target.fSubEtaGapTPC = fSubEtaGapTPC;
480 target.fEtaGap = fEtaGap;
481 target.fPtBinning = fPtBinning;
482 target.fNbBinsCentralityQCumulant = fNbBinsCentralityQCumulant;
483 target.fNbBinsPtQCumulant = fNbBinsPtQCumulant;
484 target.fMinPtQCumulant = fMinPtQCumulant;
485 target.fMaxPtQCumulant = fMaxPtQCumulant;
486 target.fAfterBurnerOn = fAfterBurnerOn;
487 target.fNonFlowNumberOfTrackClones = fNonFlowNumberOfTrackClones;
493 target.fMaxNumberOfIterations = fMaxNumberOfIterations;
494 target.fPrecisionPhi = fPrecisionPhi;
495 target.fUseMCReactionPlane = fUseMCReactionPlane;
497 target.fMCPID = fMCPID;
498 target.fNoPID = fNoPID;
499 target.fChi2OverNDFCut = fChi2OverNDFCut;
500 target.fMaxdca = fMaxdca;
501 target.fMaxopeningtheta = fMaxopeningtheta;
502 target.fMaxopeningphi = fMaxopeningphi;
503 target.fMaxopening3D = fMaxopening3D;
504 target.fMaxInvmass = fMaxInvmass;
505 target.fSetMassConstraint = fSetMassConstraint;
506 target.fMonitorEventPlane = fMonitorEventPlane;
507 target.fMonitorContamination = fMonitorContamination;
508 target.fMonitorPhotonic = fMonitorPhotonic;
509 target.fMonitorWithoutPID = fMonitorWithoutPID;
510 target.fMonitorTrackCuts = fMonitorTrackCuts;
511 target.fMonitorQCumulant = fMonitorQCumulant;
512 target.fcutsRP = fcutsRP;
513 target.fcutsPOI = fcutsPOI;
514 target.fHFECuts = fHFECuts;
516 target.fPIDTOFOnly = fPIDTOFOnly;
517 target.fPIDqa = fPIDqa;
518 target.fflowEvent = fflowEvent;
519 target.fAsFunctionOfP = fAsFunctionOfP;
520 target.fHFEBackgroundCuts = fHFEBackgroundCuts;
521 target.fPIDBackground = fPIDBackground;
522 target.fPIDBackgroundqa = fPIDBackgroundqa;
523 target.fAlgorithmMA = fAlgorithmMA;
524 target.fArraytrack = fArraytrack;
525 target.fCounterPoolBackground = fCounterPoolBackground;
526 target.fHFEVZEROEventPlane = fHFEVZEROEventPlane;
527 target.fHistEV=fHistEV;
528 target.fHistPileUp=fHistPileUp;
529 target.fPileUpCut=fPileUpCut;
530 target.fEventPlane=fEventPlane;
531 target.fEventPlaneaftersubtraction=fEventPlaneaftersubtraction;
532 target.fFractionContamination=fFractionContamination;
533 target.fContaminationv2=fContaminationv2;
534 target.fContaminationmeanpt=fContaminationmeanpt;
535 target.fCosSin2phiep=fCosSin2phiep;
536 target.fCos2phie=fCos2phie;
537 target.fSin2phie=fSin2phie;
538 target.fCos2phiep=fCos2phiep;
539 target.fSin2phiep=fSin2phiep;
540 target.fSin2phiephiep=fSin2phiephiep;
541 target.fCosResabc=fCosResabc;
542 target.fSinResabc=fSinResabc;
543 target.fProfileCosResab=fProfileCosResab;
544 target.fProfileCosResac=fProfileCosResac;
545 target.fProfileCosResbc=fProfileCosResbc;
546 target.fCosRes=fCosRes;
547 target.fSinRes=fSinRes;
548 target.fProfileCosRes=fProfileCosRes;
549 target.fTrackingCuts=fTrackingCuts;
550 target.fDeltaPhiMapsBeforePID=fDeltaPhiMapsBeforePID;
551 target.fCosPhiMapsBeforePID=fCosPhiMapsBeforePID;
552 target.fDeltaPhiMaps=fDeltaPhiMaps;
553 target.fDeltaPhiMapsContamination=fDeltaPhiMapsContamination;
554 target.fCosPhiMaps=fCosPhiMaps;
555 target.fProfileCosPhiMaps=fProfileCosPhiMaps;
556 target.fDeltaPhiMapsTaggedPhotonic=fDeltaPhiMapsTaggedPhotonic;
557 target.fDeltaPhiMapsTaggedNonPhotonic=fDeltaPhiMapsTaggedNonPhotonic;
558 target.fDeltaPhiMapsTaggedPhotonicLS=fDeltaPhiMapsTaggedPhotonicLS;
559 target.fMCSourceDeltaPhiMaps=fMCSourceDeltaPhiMaps;
560 target.fOppSignDeltaPhiMaps=fOppSignDeltaPhiMaps;
561 target.fSameSignDeltaPhiMaps=fSameSignDeltaPhiMaps;
562 target.fOppSignAngle=fOppSignAngle;
563 target.fSameSignAngle=fSameSignAngle;
566 for(Int_t k = 0; k < 10; k++) {
567 target.fBinCentralityLess[k] = fBinCentralityLess[k];
569 for(Int_t k = 0; k < 11; k++) {
570 target.fContamination[k] = fContamination[k];
571 target.fv2contamination[k] = fv2contamination[k];
573 target.fDebugStreamer=fDebugStreamer;
575 //____________________________________________________________
576 AliAnalysisTaskFlowTPCTOFEPSP::~AliAnalysisTaskFlowTPCTOFEPSP(){
582 if(fArraytrack) delete fArraytrack;
583 if(fListHist) delete fListHist;
584 if(fcutsRP) delete fcutsRP;
585 if(fcutsPOI) delete fcutsPOI;
586 if(fHFECuts) delete fHFECuts;
587 if(fPID) delete fPID;
588 if(fPIDTOFOnly) delete fPIDTOFOnly;
589 //if(fPIDqa) delete fPIDqa;
590 if(fflowEvent) delete fflowEvent;
591 if(fHFEBackgroundCuts) delete fHFEBackgroundCuts;
592 if(fPIDBackground) delete fPIDBackground;
593 //if(fBackgroundSubtraction) delete fBackgroundSubtraction;
594 //if(fPIDBackgroundqa) delete fPIDBackgroundqa;
595 //if(fHFEVZEROEventPlane) delete fHFEVZEROEventPlane;
596 if ( fDebugStreamer ) delete fDebugStreamer;
600 //________________________________________________________________________
601 void AliAnalysisTaskFlowTPCTOFEPSP::UserCreateOutputObjects()
604 //********************
606 //********************
612 //---------Data selection----------
613 //kMC, kGlobal, kTPCstandalone, kSPDtracklet, kPMD
614 //AliFlowTrackCuts::trackParameterType rptype = AliFlowTrackCuts::kGlobal;
615 //AliFlowTrackCuts::trackParameterType poitype = AliFlowTrackCuts::kGlobal;
617 //---------Parameter mixing--------
618 //kPure - no mixing, kTrackWithMCkine, kTrackWithMCPID, kTrackWithMCpt
619 //AliFlowTrackCuts::trackParameterMix rpmix = AliFlowTrackCuts::kPure;
620 //AliFlowTrackCuts::trackParameterMix poimix = AliFlowTrackCuts::kPure;
622 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: User create output objects");
626 AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
627 if(!TString(inputHandler->IsA()->GetName()).CompareTo("AliAODInputHandler")){
628 SetAODAnalysis(kTRUE);
629 AliDebug(2,"Put AOD analysis on");
631 SetAODAnalysis(kFALSE);
634 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: AOD ESD");
637 fcutsRP = AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010();
638 fcutsRP->SetName("StandartTPC");
639 fcutsRP->SetEtaRange(-0.9,0.9);
640 fcutsRP->SetQA(kTRUE);
641 //TList *qaCutsRP = fcutsRP->GetQA();
642 //qaCutsRP->SetName("QA_StandartTPC_RP");
644 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cutsRP");
647 fcutsPOI = new AliFlowTrackCuts("dummy");
648 fcutsPOI->SetParamType(AliFlowTrackCuts::kGlobal);
649 fcutsPOI->SetPtRange(+1,-1); // select nothing QUICK
650 fcutsPOI->SetEtaRange(+1,-1); // select nothing VZERO
653 fflowEvent->~AliFlowEvent();
654 new(fflowEvent) AliFlowEvent(fcutsRP,fcutsPOI);
656 else fflowEvent = new AliFlowEvent(fcutsRP,fcutsPOI);
658 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cutsPOI");
661 AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();
662 cc->SetNbinsMult(10000);
664 cc->SetMultMax(10000.);
665 cc->SetNbinsPt(fNbBinsPtQCumulant);
666 cc->SetPtMin(fMinPtQCumulant);
667 cc->SetPtMax(fMaxPtQCumulant);
668 cc->SetNbinsPhi(180);
670 cc->SetPhiMax(TMath::TwoPi());
671 cc->SetNbinsEta(200);
678 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: common constants");
683 fHFECuts = new AliHFEcuts;
684 fHFECuts->CreateStandardCuts();
686 fHFECuts->Initialize();
687 if(fAODAnalysis) fHFECuts->SetAOD();
689 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: HFE cuts");
693 //fPID->SetHasMCData(HasMCData());
695 fPID =new AliHFEpid("hfePid");
696 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: pid init 0");
698 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: pid init 1");
699 if(!fPID->GetNumberOfPIDdetectors()) fPID->AddDetector("TPC", 0);
700 AliDebug(2,Form("AliAnalysisTaskFlowTPCTOFEPSP: GetNumber of PID detectors %d",fPID->GetNumberOfPIDdetectors()));
701 fPID->InitializePID();
702 fPIDqa->Initialize(fPID);
703 fPID->SortDetectors();
705 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: pid and pidqa");
707 if(!fPIDTOFOnly->GetNumberOfPIDdetectors()) {
708 fPIDTOFOnly->AddDetector("TOF", 0);
709 fPIDTOFOnly->ConfigureTOF(3.);
711 fPIDTOFOnly->InitializePID();
712 fPIDTOFOnly->SortDetectors();
714 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: pidtof");
716 // HFE Background cuts
718 if(!fHFEBackgroundCuts){
719 fHFEBackgroundCuts = new AliESDtrackCuts();
720 fHFEBackgroundCuts->SetName("nackgroundcuts");
721 //Configure Default Track Cuts
722 fHFEBackgroundCuts->SetAcceptKinkDaughters(kFALSE);
723 fHFEBackgroundCuts->SetRequireTPCRefit(kTRUE);
724 fHFEBackgroundCuts->SetEtaRange(-0.9,0.9);
725 fHFEBackgroundCuts->SetRequireSigmaToVertex(kTRUE);
726 fHFEBackgroundCuts->SetMaxChi2PerClusterTPC(4.0);
727 fHFEBackgroundCuts->SetMinNClustersTPC(50);
728 fHFEBackgroundCuts->SetPtRange(0.3,1e10);
731 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: hfe background");
733 // PID background HFE
734 if(!fPIDBackground->GetNumberOfPIDdetectors()) fPIDBackground->AddDetector("TPC", 0);
735 fPIDBackground->InitializePID();
736 fPIDBackgroundqa->Initialize(fPIDBackground);
737 fPIDBackground->SortDetectors();
739 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: pid background");
741 if (fMonitorPhotonic) {
742 if(!fBackgroundSubtraction) fBackgroundSubtraction = new AliHFENonPhotonicElectron();
743 if(fAODAnalysis) fBackgroundSubtraction->SetAOD(kTRUE);
744 fBackgroundSubtraction->Init();
749 //**************************
750 // Bins for the THnSparse
751 //**************************
755 Double_t minPt = 0.1;
756 Double_t maxPt = 20.0;
757 Double_t binLimLogPt[nBinsPt+1];
758 Double_t binLimPt[nBinsPt+1];
759 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 ;
760 for(Int_t i=0; i<=nBinsPt; i++) binLimPt[i]=(Double_t)TMath::Power(10,binLimLogPt[i]);
764 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,
765 1.3, 1.4, 1.5, 2., 2.5, 3., 4., 6.};
767 if(!fPtBinning.GetSize()) fPtBinning.Set(nBinsPt+1, binLimPt);
769 Int_t nBinsPtPlus = fNbBinsPtQCumulant;
770 Double_t minPtPlus = fMinPtQCumulant;
771 Double_t maxPtPlus = fMaxPtQCumulant;
772 Double_t binLimPtPlus[nBinsPtPlus+1];
773 for(Int_t i=0; i<=nBinsPtPlus; i++) binLimPtPlus[i]=(Double_t)minPtPlus + (maxPtPlus-minPtPlus)/nBinsPtPlus*(Double_t)i ;
776 Double_t minEta = -0.8;
777 Double_t maxEta = 0.8;
778 Double_t binLimEta[nBinsEta+1];
779 for(Int_t i=0; i<=nBinsEta; i++) binLimEta[i]=(Double_t)minEta + (maxEta-minEta)/nBinsEta*(Double_t)i ;
782 Double_t minStep = 0.;
783 Double_t maxStep = 6.;
784 Double_t binLimStep[nBinsStep+1];
785 for(Int_t i=0; i<=nBinsStep; i++) binLimStep[i]=(Double_t)minStep + (maxStep-minStep)/nBinsStep*(Double_t)i ;
787 Int_t nBinsEtaLess = 2;
788 Double_t binLimEtaLess[nBinsEtaLess+1];
789 for(Int_t i=0; i<=nBinsEtaLess; i++) binLimEtaLess[i]=(Double_t)minEta + (maxEta-minEta)/nBinsEtaLess*(Double_t)i ;
792 Double_t minCos = -1.0;
793 Double_t maxCos = 1.0;
794 Double_t binLimCos[nBinsCos+1];
795 for(Int_t i=0; i<=nBinsCos; i++) binLimCos[i]=(Double_t)minCos + (maxCos-minCos)/nBinsCos*(Double_t)i ;
797 // Int_t nBinsCosSP = 50;
798 // Double_t minCosSP = -100.0;
799 // Double_t maxCosSP = 100.0;
800 // Double_t binLimCosSP[nBinsCosSP+1];
801 // for(Int_t i=0; i<=nBinsCosSP; i++) binLimCosSP[i]=(Double_t)minCosSP + (maxCosSP-minCosSP)/nBinsCosSP*(Double_t)i ;
805 Double_t maxC = 11.0;
806 Double_t binLimC[nBinsC+1];
807 for(Int_t i=0; i<=nBinsC; i++) binLimC[i]=(Double_t)minC + (maxC-minC)/nBinsC*(Double_t)i ;
809 Int_t nBinsCMore = 20;
810 Double_t minCMore = 0.0;
811 Double_t maxCMore = 20.0;
812 Double_t binLimCMore[nBinsCMore+1];
813 for(Int_t i=0; i<=nBinsCMore; i++) binLimCMore[i]=(Double_t)minCMore + (maxCMore-minCMore)/nBinsCMore*(Double_t)i ;
815 Int_t nBinsCEvenMore = 100;
816 Double_t minCEvenMore = 0.0;
817 Double_t maxCEvenMore = 100.0;
818 Double_t binLimCEvenMore[nBinsCEvenMore+1];
819 for(Int_t i=0; i<=nBinsCEvenMore; i++) binLimCEvenMore[i]=(Double_t)minCEvenMore + (maxCEvenMore-minCEvenMore)/nBinsCEvenMore*(Double_t)i ;
822 Double_t minPhi = 0.0;
823 Double_t maxPhi = TMath::Pi();
824 Double_t binLimPhi[nBinsPhi+1];
825 for(Int_t i=0; i<=nBinsPhi; i++) {
826 binLimPhi[i]=(Double_t)minPhi + (maxPhi-minPhi)/nBinsPhi*(Double_t)i ;
827 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
830 Int_t nBinsPhiLess = 2.0;
831 Double_t minPhiLess = 0.0;
832 Double_t maxPhiLess = 2.0;
833 Double_t binLimPhiLess[nBinsPhiLess+1];
834 for(Int_t i=0; i<=nBinsPhiLess; i++) {
835 binLimPhiLess[i]=(Double_t)minPhiLess + (maxPhiLess-minPhiLess)/nBinsPhiLess*(Double_t)i ;
838 Int_t nBinsTPCdEdx = 140;
839 Double_t minTPCdEdx = -12.0;
840 Double_t maxTPCdEdx = 12.0;
841 Double_t binLimTPCdEdx[nBinsTPCdEdx+1];
842 for(Int_t i=0; i<=nBinsTPCdEdx; i++) {
843 binLimTPCdEdx[i]=(Double_t)minTPCdEdx + (maxTPCdEdx-minTPCdEdx)/nBinsTPCdEdx*(Double_t)i ;
846 Int_t nBinsAngle = 40;
847 Double_t minAngle = 0.0;
848 Double_t maxAngle = 1.0;
849 Double_t binLimAngle[nBinsAngle+1];
850 for(Int_t i=0; i<=nBinsAngle; i++) {
851 binLimAngle[i]=(Double_t)minAngle + (maxAngle-minAngle)/nBinsAngle*(Double_t)i ;
852 AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
855 Int_t nBinsCharge = 2;
856 Double_t minCharge = -1.0;
857 Double_t maxCharge = 1.0;
858 Double_t binLimCharge[nBinsCharge+1];
859 for(Int_t i=0; i<=nBinsCharge; i++) binLimCharge[i]=(Double_t)minCharge + (maxCharge-minCharge)/nBinsCharge*(Double_t)i ;
861 Int_t nBinsSource = 10;
862 Double_t minSource = 0.;
863 Double_t maxSource = 10.;
864 Double_t binLimSource[nBinsSource+1];
865 for(Int_t i=0; i<=nBinsSource; i++) binLimSource[i]=(Double_t)minSource + (maxSource-minSource)/nBinsSource*(Double_t)i ;
867 Int_t nBinsInvMass = 50;
868 Double_t minInvMass = 0.;
869 Double_t maxInvMass = 0.3;
870 Double_t binLimInvMass[nBinsInvMass+1];
871 for(Int_t i=0; i<=nBinsInvMass; i++) binLimInvMass[i]=(Double_t)minInvMass + (maxInvMass-minInvMass)/nBinsInvMass*(Double_t)i ;
873 Int_t nBinsMult = 100;
874 Double_t minMult = 0.;
875 Double_t maxMult = 3000;
876 Double_t binLimMult[nBinsMult+1];
877 //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);
878 for(Int_t i=0; i<=nBinsMult; i++) binLimMult[i]=(Double_t)minMult + (maxMult-minMult)/nBinsMult*(Double_t)i;
880 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: variables");
886 fListHist = new TList();
887 fListHist->SetOwner();
889 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: list");
894 fHistEV = new TH2D("fHistEV", "events", 3, 0, 3, 3, 0,3);
896 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: histev");
898 // V0 multiplicity vs # of tracks vs centraliy
899 const Int_t nDimPU=4;
900 Int_t nBinPU[nDimPU] = {nBinsCEvenMore,nBinsCEvenMore,nBinsMult,nBinsMult};
901 fHistPileUp = new THnSparseF("PileUp","PileUp",nDimPU,nBinPU);
902 fHistPileUp->SetBinEdges(0,binLimCEvenMore);
903 fHistPileUp->SetBinEdges(1,binLimCEvenMore);
904 fHistPileUp->SetBinEdges(2,binLimMult);
905 fHistPileUp->SetBinEdges(3,binLimMult);
906 fHistPileUp->Sumw2();
908 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: eventplane");
910 // Event plane as function of phiep, centrality
912 Int_t nBina[nDima] = {nBinsPhi,nBinsPhi,nBinsPhi,nBinsC};
913 fEventPlane = new THnSparseF("EventPlane","EventPlane",nDima,nBina);
914 fEventPlane->SetBinEdges(0,binLimPhi);
915 fEventPlane->SetBinEdges(1,binLimPhi);
916 fEventPlane->SetBinEdges(2,binLimPhi);
917 fEventPlane->SetBinEdges(3,binLimC);
918 fEventPlane->Sumw2();
920 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: eventplane");
922 // Fraction of contamination, centrality
923 const Int_t nDimcont=2;
924 Int_t nBincont[nDimcont] = {fPtBinning.GetSize()-1,nBinsC};
925 fFractionContamination = new THnSparseF("Contamination","Contamination",nDimcont,nBincont);
926 fFractionContamination->SetBinEdges(0,fPtBinning.GetArray());
927 fFractionContamination->SetBinEdges(1,binLimC);
928 fFractionContamination->Sumw2();
930 fContaminationv2 = new TProfile2D("Contaminationv2","",nBinsC,binLimC,fPtBinning.GetSize()-1,fPtBinning.GetArray());
931 fContaminationv2->Sumw2();
933 fContaminationmeanpt = new TProfile2D("Contaminationmeanpt","",nBinsC,binLimC,fPtBinning.GetSize()-1,fPtBinning.GetArray());
934 fContaminationmeanpt->Sumw2();
936 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: fraction of contamination");
938 // Resolution cosres_abc centrality
939 const Int_t nDimfbis=4;
940 Int_t nBinfbis[nDimfbis] = {nBinsCos,nBinsCos,nBinsCos,nBinsCMore};
941 fCosResabc = new THnSparseF("CosRes_abc","CosRes_abc",nDimfbis,nBinfbis);
942 fCosResabc->SetBinEdges(0,binLimCos);
943 fCosResabc->SetBinEdges(1,binLimCos);
944 fCosResabc->SetBinEdges(2,binLimCos);
945 fCosResabc->SetBinEdges(3,binLimCMore);
948 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cosresabc");
950 // Resolution cosres centrality
952 Int_t nBinf[nDimf] = {nBinsCos, nBinsCMore};
953 fCosRes = new THnSparseF("CosRes","CosRes",nDimf,nBinf);
954 fCosRes->SetBinEdges(0,binLimCos);
955 fCosRes->SetBinEdges(1,binLimCMore);
958 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cosres");
962 Int_t nBing[nDimg] = {nBinsPhi,nBinsC,fPtBinning.GetSize()-1, nBinsCharge,nBinsEtaLess};
963 fDeltaPhiMaps = new THnSparseF("DeltaPhiMaps","DeltaPhiMaps",nDimg,nBing);
964 fDeltaPhiMaps->SetBinEdges(0,binLimPhi);
965 fDeltaPhiMaps->SetBinEdges(1,binLimC);
966 fDeltaPhiMaps->SetBinEdges(2,fPtBinning.GetArray());
967 fDeltaPhiMaps->SetBinEdges(3,binLimCharge);
968 fDeltaPhiMaps->SetBinEdges(4,binLimEtaLess);
969 fDeltaPhiMaps->Sumw2();
971 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: deltaphimaps");
975 Int_t nBinh[nDimh] = {nBinsCos,nBinsC,fPtBinning.GetSize()-1,nBinsCharge,nBinsEtaLess};
976 fCosPhiMaps = new THnSparseF("CosPhiMaps","CosPhiMaps",nDimh,nBinh);
977 fCosPhiMaps->SetBinEdges(0,binLimCos);
978 fCosPhiMaps->SetBinEdges(1,binLimC);
979 fCosPhiMaps->SetBinEdges(2,fPtBinning.GetArray());
980 fCosPhiMaps->SetBinEdges(3,binLimCharge);
981 fCosPhiMaps->SetBinEdges(4,binLimEtaLess);
982 fCosPhiMaps->Sumw2();
984 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cosphimaps");
987 // fMonitorEventPlane
991 if(fMonitorEventPlane) {
992 // Event Plane after subtraction as function of phiep, centrality, pt, eta
994 Int_t nBinb[nDimb] = {nBinsPhi, nBinsC};
995 fEventPlaneaftersubtraction = new THnSparseF("EventPlane_aftersubtraction","EventPlane_aftersubtraction",nDimb,nBinb);
996 fEventPlaneaftersubtraction->SetBinEdges(0,binLimPhi);
997 fEventPlaneaftersubtraction->SetBinEdges(1,binLimC);
998 fEventPlaneaftersubtraction->Sumw2();
1000 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: eventplane after sub");
1002 // Monitoring of the event Plane cos(2phi) sin(2phi) centrality
1003 const Int_t nDimi=3;
1004 Int_t nBini[nDimi] = {nBinsCos, nBinsCos, nBinsCMore};
1005 fCosSin2phiep = new THnSparseF("CosSin2phiep","CosSin2phiep",nDimi,nBini);
1006 fCosSin2phiep->SetBinEdges(0,binLimCos);
1007 fCosSin2phiep->SetBinEdges(1,binLimCos);
1008 fCosSin2phiep->SetBinEdges(2,binLimCMore);
1009 fCosSin2phiep->Sumw2();
1011 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cossin2phiep");
1013 // Monitoring Event plane after subtraction of the track
1014 const Int_t nDime=4;
1015 Int_t nBine[nDime] = {nBinsCos, nBinsC, fPtBinning.GetSize()-1, nBinsEta};
1016 fCos2phie = new THnSparseF("cos2phie","cos2phie",nDime,nBine);
1017 fCos2phie->SetBinEdges(2,fPtBinning.GetArray());
1018 fCos2phie->SetBinEdges(3,binLimEta);
1019 fCos2phie->SetBinEdges(0,binLimCos);
1020 fCos2phie->SetBinEdges(1,binLimC);
1022 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cos2phie");
1023 fSin2phie = new THnSparseF("sin2phie","sin2phie",nDime,nBine);
1024 fSin2phie->SetBinEdges(2,fPtBinning.GetArray());
1025 fSin2phie->SetBinEdges(3,binLimEta);
1026 fSin2phie->SetBinEdges(0,binLimCos);
1027 fSin2phie->SetBinEdges(1,binLimC);
1029 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: sin2phie");
1030 fCos2phiep = new THnSparseF("cos2phiep","cos2phiep",nDime,nBine);
1031 fCos2phiep->SetBinEdges(2,fPtBinning.GetArray());
1032 fCos2phiep->SetBinEdges(3,binLimEta);
1033 fCos2phiep->SetBinEdges(0,binLimCos);
1034 fCos2phiep->SetBinEdges(1,binLimC);
1035 fCos2phiep->Sumw2();
1036 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cos2phiep");
1037 fSin2phiep = new THnSparseF("sin2phiep","sin2phiep",nDime,nBine);
1038 fSin2phiep->SetBinEdges(2,fPtBinning.GetArray());
1039 fSin2phiep->SetBinEdges(3,binLimEta);
1040 fSin2phiep->SetBinEdges(0,binLimCos);
1041 fSin2phiep->SetBinEdges(1,binLimC);
1042 fSin2phiep->Sumw2();
1043 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: sin2phiep");
1044 fSin2phiephiep = new THnSparseF("sin2phie_phiep","sin2phie_phiep",nDime,nBine);
1045 fSin2phiephiep->SetBinEdges(2,fPtBinning.GetArray());
1046 fSin2phiephiep->SetBinEdges(3,binLimEta);
1047 fSin2phiephiep->SetBinEdges(0,binLimCos);
1048 fSin2phiephiep->SetBinEdges(1,binLimC);
1049 fSin2phiephiep->Sumw2();
1050 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: sin2phiephiep");
1052 const Int_t nDimfbiss=4;
1053 Int_t nBinfbiss[nDimfbiss] = {nBinsCos,nBinsCos,nBinsCos,nBinsC};
1054 fSinResabc = new THnSparseF("SinRes_abc","SinRes_abc",nDimfbiss,nBinfbiss);
1055 fSinResabc->SetBinEdges(0,binLimCos);
1056 fSinResabc->SetBinEdges(1,binLimCos);
1057 fSinResabc->SetBinEdges(2,binLimCos);
1058 fSinResabc->SetBinEdges(3,binLimC);
1059 fSinResabc->Sumw2();
1060 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: sinresabc");
1062 // Profile cosres centrality with 3 subevents
1063 fProfileCosResab = new TProfile("ProfileCosRes_a_b","ProfileCosRes_a_b",nBinsCMore,binLimCMore);
1064 fProfileCosResab->Sumw2();
1065 fProfileCosResac = new TProfile("ProfileCosRes_a_c","ProfileCosRes_a_c",nBinsCMore,binLimCMore);
1066 fProfileCosResac->Sumw2();
1067 fProfileCosResbc = new TProfile("ProfileCosRes_b_c","ProfileCosRes_b_c",nBinsCMore,binLimCMore);
1068 fProfileCosResbc->Sumw2();
1069 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: profilecosresbc");
1072 const Int_t nDimff=2;
1073 Int_t nBinff[nDimff] = {nBinsCos, nBinsC};
1074 fSinRes = new THnSparseF("SinRes","SinRes",nDimff,nBinff);
1075 fSinRes->SetBinEdges(0,binLimCos);
1076 fSinRes->SetBinEdges(1,binLimC);
1078 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: sinres");
1080 // Profile cosres centrality
1081 fProfileCosRes = new TProfile("ProfileCosRes","ProfileCosRes",nBinsCMore,binLimCMore);
1082 fProfileCosRes->Sumw2();
1083 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: profilecosres");
1085 // Profile Maps cos phi
1086 fProfileCosPhiMaps = new TProfile2D("ProfileCosPhiMaps","ProfileCosPhiMaps",nBinsC,binLimC,fPtBinning.GetSize()-1,fPtBinning.GetArray());
1087 fProfileCosPhiMaps->Sumw2();
1088 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: profilecosphimaps");
1092 // fMonitorTrackCuts
1095 if(fMonitorTrackCuts) {
1096 // Debugging tracking steps
1097 const Int_t nDimTrStep=2;
1098 Int_t nBinTrStep[nDimTrStep] = {fPtBinning.GetSize()-1,nBinsStep};
1099 fTrackingCuts = new THnSparseF("TrackingCuts","TrackingCuts",nDimTrStep,nBinTrStep);
1100 fTrackingCuts->SetBinEdges(0,fPtBinning.GetArray());
1101 fTrackingCuts->SetBinEdges(1,binLimStep);
1102 fTrackingCuts->Sumw2();
1103 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: trackingcuts");
1107 // fMonitorContamination
1110 if(fMonitorContamination) {
1111 // Maps delta phi contamination
1112 const Int_t nDimgcont=4;
1113 Int_t nBingcont[nDimgcont] = {nBinsPhiLess,nBinsC,fPtBinning.GetSize()-1, nBinsTPCdEdx};
1114 fDeltaPhiMapsContamination = new THnSparseF("DeltaPhiMapsContamination","DeltaPhiMapsContamination",nDimgcont,nBingcont);
1115 fDeltaPhiMapsContamination->SetBinEdges(0,binLimPhiLess);
1116 fDeltaPhiMapsContamination->SetBinEdges(1,binLimC);
1117 fDeltaPhiMapsContamination->SetBinEdges(2,fPtBinning.GetArray());
1118 fDeltaPhiMapsContamination->SetBinEdges(3,binLimTPCdEdx);
1119 fDeltaPhiMapsContamination->Sumw2();
1120 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: deltaphimapscontamination");
1124 // fMonitorWithoutPID
1127 if(fMonitorWithoutPID) {
1129 const Int_t nDimgb=3;
1130 Int_t nBingb[nDimgb] = {nBinsPhi,nBinsC,fPtBinning.GetSize()-1};
1132 fDeltaPhiMapsBeforePID = new THnSparseF("DeltaPhiMapsBeforePID","DeltaPhiMapsBeforePID",nDimgb,nBingb);
1133 fDeltaPhiMapsBeforePID->SetBinEdges(0,binLimPhi);
1134 fDeltaPhiMapsBeforePID->SetBinEdges(1,binLimC);
1135 fDeltaPhiMapsBeforePID->SetBinEdges(2,fPtBinning.GetArray());
1136 fDeltaPhiMapsBeforePID->Sumw2();
1137 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: deltaphimapsbeforepid");
1139 const Int_t nDimhb=3;
1140 Int_t nBinhb[nDimhb] = {nBinsCos,nBinsC,fPtBinning.GetSize()-1};
1142 fCosPhiMapsBeforePID = new THnSparseF("CosPhiMapsBeforePID","CosPhiMapsBeforePID",nDimhb,nBinhb);
1143 fCosPhiMapsBeforePID->SetBinEdges(0,binLimCos);
1144 fCosPhiMapsBeforePID->SetBinEdges(1,binLimC);
1145 fCosPhiMapsBeforePID->SetBinEdges(2,fPtBinning.GetArray());
1146 fCosPhiMapsBeforePID->Sumw2();
1147 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: cosphimapsbeforepid");
1153 if(fMonitorPhotonic) {
1155 const Int_t nDimgbp=3;
1156 Int_t nBingbp[nDimgbp] = {nBinsPhi,nBinsC,fPtBinning.GetSize()-1};
1158 fDeltaPhiMapsTaggedPhotonic = new THnSparseF("DeltaPhiMapsTaggedPhotonic","DeltaPhiMapsTaggedPhotonic",nDimgbp,nBingbp);
1159 fDeltaPhiMapsTaggedPhotonic->SetBinEdges(0,binLimPhi);
1160 fDeltaPhiMapsTaggedPhotonic->SetBinEdges(1,binLimC);
1161 fDeltaPhiMapsTaggedPhotonic->SetBinEdges(2,fPtBinning.GetArray());
1162 fDeltaPhiMapsTaggedPhotonic->Sumw2();
1163 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: deltaphimapstaggedphotonic");
1165 fDeltaPhiMapsTaggedNonPhotonic = new THnSparseF("DeltaPhiMapsTaggedNonPhotonic","DeltaPhiMapsTaggedNonPhotonic",nDimgbp,nBingbp);
1166 fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(0,binLimPhi);
1167 fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(1,binLimC);
1168 fDeltaPhiMapsTaggedNonPhotonic->SetBinEdges(2,fPtBinning.GetArray());
1169 fDeltaPhiMapsTaggedNonPhotonic->Sumw2();
1170 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: deltaphimapstaggednonphotonic");
1172 fDeltaPhiMapsTaggedPhotonicLS = new THnSparseF("DeltaPhiMapsTaggedPhotonicLS","DeltaPhiMapsTaggedPhotonicLS",nDimgbp,nBingbp);
1173 fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(0,binLimPhi);
1174 fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(1,binLimC);
1175 fDeltaPhiMapsTaggedPhotonicLS->SetBinEdges(2,fPtBinning.GetArray());
1176 fDeltaPhiMapsTaggedPhotonicLS->Sumw2();
1177 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: deltaphimapstaggedphotonicls");
1179 const Int_t nDimMCSource=3;
1180 Int_t nBinMCSource[nDimMCSource] = {nBinsC,fPtBinning.GetSize()-1,nBinsSource};
1181 fMCSourceDeltaPhiMaps = new THnSparseF("MCSourceDeltaPhiMaps","MCSourceDeltaPhiMaps",nDimMCSource,nBinMCSource);
1182 fMCSourceDeltaPhiMaps->SetBinEdges(0,binLimC);
1183 fMCSourceDeltaPhiMaps->SetBinEdges(1,fPtBinning.GetArray());
1184 fMCSourceDeltaPhiMaps->SetBinEdges(2,binLimSource);
1185 fMCSourceDeltaPhiMaps->Sumw2();
1186 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: mcsourcedeltaphimaps");
1188 // Maps invmass opposite
1189 const Int_t nDimOppSign=5;
1190 Int_t nBinOppSign[nDimOppSign] = {nBinsPhi,nBinsC,fPtBinning.GetSize()-1,nBinsInvMass,nBinsSource};
1191 fOppSignDeltaPhiMaps = new THnSparseF("OppSignDeltaPhiMaps","OppSignDeltaPhiMaps",nDimOppSign,nBinOppSign);
1192 fOppSignDeltaPhiMaps->SetBinEdges(0,binLimPhi);
1193 fOppSignDeltaPhiMaps->SetBinEdges(1,binLimC);
1194 fOppSignDeltaPhiMaps->SetBinEdges(2,fPtBinning.GetArray());
1195 fOppSignDeltaPhiMaps->SetBinEdges(3,binLimInvMass);
1196 fOppSignDeltaPhiMaps->SetBinEdges(4,binLimSource);
1197 fOppSignDeltaPhiMaps->Sumw2();
1198 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: oppsigndeltaphimaps");
1200 // Maps invmass same sign
1201 const Int_t nDimSameSign=5;
1202 Int_t nBinSameSign[nDimSameSign] = {nBinsPhi,nBinsC,fPtBinning.GetSize()-1,nBinsInvMass,nBinsSource};
1203 fSameSignDeltaPhiMaps = new THnSparseF("SameSignDeltaPhiMaps","SameSignDeltaPhiMaps",nDimSameSign,nBinSameSign);
1204 fSameSignDeltaPhiMaps->SetBinEdges(0,binLimPhi);
1205 fSameSignDeltaPhiMaps->SetBinEdges(1,binLimC);
1206 fSameSignDeltaPhiMaps->SetBinEdges(2,fPtBinning.GetArray());
1207 fSameSignDeltaPhiMaps->SetBinEdges(3,binLimInvMass);
1208 fSameSignDeltaPhiMaps->SetBinEdges(4,binLimSource);
1209 fSameSignDeltaPhiMaps->Sumw2();
1210 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: samesigndeltaphimaps");
1212 // Maps angle same sign
1213 const Int_t nDimAngleSameSign=3;
1214 Int_t nBinAngleSameSign[nDimAngleSameSign] = {nBinsAngle,nBinsC,nBinsSource};
1215 fSameSignAngle = new THnSparseF("SameSignAngleMaps","SameSignAngleMaps",nDimAngleSameSign,nBinAngleSameSign);
1216 fSameSignAngle->SetBinEdges(0,binLimAngle);
1217 fSameSignAngle->SetBinEdges(1,binLimC);
1218 fSameSignAngle->SetBinEdges(2,binLimSource);
1219 fSameSignAngle->Sumw2();
1220 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: samesignangle");
1222 // Maps angle opp sign
1223 const Int_t nDimAngleOppSign=3;
1224 Int_t nBinAngleOppSign[nDimAngleOppSign] = {nBinsAngle,nBinsC,nBinsSource};
1225 fOppSignAngle = new THnSparseF("OppSignAngleMaps","OppSignAngleMaps",nDimAngleOppSign,nBinAngleOppSign);
1226 fOppSignAngle->SetBinEdges(0,binLimAngle);
1227 fOppSignAngle->SetBinEdges(1,binLimC);
1228 fOppSignAngle->SetBinEdges(2,binLimSource);
1229 fOppSignAngle->Sumw2();
1230 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: oppsignangle");
1234 //**************************
1236 //******************************
1238 fListHist->Add(fHistEV);
1239 fListHist->Add(fHistPileUp);
1240 fListHist->Add(fEventPlane);
1241 fListHist->Add(fFractionContamination);
1242 fListHist->Add(fCosRes);
1243 fListHist->Add(fCosResabc);
1244 fListHist->Add(fCosPhiMaps);
1245 fListHist->Add(fDeltaPhiMaps);
1246 fListHist->Add(fPIDqa->MakeList("HFEpidQA"));
1247 fListHist->Add(fContaminationv2);
1248 fListHist->Add(fContaminationmeanpt);
1249 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add default");
1251 if(fMonitorEventPlane) {
1252 fListHist->Add(fProfileCosRes);
1253 fListHist->Add(fProfileCosResab);
1254 fListHist->Add(fProfileCosResac);
1255 fListHist->Add(fProfileCosResbc);
1256 fListHist->Add(fCosSin2phiep);
1257 fListHist->Add(fCos2phie);
1258 fListHist->Add(fSin2phie);
1259 fListHist->Add(fCos2phiep);
1260 fListHist->Add(fSin2phiep);
1261 fListHist->Add(fSin2phiephiep);
1262 fListHist->Add(fSinRes);
1263 fListHist->Add(fSinResabc);
1264 fListHist->Add(fProfileCosPhiMaps);
1266 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add monitor");
1268 if(fMonitorTrackCuts) fListHist->Add(fTrackingCuts);
1270 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add monitortrackcuts");
1272 if(fMonitorContamination) {
1273 fListHist->Add(fDeltaPhiMapsContamination);
1276 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add deltaphimapscontamination");
1278 if(fMonitorWithoutPID) {
1279 fListHist->Add(fDeltaPhiMapsBeforePID);
1280 fListHist->Add(fCosPhiMapsBeforePID);
1283 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add without pid");
1285 if(fMonitorPhotonic) {
1286 fListHist->Add(fPIDBackgroundqa->MakeList("HFEpidBackgroundQA"));
1287 fListHist->Add(fDeltaPhiMapsTaggedPhotonic);
1288 fListHist->Add(fDeltaPhiMapsTaggedNonPhotonic);
1289 fListHist->Add(fDeltaPhiMapsTaggedPhotonicLS);
1290 fListHist->Add(fMCSourceDeltaPhiMaps);
1291 fListHist->Add(fOppSignDeltaPhiMaps);
1292 fListHist->Add(fSameSignDeltaPhiMaps);
1293 fListHist->Add(fSameSignAngle);
1294 fListHist->Add(fOppSignAngle);
1295 fListHist->Add(fBackgroundSubtraction->GetListOutput());
1298 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add photonic");
1300 if(fHFEVZEROEventPlane && fMonitorEventPlane) fListHist->Add(fHFEVZEROEventPlane->GetOutputList());
1302 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: add event plane");
1306 PostData(1, fListHist);
1307 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {
1308 // PostData(bincless+2,fflowEvent);
1311 AliDebug(2,"AliAnalysisTaskFlowTPCTOFEPSP: post");
1316 //________________________________________________________________________
1317 void AliAnalysisTaskFlowTPCTOFEPSP::UserExec(Option_t */*option*/)
1324 Double_t massElectron = 0.000511;
1325 Double_t mcReactionPlane = 0.0;
1328 Double_t binct = 11.5;
1329 Double_t binctMore = 20.5;
1330 Double_t binctLess = -0.5;
1331 Float_t binctt = -1.0;
1333 Double_t valuecossinephiep[3];
1334 Double_t valuensparsea[4];
1335 Double_t valuensparseabis[5];
1336 Double_t valuensparsee[4];
1337 Double_t valuensparsef[2];
1338 Double_t valuensparsefsin[2];
1339 Double_t valuensparsefbis[4];
1340 Double_t valuensparsefbissin[4];
1341 Double_t valuensparseg[5];
1342 Double_t valuensparseh[5];
1343 Double_t valuensparsehprofile[3];
1344 Double_t valuensparseMCSourceDeltaPhiMaps[3];
1345 Double_t valuetrackingcuts[2];
1346 Double_t valuedeltaphicontamination[4];
1347 Double_t valuefractioncont[2];
1349 AliMCEvent *mcEvent = MCEvent();
1350 AliMCParticle *mctrack = NULL;
1353 Bool_t mcthere = kTRUE;
1355 AliAODEvent *aodE = dynamic_cast<AliAODEvent *>(fInputEvent);
1357 // printf("testd\n");
1358 AliError("No AOD Event");
1361 fAODMCHeader = dynamic_cast<AliAODMCHeader *>(fInputEvent->FindListObject(AliAODMCHeader::StdBranchName()));
1365 fAODArrayMCInfo = dynamic_cast<TClonesArray *>(fInputEvent->FindListObject(AliAODMCParticle::StdBranchName()));
1366 if(!fAODArrayMCInfo){
1370 fHFECuts->SetMCEvent(aodE);
1371 if(fMonitorPhotonic) fBackgroundSubtraction->SetAODArrayMCInfo(fAODArrayMCInfo);
1375 AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1376 if(!mcH) mcthere = kFALSE;
1378 if(fMonitorPhotonic) fBackgroundSubtraction->SetMCEvent(fMCEvent);
1387 //AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
1389 AliCentrality *centrality = fInputEvent->GetCentrality();
1390 //AliDebug(2,"Got the centrality");
1391 if(!centrality) return;
1392 cntr = centrality->GetCentralityPercentile("V0M");
1393 if((0.0< cntr) && (cntr<5.0)) binct = 0.5;
1394 if((5.0< cntr) && (cntr<10.0)) binct = 1.5;
1395 if((10.0< cntr) && (cntr<20.0)) binct = 2.5;
1396 if((20.0< cntr) && (cntr<30.0)) binct = 3.5;
1397 if((30.0< cntr) && (cntr<40.0)) binct = 4.5;
1398 if((40.0< cntr) && (cntr<50.0)) binct = 5.5;
1399 if((50.0< cntr) && (cntr<60.0)) binct = 6.5;
1400 if((60.0< cntr) && (cntr<70.0)) binct = 7.5;
1401 if((70.0< cntr) && (cntr<80.0)) binct = 8.5;
1402 if((80.0< cntr) && (cntr<90.0)) binct = 9.5;
1403 if((90.0< cntr) && (cntr<100.0)) binct = 10.5;
1405 if((0.< cntr) && (cntr < 20.)) binctt = 0.5;
1406 if((20.< cntr) && (cntr < 40.)) binctt = 1.5;
1407 if((40.< cntr) && (cntr < 80.)) binctt = 2.5;
1409 if((0.0< cntr) && (cntr<5.0)) binctMore = 0.5;
1410 if((5.0< cntr) && (cntr<10.0)) binctMore = 1.5;
1411 if((10.0< cntr) && (cntr<15.0)) binctMore = 2.5;
1412 if((15.0< cntr) && (cntr<20.0)) binctMore = 3.5;
1413 if((20.0< cntr) && (cntr<25.0)) binctMore = 4.5;
1414 if((25.0< cntr) && (cntr<30.0)) binctMore = 5.5;
1415 if((30.0< cntr) && (cntr<35.0)) binctMore = 6.5;
1416 if((35.0< cntr) && (cntr<40.0)) binctMore = 7.5;
1417 if((40.0< cntr) && (cntr<45.0)) binctMore = 8.5;
1418 if((45.0< cntr) && (cntr<50.0)) binctMore = 9.5;
1419 if((50.0< cntr) && (cntr<55.0)) binctMore = 10.5;
1420 if((55.0< cntr) && (cntr<60.0)) binctMore = 11.5;
1421 if((60.0< cntr) && (cntr<65.0)) binctMore = 12.5;
1422 if((65.0< cntr) && (cntr<70.0)) binctMore = 13.5;
1423 if((70.0< cntr) && (cntr<75.0)) binctMore = 14.5;
1424 if((75.0< cntr) && (cntr<80.0)) binctMore = 15.5;
1425 if((80.0< cntr) && (cntr<85.0)) binctMore = 16.5;
1426 if((85.0< cntr) && (cntr<90.0)) binctMore = 17.5;
1427 if((90.0< cntr) && (cntr<95.0)) binctMore = 18.5;
1428 if((95.0< cntr) && (cntr<100.0)) binctMore = 19.5;
1433 if(binct > 11.0) return;
1436 valuensparsea[3] = binct;
1437 valuensparseabis[1] = binct;
1438 valuensparsee[1] = binct;
1439 valuensparsef[1] = binctMore;
1440 valuensparsefsin[1] = binct;
1441 valuensparsefbis[3] = binctMore;
1442 valuensparsefbissin[3] = binct;
1443 valuensparseg[1] = binct;
1444 valuensparseh[1] = binct;
1445 valuefractioncont[1] = binct;
1446 valuensparsehprofile[1] = binct;
1447 valuecossinephiep[2] = binctMore;
1448 valuensparseMCSourceDeltaPhiMaps[0] = binct;
1449 valuedeltaphicontamination[1] = binct;
1451 //////////////////////
1453 //////////////////////
1455 Int_t runnumber = fInputEvent->GetRunNumber();
1456 AliDebug(2,Form("Run number %d",runnumber));
1458 if(!fPID->IsInitialized()){
1459 // Initialize PID with the given run number
1460 fPID->InitializePID(runnumber);
1462 if(!fPIDTOFOnly->IsInitialized()){
1463 // Initialize PID with the given run number
1464 fPIDTOFOnly->InitializePID(runnumber);
1468 if(!fPIDBackground->IsInitialized()){
1469 // Initialize PID with the given run number
1470 fPIDBackground->InitializePID(runnumber);
1473 fHFECuts->SetRecEvent(fInputEvent);
1474 if(mcEvent) fHFECuts->SetMCEvent(mcEvent);
1481 AliPIDResponse *pidResponse = fInputHandler->GetPIDResponse();
1483 AliDebug(2,"No PID response set");
1486 fPID->SetPIDResponse(pidResponse);
1487 fPIDTOFOnly->SetPIDResponse(pidResponse);
1488 fPIDBackground->SetPIDResponse(pidResponse);
1489 if(fMonitorPhotonic) fBackgroundSubtraction->InitRun(fInputEvent,pidResponse);
1491 fHistEV->Fill(binctt,0.0);
1496 if(!fHFECuts->CheckEventCuts("fEvRecCuts", fInputEvent)) {
1497 AliDebug(2,"Does not pass the event cut");
1498 PostData(1, fListHist);
1502 fHistEV->Fill(binctt,1.0);
1505 ///////////////////////////////////////////////////////////
1507 ///////////////////////////////////////////////////////////
1509 Float_t multTPC(0.); // tpc mult estimate
1510 Float_t multGlob(0.); // global multiplicity
1511 const Int_t nGoodTracks = fInputEvent->GetNumberOfTracks();
1512 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill tpc mult
1513 AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*>(fInputEvent->GetTrack(iTracks));
1514 if (!trackAOD) continue;
1515 if (!(trackAOD->TestFilterBit(1))) continue;
1516 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;
1519 for(Int_t iTracks = 0; iTracks < nGoodTracks; iTracks++) { // fill global mult
1520 AliAODTrack* trackAOD = dynamic_cast<AliAODTrack*>(fInputEvent->GetTrack(iTracks));
1521 if (!trackAOD) continue;
1522 if (!(trackAOD->TestFilterBit(16))) continue;
1523 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;
1524 Double_t b[2] = {-99., -99.};
1525 Double_t bCov[3] = {-99., -99., -99.};
1526 if (!(trackAOD->PropagateToDCA(fInputEvent->GetPrimaryVertex(), fInputEvent->GetMagneticField(), 100., b, bCov))) continue;
1527 if ((TMath::Abs(b[0]) > 0.3) || (TMath::Abs(b[1]) > 0.3)) continue;
1532 pileup[0]=fInputEvent->GetCentrality()->GetCentralityPercentile("V0M");
1533 pileup[1]=fInputEvent->GetCentrality()->GetCentralityPercentile("TRK");
1536 fHistPileUp->Fill(pileup);
1539 if (TMath::Abs(pileup[0]-pileup[1]) > 5) {
1540 AliDebug(2,"Does not pass the centrality correlation cut");
1543 if(multTPC < (-36.81+1.48*multGlob) && multTPC > (63.03+1.78*multGlob)){
1544 AliDebug(2,"Does not pass the multiplicity correlation cut");
1549 // AliVVZERO* vzeroData=fInputEvent->GetVZEROData();
1550 // Double_t mult[3],multV0A(0),multV0C(0);
1551 // for(Int_t i=0; i<32; ++i) {
1552 // multV0A += vzeroData->GetMultiplicityV0A(i);
1553 // multV0C += vzeroData->GetMultiplicityV0C(i);
1557 // for(Int_t k = 0; k < fInputEvent->GetNumberOfTracks(); k++){
1558 // AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);
1559 // if(!track) continue;
1560 // if(!(track->GetStatus()&AliVTrack::kITSrefit)) continue;
1561 // if(!(track->GetStatus()&AliVTrack::kTPCrefit)) continue;
1565 // mult[0]=fInputEvent->GetNumberOfTracks();
1566 // mult[1]=multV0A+multV0C;
1567 // mult[2]=binctMore;
1568 // fHistPileUp->Fill(mult);
1570 // if(fUpperPileUpCut&&fLowerPileUpCut){
1571 // if((mult[0]<fLowerPileUpCut->Eval(mult[1])) ||
1572 // (mult[0]>fUpperPileUpCut->Eval(mult[1]))){
1573 // AliDebug(2,"Does not pass the pileup cut");
1574 // PostData(1, fListHist);
1579 ////////////////////////////////////
1580 // First method event plane
1581 ////////////////////////////////////
1583 AliEventplane* vEPa = fInputEvent->GetEventplane();
1584 Float_t eventPlanea = 0.0;
1585 Float_t eventPlaneTPC = 0.0;
1586 Float_t eventPlaneV0A = 0.0;
1587 Float_t eventPlaneV0C = 0.0;
1588 Float_t eventPlaneV0 = 0.0;
1589 TVector2 *qTPC = 0x0;
1590 TVector2 *qsub1a = 0x0;
1591 TVector2 *qsub2a = 0x0;
1592 TVector2 qV0A,qV0C,qV0,*qAna;
1596 if(fHFEVZEROEventPlane && (!fAODAnalysis)){
1598 //AliESDEvent *esd = dynamic_cast<AliESDEvent*>(InputEvent());
1601 fHFEVZEROEventPlane->ProcessEvent(fInputEvent);
1603 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0A()+100) < 0.0000001) eventPlaneV0A = -100.0;
1605 eventPlaneV0A = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0A());
1606 if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();
1609 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0C()+100) < 0.0000001) eventPlaneV0C = -100.0;
1611 eventPlaneV0C = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0C());
1612 if(eventPlaneV0C > TMath::Pi()) eventPlaneV0C = eventPlaneV0C - TMath::Pi();
1615 if(TMath::Abs(fHFEVZEROEventPlane->GetEventPlaneV0()+100) < 0.0000001) eventPlaneV0 = -100.0;
1617 eventPlaneV0 = TVector2::Phi_0_2pi(fHFEVZEROEventPlane->GetEventPlaneV0());
1618 if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();
1624 Double_t qVx, qVy; //TR: info
1625 eventPlaneV0 = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,10,2,qVx,qVy));
1626 if(eventPlaneV0 > TMath::Pi()) eventPlaneV0 = eventPlaneV0 - TMath::Pi();
1628 eventPlaneV0A = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,8,2,qVx,qVy));
1629 if(eventPlaneV0A > TMath::Pi()) eventPlaneV0A = eventPlaneV0A - TMath::Pi();
1631 eventPlaneV0C = TVector2::Phi_0_2pi(vEPa->CalculateVZEROEventPlane(fInputEvent,9,2,qVx,qVy));
1632 if(eventPlaneV0C > TMath::Pi()) eventPlaneV0C = eventPlaneV0C - TMath::Pi();
1635 if(eventPlaneV0<-900) return;
1636 if(eventPlaneV0A<-900) return;
1637 if(eventPlaneV0C<-900) return;
1640 eventPlaneV0=TVector2::Phi_0_2pi(eventPlaneV0);
1641 eventPlaneV0A=TVector2::Phi_0_2pi(eventPlaneV0A);
1642 eventPlaneV0C=TVector2::Phi_0_2pi(eventPlaneV0C);
1648 qTPC = vEPa->GetQVector();
1655 TVector2 qVectorfortrack;
1656 qVectorfortrack.Set(qx,qy);
1657 eventPlaneTPC = TVector2::Phi_0_2pi(qVectorfortrack.Phi())/2.;
1659 // Choose the one used for v2
1661 if(fVZEROEventPlane){ //TR: info
1662 eventPlanea = eventPlaneV0;
1665 if(fVZEROEventPlaneA){
1666 eventPlanea = eventPlaneV0A;
1669 if(fVZEROEventPlaneC){
1670 eventPlanea = eventPlaneV0C;
1673 if(!fVZEROEventPlane){
1674 eventPlanea = eventPlaneTPC;
1678 valuecossinephiep[0] = TMath::Cos(2*eventPlanea);
1679 valuecossinephiep[1] = TMath::Sin(2*eventPlanea);
1681 Float_t eventPlanesub1a = -100.0;
1682 Float_t eventPlanesub2a = -100.0;
1683 Double_t diffsub1sub2a = -100.0;
1684 Double_t diffsub1sub2asin = -100.0;
1685 Double_t diffsubasubb = -100.0;
1686 Double_t diffsubasubc = -100.0;
1687 Double_t diffsubbsubc = -100.0;
1688 Double_t diffsubasubbsin = -100.0;
1689 Double_t diffsubasubcsin = -100.0;
1690 Double_t diffsubbsubcsin = -100.0;
1692 // two sub event TPC
1693 qsub1a = vEPa->GetQsub1();
1694 qsub2a = vEPa->GetQsub2();
1696 /////////////////////////////////////////////////////////
1697 // Cut for event with event plane reconstructed by all
1698 ////////////////////////////////////////////////////////
1700 if((!qTPC) || (!qsub1a) || (!qsub2a)) {
1701 AliDebug(2,"No event plane");
1705 eventPlanesub1a = TVector2::Phi_0_2pi(qsub1a->Phi())/2.;
1706 eventPlanesub2a = TVector2::Phi_0_2pi(qsub2a->Phi())/2.;
1707 diffsub1sub2a = TMath::Cos(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));
1708 diffsub1sub2asin = TMath::Sin(2.*TVector2::Phi_0_2pi(qsub1a->Phi()/2.- qsub2a->Phi()/2.));
1711 // if ( !fDebugStreamer ) {
1713 // TDirectory *backup = gDirectory;
1714 // fDebugStreamer = new TTreeSRedirector("TaskFlowTPCTOFEPSPdebug.root");
1715 // if ( backup ) backup->cd(); //we don't want to be cd'd to the debug streamer
1720 // double v0nrom = TMath::Sqrt(qV0.X()*qV0.X()+qV0.Y()*qV0.Y());
1721 // double v0Anrom = TMath::Sqrt(qV0A.X()*qV0A.X()+qV0A.Y()*qV0A.Y());
1722 // double v0Cnrom = TMath::Sqrt(qV0C.X()*qV0C.X()+qV0C.Y()*qV0C.Y());
1723 // double sub1nrom = TMath::Sqrt(qsub1a->X()*qsub1a->X()+qsub1a->Y()*qsub1a->Y());
1724 // double sub2nrom = TMath::Sqrt(qsub2a->X()*qsub2a->X()+qsub2a->Y()*qsub2a->Y());
1726 // (* fDebugStreamer) << "UserExec" <<
1727 // "binct="<<binct<<
1729 // "qV0A="<<v0Anrom<<
1730 // "qV0C="<<v0Cnrom<<
1731 // "qsub1a="<<sub1nrom<<
1732 // "qsub2a="<<sub2nrom<<
1736 // three sub events in case of VZEROA and VZEROC
1738 diffsubasubb = TMath::Cos(2.*(eventPlaneV0A - eventPlaneV0C)); //TR:
1739 diffsubasubc = TMath::Cos(2.*(eventPlaneV0A - eventPlaneTPC)); //TR:
1740 diffsubbsubc = TMath::Cos(2.*(eventPlaneV0C - eventPlaneTPC)); //TR:
1743 if(fVZEROEventPlaneA){
1744 diffsubasubb = qV0A.X()*qV0C.X()+qV0A.Y()*qV0C.Y();
1745 diffsubasubc = qV0A.X()*qTPC->X()+qV0A.Y()*qTPC->Y();
1746 diffsubbsubc = qV0C.X()*qTPC->X()+qV0C.Y()*qTPC->Y();
1748 else if(fVZEROEventPlaneC){
1749 diffsubasubb = qV0C.X()*qV0A.X()+qV0C.Y()*qV0A.Y();
1750 diffsubasubc = qV0C.X()*qTPC->X()+qV0C.Y()*qTPC->Y();
1751 diffsubbsubc = qV0A.X()*qTPC->X()+qV0A.Y()*qTPC->Y();
1755 diffsubasubbsin = TMath::Sin(2.*(eventPlaneV0A - eventPlaneV0C));
1756 diffsubasubcsin = TMath::Sin(2.*(eventPlaneV0A - eventPlaneTPC));
1757 diffsubbsubcsin = TMath::Sin(2.*(eventPlaneV0C - eventPlaneTPC));
1758 // three sub events in case of VZERO all
1759 if(fVZEROEventPlane && (!fVZEROEventPlaneA) && (!fVZEROEventPlaneC)) {
1761 diffsubasubb = TMath::Cos(2.*(eventPlaneV0 - eventPlanesub1a)); //TR:
1762 diffsubasubc = TMath::Cos(2.*(eventPlaneV0 - eventPlanesub2a)); //TR:
1763 diffsubbsubc = TMath::Cos(2.*(eventPlanesub1a - eventPlanesub2a)); //TR:
1766 diffsubasubb = qV0.X()*qsub1a->X()+qV0.Y()*qsub1a->Y();
1767 diffsubasubc = qV0.X()*qsub2a->X()+qV0.Y()*qsub2a->Y();
1768 diffsubbsubc = qsub1a->X()*qsub2a->X()+qsub1a->Y()*qsub2a->Y();
1771 diffsubasubbsin = TMath::Sin(2.*(eventPlaneV0 - eventPlanesub1a));
1772 diffsubasubcsin = TMath::Sin(2.*(eventPlaneV0 - eventPlanesub2a));
1773 diffsubbsubcsin = TMath::Sin(2.*(eventPlanesub1a - eventPlanesub2a));
1776 //////////////////////////////////////
1777 // AliFlowEvent and MC event plane
1778 /////////////////////////////////////
1780 Int_t nbtracks = fInputEvent->GetNumberOfTracks();
1781 AliDebug(2,Form("Number of tracks %d",nbtracks));
1783 if(fMonitorQCumulant) {
1785 fcutsRP->SetEvent( InputEvent(), MCEvent());
1786 fcutsPOI->SetEvent( InputEvent(), MCEvent());
1788 fflowEvent->~AliFlowEvent();
1789 new(fflowEvent) AliFlowEvent(fcutsRP,fcutsPOI);
1790 }else fflowEvent = new AliFlowEvent(fcutsRP,fcutsPOI);
1791 if(mcEvent && mcEvent->GenEventHeader()) {
1792 fflowEvent->SetMCReactionPlaneAngle(mcEvent);
1793 //if reaction plane not set from elsewhere randomize it before adding flow
1794 //if (!fflowEvent->IsSetMCReactionPlaneAngle()) fflowEvent->SetMCReactionPlaneAngle(gRandom->Uniform(0.0,TMath::TwoPi()));
1795 mcReactionPlane = TVector2::Phi_0_2pi(fflowEvent->GetMCReactionPlaneAngle());
1796 if(mcReactionPlane > TMath::Pi()) mcReactionPlane = mcReactionPlane - TMath::Pi();
1797 AliDebug(2,Form("MC reaction plane %f",mcReactionPlane));
1799 fflowEvent->SetReferenceMultiplicity( nbtracks );
1800 fflowEvent->DefineDeadZone(0,0,0,0);
1801 //fflowEvent.TagSubeventsInEta(-0.8,-0.1,0.1,0.8);
1806 if(fUseMCReactionPlane) {
1807 eventPlanea = mcReactionPlane;
1808 diffsub1sub2a = 0.0;
1813 //////////////////////
1815 //////////////////////
1817 fHistEV->Fill(binctt,2.0);
1820 valuensparsea[0] = eventPlaneV0A;
1821 valuensparsea[1] = eventPlaneV0C;
1822 valuensparsea[2] = eventPlaneTPC;
1823 if(fVZEROEventPlane && (!fVZEROEventPlaneA) && (!fVZEROEventPlaneC)) {
1825 valuensparsea[0] = eventPlaneV0;
1826 valuensparsea[1] = eventPlanesub1a;
1827 valuensparsea[2] = eventPlanesub2a;
1829 fEventPlane->Fill(&valuensparsea[0]);
1832 if(fMonitorEventPlane) fCosSin2phiep->Fill(&valuecossinephiep[0]);
1834 if(!fVZEROEventPlane) {
1835 valuensparsef[0] = diffsub1sub2a;
1836 fCosRes->Fill(&valuensparsef[0]);
1837 valuensparsefsin[0] = diffsub1sub2asin;
1838 if(fMonitorEventPlane) fSinRes->Fill(&valuensparsefsin[0]);
1839 if(fMonitorEventPlane) {
1840 fProfileCosRes->Fill(valuensparsef[1],valuensparsef[0]);
1844 valuensparsefbis[0] = diffsubasubb;
1845 valuensparsefbis[1] = diffsubasubc;
1846 valuensparsefbis[2] = diffsubbsubc;
1847 fCosResabc->Fill(&valuensparsefbis[0]); //TR: info
1848 valuensparsefbissin[0] = diffsubasubbsin;
1849 valuensparsefbissin[1] = diffsubbsubcsin;
1850 valuensparsefbissin[2] = diffsubasubcsin;
1851 if(fMonitorEventPlane) fSinResabc->Fill(&valuensparsefbissin[0]);
1852 if(fMonitorEventPlane) {
1853 fProfileCosResab->Fill(valuensparsefbis[3],valuensparsefbis[0]);
1854 fProfileCosResac->Fill(valuensparsefbis[3],valuensparsefbis[1]);
1855 fProfileCosResbc->Fill(valuensparsefbis[3],valuensparsefbis[2]);
1859 ////////////////////////////////////////
1860 // Loop to determine pool background
1861 /////////////////////////////////////////
1862 if(fMonitorPhotonic) {
1864 fBackgroundSubtraction->FillPoolAssociatedTracks(fInputEvent,binct);
1867 fArraytrack->~TArrayI();
1868 new(fArraytrack) TArrayI(nbtracks);
1871 fArraytrack = new TArrayI(nbtracks);
1873 fCounterPoolBackground = 0;
1875 for(Int_t k = 0; k < nbtracks; k++){
1877 AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);
1878 if(!track) continue;
1881 Bool_t survivedbackground = kTRUE;
1883 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
1885 AliESDtrack esdTrack(aodtrack);
1886 // set the TPC cluster info
1887 esdTrack.SetTPCClusterMap(aodtrack->GetTPCClusterMap());
1888 esdTrack.SetTPCSharedMap(aodtrack->GetTPCSharedMap());
1889 esdTrack.SetTPCPointsF(aodtrack->GetTPCNclsF());
1890 // needed to calculate the impact parameters
1891 AliAODEvent *aodeventu = dynamic_cast<AliAODEvent *>(fInputEvent);
1893 AliAODVertex *vAOD = aodeventu->GetPrimaryVertex();
1894 Double_t bfield = aodeventu->GetMagneticField();
1895 Double_t pos[3],cov[6];
1897 vAOD->GetCovarianceMatrix(cov);
1898 const AliESDVertex vESD(pos,cov,100.,100);
1899 esdTrack.RelateToVertex(&vESD,bfield,3.);
1901 if(!fHFEBackgroundCuts->IsSelected(&esdTrack)) {
1902 survivedbackground = kFALSE;
1907 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
1909 if(!fHFEBackgroundCuts->IsSelected(esdtrack)) survivedbackground = kFALSE;
1913 if(survivedbackground) {
1915 AliHFEpidObject hfetrack2;
1916 if(!fAODAnalysis) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);
1917 else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);
1918 hfetrack2.SetRecTrack(track);
1919 hfetrack2.SetCentrality((Int_t)binct);
1920 AliDebug(2,Form("centrality %f and %d",binct,hfetrack2.GetCentrality()));
1921 hfetrack2.SetPbPb();
1922 if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundqa)) {
1923 fArraytrack->AddAt(k,fCounterPoolBackground);
1924 fCounterPoolBackground++;
1925 AliDebug(2,Form("fCounterPoolBackground %d, track %d",fCounterPoolBackground,k));
1932 //////////////////////////
1934 //////////////////////////
1935 for(Int_t k = 0; k < nbtracks; k++){
1937 AliVTrack *track = (AliVTrack *) fInputEvent->GetTrack(k);
1938 if(!track) continue;
1941 AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
1943 AliError("AOD track is not there");
1946 AliDebug(2,"Find AOD track on");
1947 if(!(aodtrack->TestFilterBit(fFilter))) continue; // Only process AOD tracks where the HFE is set
1951 valuetrackingcuts[0] = track->Pt();
1952 valuetrackingcuts[1] = 0;
1954 // RecKine: ITSTPC cuts
1955 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecKineITSTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
1956 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
1959 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
1960 valuetrackingcuts[1] = 1;
1961 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
1963 // HFEcuts: ITS layers cuts
1964 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsITS + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
1965 valuetrackingcuts[1] = 2;
1966 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
1968 // HFE cuts: TOF and mismatch flag
1969 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTOF + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
1970 valuetrackingcuts[1] = 3;
1971 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
1973 // HFE cuts: TPC PID cleanup
1974 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTPC + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
1975 valuetrackingcuts[1] = 4;
1976 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
1978 // HFEcuts: Nb of tracklets TRD0
1979 if(!fHFECuts->CheckParticleCuts(AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack, (TObject *)track)) continue;
1980 valuetrackingcuts[1] = 5;
1981 if(fMonitorTrackCuts) fTrackingCuts->Fill(&valuetrackingcuts[0]);
1983 AliDebug(2,"Survived");
1985 /////////////////////////////////////////////////////////
1986 // Subtract candidate from TPC event plane
1987 ////////////////////////////////////////////////////////
1988 Float_t eventplanesubtracted = 0.0;
1990 if(!fVZEROEventPlane) {
1991 // Subtract the tracks from the event plane
1992 Double_t qX = qTPC->X() - vEPa->GetQContributionX(track); //Modify the components: subtract the track you want to look at with your analysis
1993 Double_t qY = qTPC->Y() - vEPa->GetQContributionY(track); //Modify the components: subtract the track you want to look at with your analysis
1994 TVector2 newQVectorfortrack;
1995 newQVectorfortrack.Set(qX,qY);
1996 eventplanesubtracted = TVector2::Phi_0_2pi(newQVectorfortrack.Phi())/2;
1998 else eventplanesubtracted = eventPlanea;
2000 ///////////////////////////////////////////
2002 //////////////////////////////////////////
2003 Bool_t fillEventPlane = kTRUE;
2004 if(!fVZEROEventPlane){
2005 if((!qsub1a) || (!qsub2a)) fillEventPlane = kFALSE;
2007 if(track->Eta() < (- fEtaGap/2.)) eventplanesubtracted = eventPlanesub1a;
2008 else if(track->Eta() > (fEtaGap/2.)) eventplanesubtracted = eventPlanesub2a;
2009 else fillEventPlane = kFALSE;
2016 if(fUseMCReactionPlane) {
2017 eventplanesubtracted = mcReactionPlane;
2018 fillEventPlane = kTRUE;
2021 //////////////////////////////////////////////////////////////////////////////
2022 ///////////////////////////AFTERBURNER
2023 Double_t phitrack = track->Phi();
2026 phitrack = GetPhiAfterAddV2(track->Phi(),mcReactionPlane);
2028 //////////////////////////////////////////////////////////////////////////////
2031 ///////////////////////
2032 // Calculate deltaphi
2033 ///////////////////////
2035 // Suppose phi track is between 0.0 and phi
2036 Double_t deltaphi = TVector2::Phi_0_2pi(phitrack - eventplanesubtracted);
2037 if(deltaphi > TMath::Pi()) deltaphi = deltaphi - TMath::Pi();
2039 ////////////////////////////////
2040 // Determine the deltaphi bin
2041 ///////////////////////////////
2044 if(((deltaphi<(TMath::Pi()/4.)) && (deltaphi>0.0)) || ((deltaphi>(3*TMath::Pi()/4.)) && (deltaphi<TMath::Pi()))) valuedeltaphicontamination[0] = 0.5;
2046 if((deltaphi>(TMath::Pi()/4.)) && (deltaphi<(3*TMath::Pi()/4.))) valuedeltaphicontamination[0] = 1.5;
2048 ////////////////////////////////////////
2050 ///////////////////////////////////////
2053 valuedeltaphicontamination[2] = track->Pt();
2054 valuensparsee[2] = track->Pt();
2055 valuensparsee[3] = track->Eta();
2056 valuensparseg[2] = track->Pt();
2057 valuensparseh[2] = track->Pt();
2058 valuefractioncont[0] = track->Pt();
2059 valuensparsehprofile[2] = track->Pt();
2060 valuensparseMCSourceDeltaPhiMaps[1] = track->Pt();
2061 if(track->Charge() > 0.0) {
2062 valuensparseg[3] = 0.2;
2063 valuensparseh[3] = 0.2;
2066 valuensparseg[3] = -0.2;
2067 valuensparseh[3] = -0.2;
2069 valuensparseh[4] = track->Eta();
2070 valuensparseg[4] = track->Eta();
2072 AliDebug(2,Form("charge %d",(Int_t)track->Charge()));
2074 ////////////////////////
2076 ///////////////////////
2078 if(fMonitorWithoutPID) {
2080 valuensparseg[0] = deltaphi;
2081 if(fillEventPlane) fDeltaPhiMapsBeforePID->Fill(&valuensparseg[0]);
2084 valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));
2085 if(fillEventPlane) {
2086 fCosPhiMapsBeforePID->Fill(&valuensparseh[0]);
2090 ////////////////////////
2092 ////////////////////////
2094 // Apply PID for Data
2097 AliHFEpidObject hfetrack;
2098 if(!fAODAnalysis) hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
2099 else hfetrack.SetAnalysisType(AliHFEpidObject::kAODanalysis);
2100 hfetrack.SetRecTrack(track);
2101 hfetrack.SetCentrality((Int_t)binct);
2102 AliDebug(2,Form("centrality %f and %d",binct,hfetrack.GetCentrality()));
2106 if(fMonitorContamination) {
2107 if(fPIDTOFOnly->IsSelected(&hfetrack,0x0,"recTrackCont",0x0)) {
2108 Float_t nsigma = pidResponse->NumberOfSigmasTPC(track, AliPID::kElectron);
2109 valuedeltaphicontamination[3] = nsigma;
2110 fDeltaPhiMapsContamination->Fill(&valuedeltaphicontamination[0]);
2114 // Complete PID TOF+TPC
2115 if(!fPID->IsSelected(&hfetrack,0x0,"recTrackCont",fPIDqa)) {
2120 if(!mcEvent) continue;
2121 if(!(mctrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(track->GetLabel()))))) continue;
2122 AliDebug(2,Form("PdgCode %d",TMath::Abs(mctrack->Particle()->GetPdgCode())));
2123 if(TMath::Abs(mctrack->Particle()->GetPdgCode())!=11) continue;
2128 /////////////////////////////////////////////////////////////////////////////
2129 // Add candidate to AliFlowEvent for POI and subtract from RP if needed
2130 ////////////////////////////////////////////////////////////////////////////
2131 if(fMonitorQCumulant) {
2132 Int_t idtrack = static_cast<AliVTrack*>(track)->GetID();
2133 Bool_t found = kFALSE;
2134 Int_t numberoffound = 0;
2135 AliDebug(2,Form("A: Number of tracks %d",fflowEvent->NumberOfTracks()));
2136 for(Int_t iRPs=0; iRPs< fflowEvent->NumberOfTracks(); iRPs++) {
2137 AliFlowTrack *iRP = (AliFlowTrack*) (fflowEvent->GetTrack(iRPs));
2139 //if(!iRP->InRPSelection()) continue;
2140 if( TMath::Abs(idtrack) == TMath::Abs(iRP->GetID()) ) {
2141 iRP->SetForPOISelection(kTRUE);
2146 AliDebug(2,Form("Found %d mal",numberoffound));
2148 AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*) MakeTrack(massElectron,track->Pt(),track->Phi(), track->Eta());
2149 sTrack->SetID(idtrack);
2150 fflowEvent->AddTrack(sTrack);
2151 AliDebug(2,"Add the track");
2153 AliDebug(2,Form("B: Number of tracks %d",fflowEvent->NumberOfTracks()));
2157 /////////////////////
2159 /////////////////////
2162 valuensparseabis[0] = eventplanesubtracted;
2163 if((fillEventPlane) && (fMonitorEventPlane)) fEventPlaneaftersubtraction->Fill(&valuensparseabis[0]);
2166 if(fMonitorEventPlane)
2169 valuensparsee[0] = TMath::Cos(2*phitrack);
2170 fCos2phie->Fill(&valuensparsee[0]);
2171 valuensparsee[0] = TMath::Sin(2*phitrack);
2172 fSin2phie->Fill(&valuensparsee[0]);
2174 valuensparsee[0] = TMath::Cos(2*eventplanesubtracted);
2175 if(fillEventPlane) fCos2phiep->Fill(&valuensparsee[0]);
2176 valuensparsee[0] = TMath::Sin(2*eventplanesubtracted);
2177 if(fillEventPlane) fSin2phiep->Fill(&valuensparsee[0]);
2178 valuensparsee[0] = TMath::Sin(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));
2179 if(fillEventPlane) fSin2phiephiep->Fill(&valuensparsee[0]);
2184 valuensparseg[0] = deltaphi;
2185 if(fillEventPlane) fDeltaPhiMaps->Fill(&valuensparseg[0]);
2188 valuensparseh[0] = TMath::Cos(2*TVector2::Phi_mpi_pi(phitrack-eventplanesubtracted));
2189 if(fillEventPlane) {
2190 fCosPhiMaps->Fill(&valuensparseh[0]); //TR: fCosPhiQSum+=valuensparseh[0]*TMath:Sqrt(qAna->X()*qAna->X()+qAna->Y()*qAna->Y()); fCosPhiQN++;
2191 if((valuefractioncont[1] >=0) && (valuefractioncont[1] < 11)){
2192 if(fContamination[((Int_t)valuefractioncont[1])]){
2193 Double_t weight = 1.;
2194 if(fAsFunctionOfP) weight = fContamination[((Int_t)valuefractioncont[1])]->Eval(track->P());
2195 else weight = fContamination[((Int_t)valuefractioncont[1])]->Eval(track->Pt());
2196 if(weight<0.0) weight=0.0;
2197 if(weight>1.0) weight=1.0;
2198 fFractionContamination->Fill(&valuefractioncont[0],weight);
2199 if(fv2contamination[((Int_t)valuefractioncont[1])]){
2200 Double_t v2 = fv2contamination[((Int_t)valuefractioncont[1])]->Eval(track->Pt());
2201 AliDebug(2,Form("value v2 %f, contamination %f and pt %f centrality %d\n",v2,weight,track->Pt(),(Int_t)valuefractioncont[1]));
2202 AliDebug(2,Form("Check for centrality 3: value v2 %f, contamination %f\n",fv2contamination[3]->Eval(track->Pt()),fContamination[3]->Eval(track->P())));
2203 AliDebug(2,Form("Check for centrality 4: value v2 %f, contamination %f\n",fv2contamination[4]->Eval(track->Pt()),fContamination[4]->Eval(track->P())));
2204 AliDebug(2,Form("Check for centrality 5: value v2 %f, contamination %f\n",fv2contamination[5]->Eval(track->Pt()),fContamination[5]->Eval(track->P())));
2205 fContaminationv2->Fill(valuefractioncont[1],valuefractioncont[0],v2,weight);
2207 fContaminationmeanpt->Fill(valuefractioncont[1],valuefractioncont[0],TMath::Abs(track->Pt()));
2210 if(fMonitorEventPlane) {
2212 valuensparseh[0] *= TMath::Sqrt(qAna->X()*qAna->X()+qAna->Y()*qAna->Y());
2213 fProfileCosPhiMaps->Fill(valuensparsehprofile[1],valuensparsehprofile[2],valuensparseh[0]); //TR: info
2217 if(fMonitorPhotonic) {
2218 Int_t indexmother = -1;
2220 if(mcthere) source = fBackgroundSubtraction->FindMother(mctrack->GetLabel(),indexmother);
2221 fBackgroundSubtraction->LookAtNonHFE(k, track, fInputEvent, 1, binct, deltaphi, source, indexmother);
2223 if((!fAODAnalysis && mcthere) || !mcthere) {
2227 source = FindMother(TMath::Abs(track->GetLabel()),mcEvent, indexmother);
2228 valuensparseMCSourceDeltaPhiMaps[2] = source;
2229 if(mcEvent) fMCSourceDeltaPhiMaps->Fill(&valuensparseMCSourceDeltaPhiMaps[0]);
2230 //LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);
2231 Int_t taggedvalue = LookAtNonHFE(k,track,fInputEvent,mcEvent,binct,deltaphi,source,indexmother);
2232 if(fMonitorPhotonic) {
2233 // No opposite charge partner found in the invariant mass choosen
2234 if((taggedvalue!=2) && (taggedvalue!=6)) {
2235 //fDeltaPhiMapsTaggedNonPhotonic->Fill(&valuensparseg[0]);
2236 //fCosPhiMapsTaggedNonPhotonic->Fill(&valuensparseh[0]);
2238 // One opposite charge partner found in the invariant mass choosen
2239 if((taggedvalue==2) || (taggedvalue==6)) {
2240 fDeltaPhiMapsTaggedPhotonic->Fill(&valuensparseg[0]);
2241 //fCosPhiMapsTaggedPhotonic->Fill(&valuensparseh[0]);
2243 // One same charge partner found in the invariant mass choosen
2244 if((taggedvalue==4) || (taggedvalue==6)) {
2245 fDeltaPhiMapsTaggedPhotonicLS->Fill(&valuensparseg[0]);
2246 //fCosPhiMapsTaggedPhotonicLS->Fill(&valuensparseh[0]);
2254 //////////////////////////////////////////////////////////////////////////////
2255 ///////////////////////////AFTERBURNER
2256 if (fAfterBurnerOn & fMonitorQCumulant)
2258 fflowEvent->AddFlow(fV1,fV2,fV3,fV4,fV5); //add flow
2259 fflowEvent->CloneTracks(fNonFlowNumberOfTrackClones); //add nonflow by cloning tracks
2261 //////////////////////////////////////////////////////////////////////////////
2265 //for(Int_t bincless = 0; bincless < fNbBinsCentralityQCumulant; bincless++) {
2266 // if((fBinCentralityLess[bincless]< cntr) && (cntr < fBinCentralityLess[bincless+1])) PostData(bincless+2,fflowEvent);
2269 if(fMonitorPhotonic) {
2276 if(fMonitorPhotonic) fBackgroundSubtraction->CountPoolAssociated(fInputEvent,binct);
2278 PostData(1, fListHist);
2283 //______________________________________________________________________________
2284 AliFlowCandidateTrack *AliAnalysisTaskFlowTPCTOFEPSP::MakeTrack( Double_t mass,
2285 Double_t pt, Double_t phi, Double_t eta) {
2287 // Make Track (Not needed actually)
2290 AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
2291 sTrack->SetMass(mass);
2293 sTrack->SetPhi(phi);
2294 sTrack->SetEta(eta);
2295 sTrack->SetForPOISelection(kTRUE);
2296 sTrack->SetForRPSelection(kFALSE);
2299 //_________________________________________________________________________________
2300 Double_t AliAnalysisTaskFlowTPCTOFEPSP::GetPhiAfterAddV2(Double_t phi,Double_t reactionPlaneAngle) const
2303 // Adds v2, uses Newton-Raphson iteration
2305 Double_t phiend=phi;
2309 Double_t phiprev=0.;
2311 for (Int_t i=0; i<fMaxNumberOfIterations; i++)
2313 phiprev=phiend; //store last value for comparison
2314 f = phiend-phi0+fV2*TMath::Sin(2.*(phiend-reactionPlaneAngle));
2315 fp = 1.0+2.0*fV2*TMath::Cos(2.*(phiend-reactionPlaneAngle)); //first derivative
2317 if (TMath::AreEqualAbs(phiprev,phiend,fPrecisionPhi)) break;
2321 //_____________________________________________________________________________________________
2322 Int_t AliAnalysisTaskFlowTPCTOFEPSP::LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *vEvent, AliMCEvent *mcEvent,Int_t binct,Double_t deltaphi,Int_t source,Int_t indexmother)
2328 // return -1 if nothing
2329 // return 2 if opposite charge within the mass range found
2330 // return 4 if like charge within the mass range found
2331 // return 6 if opposite charge and like charge within the mass range found
2334 Int_t taggedphotonic = -1;
2336 Bool_t oppositetaggedphotonic = kFALSE;
2337 Bool_t sametaggedphotonic = kFALSE;
2339 AliDebug(2,Form("fCounterPoolBackground %d in LookAtNonHFE!!!",fCounterPoolBackground));
2340 if(!fArraytrack) return taggedphotonic;
2341 AliDebug(2,Form("process track %d",iTrack1));
2346 Double_t valuensparseDeltaPhiMaps[5];
2347 Double_t valueangle[3];
2349 valuensparseDeltaPhiMaps[1] = binct;
2350 valuensparseDeltaPhiMaps[2] = track1->Pt();
2351 valuensparseDeltaPhiMaps[0] = deltaphi;
2352 valuensparseDeltaPhiMaps[4] = source;
2354 valueangle[2] = source;
2355 valueangle[1] = binct;
2358 Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()),mcEvent);
2359 Int_t numberfound = 0;
2362 Double_t bfield = vEvent->GetMagneticField();
2364 // Get Primary vertex
2365 const AliVVertex *pVtx = vEvent->GetPrimaryVertex();
2367 for(Int_t idex = 0; idex < fCounterPoolBackground; idex++)
2370 Int_t iTrack2 = fArraytrack->At(idex);
2371 AliDebug(2,Form("track %d",iTrack2));
2372 AliVTrack* track2 = (AliVTrack *) vEvent->GetTrack(iTrack2);
2375 printf("ERROR: Could not receive track %d", iTrack2);
2378 if(iTrack2==iTrack1) continue;
2379 AliDebug(2,"Different");
2381 // Reset the MC info
2382 valueangle[2] = source;
2383 valuensparseDeltaPhiMaps[4] = source;
2385 // track cuts and PID already done
2391 Int_t indexmother2 = -1;
2392 source2 = FindMother(TMath::Abs(track2->GetLabel()),mcEvent, indexmother2);
2393 pdg2 = CheckPdg(TMath::Abs(track2->GetLabel()),mcEvent);
2395 if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0)) {
2396 if(source == kElectronfromconversion) {
2397 valueangle[2] = kElectronfromconversionboth;
2398 valuensparseDeltaPhiMaps[4] = kElectronfromconversionboth;
2401 if(source == kElectronfrompi0) {
2402 valueangle[2] = kElectronfrompi0both;
2403 valuensparseDeltaPhiMaps[4] = kElectronfrompi0both;
2405 if(source == kElectronfrometa) {
2406 valueangle[2] = kElectronfrometaboth;
2407 valuensparseDeltaPhiMaps[4] = kElectronfrometaboth;
2413 if(fAlgorithmMA && (!fAODAnalysis))
2416 AliESDtrack *esdtrack2 = dynamic_cast<AliESDtrack *>(track2);
2417 AliESDtrack *esdtrack1 = dynamic_cast<AliESDtrack *>(track1);
2418 if((!esdtrack2) || (!esdtrack1)) continue;
2423 Double_t xt1; //radial position track 1 at the DCA point
2424 Double_t xt2; //radial position track 2 at the DCA point
2426 Double_t dca12 = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1);
2429 if(dca12 > fMaxdca) continue;
2431 //Momento of the track extrapolated to DCA track-track
2433 Bool_t hasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1);
2435 Bool_t hasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2);
2437 if(!hasdcaT1 || !hasdcaT2) AliWarning("It could be a problem in the extrapolation");
2439 //track1-track2 Invariant Mass
2440 Double_t eMass = 0.000510998910; //Electron mass in GeV
2441 Double_t pP1 = sqrt(p1[0]*p1[0]+p1[1]*p1[1]+p1[2]*p1[2]); //Track 1 momentum
2442 Double_t pP2 = sqrt(p2[0]*p2[0]+p2[1]*p2[1]+p2[2]*p2[2]); //Track 2 momentum
2443 Double_t eE1 = TMath::Sqrt(pP1*pP1+eMass*eMass);
2444 Double_t eE2 = TMath::Sqrt(pP2*pP2+eMass*eMass);
2446 //TLorentzVector v1(p1[0],p1[1],p1[2],sqrt(eMass*eMass+pP1*pP1));
2447 //TLorentzVector v2(p2[0],p2[1],p2[2],sqrt(eMass*eMass+pP2*pP2));
2448 //Double_t imass = (v1+v2).M(); //Invariant Mass
2449 //Double_t angle3D = v1.Angle(v2.Vect()); //Opening Angle (Total Angle)
2452 v3D1.SetXYZ(p1[0],p1[1],p1[2]);
2453 v3D2.SetXYZ(p2[0],p2[1],p2[2]);
2454 Double_t openingangle = TVector2::Phi_0_2pi(v3D2.Angle(v3D1));
2457 TVector3 motherrec = v3D1 + v3D2;
2458 Double_t invmass = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(motherrec.Px()*motherrec.Px()+motherrec.Py()*motherrec.Py()+motherrec.Pz()*motherrec.Pz()));
2461 //TVector3 vectordiff = v3D1 - v3D2;
2462 //Double_t diffphi = TVector2::Phi_0_2pi(vectordiff.Phi());
2463 //Double_t massxy = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(diffphi)));
2466 //Double_t difftheta = TVector2::Phi_0_2pi(vectordiff.Eta());
2467 //Double_t massrz = TMath::Sqrt((eE1+eE2)*(eE1+eE2)-(pP1*pP1+pP2*pP2+2*pP1*pP2*TMath::Cos(difftheta)));
2470 Float_t fCharge1 = track1->Charge();
2471 Float_t fCharge2 = track2->Charge();
2474 //valueangle[0] = diffphi;
2475 //valueangle[1] = difftheta;
2476 valueangle[0] = openingangle;
2477 if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);
2478 else fOppSignAngle->Fill(&valueangle[0]);
2481 if(openingangle > fMaxopening3D) continue;
2482 //if(difftheta > fMaxopeningtheta) continue;
2483 //if(diffphi > fMaxopeningphi) continue;
2486 valuensparseDeltaPhiMaps[3] = invmass;
2487 if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2488 else fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2491 if(invmass < fMaxInvmass) {
2492 if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;
2493 if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;
2500 Int_t fPDGtrack1 = 11;
2501 Int_t fPDGtrack2 = 11;
2503 Float_t fCharge1 = track1->Charge();
2504 Float_t fCharge2 = track2->Charge();
2506 if(fCharge1>0) fPDGtrack1 = -11;
2507 if(fCharge2>0) fPDGtrack2 = -11;
2509 AliKFParticle ktrack1(*track1, fPDGtrack1);
2510 AliKFParticle ktrack2(*track2, fPDGtrack2);
2511 AliKFParticle recoGamma(ktrack1, ktrack2);
2513 //Reconstruction Cuts
2514 if(recoGamma.GetNDF()<1) continue;
2515 Double_t chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();
2516 if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;
2519 //Double_t dca12 = ktrack1.GetDistanceFromParticle(ktrack2);
2520 //if(dca12 > fMaxdca) continue;
2522 // if set mass constraint
2523 if(fSetMassConstraint && pVtx) {
2524 AliKFVertex primV(*pVtx);
2528 recoGamma.SetProductionVertex(primV);
2529 recoGamma.SetMassConstraint(0,0.0001);
2535 recoGamma.GetMass(imass,width);
2537 //Opening Angle (Total Angle)
2538 Double_t angle = ktrack1.GetAngle(ktrack2);
2539 valueangle[0] = angle;
2540 if((fCharge1*fCharge2)>0.0) fSameSignAngle->Fill(&valueangle[0]);
2541 else fOppSignAngle->Fill(&valueangle[0]);
2544 if(angle > fMaxopening3D) continue;
2547 valuensparseDeltaPhiMaps[3] = imass;
2548 if((fCharge1*fCharge2)>0.0) fSameSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2550 fOppSignDeltaPhiMaps->Fill(&valuensparseDeltaPhiMaps[0]);
2552 if(valueangle[2] == kElectronfromconversionboth) {
2553 printf("Reconstructed charge1 %f, charge 2 %f and invmass %f",fCharge1,fCharge2,imass);
2554 printf("MC charge1 %d, charge 2 %d",pdg1,pdg2);
2555 printf("DCA %f",dca12);
2556 printf("Number of found %d",numberfound);
2562 if(imass < fMaxInvmass) {
2563 if((fCharge1*fCharge2)<0.0) oppositetaggedphotonic=kTRUE;
2564 if((fCharge1*fCharge2)>0.0) sametaggedphotonic=kTRUE;
2569 if(oppositetaggedphotonic && sametaggedphotonic){
2573 if(!oppositetaggedphotonic && sametaggedphotonic){
2577 if(oppositetaggedphotonic && !sametaggedphotonic){
2582 return taggedphotonic;
2584 //_________________________________________________________________________
2585 Int_t AliAnalysisTaskFlowTPCTOFEPSP::FindMother(Int_t tr, AliMCEvent *mcEvent, Int_t &indexmother){
2587 // Find the mother if MC
2590 if(!mcEvent) return 0;
2592 Int_t pdg = CheckPdg(tr,mcEvent);
2593 if(TMath::Abs(pdg)!= 11) {
2598 indexmother = IsMotherGamma(tr,mcEvent);
2599 if(indexmother > 0) return kElectronfromconversion;
2600 indexmother = IsMotherPi0(tr,mcEvent);
2601 if(indexmother > 0) return kElectronfrompi0;
2602 indexmother = IsMotherC(tr,mcEvent);
2603 if(indexmother > 0) return kElectronfromC;
2604 indexmother = IsMotherB(tr,mcEvent);
2605 if(indexmother > 0) return kElectronfromB;
2606 indexmother = IsMotherEta(tr,mcEvent);
2607 if(indexmother > 0) return kElectronfrometa;
2609 return kElectronfromother;
2613 //____________________________________________________________________________________________________________
2614 Int_t AliAnalysisTaskFlowTPCTOFEPSP::CheckPdg(Int_t tr, AliMCEvent* mcEvent) {
2617 // Return the pdg of the particle
2622 if(tr < 0) return pdgcode;
2624 if(!mcEvent) return pdgcode;
2626 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2629 if(mctrack->IsA() == AliMCParticle::Class()) {
2630 AliMCParticle *mctrackesd = NULL;
2631 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;
2632 pdgcode = mctrackesd->PdgCode();
2635 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2636 AliAODMCParticle *mctrackaod = NULL;
2637 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;
2638 pdgcode = mctrackaod->GetPdgCode();
2645 //____________________________________________________________________________________________________________
2646 Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherGamma(Int_t tr, AliMCEvent* mcEvent) {
2649 // Return the lab of gamma mother or -1 if not gamma
2652 if(tr < 0) return -1;
2653 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2655 if(mctrack->IsA() == AliMCParticle::Class()) {
2656 AliMCParticle *mctrackesd = NULL;
2657 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2658 TParticle *particle = 0x0;
2659 particle = mctrackesd->Particle();
2661 if(!particle) return -1;
2662 Int_t imother = particle->GetFirstMother();
2663 if(imother < 0) return -1;
2664 AliMCParticle *mothertrack = NULL;
2665 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2666 TParticle * mother = mothertrack->Particle();
2667 if(!mother) return -1;
2669 Int_t pdg = mother->GetPdgCode();
2670 if(TMath::Abs(pdg) == 22) return imother;
2671 if(TMath::Abs(pdg) == 11) {
2672 return IsMotherGamma(imother,mcEvent);
2677 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2678 AliAODMCParticle *mctrackaod = NULL;
2679 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2681 Int_t imother = mctrackaod->GetMother();
2682 if(imother < 0) return -1;
2683 AliAODMCParticle *mothertrack = NULL;
2684 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2686 Int_t pdg = mothertrack->GetPdgCode();
2687 if(TMath::Abs(pdg) == 22) return imother;
2688 if(TMath::Abs(pdg) == 11) {
2689 return IsMotherGamma(imother,mcEvent);
2700 //____________________________________________________________________________________________________________
2701 Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherPi0(Int_t tr, AliMCEvent* mcEvent) {
2704 // Return the lab of pi0 mother or -1 if not pi0
2707 if(tr < 0) return -1;
2708 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2710 if(mctrack->IsA() == AliMCParticle::Class()) {
2711 AliMCParticle *mctrackesd = NULL;
2712 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2713 TParticle *particle = 0x0;
2714 particle = mctrackesd->Particle();
2716 if(!particle) return -1;
2717 Int_t imother = particle->GetFirstMother();
2718 if(imother < 0) return -1;
2719 AliMCParticle *mothertrack = NULL;
2720 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2721 TParticle * mother = mothertrack->Particle();
2722 if(!mother) return -1;
2724 Int_t pdg = mother->GetPdgCode();
2725 if(TMath::Abs(pdg) == 111) return imother;
2726 if(TMath::Abs(pdg) == 11) {
2727 return IsMotherPi0(imother,mcEvent);
2732 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2733 AliAODMCParticle *mctrackaod = NULL;
2734 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2736 Int_t imother = mctrackaod->GetMother();
2737 if(imother < 0) return -1;
2738 AliAODMCParticle *mothertrack = NULL;
2739 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2741 Int_t pdg = mothertrack->GetPdgCode();
2742 if(TMath::Abs(pdg) == 111) return imother;
2743 if(TMath::Abs(pdg) == 11) {
2744 return IsMotherPi0(imother,mcEvent);
2752 //____________________________________________________________________________________________________________
2753 Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherC(Int_t tr, AliMCEvent* mcEvent) {
2756 // Return the lab of signal mother or -1 if not signal
2759 if(tr < 0) return -1;
2760 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2762 if(mctrack->IsA() == AliMCParticle::Class()) {
2763 AliMCParticle *mctrackesd = NULL;
2764 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2765 TParticle *particle = 0x0;
2766 particle = mctrackesd->Particle();
2768 if(!particle) return -1;
2769 Int_t imother = particle->GetFirstMother();
2770 if(imother < 0) return -1;
2771 AliMCParticle *mothertrack = NULL;
2772 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2773 TParticle * mother = mothertrack->Particle();
2774 if(!mother) return -1;
2776 Int_t pdg = mother->GetPdgCode();
2777 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;
2778 if(TMath::Abs(pdg) == 11) {
2779 return IsMotherC(imother,mcEvent);
2784 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2785 AliAODMCParticle *mctrackaod = NULL;
2786 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2788 Int_t imother = mctrackaod->GetMother();
2789 if(imother < 0) return -1;
2790 AliAODMCParticle *mothertrack = NULL;
2791 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2793 Int_t pdg = mothertrack->GetPdgCode();
2794 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;
2795 if(TMath::Abs(pdg) == 11) {
2796 return IsMotherC(imother,mcEvent);
2804 //____________________________________________________________________________________________________________
2805 Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherB(Int_t tr, AliMCEvent* mcEvent) {
2808 // Return the lab of signal mother or -1 if not signal
2811 if(tr < 0) return -1;
2812 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2814 if(mctrack->IsA() == AliMCParticle::Class()) {
2815 AliMCParticle *mctrackesd = NULL;
2816 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2817 TParticle *particle = 0x0;
2818 particle = mctrackesd->Particle();
2820 if(!particle) return -1;
2821 Int_t imother = particle->GetFirstMother();
2822 if(imother < 0) return -1;
2823 AliMCParticle *mothertrack = NULL;
2824 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2825 TParticle * mother = mothertrack->Particle();
2826 if(!mother) return -1;
2828 Int_t pdg = mother->GetPdgCode();
2829 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;
2830 if(TMath::Abs(pdg) == 11) {
2831 return IsMotherB(imother,mcEvent);
2836 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2837 AliAODMCParticle *mctrackaod = NULL;
2838 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2840 Int_t imother = mctrackaod->GetMother();
2841 if(imother < 0) return -1;
2842 AliAODMCParticle *mothertrack = NULL;
2843 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2845 Int_t pdg = mothertrack->GetPdgCode();
2846 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;
2847 if(TMath::Abs(pdg) == 11) {
2848 return IsMotherB(imother,mcEvent);
2856 //____________________________________________________________________________________________________________
2857 Int_t AliAnalysisTaskFlowTPCTOFEPSP::IsMotherEta(Int_t tr, AliMCEvent* mcEvent) {
2860 // Return the lab of pi0 mother or -1 if not pi0
2863 if(tr < 0) return -1;
2864 AliVParticle *mctrack = mcEvent->GetTrack(tr);
2866 if(mctrack->IsA() == AliMCParticle::Class()) {
2867 AliMCParticle *mctrackesd = NULL;
2868 if(!(mctrackesd = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2869 TParticle *particle = 0x0;
2870 particle = mctrackesd->Particle();
2872 if(!particle) return -1;
2873 Int_t imother = particle->GetFirstMother();
2874 if(imother < 0) return -1;
2875 AliMCParticle *mothertrack = NULL;
2876 if(!(mothertrack = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2877 TParticle * mother = mothertrack->Particle();
2878 if(!mother) return -1;
2880 Int_t pdg = mother->GetPdgCode();
2881 if(TMath::Abs(pdg) == 221) return imother;
2882 if(TMath::Abs(pdg) == 11) {
2883 return IsMotherEta(imother,mcEvent);
2888 if(mctrack->IsA() == AliAODMCParticle::Class()) {
2889 AliAODMCParticle *mctrackaod = NULL;
2890 if(!(mctrackaod = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(tr))))) return -1;
2892 Int_t imother = mctrackaod->GetMother();
2893 if(imother < 0) return -1;
2894 AliAODMCParticle *mothertrack = NULL;
2895 if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(mcEvent->GetTrack(TMath::Abs(imother))))) return -1;
2897 Int_t pdg = mothertrack->GetPdgCode();
2898 if(TMath::Abs(pdg) == 221) return imother;
2899 if(TMath::Abs(pdg) == 11) {
2900 return IsMotherEta(imother,mcEvent);