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