PWGGA/PHOSTasks/PHOS_Correlation/AliPHOSCorrelations.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2014, 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 // Analysis task for identified PHOS cluster from pi0 and extracting pi0-hadron correlation.
  17 // Authors:     Daniil Ponomarenko <Daniil.Ponomarenko@cern.ch>
  18 //              Dmitry Blau <Dmitry.Blau@cern.ch>
  19 // 09-Jul-2014
  20
  21 #include <Riostream.h>
  22 #include "THashList.h"
  23 #include "TObjArray.h"
  24 #include "TClonesArray.h"
  25 #include "TH1F.h"
  26 #include "TH2F.h"
  27 #include "TH2I.h"
  28 #include "TH3F.h"
  29 #include "TH3D.h"
  30 #include "TMath.h"
  31 #include "TVector3.h"
  32 #include "TProfile.h"
  33
  34 #include "AliAnalysisManager.h"
  35 #include "AliAnalysisTaskSE.h"
  36 #include "AliPHOSCorrelations.h"
  37 #include "AliPHOSGeometry.h"
  38 #include "AliESDEvent.h"
  39 #include "AliESDCaloCluster.h"
  40 #include "AliESDVertex.h"
  41 #include "AliESDtrackCuts.h"
  42 #include "AliESDtrack.h"
  43 #include "AliAODTrack.h"
  44 #include "AliVTrack.h"
  45 #include "AliPID.h"
  46 #include "AliTriggerAnalysis.h"
  47 #include "AliPIDResponse.h"
  48 #include "AliPHOSEsdCluster.h"
  49 #include "AliCDBManager.h"
  50 #include "AliPHOSCalibData.h"
  51 #include "AliCentrality.h"
  52 #include "AliAnalysisTaskSE.h"
  53 #include "AliEventplane.h"
  54 #include "AliOADBContainer.h"
  55 #include "AliAODEvent.h"
  56 #include "AliAODCaloCells.h"
  57 #include "AliAODCaloCluster.h"
  58 #include "AliCaloPhoton.h"
  59 #include "AliAODVertex.h"
  60 #include "AliInputEventHandler.h"
  61
  62 using std::cout;
  63 using std::endl;
  64
  65 ClassImp(AliPHOSCorrelations)
  66
  67 //_______________________________________________________________________________
  68 AliPHOSCorrelations::AliPHOSCorrelations()
  69 :AliAnalysisTaskSE(),
  70         fPHOSGeo(0x0),
  71         fOutputContainer(0x0),
  72         fMinClusterEnergy(0.3),
  73         fMinBCDistance(0),
  74         fMinNCells(3),
  75         fMinM02(0.2),
  76         fTOFCutEnabled(1),
  77         fTOFCut(100.e-9),
  78         fNVtxZBins(1),
  79         fCentEdges(10),
  80         fCentNMixed(),
  81         fNEMRPBins(9),
  82         fAssocBins(),
  83     fUseMEAlgoritmForReal(true),
  84     fUseMEAlgoritmForMix(true),
  85         fCheckHibridGlobal(kOnlyHibridTracks),
  86         fPHOSEvent(false),
  87         fMBEvent(false),
  88         fPeriod(kUndefinedPeriod),
  89         fMaxAbsVertexZ(10.),
  90         fManualV0EPCalc(false),
  91         fCentCutoffDown(0.),
  92         fCentCutoffUp(90),
  93         fMassInvMean(0.135),
  94         fMassInvSigma(0.05),
  95         fSigmaWidth(0.),
  96         fMassMeanP0(-20.9476),
  97         fMassMeanP1(0.1300),
  98         fMassSigmaP0(0.005),
  99         fMassSigmaP1(-0.0001),
 100         fEvent(0x0),
 101         fEventESD(0x0),
 102         fEventAOD(0x0),
 103         fEventHandler(0),
 104         fESDtrackCuts(0x0),
 105         fRunNumber(-999),
 106         fInternalRunNumber(0),
 107         fMultV0(0x0),
 108         fV0Cpol(0.),fV0Apol(0.),
 109         fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
 110         fVertexVector(),
 111         fVtxBin(0),
 112         fCentralityEstimator("V0M"),
 113         fCentrality(0.),
 114         fCentBin(0),
 115         fHaveTPCRP(0),
 116         fRP(0.),
 117         fEMRPBin(0),
 118         fCaloPhotonsPHOS(0x0),
 119         fTracksTPC(0x0),
 120         fCaloPhotonsPHOSLists(0x0),
 121         fTracksTPCLists(0x0)
 122 {
 123   //Deafult constructor, no memory allocations here
 124 }
 125
 126 //_______________________________________________________________________________
 127 AliPHOSCorrelations::AliPHOSCorrelations(const char *name)
 128 :AliAnalysisTaskSE(name),
 129     fPHOSGeo(0x0),
 130     fOutputContainer(0x0),
 131     fMinClusterEnergy(0.3),
 132     fMinBCDistance(0),
 133     fMinNCells(3),
 134     fMinM02(0.2),
 135     fTOFCutEnabled(1),
 136     fTOFCut(100.e-9),
 137     fNVtxZBins(1),
 138     fCentEdges(10),
 139     fCentNMixed(),
 140     fNEMRPBins(9),
 141     fAssocBins(),
 142     fUseMEAlgoritmForReal(true),
 143     fUseMEAlgoritmForMix(true),
 144     fCheckHibridGlobal(kOnlyHibridTracks),
 145     fPHOSEvent(false),
 146     fMBEvent(false),
 147     fPeriod(kUndefinedPeriod),
 148     fMaxAbsVertexZ(10.),
 149     fManualV0EPCalc(false),
 150     fCentCutoffDown(0.),
 151     fCentCutoffUp(90),
 152     fMassInvMean(0.135),
 153     fMassInvSigma(0.05),
 154     fSigmaWidth(0.),
 155     fMassMeanP0(-20.9476),
 156     fMassMeanP1(0.1300),
 157     fMassSigmaP0(0.005),
 158     fMassSigmaP1(-0.0001),
 159     fEvent(0x0),
 160     fEventESD(0x0),
 161     fEventAOD(0x0),
 162     fEventHandler(0),
 163     fESDtrackCuts(0x0),
 164     fRunNumber(-999),
 165     fInternalRunNumber(0),
 166     fMultV0(0x0),
 167     fV0Cpol(0.),fV0Apol(0.),
 168     fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
 169     fVertexVector(),
 170     fVtxBin(0),
 171     fCentralityEstimator("V0M"),
 172     fCentrality(0.),
 173     fCentBin(0),
 174     fHaveTPCRP(0),
 175     fRP(0.),
 176     fEMRPBin(0),
 177     fCaloPhotonsPHOS(0x0),
 178     fTracksTPC(0x0),
 179     fCaloPhotonsPHOSLists(0x0),
 180     fTracksTPCLists(0x0)
 181 {
 182     // Constructor
 183     // Output slots #0 write into a TH1 container
 184     DefineOutput(1,THashList::Class());
 185
 186     const Int_t nPtAssoc=10 ;
 187     Double_t ptAssocBins[nPtAssoc]={0.,0.5,1.0,1.5,2.0,3.,5.,7.,10.,16} ;
 188     fAssocBins.Set(nPtAssoc,ptAssocBins) ;
 189
 190     const int nbins = 9;
 191     Double_t edges[nbins+1] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80.};
 192     TArrayD centEdges(nbins+1, edges);
 193     Int_t nMixed[nbins] = {4,4,6,10,20,30,50,100,100};
 194     //Int_t nMixed[nbins] = {100,100,100,100,100,100,100,100,100};
 195     TArrayI centNMixed(nbins, nMixed);
 196     SetCentralityBinning(centEdges, centNMixed);
 197
 198     fVertex[0]=0; fVertex[1]=0; fVertex[2]=0;
 199
 200     SetGeometry();
 201
 202     ZeroingVariables();
 203 }
 204
 205 //_______________________________________________________________________________
 206 AliPHOSCorrelations::AliPHOSCorrelations(const char *name, Period period)
 207 :AliAnalysisTaskSE(name),
 208         fPHOSGeo(0x0),
 209         fOutputContainer(0x0),
 210         fMinClusterEnergy(0.3),
 211         fMinBCDistance(0),
 212         fMinNCells(3),
 213         fMinM02(0.2),
 214         fTOFCutEnabled(1),
 215         fTOFCut(100.e-9),
 216         fNVtxZBins(1),
 217         fCentEdges(10),
 218         fCentNMixed(),
 219         fNEMRPBins(9),
 220         fAssocBins(),
 221     fUseMEAlgoritmForReal(true),
 222     fUseMEAlgoritmForMix(true),
 223         fCheckHibridGlobal(kOnlyHibridTracks),
 224         fPHOSEvent(false),
 225         fMBEvent(false),
 226         fPeriod(period),
 227         fMaxAbsVertexZ(10.),
 228         fManualV0EPCalc(false),
 229         fCentCutoffDown(0.),
 230         fCentCutoffUp(90),
 231         fMassInvMean(0.135),
 232         fMassInvSigma(0.05),
 233         fSigmaWidth(0.),
 234         fMassMeanP0(-20.9476),
 235         fMassMeanP1(0.1300),
 236         fMassSigmaP0(0.005),
 237         fMassSigmaP1(-0.0001),
 238         fEvent(0x0),
 239         fEventESD(0x0),
 240         fEventAOD(0x0),
 241         fEventHandler(0),
 242         fESDtrackCuts(0x0),
 243         fRunNumber(-999),
 244         fInternalRunNumber(0),
 245         fMultV0(0x0),
 246         fV0Cpol(0.),fV0Apol(0.),
 247         fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
 248         fVertexVector(),
 249         fVtxBin(0),
 250         fCentralityEstimator("V0M"),
 251         fCentrality(0.),
 252         fCentBin(0),
 253         fHaveTPCRP(0),
 254         fRP(0.),
 255         fEMRPBin(0),
 256         fCaloPhotonsPHOS(0x0),
 257         fTracksTPC(0x0),
 258         fCaloPhotonsPHOSLists(0x0),
 259         fTracksTPCLists(0x0)
 260 {
 261         // Constructor
 262         // Output slots #0 write into a TH1 container
 263         DefineOutput(1,THashList::Class());
 264
 265         const Int_t nPtAssoc=10 ;
 266         Double_t ptAssocBins[nPtAssoc]={0.,0.5,1.0,1.5,2.0,3.,5.,7.,10.,16} ;
 267         fAssocBins.Set(nPtAssoc,ptAssocBins) ;
 268
 269         const int nbins = 9;
 270         Double_t edges[nbins+1] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80.};
 271         TArrayD centEdges(nbins+1, edges);
 272         Int_t nMixed[nbins] = {4,4,6,10,20,30,50,100,100};
 273         //Int_t nMixed[nbins] = {100,100,100,100,100,100,100,100,100};
 274         TArrayI centNMixed(nbins, nMixed);
 275         SetCentralityBinning(centEdges, centNMixed);
 276
 277         fVertex[0]=0; fVertex[1]=0; fVertex[2]=0;
 278
 279         SetGeometry();
 280
 281         ZeroingVariables();
 282 }
 283 //_______________________________________________________________________________
 284 AliPHOSCorrelations::~AliPHOSCorrelations()
 285 {
 286         if(fCaloPhotonsPHOS){
 287           delete fCaloPhotonsPHOS;
 288           fCaloPhotonsPHOS=0x0;
 289         }
 290
 291         if(fTracksTPC){
 292           delete fTracksTPC;
 293           fTracksTPC=0x0;
 294         }
 295
 296         if(fCaloPhotonsPHOSLists){
 297           fCaloPhotonsPHOSLists->SetOwner() ;
 298           delete fCaloPhotonsPHOSLists;
 299           fCaloPhotonsPHOSLists=0x0;
 300         }
 301
 302         if(fTracksTPCLists){
 303           fTracksTPCLists->SetOwner() ;
 304           delete fTracksTPCLists;
 305           fTracksTPCLists=0x0 ;
 306         }
 307
 308         if( fESDtrackCuts){
 309           delete fESDtrackCuts;
 310           fESDtrackCuts=0x0 ;
 311         }
 312
 313         if(fOutputContainer){
 314           delete fOutputContainer;
 315           fOutputContainer=0x0;
 316         }
 317 }
 318 //_______________________________________________________________________________
 319 void AliPHOSCorrelations::UserCreateOutputObjects()
 320 {
 321         // Create histograms
 322         // Called once
 323         const Int_t nRuns=200 ;
 324         const Int_t  ptMult = 200;
 325         const Double_t ptMin = 0.;
 326         const Double_t ptMax = 20.;
 327
 328         // Create histograms
 329         if(fOutputContainer != NULL) { delete fOutputContainer; }
 330         fOutputContainer = new THashList();
 331         fOutputContainer->SetOwner(kTRUE);
 332
 333     // Event selection
 334     fOutputContainer->Add(new TH1F("hTriggerPassedEvents","Event selection passed Cuts", 20, 0., 20.) );
 335     // Analysis event's progress
 336     fOutputContainer->Add(new TH1F("hTotSelEvents","Event selection", 15, 0., 15)) ;
 337         fOutputContainer->Add(new TH2F("hSelEvents","Event selection", kTotalSelected+1, 0., double(kTotalSelected+1), nRuns,0.,float(nRuns))) ;
 338     // Centrality, Reaction plane selection
 339     fOutputContainer->Add(new TH2F("hCentrality","Event centrality", 100,0.,100.,nRuns,0.,float(nRuns))) ;
 340         fOutputContainer->Add(new TH2F("phiRPflat","RP distribution with TPC flat", 100, 0., 2.*TMath::Pi(),20,0.,100.)) ;
 341     // Mass selection
 342         fOutputContainer->Add(new TH2F("massWindow","mean & sigma", 100,0.095,0.185,500,0.,0.05));
 343     fOutputContainer->Add(new TH1F("massWindowPass","Mass selection", 10, 0., 10.)) ;
 344     // Cluster multiplisity
 345         fOutputContainer->Add(new TH2F("hCluEvsClu","ClusterMult vs E",200,0.,10.,100,0.,100.)) ;
 346
 347
 348         // Set hists, with track's and cluster's angle distributions.
 349         SetHistPtNumTrigger(ptMult, ptMin, ptMax);
 350         SetHistEtaPhi();
 351         SetHistPHOSClusterMap();
 352         SetHistMass(ptMult, ptMin, ptMax);
 353         SetHistPtAssoc(ptMult, ptMin, ptMax);
 354
 355         // Setup photon lists
 356         Int_t kapacity = fNVtxZBins * GetNumberOfCentralityBins() * fNEMRPBins;
 357         fCaloPhotonsPHOSLists = new TObjArray(kapacity);
 358         fCaloPhotonsPHOSLists->SetOwner();
 359
 360         fTracksTPCLists = new TObjArray(kapacity);
 361         fTracksTPCLists->SetOwner();
 362
 363         PostData(1, fOutputContainer);
 364 }
 365 //_______________________________________________________________________________
 366 void AliPHOSCorrelations::SetHistPtNumTrigger(Int_t  ptMult, Double_t ptMin, Double_t ptMax)
 367 {
 368         TString spid[4]={"all","cpv","disp","both"} ;
 369         for(Int_t ipid=0; ipid<4; ipid++)
 370         {
 371                 fOutputContainer->Add(new TH1F(Form("nTrigger_%s", spid[ipid].Data()), Form("Num of trigger particle %s", spid[ipid].Data()), ptMult+300, ptMin, ptMax ) );
 372                 TH1F *h = static_cast<TH1F*>(fOutputContainer->Last()) ;
 373                 h->GetXaxis()->SetTitle("Pt [GEV]");
 374                 //h->GetYaxis()->SetTitle("#varepsilon"); // 1/efficiensy
 375         }
 376 }
 377 //_______________________________________________________________________________
 378 void AliPHOSCorrelations::SetHistEtaPhi()
 379 {
 380         // Set hists, with track's and cluster's angle distributions.
 381
 382         Float_t pi = TMath::Pi();
 383
 384         //===
 385         fOutputContainer->Add(new TH2F("clu_phieta","Cluster's #phi & #eta distribution", 300, double(-1.8), double(-0.6), 300, double(-0.2), double(0.2) ) );
 386         TH2F * h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 387         h->GetXaxis()->SetTitle("#phi [rad]");
 388         h->GetYaxis()->SetTitle("#eta");
 389
 390         //===
 391         fOutputContainer->Add(new TH2F("clusingle_phieta","Cluster's  #phi & #eta distribution", 300, double(-1.8), double(-0.6), 300, double(-0.2), double(0.2) ) );
 392         h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 393         h->GetXaxis()->SetTitle("#phi [rad]");
 394         h->GetYaxis()->SetTitle("#eta");
 395
 396         //===
 397         fOutputContainer->Add(new TH2F("track_phieta","TPC track's  #phi & #eta distribution", 200, double(-pi-0.3), double(pi+0.3), 200, double(-0.9), double(0.9) ) );
 398         h = static_cast<TH2F*>(fOutputContainer->FindObject("track_phieta")) ;
 399         h->GetXaxis()->SetTitle("#phi [rad]");
 400         h->GetYaxis()->SetTitle("#eta");
 401 }
 402 //_______________________________________________________________________________
 403 void AliPHOSCorrelations::SetHistMass(Int_t  ptMult, Double_t ptMin, Double_t ptMax)
 404 {
 405         // Set other histograms.
 406         // cout<<"\nSetting output SetHist_CutDistribution...";
 407
 408         Double_t massMin = fMassInvMean-fMassInvSigma;
 409         Double_t massMax = fMassInvMean+fMassInvSigma;
 410
 411         TString spid[4]={"all","cpv","disp","both"} ;
 412
 413         TH2F * h;
 414
 415         for(Int_t ipid=0; ipid<4; ipid++)
 416         {
 417                 // Real ++++++++++++++++++++++++++++++
 418
 419                 fOutputContainer->Add(new TH2F(Form("%s_mpt",spid[ipid].Data() )," real ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 420                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 421                 h->GetXaxis()->SetTitle("Mass [GeV]");
 422                 h->GetYaxis()->SetTitle("Pt [GEV]");
 423
 424                 // MIX +++++++++++++++++++++++++
 425
 426                 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt",spid[ipid].Data() )," mix ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 427                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 428                 h->GetXaxis()->SetTitle("Mass [GeV]");
 429                 h->GetYaxis()->SetTitle("Pt [GEV]");
 430
 431                 // Real ++++++++++++++++++++++++++++++
 432
 433                 fOutputContainer->Add(new TH2F(Form("%s_mpt_left",spid[ipid].Data() )," real ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
 434                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 435                 h->GetXaxis()->SetTitle("Mass [GeV]");
 436                 h->GetYaxis()->SetTitle("Pt [GEV]");
 437
 438                 fOutputContainer->Add(new TH2F(Form("%s_mpt_right",spid[ipid].Data() )," real ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
 439                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 440                 h->GetXaxis()->SetTitle("Mass [GeV]");
 441                 h->GetYaxis()->SetTitle("Pt [GEV]");
 442
 443                 // MIX +++++++++++++++++++++++++
 444
 445                 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_left",spid[ipid].Data() )," mix ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
 446                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 447                 h->GetXaxis()->SetTitle("Mass [GeV]");
 448                 h->GetYaxis()->SetTitle("Pt [GEV]");
 449
 450                 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_right",spid[ipid].Data() )," mix ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
 451                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 452                 h->GetXaxis()->SetTitle("Mass [GeV]");
 453                 h->GetYaxis()->SetTitle("Pt [GEV]");
 454         }
 455
 456         // Calibration PHOS Module Pi0peak {REAL}
 457         for(Int_t mod=1; mod<4; mod++){
 458           fOutputContainer->Add(new TH2F(Form("both%d_mpt",mod),Form("Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 459           h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 460           h->GetXaxis()->SetTitle("Mass [GeV]");
 461           h->GetYaxis()->SetTitle("Pt [GEV]");
 462
 463           // Calibration PHOS Module Pi0peak {MIX}
 464           fOutputContainer->Add(new TH2F(Form("mix_both%d_mpt",mod),Form(" Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 465           h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 466           h->GetXaxis()->SetTitle("Mass [GeV]");
 467           h->GetYaxis()->SetTitle("Pt [GEV]");
 468
 469         }
 470
 471         // For efficiensy.
 472         for(Int_t ipid=0; ipid<4; ipid++)
 473         {
 474                 // Real ++++++++++++++++++++++++++++++
 475
 476                 fOutputContainer->Add(new TH2F(Form("%s_mpt_eff",spid[ipid].Data() )," real ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 477                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 478                 h->GetXaxis()->SetTitle("Mass [GeV]");
 479                 h->GetYaxis()->SetTitle("Pt [GEV]");
 480
 481                 // MIX +++++++++++++++++++++++++
 482
 483                 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_eff",spid[ipid].Data() )," mix ", 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 484                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 485                 h->GetXaxis()->SetTitle("Mass [GeV]");
 486                 h->GetYaxis()->SetTitle("Pt [GEV]");
 487
 488                 // Real ++++++++++++++++++++++++++++++
 489
 490                 fOutputContainer->Add(new TH2F(Form("%s_mpt_left_eff",spid[ipid].Data() )," real ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
 491                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 492                 h->GetXaxis()->SetTitle("Mass [GeV]");
 493                 h->GetYaxis()->SetTitle("Pt [GEV]");
 494
 495                 fOutputContainer->Add(new TH2F(Form("%s_mpt_right_eff",spid[ipid].Data() )," real ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
 496                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 497                 h->GetXaxis()->SetTitle("Mass [GeV]");
 498                 h->GetYaxis()->SetTitle("Pt [GEV]");
 499
 500                 // MIX +++++++++++++++++++++++++
 501
 502                 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_left_eff",spid[ipid].Data() )," mix ", 100, 0.05, 0.1, ptMult, ptMin, ptMax ) );
 503                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 504                 h->GetXaxis()->SetTitle("Mass [GeV]");
 505                 h->GetYaxis()->SetTitle("Pt [GEV]");
 506
 507                 fOutputContainer->Add(new TH2F(Form("mix_%s_mpt_right_eff",spid[ipid].Data() )," mix ", 100, 0.2, 0.4, ptMult, ptMin, ptMax ) );
 508                 h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 509                 h->GetXaxis()->SetTitle("Mass [GeV]");
 510                 h->GetYaxis()->SetTitle("Pt [GEV]");
 511         }
 512
 513         // Calibration PHOS Module Pi0peak {REAL}
 514         for(Int_t mod=1; mod<4; mod++){
 515           fOutputContainer->Add(new TH2F(Form("both%d_mpt_eff",mod),Form("Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 516           h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 517           h->GetXaxis()->SetTitle("Mass [GeV]");
 518           h->GetYaxis()->SetTitle("Pt [GEV]");
 519
 520           // Calibration PHOS Module Pi0peak {MIX}
 521           fOutputContainer->Add(new TH2F(Form("mix_both%d_mpt_eff",mod),Form(" Both cuts (CPV + Disp) mod[%d]",mod), 100, massMin, massMax, ptMult, ptMin, ptMax ) );
 522            h = static_cast<TH2F*>(fOutputContainer->Last()) ;
 523           h->GetXaxis()->SetTitle("Mass [GeV]");
 524           h->GetYaxis()->SetTitle("Pt [GEV]");
 525
 526         }
 527
 528         // cout<<"  OK!"<<endl;
 529 }
 530 //_______________________________________________________________________________
 531 void AliPHOSCorrelations::SetHistPtAssoc(Int_t  ptMult, Double_t ptMin, Double_t ptMax)
 532 {
 533         Double_t pi = TMath::Pi();
 534
 535         Int_t PhiMult  =  100;
 536         Float_t PhiMin =  -0.5*pi;
 537         Float_t PhiMax =  1.5*pi;
 538         Int_t EtaMult  =  20;
 539         Float_t EtaMin = -1.;
 540         Float_t EtaMax =  1.;
 541
 542         TString spid[4]={"all","cpv","disp","both"} ;
 543         Int_t PhotonsInMod[6] = {1, 2, 3, 12, 13, 23};
 544
 545         for (int i = 0; i<fAssocBins.GetSize()-1; i++){
 546           for(Int_t ipid=0; ipid<4; ipid++){
 547                 fOutputContainer->Add(new TH3F(Form("%s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
 548                                                Form("%s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
 549                                                ptMult, ptMin, ptMax,  PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
 550                 TH3F * h = static_cast<TH3F*>(fOutputContainer->Last()) ;
 551                 h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
 552                 h->GetYaxis()->SetTitle("#phi [rad]");
 553                 h->GetZaxis()->SetTitle("#eta");
 554
 555                 fOutputContainer->Add(new TH3F(Form("mix_%s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
 556                                                Form("Mixed %s_ptphieta_ptAssoc_%3.1f",spid[ipid].Data(),fAssocBins.At(i+1)),
 557                                                ptMult, ptMin, ptMax,  PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
 558                 h = static_cast<TH3F*>(fOutputContainer->Last()) ;
 559                 h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
 560                 h->GetYaxis()->SetTitle("#phi [rad]");
 561                 h->GetZaxis()->SetTitle("#eta");
 562
 563
 564                 for(Int_t m=0; m<6; m++)
 565                 {
 566                         fOutputContainer->Add(new TH3F(Form("mix_%s_ptphieta_ptAssoc_%3.1f_mod%i",spid[ipid].Data(),fAssocBins.At(i+1), PhotonsInMod[m]),
 567                                                Form("Mixed %s_ptphieta_ptAssoc_%3.1f_mod%i",spid[ipid].Data(),fAssocBins.At(i+1), PhotonsInMod[m]),
 568                                                ptMult, ptMin, ptMax,  PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
 569                         h = static_cast<TH3F*>(fOutputContainer->Last()) ;
 570                         h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
 571                         h->GetYaxis()->SetTitle("#phi [rad]");
 572                         h->GetZaxis()->SetTitle("#eta");
 573                 }
 574
 575                 for(Int_t itpc=1; itpc<3; itpc++)
 576                 {
 577                         fOutputContainer->Add(new TH3F(Form("mix_%s_ptphieta_ptAssoc_%3.1f_tpc%i",spid[ipid].Data(),fAssocBins.At(i+1), itpc),
 578                                                Form("Mixed %s_ptphieta_ptAssoc_%3.1f_tpc%i",spid[ipid].Data(),fAssocBins.At(i+1), itpc),
 579                                                ptMult, ptMin, ptMax,  PhiMult, PhiMin, PhiMax, EtaMult, EtaMin, EtaMax ) );
 580                         h = static_cast<TH3F*>(fOutputContainer->Last()) ;
 581                         h->GetXaxis()->SetTitle("Pt_{triger} [GEV]");
 582                         h->GetYaxis()->SetTitle("#phi [rad]");
 583                         h->GetZaxis()->SetTitle("#eta");
 584                 }
 585           }
 586         }
 587 }
 588
 589 void AliPHOSCorrelations::SetHistPHOSClusterMap()
 590 {
 591         for(int i =  0; i<5; i++)
 592         {
 593                 //  Cluster X/Z/E distribution.
 594                 fOutputContainer->Add(new TH3F(Form("QA_cluXZE_mod%i", i),Form("PHOS Clusters XZE distribution of module %i", i), 70, 0, 70, 60, 0, 60, 200, 0, 20 ) );
 595                 TH3F *h = static_cast<TH3F*>(fOutputContainer->Last()) ;
 596                 h->GetXaxis()->SetTitle("X");
 597                 h->GetYaxis()->SetTitle("Z");
 598                 h->GetZaxis()->SetTitle("E");
 599         }
 600 }
 601 //_______________________________________________________________________________
 602 void AliPHOSCorrelations::UserExec(Option_t *)
 603 {
 604         // Main loop, called for each event analyze ESD/AOD
 605         // Step 0: Event Objects
 606
 607     LogProgress(0);
 608         fEvent = InputEvent();
 609         if( ! fEvent )
 610         {
 611                 AliError("Event could not be retrieved");
 612                 PostData(1, fOutputContainer);
 613                 return ;
 614         }
 615     LogProgress(1);
 616
 617     ZeroingVariables();
 618
 619         fEventESD = dynamic_cast<AliESDEvent*>(fEvent);
 620         fEventAOD = dynamic_cast<AliAODEvent*>(fEvent);
 621
 622         // Step 1(done once):
 623         if( fRunNumber != fEvent->GetRunNumber() )
 624         {
 625                 fRunNumber = fEvent->GetRunNumber();
 626                 fInternalRunNumber = ConvertToInternalRunNumber(fRunNumber);
 627                 //SetGeometry();
 628                 SetESDTrackCuts();
 629         }
 630
 631         //Data
 632         //get Event-Handler for the trigger information
 633         fEventHandler= dynamic_cast<AliInputEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 634         if (!fEventHandler)
 635         {
 636                 AliError("Could not get InputHandler");
 637                 PostData(1, fOutputContainer);
 638                 return; // Reject!
 639         }
 640     LogProgress(2);
 641
 642         if( RejectTriggerMaskSelection() )
 643         {
 644                 PostData(1, fOutputContainer);
 645                 return; // Reject!
 646         }
 647     LogProgress(3);
 648
 649         // Step 2: Vertex
 650         // fVertex, fVertexVector, fVtxBin
 651         SetVertex();
 652         if( RejectEventVertex() )
 653         {
 654                 PostData(1, fOutputContainer);
 655                 return; // Reject!
 656         }
 657     LogProgress(4);
 658
 659         // Step 3: Centrality
 660         // fCentrality, fCentBin
 661         SetCentrality();
 662         if( RejectEventCentrality() )
 663         {
 664                 PostData(1, fOutputContainer);
 665                 return; // Reject!
 666         }
 667     LogProgress(5);
 668         FillHistogram("hCentrality",fCentrality,fInternalRunNumber-0.5) ;
 669
 670         // Step 4: Reaction Plane
 671         // fHaveTPCRP, fRP, fRPV0A, fRPV0C, fRPBin
 672         EvalReactionPlane();
 673         fEMRPBin = GetRPBin();
 674
 675         // Step 5: Event Photons (PHOS Clusters) selectionMakeFlat
 676         SelectPhotonClusters();
 677         if( ! fCaloPhotonsPHOS->GetEntriesFast() )
 678                 LogSelection(kHasPHOSClusters, fInternalRunNumber);
 679
 680         // Step 6: Event Associated particles (TPC Tracks) selection
 681         SelectAccosiatedTracks();
 682         if( ! fTracksTPC->GetEntriesFast() )
 683         LogSelection(kHasTPCTracks, fInternalRunNumber);
 684         LogSelection(kTotalSelected, fInternalRunNumber);
 685
 686     // Step 7: Make TPC's mask
 687     FillTrackEtaPhi();
 688
 689     LogProgress(6);
 690
 691     // Step 8: Start correlation analysis.
 692     // Filling real histograms:
 693     if (fUseMEAlgoritmForReal) // if true than use ME algoritm
 694     {
 695         LogProgress(7);
 696         SelectTriggerPi0ME();   // Extract one most energetic pi0 candidate in this event.
 697
 698         if (fPHOSEvent)
 699         {
 700             ConsiderPi0sME(); // Consider the most energetic Pi0 in this event with all tracks of this event.
 701             LogProgress(8);
 702         }
 703     }
 704     else //using common algoritm
 705     {
 706         LogProgress(10);
 707         if (fPHOSEvent)
 708         {
 709             ConsiderPi0s(); // Extract all pi0 candidates and compare them with all track of this event.
 710             LogProgress(11);
 711         }
 712     }
 713
 714     // Filling mixing histograms:
 715     if(fUseMEAlgoritmForMix) // if true than use ME algoritm
 716     {
 717         // Filling mixing histograms
 718         if (fMBEvent)
 719         {
 720             if (!fUseMEAlgoritmForReal)  SelectTriggerPi0ME();   // Extract one most energetic pi0 candidate in this event if it not done yet.
 721
 722             ConsiderPi0sMix();      // Make background for extracting pi0 mass.
 723             ConsiderTracksMixME();  // Compare only one most energetic pi0 candidate with all tracks from previous MB events.
 724
 725             UpdatePhotonLists();    // Updating pull of photons.
 726             UpdateTrackLists();     // Updating pull of tracks.
 727             LogProgress(9);
 728         }
 729     }
 730     else //using common algoritm
 731     {
 732         if (fMBEvent)
 733         {
 734             ConsiderPi0sMix();      // Make background for extracting pi0 mass.
 735             ConsiderTracksMix();    // Compare all pi0 candidates in current event whit all tracks from previous MB events.
 736
 737             UpdatePhotonLists();
 738             UpdateTrackLists();
 739             LogProgress(12);
 740         }
 741     }
 742
 743     LogProgress(13);
 744         // Post output data.
 745         PostData(1, fOutputContainer);
 746 }
 747 //_______________________________________________________________________________
 748 void AliPHOSCorrelations::SetESDTrackCuts()
 749 {
 750   if( fEventESD ) {
 751     // Create ESD track cut
 752     fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts() ;
 753     //fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
 754     fESDtrackCuts->SetRequireTPCRefit(kTRUE);
 755   }
 756 }
 757 //_______________________________________________________________________________
 758 Int_t AliPHOSCorrelations::ConvertToInternalRunNumber(Int_t run){
 759   if(fPeriod== kLHC11h){
 760                                 switch(run)
 761                                 {
 762                                         case  170593 : return 179 ;
 763                                         case  170572 : return 178 ;
 764                                         case  170556 : return 177 ;
 765                                         case  170552 : return 176 ;
 766                                         case  170546 : return 175 ;
 767                                         case  170390 : return 174 ;
 768                                         case  170389 : return 173 ;
 769                                         case  170388 : return 172 ;
 770                                         case  170387 : return 171 ;
 771                                         case  170315 : return 170 ;
 772                                         case  170313 : return 169 ;
 773                                         case  170312 : return 168 ;
 774                                         case  170311 : return 167 ;
 775                                         case  170309 : return 166 ;
 776                                         case  170308 : return 165 ;
 777                                         case  170306 : return 164 ;
 778                                         case  170270 : return 163 ;
 779                                         case  170269 : return 162 ;
 780                                         case  170268 : return 161 ;
 781                                         case  170267 : return 160 ;
 782                                         case  170264 : return 159 ;
 783                                         case  170230 : return 158 ;
 784                                         case  170228 : return 157 ;
 785                                         case  170208 : return 156 ;
 786                                         case  170207 : return 155 ;
 787                                         case  170205 : return 154 ;
 788                                         case  170204 : return 153 ;
 789                                         case  170203 : return 152 ;
 790                                         case  170195 : return 151 ;
 791                                         case  170193 : return 150 ;
 792                                         case  170163 : return 149 ;
 793                                         case  170162 : return 148 ;
 794                                         case  170159 : return 147 ;
 795                                         case  170155 : return 146 ;
 796                                         case  170152 : return 145 ;
 797                                         case  170091 : return 144 ;
 798                                         case  170089 : return 143 ;
 799                                         case  170088 : return 142 ;
 800                                         case  170085 : return 141 ;
 801                                         case  170084 : return 140 ;
 802                                         case  170083 : return 139 ;
 803                                         case  170081 : return 138 ;
 804                                         case  170040 : return 137 ;
 805                                         case  170038 : return 136 ;
 806                                         case  170036 : return 135 ;
 807                                         case  170027 : return 134 ;
 808                                         case  169981 : return 133 ;
 809                                         case  169975 : return 132 ;
 810                                         case  169969 : return 131 ;
 811                                         case  169965 : return 130 ;
 812                                         case  169961 : return 129 ;
 813                                         case  169956 : return 128 ;
 814                                         case  169926 : return 127 ;
 815                                         case  169924 : return 126 ;
 816                                         case  169923 : return 125 ;
 817                                         case  169922 : return 124 ;
 818                                         case  169919 : return 123 ;
 819                                         case  169918 : return 122 ;
 820                                         case  169914 : return 121 ;
 821                                         case  169859 : return 120 ;
 822                                         case  169858 : return 119 ;
 823                                         case  169855 : return 118 ;
 824                                         case  169846 : return 117 ;
 825                                         case  169838 : return 116 ;
 826                                         case  169837 : return 115 ;
 827                                         case  169835 : return 114 ;
 828                                         case  169683 : return 113 ;
 829                                         case  169628 : return 112 ;
 830                                         case  169591 : return 111 ;
 831                                         case  169590 : return 110 ;
 832                                         case  169588 : return 109 ;
 833                                         case  169587 : return 108 ;
 834                                         case  169586 : return 107 ;
 835                                         case  169584 : return 106 ;
 836                                         case  169557 : return 105 ;
 837                                         case  169555 : return 104 ;
 838                                         case  169554 : return 103 ;
 839                                         case  169553 : return 102 ;
 840                                         case  169550 : return 101 ;
 841                                         case  169515 : return 100 ;
 842                                         case  169512 : return 99 ;
 843                                         case  169506 : return 98 ;
 844                                         case  169504 : return 97 ;
 845                                         case  169498 : return 96 ;
 846                                         case  169475 : return 95 ;
 847                                         case  169420 : return 94 ;
 848                                         case  169419 : return 93 ;
 849                                         case  169418 : return 92 ;
 850                                         case  169417 : return 91 ;
 851                                         case  169415 : return 90 ;
 852                                         case  169411 : return 89 ;
 853                                         case  169238 : return 88 ;
 854                                         case  169236 : return 87 ;
 855                                         case  169167 : return 86 ;
 856                                         case  169160 : return 85 ;
 857                                         case  169156 : return 84 ;
 858                                         case  169148 : return 83 ;
 859                                         case  169145 : return 82 ;
 860                                         case  169144 : return 81 ;
 861                                         case  169143 : return 80 ;
 862                                         case  169138 : return 79 ;
 863                                         case  169099 : return 78 ;
 864                                         case  169094 : return 77 ;
 865                                         case  169091 : return 76 ;
 866                                         case  169045 : return 75 ;
 867                                         case  169044 : return 74 ;
 868                                         case  169040 : return 73 ;
 869                                         case  169035 : return 72 ;
 870                                         case  168992 : return 71 ;
 871                                         case  168988 : return 70 ;
 872                                         case  168984 : return 69 ;
 873                                         case  168826 : return 68 ;
 874                                         case  168777 : return 67 ;
 875                                         case  168514 : return 66 ;
 876                                         case  168512 : return 65 ;
 877                                         case  168511 : return 64 ;
 878                                         case  168467 : return 63 ;
 879                                         case  168464 : return 62 ;
 880                                         case  168461 : return 61 ;
 881                                         case  168460 : return 60 ;
 882                                         case  168458 : return 59 ;
 883                                         case  168362 : return 58 ;
 884                                         case  168361 : return 57 ;
 885                                         case  168356 : return 56 ;
 886                                         case  168342 : return 55 ;
 887                                         case  168341 : return 54 ;
 888                                         case  168325 : return 53 ;
 889                                         case  168322 : return 52 ;
 890                                         case  168318 : return 51 ;
 891                                         case  168311 : return 50 ;
 892                                         case  168310 : return 49 ;
 893                                         case  168213 : return 48 ;
 894                                         case  168212 : return 47 ;
 895                                         case  168208 : return 46 ;
 896                                         case  168207 : return 45 ;
 897                                         case  168206 : return 44 ;
 898                                         case  168205 : return 43 ;
 899                                         case  168204 : return 42 ;
 900                                         case  168203 : return 41 ;
 901                                         case  168181 : return 40 ;
 902                                         case  168177 : return 39 ;
 903                                         case  168175 : return 38 ;
 904                                         case  168173 : return 37 ;
 905                                         case  168172 : return 36 ;
 906                                         case  168171 : return 35 ;
 907                                         case  168115 : return 34 ;
 908                                         case  168108 : return 33 ;
 909                                         case  168107 : return 32 ;
 910                                         case  168105 : return 31 ;
 911                                         case  168104 : return 30 ;
 912                                         case  168103 : return 29 ;
 913                                         case  168076 : return 28 ;
 914                                         case  168069 : return 27 ;
 915                                         case  168068 : return 26 ;
 916                                         case  168066 : return 25 ;
 917                                         case  167988 : return 24 ;
 918                                         case  167987 : return 23 ;
 919                                         case  167986 : return 22 ;
 920                                         case  167985 : return 21 ;
 921                                         case  167921 : return 20 ;
 922                                         case  167920 : return 19 ;
 923                                         case  167915 : return 18 ;
 924                                         case  167909 : return 17 ;
 925                                         case  167903 : return 16 ;
 926                                         case  167902 : return 15 ;
 927                                         case  167818 : return 14 ;
 928                                         case  167814 : return 13 ;
 929                                         case  167813 : return 12 ;
 930                                         case  167808 : return 11 ;
 931                                         case  167807 : return 10 ;
 932                                         case  167806 : return 9 ;
 933                                         case  167713 : return 8 ;
 934                                         case  167712 : return 7 ;
 935                                         case  167711 : return 6 ;
 936                                         case  167706 : return 5 ;
 937                                         case  167693 : return 4 ;
 938                                         case  166532 : return 3 ;
 939                                         case  166530 : return 2 ;
 940                                         case  166529 : return 1 ;
 941
 942                                         default : return 199;
 943                                 }
 944                         }
 945         if(fPeriod== kLHC10h){
 946                 switch(run){
 947                                                 case  139517 : return 137;
 948                                                 case  139514 : return 136;
 949                                                 case  139513 : return 135;
 950                                                 case  139511 : return 134;
 951                                                 case  139510 : return 133;
 952                                                 case  139507 : return 132;
 953                                                 case  139505 : return 131;
 954                                                 case  139504 : return 130;
 955                                                 case  139503 : return 129;
 956                                                 case  139470 : return 128;
 957                                                 case  139467 : return 127;
 958                                                 case  139466 : return 126;
 959                                                 case  139465 : return 125;
 960                                                 case  139440 : return 124;
 961                                                 case  139439 : return 123;
 962                                                 case  139438 : return 122;
 963                                                 case  139437 : return 121;
 964                                                 case  139360 : return 120;
 965                                                 case  139329 : return 119;
 966                                                 case  139328 : return 118;
 967                                                 case  139314 : return 117;
 968                                                 case  139311 : return 116;
 969                                                 case  139310 : return 115;
 970                                                 case  139309 : return 114;
 971                                                 case  139308 : return 113;
 972                                                 case  139173 : return 112;
 973                                                 case  139172 : return 111;
 974                                                 case  139110 : return 110;
 975                                                 case  139107 : return 109;
 976                                                 case  139105 : return 108;
 977                                                 case  139104 : return 107;
 978                                                 case  139042 : return 106;
 979                                                 case  139038 : return 105;
 980                                                 case  139037 : return 104;
 981                                                 case  139036 : return 103;
 982                                                 case  139029 : return 102;
 983                                                 case  139028 : return 101;
 984                                                 case  138983 : return 100;
 985                                                 case  138982 : return 99;
 986                                                 case  138980 : return 98;
 987                                                 case  138979 : return 97;
 988                                                 case  138978 : return 96;
 989                                                 case  138977 : return 95;
 990                                                 case  138976 : return 94;
 991                                                 case  138973 : return 93;
 992                                                 case  138972 : return 92;
 993                                                 case  138965 : return 91;
 994                                                 case  138924 : return 90;
 995                                                 case  138872 : return 89;
 996                                                 case  138871 : return 88;
 997                                                 case  138870 : return 87;
 998                                                 case  138837 : return 86;
 999                                                 case  138830 : return 85;
1000                                                 case  138828 : return 84;
1001                                                 case  138826 : return 83;
1002                                                 case  138796 : return 82;
1003                                                 case  138795 : return 81;
1004                                                 case  138742 : return 80;
1005                                                 case  138732 : return 79;
1006                                                 case  138730 : return 78;
1007                                                 case  138666 : return 77;
1008                                                 case  138662 : return 76;
1009                                                 case  138653 : return 75;
1010                                                 case  138652 : return 74;
1011                                                 case  138638 : return 73;
1012                                                 case  138624 : return 72;
1013                                                 case  138621 : return 71;
1014                                                 case  138583 : return 70;
1015                                                 case  138582 : return 69;
1016                                                 case  138579 : return 68;
1017                                                 case  138578 : return 67;
1018                                                 case  138534 : return 66;
1019                                                 case  138469 : return 65;
1020                                                 case  138442 : return 64;
1021                                                 case  138439 : return 63;
1022                                                 case  138438 : return 62;
1023                                                 case  138396 : return 61;
1024                                                 case  138364 : return 60;
1025                                                 case  138359 : return 59;
1026                                                 case  138275 : return 58;
1027                                                 case  138225 : return 57;
1028                                                 case  138201 : return 56;
1029                                                 case  138200 : return 55;
1030                                                 case  138197 : return 54;
1031                                                 case  138192 : return 53;
1032                                                 case  138190 : return 52;
1033                                                 case  138154 : return 51;
1034                                                 case  138153 : return 50;
1035                                                 case  138151 : return 49;
1036                                                 case  138150 : return 48;
1037                                                 case  138126 : return 47;
1038                                                 case  138125 : return 46;
1039                                                 case  137848 : return 45;
1040                                                 case  137847 : return 44;
1041                                                 case  137844 : return 43;
1042                                                 case  137843 : return 42;
1043                                                 case  137752 : return 41;
1044                                                 case  137751 : return 40;
1045                                                 case  137748 : return 39;
1046                                                 case  137724 : return 38;
1047                                                 case  137722 : return 37;
1048                                                 case  137718 : return 36;
1049                                                 case  137704 : return 35;
1050                                                 case  137693 : return 34;
1051                                                 case  137692 : return 33;
1052                                                 case  137691 : return 32;
1053                                                 case  137689 : return 31;
1054                                                 case  137686 : return 30;
1055                                                 case  137685 : return 29;
1056                                                 case  137639 : return 28;
1057                                                 case  137638 : return 27;
1058                                                 case  137608 : return 26;
1059                                                 case  137595 : return 25;
1060                                                 case  137549 : return 24;
1061                                                 case  137546 : return 23;
1062                                                 case  137544 : return 22;
1063                                                 case  137541 : return 21;
1064                                                 case  137539 : return 20;
1065                                                 case  137531 : return 19;
1066                                                 case  137530 : return 18;
1067                                                 case  137443 : return 17;
1068                                                 case  137441 : return 16;
1069                                                 case  137440 : return 15;
1070                                                 case  137439 : return 14;
1071                                                 case  137434 : return 13;
1072                                                 case  137432 : return 12;
1073                                                 case  137431 : return 11;
1074                                                 case  137430 : return 10;
1075                                                 case  137366 : return 9;
1076                                                 case  137243 : return 8;
1077                                                 case  137236 : return 7;
1078                                                 case  137235 : return 6;
1079                                                 case  137232 : return 5;
1080                                                 case  137231 : return 4;
1081                                                 case  137165 : return 3;
1082                                                 case  137162 : return 2;
1083                                                 case  137161 : return 1;
1084                                                 default : return 199;
1085                                         }
1086                                 }
1087                                 if( kLHC13 == fPeriod )
1088                                 {
1089                                         switch(run)
1090                                         {
1091                                                 case  195344 : return 1;
1092                                                 case  195346 : return 2;
1093                                                 case  195351 : return 3;
1094                                                 case  195389 : return 4;
1095                                                 case  195390 : return 5;
1096                                                 case  195391 : return 6;
1097                                                 case  195478 : return 7;
1098                                                 case  195479 : return 8;
1099                                                 case  195480 : return 9;
1100                                                 case  195481 : return 10;
1101                                                 case  195482 : return 11;
1102                                                 case  195483 : return 12;
1103                                                 case  195529 : return 13;
1104                                                 case  195531 : return 14;
1105                                                 case  195532 : return 15;
1106                                                 case  195566 : return 16;
1107                                                 case  195567 : return 17;
1108                                                 case  195568 : return 18;
1109                                                 case  195592 : return 19;
1110                                                 case  195593 : return 20;
1111                                                 case  195596 : return 21;
1112                                                 case  195633 : return 22;
1113                                                 case  195635 : return 23;
1114                                                 case  195644 : return 24;
1115                                                 case  195673 : return 25;
1116                                                 case  195675 : return 26;
1117                                                 case  195676 : return 27;
1118                                                 case  195677 : return 28;
1119                                                 case  195681 : return 29;
1120                                                 case  195682 : return 30;
1121                                                 case  195720 : return 31;
1122                                                 case  195721 : return 32;
1123                                                 case  195722 : return 33;
1124                                                 case  195724 : return 34;
1125                                                 case  195725 : return 34;
1126                                                 case  195726 : return 35;
1127                                                 case  195727 : return 36;
1128                                                 case  195760 : return 37;
1129                                                 case  195761 : return 38;
1130                                                 case  195765 : return 39;
1131                                                 case  195767 : return 40;
1132                                                 case  195783 : return 41;
1133                                                 case  195787 : return 42;
1134                                                 case  195826 : return 43;
1135                                                 case  195827 : return 44;
1136                                                 case  195829 : return 45;
1137                                                 case  195830 : return 46;
1138                                                 case  195831 : return 47;
1139                                                 case  195867 : return 48;
1140                                                 case  195869 : return 49;
1141                                                 case  195871 : return 50;
1142                                                 case  195872 : return 51;
1143                                                 case  195873 : return 52;
1144                                                 case  195935 : return 53;
1145                                                 case  195949 : return 54;
1146                                                 case  195950 : return 55;
1147                                                 case  195954 : return 56;
1148                                                 case  195955 : return 57;
1149                                                 case  195958 : return 58;
1150                                                 case  195989 : return 59;
1151                                                 case  195994 : return 60;
1152                                                 case  195998 : return 61;
1153                                                 case  196000 : return 62;
1154                                                 case  196006 : return 63;
1155                                                 case  196085 : return 64;
1156                                                 case  196089 : return 65;
1157                                                 case  196090 : return 66;
1158                                                 case  196091 : return 67;
1159                                                 case  196099 : return 68;
1160                                                 case  196105 : return 69;
1161                                                 case  196107 : return 70;
1162                                                 case  196185 : return 71;
1163                                                 case  196187 : return 72;
1164                                                 case  196194 : return 73;
1165                                                 case  196197 : return 74;
1166                                                 case  196199 : return 75;
1167                                                 case  196200 : return 76;
1168                                                 case  196201 : return 77;
1169                                                 case  196203 : return 78;
1170                                                 case  196208 : return 79;
1171                                                 case  196214 : return 80;
1172                                                 case  196308 : return 81;
1173                                                 case  196309 : return 82;
1174                                                 case  196310 : return 83;
1175                                                 case  196311 : return 84;
1176                                                 case  196433 : return 85;
1177                                                 case  196474 : return 86;
1178                                                 case  196475 : return 87;
1179                                                 case  196477 : return 88;
1180                                                 case  196528 : return 89;
1181                                                 case  196533 : return 90;
1182                                                 case  196535 : return 91;
1183                                                 case  196563 : return 92;
1184                                                 case  196564 : return 93;
1185                                                 case  196566 : return 94;
1186                                                 case  196568 : return 95;
1187                                                 case  196601 : return 96;
1188                                                 case  196605 : return 97;
1189                                                 case  196608 : return 98;
1190                                                 case  196646 : return 99;
1191                                                 case  196648 : return 100;
1192                                                 case  196701 : return 101;
1193                                                 case  196702 : return 102;
1194                                                 case  196703 : return 103;
1195                                                 case  196706 : return 104;
1196                                                 case  196714 : return 105;
1197                                                 case  196720 : return 106;
1198                                                 case  196721 : return 107;
1199                                                 case  196722 : return 108;
1200                                                 case  196772 : return 109;
1201                                                 case  196773 : return 110;
1202                                                 case  196774 : return 111;
1203                                                 case  196869 : return 112;
1204                                                 case  196870 : return 113;
1205                                                 case  196874 : return 114;
1206                                                 case  196876 : return 115;
1207                                                 case  196965 : return 116;
1208                                                 case  196967 : return 117;
1209                                                 case  196972 : return 118;
1210                                                 case  196973 : return 119;
1211                                                 case  196974 : return 120;
1212                                                 case  197003 : return 121;
1213                                                 case  197011 : return 122;
1214                                                 case  197012 : return 123;
1215                                                 case  197015 : return 124;
1216                                                 case  197027 : return 125;
1217                                                 case  197031 : return 126;
1218                                                 case  197089 : return 127;
1219                                                 case  197090 : return 128;
1220                                                 case  197091 : return 129;
1221                                                 case  197092 : return 130;
1222                                                 case  197094 : return 131;
1223                                                 case  197098 : return 132;
1224                                                 case  197099 : return 133;
1225                                                 case  197138 : return 134;
1226                                                 case  197139 : return 135;
1227                                                 case  197142 : return 136;
1228                                                 case  197143 : return 137;
1229                                                 case  197144 : return 138;
1230                                                 case  197145 : return 139;
1231                                                 case  197146 : return 140;
1232                                                 case  197147 : return 140;
1233                                                 case  197148 : return 141;
1234                                                 case  197149 : return 142;
1235                                                 case  197150 : return 143;
1236                                                 case  197152 : return 144;
1237                                                 case  197153 : return 145;
1238                                                 case  197184 : return 146;
1239                                                 case  197189 : return 147;
1240                                                 case  197247 : return 148;
1241                                                 case  197248 : return 149;
1242                                                 case  197254 : return 150;
1243                                                 case  197255 : return 151;
1244                                                 case  197256 : return 152;
1245                                                 case  197258 : return 153;
1246                                                 case  197260 : return 154;
1247                                                 case  197296 : return 155;
1248                                                 case  197297 : return 156;
1249                                                 case  197298 : return 157;
1250                                                 case  197299 : return 158;
1251                                                 case  197300 : return 159;
1252                                                 case  197302 : return 160;
1253                                                 case  197341 : return 161;
1254                                                 case  197342 : return 162;
1255                                                 case  197348 : return 163;
1256                                                 case  197349 : return 164;
1257                                                 case  197351 : return 165;
1258                                                 case  197386 : return 166;
1259                                                 case  197387 : return 167;
1260                                                 case  197388 : return 168;
1261                                                 default : return 199;
1262                                         }
1263                                 }
1264                                 if((fPeriod == kUndefinedPeriod) && (fDebug >= 1) )
1265                                 {
1266                                         AliWarning("Period not defined");
1267                                 }
1268                                 return 1;
1269                         }
1270
1271 //_______________________________________________________________________________
1272 Bool_t AliPHOSCorrelations::RejectTriggerMaskSelection()
1273 {
1274         // Analyse trigger event and reject it if it not intresting.
1275         const Bool_t REJECT = true;
1276         const Bool_t ACCEPT = false;
1277
1278         if( fDebug >= 2 )
1279                 AliInfo( Form("Event passed offline phos trigger test: %s ", fEvent->GetFiredTriggerClasses().Data() ) );
1280
1281         Int_t physSelMask = fEventHandler->IsEventSelected();
1282
1283         Bool_t isAny = physSelMask & AliVEvent::kAny;
1284
1285         Bool_t isPHI1 = physSelMask & AliVEvent::kPHI1;
1286         Bool_t isPHI7 = physSelMask & AliVEvent::kPHI7;
1287         Bool_t isPHI8 = physSelMask & AliVEvent::kPHI8;
1288         Bool_t isCentral = physSelMask & AliVEvent::kCentral;
1289         Bool_t isSemiCentral = physSelMask & AliVEvent::kSemiCentral;
1290         Bool_t isPHOSPb = physSelMask & AliVEvent::kPHOSPb;
1291
1292         Bool_t isMB = physSelMask & AliVEvent::kMB;
1293         Bool_t isINT7 = physSelMask & AliVEvent::kINT7;
1294         Bool_t isAnyINT = physSelMask & AliVEvent::kAnyINT;
1295
1296
1297         // All input events
1298         FillHistogram("hTriggerPassedEvents", 0);
1299         //if ( !isAny ) cout<<"Strange event"<<endl; // We loose some events O_o
1300         if ( isAny ) FillHistogram("hTriggerPassedEvents", 1.);
1301
1302         // PHOS events.
1303         if ( isPHI1 ) FillHistogram("hTriggerPassedEvents",  2.);
1304         if ( isPHI7 ) FillHistogram("hTriggerPassedEvents",  3.);
1305         if (isPHI8) FillHistogram("hTriggerPassedEvents",  4.);
1306         if ( isCentral ) FillHistogram("hTriggerPassedEvents", 5.);
1307         if (isSemiCentral) FillHistogram("hTriggerPassedEvents", 6.);
1308         if (isPHOSPb) FillHistogram("hTriggerPassedEvents", 7.);
1309
1310         // MB events.
1311         if ( isMB ) FillHistogram("hTriggerPassedEvents", 8.);
1312         if ( isINT7 ) FillHistogram("hTriggerPassedEvents", 9.);
1313         if ( isAnyINT ) FillHistogram("hTriggerPassedEvents", 10.);
1314
1315
1316         Bool_t isTriggerEvent = isPHI1 || isPHI7 || isPHI8 || isCentral || isSemiCentral || isPHOSPb;
1317         // PHOS trigger events.
1318         if ( isTriggerEvent )
1319         {
1320                 FillHistogram("hTriggerPassedEvents", 16.);
1321                 fPHOSEvent = true;
1322                 fMBEvent = false;
1323                 return ACCEPT;
1324         }
1325
1326         // MB events without trigger PHOS events.
1327     //if ( (isMB || isINT7 || isAnyINT) && !isTriggerEvent ) //Less enteries in PbPb.
1328     if ( !isTriggerEvent )
1329         {
1330                 FillHistogram("hTriggerPassedEvents", 17.);
1331                 fPHOSEvent = false;
1332                 fMBEvent = true;
1333                 return ACCEPT;
1334         }
1335
1336         // other events
1337         FillHistogram("hTriggerPassedEvents",  18.);
1338         return REJECT;
1339 }
1340 //_______________________________________________________________________________
1341 void AliPHOSCorrelations::SetVertex()
1342 {
1343         const AliVVertex *primaryVertex = fEvent->GetPrimaryVertex();
1344         if( primaryVertex )
1345         {
1346                 fVertex[0] = primaryVertex->GetX();
1347                 fVertex[1] = primaryVertex->GetY();
1348                 fVertex[2] = primaryVertex->GetZ();
1349         }
1350         else
1351         {
1352                 //AliError("Event has 0x0 Primary Vertex, defaulting to origo");
1353                 fVertex[0] = 0;
1354                 fVertex[1] = 0;
1355                 fVertex[2] = 0;
1356         }
1357         fVertexVector = TVector3(fVertex);
1358
1359         fVtxBin=0 ;// No support for vtx binning implemented.
1360 }
1361 //_______________________________________________________________________________
1362 Bool_t AliPHOSCorrelations::RejectEventVertex()
1363 {
1364   if( ! fEvent->GetPrimaryVertex() )
1365     return true; // reject
1366   LogSelection(kHasVertex, fInternalRunNumber);
1367
1368   if ( TMath::Abs(fVertexVector.z()) > fMaxAbsVertexZ )
1369     return true; // reject
1370   LogSelection(kHasAbsVertex, fInternalRunNumber);
1371
1372   return false; // accept event.
1373 }
1374 //_______________________________________________________________________________
1375 void AliPHOSCorrelations::SetCentrality()
1376 {
1377         AliCentrality *centrality = fEvent->GetCentrality();
1378         if( centrality )
1379                 fCentrality=centrality->GetCentralityPercentile(fCentralityEstimator);
1380         else
1381         {
1382                 AliError("Event has 0x0 centrality");
1383                 fCentrality = -1.;
1384         }
1385
1386         //cout<<"fCentrality: "<<fCentrality<<endl;
1387         //FillHistogram("hCentrality",fCentrality,fInternalRunNumber-0.5) ;
1388         fCentBin = GetCentralityBin(fCentrality);
1389 }
1390 //_______________________________________________________________________________
1391 Bool_t AliPHOSCorrelations::RejectEventCentrality()
1392 {
1393         if (fCentrality<fCentCutoffDown)
1394                 return true; //reject
1395         if(fCentrality>fCentCutoffUp)
1396                 return true;
1397
1398         return false;  // accept event.
1399 }
1400 //_______________________________________________________________________________
1401 void AliPHOSCorrelations::SetCentralityBinning(const TArrayD& edges, const TArrayI& nMixed){
1402 // Define centrality bins by their edges
1403   for(int i=0; i<edges.GetSize()-1; ++i)
1404     if(edges.At(i) > edges.At(i+1)) AliFatal("edges are not sorted");
1405   if( edges.GetSize() != nMixed.GetSize()+1) AliFatal("edges and nMixed don't have appropriate relative sizes");
1406
1407   fCentEdges = edges;
1408   fCentNMixed = nMixed;
1409 }
1410 //_______________________________________________________________________________
1411 Int_t AliPHOSCorrelations::GetCentralityBin(Float_t centralityV0M){
1412   int lastBinUpperIndex = fCentEdges.GetSize() -1;
1413   if( centralityV0M > fCentEdges[lastBinUpperIndex] ) {
1414     if( fDebug >= 1 )
1415       AliWarning( Form("centrality (%f) larger then upper edge of last centrality bin (%f)!", centralityV0M, fCentEdges[lastBinUpperIndex]) );
1416     return lastBinUpperIndex-1;
1417   }
1418   if( centralityV0M < fCentEdges[0] ) {
1419     if( fDebug >= 1 )
1420       AliWarning( Form("centrality (%f) smaller then lower edge of first bin (%f)!", centralityV0M, fCentEdges[0]) );
1421     return 0;
1422   }
1423
1424   fCentBin = TMath::BinarySearch<Double_t> ( GetNumberOfCentralityBins(), fCentEdges.GetArray(), centralityV0M );
1425   return fCentBin;
1426 }
1427 //_______________________________________________________________________________
1428 void AliPHOSCorrelations::SetCentralityBorders (double down , double up ){
1429   if (down < 0. || up > 100 || up<=down)
1430      AliError( Form("Warning. Bad value of centrality borders. Setting as default: fCentCutoffDown=%2.f, fCentCutoffUp=%2.f",fCentCutoffDown,fCentCutoffUp) );
1431   else{
1432     fCentCutoffDown = down;
1433     fCentCutoffUp = up;
1434     AliInfo( Form("Centrality border was set as fCentCutoffDown=%2.f, fCentCutoffUp=%2.f",fCentCutoffDown,fCentCutoffUp) );
1435   }
1436 }
1437
1438 //_______________________________________________________________________________
1439 void AliPHOSCorrelations::EvalReactionPlane()
1440 {
1441         // assigns: fHaveTPCRP and fRP
1442         // also does a few histogram fills
1443
1444         AliEventplane *eventPlane = fEvent->GetEventplane();
1445         if( ! eventPlane ) { AliError("Event has no event plane"); return; }
1446
1447         Double_t reactionPlaneQ = eventPlane->GetEventplane("Q");
1448
1449         if(reactionPlaneQ>=999 || reactionPlaneQ < 0.)
1450         {
1451                 //reaction plain was not defined
1452                 fHaveTPCRP = kFALSE;
1453         }
1454         else
1455         {
1456                 fHaveTPCRP = kTRUE;
1457         }
1458
1459         if(fHaveTPCRP)
1460                 fRP = reactionPlaneQ;
1461         else
1462                 fRP = 0.;
1463
1464         FillHistogram("phiRPflat",fRP,fCentrality) ;
1465 }
1466 //_______________________________________________________________________________
1467 Int_t AliPHOSCorrelations::GetRPBin()
1468 {
1469         Double_t averageRP;
1470         averageRP = fRP ;       // If possible, it is better to have EP bin from TPC
1471                                 // to have similar events for miximng (including jets etc)   (fRPV0A+fRPV0C+fRP) /3.;
1472
1473         fEMRPBin = Int_t(fNEMRPBins*(averageRP)/TMath::Pi());
1474
1475         if( fEMRPBin > (Int_t)fNEMRPBins-1 )
1476                 fEMRPBin = fNEMRPBins-1 ;
1477         else
1478         if(fEMRPBin < 0) fEMRPBin=0;
1479
1480         return fEMRPBin;
1481 }
1482 //_______________________________________________________________________________
1483 void AliPHOSCorrelations::SelectPhotonClusters()
1484 {
1485         //Selects PHOS clusters
1486
1487         // clear (or create) array for holding events photons/clusters
1488         if(fCaloPhotonsPHOS)
1489                 fCaloPhotonsPHOS->Clear();
1490         else
1491         {
1492                 fCaloPhotonsPHOS = new TClonesArray("AliCaloPhoton",200);
1493                 fCaloPhotonsPHOS->SetOwner();
1494         }
1495
1496         Int_t inPHOS=0 ;
1497
1498         for (Int_t i=0;  i<fEvent->GetNumberOfCaloClusters();  i++)
1499         {
1500                 AliVCluster *clu = fEvent->GetCaloCluster(i);
1501                 if (!clu->IsPHOS() || clu->E()< fMinClusterEnergy) continue; // reject cluster
1502
1503                 Float_t  position[3];
1504                 clu->GetPosition(position);
1505                 TVector3 global(position) ;
1506                 Int_t relId[4] ;
1507                 fPHOSGeo->GlobalPos2RelId(global,relId) ;
1508                 Int_t modPHOS  = relId[0] ;
1509                 Int_t cellXPHOS = relId[2];
1510                 Int_t cellZPHOS = relId[3] ;
1511
1512                 Double_t distBC=clu->GetDistanceToBadChannel();
1513                 if(distBC<fMinBCDistance)
1514                         continue ;
1515
1516                 if(clu->GetNCells() < fMinNCells) continue ;
1517                 if(clu->GetM02() < fMinM02)   continue ;
1518
1519                 if(fTOFCutEnabled){
1520                         Double_t tof = clu->GetTOF();
1521                         if(TMath::Abs(tof) > fTOFCut ) continue ;
1522                 }
1523                 TLorentzVector lorentzMomentum;
1524                 Double_t ecore = clu->GetCoreEnergy();
1525                 //Double_t ecore = clu->E();
1526
1527                 FillHistogram("hCluEvsClu", clu->E(), clu->GetNCells()) ;
1528
1529                 Double_t origo[3] = {0,0,0}; // don't rely on event vertex, assume (0,0,0) ?
1530                 //clu->GetMomentum(lorentzMomentum, fVertex);
1531                 clu->GetMomentum(lorentzMomentum, origo);
1532
1533                 if(inPHOS>=fCaloPhotonsPHOS->GetSize()){
1534                         fCaloPhotonsPHOS->Expand(inPHOS+50) ;
1535                 }
1536
1537                 AliCaloPhoton * ph =new((*fCaloPhotonsPHOS)[inPHOS]) AliCaloPhoton(lorentzMomentum.X(),lorentzMomentum.Py(),lorentzMomentum.Z(),lorentzMomentum.E());
1538                 inPHOS++ ;
1539                 ph->SetCluster(clu);
1540
1541                 /*Float_t cellId=clu->GetCellAbsId(0) ;
1542                 Int_t mod = (Int_t)TMath:: Ceil(cellId/(56*64) ) ; */
1543                 ph->SetModule(modPHOS) ;
1544
1545                 lorentzMomentum*=ecore/lorentzMomentum.E() ;
1546
1547                 //ph->SetNCells(clu->GetNCells());
1548                 ph->SetMomV2(&lorentzMomentum) ;
1549                 ph->SetDispBit(clu->GetDispersion()<2.5) ;
1550                 ph->SetCPVBit(clu->GetEmcCpvDistance()>2.) ;
1551
1552                 FillHistogram(Form("QA_cluXZE_mod%i", modPHOS), cellXPHOS, cellZPHOS, lorentzMomentum.E() ) ;
1553         }
1554 }
1555 //_______________________________________________________________________________
1556 void AliPHOSCorrelations::SelectAccosiatedTracks()
1557 {
1558         // clear (or create) array for holding events tracks
1559         if(fTracksTPC)
1560                 fTracksTPC->Clear();
1561         else
1562         {
1563                 fTracksTPC = new TClonesArray("TLorentzVector",12000);
1564         }
1565         Int_t iTracks=0 ;
1566         for (Int_t i=0; i<fEvent->GetNumberOfTracks(); i++)
1567         {
1568
1569                 AliVParticle *track = fEvent->GetTrack(i);
1570                         if(fEventESD){
1571                         if(!SelectESDTrack((AliESDtrack*)track)) continue ;
1572                 }
1573                 else{
1574                         if(!SelectAODTrack((AliAODTrack*)track)) continue ;
1575                 }
1576                 Double_t px = track->Px();
1577                 Double_t py = track->Py();
1578                 Double_t pz = track->Pz() ;
1579                 Double_t e = track->E() ;
1580
1581                 if(iTracks>=fTracksTPC->GetSize())
1582                                 fTracksTPC->Expand(iTracks+50) ;
1583
1584                 new((*fTracksTPC)[iTracks]) TLorentzVector(px, py, pz,e);
1585                 iTracks++ ;
1586         }
1587 }
1588 //_______________________________________________________________________________
1589 void AliPHOSCorrelations::SelectTriggerPi0ME()
1590 {
1591     const Int_t nPHOS=fCaloPhotonsPHOS->GetEntriesFast() ;
1592     for(Int_t i1=0; i1 < nPHOS-1; i1++)
1593     {
1594         AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1595         for (Int_t i2=i1+1; i2<nPHOS; i2++)
1596         {
1597             AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
1598             TLorentzVector p12  = *ph1  + *ph2;
1599
1600             Double_t phiTrigger=p12.Phi() ;
1601             Double_t etaTrigger=p12.Eta() ;
1602
1603             Double_t m=p12.M() ;
1604             Double_t pt=p12.Pt() ;
1605             Double_t eff = 1./GetEfficiency(pt);
1606             int mod1 = ph1->Module() ;
1607             int mod2 = ph2->Module() ;
1608
1609             FillHistogram("clu_phieta",phiTrigger,etaTrigger);
1610             FillHistogram("clusingle_phieta",ph1->Phi(), ph1->Eta());
1611             FillHistogram("clusingle_phieta",ph2->Phi(), ph2->Eta());
1612
1613
1614             FillHistogram("all_mpt",m, pt);
1615             FillHistogram("all_mpt_left",m, pt);
1616             FillHistogram("all_mpt_right",m, pt);
1617
1618             FillHistogram("all_mpt_eff",m, pt, eff);
1619             FillHistogram("all_mpt_left_eff",m, pt, eff);
1620             FillHistogram("all_mpt_right_eff",m, pt, eff);
1621
1622             if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1623             {
1624                 FillHistogram("cpv_mpt",m, pt);
1625                 FillHistogram("cpv_mpt_left",m, pt);
1626                 FillHistogram("cpv_mpt_right",m, pt);
1627
1628                 FillHistogram("cpv_mpt_eff",m, pt, eff);
1629                 FillHistogram("cpv_mpt_left_eff",m, pt, eff);
1630                 FillHistogram("cpv_mpt_right_eff",m, pt, eff);
1631             }
1632
1633             if ( ph1->IsDispOK() && ph2->IsDispOK() )
1634             {
1635                 FillHistogram("disp_mpt",m, pt);
1636                 FillHistogram("disp_mpt_left",m, pt);
1637                 FillHistogram("disp_mpt_right",m, pt);
1638
1639                 FillHistogram("disp_mpt_eff",m, pt, eff);
1640                 FillHistogram("disp_mpt_left_eff",m, pt, eff);
1641                 FillHistogram("disp_mpt_right_eff",m, pt, eff);
1642                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1643                 {
1644                     FillHistogram("both_mpt",m, pt);
1645                     FillHistogram("both_mpt_left",m, pt);
1646                     FillHistogram("both_mpt_right",m, pt);
1647
1648                     FillHistogram("both_mpt_eff",m, pt, eff);
1649                     FillHistogram("both_mpt_left_eff",m, pt, eff);
1650                     FillHistogram("both_mpt_right_eff",m, pt, eff);
1651                     if(mod1 == mod2) // for each module
1652                     {
1653                         FillHistogram(Form("both%d_mpt",mod1),m, pt);
1654                         FillHistogram(Form("both%d_mpt_eff",mod1),m, pt, eff);
1655                     }
1656                 }
1657             }
1658
1659             if(!TestMass(m,pt)) continue;
1660
1661             Int_t modCase = GetModCase(mod1, mod2);
1662
1663             //Now we choosing most energetic pi0.
1664             TestPi0ME(kPidAll, p12, modCase);
1665             if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1666                 TestPi0ME(kPidCPV, p12, modCase);
1667             if ( ph1->IsDispOK() && ph2->IsDispOK() )
1668             {
1669                 TestPi0ME(kPidDisp, p12, modCase);
1670                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1671                     TestPi0ME(kPidBoth, p12, modCase);
1672             }
1673         }
1674     }
1675 }
1676
1677 //_______________________________________________________________________________
1678 void AliPHOSCorrelations::ConsiderPi0s()
1679 {
1680     // Consider all photons from PHOS
1681         const Int_t nPHOS=fCaloPhotonsPHOS->GetEntriesFast() ;
1682         for(Int_t i1=0; i1 < nPHOS-1; i1++)
1683         {
1684                 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1685                 for (Int_t i2=i1+1; i2<nPHOS; i2++)
1686                 {
1687                         AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
1688                         TLorentzVector p12  = *ph1  + *ph2;
1689
1690                         Double_t phiTrigger=p12.Phi() ;
1691                         Double_t etaTrigger=p12.Eta() ;
1692
1693                         Double_t m=p12.M() ;
1694                         Double_t pt=p12.Pt() ;
1695                         Double_t eff = 1./GetEfficiency(pt);
1696                         int mod1 = ph1->Module() ;
1697                         int mod2 = ph2->Module() ;
1698
1699                         FillHistogram("clu_phieta",phiTrigger,etaTrigger);
1700                         FillHistogram("clusingle_phieta",ph1->Phi(), ph1->Eta());
1701                         FillHistogram("clusingle_phieta",ph2->Phi(), ph2->Eta());
1702
1703
1704                         FillHistogram("all_mpt",m, pt);
1705                         FillHistogram("all_mpt_left",m, pt);
1706                         FillHistogram("all_mpt_right",m, pt);
1707
1708                         FillHistogram("all_mpt_eff",m, pt, eff);
1709                         FillHistogram("all_mpt_left_eff",m, pt, eff);
1710                         FillHistogram("all_mpt_right_eff",m, pt, eff);
1711
1712                         if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1713                         {
1714                                 FillHistogram("cpv_mpt",m, pt);
1715                                 FillHistogram("cpv_mpt_left",m, pt);
1716                                 FillHistogram("cpv_mpt_right",m, pt);
1717
1718                                 FillHistogram("cpv_mpt_eff",m, pt, eff);
1719                                 FillHistogram("cpv_mpt_left_eff",m, pt, eff);
1720                                 FillHistogram("cpv_mpt_right_eff",m, pt, eff);
1721                         }
1722
1723                         if ( ph1->IsDispOK() && ph2->IsDispOK() )
1724                         {
1725                                 FillHistogram("disp_mpt",m, pt);
1726                                 FillHistogram("disp_mpt_left",m, pt);
1727                                 FillHistogram("disp_mpt_right",m, pt);
1728
1729                                 FillHistogram("disp_mpt_eff",m, pt, eff);
1730                                 FillHistogram("disp_mpt_left_eff",m, pt, eff);
1731                                 FillHistogram("disp_mpt_right_eff",m, pt, eff);
1732                                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1733                                 {
1734                                         FillHistogram("both_mpt",m, pt);
1735                                         FillHistogram("both_mpt_left",m, pt);
1736                                         FillHistogram("both_mpt_right",m, pt);
1737
1738                                         FillHistogram("both_mpt_eff",m, pt, eff);
1739                                         FillHistogram("both_mpt_left_eff",m, pt, eff);
1740                                         FillHistogram("both_mpt_right_eff",m, pt, eff);
1741                                         if(mod1 == mod2) // for each module
1742                                         {
1743                                                 FillHistogram(Form("both%d_mpt",mod1),m, pt);
1744                                                 FillHistogram(Form("both%d_mpt_eff",mod1),m, pt, eff);
1745                                         }
1746                                 }
1747                         }
1748
1749                         if(!TestMass(m,pt)) continue;
1750
1751             FillHistogram("nTrigger_all", pt, eff);
1752                         if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1753                                 FillHistogram("nTrigger_cpv", pt, eff);
1754                         if ( ph1->IsDispOK() && ph2->IsDispOK() )
1755                         {
1756                                 FillHistogram("nTrigger_disp", pt, eff);
1757                                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1758                                         FillHistogram("nTrigger_both", pt, eff);
1759                         }
1760
1761                         // Take track's angles and compare with cluster's angles.
1762                         for(Int_t i3=0; i3<fTracksTPC->GetEntriesFast(); i3++){
1763                                 TLorentzVector * track = (TLorentzVector*)fTracksTPC->At(i3);
1764
1765                                 Double_t phiAssoc = track->Phi();
1766                                 Double_t etaAssoc = track->Eta();
1767                                 Double_t ptAssoc = track->Pt();
1768
1769                                 Double_t dPhi = phiTrigger - phiAssoc;
1770                                 while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1771                                 while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1772
1773                                 Double_t dEta = etaTrigger - etaAssoc;
1774
1775                                 Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1776                                 FillHistogram(Form("all_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1777                                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1778                                         FillHistogram(Form("cpv_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1779
1780                                 if ( ph1->IsDispOK() && ph2->IsDispOK() ){
1781                                         FillHistogram(Form("disp_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1782                                         if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1783                                                 FillHistogram(Form("both_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1784                                 }
1785                         }
1786                 }
1787         }
1788 }
1789
1790 //_______________________________________________________________________________
1791 void AliPHOSCorrelations::ConsiderPi0sMix()
1792 {
1793         TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
1794         for(Int_t evi=0; evi<arrayList->GetEntries();evi++)
1795         {
1796                 TClonesArray * mixPHOS = static_cast<TClonesArray*>(arrayList->At(evi));
1797                 for (Int_t i1=0; i1 < fCaloPhotonsPHOS->GetEntriesFast(); i1++)
1798                 {
1799                         AliCaloPhoton * ph1 = (AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1800                         for(Int_t i2=0; i2 < mixPHOS->GetEntriesFast(); i2++)
1801                         {
1802                                 AliCaloPhoton * ph2 = (AliCaloPhoton*)mixPHOS->At(i2) ;
1803                                 TLorentzVector p12  = *ph1  + *ph2;
1804                                 Double_t m=p12.M() ;
1805                                 Double_t pt=p12.Pt() ;
1806                                 Double_t eff = 1./GetEfficiency(pt);
1807
1808                                 int mod1 = ph1->Module() ;
1809                                 int mod2 = ph2->Module() ;
1810
1811                                 FillHistogram("mix_all_mpt", m, pt);
1812                                 FillHistogram("mix_all_mpt_left",m, pt);
1813                                 FillHistogram("mix_all_mpt_right",m, pt);
1814
1815                                 FillHistogram("mix_all_mpt_eff", m, pt, eff);
1816                                 FillHistogram("mix_all_mpt_left_eff",m, pt, eff);
1817                                 FillHistogram("mix_all_mpt_right_eff",m, pt, eff);
1818
1819                                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1820                                 {
1821                                         FillHistogram("mix_cpv_mpt",m, pt);
1822                                         FillHistogram("mix_cpv_mpt_left",m, pt);
1823                                         FillHistogram("mix_cpv_mpt_right",m, pt);
1824
1825                                         FillHistogram("mix_cpv_mpt_eff",m, pt, eff);
1826                                         FillHistogram("mix_cpv_mpt_left_eff",m, pt, eff);
1827                                         FillHistogram("mix_cpv_mpt_right_eff",m, pt, eff);
1828                                 }
1829                                 if ( ph1->IsDispOK() && ph2->IsDispOK() )
1830                                 {
1831                                         FillHistogram("mix_disp_mpt",m, pt);
1832                                         FillHistogram("mix_disp_mpt_left",m, pt);
1833                                         FillHistogram("mix_disp_mpt_right",m, pt);
1834
1835                                         FillHistogram("mix_disp_mpt_eff",m, pt, eff);
1836                                         FillHistogram("mix_disp_mpt_left_eff",m, pt, eff);
1837                                         FillHistogram("mix_disp_mpt_right_eff",m, pt, eff);
1838
1839                                         if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1840                                         {
1841                                                 FillHistogram("mix_both_mpt",m, pt);
1842                                                 FillHistogram("mix_both_mpt_left",m, pt);
1843                                                 FillHistogram("mix_both_mpt_right",m, pt);
1844
1845                                                 FillHistogram("mix_both_mpt_eff",m, pt, eff);
1846                                                 FillHistogram("mix_both_mpt_left_eff",m, pt, eff);
1847                                                 FillHistogram("mix_both_mpt_right_eff",m, pt, eff);
1848
1849                                                 if (mod1 == mod2) // for each module
1850                                                 {
1851                                                         FillHistogram(Form("mix_both%d_mpt",mod1),m, pt);
1852                                                         FillHistogram(Form("mix_both%d_mpt_eff",mod1),m, pt, eff);
1853                                                 }
1854                                         }
1855                                 }
1856                         }
1857                 }
1858         }
1859 }
1860 //_______________________________________________________________________________
1861 void AliPHOSCorrelations::ConsiderTracksMix()
1862 {
1863         TList * arrayList = GetTracksTPCList(fVtxBin, fCentBin, fEMRPBin);
1864         for (Int_t i1=0; i1 < fCaloPhotonsPHOS->GetEntriesFast(); i1++) {
1865                 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
1866                 for (Int_t i2=0; i2<fCaloPhotonsPHOS->GetEntriesFast(); i2++){
1867                         AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
1868                         TLorentzVector p12  = *ph1  + *ph2;
1869                         Double_t phiTrigger=p12.Phi() ;
1870                         Double_t etaTrigger=p12.Eta() ;
1871
1872                         Double_t m=p12.M() ;
1873                         Double_t pt=p12.Pt() ;
1874                         Double_t eff = 1./GetEfficiency(pt);
1875                         Int_t mod1 = ph1->Module();
1876                         Int_t mod2 = ph2->Module();
1877
1878
1879             if(!TestMass(m,pt)) continue;
1880
1881                         for(Int_t evi=0; evi<arrayList->GetEntries();evi++){
1882                                 TClonesArray * mixTracks = static_cast<TClonesArray*>(arrayList->At(evi));
1883                                 for(Int_t i3=0; i3<mixTracks->GetEntriesFast(); i3++){
1884                                         TLorentzVector * track = (TLorentzVector*)mixTracks->At(i3);
1885
1886                                         Double_t phiAssoc = track->Phi();
1887                                         Double_t etaAssoc = track->Eta();
1888                     Double_t ptAssoc =  track->Pt();
1889
1890                                         Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1891
1892                                         Double_t dPhi = phiTrigger - phiAssoc;
1893                                         while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1894                                         while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1895
1896                     Double_t dEta = etaTrigger - etaAssoc;
1897
1898                                         FillHistogram(Form("mix_all_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1899                                         FillHistogram(Form("mix_all_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1900                     FillHistogram(Form("mix_all_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1901
1902                                         if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1903                                         {
1904                                                 FillHistogram(Form("mix_cpv_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1905                                                 FillHistogram(Form("mix_cpv_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1906                                                 FillHistogram(Form("mix_cpv_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1907                                         }
1908
1909                                         if ( ph1->IsDispOK() && ph2->IsDispOK() ) {
1910                                                 FillHistogram(Form("mix_disp_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1911                                                 FillHistogram(Form("mix_disp_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1912                                                 FillHistogram(Form("mix_disp_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1913                                                 if ( ph1->IsCPVOK() && ph2->IsCPVOK() )
1914                                                 {
1915                                                         FillHistogram(Form("mix_both_ptphieta_ptAssoc_%3.1f",ptAssocBin), pt, dPhi, dEta, eff);
1916                                                         FillHistogram(Form("mix_both_ptphieta_ptAssoc_%3.1f_mod%i",ptAssocBin, GetModCase(mod1, mod2)), pt, dPhi, dEta, eff);
1917                                                         FillHistogram(Form("mix_both_ptphieta_ptAssoc_%3.1f_tpc%i",ptAssocBin, CheckTriggerEta(etaTrigger)), pt, dPhi, dEta, eff);
1918                                                 }
1919                                         }
1920                                 }
1921                         }
1922                 }
1923         }
1924 }
1925
1926 //_______________________________________________________________________________
1927 void AliPHOSCorrelations::ConsiderPi0sME()
1928 {
1929         TString spid[4]={"all","cpv","disp","both"} ;
1930     // Counting number of trigger particles.
1931         for (int ipid = 0; ipid < 4; ipid++)
1932         {
1933                 if (fMEExists[ipid])
1934             FillHistogram(Form("nTrigger_%s", spid[ipid].Data()), GetMEPt(ipid), 1./GetEfficiency(GetMEPt(ipid)));
1935         }
1936
1937     // Take track's angles and compare with trigger's angles.
1938     for(Int_t i3=0; i3<fTracksTPC->GetEntriesFast(); i3++){
1939                 TLorentzVector * track = (TLorentzVector*)fTracksTPC->At(i3);
1940
1941                 Double_t phiAssoc = track->Phi();
1942                 Double_t etaAssoc = track->Eta();
1943                 Double_t ptAssoc = track->Pt();
1944
1945                 Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1946                 Double_t dPhi(0.), dEta(0.);
1947
1948                 for (int ipid = 0; ipid < 4; ipid++)
1949                 {
1950             if (GetMEExists(ipid))
1951                         {
1952                 dPhi = GetMEPhi(ipid) - phiAssoc;
1953                                 while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1954                                 while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1955                 dEta = GetMEEta(ipid) - etaAssoc;
1956                 FillHistogram(Form("%s_ptphieta_ptAssoc_%3.1f", spid[ipid].Data(), ptAssocBin),  GetMEPt(ipid), dPhi, dEta, 1./GetEfficiency(GetMEPt(ipid)) );
1957                         }
1958                 }
1959         }
1960 }
1961 //_______________________________________________________________________________
1962 void AliPHOSCorrelations::ConsiderTracksMixME()
1963 {
1964         TString spid[4]={"all","cpv","disp","both"} ;
1965
1966         TList * arrayList = GetTracksTPCList(fVtxBin, fCentBin, fEMRPBin);
1967
1968         for(Int_t evi=0; evi<arrayList->GetEntries();evi++){
1969                 TClonesArray * mixTracks = static_cast<TClonesArray*>(arrayList->At(evi));
1970                 for(Int_t i3=0; i3<mixTracks->GetEntriesFast(); i3++){
1971                         TLorentzVector * track = (TLorentzVector*)mixTracks->At(i3);
1972
1973                         Double_t phiAssoc = track->Phi();
1974                         Double_t etaAssoc = track->Eta();
1975                         Double_t ptAssoc =  track->Pt();
1976
1977                         Double_t ptAssocBin=GetAssocBin(ptAssoc) ;
1978
1979             Double_t ptTrigger (0.);
1980
1981                         Double_t dPhi(0.), dEta(0.);
1982
1983                         for (int ipid = 0; ipid < 4; ipid++)
1984                         {
1985                 if (GetMEExists(ipid))
1986                                 {
1987                     dPhi = GetMEPhi(ipid) - phiAssoc;
1988                                         while (dPhi > 1.5*TMath::Pi()) dPhi-=2*TMath::Pi();
1989                                         while (dPhi < -.5*TMath::Pi()) dPhi+=2*TMath::Pi();
1990                     dEta = GetMEEta(ipid) - etaAssoc;
1991                     ptTrigger = GetMEPt(ipid);
1992
1993                     FillHistogram(Form("mix_%s_ptphieta_ptAssoc_%3.1f", spid[ipid].Data(), ptAssocBin), ptTrigger, dPhi, dEta, 1./GetEfficiency(ptTrigger));
1994                     FillHistogram(Form("mix_%s_ptphieta_ptAssoc_%3.1f_mod%i", spid[ipid].Data(), ptAssocBin, GetMEModCase(ipid)), ptTrigger, dPhi, dEta, 1./GetEfficiency(ptTrigger));
1995                     FillHistogram(Form("mix_%s_ptphieta_ptAssoc_%3.1f_tpc%i", spid[ipid].Data(), ptAssocBin, CheckTriggerEta(GetMEEta(ipid))), ptTrigger, dPhi, dEta, 1./GetEfficiency(ptTrigger));
1996                                 }
1997                         }
1998                 }
1999         }
2000 }
2001
2002 //_______________________________________________________________________________
2003 TList* AliPHOSCorrelations::GetCaloPhotonsPHOSList(UInt_t vtxBin, UInt_t centBin, UInt_t rpBin){
2004
2005   int offset = vtxBin * GetNumberOfCentralityBins() * fNEMRPBins + centBin * fNEMRPBins + rpBin;
2006   if( fCaloPhotonsPHOSLists->At(offset) ) {
2007     TList* list = dynamic_cast<TList*> (fCaloPhotonsPHOSLists->At(offset));
2008     return list;
2009   }
2010   else{ // no list for this bin has been created, yet
2011     TList* list = new TList();
2012     fCaloPhotonsPHOSLists->AddAt(list, offset);
2013     return list;
2014   }
2015 }
2016 //_______________________________________________________________________________
2017 TList* AliPHOSCorrelations::GetTracksTPCList(UInt_t vtxBin, UInt_t centBin, UInt_t rpBin){
2018
2019   int offset = vtxBin * GetNumberOfCentralityBins() * fNEMRPBins + centBin * fNEMRPBins + rpBin;
2020   if( fTracksTPCLists->At(offset) ) { // list exists
2021      TList* list = dynamic_cast<TList*> (fTracksTPCLists->At(offset));
2022      return list;
2023   }
2024   else { // no list for this bin has been created, yet
2025     TList* list = new TList();
2026     fTracksTPCLists->AddAt(list, offset);
2027     return list;
2028   }
2029 }
2030 //_______________________________________________________________________________
2031 Double_t AliPHOSCorrelations::GetAssocBin(Double_t pt) const
2032 {
2033   //Calculates bin of associated particle pt.
2034   for(Int_t i=1; i<fAssocBins.GetSize(); i++){
2035     if(pt>fAssocBins.At(i-1) && pt<fAssocBins.At(i))
2036       return fAssocBins.At(i) ;
2037   }
2038   return fAssocBins.At(fAssocBins.GetSize()-1) ;
2039 }
2040 //_______________________________________________________________________________
2041 void AliPHOSCorrelations::FillTrackEtaPhi()
2042 {
2043         // Distribution TPC's tracks by angles.
2044         for (Int_t i1=0; i1<fTracksTPC->GetEntriesFast(); i1++){
2045                 TLorentzVector * track = (TLorentzVector*)fTracksTPC->At(i1);
2046                 FillHistogram( "track_phieta", track->Phi(), track->Eta() );
2047         }
2048 }
2049
2050 //_______________________________________________________________________________
2051 void AliPHOSCorrelations::UpdatePhotonLists()
2052 {
2053   //Now we either add current events to stack or remove
2054   //If no photons in current event - no need to add it to mixed
2055
2056   TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
2057   if( fDebug >= 2 )
2058     AliInfo( Form("fCentBin=%d, fCentNMixed[]=%d",fCentBin,fCentNMixed[fCentBin]) );
2059   if(fCaloPhotonsPHOS->GetEntriesFast()>0)
2060   {
2061     arrayList->AddFirst(fCaloPhotonsPHOS) ;
2062     fCaloPhotonsPHOS=0x0;
2063     if(arrayList->GetEntries() > fCentNMixed[fCentBin])
2064     { // Remove redundant events
2065       TClonesArray * tmp = static_cast<TClonesArray*>(arrayList->Last()) ;
2066       arrayList->RemoveLast() ;
2067       delete tmp;
2068     }
2069   }
2070 }
2071 //_______________________________________________________________________________
2072 void AliPHOSCorrelations::UpdateTrackLists()
2073 {
2074   //Now we either add current events to stack or remove
2075   //If no photons in current event - no need to add it to mixed
2076
2077   TList * arrayList = GetTracksTPCList(fVtxBin, fCentBin, fEMRPBin);
2078
2079   if( fDebug >= 2 )
2080     AliInfo( Form("fCentBin=%d, fCentNMixed[]=%d",fCentBin,fCentNMixed[fCentBin]) );
2081   if(fTracksTPC->GetEntriesFast()>0)
2082   {
2083
2084     arrayList->AddFirst(fTracksTPC) ;
2085     fTracksTPC=0x0;
2086     if(arrayList->GetEntries() > fCentNMixed[fCentBin])
2087     { // Remove redundant events
2088       TClonesArray * tmp = static_cast<TClonesArray*>(arrayList->Last()) ;
2089       arrayList->RemoveLast() ;
2090       delete tmp;
2091     }
2092   }
2093 }
2094 //_______________________________________________________________________________
2095 Bool_t AliPHOSCorrelations::SelectESDTrack(AliESDtrack * t) const
2096 // Estimate if this track can be used for the RP calculation. If all right - return "TRUE"
2097 {
2098         Float_t pt=t->Pt();
2099         if(pt<0.5 || pt>20.) return kFALSE ;
2100         if(fabs( t->Eta() )>0.8) return kFALSE;
2101         if(!fESDtrackCuts->AcceptTrack(t)) return kFALSE ;
2102         return kTRUE ;
2103 }
2104 //_______________________________________________________________________________
2105 Bool_t AliPHOSCorrelations::SelectAODTrack(AliAODTrack * t) const
2106 // Estimate if this track can be used for the RP calculation. If all right - return "TRUE"
2107 {
2108         Float_t pt=t->Pt();
2109         if(pt<0.5 || pt>20.) return kFALSE ;
2110         if(fabs( t->Eta() )>0.8) return kFALSE;
2111         if(fCheckHibridGlobal == kOnlyHibridTracks)
2112         {
2113                 if(!t->IsHybridGlobalConstrainedGlobal())
2114                         return kFALSE ;
2115         }
2116
2117         if (fCheckHibridGlobal == kWithOutHibridTracks)
2118         {
2119                 if(t->IsHybridGlobalConstrainedGlobal())
2120                         return kFALSE ;
2121         }
2122
2123         return kTRUE ;
2124 }
2125 //_______________________________________________________________________________
2126 void AliPHOSCorrelations::LogProgress(int step)
2127 // Fill "step by step" hist
2128 {
2129   //FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
2130   FillHistogram("hTotSelEvents", step+0.5);
2131 }
2132 //_______________________________________________________________________________
2133 void AliPHOSCorrelations::LogSelection(int step, int internalRunNumber)
2134 {
2135   // the +0.5 is not realy neccisarry, but oh well... -henrik
2136   FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
2137   //FillHistogram("hTotSelEvents", step+0.5);
2138 }
2139 //_______________________________________________________________________________
2140 Bool_t AliPHOSCorrelations::TestMass(Double_t m, Double_t pt)
2141 {
2142         //Check if mair in pi0 peak window
2143         //To make pT-dependent
2144         if (!fSigmaWidth)       // Default big window
2145         {
2146                 FillHistogram("massWindow", fMassInvMean, fMassInvSigma);
2147         if(fMassInvMean-fMassInvSigma<m && m<fMassInvMean+fMassInvSigma)
2148         {
2149             FillHistogram("massWindowPass", 1);
2150             return true;
2151         }
2152         else
2153         {
2154             FillHistogram("massWindowPass", 2);
2155             return false;
2156         }
2157         }
2158         else                    // Parametrization
2159         {
2160                 FillHistogram("massWindow", MassMeanFunktion(pt), MassSigmaFunktion(pt)*fSigmaWidth);
2161         if ( MassMeanFunktion(pt)-MassSigmaFunktion(pt)*fSigmaWidth<m && m<MassMeanFunktion(pt)+MassSigmaFunktion(pt)*fSigmaWidth )
2162         {
2163             FillHistogram("massWindowPass", 3);
2164             return true;
2165         }
2166         else
2167         {
2168             FillHistogram("massWindowPass", 4);
2169             return false;
2170         }
2171         }
2172 }
2173 //_______________________________________________________________________________
2174 Double_t AliPHOSCorrelations::MassMeanFunktion(Double_t &pt) const
2175 {
2176         // Parametrization mean of mass window
2177         return ( fMassMeanP1+TMath::Power(1.25,-pt+fMassMeanP0) );
2178 }
2179 //_______________________________________________________________________________
2180 Double_t AliPHOSCorrelations::MassSigmaFunktion(Double_t &pt) const
2181 {
2182         // Parametrization sigma of mass window
2183     //TODO:: Fix falling at large pT.
2184         return ( fabs(fMassSigmaP0 + fMassSigmaP1*pt) );
2185 }
2186 //_____________________________________________________________________________
2187 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x)const{
2188   //FillHistogram
2189   TH1 * hist = dynamic_cast<TH1*>(fOutputContainer->FindObject(key)) ;
2190   if(hist)
2191     hist->Fill(x) ;
2192   else
2193     AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2194 }
2195 //_____________________________________________________________________________
2196 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x,Double_t y)const{
2197   //FillHistogram
2198   TH1 * th1 = dynamic_cast<TH1*> (fOutputContainer->FindObject(key));
2199   if(th1)
2200     th1->Fill(x, y) ;
2201   else
2202     AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2203 }
2204
2205 //_____________________________________________________________________________
2206 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
2207   //Fills 1D histograms with key
2208   TObject * obj = fOutputContainer->FindObject(key);
2209
2210   TH2 * th2 = dynamic_cast<TH2*> (obj);
2211   if(th2) {
2212     th2->Fill(x, y, z) ;
2213     return;
2214   }
2215
2216   TH3 * th3 = dynamic_cast<TH3*> (obj);
2217   if(th3) {
2218     th3->Fill(x, y, z) ;
2219     return;
2220   }
2221
2222   AliError(Form("can not find histogram (of instance TH2) <%s> ",key)) ;
2223 }
2224 //_____________________________________________________________________________
2225 void AliPHOSCorrelations::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z, Double_t w) const{
2226   //Fills 1D histograms with key
2227   TObject * obj = fOutputContainer->FindObject(key);
2228
2229   TH3 * th3 = dynamic_cast<TH3*> (obj);
2230   if(th3) {
2231     th3->Fill(x, y, z, w) ;
2232     return;
2233   }
2234
2235   AliError(Form("can not find histogram (of instance TH3) <%s> ",key)) ;
2236 }
2237 //_____________________________________________________________________________
2238 void AliPHOSCorrelations::SetGeometry()
2239 {
2240   // Initialize the PHOS geometry
2241   //Init geometry
2242   if(!fPHOSGeo){
2243      AliOADBContainer geomContainer("phosGeo");
2244      geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
2245      TObjArray *matrixes = (TObjArray*)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
2246      fPHOSGeo =  AliPHOSGeometry::GetInstance("IHEP") ;
2247      for(Int_t mod=0; mod<5; mod++) {
2248         if(!matrixes->At(mod)) {
2249           if( fDebug )
2250             AliInfo(Form("No PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2251           continue;
2252         }
2253         else {
2254              fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod) ;
2255              if( fDebug >1 )
2256                AliInfo(Form("Adding PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2257         }
2258      }
2259   }
2260 }
2261 //_____________________________________________________________________________
2262 Double_t AliPHOSCorrelations::GetEfficiency(Double_t x) const {
2263         //Efficiency for Both2core only!
2264
2265         Double_t e =1.;
2266         // From 0 to 5 - 11h for different centrality.
2267         /*0: 0-5%
2268         1: 5-10%
2269         2: 10-20%
2270         3: 20-40%
2271         4: 40-60%
2272         5: 60-80%
2273         6: 0-20%
2274         7: 0-10%*/
2275         Double_t par0[9] = {-798863, 339.714, 6407.1, -457.778, 1283.65, -117.075, -19.3764, 0, 0};
2276         Double_t par1[9] = {-799344, -1852.1, 3326.29, -384.229, 504.046, 562.608, 130.518, 0, 0};
2277         Double_t par2[9] = {-858904, -1923.28, 5350.74, -568.946, 945.497, 419.647, 101.911, 0, 0};
2278         Double_t par3[9] = {-795652, -1495.97, 2926.46, -357.804, 478.961, 551.127, 128.86, 0, 0};
2279         Double_t par4[9] = {-891951, 279626, -123110, -5464.75, 27470.8, 283264, 15355.1, 192762, 44828.6};
2280         Double_t par5[9] = {-1.1094e+06, -986.915, 2127.71, -268.908, 375.594, 380.791, 89.4053, 0, 0};
2281         // Double_t par6[7] = {4.86106e+09, 4.47013e+08, -1.48079e+09, 1.47233e+08, -2.62356e+08, -1.00639e+08, -2.45629e+07, 0, 0};
2282         // Double_t par7[7] = {-1.36243e+06, -26011.1, 135838, -12161.3, 24956.8, 4985.4, 1285.57, 0, 0};
2283
2284         // 8 for pPb13 and 0-100%
2285          Double_t par8[9] = {6.87095e+06, 8.36553e+06, -3.29572e+06, 2.18688e+06, -739490, 521666, 106661, 0, 0};
2286
2287
2288         Double_t* pFitPoint;
2289
2290         if(fPeriod == kLHC11h)
2291         {
2292                 if(x<1.) x = 1.;
2293
2294                 if (fCentrality<=5)  pFitPoint = &par0[0];
2295                 if (fCentrality>5 && fCentrality<=10) pFitPoint = &par1[0];
2296                 if (fCentrality>10 && fCentrality<=20) pFitPoint = &par2[0];
2297                 if (fCentrality>20 && fCentrality<=40) pFitPoint = &par3[0];
2298                 if (fCentrality>40 && fCentrality<=60) pFitPoint = &par4[0];
2299                 if (fCentrality>60) pFitPoint = &par5[0];
2300
2301                 Double_t pFit[9];
2302                 for (int i = 0; i < 9; ++i)
2303                 {
2304                         pFit[i] = *(pFitPoint+i);
2305                 }
2306
2307                 if (fCentrality>40 && fCentrality<=60)
2308                         e = TMath::Exp(-(((((1.+(pFit[1]*x))+(pFit[2]*(x*x)))+(pFit[5]*(x*(x*x))))+(pFit[7]*(x*(x*(x*x)))))/((((pFit[3]*x)+(pFit[4]*(x*x)))+(pFit[6]*(x*(x*x))))+(pFit[8]*(x*(x*(x*x))))))) ;
2309                 else
2310                         e = TMath::Exp(-((((1.+(pFit[1]*x))+(pFit[2]*(x*x)))+(pFit[5]*(x*(x*x))))/(((pFit[3]*x)+(pFit[4]*(x*x)))+(pFit[6]*(x*(x*x)))))) ;
2311         }
2312         else
2313         if( fPeriod == kLHC13 )
2314         {
2315                 pFitPoint = &par8[0];
2316                 Double_t pFit[9];
2317                 for( int i = 0; i < 9; i++ )
2318                 {
2319                         pFit[i] = *(pFitPoint+i);
2320                 }
2321
2322                 e = TMath::Exp(-((((pFit[0]+(pFit[1]*x))+(pFit[2]*(x*x)))+(pFit[5]*(x*(x*x))))/(((1.+(pFit[3]*x))+(pFit[4]*(x*x)))+(pFit[6]*(x*(x*x)))))) ;
2323         }
2324         else
2325         {
2326                 // No case
2327                 AliWarning(Form("No efficiensy choise."));
2328                 e = 1.;
2329         }
2330
2331         return e;
2332         // return 1.; // For test.
2333 }
2334 //_____________________________________________________________________________
2335 Int_t AliPHOSCorrelations::GetModCase(Int_t &mod1, Int_t &mod2) const {
2336
2337         // Return modules pair namber.
2338         if(mod1 == mod2)
2339         {
2340                 if(mod1 == 1) return 1;
2341                 if(mod1 == 2) return 2;
2342                 if(mod1 == 3) return 3;
2343         }
2344         else
2345         {
2346                 if(mod1 == 1 || mod2 == 1)
2347                         if(mod1 == 2 || mod2 == 2)
2348                                 return 12;
2349
2350                 if(mod1 == 1 || mod2 == 1)
2351                         if(mod1 == 3 || mod2 == 3)
2352                                 return 13;
2353                 if(mod1 == 2 || mod2 == 2)
2354                         if(mod1 == 3 || mod2 == 3)
2355                                 return 23;
2356         }
2357
2358         AliError(Form("No choise for mod1 = %i, mod2 = %i", mod1, mod2));
2359         return 1;
2360 }
2361 //_____________________________________________________________________________
2362 Int_t AliPHOSCorrelations::CheckTriggerEta(Double_t eta){
2363         if (eta>=0.)
2364                 return 1;
2365         return 2;
2366 }
2367 //_____________________________________________________________________________
2368 void AliPHOSCorrelations::TestPi0ME(Int_t ipid, TLorentzVector p12, Int_t modCase)
2369 {
2370         Double_t phiTrigger=p12.Phi() ;
2371         Double_t etaTrigger=p12.Eta() ;
2372         Double_t pt=p12.Pt() ;
2373
2374     if ( GetMEExists(ipid) )
2375     {
2376         if ( pt>=GetMEPt(ipid) )
2377         {
2378             SetMEPt(ipid,pt);
2379             SetMEPhi(ipid, phiTrigger);
2380             SetMEEta(ipid, etaTrigger);
2381             SetMEModCase(ipid, modCase);
2382         }
2383     }
2384     else
2385     {
2386         SetMEPt(ipid,pt);
2387         SetMEPhi(ipid, phiTrigger);
2388         SetMEEta(ipid, etaTrigger);
2389         SetMEModCase(ipid, modCase);
2390         SetMEExists(ipid);
2391     }
2392 }
2393 //_____________________________________________________________________________
2394 void AliPHOSCorrelations::ZeroingVariables(){
2395         // Set Phi, Eta, pT, modNumber andtrigger variable of moust energetic trigger particle to zero.
2396         for (int i = 0; i < 4; ++i)
2397         {
2398                 fMEExists[i] = false;
2399                 fMEPhi[i] = fMEEta[i] = fMEPt[i] = -99;
2400                 fMEModCase[i] = 1;
2401         }
2402 }
2403