]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/PHOSTasks/PHOS_PbPb/AliAnalysisTaskPi0Flow.cxx
Added AlidNdPtAnalysisPbPb2011 to libPWGLFspectra. Changes by Phillip Luettig.
[u/mrichter/AliRoot.git] / PWGGA / PHOSTasks / PHOS_PbPb / AliAnalysisTaskPi0Flow.cxx
CommitLineData
b8bea077 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
2cde7444 16#include "TChain.h"
17#include "TTree.h"
18#include "TObjArray.h"
19#include "TF1.h"
20#include "TFile.h"
21#include "TH1F.h"
22#include "TH2F.h"
23#include "TH2I.h"
24#include "TH3F.h"
25#include "TParticle.h"
26#include "TCanvas.h"
27#include "TStyle.h"
28#include "TRandom.h"
9ed0694f 29#include "THashList.h"
2cde7444 30
31#include "AliAnalysisManager.h"
32#include "AliMCEventHandler.h"
33#include "AliMCEvent.h"
34#include "AliStack.h"
35#include "AliAnalysisTaskSE.h"
36#include "AliAnalysisTaskPi0Flow.h"
37#include "AliCaloPhoton.h"
38#include "AliPHOSGeometry.h"
ada6b882 39#include "TGeoManager.h"
2cde7444 40#include "AliPHOSEsdCluster.h"
41#include "AliPHOSCalibData.h"
42#include "AliESDEvent.h"
43#include "AliESDCaloCells.h"
44#include "AliESDVertex.h"
45#include "AliESDtrackCuts.h"
ada6b882 46#include "AliAODEvent.h"
2cde7444 47#include "AliLog.h"
48#include "AliPID.h"
49#include "AliCDBManager.h"
ada6b882 50#include <AliAODCaloCluster.h>
51#include "AliCentrality.h"
2cde7444 52#include "AliESDtrackCuts.h"
53#include "AliEventplane.h"
54#include "TProfile.h"
55#include "AliOADBContainer.h"
b509e3a5 56#include "AliPHOSEPFlattener.h"
2cde7444 57
ada6b882 58// Analysis task to fill histograms with PHOS ESD or AOD clusters and cells
59// Authors : Dmitri Peressounko
60// Date : 28.05.2011
61// Modified: 03.08.2012 Henrik Qvigstad
91b112bd 62/* $Id$ */
2cde7444 63
ada6b882 64ClassImp(AliAnalysisTaskPi0Flow);
2cde7444 65
66//________________________________________________________________________
67Double_t rnlin(Double_t *x, Double_t * /*par*/)
68{
69 //a = par[0], b = par[1].
70 //1+a*exp(-e/b)
71
ada6b882 72// return 0.0241+1.0504*x[0]+0.000249*x[0]*x[0] ;
73 return 1.015*(0.0241+1.0504*x[0]+0.000249*x[0]*x[0]) ;
2cde7444 74
75}
76
77//________________________________________________________________________
ada6b882 78AliAnalysisTaskPi0Flow::AliAnalysisTaskPi0Flow(const char *name, Period period)
2cde7444 79: AliAnalysisTaskSE(name),
ada6b882 80 fCentEdges(10),
107849ac 81 fCentNMixed(10),
ada6b882 82 fNEMRPBins(9),
ada6b882 83 fPeriod(period),
2cde7444 84 fOutputContainer(0x0),
2cde7444 85 fNonLinCorr(0),
ada6b882 86 fEvent(0x0),
87 fEventESD(0x0),
88 fEventAOD(0x0),
89 fMCStack(0x0),
90 fRunNumber(-999),
91 fInternalRunNumber(0),
92 fPHOSGeo(0),
93 fMultV0(0x0),
94 fV0Cpol(0.),fV0Apol(0.),
95 fESDtrackCuts(0x0),
96 fPHOSCalibData(0x0),
b509e3a5 97 fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
98 fTPCFlat(0x0),
99 fV0AFlat(0x0),
100 fV0CFlat(0x0),
ada6b882 101 fVertexVector(),
102 fVtxBin(0),
103 fCentralityV0M(0.),
104 fCentBin(0),
105 fHaveTPCRP(0),
2cde7444 106 fRP(0),
107 fRPV0A(0),
108 fRPV0C(0),
ada6b882 109 fEMRPBin(0),
110 fCaloPhotonsPHOS(0x0),
111 fCaloPhotonsPHOSLists(0x0)
2cde7444 112{
bba0c9ef 113 const int nbins = 9;
114 Double_t edges[nbins+1] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80.};
115 fCentEdges = TArrayD(nbins+1, edges);
116 Int_t nMixed[nbins] = {4,4,6,10,20,30,50,100,100};
117 fCentNMixed = TArrayI(nbins, nMixed);
ada6b882 118
119 for(Int_t i=0;i<kNCenBins;i++){
2cde7444 120 for(Int_t j=0;j<2; j++)
ada6b882 121 for(Int_t k=0; k<2; k++) {
122 fMeanQ[i][j][k]=0.;
123 fWidthQ[i][j][k]=0.;
124 }
2cde7444 125 }
126
ada6b882 127 fVertex[0]=0; fVertex[1]=0; fVertex[2]=0;
128
2cde7444 129 // Output slots #0 write into a TH1 container
130 DefineOutput(1,TList::Class());
131
ada6b882 132
133
2cde7444 134 // Set bad channel map
135 char key[55] ;
136 for(Int_t i=0; i<6; i++){
53c2a02a 137 snprintf(key,55,"PHOS_BadMap_mod%d",i) ;
2cde7444 138 fPHOSBadMap[i]=new TH2I(key,"Bad Modules map",64,0.,64.,56,0.,56.) ;
139 }
2cde7444 140
141
2cde7444 142 // Initialize non-linrarity correction
143 fNonLinCorr = new TF1("nonlib",rnlin,0.,40.,0);
144
145
146}
ada6b882 147//___________________________________________________________________________
148AliAnalysisTaskPi0Flow::~AliAnalysisTaskPi0Flow()
149{
150 delete fNonLinCorr;
151 delete fESDtrackCuts;
152 delete fPHOSCalibData;
153 delete fCaloPhotonsPHOSLists;
b509e3a5 154 if(fTPCFlat)delete fTPCFlat; fTPCFlat=0x0;
155 if(fV0AFlat)delete fV0AFlat; fV0AFlat=0x0;
156 if(fV0CFlat)delete fV0CFlat; fV0CFlat=0x0;
157
ada6b882 158}
2cde7444 159//________________________________________________________________________
160void AliAnalysisTaskPi0Flow::UserCreateOutputObjects()
161{
162 // Create histograms
163 // Called once
164 const Int_t nRuns=200 ;
ada6b882 165
166 // histograms
2cde7444 167 if(fOutputContainer != NULL){
168 delete fOutputContainer;
169 }
9ed0694f 170 fOutputContainer = new THashList();
2cde7444 171 fOutputContainer->SetOwner(kTRUE);
ada6b882 172
2cde7444 173 //========QA histograms=======
174
175 //Event selection
ada6b882 176 fOutputContainer->Add(new TH2F("hSelEvents","Event selection", 12,0.,13.,nRuns,0.,float(nRuns))) ;
177 fOutputContainer->Add(new TH1F("hTotSelEvents","Event selection", 12,0.,12.)) ;
178
2cde7444 179 //vertex distribution
180 fOutputContainer->Add(new TH2F("hZvertex","Z vertex position", 50,-25.,25.,nRuns,0.,float(nRuns))) ;
ada6b882 181
2cde7444 182 //Centrality
183 fOutputContainer->Add(new TH2F("hCentrality","Event centrality", 100,0.,100.,nRuns,0.,float(nRuns))) ;
184 fOutputContainer->Add(new TH2F("hCenPHOS","Centrality vs PHOSclusters", 100,0.,100.,200,0.,200.)) ;
185 fOutputContainer->Add(new TH2F("hCenPHOSCells","Centrality vs PHOS cells", 100,0.,100.,100,0.,1000.)) ;
ada6b882 186 fOutputContainer->Add(new TH2F("hCenTrack","Centrality vs tracks", 100,0.,100.,100,0.,15000.)) ;
2cde7444 187 fOutputContainer->Add(new TH2F("hCluEvsClu","ClusterMult vs E",200,0.,10.,100,0.,100.)) ;
ada6b882 188
189
2cde7444 190 //Reaction plane
191 fOutputContainer->Add(new TH3F("hPHOSphi","cos" ,10,0.,100.,20,0.,10.,100,-TMath::Pi(),TMath::Pi()));
ada6b882 192
193 fOutputContainer->Add(new TH2F("cos2AC","RP correlation between TPC subs", 100,-1.,1.,20,0.,100.)) ;
194 fOutputContainer->Add(new TH2F("cos2V0AC","RP correlation between VO A and C sides", 100,-1.,1.,20,0.,100.)) ;
195 fOutputContainer->Add(new TH2F("cos2V0ATPC","RP correlation between TPC and V0A", 100,-1.,1.,20,0.,100.)) ;
196 fOutputContainer->Add(new TH2F("cos2V0CTPC","RP correlation between TPC and V0C", 100,-1.,1.,20,0.,100.)) ;
197
198 fOutputContainer->Add(new TH2F("phiRP","RP distribution with TPC", 100,0.,TMath::Pi(),20,0.,100.)) ;
199 fOutputContainer->Add(new TH2F("phiRPflat","RP distribution with TPC flat", 100,0.,TMath::Pi(),20,0.,100.)) ;
200 fOutputContainer->Add(new TH2F("phiRPV0A","RP distribution with V0A", 100,0.,TMath::Pi(),20,0.,100.)) ;
201 fOutputContainer->Add(new TH2F("phiRPV0C","RP distribution with V0C", 100,0.,TMath::Pi(),20,0.,100.)) ;
f5e08ac9 202 fOutputContainer->Add(new TH3F("phiRPV0AC","RP distribution with V0A and V0C", 100,0.,TMath::Pi(),100,0.,TMath::Pi(),20,0.,100.)) ;
ada6b882 203 fOutputContainer->Add(new TH2F("phiRPV0Aflat","RP distribution with V0 flat", 100,0.,TMath::Pi(),20,0.,100.)) ;
204 fOutputContainer->Add(new TH2F("phiRPV0Cflat","RP distribution with V0 flat", 100,0.,TMath::Pi(),20,0.,100.)) ;
205 fOutputContainer->Add(new TH3F("phiRPV0ATPC","RP distribution with V0A + TPC", 100,0.,TMath::Pi(),100,0.,TMath::Pi(),20,0.,100.)) ;
206 fOutputContainer->Add(new TH3F("phiRPV0CTPC","RP distribution with V0C + TPC", 100,0.,TMath::Pi(),100,0.,TMath::Pi(),20,0.,100.)) ;
207
208
2cde7444 209 //PHOS QA
210 fOutputContainer->Add(new TH1I("hCellMultEvent" ,"PHOS cell multiplicity per event" ,2000,0,2000));
211 fOutputContainer->Add(new TH1I("hCellMultEventM1","PHOS cell multiplicity per event, M1",2000,0,2000));
212 fOutputContainer->Add(new TH1I("hCellMultEventM2","PHOS cell multiplicity per event, M2",2000,0,2000));
213 fOutputContainer->Add(new TH1I("hCellMultEventM3","PHOS cell multiplicity per event, M3",2000,0,2000));
214
215 fOutputContainer->Add(new TH1F("hCellEnergy" ,"Cell energy" ,3000,0.,30.));
216 fOutputContainer->Add(new TH1F("hCellEnergyM1","Cell energy in module 1",3000,0.,30.));
217 fOutputContainer->Add(new TH1F("hCellEnergyM2","Cell energy in module 2",3000,0.,30.));
218 fOutputContainer->Add(new TH1F("hCellEnergyM3","Cell energy in module 3",3000,0.,30.));
219
220 fOutputContainer->Add(new TH2F("hCellNXZM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
221 fOutputContainer->Add(new TH2F("hCellNXZM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
222 fOutputContainer->Add(new TH2F("hCellNXZM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
223 fOutputContainer->Add(new TH2F("hCellEXZM1","Cell E(X,Z), M1",64,0.5,64.5, 56,0.5,56.5));
224 fOutputContainer->Add(new TH2F("hCellEXZM2","Cell E(X,Z), M2",64,0.5,64.5, 56,0.5,56.5));
225 fOutputContainer->Add(new TH2F("hCellEXZM3","Cell E(X,Z), M3",64,0.5,64.5, 56,0.5,56.5));
ada6b882 226
2cde7444 227 //Bad Map
228 fOutputContainer->Add(new TH2F("hCluLowM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
229 fOutputContainer->Add(new TH2F("hCluLowM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
230 fOutputContainer->Add(new TH2F("hCluLowM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
231
232 fOutputContainer->Add(new TH2F("hCluHighM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
233 fOutputContainer->Add(new TH2F("hCluHighM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
234 fOutputContainer->Add(new TH2F("hCluHighM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
ada6b882 235
2cde7444 236 fOutputContainer->Add(new TH2F("hCluVetoM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
237 fOutputContainer->Add(new TH2F("hCluVetoM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
238 fOutputContainer->Add(new TH2F("hCluVetoM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
239
240 fOutputContainer->Add(new TH2F("hCluDispM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
241 fOutputContainer->Add(new TH2F("hCluDispM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
242 fOutputContainer->Add(new TH2F("hCluDispM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
243
ada6b882 244
2cde7444 245 //Single photon and pi0 spectrum
ada6b882 246 const Int_t nPtPhot = 300 ;
247 const Double_t ptPhotMax = 30 ;
248 const Int_t nM = 500;
249 const Double_t mMin = 0.0;
250 const Double_t mMax = 1.0;
251
2cde7444 252 //PHOS calibration QA
253 fOutputContainer->Add(new TH2F("hPi0M11","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
254 fOutputContainer->Add(new TH2F("hPi0M12","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
255 fOutputContainer->Add(new TH2F("hPi0M13","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
256 fOutputContainer->Add(new TH2F("hPi0M22","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
257 fOutputContainer->Add(new TH2F("hPi0M23","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
258 fOutputContainer->Add(new TH2F("hPi0M33","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
ada6b882 259
260 // Histograms for different centralities
2cde7444 261 char key[55] ;
ada6b882 262 for(Int_t cent=0; cent < fCentEdges.GetSize(); cent++){
2cde7444 263 snprintf(key,55,"hPhotAll_cen%d",cent) ;
264 fOutputContainer->Add(new TH1F(key,"All clusters",nPtPhot,0.,ptPhotMax));
265 snprintf(key,55,"hPhotAllcore_cen%d",cent) ;
266 fOutputContainer->Add(new TH1F(key,"All clusters",nPtPhot,0.,ptPhotMax));
267 snprintf(key,55,"hPhotAllwou_cen%d",cent) ;
268 fOutputContainer->Add(new TH1F(key,"All clusters",nPtPhot,0.,ptPhotMax));
269 snprintf(key,55,"hPhotDisp_cen%d",cent) ;
270 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
271 snprintf(key,55,"hPhotDisp2_cen%d",cent) ;
272 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
273 snprintf(key,55,"hPhotDispcore_cen%d",cent) ;
274 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
275 snprintf(key,55,"hPhotDispwou_cen%d",cent) ;
276 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
277 snprintf(key,55,"hPhotCPV_cen%d",cent) ;
278 fOutputContainer->Add(new TH1F(key,"CPV clusters",nPtPhot,0.,ptPhotMax));
279 snprintf(key,55,"hPhotCPVcore_cen%d",cent) ;
280 fOutputContainer->Add(new TH1F(key,"CPV clusters",nPtPhot,0.,ptPhotMax));
281 snprintf(key,55,"hPhotCPV2_cen%d",cent) ;
282 fOutputContainer->Add(new TH1F(key,"CPV clusters",nPtPhot,0.,ptPhotMax));
283 snprintf(key,55,"hPhotBoth_cen%d",cent) ;
284 fOutputContainer->Add(new TH1F(key,"Both clusters",nPtPhot,0.,ptPhotMax));
285 snprintf(key,55,"hPhotBothcore_cen%d",cent) ;
286 fOutputContainer->Add(new TH1F(key,"Both clusters",nPtPhot,0.,ptPhotMax));
287
288 snprintf(key,55,"hPi0All_cen%d",cent) ;
289 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
290 snprintf(key,55,"hPi0Allcore_cen%d",cent) ;
291 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
292 snprintf(key,55,"hPi0Allwou_cen%d",cent) ;
293 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
294 snprintf(key,55,"hPi0Disp_cen%d",cent) ;
295 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
296 snprintf(key,55,"hPi0Disp2_cen%d",cent) ;
297 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
298 snprintf(key,55,"hPi0Dispcore_cen%d",cent) ;
299 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
300 snprintf(key,55,"hPi0Dispwou_cen%d",cent) ;
301 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
302 snprintf(key,55,"hPi0CPV_cen%d",cent) ;
303 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
304 snprintf(key,55,"hPi0CPVcore_cen%d",cent) ;
305 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
306 snprintf(key,55,"hPi0CPV2_cen%d",cent) ;
307 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
308 snprintf(key,55,"hPi0Both_cen%d",cent) ;
309 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
310 snprintf(key,55,"hPi0Bothcore_cen%d",cent) ;
311 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
312
313 snprintf(key,55,"hPi0All_a07_cen%d",cent) ;
314 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
315 snprintf(key,55,"hPi0Disp_a07_cen%d",cent) ;
316 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
317 snprintf(key,55,"hPi0CPV_a07_cen%d",cent) ;
318 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
319 snprintf(key,55,"hPi0CPV2_a07_cen%d",cent) ;
320 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
321 snprintf(key,55,"hPi0Both_a07_cen%d",cent) ;
ada6b882 322 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
2cde7444 323
324 snprintf(key,55,"hSingleAll_cen%d",cent) ;
325 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
326 snprintf(key,55,"hSingleAllcore_cen%d",cent) ;
327 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
328 snprintf(key,55,"hSingleAllwou_cen%d",cent) ;
329 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
330 snprintf(key,55,"hSingleDisp_cen%d",cent) ;
331 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
332 snprintf(key,55,"hSingleDisp2_cen%d",cent) ;
333 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
334 snprintf(key,55,"hSingleDispcore_cen%d",cent) ;
335 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
336 snprintf(key,55,"hSingleDispwou_cen%d",cent) ;
337 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
338 snprintf(key,55,"hSingleCPV_cen%d",cent) ;
339 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
340 snprintf(key,55,"hSingleCPVcore_cen%d",cent) ;
341 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
342 snprintf(key,55,"hSingleCPV2_cen%d",cent) ;
343 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
344 snprintf(key,55,"hSingleBoth_cen%d",cent) ;
345 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
346 snprintf(key,55,"hSingleBothcore_cen%d",cent) ;
347 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
ada6b882 348
349
2cde7444 350 snprintf(key,55,"hMiPi0All_cen%d",cent) ;
351 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
352 snprintf(key,55,"hMiPi0Allcore_cen%d",cent) ;
353 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
354 snprintf(key,55,"hMiPi0Allwou_cen%d",cent) ;
355 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
356 snprintf(key,55,"hMiPi0Disp_cen%d",cent) ;
357 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
358 snprintf(key,55,"hMiPi0Disp2_cen%d",cent) ;
359 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
360 snprintf(key,55,"hMiPi0Dispwou_cen%d",cent) ;
361 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
362 snprintf(key,55,"hMiPi0Dispcore_cen%d",cent) ;
363 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
364 snprintf(key,55,"hMiPi0CPV_cen%d",cent) ;
365 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
366 snprintf(key,55,"hMiPi0CPVcore_cen%d",cent) ;
367 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
368 snprintf(key,55,"hMiPi0CPV2_cen%d",cent) ;
369 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
370 snprintf(key,55,"hMiPi0Both_cen%d",cent) ;
371 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
372 snprintf(key,55,"hMiPi0Bothcore_cen%d",cent) ;
373 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
ada6b882 374
2cde7444 375 snprintf(key,55,"hMiPi0All_a07_cen%d",cent) ;
376 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
377 snprintf(key,55,"hMiPi0Disp_a07_cen%d",cent) ;
378 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
379 snprintf(key,55,"hMiPi0CPV_a07_cen%d",cent) ;
380 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
381 snprintf(key,55,"hMiPi0CPV2_a07_cen%d",cent) ;
382 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
383 snprintf(key,55,"hMiPi0Both_a07_cen%d",cent) ;
ada6b882 384 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
2cde7444 385
386 snprintf(key,55,"hMiSingleAll_cen%d",cent) ;
387 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
388 snprintf(key,55,"hMiSingleAllwou_cen%d",cent) ;
389 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
390 snprintf(key,55,"hMiSingleAllcore_cen%d",cent) ;
391 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
392 snprintf(key,55,"hMiSingleDisp_cen%d",cent) ;
393 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
394 snprintf(key,55,"hMiSingleDisp2_cen%d",cent) ;
395 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
396 snprintf(key,55,"hMiSingleDispwou_cen%d",cent) ;
397 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
398 snprintf(key,55,"hMiSingleDispcore_cen%d",cent) ;
399 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
400 snprintf(key,55,"hMiSingleCPV_cen%d",cent) ;
401 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
402 snprintf(key,55,"hMiSingleCPVcore_cen%d",cent) ;
403 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
404 snprintf(key,55,"hMiSingleCPV2_cen%d",cent) ;
405 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
406 snprintf(key,55,"hMiSingleBoth_cen%d",cent) ;
407 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
408 snprintf(key,55,"hMiSingleBothcore_cen%d",cent) ;
409 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
410 }
411
2cde7444 412
ada6b882 413
414 const Int_t nPt = 20;
415 const Double_t xPt[21]={0.6,1.,1.5,2.,2.5,3.,3.5,4.,4.5,5.,5.5,6.,7.,8.,9.,10.,12.,14.,16.,18.,20.} ;
416 const Int_t nPhi=10 ;
417 Double_t xPhi[nPhi+1] ;
418 for(Int_t i=0;i<=nPhi;i++)
419 xPhi[i]=i*TMath::Pi() /nPhi ;
420 const Int_t nMm=200 ;
421 Double_t xM[nMm+1] ;
422 for(Int_t i=0;i<=nMm;i++)
423 xM[i]=i*0.5 /nMm;
2cde7444 424
425 char phiTitle[15] ;
426 for(Int_t iRP=0; iRP<3; iRP++){
427 if(iRP==0)
428 snprintf(phiTitle,15,"TPC") ;
429 if(iRP==1)
430 snprintf(phiTitle,15,"V0A") ;
431 if(iRP==2)
432 snprintf(phiTitle,15,"V0C") ;
ada6b882 433 for(Int_t cent=0; cent<fCentEdges.GetSize(); cent++){
434 snprintf(key,55,"hPhotPhi%sAll_cen%d",phiTitle,cent) ;
435 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
436 snprintf(key,55,"hPhotPhi%sAllcore_cen%d",phiTitle,cent) ;
437 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
438 snprintf(key,55,"hPhotPhi%sDisp_cen%d",phiTitle,cent) ;
439 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
440 snprintf(key,55,"hPhotPhi%sDispcore_cen%d",phiTitle,cent) ;
441 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
442 snprintf(key,55,"hPhotPhi%sCPV_cen%d",phiTitle,cent) ;
443 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
444 snprintf(key,55,"hPhotPhi%sCPVcore_cen%d",phiTitle,cent) ;
445 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
446 snprintf(key,55,"hPhotPhi%sBoth_cen%d",phiTitle,cent) ;
447 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
448 snprintf(key,55,"hPhotPhi%sBothcore_cen%d",phiTitle,cent) ;
449 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
450
451 //Pions
452 snprintf(key,55,"hMassPt%sAll_cen%d",phiTitle,cent) ;
453 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
454 snprintf(key,55,"hMassPt%sAllcore_cen%d",phiTitle,cent) ;
455 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
456 snprintf(key,55,"hMassPt%sCPV_cen%d",phiTitle,cent) ;
457 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
458 snprintf(key,55,"hMassPt%sCPVcore_cen%d",phiTitle,cent) ;
459 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
460 snprintf(key,55,"hMassPt%sDisp_cen%d",phiTitle,cent) ;
461 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
462 snprintf(key,55,"hMassPt%sDispcore_cen%d",phiTitle,cent) ;
463 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
464 snprintf(key,55,"hMassPt%sBoth_cen%d",phiTitle,cent) ;
465 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
466 snprintf(key,55,"hMassPt%sBothcore_cen%d",phiTitle,cent) ;
467 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
468
469 //Mixed
470 snprintf(key,55,"hMiMassPt%sAll_cen%d",phiTitle,cent) ;
471 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
472 snprintf(key,55,"hMiMassPt%sAllcore_cen%d",phiTitle,cent) ;
473 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
474 snprintf(key,55,"hMiMassPt%sCPV_cen%d",phiTitle,cent) ;
475 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
476 snprintf(key,55,"hMiMassPt%sCPVcore_cen%d",phiTitle,cent) ;
477 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
478 snprintf(key,55,"hMiMassPt%sDisp_cen%d",phiTitle,cent) ;
479 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
480 snprintf(key,55,"hMiMassPt%sDispcore_cen%d",phiTitle,cent) ;
481 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
482 snprintf(key,55,"hMiMassPt%sBoth_cen%d",phiTitle,cent) ;
483 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
484 snprintf(key,55,"hMiMassPt%sBothcore_cen%d",phiTitle,cent) ;
485 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
486 }
487 }
488
489 // Setup photon lists
490 Int_t kapacity = kNVtxZBins * GetNumberOfCentralityBins() * fNEMRPBins;
491 fCaloPhotonsPHOSLists = new TObjArray(kapacity);
492 fCaloPhotonsPHOSLists->SetOwner();
2cde7444 493
494 PostData(1, fOutputContainer);
2cde7444 495}
496
497//________________________________________________________________________
ada6b882 498void AliAnalysisTaskPi0Flow::UserExec(Option_t *)
2cde7444 499{
500 // Main loop, called for each event
501 // Analyze ESD/AOD
ada6b882 502
503
504 // Step 0: Event Objects
505 fEvent = GetEvent();
506 fEventESD = dynamic_cast<AliESDEvent*> (fEvent);
507 fEventAOD = dynamic_cast<AliAODEvent*> (fEvent);
508 fMCStack = GetMCStack();
c9fb2807 509 LogProgress(0);
ada6b882 510
511
512 // Step 1: Run Number, Misalignment Matrix, and Calibration
513 // fRunNumber, fInternalRunNumber, fMultV0, fV0Cpol, fV0Apol, fMeanQ, fWidthQ
514 if( fRunNumber != fEvent->GetRunNumber()) { // Check run number
515 // this should run only at first call of UserExec(),
516 // or if task runs over multiple runs, which should not occur in normal use.
517
518 // if run number has changed, set run variables
519 fRunNumber = fEvent->GetRunNumber();
520 fInternalRunNumber = ConvertToInternalRunNumber(fRunNumber);
521 // then set misalignment and V0 calibration
522 SetGeometry();
523 SetMisalignment();
524 SetV0Calibration();
525 SetESDTrackCuts();
526 SetPHOSCalibData();
b509e3a5 527 SetFlatteningData();
2cde7444 528 }
c9fb2807 529 LogProgress(1);
530 LogSelection(0, fInternalRunNumber);
ada6b882 531
532
533 // Step 2: Vertex
534 // fVertex, fVertexVector, fVtxBin
535 SetVertex();
536 if( RejectEventVertex() ) {
2cde7444 537 PostData(1, fOutputContainer);
ada6b882 538 return; // Reject!
2cde7444 539 }
c9fb2807 540 LogProgress(2);
2cde7444 541
ada6b882 542// Step 3:
543// if(event->IsPileupFromSPD()){
544// PostData(1, fOutputContainer);
545// return; // Reject!
546// }
c9fb2807 547 LogProgress(3);
2cde7444 548
2cde7444 549
ada6b882 550 // Step 4: Centrality
551 // fCentralityV0M, fCentBin
552 SetCentrality();
553 if( RejectCentrality() ){
2cde7444 554 PostData(1, fOutputContainer);
ada6b882 555 return; // Reject!
2cde7444 556 }
c9fb2807 557 LogProgress(4);
2cde7444 558
ada6b882 559
560 // Step 5: Reaction Plane
561 // fHaveTPCRP, fRP, fRPV0A, fRPV0C, fRPBin
f5e08ac9 562 EvalReactionPlane(); //TODO: uncomment this, or at least deal with it
563 EvalV0ReactionPlane(); //TODO: uncomment this, or at least deal with it
564 fEMRPBin = GetRPBin(); //TODO: uncomment this, or at least deal with it
c9fb2807 565 LogProgress(5);
ada6b882 566
567
568 // Step 6: MC
569 // ProcessMC() ;
c9fb2807 570 LogProgress(6);
ada6b882 571
572
573 // Step 7: QA PHOS cells
574 FillPHOSCellQAHists();
c9fb2807 575 LogProgress(7);
ada6b882 576
577
578 // Step 8: Event Photons (PHOS Clusters) selection
579 SelectPhotonClusters();
580 FillSelectedClusterHistograms();
c9fb2807 581 LogProgress(8);
ada6b882 582
583 // Step 9: Consider pi0 (photon/cluster) pairs.
584 ConsiderPi0s();
c9fb2807 585 LogProgress(9);
ada6b882 586
587 // Step 10; Mixing
588 ConsiderPi0sMix();
c9fb2807 589 LogProgress(10);
2cde7444 590
ada6b882 591 // Step 11: Update lists
592 UpdateLists();
c9fb2807 593 LogProgress(11);
ada6b882 594
2cde7444 595
ada6b882 596 // Post output data.
597 PostData(1, fOutputContainer);
598}
2cde7444 599
ada6b882 600//________________________________________________________________________
4c5399e6 601// void AliAnalysisTaskPi0Flow::Terminate(Option_t *)
602// {
603// // Draw result to the screen
604// // Called once at the end of the query
605// // new TCanvas;
606// // TH1 * hTotSelEvents = dynamic_cast<TH1*>(fOutputContainer->FindObject("hTotSelEvents"));
607// // hTotSelEvents->Draw();
608// }
ada6b882 609//________________________________________________________________________
610void AliAnalysisTaskPi0Flow::SetCentralityBinning(const TArrayD& edges)
611{
107849ac 612 // Define centrality bins by their edges
613
ada6b882 614 int last = edges.GetSize()-1;
615 if( edges.At(0) < 0.)
616 AliError("lower edge less then 0");
617 if( 90. < edges.At(last) )
618 AliError("upper edge larger then 90.");
619 for(int i=0; i<last-1; ++i)
620 if(edges.At(i) > edges.At(i+1))
621 AliError("edges are not sorted");
622
623 fCentEdges = edges;
624}
107849ac 625//________________________________________________________________________
626void AliAnalysisTaskPi0Flow::SetNMixedPerCentrality(const TArrayI& nMixed)
627{
628 // Set number of mixed events for all centrality bins
629
630 fCentNMixed = nMixed;
631}
2cde7444 632
ada6b882 633//________________________________________________________________________
634void AliAnalysisTaskPi0Flow::SetPHOSBadMap(Int_t mod, TH2I* badMapHist)
635 {
636 if(fPHOSBadMap[mod])
637 delete fPHOSBadMap[mod];
638
639 fPHOSBadMap[mod]=new TH2I(*badMapHist);
f5e08ac9 640 if(fDebug)
ada6b882 641 AliInfo(Form("Setting Bad Map Histogram %s",fPHOSBadMap[mod]->GetName()));
642 }
2cde7444 643
ada6b882 644//________________________________________________________________________
645Bool_t AliAnalysisTaskPi0Flow::IsGoodChannel(const char * det, Int_t mod, Int_t ix, Int_t iz)
646{
647 //Check if this channel belogs to the good ones
2cde7444 648
ada6b882 649 if(strcmp(det,"PHOS")==0){
650 if(mod>5 || mod<1){
651 AliError(Form("No bad map for PHOS module %d",mod)) ;
652 return kTRUE ;
653 }
654 if(!fPHOSBadMap[mod]){
655 AliError(Form("No Bad map for PHOS module %d, !fPHOSBadMap[mod]",mod)) ;
656 return kTRUE ;
657 }
658 if(fPHOSBadMap[mod]->GetBinContent(ix,iz)>0)
659 return kFALSE ;
660 else
661 return kTRUE ;
2cde7444 662 }
663 else{
ada6b882 664 AliError(Form("Can not find bad channels for detector %s ",det)) ;
2cde7444 665 }
2cde7444 666
ada6b882 667 //Remove 6 noisy channels in run 139036, LHC10h
668 if( 139036 == fRunNumber
669 && mod==1
670 && (ix==9||ix==10||ix==11)
671 && (iz==45 || iz==46))
672 return kFALSE;
2cde7444 673
ada6b882 674 return kTRUE ;
675}
676//_____________________________________________________________________________
677void AliAnalysisTaskPi0Flow::FillPHOSCellQAHists()
678{
e2139d72 679 // Fill cell occupancy per module
680
ada6b882 681 AliVCaloCells * cells = fEvent->GetPHOSCells();
2cde7444 682
ada6b882 683 FillHistogram("hCenPHOSCells",fCentralityV0M,cells->GetNumberOfCells()) ;
684 FillHistogram("hCenTrack",fCentralityV0M,fEvent->GetNumberOfTracks()) ;
2cde7444 685
2cde7444 686
2cde7444 687 Int_t nCellModule[3] = {0,0,0};
688 for (Int_t iCell=0; iCell<cells->GetNumberOfCells(); iCell++) {
689 Int_t cellAbsId = cells->GetCellNumber(iCell);
ada6b882 690 Int_t relId[4] = {0,0,0,0};
2cde7444 691 fPHOSGeo->AbsToRelNumbering(cellAbsId,relId);
692 Int_t mod1 = relId[0];
693 Int_t cellX = relId[2];
694 Int_t cellZ = relId[3] ;
695 Float_t energy = cells->GetAmplitude(iCell);
696 FillHistogram("hCellEnergy",energy);
697 if(mod1==1) {
698 nCellModule[0]++;
699 FillHistogram("hCellEnergyM1",cells->GetAmplitude(iCell));
700 FillHistogram("hCellNXZM1",cellX,cellZ,1.);
701 FillHistogram("hCellEXZM1",cellX,cellZ,energy);
702 }
703 else if (mod1==2) {
704 nCellModule[1]++;
705 FillHistogram("hCellEnergyM2",cells->GetAmplitude(iCell));
706 FillHistogram("hCellNXZM2",cellX,cellZ,1.);
707 FillHistogram("hCellEXZM2",cellX,cellZ,energy);
708 }
709 else if (mod1==3) {
710 nCellModule[2]++;
711 FillHistogram("hCellEnergyM3",cells->GetAmplitude(iCell));
712 FillHistogram("hCellNXZM3",cellX,cellZ,1.);
713 FillHistogram("hCellEXZM3",cellX,cellZ,energy);
714 }
715 }
716 FillHistogram("hCellMultEventM1",nCellModule[0]);
717 FillHistogram("hCellMultEventM2",nCellModule[1]);
718 FillHistogram("hCellMultEventM3",nCellModule[2]);
719
ada6b882 720}
721//_____________________________________________________________________________
722void AliAnalysisTaskPi0Flow::SelectPhotonClusters()
723{
724 // clear (or create) array for holding events photons/clusters
725 if(fCaloPhotonsPHOS)
726 fCaloPhotonsPHOS->Clear();
727 else {
728 fCaloPhotonsPHOS = new TObjArray(200);
729 fCaloPhotonsPHOS->SetOwner();
730 }
731
732
91b112bd 733 AliVCaloCells* cells = dynamic_cast<AliVCaloCells*> (fEvent->GetPHOSCells());
ada6b882 734 for (Int_t i=0; i<fEvent->GetNumberOfCaloClusters(); i++) {
735 AliVCluster *clu = fEvent->GetCaloCluster(i);
736
737 if ( !clu->IsPHOS() || clu->E()< kMinClusterEnergy) continue; // reject cluster
738
739 // check if cell/channel is good.
2cde7444 740 Float_t position[3];
741 clu->GetPosition(position);
742 TVector3 global(position) ;
743 Int_t relId[4] ;
744 fPHOSGeo->GlobalPos2RelId(global,relId) ;
745 Int_t mod = relId[0] ;
746 Int_t cellX = relId[2];
747 Int_t cellZ = relId[3] ;
ada6b882 748 if ( !IsGoodChannel("PHOS",mod,cellX,cellZ) )
749 continue ; // reject if not.
750
751
752 FillHistogram("hCluEvsClu", clu->E(), clu->GetNCells()) ;
753
754 if(clu->GetNCells() < kMinNCells) continue ;
755 if(clu->GetM02() < kMinM02) continue ;
756
757 TLorentzVector lorentzMomentum;
758 Double_t ecore;
91b112bd 759 ecore = CoreEnergy(clu,cells);
ada6b882 760
761 //if ESD, Apply re-Calibreation
762 Double_t origo[3] = {0,0,0}; // don't rely on event vertex, assume (0,0,0)
763 if( fEventESD ) {
764 AliPHOSEsdCluster cluPHOS1( *(AliESDCaloCluster*) (clu) );
91b112bd 765 cluPHOS1.Recalibrate(fPHOSCalibData, dynamic_cast<AliESDCaloCells*> (cells)); // modify the cell energies
ada6b882 766 Reclusterize(&cluPHOS1) ;
767 cluPHOS1.EvalAll(kLogWeight, fVertexVector); // recalculate the cluster parameters
768 cluPHOS1.SetE(fNonLinCorr->Eval(cluPHOS1.E()));// Users's nonlinearity
769
770 if(cluPHOS1.E()<0.3) continue; // check energy again
771
772 //correct misalignment
773 TVector3 localPos;
774 const Float_t shiftX[6]={0.,-2.3,-2.11,-1.53,0.,0.} ;
775 const Float_t shiftZ[6]={0.,-0.4, 0.52, 0.8,0.,0.} ;
776 fPHOSGeo->Global2Local(localPos,global,mod) ;
777 fPHOSGeo->Local2Global(mod,localPos.X()+shiftX[mod],localPos.Z()+shiftZ[mod],global);
778 position[0]=global.X() ;
779 position[1]=global.Y() ;
780 position[2]=global.Z() ;
781 cluPHOS1.SetPosition(position);
782
783 cluPHOS1.GetMomentum(lorentzMomentum ,origo);
ada6b882 784
785 //TODO: Check, this may be LHC10h specific:
786 if(mod==2) lorentzMomentum*=135.5/134.0 ;
787 if(mod==3) lorentzMomentum*=135.5/137.2 ;
788 if(mod==2) ecore*=135.5/134.0 ;
789 if(mod==3) ecore*=135.5/137.2 ;
790
2cde7444 791 }
ada6b882 792 else if (fEventAOD) { // is ! ESD, AOD.
793 AliESDCaloCluster* aodCluster = (AliESDCaloCluster*) (clu);
794 aodCluster->GetMomentum(lorentzMomentum ,origo);
2cde7444 795 }
ada6b882 796 else {
797 AliError("(Calo)Cluster is neither ESD nor AOD");
798 continue;
2cde7444 799 }
2cde7444 800
ada6b882 801
802 char skey[55];
53c2a02a 803 snprintf(skey,55,"hCluLowM%d",mod) ;
2cde7444 804 FillHistogram(skey,cellX,cellZ,1.);
ada6b882 805 if(lorentzMomentum.E()>1.5){
2cde7444 806 sprintf(skey,"hCluHighM%d",mod) ;
807 FillHistogram(skey,cellX,cellZ,1.);
808 }
ada6b882 809
810 fCaloPhotonsPHOS->Add(new AliCaloPhoton(lorentzMomentum.X(),lorentzMomentum.Py(),lorentzMomentum.Z(),lorentzMomentum.E()) );
811 AliCaloPhoton * ph = (AliCaloPhoton*) fCaloPhotonsPHOS->At( fCaloPhotonsPHOS->GetLast() );
812
2cde7444 813 ph->SetModule(mod) ;
ada6b882 814 lorentzMomentum*= ecore/lorentzMomentum.E() ;
815 ph->SetMomV2(&lorentzMomentum) ;
2cde7444 816 ph->SetNCells(clu->GetNCells());
817 ph->SetDispBit(TestLambda(clu->E(),clu->GetM20(),clu->GetM02())) ;
818 ph->SetDisp2Bit(TestLambda2(clu->E(),clu->GetM20(),clu->GetM02())) ;
819 if(ph->IsDispOK()){
820 sprintf(skey,"hCluDispM%d",mod) ;
821 FillHistogram(skey,cellX,cellZ,1.);
822 }
ada6b882 823
824 // Track Matching
2cde7444 825 Double_t dx=clu->GetTrackDx() ;
826 Double_t dz=clu->GetTrackDz() ;
827 Bool_t cpvBit=kTRUE ; //No track matched by default
828 Bool_t cpvBit2=kTRUE ; //More Strict criterion
ada6b882 829 if( fEventESD ) {
830
831 TArrayI * itracks = static_cast<AliESDCaloCluster*> (clu)->GetTracksMatched() ;
832 if(itracks->GetSize()>0){
833 Int_t iTr = itracks->At(0);
834 if(iTr>=0 && iTr<fEvent->GetNumberOfTracks()){
835 AliVParticle* track = fEvent->GetTrack(iTr);
836 Double_t pt = track->Pt() ;
837 Short_t charge = track->Charge() ;
838 Double_t r=TestCPV(dx, dz, pt, charge) ;
839 cpvBit=(r>2.) ;
840 cpvBit2=(r>4.) ;
841 }
842 }
843 }
844 else if ( fEventAOD ) {
845 int nTracksMatched = clu->GetNTracksMatched();
846 if(nTracksMatched > 0) {
847 AliVTrack* track = dynamic_cast<AliVTrack*> (clu->GetTrackMatched(0));
848 if ( track ) {
849 Double_t pt = track->Pt();
850 Short_t charge = track->Charge();
851 Double_t r = TestCPV(dx, dz, pt, charge) ;
852 cpvBit=(r>2.) ;
853 cpvBit2=(r>4.) ;
854 }
2cde7444 855 }
856 }
857 ph->SetCPVBit(cpvBit) ;
858 ph->SetCPV2Bit(cpvBit2) ;
859 if(cpvBit){
860 sprintf(skey,"hCluVetoM%d",mod) ;
861 FillHistogram(skey,cellX,cellZ,1.);
862 }
863 ph->SetEMCx(float(cellX)) ;
864 ph->SetEMCz(float(cellZ)) ;
ada6b882 865 // ph->SetLambdas(clu->GetM20(),clu->GetM02()) ;
866 ph->SetUnfolded(clu->GetNExMax()<2); // Remember, if it is unfolde
867
2cde7444 868 }
6ec40e8c 869 FillHistogram("hCenPHOS",fCentralityV0M, fCaloPhotonsPHOS->GetEntriesFast()) ;
ada6b882 870}
871//_____________________________________________________________________________
872void AliAnalysisTaskPi0Flow::FillSelectedClusterHistograms()
873{
6ec40e8c 874 for (Int_t i1=0; i1<fCaloPhotonsPHOS->GetEntriesFast(); i1++) {
ada6b882 875 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
2cde7444 876
2cde7444 877 Double_t dphiA=ph1->Phi()-fRPV0A ;
878 while(dphiA<0)dphiA+=TMath::Pi() ;
879 while(dphiA>TMath::Pi())dphiA-=TMath::Pi() ;
880
881 Double_t dphiC=ph1->Phi()-fRPV0C ;
882 while(dphiC<0)dphiC+=TMath::Pi() ;
883 while(dphiC>TMath::Pi())dphiC-=TMath::Pi() ;
884
885 Double_t dphiT=ph1->Phi()-fRP ;
886 while(dphiT<0)dphiT+=TMath::Pi() ;
887 while(dphiT>TMath::Pi())dphiT-=TMath::Pi() ;
ada6b882 888
889 FillHistogram(Form("hPhotPhiV0AAll_cen%d",fCentBin),ph1->Pt(),dphiA) ;
890 FillHistogram(Form("hPhotPhiV0CAll_cen%d",fCentBin),ph1->Pt(),dphiC) ;
2cde7444 891 if(fHaveTPCRP)
ada6b882 892 FillHistogram(Form("hPhotPhiTPCAll_cen%d",fCentBin),ph1->Pt(),dphiT) ;
893
894 FillHistogram(Form("hPhotAll_cen%d",fCentBin),ph1->Pt()) ;
895 FillHistogram(Form("hPhotAllcore_cen%d",fCentBin),ph1->GetMomV2()->Pt()) ;
2cde7444 896 if(ph1->IsntUnfolded()){
ada6b882 897 FillHistogram(Form("hPhotAllwou_cen%d",fCentBin),ph1->Pt()) ;
2cde7444 898 }
899 if(ph1->IsCPVOK()){
ada6b882 900 FillHistogram(Form("hPhotPhiV0ACPV_cen%d",fCentBin),ph1->Pt(),dphiA) ;
901 FillHistogram(Form("hPhotPhiV0CCPV_cen%d",fCentBin),ph1->Pt(),dphiC) ;
2cde7444 902 if(fHaveTPCRP)
ada6b882 903 FillHistogram(Form("hPhotPhiTPCCPV_cen%d",fCentBin),ph1->Pt(),dphiT) ;
2cde7444 904
ada6b882 905 FillHistogram(Form("hPhotPhiV0ACPVcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiA) ;
906 FillHistogram(Form("hPhotPhiV0CCPVcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiC) ;
2cde7444 907 if(fHaveTPCRP)
ada6b882 908 FillHistogram(Form("hPhotPhiTPCCPVcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiT) ;
2cde7444 909
ada6b882 910 FillHistogram(Form("hPhotCPV_cen%d",fCentBin),ph1->Pt()) ;
911 FillHistogram(Form("hPhotCPVcore_cen%d",fCentBin),ph1->GetMomV2()->Pt()) ;
2cde7444 912 }
913 if(ph1->IsCPV2OK()){
ada6b882 914 FillHistogram(Form("hPhotCPV2_cen%d",fCentBin),ph1->Pt()) ;
2cde7444 915 }
916 if(ph1->IsDispOK()){
ada6b882 917 FillHistogram(Form("hPhotPhiV0ADisp_cen%d",fCentBin),ph1->Pt(),dphiA) ;
918 FillHistogram(Form("hPhotPhiV0CDisp_cen%d",fCentBin),ph1->Pt(),dphiC) ;
2cde7444 919 if(fHaveTPCRP)
ada6b882 920 FillHistogram(Form("hPhotPhiTPCDisp_cen%d",fCentBin),ph1->Pt(),dphiT) ;
2cde7444 921
ada6b882 922 FillHistogram(Form("hPhotPhiV0ADispcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiA) ;
923 FillHistogram(Form("hPhotPhiV0CDispcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiC) ;
2cde7444 924 if(fHaveTPCRP)
ada6b882 925 FillHistogram(Form("hPhotPhiTPCDispcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiT) ;
926
927 FillHistogram(Form("hPhotDisp_cen%d",fCentBin),ph1->Pt()) ;
2cde7444 928 if(ph1->IsDisp2OK()){
ada6b882 929 FillHistogram(Form("hPhotDisp2_cen%d",fCentBin),ph1->Pt()) ;
2cde7444 930 }
ada6b882 931 FillHistogram(Form("hPhotDispcore_cen%d",fCentBin),ph1->GetMomV2()->Pt()) ;
2cde7444 932 if(ph1->IsntUnfolded()){
ada6b882 933 FillHistogram(Form("hPhotDispwou_cen%d",fCentBin),ph1->Pt()) ;
2cde7444 934 }
935 if(ph1->IsCPVOK()){
ada6b882 936 FillHistogram(Form("hPhotPhiV0ABoth_cen%d",fCentBin),ph1->Pt(),dphiA) ;
937 FillHistogram(Form("hPhotPhiV0CBoth_cen%d",fCentBin),ph1->Pt(),dphiC) ;
2cde7444 938 if(fHaveTPCRP)
ada6b882 939 FillHistogram(Form("hPhotPhiTPCBoth_cen%d",fCentBin),ph1->Pt(),dphiT) ;
2cde7444 940
ada6b882 941 FillHistogram(Form("hPhotPhiV0ABothcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiA) ;
942 FillHistogram(Form("hPhotPhiV0CBothcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiC) ;
2cde7444 943 if(fHaveTPCRP)
ada6b882 944 FillHistogram(Form("hPhotPhiTPCBothcore_cen%d",fCentBin),ph1->GetMomV2()->Pt(),dphiT) ;
2cde7444 945
ada6b882 946 FillHistogram(Form("hPhotBoth_cen%d",fCentBin),ph1->Pt()) ;
947 FillHistogram(Form("hPhotBothcore_cen%d",fCentBin),ph1->GetMomV2()->Pt()) ;
2cde7444 948 }
ada6b882 949 }
2cde7444 950 }
ada6b882 951}
952//_____________________________________________________________________________
953void AliAnalysisTaskPi0Flow::ConsiderPi0s()
954{
955 char key[55];
6ec40e8c 956 for (Int_t i1=0; i1 < fCaloPhotonsPHOS->GetEntriesFast()-1; i1++) {
ada6b882 957 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
6ec40e8c 958 for (Int_t i2=i1+1; i2<fCaloPhotonsPHOS->GetEntriesFast(); i2++) {
ada6b882 959 AliCaloPhoton * ph2=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i2) ;
2cde7444 960 TLorentzVector p12 = *ph1 + *ph2;
ada6b882 961 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
962 FillHistogram("hPHOSphi",fCentralityV0M,p12.Pt(),p12.Phi());
2cde7444 963 Double_t dphiA=p12.Phi()-fRPV0A ;
964 while(dphiA<0)dphiA+=TMath::Pi() ;
965 while(dphiA>TMath::Pi())dphiA-=TMath::Pi() ;
ada6b882 966
2cde7444 967 Double_t dphiC=p12.Phi()-fRPV0C ;
968 while(dphiC<0)dphiC+=TMath::Pi() ;
969 while(dphiC>TMath::Pi())dphiC-=TMath::Pi() ;
970
971 Double_t dphiT=p12.Phi()-fRP ;
972 while(dphiT<0)dphiT+=TMath::Pi() ;
973 while(dphiT>TMath::Pi())dphiT-=TMath::Pi() ;
ada6b882 974
2cde7444 975 Double_t a=TMath::Abs((ph1->E()-ph2->E())/(ph1->E()+ph2->E())) ;
ada6b882 976
977 FillHistogram(Form("hMassPtV0AAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiA) ;
978 FillHistogram(Form("hMassPtV0CAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiC) ;
2cde7444 979 if(fHaveTPCRP)
ada6b882 980 FillHistogram(Form("hMassPtTPCAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiT) ;
2cde7444 981
ada6b882 982 FillHistogram(Form("hMassPtV0AAllcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiA) ;
983 FillHistogram(Form("hMassPtV0CAllcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiC) ;
2cde7444 984 if(fHaveTPCRP)
ada6b882 985 FillHistogram(Form("hMassPtTPCAllcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiT) ;
986
987
988 FillHistogram(Form("hPi0All_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
989 FillHistogram(Form("hPi0Allcore_cen%d",fCentBin),pv12.M() ,pv12.Pt()) ;
2cde7444 990 if(ph1->IsntUnfolded() && ph2->IsntUnfolded()){
ada6b882 991 FillHistogram(Form("hPi0Allwou_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 992 }
ada6b882 993
994 FillHistogram(Form("hSingleAll_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
995 FillHistogram(Form("hSingleAll_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
996 FillHistogram(Form("hSingleAllcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
997 FillHistogram(Form("hSingleAllcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 998 if(ph1->IsntUnfolded())
ada6b882 999 FillHistogram(Form("hSingleAllwou_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1000 if(ph2->IsntUnfolded())
ada6b882 1001 FillHistogram(Form("hSingleAllwou_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1002 if(ph1->IsCPVOK()){
ada6b882 1003 FillHistogram(Form("hSingleCPV_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1004 FillHistogram(Form("hSingleCPVcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1005 }
1006 if(ph2->IsCPVOK()){
ada6b882 1007 FillHistogram(Form("hSingleCPV_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1008 FillHistogram(Form("hSingleCPV_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1009 }
1010 if(ph1->IsCPV2OK()){
ada6b882 1011 FillHistogram(Form("hSingleCPV2_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1012 }
1013 if(ph2->IsCPV2OK()){
ada6b882 1014 FillHistogram(Form("hSingleCPV2_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1015 }
2cde7444 1016 if(ph1->IsDispOK()){
ada6b882 1017 FillHistogram(Form("hSingleDisp_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1018 if(ph1->IsntUnfolded()){
ada6b882 1019 FillHistogram(Form("hSingleDispwou_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1020 }
ada6b882 1021 FillHistogram(Form("hSingleDispcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1022 }
1023 if(ph2->IsDispOK()){
ada6b882 1024 FillHistogram(Form("hSingleDisp_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1025 if(ph1->IsntUnfolded()){
ada6b882 1026 FillHistogram(Form("hSingleDispwou_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1027 }
ada6b882 1028 FillHistogram(Form("hSingleDispcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1029 }
1030 if(ph1->IsDisp2OK()){
ada6b882 1031 FillHistogram(Form("hSingleDisp2_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1032 }
1033 if(ph2->IsDisp2OK()){
ada6b882 1034 FillHistogram(Form("hSingleDisp2_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1035 }
1036 if(ph1->IsDispOK() && ph1->IsCPVOK()){
ada6b882 1037 FillHistogram(Form("hSingleBoth_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1038 FillHistogram(Form("hSingleBothcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1039 }
1040 if(ph2->IsDispOK() && ph2->IsCPVOK()){
ada6b882 1041 FillHistogram(Form("hSingleBoth_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1042 FillHistogram(Form("hSingleBothcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
1043 }
1044
1045
1046 if(a<kAlphaCut){
1047 FillHistogram(Form("hPi0All_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1048 }
1049
1050 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
ada6b882 1051 snprintf(key,55,"hMassPtCPV_cen%d",fCentBin) ;
1052 FillHistogram(Form("hMassPtV0ACPV_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiA) ;
1053 FillHistogram(Form("hMassPtV0CCPV_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiC) ;
2cde7444 1054 if(fHaveTPCRP)
ada6b882 1055 FillHistogram(Form("hMassPtTPCCPV_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiT) ;
2cde7444 1056
ada6b882 1057 FillHistogram(Form("hMassPtV0ACPVcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiA) ;
1058 FillHistogram(Form("hMassPtV0CCPVcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiC) ;
2cde7444 1059 if(fHaveTPCRP)
ada6b882 1060 FillHistogram(Form("hMassPtTPCCPVcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiT) ;
2cde7444 1061
ada6b882 1062 FillHistogram(Form("hPi0CPV_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1063 FillHistogram(Form("hPi0CPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt()) ;
1064
1065 if(a<kAlphaCut){
1066 FillHistogram(Form("hPi0CPV_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1067 }
ada6b882 1068 }
2cde7444 1069 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
ada6b882 1070 FillHistogram(Form("hPi0CPV2_cen%d",fCentBin),p12.M(),p12.Pt()) ;
1071 if(a<kAlphaCut){
1072 FillHistogram(Form("hPi0CPV2_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1073 }
ada6b882 1074 }
2cde7444 1075 if(ph1->IsDispOK() && ph2->IsDispOK()){
ada6b882 1076 snprintf(key,55,"hMassPtDisp_cen%d",fCentBin) ;
1077 FillHistogram(Form("hMassPtV0ADisp_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiA) ;
1078 FillHistogram(Form("hMassPtV0CDisp_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiC) ;
2cde7444 1079 if(fHaveTPCRP)
ada6b882 1080 FillHistogram(Form("hMassPtTPCDisp_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiT) ;
2cde7444 1081
ada6b882 1082 FillHistogram(Form("hMassPtV0ADispcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiA) ;
1083 FillHistogram(Form("hMassPtV0CDispcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiC) ;
2cde7444 1084 if(fHaveTPCRP)
ada6b882 1085 FillHistogram(Form("hMassPtTPCDispcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiT) ;
1086
1087 FillHistogram(Form("hPi0Disp_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1088 FillHistogram(Form("hPi0Dispcore_cen%d",fCentBin),pv12.M(), pv12.Pt()) ;
2cde7444 1089 if(ph1->IsDisp2OK() && ph2->IsDisp2OK()){
ada6b882 1090 FillHistogram(Form("hPi0Disp2_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1091 }
1092 if(ph1->IsntUnfolded() && ph2->IsntUnfolded()){
ada6b882 1093 FillHistogram(Form("hPi0Dispwou_cen%d",fCentBin),p12.M(), p12.Pt()) ;
2cde7444 1094 }
ada6b882 1095
1096 if(a<kAlphaCut){
1097 FillHistogram(Form("hPi0Disp_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1098 }
ada6b882 1099
2cde7444 1100 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
ada6b882 1101 FillHistogram(Form("hMassPtV0ABoth_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiA) ;
1102 FillHistogram(Form("hMassPtV0CBoth_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiC) ;
2cde7444 1103 if(fHaveTPCRP)
ada6b882 1104 FillHistogram(Form("hMassPtTPCBoth_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiT) ;
2cde7444 1105
ada6b882 1106 FillHistogram(Form("hMassPtV0ABothcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiA) ;
1107 FillHistogram(Form("hMassPtV0CBothcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiC) ;
2cde7444 1108 if(fHaveTPCRP)
ada6b882 1109 FillHistogram(Form("hMassPtTPCBothcore_cen%d",fCentBin),pv12.M() ,pv12.Pt(),dphiT) ;
1110
1111 FillHistogram(Form("hPi0Both_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1112 FillHistogram(Form("hPi0Bothcore_cen%d",fCentBin),pv12.M() ,pv12.Pt()) ;
1113
1114 if(a<kAlphaCut){
1115 snprintf(key,55,"hPi0Both_a07_cen%d",fCentBin) ;
1116 FillHistogram(Form("hPi0Both_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1117 }
ada6b882 1118 if(fCentralityV0M>20.){
2cde7444 1119 if(ph1->Module()==1 && ph2->Module()==1)
1120 FillHistogram("hPi0M11",p12.M(),p12.Pt() );
1121 else if(ph1->Module()==2 && ph2->Module()==2)
1122 FillHistogram("hPi0M22",p12.M(),p12.Pt() );
1123 else if(ph1->Module()==3 && ph2->Module()==3)
1124 FillHistogram("hPi0M33",p12.M(),p12.Pt() );
1125 else if(ph1->Module()==1 && ph2->Module()==2)
1126 FillHistogram("hPi0M12",p12.M(),p12.Pt() );
1127 else if(ph1->Module()==1 && ph2->Module()==3)
1128 FillHistogram("hPi0M13",p12.M(),p12.Pt() );
1129 else if(ph1->Module()==2 && ph2->Module()==3)
1130 FillHistogram("hPi0M23",p12.M(),p12.Pt() );
1131 }
1132
1133 }
1134 }
1135 } // end of loop i2
1136 } // end of loop i1
ada6b882 1137}
1138//_____________________________________________________________________________
1139void AliAnalysisTaskPi0Flow::ConsiderPi0sMix()
1140{
1141 char key[55];
1142
1143 TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
1144
6ec40e8c 1145 for (Int_t i1=0; i1<fCaloPhotonsPHOS->GetEntriesFast(); i1++) {
ada6b882 1146 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
89e3079e 1147 for(Int_t evi=0; evi<arrayList->GetEntries();evi++){
ada6b882 1148 TObjArray * mixPHOS = static_cast<TObjArray*>(arrayList->At(evi));
2cde7444 1149 for(Int_t i2=0; i2<mixPHOS->GetEntriesFast();i2++){
1150 AliCaloPhoton * ph2=(AliCaloPhoton*)mixPHOS->At(i2) ;
1151 TLorentzVector p12 = *ph1 + *ph2;
1152 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
1153
1154 Double_t dphiA=p12.Phi()-fRPV0A ;
1155 while(dphiA<0)dphiA+=TMath::Pi() ;
1156 while(dphiA>TMath::Pi())dphiA-=TMath::Pi() ;
1157
1158 Double_t dphiC=p12.Phi()-fRPV0C ;
1159 while(dphiC<0)dphiC+=TMath::Pi() ;
1160 while(dphiC>TMath::Pi())dphiC-=TMath::Pi() ;
1161
1162 Double_t dphiT=p12.Phi()-fRP ;
1163 while(dphiT<0)dphiT+=TMath::Pi() ;
1164 while(dphiT>TMath::Pi())dphiT-=TMath::Pi() ;
ada6b882 1165
1166
2cde7444 1167 Double_t a=TMath::Abs((ph1->E()-ph2->E())/(ph1->E()+ph2->E())) ;
ada6b882 1168
1169 snprintf(key,55,"hMiMassPtAll_cen%d",fCentBin) ;
1170 FillHistogram(Form("hMiMassPtV0AAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiA) ;
1171 FillHistogram(Form("hMiMassPtV0CAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiC) ;
2cde7444 1172 if(fHaveTPCRP)
ada6b882 1173 FillHistogram(Form("hMiMassPtTPCAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiT) ;
2cde7444 1174
ada6b882 1175 FillHistogram(Form("hMiMassPtV0AAllcore_cen%d",fCentBin),pv12.M(), pv12.Pt(), dphiA) ;
1176 FillHistogram(Form("hMiMassPtV0CAllcore_cen%d",fCentBin),pv12.M(), pv12.Pt(), dphiC) ;
2cde7444 1177 if(fHaveTPCRP)
ada6b882 1178 FillHistogram(Form("hMiMassPtTPCAllcore_cen%d",fCentBin),pv12.M(), pv12.Pt(), dphiT) ;
2cde7444 1179
ada6b882 1180 FillHistogram(Form("hMiPi0All_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1181 FillHistogram(Form("hMiPi0Allcore_cen%d",fCentBin),pv12.M() ,pv12.Pt()) ;
2cde7444 1182 if(ph1->IsntUnfolded() && ph2->IsntUnfolded()){
ada6b882 1183 FillHistogram(Form("hMiPi0Allwou_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1184 }
1185
ada6b882 1186 FillHistogram(Form("hMiSingleAll_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1187 FillHistogram(Form("hMiSingleAll_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1188 FillHistogram(Form("hMiSingleAllcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
1189 FillHistogram(Form("hMiSingleAllcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1190 if(ph1->IsntUnfolded())
ada6b882 1191 FillHistogram(Form("hMiSingleAllwou_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1192 if(ph2->IsntUnfolded())
ada6b882 1193 FillHistogram(Form("hMiSingleAllwou_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1194 if(ph1->IsCPVOK()){
ada6b882 1195 FillHistogram(Form("hMiSingleCPV_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1196 FillHistogram(Form("hMiSingleCPVcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1197 }
1198 if(ph2->IsCPVOK()){
ada6b882 1199 FillHistogram(Form("hMiSingleCPV_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1200 FillHistogram(Form("hMiSingleCPVcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1201 }
1202 if(ph1->IsCPV2OK()){
ada6b882 1203 FillHistogram(Form("hMiSingleCPV2_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1204 }
1205 if(ph2->IsCPV2OK()){
ada6b882 1206 FillHistogram(Form("hMiSingleCPV2_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1207 }
2cde7444 1208 if(ph1->IsDispOK()){
ada6b882 1209 FillHistogram(Form("hMiSingleDisp_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1210 if(ph1->IsntUnfolded()){
ada6b882 1211 FillHistogram(Form("hMiSingleDispwou_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1212 }
ada6b882 1213 FillHistogram(Form("hMiSingleDispcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1214 }
1215 if(ph2->IsDispOK()){
ada6b882 1216 FillHistogram(Form("hMiSingleDisp_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1217 if(ph1->IsntUnfolded()){
ada6b882 1218 FillHistogram(Form("hMiSingleDispwou_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1219 }
ada6b882 1220 FillHistogram(Form("hMiSingleDispcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1221 }
1222 if(ph1->IsDisp2OK()){
ada6b882 1223 FillHistogram(Form("hMiSingleDisp2_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1224 }
2cde7444 1225 if(ph2->IsDisp2OK()){
ada6b882 1226 FillHistogram(Form("hMiSingleDisp2_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1227 }
1228 if(ph1->IsDispOK() && ph1->IsCPVOK()){
ada6b882 1229 snprintf(key,55,"hMiSingleBoth_cen%d",fCentBin) ;
2cde7444 1230 FillHistogram(key,p12.M(),ph1->Pt()) ;
ada6b882 1231 snprintf(key,55,"hMiSingleBothcore_cen%d",fCentBin) ;
2cde7444 1232 FillHistogram(key,pv12.M(),ph1->GetMomV2()->Pt()) ;
1233 }
1234 if(ph2->IsDispOK() && ph2->IsCPVOK()){
ada6b882 1235 snprintf(key,55,"hMiSingleBoth_cen%d",fCentBin) ;
2cde7444 1236 FillHistogram(key,p12.M(),ph2->Pt()) ;
ada6b882 1237 snprintf(key,55,"hMiSingleBothcore_cen%d",fCentBin) ;
2cde7444 1238 FillHistogram(key,pv12.M(),ph2->GetMomV2()->Pt()) ;
ada6b882 1239 }
1240
1241
1242
1243 if(a<kAlphaCut){
1244 snprintf(key,55,"hMiPi0All_a07_cen%d",fCentBin) ;
2cde7444 1245 FillHistogram(key,p12.M() ,p12.Pt()) ;
1246 }
1247 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
ada6b882 1248 FillHistogram(Form("hMiMassPtV0ACPV_cen%d",fCentBin),p12.M(), p12.Pt(),dphiA) ;
1249 FillHistogram(Form("hMiMassPtV0CCPV_cen%d",fCentBin),p12.M(), p12.Pt(),dphiC) ;
2cde7444 1250 if(fHaveTPCRP)
ada6b882 1251 FillHistogram(Form("hMiMassPtTPCCPV_cen%d",fCentBin),p12.M(), p12.Pt(),dphiT) ;
2cde7444 1252
ada6b882 1253 FillHistogram(Form("hMiMassPtV0ACPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiA) ;
1254 FillHistogram(Form("hMiMassPtV0CCPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiC) ;
2cde7444 1255 if(fHaveTPCRP)
ada6b882 1256 FillHistogram(Form("hMiMassPtTPCCPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiT) ;
1257
1258 FillHistogram(Form("hMiPi0CPV_cen%d",fCentBin),p12.M(), p12.Pt()) ;
1259 FillHistogram(Form("hMiPi0CPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt()) ;
2cde7444 1260
ada6b882 1261 if(a<kAlphaCut){
1262 FillHistogram(Form("hMiPi0CPV_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1263 }
1264 }
1265 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
ada6b882 1266 FillHistogram(Form("hMiPi0CPV2_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1267
ada6b882 1268 if(a<kAlphaCut){
1269 FillHistogram(Form("hMiPi0CPV2_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1270 }
1271 }
1272 if(ph1->IsDispOK() && ph2->IsDispOK()){
ada6b882 1273 FillHistogram(Form("hMiMassPtV0ADisp_cen%d",fCentBin),p12.M(),p12.Pt(),dphiA) ;
1274 FillHistogram(Form("hMiMassPtV0CDisp_cen%d",fCentBin),p12.M(),p12.Pt(),dphiC) ;
2cde7444 1275 if(fHaveTPCRP)
ada6b882 1276 FillHistogram(Form("hMiMassPtTPCDisp_cen%d",fCentBin),p12.M(),p12.Pt(),dphiT) ;
1277
1278 FillHistogram(Form("hMiMassPtV0ADispcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiA) ;
1279 FillHistogram(Form("hMiMassPtV0CDispcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiC) ;
2cde7444 1280 if(fHaveTPCRP)
ada6b882 1281 FillHistogram(Form("hMiMassPtTPCDispcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiT) ;
1282
2cde7444 1283
ada6b882 1284 FillHistogram(Form("hMiPi0Disp_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1285 FillHistogram(Form("hMiPi0Dispcore_cen%d",fCentBin),pv12.M(),pv12.Pt()) ;
2cde7444 1286 if(ph1->IsntUnfolded() && ph2->IsntUnfolded()){
ada6b882 1287 FillHistogram(Form("hMiPi0Dispwou_cen%d",fCentBin),p12.M(),p12.Pt()) ;
2cde7444 1288 }
ada6b882 1289
1290 if(a<kAlphaCut){
1291 FillHistogram(Form("hMiPi0Disp_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1292 }
2cde7444 1293 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
ada6b882 1294 FillHistogram(Form("hMiMassPtV0ABoth_cen%d",fCentBin),p12.M(),p12.Pt(),dphiA) ;
1295 FillHistogram(Form("hMiMassPtV0CBoth_cen%d",fCentBin),p12.M(),p12.Pt(),dphiC) ;
2cde7444 1296 if(fHaveTPCRP)
ada6b882 1297 FillHistogram(Form("hMiMassPtTPCBoth_cen%d",fCentBin),p12.M(),p12.Pt(),dphiT) ;
2cde7444 1298
ada6b882 1299 FillHistogram(Form("hMiMassPtV0ABothcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiA) ;
1300 FillHistogram(Form("hMiMassPtV0CBothcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiC) ;
2cde7444 1301 if(fHaveTPCRP)
ada6b882 1302 FillHistogram(Form("hMiMassPtTPCBothcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiT) ;
2cde7444 1303
ada6b882 1304 FillHistogram(Form("hMiPi0Both_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1305 FillHistogram(Form("hMiPi0Bothcore_cen%d",fCentBin),pv12.M(),pv12.Pt()) ;
1306
1307 if(a<kAlphaCut){
1308 FillHistogram(Form("hMiPi0Both_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1309 }
1310 }
1311 }
1312 } // end of loop i2
1313 }
1314 } // end of loop i1
2cde7444 1315}
ada6b882 1316//_____________________________________________________________________________
1317void AliAnalysisTaskPi0Flow::UpdateLists()
2cde7444 1318{
ada6b882 1319 //Now we either add current events to stack or remove
1320 //If no photons in current event - no need to add it to mixed
107849ac 1321
ada6b882 1322 TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
cdd3a269 1323 if( fDebug >= 2 )
107849ac 1324 AliInfo( Form("fCentBin=%d, fCentNMixed[]=%d",fCentBin,fCentNMixed[fCentBin]) );
6ec40e8c 1325 if(fCaloPhotonsPHOS->GetEntriesFast()>0){
ada6b882 1326 arrayList->AddFirst(fCaloPhotonsPHOS) ;
1327 fCaloPhotonsPHOS=0;
107849ac 1328 if(arrayList->GetEntries() > fCentNMixed[fCentBin]){ // Remove redundant events
ada6b882 1329 TObjArray * tmp = static_cast<TObjArray*>(arrayList->Last()) ;
1330 arrayList->RemoveLast() ;
1331 delete tmp ; // TODO: may conflict with delete done by list being owner.
2cde7444 1332 }
2cde7444 1333 }
ada6b882 1334 else
1335 fCaloPhotonsPHOS->Clear(); // TODO: redundant???
2cde7444 1336}
1337//_____________________________________________________________________________
1338void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x)const{
1339 //FillHistogram
9ed0694f 1340 TH1 * hist = dynamic_cast<TH1*>(fOutputContainer->FindObject(key)) ;
1341 if(hist)
1342 hist->Fill(x) ;
1343 else
1344 AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2cde7444 1345}
1346//_____________________________________________________________________________
1347void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x,Double_t y)const{
1348 //FillHistogram
9ed0694f 1349 TH1 * th1 = dynamic_cast<TH1*> (fOutputContainer->FindObject(key));
1350 if(th1)
1351 th1->Fill(x, y) ;
1352 else
1353 AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2cde7444 1354}
1355
1356//_____________________________________________________________________________
1357void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
1358 //Fills 1D histograms with key
9ed0694f 1359 TObject * obj = fOutputContainer->FindObject(key);
1360
1361 TH2 * th2 = dynamic_cast<TH2*> (obj);
1362 if(th2) {
1363 th2->Fill(x, y, z) ;
1364 return;
2cde7444 1365 }
9ed0694f 1366
1367 TH3 * th3 = dynamic_cast<TH3*> (obj);
1368 if(th3) {
1369 th3->Fill(x, y, z) ;
1370 return;
2cde7444 1371 }
9ed0694f 1372
1373 AliError(Form("can not find histogram (of instance TH2) <%s> ",key)) ;
2cde7444 1374}
1375
ada6b882 1376//_____________________________________________________________________________
1377AliVEvent* AliAnalysisTaskPi0Flow::GetEvent()
1378{
1379 fEvent = InputEvent();
1380 if( ! fEvent ) {
1381 AliError("Event could not be retrieved");
1382 PostData(1, fOutputContainer);
1383 }
1384 return fEvent;
1385}
1386
1387
2cde7444 1388//___________________________________________________________________________
ada6b882 1389AliStack* AliAnalysisTaskPi0Flow::GetMCStack()
1390{
1391 fMCStack = 0;
1392 AliVEventHandler* eventHandler = AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler();
1393 if(eventHandler){
1394 AliMCEventHandler* mcEventHandler = dynamic_cast<AliMCEventHandler*> (eventHandler);
1395 if( mcEventHandler)
1396 fMCStack = static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent()->Stack();
1397 }
1398 return fMCStack;
1399}
2cde7444 1400
ada6b882 1401//___________________________________________________________________________
1402Int_t AliAnalysisTaskPi0Flow::GetCentralityBin(Float_t centralityV0M)
1403{
1404 /* fCentBin=1+Int_t(centralityV0M/100. *kNCenBins) ;
1405 if(centralityV0M < 5. || fCentBin < 0)
1406 fCentBin=0 ;
1407 if(fCentBin > kNCenBins-1)
1408 fCentBin = kNCenBins-1 ;
1409 */
1410 int lastBinUpperIndex = fCentEdges.GetSize() -1;
1411 if( centralityV0M > fCentEdges[lastBinUpperIndex] ) {
cdd3a269 1412 if( fDebug >= 1 )
f5e08ac9 1413 AliWarning( Form("centrality (%f) larger then upper edge of last centrality bin (%f)!", centralityV0M, fCentEdges[lastBinUpperIndex]) );
ada6b882 1414 return lastBinUpperIndex-1;
1415 }
1416 if( centralityV0M < fCentEdges[0] ) {
cdd3a269 1417 if( fDebug >= 1 )
f5e08ac9 1418 AliWarning( Form("centrality (%f) smaller then lower edge of first bin (%f)!", centralityV0M, fCentEdges[0]) );
ada6b882 1419 return 0;
1420 }
1421
1422 fCentBin = TMath::BinarySearch<Double_t> ( GetNumberOfCentralityBins(), fCentEdges.GetArray(), centralityV0M );
1423 return fCentBin;
1424}
1425
1426//___________________________________________________________________________
1427Int_t AliAnalysisTaskPi0Flow::GetRPBin()
1428{
1429 Double_t averageRP;
1430 if(fHaveTPCRP)
1431 averageRP = fRPV0A+fRPV0C+fRP /3.;
1432 else
1433 averageRP = fRPV0A+fRPV0C /2.;
1434
1435 fEMRPBin = Int_t(fNEMRPBins*(averageRP)/TMath::Pi());
1436
1437 if(fEMRPBin> (Int_t) fNEMRPBins-1)
1438 fEMRPBin=fNEMRPBins-1 ;
1439 else if(fEMRPBin<0)
1440 fEMRPBin=0;
1441
cdd3a269 1442 if ( fDebug >= 2 )
f5e08ac9 1443 AliInfo(Form("Event Mixing Reaction Plane bin is: %d", fEMRPBin));
1444
ada6b882 1445 return fEMRPBin;
1446}
1447
1448
1449//_____________________________________________________________________________
c9fb2807 1450void AliAnalysisTaskPi0Flow::LogProgress(int step)
ada6b882 1451{
cdd3a269 1452 if(fDebug >= 2) {
ada6b882 1453 AliInfo(Form("step %d completed", step));
1454 }
1455 // the +0.5 is not realy neccisarry, but oh well... -henrik
c9fb2807 1456 //FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
1457 //FillHistogram("hTotSelEvents", step+0.5);
1458}
1459
1460void AliAnalysisTaskPi0Flow::LogSelection(int step, int internalRunNumber)
1461{
1462 // if(fDebug > 1) {
1463 // AliInfo(Form("step %d completed", step));
1464 // }
1465 // the +0.5 is not realy neccisarry, but oh well... -henrik
ada6b882 1466 FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
1467 FillHistogram("hTotSelEvents", step+0.5);
1468}
1469
c9fb2807 1470
ada6b882 1471//___________________________________________________________________________
1472Int_t AliAnalysisTaskPi0Flow::ConvertToInternalRunNumber(Int_t run){
1473 if( kLHC11h == fPeriod ) {
1474 switch(run){
1475 case 170593 : return 179 ;
1476 case 170572 : return 178 ;
1477 case 170556 : return 177 ;
1478 case 170552 : return 176 ;
1479 case 170546 : return 175 ;
1480 case 170390 : return 174 ;
1481 case 170389 : return 173 ;
1482 case 170388 : return 172 ;
1483 case 170387 : return 171 ;
1484 case 170315 : return 170 ;
1485 case 170313 : return 169 ;
1486 case 170312 : return 168 ;
1487 case 170311 : return 167 ;
1488 case 170309 : return 166 ;
1489 case 170308 : return 165 ;
1490 case 170306 : return 164 ;
1491 case 170270 : return 163 ;
1492 case 170269 : return 162 ;
1493 case 170268 : return 161 ;
1494 case 170267 : return 160 ;
1495 case 170264 : return 159 ;
1496 case 170230 : return 158 ;
1497 case 170228 : return 157 ;
1498 case 170208 : return 156 ;
1499 case 170207 : return 155 ;
1500 case 170205 : return 154 ;
1501 case 170204 : return 153 ;
1502 case 170203 : return 152 ;
1503 case 170195 : return 151 ;
1504 case 170193 : return 150 ;
1505 case 170163 : return 149 ;
1506 case 170162 : return 148 ;
1507 case 170159 : return 147 ;
1508 case 170155 : return 146 ;
1509 case 170152 : return 145 ;
1510 case 170091 : return 144 ;
1511 case 170089 : return 143 ;
1512 case 170088 : return 142 ;
1513 case 170085 : return 141 ;
1514 case 170084 : return 140 ;
1515 case 170083 : return 139 ;
1516 case 170081 : return 138 ;
1517 case 170040 : return 137 ;
1518 case 170038 : return 136 ;
1519 case 170036 : return 135 ;
1520 case 170027 : return 134 ;
1521 case 169981 : return 133 ;
1522 case 169975 : return 132 ;
1523 case 169969 : return 131 ;
1524 case 169965 : return 130 ;
1525 case 169961 : return 129 ;
1526 case 169956 : return 128 ;
1527 case 169926 : return 127 ;
1528 case 169924 : return 126 ;
1529 case 169923 : return 125 ;
1530 case 169922 : return 124 ;
1531 case 169919 : return 123 ;
1532 case 169918 : return 122 ;
1533 case 169914 : return 121 ;
1534 case 169859 : return 120 ;
1535 case 169858 : return 119 ;
1536 case 169855 : return 118 ;
1537 case 169846 : return 117 ;
1538 case 169838 : return 116 ;
1539 case 169837 : return 115 ;
1540 case 169835 : return 114 ;
1541 case 169683 : return 113 ;
1542 case 169628 : return 112 ;
1543 case 169591 : return 111 ;
1544 case 169590 : return 110 ;
1545 case 169588 : return 109 ;
1546 case 169587 : return 108 ;
1547 case 169586 : return 107 ;
1548 case 169584 : return 106 ;
1549 case 169557 : return 105 ;
1550 case 169555 : return 104 ;
1551 case 169554 : return 103 ;
1552 case 169553 : return 102 ;
1553 case 169550 : return 101 ;
1554 case 169515 : return 100 ;
1555 case 169512 : return 99 ;
1556 case 169506 : return 98 ;
1557 case 169504 : return 97 ;
1558 case 169498 : return 96 ;
1559 case 169475 : return 95 ;
1560 case 169420 : return 94 ;
1561 case 169419 : return 93 ;
1562 case 169418 : return 92 ;
1563 case 169417 : return 91 ;
1564 case 169415 : return 90 ;
1565 case 169411 : return 89 ;
1566 case 169238 : return 88 ;
1567 case 169236 : return 87 ;
1568 case 169167 : return 86 ;
1569 case 169160 : return 85 ;
1570 case 169156 : return 84 ;
1571 case 169148 : return 83 ;
1572 case 169145 : return 82 ;
1573 case 169144 : return 81 ;
1574 case 169143 : return 80 ;
1575 case 169138 : return 79 ;
1576 case 169099 : return 78 ;
1577 case 169094 : return 77 ;
1578 case 169091 : return 76 ;
1579 case 169045 : return 75 ;
1580 case 169044 : return 74 ;
1581 case 169040 : return 73 ;
1582 case 169035 : return 72 ;
1583 case 168992 : return 71 ;
1584 case 168988 : return 70 ;
1585 case 168984 : return 69 ;
1586 case 168826 : return 68 ;
1587 case 168777 : return 67 ;
1588 case 168514 : return 66 ;
1589 case 168512 : return 65 ;
1590 case 168511 : return 64 ;
1591 case 168467 : return 63 ;
1592 case 168464 : return 62 ;
1593 case 168461 : return 61 ;
1594 case 168460 : return 60 ;
1595 case 168458 : return 59 ;
1596 case 168362 : return 58 ;
1597 case 168361 : return 57 ;
1598 case 168356 : return 56 ;
1599 case 168342 : return 55 ;
1600 case 168341 : return 54 ;
1601 case 168325 : return 53 ;
1602 case 168322 : return 52 ;
1603 case 168318 : return 51 ;
1604 case 168311 : return 50 ;
1605 case 168310 : return 49 ;
1606 case 168213 : return 48 ;
1607 case 168212 : return 47 ;
1608 case 168208 : return 46 ;
1609 case 168207 : return 45 ;
1610 case 168206 : return 44 ;
1611 case 168205 : return 43 ;
1612 case 168204 : return 42 ;
1613 case 168203 : return 41 ;
1614 case 168181 : return 40 ;
1615 case 168177 : return 39 ;
1616 case 168175 : return 38 ;
1617 case 168173 : return 37 ;
1618 case 168172 : return 36 ;
1619 case 168171 : return 35 ;
1620 case 168115 : return 34 ;
1621 case 168108 : return 33 ;
1622 case 168107 : return 32 ;
1623 case 168105 : return 31 ;
1624 case 168104 : return 30 ;
1625 case 168103 : return 29 ;
1626 case 168076 : return 28 ;
1627 case 168069 : return 27 ;
1628 case 168068 : return 26 ;
1629 case 168066 : return 25 ;
1630 case 167988 : return 24 ;
1631 case 167987 : return 23 ;
1632 case 167986 : return 22 ;
1633 case 167985 : return 21 ;
1634 case 167921 : return 20 ;
1635 case 167920 : return 19 ;
1636 case 167915 : return 18 ;
1637 case 167909 : return 17 ;
1638 case 167903 : return 16 ;
1639 case 167902 : return 15 ;
1640 case 167818 : return 14 ;
1641 case 167814 : return 13 ;
1642 case 167813 : return 12 ;
1643 case 167808 : return 11 ;
1644 case 167807 : return 10 ;
1645 case 167806 : return 9 ;
1646 case 167713 : return 8 ;
1647 case 167712 : return 7 ;
1648 case 167711 : return 6 ;
1649 case 167706 : return 5 ;
1650 case 167693 : return 4 ;
1651 case 166532 : return 3 ;
1652 case 166530 : return 2 ;
1653 case 166529 : return 1 ;
1654
1655 default : return 199;
1656 }
1657 }
1658 if( kLHC10h == fPeriod ) {
1659 switch(run){
1660 case 139517 : return 137;
1661 case 139514 : return 136;
1662 case 139513 : return 135;
1663 case 139511 : return 134;
1664 case 139510 : return 133;
1665 case 139507 : return 132;
1666 case 139505 : return 131;
1667 case 139504 : return 130;
1668 case 139503 : return 129;
1669 case 139470 : return 128;
1670 case 139467 : return 127;
1671 case 139466 : return 126;
1672 case 139465 : return 125;
1673 case 139440 : return 124;
1674 case 139439 : return 123;
1675 case 139438 : return 122;
1676 case 139437 : return 121;
1677 case 139360 : return 120;
1678 case 139329 : return 119;
1679 case 139328 : return 118;
1680 case 139314 : return 117;
1681 case 139311 : return 116;
1682 case 139310 : return 115;
1683 case 139309 : return 114;
1684 case 139308 : return 113;
1685 case 139173 : return 112;
1686 case 139172 : return 111;
1687 case 139110 : return 110;
1688 case 139107 : return 109;
1689 case 139105 : return 108;
1690 case 139104 : return 107;
1691 case 139042 : return 106;
1692 case 139038 : return 105;
1693 case 139037 : return 104;
1694 case 139036 : return 103;
1695 case 139029 : return 102;
1696 case 139028 : return 101;
1697 case 138983 : return 100;
1698 case 138982 : return 99;
1699 case 138980 : return 98;
1700 case 138979 : return 97;
1701 case 138978 : return 96;
1702 case 138977 : return 95;
1703 case 138976 : return 94;
1704 case 138973 : return 93;
1705 case 138972 : return 92;
1706 case 138965 : return 91;
1707 case 138924 : return 90;
1708 case 138872 : return 89;
1709 case 138871 : return 88;
1710 case 138870 : return 87;
1711 case 138837 : return 86;
1712 case 138830 : return 85;
1713 case 138828 : return 84;
1714 case 138826 : return 83;
1715 case 138796 : return 82;
1716 case 138795 : return 81;
1717 case 138742 : return 80;
1718 case 138732 : return 79;
1719 case 138730 : return 78;
1720 case 138666 : return 77;
1721 case 138662 : return 76;
1722 case 138653 : return 75;
1723 case 138652 : return 74;
1724 case 138638 : return 73;
1725 case 138624 : return 72;
1726 case 138621 : return 71;
1727 case 138583 : return 70;
1728 case 138582 : return 69;
1729 case 138579 : return 68;
1730 case 138578 : return 67;
1731 case 138534 : return 66;
1732 case 138469 : return 65;
1733 case 138442 : return 64;
1734 case 138439 : return 63;
1735 case 138438 : return 62;
1736 case 138396 : return 61;
1737 case 138364 : return 60;
1738 case 138359 : return 59;
1739 case 138275 : return 58;
1740 case 138225 : return 57;
1741 case 138201 : return 56;
1742 case 138200 : return 55;
1743 case 138197 : return 54;
1744 case 138192 : return 53;
1745 case 138190 : return 52;
1746 case 138154 : return 51;
1747 case 138153 : return 50;
1748 case 138151 : return 49;
1749 case 138150 : return 48;
1750 case 138126 : return 47;
1751 case 138125 : return 46;
1752 case 137848 : return 45;
1753 case 137847 : return 44;
1754 case 137844 : return 43;
1755 case 137843 : return 42;
1756 case 137752 : return 41;
1757 case 137751 : return 40;
1758 case 137748 : return 39;
1759 case 137724 : return 38;
1760 case 137722 : return 37;
1761 case 137718 : return 36;
1762 case 137704 : return 35;
1763 case 137693 : return 34;
1764 case 137692 : return 33;
1765 case 137691 : return 32;
1766 case 137689 : return 31;
1767 case 137686 : return 30;
1768 case 137685 : return 29;
1769 case 137639 : return 28;
1770 case 137638 : return 27;
1771 case 137608 : return 26;
1772 case 137595 : return 25;
1773 case 137549 : return 24;
1774 case 137546 : return 23;
1775 case 137544 : return 22;
1776 case 137541 : return 21;
1777 case 137539 : return 20;
1778 case 137531 : return 19;
1779 case 137530 : return 18;
1780 case 137443 : return 17;
1781 case 137441 : return 16;
1782 case 137440 : return 15;
1783 case 137439 : return 14;
1784 case 137434 : return 13;
1785 case 137432 : return 12;
1786 case 137431 : return 11;
1787 case 137430 : return 10;
1788 case 137366 : return 9;
1789 case 137243 : return 8;
1790 case 137236 : return 7;
1791 case 137235 : return 6;
1792 case 137232 : return 5;
1793 case 137231 : return 4;
1794 case 137165 : return 3;
1795 case 137162 : return 2;
1796 case 137161 : return 1;
1797 default : return 199;
1798 }
1799 }
cdd3a269 1800 if(kUndefinedPeriod && fDebug >= 1 ) {
ada6b882 1801 AliWarning("Period not defined");
1802 }
1803 return 1;
2cde7444 1804}
1805//_____________________________________________________________________________
1806Bool_t AliAnalysisTaskPi0Flow::TestLambda(Double_t pt,Double_t l1,Double_t l2){
ada6b882 1807
2cde7444 1808 Double_t l2Mean = 1.53126+9.50835e+06/(1.+1.08728e+07*pt+1.73420e+06*pt*pt) ;
1809 Double_t l1Mean = 1.12365+0.123770*TMath::Exp(-pt*0.246551)+5.30000e-03*pt ;
1810 Double_t l2Sigma = 6.48260e-02+7.60261e+10/(1.+1.53012e+11*pt+5.01265e+05*pt*pt)+9.00000e-03*pt;
1811 Double_t l1Sigma = 4.44719e-04+6.99839e-01/(1.+1.22497e+00*pt+6.78604e-07*pt*pt)+9.00000e-03*pt;
1812 Double_t c=-0.35-0.550*TMath::Exp(-0.390730*pt) ;
ada6b882 1813 Double_t R2=0.5*(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma +
1814 0.5*(l2-l2Mean)*(l2-l2Mean)/l2Sigma/l2Sigma +
1815 0.5*c*(l1-l1Mean)*(l2-l2Mean)/l1Sigma/l2Sigma ;
2cde7444 1816 return (R2<2.5*2.5) ;
ada6b882 1817
2cde7444 1818}
1819//_____________________________________________________________________________
1820Bool_t AliAnalysisTaskPi0Flow::TestLambda2(Double_t pt,Double_t l1,Double_t l2){
ada6b882 1821
2cde7444 1822 Double_t l2Mean = 1.53126+9.50835e+06/(1.+1.08728e+07*pt+1.73420e+06*pt*pt) ;
1823 Double_t l1Mean = 1.12365+0.123770*TMath::Exp(-pt*0.246551)+5.30000e-03*pt ;
1824 Double_t l2Sigma = 6.48260e-02+7.60261e+10/(1.+1.53012e+11*pt+5.01265e+05*pt*pt)+9.00000e-03*pt;
1825 Double_t l1Sigma = 4.44719e-04+6.99839e-01/(1.+1.22497e+00*pt+6.78604e-07*pt*pt)+9.00000e-03*pt;
1826 Double_t c=-0.35-0.550*TMath::Exp(-0.390730*pt) ;
ada6b882 1827 Double_t R2=0.5*(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma +
1828 0.5*(l2-l2Mean)*(l2-l2Mean)/l2Sigma/l2Sigma +
1829 0.5*c*(l1-l1Mean)*(l2-l2Mean)/l1Sigma/l2Sigma ;
2cde7444 1830 return (R2<1.5*1.5) ;
ada6b882 1831
2cde7444 1832}
ada6b882 1833//____________________________________________________________________________
1834TList* AliAnalysisTaskPi0Flow::GetCaloPhotonsPHOSList(UInt_t vtxBin, UInt_t centBin, UInt_t rpBin)
1835{
1836 int offset = vtxBin * GetNumberOfCentralityBins() * fNEMRPBins
1837 + centBin * fNEMRPBins
1838 + rpBin;
1839 if( fCaloPhotonsPHOSLists->At(offset) ) { // list exists
1840 TList* list = dynamic_cast<TList*> (fCaloPhotonsPHOSLists->At(offset));
1841 if( ! list )
1842 AliError("object in fCaloPhotonsPHOSLists at %i did not cast");
1843 return list;
1844 }
1845 else {// no list for this bin has been created, yet
1846 TList* list = new TList();
1847 list->SetOwner();
1848 fCaloPhotonsPHOSLists->AddAt(list, offset);
1849 return list;
1850 }
1851}
1852
2cde7444 1853//____________________________________________________________________________
1854Double_t AliAnalysisTaskPi0Flow::TestCPV(Double_t dx, Double_t dz, Double_t pt, Int_t charge){
1855 //Parameterization of LHC10h period
1856 //_true if neutral_
ada6b882 1857
2cde7444 1858 Double_t meanX=0;
1859 Double_t meanZ=0.;
1860 Double_t sx=TMath::Min(5.4,2.59719e+02*TMath::Exp(-pt/1.02053e-01)+
ada6b882 1861 6.58365e-01*5.91917e-01*5.91917e-01/((pt-9.61306e-01)*(pt-9.61306e-01)+5.91917e-01*5.91917e-01)+1.59219);
2cde7444 1862 Double_t sz=TMath::Min(2.75,4.90341e+02*1.91456e-02*1.91456e-02/(pt*pt+1.91456e-02*1.91456e-02)+1.60) ;
ada6b882 1863 AliVEvent *event = InputEvent();
53c2a02a 1864 Double_t mf = 0.; //
1865 if(event) mf = event->GetMagneticField(); //Positive for ++ and negative for --
2cde7444 1866
1867 if(mf<0.){ //field --
1868 meanZ = -0.468318 ;
1869 if(charge>0)
1870 meanX=TMath::Min(7.3, 3.89994*1.20679*1.20679/(pt*pt+1.20679*1.20679)+0.249029+2.49088e+07*TMath::Exp(-pt*3.33650e+01)) ;
1871 else
1872 meanX=-TMath::Min(7.7,3.86040*0.912499*0.912499/(pt*pt+0.912499*0.912499)+1.23114+4.48277e+05*TMath::Exp(-pt*2.57070e+01)) ;
1873 }
1874 else{ //Field ++
1875 meanZ= -0.468318;
1876 if(charge>0)
1877 meanX=-TMath::Min(8.0,3.86040*1.31357*1.31357/(pt*pt+1.31357*1.31357)+0.880579+7.56199e+06*TMath::Exp(-pt*3.08451e+01)) ;
1878 else
ada6b882 1879 meanX= TMath::Min(6.85, 3.89994*1.16240*1.16240/(pt*pt+1.16240*1.16240)-0.120787+2.20275e+05*TMath::Exp(-pt*2.40913e+01)) ;
2cde7444 1880 }
1881
1882 Double_t rz=(dz-meanZ)/sz ;
1883 Double_t rx=(dx-meanX)/sx ;
1884 return TMath::Sqrt(rx*rx+rz*rz) ;
1885}
1886//____________________________________________________________________________
b509e3a5 1887void AliAnalysisTaskPi0Flow::SetFlatteningData(){
1888 //Read objects with flattening parameters
1889 AliOADBContainer flatContainer("phosFlat");
1890 flatContainer.InitFromFile(fEPcalibFileName.Data(),"phosFlat");
1891 TObjArray *maps = (TObjArray*)flatContainer.GetObject(fRunNumber,"phosFlat");
1892 if(!maps){
1893 AliError(Form("Can not read Flattening for run %d. \n From file >%s<\n",fRunNumber,fEPcalibFileName.Data())) ;
2cde7444 1894 }
b509e3a5 1895 else{
1896 AliInfo(Form("Setting PHOS flattening with name %s \n",maps->GetName())) ;
1897 AliPHOSEPFlattener * h = (AliPHOSEPFlattener*)maps->At(0) ;
1898 if(fTPCFlat) delete fTPCFlat ;
1899 fTPCFlat = new AliPHOSEPFlattener() ;
1900 fTPCFlat = h ;
1901 h = (AliPHOSEPFlattener*)maps->At(1) ;
1902 if(fV0AFlat) delete fV0AFlat ;
1903 fV0AFlat = new AliPHOSEPFlattener() ;
1904 fV0AFlat = h ;
1905 h = (AliPHOSEPFlattener*)maps->At(2) ;
1906 if(fV0CFlat) delete fV0CFlat ;
1907 fV0CFlat = new AliPHOSEPFlattener() ;
1908 fV0CFlat = h ;
1909 }
1910
2cde7444 1911}
b509e3a5 1912 //____________________________________________________________________________
1913Double_t AliAnalysisTaskPi0Flow::ApplyFlattening(Double_t phi, Double_t c){
1914
1915 if(fTPCFlat)
1916 return fTPCFlat->MakeFlat(phi,c);
1917 return phi ;
2cde7444 1918
2cde7444 1919}
b509e3a5 1920//____________________________________________________________________________
1921Double_t AliAnalysisTaskPi0Flow::ApplyFlatteningV0A(Double_t phi, Double_t c){
1922
1923 if(fV0AFlat)
1924 return fV0AFlat->MakeFlat(phi,c);
1925 return phi ;
2cde7444 1926
1927}
1928//____________________________________________________________________________
b509e3a5 1929Double_t AliAnalysisTaskPi0Flow::ApplyFlatteningV0C(Double_t phi, Double_t c){
1930
1931 if(fV0CFlat)
1932 return fV0CFlat->MakeFlat(phi,c);
ada6b882 1933 return phi ;
1934
1935}
2cde7444 1936//____________________________________________________________________________
91b112bd 1937Double_t AliAnalysisTaskPi0Flow::CoreEnergy(AliVCluster * clu, AliVCaloCells * cells)
1938{
2cde7444 1939 //calculate energy of the cluster in the circle with radius distanceCut around the maximum
ada6b882 1940
2cde7444 1941 //Can not use already calculated coordinates?
1942 //They have incidence correction...
1943 const Double_t distanceCut =3.5 ;
1944 const Double_t logWeight=4.5 ;
ada6b882 1945
1946 Double32_t * elist = clu->GetCellsAmplitudeFraction() ;
2cde7444 1947// Calculates the center of gravity in the local PHOS-module coordinates
1948 Float_t wtot = 0;
1949 Double_t xc[100]={0} ;
1950 Double_t zc[100]={0} ;
1951 Double_t x = 0 ;
1952 Double_t z = 0 ;
1953 Int_t mulDigit=TMath::Min(100,clu->GetNCells()) ;
1954 for(Int_t iDigit=0; iDigit<mulDigit; iDigit++) {
1955 Int_t relid[4] ;
1956 Float_t xi ;
1957 Float_t zi ;
1958 fPHOSGeo->AbsToRelNumbering(clu->GetCellAbsId(iDigit), relid) ;
1959 fPHOSGeo->RelPosInModule(relid, xi, zi);
1960 xc[iDigit]=xi ;
1961 zc[iDigit]=zi ;
91b112bd 1962 elist[iDigit] *= cells->GetCellAmplitude(clu->GetCellsAbsId()[iDigit]);
cdd3a269 1963 if( fDebug >= 3 )
1964 printf("%f ",elist[iDigit]);
2cde7444 1965 if (clu->E()>0 && elist[iDigit]>0) {
1966 Float_t w = TMath::Max( 0., logWeight + TMath::Log( elist[iDigit] / clu->E() ) ) ;
1967 x += xc[iDigit] * w ;
1968 z += zc[iDigit] * w ;
1969 wtot += w ;
1970 }
1971 }
1972 if (wtot>0) {
1973 x /= wtot ;
1974 z /= wtot ;
1975 }
1976 Double_t coreE=0. ;
1977 for(Int_t iDigit=0; iDigit < mulDigit; iDigit++) {
1978 Double_t distance = TMath::Sqrt((xc[iDigit]-x)*(xc[iDigit]-x)+(zc[iDigit]-z)*(zc[iDigit]-z)) ;
1979 if(distance < distanceCut)
1980 coreE += elist[iDigit] ;
1981 }
1982 //Apply non-linearity correction
1983 return fNonLinCorr->Eval(coreE) ;
1984}
1985//____________________________________________________________________________
1986Bool_t AliAnalysisTaskPi0Flow::AreNeibors(Int_t id1,Int_t id2){
ada6b882 1987 // return true if absId are "Neighbors" (adjacent, including diagornaly,)
1988 // false if not.
2cde7444 1989
ada6b882 1990 Int_t relid1[4] ;
1991 fPHOSGeo->AbsToRelNumbering(id1, relid1) ;
2cde7444 1992
ada6b882 1993 Int_t relid2[4] ;
1994 fPHOSGeo->AbsToRelNumbering(id2, relid2) ;
2cde7444 1995
ada6b882 1996 // if inside the same PHOS module
1997 if ( (relid1[0] == relid2[0]) && (relid1[1]==relid2[1]) ) {
1998 const Int_t rowdiff = TMath::Abs( relid1[2] - relid2[2] ) ;
1999 const Int_t coldiff = TMath::Abs( relid1[3] - relid2[3] ) ;
2000
2001 // and if diff in both direction is 1 or less
2002 if (( coldiff <= 1 ) && ( rowdiff <= 1 ))
2003 return true; // are neighbors
2004 }
2cde7444 2005
ada6b882 2006 // else false
2007 return false;
2cde7444 2008}
2009//____________________________________________________________________________
ada6b882 2010void AliAnalysisTaskPi0Flow::Reclusterize(AliVCluster * clu){
2cde7444 2011 //Re-clusterize to make continues cluster
ada6b882 2012
2cde7444 2013 const Int_t oldMulDigit=clu->GetNCells() ;
ada6b882 2014 Double32_t * elist = clu->GetCellsAmplitudeFraction() ;
2cde7444 2015 UShort_t * dlist = clu->GetCellsAbsId();
2016
2017 Int_t index[oldMulDigit] ;
2018 Bool_t used[oldMulDigit] ;
2019 for(Int_t i=0; i<oldMulDigit; i++) used[i]=0 ;
2020 Int_t inClu=0 ;
2021 Double_t eMax=0. ;
2022 //find maximum
2023 for(Int_t iDigit=0; iDigit<oldMulDigit; iDigit++) {
2024 if(eMax<elist[iDigit]){
2025 eMax=elist[iDigit];
2026 index[0]=iDigit ;
2027 inClu=1 ;
2028 }
2029 }
2030 if(inClu==0){ //empty cluster
2031 return ;
2032 }
2033 used[index[0]]=kTRUE ; //mark as used
2034 for(Int_t i=0; i<inClu; i++){
2035 for(Int_t iDigit=0 ;iDigit<oldMulDigit; iDigit++){
2036 if(used[iDigit]) //already used
2037 continue ;
2038 if(AreNeibors(dlist[index[i]],dlist[iDigit])){
2039 index[inClu]= iDigit ;
2040 inClu++ ;
2041 used[iDigit]=kTRUE ;
2042 }
2043 }
2044 }
ada6b882 2045
2cde7444 2046 if(inClu==oldMulDigit) //no need to modify
ada6b882 2047 return ;
2cde7444 2048
2049 clu->SetNCells(inClu);
2050 //copy
2051 UShort_t tmpD[oldMulDigit] ;
2052 Double_t tmpE[oldMulDigit] ;
2053 for(Int_t i=0; i<oldMulDigit; i++){
ada6b882 2054 tmpD[i]=dlist[i] ;
2cde7444 2055 tmpE[i]=elist[i] ;
2056 }
2057 //change order of digits in list so that
2058 //first inClu cells were true ones
2059 for(Int_t i=0; i<inClu; i++){
2060 dlist[i]=tmpD[index[i]] ;
2061 elist[i]=tmpE[index[i]] ;
2062 }
ada6b882 2063
2064
2065}
2066
2067//_____________________________________________________________________________
2068void AliAnalysisTaskPi0Flow::SetMisalignment(){
2069 // sets the misalignment vertex if ESD
2070 if( fEventESD ) {
2071 for(Int_t mod=0; mod<5; mod++) {
2072 const TGeoHMatrix* modMatrix = fEvent->GetPHOSMatrix(mod);
2073 if( ! modMatrix) {
f5e08ac9 2074 if( fDebug )
ada6b882 2075 AliInfo(Form("no PHOS Geometric Misalignment Matrix for module %d", mod));
2076 continue;
2077 }
2078 else {
2079 fPHOSGeo->SetMisalMatrix(modMatrix, mod);
f5e08ac9 2080 if( fDebug )
cdd3a269 2081 AliInfo(Form("PHOS Geometric Misalignment Matrix set for module %d", mod));
ada6b882 2082 }
2083 }
2084 }
2cde7444 2085}
2086
2087//_____________________________________________________________________________
ada6b882 2088void AliAnalysisTaskPi0Flow::SetV0Calibration(){
2089 // assigns: fMultV0, fV0Cpol, fV0Apol, fMeanQ, and fWidthQ
2090
2091 int runNumber = this->fRunNumber;
2092
2cde7444 2093 TString oadbfilename = "$ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root";
2094 TFile *foadb = TFile::Open(oadbfilename.Data());
2095
2096 if(!foadb){
ada6b882 2097 AliError(Form("OADB file %s cannot be opened\n", oadbfilename.Data()));
2098 AliError("V0 Calibration not set !\n");
2cde7444 2099 return;
2100 }
2101
2102 AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr");
2103 if(!cont){
ada6b882 2104 AliError("OADB object hMultV0BefCorr is not available in the file");
2105 AliError("V0 Calibration not set!\n");
2106 return;
2cde7444 2107 }
2108
ada6b882 2109 if(!(cont->GetObject(runNumber))){
2110 AliError(Form("OADB object hMultV0BefCorr is not available for run %i, trying 137366)",runNumber));
2111 runNumber = 137366;
2112 }
2113 if(!(cont->GetObject(runNumber))){
2114 AliError(Form("OADB object hMultV0BefCorr is not available for run %i ",runNumber));
2115 AliError("V0 Calibration not set!\n");
2116 return;
2cde7444 2117 }
2cde7444 2118
f5e08ac9 2119 if( fDebug ) AliInfo("Setting V0 calibration") ;
ada6b882 2120 fMultV0 = ((TH2F *) cont->GetObject(runNumber))->ProfileX();
2121
2122 TF1 *fpol0 = new TF1("fpol0","pol0");
e2139d72 2123 fMultV0->Fit(fpol0,"Q0","",0,31);
2cde7444 2124 fV0Cpol = fpol0->GetParameter(0);
e2139d72 2125 fMultV0->Fit(fpol0,"Q0","",32,64);
2cde7444 2126 fV0Apol = fpol0->GetParameter(0);
2127
2128 for(Int_t iside=0;iside<2;iside++){
2129 for(Int_t icoord=0;icoord<2;icoord++){
ada6b882 2130 for(Int_t i=0;i < kNCenBins;i++){
2cde7444 2131 char namecont[100];
2132 if(iside==0 && icoord==0)
53c2a02a 2133 snprintf(namecont,100,"hQxc2_%i",i);
2cde7444 2134 else if(iside==1 && icoord==0)
53c2a02a 2135 snprintf(namecont,100,"hQxa2_%i",i);
2cde7444 2136 else if(iside==0 && icoord==1)
53c2a02a 2137 snprintf(namecont,100,"hQyc2_%i",i);
2cde7444 2138 else if(iside==1 && icoord==1)
53c2a02a 2139 snprintf(namecont,100,"hQya2_%i",i);
2cde7444 2140
2141 cont = (AliOADBContainer*) foadb->Get(namecont);
2142 if(!cont){
ada6b882 2143 AliError(Form("OADB object %s is not available in the file %s", namecont, oadbfilename.Data()));
2144 AliError("V0 Calibration not fully set!\n");
2145 return;
2cde7444 2146 }
2cde7444 2147
ada6b882 2148 if(!(cont->GetObject(runNumber))){
2149 AliError(Form("OADB object %s is not available for run %i, trying run 137366",namecont,runNumber));
2150 runNumber = 137366;
2cde7444 2151 }
ada6b882 2152 if(!(cont->GetObject(runNumber))){
2153 AliError(Form("OADB object %s is not available for run %i",namecont,runNumber));
2154 AliError("V0 Calibration not fully set!\n");
2155 return;
2cde7444 2156 }
ada6b882 2157 fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(runNumber))->GetMean();
2158 fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(runNumber))->GetRMS();
2159
2160 //for v3
2161// if(iside==0 && icoord==0)
2162// snprintf(namecont,100,"hQxc3_%i",i);
2163// else if(iside==1 && icoord==0)
2164// snprintf(namecont,100,"hQxa3_%i",i);
2165// else if(iside==0 && icoord==1)
2166// snprintf(namecont,100,"hQyc3_%i",i);
2167// else if(iside==1 && icoord==1)
2168// snprintf(namecont,100,"hQya3_%i",i);
2169//
2170// cont = (AliOADBContainer*) foadb->Get(namecont);
2171// if(!cont){
2172// AliError(Form("OADB object %s is not available in the file",namecont));
2173// AliError("V0 Calibration not fully set!\n");
2174// return;
2175// }
2176//
2177// if(!(cont->GetObject(runNumber))){
2178// AliError(Form("OADB object %s is not available for run %i, trying run 137366",namecont,runNumber));
2179// runNumber = 137366;
2180// }
2181// if(!(cont->GetObject(runNumber))){
2182// AliError(Form("OADB object %s is not available for run %i",namecont,runNumber));
2183// AliError("V0 Calibration not fully set!\n");
2184// return;
2185// }
2cde7444 2186// fMeanQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
2187// fWidthQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
2188
2189 }
2190 }
2191 }
ada6b882 2192
2193 delete fpol0; fpol0=0;
2cde7444 2194}
ada6b882 2195
2196//_____________________________________________________________________________
2197void AliAnalysisTaskPi0Flow::SetESDTrackCuts()
2198{
2199 if( fEventESD ) {
2200 // Create ESD track cut
2201 fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts() ;
2202 //fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
2203 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
2204 }
2205}
2206
2207//_____________________________________________________________________________
2208void AliAnalysisTaskPi0Flow::SetGeometry()
2209{
2210 // Initialize the PHOS geometry
f5e08ac9 2211 if( kLHC10h == fPeriod && fEventESD ) {
ada6b882 2212 TGeoManager::Import("geometry.root"); //TODO: should perhaps not be done
2213 fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP") ;
2214 if( ! fPHOSGeo )
2215 AliError("geometry (fPHOSGeo) not initialised");
2216 }
2217
2218 //Init geometry
2219 if(!fPHOSGeo){
2220 AliOADBContainer geomContainer("phosGeo");
2221 geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
2222 TObjArray *matrixes = (TObjArray*)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
2223 fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP") ;
2224 for(Int_t mod=0; mod<5; mod++) {
2225 if(!matrixes->At(mod)) {
2226 continue;
cdd3a269 2227 if( fDebug )
ada6b882 2228 AliInfo(Form("No PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2229 }
2230 else {
2231 fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod) ;
2232 if( fDebug >1 )
2233 AliInfo(Form("Adding PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2234 }
2235 }
2236 }
2237}
2238
2239//_____________________________________________________________________________
2240void AliAnalysisTaskPi0Flow::SetPHOSCalibData()
2241{
2242 if( fPHOSCalibData )
2243 delete fPHOSCalibData;
2244 fPHOSCalibData = 0;
2cde7444 2245
ada6b882 2246 // Calibration only needed for ESD
2247 if( fEventESD /*&& */ ) {
f5e08ac9 2248 if( kLHC10h == fPeriod && fEventESD ) {
ada6b882 2249 //We have to apply re-calibration for pass1 LCH10h
2250 // Initialize decalibration factors in the form of the OCDB object
2251 AliCDBManager * man = AliCDBManager::Instance();
2252 man->SetRun(140000) ; //TODO; revise, this should probably not b done.
2253 man->SetDefaultStorage("local://OCDB");
2254 }
2255 fPHOSCalibData = new AliPHOSCalibData();
2256 }
2257}
2258
2259//_____________________________________________________________________________
2260void AliAnalysisTaskPi0Flow::SetVertex()
2261{
2262 const AliVVertex *primaryVertex = fEvent->GetPrimaryVertex();
2263 if( primaryVertex ) {
2264 fVertex[0] = primaryVertex->GetX();
2265 fVertex[1] = primaryVertex->GetY();
2266 fVertex[2] = primaryVertex->GetZ();
2267 }
2268 else {
2269 AliError("Event has 0x0 Primary Vertex, defaulting to origo");
2270 fVertex[0] = 0;
2271 fVertex[1] = 0;
2272 fVertex[2] = 0;
2273 }
2274 fVertexVector = TVector3(fVertex);
2275 FillHistogram("hZvertex", fVertexVector.z(), fInternalRunNumber-0.5);
f5e08ac9 2276
cdd3a269 2277 if( fDebug >= 2 )
f5e08ac9 2278 AliInfo(Form("Vertex is set to (%.1f,%.1f,%.1f)", fVertex[0], fVertex[1], fVertex[2]));
ada6b882 2279
2280 fVtxBin=0 ;// No support for vtx binning implemented.
2281}
2282
2283//_____________________________________________________________________________
2284Bool_t AliAnalysisTaskPi0Flow::RejectEventVertex()
2285{
2286 if( ! fEvent->GetPrimaryVertex() )
2287 return true; // reject
c9fb2807 2288 LogSelection(1, fInternalRunNumber);
2289
ada6b882 2290 if ( TMath::Abs(fVertexVector.z()) > 10. )
2291 return true; // reject
c9fb2807 2292 LogSelection(2, fInternalRunNumber);
ada6b882 2293
2294 return false; // accept event.
2295}
2296
2297//_____________________________________________________________________________
2298void AliAnalysisTaskPi0Flow::SetCentrality()
2299{
2300 AliCentrality *centrality = fEvent->GetCentrality();
2301 if( centrality )
2302 fCentralityV0M=centrality->GetCentralityPercentile("V0M");
2303 else {
2304 AliError("Event has 0x0 centrality");
2305 fCentralityV0M = -1.;
2306 }
2307 FillHistogram("hCentrality",fCentralityV0M,fInternalRunNumber-0.5) ;
2308
2309 fCentBin = GetCentralityBin(fCentralityV0M);
f5e08ac9 2310
cdd3a269 2311 if ( fDebug >= 2 )
f5e08ac9 2312 AliInfo(Form("Centrality (bin) is: %f (%d)", fCentralityV0M, fCentBin));
ada6b882 2313}
2314
2315//_____________________________________________________________________________
2316Bool_t AliAnalysisTaskPi0Flow::RejectCentrality()
2317{
2318 if( ! fEvent->GetCentrality() )
2319 return true; // reject
c9fb2807 2320 LogSelection(3, fInternalRunNumber);
2321
ada6b882 2322// if( fCentralityV0M <= 0. || fCentralityV0M>80. )
2323// return true; // reject
2324
2325 int lastBinUpperIndex = fCentEdges.GetSize() -1;
f5e08ac9 2326 if( fCentralityV0M > fCentEdges[lastBinUpperIndex] ) {
2327 if( fDebug )
2328 AliInfo("Rejecting due to centrality outside of binning.");
ada6b882 2329 return true; // reject
f5e08ac9 2330 }
c9fb2807 2331 LogSelection(4, fInternalRunNumber);
2332
f5e08ac9 2333 if( fCentralityV0M < fCentEdges[0] ) {
2334 if( fDebug )
2335 AliInfo("Rejecting due to centrality outside of binning.");
ada6b882 2336 return true; // reject
f5e08ac9 2337 }
c9fb2807 2338 LogSelection(5, fInternalRunNumber);
ada6b882 2339
2340 return false;
2341}
2342
2343
2344//_____________________________________________________________________________
2345void AliAnalysisTaskPi0Flow::EvalReactionPlane()
2346{
2347 // assigns: fHaveTPCRP and fRP
2348 // also does a few histogram fills
2349
2350 AliEventplane *eventPlane = fEvent->GetEventplane();
f5e08ac9 2351 if( ! eventPlane ) { AliError("Event has no event plane"); return; }
2352
ada6b882 2353 Double_t reactionPlaneQ = eventPlane->GetEventplane("Q");
2354 FillHistogram("phiRP",reactionPlaneQ,fCentralityV0M) ;
2355
2356 if(reactionPlaneQ==999 || reactionPlaneQ < 0.){ //reaction plain was not defined
f5e08ac9 2357 if( fDebug ) AliInfo(Form("No Q Reaction Plane, value is %f", reactionPlaneQ));
2358 fHaveTPCRP = kFALSE;
ada6b882 2359 }
2360 else{
cdd3a269 2361 if( fDebug >= 2 ) AliInfo(Form("Q Reaction Plane is %f", reactionPlaneQ));
f5e08ac9 2362 fHaveTPCRP = kTRUE;
ada6b882 2363 }
f5e08ac9 2364
ada6b882 2365 if(fHaveTPCRP){
e945b860 2366 fRP = ApplyFlattening(reactionPlaneQ, fCentralityV0M) ;
ada6b882 2367
2368 while(fRP<0) fRP+=TMath::Pi();
2369 while(fRP>TMath::Pi()) fRP-=TMath::Pi();
2370 FillHistogram("phiRPflat",fRP,fCentralityV0M) ;
2371 Double_t dPsi = eventPlane->GetQsubRes() ;
2372 FillHistogram("cos2AC",TMath::Cos(2.*dPsi),fCentralityV0M) ;
2373 }
f5e08ac9 2374 else
2375 fRP=0.;
ada6b882 2376}
2377
2378
2cde7444 2379//____________________________________________________________________________
ada6b882 2380void AliAnalysisTaskPi0Flow::EvalV0ReactionPlane(){
2381 // set: fRPV0A and fRPV0C
2cde7444 2382
2383 //VZERO data
ada6b882 2384 AliVVZERO* v0 = fEvent->GetVZEROData();
2cde7444 2385
ada6b882 2386 //reset Q vector info
2cde7444 2387 Double_t Qxa2 = 0, Qya2 = 0;
2388 Double_t Qxc2 = 0, Qyc2 = 0;
2389
2390 for (Int_t iv0 = 0; iv0 < 64; iv0++) {
2391 Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
ada6b882 2392 Float_t multv0 = v0->GetMultiplicity(iv0);
2cde7444 2393 if (iv0 < 32){ // V0C
2394 Qxc2 += TMath::Cos(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
2395 Qyc2 += TMath::Sin(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
2396 } else { // V0A
2397 Qxa2 += TMath::Cos(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
2398 Qya2 += TMath::Sin(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
2399 }
2400 }
2401
ada6b882 2402 Int_t iC = -1;
2cde7444 2403 // centrality bins
ada6b882 2404 if(fCentralityV0M < 5) iC = 0;
2405 else if(fCentralityV0M < 10) iC = 1;
2406 else if(fCentralityV0M < 20) iC = 2;
2407 else if(fCentralityV0M < 30) iC = 3;
2408 else if(fCentralityV0M < 40) iC = 4;
2409 else if(fCentralityV0M < 50) iC = 5;
2410 else if(fCentralityV0M < 60) iC = 6;
2411 else if(fCentralityV0M < 70) iC = 7;
2cde7444 2412 else iC = 8;
2413
2414 //grab for each centrality the proper histo with the Qx and Qy to do the recentering
2415 Double_t Qxamean2 = fMeanQ[iC][1][0];
2416 Double_t Qxarms2 = fWidthQ[iC][1][0];
2417 Double_t Qyamean2 = fMeanQ[iC][1][1];
2418 Double_t Qyarms2 = fWidthQ[iC][1][1];
ada6b882 2419
2cde7444 2420 Double_t Qxcmean2 = fMeanQ[iC][0][0];
2421 Double_t Qxcrms2 = fWidthQ[iC][0][0];
2422 Double_t Qycmean2 = fMeanQ[iC][0][1];
ada6b882 2423 Double_t Qycrms2 = fWidthQ[iC][0][1];
2424
2cde7444 2425 Double_t QxaCor2 = (Qxa2 - Qxamean2)/Qxarms2;
2426 Double_t QyaCor2 = (Qya2 - Qyamean2)/Qyarms2;
2427 Double_t QxcCor2 = (Qxc2 - Qxcmean2)/Qxcrms2;
2428 Double_t QycCor2 = (Qyc2 - Qycmean2)/Qycrms2;
ada6b882 2429
2cde7444 2430 fRPV0A = TMath::ATan2(QyaCor2, QxaCor2)/2.;
2431 fRPV0C = TMath::ATan2(QycCor2, QxcCor2)/2.;
2432
6b1c1977 2433 while (fRPV0A<0 ) fRPV0A+=TMath::Pi() ;
2434 while (fRPV0A>TMath::Pi()) fRPV0A-=TMath::Pi() ;
2435 while (fRPV0C<0 ) fRPV0C+=TMath::Pi() ;
2436 while (fRPV0C>TMath::Pi()) fRPV0C-=TMath::Pi() ;
2437
2438 // Reaction plane histograms before flattening
2439 if( fDebug >= 2 )
2440 AliInfo(Form("V0 Reaction Plane before flattening: A side: %f, C side: %f", fRPV0A, fRPV0C));
2441
2442 FillHistogram("phiRPV0A" ,fRPV0A,fCentralityV0M);
2443 FillHistogram("phiRPV0C" ,fRPV0C,fCentralityV0M);
2444 FillHistogram("phiRPV0AC",fRPV0A,fRPV0C,fCentralityV0M) ;
2cde7444 2445
f5e08ac9 2446 // Flattening
e945b860 2447 fRPV0A=ApplyFlatteningV0A(fRPV0A,fCentralityV0M) ;
6b1c1977 2448 while (fRPV0A<0 ) fRPV0A+=TMath::Pi() ;
2449 while (fRPV0A>TMath::Pi()) fRPV0A-=TMath::Pi() ;
2450
e945b860 2451 fRPV0C=ApplyFlatteningV0C(fRPV0C,fCentralityV0M) ;
6b1c1977 2452 while (fRPV0C<0 ) fRPV0C+=TMath::Pi() ;
2453 while (fRPV0C>TMath::Pi()) fRPV0C-=TMath::Pi() ;
e945b860 2454
cdd3a269 2455 if( fDebug >= 2 )
6b1c1977 2456 AliInfo(Form("V0 Reaction Plane after flattening: A side: %f, C side: %f", fRPV0A, fRPV0C));
ada6b882 2457
2458 FillHistogram("phiRPV0Aflat",fRPV0A,fCentralityV0M) ;
ada6b882 2459 FillHistogram("cos2V0AC",TMath::Cos(2.*(fRPV0A-fRPV0C)),fCentralityV0M) ;
2460 if(fHaveTPCRP){
2461 FillHistogram("phiRPV0ATPC",fRP,fRPV0A,fCentralityV0M) ;
2462 FillHistogram("cos2V0ATPC",TMath::Cos(2.*(fRP-fRPV0A)),fCentralityV0M) ;
2463 }
2464
2465 FillHistogram("phiRPV0Cflat",fRPV0C,fCentralityV0M) ;
2466 if(fHaveTPCRP){
2467 FillHistogram("phiRPV0CTPC",fRP,fRPV0C,fCentralityV0M) ;
2468 FillHistogram("cos2V0CTPC",TMath::Cos(2.*(fRP-fRPV0C)),fCentralityV0M) ;
2469 }
2cde7444 2470}