]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/PHOSTasks/PHOS_PbPb/AliAnalysisTaskPi0Flow.cxx
add trigger bins to track QA
[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),
6a1a4e92 84 fManualV0EPCalc(false),
2cde7444 85 fOutputContainer(0x0),
2cde7444 86 fNonLinCorr(0),
ada6b882 87 fEvent(0x0),
88 fEventESD(0x0),
89 fEventAOD(0x0),
90 fMCStack(0x0),
91 fRunNumber(-999),
92 fInternalRunNumber(0),
93 fPHOSGeo(0),
94 fMultV0(0x0),
95 fV0Cpol(0.),fV0Apol(0.),
96 fESDtrackCuts(0x0),
97 fPHOSCalibData(0x0),
b509e3a5 98 fEPcalibFileName("$ALICE_ROOT/OADB/PHOS/PHOSflat.root"),
99 fTPCFlat(0x0),
100 fV0AFlat(0x0),
101 fV0CFlat(0x0),
ada6b882 102 fVertexVector(),
103 fVtxBin(0),
104 fCentralityV0M(0.),
105 fCentBin(0),
106 fHaveTPCRP(0),
2cde7444 107 fRP(0),
108 fRPV0A(0),
109 fRPV0C(0),
ada6b882 110 fEMRPBin(0),
111 fCaloPhotonsPHOS(0x0),
112 fCaloPhotonsPHOSLists(0x0)
2cde7444 113{
bba0c9ef 114 const int nbins = 9;
115 Double_t edges[nbins+1] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80.};
85394e40 116 TArrayD centEdges(nbins+1, edges);
bba0c9ef 117 Int_t nMixed[nbins] = {4,4,6,10,20,30,50,100,100};
85394e40 118 TArrayI centNMixed(nbins, nMixed);
119 SetCentralityBinning(centEdges, centNMixed);
ada6b882 120
121 for(Int_t i=0;i<kNCenBins;i++){
2cde7444 122 for(Int_t j=0;j<2; j++)
ada6b882 123 for(Int_t k=0; k<2; k++) {
124 fMeanQ[i][j][k]=0.;
125 fWidthQ[i][j][k]=0.;
126 }
2cde7444 127 }
128
ada6b882 129 fVertex[0]=0; fVertex[1]=0; fVertex[2]=0;
130
2cde7444 131 // Output slots #0 write into a TH1 container
132 DefineOutput(1,TList::Class());
133
ada6b882 134
135
2cde7444 136 // Set bad channel map
137 char key[55] ;
138 for(Int_t i=0; i<6; i++){
53c2a02a 139 snprintf(key,55,"PHOS_BadMap_mod%d",i) ;
2cde7444 140 fPHOSBadMap[i]=new TH2I(key,"Bad Modules map",64,0.,64.,56,0.,56.) ;
141 }
2cde7444 142
143
2cde7444 144 // Initialize non-linrarity correction
145 fNonLinCorr = new TF1("nonlib",rnlin,0.,40.,0);
146
147
148}
ada6b882 149//___________________________________________________________________________
150AliAnalysisTaskPi0Flow::~AliAnalysisTaskPi0Flow()
151{
152 delete fNonLinCorr;
153 delete fESDtrackCuts;
154 delete fPHOSCalibData;
155 delete fCaloPhotonsPHOSLists;
b509e3a5 156 if(fTPCFlat)delete fTPCFlat; fTPCFlat=0x0;
157 if(fV0AFlat)delete fV0AFlat; fV0AFlat=0x0;
158 if(fV0CFlat)delete fV0CFlat; fV0CFlat=0x0;
159
ada6b882 160}
2cde7444 161//________________________________________________________________________
162void AliAnalysisTaskPi0Flow::UserCreateOutputObjects()
163{
164 // Create histograms
165 // Called once
166 const Int_t nRuns=200 ;
ada6b882 167
168 // histograms
2cde7444 169 if(fOutputContainer != NULL){
170 delete fOutputContainer;
171 }
9ed0694f 172 fOutputContainer = new THashList();
2cde7444 173 fOutputContainer->SetOwner(kTRUE);
ada6b882 174
2cde7444 175 //========QA histograms=======
176
177 //Event selection
ada6b882 178 fOutputContainer->Add(new TH2F("hSelEvents","Event selection", 12,0.,13.,nRuns,0.,float(nRuns))) ;
179 fOutputContainer->Add(new TH1F("hTotSelEvents","Event selection", 12,0.,12.)) ;
180
2cde7444 181 //vertex distribution
182 fOutputContainer->Add(new TH2F("hZvertex","Z vertex position", 50,-25.,25.,nRuns,0.,float(nRuns))) ;
ada6b882 183
2cde7444 184 //Centrality
185 fOutputContainer->Add(new TH2F("hCentrality","Event centrality", 100,0.,100.,nRuns,0.,float(nRuns))) ;
186 fOutputContainer->Add(new TH2F("hCenPHOS","Centrality vs PHOSclusters", 100,0.,100.,200,0.,200.)) ;
187 fOutputContainer->Add(new TH2F("hCenPHOSCells","Centrality vs PHOS cells", 100,0.,100.,100,0.,1000.)) ;
ada6b882 188 fOutputContainer->Add(new TH2F("hCenTrack","Centrality vs tracks", 100,0.,100.,100,0.,15000.)) ;
2cde7444 189 fOutputContainer->Add(new TH2F("hCluEvsClu","ClusterMult vs E",200,0.,10.,100,0.,100.)) ;
ada6b882 190
191
2cde7444 192 //Reaction plane
193 fOutputContainer->Add(new TH3F("hPHOSphi","cos" ,10,0.,100.,20,0.,10.,100,-TMath::Pi(),TMath::Pi()));
ada6b882 194
195 fOutputContainer->Add(new TH2F("cos2AC","RP correlation between TPC subs", 100,-1.,1.,20,0.,100.)) ;
196 fOutputContainer->Add(new TH2F("cos2V0AC","RP correlation between VO A and C sides", 100,-1.,1.,20,0.,100.)) ;
197 fOutputContainer->Add(new TH2F("cos2V0ATPC","RP correlation between TPC and V0A", 100,-1.,1.,20,0.,100.)) ;
198 fOutputContainer->Add(new TH2F("cos2V0CTPC","RP correlation between TPC and V0C", 100,-1.,1.,20,0.,100.)) ;
199
200 fOutputContainer->Add(new TH2F("phiRP","RP distribution with TPC", 100,0.,TMath::Pi(),20,0.,100.)) ;
201 fOutputContainer->Add(new TH2F("phiRPflat","RP distribution with TPC flat", 100,0.,TMath::Pi(),20,0.,100.)) ;
202 fOutputContainer->Add(new TH2F("phiRPV0A","RP distribution with V0A", 100,0.,TMath::Pi(),20,0.,100.)) ;
203 fOutputContainer->Add(new TH2F("phiRPV0C","RP distribution with V0C", 100,0.,TMath::Pi(),20,0.,100.)) ;
f5e08ac9 204 fOutputContainer->Add(new TH3F("phiRPV0AC","RP distribution with V0A and V0C", 100,0.,TMath::Pi(),100,0.,TMath::Pi(),20,0.,100.)) ;
ada6b882 205 fOutputContainer->Add(new TH2F("phiRPV0Aflat","RP distribution with V0 flat", 100,0.,TMath::Pi(),20,0.,100.)) ;
206 fOutputContainer->Add(new TH2F("phiRPV0Cflat","RP distribution with V0 flat", 100,0.,TMath::Pi(),20,0.,100.)) ;
207 fOutputContainer->Add(new TH3F("phiRPV0ATPC","RP distribution with V0A + TPC", 100,0.,TMath::Pi(),100,0.,TMath::Pi(),20,0.,100.)) ;
208 fOutputContainer->Add(new TH3F("phiRPV0CTPC","RP distribution with V0C + TPC", 100,0.,TMath::Pi(),100,0.,TMath::Pi(),20,0.,100.)) ;
209
210
2cde7444 211 //PHOS QA
212 fOutputContainer->Add(new TH1I("hCellMultEvent" ,"PHOS cell multiplicity per event" ,2000,0,2000));
213 fOutputContainer->Add(new TH1I("hCellMultEventM1","PHOS cell multiplicity per event, M1",2000,0,2000));
214 fOutputContainer->Add(new TH1I("hCellMultEventM2","PHOS cell multiplicity per event, M2",2000,0,2000));
215 fOutputContainer->Add(new TH1I("hCellMultEventM3","PHOS cell multiplicity per event, M3",2000,0,2000));
216
217 fOutputContainer->Add(new TH1F("hCellEnergy" ,"Cell energy" ,3000,0.,30.));
218 fOutputContainer->Add(new TH1F("hCellEnergyM1","Cell energy in module 1",3000,0.,30.));
219 fOutputContainer->Add(new TH1F("hCellEnergyM2","Cell energy in module 2",3000,0.,30.));
220 fOutputContainer->Add(new TH1F("hCellEnergyM3","Cell energy in module 3",3000,0.,30.));
221
222 fOutputContainer->Add(new TH2F("hCellNXZM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
223 fOutputContainer->Add(new TH2F("hCellNXZM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
224 fOutputContainer->Add(new TH2F("hCellNXZM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
225 fOutputContainer->Add(new TH2F("hCellEXZM1","Cell E(X,Z), M1",64,0.5,64.5, 56,0.5,56.5));
226 fOutputContainer->Add(new TH2F("hCellEXZM2","Cell E(X,Z), M2",64,0.5,64.5, 56,0.5,56.5));
227 fOutputContainer->Add(new TH2F("hCellEXZM3","Cell E(X,Z), M3",64,0.5,64.5, 56,0.5,56.5));
ada6b882 228
2cde7444 229 //Bad Map
230 fOutputContainer->Add(new TH2F("hCluLowM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
231 fOutputContainer->Add(new TH2F("hCluLowM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
232 fOutputContainer->Add(new TH2F("hCluLowM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
233
234 fOutputContainer->Add(new TH2F("hCluHighM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
235 fOutputContainer->Add(new TH2F("hCluHighM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
236 fOutputContainer->Add(new TH2F("hCluHighM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
ada6b882 237
2cde7444 238 fOutputContainer->Add(new TH2F("hCluVetoM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
239 fOutputContainer->Add(new TH2F("hCluVetoM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
240 fOutputContainer->Add(new TH2F("hCluVetoM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
241
242 fOutputContainer->Add(new TH2F("hCluDispM1","Cell (X,Z), M1" ,64,0.5,64.5, 56,0.5,56.5));
243 fOutputContainer->Add(new TH2F("hCluDispM2","Cell (X,Z), M2" ,64,0.5,64.5, 56,0.5,56.5));
244 fOutputContainer->Add(new TH2F("hCluDispM3","Cell (X,Z), M3" ,64,0.5,64.5, 56,0.5,56.5));
245
ada6b882 246
2cde7444 247 //Single photon and pi0 spectrum
56f54833 248 const Int_t nPtPhot = 400 ;
249 const Double_t ptPhotMax = 40 ;
ada6b882 250 const Int_t nM = 500;
251 const Double_t mMin = 0.0;
252 const Double_t mMax = 1.0;
253
2cde7444 254 //PHOS calibration QA
255 fOutputContainer->Add(new TH2F("hPi0M11","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
256 fOutputContainer->Add(new TH2F("hPi0M12","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
257 fOutputContainer->Add(new TH2F("hPi0M13","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
258 fOutputContainer->Add(new TH2F("hPi0M22","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
259 fOutputContainer->Add(new TH2F("hPi0M23","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
260 fOutputContainer->Add(new TH2F("hPi0M33","Pairs in modules",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
ada6b882 261
262 // Histograms for different centralities
2cde7444 263 char key[55] ;
ada6b882 264 for(Int_t cent=0; cent < fCentEdges.GetSize(); cent++){
2cde7444 265 snprintf(key,55,"hPhotAll_cen%d",cent) ;
266 fOutputContainer->Add(new TH1F(key,"All clusters",nPtPhot,0.,ptPhotMax));
267 snprintf(key,55,"hPhotAllcore_cen%d",cent) ;
268 fOutputContainer->Add(new TH1F(key,"All clusters",nPtPhot,0.,ptPhotMax));
269 snprintf(key,55,"hPhotAllwou_cen%d",cent) ;
270 fOutputContainer->Add(new TH1F(key,"All clusters",nPtPhot,0.,ptPhotMax));
271 snprintf(key,55,"hPhotDisp_cen%d",cent) ;
272 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
273 snprintf(key,55,"hPhotDisp2_cen%d",cent) ;
274 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
275 snprintf(key,55,"hPhotDispcore_cen%d",cent) ;
276 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
277 snprintf(key,55,"hPhotDispwou_cen%d",cent) ;
278 fOutputContainer->Add(new TH1F(key,"Disp clusters",nPtPhot,0.,ptPhotMax));
279 snprintf(key,55,"hPhotCPV_cen%d",cent) ;
280 fOutputContainer->Add(new TH1F(key,"CPV clusters",nPtPhot,0.,ptPhotMax));
281 snprintf(key,55,"hPhotCPVcore_cen%d",cent) ;
282 fOutputContainer->Add(new TH1F(key,"CPV clusters",nPtPhot,0.,ptPhotMax));
283 snprintf(key,55,"hPhotCPV2_cen%d",cent) ;
284 fOutputContainer->Add(new TH1F(key,"CPV clusters",nPtPhot,0.,ptPhotMax));
285 snprintf(key,55,"hPhotBoth_cen%d",cent) ;
286 fOutputContainer->Add(new TH1F(key,"Both clusters",nPtPhot,0.,ptPhotMax));
287 snprintf(key,55,"hPhotBothcore_cen%d",cent) ;
288 fOutputContainer->Add(new TH1F(key,"Both clusters",nPtPhot,0.,ptPhotMax));
289
290 snprintf(key,55,"hPi0All_cen%d",cent) ;
291 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
292 snprintf(key,55,"hPi0Allcore_cen%d",cent) ;
293 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
294 snprintf(key,55,"hPi0Allwou_cen%d",cent) ;
295 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
296 snprintf(key,55,"hPi0Disp_cen%d",cent) ;
297 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
298 snprintf(key,55,"hPi0Disp2_cen%d",cent) ;
299 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
300 snprintf(key,55,"hPi0Dispcore_cen%d",cent) ;
301 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
302 snprintf(key,55,"hPi0Dispwou_cen%d",cent) ;
303 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
304 snprintf(key,55,"hPi0CPV_cen%d",cent) ;
305 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
306 snprintf(key,55,"hPi0CPVcore_cen%d",cent) ;
307 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
308 snprintf(key,55,"hPi0CPV2_cen%d",cent) ;
309 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
310 snprintf(key,55,"hPi0Both_cen%d",cent) ;
311 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
312 snprintf(key,55,"hPi0Bothcore_cen%d",cent) ;
313 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
314
315 snprintf(key,55,"hPi0All_a07_cen%d",cent) ;
316 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
317 snprintf(key,55,"hPi0Disp_a07_cen%d",cent) ;
318 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
319 snprintf(key,55,"hPi0CPV_a07_cen%d",cent) ;
320 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
321 snprintf(key,55,"hPi0CPV2_a07_cen%d",cent) ;
322 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
323 snprintf(key,55,"hPi0Both_a07_cen%d",cent) ;
ada6b882 324 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
2cde7444 325
326 snprintf(key,55,"hSingleAll_cen%d",cent) ;
327 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
328 snprintf(key,55,"hSingleAllcore_cen%d",cent) ;
329 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
330 snprintf(key,55,"hSingleAllwou_cen%d",cent) ;
331 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
332 snprintf(key,55,"hSingleDisp_cen%d",cent) ;
333 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
334 snprintf(key,55,"hSingleDisp2_cen%d",cent) ;
335 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
336 snprintf(key,55,"hSingleDispcore_cen%d",cent) ;
337 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
338 snprintf(key,55,"hSingleDispwou_cen%d",cent) ;
339 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
340 snprintf(key,55,"hSingleCPV_cen%d",cent) ;
341 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
342 snprintf(key,55,"hSingleCPVcore_cen%d",cent) ;
343 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
344 snprintf(key,55,"hSingleCPV2_cen%d",cent) ;
345 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
346 snprintf(key,55,"hSingleBoth_cen%d",cent) ;
347 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
348 snprintf(key,55,"hSingleBothcore_cen%d",cent) ;
349 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
ada6b882 350
351
2cde7444 352 snprintf(key,55,"hMiPi0All_cen%d",cent) ;
353 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
354 snprintf(key,55,"hMiPi0Allcore_cen%d",cent) ;
355 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
356 snprintf(key,55,"hMiPi0Allwou_cen%d",cent) ;
357 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
358 snprintf(key,55,"hMiPi0Disp_cen%d",cent) ;
359 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
360 snprintf(key,55,"hMiPi0Disp2_cen%d",cent) ;
361 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
362 snprintf(key,55,"hMiPi0Dispwou_cen%d",cent) ;
363 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
364 snprintf(key,55,"hMiPi0Dispcore_cen%d",cent) ;
365 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
366 snprintf(key,55,"hMiPi0CPV_cen%d",cent) ;
367 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
368 snprintf(key,55,"hMiPi0CPVcore_cen%d",cent) ;
369 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
370 snprintf(key,55,"hMiPi0CPV2_cen%d",cent) ;
371 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
372 snprintf(key,55,"hMiPi0Both_cen%d",cent) ;
373 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
374 snprintf(key,55,"hMiPi0Bothcore_cen%d",cent) ;
375 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
ada6b882 376
2cde7444 377 snprintf(key,55,"hMiPi0All_a07_cen%d",cent) ;
378 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
379 snprintf(key,55,"hMiPi0Disp_a07_cen%d",cent) ;
380 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
381 snprintf(key,55,"hMiPi0CPV_a07_cen%d",cent) ;
382 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
383 snprintf(key,55,"hMiPi0CPV2_a07_cen%d",cent) ;
384 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
385 snprintf(key,55,"hMiPi0Both_a07_cen%d",cent) ;
ada6b882 386 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
2cde7444 387
388 snprintf(key,55,"hMiSingleAll_cen%d",cent) ;
389 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
390 snprintf(key,55,"hMiSingleAllwou_cen%d",cent) ;
391 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
392 snprintf(key,55,"hMiSingleAllcore_cen%d",cent) ;
393 fOutputContainer->Add(new TH2F(key,"All clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
394 snprintf(key,55,"hMiSingleDisp_cen%d",cent) ;
395 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
396 snprintf(key,55,"hMiSingleDisp2_cen%d",cent) ;
397 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
398 snprintf(key,55,"hMiSingleDispwou_cen%d",cent) ;
399 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
400 snprintf(key,55,"hMiSingleDispcore_cen%d",cent) ;
401 fOutputContainer->Add(new TH2F(key,"Disp clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
402 snprintf(key,55,"hMiSingleCPV_cen%d",cent) ;
403 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
404 snprintf(key,55,"hMiSingleCPVcore_cen%d",cent) ;
405 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
406 snprintf(key,55,"hMiSingleCPV2_cen%d",cent) ;
407 fOutputContainer->Add(new TH2F(key,"CPV clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
408 snprintf(key,55,"hMiSingleBoth_cen%d",cent) ;
409 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
410 snprintf(key,55,"hMiSingleBothcore_cen%d",cent) ;
411 fOutputContainer->Add(new TH2F(key,"Both clusters",nM,mMin,mMax,nPtPhot,0.,ptPhotMax));
412 }
413
2cde7444 414
ada6b882 415
416 const Int_t nPt = 20;
417 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.} ;
418 const Int_t nPhi=10 ;
419 Double_t xPhi[nPhi+1] ;
420 for(Int_t i=0;i<=nPhi;i++)
421 xPhi[i]=i*TMath::Pi() /nPhi ;
422 const Int_t nMm=200 ;
423 Double_t xM[nMm+1] ;
424 for(Int_t i=0;i<=nMm;i++)
425 xM[i]=i*0.5 /nMm;
2cde7444 426
427 char phiTitle[15] ;
428 for(Int_t iRP=0; iRP<3; iRP++){
429 if(iRP==0)
430 snprintf(phiTitle,15,"TPC") ;
431 if(iRP==1)
432 snprintf(phiTitle,15,"V0A") ;
433 if(iRP==2)
434 snprintf(phiTitle,15,"V0C") ;
ada6b882 435 for(Int_t cent=0; cent<fCentEdges.GetSize(); cent++){
436 snprintf(key,55,"hPhotPhi%sAll_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%sAllcore_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%sDisp_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%sDispcore_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%sCPV_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%sCPVcore_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%sBoth_cen%d",phiTitle,cent) ;
449 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
450 snprintf(key,55,"hPhotPhi%sBothcore_cen%d",phiTitle,cent) ;
451 fOutputContainer->Add(new TH2F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nPt,xPt,nPhi,xPhi));
452
453 //Pions
454 snprintf(key,55,"hMassPt%sAll_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%sAllcore_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%sCPV_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%sCPVcore_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%sDisp_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%sDispcore_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%sBoth_cen%d",phiTitle,cent) ;
467 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
468 snprintf(key,55,"hMassPt%sBothcore_cen%d",phiTitle,cent) ;
469 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
470
471 //Mixed
472 snprintf(key,55,"hMiMassPt%sAll_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%sAllcore_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%sCPV_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%sCPVcore_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%sDisp_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%sDispcore_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%sBoth_cen%d",phiTitle,cent) ;
485 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
486 snprintf(key,55,"hMiMassPt%sBothcore_cen%d",phiTitle,cent) ;
487 fOutputContainer->Add(new TH3F(key,"(M,p_{T},d#phi)_{#gamma#gamma}" ,nMm,xM,nPt,xPt,nPhi,xPhi));
488 }
489 }
490
491 // Setup photon lists
492 Int_t kapacity = kNVtxZBins * GetNumberOfCentralityBins() * fNEMRPBins;
493 fCaloPhotonsPHOSLists = new TObjArray(kapacity);
494 fCaloPhotonsPHOSLists->SetOwner();
2cde7444 495
496 PostData(1, fOutputContainer);
2cde7444 497}
498
499//________________________________________________________________________
ada6b882 500void AliAnalysisTaskPi0Flow::UserExec(Option_t *)
2cde7444 501{
502 // Main loop, called for each event
503 // Analyze ESD/AOD
ada6b882 504
505
506 // Step 0: Event Objects
507 fEvent = GetEvent();
508 fEventESD = dynamic_cast<AliESDEvent*> (fEvent);
509 fEventAOD = dynamic_cast<AliAODEvent*> (fEvent);
510 fMCStack = GetMCStack();
c9fb2807 511 LogProgress(0);
ada6b882 512
513
514 // Step 1: Run Number, Misalignment Matrix, and Calibration
515 // fRunNumber, fInternalRunNumber, fMultV0, fV0Cpol, fV0Apol, fMeanQ, fWidthQ
516 if( fRunNumber != fEvent->GetRunNumber()) { // Check run number
517 // this should run only at first call of UserExec(),
518 // or if task runs over multiple runs, which should not occur in normal use.
519
520 // if run number has changed, set run variables
521 fRunNumber = fEvent->GetRunNumber();
522 fInternalRunNumber = ConvertToInternalRunNumber(fRunNumber);
523 // then set misalignment and V0 calibration
524 SetGeometry();
525 SetMisalignment();
526 SetV0Calibration();
527 SetESDTrackCuts();
528 SetPHOSCalibData();
b509e3a5 529 SetFlatteningData();
2cde7444 530 }
c9fb2807 531 LogProgress(1);
532 LogSelection(0, fInternalRunNumber);
ada6b882 533
534
535 // Step 2: Vertex
536 // fVertex, fVertexVector, fVtxBin
537 SetVertex();
538 if( RejectEventVertex() ) {
2cde7444 539 PostData(1, fOutputContainer);
ada6b882 540 return; // Reject!
2cde7444 541 }
c9fb2807 542 LogProgress(2);
2cde7444 543
ada6b882 544// Step 3:
545// if(event->IsPileupFromSPD()){
546// PostData(1, fOutputContainer);
547// return; // Reject!
548// }
c9fb2807 549 LogProgress(3);
2cde7444 550
2cde7444 551
ada6b882 552 // Step 4: Centrality
553 // fCentralityV0M, fCentBin
554 SetCentrality();
555 if( RejectCentrality() ){
2cde7444 556 PostData(1, fOutputContainer);
ada6b882 557 return; // Reject!
2cde7444 558 }
c9fb2807 559 LogProgress(4);
2cde7444 560
ada6b882 561
562 // Step 5: Reaction Plane
563 // fHaveTPCRP, fRP, fRPV0A, fRPV0C, fRPBin
f5e08ac9 564 EvalReactionPlane(); //TODO: uncomment this, or at least deal with it
565 EvalV0ReactionPlane(); //TODO: uncomment this, or at least deal with it
566 fEMRPBin = GetRPBin(); //TODO: uncomment this, or at least deal with it
c9fb2807 567 LogProgress(5);
ada6b882 568
569
570 // Step 6: MC
571 // ProcessMC() ;
c9fb2807 572 LogProgress(6);
ada6b882 573
574
575 // Step 7: QA PHOS cells
576 FillPHOSCellQAHists();
c9fb2807 577 LogProgress(7);
ada6b882 578
579
580 // Step 8: Event Photons (PHOS Clusters) selection
581 SelectPhotonClusters();
582 FillSelectedClusterHistograms();
c9fb2807 583 LogProgress(8);
ada6b882 584
3d7a884f 585 if( ! fCaloPhotonsPHOS->GetEntriesFast() )
586 return;
587 else
588 LogSelection(6, fInternalRunNumber);
589
590
ada6b882 591 // Step 9: Consider pi0 (photon/cluster) pairs.
592 ConsiderPi0s();
c9fb2807 593 LogProgress(9);
ada6b882 594
595 // Step 10; Mixing
596 ConsiderPi0sMix();
c9fb2807 597 LogProgress(10);
2cde7444 598
ada6b882 599 // Step 11: Update lists
600 UpdateLists();
c9fb2807 601 LogProgress(11);
ada6b882 602
2cde7444 603
ada6b882 604 // Post output data.
605 PostData(1, fOutputContainer);
606}
2cde7444 607
ada6b882 608//________________________________________________________________________
4c5399e6 609// void AliAnalysisTaskPi0Flow::Terminate(Option_t *)
610// {
611// // Draw result to the screen
612// // Called once at the end of the query
613// // new TCanvas;
614// // TH1 * hTotSelEvents = dynamic_cast<TH1*>(fOutputContainer->FindObject("hTotSelEvents"));
615// // hTotSelEvents->Draw();
616// }
ada6b882 617//________________________________________________________________________
85394e40 618void AliAnalysisTaskPi0Flow::SetCentralityBinning(const TArrayD& edges, const TArrayI& nMixed)
ada6b882 619{
107849ac 620 // Define centrality bins by their edges
85394e40 621 if( edges.At(0) < 0.) AliFatal("lower edge less then 0");
622 if( 90. < edges.At(edges.GetSize()-1) ) AliFatal("upper edge larger then 90.");
623 for(int i=0; i<edges.GetSize()-1; ++i)
624 if(edges.At(i) > edges.At(i+1)) AliFatal("edges are not sorted");
625 if( edges.GetSize() != nMixed.GetSize()+1) AliFatal("edges and nMixed don't have appropriate relative sizes");
ada6b882 626
627 fCentEdges = edges;
107849ac 628 fCentNMixed = nMixed;
629}
2cde7444 630
85394e40 631
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() ;
f369ddfc 827 Bool_t cpvBit=kTRUE ; //No track matched by default. True means: not from charged, according to veto.
261d581c 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 }
2cde7444 1118 if(ph1->Module()==1 && ph2->Module()==1)
1119 FillHistogram("hPi0M11",p12.M(),p12.Pt() );
1120 else if(ph1->Module()==2 && ph2->Module()==2)
1121 FillHistogram("hPi0M22",p12.M(),p12.Pt() );
1122 else if(ph1->Module()==3 && ph2->Module()==3)
1123 FillHistogram("hPi0M33",p12.M(),p12.Pt() );
1124 else if(ph1->Module()==1 && ph2->Module()==2)
1125 FillHistogram("hPi0M12",p12.M(),p12.Pt() );
1126 else if(ph1->Module()==1 && ph2->Module()==3)
1127 FillHistogram("hPi0M13",p12.M(),p12.Pt() );
1128 else if(ph1->Module()==2 && ph2->Module()==3)
1129 FillHistogram("hPi0M23",p12.M(),p12.Pt() );
2cde7444 1130
1131 }
1132 }
1133 } // end of loop i2
1134 } // end of loop i1
ada6b882 1135}
1136//_____________________________________________________________________________
1137void AliAnalysisTaskPi0Flow::ConsiderPi0sMix()
1138{
1139 char key[55];
1140
1141 TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
1142
6ec40e8c 1143 for (Int_t i1=0; i1<fCaloPhotonsPHOS->GetEntriesFast(); i1++) {
ada6b882 1144 AliCaloPhoton * ph1=(AliCaloPhoton*)fCaloPhotonsPHOS->At(i1) ;
89e3079e 1145 for(Int_t evi=0; evi<arrayList->GetEntries();evi++){
ada6b882 1146 TObjArray * mixPHOS = static_cast<TObjArray*>(arrayList->At(evi));
2cde7444 1147 for(Int_t i2=0; i2<mixPHOS->GetEntriesFast();i2++){
1148 AliCaloPhoton * ph2=(AliCaloPhoton*)mixPHOS->At(i2) ;
1149 TLorentzVector p12 = *ph1 + *ph2;
1150 TLorentzVector pv12 = *(ph1->GetMomV2()) + *(ph2->GetMomV2());
1151
1152 Double_t dphiA=p12.Phi()-fRPV0A ;
1153 while(dphiA<0)dphiA+=TMath::Pi() ;
1154 while(dphiA>TMath::Pi())dphiA-=TMath::Pi() ;
1155
1156 Double_t dphiC=p12.Phi()-fRPV0C ;
1157 while(dphiC<0)dphiC+=TMath::Pi() ;
1158 while(dphiC>TMath::Pi())dphiC-=TMath::Pi() ;
1159
1160 Double_t dphiT=p12.Phi()-fRP ;
1161 while(dphiT<0)dphiT+=TMath::Pi() ;
1162 while(dphiT>TMath::Pi())dphiT-=TMath::Pi() ;
ada6b882 1163
1164
2cde7444 1165 Double_t a=TMath::Abs((ph1->E()-ph2->E())/(ph1->E()+ph2->E())) ;
ada6b882 1166
1167 snprintf(key,55,"hMiMassPtAll_cen%d",fCentBin) ;
1168 FillHistogram(Form("hMiMassPtV0AAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiA) ;
1169 FillHistogram(Form("hMiMassPtV0CAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiC) ;
2cde7444 1170 if(fHaveTPCRP)
ada6b882 1171 FillHistogram(Form("hMiMassPtTPCAll_cen%d",fCentBin),p12.M() ,p12.Pt(),dphiT) ;
2cde7444 1172
ada6b882 1173 FillHistogram(Form("hMiMassPtV0AAllcore_cen%d",fCentBin),pv12.M(), pv12.Pt(), dphiA) ;
1174 FillHistogram(Form("hMiMassPtV0CAllcore_cen%d",fCentBin),pv12.M(), pv12.Pt(), dphiC) ;
2cde7444 1175 if(fHaveTPCRP)
ada6b882 1176 FillHistogram(Form("hMiMassPtTPCAllcore_cen%d",fCentBin),pv12.M(), pv12.Pt(), dphiT) ;
2cde7444 1177
ada6b882 1178 FillHistogram(Form("hMiPi0All_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1179 FillHistogram(Form("hMiPi0Allcore_cen%d",fCentBin),pv12.M() ,pv12.Pt()) ;
2cde7444 1180 if(ph1->IsntUnfolded() && ph2->IsntUnfolded()){
ada6b882 1181 FillHistogram(Form("hMiPi0Allwou_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1182 }
1183
ada6b882 1184 FillHistogram(Form("hMiSingleAll_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1185 FillHistogram(Form("hMiSingleAll_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1186 FillHistogram(Form("hMiSingleAllcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
1187 FillHistogram(Form("hMiSingleAllcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1188 if(ph1->IsntUnfolded())
ada6b882 1189 FillHistogram(Form("hMiSingleAllwou_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1190 if(ph2->IsntUnfolded())
ada6b882 1191 FillHistogram(Form("hMiSingleAllwou_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1192 if(ph1->IsCPVOK()){
ada6b882 1193 FillHistogram(Form("hMiSingleCPV_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1194 FillHistogram(Form("hMiSingleCPVcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1195 }
1196 if(ph2->IsCPVOK()){
ada6b882 1197 FillHistogram(Form("hMiSingleCPV_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1198 FillHistogram(Form("hMiSingleCPVcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1199 }
1200 if(ph1->IsCPV2OK()){
ada6b882 1201 FillHistogram(Form("hMiSingleCPV2_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1202 }
1203 if(ph2->IsCPV2OK()){
ada6b882 1204 FillHistogram(Form("hMiSingleCPV2_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
1205 }
2cde7444 1206 if(ph1->IsDispOK()){
ada6b882 1207 FillHistogram(Form("hMiSingleDisp_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1208 if(ph1->IsntUnfolded()){
ada6b882 1209 FillHistogram(Form("hMiSingleDispwou_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
2cde7444 1210 }
ada6b882 1211 FillHistogram(Form("hMiSingleDispcore_cen%d",fCentBin),pv12.M(),ph1->GetMomV2()->Pt()) ;
2cde7444 1212 }
1213 if(ph2->IsDispOK()){
ada6b882 1214 FillHistogram(Form("hMiSingleDisp_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1215 if(ph1->IsntUnfolded()){
ada6b882 1216 FillHistogram(Form("hMiSingleDispwou_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1217 }
ada6b882 1218 FillHistogram(Form("hMiSingleDispcore_cen%d",fCentBin),pv12.M(),ph2->GetMomV2()->Pt()) ;
2cde7444 1219 }
1220 if(ph1->IsDisp2OK()){
ada6b882 1221 FillHistogram(Form("hMiSingleDisp2_cen%d",fCentBin),p12.M(),ph1->Pt()) ;
1222 }
2cde7444 1223 if(ph2->IsDisp2OK()){
ada6b882 1224 FillHistogram(Form("hMiSingleDisp2_cen%d",fCentBin),p12.M(),ph2->Pt()) ;
2cde7444 1225 }
1226 if(ph1->IsDispOK() && ph1->IsCPVOK()){
ada6b882 1227 snprintf(key,55,"hMiSingleBoth_cen%d",fCentBin) ;
2cde7444 1228 FillHistogram(key,p12.M(),ph1->Pt()) ;
ada6b882 1229 snprintf(key,55,"hMiSingleBothcore_cen%d",fCentBin) ;
2cde7444 1230 FillHistogram(key,pv12.M(),ph1->GetMomV2()->Pt()) ;
1231 }
1232 if(ph2->IsDispOK() && ph2->IsCPVOK()){
ada6b882 1233 snprintf(key,55,"hMiSingleBoth_cen%d",fCentBin) ;
2cde7444 1234 FillHistogram(key,p12.M(),ph2->Pt()) ;
ada6b882 1235 snprintf(key,55,"hMiSingleBothcore_cen%d",fCentBin) ;
2cde7444 1236 FillHistogram(key,pv12.M(),ph2->GetMomV2()->Pt()) ;
ada6b882 1237 }
1238
1239
1240
1241 if(a<kAlphaCut){
1242 snprintf(key,55,"hMiPi0All_a07_cen%d",fCentBin) ;
2cde7444 1243 FillHistogram(key,p12.M() ,p12.Pt()) ;
1244 }
1245 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
ada6b882 1246 FillHistogram(Form("hMiMassPtV0ACPV_cen%d",fCentBin),p12.M(), p12.Pt(),dphiA) ;
1247 FillHistogram(Form("hMiMassPtV0CCPV_cen%d",fCentBin),p12.M(), p12.Pt(),dphiC) ;
2cde7444 1248 if(fHaveTPCRP)
ada6b882 1249 FillHistogram(Form("hMiMassPtTPCCPV_cen%d",fCentBin),p12.M(), p12.Pt(),dphiT) ;
2cde7444 1250
ada6b882 1251 FillHistogram(Form("hMiMassPtV0ACPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiA) ;
1252 FillHistogram(Form("hMiMassPtV0CCPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiC) ;
2cde7444 1253 if(fHaveTPCRP)
ada6b882 1254 FillHistogram(Form("hMiMassPtTPCCPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt(),dphiT) ;
1255
1256 FillHistogram(Form("hMiPi0CPV_cen%d",fCentBin),p12.M(), p12.Pt()) ;
1257 FillHistogram(Form("hMiPi0CPVcore_cen%d",fCentBin),pv12.M(), pv12.Pt()) ;
2cde7444 1258
ada6b882 1259 if(a<kAlphaCut){
1260 FillHistogram(Form("hMiPi0CPV_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1261 }
1262 }
1263 if(ph1->IsCPV2OK() && ph2->IsCPV2OK()){
ada6b882 1264 FillHistogram(Form("hMiPi0CPV2_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1265
ada6b882 1266 if(a<kAlphaCut){
1267 FillHistogram(Form("hMiPi0CPV2_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1268 }
1269 }
1270 if(ph1->IsDispOK() && ph2->IsDispOK()){
ada6b882 1271 FillHistogram(Form("hMiMassPtV0ADisp_cen%d",fCentBin),p12.M(),p12.Pt(),dphiA) ;
1272 FillHistogram(Form("hMiMassPtV0CDisp_cen%d",fCentBin),p12.M(),p12.Pt(),dphiC) ;
2cde7444 1273 if(fHaveTPCRP)
ada6b882 1274 FillHistogram(Form("hMiMassPtTPCDisp_cen%d",fCentBin),p12.M(),p12.Pt(),dphiT) ;
1275
1276 FillHistogram(Form("hMiMassPtV0ADispcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiA) ;
1277 FillHistogram(Form("hMiMassPtV0CDispcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiC) ;
2cde7444 1278 if(fHaveTPCRP)
ada6b882 1279 FillHistogram(Form("hMiMassPtTPCDispcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiT) ;
1280
2cde7444 1281
ada6b882 1282 FillHistogram(Form("hMiPi0Disp_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1283 FillHistogram(Form("hMiPi0Dispcore_cen%d",fCentBin),pv12.M(),pv12.Pt()) ;
2cde7444 1284 if(ph1->IsntUnfolded() && ph2->IsntUnfolded()){
ada6b882 1285 FillHistogram(Form("hMiPi0Dispwou_cen%d",fCentBin),p12.M(),p12.Pt()) ;
2cde7444 1286 }
ada6b882 1287
1288 if(a<kAlphaCut){
1289 FillHistogram(Form("hMiPi0Disp_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1290 }
2cde7444 1291 if(ph1->IsCPVOK() && ph2->IsCPVOK()){
ada6b882 1292 FillHistogram(Form("hMiMassPtV0ABoth_cen%d",fCentBin),p12.M(),p12.Pt(),dphiA) ;
1293 FillHistogram(Form("hMiMassPtV0CBoth_cen%d",fCentBin),p12.M(),p12.Pt(),dphiC) ;
2cde7444 1294 if(fHaveTPCRP)
ada6b882 1295 FillHistogram(Form("hMiMassPtTPCBoth_cen%d",fCentBin),p12.M(),p12.Pt(),dphiT) ;
2cde7444 1296
ada6b882 1297 FillHistogram(Form("hMiMassPtV0ABothcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiA) ;
1298 FillHistogram(Form("hMiMassPtV0CBothcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiC) ;
2cde7444 1299 if(fHaveTPCRP)
ada6b882 1300 FillHistogram(Form("hMiMassPtTPCBothcore_cen%d",fCentBin),pv12.M(),pv12.Pt(),dphiT) ;
2cde7444 1301
ada6b882 1302 FillHistogram(Form("hMiPi0Both_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
1303 FillHistogram(Form("hMiPi0Bothcore_cen%d",fCentBin),pv12.M(),pv12.Pt()) ;
1304
1305 if(a<kAlphaCut){
1306 FillHistogram(Form("hMiPi0Both_a07_cen%d",fCentBin),p12.M() ,p12.Pt()) ;
2cde7444 1307 }
1308 }
1309 }
1310 } // end of loop i2
1311 }
1312 } // end of loop i1
2cde7444 1313}
ada6b882 1314//_____________________________________________________________________________
1315void AliAnalysisTaskPi0Flow::UpdateLists()
2cde7444 1316{
ada6b882 1317 //Now we either add current events to stack or remove
1318 //If no photons in current event - no need to add it to mixed
107849ac 1319
ada6b882 1320 TList * arrayList = GetCaloPhotonsPHOSList(fVtxBin, fCentBin, fEMRPBin);
cdd3a269 1321 if( fDebug >= 2 )
107849ac 1322 AliInfo( Form("fCentBin=%d, fCentNMixed[]=%d",fCentBin,fCentNMixed[fCentBin]) );
6ec40e8c 1323 if(fCaloPhotonsPHOS->GetEntriesFast()>0){
ada6b882 1324 arrayList->AddFirst(fCaloPhotonsPHOS) ;
1325 fCaloPhotonsPHOS=0;
107849ac 1326 if(arrayList->GetEntries() > fCentNMixed[fCentBin]){ // Remove redundant events
ada6b882 1327 TObjArray * tmp = static_cast<TObjArray*>(arrayList->Last()) ;
1328 arrayList->RemoveLast() ;
1329 delete tmp ; // TODO: may conflict with delete done by list being owner.
2cde7444 1330 }
2cde7444 1331 }
ada6b882 1332 else
1333 fCaloPhotonsPHOS->Clear(); // TODO: redundant???
2cde7444 1334}
1335//_____________________________________________________________________________
1336void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x)const{
1337 //FillHistogram
9ed0694f 1338 TH1 * hist = dynamic_cast<TH1*>(fOutputContainer->FindObject(key)) ;
1339 if(hist)
1340 hist->Fill(x) ;
1341 else
1342 AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2cde7444 1343}
1344//_____________________________________________________________________________
1345void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x,Double_t y)const{
1346 //FillHistogram
9ed0694f 1347 TH1 * th1 = dynamic_cast<TH1*> (fOutputContainer->FindObject(key));
1348 if(th1)
1349 th1->Fill(x, y) ;
1350 else
1351 AliError(Form("can not find histogram (of instance TH1) <%s> ",key)) ;
2cde7444 1352}
1353
1354//_____________________________________________________________________________
1355void AliAnalysisTaskPi0Flow::FillHistogram(const char * key,Double_t x,Double_t y, Double_t z) const{
1356 //Fills 1D histograms with key
9ed0694f 1357 TObject * obj = fOutputContainer->FindObject(key);
1358
1359 TH2 * th2 = dynamic_cast<TH2*> (obj);
1360 if(th2) {
1361 th2->Fill(x, y, z) ;
1362 return;
2cde7444 1363 }
9ed0694f 1364
1365 TH3 * th3 = dynamic_cast<TH3*> (obj);
1366 if(th3) {
1367 th3->Fill(x, y, z) ;
1368 return;
2cde7444 1369 }
9ed0694f 1370
1371 AliError(Form("can not find histogram (of instance TH2) <%s> ",key)) ;
2cde7444 1372}
1373
ada6b882 1374//_____________________________________________________________________________
1375AliVEvent* AliAnalysisTaskPi0Flow::GetEvent()
1376{
1377 fEvent = InputEvent();
1378 if( ! fEvent ) {
1379 AliError("Event could not be retrieved");
1380 PostData(1, fOutputContainer);
1381 }
1382 return fEvent;
1383}
1384
1385
2cde7444 1386//___________________________________________________________________________
ada6b882 1387AliStack* AliAnalysisTaskPi0Flow::GetMCStack()
1388{
1389 fMCStack = 0;
1390 AliVEventHandler* eventHandler = AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler();
1391 if(eventHandler){
1392 AliMCEventHandler* mcEventHandler = dynamic_cast<AliMCEventHandler*> (eventHandler);
1393 if( mcEventHandler)
1394 fMCStack = static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler())->MCEvent()->Stack();
1395 }
1396 return fMCStack;
1397}
2cde7444 1398
ada6b882 1399//___________________________________________________________________________
1400Int_t AliAnalysisTaskPi0Flow::GetCentralityBin(Float_t centralityV0M)
1401{
1402 /* fCentBin=1+Int_t(centralityV0M/100. *kNCenBins) ;
1403 if(centralityV0M < 5. || fCentBin < 0)
1404 fCentBin=0 ;
1405 if(fCentBin > kNCenBins-1)
1406 fCentBin = kNCenBins-1 ;
1407 */
1408 int lastBinUpperIndex = fCentEdges.GetSize() -1;
1409 if( centralityV0M > fCentEdges[lastBinUpperIndex] ) {
cdd3a269 1410 if( fDebug >= 1 )
f5e08ac9 1411 AliWarning( Form("centrality (%f) larger then upper edge of last centrality bin (%f)!", centralityV0M, fCentEdges[lastBinUpperIndex]) );
ada6b882 1412 return lastBinUpperIndex-1;
1413 }
1414 if( centralityV0M < fCentEdges[0] ) {
cdd3a269 1415 if( fDebug >= 1 )
f5e08ac9 1416 AliWarning( Form("centrality (%f) smaller then lower edge of first bin (%f)!", centralityV0M, fCentEdges[0]) );
ada6b882 1417 return 0;
1418 }
1419
1420 fCentBin = TMath::BinarySearch<Double_t> ( GetNumberOfCentralityBins(), fCentEdges.GetArray(), centralityV0M );
1421 return fCentBin;
1422}
1423
1424//___________________________________________________________________________
1425Int_t AliAnalysisTaskPi0Flow::GetRPBin()
1426{
1427 Double_t averageRP;
1428 if(fHaveTPCRP)
1429 averageRP = fRPV0A+fRPV0C+fRP /3.;
1430 else
1431 averageRP = fRPV0A+fRPV0C /2.;
1432
1433 fEMRPBin = Int_t(fNEMRPBins*(averageRP)/TMath::Pi());
1434
1435 if(fEMRPBin> (Int_t) fNEMRPBins-1)
1436 fEMRPBin=fNEMRPBins-1 ;
1437 else if(fEMRPBin<0)
1438 fEMRPBin=0;
1439
cdd3a269 1440 if ( fDebug >= 2 )
f5e08ac9 1441 AliInfo(Form("Event Mixing Reaction Plane bin is: %d", fEMRPBin));
1442
ada6b882 1443 return fEMRPBin;
1444}
1445
1446
1447//_____________________________________________________________________________
c9fb2807 1448void AliAnalysisTaskPi0Flow::LogProgress(int step)
ada6b882 1449{
cdd3a269 1450 if(fDebug >= 2) {
ada6b882 1451 AliInfo(Form("step %d completed", step));
1452 }
1453 // the +0.5 is not realy neccisarry, but oh well... -henrik
c9fb2807 1454 //FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
1455 //FillHistogram("hTotSelEvents", step+0.5);
1456}
1457
1458void AliAnalysisTaskPi0Flow::LogSelection(int step, int internalRunNumber)
1459{
1460 // if(fDebug > 1) {
1461 // AliInfo(Form("step %d completed", step));
1462 // }
1463 // the +0.5 is not realy neccisarry, but oh well... -henrik
ada6b882 1464 FillHistogram("hSelEvents", step+0.5, internalRunNumber-0.5);
1465 FillHistogram("hTotSelEvents", step+0.5);
1466}
1467
c9fb2807 1468
ada6b882 1469//___________________________________________________________________________
1470Int_t AliAnalysisTaskPi0Flow::ConvertToInternalRunNumber(Int_t run){
1471 if( kLHC11h == fPeriod ) {
1472 switch(run){
1473 case 170593 : return 179 ;
1474 case 170572 : return 178 ;
1475 case 170556 : return 177 ;
1476 case 170552 : return 176 ;
1477 case 170546 : return 175 ;
1478 case 170390 : return 174 ;
1479 case 170389 : return 173 ;
1480 case 170388 : return 172 ;
1481 case 170387 : return 171 ;
1482 case 170315 : return 170 ;
1483 case 170313 : return 169 ;
1484 case 170312 : return 168 ;
1485 case 170311 : return 167 ;
1486 case 170309 : return 166 ;
1487 case 170308 : return 165 ;
1488 case 170306 : return 164 ;
1489 case 170270 : return 163 ;
1490 case 170269 : return 162 ;
1491 case 170268 : return 161 ;
1492 case 170267 : return 160 ;
1493 case 170264 : return 159 ;
1494 case 170230 : return 158 ;
1495 case 170228 : return 157 ;
1496 case 170208 : return 156 ;
1497 case 170207 : return 155 ;
1498 case 170205 : return 154 ;
1499 case 170204 : return 153 ;
1500 case 170203 : return 152 ;
1501 case 170195 : return 151 ;
1502 case 170193 : return 150 ;
1503 case 170163 : return 149 ;
1504 case 170162 : return 148 ;
1505 case 170159 : return 147 ;
1506 case 170155 : return 146 ;
1507 case 170152 : return 145 ;
1508 case 170091 : return 144 ;
1509 case 170089 : return 143 ;
1510 case 170088 : return 142 ;
1511 case 170085 : return 141 ;
1512 case 170084 : return 140 ;
1513 case 170083 : return 139 ;
1514 case 170081 : return 138 ;
1515 case 170040 : return 137 ;
1516 case 170038 : return 136 ;
1517 case 170036 : return 135 ;
1518 case 170027 : return 134 ;
1519 case 169981 : return 133 ;
1520 case 169975 : return 132 ;
1521 case 169969 : return 131 ;
1522 case 169965 : return 130 ;
1523 case 169961 : return 129 ;
1524 case 169956 : return 128 ;
1525 case 169926 : return 127 ;
1526 case 169924 : return 126 ;
1527 case 169923 : return 125 ;
1528 case 169922 : return 124 ;
1529 case 169919 : return 123 ;
1530 case 169918 : return 122 ;
1531 case 169914 : return 121 ;
1532 case 169859 : return 120 ;
1533 case 169858 : return 119 ;
1534 case 169855 : return 118 ;
1535 case 169846 : return 117 ;
1536 case 169838 : return 116 ;
1537 case 169837 : return 115 ;
1538 case 169835 : return 114 ;
1539 case 169683 : return 113 ;
1540 case 169628 : return 112 ;
1541 case 169591 : return 111 ;
1542 case 169590 : return 110 ;
1543 case 169588 : return 109 ;
1544 case 169587 : return 108 ;
1545 case 169586 : return 107 ;
1546 case 169584 : return 106 ;
1547 case 169557 : return 105 ;
1548 case 169555 : return 104 ;
1549 case 169554 : return 103 ;
1550 case 169553 : return 102 ;
1551 case 169550 : return 101 ;
1552 case 169515 : return 100 ;
1553 case 169512 : return 99 ;
1554 case 169506 : return 98 ;
1555 case 169504 : return 97 ;
1556 case 169498 : return 96 ;
1557 case 169475 : return 95 ;
1558 case 169420 : return 94 ;
1559 case 169419 : return 93 ;
1560 case 169418 : return 92 ;
1561 case 169417 : return 91 ;
1562 case 169415 : return 90 ;
1563 case 169411 : return 89 ;
1564 case 169238 : return 88 ;
1565 case 169236 : return 87 ;
1566 case 169167 : return 86 ;
1567 case 169160 : return 85 ;
1568 case 169156 : return 84 ;
1569 case 169148 : return 83 ;
1570 case 169145 : return 82 ;
1571 case 169144 : return 81 ;
1572 case 169143 : return 80 ;
1573 case 169138 : return 79 ;
1574 case 169099 : return 78 ;
1575 case 169094 : return 77 ;
1576 case 169091 : return 76 ;
1577 case 169045 : return 75 ;
1578 case 169044 : return 74 ;
1579 case 169040 : return 73 ;
1580 case 169035 : return 72 ;
1581 case 168992 : return 71 ;
1582 case 168988 : return 70 ;
1583 case 168984 : return 69 ;
1584 case 168826 : return 68 ;
1585 case 168777 : return 67 ;
1586 case 168514 : return 66 ;
1587 case 168512 : return 65 ;
1588 case 168511 : return 64 ;
1589 case 168467 : return 63 ;
1590 case 168464 : return 62 ;
1591 case 168461 : return 61 ;
1592 case 168460 : return 60 ;
1593 case 168458 : return 59 ;
1594 case 168362 : return 58 ;
1595 case 168361 : return 57 ;
1596 case 168356 : return 56 ;
1597 case 168342 : return 55 ;
1598 case 168341 : return 54 ;
1599 case 168325 : return 53 ;
1600 case 168322 : return 52 ;
1601 case 168318 : return 51 ;
1602 case 168311 : return 50 ;
1603 case 168310 : return 49 ;
1604 case 168213 : return 48 ;
1605 case 168212 : return 47 ;
1606 case 168208 : return 46 ;
1607 case 168207 : return 45 ;
1608 case 168206 : return 44 ;
1609 case 168205 : return 43 ;
1610 case 168204 : return 42 ;
1611 case 168203 : return 41 ;
1612 case 168181 : return 40 ;
1613 case 168177 : return 39 ;
1614 case 168175 : return 38 ;
1615 case 168173 : return 37 ;
1616 case 168172 : return 36 ;
1617 case 168171 : return 35 ;
1618 case 168115 : return 34 ;
1619 case 168108 : return 33 ;
1620 case 168107 : return 32 ;
1621 case 168105 : return 31 ;
1622 case 168104 : return 30 ;
1623 case 168103 : return 29 ;
1624 case 168076 : return 28 ;
1625 case 168069 : return 27 ;
1626 case 168068 : return 26 ;
1627 case 168066 : return 25 ;
1628 case 167988 : return 24 ;
1629 case 167987 : return 23 ;
1630 case 167986 : return 22 ;
1631 case 167985 : return 21 ;
1632 case 167921 : return 20 ;
1633 case 167920 : return 19 ;
1634 case 167915 : return 18 ;
1635 case 167909 : return 17 ;
1636 case 167903 : return 16 ;
1637 case 167902 : return 15 ;
1638 case 167818 : return 14 ;
1639 case 167814 : return 13 ;
1640 case 167813 : return 12 ;
1641 case 167808 : return 11 ;
1642 case 167807 : return 10 ;
1643 case 167806 : return 9 ;
1644 case 167713 : return 8 ;
1645 case 167712 : return 7 ;
1646 case 167711 : return 6 ;
1647 case 167706 : return 5 ;
1648 case 167693 : return 4 ;
1649 case 166532 : return 3 ;
1650 case 166530 : return 2 ;
1651 case 166529 : return 1 ;
1652
1653 default : return 199;
1654 }
1655 }
1656 if( kLHC10h == fPeriod ) {
1657 switch(run){
1658 case 139517 : return 137;
1659 case 139514 : return 136;
1660 case 139513 : return 135;
1661 case 139511 : return 134;
1662 case 139510 : return 133;
1663 case 139507 : return 132;
1664 case 139505 : return 131;
1665 case 139504 : return 130;
1666 case 139503 : return 129;
1667 case 139470 : return 128;
1668 case 139467 : return 127;
1669 case 139466 : return 126;
1670 case 139465 : return 125;
1671 case 139440 : return 124;
1672 case 139439 : return 123;
1673 case 139438 : return 122;
1674 case 139437 : return 121;
1675 case 139360 : return 120;
1676 case 139329 : return 119;
1677 case 139328 : return 118;
1678 case 139314 : return 117;
1679 case 139311 : return 116;
1680 case 139310 : return 115;
1681 case 139309 : return 114;
1682 case 139308 : return 113;
1683 case 139173 : return 112;
1684 case 139172 : return 111;
1685 case 139110 : return 110;
1686 case 139107 : return 109;
1687 case 139105 : return 108;
1688 case 139104 : return 107;
1689 case 139042 : return 106;
1690 case 139038 : return 105;
1691 case 139037 : return 104;
1692 case 139036 : return 103;
1693 case 139029 : return 102;
1694 case 139028 : return 101;
1695 case 138983 : return 100;
1696 case 138982 : return 99;
1697 case 138980 : return 98;
1698 case 138979 : return 97;
1699 case 138978 : return 96;
1700 case 138977 : return 95;
1701 case 138976 : return 94;
1702 case 138973 : return 93;
1703 case 138972 : return 92;
1704 case 138965 : return 91;
1705 case 138924 : return 90;
1706 case 138872 : return 89;
1707 case 138871 : return 88;
1708 case 138870 : return 87;
1709 case 138837 : return 86;
1710 case 138830 : return 85;
1711 case 138828 : return 84;
1712 case 138826 : return 83;
1713 case 138796 : return 82;
1714 case 138795 : return 81;
1715 case 138742 : return 80;
1716 case 138732 : return 79;
1717 case 138730 : return 78;
1718 case 138666 : return 77;
1719 case 138662 : return 76;
1720 case 138653 : return 75;
1721 case 138652 : return 74;
1722 case 138638 : return 73;
1723 case 138624 : return 72;
1724 case 138621 : return 71;
1725 case 138583 : return 70;
1726 case 138582 : return 69;
1727 case 138579 : return 68;
1728 case 138578 : return 67;
1729 case 138534 : return 66;
1730 case 138469 : return 65;
1731 case 138442 : return 64;
1732 case 138439 : return 63;
1733 case 138438 : return 62;
1734 case 138396 : return 61;
1735 case 138364 : return 60;
1736 case 138359 : return 59;
1737 case 138275 : return 58;
1738 case 138225 : return 57;
1739 case 138201 : return 56;
1740 case 138200 : return 55;
1741 case 138197 : return 54;
1742 case 138192 : return 53;
1743 case 138190 : return 52;
1744 case 138154 : return 51;
1745 case 138153 : return 50;
1746 case 138151 : return 49;
1747 case 138150 : return 48;
1748 case 138126 : return 47;
1749 case 138125 : return 46;
1750 case 137848 : return 45;
1751 case 137847 : return 44;
1752 case 137844 : return 43;
1753 case 137843 : return 42;
1754 case 137752 : return 41;
1755 case 137751 : return 40;
1756 case 137748 : return 39;
1757 case 137724 : return 38;
1758 case 137722 : return 37;
1759 case 137718 : return 36;
1760 case 137704 : return 35;
1761 case 137693 : return 34;
1762 case 137692 : return 33;
1763 case 137691 : return 32;
1764 case 137689 : return 31;
1765 case 137686 : return 30;
1766 case 137685 : return 29;
1767 case 137639 : return 28;
1768 case 137638 : return 27;
1769 case 137608 : return 26;
1770 case 137595 : return 25;
1771 case 137549 : return 24;
1772 case 137546 : return 23;
1773 case 137544 : return 22;
1774 case 137541 : return 21;
1775 case 137539 : return 20;
1776 case 137531 : return 19;
1777 case 137530 : return 18;
1778 case 137443 : return 17;
1779 case 137441 : return 16;
1780 case 137440 : return 15;
1781 case 137439 : return 14;
1782 case 137434 : return 13;
1783 case 137432 : return 12;
1784 case 137431 : return 11;
1785 case 137430 : return 10;
1786 case 137366 : return 9;
1787 case 137243 : return 8;
1788 case 137236 : return 7;
1789 case 137235 : return 6;
1790 case 137232 : return 5;
1791 case 137231 : return 4;
1792 case 137165 : return 3;
1793 case 137162 : return 2;
1794 case 137161 : return 1;
1795 default : return 199;
1796 }
1797 }
cdd3a269 1798 if(kUndefinedPeriod && fDebug >= 1 ) {
ada6b882 1799 AliWarning("Period not defined");
1800 }
1801 return 1;
2cde7444 1802}
1803//_____________________________________________________________________________
1804Bool_t AliAnalysisTaskPi0Flow::TestLambda(Double_t pt,Double_t l1,Double_t l2){
ada6b882 1805
2cde7444 1806 Double_t l2Mean = 1.53126+9.50835e+06/(1.+1.08728e+07*pt+1.73420e+06*pt*pt) ;
1807 Double_t l1Mean = 1.12365+0.123770*TMath::Exp(-pt*0.246551)+5.30000e-03*pt ;
1808 Double_t l2Sigma = 6.48260e-02+7.60261e+10/(1.+1.53012e+11*pt+5.01265e+05*pt*pt)+9.00000e-03*pt;
1809 Double_t l1Sigma = 4.44719e-04+6.99839e-01/(1.+1.22497e+00*pt+6.78604e-07*pt*pt)+9.00000e-03*pt;
1810 Double_t c=-0.35-0.550*TMath::Exp(-0.390730*pt) ;
ada6b882 1811 Double_t R2=0.5*(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma +
1812 0.5*(l2-l2Mean)*(l2-l2Mean)/l2Sigma/l2Sigma +
1813 0.5*c*(l1-l1Mean)*(l2-l2Mean)/l1Sigma/l2Sigma ;
2cde7444 1814 return (R2<2.5*2.5) ;
ada6b882 1815
2cde7444 1816}
1817//_____________________________________________________________________________
1818Bool_t AliAnalysisTaskPi0Flow::TestLambda2(Double_t pt,Double_t l1,Double_t l2){
ada6b882 1819
2cde7444 1820 Double_t l2Mean = 1.53126+9.50835e+06/(1.+1.08728e+07*pt+1.73420e+06*pt*pt) ;
1821 Double_t l1Mean = 1.12365+0.123770*TMath::Exp(-pt*0.246551)+5.30000e-03*pt ;
1822 Double_t l2Sigma = 6.48260e-02+7.60261e+10/(1.+1.53012e+11*pt+5.01265e+05*pt*pt)+9.00000e-03*pt;
1823 Double_t l1Sigma = 4.44719e-04+6.99839e-01/(1.+1.22497e+00*pt+6.78604e-07*pt*pt)+9.00000e-03*pt;
1824 Double_t c=-0.35-0.550*TMath::Exp(-0.390730*pt) ;
ada6b882 1825 Double_t R2=0.5*(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma +
1826 0.5*(l2-l2Mean)*(l2-l2Mean)/l2Sigma/l2Sigma +
1827 0.5*c*(l1-l1Mean)*(l2-l2Mean)/l1Sigma/l2Sigma ;
2cde7444 1828 return (R2<1.5*1.5) ;
ada6b882 1829
2cde7444 1830}
ada6b882 1831//____________________________________________________________________________
1832TList* AliAnalysisTaskPi0Flow::GetCaloPhotonsPHOSList(UInt_t vtxBin, UInt_t centBin, UInt_t rpBin)
1833{
1834 int offset = vtxBin * GetNumberOfCentralityBins() * fNEMRPBins
1835 + centBin * fNEMRPBins
1836 + rpBin;
1837 if( fCaloPhotonsPHOSLists->At(offset) ) { // list exists
1838 TList* list = dynamic_cast<TList*> (fCaloPhotonsPHOSLists->At(offset));
1839 if( ! list )
1840 AliError("object in fCaloPhotonsPHOSLists at %i did not cast");
1841 return list;
1842 }
1843 else {// no list for this bin has been created, yet
1844 TList* list = new TList();
1845 list->SetOwner();
1846 fCaloPhotonsPHOSLists->AddAt(list, offset);
1847 return list;
1848 }
1849}
1850
2cde7444 1851//____________________________________________________________________________
1852Double_t AliAnalysisTaskPi0Flow::TestCPV(Double_t dx, Double_t dz, Double_t pt, Int_t charge){
1853 //Parameterization of LHC10h period
1854 //_true if neutral_
ada6b882 1855
2cde7444 1856 Double_t meanX=0;
1857 Double_t meanZ=0.;
1858 Double_t sx=TMath::Min(5.4,2.59719e+02*TMath::Exp(-pt/1.02053e-01)+
ada6b882 1859 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 1860 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 1861 AliVEvent *event = InputEvent();
53c2a02a 1862 Double_t mf = 0.; //
1863 if(event) mf = event->GetMagneticField(); //Positive for ++ and negative for --
2cde7444 1864
1865 if(mf<0.){ //field --
1866 meanZ = -0.468318 ;
1867 if(charge>0)
1868 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)) ;
1869 else
1870 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)) ;
1871 }
1872 else{ //Field ++
1873 meanZ= -0.468318;
1874 if(charge>0)
1875 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)) ;
1876 else
ada6b882 1877 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 1878 }
1879
1880 Double_t rz=(dz-meanZ)/sz ;
1881 Double_t rx=(dx-meanX)/sx ;
1882 return TMath::Sqrt(rx*rx+rz*rz) ;
1883}
1884//____________________________________________________________________________
b509e3a5 1885void AliAnalysisTaskPi0Flow::SetFlatteningData(){
1886 //Read objects with flattening parameters
1887 AliOADBContainer flatContainer("phosFlat");
1888 flatContainer.InitFromFile(fEPcalibFileName.Data(),"phosFlat");
1889 TObjArray *maps = (TObjArray*)flatContainer.GetObject(fRunNumber,"phosFlat");
1890 if(!maps){
1891 AliError(Form("Can not read Flattening for run %d. \n From file >%s<\n",fRunNumber,fEPcalibFileName.Data())) ;
2cde7444 1892 }
b509e3a5 1893 else{
1894 AliInfo(Form("Setting PHOS flattening with name %s \n",maps->GetName())) ;
1895 AliPHOSEPFlattener * h = (AliPHOSEPFlattener*)maps->At(0) ;
1896 if(fTPCFlat) delete fTPCFlat ;
1897 fTPCFlat = new AliPHOSEPFlattener() ;
1898 fTPCFlat = h ;
1899 h = (AliPHOSEPFlattener*)maps->At(1) ;
1900 if(fV0AFlat) delete fV0AFlat ;
1901 fV0AFlat = new AliPHOSEPFlattener() ;
1902 fV0AFlat = h ;
1903 h = (AliPHOSEPFlattener*)maps->At(2) ;
1904 if(fV0CFlat) delete fV0CFlat ;
1905 fV0CFlat = new AliPHOSEPFlattener() ;
1906 fV0CFlat = h ;
1907 }
1908
2cde7444 1909}
b509e3a5 1910 //____________________________________________________________________________
1911Double_t AliAnalysisTaskPi0Flow::ApplyFlattening(Double_t phi, Double_t c){
1912
1913 if(fTPCFlat)
1914 return fTPCFlat->MakeFlat(phi,c);
1915 return phi ;
2cde7444 1916
2cde7444 1917}
b509e3a5 1918//____________________________________________________________________________
1919Double_t AliAnalysisTaskPi0Flow::ApplyFlatteningV0A(Double_t phi, Double_t c){
1920
1921 if(fV0AFlat)
1922 return fV0AFlat->MakeFlat(phi,c);
1923 return phi ;
2cde7444 1924
1925}
1926//____________________________________________________________________________
b509e3a5 1927Double_t AliAnalysisTaskPi0Flow::ApplyFlatteningV0C(Double_t phi, Double_t c){
1928
1929 if(fV0CFlat)
1930 return fV0CFlat->MakeFlat(phi,c);
ada6b882 1931 return phi ;
1932
1933}
2cde7444 1934//____________________________________________________________________________
91b112bd 1935Double_t AliAnalysisTaskPi0Flow::CoreEnergy(AliVCluster * clu, AliVCaloCells * cells)
1936{
2cde7444 1937 //calculate energy of the cluster in the circle with radius distanceCut around the maximum
ada6b882 1938
2cde7444 1939 //Can not use already calculated coordinates?
1940 //They have incidence correction...
1941 const Double_t distanceCut =3.5 ;
1942 const Double_t logWeight=4.5 ;
ada6b882 1943
1944 Double32_t * elist = clu->GetCellsAmplitudeFraction() ;
2cde7444 1945// Calculates the center of gravity in the local PHOS-module coordinates
1946 Float_t wtot = 0;
1947 Double_t xc[100]={0} ;
1948 Double_t zc[100]={0} ;
1949 Double_t x = 0 ;
1950 Double_t z = 0 ;
1951 Int_t mulDigit=TMath::Min(100,clu->GetNCells()) ;
1952 for(Int_t iDigit=0; iDigit<mulDigit; iDigit++) {
1953 Int_t relid[4] ;
1954 Float_t xi ;
1955 Float_t zi ;
1956 fPHOSGeo->AbsToRelNumbering(clu->GetCellAbsId(iDigit), relid) ;
1957 fPHOSGeo->RelPosInModule(relid, xi, zi);
1958 xc[iDigit]=xi ;
1959 zc[iDigit]=zi ;
91b112bd 1960 elist[iDigit] *= cells->GetCellAmplitude(clu->GetCellsAbsId()[iDigit]);
cdd3a269 1961 if( fDebug >= 3 )
1962 printf("%f ",elist[iDigit]);
2cde7444 1963 if (clu->E()>0 && elist[iDigit]>0) {
1964 Float_t w = TMath::Max( 0., logWeight + TMath::Log( elist[iDigit] / clu->E() ) ) ;
1965 x += xc[iDigit] * w ;
1966 z += zc[iDigit] * w ;
1967 wtot += w ;
1968 }
1969 }
1970 if (wtot>0) {
1971 x /= wtot ;
1972 z /= wtot ;
1973 }
1974 Double_t coreE=0. ;
1975 for(Int_t iDigit=0; iDigit < mulDigit; iDigit++) {
1976 Double_t distance = TMath::Sqrt((xc[iDigit]-x)*(xc[iDigit]-x)+(zc[iDigit]-z)*(zc[iDigit]-z)) ;
1977 if(distance < distanceCut)
1978 coreE += elist[iDigit] ;
1979 }
1980 //Apply non-linearity correction
1981 return fNonLinCorr->Eval(coreE) ;
1982}
1983//____________________________________________________________________________
1984Bool_t AliAnalysisTaskPi0Flow::AreNeibors(Int_t id1,Int_t id2){
ada6b882 1985 // return true if absId are "Neighbors" (adjacent, including diagornaly,)
1986 // false if not.
2cde7444 1987
ada6b882 1988 Int_t relid1[4] ;
1989 fPHOSGeo->AbsToRelNumbering(id1, relid1) ;
2cde7444 1990
ada6b882 1991 Int_t relid2[4] ;
1992 fPHOSGeo->AbsToRelNumbering(id2, relid2) ;
2cde7444 1993
ada6b882 1994 // if inside the same PHOS module
1995 if ( (relid1[0] == relid2[0]) && (relid1[1]==relid2[1]) ) {
1996 const Int_t rowdiff = TMath::Abs( relid1[2] - relid2[2] ) ;
1997 const Int_t coldiff = TMath::Abs( relid1[3] - relid2[3] ) ;
1998
1999 // and if diff in both direction is 1 or less
2000 if (( coldiff <= 1 ) && ( rowdiff <= 1 ))
2001 return true; // are neighbors
2002 }
2cde7444 2003
ada6b882 2004 // else false
2005 return false;
2cde7444 2006}
2007//____________________________________________________________________________
ada6b882 2008void AliAnalysisTaskPi0Flow::Reclusterize(AliVCluster * clu){
2cde7444 2009 //Re-clusterize to make continues cluster
ada6b882 2010
2cde7444 2011 const Int_t oldMulDigit=clu->GetNCells() ;
ada6b882 2012 Double32_t * elist = clu->GetCellsAmplitudeFraction() ;
2cde7444 2013 UShort_t * dlist = clu->GetCellsAbsId();
2014
2015 Int_t index[oldMulDigit] ;
2016 Bool_t used[oldMulDigit] ;
2017 for(Int_t i=0; i<oldMulDigit; i++) used[i]=0 ;
2018 Int_t inClu=0 ;
2019 Double_t eMax=0. ;
2020 //find maximum
2021 for(Int_t iDigit=0; iDigit<oldMulDigit; iDigit++) {
2022 if(eMax<elist[iDigit]){
2023 eMax=elist[iDigit];
2024 index[0]=iDigit ;
2025 inClu=1 ;
2026 }
2027 }
2028 if(inClu==0){ //empty cluster
2029 return ;
2030 }
2031 used[index[0]]=kTRUE ; //mark as used
2032 for(Int_t i=0; i<inClu; i++){
2033 for(Int_t iDigit=0 ;iDigit<oldMulDigit; iDigit++){
2034 if(used[iDigit]) //already used
2035 continue ;
2036 if(AreNeibors(dlist[index[i]],dlist[iDigit])){
2037 index[inClu]= iDigit ;
2038 inClu++ ;
2039 used[iDigit]=kTRUE ;
2040 }
2041 }
2042 }
ada6b882 2043
2cde7444 2044 if(inClu==oldMulDigit) //no need to modify
ada6b882 2045 return ;
2cde7444 2046
2047 clu->SetNCells(inClu);
2048 //copy
2049 UShort_t tmpD[oldMulDigit] ;
2050 Double_t tmpE[oldMulDigit] ;
2051 for(Int_t i=0; i<oldMulDigit; i++){
ada6b882 2052 tmpD[i]=dlist[i] ;
2cde7444 2053 tmpE[i]=elist[i] ;
2054 }
2055 //change order of digits in list so that
2056 //first inClu cells were true ones
2057 for(Int_t i=0; i<inClu; i++){
2058 dlist[i]=tmpD[index[i]] ;
2059 elist[i]=tmpE[index[i]] ;
2060 }
ada6b882 2061
2062
2063}
2064
2065//_____________________________________________________________________________
2066void AliAnalysisTaskPi0Flow::SetMisalignment(){
2067 // sets the misalignment vertex if ESD
2068 if( fEventESD ) {
2069 for(Int_t mod=0; mod<5; mod++) {
2070 const TGeoHMatrix* modMatrix = fEvent->GetPHOSMatrix(mod);
2071 if( ! modMatrix) {
f5e08ac9 2072 if( fDebug )
ada6b882 2073 AliInfo(Form("no PHOS Geometric Misalignment Matrix for module %d", mod));
2074 continue;
2075 }
2076 else {
2077 fPHOSGeo->SetMisalMatrix(modMatrix, mod);
f5e08ac9 2078 if( fDebug )
cdd3a269 2079 AliInfo(Form("PHOS Geometric Misalignment Matrix set for module %d", mod));
ada6b882 2080 }
2081 }
2082 }
2cde7444 2083}
2084
2085//_____________________________________________________________________________
ada6b882 2086void AliAnalysisTaskPi0Flow::SetV0Calibration(){
2087 // assigns: fMultV0, fV0Cpol, fV0Apol, fMeanQ, and fWidthQ
2088
6a1a4e92 2089 if ( ! fManualV0EPCalc ) {
2090 if( fDebug >=2 )
2091 AliInfo("Not setting V0Calibration, only needed for manual V0 EP Calculation");
2092 return;
2093 }
2094
ada6b882 2095 int runNumber = this->fRunNumber;
2096
2cde7444 2097 TString oadbfilename = "$ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root";
2098 TFile *foadb = TFile::Open(oadbfilename.Data());
2099
2100 if(!foadb){
ada6b882 2101 AliError(Form("OADB file %s cannot be opened\n", oadbfilename.Data()));
2102 AliError("V0 Calibration not set !\n");
2cde7444 2103 return;
2104 }
2105
2106 AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr");
2107 if(!cont){
ada6b882 2108 AliError("OADB object hMultV0BefCorr is not available in the file");
2109 AliError("V0 Calibration not set!\n");
2110 return;
2cde7444 2111 }
2112
ada6b882 2113 if(!(cont->GetObject(runNumber))){
2114 AliError(Form("OADB object hMultV0BefCorr is not available for run %i, trying 137366)",runNumber));
2115 runNumber = 137366;
2116 }
2117 if(!(cont->GetObject(runNumber))){
2118 AliError(Form("OADB object hMultV0BefCorr is not available for run %i ",runNumber));
2119 AliError("V0 Calibration not set!\n");
2120 return;
2cde7444 2121 }
2cde7444 2122
f5e08ac9 2123 if( fDebug ) AliInfo("Setting V0 calibration") ;
ada6b882 2124 fMultV0 = ((TH2F *) cont->GetObject(runNumber))->ProfileX();
2125
2126 TF1 *fpol0 = new TF1("fpol0","pol0");
e2139d72 2127 fMultV0->Fit(fpol0,"Q0","",0,31);
2cde7444 2128 fV0Cpol = fpol0->GetParameter(0);
e2139d72 2129 fMultV0->Fit(fpol0,"Q0","",32,64);
2cde7444 2130 fV0Apol = fpol0->GetParameter(0);
2131
2132 for(Int_t iside=0;iside<2;iside++){
2133 for(Int_t icoord=0;icoord<2;icoord++){
ada6b882 2134 for(Int_t i=0;i < kNCenBins;i++){
2cde7444 2135 char namecont[100];
2136 if(iside==0 && icoord==0)
53c2a02a 2137 snprintf(namecont,100,"hQxc2_%i",i);
2cde7444 2138 else if(iside==1 && icoord==0)
53c2a02a 2139 snprintf(namecont,100,"hQxa2_%i",i);
2cde7444 2140 else if(iside==0 && icoord==1)
53c2a02a 2141 snprintf(namecont,100,"hQyc2_%i",i);
2cde7444 2142 else if(iside==1 && icoord==1)
53c2a02a 2143 snprintf(namecont,100,"hQya2_%i",i);
2cde7444 2144
2145 cont = (AliOADBContainer*) foadb->Get(namecont);
2146 if(!cont){
ada6b882 2147 AliError(Form("OADB object %s is not available in the file %s", namecont, oadbfilename.Data()));
2148 AliError("V0 Calibration not fully set!\n");
2149 return;
2cde7444 2150 }
2cde7444 2151
ada6b882 2152 if(!(cont->GetObject(runNumber))){
2153 AliError(Form("OADB object %s is not available for run %i, trying run 137366",namecont,runNumber));
2154 runNumber = 137366;
2cde7444 2155 }
ada6b882 2156 if(!(cont->GetObject(runNumber))){
2157 AliError(Form("OADB object %s is not available for run %i",namecont,runNumber));
2158 AliError("V0 Calibration not fully set!\n");
2159 return;
2cde7444 2160 }
ada6b882 2161 fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(runNumber))->GetMean();
2162 fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(runNumber))->GetRMS();
2163
2164 //for v3
2165// if(iside==0 && icoord==0)
2166// snprintf(namecont,100,"hQxc3_%i",i);
2167// else if(iside==1 && icoord==0)
2168// snprintf(namecont,100,"hQxa3_%i",i);
2169// else if(iside==0 && icoord==1)
2170// snprintf(namecont,100,"hQyc3_%i",i);
2171// else if(iside==1 && icoord==1)
2172// snprintf(namecont,100,"hQya3_%i",i);
2173//
2174// cont = (AliOADBContainer*) foadb->Get(namecont);
2175// if(!cont){
2176// AliError(Form("OADB object %s is not available in the file",namecont));
2177// AliError("V0 Calibration not fully set!\n");
2178// return;
2179// }
2180//
2181// if(!(cont->GetObject(runNumber))){
2182// AliError(Form("OADB object %s is not available for run %i, trying run 137366",namecont,runNumber));
2183// runNumber = 137366;
2184// }
2185// if(!(cont->GetObject(runNumber))){
2186// AliError(Form("OADB object %s is not available for run %i",namecont,runNumber));
2187// AliError("V0 Calibration not fully set!\n");
2188// return;
2189// }
2cde7444 2190// fMeanQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
2191// fWidthQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
2192
2193 }
2194 }
2195 }
ada6b882 2196
2197 delete fpol0; fpol0=0;
2cde7444 2198}
ada6b882 2199
2200//_____________________________________________________________________________
2201void AliAnalysisTaskPi0Flow::SetESDTrackCuts()
2202{
2203 if( fEventESD ) {
2204 // Create ESD track cut
2205 fESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts() ;
2206 //fESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
2207 fESDtrackCuts->SetRequireTPCRefit(kTRUE);
2208 }
2209}
2210
2211//_____________________________________________________________________________
2212void AliAnalysisTaskPi0Flow::SetGeometry()
2213{
2214 // Initialize the PHOS geometry
f5e08ac9 2215 if( kLHC10h == fPeriod && fEventESD ) {
ada6b882 2216 TGeoManager::Import("geometry.root"); //TODO: should perhaps not be done
2217 fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP") ;
2218 if( ! fPHOSGeo )
2219 AliError("geometry (fPHOSGeo) not initialised");
2220 }
2221
2222 //Init geometry
2223 if(!fPHOSGeo){
2224 AliOADBContainer geomContainer("phosGeo");
2225 geomContainer.InitFromFile("$ALICE_ROOT/OADB/PHOS/PHOSGeometry.root","PHOSRotationMatrixes");
2226 TObjArray *matrixes = (TObjArray*)geomContainer.GetObject(fRunNumber,"PHOSRotationMatrixes");
2227 fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP") ;
2228 for(Int_t mod=0; mod<5; mod++) {
2229 if(!matrixes->At(mod)) {
2230 continue;
cdd3a269 2231 if( fDebug )
ada6b882 2232 AliInfo(Form("No PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2233 }
2234 else {
2235 fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod) ;
2236 if( fDebug >1 )
2237 AliInfo(Form("Adding PHOS Matrix for mod:%d, geo=%p\n", mod, fPHOSGeo));
2238 }
2239 }
2240 }
2241}
2242
2243//_____________________________________________________________________________
2244void AliAnalysisTaskPi0Flow::SetPHOSCalibData()
2245{
2246 if( fPHOSCalibData )
2247 delete fPHOSCalibData;
2248 fPHOSCalibData = 0;
2cde7444 2249
ada6b882 2250 // Calibration only needed for ESD
2251 if( fEventESD /*&& */ ) {
f5e08ac9 2252 if( kLHC10h == fPeriod && fEventESD ) {
ada6b882 2253 //We have to apply re-calibration for pass1 LCH10h
2254 // Initialize decalibration factors in the form of the OCDB object
2255 AliCDBManager * man = AliCDBManager::Instance();
2256 man->SetRun(140000) ; //TODO; revise, this should probably not b done.
2257 man->SetDefaultStorage("local://OCDB");
2258 }
2259 fPHOSCalibData = new AliPHOSCalibData();
2260 }
2261}
2262
2263//_____________________________________________________________________________
2264void AliAnalysisTaskPi0Flow::SetVertex()
2265{
2266 const AliVVertex *primaryVertex = fEvent->GetPrimaryVertex();
2267 if( primaryVertex ) {
2268 fVertex[0] = primaryVertex->GetX();
2269 fVertex[1] = primaryVertex->GetY();
2270 fVertex[2] = primaryVertex->GetZ();
2271 }
2272 else {
2273 AliError("Event has 0x0 Primary Vertex, defaulting to origo");
2274 fVertex[0] = 0;
2275 fVertex[1] = 0;
2276 fVertex[2] = 0;
2277 }
2278 fVertexVector = TVector3(fVertex);
2279 FillHistogram("hZvertex", fVertexVector.z(), fInternalRunNumber-0.5);
f5e08ac9 2280
cdd3a269 2281 if( fDebug >= 2 )
f5e08ac9 2282 AliInfo(Form("Vertex is set to (%.1f,%.1f,%.1f)", fVertex[0], fVertex[1], fVertex[2]));
ada6b882 2283
2284 fVtxBin=0 ;// No support for vtx binning implemented.
2285}
2286
2287//_____________________________________________________________________________
2288Bool_t AliAnalysisTaskPi0Flow::RejectEventVertex()
2289{
2290 if( ! fEvent->GetPrimaryVertex() )
2291 return true; // reject
c9fb2807 2292 LogSelection(1, fInternalRunNumber);
2293
ada6b882 2294 if ( TMath::Abs(fVertexVector.z()) > 10. )
2295 return true; // reject
c9fb2807 2296 LogSelection(2, fInternalRunNumber);
ada6b882 2297
2298 return false; // accept event.
2299}
2300
2301//_____________________________________________________________________________
2302void AliAnalysisTaskPi0Flow::SetCentrality()
2303{
2304 AliCentrality *centrality = fEvent->GetCentrality();
2305 if( centrality )
2306 fCentralityV0M=centrality->GetCentralityPercentile("V0M");
2307 else {
2308 AliError("Event has 0x0 centrality");
2309 fCentralityV0M = -1.;
2310 }
2311 FillHistogram("hCentrality",fCentralityV0M,fInternalRunNumber-0.5) ;
2312
2313 fCentBin = GetCentralityBin(fCentralityV0M);
f5e08ac9 2314
cdd3a269 2315 if ( fDebug >= 2 )
f5e08ac9 2316 AliInfo(Form("Centrality (bin) is: %f (%d)", fCentralityV0M, fCentBin));
ada6b882 2317}
2318
2319//_____________________________________________________________________________
2320Bool_t AliAnalysisTaskPi0Flow::RejectCentrality()
2321{
2322 if( ! fEvent->GetCentrality() )
2323 return true; // reject
c9fb2807 2324 LogSelection(3, fInternalRunNumber);
2325
ada6b882 2326// if( fCentralityV0M <= 0. || fCentralityV0M>80. )
2327// return true; // reject
2328
2329 int lastBinUpperIndex = fCentEdges.GetSize() -1;
f5e08ac9 2330 if( fCentralityV0M > fCentEdges[lastBinUpperIndex] ) {
2331 if( fDebug )
2332 AliInfo("Rejecting due to centrality outside of binning.");
ada6b882 2333 return true; // reject
f5e08ac9 2334 }
c9fb2807 2335 LogSelection(4, fInternalRunNumber);
2336
f5e08ac9 2337 if( fCentralityV0M < fCentEdges[0] ) {
2338 if( fDebug )
2339 AliInfo("Rejecting due to centrality outside of binning.");
ada6b882 2340 return true; // reject
f5e08ac9 2341 }
c9fb2807 2342 LogSelection(5, fInternalRunNumber);
ada6b882 2343
2344 return false;
2345}
2346
2347
2348//_____________________________________________________________________________
2349void AliAnalysisTaskPi0Flow::EvalReactionPlane()
2350{
2351 // assigns: fHaveTPCRP and fRP
2352 // also does a few histogram fills
2353
2354 AliEventplane *eventPlane = fEvent->GetEventplane();
f5e08ac9 2355 if( ! eventPlane ) { AliError("Event has no event plane"); return; }
2356
ada6b882 2357 Double_t reactionPlaneQ = eventPlane->GetEventplane("Q");
2358 FillHistogram("phiRP",reactionPlaneQ,fCentralityV0M) ;
2359
2360 if(reactionPlaneQ==999 || reactionPlaneQ < 0.){ //reaction plain was not defined
f5e08ac9 2361 if( fDebug ) AliInfo(Form("No Q Reaction Plane, value is %f", reactionPlaneQ));
2362 fHaveTPCRP = kFALSE;
ada6b882 2363 }
2364 else{
cdd3a269 2365 if( fDebug >= 2 ) AliInfo(Form("Q Reaction Plane is %f", reactionPlaneQ));
f5e08ac9 2366 fHaveTPCRP = kTRUE;
ada6b882 2367 }
f5e08ac9 2368
ada6b882 2369 if(fHaveTPCRP){
e945b860 2370 fRP = ApplyFlattening(reactionPlaneQ, fCentralityV0M) ;
ada6b882 2371
2372 while(fRP<0) fRP+=TMath::Pi();
2373 while(fRP>TMath::Pi()) fRP-=TMath::Pi();
2374 FillHistogram("phiRPflat",fRP,fCentralityV0M) ;
2375 Double_t dPsi = eventPlane->GetQsubRes() ;
2376 FillHistogram("cos2AC",TMath::Cos(2.*dPsi),fCentralityV0M) ;
2377 }
f5e08ac9 2378 else
2379 fRP=0.;
ada6b882 2380}
2381
2382
2cde7444 2383//____________________________________________________________________________
ada6b882 2384void AliAnalysisTaskPi0Flow::EvalV0ReactionPlane(){
2385 // set: fRPV0A and fRPV0C
2cde7444 2386
6a1a4e92 2387 // Do Manual V0 EP Calculation
2388 if ( fManualV0EPCalc )
2389 {
2390 //VZERO data
2391 AliVVZERO* v0 = fEvent->GetVZEROData();
2392
2393 //reset Q vector info
2394 Double_t Qxa2 = 0, Qya2 = 0;
2395 Double_t Qxc2 = 0, Qyc2 = 0;
2396
2397 for (Int_t iv0 = 0; iv0 < 64; iv0++) {
2398 Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
2399 Float_t multv0 = v0->GetMultiplicity(iv0);
2400 if (iv0 < 32){ // V0C
2401 Qxc2 += TMath::Cos(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
2402 Qyc2 += TMath::Sin(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
2403 } else { // V0A
2404 Qxa2 += TMath::Cos(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
2405 Qya2 += TMath::Sin(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
2406 }
2407 }
2cde7444 2408
6a1a4e92 2409 Int_t iC = -1;
2410 // centrality bins
2411 if(fCentralityV0M < 5) iC = 0;
2412 else if(fCentralityV0M < 10) iC = 1;
2413 else if(fCentralityV0M < 20) iC = 2;
2414 else if(fCentralityV0M < 30) iC = 3;
2415 else if(fCentralityV0M < 40) iC = 4;
2416 else if(fCentralityV0M < 50) iC = 5;
2417 else if(fCentralityV0M < 60) iC = 6;
2418 else if(fCentralityV0M < 70) iC = 7;
2419 else iC = 8;
2420
2421 //grab for each centrality the proper histo with the Qx and Qy to do the recentering
2422 Double_t Qxamean2 = fMeanQ[iC][1][0];
2423 Double_t Qxarms2 = fWidthQ[iC][1][0];
2424 Double_t Qyamean2 = fMeanQ[iC][1][1];
2425 Double_t Qyarms2 = fWidthQ[iC][1][1];
2426
2427 Double_t Qxcmean2 = fMeanQ[iC][0][0];
2428 Double_t Qxcrms2 = fWidthQ[iC][0][0];
2429 Double_t Qycmean2 = fMeanQ[iC][0][1];
2430 Double_t Qycrms2 = fWidthQ[iC][0][1];
2431
2432 Double_t QxaCor2 = (Qxa2 - Qxamean2)/Qxarms2;
2433 Double_t QyaCor2 = (Qya2 - Qyamean2)/Qyarms2;
2434 Double_t QxcCor2 = (Qxc2 - Qxcmean2)/Qxcrms2;
2435 Double_t QycCor2 = (Qyc2 - Qycmean2)/Qycrms2;
2436
2437 fRPV0A = TMath::ATan2(QyaCor2, QxaCor2)/2.;
2438 fRPV0C = TMath::ATan2(QycCor2, QxcCor2)/2.;
2439 }
2440 else // Use Official V0 EP Calculation.
2441 {
2442 AliEventplane *eventPlane = fEvent->GetEventplane();
2443 if( ! eventPlane ) { AliError("Event has no event plane"); return; }
2444 fRPV0A = eventPlane->GetEventplane("V0A", fEvent);
2445 fRPV0C = eventPlane->GetEventplane("V0C", fEvent);
2446 }
2447
2448 // Check that the A&C RP are within allowed range.
2449 if( fDebug >= 3 && (fRPV0A<0 || fRPV0A>TMath::Pi() ) )
2450 AliInfo(Form("RPV0A outside of permited range [0,pi]: %f, correcting", fRPV0A));
2451 if( fDebug >= 3 && (fRPV0C<0 || fRPV0C>TMath::Pi() ) )
2452 AliInfo(Form("RPV0C outside of permited range [0,pi]: %f, correcting", fRPV0C));
6b1c1977 2453 while (fRPV0A<0 ) fRPV0A+=TMath::Pi() ;
2454 while (fRPV0A>TMath::Pi()) fRPV0A-=TMath::Pi() ;
2455 while (fRPV0C<0 ) fRPV0C+=TMath::Pi() ;
2456 while (fRPV0C>TMath::Pi()) fRPV0C-=TMath::Pi() ;
2457
2458 // Reaction plane histograms before flattening
2459 if( fDebug >= 2 )
2460 AliInfo(Form("V0 Reaction Plane before flattening: A side: %f, C side: %f", fRPV0A, fRPV0C));
2461
2462 FillHistogram("phiRPV0A" ,fRPV0A,fCentralityV0M);
2463 FillHistogram("phiRPV0C" ,fRPV0C,fCentralityV0M);
2464 FillHistogram("phiRPV0AC",fRPV0A,fRPV0C,fCentralityV0M) ;
2cde7444 2465
f5e08ac9 2466 // Flattening
e945b860 2467 fRPV0A=ApplyFlatteningV0A(fRPV0A,fCentralityV0M) ;
6b1c1977 2468 while (fRPV0A<0 ) fRPV0A+=TMath::Pi() ;
2469 while (fRPV0A>TMath::Pi()) fRPV0A-=TMath::Pi() ;
2470
e945b860 2471 fRPV0C=ApplyFlatteningV0C(fRPV0C,fCentralityV0M) ;
6b1c1977 2472 while (fRPV0C<0 ) fRPV0C+=TMath::Pi() ;
2473 while (fRPV0C>TMath::Pi()) fRPV0C-=TMath::Pi() ;
e945b860 2474
cdd3a269 2475 if( fDebug >= 2 )
6b1c1977 2476 AliInfo(Form("V0 Reaction Plane after flattening: A side: %f, C side: %f", fRPV0A, fRPV0C));
ada6b882 2477
2478 FillHistogram("phiRPV0Aflat",fRPV0A,fCentralityV0M) ;
ada6b882 2479 FillHistogram("cos2V0AC",TMath::Cos(2.*(fRPV0A-fRPV0C)),fCentralityV0M) ;
2480 if(fHaveTPCRP){
2481 FillHistogram("phiRPV0ATPC",fRP,fRPV0A,fCentralityV0M) ;
2482 FillHistogram("cos2V0ATPC",TMath::Cos(2.*(fRP-fRPV0A)),fCentralityV0M) ;
2483 }
2484
2485 FillHistogram("phiRPV0Cflat",fRPV0C,fCentralityV0M) ;
2486 if(fHaveTPCRP){
2487 FillHistogram("phiRPV0CTPC",fRP,fRPV0C,fCentralityV0M) ;
2488 FillHistogram("cos2V0CTPC",TMath::Cos(2.*(fRP-fRPV0C)),fCentralityV0M) ;
2489 }
2cde7444 2490}