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