PWGHF/hfe/AliAnalysisTaskFlowTPCTOFEPSP.cxx

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