PWGCF/Correlations/DPhi/AliAnalysisTaskPhiCorrelations.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id:$ */
  17
  18 #include <TROOT.h>
  19 #include <TInterpreter.h>
  20 #include <TChain.h>
  21 #include <TFile.h>
  22 #include <TList.h>
  23 #include <TMath.h>
  24 #include <TTree.h>
  25 #include <TH2F.h>
  26 #include <TH3F.h>
  27 #include <TRandom.h>
  28 #include <TParameter.h>
  29
  30 #include "AliAnalysisTaskPhiCorrelations.h"
  31 #include "AliAnalyseLeadingTrackUE.h"
  32 #include "AliUEHistograms.h"
  33 #include "AliUEHist.h"
  34
  35 #include "AliAnalysisManager.h"
  36 #include "AliAODHandler.h"
  37 #include "AliAODInputHandler.h"
  38 #include "AliAODMCParticle.h"
  39 #include "AliInputEventHandler.h"
  40 #include "AliLog.h"
  41 #include "AliMCEventHandler.h"
  42 #include "AliVParticle.h"
  43 #include "AliCFContainer.h"
  44 #include "AliEventplane.h"
  45
  46 #include "AliESDEvent.h"
  47 #include "AliESDInputHandler.h"
  48 #include "AliMultiplicity.h"
  49 #include "AliCentrality.h"
  50 #include "AliStack.h"
  51 #include "AliAODMCHeader.h"
  52 #include "AliGenCocktailEventHeader.h"
  53 #include "AliGenEventHeader.h"
  54 #include "AliCollisionGeometry.h"
  55
  56 #include "AliEventPoolManager.h"
  57
  58 #include "AliESDZDC.h"
  59 #include "AliESDtrackCuts.h"
  60
  61 #include "AliHelperPID.h"
  62 #include "AliAnalysisUtils.h"
  63 #include "TMap.h"
  64
  65 ////////////////////////////////////////////////////////////////////////
  66 //
  67 // Analysis class for azimuthal correlation studies
  68 // Based on the UE task from Sara Vallero and Jan Fiete
  69 //
  70 // This class needs input AODs.
  71 // The output is a list of analysis-specific containers.
  72 //
  73 // The AOD can be either connected to the InputEventHandler
  74 // for a chain of AOD files
  75 // or
  76 // to the OutputEventHandler
  77 // for a chain of ESD files,
  78 // in this case the class should be in the train after the jet-finder
  79 //
  80 //    Authors:
  81 //    Jan Fiete Grosse-Oetringhaus
  82 //
  83 ////////////////////////////////////////////////////////////////////////
  84
  85
  86 ClassImp( AliAnalysisTaskPhiCorrelations )
  87 ClassImp( AliDPhiBasicParticle )
  88
  89 //____________________________________________________________________
  90 AliAnalysisTaskPhiCorrelations:: AliAnalysisTaskPhiCorrelations(const char* name):
  91 AliAnalysisTask(name,""),
  92 // general configuration
  93 fDebug(0),
  94 fMode(0),
  95 fReduceMemoryFootprint(kFALSE),
  96 fFillMixed(kTRUE),
  97 fMixingTracks(50000),
  98 fTwoTrackEfficiencyStudy(kFALSE),
  99 fTwoTrackEfficiencyCut(0),
 100 fTwoTrackCutMinRadius(0.8),
 101 fUseVtxAxis(kFALSE),
 102 fCourseCentralityBinning(kFALSE),
 103 fSkipTrigger(kFALSE),
 104 fInjectedSignals(kFALSE),
 105 fRandomizeReactionPlane(kFALSE),
 106 fHelperPID(0x0),
 107 fAnalysisUtils(0x0),
 108 fMap(0x0),
 109 // pointers to UE classes
 110 fAnalyseUE(0x0),
 111 fHistos(0x0),
 112 fHistosMixed(0),
 113 fEfficiencyCorrectionTriggers(0),
 114 fEfficiencyCorrectionAssociated(0),
 115 fCentralityWeights(0),
 116 // handlers and events
 117 fAOD(0x0),
 118 fESD(0x0),
 119 fArrayMC(0x0),
 120 fInputHandler(0x0),
 121 fMcEvent(0x0),
 122 fMcHandler(0x0),
 123 fPoolMgr(0x0),
 124 // histogram settings
 125 fListOfHistos(0x0),
 126 // event QA
 127 fnTracksVertex(1),  // QA tracks pointing to principal vertex (= 3 default)
 128 fZVertex(7.),
 129 fAcceptOnlyMuEvents(kFALSE),
 130 fCentralityMethod("V0M"),
 131 // track cuts
 132 fTrackEtaCut(0.8),
 133 fTrackEtaCutMin(-1.),
 134 fTrackPhiCutEvPlMin(0.),
 135 fTrackPhiCutEvPlMax(0.),
 136 fOnlyOneEtaSide(0),
 137 fPtMin(0.5),
 138 fDCAXYCut(0),
 139 fSharedClusterCut(-1),
 140 fCrossedRowsCut(-1),
 141 fFoundFractionCut(-1),
 142 fFilterBit(0xFF),
 143 fTrackStatus(0),
 144 fSelectBit(AliVEvent::kMB|AliVEvent::kUserDefined),
 145 fUseChargeHadrons(kFALSE),
 146 fParticleSpeciesTrigger(-1),
 147 fParticleSpeciesAssociated(-1),
 148 fCheckMotherPDG(kTRUE),
 149 fTrackletDphiCut(9999999.),
 150 fSelectCharge(0),
 151 fTriggerSelectCharge(0),
 152 fAssociatedSelectCharge(0),
 153 fTriggerRestrictEta(-1),
 154 fEtaOrdering(kFALSE),
 155 fCutConversions(kFALSE),
 156 fCutResonances(kFALSE),
 157 fRejectResonanceDaughters(-1),
 158 fFillOnlyStep0(kFALSE),
 159 fSkipStep6(kFALSE),
 160 fRejectCentralityOutliers(kFALSE),
 161 fRejectZeroTrackEvents(kFALSE),
 162 fRemoveWeakDecays(kFALSE),
 163 fRemoveDuplicates(kFALSE),
 164 fSkipFastCluster(kFALSE),
 165 fWeightPerEvent(kFALSE),
 166 fCustomBinning(),
 167 fPtOrder(kTRUE),
 168 fTriggersFromDetector(0),
 169 fAssociatedFromDetector(0),
 170 fMCUseUncheckedCentrality(kFALSE),
 171 fFillpT(kFALSE)
 172 {
 173   // Default constructor
 174   // Define input and output slots here
 175   // Input slot #0 works with a TChain
 176   DefineInput(0, TChain::Class());
 177   // Output slot #0 writes into a TList container
 178   DefineOutput(0, TList::Class());
 179 }
 180
 181 AliAnalysisTaskPhiCorrelations::~AliAnalysisTaskPhiCorrelations()
 182 {
 183   // destructor
 184
 185   if (fListOfHistos  && !AliAnalysisManager::GetAnalysisManager()->IsProofMode())
 186     delete fListOfHistos;
 187 }
 188
 189 //____________________________________________________________________
 190 void AliAnalysisTaskPhiCorrelations::ConnectInputData(Option_t* /*option*/)
 191 {
 192
 193   // Connect the input data
 194   if (fDebug > 1) AliInfo("ConnectInputData() ");
 195
 196   // Since AODs can either be connected to the InputEventHandler
 197   // or to the OutputEventHandler ( the AOD is created by a previus task in the train )
 198   // we need to get the pointer to the AODEvent correctly.
 199
 200   // Delta AODs are also accepted.
 201
 202   TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
 203
 204   if( handler && handler->InheritsFrom("AliAODInputHandler") )
 205   { // input AOD
 206     fAOD = ((AliAODInputHandler*)handler)->GetEvent();
 207     if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODInputHandler");
 208   }
 209   else
 210   {  //output AOD
 211     handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
 212     if (handler && handler->InheritsFrom("AliAODHandler") )
 213     {
 214       fAOD = ((AliAODHandler*)handler)->GetAOD();
 215       if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliAODHandler");
 216     }
 217     else
 218     {  // no AOD
 219       AliWarning("I can't get any AOD Event Handler");
 220     }
 221   }
 222
 223   if (handler && handler->InheritsFrom("AliESDInputHandler") )
 224   { // input ESD
 225     // pointer received per event in ::Exec
 226     if (fDebug > 1) AliInfo(" ==== Tracks and Jets from AliESDInputHandler");
 227   }
 228
 229   // Initialize common pointers
 230   Initialize();
 231 }
 232
 233 //____________________________________________________________________
 234 void  AliAnalysisTaskPhiCorrelations::CreateOutputObjects()
 235 {
 236   // Create the output container
 237
 238   if (fDebug > 1) AliInfo("CreateOutputObjects()");
 239
 240   // Initialize class with main algorithms, event and track selection.
 241   fAnalyseUE = new AliAnalyseLeadingTrackUE();
 242   fAnalyseUE->SetParticleSelectionCriteria(fFilterBit, fUseChargeHadrons, fTrackEtaCut, fTrackEtaCutMin, fPtMin);
 243   fAnalyseUE->SetDCAXYCut(fDCAXYCut);
 244   fAnalyseUE->SetSharedClusterCut(fSharedClusterCut);
 245   fAnalyseUE->SetCrossedRowsCut(fCrossedRowsCut);
 246   fAnalyseUE->SetFoundFractionCut(fFoundFractionCut);
 247   fAnalyseUE->SetTrackStatus(fTrackStatus);
 248   fAnalyseUE->SetCheckMotherPDG(fCheckMotherPDG);
 249   fAnalyseUE->SetDebug(fDebug);
 250   fAnalyseUE->DefineESDCuts(fFilterBit);
 251   fAnalyseUE->SetEventSelection(fSelectBit);
 252   fAnalyseUE->SetHelperPID(fHelperPID);
 253   if(fTrackPhiCutEvPlMax > 0.0001)
 254     fAnalyseUE->SetParticlePhiCutEventPlane(fTrackPhiCutEvPlMin,fTrackPhiCutEvPlMax);
 255   if ((fParticleSpeciesTrigger != -1 || fParticleSpeciesAssociated != -1) && !fHelperPID)
 256     AliFatal("HelperPID object should be set in the steering macro");
 257
 258   // Initialize output list of containers
 259   if (fListOfHistos != NULL){
 260         delete fListOfHistos;
 261         fListOfHistos = NULL;
 262         }
 263   if (!fListOfHistos){
 264         fListOfHistos = new TList();
 265         fListOfHistos->SetOwner(kTRUE);
 266         }
 267
 268   // Initialize class to handle histograms
 269   TString histType = "4R";
 270   if (fUseVtxAxis == 1)
 271     histType = "5R";
 272   else if (fUseVtxAxis == 2)
 273     histType = "6R";
 274   if (fCourseCentralityBinning)
 275     histType += "C";
 276   fHistos = new AliUEHistograms("AliUEHistogramsSame", histType, fCustomBinning);
 277   fHistosMixed = new AliUEHistograms("AliUEHistogramsMixed", histType, fCustomBinning);
 278
 279   fHistos->SetSelectCharge(fSelectCharge);
 280   fHistosMixed->SetSelectCharge(fSelectCharge);
 281
 282   fHistos->SetSelectTriggerCharge(fTriggerSelectCharge);
 283   fHistosMixed->SetSelectTriggerCharge(fTriggerSelectCharge);
 284
 285   fHistos->SetSelectAssociatedCharge(fAssociatedSelectCharge);
 286   fHistosMixed->SetSelectAssociatedCharge(fAssociatedSelectCharge);
 287
 288   fHistos->SetTriggerRestrictEta(fTriggerRestrictEta);
 289   fHistosMixed->SetTriggerRestrictEta(fTriggerRestrictEta);
 290
 291   fHistos->SetOnlyOneEtaSide(fOnlyOneEtaSide);
 292   fHistosMixed->SetOnlyOneEtaSide(fOnlyOneEtaSide);
 293
 294   fHistos->SetEtaOrdering(fEtaOrdering);
 295   fHistosMixed->SetEtaOrdering(fEtaOrdering);
 296
 297   fHistos->SetPairCuts(fCutConversions, fCutResonances);
 298   fHistosMixed->SetPairCuts(fCutConversions, fCutResonances);
 299
 300   fHistos->SetRejectResonanceDaughters(fRejectResonanceDaughters);
 301   fHistosMixed->SetRejectResonanceDaughters(fRejectResonanceDaughters);
 302
 303   fHistos->SetTrackEtaCut(fTrackEtaCut);
 304   fHistosMixed->SetTrackEtaCut(fTrackEtaCut);
 305
 306   fHistos->SetWeightPerEvent(fWeightPerEvent);
 307   fHistosMixed->SetWeightPerEvent(fWeightPerEvent);
 308
 309   fHistos->SetPtOrder(fPtOrder);
 310   fHistosMixed->SetPtOrder(fPtOrder);
 311
 312   fHistos->SetTwoTrackCutMinRadius(fTwoTrackCutMinRadius);
 313   fHistosMixed->SetTwoTrackCutMinRadius(fTwoTrackCutMinRadius);
 314
 315   if (fEfficiencyCorrectionTriggers)
 316    {
 317     fHistos->SetEfficiencyCorrectionTriggers(fEfficiencyCorrectionTriggers);
 318     fHistosMixed->SetEfficiencyCorrectionTriggers((THnF*) fEfficiencyCorrectionTriggers->Clone());
 319    }
 320   if (fEfficiencyCorrectionAssociated)
 321   {
 322     fHistos->SetEfficiencyCorrectionAssociated(fEfficiencyCorrectionAssociated);
 323     fHistosMixed->SetEfficiencyCorrectionAssociated((THnF*) fEfficiencyCorrectionAssociated->Clone());
 324   }
 325
 326   // add histograms to list
 327   fListOfHistos->Add(fHistos);
 328   fListOfHistos->Add(fHistosMixed);
 329   // add HelperPID to list
 330   if (fHelperPID)
 331     fListOfHistos->Add(fHelperPID);
 332   // add TMap to list
 333   if (fMap)
 334     fListOfHistos->Add(fMap);
 335
 336   fListOfHistos->Add(new TH2F("processIDs", ";#Delta#phi;process id", 100, -0.5 * TMath::Pi(), 1.5 * TMath::Pi(), kPNoProcess + 1, -0.5, kPNoProcess + 0.5));
 337   fListOfHistos->Add(new TH1F("eventStat", ";;events", 4, -0.5, 3.5));
 338   fListOfHistos->Add(new TH2F("mixedDist", ";centrality;tracks;events", 101, 0, 101, 200, 0, fMixingTracks * 1.5));
 339   fListOfHistos->Add(new TH2F("referenceMultiplicity", ";centrality;tracks;events", 101, 0, 101, 200, 0, 200));
 340   if (fCentralityMethod == "V0A_MANUAL")
 341   {
 342     fListOfHistos->Add(new TH2F("V0AMult", "V0A multiplicity;V0A multiplicity;V0A multiplicity (scaled)", 1000, -.5, 999.5, 1000, -.5, 999.5));
 343     fListOfHistos->Add(new TH2F("V0AMultCorrelation", "V0A multiplicity;V0A multiplicity;SPD tracklets", 1000, -.5, 999.5, 1000, -.5, 999.5));
 344   }
 345   if (fTriggersFromDetector == 1 || fTriggersFromDetector == 2 || fAssociatedFromDetector == 1 || fAssociatedFromDetector == 2)
 346     fListOfHistos->Add(new TH1F("V0SingleCells", "V0 single cell multiplicity;multiplicity;events", 100, -0.5, 99.5));
 347   if (fTriggersFromDetector == 3 || fAssociatedFromDetector == 3)
 348     fListOfHistos->Add(new TH1F("DphiTrklets", "tracklets Dphi;#Delta#phi,trklets (mrad);entries", 100, -100, 100));
 349
 350   PostData(0,fListOfHistos);
 351
 352   // Add task configuration to output list
 353   AddSettingsTree();
 354
 355   // event mixing
 356   Int_t poolsize   = 1000;  // Maximum number of events, ignored in the present implemention of AliEventPoolManager
 357
 358   Int_t nCentralityBins  = fHistos->GetUEHist(2)->GetEventHist()->GetNBins(1);
 359   Double_t* centralityBins = (Double_t*) fHistos->GetUEHist(2)->GetEventHist()->GetAxis(1, 0)->GetXbins()->GetArray();
 360
 361   const Int_t kNZvtxBins  = 10+(1+10)*4;
 362   // bins for further buffers are shifted by 100 cm
 363   Double_t vertexBins[kNZvtxBins+1] = { -10,   -8,  -6,  -4,  -2,   0,   2,   4,   6,   8,  10,
 364                                        90,  92,  94,  96,  98, 100, 102, 104, 106, 108, 110,
 365                                       190, 192, 194, 196, 198, 200, 202, 204, 206, 208, 210,
 366                                       290, 292, 294, 296, 298, 300, 302, 304, 306, 308, 310,
 367                                       390, 392, 394, 396, 398, 400, 402, 404, 406, 408, 410 };
 368
 369   Int_t nZvtxBins  = kNZvtxBins;
 370   Double_t* zvtxbin = vertexBins;
 371
 372   if (fMode == 0 && fHistos->GetUEHist(2)->GetEventHist()->GetNVar() > 2)
 373   {
 374     nZvtxBins = fHistos->GetUEHist(2)->GetEventHist()->GetNBins(2);
 375     zvtxbin = (Double_t*) fHistos->GetUEHist(2)->GetEventHist()->GetAxis(2, 0)->GetXbins()->GetArray();
 376   }
 377
 378   fPoolMgr = new AliEventPoolManager(poolsize, fMixingTracks, nCentralityBins, centralityBins, nZvtxBins, zvtxbin);
 379   fPoolMgr->SetTargetValues(fMixingTracks, 0.1, 5);
 380 }
 381
 382 //____________________________________________________________________
 383 void  AliAnalysisTaskPhiCorrelations::Exec(Option_t */*option*/)
 384 {
 385   // exec (per event)
 386   fAnalyseUE->NextEvent();
 387
 388   // receive ESD pointer if we are not running AOD analysis
 389   if (!fAOD)
 390   {
 391     AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
 392     if (handler && handler->InheritsFrom("AliESDInputHandler"))
 393       fESD = ((AliESDInputHandler*)handler)->GetEvent();
 394   }
 395
 396   if (fMode)
 397   {
 398     // correction mode
 399
 400     if (fMcHandler)
 401       fMcEvent = fMcHandler->MCEvent();
 402
 403     if (fAOD)
 404     {
 405       // array of MC particles
 406       fArrayMC = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
 407       if (!fArrayMC)
 408         AliFatal("No array of MC particles found !!!");
 409     }
 410
 411     AnalyseCorrectionMode();
 412   }
 413   else AnalyseDataMode();
 414 }
 415
 416 /******************** ANALYSIS METHODS *****************************/
 417
 418 //____________________________________________________________________
 419 void  AliAnalysisTaskPhiCorrelations::AddSettingsTree()
 420 {
 421   //Write settings to output list
 422   TTree *settingsTree   = new TTree("UEAnalysisSettings","Analysis Settings in UE estimation");
 423   settingsTree->Branch("fnTracksVertex", &fnTracksVertex,"nTracksVertex/I");
 424   settingsTree->Branch("fZVertex", &fZVertex,"ZVertex/D");
 425   settingsTree->Branch("fAcceptOnlyMuEvents", &fAcceptOnlyMuEvents,"AcceptOnlyMuEvents/O");
 426   //settingsTree->Branch("fCentralityMethod", fCentralityMethod.Data(),"CentralityMethod/C");
 427   settingsTree->Branch("fTrackEtaCut", &fTrackEtaCut, "TrackEtaCut/D");
 428   settingsTree->Branch("fTrackEtaCutMin", &fTrackEtaCutMin, "TrackEtaCutMin/D");
 429   settingsTree->Branch("fTrackPhiCutEvPlMin", &fTrackPhiCutEvPlMin, "TrackPhiCutEvPlMin/D");
 430   settingsTree->Branch("fTrackPhiCutEvPlMax", &fTrackPhiCutEvPlMax, "TrackPhiCutEvPlMax/D");
 431   settingsTree->Branch("fOnlyOneEtaSide", &fOnlyOneEtaSide,"OnlyOneEtaSide/I");
 432   settingsTree->Branch("fPtMin", &fPtMin, "PtMin/D");
 433   settingsTree->Branch("fFilterBit", &fFilterBit,"FilterBit/I");
 434   settingsTree->Branch("fSharedClusterCut", &fSharedClusterCut,"SharedClusterCut/D");
 435   settingsTree->Branch("fCrossedRowsCut", &fCrossedRowsCut,"CrossedRowsCut/I");
 436   settingsTree->Branch("fFoundFractionCut", &fFoundFractionCut,"FoundFractionCut/D");
 437   settingsTree->Branch("fTrackStatus", &fTrackStatus,"TrackStatus/I");
 438   settingsTree->Branch("fSelectBit", &fSelectBit,"EventSelectionBit/I");
 439   settingsTree->Branch("fUseChargeHadrons", &fUseChargeHadrons,"UseChHadrons/O");
 440   settingsTree->Branch("fParticleSpeciesTrigger", &fParticleSpeciesTrigger,"ParticleSpeciesTrigger/I");
 441   settingsTree->Branch("fParticleSpeciesAssociated", &fParticleSpeciesAssociated,"ParticleSpeciesAssociated/I");
 442   settingsTree->Branch("fCheckMotherPDG", &fCheckMotherPDG,"CheckMotherPDG/I");
 443   settingsTree->Branch("fSelectCharge", &fSelectCharge,"SelectCharge/I");
 444   settingsTree->Branch("fTriggerSelectCharge", &fTriggerSelectCharge,"TriggerSelectCharge/I");
 445   settingsTree->Branch("fAssociatedSelectCharge", &fAssociatedSelectCharge,"fAssociatedSelectCharge/I");
 446   settingsTree->Branch("fTriggerRestrictEta", &fTriggerRestrictEta,"TriggerRestrictEta/D");
 447   settingsTree->Branch("fEtaOrdering", &fEtaOrdering,"EtaOrdering/O");
 448   settingsTree->Branch("fCutConversions", &fCutConversions,"CutConversions/O");
 449   settingsTree->Branch("fCutResonances", &fCutResonances,"CutResonances/O");
 450   settingsTree->Branch("fRejectResonanceDaughters", &fRejectResonanceDaughters,"RejectResonanceDaughters/I");
 451   settingsTree->Branch("fFillpT", &fFillpT,"FillpT/O");
 452   settingsTree->Branch("fMixingTracks", &fMixingTracks,"MixingTracks/I");
 453   settingsTree->Branch("fSkipTrigger", &fSkipTrigger,"SkipTrigger/O");
 454   settingsTree->Branch("fInjectedSignals", &fInjectedSignals,"InjectedSignals/O");
 455   settingsTree->Branch("fRandomizeReactionPlane", &fRandomizeReactionPlane,"RandomizeReactionPlane/O");
 456   settingsTree->Branch("fRejectCentralityOutliers", &fRejectCentralityOutliers,"RejectCentralityOutliers/O");
 457   settingsTree->Branch("fRejectZeroTrackEvents", &fRejectZeroTrackEvents,"RejectZeroTrackEvents/O");
 458   settingsTree->Branch("fRemoveWeakDecays", &fRemoveWeakDecays,"RemoveWeakDecays/O");
 459   settingsTree->Branch("fRemoveDuplicates", &fRemoveDuplicates,"RemoveDuplicates/O");
 460   settingsTree->Branch("fSkipFastCluster", &fSkipFastCluster,"SkipFastCluster/O");
 461   settingsTree->Branch("fWeightPerEvent", &fWeightPerEvent,"WeightPerEvent/O");
 462   settingsTree->Branch("fPtOrder", &fPtOrder,"PtOrder/O");
 463   settingsTree->Branch("fTriggersFromDetector", &fTriggersFromDetector,"TriggersFromDetector/I");
 464   settingsTree->Branch("fAssociatedFromDetector", &fAssociatedFromDetector,"AssociatedFromDetector/I");
 465   settingsTree->Branch("fMCUseUncheckedCentrality", &fMCUseUncheckedCentrality,"MCUseUncheckedCentrality/O");
 466   settingsTree->Branch("fTwoTrackEfficiencyCut", &fTwoTrackEfficiencyCut,"TwoTrackEfficiencyCut/D");
 467   settingsTree->Branch("fTwoTrackCutMinRadius", &fTwoTrackCutMinRadius,"TwoTrackCutMinRadius/D");
 468
 469   //fCustomBinning
 470
 471   settingsTree->Fill();
 472   fListOfHistos->Add(settingsTree);
 473 }
 474
 475 //____________________________________________________________________
 476 void  AliAnalysisTaskPhiCorrelations::AnalyseCorrectionMode()
 477 {
 478   // Run the analysis on MC to get the correction maps
 479   //
 480
 481   if ( fDebug > 3 ) AliInfo( " Processing event in Corrections mode ..." );
 482
 483   Double_t centrality = 0;
 484
 485   if (fCentralityMethod.Length() > 0)
 486   {
 487     if (fCentralityMethod == "MC_b")
 488     {
 489       AliGenEventHeader* eventHeader = GetFirstHeader();
 490       if (!eventHeader)
 491       {
 492         // We avoid AliFatal here, because the AOD productions sometimes have events where the MC header is missing
 493         // (due to unreadable Kinematics) and we don't want to loose the whole job because of a few events
 494         AliError("Event header not found. Skipping this event.");
 495         fHistos->FillEvent(0, AliUEHist::kCFStepAnaTopology);
 496         return;
 497       }
 498
 499       AliCollisionGeometry* collGeometry = dynamic_cast<AliCollisionGeometry*> (eventHeader);
 500       if (!collGeometry)
 501       {
 502         eventHeader->Dump();
 503         AliFatal("Asking for MC_b centrality, but event header has no collision geometry information");
 504       }
 505
 506       centrality = collGeometry->ImpactParameter();
 507     }
 508     else if (fCentralityMethod == "nano")
 509     {
 510       centrality = (Float_t) gROOT->ProcessLine(Form("100.0 + 100.0 * ((AliNanoAODHeader*) %p)->GetCentrality(\"%s\")", fAOD->GetHeader(), fCentralityMethod.Data())) / 100 - 1.0;
 511     }
 512     else
 513     {
 514       AliCentrality *centralityObj = 0;
 515       if (fAOD)
 516         centralityObj = fAOD->GetHeader()->GetCentralityP();
 517       else if (fESD)
 518         centralityObj = fESD->GetCentrality();
 519
 520       if (centralityObj)
 521       {
 522         if (fMCUseUncheckedCentrality)
 523           centrality = centralityObj->GetCentralityPercentileUnchecked(fCentralityMethod);
 524         else
 525           centrality = centralityObj->GetCentralityPercentile(fCentralityMethod);
 526       }
 527       else
 528       {
 529         Printf("WARNING: Centrality object is 0");
 530         centrality = -1;
 531       }
 532     }
 533
 534     AliInfo(Form("Centrality is %f", centrality));
 535   }
 536
 537   // Support for ESD and AOD based analysis
 538   AliVEvent* inputEvent = fAOD;
 539   if (!inputEvent)
 540     inputEvent = fESD;
 541
 542   Float_t bSign = 0;
 543
 544   if (inputEvent)
 545   {
 546     fHistos->SetRunNumber(inputEvent->GetRunNumber());
 547     bSign = (inputEvent->GetMagneticField() > 0) ? 1 : -1;
 548   }
 549
 550   TObject* mc = fArrayMC;
 551   if (!mc)
 552     mc = fMcEvent;
 553
 554   // count all events
 555   fHistos->FillEvent(centrality, -1);
 556
 557   if (centrality < 0)
 558     return;
 559
 560   // Only consider MC events within the vtx-z region used also as cut on the reconstructed vertex
 561   TObject* vertexSupplier = fMcEvent;
 562   if (fAOD) // AOD
 563     vertexSupplier = fAOD->FindListObject(AliAODMCHeader::StdBranchName());
 564
 565   if (!fAnalyseUE->VertexSelection(vertexSupplier, 0, fZVertex))
 566     return;
 567
 568     Float_t zVtx = 0;
 569   if (fAOD)
 570     zVtx = ((AliAODMCHeader*) vertexSupplier)->GetVtxZ();
 571   else
 572     zVtx = fMcEvent->GetPrimaryVertex()->GetZ();
 573   Float_t weight = 1;
 574   if (fFillpT)
 575     weight = -1;
 576
 577   // For productions with injected signals, figure out above which label to skip particles/tracks
 578   Int_t skipParticlesAbove = 0;
 579   if (fInjectedSignals)
 580   {
 581     AliGenEventHeader* eventHeader = 0;
 582     Int_t headers = 0;
 583
 584     if (fMcEvent)
 585     {
 586       // ESD
 587       AliHeader* header = (AliHeader*) fMcEvent->Header();
 588       if (!header)
 589         AliFatal("fInjectedSignals set but no MC header found");
 590
 591       AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
 592       if (!cocktailHeader)
 593       {
 594         header->Dump();
 595         AliFatal("fInjectedSignals set but no MC cocktail header found");
 596       }
 597
 598       headers = cocktailHeader->GetHeaders()->GetEntries();
 599       eventHeader = dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->First());
 600
 601       if (fDebug > 4)
 602       {
 603         for (Int_t i=0; i<cocktailHeader->GetHeaders()->GetEntries(); i++)
 604         {
 605           AliGenEventHeader* headerTmp = dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->At(i));
 606           if (headerTmp)
 607             Printf("%d particles in header:", headerTmp->NProduced());
 608           cocktailHeader->GetHeaders()->At(i)->Dump();
 609         }
 610       }
 611     }
 612     else
 613     {
 614       // AOD
 615       AliAODMCHeader* header = (AliAODMCHeader*) fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 616       if (!header)
 617         AliFatal("fInjectedSignals set but no MC header found");
 618
 619       headers = header->GetNCocktailHeaders();
 620       eventHeader = header->GetCocktailHeader(0);
 621     }
 622
 623     if (!eventHeader)
 624     {
 625       // We avoid AliFatal here, because the AOD productions sometimes have events where the MC header is missing
 626       // (due to unreadable Kinematics) and we don't want to loose the whole job because of a few events
 627       AliError("First event header not found. Skipping this event.");
 628       fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
 629       return;
 630     }
 631
 632     skipParticlesAbove = eventHeader->NProduced();
 633     AliInfo(Form("Injected signals in this event (%d headers). Keeping particles/tracks of %s. Will skip particles/tracks above %d.", headers, eventHeader->ClassName(), skipParticlesAbove));
 634   }
 635
 636   if (fCentralityWeights && !AcceptEventCentralityWeight(centrality))
 637   {
 638     AliInfo(Form("Rejecting event because of centrality weighting: %f", centrality));
 639     fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
 640     return;
 641   }
 642
 643   // Get MC primaries
 644   // triggers
 645   TObjArray* tmpList = fAnalyseUE->GetAcceptedParticles(mc, 0, kTRUE, fParticleSpeciesTrigger, kTRUE);
 646   CleanUp(tmpList, mc, skipParticlesAbove);
 647   TObjArray* tracksMC = CloneAndReduceTrackList(tmpList);
 648   delete tmpList;
 649
 650   // associated
 651   TObjArray* tracksCorrelateMC = tracksMC;
 652   if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
 653   {
 654     tmpList = fAnalyseUE->GetAcceptedParticles(mc, 0, kTRUE, fParticleSpeciesAssociated, kTRUE);
 655     CleanUp(tmpList, mc, skipParticlesAbove);
 656     tracksCorrelateMC = CloneAndReduceTrackList(tmpList);
 657     delete tmpList;
 658   }
 659
 660   if (fRandomizeReactionPlane)
 661   {
 662     Double_t centralityDigits = centrality*1000. - (Int_t)(centrality*1000.);
 663     Double_t angle = TMath::TwoPi() * centralityDigits;
 664     AliInfo(Form("Shifting phi of all tracks by %f (digits %f)", angle, centralityDigits));
 665     ShiftTracks(tracksMC, angle);
 666     if (tracksCorrelateMC != tracksMC)
 667       ShiftTracks(tracksCorrelateMC, angle);
 668   }
 669
 670   if (fFillOnlyStep0)
 671     zVtx = 0;
 672
 673   // Event selection based on number of number of MC particles
 674   if (fRejectZeroTrackEvents && tracksMC->GetEntriesFast() == 0)
 675   {
 676     AliInfo(Form("Rejecting event due to kinematic selection: %f %d", centrality, tracksMC->GetEntriesFast()));
 677     fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
 678     if (tracksMC != tracksCorrelateMC)
 679       delete tracksCorrelateMC;
 680     delete tracksMC;
 681     return;
 682   }
 683
 684   // (MC-true all particles)
 685   // STEP 0
 686   fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepAll, tracksMC, tracksCorrelateMC, weight);
 687
 688   // mixed event
 689   if (fFillMixed)
 690   {
 691     AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx);
 692     if (fFillOnlyStep0)
 693       ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());
 694     if (pool->IsReady())
 695       for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++)
 696         fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepAll, tracksMC, pool->GetEvent(jMix), 1.0 / pool->GetCurrentNEvents(), (jMix == 0));
 697     pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelateMC));
 698   }
 699
 700 //   Printf("trigger: %d", ((AliInputEventHandler*)fInputHandler)->IsEventSelected());
 701
 702   // Trigger selection ************************************************
 703   if (!fFillOnlyStep0 && (fSkipTrigger || fAnalyseUE->TriggerSelection(fInputHandler)))
 704   {
 705     // (MC-true all particles)
 706     // STEP 1
 707     if (!fReduceMemoryFootprint)
 708       fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTriggered, tracksMC, tracksCorrelateMC, weight);
 709     else
 710       fHistos->FillEvent(centrality, AliUEHist::kCFStepTriggered);
 711
 712     if (!inputEvent) {
 713       AliFatal("UNEXPECTED: inputEvent is 0. Trigger selection should have failed");
 714       return;
 715     }
 716
 717     // Vertex selection *************************************************
 718     if (fAnalyseUE->VertexSelection(inputEvent, fnTracksVertex, fZVertex))
 719     {
 720       // fill here for tracking efficiency
 721       // loop over particle species
 722
 723       for (Int_t particleSpecies = 0; particleSpecies < 4; particleSpecies++)
 724       {
 725         TObjArray* primMCParticles = fAnalyseUE->GetAcceptedParticles(mc, 0x0, kTRUE, particleSpecies, kTRUE);
 726         TObjArray* primRecoTracksMatched = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, particleSpecies, kTRUE, kFALSE);
 727         TObjArray* allRecoTracksMatched  = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, particleSpecies, kTRUE, kFALSE);
 728         TObjArray* primRecoTracksMatchedPID = 0;
 729         TObjArray* allRecoTracksMatchedPID  = 0;
 730
 731         if (fHelperPID)
 732         {
 733           primRecoTracksMatchedPID = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, particleSpecies, kTRUE, kTRUE);
 734           allRecoTracksMatchedPID  = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, particleSpecies, kTRUE, kTRUE);
 735         }
 736
 737         CleanUp(primMCParticles, mc, skipParticlesAbove);
 738         CleanUp(primRecoTracksMatched, mc, skipParticlesAbove);
 739         CleanUp(allRecoTracksMatched, mc, skipParticlesAbove);
 740         CleanUp(primRecoTracksMatchedPID, mc, skipParticlesAbove);
 741         CleanUp(allRecoTracksMatchedPID, mc, skipParticlesAbove);
 742
 743         // select charges
 744         if (fTriggerSelectCharge != 0)
 745         {
 746           SelectCharge(primMCParticles);
 747           SelectCharge(primRecoTracksMatched);
 748           SelectCharge(allRecoTracksMatched);
 749           SelectCharge(primRecoTracksMatchedPID);
 750           SelectCharge(allRecoTracksMatchedPID);
 751         }
 752
 753         fHistos->FillTrackingEfficiency(primMCParticles, primRecoTracksMatched, allRecoTracksMatched, primRecoTracksMatchedPID, allRecoTracksMatchedPID, 0, particleSpecies, centrality, zVtx);
 754
 755 //      Printf("%d --> %d %d %d", particleSpecies, primMCParticles->GetEntries(), primRecoTracksMatched->GetEntries(), allRecoTracksMatched->GetEntries());
 756
 757         delete primMCParticles;
 758         delete primRecoTracksMatched;
 759         delete allRecoTracksMatched;
 760       }
 761
 762       TObjArray* fakeParticles = fAnalyseUE->GetFakeParticles(inputEvent, mc, kFALSE, -1, kTRUE);
 763       CleanUp((TObjArray*) fakeParticles->At(0), mc, skipParticlesAbove);
 764       CleanUp((TObjArray*) fakeParticles->At(1), mc, skipParticlesAbove);
 765
 766       fHistos->FillTrackingEfficiency(0, 0, 0, 0, 0, (TObjArray*) fakeParticles->At(2), 0, centrality, zVtx);
 767       fHistos->FillFakePt(fakeParticles, centrality);
 768 //       Printf(">>>>> %d %d %d fakes", ((TObjArray*) fakeParticles->At(0))->GetEntriesFast(), ((TObjArray*) fakeParticles->At(1))->GetEntriesFast(), ((TObjArray*) fakeParticles->At(2))->GetEntriesFast());
 769       delete fakeParticles;
 770
 771       // (MC-true all particles)
 772       // STEP 2
 773       if (!fReduceMemoryFootprint)
 774         fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepVertex, tracksMC, tracksCorrelateMC, weight);
 775       else
 776         fHistos->FillEvent(centrality, AliUEHist::kCFStepVertex);
 777
 778       // Get MC primaries that match reconstructed track
 779       // triggers
 780       tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, fParticleSpeciesTrigger, kTRUE);
 781       CleanUp(tmpList, mc, skipParticlesAbove);
 782       TObjArray* tracksRecoMatchedPrim = CloneAndReduceTrackList(tmpList);
 783       delete tmpList;
 784
 785       // associated
 786       TObjArray* tracksCorrelateRecoMatchedPrim = tracksRecoMatchedPrim;
 787       if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
 788       {
 789         tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kTRUE, fParticleSpeciesAssociated, kTRUE);
 790         CleanUp(tmpList, mc, skipParticlesAbove);
 791         tracksCorrelateRecoMatchedPrim = CloneAndReduceTrackList(tmpList);
 792         delete tmpList;
 793       }
 794
 795       // (RECO-matched (quantities from MC particle) primary particles)
 796       // STEP 4
 797       fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTrackedOnlyPrim, tracksRecoMatchedPrim, tracksCorrelateRecoMatchedPrim, weight);
 798
 799       // mixed event
 800       if (fFillMixed)
 801       {
 802         AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx + 200);
 803         if (pool->IsReady())
 804           for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++)
 805             fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTrackedOnlyPrim, tracksRecoMatchedPrim, pool->GetEvent(jMix), 1.0 / pool->GetCurrentNEvents(), (jMix == 0));
 806         pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelateRecoMatchedPrim));
 807       }
 808
 809       // Get MC primaries + secondaries that match reconstructed track
 810       // triggers
 811       tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, fParticleSpeciesTrigger, kTRUE);
 812       CleanUp(tmpList, mc, skipParticlesAbove);
 813       TObjArray* tracksRecoMatchedAll = CloneAndReduceTrackList(tmpList);
 814       delete tmpList;
 815
 816       // associated
 817       TObjArray* tracksCorrelateRecoMatchedAll = tracksRecoMatchedAll;
 818       if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
 819       {
 820         tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, mc, kFALSE, fParticleSpeciesAssociated, kTRUE);
 821         CleanUp(tmpList, mc, skipParticlesAbove);
 822         tracksCorrelateRecoMatchedAll = CloneAndReduceTrackList(tmpList);
 823         delete tmpList;
 824       }
 825
 826       // (RECO-matched (quantities from MC particle) all particles)
 827       // STEP 5
 828       fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTracked, tracksRecoMatchedAll, tracksCorrelateRecoMatchedAll, weight);
 829
 830       // mixed event
 831       if (fFillMixed)
 832       {
 833         AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx + 300);
 834         if (pool->IsReady())
 835           for (Int_t jMix=0; jMix<pool->GetCurrentNEvents(); jMix++)
 836             fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepTracked, tracksRecoMatchedAll, pool->GetEvent(jMix), 1.0 / pool->GetCurrentNEvents(), (jMix == 0));
 837         pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelateRecoMatchedAll));
 838       }
 839
 840       // Get RECO tracks
 841       // triggers
 842       tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesTrigger, kTRUE);
 843       CleanUp(tmpList, mc, skipParticlesAbove);
 844       TObjArray* tracks = CloneAndReduceTrackList(tmpList);
 845       delete tmpList;
 846
 847       // associated
 848       TObjArray* tracksCorrelate = tracks;
 849       if (fParticleSpeciesAssociated != fParticleSpeciesTrigger)
 850       {
 851         tmpList = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesAssociated, kTRUE);
 852         CleanUp(tmpList, mc, skipParticlesAbove);
 853         tracksCorrelate = CloneAndReduceTrackList(tmpList);
 854         delete tmpList;
 855       }
 856
 857       // (RECO all tracks)
 858       // STEP 6
 859       if (!fSkipStep6)
 860         fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracks, tracksCorrelate, weight);
 861
 862       // two track cut, STEP 8
 863       if (fTwoTrackEfficiencyCut > 0)
 864         fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracks, tracksCorrelate, weight, kTRUE, kTRUE, bSign, fTwoTrackEfficiencyCut);
 865
 866       // apply correction efficiency, STEP 10
 867       if (fEfficiencyCorrectionTriggers || fEfficiencyCorrectionAssociated)
 868       {
 869           // with or without two track efficiency depending on if fTwoTrackEfficiencyCut is set
 870         Bool_t twoTrackCut = (fTwoTrackEfficiencyCut > 0);
 871
 872         fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepCorrected, tracks, tracksCorrelate, weight, kTRUE, twoTrackCut, bSign, fTwoTrackEfficiencyCut, kTRUE);
 873       }
 874
 875       // mixed event
 876       if (fFillMixed)
 877       {
 878         AliEventPool* pool2 = fPoolMgr->GetEventPool(centrality, zVtx + 100);
 879         ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool2->NTracksInPool());
 880         if (pool2->IsReady())
 881         {
 882           for (Int_t jMix=0; jMix<pool2->GetCurrentNEvents(); jMix++)
 883           {
 884             // STEP 6
 885             if (!fSkipStep6)
 886               fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracks, pool2->GetEvent(jMix), 1.0 / pool2->GetCurrentNEvents(), (jMix == 0));
 887
 888             // two track cut, STEP 8
 889             if (fTwoTrackEfficiencyCut > 0)
 890               fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracks, pool2->GetEvent(jMix), 1.0 / pool2->GetCurrentNEvents(), (jMix == 0), kTRUE, bSign, fTwoTrackEfficiencyCut);
 891
 892             // apply correction efficiency, STEP 10
 893             if (fEfficiencyCorrectionTriggers || fEfficiencyCorrectionAssociated)
 894             {
 895               // with or without two track efficiency depending on if fTwoTrackEfficiencyCut is set
 896               Bool_t twoTrackCut = (fTwoTrackEfficiencyCut > 0);
 897
 898               fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepCorrected, tracks, pool2->GetEvent(jMix), 1.0 / pool2->GetCurrentNEvents(), (jMix == 0), twoTrackCut, bSign, fTwoTrackEfficiencyCut, kTRUE);
 899             }
 900           }
 901         }
 902         pool2->UpdatePool(CloneAndReduceTrackList(tracksCorrelate));
 903       }
 904
 905       if (0 && !fReduceMemoryFootprint)
 906       {
 907         // make list of secondaries (matched with MC)
 908         TObjArray* tracksRecoMatchedSecondaries = new TObjArray;
 909         for (Int_t i=0; i<tracksRecoMatchedAll->GetEntriesFast(); i++)
 910           if (((AliAODMCParticle*)tracksRecoMatchedAll->At(i))->IsPhysicalPrimary() == kFALSE)
 911             tracksRecoMatchedSecondaries->Add(tracksRecoMatchedAll->At(i));
 912
 913         // Study: Use only secondaries as trigger particles and plot the correlation vs. all particles; store in step 9
 914         fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy2, tracksRecoMatchedSecondaries, tracksRecoMatchedAll, weight);
 915
 916         // Study: Use only primaries as trigger particles and plot the correlation vs. secondaries; store in step 8
 917         fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracksRecoMatchedPrim, tracksRecoMatchedSecondaries, weight);
 918
 919         // plot delta phi vs process id of secondaries
 920         // trigger particles: primaries in 4 < pT < 10
 921         // associated particles: secondaries in 1 < pT < 10
 922
 923         for (Int_t i=0; i<tracksRecoMatchedPrim->GetEntriesFast(); i++)
 924         {
 925           AliVParticle* triggerParticle = (AliVParticle*) tracksRecoMatchedPrim->At(i);
 926
 927           if (triggerParticle->Pt() < 4 || triggerParticle->Pt() > 10)
 928             continue;
 929
 930           for (Int_t j=0; j<tracksRecoMatchedSecondaries->GetEntriesFast(); j++)
 931           {
 932             AliAODMCParticle* particle = (AliAODMCParticle*) tracksRecoMatchedSecondaries->At(j);
 933
 934             if (particle->Pt() < 1 || particle->Pt() > 10)
 935               continue;
 936
 937             if (particle->Pt() > triggerParticle->Pt())
 938               continue;
 939
 940             Double_t deltaPhi = triggerParticle->Phi() - particle->Phi();
 941             if (deltaPhi > 1.5 * TMath::Pi())
 942               deltaPhi -= TMath::TwoPi();
 943             if (deltaPhi < -0.5 * TMath::Pi())
 944               deltaPhi += TMath::TwoPi();
 945
 946             Int_t processID = fMcEvent->Stack()->Particle(particle->GetLabel())->GetUniqueID();
 947
 948             ((TH2F*) fListOfHistos->FindObject("processIDs"))->Fill(deltaPhi, processID);
 949           }
 950         }
 951
 952         delete tracksRecoMatchedSecondaries;
 953       }
 954
 955       if (tracksCorrelateRecoMatchedPrim != tracksRecoMatchedPrim)
 956         delete tracksCorrelateRecoMatchedPrim;
 957       delete tracksRecoMatchedPrim;
 958
 959       if (tracksCorrelateRecoMatchedAll != tracksRecoMatchedAll)
 960         delete tracksCorrelateRecoMatchedAll;
 961       delete tracksRecoMatchedAll;
 962
 963       if (tracksCorrelate != tracks)
 964         delete tracksCorrelate;
 965       delete tracks;
 966     }
 967   }
 968
 969   if (tracksMC != tracksCorrelateMC)
 970     delete tracksCorrelateMC;
 971   delete tracksMC;
 972 }
 973
 974 //____________________________________________________________________
 975 AliGenEventHeader* AliAnalysisTaskPhiCorrelations::GetFirstHeader()
 976 {
 977   // get first MC header from either ESD/AOD (including cocktail header if available)
 978
 979   if (fMcEvent)
 980   {
 981     // ESD
 982     AliHeader* header = (AliHeader*) fMcEvent->Header();
 983     if (!header)
 984       return 0;
 985
 986     AliGenCocktailEventHeader* cocktailHeader = dynamic_cast<AliGenCocktailEventHeader*> (header->GenEventHeader());
 987     if (cocktailHeader)
 988       return dynamic_cast<AliGenEventHeader*> (cocktailHeader->GetHeaders()->First());
 989
 990     return dynamic_cast<AliGenEventHeader*> (header->GenEventHeader());
 991   }
 992   else
 993   {
 994     // AOD
 995     AliAODMCHeader* header = (AliAODMCHeader*) fAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 996     if (!header)
 997       return 0;
 998
 999     return header->GetCocktailHeader(0);
1000   }
1001 }
1002
1003 //____________________________________________________________________
1004 void  AliAnalysisTaskPhiCorrelations::AnalyseDataMode()
1005 {
1006
1007   // Run the analysis on DATA or MC to get raw distributions
1008
1009   if ( fDebug > 3 ) AliInfo( " Processing event in Data mode ..." );
1010
1011   ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(0);
1012
1013   if (!fInputHandler)
1014     return;
1015
1016   // skip not selected events here (the AOD is not updated for those)
1017   if (!fSkipTrigger && !(fInputHandler->IsEventSelected() & fSelectBit))
1018     return;
1019
1020   // skip fast cluster events here if requested
1021   if (fSkipFastCluster && (fInputHandler->IsEventSelected() & AliVEvent::kFastOnly))
1022     return;
1023
1024   // Support for ESD and AOD based analysis
1025   AliVEvent* inputEvent = fAOD;
1026   if (!inputEvent)
1027     inputEvent = fESD;
1028
1029   Double_t centrality = 0;
1030
1031   AliCentrality *centralityObj = 0;
1032   if (fCentralityMethod.Length() > 0)
1033   {
1034     if (fCentralityMethod == "ZNA_MANUAL")
1035     {
1036       Bool_t zna = kFALSE;
1037       for(Int_t j = 0; j < 4; ++j) {
1038         if (fESD->GetZDCData()->GetZDCTDCData(12,j) != 0) {
1039           zna = kTRUE;
1040         }
1041       }
1042
1043 //       Printf("%d %f", zna, fZNAtower[0]);
1044       if (zna)
1045       {
1046         // code from Chiara O (23.10.12)
1047         const Double_t *fZNAtower = fESD->GetZDCData()->GetZN2TowerEnergy();
1048         Float_t znacut[4] = {681., 563., 413., 191.};
1049
1050         if(fZNAtower[0]>znacut[0]) centrality = 1;
1051         else if(fZNAtower[0]>znacut[1]) centrality = 21;
1052         else if(fZNAtower[0]>znacut[2]) centrality = 41;
1053         else if(fZNAtower[0]>znacut[3]) centrality = 61;
1054         else centrality = 81;
1055       }
1056       else
1057         centrality = -1;
1058     }
1059     else if (fCentralityMethod == "TRACKS_MANUAL")
1060     {
1061       // for pp
1062       TObjArray* tracks = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, -1, kTRUE);
1063       centrality = tracks->GetEntriesFast();
1064       if (centrality > 40)
1065         centrality = 41;
1066 //       Printf("%d %f", tracks->GetEntriesFast(), centrality);
1067       delete tracks;
1068     }
1069     else if (fCentralityMethod == "V0A_MANUAL")
1070     {
1071       // for pp
1072
1073       //Total multiplicity in the VZERO A detector
1074       Float_t MV0A=inputEvent->GetVZEROData()->GetMTotV0A();
1075       Float_t MV0AScaled=0.;
1076       if (fMap){
1077         TParameter<float>* sf=(TParameter<float>*)fMap->GetValue(Form("%d",inputEvent->GetRunNumber()));
1078         if(sf)MV0AScaled=MV0A*sf->GetVal();
1079       }
1080
1081       if (MV0AScaled > 0)
1082         centrality = MV0AScaled;
1083       else
1084         centrality = -1;
1085     }
1086     else if (fCentralityMethod == "nano")
1087     {
1088 //       fAOD->GetHeader()->Dump();
1089 //       Printf("%p %p %d", dynamic_cast<AliNanoAODHeader*> (fAOD->GetHeader()), dynamic_cast<AliNanoAODHeader*> ((TObject*) (fAOD->GetHeader())), fAOD->GetHeader()->InheritsFrom("AliNanoAODHeader"));
1090
1091       Int_t error = 0;
1092       centrality = (Float_t) gROOT->ProcessLine(Form("100.0 + 100.0 * ((AliNanoAODHeader*) %p)->GetCentrality(\"%s\")", fAOD->GetHeader(), fCentralityMethod.Data()), &error) / 100 - 1.0;
1093       if (error != TInterpreter::kNoError)
1094         centrality = -1;
1095     }
1096     else
1097     {
1098       if (fAOD)
1099         centralityObj = fAOD->GetHeader()->GetCentralityP();
1100       else if (fESD)
1101         centralityObj = fESD->GetCentrality();
1102
1103       if (centralityObj)
1104         centrality = centralityObj->GetCentralityPercentile(fCentralityMethod);
1105         //centrality = centralityObj->GetCentralityPercentileUnchecked(fCentralityMethod);
1106       else
1107         centrality = -1;
1108
1109       if (fAOD)
1110       {
1111         // remove outliers
1112         if (centrality == 0)
1113         {
1114           if (fAOD->GetVZEROData())
1115           {
1116             Float_t multV0 = 0;
1117             for (Int_t i=0; i<64; i++)
1118               multV0 += fAOD->GetVZEROData()->GetMultiplicity(i);
1119             if (multV0 < 19500)
1120             {
1121               centrality = -1;
1122               AliInfo("Rejecting event due to too small V0 multiplicity");
1123             }
1124           }
1125         }
1126       }
1127     }
1128
1129     AliInfo(Form("Centrality is %f", centrality));
1130   }
1131
1132   Float_t bSign = (inputEvent->GetMagneticField() > 0) ? 1 : -1;
1133
1134   fHistos->SetRunNumber(inputEvent->GetRunNumber());
1135
1136   // Fill the "event-counting-container", it is needed to get the number of events remaining after each event-selection cut
1137   fHistos->FillEvent(centrality, AliUEHist::kCFStepAll);
1138
1139   // Trigger selection ************************************************
1140   if (!fSkipTrigger && !fAnalyseUE->TriggerSelection(fInputHandler)) return;
1141
1142   // Fill the "event-counting-container", it is needed to get the number of events remaining after each event-selection cut
1143   fHistos->FillEvent(centrality, AliUEHist::kCFStepTriggered);
1144
1145   // Pileup selection ************************************************
1146   if (fAnalysisUtils && fAnalysisUtils->IsPileUpEvent(inputEvent))
1147   {
1148     // count the removed events
1149     fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
1150
1151     return;
1152   }
1153
1154   // Reject events without a muon in the muon arm ************************************************
1155   if(fAcceptOnlyMuEvents && !IsMuEvent())return;
1156
1157   // Vertex selection *************************************************
1158   if(!fAnalyseUE->VertexSelection(inputEvent, fnTracksVertex, fZVertex)) return;
1159
1160   // Fill the "event-counting-container", it is needed to get the number of events remaining after each event-selection cut
1161   fHistos->FillEvent(centrality, AliUEHist::kCFStepVertex);
1162
1163   // fill V0 control histograms
1164   if (fCentralityMethod == "V0A_MANUAL")
1165   {
1166     ((TH2F*) fListOfHistos->FindObject("V0AMult"))->Fill(inputEvent->GetVZEROData()->GetMTotV0A(), centrality);
1167     if (fAOD)
1168       ((TH2F*) fListOfHistos->FindObject("V0AMultCorrelation"))->Fill(inputEvent->GetVZEROData()->GetMTotV0A(), fAOD->GetTracklets()->GetNumberOfTracklets());
1169   }
1170
1171   // optimization
1172   if (centrality < 0)
1173     return;
1174
1175   TObjArray* tracks = 0;
1176
1177   Double_t evtPlanePhi = -999.; //A value outside [-pi/2,pi/2] will be ignored
1178   if(fTrackPhiCutEvPlMax > 0.0001) {
1179     AliEventplane* evtPlane = inputEvent->GetEventplane();
1180     Double_t qx = 0; Double_t qy = 0;
1181     if(evtPlane) evtPlanePhi = evtPlane->CalculateVZEROEventPlane(inputEvent, 10, 2, qx, qy);
1182     //Reject event if the plane is not available
1183     else return;
1184   }
1185
1186   if (fTriggersFromDetector == 0)
1187     tracks = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesTrigger, kTRUE, kTRUE, evtPlanePhi);
1188   else if (fTriggersFromDetector <= 4)
1189     tracks=GetParticlesFromDetector(inputEvent,fTriggersFromDetector);
1190   else
1191     AliFatal(Form("Invalid setting for fTriggersFromDetector: %d", fTriggersFromDetector));
1192
1193   //Printf("Accepted %d tracks", tracks->GetEntries());
1194
1195   // check for outlier in centrality vs number of tracks (rough constants extracted from correlation histgram)
1196   Bool_t reject = kFALSE;
1197   if (fRejectCentralityOutliers)
1198   {
1199     if (centrality > 40 && centrality <= 50 && tracks->GetEntriesFast() > 1160)
1200       reject = kTRUE;
1201     if (centrality > 50 && centrality <= 60 && tracks->GetEntriesFast() > 650)
1202       reject = kTRUE;
1203     if (centrality > 60 && centrality <= 70 && tracks->GetEntriesFast() > 370)
1204       reject = kTRUE;
1205     if (centrality > 70 && centrality <= 80 && tracks->GetEntriesFast() > 220)
1206       reject = kTRUE;
1207     if (centrality > 80 && centrality <= 90 && tracks->GetEntriesFast() > 130)
1208       reject = kTRUE;
1209     if (centrality > 90 && tracks->GetEntriesFast() > 75)
1210       reject = kTRUE;
1211   }
1212
1213   if (reject)
1214   {
1215     AliInfo(Form("Rejecting event due to centrality vs tracks correlation: %f %d", centrality, tracks->GetEntriesFast()));
1216     fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
1217     delete tracks;
1218     return;
1219   }
1220
1221   if (fRejectZeroTrackEvents && tracks->GetEntriesFast() == 0)
1222   {
1223     AliInfo(Form("Rejecting event because it has no tracks: %f %d", centrality, tracks->GetEntriesFast()));
1224     fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
1225     delete tracks;
1226     return;
1227   }
1228
1229   if (fCentralityWeights && !AcceptEventCentralityWeight(centrality))
1230   {
1231     AliInfo(Form("Rejecting event because of centrality weighting: %f", centrality));
1232     fHistos->FillEvent(centrality, AliUEHist::kCFStepAnaTopology);
1233     delete tracks;
1234     return;
1235   }
1236
1237   // correlate particles with...
1238   TObjArray* tracksCorrelate = 0;
1239   if(fAssociatedFromDetector==0){
1240     if (fParticleSpeciesAssociated != fParticleSpeciesTrigger || fTriggersFromDetector > 0 || fTrackPhiCutEvPlMax > 0.0001)
1241       tracksCorrelate = fAnalyseUE->GetAcceptedParticles(inputEvent, 0, kTRUE, fParticleSpeciesAssociated, kTRUE);
1242   }
1243   else if (fAssociatedFromDetector <= 4){
1244     if(fAssociatedFromDetector != fTriggersFromDetector)
1245       tracksCorrelate = GetParticlesFromDetector(inputEvent,fAssociatedFromDetector);
1246   }
1247   else
1248     AliFatal(Form("Invalid setting for fAssociatedFromDetector: %d", fAssociatedFromDetector));
1249
1250   // reference multiplicity
1251   Int_t referenceMultiplicity = -1;
1252   if (fESD)
1253     referenceMultiplicity = AliESDtrackCuts::GetReferenceMultiplicity(fESD);
1254   else if (fAOD)
1255     referenceMultiplicity = tracks->GetEntriesFast(); // TODO to be replaced by the estimator once available in the AOD
1256 //    referenceMultiplicity = fAOD->GetHeader()->GetRefMultiplicityComb05();
1257
1258   ((TH2F*) fListOfHistos->FindObject("referenceMultiplicity"))->Fill(centrality, referenceMultiplicity);
1259
1260   const AliVVertex* vertex = inputEvent->GetPrimaryVertex();
1261   Double_t zVtx = vertex->GetZ();
1262
1263   Float_t weight = 1;
1264   if (fFillpT)
1265     weight = -1;
1266
1267   // Fill containers at STEP 6 (reconstructed)
1268   if (centrality >= 0)
1269   {
1270     if (!fSkipStep6)
1271       fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracks, tracksCorrelate, weight, kTRUE, kFALSE, 0, 0.02, kTRUE);
1272
1273     ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(1);
1274
1275     if (fTwoTrackEfficiencyCut > 0)
1276       fHistos->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracks, tracksCorrelate, weight, kTRUE, kTRUE, bSign, fTwoTrackEfficiencyCut, kTRUE);
1277   }
1278
1279   // create a list of reduced objects. This speeds up processing and reduces memory consumption for the event pool
1280   TObjArray* tracksClone = CloneAndReduceTrackList(tracks);
1281   delete tracks;
1282
1283   if (fFillMixed)
1284   {
1285     // event mixing
1286
1287     // 1. First get an event pool corresponding in mult (cent) and
1288     //    zvertex to the current event. Once initialized, the pool
1289     //    should contain nMix (reduced) events. This routine does not
1290     //    pre-scan the chain. The first several events of every chain
1291     //    will be skipped until the needed pools are filled to the
1292     //    specified depth. If the pool categories are not too rare, this
1293     //    should not be a problem. If they are rare, you could lose
1294     //    statistics.
1295
1296     // 2. Collect the whole pool's content of tracks into one TObjArray
1297     //    (bgTracks), which is effectively a single background super-event.
1298
1299     // 3. The reduced and bgTracks arrays must both be passed into
1300     //    FillCorrelations(). Also nMix should be passed in, so a weight
1301     //    of 1./nMix can be applied.
1302
1303     AliEventPool* pool = fPoolMgr->GetEventPool(centrality, zVtx);
1304
1305     if (!pool)
1306       AliFatal(Form("No pool found for centrality = %f, zVtx = %f", centrality, zVtx));
1307
1308 //     pool->SetDebug(1);
1309
1310     if (pool->IsReady())
1311     {
1312
1313       Int_t nMix = pool->GetCurrentNEvents();
1314 //       cout << "nMix = " << nMix << " tracks in pool = " << pool->NTracksInPool() << endl;
1315
1316       ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(2);
1317       ((TH2F*) fListOfHistos->FindObject("mixedDist"))->Fill(centrality, pool->NTracksInPool());
1318       if (pool->IsReady())
1319         ((TH1F*) fListOfHistos->FindObject("eventStat"))->Fill(3);
1320
1321       // Fill mixed-event histos here
1322       for (Int_t jMix=0; jMix<nMix; jMix++)
1323       {
1324         TObjArray* bgTracks = pool->GetEvent(jMix);
1325
1326         if (!fSkipStep6)
1327           fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepReconstructed, tracksClone, bgTracks, 1.0 / nMix, (jMix == 0), kFALSE, 0, 0.02, kTRUE);
1328
1329         if (fTwoTrackEfficiencyCut > 0)
1330           fHistosMixed->FillCorrelations(centrality, zVtx, AliUEHist::kCFStepBiasStudy, tracksClone, bgTracks, 1.0 / nMix, (jMix == 0), kTRUE, bSign, fTwoTrackEfficiencyCut, kTRUE);
1331       }
1332     }
1333
1334     // ownership is with the pool now
1335     if (tracksCorrelate)
1336     {
1337       pool->UpdatePool(CloneAndReduceTrackList(tracksCorrelate));
1338       delete tracksClone;
1339     }
1340     else
1341       pool->UpdatePool(tracksClone);
1342     //pool->PrintInfo();
1343   }
1344   else
1345   {
1346     delete tracksClone;
1347     if (tracksCorrelate)
1348       delete tracksCorrelate;
1349   }
1350 }
1351
1352 TObjArray* AliAnalysisTaskPhiCorrelations::CloneAndReduceTrackList(TObjArray* tracks)
1353 {
1354   // clones a track list by using AliDPhiBasicParticle which uses much less memory (used for event mixing)
1355
1356   TObjArray* tracksClone = new TObjArray;
1357   tracksClone->SetOwner(kTRUE);
1358
1359   for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
1360   {
1361     AliVParticle* particle = (AliVParticle*) tracks->UncheckedAt(i);
1362     AliDPhiBasicParticle* copy = new AliDPhiBasicParticle(particle->Eta(), particle->Phi(), particle->Pt(), particle->Charge());
1363     copy->SetUniqueID(particle->GetUniqueID());
1364     tracksClone->Add(copy);
1365   }
1366
1367   return tracksClone;
1368 }
1369
1370 //____________________________________________________________________
1371 void  AliAnalysisTaskPhiCorrelations::Initialize()
1372 {
1373   // input handler
1374   fInputHandler = (AliInputEventHandler*)
1375          ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
1376   // MC handler
1377   fMcHandler = dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
1378 }
1379
1380 //____________________________________________________________________
1381 void AliAnalysisTaskPhiCorrelations::RemoveDuplicates(TObjArray* tracks)
1382 {
1383   // remove particles with the same label
1384
1385   Int_t before = tracks->GetEntriesFast();
1386
1387   for (Int_t i=0; i<before; ++i)
1388   {
1389     AliVParticle* part = (AliVParticle*) tracks->At(i);
1390
1391     for (Int_t j=i+1; j<before; ++j)
1392     {
1393       AliVParticle* part2 = (AliVParticle*) tracks->At(j);
1394
1395       if (part->GetLabel() == part2->GetLabel())
1396       {
1397         Printf("Removing %d with label %d (duplicated in %d)", i, part->GetLabel(), j); part->Dump(); part2->Dump();
1398         TObject* object = tracks->RemoveAt(i);
1399         if (tracks->IsOwner())
1400           delete object;
1401         break;
1402       }
1403     }
1404   }
1405
1406   tracks->Compress();
1407
1408   if (before > tracks->GetEntriesFast())
1409     AliInfo(Form("Reduced from %d to %d", before, tracks->GetEntriesFast()));
1410 }
1411
1412 void AliAnalysisTaskPhiCorrelations::CleanUp(TObjArray* tracks, TObject* mcObj, Int_t maxLabel)
1413 {
1414   // calls RemoveInjectedSignals, RemoveWeakDecays and RemoveDuplicates
1415
1416   if (!tracks)
1417     return;
1418
1419   if (fInjectedSignals)
1420     fAnalyseUE->RemoveInjectedSignals(tracks, mcObj, maxLabel);
1421   if (fRemoveWeakDecays)
1422     fAnalyseUE->RemoveWeakDecays(tracks, mcObj);
1423   if (fRemoveDuplicates)
1424     RemoveDuplicates(tracks);
1425 }
1426
1427 //____________________________________________________________________
1428 void AliAnalysisTaskPhiCorrelations::SelectCharge(TObjArray* tracks)
1429 {
1430   // remove particles with charge not selected (depending on fTriggerSelectCharge)
1431
1432   if (!tracks)
1433     return;
1434
1435   Int_t before = tracks->GetEntriesFast();
1436
1437   for (Int_t i=0; i<before; ++i)
1438   {
1439     AliVParticle* part = (AliVParticle*) tracks->At(i);
1440
1441     if (part->Charge() * fTriggerSelectCharge < -1)
1442     {
1443 //       Printf("Removing %d with charge %d", i, part->Charge());
1444       TObject* object = tracks->RemoveAt(i);
1445       if (tracks->IsOwner())
1446         delete object;
1447     }
1448   }
1449
1450   tracks->Compress();
1451
1452   if (before > tracks->GetEntriesFast())
1453     AliInfo(Form("Reduced from %d to %d", before, tracks->GetEntriesFast()));
1454 }
1455
1456 //____________________________________________________________________
1457 Bool_t AliAnalysisTaskPhiCorrelations::AcceptEventCentralityWeight(Double_t centrality)
1458 {
1459   // rejects "randomly" events such that the centrality gets flat
1460   // uses fCentralityWeights histogram
1461
1462   // TODO code taken and adapted from AliRDHFCuts; waiting for general class AliCentralityFlattening
1463
1464   Double_t weight = fCentralityWeights->GetBinContent(fCentralityWeights->FindBin(centrality));
1465   Double_t centralityDigits = centrality*100. - (Int_t)(centrality*100.);
1466
1467   Bool_t result = kFALSE;
1468   if (centralityDigits < weight)
1469     result = kTRUE;
1470
1471   AliInfo(Form("Centrality: %f; Digits: %f; Weight: %f; Result: %d", centrality, centralityDigits, weight, result));
1472
1473   return result;
1474 }
1475
1476 //____________________________________________________________________
1477 void AliAnalysisTaskPhiCorrelations::ShiftTracks(TObjArray* tracks, Double_t angle)
1478 {
1479   // shifts the phi angle of all tracks by angle
1480   // 0 <= angle <= 2pi
1481
1482   for (Int_t i=0; i<tracks->GetEntriesFast(); ++i)
1483   {
1484     AliDPhiBasicParticle* part = (AliDPhiBasicParticle*) tracks->At(i);
1485     Double_t newAngle = part->Phi() + angle;
1486     if (newAngle >= TMath::TwoPi())
1487       newAngle -= TMath::TwoPi();
1488
1489     part->SetPhi(newAngle);
1490   }
1491 }
1492
1493 //____________________________________________________________________
1494 TObjArray* AliAnalysisTaskPhiCorrelations::GetParticlesFromDetector(AliVEvent* inputEvent, Int_t idet)
1495 {
1496   //1 = VZERO_A; 2 = VZERO_C; 3 = SPD tracklets
1497   TObjArray* obj = new TObjArray;
1498   obj->SetOwner(kTRUE);
1499
1500   if (idet == 1 || idet == 2)
1501   {
1502     AliVVZERO* vZero = inputEvent->GetVZEROData();
1503
1504     const Int_t vZeroStart = (idet == 1) ? 32 : 0;
1505
1506     TH1F* singleCells = (TH1F*) fListOfHistos->FindObject("V0SingleCells");
1507     for (Int_t i=vZeroStart; i<vZeroStart+32; i++)
1508       {
1509         Float_t weight = vZero->GetMultiplicity(i);
1510         singleCells->Fill(weight);
1511
1512         // rough estimate of multiplicity
1513         for (Int_t j=0; j<TMath::Nint(weight); j++)
1514           {
1515             AliDPhiBasicParticle* particle = new AliDPhiBasicParticle((AliVVZERO::GetVZEROEtaMax(i) + AliVVZERO::GetVZEROEtaMin(i)) / 2, AliVVZERO::GetVZEROAvgPhi(i), 1.1, 0); // fit pT = 1.1 and charge = 0
1516             particle->SetUniqueID(fAnalyseUE->GetEventCounter()* 100000 + j + i * 1000);
1517
1518             obj->Add(particle);
1519           }
1520       }
1521   }
1522   else if (idet == 3)
1523     {
1524       if(!fAOD)
1525         AliFatal("Tracklets only available on AOD");
1526       AliAODTracklets* trklets=(AliAODTracklets*)fAOD->GetTracklets();
1527       if(!trklets)
1528         AliFatal("AliAODTracklets not found");
1529       for(Int_t itrklets=0;itrklets<trklets->GetNumberOfTracklets();itrklets++)
1530         {
1531           Float_t eta=-TMath::Log(TMath::Tan(trklets->GetTheta(itrklets)/2));
1532           if(TMath::Abs(eta)>fTrackEtaCut)continue;
1533           Float_t pT=1000*TMath::Abs(trklets->GetDeltaPhi(itrklets));//in mrad
1534           if(pT>fTrackletDphiCut)continue;
1535           TH1F* DphiTrklets = (TH1F*)fListOfHistos->FindObject("DphiTrklets");
1536           DphiTrklets->Fill(1000*trklets->GetDeltaPhi(itrklets)); //in mrad
1537           Float_t phi=trklets->GetPhi(itrklets);
1538           phi+=trklets->GetDeltaPhi(itrklets)*39./34.; //correction dphi*39./34. (Dphi in rad)
1539           if (phi<0) phi+=TMath::TwoPi();
1540           if (phi>TMath::TwoPi()) phi-=TMath::TwoPi();
1541
1542           AliDPhiBasicParticle* particle = new AliDPhiBasicParticle(eta,phi, pT, 0); // pT = TMath::Abs(trklets->GetDeltaPhi(itrklets)) in mrad and charge = 0
1543           particle->SetUniqueID(fAnalyseUE->GetEventCounter()* 100000 + itrklets);
1544
1545           obj->Add(particle);
1546         }
1547     }
1548   else if (idet == 4)
1549     {
1550       if(!fAOD)
1551         AliFatal("Muon selection only implemented on AOD");//FIXME to be implemented also for ESDs as in AliAnalyseLeadingTrackUE::GetAcceptedPArticles
1552       for (Int_t iTrack = 0; iTrack < fAOD->GetNTracks(); iTrack++) {
1553         AliAODTrack* track = fAOD->GetTrack(iTrack);
1554         if (!track->IsMuonTrack()) continue;
1555         //Float_t dca    = track->DCA();
1556         //Float_t chi2   = track->Chi2perNDF();
1557         Float_t rabs   = track->GetRAtAbsorberEnd();
1558         Float_t eta    = track->Eta();
1559         Int_t   matching   = track->GetMatchTrigger();
1560         if (rabs < 17.6 || rabs > 89.5) continue;
1561         if (eta < -4 || eta > -2.5) continue;
1562         if (matching < 2) continue;
1563
1564         AliDPhiBasicParticle* particle = new AliDPhiBasicParticle(eta,track->Phi(), track->Pt(), track->Charge());
1565         particle->SetUniqueID(fAnalyseUE->GetEventCounter() * 100000 + iTrack);
1566
1567         obj->Add(particle);
1568       }
1569     }
1570   else
1571     AliFatal(Form("GetParticlesFromDetector: Invalid idet value: %d", idet));
1572
1573   return obj;
1574 }
1575
1576 //____________________________________________________________________
1577 Bool_t AliAnalysisTaskPhiCorrelations::IsMuEvent(){
1578
1579   if(!fAOD)
1580     AliFatal("Muon selection only implemented on AOD");//FIXME to be implemented also for ESDs as in AliAnalyseLeadingTrackUE::GetAcceptedPArticles
1581   for (Int_t iTrack = 0; iTrack < fAOD->GetNTracks(); iTrack++) {
1582     AliAODTrack* track = fAOD->GetTrack(iTrack);
1583     if (!track->IsMuonTrack()) continue;
1584     //Float_t dca    = track->DCA();
1585     //Float_t chi2   = track->Chi2perNDF();
1586     Float_t rabs   = track->GetRAtAbsorberEnd();
1587     Float_t eta    = track->Eta();
1588     Int_t   matching   = track->GetMatchTrigger();
1589     if (rabs < 17.6 || rabs > 89.5) continue;
1590     if (eta < -4 || eta > -2.5) continue;
1591     if (matching < 2) continue;
1592     return kTRUE;
1593   }
1594   return kFALSE;
1595
1596 }